Bagaimana cara mengatur dan menafsirkan kontras ANOVA dengan paket mobil di R?

15

Katakanlah saya memiliki percobaan faktorial 2x2 sederhana yang ingin saya lakukan ANOVA. Seperti ini, misalnya:

d   <- data.frame(a=factor(sample(c('a1','a2'), 100, rep=T)),
                  b=factor(sample(c('b1','b2'), 100, rep=T)));
d$y <- as.numeric(d$a)*rnorm(100, mean=.75, sd=1) +
       as.numeric(d$b)*rnorm(100, mean=1.2, sd=1) +
       as.numeric(d$a)*as.numeric(d$b)*rnorm(100, mean=.5, sd=1) +
       rnorm(100);
  1. Dengan tidak adanya interaksi yang signifikan, secara default (yaitu contr.treatment) output dari Anova()adalah signifikansi keseluruhan dari asemua tingkatan bdan dari bsemua tingkatan a, apakah itu benar?

  2. Bagaimana saya harus menentukan kontras yang akan memungkinkan saya untuk menguji signifikansi efek adengan bdijaga konstan pada level b1, efek adengan bdijaga konstan pada level b2, dan interaksi a:b?

f1r3br4nd
sumber

Jawaban:

18

Contoh Anda mengarah ke ukuran sel yang tidak sama, yang berarti bahwa perbedaan "jenis jumlah kuadrat", dan tes untuk efek utama tidak sesederhana seperti yang Anda nyatakan. Anova()menggunakan jumlah kuadrat tipe II. Lihat pertanyaan ini sebagai permulaan.

Ada berbagai cara untuk menguji kontras. Perhatikan bahwa tipe SS tidak penting karena kami akhirnya menguji dalam desain satu faktor yang terkait. Saya sarankan menggunakan langkah-langkah berikut:

# turn your 2x2 design into the corresponding 4x1 design using interaction()
> d$ab <- interaction(d$a, d$b)       # creates new factor coding the 2*2 conditions
> levels(d$ab)                        # this is the order of the 4 conditions
[1] "a1.b1" "a2.b1" "a1.b2" "a2.b2"

> aovRes <- aov(y ~ ab, data=d)       # oneway ANOVA using aov() with new factor

# specify the contrasts you want to test as a matrix (see above for order of cells)
> cntrMat <- rbind("contr 01"=c(1, -1,  0,  0),  # coefficients for testing a within b1
+                  "contr 02"=c(0,  0,  1, -1),  # coefficients for testing a within b2
+                  "contr 03"=c(1, -1, -1,  1))  # coefficients for interaction

# test contrasts without adjusting alpha, two-sided hypotheses
> library(multcomp)                   # for glht()
> summary(glht(aovRes, linfct=mcp(ab=cntrMat), alternative="two.sided"),
+         test=adjusted("none"))
Simultaneous Tests for General Linear Hypotheses
Multiple Comparisons of Means: User-defined Contrasts
Fit: aov(formula = y ~ ab, data = d)

Linear Hypotheses:
              Estimate Std. Error t value Pr(>|t|)
contr 01 == 0  -0.7704     0.7875  -0.978    0.330
contr 02 == 0  -1.0463     0.9067  -1.154    0.251
contr 03 == 0   0.2759     1.2009   0.230    0.819
(Adjusted p values reported -- none method)    

Sekarang secara manual periksa hasilnya untuk kontras pertama.

> P       <- 2                             # number of levels factor a
> Q       <- 2                             # number of levels factor b
> Njk     <- table(d$ab)                   # cell sizes
> Mjk     <- tapply(d$y, d$ab, mean)       # cell means
> dfSSE   <- sum(Njk) - P*Q                # degrees of freedom error SS
> SSE     <- sum((d$y - ave(d$y, d$ab, FUN=mean))^2)    # error SS
> MSE     <- SSE / dfSSE                   # mean error SS
> (psiHat <- sum(cntrMat[1, ] * Mjk))      # contrast estimate
[1] -0.7703638

> lenSq <- sum(cntrMat[1, ]^2 / Njk)       # squared length of contrast
> (SE   <- sqrt(lenSq*MSE))                # standard error
[1] 0.7874602

> (tStat <- psiHat / SE)                   # t-statistic
[1] -0.9782893

> (pVal <- 2 * (1-pt(abs(tStat), dfSSE)))  # p-value
[1] 0.3303902
caracal
sumber
3
TERIMA KASIH!!! Anda baru saja menjawab pertanyaan yang belum dimiliki oleh dua semester statistik tingkat pascasarjana. Saya bahkan mempertimbangkan untuk menggunakan anova satu arah sebelumnya, tetapi tidak dapat menemukan konfirmasi bahwa ini adalah pendekatan yang sah.
f1r3br4nd
@ f1r3br4nd Ini sah karena kesalahan MS sama dalam desain satu arah dan asli dua arah.
caracal
Satu pertanyaan lanjutan terakhir, jika saya boleh: bagaimana interaksi dua arah menggeneralisasikan ke interaksi sejumlah besar variabel? Jika saya memiliki istilah A B C, akankah saya membangunnya dari A: B = (A | B = 1 - A | B = 2), C: B = (C | B = 1 - C | B = 2 ), A: B: C = A: B - C: B, dan seterusnya?
f1r3br4nd
2
@ f1r3br4nd Dalam desain 2x2x2, hanya ada satu kontras interaksi A B C yang unik (seperti hanya ada satu dalam wadah 2x2). Koefisien dalam kontras interaksi B C harus dijumlahkan ke nol di atas baris (A), kolom (B), dan pesawat (C) dalam "desain kubus". Jika urutan sel dalam desain satu arah yang terkait adalah a1.b1.c1, a2.b1.c1, a1.b2.c1, a2.b2.c1, a1.b1.c2, a2.b1.c2, a1.b2.c2, a2.b2.c2, maka koefisiennya adalah c(1, -1, -1, 1, -1, 1, 1, -1). Jika Anda memiliki lebih dari dua grup dalam faktor Anda, semua kontras yang mengikuti aturan jumlah-ke-nol adalah kontras interaksi 3-arah.
caracal