估計約翰森協同方法中的約翰森協整方法在R約爲

Question

我正在通過在R的「變量」庫中的小插曲工作示例工作。我明白大多數在小插曲期望表中的示例表5的小插曲here

運行下面的代碼，我看到協整矢量和加載點估計的派生，但我不明白在哪裏得到的t統計。

我跑這個代碼：

library("vars") data("Canada") Canada <- Canada[, c("prod", "e", "U", "rw")] vecm <- ca.jo(Canada[, c("rw", "prod", "e", "U")], type = "trace", ecdet = "trend", K = 3, spec = "transitory") vecm.r1 <- cajorls(vecm, r = 1)

，並獲得這些特徵向量和權重，

Eigenvectors, normalised to first column: 
(These are the cointegration relations) 

       rw.l1  prod.l1  e.l1  U.l1 trend.l1 
rw.l1  1.00000000 1.0000000 1.0000000 1.000000 1.0000000 
prod.l1 0.54487553 -3.0021508 0.7153696 -7.173608 0.4087221 
e.l1  -0.01299605 -3.8867890 -2.0625220 -30.429074 -3.3884676 
U.l1  1.72657188 -10.2183404 -5.3124427 -49.077209 -5.1326687 
trend.l1 -0.70918872 0.6913363 -0.3643533 11.424630 0.1157125 

Weights W: 
(This is the loading matrix) 

       rw.l1  prod.l1  e.l1   U.l1  trend.l1 
rw.d -0.084814510 0.048563997 -0.02368720 -0.0016583069 5.722004e-12 
prod.d -0.011994081 0.009204887 -0.09921487 0.0020567547 -7.478364e-12 
e.d -0.015606039 -0.038019447 -0.01140202 -0.0005559337 -1.229460e-11 
U.d -0.008659911 0.020499657 0.02896325 0.0009140795 1.103862e-11 

這給我表了正確的α和β值5

但我不知道表5中的t統計量是在代碼中得出的。外面的任何人都可以指出我正確的方向嗎？

問候，

詹姆斯

Answer 1

我最大的抱怨與urca：它不計算t統計量的測試向量。看起來作者可能在某個時候想要包含該特徵（因此表5列在你上面引用的小插圖中），但是並沒有真正做到這一點。

好消息：對於一個協整向量（r = 1）的情況，有示例代碼可用於獲取標準錯誤。看到這個discussion thread

與urca的長期協整關係的標準誤差時，有一個單一的協整向量

從你上面的例子中繼續，上述討論的nabble線程給出的代碼是：

alpha <- coef(vecm.r1$rlm)[1, ] # the coefficients on ecm1 
beta <- vecm.r1$beta # the point estimates of beta 

resids <- resid(vecm.r1$rlm) 
N <- nrow(resids) 
sigma <- crossprod(resids)/N 

## t-stats for beta 
beta.se <- sqrt(diag(kronecker(solve(crossprod(vecm@RK[, -1])), 
           solve(t(alpha) %*% solve(sigma) %*% alpha)))) 
beta.t <- c(NA, beta[-1]/beta.se) 
names(beta.t) <- rownames(vecm.r1$beta) 
beta.t 

這嚴重影響了小插曲給出的標準錯誤。

與格蕾特

比較，如果我們運行與格蕾特相同的模型，我們得到的測試版完全相同的點估計，但標準差大致同意：

獎勵：創建乳膠表使用urca和texreg的標準誤差的協整向量

您需要編寫一個稍微笨拙的「提取」函數，將標準誤差和點估計一起組裝到窗體中這texreg喜歡：

extract.cajo_beta <- function(cajo, orls) { 
    alpha <- coef(orls$rlm)[1, ]; 
    resids <- resid(orls$rlm); 
    N <- nrow(resids); 
    sigma <- crossprod(resids)/N; 

    # get standard errors and p-values 
    beta <- orls$beta 
    beta.se <- sqrt(diag(kronecker(solve(crossprod(cajo@RK[, -1])), 
            solve(t(alpha) %*% solve(sigma) %*% alpha)))); 
    beta.se2 <- c(NA, beta.se); 
    beta.t <- c(NA, beta[-1]/beta.se); 
    beta.pval <- dt(beta.t, df= orls$rlm$df.residual) 

    tr <- createTexreg(coef.names = as.character(rownames(beta)), coef = as.numeric(beta), se = beta.se2, 
     pvalues = beta.pval, 
     # remove this goodness of fit measure afterwards from the table 
     gof.names = c('Dummy'), gof = c(1), gof.decimal = c(FALSE) 
    ); 
    return(tr); 
} 

然後運行：

> screenreg(extract.cajo_beta(vecm, vecm.r1)) 

================= 
      Model 1 
----------------- 
rw.l1  1.00 

prod.l1 0.54 
      (0.61) 
e.l1  -0.01 
      (0.68) 
U.l1  1.73 
      (1.45) 
trend.l1 -0.71 * 
      (0.28) 
----------------- 
Dummy  1  
================= 
*** p < 0.001, ** p < 0.01, * p < 0.05