Jak poradzić sobie z nieistniejącymi lub brakującymi danymi?

Próbowałem metody prognozowania i chcę sprawdzić, czy moja metoda jest poprawna, czy nie.

Moje badanie porównuje różne rodzaje funduszy wspólnego inwestowania. Chcę użyć indeksu GCC jako punktu odniesienia dla jednego z nich, ale problem polega na tym, że indeks GCC zatrzymał się we wrześniu 2011 r., A moje badanie trwa od stycznia 2003 r. Do lipca 2014 r. Próbowałem więc użyć innego indeksu, indeksu MSCI, zrobić regresję liniową, ale problem polega na tym, że w indeksie MSCI brakuje danych z września 2010 r.

Aby obejść ten problem, wykonałem następujące czynności. Czy te kroki są prawidłowe?

1. W indeksie MSCI brakuje danych za okres od września 2010 r. Do lipca 2012 r. „Podałem go”, stosując średnie ruchome dla pięciu obserwacji. Czy to podejście jest prawidłowe? Jeśli tak, ile obserwacji powinienem zastosować?

2. Po oszacowaniu brakujących danych przeprowadziłem regresję indeksu GCC (jako zmiennej zależnej) w stosunku do indeksu MSCI (jako zmiennej niezależnej) dla wzajemnie dostępnego okresu (od stycznia 2007 r. Do września 2011 r.), A następnie poprawiłem model ze wszystkich problemów. Dla każdego miesiąca x zastępuję danymi z indeksu MSCI dla okresu odpoczynku. Czy to jest ważne?

Poniżej znajdują się dane w formacie wartości rozdzielanych przecinkami, zawierające lata według wierszy i miesiące według kolumn. Dane są również dostępne poprzez ten link .

Seria GCC:

,Jan,Feb,Mar,Apr,May,Jun,Jul,Aug,Sep,Oct,Nov,Dec
2002,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,117.709
2003,120.176,117.983,120.913,134.036,145.829,143.108,149.712,156.997,162.158,158.526,166.42,180.306
2004,185.367,185.604,200.433,218.923,226.493,230.492,249.953,262.295,275.088,295.005,328.197,336.817
2005,346.721,363.919,423.232,492.508,519.074,605.804,581.975,676.021,692.077,761.837,863.65,844.865
2006,947.402,993.004,909.894,732.646,598.877,686.258,634.835,658.295,672.233,677.234,491.163,488.911
2007,440.237,486.828,456.164,452.141,495.19,473.926,492.782,525.295,519.081,575.744,599.984,668.192
2008,626.203,681.292,616.841,676.242,657.467,654.66,635.478,603.639,527.326,396.904,338.696,308.085
2009,279.706,252.054,272.082,314.367,340.354,325.99,326.46,327.053,354.192,339.035,329.668,318.267
2010,309.847,321.98,345.594,335.045,311.363,299.555,310.802,306.523,315.496,324.153,323.256,334.802
2011,331.133,311.292,323.08,327.105,320.258,312.749,305.073,297.087,298.671,NA,NA,NA

Seria MSCI:

,Jan,Feb,Mar,Apr,May,Jun,Jul,Aug,Sep,Oct,Nov,Dec
2007,NA,NA,NA,NA,1000,958.645,1016.085,1049.468,1033.775,1118.854,1142.347,1298.223
2008,1197.656,1282.557,1164.874,1248.42,1227.061,1221.049,1161.246,1112.582,929.379,680.086,516.511,521.127
2009,487.562,450.331,478.255,560.667,605.143,598.611,609.559,615.73,662.891,655.639,628.404,602.14
2010,601.1,622.624,661.875,644.751,588.526,587.4,615.008,606.133,NA,NA,NA,NA
2011,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA
2012,NA,NA,NA,NA,NA,NA,NA,609.51,598.428,595.622,582.905,599.447
2013,627.561,619.581,636.284,632.099,651.995,651.39,687.194,676.76,694.575,704.806,727.625,739.842
2014,759.036,787.057,817.067,824.313,857.055,805.31,873.619,NA,NA,NA,NA,NA

