[R] Bca confidence intervals for Pearson coefficient, thanks.

varin sacha varinsacha at yahoo.fr
Sun Dec 21 21:02:58 CET 2014


Dear Professor FOX,
I really thank you lots for all your precisions. One last precision, now if I want tot calculate the BCa bootstrap CIs for the Cramer's V and the Eta-squared.
## Read in the data file, using headers and Tab separator> test = read.table(file.choose(), header = TRUE, sep = "\t")
## Check the variable names> names(test)[1] "gender"    "option"    "math.test" "geo.test"  "shopping"  "sports"   
## Cramer's V
> library(questionr)> tab = table(test$gender, test$option)> cramer.v(tab)                         [1] 0.1490712                           
## Eta.square
> library(lsr)> test.aov = aov(math.test ~ gender, data = test)
> etaSquared(test.aov)          eta.sq eta.sq.partgender 0.1207154   0.1207154

Best Regards, looking forward to reading you once more,
Sacha
      De : John Fox <jfox at mcmaster.ca>
 À : varin sacha <varinsacha at yahoo.fr> 
Cc : r-help help <r-help at r-project.org> 
 Envoyé le : Dimanche 21 décembre 2014 5h33
 Objet : Re: [R] Bca confidence intervals for Pearson coefficient, thanks.
   
Dear varin sacha,

I think that you misunderstand how boot() and boot.ci() work. The boot() function in the simplest case takes two arguments, for the data and indices into the data, while boot.ci() takes as its principal argument the object returned by boot(). All of this seems reasonably clear in ?boot and ?boot.ci.

Here's an example with different data (since as far as I can see you didn't supply yours):

------------- snip --------

> library(boot)
> 
> x <- longley$Year
> y <- longley$Population
> 
> cor(cbind(x, y))
          x        y
x 1.0000000 0.9939528
y 0.9939528 1.0000000
> 
> myCor <- function(data, index){
+  cor(data[index, ])[1, 2]
+ }
> 
> set.seed(12345)
> (b <- boot(data=cbind(x, y), statistic=myCor, R=200))

ORDINARY NONPARAMETRIC BOOTSTRAP


Call:
boot(data = cbind(x, y), statistic = myCor, R = 200)


Bootstrap Statistics :
    original      bias    std. error
t1* 0.9939528 0.0008263766 0.001850004

> boot.ci(b, type="bca")
BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
Based on 200 bootstrap replicates

CALL : 
boot.ci(boot.out = b, type = "bca")

Intervals : 
Level      BCa          
95%  ( 0.9895,  0.9969 )  
Calculations and Intervals on Original Scale
Warning : BCa Intervals used Extreme Quantiles
Some BCa intervals may be unstable
Warning message:
In norm.inter(t, adj.alpha) : extreme order statistics used as endpoints

--------------- snip ---------------

Note that 200 bootstrap replications are generally sufficient for bootstrap standard errors (a normal-theory CI would be a poor choice here, unless you transform the correlation coefficient), but really aren't enough for a BCa interval.

I hope this helps,
 John



------------------------------------------------
John Fox, Professor
McMaster University
Hamilton, Ontario, Canada
http://socserv.mcmaster.ca/jfox/
    
On Sat, 20 Dec 2014 20:36:57 +0000 (UTC)

> Hi to everyone,
> I am trying to calculate the Bca bootstrap confidence intervals for the Pearson coefficient.
> x=Dataset$math.testy=Dataset$geo.testcor(x,y,method="pearson")[1] 0.6983799
> boot.ci(cor, conf=0.95, type=bca)Erreur dans boot.out$t0 : objet de type 'closure' non indiçable
> 
> I have tried as well to calculate the Pearson coefficient using bootstrap and then to calculate the Bca bootstrap CIs of the Pearson. It doesn't work either. 
> boot(data = cbind(x, y), statistic = cor, R = 200)
> 
> ORDINARY NONPARAMETRIC BOOTSTRAP
> 
> 
> Call:
> boot(data = cbind(x, y), statistic = cor, R = 200)
> 
> 
> Bootstrap Statistics :
>      original    bias    std. error
> t1* -0.6243713 0.6295142  0.2506267
> t2* -0.1366533 0.1565392  0.2579134
> > boot.ci(cor, conf=0.95, type=bca)
> Erreur dans boot.out$t0 : objet de type 'closure' non indiçable
> Many thanks to tell me how to correct my R script to get the Bca CIs for my Pearson coefficient. Best regards, looking forward to reading you,
> SV
> 
>     [[alternative HTML version deleted]]
> 
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> 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.

    
    


  
	[[alternative HTML version deleted]]



More information about the R-help mailing list