Predictive regressions: sensitivity analysis

Packages & parameters

We start by loading the required packages & functions.

rm(list = ls())
library(tidyverse)
library(Rcpp)
library(mvtnorm)
library(reshape2)
library(cowplot)
library(latex2exp)            # For TeX in ggplot
library(viridis)
sourceCpp("c_functions.cpp")
source("R_functions.R")
#folder <- "./sensitivity"

Then, we set some default parameter values. Note: because T is reserved in R, we use TT instead.

Below, we create the grid of parameters that we will test.

a_x <- 0           # Constant in the AR process
a_y <- 0           # Constant in the AR process
n_sim <- 10^6      # Number of simulations

r_x <- c(0.1, 0.5, 0.9, 0.99)
r_y <- c(0.1, 0.5, 0.9)
k   <- c(1, 3, 6, 12, 24)

#r_y <- c(0.1, 0.3, 0.5, 0.7, 0.9, 0.96)
cor <- c(0, 0.5, 0.8)
TT  <- c(6, 12, 24, 36, 60, 84, 120)
pars <- expand.grid(r_x, r_y, cor, k, TT)
r_x <- pars[,1]
r_y <- pars[,2]
cor <- pars[,3]
k <- pars[,4]
TT  <- pars[,5]

Functions & calls

In the tests below, all processes are built so that their variance is equal to one.
The chunk below evaluates the MSE for 8,820 combinations of parameters, it is therefore very time-consuming (hours of CPU), notably because we require one million simulations. Reducing to 100,000 yields more reasonable times. We split the tasks along the k axis.

quad_err <- function(r_x, r_y, cor, TT, k){
    #print(c(r_x, r_y, cor, TT, k))
    quad_err_C_2d2(a_x, a_y, r_x, r_y, 1, 1, cor, TT, k, n_sim)
}

results <- pmap(list(r_x, r_y, cor, TT, k), quad_err) %>% 
    melt() %>% bind_cols(cbind(pars,k))
colnames(results) <- c("MSE", "ID", "r_x", "r_y", "cor", "TT", "k")
save(results, file = "results.RData")

results <- results %>%
    mutate(R2 = 1 - (1-r_y^2) * MSE)

Plot results of sensitivity

Below, we study the sentitivity of the MSE to changes in the parameters set and provide different graphical illustrations.

When k and r_y are independent

First, for k=1.

results %>% 
    filter(k == 1, r_x %in% c(0.1,0.5,0.9,0.99), r_y %in% c(0.1,0.5,0.9), cor %in% c(0,0.5,0.8)) %>%
    ggplot(aes(x = TT, y = R2, color = as.factor(r_x))) + geom_line() +
    facet_grid(cor ~ r_y, scales = "free") +
    # facet_wrap(cor ~ r_y , scales = "free", ncol=3) +
    # theme(strip.text.x = element_blank()) +
    labs(color = TeX("$\\textbf{\\rho_x}$")) +
    xlab(TeX("\\textbf{T}")) + ylab(TeX("\\textbf{R squared}")) +
    theme(legend.position = c(0.2, 0.7)) +
    scale_color_viridis(begin = 0.01, end = 0.91,        # Color management
                        option = "plasma", discrete = T,
                        direction = -1)
ggsave("sensi_k_01.pdf", width = 6, height = 4)    

Then, for k=12.

results %>% 
    filter(k == 12, r_x %in% c(0.1,0.5,0.9,0.99), r_y %in% c(0.1,0.5,0.9), cor %in% c(0,0.5,0.8)) %>%
    ggplot(aes(x = TT, y = R2, color = as.factor(r_x))) + geom_line() +
    facet_grid(cor ~ r_y, scales = "free") +
    # facet_wrap(cor ~ r_y , scales = "free", ncol=3) +
    # theme(strip.text.x = element_blank()) +
    labs(color = TeX("$\\textbf{\\rho_x}$")) +
    xlab(TeX("\\textbf{T}")) + ylab(TeX("\\textbf{R squared}")) +
    theme(legend.position = c(0.92, 0.85)) +
    scale_color_viridis(begin = 0.01, end = 0.91,        # Color management
                        option = "plasma", discrete = T,
                        direction = -1)
ggsave("sensi_k_12.pdf", width = 6, height = 4)    

When k and r_y are linked (memory tradeoff)

We proceeed to a new round of calls.

n_sim <- 10^6
cor <- 0.5 #c(0, 0.5, 0.8)
TT  <- c(1500) # c(6, 12, 24, 36, 60, 84, 120, 150, 200, 250, 300, 500, 1000)
r_x <- 0.9#c(0.1, 0.7, 0.99)
k   <- c(1, 3, 6, 12, 24)#c(1, 2, 3, 4, 5, 6, 12, 24, 36)
pars <- expand.grid(r_x, cor, k, TT)
r_x <- pars[,1]
cor <- pars[,2]
k <- pars[,3]
TT  <- pars[,4]
r_y <- 1-0.9/k
quad_err <- function(r_x, cor, TT, k){
    quad_err_C_2d2(a_x = 0, a_y = 0, 
                   r_x = r_x, r_y = 1 - 0.9/k, 
                   s_x = sqrt(1-r_x^2), s_y = sqrt(1-(1 - 0.9/k)^2), 
                   cor = cor, TT = TT, k = k, n_sim = n_sim)
}

res_sensi <- pmap(list(r_x, cor, TT, k), quad_err) %>% 
    melt() %>% bind_cols(cbind(pars,r_y))
colnames(res_sensi) <- c("value", "ID", "r_x", "cor", "k", "TT", "r_y")

# save(res_sensi, file = "sensi2.RData")


res_sensi <- res_sensi %>%
    mutate(R2 = 1 - (1-r_y^2) * value)

res_sensi %>% 
    filter(TT < 150, cor == 0.5, r_x == 0.9) %>%
    #filter(cor %in% c(0, 0.5, 0.8), r_x %in% c(0.1, 0.7, 0.99), k %in% c(1, 3, 6, 12, 24, 36)) %>%
    ggplot(aes(x = k, y = R2, color = as.factor(TT))) + geom_line() +
    #facet_grid(cor ~ r_x, scales = "free") +
    #facet_wrap(cor ~ r_x, scales = "free", ncol=3) +
    # theme(strip.text.x = element_blank()) +
    labs(color = TeX("$T$")) #+
    # theme(legend.position = c(0.78, 0.84)) 
ggsave("mem.pdf")

Saving 7.29 x 4.51 in image

res_sensi %>% 
    filter(k %in% c(1, 3, 6, 12, 24), TT < 150) %>%
    #filter(cor %in% c(0, 0.5, 0.8), r_x %in% c(0.1, 0.7, 0.99), k %in% c(1, 3, 6, 12, 24, 36)) %>%
    ggplot(aes(x = TT, y = R2, color = as.factor(r_x))) + geom_line() +
    facet_grid(cor ~ r_y, scales = "free") +
    scale_color_viridis(begin = 0.01, end = 0.91,        # Color management
                        discrete = T, option = "plasma", 
                        direction = -1) +
    #facet_wrap(cor ~ r_x, scales = "free", ncol=3) +
    # theme(strip.text.x = element_blank()) +
    labs(color = TeX("$r_x$")) #+
ggsave("sensi.pdf", width = 6, height = 4)