読者です 読者をやめる 読者になる 読者になる

草薙の研究ログ

英語教育関係。でも最近は統計(特にR)ネタが中心。

論文で報告されていない相関係数をサルベージ!(Rの関数つき)

情報処理 r 統計

背景

対応のあるデータの論文で,あるある。

平均値と標準偏差のペアが報告されている,よしよし。

対応のあるt検定の検定統計量が報告されている,ふむふむ。

でも相関係数や分散共分散行列が報告されていない!なんでだよー。。。

 

たとえば,標準化効果量は普通のd(標準化平均差)と,差得点に基づく標準化平均差では値が異なる。もちろん,相関係数がわかれば,このふたつは変換できるけど,相関係数が報告されていないとそれはむり。むり。メタ分析にも影響するよね。

 

もともと分散共分散行列,平均があればデータはだいたいシミュレートできる。

R,MASSパッケージのmvrnorm関数とか。でも,報告されていないとなー。

 

指導法の検証でも,そもそも共分散は教育的にも大きな意味を持つので,報告したほうがいいとおもうの。切実に。切実に。

 

でも,まあいいや。以下の条件なら,報告されていない相関係数がわかるから。

 

  1. 平均値のふたつが報告されている
  2. 標準偏差のふたつが報告されている
  3. 対応のあるt検定の検定統計量tがある
  4. (3の代わりにp値と標本サイズが報告されている)

 

関数つくっちゃえ!

対応ありのt検定にもちいる,「差得点に基づく標準偏差」は,完全データがなくても,二群の標準偏差と相関係数からももとまる。簡単にいえば,対応のあるt検定には相関係数が使われているので,そこからサルベージできるかもだ。

数式から計算もできるのだろうけど,ちょっとめんどくさいので,以下のようなアルゴリズムで答えをえちゃえ。

 

 

  1. 相関係数を0.000から,0.999まで用意する
  2. 論文で報告されている平均値と標準偏差と,用意したすべての相関係数をもちいて対応のあるt検定の検定統計量を計算する
  3. 論文で報告されているt値との差分をすべてもとめる
  4. 最小の差分に対応する相関係数が答え!

簡単だ♪尤度といえなくもない。

サルベージ!サルベージ!

これをRの関数でつくっちゃえ。

 

 

salvager<-function(n,m1,m2,sd1,sd2,t){


md<-abs(m1-m2) #平均差
v1<-sd1^2 #分散1
v2<-sd2^2 #分散2
sr<-seq(0,.999,.001) #相関係数の用意
sigmad<-sqrt(v1+v2-(2*sr*sd1*sd2)) #差得点の標準偏差
simt<-md/sigmad*sqrt(n) #対応あるt値を計算
diff<-t-simt #報告値との差分
ref<-data.frame(sr,abs(diff))
mindiff<-min(ref[,2]) #最小の差分をもとめる
rr<-ref[,1][ref[,2]==mindiff] #対応する相関係数


plot(ref,type="l",xlab="r", ylab="Difference") #図示
points(rr,mindiff,cex=2,pch=4,col=2,lwd=3)
varmatrix<-matrix(c(v1,rr*sd1*sd2,rr*sd1*sd2,v2),2,2) #分散共分散行列
list("r"=rr,"t"=t,"var"=varmatrix) #結果


}

 

使い方は,

salvager(標本サイズ, 平均値1, 平均値2, 標準偏差1, 標準偏差2, 報告されているt値)

 

使い方!

たとえば,

salvager(100,20,30,10,10,15.8114)


$r
[1] 0.8

$t
[1] 15.8114

$var
[,1] [,2]
[1,] 100 80
[2,] 80 100

 

こんな感じで,サルベージした相関係数と,報告されたt値と分散共分散行列を返してくれる。自動で図示もするよ。こんなかんじ。

 

f:id:kusanagik:20150406224940p:plain

たとえば,MASSのmvrnormで,

 

mvrnorm(100,mu=c(20,30),Sigma=matrix(c(100,30,30,100),2,2),empirical=T)->dat

 

とか作って,

 

t.test(dat[,1],dat[,2], paired=TRUE)


Paired t-test

data: dat[, 1] and dat[, 2]
t = -8.4515, df = 99, p-value = 2.563e-13
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
-12.347757 -7.652243
sample estimates:
mean of the differences
-10

t値が8.4515なので,いれてサルベージしてやる。

salvager(100,20,30,10,10,8.4515)
$`Salvaged r`
[1] 0.3

$t
[1] 8.4515

$var
[,1] [,2]
[1,] 100 30
[2,] 30 100

 

おうおう!分散共分散行列が再現出来てるね。

念のため,

cor(dat)

[,1] [,2]
[1,] 1.0 0.3
[2,] 0.3 1.0

 

よしよし。

メタ分析にも役立つかもだし,レビューでも役立つかも。