Temporary model results

#Data loading

#Data loading

source("estimating_paper.R")

location <- paste0(getwd(), "/")
locationdata <- paste0(location)

load(paste0(locationdata,"GDP_pc_data.RData"))
load(paste0(locationdata, "canada_maps.RData"))
load(paste0(locationdata, "codings.RData"))
load(paste0(locationdata,"population_ssp.RData"))
load(paste0(locationdata,"proj_climate.RData"))
load(paste0(locationdata,"weather_post_aug_v5anomaly.rda"))

#set the levels of provinces
province_codings <- province_codings[order(match(province_codings$abv, 
                                                 c("BC", "AB", "SK", "MB", "ON", "QC", 
                                                   "NB", "NS","PE", "NL", 
                                                   "YT", "NT", "NU"))),]
dat$provincecode <- factor(dat$provincecode, 
                           levels = province_codings$abv, 
                           ordered = TRUE)
province_colors <- province_codings$colorcode
names(province_colors) <- province_codings$abv
province_colors <- province_colors[province_codings$abv %in% unique(dat$provincecode)]

#weather data
season <- c("Spring", "Summer", "Fall", "Winter")
colnames(weather_anomaly)[4:9] <- c("MTemp", "Precip", "MTemp_M", "Precip_M", "MTemp_A", "Precip_A") #M means means and A means anomalies
weather_anomaly$province <- factor(weather_anomaly$province, 
                                   levels = province_codings$abv, 
                                   ordered = TRUE)
#projected climate variables by Winnie
colnames(projected_climate_A)[5] <- "Precip_proj"
colnames(projected_climate_A)[8] <- "Precip_proj_A"
projected_climate_A$Precip_proj_A <- projected_climate_A$Precip_proj_A*100

#process seasonal weather variables
weather_anomaly <- weather_anomaly %>% dplyr::select(-c(MTemp_M, Precip_M)) %>%
  pivot_wider(names_from = Season, values_from = c(MTemp,Precip,MTemp_A,Precip_A))
dat <- left_join(dat, weather_anomaly, by=c("Year"="Year", "provincecode"="province"))
dat <- dat[!dat$provincecode %in% c("NT", "YT", "NU"),]
dat$provincecode <- factor(dat$provincecode, 
                           levels = province_codings$abv[1:10], 
                           ordered = TRUE)
dat$Mining_DiffLog[(dat$provincecode=="PE" & dat$Year=="2004")] <- NA
# dat$GDP.per.capita_DiffLog[(dat$provincecode=="NL" & dat$Year=="2002")] <- NA

econ_indices <- read.csv('CMO_new.csv')
colnames(econ_indices)[1:2] <- c('Year','World.Bank')
econ_indices$World.Bank_log_diff <- c(NA, diff(log(econ_indices$World.Bank)))
econ_indices$Energy_log_diff <- c(NA, diff(log(econ_indices$Energy)))
econ_indices$Non.Energy_log_diff <- c(NA, diff(log(econ_indices$Non.energy)))
econ_indices$Precious.Metals_log_diff <- c(NA, diff(log(econ_indices$Precious.Metals)))


extreme_events <- read.csv('econ_events.csv')
extreme_events <- extreme_events[!extreme_events$Province %in% c("NT", "YT", "NU"),]
extreme_events$Province <- factor(extreme_events$Province, 
                           levels = province_codings$abv[1:10], 
                           ordered = TRUE)
extreme_events$extrem.econ <- 1

extreme_events_sum <- extreme_events %>%
  group_by(Province, Year) %>%
  summarise(
    n_events       = n(),                  
    any_econ       = max(extrem.econ),     
    .groups = "drop"
  ) %>%
  complete(
    Province, 
    Year = full_seq(1997:2017, 1),    # ensures consecutive years
    fill = list(
      n_events = 0,
      any_econ = 0
    )
  )

