Как вычислить средние значения из линейной модели в R?
Я работаю с набором данных о сохранении и его влиянии на биомассу, в котором пятьдесят участков земли, каждый по одному гектару, были отобраны случайным образом с площади в десять тысяч гектаров в Северной Англии.
Для каждого земельного участка были зарегистрированы следующие переменные:
* Биомасса: оценка биомассы растительности в кг на квадратный метр.• alt: средняя высота участка в метрах над уровнем моря.
• минусы: категориальная переменная, которая кодировался 1, Если участок был частью заповедной зоны, и 2 в противном случае.
* Почва категориальная переменная, грубо классифицирующая тип почвы как 1 для мела, 2 для глины и 3 для суглинка.В данный момент я борюсь с двумя вещами в частности:
Как вычислить среднюю разницу в биомассе между глинистыми (soil2) и суглинистыми (soil3) почвами на основе моей адаптированной модели (model1) и вычислить 95% доверительный интервал для этого среднего прогнозируемого значения.
И как это сделать вычислить среднюю прогнозируемую биомассу для участка, расположенного в заповедной зоне с преимущественно глинистой почвой на высоте 300 м?
Это краткое изложение линейной модели, с которой я работаю.Call:
lm(formula = biomass ~ alt + soil + cons, data = conservation)
Residuals:
Min 1Q Median 3Q Max
-0.183105 -0.052926 0.005593 0.061844 0.194402
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 2.2928629 0.0357850 64.073 < 2e-16 ***
alt -0.0029068 0.0001302 -22.318 < 2e-16 ***
soil2 -0.0862220 0.0342955 -2.514 0.0156 *
soil3 -0.2309939 0.0354480 -6.516 5.33e-08 ***
cons2 0.0488634 0.0292075 1.673 0.1013
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Residual standard error: 0.09428 on 45 degrees of freedom
Multiple R-squared: 0.9459, Adjusted R-squared: 0.9411
F-statistic: 196.7 on 4 and 45 DF, p-value: < 2.2e-16
И вот данные:
dput(conservation)
structure(list(biomass = c(2.01, 2.06, 1.7, 2.07, 1.88, 2.11,
0.98, 2.14, 1.75, 1.81, 2.15, 1.68, 2.23, 2.04, 1.67, 1.77, 1.74,
1.53, 1.79, 2.15, 1.39, 2.19, 2.14, 2.29, 1.91, 1.73, 2.21, 1.96,
2.07, 2.01, 2.2, 2.24, 1.33, 1.05, 1.36, 1.72, 1.44, 1.52, 2.09,
1.42, 1.64, 0.92, 1.65, 1.37, 0.77, 1.57, 2.25, 2.23, 2.03, 1.18
), alt = c(116L, 21L, 130L, 65L, 117L, 82L, 359L, 5L, 86L, 91L,
64L, 178L, 79L, 70L, 209L, 110L, 161L, 248L, 146L, 23L, 237L,
84L, 40L, 7L, 161L, 122L, 25L, 146L, 67L, 118L, 42L, 57L, 277L,
338L, 331L, 153L, 239L, 237L, 67L, 171L, 206L, 371L, 107L, 236L,
482L, 240L, 56L, 42L, 68L, 436L), cons = structure(c(2L, 2L,
1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 2L, 2L, 2L, 2L, 1L, 1L, 2L, 1L,
1L, 2L, 2L, 2L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 1L,
1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 2L, 2L, 2L, 2L
), .Label = c("1", "2"), class = "factor"), soil = structure(c(2L,
3L, 2L, 2L, 2L, 1L, 2L, 2L, 3L, 3L, 2L, 1L, 1L, 2L, 2L, 3L, 2L,
3L, 2L, 2L, 3L, 1L, 3L, 2L, 1L, 3L, 1L, 1L, 1L, 1L, 1L, 1L, 3L,
3L, 2L, 1L, 3L, 2L, 1L, 3L, 3L, 3L, 3L, 3L, 2L, 2L, 1L, 1L, 1L,
2L), .Label = c("1", "2", "3"), class = "factor"), alt.factor =
structure(c(1L,
1L, 2L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 2L, 1L, 1L, 2L, 1L, 2L,
2L, 2L, 1L, 2L, 1L, 1L, 1L, 2L, 2L, 1L, 2L, 1L, 2L, 1L, 1L, 2L,
2L, 2L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 1L, 2L, 2L, 2L, 1L, 1L, 1L,
2L), .Label = c("below median", "above median"), class = "factor")),
.Names = c("biomass",
"alt", "cons", "soil", "alt.factor"), row.names = c(NA, -50L), class =
"data.frame")
1 ответ:
Строго говоря, это частный случай того, что мы называем "проверкой линейной гипотезы". Но я думаю, что это не входит в намерения вашего назначения, поэтому я не пойду на такой подход. Если вас это интересует, прочитайте получить p-значение для группового среднего различия без переоснащение линейной модели новым эталонным уровнем . То, что я сделаю здесь, - это просто использовать другой уровень фактора в качестве контраста и переоборудовать вашу модель. В данный момент у вас есть" soil1 "в качестве уровня контрастности; я сброшу" soil2 " в качестве уровня контрастности. Посмотрите на , Как установить контрасты для моей переменной в регрессионном анализе с R? для общего лечения.Как вычислить среднюю разницу в биомассе между глинистыми (soil2) и суглинистыми (soil3) почвами на основе моей адаптированной модели (model1) и вычислить 95% доверительный интервал для этого среднего прогнозируемого значения.
fit <- lm(biomass ~ alt + soil + cons, data = conservation, contrasts = list(soil = contr.treatment(n = 3, base = 2))) #Coefficients: # Estimate Std. Error t value Pr(>|t|) #(Intercept) 2.2066409 0.0400572 55.087 < 2e-16 *** #alt -0.0029068 0.0001302 -22.318 < 2e-16 *** #soil1 0.0862220 0.0342955 2.514 0.0156 * #soil3 -0.1447719 0.0325295 -4.450 5.59e-05 *** #cons2 0.0488634 0.0292075 1.673 0.1013 #--- #Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 # #Residual standard error: 0.09428 on 45 degrees of freedom #Multiple R-squared: 0.9459, Adjusted R-squared: 0.9411 #F-statistic: 196.7 on 4 and 45 DF, p-value: < 2.2e-16
Теперь коэффициент для "soil3" дает разницу среднего значения группы для "soil3" и контраста уровень "soil2". Довольно просто получить доверительный интервал для этого коэффициента, исходя из стандартной ошибки и остаточной степени свободы модели, но опять же, это может быть слишком технически для вас. Рассмотрим использование
confint
:confint(fit, "soil3", level = 0.95) # 2.5 % 97.5 % #soil3 -0.2102896 -0.0792541
А как вычислить среднюю прогнозируемую биомассу для участка, расположенного в заповедной зоне с преимущественно глинистой почвой на высоте 300 м?Для предсказания ответа
biomass
, мы можем использоватьpredict
:predict(fit, newdata = list(alt = 300, soil = "2", cons = "1")) # 1 #1.334606
Так что среднее значение прогноза-около
1.3346
.