Powiedzmy, że chcemy oszacować premię za ryzyko amerykańskiej giełdy w roku 2025 i 2026. Określenie "premia za ryzyko" używam jako skrót na oznaczenie nadwyżki stopy zwrotu z indeksu (tu S&P 500) nad stopę zwrotu z 10-letnich obligacji skarbowych. Precyzyjnie przecież mówiąc o premii mamy na myśli oczekiwane lub wymagane wynagrodzenie za przyszłe ryzyko, czyli takie, którego okres jeszcze nie zaczął się. I może się zrealizować (wtedy albo stracimy, albo mniej zarobimy), albo nie zrealizować (zarobimy co najmniej tyle, ile wymagamy). Dodatnia premia za ryzyko wynika z tego, że inwestorzy są niechętni wobec ryzyka. To znaczy, że wyceniają niżej aktywa o wyższym ryzyku. Mechanizm ten pozwala łatwo zrozumieć, dlaczego oczekiwana (a więc i średnia) stopa zwrotu jest dodatnia, a nie zerowa, jak w hazardzie. Otóż cała tajemnica leży w tym, że to ludzie wyceniają aktywa, a nie natura. Jeżeli wielu inwestorów będzie wyceniać niżej akcje z powodu ich większego ryzyka, to przecież spółka nie przestaje zarabiać tyle co zarabiała, jej wartość księgowa się nie zmienia. Spadek ceny oznacza więc, że inwestorzy mogą więcej zarobić, sprzedając w przyszłości akcje. Ryzyko jest mimo wszystko miarą psychologiczną, a więc dość subiektywną - i to właśnie sprawia, że inwestowanie daje dużo większe możliwości niż hazard.
A więc dane są poprawne, bo idealnie korelują. Oczywiście w tym okresie średnia była jeszcze wyższa i dla S&P U.S. Equity Risk Premium Index wyniosła 13,2%, a Damodarana 13,5%.
Pytanie brzmi czy średnia się zmieniała w ciągu tych 80 lat. Przez pierwsze 40 lat wyniosła ok. 8,6%, a kolejne niecałe 6,8%. Sprawdźmy czy możemy uznać to za coś więcej niż przypadek.
Chow test
data: premia ~ 1
F = 0.18, p-value = 0.7
Jak widać nie ma podstaw do odrzucenia hipotezy, że średnia jest stała. Wykonałem dodatkowo najróżniejsze testy (ADF, KPSS, Bai-Perrona, CUSUM, MOSUM) i nie wykryłem niestacjonarności czy niestabilności.
Modelowanie
Jeżeli średnia się nie zmienia, to trzeba sprawdzić autokorelację. Najpierw ACF:
acf(premia, 15)
Do dziś nie mogę z tego, że R pokazuje bez sensu ten zerowy rząd. Głupota. W cząstkowej jest już normalnie od pierwszego rzędu:
pacf(premia, 15)
W Dodatku do tego artykułu zrobiłem analogiczną analizę dla nadwyżki nad bonami. Tam PACF 4-tego rzędu jest na granicy istotności. Tutaj natomiast nie można mówić o jakiejkolwiek istotności. Z tego z kolei wynika, że obligacje skarbowe mają swój udział w cyklu 4-letnim (bo je tam dodaliśmy).
Mimo wszystko sprawdzę model AR(4) i jego stabilność. W R kod będzie następujący:
# Pierwszy model - arima
premia_model_arma1 <- arima(
premia,
order = c(4, 0, 0), # ARMA(4, 0, 0)
include.mean = TRUE, # estymuj μ
fixed = c(0, 0, 0, NA, NA), # φ1=0, φ2=0, φ3=0, φ4=estymuj, μ=estymuj
transform.pars = FALSE # aby fixed=0 było literalne
)
summary(premia_model_arma1)
Coefficients:
ar1 ar2 ar3 ar4 intercept
0 0 0 0.205 0.079
s.e. 0 0 0 0.111 0.026
sigma^2 estimated as 0.0344: log likelihood = 21.16, aic = -36.33
Gdyby zastosować taki model, dostalibyśmy rzeczywiście średnią 7,9% (+/- 2,6%), ale dodatkowo współczynnik 4-tego rzędu 0,205 (+/- 0,11).
Aby sprawdzić stabilność tego modelu, trzeba niestety zamienić na formę typu:
premia4 <- tail(zmienna, -4)
premia_4 <- head(zmienna, -4)
premia.model <- premia4 ~ premia_4
i przykładowo, żeby użyć testu CUSUM-Score, wpisać:
par(mfrow=c(3,1), mar=c(2.5,5,1.5,3))
cus = efp(premia.model, type="Score-CUSUM")
plot(cus, functional=NULL)
Za jednym zamachem mamy tu sprawdzenie stabilności wszystkich współczynników, w tym wariancji. Brak jakichkolwiek wychyleń (byłyby zaznaczenia) świadczy o stabilności.
To jednak tylko punkt wyjścia, bo znacznie lepszym sposobem jest użycie autoarfima() w rugarch lub auto.arima() w forecast. Wybrałem ten pierwszy, bo dostarcza wszystkie potrzebne statystyki. Ten argument omówiłem w tym artykule. Sprawdziłem więc w autoarfima max ARMA(7, 7) w metodzie full. Wybrałem rozkład "std", bo dawał najlepsze rezultaty, jeśli chodzi o statystyki dopasowania. Wadą rugarch jest czasochłonność obliczeń, dlatego koniecznie trzeba używać argumentu cluster tam gdzie się da. Ponadto wybrałem solver "nlminb", ponieważ wydaje się dość efektywny. Użyłem funkcji dwukrotnie - wg kryterium BIC i AIC (HQIC nie używam, bo danych jest za mało, zob. tu). Okazało się, że oba dały ten sam model. Zwracam też uwagę, że nie stosuję tutaj modelu GARCH ze względu na małą liczbę danych. Efekt:
# Drugi model - autoarfima (rugarch)
if (require("rugarch")==FALSE) {
install.packages("rugarch")
}
library("rugarch")
# Kryterium BIC
klaster <- makeCluster(detectCores() - 1)
premia_model_arma_bic <- autoarfima(data = premia, ar.max = 7, ma.max = 7, criterion = "BIC",
method = "full", arfima = FALSE, include.mean = NULL,
distribution.model = "std", cluster = klaster, external.regressors = NULL,
solver = "nlminb")
stopCluster(klaster)
premia_model_arma_bic
# Kryterium AIC
klaster <- makeCluster(detectCores() - 1)
premia_model_arma_aic <- autoarfima(data = premia, ar.max = 7, ma.max = 7, criterion = "AIC",
method = "full", arfima = FALSE, include.mean = NULL,
distribution.model = "std", cluster = klaster, external.regressors = NULL,
solver = "nlminb")
stopCluster(klaster)
premia_model_arma_aic
premia_model_arma2 <- premia_model_arma_aic
*----------------------------------*
* ARFIMA Model Fit *
*----------------------------------*
Mean Model : ARFIMA(7,0,4)
Distribution : std
Optimal Parameters
------------------------------------
Estimate Std. Error t value Pr(>|t|)
ar1 0.44901 0.000926 484.974 0
ar2 0.00000 NA NA NA
ar3 0.00000 NA NA NA
ar4 0.00000 NA NA NA
ar5 -0.22529 0.000465 -484.974 0
ar6 0.00000 NA NA NA
ar7 0.34287 0.000027 12785.413 0
ma1 -1.28980 0.000061 -21191.451 0
ma2 0.00000 NA NA NA
ma3 0.61606 0.000029 21191.451 0
ma4 0.44741 0.000021 21191.451 0
sigma 0.13607 0.005403 25.184 0
shape 100.00000 0.004585 21809.383 0
Robust Standard Errors:
Estimate Std. Error t value Pr(>|t|)
ar1 0.44901 2.324983 0.19312 0.84686
ar2 0.00000 NA NA NA
ar3 0.00000 NA NA NA
ar4 0.00000 NA NA NA
ar5 -0.22529 1.166585 -0.19312 0.84686
ar6 0.00000 NA NA NA
ar7 0.34287 0.000385 891.10681 0.00000
ma1 -1.28980 0.005203 -247.88844 0.00000
ma2 0.00000 NA NA NA
ma3 0.61606 0.000766 804.54128 0.00000
ma4 0.44741 0.001619 276.38488 0.00000
sigma 0.13607 0.619235 0.21974 0.82608
shape 100.00000 1.325372 75.45050 0.00000
LogLikelihood : 46.12
Information Criteria
------------------------------------
Akaike -0.95299
Bayes -0.71479
Shibata -0.97067
Hannan-Quinn -0.85749
Weighted Ljung-Box Test on Standardized Residuals
------------------------------------
statistic p-value
Lag[1] 0.000611 0.9803
Lag[2*(p+q)+(p+q)-1][32] 11.645908 1.0000
Lag[4*(p+q)+(p+q)-1][54] 17.816677 0.9976
H0 : No serial correlation
Weighted Ljung-Box Test on Standardized Squared Residuals
------------------------------------
statistic p-value
Lag[1] 0.4785 0.4891
Lag[2*(p+q)+(p+q)-1][2] 0.7811 0.5750
Lag[4*(p+q)+(p+q)-1][5] 2.3214 0.5447
ARCH LM Tests
------------------------------------
Statistic DoF P-Value
ARCH Lag[2] 1.519 2 0.4678
ARCH Lag[5] 2.647 5 0.7542
ARCH Lag[10] 9.429 10 0.4920
Nyblom stability test
------------------------------------
Joint Statistic: 231429
Individual Statistics:
ar1 0.09969
ar5 0.05997
ar7 0.01355
ma1 0.03352
ma3 0.11951
ma4 0.09983
sigma 0.13154
shape 0.42456
Asymptotic Critical Values (10% 5% 1%)
Joint Statistic: 1.89 2.11 2.59
Individual Statistic: 0.35 0.47 0.75
BIC i AIC wskazały ten sam model: AR o rzędach 1, 5 i 7; MA o rzędach 1, 3 i 4. W nieodpornych statystykach (tj. przy założeniu danego rozkładu, tu t-Studenta) wszystkie są istotne. Jednak statystyki odporne (na odchylenia od założonego rozkładu reszt) wskazują jedynie istotność rzędu 7 dla AR, a dla MA wszystkie są istotne. Pozostałe statystyki są okej, wszystkie współczynniki są stabilne, chociaż znowu cały model mocno niestabilny (pewnie niezbyt idealny rozkład).
Test MCS
Aby się upewnić, że na pewno drugi model jest lepszy od pierwszego, użyjemy testu MCS, który omówiłem w poprzednim artykule. Tak jak wtedy za stratę uznam sytuację, gdy znak prognozy in-sample nie zgadza się ze znakiem rzeczywistej wartości.
#MCS test
Funkcja_strat <- function(rzeczywiste, prognoza) {
roznica <- rzeczywiste - prognoza
strata <- ifelse(sign(rzeczywiste) != sign(prognoza), abs(roznica), 0)
return(strata)
}
# Obliczenie funkcji strat dla każdego modelu (in-sample)
straty1 <- Funkcja_strat(premia, as.numeric(fitted(premia_model_arma1)))
straty2 <- Funkcja_strat(premia, premia_model_arma2$fit@fit$fitted.values)
straty <- cbind(straty1, straty2)
mcs_stationary <- mcsTest(straty, alpha = 0.05, nboot = 500, nblock = 10, boot = "stationary")
print(mcs_stationary)
$includedR
[1] 1 2
$pvalsR
[,1]
[1,] 0.23
[2,] 1.00
$excludedR
NULL
$includedSQ
[1] 1 2
$pvalsSQ
[,1]
[1,] 0.23
[2,] 1.00
$excludedSQ
NULL
Test MCS nie pozwolił wykluczyć żadnego z dwóch modeli, ale wskazuje dużo bardziej na zachowanie drugiego, tj. z autoarfima, więc go wybierzemy.
Wykres modelu
Wobec tego sprawdźmy jak te dopasowania dla najlepszego modelu wyglądają:
# Wartości dopasowane
dopas <- premia_model_arma2$fit@fit$fitted.values
okres <- c(1945:2024)
# Wykres najlepszego modelu z prawdziwymi obserwacjami
plot(x=okres, y=premia, type="l", ylab = "")
lines(x=okres, y=dopas, col="red")
abline(h=0, col="grey")
abline(h=mean(premia), col="blue", lwd=2)
legend("topright", legend=c("premia za ryzyko S&P 500", "prognoza in-sample", "średnia premia za ryzyko"),
col=c("black","red", "blue"), lwd=2, cex=0.6)
Prognoza na 2025-2026
Niestety nie mam wiedzy, czy można od razu w rugarch robić prognozy mając już najlepszy model z autoarfima. jestem niemal pewny, że jest sposób, ale szybciej jest mi tak zrobić. Dlatego robię to "na piechotę". Obliczam arfimaspec(), z tego arfimafit i dopiero z tego arfimaforecast(). W arfimaspec() wpisujemy normalnie maksymalne rzędy ARMA, ale za to dodajemy argument fixed.pars = list(ar2=0, ar3=0), który właśnie kontroluje ustawienie konkretnych współczynników.
rzadAR <- premia_model_arma2$fit@model$modelinc[2]
rzadMA <- premia_model_arma2$fit@model$modelinc[3]
premia_spec <- arfimaspec(mean.model = list(armaOrder = c(rzadAR, rzadMA), include.mean = TRUE,
arfima = FALSE, external.regressors = NULL), distribution.model = "std",
start.pars = list()
, fixed.pars = list(ar2=0,
ar3=0,
ar4=0,
ar6=0,
ma2=0)
)
premia_fit <- arfimafit(spec = premia_spec, data = premia, solver = "nlminb")
prognoza <- arfimaforecast(premia_fit, n.ahead = 2)
prognoza
*----------------------------------*
* ARFIMA Model Forecast *
*----------------------------------*
Horizon: 2
Roll Steps: 0
Out of Sample: 0
0-roll forecast:
T+1 T+2
0.05436 0.07451
Na rok 2025 dostajemy premię niecałe 5,4% i na 2026 7,5%. Dodajemy nowe prognozy do wykresu:
nowa_prognoza <- as.numeric(prognoza@forecast$seriesFor)
# Wykres najlepszego modelu z prawdziwymi obserwacjami
plot(x=c(okres, 2025, 2026), y=c(premia, NA, NA), type="l", ylab = "")
lines(x=c(okres, 2025, 2026), y=c(dopas, NA, NA), col="red")
lines(x=c(okres, 2025, 2026), y=c(rep(NA, length(premia)-1), tail(dopas, 1), nowa_prognoza), col="green", lwd=2, type="o", cex=0.4)
abline(h=0, col="grey")
abline(h=mean(premia), col="blue", lwd=2)
legend("topright", legend=c("premia za ryzyko S&P 500", "prognoza in-sample", "prognoza poza próbą", "średnia premia za ryzyko"),
col=c("black","red", "green", "blue"), lwd=2, cex = 0.6)
Prognoza krocząca
Czerwone linie są tak naprawdę jedynie dopasowaniem do danych, a nie prognozą. Aby uzyskać bardziej realny obraz dobroci modelu, wykonam prognozy kroczące za pomocą arfimaroll(), czyli prognozy out-of-sample. Jest to odpowiednik ugarchroll bez procesu GARCH (zob.
ten art.). Zobaczmy 10 ostatnich kroczących-nawijanych prognoz, tj. od roku 2015:
# Prognoza krocząca-nawijana
klaster <- makeCluster(detectCores() - 1)
premia_prognoza_roll <- arfimaroll(spec = premia_spec, data = premia, forecast.length = 10,
refit.every = 1, refit.window = "recursive",
window.size = NULL, solver = "nlminb", fit.control = list(ar2=0, ar3=0),
cluster=klaster)
premia_prognoza_roll <- resume(premia_prognoza_roll)
stopCluster(klaster)
prognoza_roll <- premia_prognoza_roll@forecast$density$Mu
okres_prognozy_roll <- c(2015:2024)
names(prognoza_roll) <- okres_prognozy_roll
prognoza_roll
# Wykres z prognozą kroczącą
plot(x=c(okres, 2025, 2026), y=c(premia, NA, NA), type="l", ylab = "", xlim=c(2015, 2026))
lines(x=c(okres, 2025, 2026), y=c(dopas, NA, NA), col="red")
lines(x=c(okres, 2025, 2026), y=c(rep(NA, length(premia)-1), tail(prognoza_roll, 1), nowa_prognoza), col="green", lwd=2, type="o", cex=0.4)
lines(x=c(NA, okres_prognozy_roll), y=c(tail(dopas, 1), prognoza_roll), col="green", lwd=2, type="o", cex=0.4)
abline(h=0, col="grey")
abline(h=mean(premia), col="blue", lwd=2)
legend("topright", legend=c("premia za ryzyko S&P 500", "prognoza in-sample", "prognoza poza próbą", "średnia premia za ryzyko"),
col=c("black","red", "green", "blue"), lwd=2, cex = 0.6)
Prognoza sugeruje spadek premii w 2025 i może wzrost 2026. Zauważmy tylko, że spadek nie poniżej zera, tylko poniżej średniej. No i oczywiście nie mówi to nic na temat zmian w obligacjach.
Rozkład częstości i gęstości
Skupiając się na prognozowaniu łatwo jednak zgubić szeroki obraz, tj. ekonomiczne znaczenie premii za ryzyko. Dlatego musimy mieć w pamięci, że prawdopodobnie tak będzie, ale faktycznie może się zdarzyć niespodzianka. Ostatnie lata przyzwyczajają nie tylko zwykłych ludzi, ale i ekspertów do tego, że co roku z indeksów USA jest dodatnia stopa zwrotu. Od 2009 r do 2024, czyli przez 16 lat zdarzyły się zaledwie 2 ujemne roczne (całkowite) stopy zwrotu! Tworzy to złudzenie, że ryzyka nie ma, bo wystarczy przeczekać parę lat. Jednak ekonomia robi swoje i musi to jakoś wyrównać. Dlatego muszą nie tylko przydarzać się szokowe, bardzo silne i szybkie spadki (prędzej czy później). Muszą występować dłuższe okresy, w których obligacje skarbowe przyniosą większe zwroty od akcji. Gdyby nie ich nie było, to rzeczywiście nie byłoby ryzyka.
Nawet historyczna ujemna skośność (średnia < mediana) i kurtoza w rozkładzie premii o tym świadczy, zobaczmy:
# Zapisujemy stare ustawienia graficzne
old_par <- par(no.readonly = TRUE)
# Dajemy więcej miejsca po prawej (dół, lewo, góra, prawo)
par(mar = c(3.5, 5, 4, 5) - 0.1, mgp = c(2, 1, 0))
# 1. Obliczamy histogram, ale nie rysujemy
h <- hist(premia,
breaks = 10,
plot = FALSE)
# 2. Estymujemy KDE (Kernel Density Estimation)
d <- density(premia,
kernel = "gaussian",
bw = "nrd0")
# 3. Rysujemy histogram na lewą oś Y (counts)
plot(h,
freq = TRUE,
col = "lightgray",
border = "white",
xlab = "Premia",
ylab = "Częstość",
xlim = range(h$breaks, d$x),
ylim = c(0, max(h$counts) * 1.1),
main = "Histogram i gęstość rozkładu – lewa oś: częstość, prawa oś: gęstość")
# 4. Nałóż KDE na ten sam wykres
par(new = TRUE)
plot(d$x, d$y,
type = "l",
lwd = 2,
axes = FALSE,
xlab = "",
ylab = "",
ylim = c(0, max(d$y) * 1.1))
# 5. Dodajemy prawą oś Y (density)
axis(side = 4,
at = pretty(c(0, max(d$y) * 1.1)))
mtext("Gęstość (KDE)", side = 4, line = 3)
# 7. Przywracamy poprzednie ustawienia
par(old_par)
Jeżeli rynek jest efektywny w długim okresie, to powinniśmy się spodziewać, że podobny rozkład wystąpi dla okresów wieloletnich. Zobaczmy dla 4-letnich. Zagregujmy stopy zwrotu:
grupy <- ceiling(1:length(premia)/4)
premia_4lata <- exp(tapply(X = log(1+premia), INDEX = grupy, FUN = sum)) - 1
I powtarzamy kod dla wykresu rozkładu częstości i gęstości, zastępując zmienną premia przez premia_4lata. Dostajemy:
Rozkład gęstości 4-letnich premii zawiera niewiele danych, ale też występuje na nim skośność.
Dodatek dla premii za ryzyko nad bonami skarbowymi 3-miesięcznymi
Powtórzmy całą procedurę dla nadwyżki ponad bony skarbowe. To znaczy premia zawiera również ryzyko długoterminowych obligacji skarbowych w postaci zmienności stopy procentowej.
Stabilność średniej
W całym okresie średnia wyniosła aż 9%. Odwrotnie niż dla obligacji pierwsze 40 lat dały mniej, bo 8%, a kolejne prawie 10%. Czy więc tym razem to istotna zmiana średniej? Test Chowa wskazuje na stabilność:
Chow test
data: premia ~ 1
F = 0.22, p-value = 0.6
ACF i PACF
acf(premia, 15)
pacf(premia, 15)
Okazuje się, że jedynie ta 4-tego rzędu jest na granicy istotności (przekracza 0,2). To potwierdzałoby istnienie 4-letnich cykli, chociaż wykorzystanie tego fenomenu na pewno nie jest proste.
Modelowanie
Sprawdzamy zatem AR(4). Efekt:
Coefficients:
ar1 ar2 ar3 ar4 intercept
0 0 0 0.236 0.090
s.e. 0 0 0 0.112 0.024
sigma^2 estimated as 0.0282: log likelihood = 29.05, aic = -52.11
Training set error measures:
ME RMSE MAE MPE MAPE MASE ACF1
Training set 0.000006594 0.168 0.1362 -266.2 433.2 0.6704 -0.06222 |
Gdyby zastosować taki model, dostalibyśmy średnią 9% (+/- 2,4%) - faktycznie, tyle co średnia arytmetyczna. Dodatkowo współczynnik 4-tego rzędu 0,24 (+/- 0,11).
To teraz autoarfima z tymi samymi ustawieniami.
*----------------------------------*
* ARFIMA Model Fit *
*----------------------------------*
Mean Model : ARFIMA(4,0,7)
Distribution : std
Optimal Parameters
------------------------------------
Estimate Std. Error t value Pr(>|t|)
mu 0.098811 0.000012 8073.8245 0.00000
ar1 0.000000 NA NA NA
ar2 0.000000 NA NA NA
ar3 0.945906 0.000100 9503.1496 0.00000
ar4 0.054088 0.000008 6695.8385 0.00000
ma1 0.279539 0.000033 8524.8266 0.00000
ma2 0.000000 NA NA NA
ma3 -1.685603 0.000062 -27299.6027 0.00000
ma4 -0.135260 0.000014 -9448.6750 0.00000
ma5 -0.118209 0.000014 -8582.6486 0.00000
ma6 0.000000 NA NA NA
ma7 -0.036416 0.000006 -5762.9626 0.00000
sigma 0.111497 0.009553 11.6718 0.00000
shape 16.523609 12.901423 1.2808 0.20028
Robust Standard Errors:
Estimate Std. Error t value Pr(>|t|)
mu 0.098811 0.000566 174.510530 0.00000
ar1 0.000000 NA NA NA
ar2 0.000000 NA NA NA
ar3 0.945906 0.010638 88.917365 0.00000
ar4 0.054088 0.000072 750.823674 0.00000
ma1 0.279539 0.002076 134.669756 0.00000
ma2 0.000000 NA NA NA
ma3 -1.685603 0.002733 -616.737764 0.00000
ma4 -0.135260 0.000331 -408.478534 0.00000
ma5 -0.118209 0.000610 -193.739688 0.00000
ma6 0.000000 NA NA NA
ma7 -0.036416 0.000020 -1799.605712 0.00000
sigma 0.111497 0.093611 1.191060 0.23363
shape 16.523609 1287.342190 0.012835 0.98976
LogLikelihood : 62.2351
Information Criteria
------------------------------------
Akaike -1.3059
Bayes -1.0081
Shibata -1.3327
Hannan-Quinn -1.1865
Weighted Ljung-Box Test on Standardized Residuals
------------------------------------
statistic p-value
Lag[1] 3.101 0.07826
Lag[2*(p+q)+(p+q)-1][32] 17.518 0.04329
Lag[4*(p+q)+(p+q)-1][54] 27.653 0.46556
H0 : No serial correlation
Weighted Ljung-Box Test on Standardized Squared Residuals
------------------------------------
statistic p-value
Lag[1] 0.3235 0.5695
Lag[2*(p+q)+(p+q)-1][2] 0.6261 0.6373
Lag[4*(p+q)+(p+q)-1][5] 1.6151 0.7116
ARCH LM Tests
------------------------------------
Statistic DoF P-Value
ARCH Lag[2] 0.6869 2 0.7093
ARCH Lag[5] 2.5554 5 0.7681
ARCH Lag[10] 5.7786 10 0.8335
Nyblom stability test
------------------------------------
Joint Statistic: 411.572
Individual Statistics:
mu 0.03948
ar3 0.04534
ar4 0.02661
ma1 0.04062
ma3 0.15471
ma4 0.02460
ma5 0.03678
ma7 0.02549
sigma 0.13303
shape 0.08077
Asymptotic Critical Values (10% 5% 1%)
Joint Statistic: 2.29 2.54 3.05
Individual Statistic: 0.35 0.47 0.75
Dostałem model ARMA(4, 7), ale przy usuniętych ar1, ar2, ma2 i ma6.
Jednak całkowity model okazuje się niestabilny, bo Joint Statistic = 412. Za to poszczególne współczynniki są stabilne shape, czyli kształt rozkładu reszt, jest niestabilny. Wynika to po prostu z tego, że rozkład t-Studenta nie jest idealnym rozkładem tutaj.
Test MCS
$includedR
[1] 2
$pvalsR
[,1]
[1,] 0
[2,] 1
$excludedR
[1] 1
$includedSQ
[1] 2
$pvalsSQ
[,1]
[1,] 0
[2,] 1
$excludedSQ
[1] 1
Test MCS potwierdza jednoznacznie, że drugi model (uzyskany z autoarfima) jest lepszy, albowiem każe odrzucić pierwszy model.
Porównując ze sobą obydwa rodzaje premii (nad bonami i nad obligacjami), to co się rzuca w oczy, to wpływ obligacji skarbowych na 4-letnią cykliczność stopy zwrotu. Po odjęciu bonów premia wykazuje pewną autokorelację 4-tego rzędu. Po odjęciu obligacji skarbowych premia ta autokorelacja mieści się w granicach nieistotności, a więc przypadkowości (jednak zwracam uwagę, że to PACF to taka klasyczna metoda szukania cykli, raczej wstępna). Wydaje się, że taka dysproporcja wynika z tego, że długoterminowe obligacje odzwierciedlają oczekiwania koniunktury, a ta przecież porusza się właśnie w takich kilkuletnich cyklach. Dostajemy więc pewnego rodzaju "dowód", że indeks S&P 500 podąża w jakimś zakresie w kierunku zmian gospodarczych.