Moja propozycja jest podobna do tego, co proponujecie, z tym wyjątkiem, że zamiast średnich ruchomych użyłbym modelu szeregów czasowych. Ramy modeli ARIMA są również odpowiednie do uzyskiwania prognozy obejmującej nie tylko serię MSCI jako regresor, ale także opóźnienia serii GCC, które mogą również wychwytywać dynamikę danych.

Po pierwsze, możesz dopasować model ARIMA do serii MSCI i interpolować brakujące obserwacje w tej serii. Następnie możesz dopasować model ARIMA dla serii GCC, używając MSCI jako egzogennych regresorów i uzyskać prognozy dla GCC na podstawie tego modelu. Robiąc to, musisz ostrożnie radzić sobie z przerwami, które są graficznie obserwowane w serii i które mogą zniekształcić wybór i dopasowanie modelu ARIMA.

Oto, co robię w tej analizie R. Korzystam z tej funkcji, forecast::auto.arimaaby dokonać wyboru modelu ARIMA i tsoutliers::tsowykryć możliwe przesunięcia poziomu (LS), zmiany tymczasowe (TC) lub wartości odstające addytywne (AO).

Są to dane po załadowaniu:

gcc <- structure(c(117.709, 120.176, 117.983, 120.913, 134.036, 145.829, 143.108, 149.712, 156.997, 162.158, 158.526, 166.42, 180.306, 185.367, 185.604, 200.433, 218.923, 226.493, 230.492, 249.953, 262.295, 275.088, 295.005, 328.197, 336.817, 346.721, 363.919, 423.232, 492.508, 519.074, 605.804, 581.975, 676.021, 692.077, 761.837, 863.65, 844.865, 947.402, 993.004, 909.894, 732.646, 598.877, 686.258, 634.835, 658.295, 672.233, 677.234, 491.163, 488.911, 440.237, 486.828, 456.164, 452.141, 495.19, 473.926, 
492.782, 525.295, 519.081, 575.744, 599.984, 668.192, 626.203, 681.292, 616.841, 676.242, 657.467, 654.66, 635.478, 603.639, 527.326, 396.904, 338.696, 308.085, 279.706, 252.054, 272.082, 314.367, 340.354, 325.99, 326.46, 327.053, 354.192, 339.035, 329.668, 318.267, 309.847, 321.98, 345.594, 335.045, 311.363, 
299.555, 310.802, 306.523, 315.496, 324.153, 323.256, 334.802, 331.133, 311.292, 323.08, 327.105, 320.258, 312.749, 305.073, 297.087, 298.671), .Tsp = c(2002.91666666667, 2011.66666666667, 12), class = "ts")
msci <- structure(c(1000, 958.645, 1016.085, 1049.468, 1033.775, 1118.854, 1142.347, 1298.223, 1197.656, 1282.557, 1164.874, 1248.42, 1227.061, 1221.049, 1161.246, 1112.582, 929.379, 680.086, 516.511, 521.127, 487.562, 450.331, 478.255, 560.667, 605.143, 598.611, 609.559, 615.73, 662.891, 655.639, 628.404, 602.14, 601.1, 622.624, 661.875, 644.751, 588.526, 587.4, 615.008, 606.133, 
NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, 609.51, 598.428, 595.622, 582.905, 599.447, 627.561, 619.581, 636.284, 632.099, 651.995, 651.39, 687.194, 676.76, 694.575, 704.806, 727.625, 739.842, 759.036, 787.057, 817.067, 824.313, 857.055, 805.31, 873.619), .Tsp = c(2007.33333333333, 2014.5, 12), class = "ts")

Krok 1: Dopasuj model ARIMA do serii MSCI

