Paquetes a utilizar.

library(readxl)
library(ggplot2)
require(tseries)
require(forecast)
require(timsac)
require(changepoint)

Descripción de la base de datos

El dióxido de carbono es un gas sin color ni olor conformado por un átomo de carbono y dos de oxígeno que se desprende de la respiración de seres vivos, de la combustión de carbono, e instalaciones eléctricas. Las emisiones sin control del dióxido de carbono (CO2) es uno de los principales causantes del calentamiento global y para evitar mas daños es importante reducir sus emisiones, y que esto sea responsabilidad de todos los paises y habitantes del mundo ya que es algo, que nos beneficia e impacta a todos.

Por lo anterior, trabajamos con una base de datos de https://ourworldindata.org/explorers/co2 la cual contiene las emisiones anuales de CO2 (per capita) por país en diferentes años. La base de datos contiene variables como: entidad (país), código del pais, años y emision anual de CO2. En nuestro caso, seleccionamos el país de Colombia y se tienen 100 datos desde el año 1921 al 2020, lo que se busca es identificar a lo largo de los años como ha sido las emisiones de dióxido de carbono y si estas han ido disminuyendo como tambien, predecir a futuro como seran las emisiones para tener un control sobre ellas.

library(readxl)
CO2 <- read_excel("C:/Users/Lenovo/Downloads/SEPTIMOSEMESTRE/MMF/Co2Colombia.xlsx", 
    col_types = c("text", "text", "text", 
        "numeric"))

head(CO2)
## # A tibble: 6 x 4
##   Entity   Code  Year  `Annual CO2`
##   <chr>    <chr> <chr>        <dbl>
## 1 Colombia COL   1921        0.0045
## 2 Colombia COL   1922        0.0213
## 3 Colombia COL   1923        0.0272
## 4 Colombia COL   1924        0.0281
## 5 Colombia COL   1925        0.0615
## 6 Colombia COL   1926        0.115

Seleccionamos la columna de año y CO2 Anual (per capita) las cuales seran esenciales para el trabajo a realizar.

CO2 <- CO2[,c(3:4)]

A continuación, convertimos la columna de Año en tipo Date ya que de esta manera, se facilita el manejo de fechas.

str(CO2)
## tibble [100 x 2] (S3: tbl_df/tbl/data.frame)
##  $ Year      : chr [1:100] "1921" "1922" "1923" "1924" ...
##  $ Annual CO2: num [1:100] 0.0045 0.0213 0.0272 0.0281 0.0615 ...
CO2$Year <- as.Date(CO2$Year,format="%Y")
str(CO2)
## tibble [100 x 2] (S3: tbl_df/tbl/data.frame)
##  $ Year      : Date[1:100], format: "1921-10-11" "1922-10-11" ...
##  $ Annual CO2: num [1:100] 0.0045 0.0213 0.0272 0.0281 0.0615 ...

Teniendo en cuenta que los datos utilizados estan distribuidos regularmente, es decir que tenemos un dato cada año desde 1921 hasta el 2020 entonces, en el siguiente fragmento de código definimos el objeto de la serie de tiempo ts, teniendo en cuenta que los datos son anuales e inician en el año 1921.

CO2.ts <- ts(CO2$`Annual CO2`, start = c(1921,1), frequency=1)
CO2.ts
## Time Series:
## Start = 1921 
## End = 2020 
## Frequency = 1 
##   [1]     0.0045     0.0213     0.0272     0.0281     0.0615     0.1146
##   [7]     0.1089     0.1739     0.1627     0.0806     0.0400     0.0096
##  [13]     0.0629     0.0853     0.2744     0.3133     0.3087     0.3281
##  [19]     0.4983     0.4394     0.3860     0.4214     0.3689     0.4577
##  [25]     0.4225     0.4570     0.4556     0.3942     0.4403     0.6258
##  [31]     0.6640     0.6833     0.8282     0.7194     0.8671     0.8876
##  [37]     0.9329     0.9263 10335.0000 10207.0000 10983.0000 11361.0000
##  [43] 12044.0000  1193.0000 12207.0000 12176.0000 12445.0000  1303.0000
##  [49] 13378.0000 13208.0000 13762.0000 13957.0000 14607.0000 15462.0000
##  [55]  1488.0000 15434.0000 15594.0000  1614.0000 16823.0000 16421.0000
##  [61] 16097.0000 16255.0000 17138.0000 16604.0000 16053.0000 15947.0000
##  [67] 16068.0000 16346.0000 16265.0000 17188.0000 16782.0000 17865.0000
##  [73] 18049.0000 18678.0000 16138.0000 16121.0000 17097.0000 16984.0000
##  [79] 14042.0000 14196.0000 13896.0000 13561.0000 13763.0000 12969.0000
##  [85] 14083.0000 14353.0000 13759.0000  1506.0000 16214.0000 16867.0000
##  [91] 16661.0000 17304.0000 18804.0000 18885.0000 19353.0000 19355.0000
##  [97] 17644.0000 17092.0000 17962.0000 17512.0000

Grafico de la serie de tiempo

autoplot(CO2.ts)+
  labs(title="Serie de tiempo",
       x="Año",
       y="CO2",
       color="#00a0dc")+
    theme_grey()

