Design of an optimal combination therapy with broadly neutralizing antibodies to suppress HIV-1

  1. Colin LaMont
  2. Jakub Otwinowski
  3. Kanika Vanshylla
  4. Henning Gruell
  5. Florian Klein
  6. Armita Nourmohammad  Is a corresponding author
  1. Max Planck Institute for Dynamics and Self-Organization, Germany
  2. Laboratory of Experimental Immunology, Institute of Virology Faculty of Medicine and University Hospital Cologne, University of Cologne, Germany
  3. Department of Physics, University of Washington, United States
  4. Fred Hutchinson Cancer Research Center, United States
7 figures, 2 tables and 1 additional file


Schematics for the evolutionary dynamics of viral rebound.

(A) The viral dynamics after the initiation of a treatment with bNAb infusion (t=0) is determined by two competing processes. Susceptible strains (sus) undergo exponential decay (red line) with decay rate given by r, while the resistant mutants (mut) undergo logistic growth back up to the carrying capacity (Nk) of the patient. In the deterministic limit (Equation 1), the rebound time is linearly related to the log-frequency of the mutant fraction. (B) The schematic shows the four stochastic processes of birth, death, mutation, and neutralization with their respective rates for susceptible (purple) and resistant (blue) variants. These processes define the evolution of a viral population. Note that both the susceptible and the resistant variants are subject to birth and death with their respective rates. (C) The birth and death rates can be visualized as a region of size λ=β+δ which is partitioned into birth and death events. In the absence of antibodies, the susceptible population has balanced birth and death rates, βs=δs, while the resistant population has a negative net birth rate equal to the fitness difference Δ=δm-βm. After introduction of the antibody, the susceptible population decays at rate r, and without competition from the susceptible population, the resistant population grows at the free growth-rate γ. (D) Mutational target size is inferred a priori from the genotype-phenotype mapping, which can be visualized as a bipartite graph. The nodes correspond to codons, while the edges are the mutations which link one codon to another, weighted according to the respective mutation rates. The average edge weight from codons of susceptible variants to the escape mutants determines the rate of escape mutations μ. Mutations can be divided into two types: transitions (black) are within-class, and transversions (red) are out of class nucleotide changes. Transitions occur at about 8 times higher rate than transversions (Figure 2). (E) A coarse grained fitness and mutation model for two of the escape sites (281 and 282) against antibody 3BNC117 are shown. Left: At each escape mediating site, amino acids fall into one of three groups: (i) susceptible (wild-type), (ii) escape mutant, and (iii) fatal. For an escape-class amino acid at site i the virus incurs a fitness cost Δi, and these costs are additive across sites. Right: Mutations at a given site i occur with (independent) forward μi and backward μi rates which govern the substitution events between amino acid classes.

Figure 2 with 1 supplement
Statistics of viral genome sequences from bNAb-naive HIV-1 patients.

(A) Statistics of the high-throughput longitudinal data collected from HIV-1 populations in 11 ART-naive patients from Zanini et al., 2015, is shown. Some of the have low diversity (vertical red lines) and were not usable for our study. Usable samples (vertical black lines) amount to 4–10 samples per patient, collected over 5–10 years of infection. (B) Lower panels show the relative frequencies (cube-root transformed for legibility) of different amino acids in four patients at the 3 escape sites against the 10–1074 bNAb, estimated from the polymorphism data at the nucleotide level in each patient over time. Despite 10–1074 being a broadly neutralizing antibody, mutations associated with escape (indicated by a *) are commonly observed in untreated patients. The upper right panel shows the number of individuals (out of the cohort of 11 bNAb-naive HIV-1 patients) that carry mutations associated with escape against the indicated bNAbs. (C) The nucleotide diversity associated with transversion θtv=2Neμtv is shown against the transition diversity θts=2Neμts for all patients (colors) and all time points. The covariance of these two diversity measurements yields an estimate for the transition/transversion ratio θts/θtv=7.6. (D) Left: The transition diversity is shown to grow as a function of time since infection in all the 11 patients (colors according to (A)). Right: The neutral diversity of viral populations in patients (points) from the three different clinical trials (Caskey et al., 2015; Caskey et al., 2017; Baron et al., 2018) analyzed in this study resemble the larger diversities of long-established viral populations in untreated patients. (E) The inferred forward and backward mutation rates (μ,μ), relative to the transition rate, and the median selection strength σ=2Ne(fsus-fres) at each escape site against the two bNAbs (10–1074, and 3BNC117) from the trial data used in this study are shown. Compared to the 10–1074 bNAb, escape from the 3BNC117 bNAb appears to be less costly, and is associated with a smaller mutational target.

