Регулярная регрессия по группам
Привет, у меня есть набор данных панели. Я хотел бы сделать регрессию скользящего окна для каждой фирмы и извлечь коэффициент независимой переменной. у - зависимая переменная, а х - независимая переменная. Скользящее окно равно 12. То есть, первая регрессия использует данные строки 1 для строки 12, вторая регрессия использует данные строки 2 для строки 13 и т. Д. Используется Rollapply.
Вот вопрос с той же самой ошибкой, с которой я столкнулся: Переход по группам в data.table R К счастью, этот вопрос занимает всего один столбец, а мой - два регресса, поэтому я не могу внести изменения в соответствии с рекомендуемым ответом в этом посте. Вот еще один пост, который использует цикл for. Мои реальные данные содержат более 2 миллионов наблюдений, поэтому они слишком медленные: скользящая регрессия с помощью dplyr Может ли кто-нибудь помочь?
Мой поддельный набор данных выглядит следующим образом:
dt<-rep(c("AAA","BBB","CCC"),each=24)
dt<-as.data.frame(dt)
names(dt)[names(dt)=="dt"] <- "firm"
a<-c(20100131,20100228,20100331,20100430,20100531,20100630,20100731,20100831,20100930,20101031,20101130,20101231,20110131,20110228,20110331,20110430,20110531,20110630,20110731,20110831,20110930,20111031,20111130,20111231)
dt$time<-rep(a,3)
dt<-dt%>% group_by(firm)%>%
mutate(y=rnorm(24,10,5))
dt<-dt%>% group_by(firm)%>%
mutate(x=rnorm(24,5,2))
dt<-as.data.table(dt)
Я попробовал этот код:
# create rolling regression function
roll <- function(Z)
{
t = lm(formula=y~x, data = as.data.frame(Z), na.rm=T);
return(t$coef[2])
}
dt[,beta := rollapply(dt, width=12, roll, fill=NA, by.column=FALSE, align="right") , by=firm]
Я пытаюсь создать столбец с именем "бета", который показывает коэффициент var x. Таким образом, для каждой фирмы первые данные должны быть получены из 12-го наблюдения.
Похоже, что регрессия берет x и y из 1-й строки для разных групп, а коэффициенты, кажется, немного смещены по сравнению с результатом, который я получил от EXCEL.
Второй метод, который я попробовал, это версия dplyr:
dt %>%
group_by(firm) %>%
mutate(dt,beta = rollapply(dt,12,function(x) coef(lm(y~x,data=as.data.frame(x)))[2],by.column= FALSE, fill = NA, align = "right"))
Это дает мне ту же проблему. каждая группа имеет одинаковый номер. Похоже, для каждой фирмы регрессия берет y и x из 1-го ряда.
Какие-нибудь мысли? Огромное спасибо.
2 ответа
Вот решение, которое использует rollRegres
пакет и data.table
пакет. Я также добавил модифицированную версию решения ОП, которая работает (см. Комментарий Эдди), и использовал пример с 2 миллионами наблюдений в качестве упоминания ОП.
#####
# setup data
library(rollRegres)
library(data.table)
library(dplyr)
set.seed(33700919)
n_firms <- 83334 # yields ~ the 2M firm as the OP mentions
dt <- rep(1:n_firms, each = 24)
dt <- data.frame(firm = dt)
a <-c(20100131,20100228,20100331,20100430,20100531,20100630,20100731,20100831,20100930,20101031,20101130,20101231,20110131,20110228,20110331,20110430,20110531,20110630,20110731,20110831,20110930,20111031,20111130,20111231)
dt$time <- rep(a, n_firms)
dt <- dt %>% group_by(firm) %>% mutate(y=rnorm(24,10,5))
dt <- dt %>% group_by(firm) %>% mutate(x=rnorm(24,5,2))
dt <- as.data.table(dt)
nrow(dt) # roughly the 2M rows that the OP mentions
#R [1] 2000016
#####
# fit models
setkey(dt, firm, time) # make sure data is sorted correctly
start_time <- Sys.time() # to show computation time
dt[
, beta :=
roll_regres.fit(x = cbind(1, .SD[["x"]]), y = .SD[["y"]],
width = 12L)$coefs[, 2],
by = firm]
Sys.time() - start_time
#R Time difference of 6.526595 secs
# gives the same as OP's solution with minor corrections
library(zoo)
start_time <- Sys.time()
roll <- function(Z)
lm.fit(x = cbind(1, Z[, "x"]), y = Z[, "y"])$coef[2]
dt[
, beta_zoo :=
rollapply(.SD, width=12, roll, fill=NA, by.column=FALSE, align="right"),
by=firm]
Sys.time() - start_time # much slower
#R Time difference of 1.87341 mins
# gives the same
all.equal(dt$beta, dt$beta_zoo)
#R [1] TRUE
Может быть, вы можете попробовать изменить первый аргумент в rollapply, заменить dt
в столбец, dt[, c("y","x")]
, Посмотрите, работает ли это