Correlacion.Rmd

---
title: "Pruebas"
output: html_notebook
---

# Análisis de series de tiempo económicas

```{r, echo=FALSE, message=FALSE}
library(tidyverse)
library(lubridate)
library(forecast)
library(ggrepel)
library(rvest)
library(plotly)
library(tseries)
library(fpp)

set.seed(1)
```


```{r, echo=FALSE, message=FALSE}
dja <- read_csv("data/DJA.csv",skip = 4)
gold <- read_csv("data/GOLD_1791-2018.csv",skip = 3)
interest_rate <- read_csv("data/INTERESTRATE_1857-2018.csv",skip = 1)
sap <- read_csv("data/SAP_1871-2018.csv", skip=1)
cpi <- read_csv("data/USCPI_1774-2018.csv", skip=4)
gdp <- read_csv("data/USGDP_1790-2018.csv", skip=2)
wage <- read_csv("data/USWAGE_1774-2018.csv", skip=3)
```

```{r, echo=FALSE, message=FALSE}
# Corrijo los encabezados de los csv, que no quedan bien en R. Paso los nombres de los encabezados originales a labels de las columnas (metadata)
change_header = function(tabla, nuevos_enc) {
  if (is.data.frame(tabla)) {
    for (i in names(tabla)) {
      attr(tabla[[i]], "label") = i
    }
    names(tabla) = nuevos_enc
  } else {
    stop("Se tiene que pasar un data frame o similar objeto al campo tabla")
  }
  return(tabla)
}
cpi           = change_header(cpi, c("year", "cpi"))
dja           = change_header(dja, c("date", "value"))
gold          = change_header(gold, c("year", "value"))
interest_rate = change_header(interest_rate, c("year", "short_term_ord", "short_term_surp", "long_term"))
sap           = change_header(sap, c("year", "avg_jan", "annual_yield", "accum_avg_jan"))
gdp           = change_header(gdp, c("year", "nominal", "real_2012_base", "gdp_deflator", "pop", "nominal_per_cap", "real_per_cap_2012_base"))
wage          = change_header(wage, c("year", "cost_unsk", "prod_work_hourly_comp"))
```

Función de autocorrelación sobre ruido blanco. Se puede ver que la serie está centrada en cero y que su autocorrelación sólo tiene una correlación alta para el _lag_ 0, mientras que es baja para el resto de los puntos.
```{r, echo=FALSE}
y <- ts(rnorm(161)) #le pongo la misma cantidad de observaciones que nuestra serie
autoplot(ts(y))
acf(y)
```


Autocorrelación sobre una función periodica (seno). Puede verse la fuerte autocorrelación entre los ciclos.

```{r, echo=FALSE}
seno <- ts(sin(pi*seq(-7, 7, length.out = 161))) 
autoplot(seno)
acf(seno)
```

## Autocorrelación sobre la tasa de largo plazo
El siguiente grafico muestra la serie de la tasa de interés de largo plazo de los EE.UU. desde el año 1857 como una serie anual.\centering
	\includegraphics[width=0.8\linewidth]{oro.png}
	\caption{Precios del oro por onza en el mercado de Nueva York} \label{fig:oro}
```{r, echo=FALSE}
ts(interest_rate$long_term, start = 1857) %>%
  autoplot() +
  labs(x = "Año", y = "Tasa de interés")
```
Hay que sacarle la tendencia a la serie

```{r, echo=FALSE}
plot(ts(interest_rate$long_term, start = 1857) %>%
  diff(.) %>% 
  acf(.))
```



```{r, echo=FALSE}
ma(ts(interest_rate$long_term, start = 1857),order = 3) %>% 
  autoplot()
ma(ts(interest_rate$long_term, start = 1857),order = 10) %>% 
  autoplot()
ma(ts(interest_rate$long_term, start = 1857),order = 50) %>% 
  autoplot()

```

En la siguente transformación de fourier de la serie de la tasa de interés de largo plazo, se puede observar una componente de largo plazo de gran amplitud.

```{r, echo=FALSE}
fft_ir = fft(interest_rate$long_term)
plot(Mod(fft_ir), type = "l", ylab = "Amplitud", xlab = "Frecuencia")
```

