Phylodynamic analysis in the understanding of the current COVID-19 pandemic and its utility in vaccine and antiviral design and assessment.

Over the last decades, the use of phylogenetic methods in the study of emerging infectious diseases has gained considerable traction in public health. Particularly, the integration of phylogenetic analyses with the understanding of the pathogen dynamics at the population level has provided powerful tools for epidemiological surveillance systems. In the same way, the development of statistical methods and theory, as well as improvement of computational efficiency for evolutionary analysis, has expanded the use of these tools for vaccine and antiviral development. Today with the Coronavirus Disease 2019 (COVID-19), this seems to be critical. In this article, we discuss how the application of phylodynamic analysis can improve the understanding of current pandemic dynamics as well as the design, selection, and evaluation of vaccine candidates and antivirals.

Introduction

Recent developments in phylogenetic methods have made it possible to detect and track evolutionary changes in pathogens in time and space, constituting an evolving branch of phylogenetics called phylodynamics.[1] The development and application of coalescent theory in this field have allowed the integration of phylogenetic analysis with traditional epidemiological analysis in the understanding of the geospatial and temporal pathogen dynamics at the population level.[1,2] This integration allows one to infer geospatial transmission patterns and networks, as well as epidemiological parameters of infectivity, and spread of disease when used at the “between host” level. Additionally, the construction of a viral demographic model based on viral diversity and its coalescence could be used also in assessing control measures like vaccines or therapeutics.[3]

The coalescent theory is based on a mathematical model for lineage separation and genetic drift that allows one to analyze genetic variation as a stochastic process, describing population events retrospectively.[4] In viral evolution, this theory permits analysis of genetic diversity by testing several models based on evolutionary processes (genetic drift, mutation, natural selection, etc.) and makes inferences about changes in viral population size, immune selection, and spatial dynamics.[5-7] Phylodynamics applies the coalescent theory to describe the relationship between a population’s genetic history and the shared ancestry of randomly sampled individuals and uses the molecular clock concept to infer the timing of intersecting phylogenic events.[8]

Therefore, phylodynamics details certain aspects of viral demography like rates of viral population growth and decline, the distribution of branching events in phylogenetic trees, population structure, and immune selection. It also enables the estimation of the basic reproductive number (R0) directly from gene sequence data, providing a link between evolutionary analysis of genome sequences and the epidemiology of an infectious disease.[2] In the context of the current pandemic, the analysis of the evolutionary relationship of genes or genomes with spatiotemporal variables has allowed identifying the possible origin of the cross-species transmission of the Severe Acute Respiratory Syndrome Coronavirus 2 (SARS-CoV-2), causing the Coronavirus Disease 2019 (COVID-19), and the geographical patterns of infection migration.[9]

On the other hand, patterns inferred from phylodynamic analysis allow making inferences about the immune response of the host. When the viral diversity is analyzed within an individual (within host level) and at the population level (between host level), and contrasted with the level of conferred immunity, one can infer the level of cross-immunity between different genotypes and the level of immune selective pressure.[10] In the context of the current SARS-CoV-2 pandemic, these data could be used for making inferences and predictions about the behavior of this emerging infection and the effectiveness of vaccination strategies in the short, mid, and long term, as has been used for influenza.[11,12] Moreover, the antiviral design could be improved and the effectiveness of viral control measures over time could be assessed.[8]

Although historically the study of the phylodynamic of pathogens has been limited by the cost and timescale of the generation of molecular data, during this pandemic, genomic data about SARS-CoV-2 have been widely and freely available (Figure 1).[15,16] These, along with the possibility of performing next-generation sequencing (NGS) to reveal the genetic variation present among viruses in a single sample, can provide information about the viral diversity at the individual level for vaccine and therapeutics design. In this article, we discuss how the application of phylodynamic analysis has helped not only in the understanding of current pandemic dynamics but also how it could help in the design, selection, and evaluation of vaccine candidates and antivirals.

Figure 1.

