How long someone has detectable neutralizing antibodies after vaccination is important in understanding the impact of vaccination on disease transmission. In this post I step through several different models and re-examine data from a prior paper.
A prominant scientist has been posting his highlights of pre-prints on twitter during the course of the pandemic and occasionally, I come across one that makes me scratch my head. I think that the move to open and rapid science through pre-prints has been good, but not without issues. Review (and open review especially) makes work better and is an important part of science and progress.
One paper that stood out to me was by Suthar et al. (2021) in which the authors look at the neutralizing antibody tires (nAbs) amongst T-cell response in those persons who were vaccinated with the Pfizer SARS-CoV-2 Vaccines (BNT162b2 vaccine). The authors found: “substantial waning of antibody response and T-cell Immunity” in those who were vaccinated. Naturally, this could be concerning because antibodies and T-cells are important in fighting infection and could help us understand the endgame of this first epidemic of SARS-CoV-2 (e.g., vaccine induced reduction of transmission, symptomatic disease, severe outcomes, and associated periodicity given time since last vaccine). However, there is the open question of what is the best correlate of protection (e.g., we know that antibodies wane with time, but the ability to produce antibodies in the face of a new infection from memory B-cells is also very important) (Corbett et al. 2021)
So Look at the Data
As mentioned, along with the bad things about pre-prints (media jumping on sensationalist findings, people seeking acclaim by making bold statements that are not rigorously founded in fact, general review help, etc) is that much of the data are available for secondary analysis. With the article in question, I could find some of the original data analyzed by the authors from a prior publication. We can start to look at these data ourselves:
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
x dplyr::between() masks data.table::between()
x dplyr::filter() masks stats::filter()
x dplyr::first() masks data.table::first()
x dplyr::lag() masks stats::lag()
x dplyr::last() masks data.table::last()
x purrr::transpose() masks data.table::transpose()
One things that stands out to me is participant 2015. This individual has nAbs that are much higher than the other participants.
dat_long %>%mutate(tst =ifelse(PID==2015,"PID: 2015","All Others")) %>%group_by(tst, day) %>%summarise(mu =mean(value, na.rm=TRUE)) %>%spread(day,mu) %>% knitr::kable(digits =1, caption ="Mean Value of nAbs by Group")
`summarise()` has grouped output by 'tst'. You can override using the `.groups`
argument.
Mean Value of nAbs by Group
tst
0
21
42
90
All Others
11.5
104.0
499.3
251.1
PID: 2015
10.0
2876.5
4379.5
2492.0
If this were me, I would go back and check those data just in case someone got a little too happy when titrating (been there). Regardless, accepting this level of variability in the individual, these kind of data lend themselves well to hierarchical modeling.
So What’s the Half Life?
We can model the half life of the nAbs in the participants using the multilevel model framework with lme4. I will subset the data to include only the levels after the peak antibody level and include a random effect for the participant (since we know that there is within-person variability and serial dependence).
library(lme4)
Loading required package: Matrix
Attaching package: 'Matrix'
The following objects are masked from 'package:tidyr':
expand, pack, unpack
fit_lme <-lmer(log(value) ~ day + (1|PID), data =subset(dat_long, day >21))(half_life <--log(2)/confint(fit_lme)[4,])
Computing profile confidence intervals ...
2.5 % 97.5 %
46.31701 84.74999
So based on this model, the half-life for nAbs is between 46.3 and 84.7 days.
nab_out <-function(t,hl){2^(-t/hl)}sim_dat <-tibble(day =1:365) %>%mutate(lower_bound =nab_out(day, half_life[1]),upper_bound =nab_out(day, half_life[2]))sim_dat %>%ggplot(aes(x = day))+geom_ribbon(aes(ymin = lower_bound,ymax = upper_bound), fill ="orange", alpha = .5)+geom_hline(yintercept = .25, col ="red", lty ="dashed")+labs(title ="Simulated Decay of nAbs",y ="Percent of Maximum nAbs",x ="Day After Peak nAbs" )
So it appears that in roughly 3-6 months after maximum antibody levels a given vaccinated person would have 1/4 of nAb titres.
The big question of course is first, does this matter and in what context (protection vs infection where high levels of nAbs make sense vs symptomatic disease where a trained immune system is all that is required (thinking memory B and T cells)) and what concentration is enough at each level.
Is 5% of peak is enough then we have between 6 months and a year, perhaps in time for a teaser from infection to boost nAbs.
round(days_until_perc(.05, half_life)/30,1)
2.5 % 97.5 %
6.7 12.2
What about compartmental models?
This is not my forte, but it would seem that compartmental models would be a good way to approach this problem. We would have to fix some parameters due to the data available but one could imagine that a decent data generating process could look like (at least as far as I understand how mRNA vaccines work):
Inoculum is added to the blood stream via vaccination
A latent phase for vaccinated individual’s cells to produce the spike protein which would be a function of the volume of inoculum and the number of cells encountered. The inoculum also has some decay function.
A phase where the encoded protein is in the bloodstream and recognized by the immune system as a foreign body which ramps up the innate and humoral immune response.
Eventually the spike proteins would be cleared at some rate as the immune response increases
I started down this path and then realized that I would need to make some assumptions for some of the parameters (transition rates between compartments). Being that I do not have reasonable priors for these rates, I decided not to code up this mdoel.
Other Method
I read another paper on waning immunity on coronaviruses that had much more longitudinal data (Townsend et al. 2021). The key model in the paper was locked away in a proprietary Mathematica Notebook that couldn’t open without buying the software. I coded the below model based on their methods, but we will find is that we do not have enough data to fit a very good model.
rw("decay.stan")
data {
int<lower=1> N; // Number of Samples
int<lower=1> J; //Number of participants
vector [N] nAbs; //The concentration of the antibodies
vector [N] time; //The concentration of the antibodies
}
parameters{
real omega;
real hl;
real<lower=0> sigma;
}
transformed parameters {
vector[N] theta; //vac effect
theta[1:N] = omega + (1- omega) * exp(-hl*time[1:N]);
}
model {
hl ~ normal(0,2);
omega ~ normal(0,1);
nAbs ~ normal(theta, sigma);
}
library(cmdstanr)
This is cmdstanr version 0.4.0
- Online documentation and vignettes at mc-stan.org/cmdstanr
- CmdStan path set to: /Users/michael/.cmdstanr/cmdstan-2.28.2
- Use set_cmdstan_path() to change the path
A newer version of CmdStan is available. See ?install_cmdstan() to install it.
To disable this check set option or environment variable CMDSTANR_NO_VER_CHECK=TRUE.
Chain 3 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 3 Exception: normal_lpdf: Location parameter[55] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 3 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 3 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 3
Chain 2 finished in 0.2 seconds.
Chain 4 finished in 0.1 seconds.
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[55] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: normal_lpdf: Location parameter[1] is -inf, but must be finite! (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/Rtmp8xeEVQ/model-4cb07f1e8483.stan', line 22, column 1 to column 29)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 3 finished in 0.3 seconds.
Chain 1 finished in 42.3 seconds.
All 4 chains finished successfully.
Mean chain execution time: 10.7 seconds.
Total execution time: 42.7 seconds.
Warning in seq.default(from = 1, len = along - 1): partial argument match of
'len' to 'length.out'
Warning in seq.default(to = N - 1, len = N - along): partial argument match of
'len' to 'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.list)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = ncol(arg.dim)): partial argument match of 'len' to
'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.names)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning: 232 of 2000 (12.0%) transitions ended with a divergence.
This may indicate insufficient exploration of the posterior distribution.
Possible remedies include:
* Increasing adapt_delta closer to 1 (default is 0.8)
* Reparameterizing the model (e.g. using a non-centered parameterization)
* Using informative or weakly informative prior distributions
187 of 2000 (9.0%) transitions hit the maximum treedepth limit of 15 or 2^15-1 leapfrog steps.
Trajectories that are prematurely terminated due to this limit will result in slow exploration.
Increasing the max_treedepth limit can avoid this at the expense of more computation.
If increasing max_treedepth does not remove warnings, try to reparameterize the model.
Not a great fit as we are trying to draw a complicated curve through basically two points and there are (infinitely) many ways to do this.
Warning in seq.default(from = 1, len = along - 1): partial argument match of
'len' to 'length.out'
Warning in seq.default(to = N - 1, len = N - along): partial argument match of
'len' to 'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.list)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = ncol(arg.dim)): partial argument match of 'len' to
'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.names)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
A model that is somewhat of a compromise between the compartmental model which captures more of the underlying biological process and the hierarchical model which fit the within-host variability is to model the viral decay with a direct functional form. This is similar to the approach performed by Singanayagam et al. (2021) who use this functional form for estimated RNA copies from Cycle Threshold values. The issues, as before, is that we only have a few observations per person and I don’t have strong priors, but this model is closer to what we see in the graphs where the titre increases to some point and then decay afterwards. Ideally we would have some additional measures for the participants.
Warning in readLines(x): incomplete final line found on 'rise.stan'
functions{
real viral_kinetics (real v_max, real a, real b, real t, real t_max){
return(v_max * (a + b) / (b * exp(-a*(t - t_max)) + a * exp(b*(t-t_max))));
}
}
data {
int<lower=1> N; // Number of records
int<lower=0> id[N]; //Participant IDs
vector<lower=0> [N] nAbs_max; // Maximum Ab concentration
vector<lower=0> [N] nAbs; //The concentration of the antibodies
vector [N] time; //The concentration of the antibodies
vector [N] t_max; //The concentration of the antibodies
real prior_a_mu;
real prior_b_mu;
real prior_a_sd;
real prior_b_sd;
int pred_window;
}
parameters {
real<lower=0> a;
real<lower=0> b;
real<lower=0> sigma;
}
transformed parameters{
vector<lower=0>[N] mu;
for( i in 1:N){
mu[i] = viral_kinetics(nAbs_max[id[i]], a, b, time[id[i]], t_max[id[i]]);
}
}
model {
a ~ normal(prior_a_mu, prior_a_sd);
b ~ normal(prior_b_mu, prior_b_sd);
sigma ~ normal(0,1);
nAbs ~ normal(mu, sigma);
}
generated quantities{
vector[pred_window] yhat;
vector[pred_window] mu_pred;
real avg_nabs_max = mean(nAbs_max);
real avg_tmax = mean(t_max);
for(i in 1:pred_window){
mu_pred[i] = viral_kinetics(avg_nabs_max, a, b, i, avg_tmax);
yhat[i] = normal_rng(mu_pred[i], sigma);
}
}
Chain 1 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 1 Exception: rise_model_namespace::log_prob: mu[sym1__] is nan, but must be greater than or equal to 0 (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/RtmpJ4OPS5/model-dbcf188d05cc.stan', line 28, column 1 to column 23)
Chain 1 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 1 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 1
Chain 2 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 2 Exception: rise_model_namespace::log_prob: mu[sym1__] is nan, but must be greater than or equal to 0 (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/RtmpJ4OPS5/model-dbcf188d05cc.stan', line 28, column 1 to column 23)
Chain 2 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 2 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 2
Chain 4 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 4 Exception: rise_model_namespace::log_prob: mu[sym1__] is nan, but must be greater than or equal to 0 (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/RtmpJ4OPS5/model-dbcf188d05cc.stan', line 28, column 1 to column 23)
Chain 4 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 4 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 4
Chain 4 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 4 Exception: rise_model_namespace::log_prob: mu[sym1__] is nan, but must be greater than or equal to 0 (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/RtmpJ4OPS5/model-dbcf188d05cc.stan', line 28, column 1 to column 23)
Chain 4 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 4 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 4
Chain 4 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 4 Exception: rise_model_namespace::log_prob: mu[sym1__] is nan, but must be greater than or equal to 0 (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/RtmpJ4OPS5/model-dbcf188d05cc.stan', line 28, column 1 to column 23)
Chain 4 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 4 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 4
Chain 4 Informational Message: The current Metropolis proposal is about to be rejected because of the following issue:
Chain 4 Exception: rise_model_namespace::log_prob: mu[sym1__] is nan, but must be greater than or equal to 0 (in '/var/folders/0x/6bnjy4n15kz8nbfbbwbyk3pr0000gn/T/RtmpJ4OPS5/model-dbcf188d05cc.stan', line 28, column 1 to column 23)
Chain 4 If this warning occurs sporadically, such as for highly constrained variable types like covariance matrices, then the sampler is fine,
Chain 4 but if this warning occurs often then your model may be either severely ill-conditioned or misspecified.
Chain 4
Chain 1 finished in 2.2 seconds.
Chain 2 finished in 2.3 seconds.
Chain 3 finished in 2.3 seconds.
Chain 4 finished in 2.2 seconds.
All 4 chains finished successfully.
Mean chain execution time: 2.2 seconds.
Total execution time: 2.6 seconds.
Warning in seq.default(from = 1, len = along - 1): partial argument match of
'len' to 'length.out'
Warning in seq.default(to = N - 1, len = N - along): partial argument match of
'len' to 'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.list)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = ncol(arg.dim)): partial argument match of 'len' to
'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.names)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(from = 1, len = along - 1): partial argument match of
'len' to 'length.out'
Warning in seq.default(to = N - 1, len = N - along): partial argument match of
'len' to 'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.list)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = ncol(arg.dim)): partial argument match of 'len' to
'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.names)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Note that there is a fair bit of variance in our model represented by zero. Additionally, we have some values of the predicted titre below zero. We could resolve this by changing the distribution to use a distribution that is always positive like a lognormal distribution. This would require some careful transformation of our data and associated parameters. We’ll stick with this formulation for now.
We can now visualize our results from this model.
library(posterior)
This is posterior version 1.2.1
Attaching package: 'posterior'
The following objects are masked from 'package:stats':
mad, sd, var
Warning in seq.default(from = 1, len = along - 1): partial argument match of
'len' to 'length.out'
Warning in seq.default(to = N - 1, len = N - along): partial argument match of
'len' to 'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.list)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = ncol(arg.dim)): partial argument match of 'len' to
'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.names)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(from = 1, len = along - 1): partial argument match of
'len' to 'length.out'
Warning in seq.default(to = N - 1, len = N - along): partial argument match of
'len' to 'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.list)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = ncol(arg.dim)): partial argument match of 'len' to
'length.out'
Warning in seq.default(len = N): partial argument match of 'len' to 'length.out'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(along = arg.names): partial argument match of 'along' to
'along.with'
Warning in seq.default(len = length(arg.names)): partial argument match of 'len'
to 'length.out'
Warning in seq.default(along = perm): partial argument match of 'along' to
'along.with'
pred_nab %>%ggplot(aes(i, .value))+geom_lineribbon(aes(ymin = .lower, ymax = .upper))+scale_fill_brewer()+labs(title ="Estimated Quantity of nAbs with Viral Evolution Model",y ="Quantity",x ="Days",fill ="Proportion of Scenarios" )+geom_hline(yintercept =0, col ="orange", lty ="dashed")+theme(legend.position ="top")
Again, we see evidence of detectable tire above 180 days, but there is a decent bit of variability.
Conclusions
If we had additional measurements per person and some additional information on the kinetics, In the time it took me to write this post several other papers have come out with much more information that would be fun to analyze that find that immunity is likely long lasting (Mateus et al., n.d.; Grifoni et al. 2021).
References
Corbett, Kizzmekia S., Martha C. Nason, Britta Flach, Matthew Gagne, Sarah O’Connell, Timothy S. Johnston, Shruti N. Shah, et al. 2021. “Immune Correlates of Protection by mRNA-1273 Vaccine Against SARS-CoV-2 in Nonhuman Primates.”Science, July. https://doi.org/10.1126/science.abj0299.
Grifoni, Alba, John Sidney, Randi Vita, Bjoern Peters, Shane Crotty, Daniela Weiskopf, and Alessandro Sette. 2021. “SARS-CoV-2 Human t Cell Epitopes: Adaptive Immune Response Against COVID-19.”Cell Host & Microbe 29 (7): 1076–92. https://doi.org/10.1016/j.chom.2021.05.010.
Mateus, Jose, Jennifer M. Dan, Zeli Zhang, Carolyn Rydyznski Moderbacher, Marshall Lammers, Benjamin Goodwin, Alessandro Sette, Shane Crotty, and Daniela Weiskopf. n.d. “Low-Dose mRNA-1273 COVID-19 Vaccine Generates Durable Memory Enhanced by Cross-Reactive t Cells.”Science 0 (0): eabj9853. https://doi.org/10.1126/science.abj9853.
Singanayagam, Anika, Seran Hakki, Jake Dunning, Kieran J. Madon, Michael A. Crone, Aleksandra Koycheva, Nieves Derqui-Fernandez, et al. 2021. “Community Transmission and Viral Load Kinetics of the SARS-CoV-2 Delta (b.1.617.2) Variant in Vaccinated and Unvaccinated Individuals in the UK: A Prospective, Longitudinal, Cohort Study.”The Lancet Infectious Diseases 0 (0). https://doi.org/10.1016/S1473-3099(21)00648-4.
Suthar, Mehul S., Prabhu S. Arunachalam, Mengyun Hu, Noah Reis, Meera Trisal, Olivia Raeber, Sharon Chinthrajah, et al. 2021. “Durability of Immune Responses to the BNT162b2mRNA Vaccine.”https://doi.org/10.1101/2021.09.30.462488.
Townsend, Jeffrey P., Hayley B. Hassler, Zheng Wang, Sayaka Miura, Jaiveer Singh, Sudhir Kumar, Nancy H. Ruddle, Alison P. Galvani, and Alex Dornburg. 2021. “The Durability of Immunity Against Reinfection by SARS-CoV-2: A Comparative Evolutionary Study.”The Lancet Microbe 0 (0). https://doi.org/10.1016/S2666-5247(21)00219-6.