unemploy_rate <- read.csv("UnemploymentRate.csv")
unemploy_rate$unemploy_rate_log_diff <- c(NA, diff(log(unemploy_rate$unemploy_rate)))
# ggplot(unemploy_rate) + geom_line(aes(x=year, y=unemploy_rate, col=provincecode))

Modeling

Try several models

Linear mixed effect model with robust clustered SE for the confidence interval

#data preprocessing
Difflog_names <- colnames(dat)[grepl("_DiffLog", colnames(dat))]
dat <- dat %>%
  filter(!is.na(provincecode), Year <= 2017) %>%
  arrange(provincecode, Year) %>%            # ensure order before lag
  group_by(provincecode) %>%
  mutate(
    across(all_of(Difflog_names), ~ dplyr::lag(.x, n = 1), .names = "{.col}_lag"),
    across(all_of(Difflog_names), ~ dplyr::lag(.x, n = 2), .names = "{.col}_lag2")
  ) %>%
  ungroup()

dat <- dat %>%
  left_join(
    econ_indices %>% select(Year, World.Bank_log_diff, Energy_log_diff, Non.Energy_log_diff, Precious.Metals_log_diff),
    by = c("Year" = "Year")
  )

dat$provincecode <- factor(dat$provincecode, 
                           levels = province_codings$abv[1:10], 
                           ordered = TRUE)

dat <- dat %>%
  left_join(
    extreme_events_sum,
    by = c("Year" = "Year", "provincecode"="Province")
  )

dat <- dat %>%
  left_join(
    unemploy_rate %>% select('year', 'unemploy_rate_log_diff', 'provincecode'),
    by = c("Year" = "year", "provincecode"="provincecode")
  ) 

dat$Year_num <- dat$Year
dat$Year <- factor(dat$Year)

#Model 1: linear climate anomaly + gdp_lag + province_fixed + year_fixed
tmp1 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + Year + provincecode, 
                clusteredSE = TRUE,
                note = "linear climate anomaly + year fixed")

#Model 2: quadratic climate anomaly + gdp_lag + province_fixed + year_fixed
tmp2 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                   I(MTemp_A_Spring^2) + I(MTemp_A_Summer^2) + 
                  I(MTemp_A_Fall^2) + I(MTemp_A_Winter^2) + 
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  I(Precip_A_Spring^2) + I(Precip_A_Summer^2) + 
                  I(Precip_A_Fall^2) + I(Precip_A_Winter^2) +  
                  GDP.per.capita_DiffLog_lag + Year + provincecode, 
                clusteredSE = TRUE,
                note = "quadratic climate anomaly + year fixed")

#Model 3: linear climate anomaly + gdp_lag + province_fixed + year_fixed + interactions
tmp3 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring * Precip_A_Spring +
                  MTemp_A_Summer * Precip_A_Summer +
                  MTemp_A_Fall * Precip_A_Fall +
                  MTemp_A_Winter * Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + Year + provincecode,
                clusteredSE = TRUE,
                note = "linear climate anomaly + year fixed + interactions")

#Model 4: linear climate anomaly + gdp_lag + province_fixed + year_fixed + province trend
tmp4 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + Year + provincecode + provincecode*Year_num, 
                clusteredSE = TRUE,
                note = "linear climate anomaly + year fixed + province trend")

#Model 5: linear climate anomaly + gdp_lag + province_fixed + year_fixed + economic indices

tmp5 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + provincecode +
                  Energy_log_diff + Non.Energy_log_diff + unemploy_rate_log_diff,
                clusteredSE = FALSE,
                note = "linear climate anomaly + economic indices")

#Model 6: linear climate anomaly + gdp_lag + province_fixed + year_fixed + province trend + economic indices

tmp6 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + provincecode + provincecode*Year_num + 
                  Energy_log_diff + Non.Energy_log_diff + unemploy_rate_log_diff,
                clusteredSE = TRUE,
                note = "linear climate anomaly + province trend + economic indices")