Para que la transformada de Fourier no quede distorsionada por la tendencia, se debe centrar la serie primero. Para ello, luego de varias pruebas, usamos una media movil de 10 lags, por capturar las tendencias de un plazo no tan largo en el periodo considerado y que obtiene un suavizado aceptable. Vemos cómo queda a continuación las diferentes medias móviles consideradas, todas sin centrar. Los órdenes más grandes de media móvil suavizan mucho la serie, aunque esto es subjetivo.

```{r, echo=FALSE}
ir_ma = interest_rate %>% 
  select(long_term) %>% 
  mutate(
    ma2 = ma(x = long_term, order = 2, centre = FALSE),
    ma5 = ma(x = long_term, order = 5, centre = FALthatSE),
    ma10 = ma(x = long_term, order = 10, centre = FALSE),
    ma15 = ma(x = long_term, order = 15, centre = FALSE),
    ma25 = ma(x = long_term, order = 25, centre = FALSE),
    ma50 = ma(x = long_term, order = 50, centre = FALSE),
    tiempo = 1:nrow(.)
    ) %>%
  select(-long_term) %>% 
  gather(key = "orden", value = "media_movil", -tiempo) %>% 
  mutate(orden = factor(x = orden,
                        levels = c("ma2", "ma5", "ma10", "ma15", "ma25", "ma50"),
                        ordered = TRUE)
         ) %>% 
ggplot(aes(x = tiempo, y = media_movil)) +
geom_line() +
facet_wrap(~orden)

ir_ma
```

Centrando la serie:

```{r, echo=FALSE}
ir_cntr = interest_rate %>% 
  select(year, long_term) %>% 
  mutate(ma10 = ma(long_term, order = 10, centre = FALSE),
         ir_cntr = long_term - ma10)

ir_cntr %>%
  ggplot(aes(x = year, y = ir_cntr)) +
  geom_line() +
  geom_abline(intercept = 0, slope = 0)
```

La serie tiene distinta varianza según el período, por lo que aplicamos logaritmo en base 10 para poder corregir la heterocedasticidad. Vemos a continuación el resultado de los tests de estacionariedad de la serie resultante (Dickey-Fuller aumentado y test KPSS), ya que es requisito para el correcto resultado de las transformaciones de Fourier.

```{r, include=FALSE}
kpss.test(na.omit(ir_cntr[["ir_cntr"]])) # Dice que el p valor es mayor al mostrado
adf.test(na.omit(ir_cntr[["ir_cntr"]])) # Dice que el p valor es menor al mostrado
```


Ambos tests resultan que la serie es estacionaria, por lo que se pueden hacer las transformaciones de Fourier pertinentes


Veamos cómo queda la descomposición de la serie de acuerdo con las transformaciones de Fourier.

```{r, echo=FALSE}
plot(ir_cntr[["ir_cntr"]],
     type = "l",
     ylab = "Valores centrados tasa de interés de largo plazo",
     xlab = "Tiempo",
     main = "Serie sin tendencia de la tasa de interés de largo plazo de los EE.UU.")
abline(0,0)

ir_cntr_fft = fft(na.omit(ir_cntr[["ir_cntr"]]))
plot(Mod(ir_cntr_fft),
     type = "l",
     xaxt = "n") # Parece haber algo de aliasing dentro de la transformada.
axis(side = 1, at=1:160)
```

Al descomponer las tendencias largas podemos apreciar las oscilaciones de menor amplitud dentro de la serie de tasas de interés. Se puede observar que las oscilaciones de mayor frecuencia tienen poca amplitud, por lo que son oscilaciones ruidosas, que quitaremos de la serie, quedándonos sólo con las oscilaciones de menor frecuencia.

```{r, echo=FALSE}
plot(as.vector(na.omit(ir_cntr[["ir_cntr"]])),
     type = "l",
     ylab = "Tasa de interés centrada",
     xlab = "Tiempo")
abline(0,0)
ir_antifft20 = Re(
    fft(
      c(rep(0,10), ir_cntr_fft[11:40], rep(0, length(ir_cntr_fft) - 80), ir_cntr_fft[40:11], rep(0,10)),
      inverse = TRUE
      ) / length(ir_cntr_fft)
    )
lines(
  ir_antifft20,
      col = "red"
  )

ir_rsq_antifft20 = 1 - (sum((na.omit(ir_cntr[["ir_cntr"]]) - ir_antifft20) ^ 2, na.rm = TRUE) /
                          sum((na.omit(ir_cntr[["ir_cntr"]]) - mean(ir_cntr[["ir_cntr"]], na.rm = TRUE)) ^ 2,
                              na.rm = TRUE))

# El R squared que se puede explicar por movimientos cíclicos es del orden del 42% de la varianza total de los datos centrados
```

