やる方法:「ブロック」(または - 「反復測定」?!)との相関を?
-
22-09-2019 - |
質問
私は分析するには、以下のセットアップを持っています: 私たちは約150科目を持ち、かつ各被験者のために我々は(異なる条件の下で)テストのペアを18回行いました。 試験の18点の異なる条件は、我々がここで平均に(各被験者について)試験超える場合、我々は(被験者間)テストの間に相関関係を取得しないであろうようになるように、相補的です。 私たちが知りたいことは被験者の中には、テストの間の相関(およびP値)であるが、すべての上の被験者ます。
私は今ではこれをしなかった方法では、各被験者の相関を実行し、それがその後、0異なる平均されていた場合、次に確認するためにその受信相関の分布を見ていました。 しかし、私は同じ質問を(誰かが私に「地理的相関」について何かを言ったが、浅い検索は助けにはならなかった)答えるためのより良い方法があるかもしれない疑います。
PS:私はそこに混合モデルのいくつかの並べ替えを行うには、ここの場所かもしれないが、私は「相関」を提示することを好む、および混合モデルから、このような出力を抽出する方法を確認していないでしょう理解
また、ここで私が何を言っているかのアイデアを与えるために、短いダミーコードがあります:
attach(longley)
N <- length(Unemployed)
block <- c(
rep( "a", N),
rep( "b", N),
rep( "c", N)
)
Unemployed.3 <- c(Unemployed + rnorm(1),
Unemployed + rnorm(1),
Unemployed + rnorm(1))
GNP.deflator.3 <- c(GNP.deflator + rnorm(1),
GNP.deflator + rnorm(1),
GNP.deflator + rnorm(1))
cor(Unemployed, GNP.deflator)
cor(Unemployed.3, GNP.deflator.3)
cor(Unemployed.3[block == "a"], GNP.deflator.3[block == "a"])
cor(Unemployed.3[block == "b"], GNP.deflator.3[block == "b"])
cor(Unemployed.3[block == "c"], GNP.deflator.3[block == "c"])
(I would like to somehow combine the last three correlations...)
すべてのアイデアは歓迎されます。
ベスト、 タル
解決
私はトリスタンに同意する - あなたはICCを探しています。標準実装の唯一の違いは、二つの評価者(テスト)を繰り返し、各被験者を評価することです。それを可能にする実装があるかもしれません。一方、ここでは相関関係を取得するための別のアプローチである。
あなたは明示的に残差との相関関係を許可する線形モデルの一般化されている「一般線形モデル」を、使用することができます。 gls
パッケージのnlme
機能を使用して実装し、この下のコード。私はだけでなく、他の方法があることを確認しています。 「長い」形式に最初のリシェイプにデータを我々が持っているこの機能を使用するには。私はまた、簡単にするためx
とy
に変数名を変更しました。それは私はあなたが意味考えるものだから、私はまた、代わりにあなたのコード内+rnorm(N)
の+rnorm(1)
を使用します。
library(reshape)
library(nlme)
dd <- data.frame(x=Unemployed.3, y=GNP.deflator.3, block=factor(block))
dd$occasion <- factor(rep(1:N, 3)) # variable denoting measurement occasions
dd2 <- melt(dd, id=c("block","occasion")) # reshape
# fit model with the values within a measurement occasion correlated
# and different variances allowed for the two variables
mod <- gls(value ~ variable + block, data=dd2,
cor=corSymm(form=~1|block/occasion),
weights=varIdent(form=~1|variable))
# extract correlation
mod$modelStruct$corStruct
モデリングフレームワークでは、p値を取得するために尤度比検定を使用することができます。 nlme
はまた、あなたの信頼区間を与えることができます:
mod2 <- gls(value ~ variable + block, data=dd2,
weights=varIdent(form=~1|variable))
anova(mod, mod2) # likelihood-ratio test for corr=0
intervals(mod)$corStruct # confidence interval for the correlation
他のヒント
私は専門家だが、あなたが望むもののように私には、このルックス。これは、短いコードに、上記のあなたの例と同様の相関を与え、p値を生成し、自動化されています。
> df = data.frame(block=block, Unemployed=Unemployed.3,
+ GNP.deflator=GNP.deflator.3)
> require(plyr)
Loading required package: plyr
> ddply(df, "block", function(x){
+ as.data.frame(
+ with(x,cor.test(Unemployed, GNP.deflator))[c("p.value","estimate")]
+ )})
block p.value estimate
1 a 0.01030636 0.6206334
2 b 0.01030636 0.6206334
3 c 0.01030636 0.6206334
すべての詳細を参照するには、次の操作を行います。
> dlply(df, "block", function(x){with(x,cor.test(Unemployed, GNP.deflator))})
$a
Pearson's product-moment correlation
data: Unemployed and GNP.deflator
t = 2.9616, df = 14, p-value = 0.01031
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
0.1804410 0.8536976
sample estimates:
cor
0.6206334
$b
Pearson's product-moment correlation
data: Unemployed and GNP.deflator
t = 2.9616, df = 14, p-value = 0.01031
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
0.1804410 0.8536976
sample estimates:
cor
0.6206334
$c
Pearson's product-moment correlation
data: Unemployed and GNP.deflator
t = 2.9616, df = 14, p-value = 0.01031
alternative hypothesis: true correlation is not equal to 0
95 percent confidence interval:
0.1804410 0.8536976
sample estimates:
cor
0.6206334
attr(,"split_type")
[1] "data.frame"
attr(,"split_labels")
block
1 a
2 b
3 c