Опция весов в GAM

Question

В моем наборе данных есть много избыточных наблюдений (но каждое наблюдение должно учитываться). Поэтому я считаю, что в GAM нужно использовать опцию "веса", потому что она значительно сокращает время вычислений.

gam функция (в mgcv пакет) объясняет, что они "эквивалентны" (из ?gam на аргументы weights):

"Обратите внимание, что вес 2, например, эквивалентен тому, чтобы сделать одно и то же наблюдение дважды".

Но это не кажется правильным.

yy = c(5,2,8,9)
xx = 1:4
wgts = c(3,2,4,1)
yy2 = rep(yy, wgts)
xx2 = rep(xx, wgts)
mod1 = gam(yy2 ~ xx2)
mod2 = gam(yy ~ xx, weights = wgts)
mod3 = gam(yy ~ xx, weights = wgts / mean(wgts))

predict(mod1,data.frame(xx2=1:4))
predict(mod2,data.frame(xx=1:4))
predict(mod3,data.frame(xx=1:4))

Оценки идентичны во всех трех моделях. Стандартная ошибка одинакова в моделях 2 и 3, но различна в модели 1. GCV отличается во всех трех моделях.

Я понимаю, что GCV могут быть разными. Но как мы можем сказать, что модели идентичны, если стандартные ошибки различны? Это ошибка или есть хорошее объяснение этому?

Answer 1

Проблемы, которые вы видели, не о GAM. Вы использовали gam чтобы соответствовать параметрической модели, в этом случае gam ведет себя почти так же, как lm, Чтобы ответить на ваши вопросы, достаточно сосредоточиться на случае линейной регрессии. То, что происходит с линейной моделью, будет происходить и с GLM, и с GAM. Вот как мы можем воспроизвести проблему с lm:

yy <- c(5,2,8,9)
xx <- 1:4
wgts <- c(3,2,4,1)
yy2 <- rep(yy,wgts)
xx2 <- rep(xx,wgts)
fit1 <- lm(yy2 ~ xx2)
fit2 <- lm(yy ~ xx, weights = wgts)
fit3 <- lm(yy ~ xx, weights = wgts/mean(wgts))
summary1 <- summary(fit1)
summary2 <- summary(fit2)
summary3 <- summary(fit3)
pred1 <- predict(fit1, list(xx2 = xx), interval = "confidence", se.fit = TRUE)
pred2 <- predict(fit2, list(xx = xx), interval = "confidence", se.fit = TRUE)
pred3 <- predict(fit3, list(xx = xx), interval = "confidence", se.fit = TRUE)

Все модели имеют одинаковые коэффициенты регрессии, но другие результаты могут отличаться. Вы спрашивали:

1. Для взвешенной регрессии fit2 а также fit3 почему почти все то же самое, кроме остаточной стандартной ошибки?
2. Почему взвешенная регрессия (fit2 или же fit3) не эквивалентно обычной регрессии со связями?

---

Ваш первый вопрос о масштабной инвариантности весовых наименьших квадратов к весам. Вот краткое резюме, которое я сделал:

Если мы изменим масштаб W на произвольное положительное значение, изменится только остаточная стандартная ошибка и немасштабированная ковариация. Такое изменение не подразумевает другую, неэквивалентную модель. На самом деле все, что связано с предсказанием, не затронуто. При взвешенной регрессии не просто смотрите на sigma2; это просто предельная дисперсия. Что действительно представляет интерес, так это грубая дисперсия после умножения весов. Если вы разделите свой вес на 2, вы найдете sigma2 удваивается, но вы все равно получаете тот же результат при умножении их вместе.

summary2$coef
summary3$coef

#            Estimate Std. Error   t value  Pr(>|t|)
#(Intercept) 2.128713   3.128697 0.6803832 0.5664609
#xx          1.683168   1.246503 1.3503125 0.3094222

pred2
pred3

#$fit
#       fit        lwr      upr
#1 3.811881 -5.0008685 12.62463
#2 5.495050 -0.1299942 11.12009
#3 7.178218  0.6095820 13.74685
#4 8.861386 -1.7302209 19.45299
#
#$se.fit
#       1        2        3        4 
#2.048213 1.307343 1.526648 2.461646 
#
#$df
#[1] 2
#
#$residual.scale  ## for `pred2`
#[1] 3.961448
#
#$residual.scale  ## for `pred3`
#[1] 2.50544

---

Ваш второй вопрос о значении весов. Весовые коэффициенты используются для моделирования гетероскедастического ответа для преодоления эффекта рычагов в обычной регрессии наименьших квадратов. Весовые коэффициенты пропорциональны обратной дисперсии: вы даете больший вес данным с меньшими ожидаемыми ошибками. Веса могут быть нецелыми, поэтому у них нет естественного объяснения с точки зрения повторяющихся данных. Итак, что написано в mgcv пакет не является строго правильным.

Реальная разница между fit1 а также fit2 ? это степень свободы. Проверьте приведенную выше таблицу на (n - p), n количество данных у вас есть, в то время как p это число не NA коэффициенты, так n - p это остаточная степень свободы. Для обеих моделей у нас есть p = 2 (перехват и наклон), но для fit1 у нас есть n = 10 в то время как для fit2 у нас есть n = 4, Это оказывает существенное влияние на вывод, так как теперь стандартные ошибки для коэффициентов и прогнозов (следовательно, доверительные интервалы) будут отличаться. Эти две модели далеко не эквивалентны.

summary1$coef
#            Estimate Std. Error  t value   Pr(>|t|)
#(Intercept) 2.128713  1.5643486 1.360766 0.21068210
#xx2         1.683168  0.6232514 2.700625 0.02704784

summary2$coef

#            Estimate Std. Error   t value  Pr(>|t|)
#(Intercept) 2.128713   3.128697 0.6803832 0.5664609
#xx          1.683168   1.246503 1.3503125 0.3094222

pred1

#$fit
#       fit      lwr       upr
#1 3.811881 1.450287  6.173475
#2 5.495050 3.987680  7.002419
#3 7.178218 5.417990  8.938446
#4 8.861386 6.023103 11.699669
#
#$se.fit
#        1         2         3         4 
#1.0241066 0.6536716 0.7633240 1.2308229 
#
#$df    # note, this is `10 - 2 = 8`
#[1] 8
#
#$residual.scale
#[1] 1.980724

pred2

#$fit
#       fit        lwr      upr
#1 3.811881 -5.0008685 12.62463
#2 5.495050 -0.1299942 11.12009
#3 7.178218  0.6095820 13.74685
#4 8.861386 -1.7302209 19.45299
#
#$se.fit
#       1        2        3        4 
#2.048213 1.307343 1.526648 2.461646 
#
#$df    # note, this is `4 - 2 = 2`
#[1] 2
#
#$residual.scale  ## for `pred2`
#[1] 3.961448