2

アプリオリなコントラストの設定に問題があり、助けを求めたいと思います。次のコードは、因子水準「d」に 2 つの直交対比を与える必要があります。

Response <- c(1,3,2,2,2,2,2,2,4,6,5,5,5,5,5,5,4,6,5,5,5,5,5,5)
A <- factor(c(rep("c",8),rep("d",8),rep("h",8)))
contrasts(A) <- cbind("d vs h"=c(0,1,-1),"d vs c"=c(-1,1,0))
summary.lm(aov(Response~A))

私が得るものは次のとおりです。

Call:
aov(formula = Response ~ A)

Residuals:
   Min         1Q     Median         3Q        Max 
-1.000e+00 -3.136e-16 -8.281e-18 -8.281e-18  1.000e+00 

Coefficients:
        Estimate Std. Error t value Pr(>|t|)    
(Intercept)   4.0000     0.1091  36.661  < 2e-16 ***
Ad vs h      -1.0000     0.1543  -6.481 2.02e-06 ***
Ad vs c       2.0000     0.1543  12.961 1.74e-11 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 0.5345 on 21 degrees of freedom
Multiple R-squared: 0.8889,     Adjusted R-squared: 0.8783 
F-statistic:    84 on 2 and 21 DF,  p-value: 9.56e-11

しかし、(Intercept) の推定値は 5.00 になると予想しています。これはレベル d に等しいはずですよね? また、他の見積もりは奇妙に見えます...

relevel(A, ref="d") (正しく表示される) を使用すると正しい値を簡単に取得できることはわかっていますが、自分の仮説をテストするための正しい定式化を学ぶことに興味があります。次のコード (Web サイトから) を使用して同様の例を実行すると、期待どおりに動作します。

irrigation<-factor(c(rep("Control",10),rep("Irrigated 10mm",10),rep("Irrigated20mm",10))) 
biomass<-1:30 
contrastmatrix<-cbind("10 vs 20"=c(0,1,-1),"c vs 10"=c(-1,1,0))
contrasts(irrigation)<-contrastmatrix 
summary.lm(aov(biomass~irrigation))


Call:
aov(formula = biomass ~ irrigation)

Residuals:
       Min         1Q     Median         3Q        Max 
-4.500e+00 -2.500e+00  3.608e-16  2.500e+00  4.500e+00 

Coefficients:
                   Estimate Std. Error t value Pr(>|t|)    
(Intercept)         15.5000     0.5528   28.04  < 2e-16 ***
irrigation10 vs 20 -10.0000     0.7817  -12.79 5.67e-13 ***
irrigationc vs 10   10.0000     0.7817   12.79 5.67e-13 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 3.028 on 27 degrees of freedom
Multiple R-squared: 0.8899,     Adjusted R-squared: 0.8817 
F-statistic: 109.1 on 2 and 27 DF,  p-value: 1.162e-13

これについての説明をいただければ幸いです。

ありがとう、ジェレミアス

4

1 に答える 1

6

問題は（詳細についてcontrastsはあなたかもしれません）の理解にあると思います。?contrasts詳しく説明しましょう：

のデフォルトの方法を使用する場合factor A、

A <- factor(c(rep("c",8),rep("d",8),rep("h",8)))
> contrasts(A)
  d h
c 0 0
d 1 0
h 0 1

したがって、モデルlmはあなたを与えます

Mean(Response) = Intercept + beta_1 * I(d = 1) + beta_2 * I(h = 1)

summary.lm(aov(Response~A))
Coefficients:
    Estimate Std. Error t value Pr(>|t|)    
(Intercept)    2.000      0.189    10.6  7.1e-10 ***
Ad             3.000      0.267    11.2  2.5e-10 ***
Ah             3.000      0.267    11.2  2.5e-10 ***

したがって、グループcの場合、平均は単に切片2であり、グループdの場合、平均は2 + 3 = 5グループと同じですh。

独自のコントラストを使用するとどうなりますか:

contrasts(A) <- cbind("d vs h"=c(0,1,-1),"d vs c"=c(-1,1,0))

A
[1] c c c c c c c c d d d d d d d d h h h h h h h h
attr(,"contrasts")
    d vs h d vs c
c      0     -1
d      1      1
h     -1      0

あなたがフィットするモデルは

Mean(Response) = Intercept + beta_1 * (I(d = 1) - I(h = 1)) + beta_2 * (I(d = 1) - I(c = 1))
     = Intercept + (beta_1 + beta_2) * I(d = 1) - beta_2 * I(c = 1) - beta_1 * I(h = 1)

Coefficients:
Estimate Std. Error t value Pr(>|t|)    
(Intercept)    4.000      0.109   36.66  < 2e-16 ***
Ad vs h       -1.000      0.154   -6.48  2.0e-06 ***
Ad vs c        2.000      0.154   12.96  1.7e-11 ***

したがって、グループcの場合、平均は4 - 2 = 2、グループdの場合、平均は4 - 1 + 2 = 5、グループhの場合、平均は4 - (-1) = 5です。

==========================

アップデート：

コントラストを行う最も簡単な方法は、ベース (参照) レベルを に設定することdです。

contrasts(A) <- contr.treatment(3, base = 2)
Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept)  5.00e+00   1.89e-01    26.5  < 2e-16 ***
A1          -3.00e+00   2.67e-01   -11.2  2.5e-10 ***
A3          -4.86e-17   2.67e-01     0.0        1

コントラストを使用する場合:

Response <- c(1,3,2,2,2,2,2,2,4,6,5,5,5,5,5,5,4,6,5,5,5,5,5,5)
A <- factor(c(rep("c",8),rep("d",8),rep("h",8)))
mat<- cbind(rep(1/3, 3), "d vs h"=c(0,1,-1),"d vs c"=c(-1,1,0))
mymat <- solve(t(mat))
my.contrast <- mymat[,2:3]
contrasts(A) <- my.contrast
summary.lm(aov(Response~A))

Coefficients:
Estimate Std. Error t value Pr(>|t|)    
(Intercept) 4.00e+00   1.09e-01    36.7  < 2e-16 ***
Ad vs h     7.69e-16   2.67e-01     0.0        1    
Ad vs c     3.00e+00   2.67e-01    11.2  2.5e-10 ***

参照: http://www.ats.ucla.edu/stat/r/library/contrast_coding.htm

于 2012-11-20T15:35:06.397 に答える