
5. Регрессионный анализ.pptx
- Количество слайдов: 34
Регрессионный анализ Финансовая эконометрика
Содержание • параметрическая регрессия • непараметрическая регрессия
Параметрическая регрессия
Линейная параметрическая регрессия •
Параметрическая регрессия в R # исходные данные library(datasets) ozone <- airquality$Ozone rad <- airquality$Solar. R rem <- is. na(ozone) | is. na(rad) ozone <- ozone[!rem]; rad <- rad[!rem] # разделим выборку на обучающую и экзаменующую N <- length(ozone); E <- 20; T <- N-E train. obs <- (1: T) eval. obs <- ((T+1): N) t. rad <- rad[train. obs]; t. ozone <- ozone[train. obs] e. rad <- rad[eval. obs]; e. ozone <- ozone[eval. obs]
Параметрическая регрессия в R # регрессионная модель fit. par <- lm(ozone ~ radiation, data=data. frame(radiation=t. rad, ozone=t. ozone), weights=NULL) # другой вариант fit. par <- lm(t. ozone ~ t. rad)
Анализ качества модели summary(fit. par) Estimate Std. Error t value Pr(>|t|) (Intercept) 22. 18688 8. 02140 2. 766 0. 00690 ** radiation 0. 12879 0. 03816 3. 375 0. 00109 ** Multiple R-squared: 0. 1135, Adjusted R-squared: 0. 1035 F-statistic: 11. 39 on 1 and 89 DF, p-value: 0. 001094 fit. par$coefficients fit. par$residuals fit. par$fitted. values plot(t. rad, t. ozone, pch=16, xlab="radiation", ylab="ozone") z <- order(t. rad) lines(t. rad[z], fit. par$fitted. values[z], col="blue", lwd=3)
Анализ остатков модели res <- fit. par$residuals hist(res) plot(res, type="l") # тесты на нормальность library(f. Basics) shapiro. test(res) jarquebera. Test(res) Тест p. value Shapiro – Wilks 0. 000 Jarque – Bera 0. 001
Переформулировка модели fit. par <- lm(log(ozone) ~ rad + rad 2, data=data. frame(rad=t. rad, rad 2=t. rad^2, ozone=t. ozone)) Estimate Std. Error t value Pr(>|t|) (Intercept) 1. 585 e+00 2. 514 e-01 6. 304 1. 12 e-08 *** rad 2. 337 e-02 3. 313 e-03 7. 055 3. 77 e-10 *** rad 2 -5. 660 e-05 9. 564 e-06 -5. 919 6. 11 e-08 *** Multiple R-squared: 0. 4234, Adjusted R-squared: 0. 4103 plot(t. rad, t. ozone, pch=16, xlab="radiation", ylab="ozone") z <- order(t. rad) lines(t. rad[z], exp(fit. par$fitted. values[z]), col="blue", lwd=3)
Анализ остатков модели res <- fit. par$residuals hist(res) plot(res, type="l") # тесты на нормальность shapiro. test(res) jarquebera. Test(res) Тест p. value Shapiro – Wilks 0. 288 Jarque – Bera 0. 279
Тесты на гетероскедастичность •
Тесты на гетероскедастичность •
Тесты на гетероскедастичность в R library(lmtest) # тест Бреуша–Пагана bptest(fit. par, varformula=NULL, data=NULL, studentize=FALSE) Breusch-Pagan test data: fit. par BP = 3. 2609, df = 2, p-value = 0. 1958 # тест Голдфелда–Куандта gqtest(fit. par, fraction=25, alternative="two. sided") Goldfeld-Quandt test data: fit. par GQ = 0. 6948, df 1 = 30, df 2 = 30, p-value = 0. 324
Учёт гетероскедастичности •
Тест на автокорреляцию •
Переформулировка модели ar 1 <- log(t. ozone[1: (T-1)]) fit. par <- lm(log(ozone) ~ rad + rad 2 + ar 1, data=data. frame(ozone=t. ozone[2: T], rad=t. rad[2: T], rad 2=t. rad[2: T]^2)) Estimate Std. Error t value Pr(>|t|) (Intercept) 7. 296 e-01 3. 201 e-01 2. 279 0. 025123 rad 1. 901 e-02 3. 313 e-03 5. 739 1. 40 e-07 rad 2 -4. 365 e-05 9. 601 e-06 -4. 546 1. 78 e-05 ar 1 3. 156 e-01 8. 086 e-02 3. 904 0. 000188 Тест p. value Shapiro – Wilks 0. 613 Jarque – Bera 0. 582 Breusch–Pagan 0. 213 Goldfeld–Quandt 0. 341 Durbin–Watson 0. 577 * *** ***
Построение прогноза •
Построение прогноза в R frc. par <- predict(fit. par, newdata=data. frame(rad=e. rad[2: E], rad 2=e. rad[2: E]^2, ar 1=log(e. ozone[1: (E-1)])), se. fit=TRUE, interval="prediction", level=0. 90) frc. par$fit frc. par$se. fit # прогнозные значения и интервалы # стандартные ошибки
Непараметрическая регрессия
Непараметрическая регрессия •
Непараметрическая регрессия в R # расчёт величины h library(np) bw <- npregbw(ozone ~ rad, ckertype="gaussian", bwtype="fixed", data=data. frame(ozone=t. ozone, rad=t. rad)) h <- bw$bw # ядро и функция Надарая – Уотсона kern <- function(x) exp(-(x^2/2))/sqrt(2*pi) NW <- function(x, x. dat, y. dat, h) { K 1 <- K 2 <- 0 N <- length(y. dat) for (i in 1: N) { K 1 <- K 1 + kern((x-x. dat[i])/h)*y. dat[i] K 2 <- K 2 + kern((x-x. dat[i])/h) } K 1 / K 2 }
График оценки plot(t. rad, t. ozone, pch=16) z <- order(t. rad) lines(t. rad[z], NW(t. rad, t. ozone, h)[z], col="blue")
График оценки (другой вариант) plot(t. rad, t. ozone, pch=16) fit. npar <- npreg(bw) ozone. hat <- predict(fit. npar) lines(t. rad[z], ozone. hat[z], col="blue")
•
Построение прогноза в R e <- t. ozone - ozone. hat s 2 <- var(e) g <- (4/(3*T))^(1/5)*sqrt(s 2) # симулированные значения ошибок b <- 10^4 e. star <- e[sample(1: T, size=b, replace=TRUE)]+g*rnorm(b) e. star <- sort(e. star) # прогноз и доверительные границы alpha <- 0. 1 y <- predict(fit. npar, newdata=data. frame(rad=e. rad)) bottom <- y + e. star[alpha/2*b] top <- y + e. star[(1 -alpha/2)*b]
Построение прогноза в R z <- order(e. rad) plot(e. rad, e. ozone, ylim=range(c(top, bottom))) lines(e. rad[z], y[z], col="blue", lwd=2) lines(e. rad[z], top[z], col="red", lty="dashed") lines(e. rad[z], bottom[z], col="red", lty="dashed")
•
Построение прогноза в R # моделирование условной дисперсии s 2 <- (e-mean(e))^2 h. res <- npregbw(s 2 ~ rad, ckertype="gaussian", bwtype="fixed", data=data. frame(rad=t. rad)) res. npar <- npreg(h. res) # метод бутстрапа b <- 10^4; alpha <- 0. 1 top <- bottom <- numeric(E) y <- predict(fit. npar, newdata=data. frame(rad=e. rad)) for (i in 1: E) { s 2. hat <- predict(res. npar, newdata=data. frame(rad=e. rad[i])) g <- (4/(3*T))^(1/5)*sqrt(s 2. hat) e. star <- e[sample(1: T, size=b, replace=TRUE)]+g*rnorm(b) e. star <- sort(e. star) bottom[i] <- y[i] + e. star[alpha/2*b] top[i] <- y[i] + e. star[(1 -alpha/2)*b] }
Построение прогноза в R z <- order(e. rad) plot(e. rad, e. ozone, ylim=range(c(top, bottom))) lines(e. rad[z], y[z], col="blue", lwd=2) lines(e. rad[z], top[z], col="red", lty="dashed") lines(e. rad[z], bottom[z], col="red", lty="dashed")
Непараметрическая регрессия, двумерный случай •
Множественная регрессия в R # добавочная объясняющая переменная temp <- airquality$Temp temp <- temp[!rem] t. temp <- temp[train. obs]; e. temp <- temp[eval. obs] # регрессионная модель bw 2 <- npregbw(ozone ~ rad + temp, ckertype="gaussian", bwtype="fixed", data=data. frame(ozone=t. ozone, rad=t. rad, temp=t. temp)) rad. temp <- cbind(t. rad, t. temp) ozone. hat <- NW 2(rad. temp, t. ozone, bw 2$bw, 2) # альтернативный вариант ozone. reg <- npreg(bw 2) ozone. hat <- predict(ozone. reg)
Множественная регрессия в R zt <- order(t. rad) plot(t. ozone[zt], pch=16) lines(ozone. hat[zt], col="blue", lwd=3)
Построение прогноза
Домашнее задание • рассмотреть данные о величине спроса на деньги в пакете lmtest: moneydemand • разделить выборку на обучающую и экзаменующую части • на пространстве обучающей выборки построить параметрическую и непараметрическую регрессионные модели с эндогенной переменной log. M • проверить качество параметрической модели с помощью тестов на нормальность, гетероскедастичность и автокорреляцию • построить прогноз и доверительные интервалы для эндогенной переменной на экзаменующей выборке • написать комментарии