En este tipo de series no parece haber series de muy largo plazo que abarquen los siglos XIX y XX, sobre todo se introduce una distorsión en el siglo XX en el período de mayor inflación de los EE.UU.

```{r, echo=FALSE}
plot(interest_rate$long_term, type = "l", ylim = c(-3, 14))
# lines(predict(lm(data = interest_rate, long_term ~ year), newdata = interest_rate[, "year"]), type = "l")

# Centrando la serie restando por la media (una regresion lineal tiene pendiente casi cero)
interest_rate$long_term_cntr = interest_rate$long_term - mean(interest_rate$long_term)

fft_ir2 = fft(interest_rate$long_term_cntr)

plot(Mod(fft_ir2), type = "l")
# Antitransformada de la primera, gran frecuencia.
# lines(Re(fft(fft_ir[1] * c(1, rep(0, length(fft_ir) - 2),1), inverse = TRUE) / nrow(interest_rate)),
#       col = "red")
# Antitransformada de la segunda frecuencia de mayor amplitud.
# lines(Re(fft(fft_ir[2] * c(0, 1, rep(0, length(fft_ir) - 3), 1), inverse = TRUE) / nrow(interest_rate)),
#       col = "blue")

# Antitransformada de las primeras frecuencias más altas
plot(interest_rate$long_term_cntr, type = "l")
fft_ir_frec_bajas = c(0,fft_ir2[2], rep(0, length(fft_ir2) - 3), fft_ir[2])
abline(h=0, lty=3)
lines(Re(fft(fft_ir_frec_bajas, inverse = TRUE) / nrow(interest_rate)),
      col = "red")
lines(Re(fft(
  c(rep(0, 2), fft_ir2[3], rep(0, length(fft_ir) - 5), 
    fft_ir[3], rep(0, 1)),
  inverse = TRUE
  ) / nrow(interest_rate)),
      col = "blue")
lines(Re(fft(
  c(rep(0, 3), fft_ir2[4], rep(0, length(fft_ir) - 7), 
    fft_ir[4], rep(0, 2)),
  inverse = TRUE
  ) / nrow(interest_rate)),
      col = "green")

# Antitransformada del resto de las frecuencias
lines(Re(fft(fft_ir2 * c(rep(0,4), rep(1, length(fft_ir) - 7), rep(0,3)), inverse = TRUE) / nrow(interest_rate)),
      col = "orange")
abline(h=0, lty=3)

# Guardo esta antitransformada del resto de las frecuencias
ir_short_freq_time = Re(fft(fft_ir2 * c(rep(0,4), rep(1, length(fft_ir) - 7), rep(0,3)), inverse = TRUE) / nrow(interest_rate))
fft_ir_short_freq_time = fft_ir2 * c(rep(0,4), rep(1, length(fft_ir) - 7), rep(0,3))

lines(Re(fft(c(fft_ir_short_freq_time[1:10],
                        rep(0, length(fft_ir_short_freq_time) - 20),
                        fft_ir_short_freq_time [c(10:1)]),
             inverse = TRUE) / nrow(interest_rate)),
      col = "blue", type = "l")
```