Phylogenetic clades classification of SARS-CoV-2 based on full-genome sequences by GISAID.[13,14]

Using phylodynamics for infectious disease modeling and understanding the SARS-CoV-2 pandemic

Under the family of SIR models, populations are usually divided into at least three compartments.[17] These three compartments are 1. Susceptible (S) that corresponds to nonimmune non-exposed hosts that can get infected, 2. Infected (I) that corresponds to individuals infected that can spread the disease, and 3. Recovered (R) that corresponds to individuals that recover or are “removed” from the population, this is, they get longstanding immunity for a timeframe or die. Some models under this family, add compartments, like the E compartment, under the SEIR model, that corresponds to Exposed individuals that do not spread the infection, since they have a latent asymptomatic infection but are not infectious, which is not the case for SARS-CoV-2 given that asymptomatic patients can also spread the disease (Figure 2a), and also, can consider the rate of becoming susceptible again after infection, as would be important in the case of reinfection (Figure 2a), a plausible but uncertain phenomenon in the case of SARS-CoV-2,[17] and the rate of infection from different human to human sources, like contaminated surfaces (sparking rate, Figure 2a), as has been documented for SARS-CoV-2 and other viruses.[18]

Figure 2.

Schematic diagram of between-host (a) and intra-host virus (b) dynamics. Arrows correspond to different rates, while pictures correspond to different states of the viral SIR model created with BioRender.com. *In this case, birth rate corresponds to the rate at which new individuals born, **Sparking rate corresponds to transmission from outside non-human to human sources, as it has been documented for SARS-CoV-2 and other viruses

Under these models, phylodynamic inference using Bayesian statistics allows us to estimate the size of each compartment based on the analysis of viral diversity and its coalescence back in time.[19] Under the assumption of strict or relaxed molecular clock models, different approaches could be implemented for estimation of the number of susceptible, infected, and recovered cases using open resources like BEAST or R.[20-23] These kinds of analysis require as input nucleotide sequences alignment of genic or whole-genome data sampled serially over the time or sampled at one time point. For each sequence, sampling date is specified, and this information with a scheme of substitutions and a molecular clock model helps us to estimate the timescale of the tree and the rate of evolution as it connects the temporal information contained in the sampling times to the genetic similarities embedded in the sequences.[24] Additional to these parameters, we also choose a distribution model that describes how the population size is expected to change over time.

This method provides estimations of the infection rate (births), the recovery rate (deaths), and the population size (N).[20] Under the birth–death SIR model, a constant birth-death rate is used to model the spread of an epidemic and to model speciation and extinction. Under the epidemiology view, a birth event is a transmission and a death event is a recovery or death.[19] The information provided by the branching times and the molecular clock rate can be used for estimation of the time to the most recent ancestor (TMRA), which provides information regarding the date of introduction of the virus to a specific geographic area.[22]

The outputs of this analysis are used for providing important estimations about epidemiological parameters like the basic reproductive number (R0) which corresponds to the total number of secondary cases caused by an infectious individual when introduced into a fully susceptible population, the critical community size (CCS) which is the minimum number of individuals in a population required for a pathogen to persist in a population, and the critical vaccination threshold (pc), which corresponds to the proportion of the population that must have immunity to make the epidemic disappear.[2,17] Since the control of an epidemic relies on the reduction of the R0 to <1, the pc can be calculated based on the calculation of the R0. and the estimates of the real burden of the disease could also be assessed through the estimation of the epidemic size.[25,26]

In the context of the current pandemic, phylodynamic analysis has not only analyzed and estimated those parameters but also has been in good agreement with epidemiological observations.[9,27] The origin of the outbreak has been tracked by epidemiological and phylodynamic analyses to Wuhan, China, between November and December 2019 (95% CI: November 21-December 20, 2019).[9,27] In agreement with epidemiological data, the estimated R0 has been 2.15 (95%CI: 1.79–2.75), with a doubling time of 7.1 days (95% CI: 3.0–20.5).[28-30] In Latin America for instance, the contrast of epidemiologically observed and phylodynamically inferred R0 has shown agreement and has been used for temporal analysis of the introduction of SARS-CoV-2 suggesting a detection lag of at least 21 days in some countries, supporting the complementary approaches of epidemiological and genomic surveillance.[31] In this sense, the COVID-19 pandemic and the rapid availability of genomic information of the SARS-CoV-2 virus have represented a unique opportunity for contrasting observations.

