rm(list = ls())
suppressWarnings(suppressMessages( {
library(brms)
library(beeswarm)
library(vioplot)
library(glmmTMB)
library(car)
library(cowplot)
library(ggplot2)
library(tidyverse)
library(bayesplot)
library(gtsummary)
library(vegan)
library(knitr)
library(ggplot2)
} ) )Experimental treatment with edaravone in a mouse model of spinocerebellar ataxia 1
Citrate synthase activity
This page shows R code for the study Sucha et al. (2023, International Journal of Molecular Science).
Citation:
Sucha M, Benediktova S, Tichanek F, Jedlicka J, Kapl S, Jelinkova D, Purkartova Z, Tuma J, Kuncova J, Cendelin J. Experimental Treatment with Edaravone in a Mouse Model of Spinocerebellar Ataxia 1. International Journal of Molecular Sciences. 2023; 24(13):10689. https://doi.org/10.3390/ijms241310689
GitHub page: https://github.com/filip-tichanek/edaravonSCA1
Citrate synthase activity is related to mitochondria abundance and was assessed in cerebellar (cb) and hippocampal (hp) tissue. There are 2 technical replicates per animal. The factor of animal (id) will thus represent random-intercept.
1 Upload of packages
2 Data upload and wrangling
## Data upload
load("source_data/myEnvironment.RData")
## defining new variable
csyn$sex_num <- ifelse(csyn$sex=="F",-0.5,0.5)
csyn$group_sex <- interaction(csyn$group, csyn$sex)
csyn$group_sex <- factor(csyn$group_sex,
levels=c('ctrl.wt.F','ctrl.wt.M',
'eda.wt.F','eda.wt.M',
'ctrl.sca.F','ctrl.sca.M',
'eda.sca.F','eda.sca.M'))
## subsetting cerebellar and hippocampal datasets
csyn_hp <- subset(csyn, csyn$tissue == 'hippocampus')
csyn_cb <- subset(csyn, csyn$tissue == 'cerebellum')
## data summary
summary(csyn) id sex genotype treatment tissue cs
1 : 2 F:71 sca:79 ctrl:80 cerebellum :77 Min. :0.006127
2 : 2 M:83 wt :75 eda :74 hippocampus:77 1st Qu.:0.009835
3 : 2 Median :0.011520
4 : 2 Mean :0.012113
5 : 2 3rd Qu.:0.014216
6 : 2 Max. :0.023529
(Other):142
mereni group sex_num group_sex
a:77 ctrl.wt :40 Min. :-0.50000 eda.wt.M :24
b:77 eda.wt :35 1st Qu.:-0.50000 ctrl.wt.F :20
ctrl.sca:40 Median : 0.50000 ctrl.wt.M :20
eda.sca :39 Mean : 0.03896 ctrl.sca.F:20
3rd Qu.: 0.50000 ctrl.sca.M:20
Max. : 0.50000 eda.sca.F :20
(Other) :30 3 Data exploration
## setting spacing
par(mfrow=c(3,2))
par(mar=c(4,3.5,1,1))
par(mgp=c(4,0.6,0))
## both sex together
plot(csyn_hp$cs~csyn_hp$group,
col=cola, cex.axis=0.8)
plot(csyn_cb$cs~csyn_cb$group,
col=cola, cex.axis=0.8)
## by both group and sex
plot(csyn_hp$cs~csyn_hp$group_sex,
col=cola[c(1,1,2,2,3,3,4,4)],cex.axis=0.8, las=2)
plot(csyn_cb$cs~csyn_cb$group_sex,
col=cola[c(1,1,2,2,3,3,4,4)],cex.axis=0.8, las=2)
## data distribution
hist(csyn_hp$cs, 12)
hist(csyn_cb$cs, 12)
There is no apparent sex effect. As the distribution seems slightly asymmetric, we will use gamma likelihood.
4 Hippocampus
4.1 Modelling
4.1.1 Prior specification
Data summary
## mean
mean(log(csyn_hp$cs))[1] -4.342288## SDs
sd(log(csyn_hp$cs))*5[1] 1.249979sd(log(csyn_hp[csyn_hp$genotype == 'sca',]$cs))*5[1] 1.054387Setting priors
prior01 <- c(set_prior("student_t(3, -4.34, 1.24)", class = "b", coef = "Intercept"),
set_prior("normal(0, 2)", class = "b", coef = "groupctrl.wt"),
set_prior("normal(0, 2)", class = "b", coef = "groupeda.wt"),
set_prior("normal(0, 1.2)", class = "b", coef = "groupeda.sca"))4.1.2 Model fit and diagnostics
## model fit
m_csyn_hp <- run_model(
brm(cs~0+Intercept+group+ (1|id),
data=csyn_hp, prior = prior01,
family=Gamma(link='log'),
save_pars = save_pars(all = TRUE),
iter=8000, warmup=2000,chains=2,cores=2,seed=17,
control = list(adapt_delta = 0.99)),
'output/brm/m_csyn_hp', reuse = TRUE)
## chains mixing and convergence
mcmc_trace(m_csyn_hp, pars = c('b_Intercept',
'b_groupctrl.wt',
'b_groupeda.wt',
'b_groupeda.sca'))
## model summary
summary(m_csyn_hp) Family: gamma
Links: mu = log; shape = identity
Formula: cs ~ 0 + Intercept + group + (1 | id)
Data: csyn_hp (Number of observations: 77)
Draws: 2 chains, each with iter = 8000; warmup = 2000; thin = 1;
total post-warmup draws = 12000
Group-Level Effects:
~id (Number of levels: 39)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.23 0.04 0.17 0.31 1.00 2483 4752
Population-Level Effects:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept -4.38 0.08 -4.54 -4.22 1.00 2198 3492
groupctrl.wt 0.06 0.11 -0.16 0.29 1.00 2358 3632
groupeda.wt 0.04 0.12 -0.20 0.26 1.00 2633 4567
groupeda.sca 0.07 0.11 -0.15 0.30 1.00 2467 4145
Family Specific Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
shape 58.44 13.67 34.59 88.55 1.00 4202 6452
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).Convergence is fine, chains are well mixed, Rhat values are 1, ESS sufficient
4.1.3 PPC
pp_check(m_csyn_hp,type='dens_overlay',ndraws = 50)
pp_check(m_csyn_hp,type='dens_overlay_grouped',ndraws = 50,group='group')
pp_check(m_csyn_hp,type='scatter_avg', ndraws = 50)
pp_check(m_csyn_hp,type="stat_2d", stat = c("max", "min"),ndraws=100)
pp_check(m_csyn_hp,type="stat_2d", stat = c("mean", "sd"),ndraws=100)
PPC seems fine
4.1.4 Posterior draws extraction
post_fix_hp<-as.data.frame(m_csyn_hp, variable = c("b_Intercept","b_groupctrl.wt",
"b_groupeda.wt","b_groupeda.sca"))
names(post_fix_hp)[1]<-"sca_ctrl"
post_fix_hp$wt_ctrl<-post_fix_hp$sca_ctrl+post_fix_hp$b_groupctrl.wt
post_fix_hp$wt_eda<-post_fix_hp$sca_ctrl+post_fix_hp$b_groupeda.wt
post_fix_hp$sca_eda<-post_fix_hp$sca_ctrl+post_fix_hp$b_groupeda.sca
post_fix_hp$wt_contrast<-post_fix_hp$wt_eda-post_fix_hp$wt_ctrl
summary(post_fix_hp) sca_ctrl b_groupctrl.wt b_groupeda.wt b_groupeda.sca
Min. :-4.733 Min. :-0.46234 Min. :-0.38644 Min. :-0.368559
1st Qu.:-4.435 1st Qu.:-0.01080 1st Qu.:-0.04078 1st Qu.:-0.001021
Median :-4.383 Median : 0.06427 Median : 0.03715 Median : 0.073324
Mean :-4.382 Mean : 0.06402 Mean : 0.03689 Mean : 0.073787
3rd Qu.:-4.329 3rd Qu.: 0.13976 3rd Qu.: 0.11463 3rd Qu.: 0.147812
Max. :-4.056 Max. : 0.53849 Max. : 0.48139 Max. : 0.546332
wt_ctrl wt_eda sca_eda wt_contrast
Min. :-4.673 Min. :-4.712 Min. :-4.616 Min. :-0.54445
1st Qu.:-4.372 1st Qu.:-4.401 1st Qu.:-4.361 1st Qu.:-0.10599
Median :-4.318 Median :-4.345 Median :-4.310 Median :-0.02743
Mean :-4.318 Mean :-4.345 Mean :-4.308 Mean :-0.02713
3rd Qu.:-4.264 3rd Qu.:-4.289 3rd Qu.:-4.255 3rd Qu.: 0.05123
Max. :-3.973 Max. :-4.002 Max. :-3.978 Max. : 0.48769 4.2 Data and model visualization
In data visualization, individual points represent average per subject. Violin plot (area indicating data distribution) was, in contrast, constructed on the basis of individual data points (non-averaged, 1 animals thus generated two data-points)
suppressWarnings(suppressMessages( {
range<-c(0,0.04);scal<-range[2]-range[1]
groupquant<-sapply(post_fix_hp[,c(5,6,1,7)],function(p) quantile(p, probs = c(0.025,0.975,0.5)))
groupquant<-exp(groupquant)
par(mfrow=c(1,2))
par(mar=c(2,3,0.2,0.2))
par(mgp=c(2,0.6,0))
csyn_hp$group<-factor(csyn_hp$group,levels=c('ctrl.wt', 'eda.wt','ctrl.sca','eda.sca'))
xrange<-c(0.5,4.5);xscal=xrange[2]-xrange[1]
plot(NULL,xlim=c(xrange[1],xrange[2]),ylim=c(range[1],range[2]),xlab="",
ylab="Citrate synthase in hippocampus",las=1, axes=FALSE,cex.lab=1.1)
rect(xrange[1],range[2],xrange[2],range[1],col="grey92",border=NA)
x<-range[1]
repeat{
lines(c(xrange[1],xrange[2]),c(x,x),col="white",lwd=0.7)
x=x+0.01;if(x>range[2]){break}}
x<-0
repeat{
lines(c(x,x),c(range[1],range[2]),col="white",lwd=0.7)
x=x+1;if(x>3.5){break}}
vioplot(csyn_hp$cs~csyn_hp$group,col=colb,at=c(1:4),add=T,
border=cola,axes=FALSE,drawRect=F,lwd=1,wex=0.8)
csyn_hp_sum <- csyn_hp %>% group_by(id, group, sex_num) %>% summarize(cbcsl = mean(cs))
beeswarm(csyn_hp_sum$cbcsl~csyn_hp_sum$group,col=cola, at=c(1:4),
add=T,pwpch=(16.5+csyn_hp_sum$sex_num), cex=1.2 )
tp<-(groupquant[3,])
xx<-1;wid=0.25
repeat{
lines(c(xx-wid,xx+wid),c(tp[xx],tp[xx]),lwd=2,col="black");
lines(c(xx,xx),c(groupquant[1,xx],groupquant[2,xx]),lwd=1.7,col="black")
xx<-xx+1
if(xx>4){break}}
tckk=-0.02
axis(2,las=2,cex.axis=1.1,at = seq(range[1],range[2],by=0.01),
pos=xrange[1], tck=tckk)
axis(side=1,las=1,cex.axis=1.1,at=c(seq(1,4,by=1)),
labels=c(rep("",length(seq(1,4,by=1)))),pos=range[1],tck=tckk)
lines(c(xrange[1],xrange[2]),c(range[1],range[1]))
text(c(1:4),c(rep(range[1]-0.046*scal,3)),xpd=T,cex=1.1,col=cola,
c("Wt","Wt", "SCA1", "SCA1"))
ypo<-c(range[1]-0.1*scal,range[1]-0.1*scal)
text(c(1:4),c(ypo[1],ypo[2],ypo[1],ypo[2]),xpd=T,cex=1.1,col=cola,
c("0", "0", "edv"))
ypos<-0.038
xpos=1
points(xpos,ypos,pch=16,cex=1.4,col=rgb(0.4,0.4,0.4,alpha=0.6))
text(xpos+0.24*xscal,ypos,"Females",col="grey40",cex=1.2)
points(xpos,ypos-0.08*scal,pch=17,cex=1.4,col=rgb(0.4,0.4,0.4,alpha=0.6))
text(xpos+0.24*xscal,ypos-0.08*scal,"Males",col="grey40",cex=1.2)
### probability distributions --------
ypos<-0.6;xpos=1
range=c(0, 24)
scal = 12
par(mar=c(2.2,1,0,0))
dif<-data.frame(post_fix_hp$wt_contrast,post_fix_hp$b_groupeda.sca,
post_fix_hp$b_groupctrl.wt)
dif<-exp(dif)
yla<- "Fold-difference in CS in HP"
tckk=-0.018
ste<-seq(0.4,1.7,by=0.2)
mons_poste(dif,1)
zpos=seq(0,1,1/3)
xx=1;ind=0.4;xpol<-0.65
text(xpol,zpos[xx]+zpos[2]*ind,
"Wt: ed x 0 ",cex=1)
xx=xx+1
text(xpol,zpos[xx]+zpos[2]*ind,
"SCA1: ed x 0 ",cex=1)
xx=xx+1
text(xpol,zpos[xx]+zpos[2]*ind,
"0: Wt x SCA1",cex=1)
} ) )
5 Cerebellum
5.1 Modelling
5.1.1 Prior specification
Data summary
## mean
mean(log(csyn_cb$cs))[1] -4.555682## SDs
sd(log(csyn_cb$cs))*5[1] 1.206754sd(log(csyn_cb[csyn_cb$genotype == 'sca',]$cs))*5[1] 1.308766Setting priors
prior01 <- c(set_prior("student_t(3, -4.56, 1.2)", class = "b", coef = "Intercept"),
set_prior("normal(0, 2)", class = "b", coef = "groupctrl.wt"),
set_prior("normal(0, 2)", class = "b", coef = "groupeda.wt"),
set_prior("normal(0, 1.2)" , class = "b", coef = "groupeda.sca"))5.1.2 Model fit and diagnostics
## model fit
m_csyn_cb <- run_model(
brm(cs~0+Intercept+group+ (1|id),
data=csyn_cb, prior = prior01,
family=Gamma(link='log'),
save_pars = save_pars(all = TRUE),
iter=8000, warmup=2000,chains=2,cores=2,seed=17,
control = list(adapt_delta = 0.99)),
'output/brm/m_csyn_cb', reuse = TRUE)
## chains mixing and convergence
mcmc_trace(m_csyn_cb, pars = c('b_Intercept',
'b_groupctrl.wt',
'b_groupeda.wt',
'b_groupeda.sca'))
## model summary
summary(m_csyn_cb) Family: gamma
Links: mu = log; shape = identity
Formula: cs ~ 0 + Intercept + group + (1 | id)
Data: csyn_cb (Number of observations: 77)
Draws: 2 chains, each with iter = 8000; warmup = 2000; thin = 1;
total post-warmup draws = 12000
Group-Level Effects:
~id (Number of levels: 39)
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept) 0.19 0.04 0.12 0.26 1.00 3021 3718
Population-Level Effects:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
Intercept -4.51 0.07 -4.65 -4.38 1.00 2995 4769
groupctrl.wt 0.02 0.10 -0.17 0.22 1.00 3003 5358
groupeda.wt -0.08 0.10 -0.28 0.11 1.00 3483 5663
groupeda.sca -0.06 0.10 -0.26 0.14 1.00 3505 5502
Family Specific Parameters:
Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
shape 36.62 8.41 22.07 55.16 1.00 4232 5650
Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).Convergence is fine, chains are well mixed, Rhat values are 1, ESS sufficient
5.1.3 PPC
pp_check(m_csyn_cb,type='dens_overlay',ndraws = 50)
pp_check(m_csyn_cb,type='dens_overlay_grouped',ndraws = 50,group='group')
pp_check(m_csyn_cb,type='scatter_avg', ndraws = 50)
pp_check(m_csyn_cb,type="stat_2d", stat = c("max", "min"),ndraws=100)
pp_check(m_csyn_cb,type="stat_2d", stat = c("mean", "sd"),ndraws=100)
PPC seems fine
5.1.4 Posterior draws extraction
post_fix_cb<-as.data.frame(m_csyn_cb, variable = c("b_Intercept","b_groupctrl.wt",
"b_groupeda.wt","b_groupeda.sca"))
names(post_fix_cb)[1]<-"sca_ctrl"
post_fix_cb$wt_ctrl<-post_fix_cb$sca_ctrl+post_fix_cb$b_groupctrl.wt
post_fix_cb$wt_eda<-post_fix_cb$sca_ctrl+post_fix_cb$b_groupeda.wt
post_fix_cb$sca_eda<-post_fix_cb$sca_ctrl+post_fix_cb$b_groupeda.sca
post_fix_cb$wt_contrast<-post_fix_cb$wt_eda-post_fix_cb$wt_ctrl
summary(post_fix_cb) sca_ctrl b_groupctrl.wt b_groupeda.wt b_groupeda.sca
Min. :-4.824 Min. :-0.36183 Min. :-0.53782 Min. :-0.46811
1st Qu.:-4.559 1st Qu.:-0.04349 1st Qu.:-0.15215 1st Qu.:-0.12820
Median :-4.514 Median : 0.02231 Median :-0.08505 Median :-0.06249
Mean :-4.514 Mean : 0.02327 Mean :-0.08461 Mean :-0.06160
3rd Qu.:-4.467 3rd Qu.: 0.08963 3rd Qu.:-0.01867 3rd Qu.: 0.00511
Max. :-4.198 Max. : 0.44294 Max. : 0.30466 Max. : 0.40804
wt_ctrl wt_eda sca_eda wt_contrast
Min. :-4.780 Min. :-4.931 Min. :-4.938 Min. :-0.56282
1st Qu.:-4.539 1st Qu.:-4.648 1st Qu.:-4.623 1st Qu.:-0.17759
Median :-4.490 Median :-4.599 Median :-4.575 Median :-0.10795
Mean :-4.490 Mean :-4.598 Mean :-4.575 Mean :-0.10788
3rd Qu.:-4.442 3rd Qu.:-4.549 3rd Qu.:-4.527 3rd Qu.:-0.04017
Max. :-4.205 Max. :-4.286 Max. :-4.319 Max. : 0.37646 5.2 Data and model visualization
In data visualization, individual points represent average per subject. Violin plot (area indicating data distribution) was, in contrast, constructed on the basis of individual data points (non-averaged, 1 animals thus generated two data-points)
suppressWarnings(suppressMessages( {
range<-c(0,0.04);scal<-range[2]-range[1]
groupquant<-sapply(post_fix_cb[,c(5,6,1,7)],function(p) quantile(p, probs = c(0.025,0.975,0.5)))
groupquant<-exp(groupquant)
par(mfrow=c(1,2))
par(mar=c(2,3,0.2,0.2))
par(mgp=c(2,0.6,0))
csyn_cb$group<-factor(csyn_cb$group,levels=c('ctrl.wt', 'eda.wt','ctrl.sca','eda.sca'))
xrange<-c(0.5,4.5);xscal=xrange[2]-xrange[1]
plot(NULL,xlim=c(xrange[1],xrange[2]),ylim=c(range[1],range[2]),xlab="",
ylab="Citrate synthase in cerebellum",las=1, axes=FALSE,cex.lab=1.1)
rect(xrange[1],range[2],xrange[2],range[1],col="grey92",border=NA)
x<-range[1]
repeat{
lines(c(xrange[1],xrange[2]),c(x,x),col="white",lwd=0.7)
x=x+0.01;if(x>range[2]){break}}
x<-0
repeat{
lines(c(x,x),c(range[1],range[2]),col="white",lwd=0.7)
x=x+1;if(x>3.5){break}}
vioplot(csyn_cb$cs~csyn_cb$group,col=colb,at=c(1:4),add=T,
border=cola,axes=FALSE,drawRect=F,lwd=1,wex=0.8)
csyn_cb_sum <- csyn_cb %>% group_by(id, group, sex_num) %>% summarize(cbcsl = mean(cs))
beeswarm(csyn_cb_sum$cbcsl~csyn_cb_sum$group,col=cola, at=c(1:4),
add=T,pwpch=(16.5+csyn_cb_sum$sex_num), cex=1.2 )
tp<-(groupquant[3,])
xx<-1;wid=0.25
repeat{
lines(c(xx-wid,xx+wid),c(tp[xx],tp[xx]),lwd=2,col="black");
lines(c(xx,xx),c(groupquant[1,xx],groupquant[2,xx]),lwd=1.7,col="black")
xx<-xx+1
if(xx>4){break}}
tckk=-0.02
axis(2,las=2,cex.axis=1.1,at = seq(range[1],range[2],by=0.01),
pos=xrange[1], tck=tckk)
axis(side=1,las=1,cex.axis=1.1,at=c(seq(1,4,by=1)),
labels=c(rep("",length(seq(1,4,by=1)))),pos=range[1],tck=tckk)
lines(c(xrange[1],xrange[2]),c(range[1],range[1]))
text(c(1:4),c(rep(range[1]-0.046*scal,3)),xpd=T,cex=1.1,col=cola,
c("Wt","Wt", "SCA1", "SCA1"))
ypo<-c(range[1]-0.1*scal,range[1]-0.1*scal)
text(c(1:4),c(ypo[1],ypo[2],ypo[1],ypo[2]),xpd=T,cex=1.1,col=cola,
c("0", "0", "edv"))
ypos<-0.038
xpos=1
points(xpos,ypos,pch=16,cex=1.4,col=rgb(0.4,0.4,0.4,alpha=0.6))
text(xpos+0.24*xscal,ypos,"Females",col="grey40",cex=1.2)
points(xpos,ypos-0.08*scal,pch=17,cex=1.4,col=rgb(0.4,0.4,0.4,alpha=0.6))
text(xpos+0.24*xscal,ypos-0.08*scal,"Males",col="grey40",cex=1.2)
### probability distributions --------
ypos<-0.6;xpos=1
range=c(0, 24)
scal = 12
par(mar=c(2.2,1,0,0))
dif<-data.frame(post_fix_cb$wt_contrast,post_fix_cb$b_groupeda.sca,
post_fix_cb$b_groupctrl.wt)
dif<-exp(dif)
yla<- "Fold-difference in CS in cb"
tckk=-0.018
ste<-seq(0.4,1.7,by=0.2)
mons_poste(dif,1)
zpos=seq(0,1,1/3)
xx=1;ind=0.4;xpol<-0.62
text(xpol,zpos[xx]+zpos[2]*ind,
"Wt: ed x 0 ",cex=1)
xx=xx+1
text(xpol,zpos[xx]+zpos[2]*ind,
"SCA1: ed x 0 ",cex=1)
xx=xx+1
text(xpol,zpos[xx]+zpos[2]*ind,
"0: Wt x SCA1",cex=1)
} ) )
In conclusion, there is neither evidence for the effect of genotype nor edaravone in terms of citrate synthase activity.
6 Session info
sessionInfo()R version 4.2.2 (2022-10-31 ucrt)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 10 x64 (build 19045)
Matrix products: default
locale:
[1] LC_COLLATE=English_United States.utf8
[2] LC_CTYPE=English_United States.utf8
[3] LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C
[5] LC_TIME=English_United States.utf8
attached base packages:
[1] stats graphics grDevices utils datasets methods base
other attached packages:
[1] knitr_1.42 vegan_2.6-4 lattice_0.20-45 permute_0.9-7
[5] gtsummary_1.7.0 bayesplot_1.10.0 lubridate_1.9.2 forcats_1.0.0
[9] stringr_1.5.0 dplyr_1.1.0 purrr_1.0.1 readr_2.1.4
[13] tidyr_1.3.0 tibble_3.1.8 tidyverse_2.0.0 ggplot2_3.4.1
[17] cowplot_1.1.1 car_3.1-1 carData_3.0-5 glmmTMB_1.1.5
[21] vioplot_0.4.0 zoo_1.8-11 sm_2.2-5.7.1 beeswarm_0.4.0
[25] brms_2.18.0 Rcpp_1.0.10
loaded via a namespace (and not attached):
[1] backports_1.4.1 plyr_1.8.8 igraph_1.4.1
[4] TMB_1.9.2 splines_4.2.2 crosstalk_1.2.0
[7] TH.data_1.1-2 rstantools_2.2.0 inline_0.3.19
[10] digest_0.6.31 htmltools_0.5.4 fansi_1.0.4
[13] magrittr_2.0.3 checkmate_2.1.0 cluster_2.1.4
[16] tzdb_0.3.0 RcppParallel_5.1.7 matrixStats_0.63.0
[19] xts_0.13.0 sandwich_3.0-2 timechange_0.2.0
[22] prettyunits_1.1.1 colorspace_2.1-0 xfun_0.37
[25] callr_3.7.3 crayon_1.5.2 jsonlite_1.8.4
[28] lme4_1.1-31 survival_3.5-3 glue_1.6.2
[31] gtable_0.3.1 emmeans_1.8.6 V8_4.2.2
[34] distributional_0.3.1 pkgbuild_1.4.0 rstan_2.26.16
[37] abind_1.4-5 scales_1.2.1 mvtnorm_1.1-3
[40] miniUI_0.1.1.1 xtable_1.8-4 stats4_4.2.2
[43] StanHeaders_2.26.16 DT_0.27 htmlwidgets_1.6.1
[46] threejs_0.3.3 posterior_1.4.0 ellipsis_0.3.2
[49] pkgconfig_2.0.3 loo_2.5.1 farver_2.1.1
[52] utf8_1.2.3 labeling_0.4.2 tidyselect_1.2.0
[55] rlang_1.1.1 reshape2_1.4.4 later_1.3.0
[58] munsell_0.5.0 tools_4.2.2 cli_3.6.0
[61] generics_0.1.3 evaluate_0.20 fastmap_1.1.1
[64] yaml_2.3.7 processx_3.8.0 nlme_3.1-162
[67] mime_0.12 projpred_2.6.0 xml2_1.3.3
[70] compiler_4.2.2 shinythemes_1.2.0 rstudioapi_0.14
[73] curl_5.0.0 gamm4_0.2-6 gt_0.9.0
[76] broom.helpers_1.13.0 stringi_1.7.12 ps_1.7.2
[79] Brobdingnag_1.2-9 Matrix_1.5-3 nloptr_2.0.3
[82] markdown_1.5 shinyjs_2.1.0 tensorA_0.36.2
[85] vctrs_0.5.2 pillar_1.8.1 lifecycle_1.0.3
[88] bridgesampling_1.1-2 estimability_1.4.1 httpuv_1.6.9
[91] R6_2.5.1 promises_1.2.0.1 gridExtra_2.3
[94] codetools_0.2-19 boot_1.3-28.1 colourpicker_1.2.0
[97] MASS_7.3-58.2 gtools_3.9.4 withr_2.5.0
[100] shinystan_2.6.0 multcomp_1.4-25 mgcv_1.8-41
[103] parallel_4.2.2 hms_1.1.2 grid_4.2.2
[106] coda_0.19-4 minqa_1.2.5 rmarkdown_2.23
[109] numDeriv_2016.8-1.1 shiny_1.7.4 base64enc_0.1-3
[112] dygraphs_1.1.1.6