按組擬合線性模型/方差分析

Question

我試圖在R中運行anova()並遇到一些困難。這是我迄今爲止幫助解決我的問題所做的工作。

這裏是我的數據到這一點的str()。

str(mhw) 
'data.frame': 500 obs. of 5 variables: 
$ r : int 1 2 3 4 5 6 7 8 9 10 ... 
$ c : int 1 1 1 1 1 1 1 1 1 1 ... 
$ grain: num 3.63 4.07 4.51 3.9 3.63 3.16 3.18 3.42 3.97 3.4 ... 
$ straw: num 6.37 6.24 7.05 6.91 5.93 5.59 5.32 5.52 6.03 5.66 ... 
$ Quad : Factor w/ 4 levels "NE","NW","SE",..: 2 2 2 2 2 2 2 2 2 2 ... 

列r是表示一個數值，其行中的字段中的各個描繪駐留 c欄是表示數值的列的單個情節駐留
柱四對應於地理位置在現場其中每個小區駐留

Quad <- ifelse(mhw$c > 13 & mhw$r < 11, "NE",ifelse(mhw$c < 13 & mhw$r < 11,"NW", ifelse(mhw$c < 13 & mhw$r >= 11, "SW","SE"))) 
mhw <- cbind(mhw, Quad) 

我有適合lm()如下

nov.model <-lm(mhw$grain ~ mhw$straw) 
anova(nov.model) 

這是整個田間的anova()，它測試數據集中每個地塊的穀物產量與秸稈產量。

我的麻煩是，我想運行一個單獨的anova()我的數據的四列來測試每個象限的糧食產量和秸稈產量。

也許一個with()可能會解決這個問題。我以前從未使用過，現在正在學習R的過程中。任何幫助將不勝感激。

Answer 1

我認爲您R.

fit <- with(mhw, by(mhw, Quad, function (dat) lm(grain ~ straw, data = dat))) 

尋找by設施，因爲你在Quad有4個級別，最終在fit 4個線性模型，即fit是「由」類對象長度的（一種類型的「列表」）4.

要獲得係數爲每個模型，您可以使用

sapply(fit, coef) 

要機生產線CE模型綜上所述，使用

lapply(fit, summary) 

要導出方差分析表，使用

lapply(fit, anova) 

---

由於重複的例子，我拿着從?by的例子：

tmp <- with(warpbreaks, 
      by(warpbreaks, tension, 
       function(x) lm(breaks ~ wool, data = x))) 

class(tmp) 
# [1] "by" 

mode(tmp) 
# [1] "list" 

sapply(tmp, coef) 

#     L   M   H 
#(Intercept) 44.55556 24.000000 24.555556 
#woolB  -16.33333 4.777778 -5.777778 

lapply(tmp, anova) 

#$L 
#Analysis of Variance Table 
# 
#Response: breaks 
#   Df Sum Sq Mean Sq F value Pr(>F) 
#wool  1 1200.5 1200.50 5.6531 0.03023 * 
#Residuals 16 3397.8 212.36     
#--- 
#Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 
# 
#$M 
#Analysis of Variance Table 
# 
#Response: breaks 
#   Df Sum Sq Mean Sq F value Pr(>F) 
#wool  1 102.72 102.722 1.2531 0.2795 
#Residuals 16 1311.56 81.972    
# 
#$H 
#Analysis of Variance Table 
# 
#Response: breaks 
#   Df Sum Sq Mean Sq F value Pr(>F) 
#wool  1 150.22 150.222 2.3205 0.1472 
#Residuals 16 1035.78 64.736 

---

我知道這個選項，但不是fa熟悉它。由於@Roland針對上述重複的例子，提供代碼：

library(nlme) 
lapply(lmList(breaks ~ wool | tension, data = warpbreaks), anova) 

爲您的數據，我認爲這將是

fit <- lmList(grain ~ straw | Quad, data = mhw) 
lapply(fit, anova) 

你不需要安裝nlme;它附帶R作爲推薦軟件包之一。