stochastyczny vs deterministyczny trend / sezonowość w prognozowaniu szeregów czasowych

16

Mam umiarkowane doświadczenie w prognozowaniu szeregów czasowych. Przejrzałem kilka książek o prognozowaniu i nie widzę w nich odpowiedzi na następujące pytania.

Mam dwa pytania:

  1. Jak określiłbym obiektywnie (za pomocą testu statystycznego), czy dany szereg czasowy ma:

    • Sezonowość stochastyczna lub sezonowość deterministyczna
    • Trend stochastyczny lub trend deterministyczny
  2. Co by się stało, gdybym modelował moje szeregi czasowe jako deterministyczny trend / sezonowość, gdy szereg ma wyraźnie stochastyczny komponent?

Będziemy wdzięczni za wszelką pomoc w rozwiązaniu tych pytań.

Przykładowe dane trendu:

7,657
5,451
10,883
9,554
9,519
10,047
10,663
10,864
11,447
12,710
15,169
16,205
14,507
15,400
16,800
19,000
20,198
18,573
19,375
21,032
23,250
25,219
28,549
29,759
28,262
28,506
33,885
34,776
35,347
34,628
33,043
30,214
31,013
31,496
34,115
33,433
34,198
35,863
37,789
34,561
36,434
34,371
33,307
33,295
36,514
36,593
38,311
42,773
45,000
46,000
42,000
47,000
47,500
48,000
48,500
47,000
48,900
Synoptyk
źródło
Istnieją 4 możliwe stany natury. Nie ma analitycznego rozwiązania tego pytania, ponieważ przestrzeń próbki modelu jest względnie nieograniczona. Aby empirycznie odpowiedzieć na to dokuczliwe pytanie, pomogłem opracować AUTOBOX autobox.com/cms . AUTOBOX przeprowadza turniej, aby zbadać wszystkie 4 te przypadki i ocenić jakość 4 wynikowych modeli pod kątem konieczności i wystarczalności. Dlaczego nie opublikujesz przykładowego szeregu czasowego, a ja opublikuję 4 wyniki pokazujące, w jaki sposób problem został rozwiązany.
IrishStat

Odpowiedzi:

15

1) Jeśli chodzi o twoje pierwsze pytanie, w literaturze opracowano i omówiono niektóre statystyki testów w celu przetestowania zerowej stacjonarności i zerowej wartości pierwiastka jednostkowego. Niektóre z wielu artykułów napisanych na ten temat to:

Związane z trendem:

  • Dickey, D. y Fuller, W. (1979a), Rozkład estymatorów dla autoregresywnych szeregów czasowych z pierwiastkiem jednostkowym, Journal of American Statistics Association 74, 427-31.
  • Dickey, D. y Fuller, W. (1981), Statystyka ilorazu wiarygodności dla szeregów czasowych autoregresji z pierwiastkiem jednostkowym, Econometrica 49, 1057-1071.
  • Kwiatkowski, D., Phillips, P., Schmidt, P. y Shin, Y. (1992), Testowanie zerowej hipotezy stacjonarności wobec alternatywy pierwiastka: Czy jesteśmy pewni, że ekonomiczne szeregi czasowe mają pierwiastek? , Journal of Econometrics 54, 159-178.
  • Phillips, P. y Perron, P. (1988), Testowanie pierwiastka jednostkowego w regresji szeregów czasowych, Biometrika 75, 335-46.
  • Durlauf, S. y Phillips, P. (1988), Trendy a losowe spacery w analizie szeregów czasowych, Econometrica 56, 1333-54.

Związane ze składnikiem sezonowym:

  • Hylleberg, S., Engle, R., Granger, C. y Yoo, B. (1990), Sezonowa integracja i kointegracja, Journal of Econometrics 44, 215-38.
  • Canova, F. y Hansen, BE (1995), Czy wzorce sezonowe są stałe w czasie? test stabilności sezonowej, Journal of Business and Economic Statistics 13, 237-252.
  • Franses, P. (1990), Testowanie sezonowych korzeni jednostkowych w danych miesięcznych, Raport techniczny 9032, Instytut Ekonometryczny.
  • Ghysels, E., Lee, H. y Noh, J. (1994), Testowanie pierwiastków jednostkowych w sezonowych szeregach czasowych. niektóre rozszerzenia teoretyczne i dochodzenie Monte Carlo, Journal of Econometrics 62, 415-442.

Podręcznik Banerjee, A., Dolado, J., Galbraith, J. y Hendry, D. (1993), Kointegracja, korekcja błędów oraz analiza ekonometryczna danych niestacjonarnych, Advanced Texts in Econometrics. Oxford University Press jest również dobrym punktem odniesienia.