Verificamos que la serie sea ESTACIONARIA por medio de la pruebaDicker-Fuller ‘adf.test’.

\(H_0\) = La serie es no estacionaria.

\(H_1\) = La serie es estacionaria

adf.test(CO2.ts, alternative = c("stationary", "explosive"))
## 
##  Augmented Dickey-Fuller Test
## 
## data:  CO2.ts
## Dickey-Fuller = -2.3584, Lag order = 4, p-value = 0.4282
## alternative hypothesis: stationary

Dado que p-value es mayor a 0.05 entonces, la serie no es estacionaria por tanto, identificaremos cuantas veces se debe realizar diferenciación para hacer que la serie sea estacionaria.

ndiffs(CO2.ts)
## [1] 1

Nos dice que se debe diferenciar una vez esta, la llamaremos dif.CO2.ts

dif.CO2.ts <- diff(CO2.ts)

Con la diferenciación realizada anteriormente entonces, realizaremos nuevamente la prueba Dicker-Fuller para mirar si la serie es estacionaria:

adf.test(dif.CO2.ts, alternative = c("stationary", "explosive"))
## 
##  Augmented Dickey-Fuller Test
## 
## data:  dif.CO2.ts
## Dickey-Fuller = -5.393, Lag order = 4, p-value = 0.01
## alternative hypothesis: stationary

Dado que p-value es menor a 0.05 se rechaza \(H_0\) entonces, ahora la serie es estacionaria. Graficamente se veria asi:

autoplot(dif.CO2.ts)+
  labs(title="Serie de tiempo con diferenciación",
       x="Año",
       y="CO2",
       color="#00a0dc")+
    theme_grey()

Funciones ACF Y PACF

par(mfrow=c(1,2))
acf(dif.CO2.ts)
pacf(dif.CO2.ts)

En cuanto al ACF, se puede mencionar que la mayoría de sus valores no son significativos, pues se encuentran dentro de las lineas y azules, pero el segundo valor, muestra una correlación negativa. En cuanto al PACF, se puede decir que sus valores no son muy significativos, pero, los primeros 2 rezagos si son más significativas, y presentan una correlación negativa

Puntos de cambio

Se determinan los puntos de cambio de la serie de tiempo los cuales permiten saber en qué momentos hubo cambios drásticos entre los años, para luego investigar cuales son sus causas.En este caso, el punto de cambio en la media es a partir de 2007.

Pc<-cpt.mean(dif.CO2.ts,method = "AMOC") 
cpts(Pc)
## [1] 87
plot(Pc)

Modelo ARIMA y Residuos

A continuación, por medio de la función auto.arima() se realizara un ajuste automático de un modelo ARIMA a la serie de tiempo que estamos trabajando, en este caso se obtiene ARIMA(0,0,1).

modelo<-auto.arima(dif.CO2.ts)
modelo
## Series: dif.CO2.ts 
## ARIMA(0,0,1) with non-zero mean 
## 
## Coefficients:
##           ma1      mean
##       -0.7098  188.3713
## s.e.   0.0691  102.8812
## 
## sigma^2 = 12063055:  log likelihood = -946.95
## AIC=1899.89   AICc=1900.14   BIC=1907.68
residuo <- residuals(modelo)
acf(residuo)

Como el 100% de las lineas están en el rango entre \(\frac{\pm2}{\sqrt(100)}\) (-0.2 a 0.2) se puede concluir que los residuos tienen un comportamiento de ruido blanco, por lo que el modelo arima es el adecuado para realizar predicciones.

Predicciones

Seguidamente, se realizan las predicciones de emisiones de CO2 para los próximos 5 años y se visualizan graficamente teniendo de teniendo en cuenta el modelo generado anteriormente. Se puede observar que se estima un solo valor para los 5 años, lo cual puede no ser cierto por las diferentes condiciones que se puedan presentar. in

pred<-forecast(modelo,h=5)
pred
##      Point Forecast     Lo 80    Hi 80     Lo 95    Hi 95
## 2021       864.4666 -3586.607 5315.540 -5942.863 7671.796
## 2022       188.3713 -5269.945 5646.688 -8159.404 8536.146
## 2023       188.3713 -5269.945 5646.688 -8159.404 8536.146
## 2024       188.3713 -5269.945 5646.688 -8159.404 8536.146
## 2025       188.3713 -5269.945 5646.688 -8159.404 8536.146
plot(forecast(modelo))

Descomposición de la serie de tiempo

La descomposición de la serie de tiempo mediante el comando decompose() no se pudo realizar puesto que, para usar esta función la serie de tiempo debe tener por lo menos 2 periodos. Sin embargo, esta función permite descomponer una serie temporal en sus componentes estacional, de tendencia y de ruido.

#des <- decompose(dif.CO2.ts)

Conclusiones

Finalmente, después del estudio del dióxido de carbono en Colombia, se puede concluir que para realizar un pronóstico aproximado de los valores en el tiempo, es necesario un ARIMA (0,0,1). El análisis de la serie de tiempo y sus residuos permitió elegir este modelo, y pronosticar los valores de los próximos 5 años, pero se debe tener en cuentra que los modelos nunca podrán predecir con certeza lo que pasará en un futuro, pues están sujetos a muchos eventos aleatorios.