對R中的XTS系列應用滾動窗口迴歸

我有5個貨幣對的1033個日收益點的xt，我想在其上運行滾動窗口迴歸，但rollapply不適用於我定義的使用lm的函數（）。這裏是我的數據：對R中的XTS系列應用滾動窗口迴歸

> head(fxr) 
       USDZAR  USDEUR  USDGBP  USDCHF  USDCAD 
2007-10-18 -0.005028709 -0.0064079963 -0.003878743 -0.0099537170 -0.0006153215 
2007-10-19 -0.001544470 0.0014275520 -0.001842564 0.0023058211 -0.0111410271 
2007-10-22 0.010878027 0.0086642116 0.010599365 0.0051899551 0.0173792230 
2007-10-23 -0.022783987 -0.0075236355 -0.010804304 -0.0041668499 -0.0144788687 
2007-10-24 -0.006561223 0.0008545792 0.001024275 -0.0004261666 0.0049525483 
2007-10-25 -0.014788901 -0.0048523001 -0.001434280 -0.0050425302 -0.0046422944 

> tail(fxr) 
       USDZAR  USDEUR  USDGBP  USDCHF  USDCAD 
2012-02-10 0.018619309 0.007548205 0.005526184 0.006348533 0.0067151342 
2012-02-13 -0.006449463 -0.001055966 -0.002206810 -0.001638002 -0.0016995755 
2012-02-14 0.006320364 0.006843933 0.006605875 0.005992935 0.0007001751 
2012-02-15 -0.001666872 0.004319096 -0.001568874 0.003686840 -0.0015009759 
2012-02-16 0.006419616 -0.003401364 -0.005194817 -0.002709588 -0.0019044761 
2012-02-17 -0.004339687 -0.003675992 -0.003319899 -0.003043481 0.0000000000

我可以很容易地在其上運行的LM整個數據集USDZAR對其他對型號：

> lm(USDZAR ~ ., data = fxr)$coefficients 
    (Intercept)  USDEUR  USDGBP  USDCHF  USDCAD 
-1.309268e-05 5.575627e-01 1.664283e-01 -1.657206e-01 6.350490e-01

但是我想運行一個滾動62天窗口讓這些係數隨着時間的演變，所以我創建了一個功能dolm它做到這一點：當我在這個rollapply運行

> dolm 
function(x) { 
    return(lm(USDZAR ~ ., data = x)$coefficients) 
}

但是我得到以下幾點：

> rollapply(fxr, 62, FUN = dolm) 
Error in terms.formula(formula, data = data) : 
    '.' in formula and no 'data' argument

那就是即使dolm（FXR）對自己的作品罰款：

> dolm(fxr) 
    (Intercept)  USDEUR  USDGBP  USDCHF  USDCAD 
-1.309268e-05 5.575627e-01 1.664283e-01 -1.657206e-01 6.350490e-01

這是怎麼回事嗎？它似乎工作正常，如果dolm是一個更簡單的功能，例如意思是：

> dolm <- edit(dolm) 
> dolm 
function(x) { 
    return(mean(x)) 
} 
> rollapply(fxr, 62, FUN = dolm) 
        USDZAR  USDEUR  USDGBP  USDCHF  USDCAD 
2007-11-29 -1.766901e-04 -6.899297e-04 6.252596e-04 -1.155952e-03 7.021468e-04 
2007-11-30 -1.266130e-04 -6.512204e-04 7.067767e-04 -1.098413e-03 7.247315e-04 
2007-12-03 8.949942e-05 -6.406932e-04 6.637066e-04 -1.154806e-03 8.727564e-04 
2007-12-04 2.042046e-04 -5.758493e-04 5.497422e-04 -1.116308e-03 7.124593e-04 
2007-12-05 7.343586e-04 -4.899982e-04 6.161819e-04 -1.057904e-03 9.915495e-04

任何幫助非常感謝。基本上我想要的是在62天滾動窗口中獲得USDZAR〜USDEUR + USDGBP + USDCHF + USDCAD迴歸的權重。

來源

2012-02-19 Thomas Browne

這裏有幾個問題：

rollapply通過矩陣而是lm需要data.frame。
rollapply分別對每列應用函數，除非我們指定by.column=FALSE。
你可能會或可能不希望得到的結果是與日期排列正確的，但如果你使用rollapplyr：

1）將上面我們有：

dolm <- function(x) coef(lm(USDZAR ~ ., data = as.data.frame(x)))) 
rollapplyr(fxr, 62, dolm, by.column = FALSE)

2）上述dolm中lm的替代方案是使用lm.fit，它可直接與矩陣一起使用，速度也更快：

