A computational method for predicting the most likely evolutionary trajectories in the stepwise accumulation of resistance mutations
Curation statements for this article:-
Curated by eLife
eLife assessment
This report is a useful demonstration of how to predict the mutational pathways to antimicrobial resistance (AMR) emergence, particularly in the enzyme DHFR (dihydrofolate reductase). The methodology is overall solid but some of the claims are only partially supported. The work will be of interest to microbiologists and evolutionary biologists interested in antimicrobial resistance and its population genetics dynamic.
This article has been Reviewed by the following groups
Listed in
- Evaluated articles (eLife)
- Evolutionary Biology (eLife)
Abstract
Pathogen evolution of drug resistance often occurs in a stepwise manner via the accumulation of multiple mutations that in combination have a non-additive impact on fitness, a phenomenon known as epistasis. The evolution of resistance via the accumulation of point mutations in the DHFR genes of Plasmodium falciparum ( Pf ) and Plasmodium vivax ( Pv ) has been studied extensively and multiple studies have shown epistatic interactions between these mutations determine the accessible evolutionary trajectories to highly resistant multiple mutations. Here, we simulated these evolutionary trajectories using a model of molecular evolution, parameterised using Rosetta Flex ddG predictions, where selection acts to reduce the target-drug binding affinity. We observe strong agreement with pathways determined using experimentally measured IC50 values of pyrimethamine binding, which suggests binding affinity is strongly predictive of resistance and epistasis in binding affinity strongly influences the order of fixation of resistance mutations. We also infer pathways directly from the frequency of mutations found in isolate data, and observe remarkable agreement with the most likely pathways predicted by our mechanistic model, as well as those determined experimentally. This suggests mutation frequency data can be used to intuitively infer evolutionary pathways, provided sufficient sampling of the population.
Article activity feed
-
-
Author Response
Reviewer #1 (Public Review):
This work reports an important demonstration of how to predict the mutational pathways to antimicrobial resistance (AMR) emergence, particularly in the enzyme DHFR (dihydrofolate reductase). Epistasis, or non-additive effects of mutations due to their background dependence, is a major confounding factor in the predictability of protein evolution, including proteins that confer antimicrobial resistance. In the first approach, they used the Rosetta to predict the mutant DHFRdrug binding affinity and the resulting selection coefficient, which then became inputs to a population genetics model. In the second approach, they use the observed clinical/environmental frequency of the variants to estimate the selection coefficient. Overall, this work is a compelling demonstration that a mechanistic …
Author Response
Reviewer #1 (Public Review):
This work reports an important demonstration of how to predict the mutational pathways to antimicrobial resistance (AMR) emergence, particularly in the enzyme DHFR (dihydrofolate reductase). Epistasis, or non-additive effects of mutations due to their background dependence, is a major confounding factor in the predictability of protein evolution, including proteins that confer antimicrobial resistance. In the first approach, they used the Rosetta to predict the mutant DHFRdrug binding affinity and the resulting selection coefficient, which then became inputs to a population genetics model. In the second approach, they use the observed clinical/environmental frequency of the variants to estimate the selection coefficient. Overall, this work is a compelling demonstration that a mechanistic model of the fitness landscape could recapitulate AMR evolution; however, considering that the number of mutations and pathways is small, a more compelling description of the robustness of the results and/or limitations of the model is needed.
Major strengths:
- This is a compelling multi-disciplinary work that combines a mechanistic fitness landscape of DHFR (previously articulated in literature and cited by the authors), Rosetta to determine the biophysical effects of mutations, and a population genetics model.
- The study takes advantage of extensive data on the clinical/environmental prevalence of DHFR mutations.
- Provides a careful review of the surrounding literature.
Major weakness:
- Considering that the number of mutations and pathways being recapitulated is rather small, I would suggest a more detailed description of the robustness of the results. For example:
a) Please report the P-value for the correlation of the predicted DDG_{binding, theory} and DDG_{binding, experimental}.
We thank the reviewer for the suggestion. We agree the available experimental data is small, limiting the statistical power of the Pearsons correlation test to determine how well Flex ddG predicts binding free energy change. However, as highlighted in the manuscript, two earlier studies by Aldeghi et al. 2018 & 2019 considered much larger datasets and found a correlation in a similar range to the one we found here. Furthermore, as suggested by the Reviewer, we carried out a onesided T-test with alternative hypothesis that the correlation is greater than 0 and found a p-value of 0.040, suggesting the correlation we observed is significant. We have included this test and p-value to the Results section.
If interested in showing the correct assignment of mutational effects, perhaps use a contingency matrix to derive a P-value.
As suggested by the Reviewer, we used a contingency matrix known as a confusion matrix to determine how accurate Flex ddG is at classifying mutations as stabilising or destabilising. This gave an accuracy of 0.89, sensitivity of 0.83 and a specificity of 1. The p-value associated with this continency table was 0.14, despite the high accuracy, sensitivity and specificity. This is likely due to the small sample size making it difficult to determine significance. This analysis has been included in the Results section.
b) Although the DDG_binding calculation in Rosetta seems to converge (Appendix figures 3 and 4), I do not think the DDG values before equilibration should be included in the final DDG estimate. In practice, there is a "burn in" number of runs where the force field optimizes the calculation to account for potential clashes in the structure, etc. This is particularly important since the starting structures are modeled from homology. Consequently, the distributions of DDG that include the equilibration runs are multimodal (Appendix figure 2), which means that calculating an average may be inappropriate.
Each Flex ddG prediction is independent (see Figure 1 of Barlow et al. 2018 for a summary of the Flex ddG method), i.e. the distribution of values does not represent a MCMC process in which there is a burn-in in order to equilibrate. The structures of both the wild-type and mutant are equilibrated in each run using the backrub algorithm. The reason so many runs are required is because each prediction is from a distribution of possible ddG values associated with that specific mutation and the authors of Flex ddG suggest running 35 runs or more and taking the average of the distribution. Therefore, in order to get an accurate prediction, enough simulations must be run per mutation to adequately characterise the distribution so that the average converges to a constant value.
- The geographical areas over which the mutational pathways are independently estimated are not isolated, allowing for the potential that an AMR variant in one region arose due to "migration" from another area. For example, the S58R-S117N is the most frequent double mutant of PvDHFR in geographically proximate Southern/Southeastern Asia (Fig. 4). To a certain extent, similar mutational patterns occur for PfDHFR in Southern/Southeastern Asia (Fig. 3). Although accounting for mutant migration in the model may be beyond the scope of the study, a clear argument for the validity of the "isolated island" assumption is needed.
The Reviewer is correct that some variants in one region may have arisen due to “migration” from another area. This would impact the method for inferring mutational pathways from regional isolate frequency data but not when considering the worldwide population. If this occurred, we would expect to see a multiple mutant appearing in a region without the precursor (single, double etc) mutations, even in the case of large sample size. However, this does not seem to have been an issue for the pathways we have been predicting here. If it were the case that a variant migrated, and the precursor mutations could not be found in that region, we could look to mutations from neighbouring regions to infer the pathway, under the assumption of migration.
We have added some discussion on this between lines 517-523:
“When inferring pathways at a regional level, it is possible we may encounter instances where genotypes with multiple mutations are observed in a specific region, but the precursor mutations in the pathway are absent. This could happen either due to insufficient sampling of the region or due to "migration" of the variant from a neighbouring region. To infer pathways in the former case more samples would be required, whereas in the latter case we can look to the data from neighbouring regions where the variant is present and use the frequency data of the precursor mutations.”
-
eLife assessment
This report is a useful demonstration of how to predict the mutational pathways to antimicrobial resistance (AMR) emergence, particularly in the enzyme DHFR (dihydrofolate reductase). The methodology is overall solid but some of the claims are only partially supported. The work will be of interest to microbiologists and evolutionary biologists interested in antimicrobial resistance and its population genetics dynamic.
-
Reviewer #1 (Public Review):
This work reports an important demonstration of how to predict the mutational pathways to antimicrobial resistance (AMR) emergence, particularly in the enzyme DHFR (dihydrofolate reductase). Epistasis, or non-additive effects of mutations due to their background dependence, is a major confounding factor in the predictability of protein evolution, including proteins that confer antimicrobial resistance. In the first approach, they used the Rosetta to predict the mutant DHFR-drug binding affinity and the resulting selection coefficient, which then became inputs to a population genetics model. In the second approach, they use the observed clinical/environmental frequency of the variants to estimate the selection coefficient. Overall, this work is a compelling demonstration that a mechanistic model of the …
Reviewer #1 (Public Review):
This work reports an important demonstration of how to predict the mutational pathways to antimicrobial resistance (AMR) emergence, particularly in the enzyme DHFR (dihydrofolate reductase). Epistasis, or non-additive effects of mutations due to their background dependence, is a major confounding factor in the predictability of protein evolution, including proteins that confer antimicrobial resistance. In the first approach, they used the Rosetta to predict the mutant DHFR-drug binding affinity and the resulting selection coefficient, which then became inputs to a population genetics model. In the second approach, they use the observed clinical/environmental frequency of the variants to estimate the selection coefficient. Overall, this work is a compelling demonstration that a mechanistic model of the fitness landscape could recapitulate AMR evolution; however, considering that the number of mutations and pathways is small, a more compelling description of the robustness of the results and/or limitations of the model is needed.
Major strengths:
1. This is a compelling multi-disciplinary work that combines a mechanistic fitness landscape of DHFR (previously articulated in literature and cited by the authors), Rosetta to determine the biophysical effects of mutations, and a population genetics model.
2. The study takes advantage of extensive data on the clinical/environmental prevalence of DHFR mutations.
3. Provides a careful review of the surrounding literature.Major weakness:
1. Considering that the number of mutations and pathways being recapitulated is rather small, I would suggest a more detailed description of the robustness of the results. For example:
a. Please report the P-value for the correlation of the predicted DDG_{binding, theory} and DDG_{binding, experimental}. If interested in showing the correct assignment of mutational effects, perhaps use a contingency matrix to derive a P-value.
b. Although the DDG_binding calculation in Rosetta seems to converge (Appendix figures 3 and 4), I do not think the DDG values before equilibration should be included in the final DDG estimate. In practice, there is a "burn in" number of runs where the force field optimizes the calculation to account for potential clashes in the structure, etc. This is particularly important since the starting structures are modeled from homology. Consequently, the distributions of DDG that include the equilibration runs are multimodal (Appendix figure 2), which means that calculating an average may be inappropriate.
2. The geographical areas over which the mutational pathways are independently estimated are not isolated, allowing for the potential that an AMR variant in one region arose due to "migration" from another area. For example, the S58R-S117N is the most frequent double mutant of PvDHFR in geographically proximate Southern/Southeastern Asia (Fig. 4). To a certain extent, similar mutational patterns occur for PfDHFR in Southern/Southeastern Asia (Fig. 3). Although accounting for mutant migration in the model may be beyond the scope of the study, a clear argument for the validity of the "isolated island" assumption is needed. -
Reviewer #2 (Public Review):
In this work the authors use a simple biophysical model to predict evolutionary trajectories of resistance to pyrimethamine - inhibitor of PfDHFR from P. falciparum and PvDHFR from P. vivax - pathogens causing malaria which presents a worldwide health concern. The authors use a simple fitness model that posits that selection coefficient -relative change in fitness between WT and mutant strains is determined by the fraction of unbound (to antibiotic inhibitor) DHFR. The population genetics simulations use the Kimura formula which is applicable to low mutation high selection regime where populations are monoclonal. The authors use computational tool Rosetta Flex ddG to assess binding of the antibiotic ligand to WT and mutant protein and compare their predicted evolutionary trajectories with lab evolution and …
Reviewer #2 (Public Review):
In this work the authors use a simple biophysical model to predict evolutionary trajectories of resistance to pyrimethamine - inhibitor of PfDHFR from P. falciparum and PvDHFR from P. vivax - pathogens causing malaria which presents a worldwide health concern. The authors use a simple fitness model that posits that selection coefficient -relative change in fitness between WT and mutant strains is determined by the fraction of unbound (to antibiotic inhibitor) DHFR. The population genetics simulations use the Kimura formula which is applicable to low mutation high selection regime where populations are monoclonal. The authors use computational tool Rosetta Flex ddG to assess binding of the antibiotic ligand to WT and mutant protein and compare their predicted evolutionary trajectories with lab evolution and data on naturally evolved variants worldwide and find semi-quantitative agreement, albeit sith significant variation in detail.
The paper is of potential interest as it presents one of the first (but not the first) attempts to compare evolutionary dynamics based on biophysics inspired fitness model with laboratory evolution and natural data for very important problem of emergence and fixation of antibiotic resistant alleles. As such it can be a useful starting point for more detailed and biophysical realistic models of evolution of resistance against anti-DHFR drugs.
-
Reviewer #3 (Public Review):
In this work, Eccleston et. al. use a computational method involving the Rosetta (Flex ddG) suite to infer epistasis in binding free energy changes for combinatorial sets of mutations in the DHFR gene and the drug pyrimethamine. They use this to estimate the most likely path of stepwise mutation accumulation in the evolution of antimalarial drug resistance. The authors also infer likely pathways from different geographical regions from isolated data using a method based on mutation frequencies. They report that these results are broadly consistent with their computational predictions as well.
In contrast to machine learning approaches, the Rosetta Flex ddG method uses physical models at the atomic scale to compute various macromolecular properties. The present paper, therefore, uses atomic-scale molecular …
Reviewer #3 (Public Review):
In this work, Eccleston et. al. use a computational method involving the Rosetta (Flex ddG) suite to infer epistasis in binding free energy changes for combinatorial sets of mutations in the DHFR gene and the drug pyrimethamine. They use this to estimate the most likely path of stepwise mutation accumulation in the evolution of antimalarial drug resistance. The authors also infer likely pathways from different geographical regions from isolated data using a method based on mutation frequencies. They report that these results are broadly consistent with their computational predictions as well.
In contrast to machine learning approaches, the Rosetta Flex ddG method uses physical models at the atomic scale to compute various macromolecular properties. The present paper, therefore, uses atomic-scale molecular properties to make predictions at the population level. As acknowledged by the authors, their method has the limitation that chemical factors other than the free energy changes are largely ignored, as are complications arising from complex population dynamics. Nonetheless, there is reasonable agreement between their predictions and the experimental data, especially at high drug concentrations.
The authors also infer likely trajectories of mutation acquisition from isolate data from various parts of the world. The inference method is based on a simple ranking scheme of mutation frequencies. It is difficult to gauge the reliability of this method, given the complexity of infectious disease dynamics, including confounding factors introduced by varied drug treatment regimens. However, predictions from the computational method are still able to capture some of the general trends in the inferred pathways from isolates, inspiring some confidence in both approaches. The authors emphasize the importance of geographic variation in evolutionary pathways, but their computational method is limited in its ability to provide quantitative insights into the origins of such variation.
A few limitations of the work should be mentioned. It suffers from a lack of summary metrics that quantify the performance of its computational method, which is important for a clearer understanding of its accuracy. While the work is a useful indicator of the potential usefulness of the Rosetta Flex ddG method in enabling evolutionary predictions through macromolecular modeling, the method is applied to a well-studied system and the work remains limited in the novelty of the insights it generates into the dynamics of the evolution of antimalarial drug resistance.
-
-
-
-