【統計の勉強】RでSteel-Dwass検定〜NSM3でpSDCFligがうまく使えなかった話

こんにちは。絶賛R言語勉強中の生命科学系研究者のえいこです。

前回は、Dunnett検定のノンパラメトリック版であるSteel検定用の関数を9割写経で自作しました。

今回は、Steel-Dwass検定のコードを書いていきます。

今回の結論は、「Steel-Dwass検定のコードは自作してみよう！」ということになるのですが、なぜ、世に溢れるパッケージを使わないで、自作の関数を作ることにしたのか...

今回は、その経緯を書き留めておこうと思います。

正確検定と漸近検定の間

Steel検定の時と同じように、最初はすでにあるパッケージなどを使ってサクッと計算したい！と思っていました。

Steel-Dwass検定について色々調べてみると（日本語検索で引っかかってきたもので）、

・群馬大学の青木先生のコード
・生物科学研究所の井口先生のNSM3を使った方法

という大きく分けて2つの方法があるようです。

google検索で表示される上位10位の中で、井口先生が正確検定の重要性を訴えていらっしゃるものが4〜5個ありました。

そして、青木先生のコードやEZRの検定は漸近検定であって、それを知らずに使っている人が多いとのこと。

正確検定と漸近検定は何が違うのか？

Steel-Dwass検定の場合、Mann-WhitneyのU検定をベースにしているのですが、検定統計量の計算を正規近似して算出する方法を漸近検定と呼んでいます。

これは、サンプル数が多くなれば必然的に検定統計量も正規分布に近づいてくるので、漸近検定でも問題ないのですが、サンプル数が少ない時には近似すると本来の値とのズレが大きくなってしまうようです。

要は、各サンプル数が10よりも小さい場合は正確検定をしっかりしたほうが良いということ。

そこで、井口先生がお勧めしていらっしゃったNSM3パッケージをインストールして、pSDCFlig関数で検定方法を”exact”にして検定してみることにしました。

すると...

実行するとRStudioがスタックして動かなくなってしまいました。

2回目、3回目、挑戦してみたものの、エラーは出ないもののスタックして途中で落ちてしまうことの繰り返し。

これはあまり使えないのでは？

他の人のページを見てみると、”excat”でうまくいかなかったので”Monte Carlo”にしてみたらうまくいった。という記述がありました。

Monte Carlo法は乱数を発生させてシミュレーションを行う方法です。

NSM3ではシミュレーション試行回数を指定できます。デフォルトでは1000回となっているので、1000回シミュレーションでやってみました。

1、2分止まってから、結果が出力されました。

うーん、微妙。

NSM3のコードの記録は、とっていません。
Rがスタックしたりしてべつに描いていたコードに影響があったら困るので、インストールしたパッケージは全部アンインストールしてしまいました。
とにかく、あまりうまく動いてくれないという印象です。

結論：Steel-Dwass検定をする前提を変えた方が良い

Steel-Dwass検定は漸近検定だけに限定する

少しだけ視点を変えて。

群馬大学の青木先生も、EZRのSteel-Dwass検定もなぜ、漸近検定だけを採用しているのか？

統計に精通している人たちが作っているので、正確検定があることも頭に入っているはずです。でも、漸近検定だけにしているのは、理由があるはず。

Steel-Dwass検定はサンプル数が多い（1グループあたり10サンプル以上）時にするもので、サンプルサイズが小さい時には別の検定をするという風に割り切っているのでは？

ということで調べてみると...

私のための統計処理　ー多重比較検定

Kruskal-Wallis検定の後の多重比較の手法

などを見てみると、各グループで二群比較をしてからp値を補正する、Dunn-Bonferroniを使うこともあるようです。

サンプルサイズが小さい時は、Dunn-Bonferroniを使って計算してみて、サンプルサイズが大きい時にはSteel-Dwass検定をすることしてみます。

統計学では色々な意見の人がいるようで、「この場合はこれ！」みたいな正解がないのが現状のようなので...

何か間違っているよという場合はご意見いただけますと嬉しいです。

大まかな方針は決まったので、今どんな感じでコードを書いているのかまとめてみます。

次回は、Steel-Dwass検定のコードとDunn-Bonferroniのコードを順番に確認していきたいと思います。

それでは、また！