Figure 2—figure supplement 1
Inference of neutral diversity from genetic data.

(A) The inter-decile ranges of neutral diversity estimates for 103 simulations of genomes of 102 independent neutral sites are shown. The true values of the diversity used for simulations are shown in black. Inter-decile maximum likelihood estimate (MLE) for the diversity θ from Equation 22 are shown in red, and the observed variance of allele frequency x as a measure of diversity θ=2x(1-x) is shown in blue. The inferred diversity associated with transversions θtv is shown against that of the transition θts for patients (colors) enrolled in (B) the longitudinal study with high throughput HIV-1 sequence data (Zanini et al., 2015), and in the bNAB trials with (C) 10–1074 (Caskey et al., 2017), (D) 3BNC117 (Caskey et al., 2015), and (E) combination of 10–1074 and 3BNC117 (Baron et al., 2018). We find that all four diversity distributions share similar ratios θts/θtv, indicated in each panel.

Figure 3 with 3 supplements
Statistics of viral rebound in clinical trials with bNAbs.

(A) Panels show viremia of three patients from the 3BNC117 trial over time (black circles) and the fitted model of the viral decay and rebound processes from Equation 1 (orange line). The viral rebound time T and the fitted carrying capacity Nk is shown in each panel. Shown are examples of a non-responder (NR; left), a rebound occurring during the trial window (0<T<56 days; center), and a late rebound (T>56 days; right). (B) We compare the distribution of rebound times in patients from the three clinical trials with 10–1074 Caskey et al., 2017, 3BNC117 Scheid et al., 2016, and the combination of the two bNAbs Baron et al., 2018 to the predictions from the simulations based on our evolutionary model (Figure 1, and Methods). The error bars show the inter decile range (0.1–0.9 quantiles) generated by the simulations for the corresponding trial. (C) The summary table shows the number of patients for whom the infecting HIV-1 population shows no response (NR), rebound during the trial window 0<T<56, and a late rebound (T>56 days) in each trial. Note that three patients were excluded from the 3BNC117 trial (*) because of insufficient dosage leading to weak viral response: 1mg/kg compared to the 3-30mg/kg in the other treatment groups. (D) Plotted are 1,200 trajectories of the mutant viral population simulated using our individual based model. Due to the individual birth-death events, fluctuations are larger when the population size is smaller. At a critical threshold, xext, fluctuations are large enough to lead to almost certain extinction in the existing viral population. The critical threshold (yellow line) is an order of magnitude larger than the post-treatment spontaneously-generated mutant fraction (red line). (E) The predicted fraction of escape events associated with post-treatment spontaneous mutations (red) and the pre-treatment standing variation (yellow) are shown for the three trials. Late rebound events are indicated in blue. Because the spontaneously-generated mutant fraction is smaller than the extinction threshold, these mutations contribute to less than 4% of escape events (red), and escape is likely primarily driven driven by standing variation (yellow), that is, pre-existing escape variants.

Figure 3—figure supplement 1
Dynamics of viremia in patients from the 10–1074 trial.

Panels shows the measured viremia over time (black dots) and the fitted deterministic dynamics from Equation 5 (red line) for all patients off ART in the 10–1074 trial Caskey et al., 2017. Patient-specific fitted carrying capacity Nk and the estimated rebound times are reported in each panel. A common decay rate r=0.36 day-1 is fitted to all patient data in this trial, and growth rate is set to 0.33 day-1.

Figure 3—figure supplement 2
Dynamics of viremia in patients from the 3BNC117 trial.

Similar to Figure 3—figure supplement 1 but for the 3BNC117 trial Caskey et al., 2015. The cohort treated with a low dosage of bNAb (1 mg/kg as opposed to 3-30 mg/kg) is shown in yellow; these patients exhibited a very weak response to there treatment and were excluded from our analysis. A common decay rate r=0.23 day-1 is fitted to all patient data in this trial, and growth rate is set to 0.33 day-1.

Figure 3—figure supplement 3
Dynamics of viremia in patients from the combination therapy trial.

Similar to Figure 3—figure supplement 1 but for the combination therapy with 10–1074 and 3BNC117 Baron et al., 2018. A common decay rate r=0.33 day-1 is fitted to all patients off ART in this trial, and growth rate is set to 0.33 day-1.

Figure 4 with 4 supplements
Statistics of viral escape for optimal combination therapy with bNAbs.

