laboratorio 7

7.1 Importacion de datos- ejemplo 2

library(tseries)

## Warning: package 'tseries' was built under R version 4.5.1

## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo

library(forecast)

## Warning: package 'forecast' was built under R version 4.5.1

library(ggplot2)

## Warning: package 'ggplot2' was built under R version 4.5.1

library(writexl)

## Warning: package 'writexl' was built under R version 4.5.1

uu <- "https://raw.githubusercontent.com/vmoprojs/DataLectures/master/CAEMP.DAT"
datos = read.csv(url(uu),sep=",",header=T)
d2s = ts(datos,start=c(1962,1),frequency=4)
autoplot(d2s, main="Índice de desempleo trimestral en Canada")

knitr::kable(summary(d2s)) 

	caemp
	Min. : 82.80
	1st Qu.: 93.38
	Median :102.81
	Mean :100.22
	3rd Qu.:106.71
	Max. :113.06

La serie no es estacionaria,por que tiene tendencia y ciclos que hacen que su media a lo largo del tiempo. La serie presenta un incremento exponencial del índice del desempleo el primer año y se mantiene este valor dirante los siguientes 3 años, llegando a un maximo de 113.06 y una media de 100.22.

7.2. Estacionariedad

acf(d2s,main="Función de autocorrelación")

adf.test(d2s)

## 
##  Augmented Dickey-Fuller Test
## 
## data:  d2s
## Dickey-Fuller = -2.6391, Lag order = 5, p-value = 0.3106
## alternative hypothesis: stationary

Basado en la prueba Dickey Fuller la serie no es estacionario, ademas la funcion de autocorrelacion muestra un desenso suave en los 12 primeros rezagos lo cual confirma la prueba. Por tanto se procede a eliminar la tendencia mediante una diferencia el la parte regular de la serie.

d2sd=diff(d2s)
acf(d2sd,main="función de autocorrelación (primera diferencia)")

adf.test(d2sd)

## Warning in adf.test(d2sd): p-value smaller than printed p-value

## 
##  Augmented Dickey-Fuller Test
## 
## data:  d2sd
## Dickey-Fuller = -4.0972, Lag order = 5, p-value = 0.01
## alternative hypothesis: stationary

El p-valor de la prueba de Dickey Fuller de la serie de la primera diferencia indica que la serie es estacionario, dado que este es menor que alfa=5% (MA1)

7.3. Seleción del modelo ARIMA

pacf(d2sd, main="función de autocorrelación parcial (primera diff)")

m1=arima(d2s,order = c(0,1,2))
m1$aic

## [1] 488.5535

m2=arima(d2s,order = c(1,1,0))
m2

## 
## Call:
## arima(x = d2s, order = c(1, 1, 0))
## 
## Coefficients:
##          ar1
##       0.4598
## s.e.  0.0758
## 
## sigma^2 estimated as 2.068:  log likelihood = -240.71,  aic = 485.41

m3=arima(d2s,order = c(1,1,2))
m3$aic

## [1] 489.1799

De los tres modelos presentes se escoje aquel que presente el valor mas bajo del indicador de AIC, en este caso el modelo 2 con un aic=485.41.

7.4. Pronóstico y evaluación del modelo

d2sp=predict(m2,3)
d2sp

## $pred
##          Qtr1     Qtr2     Qtr3
## 1996 92.17202 92.24426 92.27748
## 
## $se
##          Qtr1     Qtr2     Qtr3
## 1996 1.437907 2.544385 3.499883

ts.plot(d2s,d2sp$pred,col=c("pink","black"))

Box.test(m2$residuals)

## 
##  Box-Pierce test
## 
## data:  m2$residuals
## X-squared = 0.037787, df = 1, p-value = 0.8459

El pronóstico muestra un incremento en el índice de desempleo, ademas la prueba de ruido blanco muestra que el p-valor es menor a alfa, por lo que podemos concluir que los residuales del modelo 2 son ruido blanco.

7.5. Exportacion de datos

inicio=1962.
final=1996.5
fecha=seq(inicio,final,by=0.25)
length(fecha)

## [1] 139

desempleo=c(d2s,d2sp$pred)
length(desempleo)

## [1] 139

data=c(rep("real",136),rep("pronostico",3))
length(data)

## [1] 139

datosd=data.frame(fecha,desempleo,data)
write_xlsx(datosd,"datos_desempleo.xlsx")

ggplot(datosd,aes(x=fecha,y=desempleo,col=data))+
  geom_line()+theme_bw()