[R] Bootstraping GAMs for confidence intervales calculation

From: Luis Reino <luisreino_at_isa.utl.pt>
Date: Tue, 29 Jul 2008 17:11:19 +0100


Dear R-Users,

I am resending this message just to reminder my question regarding the calculation of a bootstrap confidence intervals for a GAM plot.

I am trying to apply a bootstrap to a GAM in order to calculate the 95% confidence intervals for a smooth curve obtained by the “plot.gam” function of the mgcv package. Nonetheless, I am getting some difficulties in transposing the results for the graphs.

I used the following commands in R, “mgcv” and “boot” packages:

*> attach(bbvc_11Jul08)*

*> model.boot<-function(data,indices){*

*+ sub.data<-data[indices,]*

*+
model<-gam(asin(sqrt(Cov0_30))~s(age,k=4,bs="cr"),family=gaussian,data=sub.data)*

*+ coef(model)}*

*> gam.boot<-boot(bbvc_11Jul08,model.boot,R=1000)*

*> gam.boot*

ORDINARY NONPARAMETRIC BOOTSTRAP Call:

boot(data = bbvc_11Jul08, statistic = model.boot, R = 1000)

Bootstrap Statistics :

original bias std. error

t1* 0.40370112 0.0002762674 0.02017378

t2* 0.04715501 -0.0144132280 0.07648348

t3* 0.14590936 -0.0135820501 0.05161244

t4* 0.13520098 -0.0096405733 0.04470793

I obtained 4 indexes in the “bootstrap” output. I know that for a lm function, with an independent variable, t1 is the index of the intercept and t2 is the index of the variable, but for GAMs I don't find what each indexes really mean!

Moreover, to calculate the confidence intervals of 95% I proceed as follows:

*>boot.ci(gam.boot,conf=0.95) *

BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS Based on 1000 bootstrap replicates

CALL : boot.ci(boot.out = gam.boot, conf = 0.95)

Intervals :

Level Normal Basic Studentized

95% ( 0.3639, 0.4430 ) ( 0.3635, 0.4434 ) ( 0.3583, 0.4507 )

Level Percentile BCa

95% ( 0.3640, 0.4439 ) ( 0.3623, 0.4434 )

Calculations and Intervals on Original Scale

Warning message:

In sqrt(tv[, 2]) : NaNs produced

My first doubt is when I make "boot.ci", which of the four should I choose? By default, R calculates "index=1:min(2,length(boot.out$t0)" but I don't know if this is the correct form to do it.

Finally, I don’t know how to introduce the 95% confidence intervals in the "plot.gam".

I hope you can help me. Thanks in advance.

Luís Reino

-- 
Luis Reino, PhD Student
Centro de Estudos Florestais
Departamento de Engenharia Florestal
Instituto Superior de Agronomia
Universidade Técnica de Lisboa
1349-017 Lisboa
Portugal

______________________________________________
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.
Received on Tue 29 Jul 2008 - 16:14:25 GMT

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.2.0, at Tue 29 Jul 2008 - 16:32:47 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.

list of date sections of archive