(A) The posterior distribution for inferred selection strength on the escape-mediating sites associated with each of the 9 bNAbs in this study is shown (right); white line: median, box: 50% around the median, bar: 80% around median. Each escape site is color coded by its location on the env gene (left) and each antibody by its associated epitope location. (B) The harmonic mean of the selection strength σ associated with cost of escape (scaled by transition diversity θts) is shown against the mutational target size for each bNAb; error bars indicate 50% around the median. For antibodies to be broadly neutralizing, it is sufficient that viral escape from them to be associated with a small mutational target size or a large fitness cost. The mutational target size is found to be weakly correlated with the average cost of escape from a given bNAb. We identify two distinct strategies for antibody breadth—selection limited and mutational-target-size limited escape pathways each highlighted in gray. (C) bNAb therapies with 1, 2, and 3 antibodies are ranked based on the predicted probability of early viral rebound, and in each case, six therapies with highest efficacies are shown; best ranked therapy is associated with the lowest probability of early rebound (lt56 days); indicate 50% around the median. Also for reference, the probability of early viral rebound two therapies from the trials in this study (10–1074 and 10–1074+3BNC117) are shown.

Figure 4—figure supplement 1
Robustness of selection inference to intra-patient temporal correlations of HIV-1 alleles and clade-specific sampling.

Fitness cost for escape mediating sites of different bNAbs is shown based on the selection inference using (A) the time-averaged likelihood (Equation 40 t-averaged) with data from all patients, as a reference (similar to Figure 4A), (B) the time-averaged likelihood (Equation 40 t-averaged) with data from patients infected with HIV-1 clade B only, and (C) independent time likelihood (Equation 40 t-independent) with data from all patients. Limiting the inference to patients infected with clade B virus in (B) result in minor changes in the inferred selection strength and a slight increase in uncertainties due to a smaller data. However, the general structure of the posteriors remain similar to (A). Treating all time points as independent in (C) increases the effective sample-size of the data and reduces model uncertainty. Nonetheless, except for slight changes on the selection ordering of escape sites against 3BNC117, the median of the selection posteriors remain similar to (A). The color code is similar to Figure 4A and is determined by the env region associated with each site.

Figure 4—figure supplement 2
Robustness of selection inference to genetic linkage and hitchhiking.

Results of selection inference for simulated populations that carry escape mediating sites linked with other selected sites in the genome are shown. The blue shading shows the interdecile ranges (0.1–0.9 quantiles) for the inferred maximum likelihood estimates (MLE) of scaled selection σ^=σ/θts on the escape mediating sites versus the true values used in the simulations, for different genetic linkage effects (columns) and different values of population diversity (rows). Each inferred value is based on 10 realizations of simulations (in-silico patients) for populations evolving over time, with samples taken at 10 time points spaced at 100 days. Simulations are done on whole genomes of length 512, in which 32 sites under strong selection 0.02f0 are equally spaced throughout the genome; f0 is the base growth rate of 1.0 day-1. Each day, the preferred allele in one of the 32 sites is swapped (selection changes sign) with probability pflip. The rate of selection fluctuation pflip and the recombination rate χ are given in each panel. The carrying capacity is set to Nk=5,000, and the event rate is λ=2 day-1; see Algorithm 3 (Appendix 5). To initialize, we let the populations equilibrate for 4Ne generations, and then gather data for maximum likelihood inference of selection at the escape mediating sites. In all panels, the median MLE for selection (black line) closely matches the true selection strength (redline). Quantile lines are smoothened by a Gaussian kernel. Variation is largest when diversity is low and when mutants are rare, which concords with the results from the Bayesian inference procedure.

Figure 4—figure supplement 3
Minimum disparity estimation for adjustment of diversity.

(A) The predicted rebound times without including an increase in diversity due to reservoirs (i.e. ξ0=1 in Equation 42) qualitatively matches the data from the three trials (histograms in each panel), with an overestimation in the number of non-responders and patients with late rebound. For comparison, the predictions with reservoir-corrected diversity is shown in Figure 3A. (B) The disparity R(ξ|t1:p) (Equation 35) is shown as a function of the reservoir factor ξ for the three trials (colors; left panel). The total disparity over all trials is shown in the right panel, where ΔR indicates the reduction in disparity by using the optimal reservoir factor ξ*=2.1. (C) The distribution for the reduction in disparity ΔR (Equation 41) is shown for 1000 realizations of simulated data under the null hypothesis with ξ0=1. The large reduction in disparity indicates that the null hypothesis can be rejected with p-value=0.004 (Equation 42).

Figure 4—figure supplement 4
Disparity analysis for robustness of selection inference.