dolm <- function(x) coef(lm.fit(cbind(Intercept = 1, x[,-1]), x[,1]))

來源

2012-02-19 17:24:56

非常感謝。是的，我在玩過很多遊戲後，也完成了它。傻我。當然，by.column = FALSE！非常感謝。只是閱讀你的動物園doc btw。好東西。我猜rollapply有點混亂的地方在於，雖然lm（）在整個xts上工作，但它並不在由rollapply（）返回的部分上。人們可以合理地預期rollapply返回將在lm（）下仍然工作的另一個xt，或者我錯過了什麼？ Mea在by.column上犯了錯誤，儘管如此。 – 2012-02-19 17:29:35

錯過了'rollapply'不是xts的一部分，但它是動物園及其調度'rollapply.zoo'的一部分。 – 2012-02-19 18:03:24

謝謝你澄清這一點。（fxr） > class（fxr） [1]「zoo」 > rollapply（fxr，62，function（x）coef（lm（USDZAR_x，data = x））， by.column = FALSE） model.frame.default中的錯誤（公式= USDZAR〜x，data = x，drop.unused.levels = TRUE）： 'data'必須是data.frame，而不是矩陣或array 所以我們仍然有這個問題。我明白...... R有很多這類問題，但仍然存在。我們在這裏所做的是lm適用於整個zoo對象，但不適用於它的rollapply子集。 – 2012-02-19 19:38:05

第三種選擇是更新QR分解中的R矩陣，如下面的代碼所做的那樣。你可以用C做加快這++，但比你需要從LINPACK的dchud和dchdd子程序（或其他功能更新R）

library(SamplerCompare) # for LINPACK `chdd` and `chud` 
roll_coef <- function(X, y, width){ 
    n <- nrow(X) 
    p <- ncol(X) 
    out <- matrix(NA_real_, n, p) 

    is_first <- TRUE 
    i <- width 
    while(i <= n){ 
    if(is_first){ 
     is_first <- FALSE 
     qr. <- qr(X[1:width, ]) 
     R <- qr.R(qr.) 

     # Use X^T for the rest 
     X <- t(X) 

     XtY <- drop(tcrossprod(y[1:width], X[, 1:width])) 
    } else { 
     x_new <- X[, i] 
     x_old <- X[, i - width] 

     # update R 
     R <- .Fortran(
     "dchud", R, p, p, x_new, 0., 0L, 0L, 
     0., 0., numeric(p), numeric(p), 
     PACKAGE = "SamplerCompare")[[1]] 

     # downdate R 
     R <- .Fortran(
     "dchdd", R, p, p, x_old, 0., 0L, 0L, 
     0., 0., numeric(p), numeric(p), integer(1), 
     PACKAGE = "SamplerCompare")[[1]] 

     # update XtY 
     XtY <- XtY + y[i] * x_new - y[i - width] * x_old 
    } 

    coef. <- .Internal(backsolve(R, XtY, p, TRUE, TRUE)) 
    out[i, ] <- .Internal(backsolve(R, coef., p, TRUE, FALSE)) 

    i <- i + 1 
    } 

    out 
} 

# simulate data 
set.seed(101) 
n <- 1000 
wdth = 100 
X <- matrix(rnorm(10 * n), n, 10) 
y <- drop(X %*% runif(10)) + rnorm(n) 
Z <- cbind(y, X) 

# assign other function 
dolm <- function(x) 
    coef(lm.fit(x[, -1], x[, 1])) 

# show that they yield the same 
library(zoo) 
all.equal(
    rollapply(Z, wdth, FUN = dolm, 
      by.column = FALSE, align = "right", fill = NA_real_), 
    roll_coef(X, y, wdth), 
    check.attributes = FALSE) 
#R> [1] TRUE 

# benchmark 
library(compiler) 
roll_coef <- cmpfun(roll_coef) 
dolm <- cmpfun(dolm) 
microbenchmark::microbenchmark(
    new = roll_coef(X, y, wdth), 
    prev = rollapply(Z, wdth, FUN = dolm, 
        by.column = FALSE, align = "right", fill = NA_real_), 
    times = 10) 
#R> Unit: milliseconds 
#R> expr  min   lq  mean  median   uq  max neval cld 
#R> new 8.631319 9.010579 9.808525 9.659665 9.973741 11.87083 10 a 
#R> prev 118.257128 121.734860 124.489826 122.882318 127.195410 135.21280 10 b

上面的解決方案要求形成model.matrix和model.response第一，但在致電roll_coef之前，這只是三個電話（model.frame的一個額外電話）。

來源

2018-02-18 23:01:11

對R中的XTS系列應用滾動窗口迴歸

回答

相關問題