Pomimo tego, że grafika ujawnia obecność niektórych przerw, nie wykryto żadnych wartości odstających tso. Może to wynikać z faktu, że w środku próby brakuje kilku obserwacji. Możemy sobie z tym poradzić w dwóch krokach. Najpierw dopasuj model ARIMA i użyj go do interpolacji brakujących obserwacji; po drugie, dopasuj model ARIMA do interpolowanej serii, sprawdzając możliwe LS, TC, AO i popraw interpolowane wartości, jeśli zostaną znalezione zmiany.

Wybierz model ARIMA dla serii MSCI:

require("forecast")
fit1 <- auto.arima(msci)
fit1
# ARIMA(1,1,2) with drift         
# Coefficients:
#           ar1     ma1     ma2    drift
#       -0.6935  1.1286  0.7906  -1.4606
# s.e.   0.1204  0.1040  0.1059   9.2071
# sigma^2 estimated as 2482:  log likelihood=-328.05
# AIC=666.11   AICc=666.86   BIC=678.38

Uzupełnij brakujące uwagi zgodnie z podejściem omówionym w mojej odpowiedzi na ten post :

kr <- KalmanSmooth(msci, fit1$model)
tmp <- which(fit1$model$Z == 1)
id <- ifelse (length(tmp) == 1, tmp[1], tmp[2])
id.na <- which(is.na(msci))
msci.filled <- msci
msci.filled[id.na] <- kr$smooth[id.na,id]

Dopasuj model ARIMA do wypełnionej serii msci.filled. Teraz znaleziono niektóre wartości odstające. Niemniej jednak, stosując alternatywne opcje, wykryto różne wartości odstające. Zatrzymam ten, który został znaleziony w większości przypadków, zmianę poziomu w październiku 2008 r. (Obserwacja 18). Możesz spróbować na przykład tych i innych opcji.

require("tsoutliers")
tso(msci.filled, remove.method = "bottom-up", tsmethod = "arima", 
  args.tsmethod = list(order = c(1,1,1)))
tso(msci.filled, remove.method = "bottom-up", args.tsmethod = list(ic = "bic"))

Wybrany model to teraz:

mo <- outliers("LS", 18)
ls <- outliers.effects(mo, length(msci))
fit2 <- auto.arima(msci, xreg = ls)
fit2
# ARIMA(2,1,0)                    
# Coefficients:
#           ar1     ar2       LS18
#       -0.1006  0.4857  -246.5287
# s.e.   0.1139  0.1093    45.3951
# sigma^2 estimated as 2127:  log likelihood=-321.78
# AIC=651.57   AICc=652.06   BIC=661.39

Użyj poprzedniego modelu, aby poprawić interpolację brakujących obserwacji:

kr <- KalmanSmooth(msci, fit2$model)
tmp <- which(fit2$model$Z == 1)
id <- ifelse (length(tmp) == 1, tmp[1], tmp[2])
msci.filled2 <- msci
msci.filled2[id.na] <- kr$smooth[id.na,id]

Początkowe i końcowe interpolacje można porównać na wykresie (nie pokazano tutaj, aby zaoszczędzić miejsce):

plot(msci.filled, col = "gray")
lines(msci.filled2)

Krok 2: Dopasuj model ARIMA do GCC, używając msci.filled2 jako egzogennego regresora

Ignoruję brakujące obserwacje na początku msci.filled2. W tym momencie znalazłem pewne trudności w użyciu auto.arimawraz z tso, więc próbowałem ręcznie kilka modeli ARIMA w tsoi ostatecznie wybrał ARIMA (1,1,0).

xreg <- window(cbind(gcc, msci.filled2)[,2], end = end(gcc))
fit3 <- tso(gcc, remove.method = "bottom-up", tsmethod = "arima",  
  args.tsmethod = list(order = c(1,1,0), xreg = data.frame(msci=xreg)))
fit3
# ARIMA(1,1,0)                    
# Coefficients:
#           ar1    msci     AO72
#       -0.1701  0.5131  30.2092
# s.e.   0.1377  0.0173   6.7387
# sigma^2 estimated as 71.1:  log likelihood=-180.62
# AIC=369.24   AICc=369.64   BIC=379.85
# Outliers:
#   type ind    time coefhat tstat
# 1   AO  72 2008:11   30.21 4.483