In this regard, other studies have shown how phylodynamic inference can be used in tracking transmission chains[32-34] and evaluating the effectiveness of public health measures through the assessment of the evolution of R0 in time.[32,33] Evaluation of transmission chains in Australia through cluster analysis showed good agreement between epidemiological and genomic clusters,[32] and phylodynamic analysis at the country level in Brazil has shown the role of large and highly connected populated cities in the establishment of SARS-CoV-2.[33] Therefore, phylodynamic analysis can be used as a complement to traditional surveillance methods. Additionally, this analysis can serve as a tool for evaluating transmission chains and public health measures for disease control.

Use of phylodynamics in vaccine design for preventing SARS-CoV-2 infection

Viral virulence, diversity, and fitness represent a challenge for viral vaccine development.[3] These are the most important factors in viral evolution, and phylogenetic analysis is a key tool to simultaneously evaluate these factors.[35] Diversity of RNA viruses can appear over a short time since they are prone to replication errors, generating heterogeneous viral populations.[3] Moreover, antigenic drift and shift, and the bottleneck at the moment of transmission, have an impact on their evolutionary dynamics and determine therapeutic effectiveness, resistance to antivirals, and evasion of the immune response.[36]

However, compared to other RNA viruses, coronaviruses exhibit relatively higher replication fidelity. For example, the SARS-CoV mutation rate was estimated at 9.0 × 10−7 substitutions per nucleotide per replication cycle,[37,38] while most RNA viruses have a rate of 1 × 10−3 to 1 × 10−5 substitutions per nucleotide per replication cycle.[39] Therefore, it seems that SARS-CoV-2 is not mutating significantly as it spreads,[40] which could have important implications for vaccine design.

Although evidence for antigenic drift or shift for SARS-CoV-2 is limited, the sustained human-to-human transmission could permit the acquisition of mutations with fitness advantages and immunological resistance.[36] Analysis of SARS-CoV-2 genomics data have shown evidence of human clonal evolution. In this sense, SARS-CoV-2 genomes analysis has shown mutations that have emerged independently in multiple times.[41] Interestingly, three sites in Orf1ab and one in Spike protein have experienced numerous recurrent mutations (>15 events) that are of particular interest in the context of adaptation of SARS-CoV-2 to the human host.[41] Consistent with that, characterizing frequently mutated residues by aligning ~660 SARS-CoV-2 genomes and validated in 10,000 datasets has shown that Orf1ab and Spike protein contained numerous nonsynonymous mutations, after available in GISAID Nextstrain.[42] However, the degree to which these observed changes are related to antigenic escape or if SARS-CoV-2 may undergo antigenic drift is unclear and has not been documented beyond in vitro experiments,[43] and inferred by the analysis of variants and glycosylation sites.[44,45]

Moreover, once a vaccine for SARS-CoV-2 is developed, it will be important to continue monitoring the genetic diversity and population dynamics of circulating strains to assess their adaptive evolution. At this point, phylodynamics provides a powerful approach for updating vaccines by analyzing evolutionary dynamics and predicting evolutionary changes in populations. Similar approaches have used with the influenza vaccine, in which the assessment of the distribution of variants and its departure from neutrality has helped in predicting emerging lineages and escape variants,[12,46] becoming an integral part of the yearly vaccine design cycle.[35] As another example, phylodynamic approaches have been used as a continuous surveillance tool to evaluate the efficacy of the vaccine for group A rotaviruses, assessing the genotypic diversity and mutation of antigenic epitopes, compared to the prevalent genotypes of the vaccine strains.[47]

