Тензорное произведение ti () в пакете GAM дает неверные результаты - PullRequest
Новая
       18

Тензорное произведение ti () в пакете GAM дает неверные результаты

1 голос
/ 24 октября 2019

Я с удивлением замечаю, что как-то трудно получить правильную подгонку функции взаимодействия от gam().

Чтобы быть более точным, я хочу оценить аддитивную функцию:

y=m_1(x)+m_2(z)+m_{12}(x,z)+u,

, где m_1(x)=x^2, m_2(z)=z^2, m_{12}(x,z)=xz. Следующий код генерирует эту модель: 

test1 <- function(x,z,sx=1,sz=1) { 

  #--m1(x) function
  m.x<-x^2
  m.x<-m.x-mean(m.x)

  #--m2(z) function
  m.z<-z^2
  m.z<-m.z-mean(m.z)

  #--m12(x,z) function
  m.xz<-x*z
  m.xz<-m.xz-mean(m.xz)
  m<-m.x+m.z+m.xz

  return(list(m=m,m.x=m.x,m.z=m.z,m.xz=m.xz))
}

n <- 1000
a=0
b=2

x <- runif(n,a,b)/20
z <- runif(n,a,b)
u <- rnorm(n,0,0.5)
model<-test1(x,z)

y <- model$m + u

Поэтому я использую gam(), подгоняя модель как 

b3 <- gam(y~ ti(x) + ti(z) + ti(x,z))
vis.gam(b3);title("tensor anova")

#---extracting basis matrix
B.f3<-model.matrix.gam(b3)

#---extracting series estimator
b3.hat<-b3$coefficients

Вопрос: при построении оценочной функциина gam() выше относительно ее истинной функции, я получаю 

par(mfrow=c(1,3))

#---m1(x)
B.x<-B.f3[,c(2:5)]
b.x.hat<-b3.hat[c(2:5)]
plot(x,B.x%*%b.x.hat)
points(x,model$m.x,col='red')
legend('topleft',c('Estimate','True'),lty=c(1,1),col=c('black','red'))

#---m2(z)
B.z<-B.f3[,c(6:9)]
b.z.hat<-b3.hat[c(6:9)]
plot(z,B.z%*%b.z.hat)
points(z,model$m.z,col='red')
legend('topleft',c('Estimate','True'),lty=c(1,1),col=c('black','red'))

#---m12(x,z)
B.xz<-B.f3[,-c(1:9)]
b.xz.hat<-b3.hat[-c(1:9)]
plot(x,B.xz%*%b.xz.hat)
points(x,model$m.xz,col='red')
legend('topleft',c('Estimate','True'),lty=c(1,1),col=c('black','red'))

Однако оценка функции m_1(x) в значительной степени отличается от x^2, и оценка функции взаимодействия m_{12}(x,z) такжев значительной степени отличается от xz, определенного в test1 выше. Результаты одинаковы, если я использую predict(b3).

Я действительно не могу понять это. Кто-нибудь может мне помочь, объяснив, почему результаты заканчиваются этим? Очень ценю это!

1 Ответ

0 голосов
/ 29 октября 2019

Во-первых, проблема вышеуказанного вопроса, конечно, не из-за пакета. Он тесно связан с условиями идентификации гладких функций. Одной из распространенных практик является навязывание предположений, что E(mj(.))=0 для всех отдельных функций j=1,...,d, и E(m_ij(x_i,x_j)|x_i)=E(m_ij(x_i,x_j)|x_j)=0 для i, не равных j. Эти условия требуют использования центрированной базисной функции в оценщике серий, что уже было сделано в пакете GAM. Однако, в моем случае выше, функция m(x,z)=x*z, определенная в test1, не удовлетворяет вышеприведенным предположениям идентификации, поскольку интеграл x*z относительно x или z не равен нулю, когда x и z имеютдиапазон от нуля до двух.

Кроме того, последовательный оценщик позволяет идентифицировать индивида и функцию взаимодействия, если кто-то навязывает m(0)=0 или m(0,x_j)=m(x_i,0)=0. Это может быть легко достигнуто, если мы сосредоточим базисную функцию вокруг нуля. Я пробовал оба случая, и они хорошо работают, когда DGP удовлетворяет условиям идентификации.

...