时间戳

首页 > 解决方案 > R 上相关系数的自举 p 值

r - R 上相关系数的自举 p 值

问题描述

在 上R，我使用 boostrap 方法来获得相关系数估计和置信区间。为了得到 p 值，我想，我可以计算不包含零的置信区间的比例。但这不是解决方案。

在这种情况下如何获得 p 值？

我cor.test用来获取系数估计。cor.test还可以给我每次测试的 p 值。但是我怎样才能得到引导的 p 值呢？

非常感谢 ！

下面是一个例子：

n=30
data = matrix (data = c (rnorm (n), rnorm (n),rnorm (n), rpois(n,1), 
rbinom(n,1,0.6)), nrow =  n, byrow = F)
data= as.data.frame(data)
z1  = replicate( Brep, sample(1:dim(data)[1], dim(data)[1], replace = T))
res = do.call  ( rbind, apply(z1, 2, function(x){ res=cor.test(data$V1[x], data$V2[x]) ; return ((list(res$p.value,res$estimate))) }))

 coeffcorr  = mean(unlist(res[,2]), na.rm = T) #bootstrapped coefficient
 confInter1 = quantile(unlist(res[,2]), c(0.025, 0.975), na.rm = T)[1] #confidence interval 1
 confInter2 = quantile(unlist(res[,2]), c(0.025, 0.975), na.rm = T)[2] #confidence interval 2  
 p.value    =  mean    (unlist(res[,1]), na.rm = T )  # pvalue

标签： rcorrelationconfidence-intervalp-valuestatistics-bootstrap

解决方案

R 中引导的标准方法是使用 base package boot。首先定义 bootstrap 函数，该函数接受两个参数，即数据集和数据集的索引。这是bootCorTest下面的功能。在函数中，您对数据集进行子集化，仅选择索引定义的行。

其余的很简单。

library(boot)

bootCorTest <- function(data, i){
    d <- data[i, ]
    cor.test(d$x, d$y)$p.value
}


# First dataset in help("cor.test")
x <- c(44.4, 45.9, 41.9, 53.3, 44.7, 44.1, 50.7, 45.2, 60.1)
y <- c( 2.6,  3.1,  2.5,  5.0,  3.6,  4.0,  5.2,  2.8,  3.8)
dat <- data.frame(x, y)

b <- boot(dat, bootCorTest, R = 1000)

b$t0
#[1] 0.10817

mean(b$t)
#[1] 0.134634

boot.ci(b)

有关函数结果的更多信息boot，boot.ci请参阅它们各自的帮助页面。

编辑。

如果要从引导统计函数返回多个值bootCorTest，则应返回一个向量。在以下情况下，它返回一个具有所需值的命名向量。

请注意，我设置了 RNG 种子，以使结果可重现。我应该已经在上面做了。

set.seed(7612)    # Make the results reproducible

bootCorTest2 <- function(data, i){
    d <- data[i, ]
    res <- cor.test(d$x, d$y)
    c(stat = res$statistic, p.value = res$p.value)
}


b2 <- boot(dat, bootCorTest, R = 1000)

b2$t0
#  stat.t  p.value 
#1.841083 0.108173


colMeans(b2$t)
#[1] 2.869479 0.133857

推荐阅读