possible bug in line+group/color+facet_wrap used for plotting ECDFs

[Szilard --]

Hi:

I'd like to plot empirical cummulative distribution functions for a number of variables (20)

and groups (2), but the steps/lines I get in the plots are not always monotonic. Here is a

simplified example to illustrate the behavior:

library(ggplot2)

d=structure(list(gr = c("a", "b", "b", "b", "b", "a", "a", "b", 
"b", "a", "b", "a", "b", "b", "a", "a", "b", "b", "a", "a", "b", 
"a", "b", "a", "b", "b", "b", "a", "b", "a", "a", "a", "b", "b", 
"b", "a", "a", "b", "b", "b", "b", "b", "b", "b", "b", "b", "b", 
"a", "b", "a", "a", "a", "a", "b", "b", "a", "b", "a", "a", "a", 
"a", "a", "b", "a", "a", "a", "a", "a", "b", "b", "b", "b", "b", 
"a", "a", "a", "a", "b", "a", "b", "b", "b", "a", "b", "b", "a", 
"a", "a", "a", "a", "b", "a", "a", "a", "a", "b", "b", "b", "a", 
"b"), variable = c("V2", "V1", "V2", "V2", "V2", "V1", "V2", 
"V1", "V2", "V1", "V1", "V1", "V2", "V1", "V2", "V2", "V2", "V1", 
"V1", "V1", "V1", "V1", "V1", "V2", "V1", "V2", "V1", "V2", "V1", 
"V2", "V1", "V2", "V1", "V2", "V1", "V1", "V2", "V1", "V1", "V2", 
"V1", "V1", "V2", "V2", "V1", "V1", "V1", "V2", "V2", "V1", "V2", 
"V2", "V2", "V1", "V2", "V2", "V1", "V1", "V2", "V1", "V1", "V2", 
"V2", "V2", "V2", "V2", "V2", "V2", "V2", "V2", "V1", "V2", "V2", 
"V2", "V2", "V1", "V2", "V2", "V2", "V1", "V1", "V2", "V2", "V2", 
"V1", "V1", "V2", "V2", "V2", "V1", "V1", "V1", "V1", "V1", "V1", 
"V1", "V1", "V1", "V1", "V1"), value = c(0, 0.5, 0, 0, 0, 1, 
0, 1, 0, 1, 1, 0.001, 0, 0.5, 0, 0, 0, 0.6, 0.2, 0.75, 0.3, 0.5, 
0.429, 0.333, 0.019, 0, 0.333, 0, 0.75, 0, 0.6, 0.062, 0.15, 
0.064, 0.333, 0.6, 0, 0.3, 0.231, 0, 0.75, 0.188, 0, 0, 0.15, 
0.5, 0.6, 0, 0, 0.429, 0, 0, 0, 0.5, 1, 0, 1, 0.75, 0, 1, 0.273, 
0, 0, 0, 0, 0, 0.05, 0, 0, 0, 0.111, 1, 0.429, 0, 0, 0.15, 0, 
0, 0, 0.75, 0.75, 0, 0, 0.111, 0.231, 1, 0, 0, 0, 1, 0.5, 1, 
0.6, 1, 0.3, 1, 0.75, 0.176, 0.094, 0.035)), .Names = c("gr", 
"variable", "value"), row.names = c(6084L, 1706L, 7789L, 7830L, 
4624L, 3313L, 4877L, 562L, 4856L, 229L, 887L, 846L, 4939L, 352L, 
6860L, 5809L, 5242L, 2888L, 2395L, 2973L, 1650L, 3085L, 3955L, 
7664L, 3485L, 7098L, 935L, 6447L, 2298L, 5953L, 762L, 6702L, 
1170L, 6203L, 957L, 1106L, 5952L, 3484L, 1473L, 4208L, 717L, 
3021L, 5760L, 5026L, 3968L, 3002L, 2282L, 6465L, 4244L, 1526L, 
4325L, 7060L, 6861L, 2852L, 8014L, 4525L, 2859L, 459L, 5092L, 
3802L, 250L, 5451L, 7752L, 5942L, 7031L, 4910L, 5984L, 5419L, 
5762L, 6171L, 581L, 7742L, 6624L, 5827L, 6706L, 2391L, 7473L, 
5497L, 5402L, 2281L, 3204L, 4243L, 5168L, 5296L, 1676L, 1909L, 
6038L, 6812L, 7454L, 1749L, 2584L, 2015L, 1358L, 1371L, 1112L, 
574L, 3201L, 2825L, 1509L, 1452L), class = "data.frame")

