Model reproduction checklist
Breck Baldwin
March 21, 2021
knitr::opts_chunk$set(
echo = TRUE,
message = FALSE,
warning = FALSE,
results = "hide",
error = FALSE,
comment = ''
)Model overview
This document exists in the repository at https://github.com/codatmo/Liverpool as model_reproduction_checklist.Rmd and rendered in html as https://github.com/codatmo/Liverpool/model_reproduction_checklist.html.
The sister page https://codatmo.github.io/Liverpool/index.html describes the model more fully.
This model is part of the CoDatMo (Co)vid (Dat)a (Mo)deling site (https://codatmo.github.io/) which is intended to replicate and make available important COVID models written in Bayesian modeling languages like Stan or PyMC3.
Validation checklist items
- Model released on github: https://github.com/codatmo/Liverpool
- Problem description
- Model description: Research goals, references, supplementary description as necessary.
- Data
- Data generating process
- Data munging and examination
- Stan program
- Running model
- Small data set to validate model execution (not done)
- Run on availble data
- Examine model output
- Model
- Model validation
- Report posterior diagnostics
- Prior predictive check parameters/predictions
- Parameter recovery with simulated data
- Posterior predictive check
- Parameter interpretation (not done)
- Cross validation (not done)
- Simulation based calibration (SBC) (not done)
- Compute environment
Problem description
Model goal is to estimate COVID hospitalizations based on reports from 111 calls to government services and/or website reports.
Model description
See https://codatmo.github.io/Liverpool/index.html for detailed description.
Data
Actual data is available as described at https://codatmo.github.io/Data/index.html. An .RDS file exists locally at nhs_regions_111_calls.RDS.
calls_111.df <- readRDS("data/nhs_regions_111_calls.Rds")
deaths.df <- readRDS("data/nhs_regions_deaths.Rds")
calls_111 <- as.numeric(calls_111.df$London[1:71])
deaths <- as.numeric(deaths.df$London[1:14])There is a script to generate simulated data.
source("src/scripts/generateSimulatedData.R")
generatedData <- generateSimulatedData(
seed = 0,
maxTime = 100,
population = 1000000,
n_beta_pieces = 20,
n_rho_calls_111_pieces = 10,
calls_111_start = 30
)
stan_data_sim <- generatedData$stan_data
stan_data_sim$compute_likelihood <- 1 # controls whether likelihood is computed
ground_truth_sim <- generatedData$ground_truthThe script does a fair amount of data munging, see below, but the data analogous to production data are two data sources:
#deaths <- stan_data$deaths
#calls_111 <- stan_data$calls_111
cat(paste(sprintf("Weekly COVID death count=%d, Days of 111_calls (7 days per weekly count)=%d",
length(deaths),
length(calls_111))))Weekly COVID death count=14, Days of 111_calls (7 days per weekly count)=71library(ggplot2)
library(tidyr)
library(gridExtra)
day <- 1:length(calls_111)
data.df <- data.frame(calls_111,day)
data.df$deaths <- ifelse(data.df$day %% 7 == 0, deaths[data.df$day %/% 7], NA)
ggplot(data.df) +
aes(x=1:length(calls_111), y=calls_111) +
geom_line() -> p1
ggplot(data.df) +
aes(x=1:length(calls_111), y=deaths) +
geom_point() -> p2
grid.arrange(p1, p2, ncol=2)
Data munging
The data undergo considerable processing before being ready to fit with Stan. There are 22 data parameters of which two are the deaths and calls_111 shown above. The data elaboration proceeds as follows:
# deaths <- stan_data$deaths
# calls_111 <- stan_data$calls_111
# below are constants
T <- 100
maxTime <- T
initial_time <- 0
n_disease_states <- 8
n_beta_pieces <- 20
calls_111_start <- 30
n_rho_calls_111_pieces <- 10
population <- 1e+06
beta_pieceLength <- ceiling(maxTime / n_beta_pieces)
beta_left_t <- seq(0, maxTime - 1, by = beta_pieceLength)
beta_right_t <- c(beta_left_t[2:n_beta_pieces], maxTime + 1)
real_data_length = length(beta_left_t) + length(beta_right_t) + 1
real_data = c(beta_left_t, beta_right_t, population)
integer_data = c(
maxTime,
length(beta_left_t),
length(beta_right_t),
length(beta_left_t),
n_disease_states
)
integer_data_length <- length(integer_data)
rho_calls_pieceLength <-
ceiling((maxTime - calls_111_start) / n_rho_calls_111_pieces)
rho_calls_111_left_t <-
seq(calls_111_start, maxTime - 1, by = rho_calls_pieceLength)
rho_calls_111_right_t <-
c(rho_calls_111_left_t[2:n_rho_calls_111_pieces], maxTime + 1)
times <- 1:T
deaths_length <- length(deaths)
deaths_starts <- seq(from = 1, to = deaths_length * 7, by = 7)
deaths_stops <- seq(from = 7, to = deaths_length * 7, by = 7)
calls_111_length <- length(calls_111)
stan_data <- list(
initial_time = initial_time,
n_beta_pieces = n_beta_pieces,
beta_left_t = beta_left_t,
beta_right_t = beta_right_t,
T = maxTime,
times = times,
n_disease_states = n_disease_states,
population = population,
deaths_length = length(deaths),
deaths_starts = deaths_starts,
deaths_stops = deaths_stops,
deaths = deaths,
real_data_length = length(beta_left_t) + length(beta_right_t) + 1,
real_data = c(beta_left_t, beta_right_t, population),
integer_data_length = 5,
integer_data = c(
maxTime,
length(beta_left_t),
length(beta_right_t),
length(beta_left_t),
n_disease_states
),
n_rho_calls_111_pieces = n_rho_calls_111_pieces,
rho_calls_111_left_t = rho_calls_111_left_t,
rho_calls_111_right_t = rho_calls_111_right_t,
calls_111_start = calls_111_start,
calls_111 = calls_111,
calls_111_length = calls_111_length,
compute_likelihood = 1
)
# test that data has been properly munged.
# for (varname in names(stan_data)) {
# # print(paste("trying:",varname))
# if (!varname %in% names(stan_data2)) {
# print(paste("NULL for ", varname))
# }
# if (!all(stan_data[[varname]] == stan_data2[[varname]])) {
# print(paste("mismatch ", varname))
# }
# }Stan program
The Stan model is located at src/models/deaths_and_111_calls.stan:
stan_file <- "src/models/deaths_and_111_calls.stan"
lines <- readLines(stan_file)
cat(paste(lines, "\n", sep=""), sep="")functions {
real[] seeiittd(real time,
real[] state,
real[] params,
real[] real_data,
int[] integer_data) {
// Unpack integer data values
int T = integer_data[1];
int n_beta_pieces = integer_data[2];
int n_disease_states = integer_data[3];
// Unpack real data values
real beta_left_t[n_beta_pieces] = real_data[1:n_beta_pieces];
real beta_right_t[n_beta_pieces] = real_data[n_beta_pieces+1:2*n_beta_pieces];
real population = real_data[2*n_beta_pieces+1];
// Unpack parameter values
real beta_left[n_beta_pieces] = params[1:n_beta_pieces];
real grad_beta[n_beta_pieces] = params[n_beta_pieces+1:2*n_beta_pieces];
real nu = params[2*n_beta_pieces+1];
real gamma = params[2*n_beta_pieces+2];
real kappa = params[2*n_beta_pieces+3];
real omega = params[2*n_beta_pieces+4];
// Unpack state
real S = state[1];
real E1 = state[2];
real E2 = state[3];
real I1 = state[4];
real I2 = state[5];
real T1 = state[6];
real T2 = state[7];
real D = state[8];
real infection_rate;
real nuE1 = nu * E1;
real nuE2 = nu * E2;
real gammaI1 = gamma * I1;
real gammaI2 = gamma * I2;
real kappaT1 = kappa * T1;
real kappaT2 = kappa * T2;
real dS_dt;
real dE1_dt;
real dE2_dt;
real dI1_dt;
real dI2_dt;
real dT1_dt;
real dT2_dt;
real dD_dt;
for (i in 1:n_beta_pieces) {
if(time >= beta_left_t[i] && time < beta_right_t[i]) {
real beta = grad_beta[i] * (time - beta_left_t[i]) + beta_left[i];
infection_rate = beta * (I1 + I2) * S / population;
}
}
dS_dt = -infection_rate;
dE1_dt = infection_rate - nuE1;
dE2_dt = nuE1 - nuE2;
dI1_dt = nuE2 - gammaI1;
dI2_dt = gammaI1 - gammaI2;
dT1_dt = gammaI2 * omega - kappaT1;
dT2_dt = kappaT1 - kappaT2;
dD_dt = kappaT2;
return {dS_dt, dE1_dt, dE2_dt, dI1_dt, dI2_dt, dT1_dt, dT2_dt, dD_dt};
}
real[ , ] integrate_ode_explicit_trapezoidal(real[] initial_state, real initial_time, real[] times, real[] params, real[] real_data, int[] integer_data) {
real h;
vector[size(initial_state)] dstate_dt_initial_time;
vector[size(initial_state)] dstate_dt_tidx;
vector[size(initial_state)] k;
real state_estimate[size(times),size(initial_state)];
h = times[1] - initial_time;
dstate_dt_initial_time = to_vector(seeiittd(initial_time, initial_state, params, real_data, integer_data));
k = h*dstate_dt_initial_time;
state_estimate[1,] = to_array_1d(to_vector(initial_state) + h*(dstate_dt_initial_time + to_vector(seeiittd(times[1], to_array_1d(to_vector(initial_state)+k), params, real_data, integer_data)))/2);
for (tidx in 1:size(times)-1) {
h = (times[tidx+1] - times[tidx]);
dstate_dt_tidx = to_vector(seeiittd(times[tidx], state_estimate[tidx], params, real_data, integer_data));
k = h*dstate_dt_tidx;
state_estimate[tidx+1,] = to_array_1d(to_vector(state_estimate[tidx,]) + h*(dstate_dt_tidx + to_vector(seeiittd(times[tidx+1], to_array_1d(to_vector(state_estimate[tidx,])+k), params, real_data, integer_data)))/2);
}
return state_estimate;
}
}
data {
real initial_time;
int<lower=1> n_beta_pieces;
real<lower=0> beta_left_t[n_beta_pieces];
real<lower=0> beta_right_t[n_beta_pieces];
int<lower=1> n_rho_calls_111_pieces;
int<lower=0> rho_calls_111_left_t[n_rho_calls_111_pieces];
int<lower=0> rho_calls_111_right_t[n_rho_calls_111_pieces];
int<lower=1> T;
real times[T];
int<lower=1> n_disease_states;
real<lower=0> population;
int<lower=1> deaths_length;
int<lower=1> deaths_starts[deaths_length];
int<lower=1> deaths_stops[deaths_length];
int<lower=0> deaths[deaths_length];
int<lower=1> calls_111_length;
int<lower=1> calls_111_start;
int<lower=0> calls_111[calls_111_length];
int real_data_length;
real real_data[real_data_length];
int integer_data_length;
int integer_data[integer_data_length];
int<lower=0, upper=1> compute_likelihood;
}
transformed data {
real mu_dL = 4.00;
real sigma_dL = 0.20;
real mu_dI = 3.06;
real sigma_dI = 0.21;
real mu_dT = 16.00;
real sigma_dT = 0.71;
int max_lag = 13;
}
parameters {
real<lower=0, upper=1> initial_state_raw[2];
real<lower=0> beta_left[n_beta_pieces];
real<lower=0> beta_right[n_beta_pieces];
real<lower=0> dL;
real<lower=0> dI;
real<lower=0> dT;
real<lower=0, upper=1> omega;
real<lower=0> reciprocal_phi_deaths;
real<lower=0> reciprocal_phi_calls_111;
real<lower=0> rho_calls_111[n_rho_calls_111_pieces];
simplex[max_lag+1] lag_weights_calls_111;
}
transformed parameters {
real initial_state[n_disease_states];
real grad_beta[n_beta_pieces];
real nu;
real gamma;
real kappa;
real phi_deaths;
real phi_calls_111;
real state_estimate[T,n_disease_states];
vector[T+1] S;
vector[T+1] E1;
vector[T+1] E2;
vector[T+1] I1;
vector[T+1] I2;
vector[T+1] T1;
vector[T+1] T2;
vector[T+1] D;
vector[T] daily_infections;
vector[T] daily_deaths;
vector[T] calls_111_lagged_daily_infections;
vector[T] daily_calls_111;
initial_state[1] = (population-5.0)*initial_state_raw[1] + 1.0;
initial_state[2] = (population-5.0)*(1.0-initial_state_raw[1])*initial_state_raw[2]/2.0 + 1.0;
initial_state[3] = (population-5.0)*(1.0-initial_state_raw[1])*initial_state_raw[2]/2.0 + 1.0;
initial_state[4] = (population-5.0)*(1.0-initial_state_raw[1])*(1.0-initial_state_raw[2])/2.0 + 1.0;
initial_state[5] = (population-5.0)*(1.0-initial_state_raw[1])*(1.0-initial_state_raw[2])/2.0 + 1.0;
initial_state[6] = 0.0;
initial_state[7] = 0.0;
initial_state[8] = 0.0;
grad_beta = to_array_1d((to_vector(beta_right) - to_vector(beta_left))./(to_vector(beta_right_t) -
to_vector(beta_left_t)));
nu = 2.0/dL;
gamma = 2.0/dI;
kappa = 2.0/dT;
phi_deaths = 1.0 / reciprocal_phi_deaths;
phi_calls_111 = 1.0 / reciprocal_phi_calls_111;
{
real params[2*n_beta_pieces+4];
params[1:n_beta_pieces] = beta_left;
params[n_beta_pieces+1:2*n_beta_pieces] = grad_beta;
params[2*n_beta_pieces+1] = nu;
params[2*n_beta_pieces+2] = gamma;
params[2*n_beta_pieces+3] = kappa;
params[2*n_beta_pieces+4] = omega;
state_estimate = integrate_ode_explicit_trapezoidal(initial_state, initial_time, times, params, real_data, integer_data);
}
S = append_row(initial_state[1], to_vector(state_estimate[, 1]));
E1 = append_row(initial_state[2], to_vector(state_estimate[, 2]));
E2 = append_row(initial_state[3], to_vector(state_estimate[, 3]));
I1 = append_row(initial_state[4], to_vector(state_estimate[, 4]));
I2 = append_row(initial_state[5], to_vector(state_estimate[, 5]));
T1 = append_row(initial_state[6], to_vector(state_estimate[, 6]));
T2 = append_row(initial_state[7], to_vector(state_estimate[, 7]));
D = append_row(initial_state[8], to_vector(state_estimate[, 8]));
daily_infections = S[:T] - S[2:] + machine_precision();
daily_deaths = D[2:] - D[:T];
calls_111_lagged_daily_infections = lag_weights_calls_111[1]*daily_infections;
for (i in 1:max_lag) {
calls_111_lagged_daily_infections += lag_weights_calls_111[i+1]*
append_row(rep_vector(0.0, i), daily_infections[:T-i]);
}
daily_calls_111 = rep_vector(0.0, T);
for (i in 1:n_rho_calls_111_pieces) {
daily_calls_111[rho_calls_111_left_t[i]:rho_calls_111_right_t[i]-1] =
calls_111_lagged_daily_infections[rho_calls_111_left_t[i]:rho_calls_111_right_t[i]-1] *
rho_calls_111[i];
}
}
model {
initial_state_raw[1] ~ beta(5.0, 0.5);
initial_state_raw[2] ~ beta(1.1, 1.1);
beta_left ~ normal(0, 0.5);
beta_right ~ normal(0, 0.5);
dL ~ normal(mu_dL, sigma_dL);
dI ~ normal(mu_dI, sigma_dI);
dT ~ normal(mu_dT, sigma_dT);
omega ~ beta(100, 9803);
reciprocal_phi_deaths ~ exponential(5);
reciprocal_phi_calls_111 ~ exponential(5);
rho_calls_111 ~ normal(0, 0.5);
lag_weights_calls_111 ~ dirichlet(rep_vector(0.1, max_lag+1));
if (compute_likelihood == 1) {
for (i in 1:deaths_length) {
target += neg_binomial_2_lpmf(deaths[i] |
sum(daily_deaths[deaths_starts[i]:deaths_stops[i]]), phi_deaths);
}
target += neg_binomial_2_lpmf(calls_111 |
daily_calls_111[calls_111_start:(calls_111_start-1)+calls_111_length], phi_calls_111);
}
}
generated quantities {
vector[T-1] growth_rate = (log(daily_infections[2:]) - log(daily_infections[:T-1]))*100;
int pred_deaths[deaths_length];
int pred_calls_111[calls_111_length];
for (i in 1:deaths_length) {
pred_deaths[i] = neg_binomial_2_rng(sum(daily_deaths[deaths_starts[i]:deaths_stops[i]]),
phi_deaths);
}
for (i in 1:calls_111_length) {
pred_calls_111[i] = neg_binomial_2_rng(daily_calls_111[calls_111_start - 1 + i],
phi_calls_111);
}
}Running model
First some inits are setup:
set.seed(0)
initf <- function() {
n_beta_pieces <- stan_data$n_beta_pieces
initial_state_raw = c(runif(1, min = 0.99999, max = 1.0), runif(1, min = 0.0, max = 1.0))
beta_left = exp(runif(n_beta_pieces, min = -2, max = 0.5))
beta_right = exp(runif(n_beta_pieces, min = -2, max = 0.5))
dL = runif(1, min = 3.5, max = 4.0)
dI = runif(1, min = 2.2, max = 2.6)
dT = runif(1, min = 11.0, max = 13.0)
omega = 1/(1 + exp(-runif(1, min = -5, max = -3)))
list(
initial_state_raw = initial_state_raw,
beta_left = beta_left,
beta_right = beta_right,
dL = dL,
dI = dI,
dT = dT,
omega = omega
)
}The code is run as:
library(cmdstanr)
library(rstan)
library(xfun)
model <- cmdstan_model(file.path("src","models","deaths_and_111_calls.stan"))
num_warmup <- 1000
num_samples <- 1000
num_chains <- 4
fit_1 <-
xfun::cache_rds({model$sample(data=stan_data,
seed=999,
chains=num_chains,
init = initf,
num_warmup=num_warmup,
num_samples=num_samples,
output_dir="output")},
rerun=FALSE, file="fit1", dir="_cache/")
rs_fit_1 <- rstan::read_stan_csv(fit_1$output_files())
print(fit_1$time())$total
[1] 2287.373
$chains
chain_id warmup sampling total
1 1 256.025 185.119 441.144
2 2 275.761 173.661 449.422
3 3 472.499 221.388 693.887
4 4 442.008 260.252 702.260#running some warmup and just 100 draws with two chains.
Viewing a text based summary of the fit, this will be very long, it is commented out and not run:
# fit_1$cmdstan_summary()Model validation
Below are the diagnostics used to help validate the model.
Run posterior diagnostics
There are standard diagnostics that look for errors in the posterior.
fit_1$cmdstan_diagnose()- Treedepth warnings
- Divergence check passed
- E-BFMI satisfactory
- R-hat values satisfactory
Prior predictive check
The prior predictive check estimates the model parameters without the likelihood being used. The resulting draws are then used to predict new data via predictive application of the likelihood given the draws. Note that compute_likelihood = 0 prevents the likelihood being computed in the model fitting step.
stan_data_no_likelihood <- stan_data
stan_data_no_likelihood$compute_likelihood <- 0
fit_2 <- xfun::cache_rds({model$sample(data=stan_data_no_likelihood, seed=999,
chains=num_chains, init = initf,
num_warmup=num_warmup,
num_samples=num_samples,
output_dir="output")},
rerun=TRUE, file="fit2", dir="_cache/")
rs_fit_2 <- rstan::read_stan_csv(fit_2$output_files())
rs_ex_fit_2 <- rstan::extract(rs_fit_2)
random_draws <- sample(1:nrow(rs_ex_fit_2$pred_deaths), 10, replace=FALSE)
deaths_draws <- t(rs_ex_fit_2$pred_deaths[random_draws,])
deaths_draws.df <- data.frame(deaths_draws)
names(deaths_draws.df) <- random_draws
draw_names <- colnames(deaths_draws.df)
weeks <- 1:nrow(deaths_draws.df)
p_data2.df <- cbind(deaths,weeks,deaths_draws.df)
p_long_data <- gather(p_data2.df, draw, deaths_sim, draw_names)
p3 <- ggplot(data=p_long_data, aes(x=weeks)) +
geom_line(aes(y=deaths_sim, group=draw, color=draw), size=.5) +
geom_line(aes(y=deaths), color="black", size=.5)
p3
Above we see the posterior of the priors without seeing any data. Generally some justification for prior distributions is expected. The above are weakly informative in that they cover a broad range of plausible values.
call_draws <- t(rs_ex_fit_2$pred_calls_111[random_draws,]) #keep same draws
call_draws.df <- data.frame(call_draws)
names(call_draws.df) <- random_draws
days <- 1:nrow(call_draws.df)
p_data3.df <- cbind(calls_111,days,call_draws.df)
p_long_data3 <- gather(p_data3.df, draw, calls_sim, draw_names)
p4 <- ggplot(data=p_long_data3, aes(x=days)) +
geom_point(aes(y=calls_sim, group=draw, color=draw), size=.5) +
geom_line(aes(y=calls_111), color="black", size=.5)
p4
Parameter recovery with simulated data
Parameter recovery establishes that for some small set of values the model reasons properly. We pick a draw from the above distributions, simulate data with it and then attempt to recover the parameters that we simulated with. We can look at the above graph and pick expected outliers or close to actual data. For this we pick 42 as a middle-of-the-road example.
library(bayesplot)
#Pick one arbitrary draw from the prior distribution
draw <- 42
stan_data_3 <- stan_data
stan_data_3$calls_111 <- rs_ex_fit_2$pred_calls_111[42,]
stan_data_3$deaths <- rs_ex_fit_2$pred_deaths[42,]
stan_data_3$compute_likelihood <- 1
fit_3 <- xfun::cache_rds({model$sample(data=stan_data_3, seed=999,
chains=num_chains, init = initf,
num_warmup=num_warmup,
num_samples=num_samples,
output_dir="output")},
rerun=TRUE, file="fit3", dir="_cache/")
rs_fit_3 <- rstan::read_stan_csv(fit_3$output_files())
rs_fit_3_draws <- rstan::extract(rs_fit_3)
mcmc_intervals(fit_3$draws(variables=c("dL", "dI", "dT","omega")))
A subset of the estimated paramters shown, mean, 50% and 90% intervals indicated. Actual values below:
report <- paste(sprintf("Draw number %d", draw),
sprintf("actual dL=%.2f", rs_ex_fit_2$dL[draw]),
sprintf("actual dI=%.2f", rs_ex_fit_2$dI[draw]),
sprintf("actual dT=%.2f", rs_ex_fit_2$dT[draw]),
sprintf("actual omega=%.2f", rs_ex_fit_2$omega[draw]),
sep="\n")
cat(report)Draw number 42
actual dL=3.81
actual dI=3.36
actual dT=16.43
actual omega=0.01Posterior predictive check
Compare predicted data from original run that used actual data.
Predict deaths
rs_fit_1_draws <- rstan::extract(rs_fit_1)
summary.df <- as.data.frame(summary(rs_fit_1,pars=c("pred_deaths"), probs=c(.05,.5,.95))$summary)
weeks <- 1:length(deaths)
all.df <- cbind(summary.df, deaths, weeks)
colnames(all.df) <- make.names(colnames(all.df)) # to remove % in the col names
p4 <- ggplot(all.df, mapping = aes(x = weeks)) +
geom_ribbon(aes(ymin = X5., ymax = X95.), fill = "orange", alpha = 0.6) +
geom_line(mapping = aes(y = X50.)) +
geom_point(mapping = aes(y = deaths)) +
labs(x = "Week", y = "Deaths")
p4
Actual deaths are black dots, mean predicted deaths black line, orange ribbon spans central 95% of draws.
Predict 111 calls
rs_fit_1_draws <- rstan::extract(rs_fit_1)
summary_2.df <- as.data.frame(summary(rs_fit_1,pars=c("pred_calls_111"), probs=c(.05,.5,.95))$summary)
days <- 1:length(calls_111)
all_2.df <- cbind(summary_2.df, calls_111, days)
colnames(all_2.df) <- make.names(colnames(all_2.df)) # to remove % in the col names
p5 <- ggplot(all_2.df, mapping = aes(x = days)) +
geom_ribbon(aes(ymin = X5., ymax = X95.), fill = "green", alpha = 0.6) +
geom_line(mapping = aes(y = X50.)) +
geom_point(mapping = aes(y = calls_111)) +
labs(x = "Day", y = "Calls 111")
p5
Actual calls are black dots, mean predicted calls black line, green ribbon spans central 90% of draws.
Parameter interpretation
Impact of param changes on predictions.
Cross validation
SBC validation
Computing environment
devtools::session_info()Error in get(genname, envir = envir) : object 'testthat_print' not found─ Session info ───────────────────────────────────────────────────────────────
setting value
version R version 3.6.2 (2019-12-12)
os macOS Catalina 10.15.7
system x86_64, darwin15.6.0
ui X11
language (EN)
collate en_US.UTF-8
ctype en_US.UTF-8
tz America/New_York
date 2021-03-26
─ Packages ───────────────────────────────────────────────────────────────────
package * version date lib source
abind 1.4-5 2016-07-21 [1] CRAN (R 3.6.0)
assertthat 0.2.1 2019-03-21 [1] CRAN (R 3.6.0)
backports 1.1.5 2019-10-02 [1] CRAN (R 3.6.0)
bayesplot * 1.7.1 2019-12-01 [1] CRAN (R 3.6.0)
callr 3.4.3 2020-03-28 [1] CRAN (R 3.6.2)
checkmate 1.9.4 2019-07-04 [1] CRAN (R 3.6.0)
cli 2.0.2 2020-02-28 [1] CRAN (R 3.6.0)
cmdstanr * 0.3.0 2021-02-24 [1] local
colorspace 1.4-1 2019-03-18 [1] CRAN (R 3.6.0)
crayon 1.3.4 2017-09-16 [1] CRAN (R 3.6.0)
data.table 1.12.8 2019-12-09 [1] CRAN (R 3.6.0)
desc 1.2.0 2018-05-01 [1] CRAN (R 3.6.0)
devtools 2.3.1 2020-07-21 [1] CRAN (R 3.6.2)
digest 0.6.23 2019-11-23 [1] CRAN (R 3.6.0)
DirichletReg * 0.7-0 2020-05-29 [1] CRAN (R 3.6.2)
dplyr 0.8.3 2019-07-04 [1] CRAN (R 3.6.0)
ellipsis 0.3.0 2019-09-20 [1] CRAN (R 3.6.0)
evaluate 0.14 2019-05-28 [1] CRAN (R 3.6.0)
fansi 0.4.0 2018-10-05 [1] CRAN (R 3.6.0)
farver 2.0.1 2019-11-13 [1] CRAN (R 3.6.0)
Formula * 1.2-3 2018-05-03 [1] CRAN (R 3.6.0)
fs 1.3.1 2019-05-06 [1] CRAN (R 3.6.0)
ggplot2 * 3.2.1 2019-08-10 [1] CRAN (R 3.6.0)
ggridges 0.5.2 2020-01-12 [1] CRAN (R 3.6.0)
glue 1.3.1 2019-03-12 [1] CRAN (R 3.6.0)
gridExtra * 2.3 2017-09-09 [1] CRAN (R 3.6.0)
gtable 0.3.0 2019-03-25 [1] CRAN (R 3.6.0)
htmltools 0.4.0 2019-10-04 [1] CRAN (R 3.6.0)
inline 0.3.15 2018-05-18 [1] CRAN (R 3.6.0)
jsonlite 1.7.0 2020-06-25 [1] CRAN (R 3.6.2)
knitr 1.28 2020-02-06 [1] CRAN (R 3.6.0)
labeling 0.3 2014-08-23 [1] CRAN (R 3.6.0)
lattice 0.20-38 2018-11-04 [1] CRAN (R 3.6.2)
lazyeval 0.2.2 2019-03-15 [1] CRAN (R 3.6.0)
lifecycle 0.2.0 2020-03-06 [1] CRAN (R 3.6.0)
loo 2.2.0 2019-12-19 [1] CRAN (R 3.6.0)
magrittr 1.5 2014-11-22 [1] CRAN (R 3.6.0)
matrixStats 0.55.0 2019-09-07 [1] CRAN (R 3.6.0)
maxLik 1.4-6 2020-11-24 [1] CRAN (R 3.6.2)
memoise 1.1.0 2017-04-21 [1] CRAN (R 3.6.0)
miscTools 0.6-26 2019-12-08 [1] CRAN (R 3.6.0)
munsell 0.5.0 2018-06-12 [1] CRAN (R 3.6.0)
pillar 1.4.3 2019-12-20 [1] CRAN (R 3.6.0)
pkgbuild 1.0.6 2019-10-09 [1] CRAN (R 3.6.0)
pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 3.6.0)
pkgload 1.0.2 2018-10-29 [1] CRAN (R 3.6.0)
plyr 1.8.5 2019-12-10 [1] CRAN (R 3.6.0)
posterior 0.1.3 2021-02-24 [1] local
prettyunits 1.0.2 2015-07-13 [1] CRAN (R 3.6.0)
processx 3.4.5 2020-11-30 [1] CRAN (R 3.6.2)
ps 1.3.0 2018-12-21 [1] CRAN (R 3.6.0)
purrr 0.3.3 2019-10-18 [1] CRAN (R 3.6.0)
R6 2.4.1 2019-11-12 [1] CRAN (R 3.6.0)
Rcpp 1.0.3 2019-11-08 [1] CRAN (R 3.6.0)
remotes 2.2.0 2020-07-21 [1] CRAN (R 3.6.2)
reshape2 1.4.3 2017-12-11 [1] CRAN (R 3.6.0)
rlang 0.4.10 2020-12-30 [1] CRAN (R 3.6.2)
rmarkdown 2.0 2019-12-12 [1] CRAN (R 3.6.0)
rprojroot 1.3-2 2018-01-03 [1] CRAN (R 3.6.0)
rstan * 2.19.3 2020-02-11 [1] CRAN (R 3.6.0)
sandwich 3.0-0 2020-10-02 [1] CRAN (R 3.6.2)
scales 1.1.0 2019-11-18 [1] CRAN (R 3.6.0)
sessioninfo 1.1.1 2018-11-05 [1] CRAN (R 3.6.0)
StanHeaders * 2.21.0-1 2020-01-19 [1] CRAN (R 3.6.0)
stringi 1.4.3 2019-03-12 [1] CRAN (R 3.6.0)
stringr 1.4.0 2019-02-10 [1] CRAN (R 3.6.0)
testthat 2.3.2 2020-03-02 [1] CRAN (R 3.6.0)
tibble 3.0.3 2020-07-10 [1] CRAN (R 3.6.2)
tidyr * 1.0.0 2019-09-11 [1] CRAN (R 3.6.0)
tidyselect 1.1.0 2020-05-11 [1] CRAN (R 3.6.2)
usethis 1.6.1 2020-04-29 [1] CRAN (R 3.6.2)
vctrs 0.3.1 2020-06-05 [1] CRAN (R 3.6.2)
withr 2.1.2 2018-03-15 [1] CRAN (R 3.6.0)
xfun * 0.22 2021-03-11 [1] CRAN (R 3.6.2)
yaml 2.2.0 2018-07-25 [1] CRAN (R 3.6.0)
zoo 1.8-7 2020-01-10 [1] CRAN (R 3.6.0)
[1] /Library/Frameworks/R.framework/Versions/3.6/Resources/library