2) Twoja druga troska jest uzasadniona literaturą. Jeśli istnieje test na pierwiastek jednostkowy, tradycyjna statystyka t, którą zastosowałbyś do trendu liniowego, nie jest zgodna z rozkładem standardowym. Patrz na przykład Phillips, P. (1987), Regresja szeregów czasowych z pierwiastkiem jednostkowym, Econometrica 55 (2), 277-301.

Jeśli pierwiastek jednostkowy istnieje i jest ignorowany, wówczas zmniejsza się prawdopodobieństwo odrzucenia wartości zerowej, że współczynnik trendu liniowego wynosi zero. Oznacza to, że zbyt często skończylibyśmy na modelowaniu deterministycznego trendu liniowego dla danego poziomu istotności. W obecności korzenia jednostki powinniśmy zamiast tego przekształcić dane, biorąc regularne różnice do danych.

3) Na przykład, jeśli używasz R, możesz wykonać następującą analizę z danymi.

x <- structure(c(7657, 5451, 10883, 9554, 9519, 10047, 10663, 10864, 
  11447, 12710, 15169, 16205, 14507, 15400, 16800, 19000, 20198, 
  18573, 19375, 21032, 23250, 25219, 28549, 29759, 28262, 28506, 
  33885, 34776, 35347, 34628, 33043, 30214, 31013, 31496, 34115, 
  33433, 34198, 35863, 37789, 34561, 36434, 34371, 33307, 33295, 
  36514, 36593, 38311, 42773, 45000, 46000, 42000, 47000, 47500, 
  48000, 48500, 47000, 48900), .Tsp = c(1, 57, 1), class = "ts")

Najpierw możesz zastosować test Dickeya-Fullera dla wartości zerowej katalogu głównego:

require(tseries)
adf.test(x, alternative = "explosive")
#   Augmented Dickey-Fuller Test
#   Dickey-Fuller = -2.0685, Lag order = 3, p-value = 0.453
#   alternative hypothesis: explosive

oraz test KPSS dla hipotezy odwrotnego zerowania, stacjonarność wobec alternatywy stacjonarności wokół trendu liniowego:

kpss.test(x, null = "Trend", lshort = TRUE)
#   KPSS Test for Trend Stationarity
#   KPSS Trend = 0.2691, Truncation lag parameter = 1, p-value = 0.01

Wyniki: test ADF, przy poziomie istotności 5% pierwiastek jednostkowy nie jest odrzucany; Test KPSS, zerowa stacjonarność jest odrzucana na rzecz modelu z trendem liniowym.

Na marginesie: użycie lshort=FALSEwartości zerowej testu KPSS nie jest odrzucane na poziomie 5%, jednak wybiera 5 opóźnień; kolejna kontrola, której nie pokazano, sugeruje, że wybór 1-3 opóźnień jest odpowiedni dla danych i prowadzi do odrzucenia hipotezy zerowej.

Zasadniczo powinniśmy kierować się testem, dla którego byliśmy w stanie odrzucić hipotezę zerową (zamiast testem, dla którego nie odrzuciliśmy (przyjęliśmy) zerową). Jednak regresja oryginalnej serii według trendu liniowego okazuje się niewiarygodna. Z jednej strony kwadrat R jest wysoki (ponad 90%), co jest wskazane w literaturze jako wskaźnik fałszywej regresji.

fit <- lm(x ~ 1 + poly(c(time(x))))
summary(fit)
#Coefficients:
#                 Estimate Std. Error t value Pr(>|t|)    
#(Intercept)       28499.3      381.6   74.69   <2e-16 ***
#poly(c(time(x)))  91387.5     2880.9   31.72   <2e-16 ***
#---
#Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
#Residual standard error: 2881 on 55 degrees of freedom
#Multiple R-squared:  0.9482,   Adjusted R-squared:  0.9472 
#F-statistic:  1006 on 1 and 55 DF,  p-value: < 2.2e-16

Z drugiej strony reszty są autokorelowane:

acf(residuals(fit)) # not displayed to save space

Co więcej, zerowy rdzeń jednostkowy w resztach nie może zostać odrzucony.

adf.test(residuals(fit))
#   Augmented Dickey-Fuller Test
#Dickey-Fuller = -2.0685, Lag order = 3, p-value = 0.547
#alternative hypothesis: stationary

W tym momencie możesz wybrać model, który będzie używany do uzyskiwania prognoz. Na przykład prognozy oparte na strukturalnym modelu szeregów czasowych i na modelu ARIMA można uzyskać w następujący sposób.

