No.04766 Re: anova()を使う際の因子の順番 【波音】 2007/11/21(Wed) 21:52
私はいつも次のようにしていますが，このやり方だと説明変数の順序を入れ換えても結果（分散分析表）は変わりませ んよ。同じデータできみかさんが提示してくださった方法でやっても，有意であるかどうか（すなわちp値）は同じでしたが，順番を入れ換える と"Resid. Df"と"Resid. Dev"の値は変わってしまいますね。
これについては「R側が計算過程でうんぬんかんぬん」といったことしかいえませんね。私にはなぜそうなるのか詳しいことは分かりません(^_^;) 他の先生方，補足お願いします。> YIELD
[1] 8.819 7.826 8.341 8.818 8.527 7.762 8.530 8.990 8.492 8.644
[11] 8.950 9.815 8.166 7.331 8.187 10.154 8.487 8.846 8.715 10.273
[21] 7.045 8.128 7.875 8.152
> BLOCK
[1] 1 1 1 1 1 1 1 1 2 2 2 2 2 2 2 2 3 3 3 3 3 3 3 3
Levels: 1 2 3
> TRTMT
[1] 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8 1 2 3 4 5 6 7 8
Levels: 1 2 3 4 5 6 7 8
> result <- glm(YIELD ~ BLOCK + TRTMT)
> summary(aov(result))
Df Sum Sq Mean Sq F value Pr(>F)
BLOCK 2 0.3937 0.1968 0.5193 0.60599
TRTMT 7 8.0776 1.1539 3.0442 0.03623 *
Residuals 14 5.3069 0.3791
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> result2 <- glm(YIELD ~ TRTMT + BLOCK)
> summary(aov(result2))
Df Sum Sq Mean Sq F value Pr(>F)
TRTMT 7 8.0776 1.1539 3.0442 0.03623 *
BLOCK 2 0.3937 0.1968 0.5193 0.60599
Residuals 14 5.3069 0.3791
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
No.04767 Re: anova()を使う際の因子の順番 【きみか】 2007/11/21(Wed) 22:21
波音さん，お返事ありがとうございます！
以下に自分の結果を示しました。
こちらでは，DとBDを入れ替えると，P値も変わってしまいました・・・。
全て連続変数です。DとBDの値は似通っていますが，どちらがFに効いているのか，知りたくて，この解析を行ないました。ですが，順番を入れ替えると結果が変わってしまい，どうしたものか，と思いまして・・・。> anova(glm(F~D+BD,family=inverse.gaussian),test='F')
Analysis of Deviance Table
Model: inverse.gaussian, link: 1/mu^2
Response: F
Terms added sequentially (first to last)
Df Deviance Resid. Df Resid. Dev F Pr(>F)
NULL 92 2.87205
D 1 0.50347 91 2.36858 19.4932 2.801e-05 ***
BD 1 0.10032 90 2.26825 3.8843 0.05181 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> anova(glm(F~BD+D,family=inverse.gaussian),test='F')
Analysis of Deviance Table
Model: inverse.gaussian, link: 1/mu^2
Response: F
Terms added sequentially (first to last)
Df Deviance Resid. Df Resid. Dev F Pr(>F)
NULL 92 2.87205
BD 1 0.51068 91 2.36137 19.7723 2.484e-05 ***
D 1 0.09312 90 2.26825 3.6052 0.0608 .
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
No.04768 Re: anova()を使う際の因子の順番 【青木繁伸】 2007/11/21(Wed) 23:30
私はglmもよく分かりませんが，きみかさんがやっているのは重回帰分析で言えば変数増加法と同じでしょうか。従 属変数をよく説明できる順で独立変数を加えていくと。"Terms added sequentially (first to last)" というのが，その説明なんでしょうか。
で，最初のNONEというのは説明変数を加えないもの，次の行にあるのが，その説明変数を加えたら前のモデルと比べて有意に優れたモデルであるかの検定ではないですか?
だとすれば，最初に取り込まれる方が有意になって当たり前でしょう。その次以降に取り込まれる独立変数はたとえ十分な説明力があっても，前に取り込まれた独立変数が説明してしまった残りを説明するしか道がないので，それが加えられてもモデルへの貢献は少ない。
波音さんがやっているやり方でよいのでは?
No.04769 Re: anova()を使う際の因子の順番 【波音】 2007/11/22(Thu) 02:21
"因子"といっておられましたので，てっきり説明変数の2つはカテゴリカル型のデータかと思っていました。説明変 数の順番を変えても結果に変わりはないと思っていましたが，連続変数であるデータセットで再度，（私が提示した方法と）同じ方法で説明変数を入れ換えて やってみると分散分析表の結果が変わってしまいますね。
http://www.oup.com/uk/orc/bin/9780199252312/01resources/datasets/excel/separate/にあるDatasheet chapter 04のワークシート内にあるSchool Children's Mathのデータを使用しました。> result <- glm(AMA ~ YEARS + HGT)carパッケージにあるAnova()を利用すると，どちらのモデルでも同じ分散分析表が出力されます。これを使うとSASでいうところの，タイプ2の平方和を計算してくれます。ただし，最初から実装されてはいないのでインストールしてください。
> summary(aov(result))
Df Sum Sq Mean Sq F value Pr(>F)
YEARS 1 422.60 422.60 1702.4271 <2e-16 ***
HGT 1 0.01 0.01 0.0317 0.8599
Residuals 29 7.20 0.25
> result2 <- glm(AMA ~ HGT + YEARS)
> summary(aov(result2))
Df Sum Sq Mean Sq F value Pr(>F)
HGT 1 412.77 412.77 1662.83 < 2.2e-16 ***
YEARS 1 9.84 9.84 39.63 7.109e-07 ***
Residuals 29 7.20 0.25> library(car) #carパッケージの読み込み連 続型のデータである変数を複数用いる場合は，こちらを使った方が良いということでしょうかね。少なくとも，summary(aov(result))とし て出力される分散分析表の値は参考書にあるものと一致しているので，AMAに対する分散分析表は間違いではないでしょうが，なぜ説明変数の順番を入れ換え ると，この方法で結果が異なってしまうのかは（私には）分かりません。。。
> Anova(result, test="F")
Anova Table (Type II tests)
Response: AMA
SS Df F Pr(>F)
YEARS 9.8376 1 39.6300 7.109e-07 ***
HGT 0.0079 1 0.0317 0.8599
Residuals 7.1989 29
> Anova(result2, test="F")
Anova Table (Type II tests)
Response: AMA
SS Df F Pr(>F)
HGT 0.0079 1 0.0317 0.8599
YEARS 9.8376 1 39.6300 7.109e-07 ***
Residuals 7.1989 29
ホント，何でなんでしょう？
No.04771 Re: anova()を使う際の因子の順番 【きみか】 2007/11/22(Thu) 09:46
波音さん，青木先生，お返事ありがとうございます！
説明不足，勉強不足で申し訳ありません。
できる限り，ちゃんと伝えられるように説明したいと思います。
長文お許しください。
説明変数の順番を入れ替えると，変わってしまうのは以下のような原因らしいです（特に後半が，ちゃんと理解できていないですが・・・）。
↓R-FAQからの引用です。
---------------------------------------------------------------------------
http://www.okada.jp.org/RWiki/?R-FAQ%C6%FC%CB%DC%B8%EC%CC%F5(7%BE%CF)
A+B+A:B のようなモデルでは, R は ~1, ~A, ~A+B および ~A+B+A:B の間の 2乗和の差をレポートすることになります. もしもモデルが ~B+A+A:B ならば, ~1, ~B, ~A+B, および ~A+B+A:B の間の (2乗和の) 差をレポートすることになります. 最初の例では, A に関する 2乗和は ~1 および ~A を比較していますが, 2番目の例では, ~B および ~B+A を比較しています. 非直交型のデザイン (つまり, 不均衡デザインの大部分) では, これらの比較は (概念的にも数値的にも) 異なります.
パッケージによっては, モデル全体と各因子を 1つずつ除いたモデル群との比較に基づく 2乗和をレポートするものもあります (たとえば, SAS の有名な「タイプ III の 2乗和」がそうです). このようなパッケージはモデルにおける因子の順番に依存しません. どちらの 2乗和の取り方が正しいのかという質問は, ときどき r-help で低レヴェルの聖戦を引き起こしています.
R がレポートする特定の 2乗和について動揺する必要はありません. まったく簡単に, お気に入りの方法で 2乗和を計算することができます. anova(model1, model2) で, 任意の 2つのモデルを比較することができ, また, drop1(model1) によって 単一の項目を落とした後の 2乗和を示すことができます.
---------------------------------------------------------------------------
そのため，タイプ3の平方和の計算の仕方で，分散分析を行なったら良いのかな，と思っていました。ですがどのようにしたら，良いかよくわかりません・・・。
以下1回目，波音さんがやっていたのと同じ方法で行なった結果です。これは，タイプ1での平方和なんでしょうか？> result<-glm(F~D+BD,family=inverse.gaussian)--------------------------------------------------
>
> result2<-glm(F~BD+D,family=inverse.gaussian)
>
> summary(aov(result))
Df Sum Sq Mean Sq F value Pr(>F)
D 1 471.8 471.8 4.0397 0.04747 *
BD 1 679.9 679.9 5.8211 0.01789 *
Residuals 89 10394.8 116.8
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
94 observations deleted due to missingness
> summary(aov(result2))
Df Sum Sq Mean Sq F value Pr(>F)
BD 1 1150.8 1150.8 9.8528 0.002301 **
D 1 0.9 0.9 0.0080 0.928988
Residuals 89 10394.8 116.8
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
94 observations deleted due to missingness
教えて頂いたタイプ2でもやってみました。> library(car)----------------------------------------------------------
> Anova(result, test="F")
Anova Table (Type II tests)
Response: F
SS DF F Pr(>F)
D 0.000005 1 0.0036 0.95201
BD 0.009537 1 6.3507 0.01352 *
Residuals 0.133655 89
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> Anova(result2, test="F")
Anova Table (Type II tests)
Response: F
SS Df F Pr(>F)
BD 0.009537 1 6.3507 0.01352 *
D 0.000005 1 0.0036 0.95201
Residuals 0.133655 89
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
入れ替えると，結果が入れ替わってしまいました。
たぶん，BDとDがとても似通った値だからだと思うのですが・・・。
BDとDが似通っているので，分散分析を行なう前にBDとDのどちらを使うか，AICなどを基準にして決めた方が良いのでしょうか。
私はある動物の行動を調べています。その動物は，ある場所から出発して別の場所へ向かい，ある場所に戻ってくるという行動を繰り返しています。
上 記のデータでは，Fが'出発する時の速さ'で，BDが戻ってくる前に向かった場所までの距離，Dが次に向かう場所までの距離です。いずれにせよ，遠ければ 遠いほど速く移動します。出発する速さが，戻ってくる前の場所までの距離なのか，これから向かう場所までの距離によるのかを知りたくて，この解析を行ない ました。
とりとめなくなってしまい，申し訳ありません。
よろしくお願い致します。
No.04772 Re: anova()を使う際の因子の順番 【波音】 2007/11/22(Thu) 11:00
> タイプ3の平方和の計算の仕方で，分散分析を行なったら良いのかな，と思っていました。ですがどのようにしたら，良いかよくわかりません・・・。
先ほど私が紹介したcarパッケージのAnova()の引数で，
> Anova(result, test="III")
とすればできますよ。
> 波音さんがやっていたのと同じ方法で行なった結果です。これは，タイプ1での平方和なんでしょうか？
は い，そうです。ちなみに，タイプ1とか2とか3とかのはSASでそう書かれている（採用している名称）のようで，一般的には逐次平方和とか調整平方和など と呼ばれるようです。たしか，この掲示板の過去ログにも同じような内容があったので，検索してみると良いかもしれません。
要するに，変数の順番を入れ換えるとsummary(aov())で計算される平方和が異なるということで，「これはどうしてなの？」ということについては私が下手に説明するよりも，
Alan Grafen, Rosie Hails（野間口謙太郎，野間口眞太郎訳） 一般線形モデルによる生物科学のための現代統計学 −あなたの実験をどのように解析するか− 共立出版
のp64-68に，まさにこれについての説明があるので，それを読んだ方がよいかと思います。
No.04785 Re: anova()を使う際の因子の順番 【きみか】 2007/11/22(Thu) 18:05
波音さん，お返事ありがとうございます！！！
すごくすっきりしました！
参考文献と過去ログを読んでみます！
本当にありがとうございました(^-^)> Anova(result,type='III')
Anova Table (Type III tests)
Response: F
LR Chisq Df Pr(>Chisq)
D 0.0036 1 0.95187
BD 6.3507 1 0.01173 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
> Anova(result2,type='III')
Anova Table (Type III tests)
Response: F
LR Chisq Df Pr(>Chisq)
BD 6.3507 1 0.01173 *
D 0.0036 1 0.95187
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
● 「統計学関連なんでもあり」の過去ログ--- 041 の目次へジャンプ
● 「統計学関連なんでもあり」の目次へジャンプ
● 直前のページへ戻る