Полосы уверенности регрессии из моделей plm

В настоящее время я использую модель фиксированных эффектов с использованием функции plm из одноименного пакета. У меня есть кластеризация данных, поэтому я использую функцию vcovCR из пакета clubSandwich, чтобы получить стандартные ошибки оцененных коэффициентов регрессии.

Я хочу сделать следующее: создать график, который показывает прогнозируемое значение Y при заданном X1 (когда другие переменные сохраняются постоянными при их средних значениях), а также добавить к графику диапазоны уверенности 95%.

Что я пробовал:

  1. Использование функции ggpredict() из пакета ggeffects. Это не сработало, потому что пакет plm, очевидно, был проблематичным для ggeffects.
  2. Вручную с использованием подхода оценки фиксированных эффектов «центрирования группы» с последующим запуском обычной lm() функции для данных, центрированных по группе. Запуск ggpredict() на этих данных дает мне хорошие полосы уверенности, но значения X и Y затем центрируются по группе (сюрприз, сюрприз!), И я хочу, чтобы прогнозы и значения X были в исходной шкале, а не в центре среднего. шкала.

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


person Phil    schedule 20.03.2020    source источник


Ответы (1)


Я пробовал нечто подобное, но не смог найти прямого решения. Я оцениваю исследование событий с запаздыванием и опережением на модели fe. Для меня сработало следующее:

ES1<-plm(total_thefts ~ lag6 + lag5 + lag4 + lag3 + lag2 + E + lead1 + lead2 + lead3 + lead4 + lead5 + lead6 + factor(month) + factor(year), data = panel, model="within")

Поскольку в plm нет прямого способа оценки надежных стандартных ошибок, я делаю это вручную и заменяю исходные:

G <- length(unique(panel$id))
c <- G/(G - 1)
robust_se <- coeftest(ES1, c * vcovHC(ES1, type = 'HC1', cluster = 'group'))
ES1<-summary(ES1)
ES1$coefficients[,2:4] <- robust_se[,2:4]

После этого я разделяю интересующие меня коэффициенты и оцениваю диапазоны CI (95%), следуя этой сообщение:

coefs <- tidy(ES1[["coefficients"]])
coefs$conf.low <- coefs$Estimate+c(-1)*coefs$Std..Error*qt(0.975,42)
coefs$conf.high <- coefs$Estimate+c(1)*coefs$Std..Error*qt(0.975,42)

Наконец, после этого я сохраняю интересующие меня параметры и рисую исследование события:

interest <- c("lag6","lag5","lag4","lag3","lag2","E","lead1","lead2","lead3","lead4","lead5","lead6")
coefs <- subset(coefs,coefs$.rownames%in%interest)
coefs$time <- c(seq(-6,-2,1),seq(0,6,1))

ggplot(coefs, aes(time, Estimate))+
           geom_line() +
           geom_point()+
           geom_pointrange(aes(ymin = conf.low, ymax = conf.high))+
           geom_hline(yintercept = 0, linetype = 2) +
           ggtitle(paste0("N: ",nobs(ES1),", Sample: ",pdim(ES1)$nT$n,", Buff: 1000m"))) +
           theme(plot.title = element_text(hjust = 0.5)))

введите здесь описание изображения

person Pedro Esteban Rodriguez    schedule 25.03.2020