Mam comiesięczny szereg czasowy z interwencją i chciałbym oszacować wpływ tej interwencji na wynik. Zdaję sobie sprawę, że seria jest raczej krótka, a efekt nie jest jeszcze zakończony.
Dane
cds <- structure(c(2580L, 2263L, 3679L, 3461L, 3645L, 3716L, 3955L, 3362L,
2637L, 2524L, 2084L, 2031L, 2256L, 2401L, 3253L, 2881L,
2555L, 2585L, 3015L, 2608L, 3676L, 5763L, 4626L, 3848L,
4523L, 4186L, 4070L, 4000L, 3498L),
.Dim=c(29L, 1L),
.Dimnames=list(NULL, "CD"),
.Tsp=c(2012, 2014.33333333333, 12), class="ts")
Metodologia
1) W przypadku tej auto.arima
funkcji zastosowano serię przedinterwencyjną (do października 2013 r.) . Sugerowanym modelem był ARIMA (1,0,0) ze średnią niezerową. Fabuła ACF wyglądała dobrze.
pre <- window(cds, start=c(2012, 01), end=c(2013, 09))
mod.pre <- auto.arima(log(pre))
# Coefficients:
# ar1 intercept
# 0.5821 7.9652
# s.e. 0.1763 0.0810
#
# sigma^2 estimated as 0.02709: log likelihood=7.89
# AIC=-9.77 AICc=-8.36 BIC=-6.64
2) Biorąc pod uwagę wykres pełnej serii, odpowiedź impulsowa została wybrana poniżej, przy T = październik 2013,
które według cryer i chan można dopasować w następujący sposób do funkcji arimax:
mod.arimax <- arimax(log(cds), order=c(1, 0, 0),
seasonal=list(order=c(0, 0, 0), frequency=12),
include.mean=TRUE,
xtransf=data.frame(Oct13=1 * (seq(cds) == 22)),
transfer=list(c(1, 1)))
mod.arimax
# Series: log(cds)
# ARIMA(1,0,0) with non-zero mean
#
# Coefficients:
# ar1 intercept Oct13-AR1 Oct13-MA0 Oct13-MA1
# 0.7619 8.0345 -0.4429 0.4261 0.3567
# s.e. 0.1206 0.1090 0.3993 0.1340 0.1557
#
# sigma^2 estimated as 0.02289: log likelihood=12.71
# AIC=-15.42 AICc=-11.61 BIC=-7.22
Resztki z tego wydawały się OK:
Fabuła wyposażenia i aktualności:
plot(fitted(mod.arimax), col="red", type="b")
lines(window(log(cds), start=c(2012, 02)), type="b")
Pytania
1) Czy ta metodologia jest odpowiednia do analizy interwencji?
2) Czy mogę spojrzeć na oszacowanie / SE dla składników funkcji przeniesienia i stwierdzić, że efekt interwencji był znaczący?
3) Jak wizualizować efekt funkcji przenoszenia (wykreślić?)
4) Czy istnieje sposób oszacowania, o ile interwencja zwiększyła produkcję po miesiącach „x”? Chyba do tego (a może # 3) pytam, jak pracować z równaniem modelu - gdyby to była prosta regresja liniowa ze zmiennymi fikcyjnymi (na przykład), mogłabym uruchamiać scenariusze z interwencją i bez i mierzyć wpływ - ale nie jestem pewien, jak pracować z tego typu modelem.
DODAJ
Na żądanie przedstawiamy resztki z dwóch parametryzacji.
Najpierw od dopasowania:
fit <- arimax(log(cds), order=c(1, 0, 0),
xtransf=
data.frame(Oct13a=1 * (seq_along(cds) == 22),
Oct13b=1 * (seq_along(cds) == 22)),
transfer=list(c(0, 0), c(1, 0)))
plot(resid(fit), type="b")
Następnie z tego dopasowania
mod.arimax <- arimax(log(cds), order=c(1, 0, 0),
seasonal=list(order=c(0, 0, 0), frequency=12),
include.mean=TRUE,
xtransf=data.frame(Oct13=1 * (seq(cds) == 22)),
transfer=list(c(1, 1)))
mod.arimax
plot(resid(mod.arimax), type="b")
źródło
Odpowiedzi:
Model AR (1) z interwencją zdefiniowaną w równaniu podanym w pytaniu można dopasować, jak pokazano poniżej. Zauważ, jak
transfer
zdefiniowano argument ; potrzebujesz również jednej zmiennej wskaźnikowejxtransf
dla każdej interwencji (puls i zmiana przejściowa):coeftest
Efekt interwencji można określić ilościowo w następujący sposób:
Możesz wykreślić efekt interwencji w następujący sposób:
Liczbowo są to szacunkowe wzrosty określone ilościowo w każdym punkcie czasowym spowodowane interwencją w październiku 2013 r .:
Interwencja zwiększa wartość obserwowanej zmiennej w październiku 2013 r. O około75 %
Możemy również tworzyć interwencje ręcznie i przekazywać je0,9231
stats::arima
jako zewnętrzne regresory. Interwencje to puls plus przejściowa zmiana z parametremxreg
Interwencje te są równoważne z wartościami odstającymi dodatków (AO) i zmianami przejściowymi (TC) zdefiniowanymi w opakowaniu
tsoutliers
. Możesz użyć tego pakietu do wykrycia tych efektów, jak pokazano w odpowiedzi @forecaster lub do zbudowania używanych wcześniej regresorów. Na przykład w tym przypadku:Edytuj 1
Widziałem, że podane równanie można przepisać jako:
i można to określić tak jak za pomocą
transfer=list(c(1, 1))
.Jak pokazano poniżej, ta parametryzacja prowadzi w tym przypadku do oszacowania parametrów, które wymagają innego efektu niż poprzednia parametryzacja. Przypomina mi to efekt innowacyjnej wartości odstającej, a nie impulsu plus przejściowej zmiany.
Nie znam się na oznaczeniu opakowania,
TSA
ale myślę, że efekt interwencji można teraz określić ilościowo w następujący sposób:Efekt ten można teraz opisać jako gwałtowny wzrost w październiku 2013 r., Po którym następuje spadek w przeciwnym kierunku; wtedy efekt interwencji zanika szybko na przemian pozytywne i negatywne skutki niszczenia masy.
Ten efekt jest nieco szczególny, ale może być możliwy w rzeczywistych danych. W tym momencie przyjrzałbym się kontekstowi twoich danych i wydarzeniom, które mogły mieć wpływ na dane. Na przykład, czy nastąpiła zmiana zasad, kampania marketingowa, odkrycie ..., które mogą wyjaśnić interwencję w październiku 2013 r. Jeśli tak, czy bardziej sensowne jest, aby to wydarzenie miało wpływ na dane, jak opisano wcześniej lub jak stwierdziliśmy z początkową parametryzacją?
Edytuj 2
Prognozy można uzyskać i wyświetlić w następujący sposób:
Pierwsze prognozy odpowiadają stosunkowo dobrze obserwowanym wartościom (szara linia przerywana). Pozostałe prognozy pokazują, w jaki sposób seria będzie kontynuować ścieżkę do pierwotnej średniej. Przedziały ufności są jednak duże, co odzwierciedla niepewność. Dlatego powinniśmy być ostrożni i zmieniać model w miarę rejestrowania nowych danych.
źródło
Czasami mniej znaczy więcej. Mając 30 obserwacji w ręku, przesłałem dane do AUTOBOX, oprogramowania, które pomogłem opracować. Przesyłam następującą analizę w nadziei na zdobycie nagrody +200 (tylko żartuję!). Naszkicowałem rzeczywiste i oczyszczone wartości, sugerując wizualnie wpływ „ostatniej aktywności”. . Model, który został opracowany automatycznie, pokazano tutaj. i tutaj . Resztki tej raczej prostej serii z przesuniętym poziomem są tutaj przedstawione . Statystyki modelu są tutaj . Podsumowując, były interwencje, które można zidentyfikować empirycznie, czyniąc proces ARIMA; dwa impulsy i 1 zmiana poziomu . Wykres Rzeczywiste / Dopasowanie i Prognoza dodatkowo podkreśla analizę.
Chciałbym na przykład zobaczyć wykres reszt z wcześniej określonych i, moim zdaniem, potencjalnie nadmiernie określonych modeli.
źródło
Poniżej znajduje się kod:
Poniżej znajduje się szacunek, w październiku 2013 r. Nastąpił wzrost o ~ 236,3 jednostki przy standardowym błędzie ~ 481,8, a następnie ma on rozkład niszczący. Funkcja automatycznie identyfikuje AR (1). Musiałem wykonać kilka iteracji i wprowadzić różnicowanie zarówno sezonowe, jak i niesezonowe do 0, co znajduje odzwierciedlenie w args.tsmethod w funkcji tso.
Poniżej znajduje się fabuła, tsoutlier to jedyny znany mi pakiet, który może ładnie wydrukować tymczasowe zmiany na wykresie.
Ta analiza, miejmy nadzieję, dostarczyła odpowiedzi na twoje 2, 3 i 4 pytania, aczkolwiek przy użyciu innej metodologii. Zwłaszcza wykres i współczynniki zapewniły efekt tej interwencji oraz to, co by się stało, gdyby nie interwencja.
Mam również nadzieję, że ktoś inny może powielić ten wykres / analizę przy użyciu modelowania funkcji przenoszenia w R. Nie jestem pewien, czy można to zrobić w R, może ktoś inny może mnie sprawdzić.
źródło