Revisando la capacidad explicativa de las primeras componentes de la transformada de Fourier
```{r, echo=FALSE}
plot(interest_rate$long_term_cntr, type = "l")
abline(h=0, lty=3)

fft_freqs_explicativas = fft_ir2 * c(0, rep(1, 6), rep(0, length(fft_ir2) - 13), rep(1, 6))

ir_antifft = Re(fft(fft_freqs_explicativas, inverse = TRUE) / nrow(interest_rate))
lines(ir_antifft, col = "blue", type = "l")

ir_rsq_antifft = 1 - (sum((na.omit(interest_rate$long_term_cntr) - ir_antifft) ^ 2,
                          na.rm = TRUE) /
                          sum(
                            (na.omit(interest_rate$long_term_cntr) - mean(interest_rate$long_term_cntr, na.rm = TRUE)) ^ 2, na.rm = TRUE))
ir_rsq_antifft

# En este caso explica mucho más la serie, centrando previamente solo por la media.

# Pruebo las primeras 4 componentes, como en el apartado anterior.
plot(interest_rate$long_term_cntr, type = "l")
fft_freqs_explicativas_4 = fft_ir2 * c(0, rep(1, 4), rep(0, length(fft_ir2) - 9), rep(1, 4))
ir_antifft_4 = Re(fft(fft_freqs_explicativas_4, inverse = TRUE) / nrow(interest_rate))
lines(ir_antifft_4, col = "blue", type = "l")
ir_rsq_antifft_4 = 1 - (sum((na.omit(interest_rate$long_term_cntr) - ir_antifft_4) ^ 2,
                          na.rm = TRUE) /
                          sum(
                            (na.omit(interest_rate$long_term_cntr) - mean(interest_rate$long_term_cntr, na.rm = TRUE)) ^ 2, na.rm = TRUE))
ir_rsq_antifft_4
```

### Nivel de precios

```{r, echo=FALSE}
cpi_ts = ts(cpi$cpi, start = min(cpi$year))
infl_ts = diff(cpi_ts)
```

Pruebo la cross-correlation entre la serie del oro y la inflación de los EE.UU.
```{r, echo=FALSE}
gold_ts = ts(gold$value, start = min(gold$year))
gold_ts_diff = diff(gold_ts)
cpi_ts_diff = diff(cpi_ts)

plot(gold_ts)
plot(gold_ts_diff)
ccf(gold_ts_diff, window(cpi_ts_diff,1791,2017))
```

Veamos la cómo se compone la serie del nivel general de precios de los EE.UU., con base en 1983.

```{r, echo=FALSE}
plot(cpi[["cpi"]],
     type = "l",
     ylab = "Nivel de precios",
     xlab = "Tiempo")
```

Naturalmente, la inflación describe una curva exponencial en el nivel de precios, por lo que se debe centrar la serie. Para ello, probamos varias combinaciones de medias móviles no centradas.

```{r, echo=FALSE}
medias = c(2,5,10,25,50)
medias %>%
  set_names(nm = paste0("ma", .)) %>%
  map_dfc(., 
          ma,
          x = cpi[["cpi"]],
          centre = FALSE) %>%
  mutate(tiempo = seq_len(length.out = nrow(.))) %>% 
  gather(key = "lag", value = "valor", -tiempo) %>% 
  mutate(lag = factor(lag, levels = paste0("ma", medias), ordered = TRUE)) %>% 
  ggplot(aes(x=tiempo, y=valor)) +
  geom_line() +
  facet_wrap(~lag)
```

Centrando la serie con una media movil de 10

