Como lidar com dados inexistentes ou ausentes?

Tentei um método de previsão e quero verificar se meu método está correto ou não.

Meu estudo está comparando diferentes tipos de fundos mútuos. Quero usar o índice GCC como referência para um deles, mas o problema é que o índice GCC parou em setembro de 2011 e meu estudo é de janeiro de 2003 a julho de 2014. Portanto, tentei usar outro índice, o índice MSCI, para fazer uma regressão linear, mas o problema é que faltam dados do índice MSCI de setembro de 2010.

Para contornar isso, fiz o seguinte. Essas etapas são válidas?

1. Faltam dados para o índice MSCI de setembro de 2010 a julho de 2012. Eu o "forneci" aplicando médias móveis para cinco observações. Essa abordagem é válida? Em caso afirmativo, quantas observações devo usar?

2. Após estimar os dados ausentes, realizei uma regressão no índice GCC (como variável dependente) versus o índice MSCI (como variável independente) para o período mutuamente disponível (de janeiro de 2007 a setembro de 2011) e, em seguida, corrigi o modelo de todos os problemas. Para cada mês, substituo x pelos dados do índice MSCI para o período restante. Isso é válido?

Abaixo estão os dados no formato Valores separados por vírgulas, contendo os anos por linhas e os meses por colunas. Os dados também estão disponíveis neste link .

Série 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

Série 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

Minha sugestão é semelhante à que você propõe, exceto que eu usaria um modelo de série temporal em vez de médias móveis. A estrutura dos modelos ARIMA também é adequada para obter previsões, incluindo não apenas a série MSCI como regressor, mas também as séries GCC que também podem capturar a dinâmica dos dados.

Primeiro, você pode ajustar um modelo ARIMA para a série MSCI e interpolar as observações ausentes nesta série. Em seguida, você pode ajustar um modelo ARIMA para a série GCC usando o MSCI como regressores exógenos e obter as previsões para o GCC com base nesse modelo. Ao fazer isso, você deve ter cuidado ao lidar com as quebras que são observadas graficamente na série e que podem distorcer a seleção e o ajuste do modelo ARIMA.

Aqui está o que eu faço nessa análise R. Utilizo a função forecast::auto.arimapara fazer a seleção do modelo ARIMA e tsoutliers::tsodetectar possíveis mudanças de nível (LS), alterações temporárias (TC) ou outliers aditivos (AO).

Estes são os dados depois de carregados:

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")

Etapa 1: Ajustar um modelo ARIMA à série MSCI

Apesar do gráfico revelar a presença de algumas quebras, não foram detectados outliers tso. Isso pode ser devido ao fato de que existem várias observações ausentes no meio da amostra. Podemos lidar com isso em duas etapas. Primeiro, ajuste um modelo ARIMA e use-o para interpolar as observações ausentes; segundo, ajuste um modelo ARIMA para a série interpolada, verificando possíveis LS, TC, AO e refine os valores interpolados se forem encontradas alterações.

Escolha o modelo ARIMA para a série 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

Preencha as observações ausentes seguindo a abordagem discutida na minha resposta a este 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]

Ajuste um modelo ARIMA à série preenchida msci.filled. Agora, alguns outliers são encontrados. No entanto, usando opções alternativas, foram detectados diferentes outliers. Vou manter o que foi encontrado na maioria dos casos, uma mudança de nível em outubro de 2008 (observação 18). Você pode tentar, por exemplo, essas e outras opções.

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"))

O modelo escolhido é agora:

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

Use o modelo anterior para refinar a interpolação de observações ausentes:

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]

As interpolações inicial e final podem ser comparadas em um gráfico (não mostrado aqui para economizar espaço):

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

Etapa 2: ajuste um modelo ARIMA ao GCC usando msci.filled2 como regressor exógeno

Eu ignoro as observações ausentes no começo de msci.filled2. Neste ponto, eu encontrei algumas dificuldades para usar auto.arimajunto com tso, então eu tentei com a mão vários modelos ARIMA em tsoe, finalmente, escolheu o 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

O gráfico do GCC mostra uma mudança no início de 2008. No entanto, parece que ele já foi capturado pelo regressor MSCI e nenhum regressor adicional foi incluído, exceto um aditivo outlier em novembro de 2008.

O gráfico dos resíduos não sugeriu nenhuma estrutura de autocorrelação, mas sugeriu uma mudança de nível em novembro de 2008 e um aditivo em fevereiro de 2011. No entanto, ao adicionar as intervenções correspondentes, o diagnóstico do modelo foi pior. Mais análises podem ser necessárias neste momento. Aqui, continuarei obtendo as previsões com base no último modelo 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. A imputação (ou seja, 'fornecida' pela média móvel ') é válida se os valores estiverem ausentes aleatoriamente. Se for um período ininterrupto de duração considerável, isso se torna improvável. A segunda parte da pergunta não é clara.
2. Dependendo da pergunta, considera-se algo entre subótimo e inválido usar seu modelo para previsão além do escopo de seus dados: por exemplo, e se a relação entre os dois índices mudar em 2012-2014? Você pode usar valores estimados por regressão (mas não substituindo diretamente pelos valores brutos de outro índice) para os pontos de dados ausentes, mas isso só faria sentido se houver uma forte relação entre os dois índices e é fundamental que esses valores sejam claramente marcado como estimativas. E o que você quer dizer com "corrigiu o modelo de todos os problemas"?

2 Parece bom. Eu iria com isso.

Quanto ao 1. Eu sugiro que você treine um modelo para prever o GCC usando todos os recursos disponíveis no conjunto de dados (que não são NA durante o período de setembro de 2011 em diante) (omita as linhas que tenham qualquer valor NA antes de setembro de 2011 durante o treinamento). O modelo deve ser muito bom (use a validação cruzada K-fold). Agora preveja o CCG para o período de setembro de 2011 em diante.

Como alternativa, você pode treinar um modelo que prevê o MSCI, usá-lo para prever os valores ausentes do MSCI. Agora, treine um modelo para prever o GCC usando o MSCI e, em seguida, preveja o GCC para o período de setembro de 2011 em diante