[R] plot of Bernoulli data
Frank E Harrell Jr
fharrell at virginia.edu
Thu Oct 4 18:19:55 CEST 2001
Thanks for elaborating Greg - I see this more clearly now.
Trellis/lattice handles overlapping intervals that are
used to make multiple panels (using shingles). That
would be for a third variable. I don't think shingles
will help with what you are doing, so I see the
utility of your function more better now. Thanks -Frank
"Warnes, Gregory R" wrote:
>
> Hi Frank,
>
> Actually the 'wapply' function I have does 'overlapping' windows. (It would
> have to if it uses 1/5 of the range at each of 10 points, n'est-pas?) The
> main reason I created it was to add a 'local variance' estimate to plots in
> an attempt to more easily detect changes in variability. EG:
>
> x <- rnorm(500,sd=10)
> y <- rnorm(500,sd=sqrt(abs(x)))
> plot(x,y)
>
> sdplus <- function(x,sign=1)
> {
> if(length(x)>0)
> mean(x) + sign * sqrt(var(x))
> else
> NA
> }
> lines(wapply(x,y,fun=sdplus,width=1/5),col="red")
> lines(wapply(x,y,fun=sdplus,sign=-1,width=1/5),col="red")
>
> Not very pretty, but does what I need.
>
> I'm not familiar with the trellace/lattice functions. Is there a better or
> more elegant way to do this?
>
> -Greg
>
> > -----Original Message-----
> > From: Frank E Harrell Jr [mailto:fharrell at virginia.edu]
> > Sent: Thursday, October 04, 2001 9:06 AM
> > To: Warnes, Gregory R
> > Cc: 'Bill Simpson'; 'r-help at stat.math.ethz.ch'
> > Subject: Re: [R] plot of Bernoulli data
> >
> >
> > Greg - My model for this type of application is to
> > either use
> >
> > - shingle objects with Trellis/Lattice to allow overlapping
> > intervals (and for this application there is good reason
> > for intervals to overlap if you are not using the
> > better loess/lowess solution)
> > or
> > - the cut2 function in the soon (I hope) to be beta-released
> > Hmisc library, with tapply and all the other stratified
> > computation functions. cut2 has the same arguments as
> > cut but also arguments g (number of quantile groups) and
> > m (number of observations per interval)
> >
> > cut2 at present does not handle fixed windows as you
> > describe, but in general such windows are problematic.
> >
> > The bigger point though is that there may be a slight
> > advantage in separating the stratification code (e.g.,
> > cut2) from the computation code (e.g., tapply and
> > summarize (in Hmisc)).
> >
> > Frank Harrell
> >
> >
> > "Warnes, Gregory R" wrote:
> > >
> > > In the next release of the gregmisc library there is a
> > function "wapply"
> > > that applies a specified function over 'windows' of x-y
> > data (code below).
> > > You would use it in this case as:
> > >
> > > x<-sort(runif(100,1,20))
> > > p<-pnorm(x,10,3)
> > > y<-as.numeric(runif(x)<p)
> > > plot(x,y)
> > > lines(x,p)
> > >
> > > # draw a line with 10 points, each computed on regions
> > spanning 1/5th of
> > > the x range
> > > lines(wapply(x,y,n=10,width=1/5),col="red")
> > >
> > > -Greg
> > >
> > > # $Id: wapply.R,v 1.2 2001/08/31 23:45:45 warneg Exp $
> > > #
> > > "wapply" _ function( x, y, fun=mean, method="range",
> > > width=1/10, n=50, ...)
> > > {
> > > method <- match.arg(method,
> > c("width","range","nobs","fraction"))
> > >
> > > if(method == "width" || method == "range" )
> > > {
> > > if(method=="range")
> > > width <- width * range(x)
> > >
> > >
> > > pts _ seq(min(x),max(x),length.out=n)
> > >
> > > result _ sapply( pts, function(pts,y,width,fun,XX,...)
> > > {
> > > low _ min((pts-width/2),max(XX))
> > > high _ max((pts+width/2), min(XX))
> > > return (fun(y[(XX>= low) &
> > (XX<=high)],...))
> > > },
> > > y=y,
> > > width=width,
> > > fun=fun,
> > > XX = x,
> > > ...)
> > >
> > > return(x=pts,y=result)
> > > }
> > > else # method=="nobs" || method=="fraction"
> > > {
> > > if( method=="fraction")
> > > width <- floor(length(x) * width)
> > >
> > > ord <- order(x)
> > > x <- x[ord]
> > > y <- y[ord]
> > >
> > > n <- length(x)
> > > center <- 1:n
> > > below <- sapply(center - width/2, function(XX) max(1,XX) )
> > > above <- sapply(center + width/2, function(XX) min(n,XX) )
> > >
> > > retval <- list()
> > > retval$x <- x
> > > retval$y <- apply(cbind(below,above), 1,
> > > function(x) fun(y[x[1]:x[2]],...) )
> > >
> > > return(retval)
> > > }
> > >
> > > }
> > >
> > > > -----Original Message-----
> > > > From: Bill Simpson [mailto:wsi at gcal.ac.uk]
> > > > Sent: Tuesday, October 02, 2001 7:57 AM
> > > > To: r-help at stat.math.ethz.ch
> > > > Subject: Re: [R] plot of Bernoulli data
> > > >
> > > >
> > > > > x<-sort(runif(100,1,20))
> > > > > p<-pnorm(x,10,3)
> > > > > y<-as.numeric(runif(x)<p)
> > > > > plot(x,y)
> > > > > lines(x,p)
> > > > df<-data.frame(x,y)
> > > > aggregate(df,list(x=(x<5),(x>5)&(x<10),(x>10) &
> > > > (x<15),(x>15)), FUN=mean)
> > > > gives me what I want but if anyone has a better way
> > to collect the
> > > > observations into bins I'd like to hear it. It would
> > be nice to
> > > > pass along something like
> > > > breaks<-c(5,10,15,20)
> > > >
> > > > Thanks
> > > >
> > > > Bill Simpson
> > > >
> > > > -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.
> > > > -.-.-.-.-.-.-.-.-.-
> > > > r-help mailing list -- Read
> > > > http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
> > > > Send "info", "help", or "[un]subscribe"
> > > > (in the "body", not the subject !) To:
> > > > r-help-request at stat.math.ethz.ch
> > > > _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._.
> > > > _._._._._._._._._._
> > > >
> > >
> > > LEGAL NOTICE
> > > Unless expressly stated otherwise, this message is
> > confidential and may be privileged. It is intended for the
> > addressee(s) only. Access to this E-mail by anyone else is
> > unauthorized. If you are not an addressee, any disclosure
> > or copying of the contents of this E-mail or any action
> > taken (or not taken) in reliance on it is unauthorized and
> > may be unlawful. If you are not an addressee, please inform
> > the sender immediately.
> > >
> > -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.
> > -.-.-.-.-.-.-.-.-.-
> > > r-help mailing list -- Read
> http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
> > Send "info", "help", or "[un]subscribe"
> > (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch
> >
> _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._.
> _._
>
> --
> Frank E Harrell Jr Prof. of Biostatistics & Statistics
> Div. of Biostatistics & Epidem. Dept. of Health Evaluation Sciences
> U. Virginia School of Medicine http://hesweb1.med.virginia.edu/biostat
> -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.
> -.-
> r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
> Send "info", "help", or "[un]subscribe"
> (in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch
> _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._.
> _._
--
Frank E Harrell Jr Prof. of Biostatistics & Statistics
Div. of Biostatistics & Epidem. Dept. of Health Evaluation Sciences
U. Virginia School of Medicine http://hesweb1.med.virginia.edu/biostat
-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !) To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
More information about the R-help
mailing list