dd=ddply(d, .(gr, variable), function(df) {
  df=df[order(df$value),]
  df$y=(1:nrow(df))/nrow(df)
  df=rbind(data.frame(gr=df$gr[1], variable=df$variable[1], value=min(df$value), y=0), df)
  df
}) 

qplot(value, y, data=dd, geom='step', colour=gr) + facet_wrap(~variable)  

The plot I obtained is attached here (fwrap-problem.png). As you can see the

ECDF looks strange in some points where there are repetitions of the same value in

the sample. I wasn't able to replicate the behavior on a smaller artificial example.

The same problem arises if geom='line' or 'path' is used (looking at 'path' might give a 

hint of what's wrong).

The problem disappears if facet_grid is used (see attached fgrid-OK.png), but since

in my real world data I have 20 variables, using that is not really satisfactory.

Any ideas how to fix this problem or get around?

Thanks,

Szilard

PS: Here is my sessionInfo():

R version 2.9.1 (2009-06-26) 
i386-apple-darwin8.11.1 

locale:
en_US.UTF-8/en_US.UTF-8/C/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] grid stats graphics grDevices utils datasets methods base  

other attached packages:
[1] ggplot2_0.8.3 reshape_0.8.3 plyr_0.1.9 proto_0.3-8  

[David Kahle]

Hi Szilard -

Try this with the d variable you provided :

  df  <- d

  n   <- 201

  rng <- extendrange(df$value)

  x   <- seq(rng[1], rng[2], length = n)

  

  ecdf_wrapper <- function(df){

  data.frame(

    gr = df$gr[1],

    variable = df$variable[1],

    x = x,

    Fx = ecdf(df$value)(x)

  )

  }

  df <- ddply(df, .(gr, variable), ecdf_wrapper)

  

  qplot(x, Fx, data = df, geom = 'step', colour = gr, facets = . ~ variable)

Is that what you're looking for?

Let me know how it goes

david.

<fwrap-problem.png><fgrid-OK.png>

[Szilard]

Thanks David, this is indeed what I wanted.

It is still strange why the original code doesn't work, possible bug?

Szilard

[hadley wickham]

> It is still strange why the original code doesn't work, possible bug?

I'm pretty sure it's a bug in facet_wrap - I use merge to correctly
line up data in different facets and I forgot to take into account
that it messes with the order of the original data.

Hadley

--
http://had.co.nz/

[Szilard]

> I'm pretty sure it's a bug in facet_wrap - I use merge to correctly
> line up data in different facets and I forgot to take into account
> that it messes with the order of the original data.