#Model 7: linear climate anomaly + gdp_lag + province_fixed + year_fixed + economic indices + dummies for structural breaks
tmp7 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + provincecode + 
                  Energy_log_diff + Non.Energy_log_diff + unemploy_rate_log_diff + any_econ, 
                clusteredSE = TRUE,
                note = "linear climate anomaly + economic indices + dummies for structural breaks")


#Model 8: linear climate anomaly + gdp_lag + province_fixed + year_random + economic indices
tmp8 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + provincecode + (1|Year_num) +
                  Energy_log_diff + Non.Energy_log_diff + unemploy_rate_log_diff,
                clusteredSE = TRUE,
                note = "linear climate anomaly + year random + economic indices")

#Model 9: linear climate anomaly + gdp_lag + province_fixed + year_random + economic indices + dummies for structural breaks
tmp9 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + provincecode + (1|Year_num) +
                   Energy_log_diff + Non.Energy_log_diff + unemploy_rate_log_diff + any_econ, 
                clusteredSE = TRUE,
                note = "linear climate anomaly + year_random + economic indices + dummies for structural breaks")

library(DT)

a <- combine_results(tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7, tmp8, tmp9)
a <- as.data.frame(a)
t.a <- as.data.frame(t(a))
datatable(t.a, options = list(
  scrollX = TRUE,
  scrollY = 300,
  paging = FALSE
))

Robustness check

rb1 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + (1|Year) + provincecode, 
                clusteredSE = TRUE,
                note = "linear climate anomaly + year random effect")

ra <- combine_results(tmp1, rb1)
ra <- as.data.frame(ra)
datatable(ra, options = list(
  scrollX = TRUE,
  scrollY = 300,
  paging = FALSE
))

pdata <- pdata.frame(dat, index = c("provincecode", "Year"))
model_fe <- plm(GDP.per.capita_DiffLog ~ MTemp_A_Spring + MTemp_A_Summer +
                 MTemp_A_Fall + MTemp_A_Winter +
                 Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + 
                 GDP.per.capita_DiffLog_lag,
               data = pdata, model = "within", effect = "twoways")
model_re <- plm(GDP.per.capita_DiffLog ~ MTemp_A_Spring + MTemp_A_Summer +
                 MTemp_A_Fall + MTemp_A_Winter +
                Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + 
                 GDP.per.capita_DiffLog_lag,
               data = pdata, random.method="swar", model = "random")
phtest(model_fe, model_re)

## 
##  Hausman Test
## 
## data:  GDP.per.capita_DiffLog ~ MTemp_A_Spring + MTemp_A_Summer + MTemp_A_Fall +  ...
## chisq = 68.773, df = 8, p-value = 8.621e-12
## alternative hypothesis: one model is inconsistent

sensitivity analysis on weather data

load(paste0(locationdata,"GDP_pc_data.RData"))
load(paste0(locationdata, "canada_maps.RData"))
load(paste0(locationdata, "codings.RData"))
load(paste0(locationdata,"population_ssp.RData"))
load(paste0(locationdata,"proj_climate.RData"))
load(paste0(locationdata,"weather_post_aug_anomaly_weighted_by_population.rda"))

province_codings <- province_codings[order(match(province_codings$abv, 
                                                 c("BC", "AB", "SK", "MB", "ON", "QC", 
                                                   "NB", "NS","PE", "NL", 
                                                   "YT", "NT", "NU"))),]
dat$provincecode <- factor(dat$provincecode, 
                           levels = province_codings$abv, 
                           ordered = TRUE)
province_colors <- province_codings$colorcode
names(province_colors) <- province_codings$abv
province_colors <- province_colors[province_codings$abv %in% unique(dat$provincecode)]

season <- c("Spring", "Summer", "Fall", "Winter")
colnames(weather_anomaly)[4:9] <- c("MTemp", "Precip", "MTemp_M", "Precip_M", "MTemp_A", "Precip_A") #M means means and A means anomalies
weather_anomaly$province <- factor(weather_anomaly$province, 
                                   levels = province_codings$abv, 
                                   ordered = TRUE)