```{r, echo=TRUE}
cpi_new_vars = . %>%
  mutate(
    cpi_ma10 = ma(cpi, order = 10, centre = FALSE),
    cpi_cntr = cpi - cpi_ma10,
    cpi_log  = log10(cpi),
    cpi_fft  = c(rep(NA, 4), fft(na.omit(cpi_cntr)), rep(NA, 5))
    )
cpi = cpi %>% cpi_new_vars

# Descomponiendo la serie logaritmica de precios antes y despues de 1900
# Suponiendo frecuencia de 50 años
plot(decompose(ts(cpi$cpi_log[cpi$year < 1900], frequency = 5), type = "additive"))
abline(h=rep(0,4))

plot(cpi_cntr ~ year,
     data = cpi,
     type = "l",
     ylab = "Índice de precios centrado",
     xlab = "Año")
abline(a = 0, b = 0)
abline(v = 1971, lty = 2, col = "red")

# Hay que dividir el análisis en dos etapas. Hasta y después del 1900

cpi_log_lm_pre  = lm(cpi_log ~ year, data = cpi[cpi$year <= 1900,])
cpi_log_lm_post = lm(cpi_log ~ year, data = cpi[cpi$year > 1900,])
plot(y = cpi$cpi_log, x = cpi$year, type = "l")
# abline(h = 1, lty = 3)
abline(v = c(1971, 1945, 1929), col = "red", lty = 3)
text(x = 1910, y = 2.0, labels = "Crisis '30")
text(x = 1945, y = 2.2, labels = "Fin SGM")
text(x = 1971, y = 2.3, labels = "Fin patrón\noro")
abline(reg = cpi_log_lm_pre, col = "blue")
abline(reg = cpi_log_lm_post, col = "blue")

cpi_log_cntr_pre = cpi$cpi_log[cpi$year <= 1900] - predict(cpi_log_lm_pre,
                                                           newdata = cpi[cpi$year <= 1900, "year"])
cpi_log_cntr_post = cpi$cpi_log[cpi$year > 1900] - predict(cpi_log_lm_post,
                                                           newdata = cpi[cpi$year > 1900, "year"])

fft_cpi_log_cntr_pre  = fft(cpi_log_cntr_pre)
fft_cpi_log_cntr_post = fft(cpi_log_cntr_post)

plot(Mod(fft_cpi_log_cntr_pre), type = "l", xlab = "Frecuencias", ylab = "Modulo", main = "Pre-1900")
plot(Mod(fft_cpi_log_cntr_post), type = "l", xlab = "Frecuencias", ylab = "Modulo", main = "Post-1900")

cpi_antifft_10_pre = Re(
    fft(
      fft_cpi_log_cntr_pre * c(rep(1,10), rep(0, length(fft_cpi_log_cntr_pre) - 19), rep(1,9)),
      inverse = TRUE
    ) / length(fft_cpi_log_cntr_pre)
  )

cpi_antifft_10_post = Re(
    fft(
      fft_cpi_log_cntr_post * c(rep(1,10), rep(0, length(fft_cpi_log_cntr_post) - 19), rep(1,9)),
      inverse = TRUE
    ) / length(fft_cpi_log_cntr_post)
  )

plot(y=c(cpi_log_cntr_pre, cpi_log_cntr_post), x=cpi$year, type="l")
lines(y=c(cpi_antifft_10_pre, cpi_antifft_10_post), x=cpi$year, col="red")
abline(v=1900,col="blue",lty=3)

cpi_rsq_antifft10 = 1 - (sum((na.omit(c(cpi_log_cntr_pre, cpi_log_cntr_post)) - c(cpi_antifft_10_pre, cpi_antifft_10_post)) ^ 2,
                          na.rm = TRUE) /
                          sum((na.omit(c(cpi_log_cntr_pre, cpi_log_cntr_post)) - mean(c(cpi_log_cntr_pre, cpi_log_cntr_post),
                                                                            na.rm = TRUE)) ^ 2,
                              na.rm = TRUE))

plot(
  x = cpi$year[cpi$year <= 1900],
  y = cpi_antifft_10_pre,
  type = "l", ylab = "CPI centrado", xlab = "Años", lwd = 2
)
lines(
  x = cpi$year[cpi$year <= 1900],
  y = Re(
    fft(
      fft_cpi_log_cntr_pre * c(rep(0,10), rep(1, length(fft_cpi_log_cntr_pre) - 19), rep(0,9)),
      inverse = TRUE
    ) / nrow(cpi)
  ),
  type = "l", col = "blue", lwd = 2
)
abline(h = 0, lty = 3, lwd = 2)

acf(Re(
    fft(
      fft_cpi_log_cntr_pre * c(rep(0,10), rep(1, length(fft_cpi_log_cntr_pre) - 19), rep(0,9)),
      inverse = TRUE
    ) / length(fft_cpi_log_cntr_pre)
  ))

# Con un filtro pasaaltos todavía queda estructura de corto plazo, pero tenemos una buena aproximación de las ondas de largo plazo.
```

Puede observarse que en los años en torno a 1971, año en que EE.UU. abandona definitivamente el patrón oro, la serie cambia completamente, por lo que habría que analizar la serie antes y después de esa fecha, ya que el índice de precios, bajo la teoría monetarista de la inflación, puede estar directamente ligada a la emisión de dinero, que desde 1971 puede darse sin respaldo en oro en los EE.UU.

Midiendo el R cuadrado que podemos obtener de estos componentes.

