[R] mcr mcreg: Why are BCa and quantile se values calculated, stored, then the stored values set to NA?

[Fabian Model]

This piece of code is indeed confusing. Generally quantile and BCa 
bootstrap do not estimate a global SE, so the glob.sigma slot you want to 
access is not really meaningful for those methods.
 
The result object of a BCa bootstrap as calculated by the mc.bootstrap 
function will contain a glob.sigma slot that is set to the analytical SE 
estimates of the complete data set for regression methods where analytical 
estimates of SE are available (LinReg, WLinReg, Deming). For the other 
regression methods (WDeming, PaBa) glob.sigma will already be NA. Since the 
analytical SE estimates are not relevant when you perform quantile or BCa 
bootstrap we decided at some point to always set them NA for these 
bootstrap methods - independent of regression method. The code copying the 
glob.sigma results is pointless and we'll remove it in the next release.
 
Now, if you want to calculate the mean or the sd of the individual 
bootstrap estimates, you can directly access the MCResultResampling and 
MCResultBCa objects. They store the estimates of regression coefficients 
for each bootstrap sample in slots B0 (intercept) and B1 (slope). Example:

library(mcr)
data(creatinine)
m <- mcreg(as.matrix(creatinine), method.bootstrap.ci = "quantile", na.rm = 
TRUE)
mean.slope.est<- mean(m at B1)
sd.slope.est<- sd(m at B1)
mean.intercept.est <- mean(m at B0)
sd.intercept.est <- sd(m at B1)
 
However, note that these estimates have no direct relationship to the BCa 
confidence intervals.
Best regards,
Fabian
 

On Tuesday, February 25, 2014 3:54:11 PM UTC+1, Peter Crowther wrote:

> In mcreg, there's the following snippet for bootstrap CI method of BCa 
> (and a very similar one for quantile): 
>
> else if (method.bootstrap.ci == "BCa") { 
>             bootB0 <- mc.calc.bca(Xboot = B0, Xjack = B0jack, 
>                 xhat = glob.coef[1], alpha) 
>             bootB1 <- mc.calc.bca(Xboot = B1, Xjack = B1jack, 
>                 xhat = glob.coef[2], alpha) 
>             bootB0$se <- glob.sigma[1] 
>             bootB1$se <- glob.sigma[2] 
>             bootB0$se <- NA 
>             bootB1$se <- NA 
>         } 
> [...] 
>         mc.res <- mc.make.CIframe(b0 = glob.coef[1], b1 = glob.coef[2], 
>             se.b0 = bootB0$se, se.b1 = bootB1$se, CI.b0 = bootB0$CI, 
>             CI.b1 = bootB1$CI) 
>
> I'd quite like to use in my code the SE values that were calculated in 
> mcreg.  Can any kind soul shed light on why they're calculated, 
> stored, then the stored values set to NA? 
>
> Thanks in advance for any insight! 
>
> - Peter 
> -- 
> Peter Crowther, Director, Melandra Limited 
>
> ______________________________________________ 
> R-h... at r-project.org <javascript:> 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. 
>