#projected climate variables by Winnie
colnames(projected_climate_A)[5] <- "Precip_proj"
colnames(projected_climate_A)[8] <- "Precip_proj_A"
projected_climate_A$Precip_proj_A <- projected_climate_A$Precip_proj_A*100

#process seasonal weather variables
weather_anomaly <- weather_anomaly %>% dplyr::select(-c(MTemp_M, Precip_M)) %>%
  pivot_wider(names_from = Season, values_from = c(MTemp,Precip,MTemp_A,Precip_A))
dat <- left_join(dat, weather_anomaly, by=c("Year"="Year", "provincecode"="province"))
dat <- dat[!dat$provincecode %in% c("NT", "YT", "NU"),]
dat$provincecode <- factor(dat$provincecode, 
                           levels = province_codings$abv[1:10], 
                           ordered = TRUE)
dat$Mining_DiffLog[(dat$provincecode=="PE" & dat$Year=="2004")] <- NA

econ_indices <- read.csv('CMO_new.csv')
colnames(econ_indices)[1:2] <- c('Year','World.Bank')
econ_indices$World.Bank_log_diff <- c(NA, diff(log(econ_indices$World.Bank)))
econ_indices$Energy_log_diff <- c(NA, diff(log(econ_indices$Energy)))
econ_indices$Non.Energy_log_diff <- c(NA, diff(log(econ_indices$Non.energy)))
econ_indices$Precious.Metals_log_diff <- c(NA, diff(log(econ_indices$Precious.Metals)))


extreme_events <- read.csv('econ_events.csv')
extreme_events <- extreme_events[!extreme_events$Province %in% c("NT", "YT", "NU"),]
extreme_events$Province <- factor(extreme_events$Province, 
                           levels = province_codings$abv[1:10], 
                           ordered = TRUE)
extreme_events$extrem.econ <- 1

extreme_events_sum <- extreme_events %>%
  group_by(Province, Year) %>%
  summarise(
    n_events       = n(),                  
    any_econ       = max(extrem.econ),     
    .groups = "drop"
  ) %>%
  complete(
    Province, 
    Year = full_seq(1997:2017, 1),    # ensures consecutive years
    fill = list(
      n_events = 0,
      any_econ = 0
    )
  )

unemploy_rate <- read.csv("UnemploymentRate.csv")
unemploy_rate$unemploy_rate_log_diff <- c(NA, diff(log(unemploy_rate$unemploy_rate)))


Difflog_names <- colnames(dat)[grepl("_DiffLog", colnames(dat))]
dat <- dat %>%
  filter(!is.na(provincecode), Year <= 2017) %>%
  arrange(provincecode, Year) %>%            # ensure order before lag
  group_by(provincecode) %>%
  mutate(
    across(all_of(Difflog_names), ~ dplyr::lag(.x, n = 1), .names = "{.col}_lag"),
    across(all_of(Difflog_names), ~ dplyr::lag(.x, n = 2), .names = "{.col}_lag2")
  ) %>%
  ungroup()

dat <- dat %>%
  left_join(
    econ_indices %>% select(Year, World.Bank_log_diff, Energy_log_diff, Non.Energy_log_diff, Precious.Metals_log_diff),
    by = c("Year" = "Year")
  )

dat$provincecode <- factor(dat$provincecode, 
                           levels = province_codings$abv[1:10], 
                           ordered = TRUE)

dat <- dat %>%
  left_join(
    extreme_events_sum,
    by = c("Year" = "Year", "provincecode"="Province")
  )

dat <- dat %>%
  left_join(
    unemploy_rate %>% select('year', 'unemploy_rate_log_diff', 'provincecode'),
    by = c("Year" = "year", "provincecode"="provincecode")
  ) 

dat$Year_num <- dat$Year
dat$Year <- factor(dat$Year)

#Model 1: linear climate anomaly + gdp_lag + province_fixed + year_fixed
tmp1 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + Year + provincecode, 
                clusteredSE = TRUE,
                note = "linear climate anomaly + year fixed")