```{r, echo=FALSE}
cpi_rsq_fft_pre = 1 - (sum((cpi_log_cntr_pre - cpi_antifft_10_pre) ^ 2, na.rm = TRUE) /
                     sum((cpi_log_cntr_pre - mean(cpi_log_cntr_pre, na.rm = TRUE)) ^ 2, na.rm = TRUE))
cpi_rsq_fft_post = 1 - (sum((cpi_log_cntr_post - cpi_antifft_10_post) ^ 2, na.rm = TRUE) /
                     sum((cpi_log_cntr_post - mean(cpi_log_cntr_post, na.rm = TRUE)) ^ 2, na.rm = TRUE))
cpi_rsq_fft_pre
cpi_rsq_fft_post
```


```{r, echo=FALSE}
cpi_prior_1971 = cpi %>%
  filter(year < 1971) %>% 
  cpi_new_vars
cpi_post_1971  = cpi %>%
  filter(year >= 1971) %>% 
  cpi_new_vars

old_par = par("mfrow")
par(mfrow = c(1, 2))
plot(Mod(cpi_prior_1971$cpi_fft), type="l")
plot(Mod(cpi_post_1971$cpi_fft), type="l")
```

```{r, echo=FALSE}

cpi_fft1 = na.omit(cpi_prior_1971[["cpi_fft"]])
cpi_fft2 = na.omit(cpi_post_1971[["cpi_fft"]])
plot(x = index(cpi[["cpi_cntr"]]),
     y = cpi[["cpi_cntr"]],
     type = "l",
     ylab = "cpi",
     xlab = "tiempo")
plot(c(
  Re(fft(cpi_fft1 * c(
    rep(0, 20), rep(1, 10), rep(0, length(cpi_fft1) - 60), rep(1,10), rep(0,20)
  ), inverse = TRUE) / length(cpi_fft1)),
  Re(fft(cpi_fft2 * c(
    rep(1, 10), rep(0, length(cpi_fft2) - 20), rep(1,10)
  ), inverse = TRUE) / length(cpi_fft2)
  )
), col = "red", lwd = 3, type = "l", ylab = "Indice de precios")
```

## Serie del pbi real
Por lo que se puede apreciar en la figura, el pbi no parece ser un proceso estacionario, puesto que la varianza crece al limpiar la tendencia central (centrar la serie) luego de aplicar una media movil centrada de orden 10.

```{r, echo=FALSE}
real_gdp_ts = ts(gdp$real_2012_base, start = min(gdp$year))
real_gdp_ts_ma3 = ma(x = real_gdp_ts, order = 3, centre = TRUE)
real_gdp_ts_ma10 = ma(x = real_gdp_ts, order = 10, centre = TRUE)

plot(real_gdp_ts)
plot(log10(real_gdp_ts))

tend = lm(log10(real_gdp_ts) ~ gdp$year)
summary(tend) # El r2 ajustado es excelente
plot(tend) # pareciera haber algo de estructura en los residuos, no obstante estan muy cerca de cero.
# Esto es porque (observar el gráfico siguiente de la recta generada) hay zonas donde la recta
# está consistentemente por debajo o encima de la curva

plot(log10(real_gdp_ts))
abline(reg = tend)

plot(real_gdp_ts_ma3)
plot(real_gdp_ts_ma10)

# Centrando la serie con la regresión
real_gdp_ts_center = log10(real_gdp_ts) - predict(tend, newdata = log10(real_gdp_ts))

plot(real_gdp_ts_center)
abline(h = 0)

real_gdp_ts_ma10_center = real_gdp_ts - real_gdp_ts_ma10
plot(real_gdp_ts_ma10_center)
abline(a = 0, b = 0)
real_gdp_ts_ma3_center = real_gdp_ts - real_gdp_ts_ma3
plot(real_gdp_ts_ma3_center)
abline(a = 0, b = 0)
```

Una vez centrada la serie podemos buscar ciclos en la serie usando la transformada de Fourier.

```{r, echo=FALSE}
fft_real_gdp_ts_center = fft(real_gdp_ts_center)
plot(Mod(fft_real_gdp_ts_center))
```

