[R] 95% confidence interval of the coefficients from a bootstrap analysis

Rui Barradas ruipbarradas at sapo.pt
Tue May 1 15:41:29 CEST 2012


Hello,


baconbeach wrote
> 
> Thanks again for your swift response!!
> 
> With your last line, I get 
> 
>> rowMeans(sapply(stor.confint, colMeans)) 
>     2.5 %    97.5 % 
> 0.3256882 0.4604677 
> 
> I need the values (2.5% and 97.5%) for each variable of my model. I don't
> think this what I am getting.
> 
> This is what my script looks like now, after your help:
> 
> N = length (data_Pb[,1])
> B = 10000
> stor.r2 = rep(0,B)
> stor.coeffs <- vector("list", B) 
> stor.confint <- vector("list", B) 
> 
> for (i in 1:B){
> idx = sample(1:N, replace=T)
> newdata = data_Pb[idx,]
> L_NPRI_25k <- log(newdata$NPRI_25k+1)
> data_Pb.boot = lm(newdata$Log_Level ~   
>     newdata$Ind_5k + newdata$MineP_50k + 
>   newdata$NPRI_10k + L_NPRI_25k )
> 
> stor.r2[i] = summary(data_Pb.boot)$r.squared
> stor.coeffs [[i]] <- coef(data_Pb.boot) 
> stor.confint[[i]] <- confint(data_Pb.boot) 
> }
> 
> hist(stor.r2, xlab="R-squared",main="Distribution of R-squared - Lead
> (log)")
> summary(stor.r2)
> rowMeans(sapply(stor.confint, colMeans)) 
> rowMeans(sapply(stor.coeffs, colMeans))
> 
> 
> Thanks 
> 
> Steeve
> 

Ok, my mistake. Bootstrapped values for each coefficient, obviously.
Try, after the loop,

# Transform list into matrix, before applying 'colMeans'
stor.coeffs <- do.call(rbind, stor.coeffs)
colMeans(stor.coeffs)

# The same, but in two steps
Q2.5 <- lapply(stor.confint, function(x) x[ , 1])
Q2.5 <- do.call(rbind, Q2.5)
head(Q2.5)
Q97.5 <- lapply(stor.confint, function(x) x[ , 2])
Q97.5 <- do.call(rbind, Q97.5)
head(Q97.5)

colMeans(Q2.5)
colMeans(Q97.5)
# Maybe this is better
stor.ci <- data.frame(Q2.5=colMeans(Q2.5), Q97.5=colMeans(Q97.5))
rm(Q2.5, Q97.5) # Final clean-up

Maybe you could reuse the name 'stor.confint', it would save some memory.

Rui Barradas


--
View this message in context: http://r.789695.n4.nabble.com/95-confidence-interval-of-the-coefficients-from-a-bootstrap-analysis-tp4599692p4600737.html
Sent from the R help mailing list archive at Nabble.com.



More information about the R-help mailing list