temp-bikes/8-bikes-dynamic-aggregation.R

In ooelrich/oscbvar: Stuff for Bayesian vector autoregressions

# First create a data table with cumulative caliper width vs lpdens.
# Next, all we need to do is to select the caliper width corresponding
# to the largest value. Note that at time t, pred_abil is the sum of
# all previous predictions.

#####################################################################
### SETUP      ######################################################
#####################################################################

library(devtools)
load_all()

library(dplyr)
library(ggplot2)
library(data.table)

dfx <- readRDS("data-raw/bikesharing-tempdata/bikes_allcw.Rds")

############################################################
### CREATE DF OF "HISTORICALLY BEST CW"                  ###
############################################################

# This is the only additional thing we should need to turn the
# aggregation function into a dynamic one. The only difference between
# the static version and the dynamic version is that we need to select
# a new/different cw at each time point. So this part just creates a
# data frame that gives the "optimal" cw at each time point. It is
# based on predictions made previously to the respective time points.

# Create a data-set that gives us a caliper width for each time point
# To re-iterate, this cw is the cw that optimizes the predictions of
# all previous predictions.

temp_list <- list()
for (i in 402:730) {
    cw_data <- dfx[
        method == "caliper_propto" & t < i,
        .(pred_abil = sum(lpdens), time = i, calw = cw),
        .(cw)
    ]
    temp_list[[i - 172]] <- cw_data
}

all_things <- do.call(rbind, temp_list)
opt_cal <- all_things[
    ,
    .SD[which.max(pred_abil)],
    .(time)
][, .(time, cw)]
colnames(opt_cal) <- c("t", "optcw")

# Merge the optimal cw with the data frame containing logscores for all
# values of cw, then for each value of t select the row for which the cw
# is equal to the optimal cw.
alll <- data.table(merge.data.frame(dfx, opt_cal, by.x = "t", by.y = "t"))

optimal_dfx <- alll[, .SD[cw == optcw | method != "caliper_propto"], by = .(t)]

# Take a look
optimal_dfx[, .(meanlpd = mean(lpdens)), .(method)]
View(optimal_dfx)
table(opt_cal$optcw)







############################################################
### GAMMA PLOTKOD FÖR DYNAMISKA SETUPEN, SPARA ETT TAG   ###
############################################################

### "NYA" KODEN BÖRJAR HÄR
dff <- rbind(dfx[dfx$t > 450, ], aggpred_data)

# FOR MAKING THE CW VS LPDENS OVER TIME THINGS
# aggpred_gdp <- dff
# aggpred_gdp$outc <- "gdp"

# aggpred_tcpi <- dff
# aggpred_tcpi$outc <- "tcpi"

aggpred_fed <- dff
aggpred_fed$outc <- "fed"

# This should be done with a lookup table
dff$cat <- 1
for (i in seq_len(nrow(dff))){
    if (dff$method[i] == "equal_wt") {
        dff$cat[i] <- 2
    }
    if (dff$method[i] == "gewisano") {
        dff$cat[i] <- 3
    }
    if (dff$method[i] == "caliper_propto") {
        dff$cat[i] <- 4
    }
}
dff$cat <- factor(dff$cat)

aggpreds <- ggplot(
        dff,
        aes(y = lpdens, x = t, group = method, col = cat)) +
    geom_line(position = position_dodge(width = 0.5), size = 0.7) +
    theme_classic() + 
    labs(
        y = "Log predictive density",
        color = "Method"
    ) +
    xlim(174, 215) +
    theme(aspect.ratio = 7/16)

# Make data frame for the labels
dff <- data.table(dff)
df_lab <- dff[t == 214] # select the last point, (t, elpd) will be (x,y)
df_lab[7, 2] <- df_lab[7, 2] - 0.1  # Move caliper version down slighlty
df_lab[6, 2] <- df_lab[6, 2] + 0.05 # Move linear pool up a smidge

# Better names
df_lab[1, 3] <- "BVAR"
df_lab[2, 3] <- "SVBVAR"
df_lab[3, 3] <- "BART"
df_lab[4, 3] <- "TVPSV"
df_lab[5, 3] <- "Equal weight"
df_lab[6, 3] <- "Linear pool"
df_lab[7, 3] <- "Local pool"

final <- aggpreds + geom_text(
    df_lab,
    mapping = aes(x = t, y = lpdens, label = method, color = cat),
    nudge_x = 2.3,
    size = 2
) + theme(legend.position = "none")

my_colors <- c("grey", RColorBrewer::brewer.pal(3, "Dark2"))
names(my_colors) <- levels(factor(c(levels(dff$cat), levels(df_lab$cat)))) 
my_scale <- scale_color_manual(name = "cat", values = my_colors)
final <- final + my_scale  

# Finally, fixing the x-axis
final <- final + scale_x_continuous(
    name = "",
    breaks = c(178, 190, 202, 214),
    labels = c("2010Q1", "2013Q1", "2016Q1", "2019Q1")
)

plottit <- sprintf("temp/final_%s.pdf", outc)
#plottit <- plot_crop(plottit) # trying to automate pdfcrop with knitr
ggsave(plottit, final)

# Then run pdfcrop fedfunds_final.pdf fedfunds_final.pdf
# in the terminal to get a version that doesn't look ridiculous