2016-07-21 8 views
2

Моей текущей проблемой заключается в вычислении дисперсии объясняется различными переменными общие аддитивная модель (GAM) с R.GAM R дисперсия объясняется переменной

Я последовал за разъяснение Вуд здесь: https://stat.ethz.ch/pipermail/r-help/2007-October/142743.html

Но я хотел бы сделать это с тремя переменными. Я попытался это:

library(mgcv) 

set.seed(0) 
n<-400 
x1 <- runif(n, 0, 1) 
x2 <- runif(n, 0, 1) 
x3 <- runif(n, 0, 1) 

f1 <- function(x) exp(2 * x) - 3.75887 
f2 <- function(x) 0.2*x^11*(10*(1-x))^6+10*(10*x)^3*(1-x)^10 
f3 <- function(x) 0.008*x^2 - 1.8*x + 874 
f <- f1(x1) + f2(x2) + f3(x3) 
e <- rnorm(n, 0, 2) 
y <- f + e 

b <- gam(y ~ s(x1, k = 3)+s(x2, k = 3)+ s(x3, k = 3)) 
b3 <- gam(y ~ s(x1) + s(x2), sp = c(b$sp[1], b$sp[2])) 
b2 <- gam(y ~ s(x1) + s(x3), sp = c(b$sp[1], b$sp[3])) 
b1 <- gam(y ~ s(x2) + s(x3), sp = c(b$sp[2], b$sp[3])) 

b0 <- gam(y~1) 

(deviance(b1)-deviance(b))/deviance(b0) 
(deviance(b2)-deviance(b))/deviance(b0) 
(deviance(b3)-deviance(b))/deviance(b0) 

Но я не понимаю результаты. Например, модель с только x1 и x2 имеет отклонение, меньшее отклонения от трех объясняющих переменных.

Является ли метод, который я использовал для извлечения дисперсии, объясняемого переменной тремя переменными, является правильным?

Означает ли это, что в глобальной модели существует смешанный эффект? Или есть другое объяснение?

Большое спасибо.

ответ

4

Вы сделали что-то здесь не так:

b <- gam(y ~ s(x1, k = 3) + s(x2, k = 3) + s(x3, k = 3)) 
b3 <- gam(y ~ s(x1) + s(x2), sp = c(b$sp[1], b$sp[2])) 
b2 <- gam(y ~ s(x1) + s(x3), sp = c(b$sp[1], b$sp[3])) 
b1 <- gam(y ~ s(x2) + s(x3), sp = c(b$sp[2], b$sp[3])) 

Почему вы задали k = 3 в первой строке, а не установка k = 3 для остальные? Без указания k, s() примет значение по умолчанию k = 10. Теперь у вас возникли проблемы: b1, b2, b3 не вложены в b.

В оригинальном примере Саймона Вуда он оставил k неуказанным, так что k=10 берется за всех s(). Фактически, вы можете варьировать значения k, но вы должны гарантировать, что у вас всегда есть тот же k для того же ковариата (для обеспечения вложенности). Например, вы можете сделать:

b <- gam(y ~ s(x1, k = 4) + s(x2, k = 6) + s(x3, k = 3)) 
b3 <- gam(y ~ s(x1, k = 4) + s(x2, k = 6), sp = c(b$sp[1], b$sp[2])) ## droping s(x3) from b 
b2 <- gam(y ~ s(x1, k = 4) + s(x3, k = 3), sp = c(b$sp[1], b$sp[3])) ## droping s(x2) from b 
b1 <- gam(y ~ s(x2, k = 6) + s(x3, k = 3), sp = c(b$sp[2], b$sp[3])) ## droping s(x1) from b 

Тогда давайте сделаем:

(deviance(b1)-deviance(b))/deviance(b0) 
# [1] 0.2073421 
(deviance(b2)-deviance(b))/deviance(b0) 
# [1] 0.4323154 
(deviance(b3)-deviance(b))/deviance(b0) 
# [1] 0.02094997 

Положительные значения означают, что сбросив любой член модели будет раздувать девиантности, который разумно, как наша истинная модель есть все три условия ,

+0

Спасибо за ваш ответ, –

+0

Спасибо за ваш ответ, я согласен с вами. у него все еще есть что-то, чего я не понимаю. Почему мне нужно указать k (максимальная степень сглаживания), если оптимизация сглаживания отсутствует? С моей точки зрения, указание «sp = ...» означает, что нет необходимости в оптимизации, и поэтому не нужно устанавливать k. –