(A) The disparity R(ξs|t1:p) (Equation 35) is shown as a function of rescaling of the selection strength σ/θts by a multiplicative factor ξs. The total disparity over all trials is shown in the right panel, where ΔR indicates the reduction in disparity by using the optimal selection factor ξs=.8. In this case the null hypothesis is that the inferred selection by the Bayesian procedure in Equation 31 requires no further rescaling for the model to fit the distribution of the rebound times (i.e. ξs=1). The large p-value=0.49 shown in (B) indicates that the null hypothesis cannot be rejected and we have no statistical justification for adding an adjustment factor for selection inferred from untreated patients.

Appendix 3—figure 1
The effect of incomplete neutralization.

Viral load curves, simulated under the noise-averaged model (deterministic limit) Equation 4, are plotted under the assumption of incomplete neutralization for different values of [A]0/IC50 (colors); incomplete neutralization is defined by a changing growth rate given by Equation 53. Small values of [A]0/IC501 (red) lead to similar trajectories to the perfectly resistant variant (solid black line), while large [A]0/IC501 (green line) is indistinguishable from the fully susceptible variants (dotted black line). For all curves, the escaping virus has an initial frequency xres=10-3, with decay rate r and free growth rate γ0, r=γ0=1/3 days.

Appendix 3—figure 2
The distribution of IC50 and [A]0.

The distribution of IC50for envelope proteins derived from circulating viruses in patients after 10–1074 infusion is shown against (A) 10–1074, (B) 3BNC117 as a control. For comparison, the distribution of the maximum (initial) 10–1074 bNAb concentration used in the trial is shown in blue. The maximum measurable IC50=50μg /ml associated with the TZM-bl neutralization assays used in these experiments is indicated in each panel.

Appendix 4—figure 1
Statistics of viral rebound with using escape pathways inferred from the DMS and the trial data.

We compare the distribution of rebound times in patients from the clinical trials with 10–1074 Caskey et al., 2017 (top) and PGT121 Stephenson et al., 2021 (bottom) bNAbs, using escape pathways inferred from the respective trial data (left) and the DMS data Dingens et al., 2017 (right). The error bars show the inter decile range (0.1–0.9 quantiles) generated by the simulations for the corresponding trial. The PGT121 trial predictions relied on the neutral diversity estimates θts from the other three trials Caskey et al., 2015; Caskey et al., 2017; Baron et al., 2018 due to the relatively limited genetic data available from this study. The fitted reservoir value of rresv.=2.07 is used for all predictions.


Appendix 1—table 1
Selection and mutational target size for escape-mediating sites against each bNAb.

Shown are the sites (column 1) and the susceptible and the escape amino acids (column 2) for each bNAb. We called patterns for 10–1074 and CD4 binding site targetting antibodies VRC01 and 3BNC117 using genetic trial data and the remainder using DMS data. The inferred mutational target size of escape at each site (forward μ and backward μ mutation rates) is shown in column 3. The major quantiles (10%, 50%(median),90%) associated with the inferred site-specific Bayesian posterior of the scaled selection strength σ^=σ/θts are shown in column 4. The corresponding quantiles for the strength of selection,σ=σ^θts after convolving the posterior for the scaled selection σ^ with the reservoir-corrected intra-patient diversity of HIV-1 ξθ are shown in column 5.

σ/θ quantilesσ quantiles
bNAbsitesusceptible AAescape AAμμ0.
Appendix 4—table 1
Escape variants inferred from trial and DMS data.

Tables show the sites mediating escape of HIV-1 from the 10–1074 and the PGT121 bNAbs, using the respective trial data Caskey et al., 2017; Stephenson et al., 2021 (left) and the DMS data Dingens et al., 2017 (right). The susceptible (sus) and resistant (res) variants at each site are identified according to the procedure detailed in the Methods. Sites that are identified as escape-mediating by only one of the methods are marked with (*). While substitutions at site 325 were called as escape mediating based on the PGT121 trial data Stephenson et al., 2021 but are missed by the DMS data, this site was only inconsistently associated with escape, and its associated escape variants usually co-occur with other escape substitutions; see the Extended Data Figure 3 in Stephenson et al., 2021. Sites are numbered according to HXB2.

10–1074 trial10–1074 trial
HXB2 sitesusresHXB2 sitesusres
PGT121 trialPGT121 DMS
HXB2 sitesusresHXB2 sitesusres

Additional files

Download links

A two-part list of links to download the article, or parts of the article, in various formats.

Downloads (link to download the article as PDF)

Open citations (links to open the citations from this article in various online reference manager services)

Cite this article (links to download the citations from this article in formats compatible with various reference manager tools)

  1. Colin LaMont
  2. Jakub Otwinowski
  3. Kanika Vanshylla
  4. Henning Gruell
  5. Florian Klein
  6. Armita Nourmohammad
Design of an optimal combination therapy with broadly neutralizing antibodies to suppress HIV-1
eLife 11:e76004.