prior predictive checks
See Rethinking p.95
code:r
`{r}
set.seed(2971)
R <- rlkjcorr(1e4,K=2,eta=3)
`
`{r}
R <- rexp(1e4,1)
dens( R )
`
`{r}
N <- 1000 # 100 lines
`
`{r}
a <- rstudent(N, 3, 1.3, 1)
dens(a)
`
`{r}
b <- rstudent( N , 3 , 0, 1)
dens(b)
`
`{r}
plot( NULL , xlim=range(aggregated_col$time_10week) , ylim=c(-20,20) ,
xlab="time per 10 week" , ylab="MI" )
mtext( "b ~ dnorm(0,10)" )
xbar <- mean(aggregated_col$time_10week)
for ( i in 1:N ) curve( ai + bi*(x - xbar) , from=min(aggregated_col$time_10week) , to=max(aggregated_col$time_10week) , add=TRUE ,
col=col.alpha("black",0.2) )
`
`{r}
wk_prior <- c(prior(student_t(3,1.3,1), class = Intercept),
prior(student_t(3, 0, 0.5), class = b),
prior(exponential(1), class = sd),
prior(lkj(2), class = cor),
prior(exponential(1), class = sigma))
wk_prior_wo_cor <- c(prior(student_t(3, 1.3, 1), class = Intercept),
prior(student_t(3, 0, 0.5), class = b),
prior(exponential(1), class = sd),
prior(exponential(1), class = sigma))
`
`{r}
prior1 <- brm(MI ~ dep_type * time_10week + (1|anon_id) + (1|question_id),
data = aggregated_col,
prior = wk_prior_wo_cor,
chain = 4,
warmup = 700,
iter = 3000,
backend = 'cmdstanr',
seed = 1234,
sample_prior = "only"
)
`
`{r}
conditional_effects(prior1)
`
`{r}
prior2 <- brm(MI ~ dep_type * time_10week * placement_test_center + (1+ dep_type + time_10week|anon_id) + (1|question_id),
data = aggregated_col,
prior = wk_prior,
chain = 4,
warmup = 700,
iter = 3000,
backend = 'cmdstanr',
seed = 1234,
sample_prior = "only"
)
`
`{r}
conditional_effects(prior2, spaghetti = TRUE)
`