# StructTS
fit1 <- StructTS(x, type = "trend")
fit1
#Variances:
# level    slope  epsilon  
#2982955        0   487180 
# 
# forecasts
p1 <- predict(fit1, 10, main = "Local trend model")
p1$pred
# [1] 49466.53 50150.56 50834.59 51518.62 52202.65 52886.68 53570.70 54254.73
# [9] 54938.76 55622.79

# ARIMA
require(forecast)
fit2 <- auto.arima(x, ic="bic", allowdrift = TRUE)
fit2
#ARIMA(0,1,0) with drift         
#Coefficients:
#         drift
#      736.4821
#s.e.  267.0055
#sigma^2 estimated as 3992341:  log likelihood=-495.54
#AIC=995.09   AICc=995.31   BIC=999.14
#
# forecasts
p2 <- forecast(fit2, 10, main = "ARIMA model")
p2$mean
# [1] 49636.48 50372.96 51109.45 51845.93 52582.41 53318.89 54055.37 54791.86
# [9] 55528.34 56264.82

Fabuła prognoz:

par(mfrow = c(2, 1), mar = c(2.5,2.2,2,2))
plot((cbind(x, p1$pred)), plot.type = "single", type = "n", 
  ylim = range(c(x, p1$pred + 1.96 * p1$se)), main = "Local trend model")
grid()
lines(x)
lines(p1$pred, col = "blue")
lines(p1$pred + 1.96 * p1$se, col = "red", lty = 2)
lines(p1$pred - 1.96 * p1$se, col = "red", lty = 2)
legend("topleft", legend = c("forecasts", "95% confidence interval"), 
  lty = c(1,2), col = c("blue", "red"), bty = "n")
plot((cbind(x, p2$mean)), plot.type = "single", type = "n", 
  ylim = range(c(x, p2$upper)), main = "ARIMA (0,1,0) with drift")
grid()
lines(x)
lines(p2$mean, col = "blue")
lines(ts(p2$lower[,2], start = end(x)[1] + 1), col = "red", lty = 2)
lines(ts(p2$upper[,2], start = end(x)[1] + 1), col = "red", lty = 2)

prognozy trendów

Prognozy są podobne w obu przypadkach i wyglądają rozsądnie. Zauważ, że prognozy są zgodne ze stosunkowo deterministycznym wzorcem podobnym do trendu liniowego, ale nie modelowaliśmy wyraźnie trendu liniowego. Powód jest następujący: i) w lokalnym modelu trendu wariancja składnika nachylenia jest szacowana na zero. To zmienia składową trendu w dryf, który daje efekt trendu liniowego. ii) ARIMA (0,1,1), model z dryfem jest wybierany w modelu dla szeregu różnicowego. Wpływ terminu stałego na szereg różnicowy jest trendem liniowym. Jest to omówione w tym poście .

Możesz sprawdzić, czy jeśli wybierzesz model lokalny lub ARIMA (0,1,0) bez dryfu, prognozy będą prostą linią poziomą, a zatem nie będą podobne do obserwowanej dynamiki danych. Cóż, jest to część układanki testów na pierwiastki jednostkowe i komponentów deterministycznych.

Edycja 1 (kontrola pozostałości): autokorelacja i częściowy ACF nie sugerują struktury w pozostałościach.

resid1 <- residuals(fit1)
resid2 <- residuals(fit2)
par(mfrow = c(2, 2))
acf(resid1, lag.max = 20, main = "ACF residuals. Local trend model")
pacf(resid1, lag.max = 20, main = "PACF residuals. Local trend model")
acf(resid2, lag.max = 20, main = "ACF residuals. ARIMA(0,1,0) with drift")
pacf(resid2, lag.max = 20, main = "PACF residuals. ARIMA(0,1,0) with drift")

ACF-PACF

Jak sugeruje IrishStat, wskazane jest również sprawdzenie wartości odstających. Dwie dodatkowe wartości odstające są wykrywane za pomocą pakietu tsoutliers.

require(tsoutliers)
resol <- tsoutliers(x, types = c("AO", "LS", "TC"), 
  remove.method = "bottom-up", 
  args.tsmethod = list(ic="bic", allowdrift=TRUE))
resol
#ARIMA(0,1,0) with drift         
#Coefficients:
#         drift        AO2       AO51
#      736.4821  -3819.000  -4500.000
#s.e.  220.6171   1167.396   1167.397
#sigma^2 estimated as 2725622:  log likelihood=-485.05
#AIC=978.09   AICc=978.88   BIC=986.2
#Outliers:
#  type ind time coefhat  tstat
#1   AO   2    2   -3819 -3.271
#2   AO  51   51   -4500 -3.855

