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

統計学備忘録(R言語のメモ)

since2016 ときどきTEXのメモ

ANOVA君 で分散分析

投稿日2017.8.25, 更新日2020.8.17

ANOVA君はフリーの統計ソフトウェア「R」で動作する分散分析関数です。ダウンロード、使い方、リリース情報などは「井関龍太のページ」を必ずご確認ください.

ANOVA君 - 井関龍太のページ

このブログではanovakun version 4.8.5を使用させていただきます 。このテキストファイルにプログラムが記載してあり、マニュアルにもなります。

まずはテキストファイル(コード)を読み込みます

source("C:\\Users\\-----------\\anovakun_485.txt")

A~Zが要因、被験者間要因は「s」の左側,被験者内要因は「s」の右側に記載、数字は各要因の水準の数。入れ子の最も外側にある要因から順番に「A」~「Z」のラベルが対応。

一元配置分散分析と多重比較

n=30, 繰り返し数10の例

testA <- c(108, 105, 98, 112, 110, 109, 115, 120, 107, 114)
testB <- c(91, 93, 89, 109, 98, 88, 102, 110, 97, 109)
testC <- c(79, 90, 87, 87, 78, 75, 96, 109, 89, 90)
dat1 <- data.frame(testA , testB , testC)

#long型に変換
library(tidyr)
df <- gather(dat1, test, data) 

    test data
1  testA  108
2  testA  105
3  testA   98
4  testA  112
5  testA  110
6  testA  109
7  testA  115
8  testA  120

anovakun(df, "As", 3)
#被験者間要因は「s」の左側のA,3は各要因の水準の数(A, B, C)

多重比較:デフォルトはShafferの方法(Bonferroniの方法以外が可能)

反復測定と球面性検定とイプシロン修正

n=10の例

運動後1<- c(108, 105, 98, 112, 110, 109, 115, 120, 107, 114)
運動後2<- c(98, 98, 95, 98, 100, 100, 105, 100, 98, 109)
運動後3<- c(95, 87, 90, 84, 88, 89, 90, 89, 90, 99)
dat1 <- data.frame(運動後1,運動後2, 運動後3) #反復測定はこの型のまま使用する

   運動後1分 運動後2分 運動後31        108        98        95
2        105        98        87
3         98        95        90
4        112        98        84
5        110       100        88
6        109       100        89
7        115       105        90
8        120       100        89
9        107        98        90
10       114       109        99

anovakun(dat1, "sA" ,3 , mau = T)  
#被験者内要因は「s」の右側のA,3は各要因の水準の数(運動後1分, 2分, 3分)
#mau = T:Mauchlyの球面性検定(デフォルトは、Mendozaの多標本球面性検定)  

== Mauchly's Sphericity Test and Epsilons ==

------------------------------------------------------------------------
 Effect      W  approx.Chi  df      p         LB     GG     HF     CM 
------------------------------------------------------------------------
      A 0.4396      6.5750   2 0.0373 *   0.5000 0.6409 0.7008 0.6081 

分散分析はそれぞれの群が独立していることが前提となります。検定の結果(p値=0.0373)から球面性の仮定が棄却されるので、被験者内効果をイプシロン修正します。さらに一般的な手法としてMANOVA混合モデルがあります。

anovakun(dat1, "sA" ,3 , gg = T)

イプシロン調整済p値が算出される

== Adjusted by Greenhouse-Geisser's Epsilon ==

-------------------------------------------------------------
 Source         SS    df         MS  F-ratio  p-value      
-------------------------------------------------------------
      s   374.6667     9    41.6296                        
-------------------------------------------------------------
      A  1940.6000  1.28  1514.0480  68.0293   0.0000 ***  
  s x A   256.7333 11.54    22.2558                        
-------------------------------------------------------------     

もし球面性の仮定が成立(p>0.05)しているのであれば、独立ではないのですが一元配置分散分析の結果を利用することになります。

二元配置分散分析と交互作用

n=30の例

性別 <- c("男", "男", "男", "男", "女", "女", "女", "女", "女", "女")
A群 <- c(108, 105, 98, 112, 110, 109, 115, 120, 107, 114)
B群 <- c(91, 93, 89, 109, 98, 88, 102, 110, 97, 109)
C群 <- c(79, 90, 87, 87, 78, 75, 96, 109, 89, 90)
dat2 <- data.frame(性別, A群, B群, C群)

#long型に変換
library(tidyr)
df2 <- gather(dat2, test, data, -性別) 

   性別 test data
1    男  A群  108
2    男  A群  105
3    男  A群   98
4    男  A群  112
5    女  A群  110
6    女  A群  109
7    女  A群  115
8    女  A群  120
9    女  A群  107
10   女  A群  114

anovakun(df2, "ABs" ,2 ,3)

Rを使用して分散分析を実施するのであれば、ANOVA君をお勧めします。