寫
在前面
“
對于任意一天t，在[t - n, t]的區間内進行回歸。如果資料一共有N天，那麼就會得到N - n個資料點
這就是滾動回歸，一個非常容易了解而且在研究中常常遇見，然而實作起來卻不是那麼容易的問題。在今天的大貓課堂中，大貓教大家用10行代碼搞定它！
PS：由于微信的限制，給大貓留言的小夥伴超過48小時後大貓就不能回複你們了。是以如果想聯系大貓，可以按照文章最後的微信号加大貓微信哦。
問
題引入
假設我們現在有N天的資料，我們希望對于每一天t，用n作為視窗期，在[t - n, t]的視窗中進行資料回歸。顯然，最終的資料會有N-n天。
為此我們構造樣例資料集，我們假設一共有5個不同的組，每組都由1000天的觀測，包含自變量x與因變量y。此外，我們假設視窗期 n = 100 天。
圖：樣例資料集
構造樣例資料集代碼如下：
# 設定随機數種子
set.seed(42)
# 生成樣例資料集，一共有a,b,c,d,e五個group，每個group都有1000日的觀測
dt <- data.table(id = rep(letters[1:5], each = 1000), date = seq(as.Date("2001-01-01"), by = 'day', len = 1000 * 5), y = rnorm(1000 * 5), x = rnorm(1000 * 5)) %>% unique(by = c("id", "date"))
要實作一行代碼完成分組回歸，需要用到data.table包！
解
題思路
解決的思路并不複雜，假設我們現在要處理的是第t行，自變量和因變量分别是x和y，滾動視窗是n天，那麼我們隻要能夠取到x[t-n, t]以及y[t-n,t]兩個向量，把他放到lm函數中就可以進行回歸得到結果。需要注意的是我們需要周遊每一個符合條件的t，以及需要把最終結果輸出成一個漂亮的資料集。
步
驟分解
大貓先把代碼放上來：
# 設定滾動視窗期，這裡為100天
n <- 100
# 計算滾動回歸！
re <- dt[, {
l <- list()
for (t in (n + 1) : .N) {
l[[t]] <- as.list(c(coef(lm(y ~ x, data = .SD[(t - n) : t])), date = date[t]))
}
rbindlist(l)
},
keyby = .(id)]
re[, ":="(date = as.Date(date, origin = "1970-01-01"))]
最終的輸出資料集是這個樣子的：
現在我們逐一分析這幾行代碼。
- 此處每個id有n = 1000天觀測，由于視窗期為100天，是以最終每個id會有1000 - 100 = 900個回歸結果
- keyby語句将原資料集按照id進行分組，具體作用可以看上期的《一行代碼搞定分組回歸》
- l <- list()語句：對于每個id的值，建立一個空白list，用于存貯回歸結果。
- for (t in (n + 1) : .N)語句：每個不同的id會有.N天觀測（這裡是1000），我們需要從 第(n + 1)天（在本例中是101）開始一直循環到第.N天。對于其中的每一天t，都進行一次回歸。
- data = .SD[(t - n) : t]語句：對于每一天t，都截取原資料集的(t-n)至t天作為lm函數的輸入資料集。關于.SD的具體使用可以見上期《一行代碼搞定分組回歸》
- rbindlist()語句：上面對于每一天t我們都生成了一個回歸，rbindlist語句将這些回歸結果打包起來輸出。
總
結
是不是很簡單？其實要完成滾動回歸并不止這一種方法，stackoverflow上有很多相關的文章，但是大貓在比較幾種方法之後，發現自己寫的這個版本是代碼最短、最容易了解、并且效率最高的！但是，這個滾動回歸的代碼也不是完美的，最大的劣勢就在于我們的滾動視窗是用“期”而不是用“天”來定義的，也就是說，程式在每次滾動的時候都會固定找前面n期的觀測，而不管這n期之間可能間隔的是10天，20天還是一個月。但是，對于大多數研究中，這種情況并不常見，是以大貓給出的代碼還是能夠應對大多數情況。