Тест отношения правдоподобия для кажущихся несвязанными регрессий - PullRequest
0 голосов
/ 06 февраля 2020

Я хотел бы запустить симуляцию по методу Монте-Карло в 2000 раз по статистике правдоподобия c для следующей системы уравнений:

y1 = a1 + b1X + c1K

y2 = a2 + b2X + c2K

Я хочу применить ограничение. Гипотеза, которую я хочу проверить: Ho: a1 + 39 = a2.

Предположим, что каждая выборка имеет 50 наблюдений.

Я использовал пакет systemfit для расчета моей имитированной модели SUR. Я просто не уверен, как применить ограничение в R.

Вот как я начал моделировать модель с ограничениями:

#Set Parameters for Monte Carlo Simulation
N <- 50
nrsim <- 2000
#Equation 1
a1 <- -29
b1 <- 0.12
g1 <- 0.12
mu1 <- 0
sig2e1 <- 1700
sde1 <- sqrt(sig2e1)
F1 <- unlist(mydata %>%
  filter(eqn == 1) %>%
  select(`F`))
C1 <- unlist(mydata %>%
  filter(eqn == 1) %>%
  select(C))

va1 <- numeric(nrsim) #stores the estimates of alpha1
vb1 <- numeric(nrsim) #stores the estimates of beta1
vg1 <- numeric(nrsim) #stores the estimates of gamma1

#Equation 2
a2 <- -12
b2 <- 0.078
g2 <- 0.9
mu2 <- 0
sig2e2 <- 1900
sde2 <- sqrt(sig2e2)
F2 <- unlist(mydata %>%
  filter(eqn == 2) %>%
  select(`F`))
C2 <- unlist(mydata %>%
  filter(eqn == 2) %>%
  select(C))   

#unrestricted model
unr.vsura1 <- numeric(nrsim) #stores the estimates of alpha1
unr.vsura2 <- numeric(nrsim) #stores the estimates of alpha2
unr.vsurb1 <- numeric(nrsim) #stores the estimates of beta1
unr.vsurb2 <- numeric(nrsim) #stores the estimates of beta2
unr.vsurg1 <- numeric(nrsim) #stores the estimates of gamma1
unr.vsurg2 <- numeric(nrsim) #stores the estimates of gamma2

for (j in 1:nrsim) {

  set.seed(12345+10*j)
  y1 <-  as.numeric(unlist(a1 + b1*F1 + g1*C1 + rnorm(N, mean = mu1, sd = sde1)))
  y2 <- as.numeric(unlist(a2 + b2*F2 + g2*C2 + rnorm(N, mean = mu2, sd = sde2)))
  X1 <- cbind(F1, C1)
  X2 <- cbind(F2, C2)
  eq1 <- y1~X1
  eq2 <- y2~X2
  system <- list(eq1 = eq1, eq2 = eq2)
  #model
  sur <- systemfit(system, method = "SUR", data = mydata)
  #store each iteration of the estimate
  unr.vsura1[j] <- coef(sur)[[1]]
  unr.vsura2[j] <- coef(sur)[[4]]
  unr.vsurb1[j] <- coef(sur)[[2]]
  unr.vsurb2[j] <- coef(sur)[[5]]
  unr.vsurg1[j] <- coef(sur)[[3]]
  unr.vsurg2[j] <- coef(sur)[[6]]
}
...