Morning Administration Enhances Humoral Response to SARS-CoV-2 Vaccination in Kidney Transplant Recipients
Statistical report
Authors’ affiliations:
Ivan Zahradka1, Filip Tichanek2, Maria Magicova1, Istvan Modos2, Ondrej Viklicky1, Vojtech Petr1
1Department of Nephrology, Transplantation Center
2Department of Data Science
Institute for Clinical and Experimental Medicine, Prague, Czech Republic
Original publication
This is statistical report of for the publication Morning Administration Enhances Humoral Response to SARS-CoV-2 Vaccination in Kidney Transplant Recipients (Zahradka et al. 2024), published in the American Journal of Transplantation
When using this code, cite original publication:
Zahradka, I., Tichanek, F., Magicova, M., Modos, I., Viklicky, O., & Petr, V. (2024). Morning administration enhances humoral response to SARS-CoV-2 vaccination in kidney transplant recipients. American Journal of Transplantation. https://doi.org/10.1016/j.ajt.2024.03.004
Original GitHub repository: https://github.com/filip-tichanek/CovidTimeRTX
1 Introduction
This project aims to evaluate whether the timing of anti-COVID-19 mRNA vaccine administration (Pfizer or Moderna) affects the probability of seroconversion (\(IgG \geq 9.5 kAU/L\)) in Kidney Transplant Receivers (KTR). The project has been inspired by a study exploring the same phenomenon in hemodialysis patients by Lin and Hung in 2023 (Lin and Hung 2023).
All patients received two vaccine doses, with a time interval of 20-29 days between the doses. The impact of timing was evaluated for both the first and second doses. Seroconversion measurements were conducted upon the completion of the full vaccination regimen, specifically after the administration of the second dose.
The impact of vaccination timing was adjusted to account for other covariates known to influence the probability of anti-COVID seroconversion in the population of renal failure patients (Lin and Hung 2023) or KTRs (Magicova et al. 2022).
The study population comprises COVID-naive patients who have undergone kidney transplantation, with varying duration between the transplantation and anti-COVID vaccination. Inclusion criteria stipulated that only patients receiving two doses of mRNA vaccine in our institution, with a time interval between the two doses ranging from 20 to 29 days, were considered (a total of 553 patients).
As the time between the 2nd dose and blood sampling differed across patients, the time between the 2nd dose and blood collection (for IgG measurements) was also adjusted. Thus, it was included as a covariate with non-linear effect.
Timings of vaccinations were primarily treated as continuous numerical variable. Final models were fitted with Bayesian framework using ‘brms’ R package (Bürkner 2017), with all continuous variables standardized by 2 SD (as recommended by Gelman (Gelman 2008)).
Everything was done in R, version 4.3.2 (R Core Team 2023)
1.1 Variables recorded
We recorded predictors shown to affect the probability of seroconversion after mRNA anti-Covid vaccines in the population of patients with renal dysfunction (Lin and Hung 2023) and KTR (Magicova et al. 2022), plus few other covariates that we considered to be important for our specific context. These include:
Outcomes
SARS-CoV-2 IgG antibody level (\(AU/ml\)),
antibody_level. This does not serve as the primary outcome due to its non-normal distribution and the presence of values beyond detection limits (lower limit: 3.6 AU/ml, upper limit: 8000 AU/ml). Instances below the detection limits are denoted as “0,” while values surpassing the limits are represented as “8000.” However, it is employed as the outcome in sensitivity analysis, specifically in themodel_quantile. In this analysis, median values of log2[antibody_levels] were modeled using quantile regression, and values beyond detection limits were substituted with 2/3 and 3/2 of the minimum and maximum values, respectively.seroconversion: \(IgG > 9.5\) (0/1)seroconversion2: alternative definition ofseroconversion, showing if theantibody_levelis above detection limit, \(IgG > 3.6\) (0/1). It is used for sensitivity analysis examining robustness of the results when the definition of seroconversion changes.
Main predictors
time of the 1st dose of the antiCovid vaccine: continuous time in hours unit (
vacc_time_cont). In Bayesian models, continuous time was included as a variable standardized by 2 standard deviations (SD),vacc_time_scaltime of the 2nd dose of the antiCovid vaccine: continuous time in hours unit (
vacc_time_2nd_cont). In Bayesian models, continuous time was included as a variable standardized by 2 standard deviations (SD),vacc_time2_scal
For sensitivity analysis, we used dichotomized time of vaccine administration, afternon1 for the 1st and afternoon2 for the 2nd vaccination dose. The threshold for the afternoon was set to either 12pm or 1pm.
Other covariates
Moderna vaccine: 0: Phizer, 1: Moderna,
vaccine_modernacalcineurin inhibitors use (0/1),
calcineurin_inhibitor. his variable originated from two distinct variables, namelytacrolimusandciclosporinA, which were later merged into a single variable.steroids use (0/1),
steroidstime between renal transplantation and blood collection (for IgG measurement) in months,
months_afterTXtime between 2nd dose and blood sampling (in days), fitted as a predictor with non-linear effect,
days_after_2nd_dosesex (0: female, 1: male),
sexeGFR prior vaccination (\(mL/min/1.73m^2\))
EGFR_vaccinationmycophenolate mofetil or mycophenolic acid use (0/1),
mmf_mpadepleting therapy in the last year (0/1),
depletationTreat_yearbody mass index,
bmiserum Albumin (\(g/L\)),
Albuminelymphocytes (\(10^9/L\)), included to models in log-transformed form,
Lymphocytesdiabetes mellitus diagnosis (0/1),
DM
We primarily fitted a full model that included all the above-mentioned covariates.
Additional variables have also been recorded. While not automatically included in the final models, we ensured that their omission does not compromise the accuracy of the models. Alternatively, we utilized them to provide a broader overview of patients’ characteristics. These are: IgG_after_time (time of the blood sampling the measurement of SARS-CoV-2 IgG), sample_Lym_days (time distance between blood sampling for lymphocytes and sampling for IgG [days prior sampling for IgG]), Lymphocytes_time (time of blood sampling for lymphocytes count evaluation, in hours), tacrolimus_level (µg/l in serum, measured in patients taking tacrolimus), cyclospirin_level (µg/l in serum, measured in patients taking ciclosporin A), mTOR use (0/1), mmf_dose (dose of mycophenolate mofetil [MMF]. If mycophenolic acid [MPA, mg] was used instead, we standardize its amount to MMF so that MPA dose was multiplied accordingly \(MMF = MPA * (250/180)\)). mmf_total is as mmf_dose but NA values are replaced for ‘0’
1.2 Exploratory models
Exploratory models were fitted model in frequentist framework using ‘mgcv’ package(Wood 2011). Firstly, we evaluated the importance of a possible interaction between the two doses of vaccination. Comparisons of models via AIC and BIC did not show any evidence that the interaction improves out-of-sample predictive accuracy. Thus, the interaction (vacc_time_scal : vacc_time2_scal) was not included in the final models.
We used generalized additive models to explore possible non-linear relationship between seroconversion and vaccination times. However, as we observed straight lines, the effects of vaccination times (vacc_time_scal, vacc_time2_scal) were fitted as predictors with linear effects.
We predicted that the time elapsed between vaccination and blood collection (days_after_2nd_dose) could have a non-linear impact on seroconversion. Our initial GAM models supported this idea. To address this non-linear effect, we included it as a predictor with non-linear effect using a thin-plate spline limited to 4 knots.
We also investigated if adding more details could improve how well the model predicts outcomes. Specifically, we looked at factors like the mmf_dose, cyclosporin_level or tacrolimus_level, and their interactions with mmf_mpa, ciclosporinA or tacrolismus. We assessed the impact on predictive accuracy using Akaike (AIC) or Bayesian (BIC) information criteria. Considering mTOR’s role in circadian rhythms, we explored if there’s a significant interaction with vaccination time. Additionally, considering the circadian rhythmicity of IgG levels, we ensured that the results were not confounded by the time of blood sampling (IgG_after_time) for measuring SARS-CoV-2 IgG antibody levels.
1.3 Bayesian analysis and prior distribution
We decided to use the Bayesian approach as it enables to directly incorporate external knowledge to statistical models (in the form of so-called prior). Since the effect of morning vaccination on seroconversion has been already reported in hemodialysis patients (Lin and Hung 2023), another vulnerable group of patients, we believe that employing a methodology that takes these insights into account can contribute to more informed inferences.
The prior, or prior probability distribution, reflects our initial beliefs about model parameters (especially odds ratios in this study) before seeing any data. When using Gaussian zero-effect-centered priors, we give the highest probability to the possibility of zero effect. This choice makes our estimates more conservative, leaning towards smaller effects.
Similarly, if we incorporate external knowledge, like insights from a previous study on a similar topic, we adjust our estimates to align with already gained information. This means the results of the earlier study guide our model, incorporating past insights and modifying our estimates accordingly. It’s a way to integrate what we already know into our current analysis.
Mathematically, \(P(\beta | \mathbf{X}, Y)\) is proportional to the product of the prior distribution \(P(\beta)\) and the likelihood \(P(Y | \mathbf{X}, \beta)\):
\[ P(\beta | \mathbf{X}, Y) \propto P(Y | \mathbf{X}, \beta) \times P(\beta) \]
Before fitting Bayesian models, we standardized all numerical predictors by 2 standard deviations (SD). This ensures that all predictors are on the same scale, allowing the establishment of a single prior for all predictors and ensuring a consistent level of shrinkage for each one (Gelman 2008).
Priors for all variables except vaccination timing were established as Gaussian priors centered around zero effect, \(\mu = 0\), and with \(\sigma = 5\). For the effect of vaccination timing, we drew guidance from recent findings by (Lin and Hung 2023), who reported a 3.81-fold increase in seroconversion odds (95% CI: 1.59 to 9.15) on the 28th day post-vaccination and a 2.54-fold increase (CI: 1.15 to 5.61) on the 56th day post-vaccination for morning vaccinations, on patients with renal failure.
Considering the average and median interval between the 2nd vaccine dose and blood collection in our patient group (mean 50, median 47 days), we adopted the estimate for the 56th day post-vaccination as our prior.
As 2 SD of continuous vaccination time in our study is similar to the difference between the mean times of vaccination in the morning vs. afternoon (with the threshold set either to 12 pm or 1 pm), we expect that similar distribution may be present also in the previous study (Lin and Hung 2023) that was used for setting the prior.
To account for uncertainty effectively, we selected a sigma value 5 times the standard error of the estimate on the logit scale. Thus the prior for the 2SD-standardized vaccination times (vacc_time_scal, vacc_time2_scal) can be expressed as \(normal(\mu, \sigma)\), where
\[\mu = -\beta_{morning}\]
\[\sigma = SE[\beta_{morning}] * 5\]
with \(\beta_{morning}\) being estimated morning effect from (Lin and Hung 2023), 56th days post-vaccination, in log[Odds Ratio], and \(SE[\beta_{morning}]\) being standard error of estimated morning effect according to (Lin and Hung 2023).
In our perspective, this prior acknowledges differences in our study population and design but still navigate the posterior distribution to the results which are likely (given the previous study), balancing informed inference with inherent variability.
Specifically, for vacc_time_scal and vacc_time2_scal, prior was set to \(\mu = -0.9\) and \(\sigma = 2\).
When we used dichotomized time for a sensitivity analysis (afternoon1 and afteroon2), prior has \(\mu = 0.9\) (exactly as estimated by (Lin and Hung 2023)) and again \(\sigma = SE[\beta_{morning}] * 5\)
As we used ‘brms’ package(Bürkner 2017), employing Stan software for probabilistic sampling (Stan Development Team 2021), the model was estimated using Hamiltonian Monte Carlo sampling (this is Markov Chain Monte Carlo - MCMC - algorithm used for sampling from complex probability distributions, particularly effective in high-dimensional spaces). We use 6,000 iterations for each of 3 chains. From these, 2,000 are warm-ups, so posterior distributions are described by 12,000 posterior draws in total.
For each model, we present summary tables containing crucial details concerning the model’s convergence. The parameter rhat denotes the autocorrelation of posterior samples and ideally holds a value of 1.00. Additionally, we ensure an ample number of posterior draws, with Bulk_ESS and Tail_ESS both exceeding 1,000 to 2,000 for reliable results.
Next, we show effect size for each predictors. Take into account that continuous covariates/predictors were transformed by 2 standard deviations. Effect sizes are thus comparable.
1.4 Sensitivity analyses
Several sensitivity analyses were conducted to assess the robustness of our findings:
Model Specification Sensitivity Analysis (Model:
model_reduced): The original model was refitted with reduced covariate adjustments, incorporating only ‘important’ covariates identified in Table 1. Specifically, variables differing between seropositive and seronegative patients were included.Prior Sensitivity Analysis (
model_uniprior): To evaluate the impact of prior specifications, the original model was refitted using zero-effect-centered non-informative priors (\(\mu = 0\), \(\sigma = 5\)) for all predictors including vaccination times.Seroconversion Threshold Sensitivity Analysis (
model_seroconversion_min): The used an alternative definition of seroconversion (seroconversion2), defined as IgG levels exceeding the minimum detectable level of SARS-CoV-2 IgG (\(\text{IgG} > 3.6 , \text{AU/ml}\)).Dichotomized Vaccination Time Analysis:
model_dichotomized_12: Dichotomized vaccination times were used using a threshold of 12 pm.model_dichotomized_13: Another analysis was performed with a threshold of 1 pm.
Quantile Regression Analysis (
model_quantile): Median IgG levels (log2-transformed) were modeled using quantile regression (within frequentist framework). Values ofantibody_levelbelow or above detection limits were replaced with 2/3 and 3/2 of the lower and upper bound of detectable range respectively.
2 Analysis
2.1 Packages and functions
Open code
if (T) {rm(list = ls() )}
if (T) {
suppressWarnings(suppressMessages({
library(tidyverse)
library(stringr)
library(ggpubr)
library(emmeans)
library(gtsummary)
library(car)
library(sjPlot)
library(flextable)
library(openxlsx)
library(mgcv)
library(pROC)
library(cowplot)
library(boot)
library(glmnet)
library(brms)
library(projpred)
library(RJDBC)
library(janitor)
library(arm)
library(corrplot)
library(lubridate)
library(kableExtra)
library(ggdist)
library(bayesplot)
library(coxed)
library(quantreg)
library(ggbeeswarm)
# Functions clashes
select <- dplyr::select
rename <- dplyr::rename
mutate <- dplyr::mutate
recode <- dplyr::recode
summarize <- dplyr::summarize
count <- dplyr::count
# Simple math functions
logit <-function(x){log(x/(1-x))}
inv.logit <- function(x){exp(x)/(1+exp(x))}
}))
}2.1.1 run_model function
The function to load or run and save (if not done previously) models. It was copied from online proposal of Paul-Christian Bürkner.
Open code
run_model <- function(expr, path, reuse = TRUE) {
path <- paste0(path, ".Rds")
if (reuse) {
fit <- suppressWarnings(try(readRDS(path), silent = TRUE))
}
if (is(fit, "try-error")) {
fit <- eval(expr)
saveRDS(fit, file = path)
}
fit
}2.1.2 p_dir function
The function calculates proportion of posterior probability behind specified threshold (argument dir set to > or <) or probability of direction (dir = max, as defined in Makowski et al., 2019 (Makowski et al. 2019) ). data: data frame containing posterior draws (in columns). dir: probability under/above threshold (value > or <) or probability of direction (max). tres: threshold (e.g. value of zero effect, i.e. “1” for odds ratio)
Open code
p_dir <- function(data, dir, tres){
if(dir == 'max'){
return(
max(length(data[data>tres]),length(data[data<tres]))/length(data)
)
} else if(dir == '<'){
return(
length(data[data<tres])/length(data)
)
} else if(dir == '>'){
return(
length(data[data>tres])/length(data)
)
} else {print('ERROR')}
}2.1.3 repor function
The function takes brms model with logit-link and reports odds ratios and their 95% credible intervals. If the variables were transformed with rescale function of the ‘arm’ package (Gelman and Su 2022), the function clean the name of variables from redundant strings. If labels and scals are specified, the function will re-name all rows and transform the estimated effects and CI accordingly. The function prints kableExtra html table in default. If nhtm = TRUE, it is printed as data.frame
Open code
repor <- function(model, labels, scals, nhtm = FALSE){
df <- data.frame(
fixef(
model, probs = c(0.025, 0.975) )
) %>% select(
-Est.Error)
df <- df %>% filter(
!str_detect(rownames(df), "srescale"),
!str_detect(rownames(df), "Intercept"))
if(hasArg(labels) == TRUE){
labels = labels
for (i in 1:dim(df)[1] ){
rownames(df)[i] = labels[i]}
} else{
rownames(df) <- gsub(
'binary.inputsEQM0.50.5', '',
gsub(
'rescale', '', rownames(df)
)
)
}
if(hasArg(scals) == TRUE){
scals = scals
for(i in 1:dim(df)[1]){
df[i,] <- df[i,] / scals[i]
}
}
df <- round(exp(df), 3)
draw <- as.data.frame(model) %>% select(
matches('b_'), -b_Intercept)
df$PD <- rep(NA, dim(df)[1] )
for (i in 1:dim(draw)[2])
df$PD[i] <- round(p_dir(draw[,i], 'max', 0), 4)
colnames(df)[1] = 'OR'
if(nhtm == TRUE){
df} else{kableExtra::kable(df)
}
}2.1.4 sca and inv.sca functions
Functions either standardize vector (dat) by specified sd and mean (the function sca), or provide its reverse (inv.sca)
Open code
sca <- function(dat, mean, sd) {
(dat - mean) / (2*sd)}
inv.sca <- function(dat, mean, sd) {
(dat * (2*sd)) + mean }2.2 Data exploration
2.2.1 Import data
Note that antiboody_level (SARS-CoV-2 IgG antibody level [AU/ml]) measurement has a detection limits, restricted within the interval (3.6, 8000). If the measured value is below the detection limit, we wrote 0, if it is above detection limit, we write 8000
Open code
findat2 <- read.table("data/table2.txt")
## changing unit of EGFR (from ml/s/1.73 m2 to ml/m/1.73 m2)
findat2 <- findat2 %>% mutate(
EGFR_vaccination = EGFR_vaccination*60
)
summary(findat2)
## patient months_afterTX seroconversion antibody_level
## Min. : 38 Min. : 4.267 Min. :0.0000 Min. : 0.0
## 1st Qu.: 114048 1st Qu.: 26.533 1st Qu.:0.0000 1st Qu.: 0.0
## Median :1175387 Median : 74.267 Median :0.0000 Median : 0.0
## Mean : 724498 Mean : 94.561 Mean :0.3834 Mean : 155.1
## 3rd Qu.:1262352 3rd Qu.:140.267 3rd Qu.:1.0000 3rd Qu.: 31.5
## Max. :1345473 Max. :420.700 Max. :1.0000 Max. :8000.0
##
## vacc_time vacc_time_2nd days_after_2nd_dose age
## Length:553 Length:553 Min. : 14.00 Min. :25.89
## Class :character Class :character 1st Qu.: 29.00 1st Qu.:57.02
## Mode :character Mode :character Median : 47.00 Median :66.00
## Mean : 50.11 Mean :64.15
## 3rd Qu.: 67.00 3rd Qu.:71.19
## Max. :113.00 Max. :87.89
##
## DM BMI sex EGFR_vaccination
## Min. :0.0000 Min. :17.9 Min. :0.0000 Min. : 5.40
## 1st Qu.:0.0000 1st Qu.:24.8 1st Qu.:0.0000 1st Qu.: 34.80
## Median :0.0000 Median :28.0 Median :1.0000 Median : 47.40
## Mean :0.3002 Mean :28.3 Mean :0.6401 Mean : 49.23
## 3rd Qu.:1.0000 3rd Qu.:31.3 3rd Qu.:1.0000 3rd Qu.: 61.20
## Max. :1.0000 Max. :42.9 Max. :1.0000 Max. :109.80
##
## mmf_mpa depletationTreat_year Albumine Lymphocytes
## Min. :0.0000 Min. :0.0000 Min. :23.20 Min. : 0.200
## 1st Qu.:1.0000 1st Qu.:0.0000 1st Qu.:40.10 1st Qu.: 1.590
## Median :1.0000 Median :0.0000 Median :42.40 Median : 2.120
## Mean :0.7703 Mean :0.0434 Mean :42.07 Mean : 2.239
## 3rd Qu.:1.0000 3rd Qu.:0.0000 3rd Qu.:44.40 3rd Qu.: 2.740
## Max. :1.0000 Max. :1.0000 Max. :52.50 Max. :12.490
##
## tacrolimus ciclosporinA steroids vaccine_moderna
## Min. :0.0000 Min. :0.0000 Min. :0.0000 Min. :0.00000
## 1st Qu.:1.0000 1st Qu.:0.0000 1st Qu.:1.0000 1st Qu.:0.00000
## Median :1.0000 Median :0.0000 Median :1.0000 Median :0.00000
## Mean :0.7957 Mean :0.1013 Mean :0.9096 Mean :0.03074
## 3rd Qu.:1.0000 3rd Qu.:0.0000 3rd Qu.:1.0000 3rd Qu.:0.00000
## Max. :1.0000 Max. :1.0000 Max. :1.0000 Max. :1.00000
##
## IgG_after_time sample_Lym_days Lymphocytes_time tacrolimus_level
## Min. : 6.000 Min. : 0.00 Min. : 6.000 Min. : 2.100
## 1st Qu.: 6.850 1st Qu.: 47.00 1st Qu.: 6.900 1st Qu.: 5.600
## Median : 7.350 Median : 84.00 Median : 7.450 Median : 6.500
## Mean : 7.395 Mean : 74.35 Mean : 7.471 Mean : 6.646
## 3rd Qu.: 7.833 3rd Qu.: 91.00 3rd Qu.: 7.917 3rd Qu.: 7.600
## Max. :12.083 Max. :1172.00 Max. :13.617 Max. :17.100
## NA's :113
## cyclospirin_level mTOR mmf_dose
## Min. : 52.00 Min. :0.00000 Min. : 500
## 1st Qu.: 91.88 1st Qu.:0.00000 1st Qu.: 500
## Median :106.00 Median :0.00000 Median :1000
## Mean :112.76 Mean :0.07776 Mean :1040
## 3rd Qu.:137.00 3rd Qu.:0.00000 3rd Qu.:1500
## Max. :187.00 Max. :1.00000 Max. :2000
## NA's :497 NA's :127
## mmf dose
findat2$mmf_dose %>% as.factor() %>% table()
## .
## 500 750 1000 1500 2000
## 108 3 196 94 25There are several variables that may have strongly right-tailed distribution, especially antibody_level
Merging tacrolimus and ciclosporin1 to a single variable, calcineurin_inhibitor
Open code
findat2$calcineurin_inhibitor <- findat2$tacrolimus + findat2$ciclosporinA2.2.2 Preparation of table for models
Making vacc_times to be continuous variable
Open code
findat2 <- findat2 %>% mutate(
vacc_time_cont = as.numeric(seconds(hms(vacc_time)))/3600,
vacc_time_2nd_cont = as.numeric(seconds(hms(vacc_time_2nd)))/3600)Defining afternoon as the time before 1 pm. There will be 2 variables, afternoon1 indicating afternoon time of the 1st vaccine dose and afternoon2 for the 2nd dose
Open code
# defining 'morning' time as the time before median vaccination time
findat2 <- findat2 %>% mutate(
afternoon1 = if_else(vacc_time_cont > 13, 1, 0),
afternoon2 = if_else(vacc_time_2nd_cont > 13, 1, 0)) %>%
mutate(afternoon_comb = interaction(afternoon1, afternoon2)) %>%
mutate(afternoon_comb = recode(afternoon_comb,
'1.1' = 'afternoon-afternoon',
'1.0' = 'afternoon-morning',
'0.1' = 'morning-afternoon',
'0.0' = 'morning-morning'),
mmf_total = as.numeric(if_else(mmf_mpa == 0, 0, mmf_dose)))2.2.3 Characteristics according to seroconversion, table1
Open code
labe <- list(
vacc_time_cont ~ '1st vaccination time [hours]',
vacc_time_2nd_cont~ '2nd vaccination time [hours]',
antibody_level ~ 'SARS-CoV-2 IgG antibody level [AU/ml]',
days_after_2nd_dose~ 'Time after 2nd dose [days]',
months_afterTX~ 'Time after TX [months]',
age ~ 'Age [years]',
sex ~ 'Male sex [0/1]',
EGFR_vaccination ~ 'eGFR [mL/min/1.73 m2]',
mmf_mpa ~ 'MMF/MPA [0/1]',
depletationTreat_year ~ 'Depleting therapy [0/1]',
DM ~ 'Diabetes mellitus [0/1]',
BMI ~ 'BMI',
Albumine ~ 'Albumine [g/L]',
Lymphocytes ~ 'Lymphocyte counts [10^9/L]',
vaccine_moderna ~ 'Vaccine',
tacrolimus ~ 'Tacrolimus [0/1]',
ciclosporinA ~ 'CicloslorinA [0/1]',
calcineurin_inhibitor ~ 'Calcineurin inhibitor [0/1]',
steroids ~ 'Steroids [0/1]',
IgG_after_time ~ 'Blood sampling time for IgG [hours]',
#Lymphocytes_time ~ 'Blood sampling time of lymphocytes [hours]',
mmf_dose ~ 'MMF dose [mg]',
tacrolimus_level ~ 'Tacrolimus level [µg/l]',
cyclospirin_level ~ 'Ciclosporin level [µg/l]'
)
table1 <- findat2 %>% select(vacc_time_cont,
vacc_time_2nd_cont,
antibody_level,
days_after_2nd_dose,
months_afterTX,
age,
sex,
EGFR_vaccination,
mmf_mpa,
depletationTreat_year,
DM,
BMI,
Albumine,
Lymphocytes,
seroconversion,
vaccine_moderna,
tacrolimus,
ciclosporinA,
calcineurin_inhibitor,
steroids,
IgG_after_time,
#Lymphocytes_time,
mmf_dose,
tacrolimus_level,
cyclospirin_level
) %>% mutate(
seroconversion = as.factor(seroconversion),
vaccine_moderna = if_else( vaccine_moderna == 1, 'Moderna', 'Phizer')
) %>% tbl_summary(
label = labe,
by='seroconversion',
type = list(mmf_dose ~ "continuous")) %>%
modify_caption("Table 1. Patient characteristics according to seroconversion") %>%
add_p()
## Warning for variable 'cyclospirin_level':
## simpleWarning in wilcox.test.default(x = DATA[[1L]], y = DATA[[2L]], ...): cannot compute exact p-value with ties
table1| Characteristic | 0, N = 3411 | 1, N = 2121 | p-value2 |
|---|---|---|---|
| 1st vaccination time [hours] | 12.15 (11.15, 13.37) | 12.28 (11.34, 13.43) | 0.4 |
| 2nd vaccination time [hours] | 12.42 (11.40, 13.33) | 12.31 (11.16, 13.18) | 0.2 |
| SARS-CoV-2 IgG antibody level [AU/ml] | 0 (0, 0) | 57 (23, 171) | <0.001 |
| Time after 2nd dose [days] | 43 (28, 67) | 49 (31, 66) | 0.2 |
| Time after TX [months] | 55 (20, 115) | 100 (42, 170) | <0.001 |
| Age [years] | 67 (58, 72) | 63 (56, 71) | 0.017 |
| Male sex [0/1] | 202 (59%) | 152 (72%) | 0.003 |
| eGFR [mL/min/1.73 m2] | 44 (33, 58) | 51 (41, 68) | <0.001 |
| MMF/MPA [0/1] | 295 (87%) | 131 (62%) | <0.001 |
| Depleting therapy [0/1] | 22 (6.5%) | 2 (0.9%) | 0.002 |
| Diabetes mellitus [0/1] | 104 (30%) | 62 (29%) | 0.8 |
| BMI | 28.3 (24.8, 31.3) | 27.7 (24.9, 31.3) | 0.5 |
| Albumine [g/L] | 42.4 (40.1, 44.4) | 42.6 (40.2, 44.3) | 0.5 |
| Lymphocyte counts [10^9/L] | 1.99 (1.46, 2.47) | 2.37 (1.76, 3.13) | <0.001 |
| Vaccine | 0.2 | ||
| Moderna | 8 (2.3%) | 9 (4.2%) | |
| Phizer | 333 (98%) | 203 (96%) | |
| Tacrolimus [0/1] | 279 (82%) | 161 (76%) | 0.10 |
| CicloslorinA [0/1] | 27 (7.9%) | 29 (14%) | 0.029 |
| Calcineurin inhibitor [0/1] | 306 (90%) | 190 (90%) | >0.9 |
| Steroids [0/1] | 316 (93%) | 187 (88%) | 0.075 |
| Blood sampling time for IgG [hours] | 7.33 (6.83, 7.87) | 7.36 (6.87, 7.79) | 0.8 |
| MMF dose [mg] | 1,000 (875, 1,500) | 1,000 (500, 1,500) | 0.9 |
| Unknown | 46 | 81 | |
| Tacrolimus level [µg/l] | 6.50 (5.50, 7.60) | 6.30 (5.70, 7.40) | >0.9 |
| Unknown | 62 | 51 | |
| Ciclosporin level [µg/l] | 109 (90, 138) | 105 (95, 129) | 0.8 |
| Unknown | 314 | 183 | |
| 1 Median (IQR); n (%) | |||
| 2 Wilcoxon rank sum test; Pearson’s Chi-squared test | |||
2.2.4 Characteristics by dichotomized vaccination time, suppl_table2
Open code
labe <- list(
vacc_time_cont ~ '1st vaccination time [hours]',
vacc_time_2nd_cont~ '2nd vaccination time [hours]',
antibody_level ~ 'SARS-CoV-2 IgG antibody level [AU/ml]',
days_after_2nd_dose~ 'Time after 2nd dose [days]',
months_afterTX~ 'Time after TX [months]',
age ~ 'Age [years]',
sex ~ 'Male sex [0/1]',
EGFR_vaccination ~ 'eGFR [mL/min/1.73 m2]',
mmf_mpa ~ 'MMF/MPA [0/1]',
depletationTreat_year ~ 'Depleting therapy [0/1]',
DM ~ 'Diabetes mellitus [0/1]',
BMI ~ 'BMI',
Albumine ~ 'Albumine [g/L]',
Lymphocytes ~ 'Lymphocyte counts [10^9/L]',
seroconversion ~ 'Seroconversion [0/1]',
vaccine_moderna ~ 'Vaccine',
tacrolimus ~ 'Tacrolimus [0/1]',
ciclosporinA ~ 'CicloslorinA [0/1]',
calcineurin_inhibitor ~ 'Calcineurin inhibitor [0/1]',
steroids ~ 'Steroids [0/1]',
IgG_after_time ~ 'Blood sampling time for IgG [hours]',
#Lymphocytes_time ~ 'Blood sampling time of lymphocytes [hours]',
mmf_dose~ 'MMF dose [mg]',
tacrolimus_level ~ 'Tacrolimus level [µg/l]',
cyclospirin_level ~ 'Ciclosporin level [µg/l]'
)
suppl_table2 <- findat2 %>% select(vacc_time_cont,
vacc_time_2nd_cont,
antibody_level,
days_after_2nd_dose,
months_afterTX,
age,
sex,
EGFR_vaccination,
mmf_mpa,
depletationTreat_year,
DM,
BMI,
Albumine,
Lymphocytes,
seroconversion,
vaccine_moderna,
tacrolimus,
ciclosporinA,
calcineurin_inhibitor,
steroids,
IgG_after_time,
mmf_dose,
tacrolimus_level,
cyclospirin_level,
afternoon_comb
) %>% mutate(
vaccine_moderna = if_else( vaccine_moderna == 1, 'Moderna', 'Phizer')
) %>% tbl_summary(label = labe,
type = list(mmf_dose ~ "continuous"),
by='afternoon_comb') %>%
modify_caption("Supplemetary Table 2. Patient characteristics according to dichotomized time of vaccination (1st dose - 2nd dose). Threshold for 'afternoon' was set to 1 pm") %>%
add_p()
suppl_table2| Characteristic | morning-morning, N = 3011 | afternoon-morning, N = 721 | morning-afternoon, N = 791 | afternoon-afternoon, N = 1011 | p-value2 |
|---|---|---|---|---|---|
| 1st vaccination time [hours] | 11.63 (10.78, 12.20) | 13.74 (13.43, 14.08) | 11.73 (11.05, 12.38) | 13.95 (13.47, 14.33) | <0.001 |
| 2nd vaccination time [hours] | 11.67 (10.63, 12.35) | 11.97 (11.21, 12.49) | 13.48 (13.23, 13.98) | 13.78 (13.40, 14.12) | <0.001 |
| SARS-CoV-2 IgG antibody level [AU/ml] | 0 (0, 36) | 6 (0, 40) | 0 (0, 13) | 0 (0, 31) | 0.2 |
| Time after 2nd dose [days] | 47 (28, 65) | 40 (26, 62) | 69 (51, 92) | 36 (27, 55) | <0.001 |
| Time after TX [months] | 64 (21, 127) | 81 (27, 150) | 103 (47, 182) | 73 (28, 131) | 0.003 |
| Age [years] | 66 (57, 71) | 62 (55, 69) | 73 (71, 76) | 62 (54, 67) | <0.001 |
| Male sex [0/1] | 197 (65%) | 47 (65%) | 46 (58%) | 64 (63%) | 0.7 |
| eGFR [mL/min/1.73 m2] | 47 (36, 61) | 44 (33, 59) | 45 (32, 57) | 52 (38, 67) | 0.087 |
| MMF/MPA [0/1] | 233 (77%) | 53 (74%) | 56 (71%) | 84 (83%) | 0.2 |
| Depleting therapy [0/1] | 16 (5.3%) | 3 (4.2%) | 0 (0%) | 5 (5.0%) | 0.2 |
| Diabetes mellitus [0/1] | 92 (31%) | 19 (26%) | 31 (39%) | 24 (24%) | 0.13 |
| BMI | 28.4 (25.0, 31.8) | 27.4 (23.7, 30.8) | 28.1 (25.3, 30.8) | 27.8 (24.5, 30.8) | 0.2 |
| Albumine [g/L] | 42.4 (40.0, 44.3) | 42.9 (40.4, 44.6) | 41.3 (40.0, 43.6) | 43.0 (40.9, 44.6) | 0.038 |
| Lymphocyte counts [10^9/L] | 2.08 (1.56, 2.69) | 2.03 (1.63, 2.71) | 2.18 (1.57, 2.84) | 2.28 (1.72, 2.82) | 0.5 |
| Seroconversion [0/1] | 120 (40%) | 32 (44%) | 25 (32%) | 35 (35%) | 0.3 |
| Vaccine | 0.12 | ||||
| Moderna | 7 (2.3%) | 2 (2.8%) | 1 (1.3%) | 7 (6.9%) | |
| Phizer | 294 (98%) | 70 (97%) | 78 (99%) | 94 (93%) | |
| Tacrolimus [0/1] | 245 (81%) | 59 (82%) | 52 (66%) | 84 (83%) | 0.013 |
| CicloslorinA [0/1] | 28 (9.3%) | 5 (6.9%) | 19 (24%) | 4 (4.0%) | <0.001 |
| Calcineurin inhibitor [0/1] | 273 (91%) | 64 (89%) | 71 (90%) | 88 (87%) | 0.8 |
| Steroids [0/1] | 271 (90%) | 69 (96%) | 70 (89%) | 93 (92%) | 0.4 |
| Blood sampling time for IgG [hours] | 7.28 (6.78, 7.80) | 7.53 (6.92, 7.95) | 7.52 (7.07, 7.79) | 7.33 (6.88, 7.93) | 0.12 |
| MMF dose [mg] | 1,000 (500, 1,500) | 1,000 (500, 1,500) | 1,000 (875, 1,000) | 1,000 (1,000, 1,500) | 0.5 |
| Unknown | 68 | 19 | 23 | 17 | |
| Tacrolimus level [µg/l] | 6.50 (5.60, 7.30) | 6.40 (5.65, 7.65) | 6.30 (5.33, 7.90) | 6.50 (5.68, 7.73) | >0.9 |
| Unknown | 56 | 13 | 27 | 17 | |
| Ciclosporin level [µg/l] | 110 (91, 138) | 105 (85, 106) | 106 (92, 135) | 135 (114, 155) | 0.3 |
| Unknown | 273 | 67 | 60 | 97 | |
| 1 Median (IQR); n (%) | |||||
| 2 Kruskal-Wallis rank sum test; Pearson’s Chi-squared test; Fisher’s exact test | |||||
2.2.5 Characteristics of continuous variables
Open code
sumtab <- findat2 %>% mutate(
log_Lymphocytes = log(Lymphocytes),
# Making *vacc_times* to be continuous variable
# ! Note that the time will be in hours
vacc_time_cont = as.numeric(seconds(hms(vacc_time)))/3600,
vacc_time_2nd_cont = as.numeric(seconds(hms(vacc_time_2nd)))/3600
) %>% select(
c(months_afterTX, antibody_level,
days_after_2nd_dose, age, BMI,
EGFR_vaccination, Albumine, Lymphocytes, mmf_total,
log_Lymphocytes,
vacc_time_cont,
vacc_time_2nd_cont,
IgG_after_time) ) %>%
pivot_longer(cols = everything()) %>%
group_by(name) %>%
summarise(mean = mean(value) %>% round(2),
std = sd(value) %>% round(2),
std_2x = 2*sd(value) %>% round(2),
median = median(value) %>% round(2),
Q1 = quantile(value, probs = c(0.25)) %>% round(2), # rozdil mezi 1. a 3. quartilem
Q3 = quantile(value, probs = c(0.75)) %>% round(2),
min = min(value) %>% round(2),
max = max(value) %>% round(2))
kableExtra::kable(sumtab) | name | mean | std | std_2x | median | Q1 | Q3 | min | max |
|---|---|---|---|---|---|---|---|---|
| Albumine | 42.07 | 3.69 | 7.38 | 42.40 | 40.10 | 44.40 | 23.20 | 52.50 |
| BMI | 28.30 | 4.84 | 9.68 | 28.00 | 24.80 | 31.30 | 17.90 | 42.90 |
| EGFR_vaccination | 49.23 | 19.87 | 39.74 | 47.40 | 34.80 | 61.20 | 5.40 | 109.80 |
| IgG_after_time | 7.40 | 0.77 | 1.54 | 7.35 | 6.85 | 7.83 | 6.00 | 12.08 |
| Lymphocytes | 2.24 | 1.01 | 2.02 | 2.12 | 1.59 | 2.74 | 0.20 | 12.49 |
| age | 64.15 | 10.16 | 20.32 | 66.00 | 57.02 | 71.19 | 25.89 | 87.89 |
| antibody_level | 155.11 | 727.22 | 1454.44 | 0.00 | 0.00 | 31.50 | 0.00 | 8000.00 |
| days_after_2nd_dose | 50.11 | 24.79 | 49.58 | 47.00 | 29.00 | 67.00 | 14.00 | 113.00 |
| log_Lymphocytes | 0.71 | 0.45 | 0.90 | 0.75 | 0.46 | 1.01 | -1.61 | 2.52 |
| mmf_total | 801.54 | 572.48 | 1144.96 | 1000.00 | 500.00 | 1000.00 | 0.00 | 2000.00 |
| months_afterTX | 94.56 | 82.10 | 164.20 | 74.27 | 26.53 | 140.27 | 4.27 | 420.70 |
| vacc_time_2nd_cont | 12.21 | 1.48 | 2.96 | 12.38 | 11.30 | 13.32 | 8.75 | 17.17 |
| vacc_time_cont | 12.28 | 1.40 | 2.80 | 12.20 | 11.30 | 13.38 | 7.17 | 17.37 |
2.2.6 Correlations between predictors
Open code
findat2 <- findat2 %>% mutate(
log_Lymphocytes = log(Lymphocytes)
)
cordat <- findat2 %>% select(vacc_time_cont,
vacc_time_2nd_cont,
antibody_level,
days_after_2nd_dose,
months_afterTX,
age,
sex,
EGFR_vaccination,
mmf_mpa,
depletationTreat_year,
DM,
BMI,
Albumine,
Lymphocytes,
tacrolimus,
ciclosporinA,
steroids,
IgG_after_time
)
sds <- findat2 %>% select(vacc_time_cont,
vacc_time_2nd_cont,
days_after_2nd_dose,
months_afterTX,
age,
sex,
EGFR_vaccination,
mmf_mpa,
depletationTreat_year,
DM,
BMI,
Albumine,
Lymphocytes,
log_Lymphocytes) %>% summarise_all(sd)*2
# correlations visually
corrplot(cor(cordat, method='spearman'))
Open code
# correlations in numbers
cor(cordat, method='spearman') %>% round(2)
## vacc_time_cont vacc_time_2nd_cont antibody_level
## vacc_time_cont 1.00 0.41 0.05
## vacc_time_2nd_cont 0.41 1.00 -0.06
## antibody_level 0.05 -0.06 1.00
## days_after_2nd_dose -0.24 0.09 0.04
## months_afterTX -0.01 0.13 0.25
## age -0.35 0.14 -0.15
## sex 0.01 0.01 0.16
## EGFR_vaccination 0.04 -0.01 0.20
## mmf_mpa 0.02 0.02 -0.27
## depletationTreat_year 0.03 -0.08 -0.15
## DM -0.06 0.03 -0.03
## BMI -0.03 -0.05 -0.01
## Albumine 0.07 -0.02 0.06
## Lymphocytes 0.03 0.10 0.24
## tacrolimus 0.04 -0.09 -0.05
## ciclosporinA -0.07 0.08 0.08
## steroids 0.02 -0.03 -0.11
## IgG_after_time 0.10 0.05 -0.03
## days_after_2nd_dose months_afterTX age sex
## vacc_time_cont -0.24 -0.01 -0.35 0.01
## vacc_time_2nd_cont 0.09 0.13 0.14 0.01
## antibody_level 0.04 0.25 -0.15 0.16
## days_after_2nd_dose 1.00 0.25 0.47 -0.01
## months_afterTX 0.25 1.00 0.21 0.04
## age 0.47 0.21 1.00 -0.06
## sex -0.01 0.04 -0.06 1.00
## EGFR_vaccination 0.01 -0.15 -0.13 0.08
## mmf_mpa -0.05 -0.16 -0.15 0.07
## depletationTreat_year -0.14 -0.32 -0.07 -0.03
## DM 0.15 -0.05 0.23 0.08
## BMI 0.01 -0.05 -0.07 0.03
## Albumine -0.08 -0.17 -0.22 0.04
## Lymphocytes 0.00 0.21 -0.01 0.01
## tacrolimus -0.10 -0.33 -0.26 -0.03
## ciclosporinA 0.08 0.26 0.20 0.00
## steroids -0.10 -0.32 -0.15 -0.03
## IgG_after_time -0.06 -0.07 -0.04 -0.08
## EGFR_vaccination mmf_mpa depletationTreat_year DM
## vacc_time_cont 0.04 0.02 0.03 -0.06
## vacc_time_2nd_cont -0.01 0.02 -0.08 0.03
## antibody_level 0.20 -0.27 -0.15 -0.03
## days_after_2nd_dose 0.01 -0.05 -0.14 0.15
## months_afterTX -0.15 -0.16 -0.32 -0.05
## age -0.13 -0.15 -0.07 0.23
## sex 0.08 0.07 -0.03 0.08
## EGFR_vaccination 1.00 0.19 0.03 0.00
## mmf_mpa 0.19 1.00 0.05 -0.05
## depletationTreat_year 0.03 0.05 1.00 0.02
## DM 0.00 -0.05 0.02 1.00
## BMI 0.01 0.10 0.09 0.25
## Albumine 0.21 0.23 -0.03 -0.07
## Lymphocytes 0.16 -0.06 -0.22 0.04
## tacrolimus 0.20 0.18 0.04 0.02
## ciclosporinA -0.17 -0.09 -0.04 0.00
## steroids -0.01 -0.02 0.07 0.00
## IgG_after_time 0.02 0.03 0.00 -0.04
## BMI Albumine Lymphocytes tacrolimus ciclosporinA
## vacc_time_cont -0.03 0.07 0.03 0.04 -0.07
## vacc_time_2nd_cont -0.05 -0.02 0.10 -0.09 0.08
## antibody_level -0.01 0.06 0.24 -0.05 0.08
## days_after_2nd_dose 0.01 -0.08 0.00 -0.10 0.08
## months_afterTX -0.05 -0.17 0.21 -0.33 0.26
## age -0.07 -0.22 -0.01 -0.26 0.20
## sex 0.03 0.04 0.01 -0.03 0.00
## EGFR_vaccination 0.01 0.21 0.16 0.20 -0.17
## mmf_mpa 0.10 0.23 -0.06 0.18 -0.09
## depletationTreat_year 0.09 -0.03 -0.22 0.04 -0.04
## DM 0.25 -0.07 0.04 0.02 0.00
## BMI 1.00 0.00 0.07 0.08 0.00
## Albumine 0.00 1.00 0.12 0.18 -0.08
## Lymphocytes 0.07 0.12 1.00 -0.01 0.08
## tacrolimus 0.08 0.18 -0.01 1.00 -0.66
## ciclosporinA 0.00 -0.08 0.08 -0.66 1.00
## steroids -0.02 0.04 0.03 0.14 -0.14
## IgG_after_time -0.09 0.04 0.01 0.00 0.03
## steroids IgG_after_time
## vacc_time_cont 0.02 0.10
## vacc_time_2nd_cont -0.03 0.05
## antibody_level -0.11 -0.03
## days_after_2nd_dose -0.10 -0.06
## months_afterTX -0.32 -0.07
## age -0.15 -0.04
## sex -0.03 -0.08
## EGFR_vaccination -0.01 0.02
## mmf_mpa -0.02 0.03
## depletationTreat_year 0.07 0.00
## DM 0.00 -0.04
## BMI -0.02 -0.09
## Albumine 0.04 0.04
## Lymphocytes 0.03 0.01
## tacrolimus 0.14 0.00
## ciclosporinA -0.14 0.03
## steroids 1.00 0.02
## IgG_after_time 0.02 1.002.2.7 Distribution of continuous variables
Open code
suppressMessages({
p1 <- findat2 %>% ggplot(aes(x = months_afterTX)) + geom_histogram()
p2 <- findat2 %>% ggplot(aes(x = antibody_level)) + geom_histogram()
p25 <- findat2 %>% ggplot(aes(x = log2(antibody_level+.1))) + geom_histogram()
p3 <- findat2 %>% ggplot(aes(x = days_after_2nd_dose)) + geom_histogram()
p4 <- findat2 %>% ggplot(aes(x = age)) + geom_histogram()
p5 <- findat2 %>% ggplot(aes(x = BMI)) + geom_histogram()
p6 <- findat2 %>% ggplot(aes(x = EGFR_vaccination ) ) + geom_histogram()
p7 <- findat2 %>% ggplot(aes(x = Albumine ) ) + geom_histogram()
p8 <- findat2 %>% ggplot(aes(x = Lymphocytes ) ) + geom_histogram()
p8b <- findat2 %>% ggplot(aes(x = log_Lymphocytes ) ) + geom_histogram()
p9 <- findat2 %>% ggplot(aes(x = vacc_time_cont) ) + geom_histogram()
p11 <- findat2 %>% ggplot(aes(x = vacc_time_2nd_cont) ) + geom_histogram()
p12 <- findat2 %>% ggplot(aes(x = IgG_after_time) ) + geom_histogram()
p13 <- findat2 %>% ggplot(aes(x = Lymphocytes_time) ) + geom_histogram()
p14 <- findat2 %>% ggplot(aes(x = tacrolimus_level) ) + geom_histogram()
p15 <- findat2 %>% ggplot(aes(x = cyclospirin_level) ) + geom_histogram()
ggpubr::ggarrange(p1, p2, p25, p3, p4, p5, p6,
p7, p8, p8b, p9, p11,
p12, p13, p14, p15)
})
## Warning: Removed 113 rows containing non-finite values (`stat_bin()`).
## Warning: Removed 497 rows containing non-finite values (`stat_bin()`).
2.2.8 Histogram of vaccination times (`suppl_fig2’)
Histograms show counts
Open code
cole <- c('#CD7006', '#0028F0')
# Manually define the bin edges
bin_edges <- c(seq(6.999, 17.499, by = 0.5))
suppl_fig2count <- findat2 %>%
select(seroconversion, vacc_time_cont, vacc_time_2nd_cont) %>%
pivot_longer(cols = c(vacc_time_cont, vacc_time_2nd_cont),
names_to = 'dose',
values_to = 'time') %>%
mutate(dose = if_else(dose == 'vacc_time_cont', '1st dose', '2nd dose'),
seroconversion = if_else(seroconversion == 1, 'Successful seroconversion',
'Unsuccessful seroconversion')) %>%
ggplot(aes(x = time, fill = dose)) +
geom_histogram(breaks = bin_edges, color = "grey50") +
facet_wrap(~ seroconversion + dose, nrow = 2) +
xlab("Vaccination time [hours]") +
ylab("Count") +
scale_fill_manual(values = cole)
suppl_fig2count
An alternative - histograms show percentages of patients
Open code
cole <- c('#CD7006', '#0028F0')
# Manually define the bin edges
bin_edges <- c(seq(6.999, 17.499, by = 0.5))
suppl_fig2a <- findat2 %>%
select(seroconversion, vacc_time_cont, vacc_time_2nd_cont) %>%
pivot_longer(cols = c(vacc_time_cont, vacc_time_2nd_cont),
names_to = 'dose',
values_to = 'time') %>%
mutate(dose = if_else(dose == 'vacc_time_cont', '1st dose', '2nd dose'),
seroconversion = if_else(seroconversion == 1, 'Successful seroconversion',
'Unsuccessful seroconversion')) %>%
mutate(count2 = if_else(seroconversion == 'Successful seroconversion',
212, 341)) %>%
filter(seroconversion == 'Successful seroconversion') %>%
ggplot(aes(x = time, fill = dose)) +
geom_histogram(breaks = bin_edges, color = "grey50",
aes(y=stat(count)/2.12)) +
facet_wrap(~ dose+ seroconversion , nrow = 1) +
xlab("") +
ylab("% of patients") +
scale_fill_manual(values = cole)+
coord_cartesian(xlim = c(6.99, 17.5), ylim = c(0, 18))
suppl_fig2b <- findat2 %>%
select(seroconversion, vacc_time_cont, vacc_time_2nd_cont) %>%
pivot_longer(cols = c(vacc_time_cont, vacc_time_2nd_cont),
names_to = 'dose',
values_to = 'time') %>%
mutate(dose = if_else(dose == 'vacc_time_cont', '1st dose', '2nd dose'),
seroconversion = if_else(seroconversion == 1, 'Successful seroconversion',
'Unsuccessful seroconversion')) %>%
mutate(count2 = if_else(seroconversion == 'Successful seroconversion',
212, 341)) %>%
filter(seroconversion == 'Unsuccessful seroconversion') %>%
ggplot(aes(x = time, fill = dose)) +
geom_histogram(breaks = bin_edges, color = "grey50",
aes(y=stat(count)/3.41)) +
facet_wrap(~ dose + seroconversion, nrow = 1) +
xlab("Vaccination time [hours]") +
ylab("% of patients") +
scale_fill_manual(values = cole)+
coord_cartesian(xlim = c(6.99, 17.5), ylim = c(0, 18))
prow <- plot_grid(
suppl_fig2a + theme(legend.position="none"),
suppl_fig2b + theme(legend.position="none"), ncol=1)
legend_b <- get_legend(
suppl_fig2b +
guides(color = guide_legend(nrow = 1)) +
theme(legend.position = "right")
)
suppl_fig2 <- plot_grid(prow, legend_b, ncol = 2, rel_widths = c(1, .2))
suppl_fig2
2.2.9 Histogram of antibody levels (`suppl_fig3’)
Open code
suppl_fig3 <- findat2 %>%
select(antibody_level, afternoon1, afternoon2, afternoon_comb, seroconversion) %>%
filter(antibody_level>9.5) %>%
mutate(antibody_level = if_else(antibody_level == 0, 2.4, antibody_level),
seroconversion = factor(seroconversion)) %>%
ggplot(aes(y=antibody_level)) +
geom_boxplot(aes(x=seroconversion), outlier.shape = NA) +
geom_beeswarm(aes(x=seroconversion,dodge.width = 2),
pch=1, cex= 3, size=2, alpha=0.8) +
scale_y_continuous(trans = 'log10') +
ylab("SARS-CoV-2 IgG antibody level [AU/ml], log10-scale") +
ggtitle("IgG titers in patients with successful seroconversion") +
coord_cartesian(xlim=c(0.5, 1.5)) +
theme(plot.title = element_text(size = 11))
suppl_fig3
Alternative
Open code
suppl_fig3_alt <- findat2 %>%
select(antibody_level, afternoon1, afternoon2, afternoon_comb, seroconversion) %>%
filter(antibody_level>3) %>%
mutate(antibody_level = if_else(antibody_level == 0, 2.4, antibody_level),
seroconversion = factor(seroconversion)) %>%
ggplot(aes(y=antibody_level, color=seroconversion)) +
geom_boxplot(aes(x=seroconversion), outlier.shape = NA) +
geom_beeswarm(aes(x=seroconversion), alpha=0.99, pch=1) +
scale_y_continuous(trans = 'log10') +
ylab("SARS-CoV-2 IgG antibody level [AU/ml], log10-scale") +
ggtitle("IgG titers, exluding values under detection limit") +
geom_hline(yintercept = 3.6, color = "grey20", linetype = "dashed") +
annotate("text", x = 1, y = 3, label = "Detection limit",
color = "grey20")
suppl_fig3_alt
2.3 Explorative (non-Bayesian) models
2.3.1 Does interaction between timings of two doses matter? suppl_table1
Open code
# no interaction
model_all <- gam(seroconversion ~
vacc_time_cont +
vacc_time_2nd_cont +
vaccine_moderna+calcineurin_inhibitor+steroids+
months_afterTX +
s(days_after_2nd_dose, k=4) +
age +
sex +
EGFR_vaccination +
mmf_mpa +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes),
data = findat2,
family = binomial, method='ML')
# with interaction
model_all_int <- gam(seroconversion ~
vacc_time_cont *
vacc_time_2nd_cont +
vaccine_moderna+calcineurin_inhibitor+steroids+
months_afterTX +
s(days_after_2nd_dose, k=4) +
age +
sex +
EGFR_vaccination +
mmf_mpa +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes),
data = findat2,
family = binomial, method='ML')
tab_model(model_all,
model_all_int,
title = "Seroconversion in response to anti-COVID-19 mRNA vaccines at patients after kidney TX; continuous daytime",
dv.labels = "" ,show.intercept = FALSE)| Predictors | Odds Ratios | CI | p | Odds Ratios | CI | p |
| vacc time cont | 1.12 | 0.94 – 1.35 | 0.215 | 1.12 | 0.40 – 3.11 | 0.829 |
| vacc time 2nd cont | 0.83 | 0.71 – 0.99 | 0.035 | 0.83 | 0.29 – 2.40 | 0.734 |
| vaccine moderna | 1.39 | 0.42 – 4.63 | 0.592 | 1.39 | 0.39 – 4.97 | 0.614 |
| calcineurin inhibitor | 1.76 | 0.80 – 3.91 | 0.161 | 1.76 | 0.80 – 3.91 | 0.161 |
| steroids | 0.74 | 0.36 – 1.55 | 0.427 | 0.74 | 0.36 – 1.55 | 0.427 |
| months afterTX | 1.01 | 1.00 – 1.01 | <0.001 | 1.01 | 1.00 – 1.01 | <0.001 |
| age | 0.97 | 0.94 – 1.00 | 0.020 | 0.97 | 0.94 – 1.00 | 0.021 |
| sex | 1.97 | 1.27 – 3.07 | 0.003 | 1.97 | 1.27 – 3.07 | 0.003 |
| EGFR vaccination | 1.03 | 1.02 – 1.04 | <0.001 | 1.03 | 1.02 – 1.04 | <0.001 |
| mmf mpa | 0.13 | 0.08 – 0.22 | <0.001 | 0.13 | 0.08 – 0.22 | <0.001 |
| depletationTreat year | 0.28 | 0.05 – 1.46 | 0.131 | 0.28 | 0.05 – 1.46 | 0.132 |
| DM | 0.92 | 0.57 – 1.48 | 0.720 | 0.92 | 0.57 – 1.48 | 0.721 |
| BMI | 1.00 | 0.95 – 1.04 | 0.867 | 1.00 | 0.95 – 1.04 | 0.867 |
| Albumine | 1.03 | 0.97 – 1.10 | 0.312 | 1.03 | 0.97 – 1.10 | 0.312 |
| Lymphocytes [log] | 2.47 | 1.46 – 4.18 | 0.001 | 2.47 | 1.46 – 4.18 | 0.001 |
| Smooth term (days after 2nd dose) |
0.179 | 0.181 | ||||
| vacc time cont × vacc time 2nd cont |
1.00 | 0.92 – 1.09 | 0.996 | |||
| Observations | 553 | 553 | ||||
| R2 | 0.259 | 0.257 | ||||
Comparison of models
Open code
ai <- AIC( model_all, model_all_int)
bi <- BIC(model_all, model_all_int)
suppl_table1 <- round(cbind(ai, BIC=bi[,2]), 1)
suppl_table1 <- kableExtra::kable(suppl_table1, caption =
"Supplementary Table 1. Comparison of models including ('model_all_int') vs. ignoring ('model_all') interaction between the timings of the 1st and the 2nd doses of anti-Covid mRNA vaccination via Akaike information criterion (AIC) and Bayesian Information Criterion (BIC). Models were fitted with continuous (non-binarized) times of vaccination. Both AIC and BIC indicates out-of-sample predictive accuracy, with lower value indicating better predictive accuracy. As models with interaction have larger AIC and BIC values, interaction term seems to be redundant and will not be included to final models. See https://github.com/filip-tichanek/CovidTimeRTX for details")
suppl_table1| df | AIC | BIC | |
|---|---|---|---|
| model_all | 17.7 | 604 | 680.5 |
| model_all_int | 18.7 | 606 | 686.9 |
Inclusion of interaction likely does not improve out-of-sample predictive accuracy and can be avoided
2.3.2 Is the effect of vaccination timing non-linear? suppl_figure1
Lets look if inclusion of non-linear morning vaccination effect improve the fit and out-of-sample predictive accuracy
Open code
model_all_nl <- gam(seroconversion ~
s(vacc_time_cont, k=4) +
s(vacc_time_2nd_cont, k=4) +
vaccine_moderna+calcineurin_inhibitor+steroids+
months_afterTX +
s(days_after_2nd_dose, k=4) +
age +
sex +
EGFR_vaccination +
mmf_mpa +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes),
data = findat2,
family = binomial, method='ML')
summary(model_all_nl)
##
## Family: binomial
## Link function: logit
##
## Formula:
## seroconversion ~ s(vacc_time_cont, k = 4) + s(vacc_time_2nd_cont,
## k = 4) + vaccine_moderna + calcineurin_inhibitor + steroids +
## months_afterTX + s(days_after_2nd_dose, k = 4) + age + sex +
## EGFR_vaccination + mmf_mpa + depletationTreat_year + DM +
## BMI + Albumine + log(Lymphocytes)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.648943 1.868741 -0.882 0.377570
## vaccine_moderna 0.329072 0.613956 0.536 0.591968
## calcineurin_inhibitor 0.568031 0.405410 1.401 0.161176
## steroids -0.298201 0.375071 -0.795 0.426583
## months_afterTX 0.006058 0.001566 3.869 0.000109 ***
## age -0.030845 0.013248 -2.328 0.019895 *
## sex 0.679310 0.225519 3.012 0.002594 **
## EGFR_vaccination 0.030996 0.005821 5.325 1.01e-07 ***
## mmf_mpa -2.037649 0.273102 -7.461 8.58e-14 ***
## depletationTreat_year -1.263213 0.836283 -1.511 0.130914
## DM -0.087812 0.245060 -0.358 0.720096
## BMI -0.003811 0.022690 -0.168 0.866626
## Albumine 0.031308 0.030975 1.011 0.312138
## log(Lymphocytes) 0.904007 0.268759 3.364 0.000769 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(vacc_time_cont) 1.000 1.000 1.535 0.2153
## s(vacc_time_2nd_cont) 1.000 1.000 4.451 0.0349 *
## s(days_after_2nd_dose) 1.443 1.743 2.207 0.1786
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.259 Deviance explained = 22.8%
## -ML = 285.47 Scale est. = 1 n = 553No, penalized maximum likelihood shrunk the non-linear effects of vaccine administration times to a straight line. See plot
Open code
par(mfrow=c(1,3))
par(mar=c(3,3,1,1))
par(mgp=c(1.5,0.5,0))
plot.gam(model_all_nl, ylab=c('Partial residuals'))
Non-linear effect was shrunk by penalized spline regression to straight line. Models with non-linear effect for vacc_time_cont and vacc_time_2nd_cont are thus exactly the same as when fitted as linear predictors. There is no reason to include timings of vaccine administration as non-linear effect.
On the other hand, the effect of days_after_2nd_dose may be non-linear, although the evidence is not very strong
2.3.3 Does the time of blood sampling for the IgG levels measurements matters?
Even the IgG levels may show circadian rhythms and the timing of the blood collection thus may matter.
Open code
model_all_st <- gam(seroconversion ~
vacc_time_cont +
vacc_time_2nd_cont +
vaccine_moderna+calcineurin_inhibitor+steroids+
months_afterTX +
s(days_after_2nd_dose, k=4) +
age +
sex +
EGFR_vaccination +
mmf_mpa +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes)+
IgG_after_time,
data = findat2,
family = binomial, method = "ML")
BIC(model_all, model_all_st)
## df BIC
## model_all 17.74275 680.5494
## model_all_st 18.76908 685.8333
AIC(model_all, model_all_st)
## df AIC
## model_all 17.74275 603.9831
## model_all_st 18.76908 604.8380There is no sign that time of blood sampling affects the probability of seroconversion.
2.3.4 Should MMF dose or serum levels of immunosuppresants be included?
mmf_total, tacrolimus_level and cyclospirin_level will be included as interaction with (not)taking a given immunosupressant
Open code
findat2 <- findat2 %>% mutate(
cyclospirin_level = if_else(is.na(cyclospirin_level), 0, cyclospirin_level),
tacrolimus_level = if_else(is.na(tacrolimus_level), 0, tacrolimus_level)
)
model_all_imunLev <- gam(seroconversion ~
vacc_time_cont +
vacc_time_2nd_cont +
vaccine_moderna+
tacrolimus +
tacrolimus:tacrolimus_level +
ciclosporinA +
ciclosporinA:cyclospirin_level +
steroids +
months_afterTX +
s(days_after_2nd_dose, k=4) +
age +
sex +
EGFR_vaccination +
mmf_mpa +
mmf_mpa:mmf_total +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes),
data = findat2,
family = binomial, method = "ML")
BIC(model_all, model_all_imunLev)
## df BIC
## model_all 17.74275 680.5494
## model_all_imunLev 21.89357 699.9536
AIC(model_all, model_all_imunLev)
## df AIC
## model_all 17.74275 603.9831
## model_all_imunLev 21.89357 605.4750
summary(model_all_imunLev)
##
## Family: binomial
## Link function: logit
##
## Formula:
## seroconversion ~ vacc_time_cont + vacc_time_2nd_cont + vaccine_moderna +
## tacrolimus + tacrolimus:tacrolimus_level + ciclosporinA +
## ciclosporinA:cyclospirin_level + steroids + months_afterTX +
## s(days_after_2nd_dose, k = 4) + age + sex + EGFR_vaccination +
## mmf_mpa + mmf_mpa:mmf_total + depletationTreat_year + DM +
## BMI + Albumine + log(Lymphocytes)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.7605555 2.2244733 -0.342 0.732423
## vacc_time_cont 0.1202698 0.0934270 1.287 0.197985
## vacc_time_2nd_cont -0.1760747 0.0862591 -2.041 0.041228 *
## vaccine_moderna 0.1487850 0.6254091 0.238 0.811958
## tacrolimus 0.0115972 0.6197150 0.019 0.985069
## ciclosporinA 1.5772051 1.2781541 1.234 0.217214
## steroids -0.2566103 0.3798467 -0.676 0.499318
## months_afterTX 0.0055212 0.0016364 3.374 0.000741 ***
## age -0.0360585 0.0136087 -2.650 0.008057 **
## sex 0.7021383 0.2281872 3.077 0.002091 **
## EGFR_vaccination 0.0336903 0.0060106 5.605 2.08e-08 ***
## mmf_mpa -1.6295212 0.3978698 -4.096 4.21e-05 ***
## depletationTreat_year -1.4617780 0.8462762 -1.727 0.084113 .
## DM -0.1354842 0.2477470 -0.547 0.584471
## BMI -0.0004730 0.0228948 -0.021 0.983517
## Albumine 0.0308921 0.0313723 0.985 0.324774
## log(Lymphocytes) 0.8878499 0.2715069 3.270 0.001075 **
## tacrolimus:tacrolimus_level 0.0653299 0.0704672 0.927 0.353877
## ciclosporinA:cyclospirin_level -0.0041117 0.0104820 -0.392 0.694867
## mmf_mpa:mmf_total -0.0004139 0.0002973 -1.392 0.163924
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(days_after_2nd_dose) 1.557 1.894 3.443 0.11
##
## R-sq.(adj) = 0.267 Deviance explained = 23.7%
## -ML = 282.35 Scale est. = 1 n = 553There is not strong evidence that the predictive accuracy has improved (note that the lower is AIC and BIC, the higher is estimated out-of-sample predictive accuracy of the model)
What about if only MFA doses were added
Open code
model_all_mmfDose <- gam(seroconversion ~
vacc_time_cont +
vacc_time_2nd_cont +
vaccine_moderna+
calcineurin_inhibitor +
steroids +
months_afterTX +
s(days_after_2nd_dose, k=4) +
age +
sex +
EGFR_vaccination +
mmf_mpa +
mmf_mpa:mmf_total +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes),
data = findat2,
family = binomial, method = "ML")
BIC(model_all, model_all_mmfDose)
## df BIC
## model_all 17.74275 680.5494
## model_all_mmfDose 18.63704 684.6035
AIC(model_all, model_all_mmfDose)
## df AIC
## model_all 17.74275 603.9831
## model_all_mmfDose 18.63704 604.1780
summary(model_all_mmfDose)
##
## Family: binomial
## Link function: logit
##
## Formula:
## seroconversion ~ vacc_time_cont + vacc_time_2nd_cont + vaccine_moderna +
## calcineurin_inhibitor + steroids + months_afterTX + s(days_after_2nd_dose,
## k = 4) + age + sex + EGFR_vaccination + mmf_mpa + mmf_mpa:mmf_total +
## depletationTreat_year + DM + BMI + Albumine + log(Lymphocytes)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.0473754 2.2075737 -0.474 0.635182
## vacc_time_cont 0.1199191 0.0929135 1.291 0.196824
## vacc_time_2nd_cont -0.1737961 0.0859166 -2.023 0.043089 *
## vaccine_moderna 0.2339928 0.6189771 0.378 0.705407
## calcineurin_inhibitor 0.6041694 0.4076761 1.482 0.138345
## steroids -0.2685242 0.3766829 -0.713 0.475929
## months_afterTX 0.0061564 0.0015654 3.933 8.40e-05 ***
## age -0.0327195 0.0134015 -2.441 0.014627 *
## sex 0.7092134 0.2270523 3.124 0.001787 **
## EGFR_vaccination 0.0314053 0.0058386 5.379 7.50e-08 ***
## mmf_mpa -1.6362042 0.3937323 -4.156 3.24e-05 ***
## depletationTreat_year -1.3498981 0.8345091 -1.618 0.105750
## DM -0.1186121 0.2467018 -0.481 0.630665
## BMI -0.0018797 0.0227846 -0.082 0.934250
## Albumine 0.0321751 0.0310433 1.036 0.299988
## log(Lymphocytes) 0.9011235 0.2691896 3.348 0.000815 ***
## mmf_mpa:mmf_total -0.0004048 0.0002912 -1.390 0.164417
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(days_after_2nd_dose) 1.369 1.637 2.015 0.18
##
## R-sq.(adj) = 0.262 Deviance explained = 23%
## -ML = 284.49 Scale est. = 1 n = 553There is relatively little evidence that inclusion of the MMF dose (mmf_total) improves accuracy. Moreover, its inclusion does not modify the estimated effect of the timing of vaccine administration.
What if we include only MMF dose (not in interaction with mmf_mpa) as non-linear effect?
Open code
model_all_mmfDose_nl <- gam(seroconversion ~
vacc_time_cont +
vacc_time_2nd_cont +
vaccine_moderna+
calcineurin_inhibitor +
steroids +
months_afterTX +
s(days_after_2nd_dose, k=4) +
age +
sex +
EGFR_vaccination +
s(mmf_total, k=4) +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes),
data = findat2,
family = binomial, method = "ML")
BIC(model_all, model_all_mmfDose_nl)
## df BIC
## model_all 17.74275 680.5494
## model_all_mmfDose_nl 19.49166 687.3627
AIC(model_all, model_all_mmfDose_nl)
## df AIC
## model_all 17.74275 603.9831
## model_all_mmfDose_nl 19.49166 603.2493
summary(model_all_mmfDose_nl)
##
## Family: binomial
## Link function: logit
##
## Formula:
## seroconversion ~ vacc_time_cont + vacc_time_2nd_cont + vaccine_moderna +
## calcineurin_inhibitor + steroids + months_afterTX + s(days_after_2nd_dose,
## k = 4) + age + sex + EGFR_vaccination + s(mmf_total, k = 4) +
## depletationTreat_year + DM + BMI + Albumine + log(Lymphocytes)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -2.638786 2.227641 -1.185 0.236189
## vacc_time_cont 0.117320 0.093063 1.261 0.207434
## vacc_time_2nd_cont -0.168696 0.085803 -1.966 0.049288 *
## vaccine_moderna 0.191649 0.617920 0.310 0.756446
## calcineurin_inhibitor 0.625665 0.405119 1.544 0.122492
## steroids -0.271591 0.377018 -0.720 0.471300
## months_afterTX 0.006229 0.001563 3.986 6.72e-05 ***
## age -0.033435 0.013429 -2.490 0.012786 *
## sex 0.718278 0.226804 3.167 0.001540 **
## EGFR_vaccination 0.031508 0.005854 5.382 7.35e-08 ***
## depletationTreat_year -1.350325 0.833988 -1.619 0.105422
## DM -0.125602 0.246217 -0.510 0.609962
## BMI -0.002020 0.022748 -0.089 0.929234
## Albumine 0.032273 0.030991 1.041 0.297703
## log(Lymphocytes) 0.887962 0.268128 3.312 0.000927 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(days_after_2nd_dose) 1.337 1.589 2.075 0.166
## s(mmf_total) 2.640 2.902 58.528 <2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## R-sq.(adj) = 0.262 Deviance explained = 23.4%
## -ML = 287.51 Scale est. = 1 n = 553Again, there is not substantial improvement of predictive accuracy when mmf_mpa is replaced with mmf_total. Similarly, estimated effects of the timing of vaccine administration are also very similar.
2.3.5 Does mTOR modulate the effect of vaccine timing?
mTOR may affect circadian rhythms. Does intake of mTOR interact with the timings of vaccine administration?
Open code
model_all_mTOR <- gam(seroconversion ~
vacc_time_cont*mTOR +
vacc_time_2nd_cont*mTOR +
vaccine_moderna+
calcineurin_inhibitor +
steroids +
months_afterTX +
s(days_after_2nd_dose, k=4) +
age +
sex +
EGFR_vaccination +
mmf_mpa +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes),
data = findat2,
family = binomial, method = "ML")
BIC(model_all, model_all_mTOR)
## df BIC
## model_all 17.74275 680.5494
## model_all_mTOR 20.75451 697.9639
AIC(model_all, model_all_mTOR)
## df AIC
## model_all 17.74275 603.9831
## model_all_mTOR 20.75451 608.4008
summary(model_all_mTOR)
##
## Family: binomial
## Link function: logit
##
## Formula:
## seroconversion ~ vacc_time_cont * mTOR + vacc_time_2nd_cont *
## mTOR + vaccine_moderna + calcineurin_inhibitor + steroids +
## months_afterTX + s(days_after_2nd_dose, k = 4) + age + sex +
## EGFR_vaccination + mmf_mpa + depletationTreat_year + DM +
## BMI + Albumine + log(Lymphocytes)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.408433 2.265082 -0.622 0.534072
## vacc_time_cont 0.123155 0.095175 1.294 0.195672
## mTOR 0.201141 3.836161 0.052 0.958184
## vacc_time_2nd_cont -0.186546 0.088804 -2.101 0.035671 *
## vaccine_moderna 0.308443 0.621063 0.497 0.619445
## calcineurin_inhibitor 1.061501 0.603503 1.759 0.078595 .
## steroids -0.283732 0.374536 -0.758 0.448717
## months_afterTX 0.006012 0.001577 3.812 0.000138 ***
## age -0.029860 0.013357 -2.236 0.025379 *
## sex 0.692064 0.225915 3.063 0.002188 **
## EGFR_vaccination 0.030460 0.005871 5.188 2.12e-07 ***
## mmf_mpa -1.992365 0.281160 -7.086 1.38e-12 ***
## depletationTreat_year -1.248403 0.837171 -1.491 0.135905
## DM -0.099049 0.245650 -0.403 0.686791
## BMI -0.003097 0.022764 -0.136 0.891799
## Albumine 0.029070 0.031020 0.937 0.348694
## log(Lymphocytes) 0.919315 0.270141 3.403 0.000666 ***
## vacc_time_cont:mTOR -0.131061 0.363300 -0.361 0.718286
## mTOR:vacc_time_2nd_cont 0.174949 0.333567 0.524 0.599946
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(days_after_2nd_dose) 1.451 1.755 2.223 0.179
##
## R-sq.(adj) = 0.258 Deviance explained = 23%
## -ML = 284.69 Scale est. = 1 n = 553The inclusion of the mTOR factor and its interactions with vaccine administration time did not improve predictive accuracy of the model. However, the effect of vaccine timing for the 2nd dose was practically diminished with taking mTOR.
What if we worked only with mTOR non-takers?
Open code
findat2_no_mTOR <- findat2 %>% filter(mTOR == 0)
dim(findat2_no_mTOR)
## [1] 510 35
model_all_no_mTOR <- gam(seroconversion ~
vacc_time_cont +
vacc_time_2nd_cont +
vaccine_moderna+
calcineurin_inhibitor +
steroids +
months_afterTX +
s(days_after_2nd_dose, k=4) +
age +
sex +
EGFR_vaccination +
mmf_mpa +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes),
data = findat2_no_mTOR,
family = binomial)
summary(model_all)
##
## Family: binomial
## Link function: logit
##
## Formula:
## seroconversion ~ vacc_time_cont + vacc_time_2nd_cont + vaccine_moderna +
## calcineurin_inhibitor + steroids + months_afterTX + s(days_after_2nd_dose,
## k = 4) + age + sex + EGFR_vaccination + mmf_mpa + depletationTreat_year +
## DM + BMI + Albumine + log(Lymphocytes)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -0.851378 2.195512 -0.388 0.698178
## vacc_time_cont 0.114933 0.092752 1.239 0.215291
## vacc_time_2nd_cont -0.180905 0.085754 -2.110 0.034895 *
## vaccine_moderna 0.329081 0.613955 0.536 0.591957
## calcineurin_inhibitor 0.568021 0.405410 1.401 0.161184
## steroids -0.298206 0.375072 -0.795 0.426577
## months_afterTX 0.006058 0.001566 3.869 0.000109 ***
## age -0.030845 0.013248 -2.328 0.019894 *
## sex 0.679311 0.225520 3.012 0.002594 **
## EGFR_vaccination 0.030996 0.005821 5.325 1.01e-07 ***
## mmf_mpa -2.037649 0.273102 -7.461 8.58e-14 ***
## depletationTreat_year -1.263212 0.836284 -1.511 0.130915
## DM -0.087813 0.245060 -0.358 0.720094
## BMI -0.003811 0.022690 -0.168 0.866601
## Albumine 0.031307 0.030975 1.011 0.312140
## log(Lymphocytes) 0.904000 0.268759 3.364 0.000769 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(days_after_2nd_dose) 1.443 1.743 2.208 0.179
##
## R-sq.(adj) = 0.259 Deviance explained = 22.8%
## -ML = 285.47 Scale est. = 1 n = 553
summary(model_all_no_mTOR)
##
## Family: binomial
## Link function: logit
##
## Formula:
## seroconversion ~ vacc_time_cont + vacc_time_2nd_cont + vaccine_moderna +
## calcineurin_inhibitor + steroids + months_afterTX + s(days_after_2nd_dose,
## k = 4) + age + sex + EGFR_vaccination + mmf_mpa + depletationTreat_year +
## DM + BMI + Albumine + log(Lymphocytes)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.778196 2.424478 -0.733 0.463293
## vacc_time_cont 0.107339 0.097110 1.105 0.269016
## vacc_time_2nd_cont -0.172046 0.090915 -1.892 0.058441 .
## vaccine_moderna 0.460370 0.668486 0.689 0.491027
## calcineurin_inhibitor 1.671210 0.848267 1.970 0.048822 *
## steroids -0.065391 0.396339 -0.165 0.868953
## months_afterTX 0.006337 0.001651 3.838 0.000124 ***
## age -0.034185 0.014065 -2.430 0.015080 *
## sex 0.845308 0.238973 3.537 0.000404 ***
## EGFR_vaccination 0.028270 0.006110 4.627 3.71e-06 ***
## mmf_mpa -2.033396 0.301103 -6.753 1.45e-11 ***
## depletationTreat_year -1.199238 0.850261 -1.410 0.158411
## DM -0.116193 0.259392 -0.448 0.654194
## BMI -0.008454 0.024776 -0.341 0.732935
## Albumine 0.028980 0.032327 0.896 0.370000
## log(Lymphocytes) 0.947763 0.284674 3.329 0.000871 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(days_after_2nd_dose) 2.358 2.717 4.72 0.107
##
## R-sq.(adj) = 0.253 Deviance explained = 22.8%
## UBRE = 0.087013 Scale est. = 1 n = 510The estimated effects did not change substantially
2.3.6 Does definition of afternoon affect model with dichotomized time?
Open code
findat2 <- findat2 %>% mutate(
afternoon1 = if_else(vacc_time_cont > 13, 1, 0),
afternoon2 = if_else(vacc_time_2nd_cont > 13, 1, 0))
model_13 <- gam(seroconversion ~
afternoon1 +
afternoon2 +
vaccine_moderna +
calcineurin_inhibitor +
steroids+
s(days_after_2nd_dose, k=4) +
age +
sex +
months_afterTX +
EGFR_vaccination +
mmf_mpa +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes),
data = findat2,
method='ML',
family = binomial)
findat2 <- findat2 %>% mutate(
afternoon1 = if_else(vacc_time_cont > 12, 1, 0),
afternoon2 = if_else(vacc_time_2nd_cont > 12, 1, 0))
model_12 <- gam(seroconversion ~
afternoon1 +
afternoon2 +
vaccine_moderna +
calcineurin_inhibitor +
steroids+
s(days_after_2nd_dose, k=4) +
age +
sex +
months_afterTX +
EGFR_vaccination +
mmf_mpa +
depletationTreat_year +
DM +
BMI +
Albumine +
log(Lymphocytes),
data = findat2,
method='ML',
family = binomial)
summary(model_13)
##
## Family: binomial
## Link function: logit
##
## Formula:
## seroconversion ~ afternoon1 + afternoon2 + vaccine_moderna +
## calcineurin_inhibitor + steroids + s(days_after_2nd_dose,
## k = 4) + age + sex + months_afterTX + EGFR_vaccination +
## mmf_mpa + depletationTreat_year + DM + BMI + Albumine + log(Lymphocytes)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.153823 1.864648 -0.619 0.53606
## afternoon1 0.160436 0.258221 0.621 0.53439
## afternoon2 -0.582066 0.252555 -2.305 0.02118 *
## vaccine_moderna 0.308227 0.610848 0.505 0.61385
## calcineurin_inhibitor 0.559768 0.410081 1.365 0.17225
## steroids -0.325726 0.377071 -0.864 0.38768
## age -0.035886 0.012703 -2.825 0.00473 **
## sex 0.643139 0.224916 2.859 0.00424 **
## months_afterTX 0.006154 0.001578 3.899 9.66e-05 ***
## EGFR_vaccination 0.031480 0.005784 5.442 5.26e-08 ***
## mmf_mpa -2.054721 0.273396 -7.516 5.67e-14 ***
## depletationTreat_year -1.204303 0.825471 -1.459 0.14458
## DM -0.069458 0.244344 -0.284 0.77621
## BMI -0.004291 0.022864 -0.188 0.85112
## Albumine 0.031941 0.031145 1.026 0.30511
## log(Lymphocytes) 0.872681 0.269093 3.243 0.00118 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(days_after_2nd_dose) 1.639 1.995 4.184 0.111
##
## R-sq.(adj) = 0.261 Deviance explained = 23.1%
## -ML = 284.83 Scale est. = 1 n = 553
summary(model_12)
##
## Family: binomial
## Link function: logit
##
## Formula:
## seroconversion ~ afternoon1 + afternoon2 + vaccine_moderna +
## calcineurin_inhibitor + steroids + s(days_after_2nd_dose,
## k = 4) + age + sex + months_afterTX + EGFR_vaccination +
## mmf_mpa + depletationTreat_year + DM + BMI + Albumine + log(Lymphocytes)
##
## Parametric coefficients:
## Estimate Std. Error z value Pr(>|z|)
## (Intercept) -1.224324 1.863816 -0.657 0.511251
## afternoon1 0.110185 0.236891 0.465 0.641839
## afternoon2 -0.235233 0.237352 -0.991 0.321651
## vaccine_moderna 0.202526 0.608818 0.333 0.739396
## calcineurin_inhibitor 0.536479 0.404042 1.328 0.184251
## steroids -0.323205 0.375210 -0.861 0.389018
## age -0.035326 0.013040 -2.709 0.006748 **
## sex 0.664483 0.224250 2.963 0.003045 **
## months_afterTX 0.005900 0.001555 3.794 0.000148 ***
## EGFR_vaccination 0.030550 0.005814 5.255 1.48e-07 ***
## mmf_mpa -2.023115 0.274014 -7.383 1.54e-13 ***
## depletationTreat_year -1.205315 0.831979 -1.449 0.147412
## DM -0.072789 0.243921 -0.298 0.765389
## BMI -0.005207 0.022572 -0.231 0.817564
## Albumine 0.033504 0.030989 1.081 0.279620
## log(Lymphocytes) 0.876223 0.268556 3.263 0.001103 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Approximate significance of smooth terms:
## edf Ref.df Chi.sq p-value
## s(days_after_2nd_dose) 1.578 1.921 2.888 0.164
##
## R-sq.(adj) = 0.253 Deviance explained = 22.4%
## -ML = 287.16 Scale est. = 1 n = 553The threshold for the “afternoon” definition has a relatively large effect in terms of estimated effect size and the clarity of the vaccination timing effect.
In summary, final models will be fitted without interaction between the timings of two vaccination doses. Vaccination times will be included as either continuous variable with linear effect or in dichotomized form. Interaction between timings of two doses will not be included to the final models.
2.4 Main (Bayesian) model
As stated above, we primarily use model with (weakly) informative priors for the effect of vaccination timing. The informative priors use previously acquired information, e.g. results from previous studies (Lin and Hung 2023), and use it for more informed inference.
As 2 SD of continuous vaccination time is only 1.25 times larger than difference between the mean times of vaccination in the morning vs. afternoon, and we expect similar distribution also in the previous study (Lin and Hung 2023), we will adopt our prior distribution directly on the basis of the results of Lin et al. (Lin and Hung 2023). The prior will be relatively very wide, accounting for limited comparability of study designs and patient population between our study vs. Lin et al. (Lin and Hung 2023).
Open code
# 1st dose
sd(findat2$vacc_time_cont)*2
## [1] 2.803969
mean(findat2[findat2$afternoon1==0,]$vacc_time_cont)-
mean(findat2[findat2$afternoon1==1,]$vacc_time_cont)
## [1] -2.270227
# 2nd dose
sd(findat2$vacc_time_2nd_cont)*2
## [1] 2.952956
mean(findat2[findat2$afternoon2==0,]$vacc_time_2nd_cont)-
mean(findat2[findat2$afternoon2==1,]$vacc_time_2nd_cont)
## [1] -2.464559Setting priors
Open code
prior1 <- c(prior(normal(0, 5), class = b),
prior(normal(-0.9, 2), class = b, coef= vacc_time_scal),
prior(normal(-0.9, 2), class = b, coef= vacc_time2_scal))
prior2 <- c(prior(normal(0, 5), class = b),
prior(normal(-0.9, 2), class = b, coef= afternoon1),
prior(normal(-0.9, 2), class = b, coef= afternoon2))2.4.1 Fitting main model
Open code
# scaling the continuous time by 2SD and centering to have zero means
findat2 <- findat2 %>% mutate(vacc_time_scal = rescale(vacc_time_cont),
vacc_time2_scal = rescale(vacc_time_2nd))
model_main <- run_model(brm(seroconversion ~
vacc_time_scal +
vacc_time2_scal +
rescale(vaccine_moderna, binary.inputs = '-0.5,0.5') +
rescale(calcineurin_inhibitor, binary.inputs = '-0.5,0.5') +
rescale(steroids, binary.inputs = '-0.5,0.5') +
rescale(months_afterTX) +
s(rescale(days_after_2nd_dose), k=4) +
rescale(DM, binary.inputs = '-0.5,0.5') +
rescale(age) +
rescale(BMI) +
rescale(Albumine) +
rescale(log(Lymphocytes))+
rescale(sex, binary.inputs = '-0.5,0.5') +
rescale(EGFR_vaccination) +
rescale(mmf_mpa, binary.inputs = '-0.5,0.5') +
rescale(depletationTreat_year, binary.inputs = '-0.5,0.5'),
data = findat2,
prior = prior1,
chains = 3, iter = 6000, warmup = 2000, seed = 17,
control = list(adapt_delta = 0.99),
cores = 1,
family = bernoulli(link='logit')),
'gitignore/brm3/model_main', reuse = TRUE)
summary(model_main)
## Family: bernoulli
## Links: mu = logit
## Formula: seroconversion ~ vacc_time_scal + vacc_time2_scal + rescale(vaccine_moderna, binary.inputs = "-0.5,0.5") + rescale(calcineurin_inhibitor, binary.inputs = "-0.5,0.5") + rescale(steroids, binary.inputs = "-0.5,0.5") + rescale(months_afterTX) + s(rescale(days_after_2nd_dose), k = 4) + rescale(DM, binary.inputs = "-0.5,0.5") + rescale(age) + rescale(BMI) + rescale(Albumine) + rescale(log(Lymphocytes)) + rescale(sex, binary.inputs = "-0.5,0.5") + rescale(EGFR_vaccination) + rescale(mmf_mpa, binary.inputs = "-0.5,0.5") + rescale(depletationTreat_year, binary.inputs = "-0.5,0.5")
## Data: findat2 (Number of observations: 553)
## Draws: 3 chains, each with iter = 6000; warmup = 2000; thin = 1;
## total post-warmup draws = 12000
##
## Smooth Terms:
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sds(srescaledays_after_2nd_dose_1) 2.17 1.65 0.13 6.29 1.00
## Bulk_ESS Tail_ESS
## sds(srescaledays_after_2nd_dose_1) 4110 3781
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI
## Intercept -0.82 0.60 -2.04
## vacc_time_scal 0.29 0.26 -0.22
## vacc_time2_scal -0.50 0.25 -1.00
## rescalevaccine_modernabinary.inputsEQM0.50.5 0.33 0.63 -0.92
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 0.58 0.42 -0.24
## rescalesteroidsbinary.inputsEQM0.50.5 -0.32 0.38 -1.06
## rescalemonths_afterTX 1.03 0.26 0.52
## rescaleDMbinary.inputsEQM0.50.5 -0.09 0.25 -0.59
## rescaleage -0.66 0.27 -1.20
## rescaleBMI -0.06 0.23 -0.51
## rescaleAlbumine 0.24 0.23 -0.21
## rescalelogLymphocytes 0.84 0.25 0.36
## rescalesexbinary.inputsEQM0.50.5 0.70 0.23 0.26
## rescaleEGFR_vaccination 1.27 0.23 0.82
## rescalemmf_mpabinary.inputsEQM0.50.5 -2.12 0.28 -2.67
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 -1.44 0.90 -3.36
## srescaledays_after_2nd_dose_1 0.04 2.16 -4.77
## u-95% CI Rhat Bulk_ESS
## Intercept 0.29 1.00 13053
## vacc_time_scal 0.80 1.00 10219
## vacc_time2_scal -0.01 1.00 11193
## rescalevaccine_modernabinary.inputsEQM0.50.5 1.57 1.00 13459
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 1.41 1.00 15062
## rescalesteroidsbinary.inputsEQM0.50.5 0.43 1.00 13466
## rescalemonths_afterTX 1.55 1.00 12351
## rescaleDMbinary.inputsEQM0.50.5 0.40 1.00 13042
## rescaleage -0.13 1.00 11398
## rescaleBMI 0.39 1.00 13810
## rescaleAlbumine 0.71 1.00 13204
## rescalelogLymphocytes 1.33 1.00 16104
## rescalesexbinary.inputsEQM0.50.5 1.16 1.00 15488
## rescaleEGFR_vaccination 1.74 1.00 12916
## rescalemmf_mpabinary.inputsEQM0.50.5 -1.58 1.00 12141
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 0.13 1.00 14230
## srescaledays_after_2nd_dose_1 3.71 1.00 6416
## Tail_ESS
## Intercept 7833
## vacc_time_scal 8714
## vacc_time2_scal 9604
## rescalevaccine_modernabinary.inputsEQM0.50.5 8037
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 9066
## rescalesteroidsbinary.inputsEQM0.50.5 8782
## rescalemonths_afterTX 9368
## rescaleDMbinary.inputsEQM0.50.5 8410
## rescaleage 8804
## rescaleBMI 8189
## rescaleAlbumine 9044
## rescalelogLymphocytes 9864
## rescalesexbinary.inputsEQM0.50.5 9159
## rescaleEGFR_vaccination 8310
## rescalemmf_mpabinary.inputsEQM0.50.5 8517
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 7578
## srescaledays_after_2nd_dose_1 8192
##
## 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).
repor(model_main)| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| vacc_time_scal | 1.339 | 0.805 | 2.222 | 0.8661 |
| vacc_time2_scal | 0.605 | 0.368 | 0.993 | 0.9772 |
| vaccine_moderna | 1.386 | 0.398 | 4.802 | 0.7033 |
| calcineurin_inhibitor | 1.789 | 0.786 | 4.111 | 0.9225 |
| steroids | 0.728 | 0.346 | 1.541 | 0.7910 |
| months_afterTX | 2.791 | 1.685 | 4.704 | 1.0000 |
| DM | 0.911 | 0.556 | 1.498 | 0.6451 |
| age | 0.517 | 0.301 | 0.874 | 0.9925 |
| BMI | 0.939 | 0.603 | 1.478 | 0.6114 |
| Albumine | 1.276 | 0.814 | 2.034 | 0.8549 |
| logLymphocytes | 2.320 | 1.440 | 3.786 | 0.9998 |
| sex | 2.015 | 1.291 | 3.201 | 0.9996 |
| EGFR_vaccination | 3.576 | 2.279 | 5.704 | 1.0000 |
| mmf_mpa | 0.120 | 0.069 | 0.205 | 1.0000 |
| depletationTreat_year | 0.237 | 0.035 | 1.137 | 0.9608 |
2.4.2 Reporting the main model
Preparation of posterior draws
Open code
sd1 <- sd(findat2$vacc_time_cont)
sd2 <- sd(findat2$vacc_time_2nd_cont)
mean1 <- mean(findat2$vacc_time_cont)
mean2 <- mean(findat2$vacc_time_2nd_cont)
time_range <- c(7.1, 17.4)
time <-seq(
time_range[1],
time_range[2], length=100)
time_z1 <- seq(
sca(time_range, mean1, sd1)[1],
sca(time_range, mean1, sd1)[2], length=100)
time_z2 <- seq(
sca(time_range, mean2, sd2)[1],
sca(time_range, mean2, sd2)[2], length=100)
post_fix <- as.data.frame(model_main) %>% select(
matches('b_')) Extraction of estimates
Open code
# estimation of seroconversion log-odds in relation to 1st dose time
vacc_1est <- data.frame()
for(i in 1:length(time_z1)){
vacc_1est[1:dim(post_fix)[1], i] <- post_fix$b_Intercept + post_fix$b_vacc_time_scal*time_z1[i]}
vacc1_draw <- t(rbind(time, inv.logit(vacc_1est)))
vacc_1est <- sapply(vacc_1est, function(p) quantile(p, probs = c(1/40, 39/40, 0.5, 1/5, 4/5)))
# estimation of seroconversion log-odds in relation to 2nd dose time
vacc_2est <- data.frame()
for(i in 1:length(time_z2)){
vacc_2est[1:dim(post_fix)[1],i] <- post_fix$b_Intercept + post_fix$b_vacc_time2_scal*time_z2[i]}
vacc2_draw <- t(rbind(time, inv.logit(vacc_2est)))
vacc_2est <- sapply(vacc_2est, function(p) quantile(p, probs = c(1/40, 39/40, 0.5, 1/5, 4/5)))
predict <- data.frame(
prediction =
unlist(c(
inv.logit(vacc_1est[3,]), inv.logit(vacc_2est[3,])
)),
cil =
unlist(c(
inv.logit(vacc_1est[1,]), inv.logit(vacc_2est[1,])
)),
ciu =
unlist(c(
inv.logit(vacc_1est[2,]), inv.logit(vacc_2est[2,])
)),
cilw =
unlist(c(
inv.logit(vacc_1est[4,]), inv.logit(vacc_2est[4,])
)),
ciuw =
unlist(c(
inv.logit(vacc_1est[5,]), inv.logit(vacc_2est[5,])
)),
time = rep(time, 2), group = c(
rep("1st dose",100), rep("2nd dose",100)
))2.4.2.1 Fig1a: Plot for estimates
Open code
cole <- c('#CD7006', '#0028F0')
range <- c(0, 1)
xpo <- 14
ypo <- 0.8
pd1 <- round(p_dir(post_fix$b_vacc_time_scal, '<', 0),2)
pd2 <- round(p_dir(post_fix$b_vacc_time2_scal, '<', 0), 3)
findat2_tra <- data.frame(time=c(findat2$vacc_time_cont,
findat2$vacc_time_2nd_cont),
group=as.factor(c(rep('1st dose', dim(findat2)[1]),
rep('2nd dose', dim(findat2)[1]))),
seroconversion = c(findat2$seroconversion, findat2$seroconversion))
fig1a <- predict %>%
ggplot(aes(x = time, y = prediction, by = group, col=group, fill=group)) +
geom_line(aes(y=prediction), linewidth = 1) +
scale_x_continuous(breaks=c(seq(7, 18, by = 2))) +
scale_y_continuous(limits = c(range[1], range[2]),
breaks = c(seq(0, 1, by=.2))) +
geom_ribbon(aes(ymin = cil, ymax = ciu),
alpha = 0.35, colour = NA) +
labs(x = "Time of vaccination", y = 'Probability of seroconversion') +
scale_color_manual(values=cole,
name="Dose",
breaks=c('1st dose', '2nd dose'),
labels=c('1st', '2nd')) +
scale_fill_manual(values=cole,
name="Dose",
breaks=c('1st dose', '2nd dose'),
labels=c('1st', '2nd')) +
facet_grid(~group) +
# stat_dots(data = findat2_tra,
# aes(y = seroconversion,
# side = ifelse(seroconversion == 0, "top", "bottom")),
# scale = 1/8) +
theme(axis.text=element_text(size=10),
axis.title=element_text(size=12),
strip.text.x = element_text(size = 12),
legend.position = "none")
fig1a
2.4.2.2 Spaghetti cooking (Fig 1a alternative)
Data
Open code
vacc_draw <- rbind(vacc1_draw, vacc2_draw)
vacc_draw <- data.frame(time = vacc_draw[,1],
group = as.factor(c(rep('1st dose', dim(vacc1_draw)[1]),
rep('2nd dose', dim(vacc1_draw)[1]))),
vacc_draw[,-1])
set.seed(10)
ndraws = 100
vacc_draw_long <- vacc_draw %>%
pivot_longer(
cols = c(sample(3:12000, ndraws)),
names_to = "variable",
values_to = "value"
) %>% select(time, group, variable, value)Print cooked spaghetti
Open code
cole <- c('#CD7006', '#0028F0')
range <- c(0, 1)
xpo <- 14
ypo <- 0.8
pd1 <- round(p_dir(post_fix$b_vacc_time_scal, '<', 0),2)
pd2 <- round(p_dir(post_fix$b_vacc_time2_scal, '<', 0), 3)
findat2_tra <- data.frame(time=c(findat2$vacc_time_cont,
findat2$vacc_time_2nd_cont),
group=as.factor(c(rep('1st dose', dim(findat2)[1]),
rep('2nd dose', dim(findat2)[1]))),
seroconversion = c(findat2$seroconversion, findat2$seroconversion))
fig1as <- predict %>%
ggplot(aes(x = time, y = prediction, by = group, col=group, fill=group)) +
geom_line(aes(y=prediction), linewidth = 1) +
geom_line(data = vacc_draw_long, aes(x = time, y=value, by=variable), linewidth = 0.2, alpha=0.3) +
scale_x_continuous(breaks=c(seq(7, 18, by = 2))) +
scale_y_continuous(limits = c(range[1], range[2]),
breaks = c(seq(0, 1, by=.2))) +
labs(x = "Time of vaccination", y = 'Probability of seroconversion)') +
scale_color_manual(values=cole,
name="Dose",
breaks=c('1st dose', '2nd dose'),
labels=c('1st', '2nd')) +
scale_fill_manual(values=cole,
name="Dose",
breaks=c('1st dose', '2nd dose'),
labels=c('1st', '2nd')) +
facet_grid(rows = vars(group)) +
stat_dots(data = findat2_tra,
aes(y = seroconversion, side = ifelse(seroconversion == 0, "top", "bottom")),
scale = 1/8) +
theme(axis.text=element_text(size=10),
axis.title=element_text(size=12),
strip.text.y = element_text(size = 12)) +
theme(legend.position = "none")
fig1as
What is the probability of seroconversion if the 2nd dose was administered at different times?
Open code
predict %>% filter(group == '2nd dose') %>% select(time, prediction) %>%
filter(row_number() %% 10 == 0)
## time prediction
## 1 8.036364 0.4757788
## 2 9.076768 0.4328700
## 3 10.117172 0.3910524
## 4 11.157576 0.3502532
## 5 12.197980 0.3106888
## 6 13.238384 0.2737142
## 7 14.278788 0.2396603
## 8 15.319192 0.2082087
## 9 16.359596 0.1800016
## 10 17.400000 0.15483712.4.2.3 Fig1b: Posterior of vaccination time effects
Open code
sd1 <- sd(findat2$vacc_time_cont)*2
sd2 <- sd(findat2$vacc_time_2nd_cont)*2
pd <- c(round(p_dir(post_fix$b_vacc_time_scal, '<', 0),2),
round(p_dir(post_fix$b_vacc_time2_scal, '<', 0), 3)
)
tr <- as_tibble(model_main) %>%
mutate(dose_1 = exp(b_vacc_time_scal/sd1),
dose_2 = exp(b_vacc_time2_scal/sd2)) %>%
select(dose_1, dose_2) %>% data.frame()
CIS <- sapply(tr, function(p) quantile(p, probs = c(1/40, 39/40, 1/2, 0.055, 0.945)))
CIS
## dose_1 dose_2
## 2.5% 0.9254916 0.7128711
## 97.5% 1.3293382 0.9975222
## 50% 1.1091260 0.8431931
## 5.5% 0.9561369 0.7355519
## 94.5% 1.2862285 0.9689681
CIS[1, ] <- round(CIS[1, ], 2)
CIS[3, ] <- round(CIS[3, ], 2)
CIS[2, 2] <- round(CIS[2, 2], 3)
CIS[2, 1] <- round(CIS[2, 1], 2)
fig1b <- as_tibble(model_main) %>%
mutate(dose_1 = exp(b_vacc_time_scal/sd1),
dose_2 = exp(b_vacc_time2_scal/sd2)) %>%
select(dose_2, dose_1) %>%
gather(key, tau) %>%
mutate(key = factor(key, levels = c("dose_2", "dose_1"))) %>%
ggplot(aes(x = tau, y = key, fill = key)) +
stat_halfeye(.width = c(0.95), slab_alpha=0.5,
linewidth = 5,
shape = 18,
point_size=5,
normalize = "groups") +
labs(x = "Odds ratio (1-hour change)",
y = 'Dose of anti-Covid vaccination') +
scale_fill_manual(values = cole,
name="Dose",
breaks=c('dose_1', 'dose_2'),
labels=c('1st', '2nd')) +
scale_y_discrete(expand = expansion(add = 0.1),
breaks=c('dose_2', 'dose_1'),
labels=c('2nd', '1st')) +
geom_vline(xintercept=1, linetype=2,
color = "red", size=0.6) +
theme(axis.text=element_text(size=12),
axis.title=element_text(size=12)) +
theme(legend.position = "none") +
ggtitle('Posterior distributions, main model') +
annotate("text", x=1.4, y=1.9 ,
label = paste0("OR: ", CIS[3,2]),
color = cole[2] ) +
annotate("text", x=1.4, y=2.9 ,
label = paste0("OR: ", CIS[3,1]),
color = cole[1] ) +
annotate("text", x=1.4, y=1.7 ,
label = paste0("95% CI: [", CIS[1,2], ", ", CIS[2,2], "]"),
color = cole[2] ) +
annotate("text", x=1.4, y=2.7 ,
label = paste0("95% CI: [", CIS[1,1], ", ", CIS[2,1], "]"),
color = cole[1] ) +
annotate("text", x=1.4, y=1.5 , label = paste0("P [OR<1]: ", pd[2]),
color = cole[2] ) +
annotate("text", x=1.4, y=2.5 , label = paste0("P [OR<1]: ", pd[1]),
color = cole[1] )
## Warning: Using `size` aesthetic for lines was deprecated in ggplot2 3.4.0.
## ℹ Please use `linewidth` instead.
fig1b
2.4.2.4 Fig2: Posterior of other covariates
Open code
sds2 <- c(
sd(findat2$vacc_time_cont)*2,
sd(findat2$vacc_time_2nd_cont)*2,
1, # moderna
1, # caclcineurin inhibitors
1, # steroids
sd(findat2$months_afterTX)*2,
1, # DM
(sd(findat2$age)*2),
(sd(findat2$BMI)*2),
sd(findat2$Albumine)*2,
sd(log(findat2$Lymphocytes))*2,
1, # male sex
sd(findat2$EGFR_vaccination)*2,
1, ## MMF/MPA
1 # deplating
)
fig2 <- as_tibble(model_main) %>%
mutate('1st vaccination time [hour]' = exp(b_vacc_time_scal/sds2[1]),
'2nd vaccination time [hour]'= exp(b_vacc_time2_scal/sds2[2]),
'Moderna [0/1]' = exp(b_rescalevaccine_modernabinary.inputsEQM0.50.5),
'Calcineurin inhibitor [0/1]' = exp(b_rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5),
'Steroids [0/1]' = exp(b_rescalesteroidsbinary.inputsEQM0.50.5),
'Time after TX [10 years]' = exp((120*b_rescalemonths_afterTX)/sds2[6]),
'Diabetes Mellitus [0/1]' = exp(b_rescaleDMbinary.inputsEQM0.50.5),
'Age [10 years]' = exp((b_rescaleage*10)/sds2[8]),
'BMI [5 units]'= exp((b_rescaleBMI *5)/sds2[9]),
'Albumine [10 g/L]' = exp((b_rescaleAlbumine*10)/sds2[10]),
'Lymphocyte counts [log(10^9/L)]' = exp(b_rescalelogLymphocytes/sds2[11]),
'Male sex [0/1]' = exp(b_rescalesexbinary.inputsEQM0.50.5),
'eGFR [10 mL/min/1.73 m2]' = exp((b_rescaleEGFR_vaccination*10)/sds2[13]),
'MMF/MPA [0/1]' = exp(b_rescalemmf_mpabinary.inputsEQM0.50.5),
'Depleting therapy [0/1]' = exp(b_rescaledepletationTreat_yearbinary.inputsEQM0.50.5) ) %>%
select( 'Moderna [0/1]',
'Calcineurin inhibitor [0/1]',
'Steroids [0/1]',
'Time after TX [10 years]',
'Diabetes Mellitus [0/1]',
'Age [10 years]',
'BMI [5 units]',
'Albumine [10 g/L]',
'Lymphocyte counts [log(10^9/L)]',
'Male sex [0/1]',
'eGFR [10 mL/min/1.73 m2]',
'MMF/MPA [0/1]',
'Depleting therapy [0/1]'
) %>%
gather(key, tau) %>%
ggplot(aes(x = tau, y = key)) +
stat_halfeye(.width = c(0.95), slab_alpha=0.8,
linewidth = 5,
shape = 18,
point_size=5,
normalize = "groups",
fill = 'grey60') +
labs(x = "Odds ratio (log-scaled axis)", y = NULL) +
#scale_y_discrete(expand = expansion(add = 0.1)) +
scale_x_continuous(trans='log2', limits=c(0.018, 30),
breaks=c(1/64, 1/16, 1/4, 1, 4, 16),
labels=c('', '1/16', '1/4', 1, 4, 16)) +
geom_vline(xintercept=1, linetype=2,
color = "red", size=0.6) +
theme(axis.text=element_text(size=12),
axis.title=element_text(size=12)) +
theme(legend.position = "none")
fig2
2.4.2.5 Table 2: main model results
Open code
sds1 <- c(
sd(findat2$vacc_time_cont)*2,
sd(findat2$vacc_time_2nd_cont)*2,
1, # moderna
1, # caclcineurin inhibitors
1, # steroids
(sd(findat2$months_afterTX)*2)/120,
1, # DM
(sd(findat2$age)*2)/10,
(sd(findat2$BMI)*2)/5,
sd(findat2$Albumine)*2,
sd(log(findat2$Lymphocytes))*2,
1, # male sex
(sd(findat2$EGFR_vaccination)*2)/10,
1, ## MMF/MPA
1 # deplating
)
labe1 <- c('1st vaccination time [hour]',
'2nd vaccination time [hour]',
'Moderna [0/1]',
'Calcineurin inhibitor[0/1]',
'Steroids [0/1]',
'Time after TX [10 years]',
'Diabetes Mellitus [0/1]',
'Age [10 years]',
'BMI [5 units]',
'Albumine [g/L]',
'Lymphocyte counts [log(10^9/L)]',
'Male sex [0/1]',
'eGFR [10 mL/min/1.73 m2]',
'MMF/MPA [0/1]',
'Depleting therapy [0/1]')
tr1 <- repor(model_main, labels=labe1, scals = sds1, nhtm=TRUE)
# full model
table2 <- kableExtra::kable(tr1, caption =
"Table 2. Results of Bayesian logistic regression modeling the probability of seroconversion following mRNA anti-Covid vaccination. 'OR' represents the odds ratio, indicating the expected fold-change in odds when the predictor changes by the unit defined in the square brackets. 'Q2.5' and 'Q97.5' denote the lower and upper bounds of the 95% credible interval, while 'PD' shows the probability of direction")
table2| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| 1st vaccination time [hour] | 1.110 | 0.925 | 1.329 | 0.8661 |
| 2nd vaccination time [hour] | 0.843 | 0.713 | 0.998 | 0.9772 |
| Moderna [0/1] | 1.386 | 0.398 | 4.802 | 0.7033 |
| Calcineurin inhibitor[0/1] | 1.789 | 0.786 | 4.111 | 0.9225 |
| Steroids [0/1] | 0.728 | 0.346 | 1.541 | 0.7910 |
| Time after TX [10 years] | 2.117 | 1.464 | 3.101 | 1.0000 |
| Diabetes Mellitus [0/1] | 0.911 | 0.556 | 1.498 | 0.6451 |
| Age [10 years] | 0.722 | 0.553 | 0.936 | 0.9925 |
| BMI [5 units] | 0.968 | 0.770 | 1.223 | 0.6114 |
| Albumine [g/L] | 1.034 | 0.973 | 1.101 | 0.8549 |
| Lymphocyte counts [log(10^9/L)] | 2.533 | 1.495 | 4.350 | 0.9998 |
| Male sex [0/1] | 2.015 | 1.291 | 3.201 | 0.9996 |
| eGFR [10 mL/min/1.73 m2] | 1.378 | 1.230 | 1.550 | 1.0000 |
| MMF/MPA [0/1] | 0.120 | 0.069 | 0.205 | 1.0000 |
| Depleting therapy [0/1] | 0.237 | 0.035 | 1.137 | 0.9608 |
2.5 Sensitivity analyses
2.5.1 Model specification sensitivity, reduced_model
Fitting the model
Open code
model_reduced <- run_model(brm(seroconversion ~
vacc_time_scal +
vacc_time2_scal +
rescale(months_afterTX) +
rescale(age) +
rescale(log(Lymphocytes))+
rescale(sex, binary.inputs = '-0.5,0.5') +
rescale(EGFR_vaccination) +
rescale(mmf_mpa, binary.inputs = '-0.5,0.5') +
rescale(depletationTreat_year, binary.inputs = '-0.5,0.5'),
data = findat2,
prior = prior1,
chains = 3, iter = 6000, warmup = 2000, seed = 17,
control = list(adapt_delta = 0.99),
cores = 1,
family = bernoulli(link='logit')),
'gitignore/brm3/model_reduced', reuse = TRUE)
summary(model_reduced)
## Family: bernoulli
## Links: mu = logit
## Formula: seroconversion ~ vacc_time_scal + vacc_time2_scal + rescale(months_afterTX) + rescale(age) + rescale(log(Lymphocytes)) + rescale(sex, binary.inputs = "-0.5,0.5") + rescale(EGFR_vaccination) + rescale(mmf_mpa, binary.inputs = "-0.5,0.5") + rescale(depletationTreat_year, binary.inputs = "-0.5,0.5")
## Data: findat2 (Number of observations: 553)
## Draws: 3 chains, each with iter = 6000; warmup = 2000; thin = 1;
## total post-warmup draws = 12000
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI
## Intercept -0.91 0.44 -1.86
## vacc_time_scal 0.24 0.25 -0.24
## vacc_time2_scal -0.48 0.25 -0.96
## rescalemonths_afterTX 1.03 0.23 0.57
## rescaleage -0.58 0.23 -1.03
## rescalelogLymphocytes 0.84 0.24 0.37
## rescalesexbinary.inputsEQM0.50.5 0.68 0.23 0.24
## rescaleEGFR_vaccination 1.33 0.23 0.90
## rescalemmf_mpabinary.inputsEQM0.50.5 -1.92 0.26 -2.45
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 -1.53 0.87 -3.42
## u-95% CI Rhat Bulk_ESS
## Intercept -0.14 1.00 13277
## vacc_time_scal 0.74 1.00 9819
## vacc_time2_scal -0.00 1.00 9938
## rescalemonths_afterTX 1.49 1.00 13511
## rescaleage -0.12 1.00 11079
## rescalelogLymphocytes 1.31 1.00 14785
## rescalesexbinary.inputsEQM0.50.5 1.13 1.00 14688
## rescaleEGFR_vaccination 1.79 1.00 12876
## rescalemmf_mpabinary.inputsEQM0.50.5 -1.41 1.00 13406
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 -0.02 1.00 13845
## Tail_ESS
## Intercept 8267
## vacc_time_scal 8931
## vacc_time2_scal 8154
## rescalemonths_afterTX 8991
## rescaleage 9462
## rescalelogLymphocytes 8627
## rescalesexbinary.inputsEQM0.50.5 8781
## rescaleEGFR_vaccination 8766
## rescalemmf_mpabinary.inputsEQM0.50.5 8231
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 8225
##
## 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).
repor(model_reduced)| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| vacc_time_scal | 1.271 | 0.789 | 2.086 | 0.8339 |
| vacc_time2_scal | 0.618 | 0.384 | 0.999 | 0.9752 |
| months_afterTX | 2.795 | 1.769 | 4.437 | 1.0000 |
| age | 0.563 | 0.357 | 0.884 | 0.9944 |
| logLymphocytes | 2.306 | 1.452 | 3.719 | 0.9998 |
| sex | 1.970 | 1.268 | 3.093 | 0.9986 |
| EGFR_vaccination | 3.783 | 2.459 | 5.961 | 1.0000 |
| mmf_mpa | 0.146 | 0.086 | 0.243 | 1.0000 |
| depletationTreat_year | 0.217 | 0.033 | 0.979 | 0.9773 |
Reporting: suppl_table3
Open code
sds2 <- c(
sd(findat2$vacc_time_cont)*2,
sd(findat2$vacc_time_2nd_cont)*2,
(sd(findat2$months_afterTX)*2)/120,
(sd(findat2$age)*2)/10,
sd(log(findat2$Lymphocytes))*2,
1, # male sex
(sd(findat2$EGFR_vaccination)*2)/10,
1, ## MMF/MPA
1 # deplating
)
labe2 = c('1st vaccination time [hour]',
'2nd vaccination time [hour]',
'Time after TX [10 years]',
'Age [10 years]',
'Lymphocyte counts [log(10^9/L)]',
'Male sex [0/1]',
'eGFR [10 mL/min/1.73 m2]',
'MMF/MPA [0/1]',
'Depleting therapy [0/1]')
tr2 <- repor(model_reduced, labels=labe2, scals = sds2, nhtm=TRUE)
suppl_table3 <- kableExtra::kable(tr2, caption =
"Supplementary Table 3. Model specification sensitivity analysis, conducted through Bayesian logistic regression, examining the probability of seroconversion following mRNA anti-Covid vaccination. In contrast to the main model, only selected covariates were included here. 'OR' represents the odds ratio, indicating the expected fold-change in odds when the predictor changes by the unit defined in square brackets. 'Q2.5' and 'Q97.5' denote the lower and upper bounds of the 95% credible interval, respectively, while 'PD' shows the probability of direction. Refer to the methods section for details."
)
suppl_table3| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| 1st vaccination time [hour] | 1.089 | 0.919 | 1.300 | 0.8339 |
| 2nd vaccination time [hour] | 0.849 | 0.723 | 1.000 | 0.9752 |
| Time after TX [10 years] | 2.120 | 1.517 | 2.971 | 1.0000 |
| Age [10 years] | 0.753 | 0.602 | 0.941 | 0.9944 |
| Lymphocyte counts [log(10^9/L)] | 2.516 | 1.510 | 4.264 | 0.9998 |
| Male sex [0/1] | 1.970 | 1.268 | 3.093 | 0.9986 |
| eGFR [10 mL/min/1.73 m2] | 1.398 | 1.254 | 1.567 | 1.0000 |
| MMF/MPA [0/1] | 0.146 | 0.086 | 0.243 | 1.0000 |
| Depleting therapy [0/1] | 0.217 | 0.033 | 0.979 | 0.9773 |
2.5.2 Prior sensitivity, model_uniprior
Next we will fit the same models as above, but we will ignore known information about the effect of morning vaccination on the seroconversion. Thus, we will apply the same models as above but with a zero-effect-centered prior (\(\mu = 0\), \(\sigma = 5\)) for all predictors including vaccination time.
Fitting model
Open code
findat2 <- findat2 %>% mutate(vacc_time_scal = rescale(vacc_time_cont),
vacc_time2_scal = rescale(vacc_time_2nd))
model_uniprior <- run_model(brm(seroconversion ~
vacc_time_scal +
vacc_time2_scal +
rescale(vaccine_moderna, binary.inputs = '-0.5,0.5') +
rescale(calcineurin_inhibitor, binary.inputs = '-0.5,0.5') +
rescale(steroids, binary.inputs = '-0.5,0.5') +
rescale(months_afterTX) +
s(rescale(days_after_2nd_dose), k=4) +
rescale(DM, binary.inputs = '-0.5,0.5') +
rescale(age) +
rescale(BMI) +
rescale(Albumine) +
rescale(log(Lymphocytes))+
rescale(sex, binary.inputs = '-0.5,0.5') +
rescale(EGFR_vaccination) +
rescale(mmf_mpa, binary.inputs = '-0.5,0.5') +
rescale(depletationTreat_year, binary.inputs = '-0.5,0.5'),
data = findat2,
prior = prior(normal(0, 5), class = "b"),
chains = 3, iter = 6000, warmup = 2000, seed = 17,
control = list(adapt_delta = 0.99),
cores = 1,
family = bernoulli(link='logit')),
'gitignore/brm3/model_uniprior', reuse = TRUE)
summary(model_uniprior)
## Family: bernoulli
## Links: mu = logit
## Formula: seroconversion ~ vacc_time_scal + vacc_time2_scal + rescale(vaccine_moderna, binary.inputs = "-0.5,0.5") + rescale(calcineurin_inhibitor, binary.inputs = "-0.5,0.5") + rescale(steroids, binary.inputs = "-0.5,0.5") + rescale(months_afterTX) + s(rescale(days_after_2nd_dose), k = 4) + rescale(DM, binary.inputs = "-0.5,0.5") + rescale(age) + rescale(BMI) + rescale(Albumine) + rescale(log(Lymphocytes)) + rescale(sex, binary.inputs = "-0.5,0.5") + rescale(EGFR_vaccination) + rescale(mmf_mpa, binary.inputs = "-0.5,0.5") + rescale(depletationTreat_year, binary.inputs = "-0.5,0.5")
## Data: findat2 (Number of observations: 553)
## Draws: 3 chains, each with iter = 6000; warmup = 2000; thin = 1;
## total post-warmup draws = 12000
##
## Smooth Terms:
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sds(srescaledays_after_2nd_dose_1) 2.21 1.72 0.16 6.49 1.00
## Bulk_ESS Tail_ESS
## sds(srescaledays_after_2nd_dose_1) 4401 3865
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI
## Intercept -0.82 0.60 -2.03
## vacc_time_scal 0.31 0.26 -0.20
## vacc_time2_scal -0.51 0.26 -1.01
## rescalevaccine_modernabinary.inputsEQM0.50.5 0.33 0.62 -0.88
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 0.58 0.41 -0.22
## rescalesteroidsbinary.inputsEQM0.50.5 -0.32 0.39 -1.08
## rescalemonths_afterTX 1.02 0.26 0.51
## rescaleDMbinary.inputsEQM0.50.5 -0.09 0.25 -0.58
## rescaleage -0.65 0.27 -1.18
## rescaleBMI -0.06 0.23 -0.51
## rescaleAlbumine 0.24 0.23 -0.20
## rescalelogLymphocytes 0.84 0.25 0.36
## rescalesexbinary.inputsEQM0.50.5 0.70 0.23 0.26
## rescaleEGFR_vaccination 1.27 0.24 0.82
## rescalemmf_mpabinary.inputsEQM0.50.5 -2.12 0.28 -2.66
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 -1.45 0.90 -3.42
## srescaledays_after_2nd_dose_1 -0.02 2.20 -4.98
## u-95% CI Rhat Bulk_ESS
## Intercept 0.28 1.00 9942
## vacc_time_scal 0.83 1.00 8963
## vacc_time2_scal -0.01 1.00 8714
## rescalevaccine_modernabinary.inputsEQM0.50.5 1.55 1.00 11686
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 1.41 1.00 11187
## rescalesteroidsbinary.inputsEQM0.50.5 0.44 1.00 9680
## rescalemonths_afterTX 1.53 1.00 9189
## rescaleDMbinary.inputsEQM0.50.5 0.39 1.00 11785
## rescaleage -0.12 1.00 9060
## rescaleBMI 0.39 1.00 12504
## rescaleAlbumine 0.69 1.00 11895
## rescalelogLymphocytes 1.34 1.00 10932
## rescalesexbinary.inputsEQM0.50.5 1.15 1.00 13100
## rescaleEGFR_vaccination 1.74 1.00 10800
## rescalemmf_mpabinary.inputsEQM0.50.5 -1.58 1.00 10861
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 0.14 1.00 11259
## srescaledays_after_2nd_dose_1 3.70 1.00 5251
## Tail_ESS
## Intercept 8316
## vacc_time_scal 8560
## vacc_time2_scal 8753
## rescalevaccine_modernabinary.inputsEQM0.50.5 9285
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 9188
## rescalesteroidsbinary.inputsEQM0.50.5 8871
## rescalemonths_afterTX 8834
## rescaleDMbinary.inputsEQM0.50.5 9053
## rescaleage 8964
## rescaleBMI 8590
## rescaleAlbumine 9089
## rescalelogLymphocytes 8912
## rescalesexbinary.inputsEQM0.50.5 9481
## rescaleEGFR_vaccination 9120
## rescalemmf_mpabinary.inputsEQM0.50.5 9459
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 7531
## srescaledays_after_2nd_dose_1 7389
##
## 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).
repor(model_uniprior)| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| vacc_time_scal | 1.364 | 0.820 | 2.292 | 0.8831 |
| vacc_time2_scal | 0.601 | 0.363 | 0.986 | 0.9781 |
| vaccine_moderna | 1.395 | 0.414 | 4.714 | 0.7036 |
| calcineurin_inhibitor | 1.782 | 0.801 | 4.078 | 0.9212 |
| steroids | 0.727 | 0.339 | 1.545 | 0.7962 |
| months_afterTX | 2.776 | 1.670 | 4.637 | 1.0000 |
| DM | 0.909 | 0.558 | 1.480 | 0.6517 |
| age | 0.520 | 0.307 | 0.885 | 0.9933 |
| BMI | 0.942 | 0.601 | 1.472 | 0.6039 |
| Albumine | 1.274 | 0.818 | 2.003 | 0.8538 |
| logLymphocytes | 2.318 | 1.427 | 3.827 | 0.9996 |
| sex | 2.010 | 1.294 | 3.169 | 0.9991 |
| EGFR_vaccination | 3.563 | 2.269 | 5.713 | 1.0000 |
| mmf_mpa | 0.121 | 0.070 | 0.206 | 1.0000 |
| depletationTreat_year | 0.235 | 0.033 | 1.149 | 0.9606 |
Reporting: suppl_table4
Open code
tr1n <- repor(model_uniprior, labels=labe1, scals = sds1, nhtm=TRUE)
# full model
suppl_table4 <- kableExtra::kable(tr1n, caption =
"Supplementary table 4. Prior sensitivity analysis, conducted using Bayesian logistic regression to model the probability of seroconversion following mRNA anti-Covid vaccination. In contrast to the main model, here we applied a zero-centered Gaussian prior also for the effect of vaccination timing. 'OR' represents the odds ratio, indicating the expected fold-change in odds when the predictor changes by the unit defined in square brackets. 'Q2.5' and 'Q97.5' denote the lower and upper bounds of the 95% credible interval, respectively, while 'PD' shows the probability of direction. Please refer to the methods section for additional details.")
suppl_table4| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| 1st vaccination time [hour] | 1.117 | 0.932 | 1.344 | 0.8831 |
| 2nd vaccination time [hour] | 0.842 | 0.709 | 0.995 | 0.9781 |
| Moderna [0/1] | 1.395 | 0.414 | 4.714 | 0.7036 |
| Calcineurin inhibitor[0/1] | 1.782 | 0.801 | 4.078 | 0.9212 |
| Steroids [0/1] | 0.727 | 0.339 | 1.545 | 0.7962 |
| Time after TX [10 years] | 2.109 | 1.454 | 3.068 | 1.0000 |
| Diabetes Mellitus [0/1] | 0.909 | 0.558 | 1.480 | 0.6517 |
| Age [10 years] | 0.724 | 0.559 | 0.942 | 0.9933 |
| BMI [5 units] | 0.970 | 0.769 | 1.221 | 0.6039 |
| Albumine [g/L] | 1.033 | 0.973 | 1.099 | 0.8538 |
| Lymphocyte counts [log(10^9/L)] | 2.531 | 1.480 | 4.402 | 0.9996 |
| Male sex [0/1] | 2.010 | 1.294 | 3.169 | 0.9991 |
| eGFR [10 mL/min/1.73 m2] | 1.377 | 1.229 | 1.551 | 1.0000 |
| MMF/MPA [0/1] | 0.121 | 0.070 | 0.206 | 1.0000 |
| Depleting therapy [0/1] | 0.235 | 0.033 | 1.149 | 0.9606 |
2.5.3 Seroconnversion threshold sensitivity, model_seroconversion_min
Fitting model
Open code
# scaling the continuous time by 2SD and centering to have zero means
findat2 <- findat2 %>% mutate(vacc_time_scal = rescale(vacc_time_cont),
vacc_time2_scal = rescale(vacc_time_2nd))
### Continuous times models
findat2$seroconversion2 <- ifelse(findat2$antibody_level>3.6, 1, 0)
# full model
model_seroconversion_min <- run_model(brm(seroconversion2 ~
vacc_time_scal +
vacc_time2_scal +
rescale(vaccine_moderna, binary.inputs = '-0.5,0.5') +
rescale(calcineurin_inhibitor, binary.inputs = '-0.5,0.5') +
rescale(steroids, binary.inputs = '-0.5,0.5') +
rescale(months_afterTX) +
s(rescale(days_after_2nd_dose), k=4) +
rescale(DM, binary.inputs = '-0.5,0.5') +
rescale(age) +
rescale(BMI) +
rescale(Albumine) +
rescale(log(Lymphocytes))+
rescale(sex, binary.inputs = '-0.5,0.5') +
rescale(EGFR_vaccination) +
rescale(mmf_mpa, binary.inputs = '-0.5,0.5') +
rescale(depletationTreat_year, binary.inputs = '-0.5,0.5'),
data = findat2,
prior = prior1,
chains = 3, iter = 6000, warmup = 2000, seed = 17,
control = list(adapt_delta = 0.99),
cores = 1,
family = bernoulli(link='logit')),
'gitignore/brm3/model_seroconversion_min', reuse = TRUE)
summary(model_seroconversion_min)
## Family: bernoulli
## Links: mu = logit
## Formula: seroconversion2 ~ vacc_time_scal + vacc_time2_scal + rescale(vaccine_moderna, binary.inputs = "-0.5,0.5") + rescale(calcineurin_inhibitor, binary.inputs = "-0.5,0.5") + rescale(steroids, binary.inputs = "-0.5,0.5") + rescale(months_afterTX) + s(rescale(days_after_2nd_dose), k = 4) + rescale(DM, binary.inputs = "-0.5,0.5") + rescale(age) + rescale(BMI) + rescale(Albumine) + rescale(log(Lymphocytes)) + rescale(sex, binary.inputs = "-0.5,0.5") + rescale(EGFR_vaccination) + rescale(mmf_mpa, binary.inputs = "-0.5,0.5") + rescale(depletationTreat_year, binary.inputs = "-0.5,0.5")
## Data: findat2 (Number of observations: 553)
## Draws: 3 chains, each with iter = 6000; warmup = 2000; thin = 1;
## total post-warmup draws = 12000
##
## Smooth Terms:
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sds(srescaledays_after_2nd_dose_1) 2.50 1.68 0.26 6.61 1.00
## Bulk_ESS Tail_ESS
## sds(srescaledays_after_2nd_dose_1) 5548 4031
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI
## Intercept 0.29 0.57 -0.84
## vacc_time_scal 0.23 0.26 -0.27
## vacc_time2_scal -0.48 0.25 -0.96
## rescalevaccine_modernabinary.inputsEQM0.50.5 0.88 0.68 -0.44
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 0.44 0.39 -0.31
## rescalesteroidsbinary.inputsEQM0.50.5 -0.77 0.39 -1.55
## rescalemonths_afterTX 0.81 0.26 0.32
## rescaleDMbinary.inputsEQM0.50.5 -0.16 0.24 -0.63
## rescaleage -0.67 0.27 -1.22
## rescaleBMI -0.07 0.22 -0.50
## rescaleAlbumine 0.26 0.22 -0.17
## rescalelogLymphocytes 1.05 0.24 0.59
## rescalesexbinary.inputsEQM0.50.5 0.89 0.22 0.46
## rescaleEGFR_vaccination 1.01 0.23 0.56
## rescalemmf_mpabinary.inputsEQM0.50.5 -2.02 0.28 -2.60
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 -1.26 0.76 -2.87
## srescaledays_after_2nd_dose_1 0.40 2.16 -4.20
## u-95% CI Rhat Bulk_ESS
## Intercept 1.39 1.00 15154
## vacc_time_scal 0.73 1.00 12815
## vacc_time2_scal -0.00 1.00 13158
## rescalevaccine_modernabinary.inputsEQM0.50.5 2.27 1.00 15963
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 1.20 1.00 14890
## rescalesteroidsbinary.inputsEQM0.50.5 -0.03 1.00 14911
## rescalemonths_afterTX 1.32 1.00 12947
## rescaleDMbinary.inputsEQM0.50.5 0.31 1.00 16642
## rescaleage -0.15 1.00 13431
## rescaleBMI 0.35 1.00 14262
## rescaleAlbumine 0.71 1.00 15466
## rescalelogLymphocytes 1.53 1.00 15000
## rescalesexbinary.inputsEQM0.50.5 1.33 1.00 15121
## rescaleEGFR_vaccination 1.46 1.00 14675
## rescalemmf_mpabinary.inputsEQM0.50.5 -1.48 1.00 13860
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 0.13 1.00 15596
## srescaledays_after_2nd_dose_1 4.39 1.00 8804
## Tail_ESS
## Intercept 8646
## vacc_time_scal 9576
## vacc_time2_scal 9386
## rescalevaccine_modernabinary.inputsEQM0.50.5 8822
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 8936
## rescalesteroidsbinary.inputsEQM0.50.5 9406
## rescalemonths_afterTX 9041
## rescaleDMbinary.inputsEQM0.50.5 9390
## rescaleage 8798
## rescaleBMI 9031
## rescaleAlbumine 8625
## rescalelogLymphocytes 8817
## rescalesexbinary.inputsEQM0.50.5 7632
## rescaleEGFR_vaccination 8994
## rescalemmf_mpabinary.inputsEQM0.50.5 9507
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 7933
## srescaledays_after_2nd_dose_1 7800
##
## 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).
repor(model_seroconversion_min)| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| vacc_time_scal | 1.262 | 0.760 | 2.081 | 0.8165 |
| vacc_time2_scal | 0.619 | 0.384 | 0.998 | 0.9757 |
| vaccine_moderna | 2.417 | 0.643 | 9.665 | 0.9058 |
| calcineurin_inhibitor | 1.546 | 0.737 | 3.320 | 0.8718 |
| steroids | 0.461 | 0.213 | 0.966 | 0.9803 |
| months_afterTX | 2.253 | 1.375 | 3.757 | 0.9996 |
| DM | 0.852 | 0.530 | 1.367 | 0.7454 |
| age | 0.511 | 0.297 | 0.864 | 0.9939 |
| BMI | 0.934 | 0.606 | 1.421 | 0.6203 |
| Albumine | 1.301 | 0.845 | 2.041 | 0.8825 |
| logLymphocytes | 2.865 | 1.796 | 4.635 | 1.0000 |
| sex | 2.429 | 1.580 | 3.786 | 0.9998 |
| EGFR_vaccination | 2.733 | 1.756 | 4.308 | 1.0000 |
| mmf_mpa | 0.133 | 0.074 | 0.228 | 1.0000 |
| depletationTreat_year | 0.284 | 0.057 | 1.138 | 0.9623 |
Reporting: suppl_table5
Open code
tr1ss <- repor(model_seroconversion_min, labels=labe1, scals = sds1, nhtm=TRUE)
# full model
suppl_table5 <- kableExtra::kable(tr1ss, caption =
"Supplementary table 5. Seroconversion threshold sensitivity analysis, conducted using Bayesian logistic regression to examine the probability of seroconversion following mRNA anti-Covid vaccination. In contrast to the main model, seroconversion is defined here as the minimum detectable level of SARS-CoV-2 IgG antibody (>3.6). 'OR' represents the odds ratio, indicating the expected fold-change in odds when the predictor changes by the unit defined in square brackets. 'Q2.5' and 'Q97.5' denote the lower and upper bounds of the 95% credible interval, respectively, while 'PD' shows the probability of direction. Please refer to the methods section for additional details.")
suppl_table5| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| 1st vaccination time [hour] | 1.086 | 0.907 | 1.299 | 0.8165 |
| 2nd vaccination time [hour] | 0.850 | 0.723 | 0.999 | 0.9757 |
| Moderna [0/1] | 2.417 | 0.643 | 9.665 | 0.9058 |
| Calcineurin inhibitor[0/1] | 1.546 | 0.737 | 3.320 | 0.8718 |
| Steroids [0/1] | 0.461 | 0.213 | 0.966 | 0.9803 |
| Time after TX [10 years] | 1.810 | 1.262 | 2.631 | 0.9996 |
| Diabetes Mellitus [0/1] | 0.852 | 0.530 | 1.367 | 0.7454 |
| Age [10 years] | 0.718 | 0.550 | 0.930 | 0.9939 |
| BMI [5 units] | 0.965 | 0.772 | 1.199 | 0.6203 |
| Albumine [g/L] | 1.036 | 0.977 | 1.101 | 0.8825 |
| Lymphocyte counts [log(10^9/L)] | 3.197 | 1.909 | 5.439 | 1.0000 |
| Male sex [0/1] | 2.429 | 1.580 | 3.786 | 0.9998 |
| eGFR [10 mL/min/1.73 m2] | 1.288 | 1.152 | 1.444 | 1.0000 |
| MMF/MPA [0/1] | 0.133 | 0.074 | 0.228 | 1.0000 |
| Depleting therapy [0/1] | 0.284 | 0.057 | 1.138 | 0.9623 |
2.5.4 Models with dichotomized time
(“morning” vs. “afternoon”, with two different definitions of afternoon)
Open code
findat2 <- findat2 %>%
mutate(
afternoon1 = if_else(vacc_time_cont>13, 1, 0),
afternoon2 = if_else(vacc_time_2nd_cont>13, 1, 0)
)
model_dichotomized_13 <- run_model(brm(seroconversion ~
afternoon1 +
afternoon2 +
rescale(vaccine_moderna, binary.inputs = '-0.5,0.5') +
rescale(calcineurin_inhibitor, binary.inputs = '-0.5,0.5') +
rescale(steroids, binary.inputs = '-0.5,0.5') +
rescale(months_afterTX) +
s(rescale(days_after_2nd_dose), k=4) +
rescale(DM, binary.inputs = '-0.5,0.5') +
rescale(age) +
rescale(BMI) +
rescale(Albumine) +
rescale(log(Lymphocytes))+
rescale(sex, binary.inputs = '-0.5,0.5') +
rescale(EGFR_vaccination) +
rescale(mmf_mpa, binary.inputs = '-0.5,0.5') +
rescale(depletationTreat_year, binary.inputs = '-0.5,0.5'),
data = findat2,
prior = prior2,
chains = 3, iter = 6000, warmup = 2000, seed = 17,
control = list(adapt_delta = 0.99),
cores = 1,
family = bernoulli(link='logit')),
'gitignore/brm3/model_dichotomized_13', reuse = TRUE)
summary(model_dichotomized_13)
## Warning: There were 2 divergent transitions after
## warmup. Increasing adapt_delta above 0.99 may help. See
## http://mc-stan.org/misc/warnings.html#divergent-transitions-after-warmup
## Family: bernoulli
## Links: mu = logit
## Formula: seroconversion ~ afternoon1 + afternoon2 + rescale(vaccine_moderna, binary.inputs = "-0.5,0.5") + rescale(calcineurin_inhibitor, binary.inputs = "-0.5,0.5") + rescale(steroids, binary.inputs = "-0.5,0.5") + rescale(months_afterTX) + s(rescale(days_after_2nd_dose), k = 4) + rescale(DM, binary.inputs = "-0.5,0.5") + rescale(age) + rescale(BMI) + rescale(Albumine) + rescale(log(Lymphocytes)) + rescale(sex, binary.inputs = "-0.5,0.5") + rescale(EGFR_vaccination) + rescale(mmf_mpa, binary.inputs = "-0.5,0.5") + rescale(depletationTreat_year, binary.inputs = "-0.5,0.5")
## Data: findat2 (Number of observations: 553)
## Draws: 3 chains, each with iter = 6000; warmup = 2000; thin = 1;
## total post-warmup draws = 12000
##
## Smooth Terms:
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sds(srescaledays_after_2nd_dose_1) 2.31 1.75 0.17 6.50 1.00
## Bulk_ESS Tail_ESS
## sds(srescaledays_after_2nd_dose_1) 4722 3487
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI
## Intercept -0.63 0.60 -1.87
## afternoon1 0.15 0.26 -0.37
## afternoon2 -0.60 0.25 -1.10
## rescalevaccine_modernabinary.inputsEQM0.50.5 0.33 0.63 -0.93
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 0.57 0.42 -0.23
## rescalesteroidsbinary.inputsEQM0.50.5 -0.33 0.38 -1.09
## rescalemonths_afterTX 1.05 0.27 0.52
## rescaleDMbinary.inputsEQM0.50.5 -0.07 0.25 -0.56
## rescaleage -0.76 0.26 -1.27
## rescaleBMI -0.06 0.22 -0.49
## rescaleAlbumine 0.24 0.23 -0.21
## rescalelogLymphocytes 0.82 0.25 0.34
## rescalesexbinary.inputsEQM0.50.5 0.67 0.23 0.23
## rescaleEGFR_vaccination 1.30 0.24 0.84
## rescalemmf_mpabinary.inputsEQM0.50.5 -2.13 0.28 -2.71
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 -1.38 0.88 -3.29
## srescaledays_after_2nd_dose_1 -0.01 2.18 -4.94
## u-95% CI Rhat Bulk_ESS
## Intercept 0.50 1.00 13475
## afternoon1 0.66 1.00 13444
## afternoon2 -0.11 1.00 13389
## rescalevaccine_modernabinary.inputsEQM0.50.5 1.58 1.00 15161
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 1.41 1.00 15389
## rescalesteroidsbinary.inputsEQM0.50.5 0.43 1.00 15027
## rescalemonths_afterTX 1.58 1.00 12808
## rescaleDMbinary.inputsEQM0.50.5 0.42 1.00 14936
## rescaleage -0.24 1.00 13344
## rescaleBMI 0.38 1.00 15149
## rescaleAlbumine 0.71 1.00 15264
## rescalelogLymphocytes 1.32 1.00 16397
## rescalesexbinary.inputsEQM0.50.5 1.12 1.00 16822
## rescaleEGFR_vaccination 1.76 1.00 14805
## rescalemmf_mpabinary.inputsEQM0.50.5 -1.59 1.00 12815
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 0.16 1.00 15691
## srescaledays_after_2nd_dose_1 3.73 1.00 7238
## Tail_ESS
## Intercept 8396
## afternoon1 9348
## afternoon2 9740
## rescalevaccine_modernabinary.inputsEQM0.50.5 8917
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 9384
## rescalesteroidsbinary.inputsEQM0.50.5 9804
## rescalemonths_afterTX 9727
## rescaleDMbinary.inputsEQM0.50.5 8954
## rescaleage 8805
## rescaleBMI 9321
## rescaleAlbumine 9753
## rescalelogLymphocytes 9273
## rescalesexbinary.inputsEQM0.50.5 8883
## rescaleEGFR_vaccination 8935
## rescalemmf_mpabinary.inputsEQM0.50.5 8659
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 7412
## srescaledays_after_2nd_dose_1 8385
##
## 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).
repor(model_dichotomized_13)| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| afternoon1 | 1.158 | 0.689 | 1.928 | 0.7152 |
| afternoon2 | 0.550 | 0.333 | 0.898 | 0.9922 |
| vaccine_moderna | 1.393 | 0.395 | 4.855 | 0.7047 |
| calcineurin_inhibitor | 1.774 | 0.793 | 4.080 | 0.9152 |
| steroids | 0.718 | 0.337 | 1.545 | 0.8105 |
| months_afterTX | 2.855 | 1.688 | 4.847 | 1.0000 |
| DM | 0.929 | 0.573 | 1.521 | 0.6148 |
| age | 0.467 | 0.280 | 0.785 | 0.9982 |
| BMI | 0.945 | 0.610 | 1.462 | 0.5973 |
| Albumine | 1.278 | 0.814 | 2.043 | 0.8598 |
| logLymphocytes | 2.274 | 1.401 | 3.729 | 0.9998 |
| sex | 1.956 | 1.256 | 3.066 | 0.9989 |
| EGFR_vaccination | 3.659 | 2.319 | 5.807 | 1.0000 |
| mmf_mpa | 0.119 | 0.066 | 0.204 | 1.0000 |
| depletationTreat_year | 0.253 | 0.037 | 1.169 | 0.9572 |
Open code
findat2 <- findat2 %>%
mutate(
afternoon1 = if_else(vacc_time_cont>12, 1, 0),
afternoon2 = if_else(vacc_time_2nd_cont>12, 1, 0)
)
model_dichotomized_12 <- run_model(brm(seroconversion ~
afternoon1 +
afternoon2 +
rescale(vaccine_moderna, binary.inputs = '-0.5,0.5') +
rescale(calcineurin_inhibitor, binary.inputs = '-0.5,0.5') +
rescale(steroids, binary.inputs = '-0.5,0.5') +
rescale(months_afterTX) +
s(rescale(days_after_2nd_dose), k=4) +
rescale(DM, binary.inputs = '-0.5,0.5') +
rescale(age) +
rescale(BMI) +
rescale(Albumine) +
rescale(log(Lymphocytes))+
rescale(sex, binary.inputs = '-0.5,0.5') +
rescale(EGFR_vaccination) +
rescale(mmf_mpa, binary.inputs = '-0.5,0.5') +
rescale(depletationTreat_year, binary.inputs = '-0.5,0.5'),
data = findat2,
prior = prior2,
chains = 3, iter = 6000, warmup = 2000, seed = 17,
control = list(adapt_delta = 0.99),
cores = 1,
family = bernoulli(link='logit')),
'gitignore/brm3/model_dichotomized_12', reuse = TRUE)
summary(model_dichotomized_12)
## Family: bernoulli
## Links: mu = logit
## Formula: seroconversion ~ afternoon1 + afternoon2 + rescale(vaccine_moderna, binary.inputs = "-0.5,0.5") + rescale(calcineurin_inhibitor, binary.inputs = "-0.5,0.5") + rescale(steroids, binary.inputs = "-0.5,0.5") + rescale(months_afterTX) + s(rescale(days_after_2nd_dose), k = 4) + rescale(DM, binary.inputs = "-0.5,0.5") + rescale(age) + rescale(BMI) + rescale(Albumine) + rescale(log(Lymphocytes)) + rescale(sex, binary.inputs = "-0.5,0.5") + rescale(EGFR_vaccination) + rescale(mmf_mpa, binary.inputs = "-0.5,0.5") + rescale(depletationTreat_year, binary.inputs = "-0.5,0.5")
## Data: findat2 (Number of observations: 553)
## Draws: 3 chains, each with iter = 6000; warmup = 2000; thin = 1;
## total post-warmup draws = 12000
##
## Smooth Terms:
## Estimate Est.Error l-95% CI u-95% CI Rhat
## sds(srescaledays_after_2nd_dose_1) 2.24 1.64 0.17 6.20 1.00
## Bulk_ESS Tail_ESS
## sds(srescaledays_after_2nd_dose_1) 5711 5053
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI
## Intercept -0.73 0.62 -1.98
## afternoon1 0.10 0.24 -0.36
## afternoon2 -0.23 0.24 -0.71
## rescalevaccine_modernabinary.inputsEQM0.50.5 0.24 0.63 -1.01
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 0.55 0.42 -0.25
## rescalesteroidsbinary.inputsEQM0.50.5 -0.34 0.38 -1.10
## rescalemonths_afterTX 1.00 0.26 0.50
## rescaleDMbinary.inputsEQM0.50.5 -0.08 0.25 -0.57
## rescaleage -0.74 0.27 -1.28
## rescaleBMI -0.07 0.23 -0.52
## rescaleAlbumine 0.26 0.23 -0.19
## rescalelogLymphocytes 0.82 0.25 0.32
## rescalesexbinary.inputsEQM0.50.5 0.69 0.23 0.24
## rescaleEGFR_vaccination 1.26 0.23 0.81
## rescalemmf_mpabinary.inputsEQM0.50.5 -2.10 0.28 -2.65
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 -1.39 0.88 -3.27
## srescaledays_after_2nd_dose_1 -0.19 2.18 -5.05
## u-95% CI Rhat Bulk_ESS
## Intercept 0.40 1.00 15863
## afternoon1 0.57 1.00 14797
## afternoon2 0.24 1.00 14699
## rescalevaccine_modernabinary.inputsEQM0.50.5 1.48 1.00 16060
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 1.39 1.00 18358
## rescalesteroidsbinary.inputsEQM0.50.5 0.42 1.00 14455
## rescalemonths_afterTX 1.52 1.00 14042
## rescaleDMbinary.inputsEQM0.50.5 0.40 1.00 17431
## rescaleage -0.22 1.00 13410
## rescaleBMI 0.38 1.00 18016
## rescaleAlbumine 0.72 1.00 15888
## rescalelogLymphocytes 1.32 1.00 15287
## rescalesexbinary.inputsEQM0.50.5 1.16 1.00 18583
## rescaleEGFR_vaccination 1.72 1.00 14882
## rescalemmf_mpabinary.inputsEQM0.50.5 -1.57 1.00 13795
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 0.14 1.00 17402
## srescaledays_after_2nd_dose_1 3.58 1.00 8115
## Tail_ESS
## Intercept 8737
## afternoon1 9378
## afternoon2 9373
## rescalevaccine_modernabinary.inputsEQM0.50.5 9035
## rescalecalcineurin_inhibitorbinary.inputsEQM0.50.5 9672
## rescalesteroidsbinary.inputsEQM0.50.5 8744
## rescalemonths_afterTX 9521
## rescaleDMbinary.inputsEQM0.50.5 8313
## rescaleage 9697
## rescaleBMI 9027
## rescaleAlbumine 8305
## rescalelogLymphocytes 8601
## rescalesexbinary.inputsEQM0.50.5 8098
## rescaleEGFR_vaccination 9090
## rescalemmf_mpabinary.inputsEQM0.50.5 9387
## rescaledepletationTreat_yearbinary.inputsEQM0.50.5 7666
## srescaledays_after_2nd_dose_1 8292
##
## 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).
repor(model_dichotomized_12)| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| afternoon1 | 1.110 | 0.696 | 1.768 | 0.6705 |
| afternoon2 | 0.793 | 0.493 | 1.274 | 0.8290 |
| vaccine_moderna | 1.272 | 0.365 | 4.393 | 0.6509 |
| calcineurin_inhibitor | 1.730 | 0.777 | 4.008 | 0.9050 |
| steroids | 0.714 | 0.334 | 1.526 | 0.8122 |
| months_afterTX | 2.722 | 1.654 | 4.573 | 1.0000 |
| DM | 0.924 | 0.567 | 1.494 | 0.6289 |
| age | 0.475 | 0.278 | 0.805 | 0.9971 |
| BMI | 0.934 | 0.595 | 1.457 | 0.6181 |
| Albumine | 1.292 | 0.828 | 2.051 | 0.8661 |
| logLymphocytes | 2.270 | 1.380 | 3.751 | 0.9997 |
| sex | 1.995 | 1.271 | 3.179 | 0.9987 |
| EGFR_vaccination | 3.523 | 2.249 | 5.600 | 1.0000 |
| mmf_mpa | 0.122 | 0.070 | 0.209 | 1.0000 |
| depletationTreat_year | 0.248 | 0.038 | 1.151 | 0.9592 |
Report: suppl_table6
Open code
sds2 <- c(
1,
1,
1, # moderna
1, # caclcineurin inhibitors
1, # steroids
(sd(findat2$months_afterTX)*2)/120,
1, # DM
(sd(findat2$age)*2)/10,
(sd(findat2$BMI)*2)/5,
sd(findat2$Albumine)*2,
sd(log(findat2$Lymphocytes))*2,
1, # male sex
(sd(findat2$EGFR_vaccination)*2)/10,
1, ## MMF/MPA
1 # deplating
)
labe2 <- c('1st dose afternoon[0/1]',
'2nd dose afternoon [0/1]',
'Moderna [0/1]',
'Calcineurin inhibitor[0/1]',
'Steroids [0/1]',
'Time after TX [10 years]',
'Diabetes Mellitus [0/1]',
'Age [10 years]',
'BMI [5 units]',
'Albumine [g/L]',
'Lymphocyte counts [log(10^9/L)]',
'Male sex [0/1]',
'eGFR [10 mL/min/1.73 m2]',
'MMF/MPA [0/1]',
'Depleting therapy [0/1]')
tr4 <- repor(model_dichotomized_12, labels=labe2, scals = sds2, nhtm=TRUE)
tr6 <- repor(model_dichotomized_13, labels=labe2, scals = sds2, nhtm=TRUE)
tr46 <- cbind(tr4, tr6)
# full model
suppl_table6 <- kableExtra::kable(tr46, caption =
"Supplementary Table 6. Results of Bayesian logistic regression modeling he probability of seroconversion following mRNA anti-Covid vaccination. In contrast to the main model, here we used dichotomized vaccination times ('afternoon' vs. 'morning') with two different thresholds (12 pm and 1pm, left and right section respectively). 'OR' implies odds ratio (expected fold-change in odds when predictor changes by the unit defined in the square brackets). 'Q2.5' and 'Q97.5' are lower and upper bounds of 95% credible interval. 'PD' implies proability of direction") %>%
add_header_above(c(" " = 1, "Threshold: 12 pm" = 4, "Threshold: 1 pm" = 4)) %>%
column_spec(1, width_min = '2.5in')
suppl_table6| OR | Q2.5 | Q97.5 | PD | OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|---|---|---|---|
| 1st dose afternoon[0/1] | 1.110 | 0.696 | 1.768 | 0.6705 | 1.158 | 0.689 | 1.928 | 0.7152 |
| 2nd dose afternoon [0/1] | 0.793 | 0.493 | 1.274 | 0.8290 | 0.550 | 0.333 | 0.898 | 0.9922 |
| Moderna [0/1] | 1.272 | 0.365 | 4.393 | 0.6509 | 1.393 | 0.395 | 4.855 | 0.7047 |
| Calcineurin inhibitor[0/1] | 1.730 | 0.777 | 4.008 | 0.9050 | 1.774 | 0.793 | 4.080 | 0.9152 |
| Steroids [0/1] | 0.714 | 0.334 | 1.526 | 0.8122 | 0.718 | 0.337 | 1.545 | 0.8105 |
| Time after TX [10 years] | 2.079 | 1.444 | 3.037 | 1.0000 | 2.152 | 1.466 | 3.169 | 1.0000 |
| Diabetes Mellitus [0/1] | 0.924 | 0.567 | 1.494 | 0.6289 | 0.929 | 0.573 | 1.521 | 0.6148 |
| Age [10 years] | 0.693 | 0.532 | 0.899 | 0.9971 | 0.688 | 0.534 | 0.888 | 0.9982 |
| BMI [5 units] | 0.965 | 0.765 | 1.214 | 0.6181 | 0.971 | 0.775 | 1.216 | 0.5973 |
| Albumine [g/L] | 1.035 | 0.975 | 1.102 | 0.8661 | 1.034 | 0.972 | 1.102 | 0.8598 |
| Lymphocyte counts [log(10^9/L)] | 2.473 | 1.427 | 4.306 | 0.9997 | 2.478 | 1.452 | 4.278 | 0.9998 |
| Male sex [0/1] | 1.995 | 1.271 | 3.179 | 0.9987 | 1.956 | 1.256 | 3.066 | 0.9989 |
| eGFR [10 mL/min/1.73 m2] | 1.373 | 1.226 | 1.543 | 1.0000 | 1.386 | 1.236 | 1.557 | 1.0000 |
| MMF/MPA [0/1] | 0.122 | 0.070 | 0.209 | 1.0000 | 0.119 | 0.066 | 0.204 | 1.0000 |
| Depleting therapy [0/1] | 0.248 | 0.038 | 1.151 | 0.9592 | 0.253 | 0.037 | 1.169 | 0.9572 |
2.5.5 Quantile regression analysis, model_quantile
This is another sensitivity analysis with respect to the seroconversion definition. Here are the IgG levels modelled directly.
Due to detection limit (3.6 on left and 8000 on right), it is difficult to model the antibody levels with common regression methods. We will write the censored value as 2/3 of the detection limit on the left and 3/2 of the detection limit on the right. Next, we will provide quantile regression which is not affected by specific values and is more robust toward assumptions of parametric methods. We will use log2-transformed antibody level as the outcome as it has strongly right-tailed distribution. As rq function does not support splines, days_after_2nd_dose will be fitted as a predictor with linear function.
Open code
set.seed(17)
findat2 <- findat2 %>% mutate(
antibody_level_cens = if_else(antibody_level <= 3.6, 2.4,
antibody_level)) %>%
mutate(antibody_level_cens = if_else(antibody_level_cens >= 8000, 12000,
antibody_level_cens))
model_quant <- rq(log2(antibody_level_cens) ~
vacc_time_cont +
vacc_time_2nd_cont +
vaccine_moderna +
calcineurin_inhibitor +
steroids +
I(months_afterTX/120) +
DM +
I(age/10) +
I(BMI/5) +
Albumine +
log(Lymphocytes) +
sex +
I(EGFR_vaccination/10) +
mmf_mpa +
depletationTreat_year +
I(days_after_2nd_dose/7),
data = findat2)
kableExtra::kable(summary(model_quant)$coefficients[,])| coefficients | lower bd | upper bd | |
|---|---|---|---|
| (Intercept) | 3.5531088 | 0.9903022 | 6.6424111 |
| vacc_time_cont | -0.0447014 | -0.1161227 | 0.0689589 |
| vacc_time_2nd_cont | -0.0683017 | -0.1706167 | -0.0113548 |
| vaccine_moderna | 1.2152475 | -0.7236033 | 2.8099856 |
| calcineurin_inhibitor | -0.0360449 | -0.2570249 | 0.4818535 |
| steroids | -0.3691380 | -1.3110228 | -0.0708419 |
| I(months_afterTX/120) | 0.7285906 | 0.5227586 | 1.0693897 |
| DM | 0.1656071 | -0.1358856 | 0.4255623 |
| I(age/10) | -0.3521444 | -0.5755477 | -0.1525002 |
| I(BMI/5) | -0.0474585 | -0.1681723 | 0.0747839 |
| Albumine | 0.0603265 | 0.0139782 | 0.0865520 |
| log(Lymphocytes) | 0.4647274 | 0.1132845 | 0.9065730 |
| sex | 0.6639333 | 0.4261897 | 0.9147813 |
| I(EGFR_vaccination/10) | 0.2491023 | 0.1715190 | 0.3328568 |
| mmf_mpa | -2.4599550 | -2.8039815 | -1.9694721 |
| depletationTreat_year | -0.2040382 | -0.4739763 | 0.5367153 |
| I(days_after_2nd_dose/7) | 0.0196654 | -0.0145490 | 0.0595691 |
Report: suppl_table7
Open code
tr <- data.frame(estimate = summary(model_quant)$coefficients[-1, 1],
Q2.5 = summary(model_quant)$coefficients[-1,2],
Q97.5 = summary(model_quant)$coefficients[-1,3])
tr <- round(tr, 2)
rownames(tr) <- c(labe1, 'Weeks after 2nd dose')
suppl_table7 <- kableExtra::kable(tr, caption =
"Supplementary Table 7. Results of quantile regression modelling for the median level of log2-transformed SARS-CoV-2 IgG. 'Estimate' indicates the expected changes in log2(IgG) when the predictor changes by the unit defined in the square brackets. 'Q2.5' and 'Q97.5' represent the lower and upper bounds of the 95% confidence interval, respectively.")
suppl_table7| estimate | Q2.5 | Q97.5 | |
|---|---|---|---|
| 1st vaccination time [hour] | -0.04 | -0.12 | 0.07 |
| 2nd vaccination time [hour] | -0.07 | -0.17 | -0.01 |
| Moderna [0/1] | 1.22 | -0.72 | 2.81 |
| Calcineurin inhibitor[0/1] | -0.04 | -0.26 | 0.48 |
| Steroids [0/1] | -0.37 | -1.31 | -0.07 |
| Time after TX [10 years] | 0.73 | 0.52 | 1.07 |
| Diabetes Mellitus [0/1] | 0.17 | -0.14 | 0.43 |
| Age [10 years] | -0.35 | -0.58 | -0.15 |
| BMI [5 units] | -0.05 | -0.17 | 0.07 |
| Albumine [g/L] | 0.06 | 0.01 | 0.09 |
| Lymphocyte counts [log(10^9/L)] | 0.46 | 0.11 | 0.91 |
| Male sex [0/1] | 0.66 | 0.43 | 0.91 |
| eGFR [10 mL/min/1.73 m2] | 0.25 | 0.17 | 0.33 |
| MMF/MPA [0/1] | -2.46 | -2.80 | -1.97 |
| Depleting therapy [0/1] | -0.20 | -0.47 | 0.54 |
| Weeks after 2nd dose | 0.02 | -0.01 | 0.06 |
Estimation of median IgG when the 2nd dose was administrated at 9 am vs. 3 pm
Open code
model_quant_cent <- rq((antibody_level_cens) ~
rescale(vacc_time_cont) +
vacc_time_2nd_cont +
rescale(vaccine_moderna) +
rescale(calcineurin_inhibitor) +
rescale(steroids) +
rescale(months_afterTX) +
rescale(DM) +
rescale(age) +
rescale(BMI) +
rescale(Albumine) +
rescale(log(Lymphocytes)) +
rescale(sex) +
rescale(EGFR_vaccination) +
rescale(mmf_mpa) +
rescale(depletationTreat_year) +
rescale(days_after_2nd_dose),
data = findat2)
## predicted median IgG at 9 am
pred9 <- 2^(coef(model_quant_cent)[1] + (9*coef(model_quant_cent)[3]))
pred9
## (Intercept)
## 1949.066
## predicted median IgG at 15 pm
pred15 <- 2^(coef(model_quant_cent)[1] + (15*coef(model_quant_cent)[3]))
pred15
## (Intercept)
## 549.356
## predicted percentage change of IgG from 9 am to 3 pm
## (both calculation should show the same results)
(pred15/pred9)
## (Intercept)
## 0.281856
(2^(coef(model_quant_cent)[3]))**6
## vacc_time_2nd_cont
## 0.281856
## predicted percentage change of median IgG by 1 hour
## (both calculation should lead to the same result)
2^(coef(model_quant_cent)[3])
## vacc_time_2nd_cont
## 0.8097256
(2^(coef(model_quant_cent)[1] + (10*coef(model_quant_cent)[3])))/
(2^(coef(model_quant_cent)[1] + (9*coef(model_quant_cent)[3])))
## (Intercept)
## 0.80972562.5.6 Fig1c: joint visualisation of sensitivity analyses
Data preparation
Open code
cole <- c('#CD7006', '#0028F0')
sd1 <- sd(findat2$vacc_time_cont)*2
sd2 <- sd(findat2$vacc_time_2nd_cont)*2
## model specification sensitivity
tr_reduced <- as_tibble(model_reduced) %>%
mutate(dose_1_reduced = exp(b_vacc_time_scal/sd1),
dose_2_reduced = exp(b_vacc_time2_scal/sd2)) %>%
select(dose_1_reduced, dose_2_reduced) %>% data.frame()
CI_reduced <- sapply(tr_reduced,
function(p) quantile(p, probs = c(1/40, 39/40)))
## prior sensitivity
tr_unprior <- as_tibble(model_uniprior) %>%
mutate(dose_1_unprior = exp(b_vacc_time_scal/sd1),
dose_2_unprior = exp(b_vacc_time2_scal/sd2)) %>%
select(dose_1_unprior, dose_2_unprior) %>% data.frame()
## seroconversion definition sensitivity
tr_seroconversion_min <- as_tibble(model_seroconversion_min) %>%
mutate(dose_1_seroconversion_min = exp(b_vacc_time_scal/sd1),
dose_2_seroconversion_min = exp(b_vacc_time2_scal/sd2)) %>%
select(dose_1_seroconversion_min,
dose_2_seroconversion_min) %>% data.frame()
## dichotomized time 12 pm
tr_dichotomized_12 <- as_tibble(model_dichotomized_12) %>%
mutate(dose_1_dichotomized_12 = exp(b_afternoon1),
dose_2_dichotomized_12 = exp(b_afternoon2)) %>%
select(dose_1_dichotomized_12,
dose_2_dichotomized_12) %>% data.frame()
## dichotomized time 1 pm
tr_dichotomized_13 <- as_tibble(model_dichotomized_13) %>%
mutate(dose_1_dichotomized_13 = exp(b_afternoon1),
dose_2_dichotomized_13 = exp(b_afternoon2)) %>%
select(dose_1_dichotomized_13,
dose_2_dichotomized_13) %>% data.frame()
tr_joint <- cbind(tr_reduced,
tr_unprior,
tr_seroconversion_min,
tr_dichotomized_12,
tr_dichotomized_13)
CI_joint <- t(sapply(tr_joint,
function(p) quantile(p, probs = c(1/40, 39/40, 0.055, 0.945, 0.5)))) %>%
data.frame() %>%
rownames_to_column(var = "cat") %>%
mutate(dose = if_else(grepl("dose_1", cat), "dose1", "dose2")) %>%
mutate(X97.5. = if_else(dose == 'dose2',
round(X97.5., 4),
round(X97.5., 2)),
X2.5. = round(X2.5., 2),
X50. = round(X50., 2),
key = c("Reduced", "Reduced",
"non-inf. pr.", "non-inf. pr.",
"Sero+: 3.6", "Sero+: 3.6",
"By 12 pm", "By 12 pm",
"By 1 pm", "By 1 pm"))
data.frame(CI_joint)
## cat X2.5. X97.5. X5.5. X94.5. X50. dose
## 1 dose_1_reduced 0.92 1.3000 0.9471971 1.2572603 1.09 dose1
## 2 dose_2_reduced 0.72 0.9995 0.7435302 0.9712552 0.85 dose2
## 3 dose_1_unprior 0.93 1.3400 0.9631583 1.2995047 1.12 dose1
## 4 dose_2_unprior 0.71 0.9954 0.7340939 0.9652361 0.84 dose2
## 5 dose_1_seroconversion_min 0.91 1.3000 0.9377508 1.2582947 1.09 dose1
## 6 dose_2_seroconversion_min 0.72 0.9994 0.7439605 0.9701603 0.85 dose2
## 7 dose_1_dichotomized_12 0.70 1.7700 0.7609404 1.6238518 1.11 dose1
## 8 dose_2_dichotomized_12 0.49 1.2740 0.5409587 1.1671746 0.79 dose2
## 9 dose_1_dichotomized_13 0.69 1.9300 0.7597255 1.7548949 1.16 dose1
## 10 dose_2_dichotomized_13 0.33 0.8980 0.3669853 0.8228976 0.55 dose2
## key
## 1 Reduced
## 2 Reduced
## 3 non-inf. pr.
## 4 non-inf. pr.
## 5 Sero+: 3.6
## 6 Sero+: 3.6
## 7 By 12 pm
## 8 By 12 pm
## 9 By 1 pm
## 10 By 1 pm
tr_d1 <- tr_joint %>%
select(matches("dose_1")) %>%
gather(key, tau) %>%
mutate(key = factor(key, levels = c("dose_1_dichotomized_13",
"dose_1_dichotomized_12",
"dose_1_seroconversion_min",
"dose_1_unprior",
"dose_1_reduced"),
labels = rev(c("Reduced",
"non-inf. pr.",
"Sero+: 3.6",
"By 12 pm",
"By 1 pm"))),
dose = 'dose1')
tr_d2 <- tr_joint %>%
select(matches("dose_2")) %>%
gather(key, tau) %>%
mutate(key = factor(key, levels = c("dose_2_dichotomized_13",
"dose_2_dichotomized_12",
"dose_2_seroconversion_min",
"dose_2_unprior",
"dose_2_reduced"),
labels = rev(c("Reduced",
"non-inf. pr.",
"Sero+: 3.6",
"By 12 pm",
"By 1 pm"))),
dose = 'dose2')
tr_d <- rbind(tr_d1, tr_d2)Printing
Open code
fig1c <- tr_d %>% ggplot(
aes(x = tau, y = key, fill = dose)) +
stat_halfeye(.width = c(0.95), slab_alpha=0.5,
linewidth = 5,
shape = 18,
point_size = 5,
normalize = "panels") +
labs(x = "Odds ratio (1-hour change or afternoon vs. morning)",
y = 'Model') +
geom_vline(xintercept=1, linetype=2,
color = "red", size=0.6) +
scale_x_continuous(limits=c(0.3, 2),
breaks=c(0.5, 1, 1.5, 2)) +
scale_fill_manual(values = cole) +
facet_wrap(~dose,
labeller = labeller(dose =
c("dose1" = "1st dose",
"dose2" = "2nd dose"))) +
theme(legend.position = "none",
strip.text = element_text(size = 14),
axis.text.y = element_text(size=11),
axis.text.x = element_text(size=11),
axis.title = element_text(size=12)) +
ggtitle("Posterior distributions, sensitivity analyses") +
geom_text(data = CI_joint,
label = paste0("95% CI: [", CI_joint$X2.5., ", ", CI_joint$X97.5., "]"),
x = 1.7, vjust = -1,
color = rep(c(cole[1], cole[2]),5))+
geom_text(data = CI_joint,
label = paste0("OR: ", CI_joint$X50),
x = 1.7, vjust = -3,
color = rep(c(cole[1], cole[2]),5))
fig1c
3 Published results
3.1 Tables
3.1.1 Table 1
Open code
table1| Characteristic | 0, N = 3411 | 1, N = 2121 | p-value2 |
|---|---|---|---|
| 1st vaccination time [hours] | 12.15 (11.15, 13.37) | 12.28 (11.34, 13.43) | 0.4 |
| 2nd vaccination time [hours] | 12.42 (11.40, 13.33) | 12.31 (11.16, 13.18) | 0.2 |
| SARS-CoV-2 IgG antibody level [AU/ml] | 0 (0, 0) | 57 (23, 171) | <0.001 |
| Time after 2nd dose [days] | 43 (28, 67) | 49 (31, 66) | 0.2 |
| Time after TX [months] | 55 (20, 115) | 100 (42, 170) | <0.001 |
| Age [years] | 67 (58, 72) | 63 (56, 71) | 0.017 |
| Male sex [0/1] | 202 (59%) | 152 (72%) | 0.003 |
| eGFR [mL/min/1.73 m2] | 44 (33, 58) | 51 (41, 68) | <0.001 |
| MMF/MPA [0/1] | 295 (87%) | 131 (62%) | <0.001 |
| Depleting therapy [0/1] | 22 (6.5%) | 2 (0.9%) | 0.002 |
| Diabetes mellitus [0/1] | 104 (30%) | 62 (29%) | 0.8 |
| BMI | 28.3 (24.8, 31.3) | 27.7 (24.9, 31.3) | 0.5 |
| Albumine [g/L] | 42.4 (40.1, 44.4) | 42.6 (40.2, 44.3) | 0.5 |
| Lymphocyte counts [10^9/L] | 1.99 (1.46, 2.47) | 2.37 (1.76, 3.13) | <0.001 |
| Vaccine | 0.2 | ||
| Moderna | 8 (2.3%) | 9 (4.2%) | |
| Phizer | 333 (98%) | 203 (96%) | |
| Tacrolimus [0/1] | 279 (82%) | 161 (76%) | 0.10 |
| CicloslorinA [0/1] | 27 (7.9%) | 29 (14%) | 0.029 |
| Calcineurin inhibitor [0/1] | 306 (90%) | 190 (90%) | >0.9 |
| Steroids [0/1] | 316 (93%) | 187 (88%) | 0.075 |
| Blood sampling time for IgG [hours] | 7.33 (6.83, 7.87) | 7.36 (6.87, 7.79) | 0.8 |
| MMF dose [mg] | 1,000 (875, 1,500) | 1,000 (500, 1,500) | 0.9 |
| Unknown | 46 | 81 | |
| Tacrolimus level [µg/l] | 6.50 (5.50, 7.60) | 6.30 (5.70, 7.40) | >0.9 |
| Unknown | 62 | 51 | |
| Ciclosporin level [µg/l] | 109 (90, 138) | 105 (95, 129) | 0.8 |
| Unknown | 314 | 183 | |
| 1 Median (IQR); n (%) | |||
| 2 Wilcoxon rank sum test; Pearson’s Chi-squared test | |||
3.1.2 Table 2
Open code
table2| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| 1st vaccination time [hour] | 1.110 | 0.925 | 1.329 | 0.8661 |
| 2nd vaccination time [hour] | 0.843 | 0.713 | 0.998 | 0.9772 |
| Moderna [0/1] | 1.386 | 0.398 | 4.802 | 0.7033 |
| Calcineurin inhibitor[0/1] | 1.789 | 0.786 | 4.111 | 0.9225 |
| Steroids [0/1] | 0.728 | 0.346 | 1.541 | 0.7910 |
| Time after TX [10 years] | 2.117 | 1.464 | 3.101 | 1.0000 |
| Diabetes Mellitus [0/1] | 0.911 | 0.556 | 1.498 | 0.6451 |
| Age [10 years] | 0.722 | 0.553 | 0.936 | 0.9925 |
| BMI [5 units] | 0.968 | 0.770 | 1.223 | 0.6114 |
| Albumine [g/L] | 1.034 | 0.973 | 1.101 | 0.8549 |
| Lymphocyte counts [log(10^9/L)] | 2.533 | 1.495 | 4.350 | 0.9998 |
| Male sex [0/1] | 2.015 | 1.291 | 3.201 | 0.9996 |
| eGFR [10 mL/min/1.73 m2] | 1.378 | 1.230 | 1.550 | 1.0000 |
| MMF/MPA [0/1] | 0.120 | 0.069 | 0.205 | 1.0000 |
| Depleting therapy [0/1] | 0.237 | 0.035 | 1.137 | 0.9608 |
3.1.3 Suppl. Table 1
Open code
suppl_table1| df | AIC | BIC | |
|---|---|---|---|
| model_all | 17.7 | 604 | 680.5 |
| model_all_int | 18.7 | 606 | 686.9 |
3.1.4 Suppl. Table 2
Open code
suppl_table2| Characteristic | morning-morning, N = 3011 | afternoon-morning, N = 721 | morning-afternoon, N = 791 | afternoon-afternoon, N = 1011 | p-value2 |
|---|---|---|---|---|---|
| 1st vaccination time [hours] | 11.63 (10.78, 12.20) | 13.74 (13.43, 14.08) | 11.73 (11.05, 12.38) | 13.95 (13.47, 14.33) | <0.001 |
| 2nd vaccination time [hours] | 11.67 (10.63, 12.35) | 11.97 (11.21, 12.49) | 13.48 (13.23, 13.98) | 13.78 (13.40, 14.12) | <0.001 |
| SARS-CoV-2 IgG antibody level [AU/ml] | 0 (0, 36) | 6 (0, 40) | 0 (0, 13) | 0 (0, 31) | 0.2 |
| Time after 2nd dose [days] | 47 (28, 65) | 40 (26, 62) | 69 (51, 92) | 36 (27, 55) | <0.001 |
| Time after TX [months] | 64 (21, 127) | 81 (27, 150) | 103 (47, 182) | 73 (28, 131) | 0.003 |
| Age [years] | 66 (57, 71) | 62 (55, 69) | 73 (71, 76) | 62 (54, 67) | <0.001 |
| Male sex [0/1] | 197 (65%) | 47 (65%) | 46 (58%) | 64 (63%) | 0.7 |
| eGFR [mL/min/1.73 m2] | 47 (36, 61) | 44 (33, 59) | 45 (32, 57) | 52 (38, 67) | 0.087 |
| MMF/MPA [0/1] | 233 (77%) | 53 (74%) | 56 (71%) | 84 (83%) | 0.2 |
| Depleting therapy [0/1] | 16 (5.3%) | 3 (4.2%) | 0 (0%) | 5 (5.0%) | 0.2 |
| Diabetes mellitus [0/1] | 92 (31%) | 19 (26%) | 31 (39%) | 24 (24%) | 0.13 |
| BMI | 28.4 (25.0, 31.8) | 27.4 (23.7, 30.8) | 28.1 (25.3, 30.8) | 27.8 (24.5, 30.8) | 0.2 |
| Albumine [g/L] | 42.4 (40.0, 44.3) | 42.9 (40.4, 44.6) | 41.3 (40.0, 43.6) | 43.0 (40.9, 44.6) | 0.038 |
| Lymphocyte counts [10^9/L] | 2.08 (1.56, 2.69) | 2.03 (1.63, 2.71) | 2.18 (1.57, 2.84) | 2.28 (1.72, 2.82) | 0.5 |
| Seroconversion [0/1] | 120 (40%) | 32 (44%) | 25 (32%) | 35 (35%) | 0.3 |
| Vaccine | 0.12 | ||||
| Moderna | 7 (2.3%) | 2 (2.8%) | 1 (1.3%) | 7 (6.9%) | |
| Phizer | 294 (98%) | 70 (97%) | 78 (99%) | 94 (93%) | |
| Tacrolimus [0/1] | 245 (81%) | 59 (82%) | 52 (66%) | 84 (83%) | 0.013 |
| CicloslorinA [0/1] | 28 (9.3%) | 5 (6.9%) | 19 (24%) | 4 (4.0%) | <0.001 |
| Calcineurin inhibitor [0/1] | 273 (91%) | 64 (89%) | 71 (90%) | 88 (87%) | 0.8 |
| Steroids [0/1] | 271 (90%) | 69 (96%) | 70 (89%) | 93 (92%) | 0.4 |
| Blood sampling time for IgG [hours] | 7.28 (6.78, 7.80) | 7.53 (6.92, 7.95) | 7.52 (7.07, 7.79) | 7.33 (6.88, 7.93) | 0.12 |
| MMF dose [mg] | 1,000 (500, 1,500) | 1,000 (500, 1,500) | 1,000 (875, 1,000) | 1,000 (1,000, 1,500) | 0.5 |
| Unknown | 68 | 19 | 23 | 17 | |
| Tacrolimus level [µg/l] | 6.50 (5.60, 7.30) | 6.40 (5.65, 7.65) | 6.30 (5.33, 7.90) | 6.50 (5.68, 7.73) | >0.9 |
| Unknown | 56 | 13 | 27 | 17 | |
| Ciclosporin level [µg/l] | 110 (91, 138) | 105 (85, 106) | 106 (92, 135) | 135 (114, 155) | 0.3 |
| Unknown | 273 | 67 | 60 | 97 | |
| 1 Median (IQR); n (%) | |||||
| 2 Kruskal-Wallis rank sum test; Pearson’s Chi-squared test; Fisher’s exact test | |||||
3.1.5 Suppl. Table 3
Open code
suppl_table3| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| 1st vaccination time [hour] | 1.089 | 0.919 | 1.300 | 0.8339 |
| 2nd vaccination time [hour] | 0.849 | 0.723 | 1.000 | 0.9752 |
| Time after TX [10 years] | 2.120 | 1.517 | 2.971 | 1.0000 |
| Age [10 years] | 0.753 | 0.602 | 0.941 | 0.9944 |
| Lymphocyte counts [log(10^9/L)] | 2.516 | 1.510 | 4.264 | 0.9998 |
| Male sex [0/1] | 1.970 | 1.268 | 3.093 | 0.9986 |
| eGFR [10 mL/min/1.73 m2] | 1.398 | 1.254 | 1.567 | 1.0000 |
| MMF/MPA [0/1] | 0.146 | 0.086 | 0.243 | 1.0000 |
| Depleting therapy [0/1] | 0.217 | 0.033 | 0.979 | 0.9773 |
3.1.6 Suppl. Table 4
Open code
suppl_table4| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| 1st vaccination time [hour] | 1.117 | 0.932 | 1.344 | 0.8831 |
| 2nd vaccination time [hour] | 0.842 | 0.709 | 0.995 | 0.9781 |
| Moderna [0/1] | 1.395 | 0.414 | 4.714 | 0.7036 |
| Calcineurin inhibitor[0/1] | 1.782 | 0.801 | 4.078 | 0.9212 |
| Steroids [0/1] | 0.727 | 0.339 | 1.545 | 0.7962 |
| Time after TX [10 years] | 2.109 | 1.454 | 3.068 | 1.0000 |
| Diabetes Mellitus [0/1] | 0.909 | 0.558 | 1.480 | 0.6517 |
| Age [10 years] | 0.724 | 0.559 | 0.942 | 0.9933 |
| BMI [5 units] | 0.970 | 0.769 | 1.221 | 0.6039 |
| Albumine [g/L] | 1.033 | 0.973 | 1.099 | 0.8538 |
| Lymphocyte counts [log(10^9/L)] | 2.531 | 1.480 | 4.402 | 0.9996 |
| Male sex [0/1] | 2.010 | 1.294 | 3.169 | 0.9991 |
| eGFR [10 mL/min/1.73 m2] | 1.377 | 1.229 | 1.551 | 1.0000 |
| MMF/MPA [0/1] | 0.121 | 0.070 | 0.206 | 1.0000 |
| Depleting therapy [0/1] | 0.235 | 0.033 | 1.149 | 0.9606 |
3.1.7 Suppl. Table 5
Open code
suppl_table5| OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|
| 1st vaccination time [hour] | 1.086 | 0.907 | 1.299 | 0.8165 |
| 2nd vaccination time [hour] | 0.850 | 0.723 | 0.999 | 0.9757 |
| Moderna [0/1] | 2.417 | 0.643 | 9.665 | 0.9058 |
| Calcineurin inhibitor[0/1] | 1.546 | 0.737 | 3.320 | 0.8718 |
| Steroids [0/1] | 0.461 | 0.213 | 0.966 | 0.9803 |
| Time after TX [10 years] | 1.810 | 1.262 | 2.631 | 0.9996 |
| Diabetes Mellitus [0/1] | 0.852 | 0.530 | 1.367 | 0.7454 |
| Age [10 years] | 0.718 | 0.550 | 0.930 | 0.9939 |
| BMI [5 units] | 0.965 | 0.772 | 1.199 | 0.6203 |
| Albumine [g/L] | 1.036 | 0.977 | 1.101 | 0.8825 |
| Lymphocyte counts [log(10^9/L)] | 3.197 | 1.909 | 5.439 | 1.0000 |
| Male sex [0/1] | 2.429 | 1.580 | 3.786 | 0.9998 |
| eGFR [10 mL/min/1.73 m2] | 1.288 | 1.152 | 1.444 | 1.0000 |
| MMF/MPA [0/1] | 0.133 | 0.074 | 0.228 | 1.0000 |
| Depleting therapy [0/1] | 0.284 | 0.057 | 1.138 | 0.9623 |
3.1.8 Suppl. Table 6
Open code
suppl_table6| OR | Q2.5 | Q97.5 | PD | OR | Q2.5 | Q97.5 | PD | |
|---|---|---|---|---|---|---|---|---|
| 1st dose afternoon[0/1] | 1.110 | 0.696 | 1.768 | 0.6705 | 1.158 | 0.689 | 1.928 | 0.7152 |
| 2nd dose afternoon [0/1] | 0.793 | 0.493 | 1.274 | 0.8290 | 0.550 | 0.333 | 0.898 | 0.9922 |
| Moderna [0/1] | 1.272 | 0.365 | 4.393 | 0.6509 | 1.393 | 0.395 | 4.855 | 0.7047 |
| Calcineurin inhibitor[0/1] | 1.730 | 0.777 | 4.008 | 0.9050 | 1.774 | 0.793 | 4.080 | 0.9152 |
| Steroids [0/1] | 0.714 | 0.334 | 1.526 | 0.8122 | 0.718 | 0.337 | 1.545 | 0.8105 |
| Time after TX [10 years] | 2.079 | 1.444 | 3.037 | 1.0000 | 2.152 | 1.466 | 3.169 | 1.0000 |
| Diabetes Mellitus [0/1] | 0.924 | 0.567 | 1.494 | 0.6289 | 0.929 | 0.573 | 1.521 | 0.6148 |
| Age [10 years] | 0.693 | 0.532 | 0.899 | 0.9971 | 0.688 | 0.534 | 0.888 | 0.9982 |
| BMI [5 units] | 0.965 | 0.765 | 1.214 | 0.6181 | 0.971 | 0.775 | 1.216 | 0.5973 |
| Albumine [g/L] | 1.035 | 0.975 | 1.102 | 0.8661 | 1.034 | 0.972 | 1.102 | 0.8598 |
| Lymphocyte counts [log(10^9/L)] | 2.473 | 1.427 | 4.306 | 0.9997 | 2.478 | 1.452 | 4.278 | 0.9998 |
| Male sex [0/1] | 1.995 | 1.271 | 3.179 | 0.9987 | 1.956 | 1.256 | 3.066 | 0.9989 |
| eGFR [10 mL/min/1.73 m2] | 1.373 | 1.226 | 1.543 | 1.0000 | 1.386 | 1.236 | 1.557 | 1.0000 |
| MMF/MPA [0/1] | 0.122 | 0.070 | 0.209 | 1.0000 | 0.119 | 0.066 | 0.204 | 1.0000 |
| Depleting therapy [0/1] | 0.248 | 0.038 | 1.151 | 0.9592 | 0.253 | 0.037 | 1.169 | 0.9572 |
3.1.9 Suppl. Table 7
Open code
suppl_table7| estimate | Q2.5 | Q97.5 | |
|---|---|---|---|
| 1st vaccination time [hour] | -0.04 | -0.12 | 0.07 |
| 2nd vaccination time [hour] | -0.07 | -0.17 | -0.01 |
| Moderna [0/1] | 1.22 | -0.72 | 2.81 |
| Calcineurin inhibitor[0/1] | -0.04 | -0.26 | 0.48 |
| Steroids [0/1] | -0.37 | -1.31 | -0.07 |
| Time after TX [10 years] | 0.73 | 0.52 | 1.07 |
| Diabetes Mellitus [0/1] | 0.17 | -0.14 | 0.43 |
| Age [10 years] | -0.35 | -0.58 | -0.15 |
| BMI [5 units] | -0.05 | -0.17 | 0.07 |
| Albumine [g/L] | 0.06 | 0.01 | 0.09 |
| Lymphocyte counts [log(10^9/L)] | 0.46 | 0.11 | 0.91 |
| Male sex [0/1] | 0.66 | 0.43 | 0.91 |
| eGFR [10 mL/min/1.73 m2] | 0.25 | 0.17 | 0.33 |
| MMF/MPA [0/1] | -2.46 | -2.80 | -1.97 |
| Depleting therapy [0/1] | -0.20 | -0.47 | 0.54 |
| Weeks after 2nd dose | 0.02 | -0.01 | 0.06 |
3.2 Figures
3.2.1 Figure 1
Open code
bot <- plot_grid(fig1a, fig1b, rel_widths = c(0.5, 0.5), labels=c("A", "B"), ncol=2)
plot_grid(bot, fig1c, nrow=2, rel_heights = c(0.45, 0.55), labels=c("", "C"))
3.2.2 Figure 2
Open code
fig2
3.2.3 Suppl. Figure 1
Open code
par(mfrow=c(1,3))
par(mar=c(3,3,1,1))
par(mgp=c(1.5,0.5,0))
plot.gam(model_all_nl, ylab=c('Partial residuals'))
3.2.4 Suppl. Figure 2
Open code
suppl_fig2
3.2.5 Suppl. Figure 3
Open code
suppl_fig3
4 Reproducibility
Open code
sessionInfo()
## R version 4.3.2 (2023-10-31)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 22.04.4 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.10.0
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.10.0
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=cs_CZ.UTF-8 LC_COLLATE=en_US.UTF-8
## [5] LC_MONETARY=cs_CZ.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=cs_CZ.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=cs_CZ.UTF-8 LC_IDENTIFICATION=C
##
## time zone: Europe/Prague
## tzcode source: system (glibc)
##
## attached base packages:
## [1] stats graphics grDevices utils datasets methods base
##
## other attached packages:
## [1] ggbeeswarm_0.7.2 quantreg_5.88 SparseM_1.81 coxed_0.3.3
## [5] survival_3.5-7 rms_6.7-1 Hmisc_5.1-1 bayesplot_1.8.1
## [9] ggdist_3.3.0 kableExtra_1.3.4 lubridate_1.8.0 corrplot_0.92
## [13] arm_1.13-1 lme4_1.1-34 MASS_7.3-60 janitor_2.2.0
## [17] RJDBC_0.2-10 rJava_1.0-6 DBI_1.1.2 projpred_2.6.0
## [21] brms_2.16.3 Rcpp_1.0.11 glmnet_4.1-8 Matrix_1.6-3
## [25] boot_1.3-28 cowplot_1.1.1 pROC_1.18.0 mgcv_1.9-1
## [29] nlme_3.1-163 openxlsx_4.2.5 flextable_0.9.3 sjPlot_2.8.15
## [33] car_3.1-2 carData_3.0-5 gtsummary_1.7.2 emmeans_1.7.2
## [37] ggpubr_0.4.0 forcats_1.0.0 stringr_1.5.0 dplyr_1.1.3
## [41] purrr_1.0.2 readr_2.1.2 tidyr_1.2.0 tibble_3.2.1
## [45] ggplot2_3.4.3 tidyverse_1.3.1
##
## loaded via a namespace (and not attached):
## [1] fs_1.6.3 matrixStats_1.0.0 httr_1.4.2
## [4] webshot_0.5.5 insight_0.19.5 tools_4.3.2
## [7] backports_1.4.1 sjlabelled_1.2.0 utf8_1.2.3
## [10] R6_2.5.1 DT_0.17 withr_2.5.0
## [13] Brobdingnag_1.2-7 prettyunits_1.1.1 gridExtra_2.3
## [16] cli_3.6.1 textshaping_0.3.6 performance_0.10.5
## [19] gt_0.9.0 shinyjs_2.1.0 officer_0.6.2
## [22] sandwich_3.0-1 labeling_0.4.2 sass_0.4.7
## [25] mvtnorm_1.1-3 polspline_1.1.24 ggridges_0.5.3
## [28] askpass_1.1 commonmark_1.9.0 systemfonts_1.0.4
## [31] StanHeaders_2.21.0-7 foreign_0.8-86 gfonts_0.2.0
## [34] svglite_2.1.0 readxl_1.3.1 rstudioapi_0.13
## [37] httpcode_0.3.0 generics_0.1.2 shape_1.4.6
## [40] gtools_3.9.2 crosstalk_1.2.0 distributional_0.3.2
## [43] zip_2.2.0 inline_0.3.19 loo_2.4.1
## [46] fansi_1.0.4 abind_1.4-5 lifecycle_1.0.3
## [49] multcomp_1.4-18 yaml_2.3.5 snakecase_0.11.1
## [52] grid_4.3.2 promises_1.2.0.1 crayon_1.5.2
## [55] miniUI_0.1.1.1 lattice_0.22-5 haven_2.4.3
## [58] pillar_1.9.0 knitr_1.37 estimability_1.3
## [61] shinystan_2.5.0 codetools_0.2-19 glue_1.6.2
## [64] fontLiberation_0.1.0 data.table_1.14.8 broom.helpers_1.14.0
## [67] vctrs_0.6.3 cellranger_1.1.0 gtable_0.3.4
## [70] datawizard_0.9.0 assertthat_0.2.1 xfun_0.40
## [73] mime_0.12 coda_0.19-4 rsconnect_0.8.25
## [76] iterators_1.0.14 shinythemes_1.2.0 ellipsis_0.3.2
## [79] TH.data_1.1-0 xts_0.12.1 fontquiver_0.2.1
## [82] threejs_0.3.3 rstan_2.21.3 tensorA_0.36.2
## [85] vipor_0.4.5 rpart_4.1.23 colorspace_2.1-0
## [88] nnet_7.3-19 tidyselect_1.2.0 processx_3.8.2
## [91] compiler_4.3.2 curl_4.3.2 rvest_1.0.2
## [94] htmlTable_2.4.0 xml2_1.3.3 fontBitstreamVera_0.1.1
## [97] bayestestR_0.13.1 colourpicker_1.1.1 posterior_1.2.0
## [100] checkmate_2.0.0 scales_1.2.1 dygraphs_1.1.1.6
## [103] quadprog_1.5-8 callr_3.7.3 digest_0.6.29
## [106] minqa_1.2.4 rmarkdown_2.25 htmltools_0.5.6
## [109] pkgconfig_2.0.3 base64enc_0.1-3 highr_0.9
## [112] dbplyr_2.1.1 fastmap_1.1.1 rlang_1.1.1
## [115] htmlwidgets_1.6.2 shiny_1.8.0 farver_2.1.1
## [118] zoo_1.8-9 jsonlite_1.7.3 magrittr_2.0.3
## [121] Formula_1.2-4 parameters_0.21.2 munsell_0.5.0
## [124] gdtools_0.3.3 stringi_1.7.12 plyr_1.8.8
## [127] pkgbuild_1.3.1 parallel_4.3.2 sjmisc_2.8.9
## [130] ggeffects_1.3.1 splines_4.3.2 hms_1.1.1
## [133] sjstats_0.18.2 ps_1.6.0 igraph_1.5.1
## [136] uuid_1.0-3 ggsignif_0.6.3 markdown_1.8
## [139] effectsize_0.8.6 reshape2_1.4.4 stats4_4.3.2
## [142] rstantools_2.1.1 crul_1.4.0 reprex_2.0.1
## [145] evaluate_0.15 RcppParallel_5.1.7 modelr_0.1.8
## [148] nloptr_2.0.0 tzdb_0.4.0 foreach_1.5.2
## [151] httpuv_1.6.5 MatrixModels_0.5-0 openssl_1.4.6
## [154] broom_1.0.5 xtable_1.8-4 rstatix_0.7.0
## [157] later_1.3.0 viridisLite_0.4.2 ragg_1.2.1
## [160] beeswarm_0.4.0 cluster_2.1.6 gamm4_0.2-6
## [163] bridgesampling_1.1-2References
Bürkner, Paul-Christian. 2017. “Brms: An r Package for Bayesian Multilevel Models Using Stan” 80. https://doi.org/10.18637/jss.v080.i01.
Gelman, Andrew. 2008. “Scaling Regression Inputs by Dividing by Two Standard Deviations.” Statistics in Medicine 27 (15): 2865–73. https://doi.org/10.1002/sim.3107.
Gelman, Andrew, and Yu-Sung Su. 2022. “Arm: Data Analysis Using Regression and Multilevel/Hierarchical Models.” https://CRAN.R-project.org/package=arm.
Lin, Ting-Yun, and Szu-Chun Hung. 2023. “Morning Vaccination Compared with Afternoon/Evening Vaccination on Humoral Response to SARS-CoV-2 Adenovirus Vector-Based Vaccine in Hemodialysis Patients.” Clinical Journal of the American Society of Nephrology, May. https://doi.org/10.2215/cjn.0000000000000207.
Magicova, Maria, Ivan Zahradka, Martina Fialova, Tomas Neskudla, Jiri Gurka, Istvan Modos, Michal Hojny, et al. 2022. “Determinants of Immune Response to AntiSARS-CoV-2 mRNA Vaccines in Kidney Transplant Recipients: A Prospective Cohort Study.” Transplantation 106 (4): 842–52. https://doi.org/10.1097/tp.0000000000004044.
Makowski, Dominique, Mattan S. Ben-Shachar, S. H. Annabel Chen, and Daniel Lüdecke. 2019. “Indices of Effect Existence and Significance in the Bayesian Framework.” Frontiers in Psychology 10 (December). https://doi.org/10.3389/fpsyg.2019.02767.
R Core Team. 2023. R: A Language and Environment for Statistical Computing. Vienna, Austria: R Foundation for Statistical Computing. https://www.R-project.org/.
Stan Development Team. 2021. “RStan: The R Interface to Stan.” https://mc-stan.org/.
Wood, S. N. 2011. “Fast Stable Restricted Maximum Likelihood and Marginal Likelihood Estimation of Semiparametric Generalized Linear Models” 73: 3–36.
Zahradka, Ivan, Filip Tichanek, Maria Magicova, Istvan Modos, Ondrej Viklicky, and Vojtech Petr. 2024. “Morning Administration Enhances Humoral Response to SARS-CoV-2 Vaccination in Kidney Transplant Recipients.” American Journal of Transplantation. https://doi.org/10.1016/j.ajt.2024.03.004.