{ "cells": [ { "cell_type": "code", "execution_count": 1, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "# https://cran.r-project.org/web/views/TimeSeries.html\n", "# https://www.datacamp.com/tracks/time-series-with-r\n", "# https://www.datascience.com/blog/introduction-to-forecasting-with-arima-in-r-learn-data-science-tutorials\n", "# https://www.otexts.org/fpp/8/7" ] }, { "cell_type": "code", "execution_count": 2, "metadata": { "collapsed": false }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Loading required package: fit.models\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "\n", "Attaching package: ‘tseries’\n", "\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "The following object is masked from ‘package:imputeTS’:\n", "\n", " na.remove\n", "\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "Loading required package: timeDate\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "Loading required package: timeSeries\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "Loading required package: fBasics\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "Loading required package: parallel\n" ] } ], "source": [ "# install.packages('robust')\n", "# install.packages('TSPred')\n", "\n", "library(TSPred)\n", "library(hexbin)\n", "library(robust)\n", "library(imputeTS)\n", "library(tseries)\n", "library(forecast)\n", "library(fUnitRoots)\n", "library(portes)\n", "library(nortest)\n", "library(tsoutliers)\n", "library(boot)" ] }, { "cell_type": "code", "execution_count": 3, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "# T <- read.xls(\"Dados das estacoes com 2018.xlsx\", sheet = 1, header = TRUE, stringsAsFactors = FALSE)\n", "T <- read.csv(\"Dados das estacoes com 2018.csv\")" ] }, { "cell_type": "code", "execution_count": 4, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\n", "\n", "\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\n", "
AnoMesAnoFT_INMETT_ICEAUR_INMETUR_ICEA
1961 1 15/01/196126.38193 NA 80.32258 NA
1961 2 15/02/196126.53500 NA 79.60714 NA
1961 3 15/03/196126.70000 NA 78.28226 NA
1961 4 15/04/196126.64600 NA 77.70000 NA
1961 5 15/05/196125.50000 NA 73.04839 NA
1961 6 15/06/196122.78000 NA 67.89167 NA
\n" ] }, "execution_count": 4, "metadata": { }, "output_type": "execute_result" } ], "source": [ "head(T)" ] }, { "cell_type": "code", "execution_count": 5, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "
\n", "\t
Ano
\n", "\t\t
'integer'
\n", "\t
Mes
\n", "\t\t
'integer'
\n", "\t
AnoF
\n", "\t\t
'factor'
\n", "\t
T_INMET
\n", "\t\t
'numeric'
\n", "\t
T_ICEA
\n", "\t\t
'numeric'
\n", "\t
UR_INMET
\n", "\t\t
'numeric'
\n", "\t
UR_ICEA
\n", "\t\t
'numeric'
\n", "
\n" ] }, "execution_count": 5, "metadata": { }, "output_type": "execute_result" } ], "source": [ "sapply(T, class) # Verificando os tipos" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "# Análise da Temperatura" ] }, { "cell_type": "code", "execution_count": 6, "metadata": { "collapsed": false, "scrolled": true }, "outputs": [ { "data": { "text/html": [ "png: 2" ] }, "execution_count": 6, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "9db8e7ad30ed62618b843dd9a3349e53ff92ac9f" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "options(repr.plot.width=16, repr.plot.height=8)\n", "\n", "T_INMET <- ts(T$T_INMET,start=c(1961,1),f=12)\n", "\n", "jpeg(\"g1.jpg\", width = 1000, height = 700)\n", "plot(T_INMET)\n", "dev.off()\n", "\n", "plot(T_INMET)" ] }, { "cell_type": "code", "execution_count": 7, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "png: 2" ] }, "execution_count": 7, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "29ed61ae141eec07386d3eab26c8ca1e9b108bd6" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "T_ICEA <- ts(T$T_ICEA,start=c(1961,1),f=12)\n", "\n", "jpeg(\"g2.jpg\", width = 1000, height = 700)\n", "plot(T_ICEA)\n", "dev.off()\n", "\n", "plot(T_ICEA)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Observe que as séries possuem falhas. Inicialmente as preencheremos com uma regressão robusta entre INMET e ICEA." ] }, { "cell_type": "code", "execution_count": 8, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "Call:\n", "lmRob(formula = T_INMET ~ T_ICEA, data = T)\n", "\n", "Residuals:\n", " Min 1Q Median 3Q Max \n", "-3.8137 -0.7333 -0.1029 0.6328 5.7032 \n", "\n", "Coefficients:\n", " Estimate Std. Error t value Pr(>|t|) \n", "(Intercept) 4.39915 0.77383 5.685 2.24e-08 ***\n", "T_ICEA 0.80588 0.02864 28.142 < 2e-16 ***\n", "---\n", "Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n", "\n", "Residual standard error: 0.9937 on 495 degrees of freedom\n", "Multiple R-Squared: 0.4702 \n", "\n", "Test for Bias:\n", " statistic p-value\n", "M-estimate 9.596 8.246e-03\n", "LS-estimate 45.666 1.213e-10\n", "194 observations deleted due to missingness " ] }, "execution_count": 8, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "png: 2" ] }, "execution_count": 8, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "6ec580990940fe3e20f15d2131e29523a353c03a" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "options(repr.plot.width=10, repr.plot.height=10)\n", "\n", "mod1 = lmRob(T_INMET ~ T_ICEA, data = T)\n", "summary(mod1)\n", "\n", "bin <- hexbin(T_ICEA, T_INMET, xbins=30)\n", "\n", "jpeg(\"g3.jpg\", width = 1000, height = 700)\n", "plot(bin)\n", "dev.off()\n", "\n", "plot(bin)" ] }, { "cell_type": "code", "execution_count": 9, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\n", "\n", "\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\n", "
AnoMesAnoFT_INMETT_ICEAUR_INMETUR_ICEA
3491990 1 15/01/1990NA 27.52 NA 75.2
3501990 2 15/02/1990NA 27.01 NA 78.4
3511990 3 15/03/1990NA 27.93 NA 74.8
3521990 4 15/04/1990NA 27.41 NA 74.0
3531990 5 15/05/1990NA 25.06 NA 75.6
3541990 6 15/06/1990NA 23.97 NA 77.7
\n" ] }, "execution_count": 9, "metadata": { }, "output_type": "execute_result" } ], "source": [ "aux <- T[is.na(T$T_INMET)&!is.na(T$T_ICEA),]\n", "head(aux)" ] }, { "cell_type": "code", "execution_count": 10, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "aux <- predict(mod1, newdata=T[is.na(T$T_INMET)&!is.na(T$T_ICEA),]) # Estima INMET com base no ICEA" ] }, { "cell_type": "code", "execution_count": 11, "metadata": { "collapsed": false, "scrolled": true }, "outputs": [ { "data": { "text/html": [ "
\n", "\t
349
\n", "\t\t
26.5770507637899
\n", "\t
350
\n", "\t\t
26.1660503977684
\n", "\t
351
\n", "\t\t
26.9074628227484
\n", "\t
352
\n", "\t\t
26.4884036260206
\n", "\t
353
\n", "\t\t
24.594578410039
\n", "\t
354
\n", "\t\t
23.7161658630518
\n", "
\n" ] }, "execution_count": 11, "metadata": { }, "output_type": "execute_result" } ], "source": [ "head(aux)" ] }, { "cell_type": "code", "execution_count": 12, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "# T[is.na(T$T_INMET)&!is.na(T$T_ICEA),] <- aux" ] }, { "cell_type": "code", "execution_count": 13, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "T[is.na(T$T_INMET)&!is.na(T$T_ICEA),'T_INMET'] <- aux" ] }, { "cell_type": "code", "execution_count": 14, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/svg+xml": "5655613562eb2714cb3322864133a4c77bcda5ea" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "T_INMET <- ts(T$T_INMET,start=c(1961,1),f=12)\n", "plot(T_INMET)" ] }, { "cell_type": "code", "execution_count": 15, "metadata": { "collapsed": false, "scrolled": true }, "outputs": [ { "data": { "text/html": [ "png: 2" ] }, "execution_count": 15, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "275148df6ed2a3ca9b5812b06cdc417044a9a2f3" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "# As falhas que sobraram foram preenchidas pelo estimador de Kalman.\n", "T_INMET <- na.kalman(T_INMET, model = \"auto.arima\")\n", "\n", "options(repr.plot.width=16, repr.plot.height=8)\n", "\n", "jpeg(\"g4.jpg\", width = 1000, height = 700)\n", "plot(T_INMET)\n", "dev.off()\n", "\n", "plot(T_INMET)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "# Análise da Umidade Relativa" ] }, { "cell_type": "code", "execution_count": 16, "metadata": { "collapsed": false, "scrolled": true }, "outputs": [ { "data": { "image/svg+xml": "406fa24476d3764084e635a54fcf57bdf4695836" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "UR_INMET <- ts(T$UR_INMET,start=c(1961,1),f=12)\n", "plot(UR_INMET)" ] }, { "cell_type": "code", "execution_count": 17, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/svg+xml": "792519dbfebd40f1db4cc03050f1406269bd0b6e" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "UR_ICEA <- ts(T$UR_ICEA,start=c(1961,1),f=12)\n", "plot(UR_ICEA)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Observe que as séries possuem falhas. Inicialmente as preencheremos com uma regressão entre INMET e ICEA." ] }, { "cell_type": "code", "execution_count": 18, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "Call:\n", "lmRob(formula = UR_INMET ~ UR_ICEA, data = T)\n", "\n", "Residuals:\n", " Min 1Q Median 3Q Max \n", "-14.8431 -2.7160 -0.1321 2.7989 16.2860 \n", "\n", "Coefficients:\n", " Estimate Std. Error t value Pr(>|t|) \n", "(Intercept) 14.02245 1.56939 8.935 <2e-16 ***\n", "UR_ICEA 0.86527 0.02302 37.593 <2e-16 ***\n", "---\n", "Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n", "\n", "Residual standard error: 4.127 on 432 degrees of freedom\n", "Multiple R-Squared: 0.5769 \n", "\n", "Test for Bias:\n", " statistic p-value\n", "M-estimate 7.301 2.598e-02\n", "LS-estimate 26.278 1.967e-06\n", "257 observations deleted due to missingness " ] }, "execution_count": 18, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "d90da83a05604ef9c2d9d46ee14bb46363cec0a6" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "mod2 = lmRob(UR_INMET ~ UR_ICEA, data = T)\n", "summary(mod2)\n", "\n", "bin <- hexbin(UR_INMET, UR_ICEA, xbins=30)\n", "plot(bin)" ] }, { "cell_type": "code", "execution_count": 19, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\n", "\n", "\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\n", "
AnoMesAnoFT_INMETT_ICEAUR_INMETUR_ICEA
2171979 1 15/01/197926.40 27.56 NA 78.8
2181979 2 15/02/197926.68 27.51 NA 77.8
2191979 3 15/03/197926.36 27.63 NA 78.8
2201979 4 15/04/197925.35 26.23 NA 77.7
2211979 5 15/05/197924.63 26.33 NA 70.8
2221979 6 15/06/197921.82 24.14 NA 65.4
\n" ] }, "execution_count": 19, "metadata": { }, "output_type": "execute_result" } ], "source": [ "aux <- T[is.na(T$UR_INMET)&!is.na(T$UR_ICEA),]\n", "head(aux)" ] }, { "cell_type": "code", "execution_count": 20, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "aux <- predict(mod2, newdata=T[is.na(T$UR_INMET)&!is.na(T$UR_ICEA),]) # Estima INMET com base no ICEA" ] }, { "cell_type": "code", "execution_count": 21, "metadata": { "collapsed": false, "scrolled": true }, "outputs": [ { "data": { "text/html": [ "
\n", "\t
217
\n", "\t\t
82.2057788962603
\n", "\t
218
\n", "\t\t
81.3405082703494
\n", "\t
219
\n", "\t\t
82.2057788962603
\n", "\t
220
\n", "\t\t
81.2539812077583
\n", "\t
221
\n", "\t\t
75.283613888973
\n", "\t
222
\n", "\t\t
70.611152509054
\n", "
\n" ] }, "execution_count": 21, "metadata": { }, "output_type": "execute_result" } ], "source": [ "head(aux)" ] }, { "cell_type": "code", "execution_count": 22, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "T[is.na(T$UR_INMET)&!is.na(T$UR_ICEA),'UR_INMET'] <- aux" ] }, { "cell_type": "code", "execution_count": 23, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/svg+xml": "99e4f4ad0caf3007e755ca85d7f198fdce57d0fd" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "UR_INMET <- ts(T$UR_INMET,start=c(1961,1),f=12)\n", "plot(UR_INMET)" ] }, { "cell_type": "code", "execution_count": 24, "metadata": { "collapsed": false, "scrolled": true }, "outputs": [ { "data": { "image/svg+xml": "c6297d923c42a6b8896b68f0b41448c6cb73c0ac" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "# As falhas que sobraram foram preenchidas pelo estimador de Kalman.\n", "UR_INMET <- na.kalman(UR_INMET, model = \"auto.arima\")\n", "plot(UR_INMET)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "# Cálculo do índice térmico" ] }, { "cell_type": "code", "execution_count": 25, "metadata": { "collapsed": false }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Installing package into ‘/home/user/R/x86_64-pc-linux-gnu-library/3.4’\n", "(as ‘lib’ is unspecified)\n" ] } ], "source": [ "install.packages('weathermetrics')" ] }, { "cell_type": "code", "execution_count": 26, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "library(weathermetrics)" ] }, { "cell_type": "code", "execution_count": 27, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "png: 2" ] }, "execution_count": 27, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "1ece29e7551b40fc03f0274aed88b5dd8845b6a3" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "it <- heat.index(t=T_INMET, rh=UR_INMET, temperature.metric=\"celsius\", output.metric=\"celsius\", round=2)\n", "it <- ts(it,start=c(1961,1),f=12)\n", "\n", "jpeg(\"g9.jpg\", width = 1000, height = 700)\n", "plot(it)\n", "dev.off()\n", "\n", "plot(it)" ] }, { "cell_type": "code", "execution_count": 28, "metadata": { "collapsed": false, "scrolled": true }, "outputs": [ { "data": { "text/html": [ "png: 2" ] }, "execution_count": 28, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "36f13c7e89fe6d8af6705b85609005d7d474dac5" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "# Decomposição da série temporal (aditiva, possivelmente baseada na transformáda rápida de Fourier, e apenas para fins de análise preliminar)\n", "\n", "jpeg(\"g10.jpg\", width = 1000, height = 700)\n", "plot(decompose(it))\n", "dev.off()\n", "\n", "plot(decompose(it))" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "# ARIMA" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Etapas:\n", "a) Verificar se existe a necessidade de uma transformação na série original, com objetivo de estabilizar a variância;\n", "b) Tornar a série estacionária por meio de diferenças, de modo que o processo dZt seja reduzido a um ARMA(p,q)\n", "c) Identificar o processo ARMA(p,q) resultante.\n", "d) Verificação da estacionariedade e da invertibilidade.\n", "\n", "\n", "FAC : correlação simples entre Zt e Zt – k em função da defasagem k.\n", "FACP: correlação entre Zt e Zt – k em função da defasagem k, filtrado o efeito de todas as outras defasagens sobre Zt e Zt – k.\n", "\n", "FACP -> AR\n", "\n", "D -> I\n", "\n", "FAC -> MA\n", "\n", "1- Número de AR (auto-regressivo) termos (p): termos AR são apenas defasagens da variável dependente. Por exemplo, se o símbolo p representa 5, os preditores de x (t) irá ser X (t-1) … .x (T-5).\n", "2- Número de MA (média móvel) termos (q): termos MA estão defasados erros de previsão na equação de predição. Por exemplo, se q é 5, os preditores para x (t) será E (t-1) … .e (t-5) onde e (i) é a diferença entre a média móvel ao valor imediato e real.\n", "3- Número de Diferenças (d): Estes são o número de diferenças não sazonal, ou seja, neste caso, tomamos a primeira diferença de ordem. Assim, ou nós podemos passar essa variável e colocar d = 0, ou passar a variável original e coloca -d = 1. Ambos irão gerar mesmos resultados." ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Notação: arima(p, d, q), sendo p relacionado a autocorrelação parcial, d a diferença entre os valores, q associado a autocorrelação)." ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "## Testes da estacionariedade da série" ] }, { "cell_type": "code", "execution_count": 29, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "Title:\n", " KPSS Unit Root Test\n", "\n", "Test Results:\n", " NA\n", "\n", "Description:\n", " Thu Oct 18 18:23:32 2018 by user: \n" ] }, "execution_count": 29, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "850be4912129427e0110cd3133e956f7d15dfa4f" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "urkpssTest(it, type = c(\"tau\"), lags = c(\"short\"),use.lag = NULL, doplot = TRUE)" ] }, { "cell_type": "code", "execution_count": 30, "metadata": { "collapsed": false }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Warning message in adf.test(it, alternative = \"stationary\"):\n", "“p-value smaller than printed p-value”" ] }, { "data": { "text/plain": [ "\n", "\tAugmented Dickey-Fuller Test\n", "\n", "data: it\n", "Dickey-Fuller = -14.159, Lag order = 8, p-value = 0.01\n", "alternative hypothesis: stationary\n" ] }, "execution_count": 30, "metadata": { }, "output_type": "execute_result" } ], "source": [ "adf.test(it, alternative=\"stationary\")" ] }, { "cell_type": "code", "execution_count": 31, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "png: 2" ] }, "execution_count": 31, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "png: 2" ] }, "execution_count": 31, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "fcac4d8f4a9a9e4ed1709e1a3d1a527fa6da7c90" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "4483e1210f16885737d224b276d9135c59921b86" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "jpeg(\"g11.jpg\", width = 1000, height = 700)\n", "acf(it)\n", "dev.off()\n", "acf(it)\n", "\n", "jpeg(\"g12.jpg\", width = 1000, height = 700)\n", "pacf(it)\n", "dev.off()\n", "pacf(it)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Se rejeitarmos a hipótese da estacionariedade, desta forma, partimos de d=1." ] }, { "cell_type": "code", "execution_count": 32, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "1" ] }, "execution_count": 32, "metadata": { }, "output_type": "execute_result" } ], "source": [ "ndiffs(it) # Valor de d para fazer a série estacionária." ] }, { "cell_type": "code", "execution_count": 33, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "png: 2" ] }, "execution_count": 33, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "png: 2" ] }, "execution_count": 33, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "94955b7ba0c850807b1c1b2e92e78c2de0e97d21" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "5447248712808c1e30f76996a6e2763a3bf082a8" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "# Calculando as correlações para uma diferença de ordem 1\n", "jpeg(\"g12.jpg\", width = 1000, height = 700)\n", "acf(diff(it, differences=1))\n", "dev.off()\n", "acf(diff(it, differences=1))\n", "\n", "jpeg(\"g13.jpg\", width = 1000, height = 700)\n", "pacf(diff(it, differences=1))\n", "dev.off()\n", "pacf(diff(it, differences=1))" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "Como o ACF cai depois do primeiro lag, podemos partir de p = 1-3. Para o PACF o q poderia ser igual a 2, temos, então, um ARIMA(p, d, q) = ARIMA (1, 1, 2) com os mesmos parâmetros para sazonalidade." ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "## Ajustando os modelos" ] }, { "cell_type": "code", "execution_count": 34, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "Series: it \n", "ARIMA(1,1,2) \n", "\n", "Coefficients:\n", " ar1 ma1 ma2\n", " -0.2744 0.3753 0.0765\n", "s.e. 0.3108 0.3057 0.0477\n", "\n", "sigma^2 estimated as 4.389: log likelihood=-1487.86\n", "AIC=2983.72 AICc=2983.78 BIC=3001.87" ] }, "execution_count": 34, "metadata": { }, "output_type": "execute_result" } ], "source": [ "ajuste_0 <- Arima(it, c(1, 1, 2)) # Sem sazonalidade\n", "ajuste_0" ] }, { "cell_type": "code", "execution_count": 35, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "Series: it \n", "ARIMA(1,1,2)(1,1,2)[12] \n", "\n", "Coefficients:\n", " ar1 ma1 ma2 sar1 sma1 sma2\n", " 0.5930 -1.3374 0.3538 -0.9781 -0.0084 -0.9405\n", "s.e. 0.1312 0.1503 0.1430 0.0489 0.0696 0.0716\n", "\n", "sigma^2 estimated as 1.513: log likelihood=-1119.62\n", "AIC=2253.23 AICc=2253.4 BIC=2284.87" ] }, "execution_count": 35, "metadata": { }, "output_type": "execute_result" } ], "source": [ "ajuste_1 <- Arima(it, c(1, 1, 2), seasonal=list(order = c(1, 1, 2), period = 12)) # Com sazonalidade\n", "ajuste_1" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Utilizando via auto.arima" ] }, { "cell_type": "code", "execution_count": 36, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", " ARIMA(0,0,0)(0,1,0)[12] : 2686.114\n", " ARIMA(0,0,0)(0,1,0)[12] with drift : 2688.049\n", " ARIMA(0,0,0)(0,1,1)[12] : 2324.46\n", " ARIMA(0,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,0)(0,1,2)[12] : 2326.375\n", " ARIMA(0,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,0)(1,1,0)[12] : 2477.808\n", " ARIMA(0,0,0)(1,1,0)[12] with drift : 2479.523\n", " ARIMA(0,0,0)(1,1,1)[12] : 2326.363\n", " ARIMA(0,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,0)(1,1,2)[12] : Inf\n", " ARIMA(0,0,0)(1,1,2)[12] with drift : Inf\n", " ARIMA(0,0,0)(2,1,0)[12] : 2394.978\n", " ARIMA(0,0,0)(2,1,0)[12] with drift : 2396.522\n", " ARIMA(0,0,0)(2,1,1)[12] : 2326.994\n", " ARIMA(0,0,0)(2,1,1)[12] with drift : Inf\n", " ARIMA(0,0,0)(2,1,2)[12] : 2330.21\n", " ARIMA(0,0,0)(2,1,2)[12] with drift : Inf\n", " ARIMA(0,0,1)(0,1,0)[12] : 2640.773\n", " ARIMA(0,0,1)(0,1,0)[12] with drift : 2642.737\n", " ARIMA(0,0,1)(0,1,1)[12] : 2276.591\n", " ARIMA(0,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,1)(0,1,2)[12] : 2278.55\n", " ARIMA(0,0,1)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,1)(1,1,0)[12] : 2439.248\n", " ARIMA(0,0,1)(1,1,0)[12] with drift : 2441.059\n", " ARIMA(0,0,1)(1,1,1)[12] : 2278.543\n", " ARIMA(0,0,1)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,1)(1,1,2)[12] : Inf\n", " ARIMA(0,0,1)(1,1,2)[12] with drift : Inf\n", " ARIMA(0,0,1)(2,1,0)[12] : 2366.268\n", " ARIMA(0,0,1)(2,1,0)[12] with drift : 2367.953\n", " ARIMA(0,0,1)(2,1,1)[12] : Inf\n", " ARIMA(0,0,1)(2,1,1)[12] with drift : Inf\n", " ARIMA(0,0,1)(2,1,2)[12] : 2281.766\n", " ARIMA(0,0,1)(2,1,2)[12] with drift : Inf\n", " ARIMA(0,0,2)(0,1,0)[12] : 2640.657\n", " ARIMA(0,0,2)(0,1,0)[12] with drift : 2642.63\n", " ARIMA(0,0,2)(0,1,1)[12] : 2270.698\n", " ARIMA(0,0,2)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,2)(0,1,2)[12] : 2272.719\n", " ARIMA(0,0,2)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,2)(1,1,0)[12] : 2435.692\n", " ARIMA(0,0,2)(1,1,0)[12] with drift : 2437.532\n", " ARIMA(0,0,2)(1,1,1)[12] : 2272.718\n", " ARIMA(0,0,2)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,2)(1,1,2)[12] : Inf\n", " ARIMA(0,0,2)(1,1,2)[12] with drift : Inf\n", " ARIMA(0,0,2)(2,1,0)[12] : 2360.922\n", " ARIMA(0,0,2)(2,1,0)[12] with drift : 2362.653\n", " ARIMA(0,0,2)(2,1,1)[12] : Inf\n", " ARIMA(0,0,2)(2,1,1)[12] with drift : Inf\n", " ARIMA(0,0,3)(0,1,0)[12] : 2626.19\n", " ARIMA(0,0,3)(0,1,0)[12] with drift : 2628.18\n", " ARIMA(0,0,3)(0,1,1)[12] : Inf\n", " ARIMA(0,0,3)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,3)(0,1,2)[12] : Inf\n", " ARIMA(0,0,3)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,3)(1,1,0)[12] : 2422.329\n", " ARIMA(0,0,3)(1,1,0)[12] with drift : 2424.215\n", " ARIMA(0,0,3)(1,1,1)[12] : Inf\n", " ARIMA(0,0,3)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,3)(2,1,0)[12] : 2351.283\n", " ARIMA(0,0,3)(2,1,0)[12] with drift : 2353.079\n", " ARIMA(0,0,4)(0,1,0)[12] : 2627.224\n", " ARIMA(0,0,4)(0,1,0)[12] with drift : 2629.224\n", " ARIMA(0,0,4)(0,1,1)[12] : Inf\n", " ARIMA(0,0,4)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,4)(1,1,0)[12] : 2423.745\n", " ARIMA(0,0,4)(1,1,0)[12] with drift : 2425.646\n", " ARIMA(0,0,5)(0,1,0)[12] : 2627.556\n", " ARIMA(0,0,5)(0,1,0)[12] with drift : 2629.565\n", " ARIMA(1,0,0)(0,1,0)[12] : 2634.957\n", " ARIMA(1,0,0)(0,1,0)[12] with drift : 2636.929\n", " ARIMA(1,0,0)(0,1,1)[12] : 2265.29\n", " ARIMA(1,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,0)(0,1,2)[12] : Inf\n", " ARIMA(1,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(1,0,0)(1,1,0)[12] : 2431.531\n", " ARIMA(1,0,0)(1,1,0)[12] with drift : 2433.377\n", " ARIMA(1,0,0)(1,1,1)[12] : Inf\n", " ARIMA(1,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(1,0,0)(1,1,2)[12] : Inf\n", " ARIMA(1,0,0)(1,1,2)[12] with drift : Inf\n", " ARIMA(1,0,0)(2,1,0)[12] : 2359.552\n", " ARIMA(1,0,0)(2,1,0)[12] with drift : 2361.286\n", " ARIMA(1,0,0)(2,1,1)[12] : Inf\n", " ARIMA(1,0,0)(2,1,1)[12] with drift : Inf\n", " ARIMA(1,0,0)(2,1,2)[12] : Inf\n", " ARIMA(1,0,0)(2,1,2)[12] with drift : Inf\n", " ARIMA(1,0,1)(0,1,0)[12] : 2629.704\n", " ARIMA(1,0,1)(0,1,0)[12] with drift : 2631.698\n", " ARIMA(1,0,1)(0,1,1)[12] : Inf\n", " ARIMA(1,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,1)(0,1,2)[12] : Inf\n", " ARIMA(1,0,1)(0,1,2)[12] with drift : Inf\n", " ARIMA(1,0,1)(1,1,0)[12] : 2421.485\n", " ARIMA(1,0,1)(1,1,0)[12] with drift : 2423.41\n", " ARIMA(1,0,1)(1,1,1)[12] : Inf\n", " ARIMA(1,0,1)(1,1,1)[12] with drift : Inf\n", " ARIMA(1,0,1)(1,1,2)[12] : Inf\n", " ARIMA(1,0,1)(1,1,2)[12] with drift : Inf\n", " ARIMA(1,0,1)(2,1,0)[12] : 2347.271\n", " ARIMA(1,0,1)(2,1,0)[12] with drift : 2349.139\n", " ARIMA(1,0,1)(2,1,1)[12] : Inf\n", " ARIMA(1,0,1)(2,1,1)[12] with drift : Inf\n", " ARIMA(1,0,2)(0,1,0)[12] : 2630.53\n", " ARIMA(1,0,2)(0,1,0)[12] with drift : 2632.532\n", " ARIMA(1,0,2)(0,1,1)[12] : Inf\n", " ARIMA(1,0,2)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,2)(0,1,2)[12] : Inf\n", " ARIMA(1,0,2)(0,1,2)[12] with drift : Inf\n", " ARIMA(1,0,2)(1,1,0)[12] : 2423.202\n", " ARIMA(1,0,2)(1,1,0)[12] with drift : 2425.136\n", " ARIMA(1,0,2)(1,1,1)[12] : Inf\n", " ARIMA(1,0,2)(1,1,1)[12] with drift : Inf\n", " ARIMA(1,0,2)(2,1,0)[12] : 2349.278\n", " ARIMA(1,0,2)(2,1,0)[12] with drift : 2351.154\n", " ARIMA(1,0,3)(0,1,0)[12] : 2626.242\n", " ARIMA(1,0,3)(0,1,0)[12] with drift : 2628.245\n", " ARIMA(1,0,3)(0,1,1)[12] : Inf\n", " ARIMA(1,0,3)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,3)(1,1,0)[12] : 2422.703\n", " ARIMA(1,0,3)(1,1,0)[12] with drift : 2424.623\n", " ARIMA(1,0,4)(0,1,0)[12] : 2627.554\n", " ARIMA(1,0,4)(0,1,0)[12] with drift : 2629.567\n", " ARIMA(2,0,0)(0,1,0)[12] : 2634.301\n", " ARIMA(2,0,0)(0,1,0)[12] with drift : 2636.284\n", " ARIMA(2,0,0)(0,1,1)[12] : Inf\n", " ARIMA(2,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(2,0,0)(0,1,2)[12] : Inf\n", " ARIMA(2,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(2,0,0)(1,1,0)[12] : 2427.388\n", " ARIMA(2,0,0)(1,1,0)[12] with drift : 2429.271\n", " ARIMA(2,0,0)(1,1,1)[12] : Inf\n", " ARIMA(2,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(2,0,0)(1,1,2)[12] : Inf\n", " ARIMA(2,0,0)(1,1,2)[12] with drift : Inf\n", " ARIMA(2,0,0)(2,1,0)[12] : 2353.586\n", " ARIMA(2,0,0)(2,1,0)[12] with drift : 2355.382\n", " ARIMA(2,0,0)(2,1,1)[12] : Inf\n", " ARIMA(2,0,0)(2,1,1)[12] with drift : Inf\n", " ARIMA(2,0,1)(0,1,0)[12] : 2630.941\n", " ARIMA(2,0,1)(0,1,0)[12] with drift : 2632.943\n", " ARIMA(2,0,1)(0,1,1)[12] : Inf\n", " ARIMA(2,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(2,0,1)(0,1,2)[12] : Inf\n", " ARIMA(2,0,1)(0,1,2)[12] with drift : Inf\n", " ARIMA(2,0,1)(1,1,0)[12] : 2423.271\n", " ARIMA(2,0,1)(1,1,0)[12] with drift : 2425.205\n", " ARIMA(2,0,1)(1,1,1)[12] : Inf\n", " ARIMA(2,0,1)(1,1,1)[12] with drift : Inf\n", " ARIMA(2,0,1)(2,1,0)[12] : Inf\n", " ARIMA(2,0,1)(2,1,0)[12] with drift : Inf\n", " ARIMA(2,0,2)(0,1,0)[12] : Inf\n", " ARIMA(2,0,2)(0,1,0)[12] with drift : Inf\n", " ARIMA(2,0,2)(0,1,1)[12] : Inf\n", " ARIMA(2,0,2)(0,1,1)[12] with drift : Inf\n", " ARIMA(2,0,2)(1,1,0)[12] : 2422.676\n", " ARIMA(2,0,2)(1,1,0)[12] with drift : 2424.618\n", " ARIMA(2,0,3)(0,1,0)[12] : 2627.22\n", " ARIMA(2,0,3)(0,1,0)[12] with drift : 2629.234\n", " ARIMA(3,0,0)(0,1,0)[12] : 2624.389\n", " ARIMA(3,0,0)(0,1,0)[12] with drift : 2626.391\n", " ARIMA(3,0,0)(0,1,1)[12] : Inf\n", " ARIMA(3,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(3,0,0)(0,1,2)[12] : Inf\n", " ARIMA(3,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(3,0,0)(1,1,0)[12] : 2419.92\n", " ARIMA(3,0,0)(1,1,0)[12] with drift : 2421.846\n", " ARIMA(3,0,0)(1,1,1)[12] : Inf\n", " ARIMA(3,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(3,0,0)(2,1,0)[12] : 2347.963\n", " ARIMA(3,0,0)(2,1,0)[12] with drift : 2349.818\n", " ARIMA(3,0,1)(0,1,0)[12] : 2625.213\n", " ARIMA(3,0,1)(0,1,0)[12] with drift : 2627.217\n", " ARIMA(3,0,1)(0,1,1)[12] : Inf\n", " ARIMA(3,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(3,0,1)(1,1,0)[12] : 2420.888\n", " ARIMA(3,0,1)(1,1,0)[12] with drift : 2422.81\n", " ARIMA(3,0,2)(0,1,0)[12] : 2626.973\n", " ARIMA(3,0,2)(0,1,0)[12] with drift : 2628.986\n", " ARIMA(4,0,0)(0,1,0)[12] : 2625.894\n", " ARIMA(4,0,0)(0,1,0)[12] with drift : 2627.9\n", " ARIMA(4,0,0)(0,1,1)[12] : Inf\n", " ARIMA(4,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(4,0,0)(1,1,0)[12] : 2421.564\n", " ARIMA(4,0,0)(1,1,0)[12] with drift : 2423.489\n", " ARIMA(4,0,1)(0,1,0)[12] : 2627.148\n", " ARIMA(4,0,1)(0,1,0)[12] with drift : 2629.159\n", " ARIMA(5,0,0)(0,1,0)[12] : 2626.384\n", " ARIMA(5,0,0)(0,1,0)[12] with drift : 2628.399\n", "\n", "\n", "\n", " Best model: ARIMA(1,0,0)(0,1,1)[12] \n", "\n" ] }, { "data": { "text/plain": [ "Series: it \n", "ARIMA(1,0,0)(0,1,1)[12] \n", "\n", "Coefficients:\n", " ar1 sma1\n", " 0.2988 -0.8859\n", "s.e. 0.0375 0.0245\n", "\n", "sigma^2 estimated as 1.592: log likelihood=-1129.63\n", "AIC=2265.25 AICc=2265.29 BIC=2278.82" ] }, "execution_count": 36, "metadata": { }, "output_type": "execute_result" } ], "source": [ "# Automatizando\n", "ajuste_auto <- auto.arima(it, stepwise=FALSE, approximation=FALSE, trace=TRUE)\n", "ajuste_auto" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### No caso, obteve um AIC maior, mas com SE menores. Pelo princípio da parcimônia, e sabendo que existe sazonalidade, utilizaremos, de início, o modelo (0,0,2)(0,1,1)12." ] }, { "cell_type": "code", "execution_count": 37, "metadata": { "collapsed": false }, "outputs": [ { "data": { "image/svg+xml": "988d07f26294bf6d5bd42a23c451f6a8b3cbc201" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "d7872c2ab67fddbdeb4b26268484c56453b69c5f" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "acf(ajuste_auto$residuals)\n", "pacf(ajuste_auto$residuals)" ] }, { "cell_type": "code", "execution_count": 38, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "png: 2" ] }, "execution_count": 38, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/plain": [ "\n", "\tLilliefors (Kolmogorov-Smirnov) normality test\n", "\n", "data: ajuste_auto$residuals\n", "D = 0.04776, p-value = 0.0007485\n" ] }, "execution_count": 38, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "2acc4ed404ed68d079ef1f4281255fe03ae7adfd" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "jpeg(\"g14.jpg\", width = 1000, height = 700)\n", "qqnorm(ajuste_auto$residuals)\n", "qqline(ajuste_auto$residuals)\n", "dev.off()\n", "\n", "qqnorm(ajuste_auto$residuals)\n", "qqline(ajuste_auto$residuals)\n", "\n", "lillie.test(ajuste_auto$residuals)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Não adere a Normal, mas também não apresenta comportamento anômalo. É o caso de tentar uma transformação de Box-Cox." ] }, { "cell_type": "code", "execution_count": 39, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "-0.75570192309047" ] }, "execution_count": 39, "metadata": { }, "output_type": "execute_result" } ], "source": [ "l <- BoxCox.lambda(it)\n", "l" ] }, { "cell_type": "code", "execution_count": 40, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", " ARIMA(0,0,0)(0,1,0)[12] : -5184.219\n", " ARIMA(0,0,0)(0,1,0)[12] with drift : -5182.264\n", " ARIMA(0,0,0)(0,1,1)[12] : -5533.35\n", " ARIMA(0,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,0)(0,1,2)[12] : -5531.396\n", " ARIMA(0,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,0)(1,1,0)[12] : -5382.485\n", " ARIMA(0,0,0)(1,1,0)[12] with drift : -5380.732\n", " ARIMA(0,0,0)(1,1,1)[12] : -5531.404\n", " ARIMA(0,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,0)(1,1,2)[12] : Inf\n", " ARIMA(0,0,0)(1,1,2)[12] with drift : Inf\n", " ARIMA(0,0,0)(2,1,0)[12] : -5456.414\n", " ARIMA(0,0,0)(2,1,0)[12] with drift : -5454.814\n", " ARIMA(0,0,0)(2,1,1)[12] : -5531.065\n", " ARIMA(0,0,0)(2,1,1)[12] with drift : Inf\n", " ARIMA(0,0,0)(2,1,2)[12] : Inf\n", " ARIMA(0,0,0)(2,1,2)[12] with drift : Inf\n", " ARIMA(0,0,1)(0,1,0)[12] : -5218.002\n", " ARIMA(0,0,1)(0,1,0)[12] with drift : -5216.026\n", " ARIMA(0,0,1)(0,1,1)[12] : -5567.573\n", " ARIMA(0,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,1)(0,1,2)[12] : -5566.239\n", " ARIMA(0,0,1)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,1)(1,1,0)[12] : -5408.399\n", " ARIMA(0,0,1)(1,1,0)[12] with drift : -5406.568\n", " ARIMA(0,0,1)(1,1,1)[12] : -5566.321\n", " ARIMA(0,0,1)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,1)(1,1,2)[12] : -5563.524\n", " ARIMA(0,0,1)(1,1,2)[12] with drift : Inf\n", " ARIMA(0,0,1)(2,1,0)[12] : -5474.445\n", " ARIMA(0,0,1)(2,1,0)[12] with drift : -5472.733\n", " ARIMA(0,0,1)(2,1,1)[12] : Inf\n", " ARIMA(0,0,1)(2,1,1)[12] with drift : Inf\n", " ARIMA(0,0,1)(2,1,2)[12] : Inf\n", " ARIMA(0,0,1)(2,1,2)[12] with drift : Inf\n", " ARIMA(0,0,2)(0,1,0)[12] : -5216.614\n", " ARIMA(0,0,2)(0,1,0)[12] with drift : -5214.63\n", " ARIMA(0,0,2)(0,1,1)[12] : -5571.615\n", " ARIMA(0,0,2)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,2)(0,1,2)[12] : -5570\n", " ARIMA(0,0,2)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,2)(1,1,0)[12] : -5410.267\n", " ARIMA(0,0,2)(1,1,0)[12] with drift : -5408.412\n", " ARIMA(0,0,2)(1,1,1)[12] : Inf\n", " ARIMA(0,0,2)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,2)(1,1,2)[12] : Inf\n", " ARIMA(0,0,2)(1,1,2)[12] with drift : Inf\n", " ARIMA(0,0,2)(2,1,0)[12] : -5478.05\n", " ARIMA(0,0,2)(2,1,0)[12] with drift : -5476.297\n", " ARIMA(0,0,2)(2,1,1)[12] : Inf\n", " ARIMA(0,0,2)(2,1,1)[12] with drift : Inf\n", " ARIMA(0,0,3)(0,1,0)[12] : -5226.15\n", " ARIMA(0,0,3)(0,1,0)[12] with drift : -5224.152\n", " ARIMA(0,0,3)(0,1,1)[12] : -5579.704\n", " ARIMA(0,0,3)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,3)(0,1,2)[12] : Inf\n", " ARIMA(0,0,3)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,3)(1,1,0)[12] : -5418.451\n", " ARIMA(0,0,3)(1,1,0)[12] with drift : -5416.558\n", " ARIMA(0,0,3)(1,1,1)[12] : Inf\n", " ARIMA(0,0,3)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,3)(2,1,0)[12] : -5483.142\n", " ARIMA(0,0,3)(2,1,0)[12] with drift : -5481.336\n", " ARIMA(0,0,4)(0,1,0)[12] : -5224.486\n", " ARIMA(0,0,4)(0,1,0)[12] with drift : -5222.48\n", " ARIMA(0,0,4)(0,1,1)[12] : Inf\n", " ARIMA(0,0,4)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,4)(1,1,0)[12] : -5416.737\n", " ARIMA(0,0,4)(1,1,0)[12] with drift : -5414.831\n", " ARIMA(0,0,5)(0,1,0)[12] : -5223.889\n", " ARIMA(0,0,5)(0,1,0)[12] with drift : -5221.875\n", " ARIMA(1,0,0)(0,1,0)[12] : -5220.45\n", " ARIMA(1,0,0)(0,1,0)[12] with drift : -5218.47\n", " ARIMA(1,0,0)(0,1,1)[12] : -5574.769\n", " ARIMA(1,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,0)(0,1,2)[12] : Inf\n", " ARIMA(1,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(1,0,0)(1,1,0)[12] : -5412.812\n", " ARIMA(1,0,0)(1,1,0)[12] with drift : -5410.959\n", " ARIMA(1,0,0)(1,1,1)[12] : Inf\n", " ARIMA(1,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(1,0,0)(1,1,2)[12] : -5570.727\n", " ARIMA(1,0,0)(1,1,2)[12] with drift : Inf\n", " ARIMA(1,0,0)(2,1,0)[12] : -5478.239\n", " ARIMA(1,0,0)(2,1,0)[12] with drift : -5476.494\n", " ARIMA(1,0,0)(2,1,1)[12] : Inf\n", " ARIMA(1,0,0)(2,1,1)[12] with drift : Inf\n", " ARIMA(1,0,0)(2,1,2)[12] : Inf\n", " ARIMA(1,0,0)(2,1,2)[12] with drift : Inf\n", " ARIMA(1,0,1)(0,1,0)[12] : -5222.121\n", " ARIMA(1,0,1)(0,1,0)[12] with drift : -5220.123\n", " ARIMA(1,0,1)(0,1,1)[12] : Inf\n", " ARIMA(1,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,1)(0,1,2)[12] : Inf\n", " ARIMA(1,0,1)(0,1,2)[12] with drift : Inf\n", " ARIMA(1,0,1)(1,1,0)[12] : -5420.563\n", " ARIMA(1,0,1)(1,1,0)[12] with drift : -5418.635\n", " ARIMA(1,0,1)(1,1,1)[12] : Inf\n", " ARIMA(1,0,1)(1,1,1)[12] with drift : Inf\n", " ARIMA(1,0,1)(1,1,2)[12] : Inf\n", " ARIMA(1,0,1)(1,1,2)[12] with drift : Inf\n", " ARIMA(1,0,1)(2,1,0)[12] : -5487.814\n", " ARIMA(1,0,1)(2,1,0)[12] with drift : -5485.94\n", " ARIMA(1,0,1)(2,1,1)[12] : Inf\n", " ARIMA(1,0,1)(2,1,1)[12] with drift : Inf\n", " ARIMA(1,0,2)(0,1,0)[12] : -5222.038\n", " ARIMA(1,0,2)(0,1,0)[12] with drift : -5220.031\n", " ARIMA(1,0,2)(0,1,1)[12] : Inf\n", " ARIMA(1,0,2)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,2)(0,1,2)[12] : Inf\n", " ARIMA(1,0,2)(0,1,2)[12] with drift : Inf\n", " ARIMA(1,0,2)(1,1,0)[12] : -5418.903\n", " ARIMA(1,0,2)(1,1,0)[12] with drift : -5416.966\n", " ARIMA(1,0,2)(1,1,1)[12] : Inf\n", " ARIMA(1,0,2)(1,1,1)[12] with drift : Inf\n", " ARIMA(1,0,2)(2,1,0)[12] : -5485.803\n", " ARIMA(1,0,2)(2,1,0)[12] with drift : -5483.921\n", " ARIMA(1,0,3)(0,1,0)[12] : -5225.087\n", " ARIMA(1,0,3)(0,1,0)[12] with drift : -5223.078\n", " ARIMA(1,0,3)(0,1,1)[12] : Inf\n", " ARIMA(1,0,3)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,3)(1,1,0)[12] : -5418.182\n", " ARIMA(1,0,3)(1,1,0)[12] with drift : -5416.243\n", " ARIMA(1,0,4)(0,1,0)[12] : -5223.891\n", " ARIMA(1,0,4)(0,1,0)[12] with drift : -5221.873\n", " ARIMA(2,0,0)(0,1,0)[12] : -5219.32\n", " ARIMA(2,0,0)(0,1,0)[12] with drift : -5217.331\n", " ARIMA(2,0,0)(0,1,1)[12] : Inf\n", " ARIMA(2,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(2,0,0)(0,1,2)[12] : Inf\n", " ARIMA(2,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(2,0,0)(1,1,0)[12] : -5415.145\n", " ARIMA(2,0,0)(1,1,0)[12] with drift : -5413.26\n", " ARIMA(2,0,0)(1,1,1)[12] : Inf\n", " ARIMA(2,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(2,0,0)(1,1,2)[12] : Inf\n", " ARIMA(2,0,0)(1,1,2)[12] with drift : Inf\n", " ARIMA(2,0,0)(2,1,0)[12] : -5482.265\n", " ARIMA(2,0,0)(2,1,0)[12] with drift : -5480.467\n", " ARIMA(2,0,0)(2,1,1)[12] : Inf\n", " ARIMA(2,0,0)(2,1,1)[12] with drift : Inf\n", " ARIMA(2,0,1)(0,1,0)[12] : -5221.426\n", " ARIMA(2,0,1)(0,1,0)[12] with drift : -5219.42\n", " ARIMA(2,0,1)(0,1,1)[12] : Inf\n", " ARIMA(2,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(2,0,1)(0,1,2)[12] : Inf\n", " ARIMA(2,0,1)(0,1,2)[12] with drift : Inf\n", " ARIMA(2,0,1)(1,1,0)[12] : -5418.844\n", " ARIMA(2,0,1)(1,1,0)[12] with drift : -5416.907\n", " ARIMA(2,0,1)(1,1,1)[12] : Inf\n", " ARIMA(2,0,1)(1,1,1)[12] with drift : Inf\n", " ARIMA(2,0,1)(2,1,0)[12] : Inf\n", " ARIMA(2,0,1)(2,1,0)[12] with drift : Inf\n", " ARIMA(2,0,2)(0,1,0)[12] : Inf\n", " ARIMA(2,0,2)(0,1,0)[12] with drift : Inf\n", " ARIMA(2,0,2)(0,1,1)[12] : Inf\n", " ARIMA(2,0,2)(0,1,1)[12] with drift : Inf\n", " ARIMA(2,0,2)(1,1,0)[12] : -5418.63\n", " ARIMA(2,0,2)(1,1,0)[12] with drift : -5416.687\n", " ARIMA(2,0,3)(0,1,0)[12] : -5224.318\n", " ARIMA(2,0,3)(0,1,0)[12] with drift : -5222.298\n", " ARIMA(3,0,0)(0,1,0)[12] : -5226.97\n", " ARIMA(3,0,0)(0,1,0)[12] with drift : -5224.964\n", " ARIMA(3,0,0)(0,1,1)[12] : Inf\n", " ARIMA(3,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(3,0,0)(0,1,2)[12] : Inf\n", " ARIMA(3,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(3,0,0)(1,1,0)[12] : -5420.527\n", " ARIMA(3,0,0)(1,1,0)[12] with drift : -5418.604\n", " ARIMA(3,0,0)(1,1,1)[12] : Inf\n", " ARIMA(3,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(3,0,0)(2,1,0)[12] : -5485.637\n", " ARIMA(3,0,0)(2,1,0)[12] with drift : -5483.787\n", " ARIMA(3,0,1)(0,1,0)[12] : -5226.009\n", " ARIMA(3,0,1)(0,1,0)[12] with drift : -5223.999\n", " ARIMA(3,0,1)(0,1,1)[12] : Inf\n", " ARIMA(3,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(3,0,1)(1,1,0)[12] : -5418.923\n", " ARIMA(3,0,1)(1,1,0)[12] with drift : -5417.005\n", " ARIMA(3,0,2)(0,1,0)[12] : -5224.479\n", " ARIMA(3,0,2)(0,1,0)[12] with drift : -5222.461\n", " ARIMA(4,0,0)(0,1,0)[12] : -5225.24\n", " ARIMA(4,0,0)(0,1,0)[12] with drift : -5223.23\n", " ARIMA(4,0,0)(0,1,1)[12] : Inf\n", " ARIMA(4,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(4,0,0)(1,1,0)[12] : -5418.561\n", " ARIMA(4,0,0)(1,1,0)[12] with drift : -5416.635\n", " ARIMA(4,0,1)(0,1,0)[12] : -5224.07\n", " ARIMA(4,0,1)(0,1,0)[12] with drift : -5222.054\n", " ARIMA(5,0,0)(0,1,0)[12] : -5225.74\n", " ARIMA(5,0,0)(0,1,0)[12] with drift : -5223.72\n", "\n", "\n", "\n", " Best model: ARIMA(0,0,3)(0,1,1)[12] \n", "\n" ] }, { "data": { "text/plain": [ "Series: it \n", "ARIMA(0,0,3)(0,1,1)[12] \n", "Box Cox transformation: lambda= -0.7557019 \n", "\n", "Coefficients:\n", " ma1 ma2 ma3 sma1\n", " 0.2236 0.1118 0.1139 -0.8862\n", "s.e. 0.0388 0.0381 0.0353 0.0249\n", "\n", "sigma^2 estimated as 1.526e-05: log likelihood=2794.9\n", "AIC=-5579.79 AICc=-5579.7 BIC=-5557.19" ] }, "execution_count": 40, "metadata": { }, "output_type": "execute_result" } ], "source": [ "ajuste_auto_bc <- auto.arima(it, stepwise=FALSE, approximation=FALSE, trace=TRUE, lambda=l)\n", "ajuste_auto_bc" ] }, { "cell_type": "code", "execution_count": 41, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "\tLilliefors (Kolmogorov-Smirnov) normality test\n", "\n", "data: ajuste_auto_bc$residuals\n", "D = 0.048052, p-value = 0.0006699\n" ] }, "execution_count": 41, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "aedea41b3fd33d3a949a3df0ccc6c4c4cb38bfc6" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "cafe9682e76d2dc63ebf3011c8615b5f6914e282" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "qqnorm(ajuste_auto_bc$residuals)\n", "qqline(ajuste_auto_bc$residuals)\n", "lillie.test(ajuste_auto_bc$residuals)\n", "hist(ajuste_auto_bc$residuals)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Melhorou as estatísticas, mas não a questão da aderência a Normal. Tem outliers..." ] }, { "cell_type": "code", "execution_count": 42, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "# A função não possui o parâmetro lambda\n", "# outliers <- tso(BoxCox(it,l),tsmethod = \"auto.arima\", args.tsmethod = list(lambda=l), types = c(\"AO\",\"LS\",\"TC\"))" ] }, { "cell_type": "code", "execution_count": 43, "metadata": { "collapsed": false, "scrolled": true }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Warning message in locate.outliers.oloop(y = y, fit = fit, types = types, cval = cval, :\n", "“stopped when 'maxit.oloop = 4' was reached”" ] } ], "source": [ "outliers <- tso(BoxCox(it,l),tsmethod = \"arima\", args.tsmethod = list(order=c(0, 0, 3), seasonal=list(order=c(0, 1, 1), period=12)), types = c(\"AO\",\"LS\",\"TC\"), maxit.iloop=50)" ] }, { "cell_type": "code", "execution_count": 44, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "Call:\n", "structure(list(method = NULL), .Names = \"method\")\n", "\n", "Coefficients:\n", " ma1 ma2 ma3 sma1 AO507\n", " 0.2195 0.1154 0.1241 -0.8833 -0.0170\n", "s.e. 0.0388 0.0378 0.0355 0.0252 0.0036\n", "\n", "sigma^2 estimated as 1.468e-05: log likelihood = 2805.82, aic = -5599.64\n", "\n", "Outliers:\n", " type ind time coefhat tstat\n", "1 AO 507 2003:03 -0.01704 -4.709" ] }, "execution_count": 44, "metadata": { }, "output_type": "execute_result" } ], "source": [ "outliers" ] }, { "cell_type": "code", "execution_count": 45, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "it[outliers$outliers$ind] <- NA" ] }, { "cell_type": "code", "execution_count": 46, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "# Substituindo os outliers por estimativas...\n", "it <- na.kalman(it, model = \"auto.arima\")" ] }, { "cell_type": "code", "execution_count": 47, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "-0.706569225958772" ] }, "execution_count": 47, "metadata": { }, "output_type": "execute_result" } ], "source": [ "l <- BoxCox.lambda(it)\n", "l" ] }, { "cell_type": "code", "execution_count": 48, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "Call:\n", "structure(list(method = NULL), .Names = \"method\")\n", "\n", "Coefficients:\n", " ar1 ma1 ma2 sma1 TC88\n", " 0.9910 -0.7919 -0.1484 -0.9594 -0.0157\n", "s.e. 0.0118 0.0388 0.0396 0.0317 0.0036\n", "\n", "sigma^2 estimated as 1.922e-05: log likelihood = 2709.07, aic = -5406.15\n", "\n", "Outliers:\n", " type ind time coefhat tstat\n", "1 TC 88 1968:04 -0.01565 -4.379" ] }, "execution_count": 48, "metadata": { }, "output_type": "execute_result" } ], "source": [ "outliers2 <- tso(BoxCox(it,l),tsmethod = \"arima\", args.tsmethod = list(order=c(1, 0, 2), seasonal=list(order=c(0, 1, 1), period=12)), types = c(\"AO\",\"LS\",\"TC\"), maxit.iloop=50)\n", "outliers2" ] }, { "cell_type": "code", "execution_count": 49, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "-0.83974114653526" ] }, "execution_count": 49, "metadata": { }, "output_type": "execute_result" } ], "source": [ "it[outliers2$outliers$ind] <- NA\n", "it <- na.kalman(it, model = \"auto.arima\")\n", "l <- BoxCox.lambda(it)\n", "l" ] }, { "cell_type": "code", "execution_count": 50, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "Call:\n", "structure(list(method = NULL), .Names = \"method\")\n", "\n", "Coefficients:\n", " ar1 ma1 ma2 sma1 TC89\n", " 0.9907 -0.7903 -0.1505 -0.9571 -0.0093\n", "s.e. 0.0119 0.0388 0.0395 0.0306 0.0023\n", "\n", "sigma^2 estimated as 8.025e-06: log likelihood = 3005.94, aic = -5999.88\n", "\n", "Outliers:\n", " type ind time coefhat tstat\n", "1 TC 89 1968:05 -0.00934 -4.052" ] }, "execution_count": 50, "metadata": { }, "output_type": "execute_result" } ], "source": [ "outliers3 <- tso(BoxCox(it,l),tsmethod = \"arima\", args.tsmethod = list(order=c(1, 0, 2), seasonal=list(order=c(0, 1, 1), period=12)), types = c(\"AO\",\"LS\",\"TC\"), maxit.iloop=50)\n", "outliers3" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Reajustando" ] }, { "cell_type": "code", "execution_count": 51, "metadata": { "collapsed": false }, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "\n", " ARIMA(0,0,0)(0,1,0)[12] : -5595.614\n", " ARIMA(0,0,0)(0,1,0)[12] with drift : -5593.661\n", " ARIMA(0,0,0)(0,1,1)[12] : -5935.103\n", " ARIMA(0,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,0)(0,1,2)[12] : -5933.397\n", " ARIMA(0,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,0)(1,1,0)[12] : -5784.9\n", " ARIMA(0,0,0)(1,1,0)[12] with drift : -5783.151\n", " ARIMA(0,0,0)(1,1,1)[12] : -5933.43\n", " ARIMA(0,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,0)(1,1,2)[12] : -5931.091\n", " ARIMA(0,0,0)(1,1,2)[12] with drift : Inf\n", " ARIMA(0,0,0)(2,1,0)[12] : -5857.596\n", " ARIMA(0,0,0)(2,1,0)[12] with drift : -5856.001\n", " ARIMA(0,0,0)(2,1,1)[12] : -5932.721\n", " ARIMA(0,0,0)(2,1,1)[12] with drift : Inf\n", " ARIMA(0,0,0)(2,1,2)[12] : Inf\n", " ARIMA(0,0,0)(2,1,2)[12] with drift : Inf\n", " ARIMA(0,0,1)(0,1,0)[12] : -5624.147\n", " ARIMA(0,0,1)(0,1,0)[12] with drift : -5622.173\n", " ARIMA(0,0,1)(0,1,1)[12] : -5965.234\n", " ARIMA(0,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,1)(0,1,2)[12] : -5964.415\n", " ARIMA(0,0,1)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,1)(1,1,0)[12] : -5806.584\n", " ARIMA(0,0,1)(1,1,0)[12] with drift : -5804.76\n", " ARIMA(0,0,1)(1,1,1)[12] : -5964.537\n", " ARIMA(0,0,1)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,1)(1,1,2)[12] : Inf\n", " ARIMA(0,0,1)(1,1,2)[12] with drift : Inf\n", " ARIMA(0,0,1)(2,1,0)[12] : -5872.149\n", " ARIMA(0,0,1)(2,1,0)[12] with drift : -5870.45\n", " ARIMA(0,0,1)(2,1,1)[12] : Inf\n", " ARIMA(0,0,1)(2,1,1)[12] with drift : Inf\n", " ARIMA(0,0,1)(2,1,2)[12] : Inf\n", " ARIMA(0,0,1)(2,1,2)[12] with drift : Inf\n", " ARIMA(0,0,2)(0,1,0)[12] : -5622.49\n", " ARIMA(0,0,2)(0,1,0)[12] with drift : -5620.509\n", " ARIMA(0,0,2)(0,1,1)[12] : -5969.018\n", " ARIMA(0,0,2)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,2)(0,1,2)[12] : -5967.759\n", " ARIMA(0,0,2)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,2)(1,1,0)[12] : -5808.343\n", " ARIMA(0,0,2)(1,1,0)[12] with drift : -5806.496\n", " ARIMA(0,0,2)(1,1,1)[12] : Inf\n", " ARIMA(0,0,2)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,2)(1,1,2)[12] : -5965.001\n", " ARIMA(0,0,2)(1,1,2)[12] with drift : Inf\n", " ARIMA(0,0,2)(2,1,0)[12] : -5875.818\n", " ARIMA(0,0,2)(2,1,0)[12] with drift : -5874.077\n", " ARIMA(0,0,2)(2,1,1)[12] : Inf\n", " ARIMA(0,0,2)(2,1,1)[12] with drift : Inf\n", " ARIMA(0,0,3)(0,1,0)[12] : -5634.082\n", " ARIMA(0,0,3)(0,1,0)[12] with drift : -5632.085\n", " ARIMA(0,0,3)(0,1,1)[12] : -5978.213\n", " ARIMA(0,0,3)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,3)(0,1,2)[12] : Inf\n", " ARIMA(0,0,3)(0,1,2)[12] with drift : Inf\n", " ARIMA(0,0,3)(1,1,0)[12] : -5816.817\n", " ARIMA(0,0,3)(1,1,0)[12] with drift : -5814.928\n", " ARIMA(0,0,3)(1,1,1)[12] : Inf\n", " ARIMA(0,0,3)(1,1,1)[12] with drift : Inf\n", " ARIMA(0,0,3)(2,1,0)[12] : -5880.987\n", " ARIMA(0,0,3)(2,1,0)[12] with drift : -5879.188\n", " ARIMA(0,0,4)(0,1,0)[12] : -5633.091\n", " ARIMA(0,0,4)(0,1,0)[12] with drift : -5631.085\n", " ARIMA(0,0,4)(0,1,1)[12] : -5977.304\n", " ARIMA(0,0,4)(0,1,1)[12] with drift : Inf\n", " ARIMA(0,0,4)(1,1,0)[12] : -5815.636\n", " ARIMA(0,0,4)(1,1,0)[12] with drift : -5813.731\n", " ARIMA(0,0,5)(0,1,0)[12] : -5631.492\n", " ARIMA(0,0,5)(0,1,0)[12] with drift : -5629.478\n", " ARIMA(1,0,0)(0,1,0)[12] : -5626.03\n", " ARIMA(1,0,0)(0,1,0)[12] with drift : -5624.052\n", " ARIMA(1,0,0)(0,1,1)[12] : -5971.74\n", " ARIMA(1,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,0)(0,1,2)[12] : Inf\n", " ARIMA(1,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(1,0,0)(1,1,0)[12] : -5810.389\n", " ARIMA(1,0,0)(1,1,0)[12] with drift : -5808.544\n", " ARIMA(1,0,0)(1,1,1)[12] : Inf\n", " ARIMA(1,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(1,0,0)(1,1,2)[12] : Inf\n", " ARIMA(1,0,0)(1,1,2)[12] with drift : Inf\n", " ARIMA(1,0,0)(2,1,0)[12] : -5875.374\n", " ARIMA(1,0,0)(2,1,0)[12] with drift : -5873.644\n", " ARIMA(1,0,0)(2,1,1)[12] : Inf\n", " ARIMA(1,0,0)(2,1,1)[12] with drift : Inf\n", " ARIMA(1,0,0)(2,1,2)[12] : Inf\n", " ARIMA(1,0,0)(2,1,2)[12] with drift : Inf\n", " ARIMA(1,0,1)(0,1,0)[12] : -5629.227\n", " ARIMA(1,0,1)(0,1,0)[12] with drift : -5627.228\n", " ARIMA(1,0,1)(0,1,1)[12] : Inf\n", " ARIMA(1,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,1)(0,1,2)[12] : Inf\n", " ARIMA(1,0,1)(0,1,2)[12] with drift : Inf\n", " ARIMA(1,0,1)(1,1,0)[12] : -5819.324\n", " ARIMA(1,0,1)(1,1,0)[12] with drift : -5817.396\n", " ARIMA(1,0,1)(1,1,1)[12] : Inf\n", " ARIMA(1,0,1)(1,1,1)[12] with drift : Inf\n", " ARIMA(1,0,1)(1,1,2)[12] : Inf\n", " ARIMA(1,0,1)(1,1,2)[12] with drift : Inf\n", " ARIMA(1,0,1)(2,1,0)[12] : -5885.866\n", " ARIMA(1,0,1)(2,1,0)[12] with drift : -5883.994\n", " ARIMA(1,0,1)(2,1,1)[12] : Inf\n", " ARIMA(1,0,1)(2,1,1)[12] with drift : Inf\n", " ARIMA(1,0,2)(0,1,0)[12] : -5628.545\n", " ARIMA(1,0,2)(0,1,0)[12] with drift : -5626.539\n", " ARIMA(1,0,2)(0,1,1)[12] : Inf\n", " ARIMA(1,0,2)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,2)(0,1,2)[12] : Inf\n", " ARIMA(1,0,2)(0,1,2)[12] with drift : Inf\n", " ARIMA(1,0,2)(1,1,0)[12] : -5817.44\n", " ARIMA(1,0,2)(1,1,0)[12] with drift : -5815.504\n", " ARIMA(1,0,2)(1,1,1)[12] : Inf\n", " ARIMA(1,0,2)(1,1,1)[12] with drift : Inf\n", " ARIMA(1,0,2)(2,1,0)[12] : -5883.839\n", " ARIMA(1,0,2)(2,1,0)[12] with drift : -5881.961\n", " ARIMA(1,0,3)(0,1,0)[12] : -5633.48\n", " ARIMA(1,0,3)(0,1,0)[12] with drift : -5631.472\n", " ARIMA(1,0,3)(0,1,1)[12] : Inf\n", " ARIMA(1,0,3)(0,1,1)[12] with drift : Inf\n", " ARIMA(1,0,3)(1,1,0)[12] : -5816.716\n", " ARIMA(1,0,3)(1,1,0)[12] with drift : -5814.777\n", " ARIMA(1,0,4)(0,1,0)[12] : -5632.07\n", " ARIMA(1,0,4)(0,1,0)[12] with drift : -5630.052\n", " ARIMA(2,0,0)(0,1,0)[12] : -5624.965\n", " ARIMA(2,0,0)(0,1,0)[12] with drift : -5622.978\n", " ARIMA(2,0,0)(0,1,1)[12] : -5976.243\n", " ARIMA(2,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(2,0,0)(0,1,2)[12] : Inf\n", " ARIMA(2,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(2,0,0)(1,1,0)[12] : -5813.072\n", " ARIMA(2,0,0)(1,1,0)[12] with drift : -5811.194\n", " ARIMA(2,0,0)(1,1,1)[12] : Inf\n", " ARIMA(2,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(2,0,0)(1,1,2)[12] : -5972.219\n", " ARIMA(2,0,0)(1,1,2)[12] with drift : Inf\n", " ARIMA(2,0,0)(2,1,0)[12] : -5879.922\n", " ARIMA(2,0,0)(2,1,0)[12] with drift : -5878.135\n", " ARIMA(2,0,0)(2,1,1)[12] : Inf\n", " ARIMA(2,0,0)(2,1,1)[12] with drift : Inf\n", " ARIMA(2,0,1)(0,1,0)[12] : -5628.116\n", " ARIMA(2,0,1)(0,1,0)[12] with drift : -5626.11\n", " ARIMA(2,0,1)(0,1,1)[12] : Inf\n", " ARIMA(2,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(2,0,1)(0,1,2)[12] : Inf\n", " ARIMA(2,0,1)(0,1,2)[12] with drift : Inf\n", " ARIMA(2,0,1)(1,1,0)[12] : -5817.417\n", " ARIMA(2,0,1)(1,1,0)[12] with drift : -5815.482\n", " ARIMA(2,0,1)(1,1,1)[12] : Inf\n", " ARIMA(2,0,1)(1,1,1)[12] with drift : Inf\n", " ARIMA(2,0,1)(2,1,0)[12] : Inf\n", " ARIMA(2,0,1)(2,1,0)[12] with drift : Inf\n", " ARIMA(2,0,2)(0,1,0)[12] : Inf\n", " ARIMA(2,0,2)(0,1,0)[12] with drift : Inf\n", " ARIMA(2,0,2)(0,1,1)[12] : Inf\n", " ARIMA(2,0,2)(0,1,1)[12] with drift : Inf\n", " ARIMA(2,0,2)(1,1,0)[12] : -5816.598\n", " ARIMA(2,0,2)(1,1,0)[12] with drift : -5814.659\n", " ARIMA(2,0,3)(0,1,0)[12] : -5632.052\n", " ARIMA(2,0,3)(0,1,0)[12] with drift : -5630.034\n", " ARIMA(3,0,0)(0,1,0)[12] : -5635.615\n", " ARIMA(3,0,0)(0,1,0)[12] with drift : -5633.61\n", " ARIMA(3,0,0)(0,1,1)[12] : Inf\n", " ARIMA(3,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(3,0,0)(0,1,2)[12] : Inf\n", " ARIMA(3,0,0)(0,1,2)[12] with drift : Inf\n", " ARIMA(3,0,0)(1,1,0)[12] : -5819.236\n", " ARIMA(3,0,0)(1,1,0)[12] with drift : -5817.316\n", " ARIMA(3,0,0)(1,1,1)[12] : Inf\n", " ARIMA(3,0,0)(1,1,1)[12] with drift : Inf\n", " ARIMA(3,0,0)(2,1,0)[12] : -5883.556\n", " ARIMA(3,0,0)(2,1,0)[12] with drift : -5881.712\n", " ARIMA(3,0,1)(0,1,0)[12] : -5633.852\n", " ARIMA(3,0,1)(0,1,0)[12] with drift : -5631.841\n", " ARIMA(3,0,1)(0,1,1)[12] : Inf\n", " ARIMA(3,0,1)(0,1,1)[12] with drift : Inf\n", " ARIMA(3,0,1)(1,1,0)[12] : -5817.198\n", " ARIMA(3,0,1)(1,1,0)[12] with drift : -5815.272\n", " ARIMA(3,0,2)(0,1,0)[12] : -5632.358\n", " ARIMA(3,0,2)(0,1,0)[12] with drift : -5630.341\n", " ARIMA(4,0,0)(0,1,0)[12] : -5633.71\n", " ARIMA(4,0,0)(0,1,0)[12] with drift : -5631.7\n", " ARIMA(4,0,0)(0,1,1)[12] : Inf\n", " ARIMA(4,0,0)(0,1,1)[12] with drift : Inf\n", " ARIMA(4,0,0)(1,1,0)[12] : -5817.205\n", " ARIMA(4,0,0)(1,1,0)[12] with drift : -5815.274\n", " ARIMA(4,0,1)(0,1,0)[12] : -5632.009\n", " ARIMA(4,0,1)(0,1,0)[12] with drift : -5629.992\n", " ARIMA(5,0,0)(0,1,0)[12] : -5632.617\n", " ARIMA(5,0,0)(0,1,0)[12] with drift : -5630.598\n", "\n", "\n", "\n", " Best model: ARIMA(0,0,3)(0,1,1)[12] \n", "\n" ] }, { "data": { "text/plain": [ "Series: it \n", "ARIMA(0,0,3)(0,1,1)[12] \n", "Box Cox transformation: lambda= -0.8397411 \n", "\n", "Coefficients:\n", " ma1 ma2 ma3 sma1\n", " 0.2070 0.1108 0.1212 -0.8823\n", "s.e. 0.0387 0.0377 0.0356 0.0253\n", "\n", "sigma^2 estimated as 8.5e-06: log likelihood=2994.15\n", "AIC=-5978.3 AICc=-5978.21 BIC=-5955.7" ] }, "execution_count": 51, "metadata": { }, "output_type": "execute_result" } ], "source": [ "# Reajustando\n", "ajuste_auto_bc_outliers <- auto.arima(it, stepwise=FALSE, approximation=FALSE, trace=TRUE, lambda=l)\n", "ajuste_auto_bc_outliers" ] }, { "cell_type": "code", "execution_count": 52, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "\tLilliefors (Kolmogorov-Smirnov) normality test\n", "\n", "data: ajuste_auto_bc_outliers$residuals\n", "D = 0.045034, p-value = 0.002032\n" ] }, "execution_count": 52, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "53f57e35d3f2fba0eb5a603005793da8afa8db06" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "a2e5b2ae6900f3cbb2d1535b5da35c2b1a83ebd3" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "qqnorm(ajuste_auto_bc_outliers$residuals)\n", "qqline(ajuste_auto_bc_outliers$residuals)\n", "lillie.test(ajuste_auto_bc_outliers$residuals)\n", "hist(ajuste_auto_bc_outliers$residuals)" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "# Testando o modelo" ] }, { "cell_type": "code", "execution_count": 53, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "\n", "\n", "\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\t\n", "\n", "
lagsstatisticdfp-value
5 1.062512 1 0.3026430
10 3.328105 6 0.7666811
20 12.33643616 0.7205137
25 15.10461521 0.8176662
39 22.61880435 0.9473957
35 20.70203431 0.9195489
40 23.48525836 0.9463749
45 24.75732241 0.9788522
50 27.18184646 0.9877001
\n" ] }, "execution_count": 53, "metadata": { }, "output_type": "execute_result" } ], "source": [ "LjungBox(ajuste_auto_bc_outliers, lag=c(5, 10, 20, 25, 39, 35, 40, 45, 50), season=12) # Testa se existe autocorrelação entre os dados considerando a sazonalidade" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "### Previsões para os próximos 3 meses" ] }, { "cell_type": "code", "execution_count": 54, "metadata": { "collapsed": false, "scrolled": true }, "outputs": [ { "data": { "text/plain": [ " Point Forecast Lo 80 Hi 80 Lo 95 Hi 95\n", "Aug 2018 26.33039 24.87485 27.95070 24.16290 28.88467\n", "Sep 2018 28.98351 27.21961 30.96980 26.36381 32.12565\n", "Oct 2018 30.73892 28.76825 32.97308 27.81667 34.28051" ] }, "execution_count": 54, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "\n", "\n", "\n", "\t\n", "\n", "
MERMSEMAEMPEMAPEMASEACF1
Training set0.1303623 1.229289 0.9624951 0.2885044 3.526468 0.7180662 0.04779255
\n" ] }, "execution_count": 54, "metadata": { }, "output_type": "execute_result" } ], "source": [ "previsoes_auto_bc_outliers <- forecast(ajuste_auto_bc_outliers, h=3)\n", "previsoes_auto_bc_outliers\n", "accuracy(previsoes_auto_bc_outliers)" ] }, { "cell_type": "code", "execution_count": 55, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "Call:\n", "lmRob(formula = ajuste_auto_bc_outliers$x ~ ajuste_auto_bc_outliers$fitted, \n", " data = T)\n", "\n", "Residuals:\n", " Min 1Q Median 3Q Max \n", "-3.915993 -0.759128 0.005429 0.725391 4.303154 \n", "\n", "Coefficients:\n", " Estimate Std. Error t value Pr(>|t|) \n", "(Intercept) 0.59538 0.57686 1.032 0.302 \n", "ajuste_auto_bc_outliers$fitted 0.98361 0.02091 47.032 <2e-16 ***\n", "---\n", "Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n", "\n", "Residual standard error: 1.133 on 689 degrees of freedom\n", "Multiple R-Squared: 0.5979 \n", "\n", "Test for Bias:\n", " statistic p-value\n", "M-estimate 0.1828 9.127e-01\n", "LS-estimate 39.6429 2.464e-09" ] }, "execution_count": 55, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "png: 2" ] }, "execution_count": 55, "metadata": { }, "output_type": "execute_result" }, { "data": { "text/html": [ "png: 2" ] }, "execution_count": 55, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "212d1e66231e95c8b140b4ab09a425411407e8ad" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "36e3c259fb50f14354c54642d92e84c0c722f940" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "mod2 = lmRob(ajuste_auto_bc_outliers$x ~ ajuste_auto_bc_outliers$fitted, data = T)\n", "summary(mod2)\n", "\n", "options(repr.plot.width=8, repr.plot.height=8)\n", "bin <- hexbin(ajuste_auto_bc_outliers$fitted, ajuste_auto_bc_outliers$x, xbins=30)\n", "\n", "jpeg(\"g15.jpg\", width = 1000, height = 700)\n", "plot(bin)\n", "dev.off()\n", "plot(bin)\n", "\n", "options(repr.plot.width=16, repr.plot.height=8)\n", "jpeg(\"g16.jpg\", width = 1000, height = 700)\n", "ts.plot(ajuste_auto_bc_outliers$fitted,ajuste_auto_bc_outliers$x, lty=c(2,1), gpars = list(col = c(\"blue\", \"black\")))\n", "dev.off()\n", "ts.plot(ajuste_auto_bc_outliers$fitted,ajuste_auto_bc_outliers$x, lty=c(2,1), gpars = list(col = c(\"blue\", \"black\")))" ] }, { "cell_type": "code", "execution_count": 56, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "png: 2" ] }, "execution_count": 56, "metadata": { }, "output_type": "execute_result" }, { "data": { "image/svg+xml": "656b5436f7f8a4514b9cbc16ff50685202bf6fbf" }, "metadata": { "image/svg+xml": { "isolated": true } }, "output_type": "execute_result" } ], "source": [ "jpeg(\"g17.jpg\", width = 1000, height = 700)\n", "plotarimapred(it, ajuste_auto_bc_outliers, xlim=c(2010, 2019))\n", "dev.off()\n", "plotarimapred(it, ajuste_auto_bc_outliers, xlim=c(2010, 2019))" ] }, { "cell_type": "markdown", "metadata": { "collapsed": false }, "source": [ "## Bootstrap" ] }, { "cell_type": "code", "execution_count": 57, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "AIC_Boot <- function(ts) {\n", " aux <- Arima(ts, lambda=l, order=c(1, 0, 0), seasonal=list(order=c(0, 1, 1), period=12))\n", " return(aux$aic)\n", "}" ] }, { "cell_type": "code", "execution_count": 58, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "aux <- Arima(it, lambda=l, order=c(1, 0, 0), seasonal=list(order=c(0, 1, 1), period=12))" ] }, { "cell_type": "code", "execution_count": 59, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/html": [ "-5971.77551730091" ] }, "execution_count": 59, "metadata": { }, "output_type": "execute_result" } ], "source": [ "AIC_Boot(it)" ] }, { "cell_type": "code", "execution_count": 60, "metadata": { "collapsed": false }, "outputs": [ ], "source": [ "# Bootstrap estacionário com bloco de 20\n", "boot_ajuste_auto_bc_outliers <- tsboot(it, AIC_Boot, R = 1000, l = 20, sim = \"geom\")" ] }, { "cell_type": "code", "execution_count": 61, "metadata": { "collapsed": false }, "outputs": [ { "data": { "text/plain": [ "\n", "STATIONARY BOOTSTRAP FOR TIME SERIES\n", "\n", "Average Block Length of 20 \n", "\n", "Call:\n", "tsboot(tseries = it, statistic = AIC_Boot, R = 1000, l = 20, \n", " sim = \"geom\")\n", "\n", "\n", "Bootstrap Statistics :\n", " original bias std. error\n", "t1* -5971.776 673.5559 48.29323" ] }, "execution_count": 61, "metadata": { }, "output_type": "execute_result" } ], "source": [ "boot_ajuste_auto_bc_outliers" ] }, { "cell_type": "code", "execution_count": 62, "metadata": { "collapsed": false }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "Warning message in boot.ci(boot_ajuste_auto_bc_outliers):\n", "“bootstrap variances needed for studentized intervals”" ] }, { "name": "stderr", "output_type": "stream", "text": [ "Warning message in boot.ci(boot_ajuste_auto_bc_outliers):\n", "“BCa intervals not defined for time series bootstraps”" ] }, { "data": { "text/plain": [ "BOOTSTRAP CONFIDENCE INTERVAL CALCULATIONS\n", "Based on 1000 bootstrap replicates\n", "\n", "CALL : \n", "boot.ci(boot.out = boot_ajuste_auto_bc_outliers)\n", "\n", "Intervals : \n", "Level Normal Basic Percentile \n", "95% (-6740, -6551 ) (-6731, -6548 ) (-5395, -5213 ) \n", "Calculations and Intervals on Original Scale" ] }, "execution_count": 62, "metadata": { }, "output_type": "execute_result" } ], "source": [ "boot.ci(boot_ajuste_auto_bc_outliers)" ] } ], "metadata": { "kernelspec": { "display_name": "R (R-Project)", "language": "r", "name": "ir" }, "language_info": { "codemirror_mode": "r", "file_extension": ".r", "mimetype": "text/x-r-source", "name": "R", "pygments_lexer": "r", "version": "3.4.4" } }, "nbformat": 4, "nbformat_minor": 0 }