Currently, different approaches such as recombinant proteins, viral vectors, DNA, RNA, and live attenuated vaccines (LAV) are being developed as SARS-CoV-2 vaccine candidates, with at least 172 candidates in pre-clinical evaluation, 61 in clinical trials, and 2 vaccines authorized for emergency use by the Food and Drug Administration in the USA (Pfizer-BioNTech COVID-19 Vaccine and Moderna COVID-19 Vaccine).[48] A critical goal of vaccine design is to elicit a similar immune response compared with the natural infection, but also to produce a wide neutralizing immune response against different genotypes to produce a durable and protective response. In the case of vaccines against SARS-CoV-2, the strategies aim to induce neutralizing antibodies, which are directed to the viral spike (S) protein as it mediates virus attachment and entry to host cell, and mainly to the RBD domain of this protein.[49,50]

Single nucleotide polymorphism (SNP) genotyping of SARS-CoV-2 has revealed that some of the more common SNP mutations are located in the S protein, like the 23403A > G spike glycoprotein mutation D614G, but their implications in pathogenesis and immune escape are not clear.[51] While this variant has been effectively neutralized by isolated monoclonal neutralizing antibodies, other variants with mutations in the receptor-binding domain (RBD) showed complete resistance to one of the tested antibodies, highlighting the impact of viral genetic diversity on the ability to produce a broadly neutralizing response.[52] However, little is known about the features of the immune response during natural infection with SARS-CoV-2. The kinetics of neutralizing antibodies suggest a decline after three months, and there are reports consistent with the possibility of reinfection.[53-57]

In this point, the analysis of phylodynamic patterns can be used to infer the level of cross-immunity among different viral clades or lineages, since phylogenies of the virus will be structured depending on the level of cross-immunity between different genotypes or variants. Grenfell et al. have previously discussed the role of immunity in structuring phylodynamic patterns.[10] Briefly, viruses that induce a high cross-immunity across different genotypes will have a phylogeny with a strong spatial structure. In contrast, viruses that induce partial cross-immunity across genotypes will have a phylogeny structured through time, and viruses that produce immune enhancement will have chaotic or complex dynamics.[10] Therefore, although there is no evidence of lack of cross-immunity among SARS-CoV-2, analyzing the structure of the obtained phylogeny would provide insights in this regard between genetic variants. At the same time, information can be gathered from the level of immune response and the number of viral copies through time in groups of individuals. This allows elucidating the Evolutionary Immunity Profile (EIP), which provides information about the probability for a host to be a major source of host-selected variants.[10] The analysis of genome regions and their rate of evolution would allow to identify regions under selective pressure and to suggest new antigenic determinants targeted by neutralizing antibodies.

In this context, a phylodynamic approach could be used to test if the immune response is operating as a selection pressure and diversity accumulated in this region driven by immune escape. The effect of directional selection driven by immune escape can be inferred from the phylogenetic tree balance. For example, in a phylodynamic analysis, the ladder-like shape of a viral phylogeny reveals the hallmarks of selection driven by immune escape. In contrast, a more balanced phylogeny may occur when a virus is not subject to immune selection or other sources of directional selection,[8] although this inference should be taken with caution because these phylogenetic patterns could also reflect sequential bottlenecks occurring with a rapid spatial spread of some types of virus. However, phylogenetic supported with statistical analysis could confirm viral lineages with a larger number of positively selected sites and predicting evolutionary changes, which are important for vaccine re-design.[3]

Taken together, the information obtained through phylodynamic analysis can enhance vaccine design by 1. Predicting the targets of immune response for neutralizing infection, 2. Predicting the level of cross-immunity produced during natural infection, 3. Evaluating the role of vaccine-induced immunity as a selective pressure.

Use of phylodynamic in antiviral research for treating SARS-CoV-2 infection

