「統計学関連なんでもあり」の過去ログ--- 040

No.02652　分散共分散行列の同等性の検定について　　【nakamura】　2007/02/08(Thu) 17:38

分散共分散行列の同等性の検定の際に用いられているbという係数はどのような過程で算出されるのでしょうか？
他の文献なども参考にしたのですが係数がないケースなのでわかりませんでした。
ご教授よろしくお願いいたします。

No.02653　Re: 分散共分散行列の同等性の検定について　　【青木繁伸】　2007/02/08(Thu) 17:56

http://aoki2.si.gunma-u.ac.jp/lecture/Wilks/wilks3.html
のことですか？

> bという係数はどのような過程で算出されるのでしょうか

とは，どのような答えを求めているのか私にはよくわかりませんが，b の計算式がわからないと言うことではないのでしょうね（5. のところに式が書いてありますから）。

No.02654　Re: 分散共分散行列の同等性の検定について　　【nakamura】　2007/02/08(Thu) 18:04

ご返事ありがとうございます。質問がわかりにくく申し訳ありませんでした。

記入されているURLの内容です。
他の文献を参考にしていたらbという係数が書かれていませんでした。
bという係数がどうして必要なのか？またどのように導出されたのか疑問に思い質問させていただきました

No.02656　Re: 分散共分散行列の同等性の検定について　　【青木繁伸】　2007/02/08(Thu) 19:44

どのような計算式でしたか？
その文献の情報を付けていただけると良いかも。（文献をすぐに手に入れることができないことも多いので，その計算式を是非）

No.02661　Re: 分散共分散行列の同等性の検定について　　【nakamura】　2007/02/08(Thu) 23:21

参考にした文献は数理統計学竹内啓著という本です。

2つの分布Σ1=Σ2の検定
群1の要素数がm分散共分散行列がΣ1,群2の要素数がn分散共分散行列がΣ2のとき
W={|Σ1|**(m/2)|Σ2|**(n/2)}/|Σ|**((m+n)/2)
Σ={(m-1)Σ1+(n-1)Σ2}/(m+n-2)
このとき-2logWがχ**2分布に従う

と書かれておりました。
いかがでしょうか？

No.02663　Re: 分散共分散行列の同等性の検定について　　【青木繁伸】　2007/02/09(Fri) 11:10

B.ボルチ，C.ファング原著　中村慶一訳
応用多変量解析
森北出版株式会社
の，86ページに基づいています。
b は修正項で，
Box, G.E.P.(1949): A General Distribution Theory for a Class of Likelihood Criteria, Biometrika, Vol. 36, pp. 317-346.
で提案されたもの。
ということです。

No.02666　Re: 分散共分散行列の同等性の検定について　　【nakamura】　2007/02/09(Fri) 11:52

わかりました
ぜひ文献探してみようと思います
ありがとうございました

No.02669　Re: 分散共分散行列の同等性の検定について　　【青木繁伸】　2007/02/09(Fri) 14:16

どれくらい補正されるのか，R でシミュレーションしてみました。

> eq.cov2 <- function(x, y)
+ {
+         p <- ncol(x)
+         s1 <- var(x)
+         s2 <- var(y)
+         n1 <- nrow(x)-1
+         n2 <- nrow(y)-1
+         sa <- (n1*s1+n2*s2)/(n1+n2)
+         b <-  1-(1/n1+1/n2-1/(n1+n2))*(2*p^2+3*p-1)/(6*p+6)
+         chi <-(n1+n2)*log(det(sa))-n1*log(det(s1))-n2*log(det(s2))
+         df <- p*(p+1)/2
+         p0 <- pchisq(chi, df, lower=FALSE)
+         chi1 <-b*chi
+         p1 <- pchisq(chi1, df, lower=FALSE)
+        return(c(p0, p1))
+ }
> 
> sim <- function(p=3, n=100, loop=100000)
+ {
+ 	result <- sapply(1:loop, function(i) {
+ 		x <- matrix(rnorm(n*p), n)
+ 		y <- matrix(rnorm(n*p), n)
+ 		eq.cov2(x, y)
+ 		} )
+ 	return(rowSums(result<0.05)/loop)
+ }

平均値0，分散1の正規母集団からサンプルサイズ n の p 変量データを2組取り出します。この場合，理論的にはどちらの分散・共分散行列も対角行列になるはずです。この同じ条件での検定を10万回行い，修正項ありの場合となしの場合の有意確率を記録します。帰無仮説が棄却された回数の割合（第一種の過誤を求めます）

> sim(2, 50) # 2変数，サンプルサイズ50
[1] 0.05401 0.05020 # 一番目が修正項なし，二番目が修正項ありのときの第一種の過誤
> sim(2, 100)
[1] 0.05216 0.05031
> sim(2, 500)
[1] 0.05169 0.05136
> sim(2, 1000)
[1] 0.05035 0.05013
> sim(3, 50)
[1] 0.05780 0.04923
> sim(3, 100)
[1] 0.05294 0.04895
> sim(3, 500)
[1] 0.05137 0.05049
> sim(3, 1000)
[1] 0.05022 0.04977
> sim(5, 50)
[1] 0.07114 0.04972
> sim(5, 100)
[1] 0.06001 0.05031
> sim(5, 1000)
[1] 0.05069 0.04973結論としては，「大差はないが，修正した方がましである」ということのようです。