Sponsored Links
データの読み込みと準備
サンプルデータは、「2-3d 因子分析実例《セールスマンデータ）（2-3d.xls）」[1]を使わせていただきます。もともと因子分析の練習をするためのデータで、セールスマンに対する評価項目がたくさん掲載されていて、多分それぞれ10段階順序尺で評価されています。ここでは容姿、学力、好感度の3因子を利用します。
> x <- read.table(pipe("pbpaste"), header=TRUE) # クリップボードからデータを読み込み、xに代入。ちなみにこちらはMac用のコードでWindowsは、read.table("clipboard", header=TRUE)。
> str(x)
'data.frame': 48 obs. of 3 variables:
$ youshi : int 6 9 7 5 6 7 9 9 9 4 ...
$ gakuryoku: int 7 10 8 6 8 7 9 9 9 7 ...
$ koukando : int 2 5 3 8 8 7 8 9 7 10 ...
> xx <- stack(x) # 多重比較検定のためにスタック形式のデータも作成しておく。
> colnames(xx) <- c("point", "items") # 列名をpoint（得点）とitems（評価項目）に変更。 四分位数を計算します。
> quantile(x$youshi)
0% 25% 50% 75% 100%
0 4 6 9 10
> quantile(x$gakuryoku)
0% 25% 50% 75% 100%
3.00 6.00 7.00 8.25 10.00
> quantile(x$koukando)
0% 25% 50% 75% 100%
2 6 7 8 10
フリードマン検定の実行
フリードマン検定（friedman.test）はmatrix（行列）にデータを変換しないといけないようです。
> y <- matrix(c(x$youshi, x$gakuryoku, x$koukando), ncol = 3) # ncolで列数を指定。yに代入。
> friedman.test(y)
Friedman rank sum test
data: y
Friedman chi-squared = 6.2874, df = 2, p-value = 0.04312
多重比較検定の実行
フリードマン検定の後の多重比較検定には、ウィルコクソンの符号順位検定で得られるp値をBonferroniやHolmによる補正をして用いる方法があります。計算方法は以下の通りです。
> pairwise.wilcox.test(xx$point, xx$items, p.adj="bonferroni",paired = TRUE, exact=F) # Bonferroniの方法
> pairwise.wilcox.test(xx$point, xx$items, p.adj="holm",paired = TRUE, exact=F) # Holmの方法
Post hoc検定については、私のための統計処理というページのこちらに分かりやすくまとめてくださっています。これを見ると、Post hocとしてノンパラメトリックの2群の比較検定を繰り返す方法以外にも、Steel-Dwass法やGames-Howelという方法があるのですね[1]。
ここではSteel-Dwass法を紹介します。NSM3というパッケージを使ってみます。
> library(NSM3)
> pSDCFlig(xx$point, xx$items) # pSDCFligはスタック形式のデータを指定。Steel-Dwass検定を実行。
Ties are present, so p-values are based on conditional null distribution.
Group sizes: 48 48 48
Using the Monte Carlo (with 10000 Iterations) method:
For treatments youshi - gakuryoku, the Dwass, Steel, Critchlow-Fligner W Statistic is 2.7743.
The smallest experimentwise error rate leading to rejection is 0.1234 .
For treatments youshi - koukando, the Dwass, Steel, Critchlow-Fligner W Statistic is 2.8785.
The smallest experimentwise error rate leading to rejection is 0.1058 .
For treatments gakuryoku - koukando, the Dwass, Steel, Critchlow-Fligner W Statistic is -0.2576.
The smallest experimentwise error rate leading to rejection is 0.9826 .
フリードマン検定では有意差が認められましたが、Steel-Dwass法を用いたPost hoc検定では有意差は認められませんでした。こういうことがあるから、いきなり多重比較検定を行っても良いという意見はもっともです。
Steel-Dwass法は、統計学的手法の話題 – 生物科学研究所のこちらのページを参考にさせていただきました[2]。pSDCFligではmethod=”Exact”で正確検定を指定できるようですが、今回使ったサンプルサイズが大きかったからか上手く動かなかったので、何も設定しませんでしたが、何も設定しないとmethod=”Monte Carlo”として計算されるようです。それでも計算にかなりの負荷がかかるらしく、数分固まっていました。
効果量の計算
クラスカル=ウォリス検定やフリードマン検定後の多重比較検定における効果量は、組み合わせごとに2群間のノンパラメトリック検定の効果量を用います。
今回は対応ある2群の比較なので、ウィルコクソンの符号順位検定を用います（対応のない場合は、ウィルコクソンの順位和検定）。Rにはwilcox.testという関数が用意されていますが、この関数はデータにタイ（同じ値）があった時に正確な値を計算できないらしく、パッケージexactRankTestsをインストールして、正確検定を行うことができるwilcox.exactを用います。ただし、データ数が多いと正確検定はできないようです[1]（データ科学便覧のこちらの記事を参考にしました）。
> library(exactRankTests) # パッケージexactRankTestsの読み込み
> wilcox.exact(x$youshi, x$gakuryoku, paired = TRUE) # paired=TRUEで符号順位検定、paired=FALSEで順位和検定を指示。
> wilcox.exact(x$youshi, x$koukando, paired = TRUE)
> wilcox.exact(x$koukando, x$youshi, paired = TRUE)
ウィルコクソンの符号順位検定そのものは今回の計算では関係なので、結果は割愛します。効果量の計算にはパッケージorddomを使います。たくさん結果が出ますが、その中のdeltaを参照します。
> library(orddom) # パッケージorddomをインストール。
> orddom(x$youshi, x$gakuryoku, paired = TRUE) # 対応のあるデータなのでpaired=TRUEと指定。
> orddom(x$youshi, x$koukando, paired = TRUE)
> orddom(x$gakuryoku, x$koukando, paired = TRUE)
容姿と学力：0.354166666666667
容姿と好感度：0.166666666666667
学力と好感度：-0.0625
となります。
[1] データ科学便覧 – Rによるウィルコクソンの符号順位検定
Sponsored Links
