第4章は統計を扱います。


今回も「シロート統計学」のハルさんとコラボレーションすることとなりました。


シロート統計学はEZRを使った統計分析をわかりやすく解説されています。


 


第4章はシロート統計学で使われていたEZRをRで行うとどうなるのか?といった視点で進めていきます。


今回使うデータもハルさんのサイトと同じものを使わせでいただく事になりました。それぞれ見比べることで参考にしてみてください!


今回はWilcoxon符号付順位和検定を紹介します



まずWilcoxon符号付順位和検定についてはハルさんのサイトをご参照ください。



また1.準備〜4.データの読み込みまでは【4-1】Rでt検定を行う方法と全く同じ流れになります。
もし1〜4まででわからない部分があれば確認してみてください。

 

1.準備

第4章は毎回ExcelデータをダウンロードしてRを使うのでプロジェクトで管理して行うことを勧めています。

 

ここではR練習というプロジェクトを作り、Excelファイルを入れるためのdataフォルダを作っています。
これを前提に次から進めていきます。
スクリーンショット 2019-10-20 7.54.14


2.スクリプトファイルの作成

次にRのコードを書くためのスクリプトファイルを作ります。

スクリーンショット 2019-11-04 0.47.12



3.データのダウンロード

今回もハルさんのサイトのデータを使わせていただきます。

デモデータ(wilcoxon符号付順位和検定)


この章ではRを使ってダウンロードしています。


download.file(url = “ファイルのURL”,
        destfile = “保存したい場所/ファイル名”)
urlはデモデータで右クリック → リンクのアドレスをコピー
destfileは保存場所と保存のファイル名を指定します。


実際のコードは以下になります。
前回のコードのURL(" "の中)とdestfileのdata/以降を変更するだけでOKです。
url <- "https://haru-reha.com/wp-content/uploads/2018/04/demo-wilcoxon-rank-sum-test.xlsx"
destfile = "data/data/demo-wilcoxon-rank-sum-test.xlsx"

download.file(url, destfile)
スクリーンショット 2019-11-04 7.41.18

dataフォルダにダウンロードできたことを確認します。




4.データの読み込み

データを読み込みます。
今回は【4-0】第4章を進めていく上での準備で行った方法で進めます。

View Fileでデータを確認します。

スクリーンショット 2019-11-04 7.44.03



データが入っているセルを確認します。
B2からC62までデータが入っています(B2:C22と表記)
スクリーンショット 2019-11-04 1.12.29


次にImport Datasetでデータを取り込みます。
スクリーンショット 2019-11-03 18.53.56

Import画面ではName, Sheet,Rangeを指定します。

Name:ハルさんのサイトと同じborgとします(大文字・小文字は別物とされます)
Sheet:このExcelは1つしかデータがないのでDefaultのままでOK
Range:先ほど確認したB2:C22

スクリーンショット 2019-11-04 7.41.10

Importボタンを押す前に右にあるコードをコピーしスクリプトファイルに貼り付けることも忘れずに行います。
library(readxl)
borg <- read_excel("data/demo-wilcoxon-rank-sum-test.xlsx", 
                   range = "B2:C22")
View(borg)

データが正しく入っていることを確認します。
スクリーンショット 2019-11-04 7.46.22

これでデータの取り込みは完成です。


5.データを要約する

次にデータを要約します。

データの要約は【4-1】Rでt検定を行う方法で紹介しました。
group_by関数とsummarize関数を使って要約しましたが、今回はsummary関数を使います。
summary関数は平均、中央値、最大・最小値、四分位範囲をまとめて出してくれますが、標準偏差はだしてくれません。ただ今回はノンパラメトリックなので標準偏差はいらないだろうという理由です。

そして今回はwideデータになっています。
スクリーンショット 2019-11-04 8.05.00
なのでコードはこうなります。
summary(borg$pre)
summary(borg$post)
スクリーンショット 2019-11-04 19.56.42


6.Wilcoxon符号付順位和検定を行う

次にWilcoxon符号付順位和検定を行います。

【4-3】Rで対応のあるt検定を行う方法で紹介しましたが、対応のあるとなしはt.test関数paired = TRUEをつけるかどうかの違いでした。

実はWilcoxon符号付順位和検定も同じです。

wilcox.testにpaired = FALSEをつける(もしくは何も付けない)とMann-Whitney U 検定
wilcox.testにpaired = TRUEをつけるとWilcoxon符号付順位和検定


