1

aov() 用に書かれたコードを car-package の Anova() 関数に変換したい。

anovadata3 <- within(anovadata3, {
  subject <- factor(subject)
  time <- factor(time)
  gender <- factor(gender)
  group <- factor(group)
  groupgender <- factor(groupgender)
})

anovadata3.aov <- aov(values ~ time*group*gender + Error(subject),
                      data = anovadata3)
summary(anovadata3.aov)

このコードにより、次の出力が得られます。

Error: subject
              Df Sum Sq Mean Sq F value  Pr(>F)   
group          1  32220   32220   8.632 0.00365 **
gender         1     30      30   0.008 0.92819   
group:gender   1     15      15   0.004 0.94952   
Residuals    221 824913    3733                   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Error: Within
                   Df Sum Sq Mean Sq F value   Pr(>F)    
time                3  21160    7053   9.223 5.53e-06 ***
time:group          3  18338    6113   7.993 3.06e-05 ***
time:gender         3   1916     639   0.835  0.47486    
time:group:gender   3  11679    3893   5.091  0.00172 ** 
Residuals         663 507012     765                     
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

ただし、 car-package の Anova() 関数を使用しようとすると、次のようになります。

library(car)
mlm <- lm(values ~ time*group*gender, data = anovadata3)
rfactor <- factor(c("time", "gender","group","groupgender","subject"))
anovadata3.aov <- Anova(mlm, idata = dataframe(rfactor), idesign = ~ rfactor, type ="III")
summary(anovadata3.aov)

これにより、この出力が得られます。

     Sum Sq                Df         F value              Pr(>F)       
 Min.   :     58.3   Min.   :  1   Min.   :  0.03871   Min.   :0.00000  
 1st Qu.:   1324.9   1st Qu.:  1   1st Qu.:  0.54230   1st Qu.:0.04997  
 Median :  10281.5   Median :  3   Median :  1.57697   Median :0.21357  
 Mean   : 196286.3   Mean   :100   Mean   : 34.03103   Mean   :0.31053  
 3rd Qu.:  12290.4   3rd Qu.:  3   3rd Qu.:  2.61758   3rd Qu.:0.50989  
 Max.   :1331924.5   Max.   :884   Max.   :262.67095   Max.   :0.84408  
                                   NA's   :1           NA's   :1

aov() に使用するコードを Anova() に合わせて作り直す方法を知っている人はいますか。https://gribblelab.wordpress.com/2009/03/09/repeated-measures-anova-using-r/のチュートリアルに従ってみました 。

Anova() を正しく取得しようとします。しかし、似たような出力は得られません。また、ウェブページから、Mauchlys と Greenhouse が提供されていると思われることがわかりますが、これはわかりません。また、Anovaの結果でイータ二乗を取得する方法を知っている人はいますか? または、別の関数を使用して eta (etaSquared()) を計算する必要があります。

以下のデータはテストに使用されました。時間、性別、グループ間の「値」と、因子間の交互作用効果に有意差があるかどうかをテストしようとしています。

          values  testperiod subject gender group groupgender time
1   118.82660110     Pretest       1      2     2   BSTfemale    1
2    61.07615138     Pretest       2      2     2   BSTfemale    1
3    57.51022740     Pretest       3      2     2   BSTfemale    1
4    70.73637347     Pretest       4      2     2   BSTfemale    1
5     9.86907880     Pretest       5      2     2   BSTfemale    1
6    64.51579546     Pretest       6      2     2   BSTfemale    1
7    63.25669342     Pretest       7      2     2   BSTfemale    1
8   109.09354856     Pretest       8      2     2   BSTfemale    1
9   140.69340502     Pretest       9      2     2   BSTfemale    1
10   93.94269807     Pretest      10      2     2   BSTfemale    1
11   43.76802256     Pretest      11      2     2   BSTfemale    1

...
898  60.85271722 FU_12_month     223      1     2     BSTmale    4
899  82.75598576 FU_12_month     224      1     2     BSTmale    4
900 -32.38497309 FU_12_month     225      1     2     BSTmale    4
4

1 に答える 1

1

何が望まれているのかという質問は少し曖昧に思えますが、R の多変量関数から望むものを得るのに苦労したのはこれが初めてではありません。eta-squaredこの用語が「多変量検定:」のセクションで満たされることを願っています。その結果summary.Anova.mlmでmultivariate=TRUEは、Pillai、Wilks、Hotelling-Lawley、および Roy の検定が返され、次に (multivariateパラメーターの値に関係なく) 「 Mauchly Tests for Sphericity」および「Greenhouse-Geisser and Huynh-Feldt Corrections for Departure from Sphericity」。結果をコンソール セッションから取得したいだけなのか、それとも結果をさらに処理したいのか、私にはわかりません。(そして、後者の希望により、過去に Rhelp で複数のリクエストが発生しました。)

car::Anova のヘルプ ページのコードを使用して、必要に応じてさらに議論するための例を作成できます (コピーできるデータ例が含まれていないため)。

phase <- factor(rep(c("pretest", "posttest", "followup"), c(5, 5, 5)),
    levels=c("pretest", "posttest", "followup"))
hour <- ordered(rep(1:5, 3))
idata <- data.frame(phase, hour)
idata
mod.ok <- lm(cbind(pre.1, pre.2, pre.3, pre.4, pre.5, 
                     post.1, post.2, post.3, post.4, post.5, 
                     fup.1, fup.2, fup.3, fup.4, fup.5) ~  treatment*gender, 
                data=OBrienKaiser)
class(AnoOBK <- Anova(mod.ok, idata= idata, idesign= ~phase*hour, type="III") # "Anova.mlm"

summary(AnoOBK)

> names(summary(AnoOBK))
[1] "type"               "repeated"           "multivariate.tests"
[4] "univariate.tests"   "pval.adjustments"   "sphericity.tests"  
[7] "SSPE"
于 2016-01-30T20:39:38.773 に答える