Descomponiendo la serie
```{r, echo=FALSE}
plot(real_gdp_ts_center, type = "l")
gdp_antifft_9 = Re(fft(
  fft_real_gdp_ts_center * c(rep(1,10),
                             rep(0, length(fft_real_gdp_ts_center) - 19),
                             rep(1,9)),
             inverse = TRUE) / length(fft_real_gdp_ts_center))
lines(gdp_antifft_9, col = "red", lwd = 3)

plot(Re(fft(
  fft_real_gdp_ts_center * c(rep(1,2),
                             rep(0, length(fft_real_gdp_ts_center) - 3),
                             rep(1,1)),
             inverse = TRUE) / length(fft_real_gdp_ts_center)),
  , col = "red", lwd = 3, ylab = "PBI centrado")
lines(Re(fft(
  fft_real_gdp_ts_center * c(rep(0,2), 1,
                             rep(0, length(fft_real_gdp_ts_center) - 5),
                             1, rep(0,1)),
             inverse = TRUE) / length(fft_real_gdp_ts_center)),
  , col = "blue", lwd = 3, ylab = "PBI centrado")
lines(Re(fft(
  fft_real_gdp_ts_center * c(rep(0,3), 1,
                             rep(0, length(fft_real_gdp_ts_center) - 7),
                             1, rep(0,2)),
             inverse = TRUE) / length(fft_real_gdp_ts_center)),
  , col = "green", lwd = 3, ylab = "PBI centrado")
lines(Re(fft(
  fft_real_gdp_ts_center * c(rep(0,4), 1,
                             rep(0, length(fft_real_gdp_ts_center) - 9),
                             1, rep(0,3)),
             inverse = TRUE) / length(fft_real_gdp_ts_center)),
  , lwd = 3, ylab = "PBI centrado", col = "orange")
lines(Re(fft(
  fft_real_gdp_ts_center * c(rep(0,5), 1,
                             rep(0, length(fft_real_gdp_ts_center) - 11),
                             1, rep(0,4)),
             inverse = TRUE) / length(fft_real_gdp_ts_center)),
  , lwd = 3, ylab = "PBI centrado")

# R cuadrado explicado por las primeras 9 frecuencias.
gdp_rsq_fft = 1 - (sum((real_gdp_ts_center - gdp_antifft_9) ^ 2, na.rm = TRUE) /
                     sum((real_gdp_ts_center - mean(real_gdp_ts_center, na.rm = TRUE)) ^ 2, na.rm = TRUE))
gdp_rsq_fft
```


Revisando el histograma de esta serie de tiempo
```{r, echo=FALSE}
hist(real_gdp_ts_ma10_center, breaks="FD")
```

Su autocorrelacion: no parece ser puramente aleatoria
```{r, echo=FALSE}
acf(na.omit(real_gdp_ts_ma10_center))
```

Aunque esta idea podria no ser cierta debida a la falta de estacionariedad.
```{r}
adf.test(na.omit(real_gdp_ts_ma10_center))
kpss.test(na.omit(real_gdp_ts_ma10_center))
```

Al probar los tests de Dickey-Fuller aumentado y KPSS, ambos rechazan su hipótesis nula, pero como ambos tests tienen una hipótesis nula intercambiada no se puede llegar a una conclusión sobre la estacionariedad de la serie. Esto puede deberse a la alta volatilidad de la serie en los últimos tiempos

Vemos a continuación lo que pasa con la tasa de crecimiento del PBI de los EE.UU. en lugar de su nivel absoluto, como en el apartado anterior.

```{r, echo=FALSE}
# Tasas de crecimiento del PBI
real_gdp_ts_d = round(diff(real_gdp_ts) / real_gdp_ts[1:(length(real_gdp_ts)-1)], 4)
plot(diff(real_gdp_ts))
plot(real_gdp_ts_d)
abline(a=0,b=0)
abline(a=mean(real_gdp_ts_d),b=0, lty=2)
mean(real_gdp_ts_d)

# Centrando la serie
real_gdp_ts_d_cntr = real_gdp_ts_d - mean(real_gdp_ts_d)
```

La media de crecimiento durante todo este periodo es del 3,8% con distinta varianza a lo largo de la serie.

Realizamos tests de estacionariedad.

```{r}
adf.test(real_gdp_ts_d)
kpss.test(real_gdp_ts_d)
```
Ambos rechazan, por lo que no se pueden obtener conclusiones acerca de la estacionariedad del proceso.

Vemos ahora la transformada de Fourier para la serie.

```{r,echo=FALSE}
real_gdp_ts_d_cntr_fft = fft(real_gdp_ts_d_cntr)
plot(Mod(real_gdp_ts_d_cntr_fft))
```

