[R] Confidence bands around LOESS

Bert Gunter gunter.berton at gene.com
Wed Aug 8 18:22:23 CEST 2012


Please Note:

loess() is a fitting algorithm, so no parameters are estimated and no
_confidence_ intervals nor bands can be produced.

What predict.loess() _does_  give is standard errors for the
fits/predictions, and you can add and subtract as many of these as you
like to produce standard error bands. This sort of vague "uncertainty
interval" is presumably what is wanted anyway, but I do think it is
important to make sure you get the nomenclature right --  and, in
particular, no probabilistic claims are made about true underlying
curves lying within the bands produced.

Cheers,
Bert

On Wed, Aug 8, 2012 at 8:51 AM, John Kane <jrkrideau at inbox.com> wrote:
> I may have missed something entirely but I think ggplot2 will do this for you pretty well automatically
>
> library(ggplot2)
> mydata <- read.csv("http://doylesdartden.com/smoothing.csv", sep=",")
>
> p  <-  ggplot(mydata , aes( X_Axis_Parameter, Y_Axis_Parameter  )) + geom_point() +
>              geom_smooth()
> p
>
>
> John Kane
> Kingston ON Canada
>
>
>> -----Original Message-----
>> From: kydaviddoyle at gmail.com
>> Sent: Tue, 7 Aug 2012 21:22:41 -0500
>> To: r-help at r-project.org
>> Subject: [R] Confidence bands around LOESS
>>
>> Hi Folks,
>>
>> I'm looking to do Confidence bands around LOESS smoothing curve.
>>
>> If found the older post about using the Standard error to approximate it
>> https://stat.ethz.ch/pipermail/r-help/2008-August/170011.html
>>
>> Also found this one
>> http://www.r-bloggers.com/sab-r-metrics-basics-of-loess-regression/
>>
>> But they both seem to be approximations of confidence intervals and I was
>> wonder if there was a way to get the CIs?
>>
>> Below is the code I have so far and my data is no the net.
>>
>> Any help would be greatly appreciated.
>>
>> Take Care
>> David
>> -----------------------------
>>
>> #Load your data.  Is located on the web at the address below
>>
>> mydata <- read.csv("http://doylesdartden.com/smoothing.csv", sep=",")
>>
>> mydata <- read.table("x.csv", header=TRUE, sep=",",)
>>
>>
>>
>> attach(mydata)
>>
>> reg1 <- lm(Y_Axis_Parameter~X_Axis_Parameter)
>>
>> par(cex=1)
>>
>> * *
>>
>> * *
>>
>> #Plots the data but makes nondetects a different color and type based on
>> column D_Y_Axis_Parameter being a 0 for ND and 1 for detect.
>>
>> plot(X_Axis_Parameter, Y_Axis_Parameter, col=ifelse(D_Y_Axis_Parameter,
>> "black", "red"),ylab = "Y_Axis_Parameter", pch=ifelse(D_Y_Axis_Parameter,
>> 19, 17), cex = 0.7)
>>
>>
>>
>> plx<-predict(loess(Y_Axis_Parameter ~ X_Axis_Parameter, data=mydata),
>> se=T)
>>
>>
>>
>>
>>
>> lines(mydata$X_Axis_Parameter,plx$fit+2*plx$s, lty=2) #rough & ready CI
>>
>> lines(mydata$X_Axis_Parameter,plx$fit-2*plx$s, lty=2)
>>
>>
>>
>> # Apply loess smoothing using the default span value of 0.8.  You can
>> change the curve by changing the span value.
>>
>> y.loess <- loess(y ~ x, span=0.8, data.frame(x=X_Axis_Parameter,
>> y=Y_Axis_Parameter))
>>
>>
>>
>> # Compute loess smoothed values for all points along the curve
>>
>> y.predict <- predict(y.loess, data.frame(x=X_Axis_Parameter))
>>
>>
>>
>> # Plots the curve.
>>
>> lines(X_Axis_Parameter,y.predict)
>>
>> * *
>>
>> #Add Legend to graY_Axis_Parameter.  You can change the size of the box
>> by
>> changing cex = 0.75  Large # makes it larger.
>>
>> legend("topleft", c("Smoothing Curve", "Detected", "NonDetect"), col =
>> c(1,
>> "black","red"), cex = 0.75,
>>
>>        text.col = "black", lty = c(1 ,-1, -1), pch = c(-1, 19, 17),
>>
>>        merge = TRUE, bg = 'gray90')
>>
>> * *
>>
>> #Add title
>>
>> title(main="Locally Weighted Scatterplot Smoothing Curve")
>>
>>
>>
>> # Done
>>
>>       [[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide
>> http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>
> ____________________________________________________________
> GET FREE SMILEYS FOR YOUR IM & EMAIL - Learn more at http://www.inbox.com/smileys
> Works with AIM®, MSN® Messenger, Yahoo!® Messenger, ICQ®, Google Talk™ and most webmails
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.



-- 

Bert Gunter
Genentech Nonclinical Biostatistics

Internal Contact Info:
Phone: 467-7374
Website:
http://pharmadevelopment.roche.com/index/pdb/pdb-functional-groups/pdb-biostatistics/pdb-ncb-home.htm



More information about the R-help mailing list