#Model 2: quadratic climate anomaly + gdp_lag + province_fixed + year_fixed
tmp2 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                   I(MTemp_A_Spring^2) + I(MTemp_A_Summer^2) + 
                  I(MTemp_A_Fall^2) + I(MTemp_A_Winter^2) + 
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  I(Precip_A_Spring^2) + I(Precip_A_Summer^2) + 
                  I(Precip_A_Fall^2) + I(Precip_A_Winter^2) +  
                  GDP.per.capita_DiffLog_lag + Year + provincecode, 
                clusteredSE = TRUE,
                note = "quadratic climate anomaly + year fixed")

#Model 3: linear climate anomaly + gdp_lag + province_fixed + year_fixed + interactions
tmp3 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring * Precip_A_Spring +
                  MTemp_A_Summer * Precip_A_Summer +
                  MTemp_A_Fall * Precip_A_Fall +
                  MTemp_A_Winter * Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + Year + provincecode,
                clusteredSE = TRUE,
                note = "linear climate anomaly + year fixed + interactions")

#Model 4: linear climate anomaly + gdp_lag + province_fixed + year_fixed + province trend
tmp4 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + Year + provincecode + provincecode*Year_num, 
                clusteredSE = TRUE,
                note = "linear climate anomaly + year fixed + province trend")

#Model 5: linear climate anomaly + gdp_lag + province_fixed + year_fixed + economic indices

tmp5 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + provincecode +
                  Energy_log_diff + Non.Energy_log_diff + unemploy_rate_log_diff,
                clusteredSE = FALSE,
                note = "linear climate anomaly + economic indices")

#Model 6: linear climate anomaly + gdp_lag + province_fixed + year_fixed + province trend + economic indices

tmp6 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + provincecode + provincecode*Year_num + 
                  Energy_log_diff + Non.Energy_log_diff + unemploy_rate_log_diff,
                clusteredSE = TRUE,
                note = "linear climate anomaly + province trend + economic indices")

#Model 7: linear climate anomaly + gdp_lag + province_fixed + year_fixed + economic indices + dummies for structural breaks
tmp7 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + provincecode + 
                  Energy_log_diff + Non.Energy_log_diff + unemploy_rate_log_diff + any_econ, 
                clusteredSE = TRUE,
                note = "linear climate anomaly + economic indices + dummies for structural breaks")


#Model 8: linear climate anomaly + gdp_lag + province_fixed + year_random + economic indices
tmp8 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + provincecode + (1|Year_num) +
                  Energy_log_diff + Non.Energy_log_diff + unemploy_rate_log_diff,
                clusteredSE = TRUE,
                note = "linear climate anomaly + year random + economic indices")

## Warning: Some predictor variables are on very different scales: consider
## rescaling

## Warning: Some predictor variables are on very different scales: consider
## rescaling

#Model 9: linear climate anomaly + gdp_lag + province_fixed + year_random + economic indices + dummies for structural breaks
tmp9 <- climfit(dat, GDP.per.capita_DiffLog ~ 0 + MTemp_A_Spring + MTemp_A_Summer +
                  MTemp_A_Fall + MTemp_A_Winter +
                  Precip_A_Spring + Precip_A_Summer + Precip_A_Fall + Precip_A_Winter +
                  GDP.per.capita_DiffLog_lag + provincecode + (1|Year_num) +
                   Energy_log_diff + Non.Energy_log_diff + unemploy_rate_log_diff + any_econ, 
                clusteredSE = TRUE,
                note = "linear climate anomaly + year_random + economic indices + dummies for structural breaks")

## Warning: Some predictor variables are on very different scales: consider
## rescaling

## Warning: Some predictor variables are on very different scales: consider
## rescaling

b <- combine_results(tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7, tmp8, tmp9)
b <- as.data.frame(b)
t.b <- as.data.frame(t(b))
datatable(t.b, options = list(
  scrollX = TRUE,
  scrollY = 300,
  paging = FALSE
))