理学療法士がまとめたノート

統計学備忘録 since2016

Rを使って統計学を勉強するブログです

等分散性の検定

準備としてカイ二乗分布とF分布が必要になります

カイ二乗分布
確率変数Z1,Z2,…Znが互いに独立で標準正規分布にしたがうとき
W=Z1^2+Z2^2+…+Zn^2
の従う分布を自由度nのカイ二乗分布とよび、χ^2(n)で表します.
期待値E[W]=n, 分散V[E]=2n

カイ二乗分布の性質
① 二つの確率変数W1,W2が独立にχ^2(n1), χ^2(n2)に従うとき
W1+ W2は、χ^2(n1+n2)に従います.(カイ二乗分布の再生性)
② 母集団N(μ,σ^2)からの無作為標本X(X1,X2,…X3)について
Xiを標準すると Zi= (Xi-μ)/ σとなる
Ziは互いに独立で標準正規分布に従うので
∑{(Xi-μ)/ σ}^2 はカイ二乗分布χ^2(n)に従うことになります
③ 標本不偏分散
s^2={∑{(Xi-X平均) ^2 }/(n-1) より
(n-1) ×s^2=∑(Xi-X平均) ^2
母集団N(μ,σ^2)からの無作為標本Xi=(X1,X2,…X3)は互いに独立でN(μ,σ^2)に従うので
W=∑{(Xi-X平均) ^2}/ σ^2 = {(n-1) ×s^2}/σ^2
母平均μをX平均で置き換えると自由度がn-1になります
したがって{(n-1) ×s^2}/σ^2は自由度n-1のカイ二乗分布χ^2(n-1)に従います.

f:id:yoshida931:20170222184024p:plain

http://tips-r.blogspot.jp/2015/01/r.htmlより)


2標本の母分散の比の検定
散らばりを表す分散は非負の値なので、その大きさを考える場合には差ではなく比で考える方が適切です.

例)単位はmm
A群の身長の標準偏差はσa=5 (分散25)
B群の身長の標準偏差はσb=10 (分散100)
σa/σb=0.5 A群の方がB群に比べると2倍散らばりの小さい集団である
現実的には標準偏差が便利で理解しやすいのですが、推定での取り扱いは上述のような分散を使います.

 

F分布
条件
・確率変数U:自由度mのカイ二乗分布に従う
・確率変数V:自由度nのカイ二乗分布に従う
・UとVは独立である

フィッシャーの分散比 F= (U/m)/(V/n)

Fが従う確率分布を自由度(m,n)のF分布F(m,n)という