Patrząc na ACF, możemy powiedzieć, że na poziomie 5% istotności, reszty są również losowe w tym modelu.

par(mfrow = c(2, 1))
acf(residuals(resol$fit), lag.max = 20, main = "ACF residuals. ARIMA with additive outliers")
pacf(residuals(resol$fit), lag.max = 20, main = "PACF residuals. ARIMA with additive outliers")

wprowadź opis zdjęcia tutaj

W tym przypadku obecność potencjalnych wartości odstających nie wydaje się zakłócać działania modeli. Potwierdza to test Jarque-Bera na normalność; zero normalności w resztach z początkowych modeli ( fit1, fit2) nie jest odrzucane na poziomie istotności 5%.

jarque.bera.test(resid1)[[1]]
# X-squared = 0.3221, df = 2, p-value = 0.8513
jarque.bera.test(resid2)[[1]]
#X-squared = 0.426, df = 2, p-value = 0.8082

Edycja 2 (wykres reszt i ich wartości) Tak wyglądają resztki:

pozostałości

A oto ich wartości w formacie csv:

0;6.9205
-0.9571;-2942.4821
2.6108;4695.5179
-0.5453;-2065.4821
-0.2026;-771.4821
0.1242;-208.4821
0.1909;-120.4821
-0.0179;-535.4821
0.1449;-153.4821
0.484;526.5179
1.0748;1722.5179
0.3818;299.5179
-1.061;-2434.4821
0.0996;156.5179
0.4805;663.5179
0.8969;1463.5179
0.4111;461.5179
-1.0595;-2361.4821
0.0098;65.5179
0.5605;920.5179
0.8835;1481.5179
0.7669;1232.5179
1.4024;2593.5179
0.3785;473.5179
-1.1032;-2233.4821
-0.3813;-492.4821
2.2745;4642.5179
0.2935;154.5179
-0.1138;-165.4821
-0.8035;-1455.4821
-1.2982;-2321.4821
-1.9463;-3565.4821
-0.1648;62.5179
-0.1022;-253.4821
0.9755;1882.5179
-0.5662;-1418.4821
-0.0176;28.5179
0.5;928.5179
0.6831;1189.5179
-1.8889;-3964.4821
0.3896;1136.5179
-1.3113;-2799.4821
-0.9934;-1800.4821
-0.4085;-748.4821
1.2902;2482.5179
-0.0996;-657.4821
0.5539;981.5179
2.0007;3725.5179
1.0227;1490.5179
0.27;263.5179
-2.336;-4736.4821
1.8994;4263.5179
0.1301;-236.4821
-0.0892;-236.4821
-0.1148;-236.4821
-1.1207;-2236.4821
0.4801;1163.5179
javlacalle
źródło
1
Czy zweryfikowałeś, że reszty z twoich modeli były losowe, tj. Brak wartości odstających lub struktury ARIMA, które są wymagane do testu istotności szacowanych współczynników. Zauważ, że jeśli masz wartości odstające w resztkach, ACF jest bez znaczenia, ponieważ rozdęta wariancja błędu prowadzi do niedoszacowania ACF. Czy możesz podać wykresy błędów, które dowodzą / sugerują przypadkowość, w przeciwnym razie twoje wnioski na temat nieskorelowanych reszt mogą być fałszywe.
IrishStat
Zdecydowanie pełna analiza wymaga sprawdzenia pozostałości. Swoją odpowiedź ograniczyłem do podania niektórych narzędzi, których można użyć do przeprowadzenia testów, o które pytał „Forecaster” i zilustrował ich użycie. Cieszę się, że interesuje Cię więcej szczegółów, zredagowałem swoje pytanie.
javlacalle
Poprosiłem o wykres czasowy resztek. Czy możesz je podać, a także same resztki, abym mógł przetworzyć je za pomocą AUTOBOX, aby potwierdzić, że są wolne od widocznej struktury. Test JB nie jest preferowany podczas testowania pulsów, przesunięć poziomów, sezonowych pulsów i / lub trendów czasu lokalnego w zbiorze danych, chociaż obecność tego rodzaju struktur może spowodować odrzucenie założenia o normalności. Pomysł, że jeśli wartość zerowa nie zostanie odrzucona, jest dowodem jej akceptacji, może być niebezpieczny. Proszę zobaczyć unc.edu/~jbhill/tsay.pdf
IrishStat
1
Dzięki. Złożyłem 57 odpowiedzi, a 5 z nich zostało wstępnie oznaczonych jako wyjątkowe. W kolejności ważności (51,3,27,52 i 48). Twój wykres wizualnie obsługuje te punkty. Wynikłe błędy nie wykazują naruszenia losowości, a następnie nie mają znaczącego ACF. Aby skorygować obserwowane wartości w celu dostosowania do wykrywania anomalii, użyj następujących poleceń: + [X1 (T)] [(- 4494.5)]: PULSE 51 + [X2 (T)] [(+ 4937.5)]: PULSE 3 + [X3 (T)] [(+ 4884.5)]: PULSE 27 + [X4 (T)] [(+ 4505,5)]: PULSE 52 + [X5 (T)] [(+ 3967.5)]: PULSE 48
IrishStat
1
@B_Miner Zwykle zaczniesz od spojrzenia na funkcję autokorelacji reszt. Jeśli autokorelacje są znaczące i duże dla dużych zamówień (tj. ACF nie rozkłada się wykładniczo do zera), możesz rozważyć zastosowanie jednostkowego testu pierwiastkowego dla reszt. Jeśli analiza reszt wskazuje, że istnieje pierwiastek jednostkowy, oznacza to, że prawdopodobnie powinieneś najpierw wziąć dwie różnice w stosunku do oryginalnych danych (tj. Wziąć ponownie różnice w szeregach różnicowych).
javlacalle
4

