Начальная корреляция в R - PullRequest
Новая
       101

Начальная корреляция в R

0 голосов
/ 15 октября 2019

Я пытаюсь сделать загрузочную корреляцию в R. У меня есть две переменные Var1 и Var2, и я хочу получить загрузочную p.value корреляции Пирсона. 

my variables look like this:
      x            y
1   .6080522    1.707642
2   1.4307273   1.772616
3   0.8226198   1.768537
4   1.7714221   1.265276
5   1.5986213   1.855719
6   1.0000000   1.606106
7   1.1678940   1.671457
8   0.6630012   1.608428
9   1.0842423   1.670619
10  0.5592512   1.107783
11  1.6442616   1.492832
12  0.8326965   1.643923
13  1.1696954   1.763181
14  0.7484543   1.762921
15  1.0842423   1.591566
16  0.9014748   1.718669
17  0.7604917   1.782863
18  0.8566499   1.796216
19  1.4307273   1.913675
20  1.7579695   1.903155

Пока у меня есть это: 

data = as.data.frame(data)
x = data$Var1
y = data$Var2
dat = data.frame(x,y)

library(boot)
set.seed(1)
bootCorTest3 <- function(data, i){
  d <- data[i, ]
  results  <- cor.test(d$x, d$y, method='pearson')
  c(est = results$estimate, stat = results$statistic, param = results$parameter, p.value = results$p.value, CI = results$conf.int)
}


b3 <- boot(dat, bootCorTest3, R = 1000)
b3

# Original (non-bootstrap) statistics with label
b3$t0
colMeans(b3$t)
boot.ci(b3, type = c("norm", "basic", "perc", "bca")) #bootstrapped CI.

Начальное значение p должно быть тем, которое я получаю с colMeans (b3 $ t), верно?

colMeans (b3 $ t) дает мне это: 

est.cor      stat.t    param.df     p.value         CI1         CI2
 0.28495324  2.13981008 48.00000000  0.14418623  0.01438146  0.51726022

Кажется, все работает нормально. Проблема в том, что я запустил одну и ту же статистику для другого программного обеспечения, и результаты сильно различаются. Значение p, которое я здесь получаю, намного выше, чем на другом. Я думаю, что я, возможно, сделал что-то не так, потому что я не силен в R.

Может кто-нибудь дать мне отзыв об этом коде? Я делаю что-то неправильно? Как бы вы получили начальное значение p.value для корреляции Пирсона?

Спасибо за потраченное время.

1 Ответ

0 голосов
/ 08 ноября 2019

Если вы хотите запустить тест корреляции, вам нужно только вернуть коэффициент корреляции из функции статистики начальной загрузки. Самозагрузка p-значения теста корреляции в этом случае не подходит, так как вы игнорируете направленность корреляционного теста.

Проверьте этот вопрос на CrossValidated для некоторых хороших ответов на выполнение тестов начальной загрузки: https://stats.stackexchange.com/questions/20701/computing-p-value-using-bootstrap-with-r

library("boot")
data <- read.csv("~/Documents/stack/tmp.csv", header = FALSE)
colnames(data) <- c("x", "y")

data <- as.data.frame(data)
x <- data$Var1
y <- data$Var2
dat <- data.frame(x,y)

set.seed(1)

b3 <- boot(data, 
  statistic = function(data, i) {
    cor(data[i, "x"], data[i, "y"], method='pearson')
  },
  R = 1000
)
b3
#> 
#> ORDINARY NONPARAMETRIC BOOTSTRAP
#> 
#> 
#> Call:
#> boot(data = data, statistic = function(data, i) {
#>     cor(data[i, "x"], data[i, "y"], method = "pearson")
#> }, R = 1000)
#> 
#> 
#> Bootstrap Statistics :
#>      original        bias    std. error
#> t1* 0.1279691 -0.0004316781    0.314056
boot.ci(b3, type = c("norm", "basic", "perc", "bca")) #bootstrapped CI. 
#> BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS
#> Based on 1000 bootstrap replicates
#> 
#> CALL : 
#> boot.ci(boot.out = b3, type = c("norm", "basic", "perc", "bca"))
#> 
#> Intervals : 
#> Level      Normal              Basic         
#> 95%   (-0.4871,  0.7439 )   (-0.4216,  0.7784 )  
#> 
#> Level     Percentile            BCa          
#> 95%   (-0.5225,  0.6775 )   (-0.5559,  0.6484 )  
#> Calculations and Intervals on Original Scale

plot(density(b3$t))
abline(v = 0, lty = "dashed", col = "grey60")

В этом случае без значения p вполне можно сказать, что большая часть массы распределения выборки очень близка к нулю.

...