第4章は統計を扱います。今回も「シロート統計学」のハルさんとコラボレーションすることとなりました。シロート統計学はEZRを使った統計分析をわかりやすく解説されています。第4章はシロート統計学で使われていたEZRをRで行うとどうなるのか?といった視点で進めていきます。今回使うデータもハルさんのサイトと同じものを使わせでいただく事になりました。それぞれ見比べることで参考にしてみてください!今回はWilcoxon符号付順位和検定を紹介しますまずWilcoxon符号付順位和検定についてはハルさんのサイトをご参照ください。1.準備〜4.データの読み込みまでは【4-1】Rでt検定を行う方法と全く同じ流れになります。1〜4まででわからない部分があれば確認してみてください。第4章は毎回ExcelデータをダウンロードしてRを使うのでプロジェクトで管理して行うことを勧めています。ここではR練習というプロジェクトを作り、Excelファイルを入れるためのdataフォルダを作っています。これを前提に次から進めていきます。次にRのコードを書くためのスクリプトファイルを作ります。3.データのダウンロード今回もハルさんのサイトのデータを使わせていただきます。デモデータ(wilcoxon符号付順位和検定)この章ではRを使ってダウンロードしています。リンクのアドレスをコピーdestfileは保存場所と保存のファイル名を指定します。実際のコードは以下になります。"の中)とdestfileのdata/以降を変更するだけでOKです。destfile)dataフォルダにダウンロードできたことを確認します。4.データの読み込みデータを読み込みます。今回は【4-0】第4章を進めていく上での準備で行った方法で進めます。Fileでデータを確認します。データが入っているセルを確認します。Datasetでデータを取り込みます。Sheet,Rangeを指定します。Name:ハルさんのサイトと同じborgとします(大文字・小文字は別物とされます)Sheet:このExcelは1つしかデータがないのでDefaultのままでOKRange:先ほど確認したB2:C22Importボタンを押す前に右にあるコードをコピーしスクリプトファイルに貼り付けることも忘れずに行います。"B2:C22")View(borg)データが正しく入っていることを確認します。これでデータの取り込みは完成です。5.データを要約する次にデータを要約します。データの要約は【4-1】Rでt検定を行う方法で紹介しました。group_by関数とsummarize関数を使って要約しましたが、今回はsummary関数を使います。まとめて出してくれますが、標準偏差はだしてくれません。今回はノンパラメトリックなので標準偏差はいらないだろうという理由です。今回はwideデータになっています。コードはこうなります。summary(borg$pre)summary(borg$post)6.Wilcoxon符号付順位和検定を行う次にWilcoxon符号付順位和検定を行います。TRUEをつけるかどうかの違いでした。実はWilcoxon符号付順位和検定も同じです。TRUEをつけるとWilcoxon符号付順位和検定ということで、wilcoxテストを行ってみます。t.test関数もそうでしたが、longデータとwideデータで書き方が違います。TRUE)今回はwideデータなのでこうなります。検定のときに悩ませたたアレが出てきました。今度はもう1行増えてます。値を計算することができませんちなみにEZRでも警告が出ています。zeroes7.タイのあるデータの対処法EZRでも使われているwilcox.test関数はタイ(同順位)があると正確なp値を計算できず、近似値を計算する設定になっていました。今回も同じ問題が出ています。ブログではある程度のn数があればEZR(wilcox.test)でもいいのではという話がありました。警告が気持ち悪い!正確なp値も知りたいというための方法も紹介します。タイに対しては奥村先生の記事が参考になります。検定ではタイがあっても正確なp値を計算するcoinパッケージのwilcox_test関数がありました。coinパッケージを使ってWilcoxon符号付順位和検定を行う場合はwilcoxsign_test関数を使います。coinパッケージを1度も使ったことがなければインストールします。coinパッケージのインストールinstall.packages("coin")wilcoxsign_test関数の書き方はちょっとクセがあります。zero.method="Wilcoxon")2exactRankTestsパッケージのwilcox.exact関数もう1つexactRankTestsパッケージがあります。このパッケージは開発が終わっており、インストールするとcoinパッケージ使ってねと警告が出ます。それでもcoinパッケージのwilcoxsign_test関数と同じ結果になります。exactRankTestsパッケージをインストールします。#exactRankTestsパッケージのインストールinstall.packages("exactRankTests")wilcox.exact関数はwilcox.testと似たような書き方ができるのでわかりやすいのが特徴です。exact=TRUE)どちらも同じになりました。8.まとめ今回はWilcoxon符号付順位和検定を紹介しました。検定との共通点や相違点を比較するとイメージが深まると思います。4章を順に見ていくと重複する箇所も出てきますので検索で来られた方はサイトマップを見ていただければ別の発見があるかもしれません。9.今回使ったコード今回使ったコードをまとめて置いておきます。第4章は統計を扱います。今回も「シロート統計学」のハルさんとコラボレーションすることとなりました。シロート統計学はEZRを使った統計分析をわかりやすく解説されています。第4章はシロート統計学で使われていたEZRをRで行うとどうなるのか?といった視点で進めていきます。今回使うデータもハルさんのサイトと同じものを使わせでいただく事になりました。それぞれ見比べることで参考にしてみてください!今回は対応のあるt検定を紹介しますまず対応のあるt検定についてはハルさんのサイトをご参照ください。1.準備〜4.データの読み込みまでは【4-1】Rでt検定を行う方法と全く同じ流れになります。1〜4まででわからない部分があれば確認してみてください。第4章は毎回ExcelデータをダウンロードしてRを使うのでプロジェクトで管理して行うことを勧めています。ここではR練習というプロジェクトを作り、Excelファイルを入れるためのdataフォルダを作っています。これを前提に次から進めていきます。次にRのコードを書くためのスクリプトファイルを作ります。3.データのダウンロード今回もハルさんのサイトのデータを使わせていただきます。デモデータ(対応のあるt検定)この章ではRを使ってダウンロードしています。リンクのアドレスをコピーdestfileは保存場所と保存のファイル名を指定します。実際のコードは以下になります。"の中)とdestfileのdata/以降を変更するだけでOKです。destfile)dataフォルダにダウンロードできました!4.データの読み込みデータを読み込みます。今回は【4-0】第4章を進めていく上での準備で行った方法で進めます。Fileでデータを確認します。データが入っているセルを確認します。Datasetでデータを取り込みます。Sheet,Rangeを指定します。Name:ハルさんのサイトと同じgaitとします(大文字・小文字は別物とされます)Sheet:このExcelは1つしかデータがないのでDefaultのままでOKRange:先ほど確認したB2:D32Importボタンを押す前に右にあるコードをコピーしスクリプトファイルに貼り付けることも忘れずに行います。"B2:D32")View(gait)データが正しく入っていることを確認します。これでデータの取り込みは完成です。正規性の確認は【4-1】Rでt検定を行う方法でも紹介しました今回ハルさんのサイトではdifferenceの列の正規性をQQプロットを使って確認しています。QQプロットの出し方はいくつかありますが、ここでは2つ紹介します。qqline関数です。gait$difference)もしEZRと同じグラフが出したい場合はcarパッケージのqqPlot関数を使います(Pが大文字なので注意!)Packagesにcarパッケージが入っていれなければinstall.packages関数でインストール後library関数で呼び出します。gait$difference)こうすることでEZRを開かなくても同様のことが行なえます。確認するのが目的であればEZRでいいと思いますが、1年後に再現しようと思ったらスクリプトに残しておくと再現性が高まります。ヒストグラムはgeom_histgramもしくはhist関数を使います。5)hist(gait$difference)シャピロウィルク検定は【4-1】Rでt検定を行う方法でも紹介しましたがshapiro.testでした。TRUEを足すだけです。FALSEをしていることになります。FALSEになります。t検定は2種類の書き方があります。その違いはデータの形にあります。ハルさんのサイトでEZRでは右の形しかできませんでしたがRならどちらも可能です。ちなみに左をlongデータ(縦に長い)、右をwideデータ(横に広い)といいます。TRUEを付け加えます。longデータで行う場合はAの並び順とBの並び順が違うと計算結果が異なりますので注意が必要です。今回はwideデータなので以下のようになります。TRUE)EZRと同じ結果になりました。86.86993となっており、差の平均(期待値)は64.16667となっています。pre)が0であると仮定しています。95%信頼区間に0が含まれているとpは0.05以上と判断されます。〜50といった具合です。上のsampleの場合、差の平均(期待値)は20とプラスなのでpostの方がより高い値と言えそうですが、95%信頼区間に0やマイナスが含まれています。言い換えると20だったけど-10にだってなり得るし、0にだってなり得るとなります。なればpostの方が高い値といえません。ちなみに今回は差(期待値)が0を仮定していたので0でしたが、例えばロジスティック回帰分析だとodd比が1であると仮定するので、そのときは1を挟むかどうかを確認することになります。p値だけでなく信頼区間を確認するとまた発見があるかもしれません。7.まとめ今回は対応のあるt検定を行いました。TRUEを付け加えるだけでした。4章を順に見ていくと重複する箇所も出てきますので検索で来られた方はサイトマップを見ていただければ別の発見があるかもしれません。8.今回使ったコード今回使ったコードをまとめて置いておきます。95%信頼区間のコードも載せています。第4章は統計を扱います。今回「シロート統計学」のハルさんとコラボレーションすることとなりました。ハルさん、ありがとうございます!シロート統計学はEZRを使った統計分析をわかりやすく解説されています。第4章はシロート統計学で使われていたEZRをRで行うとどうなるのか?といった視点で進めていきます。今回使うデータもハルさんのサイトと同じものを使わせでいただく事になりました。検定についてはハルさんのサイトをご参照ください。1.準備〜4.データの読み込みまでは【4-1】Rでt検定を行う方法と全く同じ流れになります。1〜4まででわからない部分があれば確認してみてください。第4章は毎回ExcelデータをダウンロードしてRを使うのでプロジェクトで管理して行うことを勧めています。ここではR練習というプロジェクトを作り、Excelファイルを入れるためのdataフォルダを作っています。これを前提に次から進めていきます。次にRのコードを書くためのスクリプトファイルを作ります。3.データのダウンロード今回もハルさんのサイトのデータを使わせていただきます。U検定)前回Rを使ってダウンロードしましたが、まだ慣れないと思いますので再度説明します。リンクのアドレスをコピーdestfileは保存場所と保存のファイル名を指定します。実際のコードは以下になります。"の中)とdestfileのdata/以降を変更するだけでOKです。destfile)dataフォルダにダウンロードできました!4.データの読み込みデータを読み込みます。今回は【4-0】第4章を進めていく上での準備で行った方法で進めます。Fileでデータを確認します。データが入っているセルを確認します。Datasetでデータを取り込みます。Sheet,Rangeを指定します。Name:ハルさんのサイトではMMTでしたがここではmmtとします(大文字・小文字は別物とされます)Sheet:このExcelは1つしかデータがないのでDefaultのままでOKRange:先ほど確認したC2:C62Importボタンを押す前に右にあるコードをコピーしスクリプトファイルに貼り付けることも忘れずに行います。"C2:D62")View(mmt)データが正しく入っていることを確認します。これでデータの取り込みは完成です。検定(マン・ホイットニーのU検定)は別な言い方もあります。検定はwilcox.test関数を使っています。mmt)EZRの結果はこうです。EZRと同じ結果がでました。警告が気になりますがひとまず置きます。p値だけがほしいときは変数名を付けて変数名$p.valueとします。mmt)mmt_wilcox$p.value6.グラフの作成EZRでは同時に箱ひげ図を作成しますがRでは自作する必要があります。箱ひげ図の作成に関しては第3章で紹介しています。MMT))もしEZRと同じグラフにしたい場合はboxplot関数を使うという選択肢もあります。検定を行う方法を紹介しました。wlcox.testを使えばEZRと同じ結果になります。警告が出てしまい心配な方はいると思います。以下からは少し応用編となりますがこの警告について進めていきます。値を計算することができませんRでノンパラメトリックの検定を行っているとこういった画面が出ることがあります。値を計算することができません正確なp値が計算できないと言われるととても心配になります。実は先程のEZRでも警告は出ていました。ties8.タイ(ties)って何?タイというのは同順(同じ数値)のことです。検定は平均ではなくデータを少ない順に並べてデータに順位を付け計算するため、同順位があるとそれをどう計算するかといった問題や様々な対処法が出てきます。9.タイに対してどう考えるか?正確なp値を求めるには膨大な計算が必要になるらしく、そのため正規近似というものを使って計算するようです。EZR(wilcox.test)はタイがなければ正確なp値を求めるオプションがありますが、タイがあるとオプションを付けても正確なp値を求めることができないようです。近似値ではなく正確なp値を求めるにはcoinパッケージのwilcox_test関数を使う必要があります。RjpWikiの記事に詳しく説明があります。10.coinパッケージを使ってみるwilcox_test関数を使うためにはcoinパッケージをインストールする必要があります。パッケージのインストールは今までinstall.packages関数を紹介していましたが、今回はRStudioのPackagesタブからインストールしてみます。PackagesタブのInstallをクリックしてパッケージ名を入れます。これでインストールは完了です。これでlibrary関数を使ってcoinパッケージを呼び出せます。wilcox.testとの違いは主に3点あります。distribution="exact")mmt_coinもしp値だけを出したいときはpvalue関数、95%信頼区間を出したいときはpvalue_interval関数を使います。pvalue(mmt_coin)pvalue_interval(mmt_coin)これで正確なp値を求めることができました。11.じゃあEZRはダメなのか?ただ必ず正確なp値を求めないといけないかというとそうではないようです。さんからアドバイスを頂きました。11:05:21ある程度n数があればEZRでも利用可能と考えることもできます。このあたりは統一見解があるわけではなさそうです。そもそもp値がより小さければより差が大きいわけではありませんので、ここのp値をどうするかよりも研究デザインが適切かどうかの方が大切になるかもと言うのが個人的なイメージです。実はwilcox.test(EZR)でタイが無いときの正確なp値の求め方などwilcox.testのオプションなど話は割愛させてもらっています。qiitaの記事ではwilcox.testとwilcox_testのオプションやタイの有無で結果がどう変わるかについての解説をしていますので、興味があればぜひこちらも併せてご覧ください。13.今回使ったコード今回使ったコードをまとめて置いておきます。coinパッケージのインストールに関してもコードにしています。検定を行うためにcoinパッケージのインストール#1度でもしていればしなくてOKinstall.packages("coin")#もしtidyverseパッケージをインストールしていなければインストール。第4章は統計を扱います。今回「シロート統計学」のハルさんとコラボレーションすることとなりました。ハルさん、ありがとうございます!シロート統計学はEZRを使った統計分析をわかりやすく解説されています。第4章はシロート統計学で使われていたEZRをRで行うとどうなるのか?といった視点で進めていきます。今回使うデータもハルさんのサイトと同じものを使わせでいただく事になりました。それぞれ見比べることで参考にしてみてください!今回はt検定を紹介しますまずt検定についてはハルさんのサイトをご参照ください。第4章は毎回ExcelデータをダウンロードしてRを使うのでプロジェクトで管理して行うことを勧めています。ここではR練習というプロジェクトを作り、Excelファイルを入れるためのdataフォルダを作っています。これを前提に次から進めていきます。次にRのコードを書くためのスクリプトファイルを作ります。3.データのダウンロード今回はハルさんのサイトのデータを使わせていただきます。デモデータ(t検定)これをダウンロードしてdataフォルダに入れればいいのですが実はRでできてしまいます。download.file関数を使います。"を忘れないようにしてください。destfile)以下説明します。リンクのアドレスをコピーdestfileは保存場所と保存のファイル名を指定します。保存場所は今回プロジェクトを使っているのでR練習フォルダになります。加えてdata/を付け足すことでR練習フォルダ内にあるdataフォルダという意味になります。ファイル名は自由に決められますが今回は元のファイルと同じにしました。拡張子も忘れないようにしましょう。プロジェクトを使っていなければ保存場所はgetwd関数で出てきたフォルダになります。getwd関数の()には何も入れません。getwd()この方法を使う最大のメリットは、次回使うExcelデータはurlの部分を変えるだけでできてしまうことです。保存したファイルを指定したところに移動させて・・・といった手作業必要ありません。いった作業もRで行っていくことでRにも早く慣れてくると思います。4.データの読み込みデータを読み込みます。今回は【4-0】第4章を進めていく上での準備で行った方法で進めます。Fileでデータを確認します。今回は握力のデータです。A群とB群を比較します。A1にデータが入っていません。実際にはB2からC62までデータが入っていることを確認します。Datasetを選びます。1データの名前(変数名)を付ける何でもいいのですが今回はハルさんのサイトと同じgripにしました。2読み込む範囲を指定する今回A1からのデータではないので先程確認したB2からC62を指定します。B2:C62のように左上と右下を:でつなげます。右下にコードが自動的に作られます。Importを押せば完了なのですが、このコードをコピーしスクリプトに貼り付けておけば1年たった後でも同じことができます。EZRでもスクリプトを保存することができないわけではないのですが、再現性(後でしても、他の人がしても同じ事ができる)を保つためにもこういったコードを残しておく習慣をつけるようにしましょう。コードの一番下にあるView関数を使うことでRStudio内でもデータの確認ができます。このタブを消してもデータに影響はありません。1回View関数を使えばまた表示できます。View関数はEZRで言う表示と同じです(下図はハルさんのサイトより。比べてみてください)これでデータの取り込みは完了です!5.データの要約ハルさんは次にデータの要約をしています。EZRでのデータの要約と全く同じ機能はないですが、第2章で紹介したtidyverseパッケージのgroup_by関数とsummarize関数が使えます。group_by関数とsummarize関数に関してはこちらで紹介しています。%>%やgroup_by関数、summarize関数はtidyverseパッケージに含まれていますのでtidyverseパッケージを呼び出します。tidyverseパッケージを全く使ったことが無い方はパッケージをインストールします。1度でも使ったことがあれば次の1行は必要ありません。install.packages("readxl")実際のコードは以下になります。イメージ図も添付します。n())これでA群とB群のデータのばらつきを確認することができます。ちなみにこのコードをコピーして色がついた箇所を変更すれば他の場面でも使えます!このままでもいいのですが、データ要約は後でグラフ作成に使うのでgrip_summaryという名前をつけます。次に正規性の確認を行います。ハルさんのサイトではヒストグラムを作成しました。ヒストグラムの作り方はこちらで紹介しています。"dodge"を使います。使うので注意してください。5)加えてシャピロ・ウィルク検定も紹介されています。シャピロ・ウィルク検定はshapriro.test関数を使います。A群とB群それぞれで行いますのでgrip$gripの列をA群とB群に分ける必要があります。ハルさんのサイトのEZRで行う場合と見比べてみてください。書いてあるところの意味が見えてきます。ちなみに変数(1つ選択)のgripはgrip$gripの色のついた部分です。結果は以下のとおりです。どちらも0.5を超えているので正規分布であるという仮説は棄却されませんでした。このあたりの解釈はハルさんのサイトをご参照ください。"B"])上記のように1つずつ計算する方法もいいのですが群の数だけ繰り返します。プログラミングであるRは繰り返しに強いという特徴があります。群ごとにデータを分割し、まとめて計算する方法として今回はsplit関数とmap関数を使います。map(~shapiro.test(.$grip))まずgripのデータをsplit関数を使ってA群とB群の2つのデータに分割し、map(~シャピロウィルク検定)でシャピロウィルク検定を繰り返します。EZRだと1つ1つ検定を繰り返す必要があります。Rを活用すると1度に作業が終わるため、繰り返すことでのやり間違えの予防や時間の短縮にも繋がります。8.等分散性を確認する先程は正規性を調べましたが、今回は等分散性を調べます。F検定と言いますがRではvar.test関数を使います。Rで行うとこうなります。使えば$が要らなくなります。0.05だったので等分散性は棄却されませんでした。9.t検定を行ういままで正規分布かどうか、等分散性かどうかを確認しました。これでt検定が使えます。t検定はt.test関数を使います。t検定は2種類の書き方があります。その違いはデータの形にあります。ハルさんのサイトでEZRでは左の形しかできませんでしたがRならどちらも可能です。ちなみに左をlongデータ(縦に長い)、右をwideデータ(横に広い)といいます。今回はlongデータなのでこのようになります。gripを使うとgrip$は外せます。どちらでも好きな方で大丈夫です。grip)EZRと同じ結果になりました!10.グラフの作成EZRでは丁寧にグラフも作ってくれますがRでは自作する必要があります。EZRでは平均と標準偏差の棒グラフを作成していました。先程のgrip_summaryに平均と標準偏差を求めていたのでそれを使います。棒グラフの作り方はここで紹介しています。次は以下を修正します。"grip")だいぶ見た目が良くなりました。3EZRのグラフにできるだけ近づけるEZRのグラフに似せるためにもうひと工夫します。c(0,50))これでEZRにだいぶ似ました。11.まとめかなり長くなりましたが今回はt検定を紹介しました。最初だったので1つ1つ説明しました。ボリュームが多すぎて慣れないところもあると思いますが、これから何度も出てくるものもありますので少しずつ慣れて貰えればと思います。12.コードをまとめました今回のコードを1つにまとめました。パッケージはlibrary()を1度だけ行えばいいので最初に出しています。第4章は統計を扱います。今回「シロート統計学」のハルさんとコラボレーションすることとなりました。ハルさん、ありがとうございます!シロート統計学はEZRを使った統計分析をわかりやすく解説されています。第4章はシロート統計学で使われていたEZRをRで行うとどうなるのか?といった視点で進めていきます。今回使うデータもハルさんのサイトと同じものを使わせでいただく事になりました。それぞれ見比べることで参考にしてみてください!今回は4章を進めるにあたって抑えておきたいRの基本的な操作について復習します。1.プロジェクトを作る第4章では色々なデータを扱うため、フォルダやRのファイルが煩雑になりやすいです。そこで今回はプロジェクトを使って管理します。RStudioでプロジェクトを作成する方法は第1章で紹介しています。今回は「R練習」というプロジェクトを作りました。データを入れるフォルダを作ります。分析に使うExcelファイルはデスクトップやマイドキュメントに置かずにここに置くようにします!そしてRのコードを書くスクリプトファイルを作ります。今回は【4-1】で使うためのスクリプトファイルを作成しました。これで完成です!2.Excelのデータを読み込むExcelを読み込む方法はいくつかありますが、ここではコードを書かない方法を紹介します。DatasetでExcelファイルをRで読み込みます。Datasetを選びます。Name:Rで使うときのデータの名前(自由に決められる)sheet:Excelに複数タブあればどれを使うか指定できるRange:データの読み込む範囲を指定できる主には上記3つで読み込む範囲を指定します。1枚目のタブでデータの読み込む範囲の左上がA1だったらNameの指定だけでも大丈夫です。この方法の良いところは右にコードを自動で書いてくれることです。これを忘れずにスクリプトにコピーします。することで2回目以降は1秒もかからずに行えるのでかなりの時短になります!まとめ今後色々な検定やデモデータを使うことになります。EZRを使うとあまり意識されない部分ですが、プロジェクトファイルを作ることでデータの管理やスクリプトファイルの管理が進めやすくなります。次回はt検定を紹介します。