La transformada de Fourier no es conclusiva acerca de ciclos de Kondratieff en la economía de los EE.UU.. Lo que se puede observar es que la economía de los EE.UU. es altamente estable en cuanto a producto, experimentando solo ciclos de baja intensidad y aleatorios alrededor de una media marcada de crecimiento. En el largo plazo (considerando la cantidad de años que estamos tomando para este análisis) las grandes crisis económicos (décadas del '30, '70, principios de los 2000, 2008) quedan diluidas en una gran regresión a la media de tasa de crecimiento de los EE.UU.

## Tasas de interés de los EE.UU.

Ahondaremos aún más en la búsqueda de ondas de Kondratieff en la serie de los EE.UU.

```{r,echo=FALSE}
ir_long_ma3  = ma(ir_long, order = 3, centre = TRUE)
ir_long_ma5  = ma(ir_long, order = 3, centre = TRUE)
ir_long_ma10 = ma(ir_long, order = 3, centre = TRUE)

ir_long_cntr = na.omit(ir_long - ir_long_ma3)
ir_long_cntr %>% autoplot() # Queda algo parecido al error de la serie
ir_long %>% autoplot()
ir_long_fft = fft(interest_rate$long_term)

plot(Mod(ir_long_fft), type = "l")
plot(Mod(ir_long_fft[1:10]), type = "l") # Filtro pasabajos

# Analizo los primeros 10 puntos
Re(fft(ir_long_fft[2:10], inverse = TRUE) / length(interest_rate$long_term)) %>% 
  plot(type = "l")
Re(fft(ir_long_fft[1:10], inverse = TRUE) / length(interest_rate$long_term)) %>% 
  plot(type = "l")
```

Las primeras 10 items de la transformación discreta de Fourier de la serie de la tasa de interés de largo plazo del tesoro de los EE.UU. muestran ciclos cortos y largos que pesan en la determinación de la tasa de interés de los EE.UU.

## Análisis de la serie del Dow Jones

El índice de Dow Jones es uno de los índices más seguidos en el mundo financiero ya que nuclea importantes empresas del mercado mundial.

La forma de construir el índice ha cambiado en el transcurso del tiempo, por lo que no se estan realmente observando las mismas empresas a lo largo de toda la serie.

```{r,echo=FALSE}
# TODO: revisar en la fuente de qué índice se está hablando cuando digo del Dow Jones, para incluir una breve descripción del índice observado en el informe.
```

```{r,echo=FALSE}
dja_start=mdy(min(dja$date))
dja_ts=ts(dja$value,start=dja_start,frequency = 365)
autoplot(dja_ts)
```

La serie del Dow Jones tiene un comportamiento exponencial en el tiempo con marcadas caídas en las épocas de crisis. Veamos las distintas medias móviles para la serie.

```{r, echo=FALSE}
dja_ts_ma3=ma(dja_ts,order=3,centre=TRUE)
dja_ts_ma5=ma(dja_ts,order=5,centre=TRUE)
dja_ts_ma10=ma(dja_ts,order=10,centre=TRUE)
dja_ts_ma30=ma(dja_ts,order=30,centre=TRUE)
dja_ts_ma60=ma(dja_ts,order=60,centre=TRUE)
dja_ts_ma180=ma(dja_ts,order=180,centre=TRUE)

plot(dja_ts_ma3)
plot(dja_ts_ma5)
plot(dja_ts_ma10)
plot(dja_ts_ma30)
plot(dja_ts_ma60)
plot(dja_ts_ma180)
```

Centrando con una media movil de 180 días

```{r,echo=FALSE}
dja_ts_cntr=dja_ts-dja_ts_ma180
autoplot(dja_ts_cntr)
```

Nótese el fuerte incremento de la varianza del índice con el progreso del tiempo.

A continuación intentamos realizar predicciones sobre el índice.

```{r,echo=FALSE}
dja_decomp = decompose(dja_ts)
autoplot(dja_decomp)
BoxCox(dja_ts_cntr, "auto") %>% head(150)
dja_ts_cntr %>% head(150)
```

```{r}
dja_modelo = auto.arima(y = dja_ts)
plot(forecast(dja_modelo,h=10), include = c(80, 95))
```