Fabuła GCC pokazuje zmianę na początku 2008 r. Wydaje się jednak, że została już uchwycona przez regresor MSCI i nie uwzględniono żadnych dodatkowych regresorów, z wyjątkiem dodatkowej wartości odstającej w listopadzie 2008 r.

Wykres reszt nie sugerował żadnej struktury autokorelacji, ale wykres sugerował przesunięcie poziomu w listopadzie 2008 r. I wartość odstającą dodatku w lutym 2011 r. Jednak po dodaniu odpowiednich interwencji diagnostyka modelu była gorsza. W tym momencie może być konieczna dalsza analiza. Tutaj będę nadal uzyskiwać prognozy na podstawie ostatniego modelu fit3.

$95\%$

newxreg <- data.frame(msci=window(msci.filled2, start = c(2011, 10)), AO72=rep(0, 34))
p <- predict(fit3$fit, n.ahead = 34, newxreg = newxreg)
head(p$pred)
# [1] 298.3544 298.2753 298.0958 298.0641 297.6829 297.7412
par(mar = c(3,3.5,2.5,2), las = 1)
plot(cbind(gcc, msci), xaxt = "n", xlab = "", ylab = "", plot.type = "single", type = "n")
grid()
lines(gcc, col = "blue", lwd = 2)
lines(msci, col = "green3", lwd = 2)
lines(window(msci.filled2, start = c(2010, 9), end = c(2012, 7)), col = "green", lwd = 2)
lines(p$pred, col = "red", lwd = 2)
lines(p$pred + 1.96 * p$se, col = "red", lty = 2)
lines(p$pred - 1.96 * p$se, col = "red", lty = 2)
xaxis1 <- seq(2003, 2014)
axis(side = 1, at = xaxis1, labels = xaxis1)
legend("topleft", col = c("blue", "green3", "green", "red", "red"), lwd = 2, bty = "n", lty = c(1,1,1,1,2), legend = c("GCC", "MSCI", "Interpolated values", "Forecasts", "95% confidence interval"))

1. Imputacja (tj. „Zapewniona” przez średnią ruchomą) jest ważna, jeśli wartości losowo brakuje. Jeśli jest to nieprzerwany okres o znacznej długości, staje się to mało prawdopodobne. Druga część pytania jest niejasna.
2. W zależności od pytania uważa się, że użycie modelu do prognozowania wykraczającego poza zakres danych uważa się za coś od optymalnego do niepoprawnego: np. Co, jeśli relacja między tymi dwoma indeksami zmieni się w latach 2012-2014? Możesz użyć wartości szacowanych metodą regresji (ale nie zastępować bezpośrednio surowymi wartościami innego indeksu) dla brakujących punktów danych, ale miałoby to sens tylko wtedy, gdy istnieje silny związek między tymi dwoma indeksami i jest bardzo ważne, aby te wartości były wyraźnie oznaczone jako dane szacunkowe. Co rozumiesz przez „poprawienie modelu ze wszystkich problemów”?

2 Wydaje się w porządku. Poszedłbym z tym.

Jeśli chodzi o 1. Sugeruję wytrenowanie modelu do przewidywania GCC z wykorzystaniem wszystkich funkcji dostępnych w zbiorze danych (które nie są NA w okresie od września 2011 r.) (Pomiń wiersze, które mają dowolną wartość NA przed sep2011 podczas treningu). Model powinien być bardzo dobry (użyj walidacji krzyżowania K-fold). Teraz przewidujemy GCC na okres od września 2011 r. I później.

Alternatywnie możesz wytrenować model, który przewiduje MSCI, użyj go do przewidywania brakujących wartości MSCI. Teraz wytrenuj model przewidywania GCC za pomocą MSCI, a następnie przewiduj GCC na okres od września 2011 r. I później