Previsão do número de homicídios no Rio de Janeiro (Sim, eu voltei a postar!!!)

Um pequeno disclaimer: Depois de muito relutar, estou voltando com as atividades do blog. Isso significa de maneira prática, que apesar do clipping ser o grande driver dos acessos, sinto que há uma lacuna em relação a posts técnicos aqui na blogosfera brasileira. Dessa forma, vou fazer mais posts técnicos com assuntos ligados à Machine Learning, Data Mining e afins.

Disclaimer feito, bora pro post original.

*******

Fala pessoal, tudo bem?

Um dos principais problemas em utilizar processadores numéricos como o Excel ou mesmo o LibreOffice é que se perde um pouco a flexibilidade de parametrizar os modelos de previsão.

Além disso, há um problema grave nesses processadores de reprodutibilidade, dado que pode haver inúmeras versões, planilhas voando pra todo o lado via e-mail, e por aí vai.

Nesse post eu vou utilizar dados reais da quantidade de homicídios que houveram no Rio de Janeiro de 1991 até 2011.

Os dados foram extraídos e consolidados pelo CESeC – Centro de Estudos de Segurança e Cidadania e o link está no final do script.

Antes de mais nada, vamos carregar o package Forecast:

[code language=”r”] # Load library if(!require(forecast)){ install.packages(“forecast”) library(forecast) } [/code]

Depois disso, vamos realizar a carga dos dados. Nesse caso eu fiz uma transformação na planilha original e passei os dados para csv.

[code language=”r”] # Dataset rj_murder <- read.csv(“https://raw.githubusercontent.com/fclesio/learning-space/master/Datasets/rj-homicidios-1991-2011.csv”) [/code] Vamos dar uma olhada nos dados para nos certificar que tudo está OK. [code language=”r”] # Check data head(rj_murder) [/code]

Com os dados OK, para incluir um conjunto de dados no package de Forecast e ARIMA, temos que antes de mais nada, pegar todos os datapoints e construir uma série através de um objeto de timeseries como no bloco de código abaixo:

[code language=”r”] # Transform in time series dataset all metrics inside the dataset rj_murder_total <- ts(rj_murder$total, start=c(1991, 1), end=c(2011, 12), frequency=12) rj_murder_capital <- ts(rj_murder$capital, start=c(1991, 1), end=c(2011, 12), frequency=12) rj_murder_baixada <- ts(rj_murder$baixada, start=c(1991, 1), end=c(2011, 12), frequency=12) rj_murder_interior <- ts(rj_murder$interior, start=c(1991, 1), end=c(2011, 12), frequency=12) [/code] Para entendermos o que aconteceu, vamos destrinchar erssa função TS: ts() função para a construção de série temporal rj_murder$total conjunto de datapoints do dataframe da coluna total start=c(1991, 1) início da série é janeiro de 1991 end=c(2011, 12) o fim da série é dezebro de 2011 frequency=12 como estamos falando de dados mensais, a frequencia da série é de 12 datapoints por ano

Agora que entendemos o que fizemos, vamos dar uma olhada nos nossos objetos ts() e ver se os mesmos estão corretos.

[code language=”r”] #Plot series par(mfrow=c(2,2)) plot(rj_murder_total ,main = “Total de Homicidios no RJ de 1991-2011” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5 ) plot(rj_murder_capital ,main = “Total de Homicidios no RJ de 1991-2011 (Capital)” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5 ) plot(rj_murder_baixada ,main = “Total de Homicidios no RJ de 1991-2011 (Baixada)” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5 ) plot(rj_murder_interior ,main = “Total de Homicidios no RJ de 1991-2011 (Interior)” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5 ) [/code]     Agora que vimos que as nossas séries estão corretas, vamos criar alguns objetos usando a função auto.arima(). [code language=”r”] # Fit ARIMA models fit_rj_murder_total <- auto.arima(rj_murder_total) fit_rj_murder_capital <- auto.arima(rj_murder_capital) fit_rj_murder_baixada <- auto.arima(rj_murder_baixada) fit_rj_murder_interior <- auto.arima(rj_murder_interior) [/code]

Essa função retorna o melhor modelo ARIMA de acordo com os valores de AIC, AICc e BIC.

Mas deixando um pouco a teoria de lado, vamos ver como ficaram os ajustes dos modelos ARIMA de acordo com cada uma das séries.

[code language=”r”] #Plot ARIMA Models par(mfrow=c(2,2)) plot(forecast(fit_rj_murder_total,h=12) ,main = “Total de Homicidios no RJ de 1991-2011 \n Previsão usando ARIMA” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5) plot(forecast(fit_rj_murder_capital,h=12) ,main = “Total de Homicidios no RJ de 91-2011 (Capital) \n Previsão usando ARIMA” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5) plot(forecast(fit_rj_murder_baixada,h=12) ,main = “Total de Homicidios no RJ de 91-2011 (Baixada) \n Previsão usando ARIMA” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5) plot(forecast(fit_rj_murder_interior,h=12) ,main = “Total de Homicidios no RJ de 91-2011 (Interior) \n Previsão usando ARIMA” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5) [/code]    

Um bom resultado, mas agora vamos ver o que conseguimos fazer usando o package forecast default, com o horizonte de previsão de 12 meses:

[code language=”r”] #Forecasting using confidence intervals in 80%, 95% and 99% with 12 months ahead ahead_rj_murder_total <- forecast(rj_murder_total, level=c(80, 95, 99), h=12) ahead_rj_murder_capital <- forecast(rj_murder_capital, level=c(80, 95, 99), h=12) ahead_rj_murder_interior <- forecast(rj_murder_interior, level=c(80, 95, 99), h=12) ahead_rj_murder_baixada <- forecast(rj_murder_baixada, level=c(80, 95, 99), h=12) [/code] Modelos ajustados, vamos plotar os resultados: [code language=”r”] #Plot forecasting par(mfrow=c(2,2)) plot(ahead_rj_murder_total ,main = “Total de Homicidios no RJ de 91-2011 (Total) \n Previsão usando Forecast package” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5 ,shadecols=”oldstyle”) plot(ahead_rj_murder_capital ,main = “Total de Homicidios no RJ de 91-2011 (Capital) \n Previsão usando Forecast package” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5 ,shadecols=”oldstyle”) plot(ahead_rj_murder_baixada ,main = “Total de Homicidios no RJ de 91-2011 (Baixada) \n Previsão usando Forecast package” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5 ,shadecols=”oldstyle”) plot(ahead_rj_murder_interior ,main = “Total de Homicidios no RJ de 91-2011 (Interior) \n Previsão usando Forecast package” ,xlab = “Ano” ,ylab = “Qtde Homicidios” ,type = “l” ,lwd = 2.5 ,shadecols=”oldstyle”) [/code]

Podemos ver que tanto usando o modelo ARIMA quanto com o uso do package forecast vimos que haverá uma queda no número de homicídios. Com isso podem ser elaboradas (usando mais dados) políticas públicas de segurança.

Refs: [1] - http://www.inside-r.org/packages/cran/forecast/docs/auto.arima [2] - http://www.ucamcesec.com.br/wordpress/wp-content/uploads/2011/04/RJ_2011_evolucao.xls [3] - http://www.statmethods.net/advstats/timeseries.html