Hadley: Thanks a lot for having a look. I have the same feeling as if
I used geom='path' I could see that the order is messed up. Kind of
strange though, as ddply orders the data by 'gr' and 'variable'.
(Also
note that the problem doesn't show up in facet_grid.)

Szilard

[learnr]

Another alternative:

dp <- ddply(d, .(gr, variable), transform, ecdf = ecdf(value)(value))
qplot(value, ecdf, data=dp, geom='step', colour=gr) + facet_wrap
(~variable)

--
http://learnr.wordpress.com

[Szilard]

Thanks LearnR for very efficient/compact solution.

However, some of my variables have a considerable part of their
distribution at the lowest value (at 0), and it would be nice therefore
for the ECDF line to start at (0,0) (actually (xmin,0) more generally).
I came up with inserting this code before the plot:
dp <- ddply(dp, .(gr, variable), function(df) {

rbind(data.frame(gr=df$gr[1],
variable=df$variable[1],
value=min(df$value),

ecdf=0), df)
})
but I wonder if there is a better way to do that.

Szilard

[hadley wickham]

On Wed, Aug 12, 2009 at 11:36 AM, Szilard<szilard.ma...@gmail.com> wrote:
>
> Thanks LearnR for very efficient/compact solution.
>
> However, some of my variables have a considerable part of their
> distribution at the lowest value (at 0), and it would be nice therefore
> for the ECDF line to start at (0,0) (actually (xmin,0) more generally).
> I came up with inserting this code before the plot:
> dp <- ddply(dp, .(gr, variable), function(df) {
>  rbind(data.frame(gr=df$gr[1],
>                   variable=df$variable[1],
>                   value=min(df$value),
>                   ecdf=0), df)
> })
> but I wonder if there is a better way to do that.

I'm not entirely sure what you want. The following code would ensure
that you always get a 0 on the x-axis (as well as being a little more
efficient), but your sample data already has zeros.

dp <- ddply(d, .(gr, variable), function(df) {
xs <- unique(c(0, df$value))
data.frame(value = xs, ecdf = ecdf(df$value)(xs))
})

Hadley

--
http://had.co.nz/

[Szilard]

Thanks Hadley. What I actually wanted is to ensure 0 on the y-axis, not
the x-axis as you assumed.

In a simpler example with no groups/facets:

d = data.frame(x=2:7)
dd = transform(d, ecdf=ecdf(x)(x))

I'm referring to the difference between this:

qplot(x, ecdf, data=dd, geom='step') + ylim(c(0,1))

and this:

dd2 = rbind( data.frame(x=min(d$x), ecdf=0), dd )
qplot(x, ecdf, data=dd2, geom='step') + ylim(c(0,1))

I would prefer the 2nd, where the steps start from y=0.

From this point of view the previous solution given by David fits better:

resolution = 1000
rng = extendrange(d$x)
xsteps = seq(rng[1], rng[2], length = resolution)
dd3 = data.frame( x = xsteps, ecdf = ecdf(d$x)(xsteps) )

qplot(x, ecdf, data=dd3, geom='step') + ylim(c(0,1))

However, in this case setting resolution very low would distort the
plot, while having it very high might slow down the plot in case of many
variables (groups/facets).

Szilard

[hadley wickham]

On Thu, Aug 13, 2009 at 5:47 PM, Szilard<szilard.ma...@gmail.com> wrote:
>
> Thanks Hadley. What I actually wanted is to ensure 0 on the y-axis, not
> the x-axis as you assumed.

Ahhh, ok. In that case I think you want this:

dp <- ddply(d, .(gr, variable), function(df) {

xs <- unique(df$value)
data.frame(value = c(min(xs), xs), ecdf = c(0, ecdf(df$value)(xs)))
})

> resolution = 1000
> rng = extendrange(d$x)
> xsteps = seq(rng[1], rng[2], length = resolution)
> dd3 = data.frame( x = xsteps, ecdf = ecdf(d$x)(xsteps) )
>
> qplot(x, ecdf, data=dd3, geom='step') + ylim(c(0,1))
>
> However, in this case setting resolution very low would distort the
> plot, while having it very high might slow down the plot in case of many
> variables (groups/facets).

Using a grid with an ecdf isn't terribly appealing since you know
where all the changes are going to be.

Hadley

--
http://had.co.nz/

[Szilard]

Thanks Hadley. Indeed, this is probably the easiest way to do it. Thank
you to all that have contributed to this thread.

Szilard

hadley wickham wrote:
> > Ahhh, ok. In that case I think you want this:
> >
> > dp <- ddply(d, .(gr, variable), function(df) {
> > xs <- unique(df$value)
> > data.frame(value = c(min(xs), xs), ecdf = c(0, ecdf(df$value)(xs)))
> > })
> >