Fit a GARCH (1,1) - model dengan kovariat di R

10

Saya memiliki beberapa pengalaman dengan pemodelan deret waktu, dalam bentuk model ARIMA sederhana dan sebagainya. Sekarang saya memiliki beberapa data yang menunjukkan pengelompokan volatilitas, dan saya ingin mencoba memulai dengan memasang model GARCH (1,1) pada data.

Saya memiliki serangkaian data dan sejumlah variabel yang saya pikir mempengaruhinya. Jadi dalam istilah regresi dasar, sepertinya:

yt=α+β1xt1+β2xt2+ϵt.

Tapi saya benar-benar bingung bagaimana menerapkan ini menjadi model GARCH (1,1)? Saya telah melihat rugarch-paket dan fGarch-paket dalam R, tapi saya belum bisa melakukan sesuatu yang berarti selain contoh yang dapat ditemukan di internet.

Thorst
sumber

Jawaban:

12

Berikut adalah contoh implementasi menggunakan rugarchpaket dan dengan beberapa data palsu. Fungsi ini ugarchfitmemungkinkan untuk dimasukkannya regressor eksternal dalam persamaan rata-rata (perhatikan penggunaan external.regressorsdi fit.specdalam kode di bawah).

Untuk memperbaiki notasi, modelnya adalah mana dan menunjukkan kovariat pada waktu , dan dengan asumsi / persyaratan "biasa" pada parameter dan proses inovasi .

yt=λ0+λ1xt,1+λ2xt,2+ϵt,ϵt=σtZt,σt2=ω+αϵt12+βσt12,
xt,1xt,2tZt

Nilai parameter yang digunakan dalam contoh adalah sebagai berikut.

## Model parameters
nb.period  <- 1000
omega      <- 0.00001
alpha      <- 0.12
beta       <- 0.87
lambda     <- c(0.001, 0.4, 0.2)

Gambar di bawah ini menunjukkan seri kovariat dan serta seri . The kode yang digunakan untuk menghasilkan mereka disediakan di bawah ini.xt,1xt,2ytR

masukkan deskripsi gambar di sini

## Dependencies
library(rugarch)

## Generate some covariates
set.seed(234)
ext.reg.1 <- 0.01 * (sin(2*pi*(1:nb.period)/nb.period))/2 + rnorm(nb.period, 0, 0.0001)
ext.reg.2 <- 0.05 * (sin(6*pi*(1:nb.period)/nb.period))/2 + rnorm(nb.period, 0, 0.001)
ext.reg   <- cbind(ext.reg.1, ext.reg.2)

## Generate some GARCH innovations
sim.spec    <- ugarchspec(variance.model     = list(model = "sGARCH", garchOrder = c(1,1)), 
                          mean.model         = list(armaOrder = c(0,0), include.mean = FALSE),
                          distribution.model = "norm", 
                          fixed.pars         = list(omega = omega, alpha1 = alpha, beta1 = beta))
path.sgarch <- ugarchpath(sim.spec, n.sim = nb.period, n.start = 1)
epsilon     <- as.vector(fitted(path.sgarch))

## Create the time series
y <- lambda[1] + lambda[2] * ext.reg[, 1] + lambda[3] * ext.reg[, 2] + epsilon

## Data visualization
par(mfrow = c(3,1))
plot(ext.reg[, 1], type = "l", xlab = "Time", ylab = "Covariate 1")
plot(ext.reg[, 2], type = "l", xlab = "Time", ylab = "Covariate 2")
plot(y, type = "h", xlab = "Time")
par(mfrow = c(1,1))

Penyesuaian dilakukan dengan ugarchfitsebagai berikut.

## Fit
fit.spec <- ugarchspec(variance.model     = list(model = "sGARCH",
                                                 garchOrder = c(1, 1)), 
                       mean.model         = list(armaOrder = c(0, 0),
                                                 include.mean = TRUE,
                                                 external.regressors = ext.reg), 
                       distribution.model = "norm")
fit      <- ugarchfit(data = y, spec = fit.spec)

Estimasi parameter adalah

## Results review
fit.val     <- coef(fit)
fit.sd      <- diag(vcov(fit))
true.val    <- c(lambda, omega, alpha, beta)
fit.conf.lb <- fit.val + qnorm(0.025) * fit.sd
fit.conf.ub <- fit.val + qnorm(0.975) * fit.sd
> print(fit.val)
#     mu       mxreg1       mxreg2        omega       alpha1        beta1 
#1.724885e-03 3.942020e-01 7.342743e-02 1.451739e-05 1.022208e-01 8.769060e-01 
> print(fit.sd)
#[1] 4.635344e-07 3.255819e-02 1.504019e-03 1.195897e-10 8.312088e-04 3.375684e-04

Dan nilai sebenarnya yang sesuai adalah

> print(true.val)
#[1] 0.00100 0.40000 0.20000 0.00001 0.12000 0.87000

Gambar berikut menunjukkan estimasi parameter dengan interval kepercayaan 95%, dan nilai sebenarnya. The Rkode yang digunakan untuk menghasilkan itu disediakan di bawah.

masukkan deskripsi gambar di sini

plot(c(lambda, omega, alpha, beta), pch = 1, col = "red",
     ylim = range(c(fit.conf.lb, fit.conf.ub, true.val)),
     xlab = "", ylab = "", axes = FALSE)
box(); axis(1, at = 1:length(fit.val), labels = names(fit.val)); axis(2)
points(coef(fit), col = "blue", pch = 4)
for (i in 1:length(fit.val)) {
    lines(c(i,i), c(fit.conf.lb[i], fit.conf.ub[i]))
}
legend( "topleft", legend = c("true value", "estimate", "confidence interval"),
        col = c("red", "blue", 1), pch = c(1, 4, NA), lty = c(NA, NA, 1), inset = 0.01)
QuantIbex
sumber
bagaimana Anda memperkirakan parameter (lambda, omega, alpha, beta) ??
chs
1
@ chs Perkiraan parameter diperoleh dengan ugarchfitfungsi.
QuantIbex