ということで、wilcoxテストを行ってみます。
t.test関数もそうでしたが、longデータとwideデータで書き方が違います。

wideデータの場合 → , を使う
wilcox.test(1列目, 2列目, paired = TRUE)

longデータの場合 → を使う
wilcox.test(数値 ~ グループ, paired = TRUE)


今回はwideデータなのでこうなります。
wilcox.test(borg$pre, borg$post, paired = TRUE)
スクリーンショット 2019-11-04 20.38.57

EZRと同じ結果になりましたが、Mann-Whitney U 検定のときに悩ませたたアレが出てきました。
今度はもう1行増えてます。
タイがあるため、正確な p 値を計算することができません 
ゼロ値のため、正確な p 値を計算することができません
ちなみにEZRでも警告が出ています。
スクリーンショット 2019-11-04 20.59.27

cannot compute exact p-value with ties
cannot compute exact p-value with zeroes


7.タイのあるデータの対処法

EZRでも使われているwilcox.test関数はタイ(同順位)があると正確なp値を計算できず、近似値を計算する設定になっていました。





今回も同じ問題が出ています。
ブログではある程度のn数があればEZR(wilcox.test)でもいいのではという話がありました。
ただ警告が気持ち悪い!正確なp値も知りたいというための方法も紹介します。

タイに対しては奥村先生の記事が参考になります。



①coinパッケージのwilcoxsign_test関数

Mann-Whitney U 検定ではタイがあっても正確なp値を計算するcoinパッケージwilcox_test関数がありました。coinパッケージを使ってWilcoxon符号付順位和検定を行う場合はwilcoxsign_test関数を使います。

まだcoinパッケージを1度も使ったことがなければインストールします。

coinパッケージのインストール
install.packages("coin")

wilcoxsign_test関数の書き方はちょっとクセがあります。。。

スクリーンショット 2019-11-04 23.11.32
#パッケージの読み込み
library(coin)

wilcoxsign_test(borg$pre ~ borg$post, distribution = "exact", zero.method="Wilcoxon")
スクリーンショット 2019-11-04 23.12.57


exactRankTestsパッケージのwilcox.exact関数

もう1つexactRankTestsパッケージがあります。
このパッケージは開発が終わっており、インストールするとcoinパッケージ使ってねと警告が出ます。それでもcoinパッケージのwilcoxsign_test関数と同じ結果になります。
スクリーンショット 2019-11-05 0.06.33

まずexactRankTestsパッケージをインストールします。
#exactRankTestsパッケージのインストール
install.packages("exactRankTests")
wilcox.exact関数はwilcox.testと似たような書き方ができるのでわかりやすいのが特徴です。
スクリーンショット 2019-11-04 23.53.16
library(exactRankTests)
wilcox.exact(borg$pre, borg$post, paired = TRUE, exact=TRUE)
スクリーンショット 2019-11-04 23.55.13

どちらも同じになりました。


8.まとめ
今回はWilcoxon符号付順位和検定を紹介しました。
Mann-Whitney U 検定との共通点や相違点を比較するとイメージが深まると思います。

4章を順に見ていくと重複する箇所も出てきますので検索で来られた方はサイトマップを見ていただければ別の発見があるかもしれません。

 


9.今回使ったコード

今回使ったコードをまとめて置いておきます。

#データのダウンロード
url <- "https://haru-reha.com/wp-content/uploads/2018/04/demo-wilcoxon-rank-sum-test.xlsx"
destfile = "data/demo-wilcoxon-rank-sum-test.xlsx"

download.file(url, destfile)


library(readxl)
borg <- read_excel("data/demo-wilcoxon-rank-sum-test.xlsx", 
                   range = "B2:C22")
View(borg)

summary(borg$pre)
summary(borg$post)

summary_pre <- summary(borg$pre)
summary_post <- summary(borg$post)


#Wilcoxon符号付順位和検定
wilcox.test(borg$pre, borg$post, paired = TRUE)


#coinパッケージのインストール
install.packages("coin")

library(coin)
wilcoxsign_test(borg$pre ~ borg$post, distribution = "exact", zero.method="Wilcoxon")

#exactRankTestsパッケージのインストール install.packages("exactRankTests") library(exactRankTests) wilcox.exact(borg$pre, borg$post, exact=TRUE, paired = TRUE)