Just as an epidemic/pandemic is understood through an SIR model, viral infection within a host can be understood in the same way. The viral infection is treated as a microepidemic among host cells. Those cells are treated as susceptible (S) individuals that get infected (I) and then are removed (R) as they die (Figure 2b).[3] Whether the virus infection will spread within the host depends on a condition very similar to the spread of an epidemic at the population level and would depend on the rate at which susceptible cells are infected, and the rate at which infected cells die or are removed by the immune response or the cytotoxic effect of the virus.[3] To propagate the infection in the host, the conditions of viral dynamics must favor R0 > 1, and hence therapeutic interventions such as antivirals should seek to reduce R0 below 1 to promote infection control by reducing the rate at which cells get infected.[58] At the same time, RNA viruses can rapidly adapt to changes in the host environment, since they exist as genetically diverse quasispecies.[59] SARS-CoV-2 is not the exemption, and sequence analysis indicated the presence of viral quasispecies, as also reported for SARS-CoV and MERS-CoV.[60]

The assessment of viral diversity in antiviral experiments, both in vitro and in vivo, can help to estimate their effect on R0 through the analysis of their effect on diversity. This diversity could be also tested looking for evidence of positive selection and directional selection, through the use of probabilistic models of accumulating mutations, assessing the genetic barrier, but also analyzing the diversity control that must be achieved to maintain an R0 < 1. Depending on the antiviral effect on diversity, the need of using multiple drugs and regimens for disease control could be analyzed, and also the probability of escape variants based on the understanding of viral divergence and the substitution rate could be determined.[3,61,62]

Different drugs have been tested against SARS-CoV-2 as antivirals, both as therapeutics and prophylactics including biological therapy, immune-therapeutics, and antiretrovirals.[63-68] The absence of effective antiviral drugs for the treatment of SARS-CoV-2 infection has led to the evaluation of drugs already FDA-approved for the treatment of other viral infections such as lopinavir, ritonavir, ribavirin, remdesivir,[69] oseltamivir, favipiravir, umifenovir,[67] penciclovir, and chloroquine.[70,71] When the introduction of an approved drug happens, phylodynamic analysis can provide tools for monitoring the effectiveness of the intervention by analyzing and quantifying how intra-host viral evolution over time in serial samples and the effect of cessation of replication with effective treatment.[8]

Therefore, the use of the coalescent model would allow the estimation of the effective population size (Ne), which in terms of phylodynamic is proportional to the effective number of infected population,[72] and the assessment of this parameter would help to identify the effect of the introduction of antiviral therapy in the epidemic curve and transmission dynamics.[73] Additionally, although antiviral resistance is more unlikely for an acute virus infection as COVID-19 than for a chronic infection like HIV, the example of anti-retroviral therapy (ART) has shown that evolutionary analysis can be used for the detection of adaptive mutations that if fixed can promote antiviral resistance.[74] This approach could also be used for timely detection of the emergence of drug-resistant strains and the need to switch to a new three-drug regimen (switching therapy).[73] Coupling these data with clinical and epidemiological metadata would help also to identify transmission networks associated with the emergence of resistance in case it happens.[74]

Conclusions

The development of computational and statistical techniques and their application for the understanding of the relationship between spatiotemporal variables (i.e., geographical area, dates) and the genomic diversity of pathogens like SARS-CoV-2 can provide important insights in the understanding of the current pandemics, but also in the design and evaluation of interventions for prevention and treatment of infection. The use of these methods has been historically constrained by the availability of genomic data of the pathogens. However, initiatives like EpiCoVTM by GISAID have made genomic information about the virus free and widely available.[13,75] Given the implications of being able to design, assess, and adapt vaccines through the use of genomic data and computational diversity analysis, local authorities should invest in sequencing technologies for genomic surveillance of current pandemics,[76] especially in developing countries where this technology is not widely available. As the pandemic has expanded, studies contrasting traditional epidemiological surveillance and phylodynamics have emerged. However, this review has been limited by the amount of evidence linking phylodynamic analysis with host immune response and the lack of studies using this approach for the assessment of antiviral drugs. Given the potential benefits of these analyses in the design of vaccines and antivirals, further research is recommended.