W odniesieniu do danych niesezonowych ... Trendy mogą mieć dwie formy y (t) = y (t − 1) + θ0 (A) Trend stochastyczny lub Y (t) = a + bx1 + cx2 (B) Deterministyczny Trend itp. Gdzie x1 = 1,2,3,4 .... ti x2 = 0,0,0,0,0,1,2,3,4 więc jeden trend dotyczy obserwacji 1 − t, a drugi trend dotyczy obserwacji 6 do t.

Twoja niesezonowa seria zawierała 29 wartości. Użyłem AUTOBOX oprogramowania, które pomogłem opracować w sposób całkowicie automatyczny. AUTOBOX jest przejrzystą procedurą, ponieważ szczegółowo opisuje każdy etap procesu modelowania. Tutaj przedstawiono wykres serii / dopasowanych wartości / prognoz wprowadź opis zdjęcia tutaj. Użycie AUTOBOX do utworzenia modelu typu A doprowadziło do następujących wniosków wprowadź opis zdjęcia tutaj. Równanie jest ponownie przedstawione tutaj wprowadź opis zdjęcia tutaj, Statystyka modelu wprowadź opis zdjęcia tutaj. Wykres reszt znajduje się tutaj, wprowadź opis zdjęcia tutaja tabela prognozowanych wartości jest tutaj wprowadź opis zdjęcia tutaj. Ograniczenie AUTOBOX do modelu typu B doprowadziło do wykrycia przez AUTOBOX podwyższonego trendu w okresie 14: wprowadź opis zdjęcia tutaj wprowadź opis zdjęcia tutajwprowadź opis zdjęcia tutaj!wprowadź opis zdjęcia tutajwprowadź opis zdjęcia tutajwprowadź opis zdjęcia tutaj

Jeśli chodzi o porównywanie modeli: Ponieważ liczba dopasowanych obserwacji jest różna (odpowiednio 26 i 29), nie jest możliwe użycie standardowych wskaźników (tj. R-kwadrat, standardowy błąd dev, AIC itp.) W celu ustalenia dominacji, chociaż w tym przypadku skinienie głową przejdź do A. Resztki z A są lepsze ze względu na strukturę AR (2). Prognozy z B są nieco agresywne, podczas gdy wzór prognoz A jest bardziej intuicyjny. Można powstrzymać się od powiedzenia 4 obserwacji i ocenić dokładność prognozy dla prognozy na 1 okres z 4 różnych źródeł (25, 26, 27 i 28).

IrishStat
źródło
2Irlandzka statystyka zapewnia doskonałą odpowiedź. Czytałem niektóre, że połączylibyśmy trendy stochastyczne i deterministyczne, czyli yt = y (t-1) + a + bt = ct? czy byłoby to pomocne
przepowiednia
Modelowa forma y (t) = B0 + B1 * t + a (t) [thetha / phi] załamuje się, jeśli phi jest powiedziane [1-B], ponieważ usuwanie ułamków zasadniczo różnicuje zmienną t dając stałą zderzającą się z B0. Innymi słowy, struktura ARIMA połączona ze wskaźnikami czasu może spowodować spustoszenie. Podany model można oszacować, ale zdecydowanie nie jest to preferowane podejście (być może brak endogenności!). Ktoś inny, kto to czyta, może komentować, może w tym pomóc. Nie jest to właściwy podzbiór funkcji przenoszenia i.imgur.com/dv4bAts.png
IrishStat