2標本の標本分散s1^2, s2^2を考えます
標本不偏分散s^2={∑{(Xi-X平均) ^2 }/(n-1) 
標本s1^2:{(m-1) ×s1^2}/σ1^2は自由度m-1のカイ二乗分布χ^2(m-1)に従います
標本s2^2:{(n-1) ×s2^2}/σ2^2は自由度n-1のカイ二乗分布χ^2(n-1)に従います
s1^2とs2^2は独立です

F={(m-1) ×s1^2}/(σ1^2)*(m-1)÷{(n-1) ×s2^2}/(σ2^2)*(n-1)
= (s1^2/σ1^2)÷(s2^2/σ2^2)
= (s1^2/s2^2)*(σ2^2/σ1^2)
の従う分布は、自由度(m-1,n-1)のF分布F(m-1,n-1)となります

 

 

それでは始めます
赤文字のみRで実行

2標本の場合

母分散の比の検定になります
二つの母集団からの標本をXi~N(μx,σx^2),Yi~N(μy,σy^2)とします.
帰無仮説:σx^2=σy^2、(σx^2)/(σy^2)=1
対立仮設:σx^2≠σy^2、(σx^2)/(σy^2)≠1

サンプルは「石村卓夫;分散分析のはなし,東京図書,1998,p69」
Xi<-c(12.2,18.8,18.2)
Yi<-c(26.4,32.6,31.3)
σx二乗<-sum((Xi-mean(Xi) )^2)/2;σx二乗
σy二乗<-sum((Yi-mean(Yi) )^2)/2;σy二乗

T<-σx二乗/σy二乗
1.25
F分布を使用します

qf(0.025,2,2,lower=F) # F 分布の上側5%点
39
qf(0.975,2,2,lower=F) # F 分布の下側5%点
0.026

T=1.25は棄却域には入らないので帰無仮説は棄却できない
標本をXiとYiの母分散は等しいと仮定できることになります

 

SPSSではlenene検定が使用されていますので以下に算出しておきます
Levene's Test もF分布を使用していますが、上記方法とは異なります
leveneTest(x,group)

Levene's Test for Homogeneity of Variance (center = median)
    Df  F value  Pr(>F)
group  1  0.0031  0.9585
    4

 

 

 

 

3標本以上の場合

 

繰り返し数が等しいとき:Hartleyの検定
各水準の分散(Si^2 )について
Si^2の最大分散÷最小分散が1に近づくことで等分散性を仮定します
最大分散/最小分散=max{Si^2}/min{Si^2}

 

f:id:yoshida931:20170222184102p:plain

http://www2.vmas.kitasato-u.ac.jp/lecture0/statistics/hartley.pdf


検定表の見方:Fmax(k , f)(k) = Fmax(水準数,繰り返し数-1)(水準数)
帰無仮説:各水準の分散は等しい

サンプルは「石村卓夫;分散分析のはなし,東京図書,1998,p108」
A<-c(12.2,18.8,18.2)
B<-c(22.2,20.5,14.6)
C<-c(20.8,19.5,26.3)
D<-c(26.4,32.6,31.3)
E<-c(24.5,21.2,22.4)

x<-c(A,B,C,D,E)

SuppDistsをインストトールします
library(SuppDists)
factor()カテゴリーを要素としたベクトル
group=factor(rep(c("A", "B", "C", "D", "E"), c(3, 3, 3, 3, 3)))
levels()でグループ化をstr()で構成を確認できます
fmax=max( by(x, group, FUN=var) ) / min( by(x, group, FUN=var) )  
5.70

by関数の引数(data,INDICES,FUN=,simplify=)
data:data.frameやmatrixもOK
INDICES:factorをsubsetに指定
FUN :引数dataに適用する関数を指定。ここでは分散varを指定。
simplify :結果のデータ形式。TRUE:scalar、FALSE:list or array。

 

積分布pmaxFratio
pmaxFratio( fmax, 繰り返し数の自由度, 水準数)
1-pmaxFratio(fmax, 2, 5)  #p値
0.8154272

帰無仮説は棄却されず、等分散性を仮定することができる

 

 

繰り返し数が異なるとき:Bartlettの検定
Rヘルプより
繰り返し数は異なりませんが同じサンプルで試行してみます
以下実施方法のみ記載
A<-c(12.2,18.8,18.2)
B<-c(22.2,20.5,14.6)
C<-c(20.8,19.5,26.3)
D<-c(26.4,32.6,31.3)
E<-c(24.5,21.2,22.4)

x<-c(A,B,C,D,E)
group=factor(rep(c("A", "B", "C", "D", "E"), c(3, 3, 3, 3, 3)))

bartlett.test(x,group)


Bartlett test of homogeneity of variances
data: x and group
Bartlett's K-squared = 1.2291, df = 4, p-value = 0.8733

 


Leveneの検定
最もよくみかける等分散性の検定だと思われます
水準の繰り返し数に縛りはありません

A<-c(12.2,18.8,18.2)
B<-c(22.2,20.5,14.6)
C<-c(20.8,19.5,26.3)
D<-c(26.4,32.6,31.3)
E<-c(24.5,21.2,22.4)

統計量(Wikipediaより)
W={(N-k)/(k-1)}×{∑Ni(Zi.-Z..)^2/∑∑(Zij-Zi.)^2}

数値はsample
k:比較する群の数水準
水準数i  1,2,3,4,5
繰り返しj  1,2,3
N:総サンプル数
Ni:第i群のサンプル数
自由度N-k=15-5=10
自由度k-1=5-1=4
Z.. 全データの平均値
Zi. i番目の水準内の平均値 
Yij:第i群のj番目の変数
Zij:| Yij – Yiの平均|、| Yij – Yiの中央値|


carパッケージをインストールします
x<-c(A,B,C,D,E)
group=factor(rep(c("A", "B", "C", "D", "E"), c(3, 3, 3, 3, 3)))

leveneTest(x,group)

Levene's Test for Homogeneity of Variance (center = median)
Df F value Pr(>F)
group 4 0.1256 0.9698
10

p値=0.9698で帰無仮説は棄却されず、等分散性を仮定した検定が可能になります

 

引用・参考
バイオインフォマティクス入門 http://bio-info.biz/
石村卓夫;分散分析のはなし,東京図書,1998