Introduction

The coronavirus disease 2019 (COVID-19) pandemic has triggered an unprecedented global response in pathogen genome sequencing, and nearly 400,000 full or partial severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) genomes were generated and shared publicly within its first year. Although phylogenetic tools have become increasingly relevant to the public health management of a range of viral epidemics1,2,3,4, the COVID-19 crisis is the first global health emergency during which large-scale, real-time genomic sequencing and analysis have underpinned public health decisions. The first 12 months of the pandemic were characterized by continual change in the global epidemiological and virological situation, and the analysis of genome sequences was essential in tracking the changing situation. Phylogenetic and phylodynamic approaches (Box 1) can unlock information contained in sampled genomes and are often analysed in conjunction with other data sources5. Such analyses have been used to quantify international virus spread, identify outbreaks and transmission chains in specific settings, estimate growth rates and reproduction numbers, account for surveillance gaps and lags, identify and track mutations of interest, discover and analyse variants of concern, and investigate intra-host virus evolution.

This Review focuses on how SARS-CoV-2 transmission, epidemiology and spatial dispersal have been measured and investigated through phylogenetic and phylodynamic analyses of SARS-CoV-2 genomes (Fig. 1). It is intended to be a retrospective overview that uses examples from the first year of the pandemic to demonstrate the contributions of phylogenetics in the context of different phases of pandemic response. We examine how such analyses have informed global efforts to understand, control and predict the pandemic, and outline arising new challenges and how they are being addressed. We do not review events that precede the widespread emergence of SARS-CoV-2 (such as the evolutionary origins of the pandemic in non-human host species) or its functional genomics (that is, how virus mutations contribute to phenotypes such as transmissibility). Given the scale of the field and the size of the literature on SARS-CoV-2 genomic epidemiology, we do not attempt to provide a systematic review. Instead we focus on studies that represent the first year of the pandemic, which saw evolutionary approaches applied to a wide variety of public health interventions worldwide, often in an ad hoc or pragmatic manner. We further highlight research that was influential in contributing to epidemiological understanding and public health decision making. The pandemic’s first year also best illustrates the potential of these methods for urgent risk assessment, prediction and control of future emerging viruses. We mostly refer to the genetic diversity of SARS-CoV-2 using the Pango dynamic nomenclature6 (Box 1), but also sometimes use the WHO ‘Greek letter’ nomenclature scheme for particular variants of concern (VOCs) and variants of interest (VOIs).

Fig. 1: Phylodynamic approaches to the investigation of SARS-CoV-2 transmission.
figure 1

Relevant clinical and public health questions are defined (top row), phylodynamic and epidemiological data and models are then combined (middle row), and used in combined or joint analyses to provide actionable insight into virus transmission (bottom row). a | Phylogenetic approaches estimate the rate of international lineage introductions and distinguish introductions from community transmission. b | Genome sequences and phylogenetics support outbreak analyses by identifying or refuting links between local cases; this can lead to identification of outbreak sources and drivers or assessment of nosocomial transmission. c | Phylodynamic techniques using epidemiological demographic models, such as the susceptible–exposed–infected–recovered (SEIR) model, allow us to compare transmission rates between lineages bearing different key genotypes (for example, variants of concern (VOCs) and pre-existing lineages). d | Relative timing of variant and lineage emergence from the global (or regional) phylogeny, and scattering of case genomes across clades can distinguish persistent from repeat infections in some scenarios. Phylogenetics is also useful in studies of lineage turnover and interactions within the host. Panel colours indicate related themes: blue, public health; green, epidemiological parameters; red, clinical parameters. TMRCA, time to the most recent common ancestor.

Tracking the global pandemic

Revealing how SARS-CoV-2 spread globally in early 2020 was important in informing public health strategies. Phylodynamic, particularly phylogeographic, methods can be used to estimate the timing and location of ancestral nodes within a molecular phylogeny7,8,9, allowing inference of the route and rate of spread of pandemic lineages, from the site of its initial detection in Wuhan, China, to the location of each sampled patient from which a virus genome was obtained.

International travel restrictions

Phylogeographic studies have investigated the impact of international travel restrictions, quantifying the number of lineage introductions from abroad and the relative contribution of local transmission. For example, a global phylogeny of the pandemic showed that earlier lineages were highly cosmopolitan, whereas later lineages tended to be continent-specific, which may reflect the rapid declines in mobility as many countries concurrently imposed restrictions on international travel10, although early sampling in some countries may have been biased towards cases in international travellers.

At the national scale, studies have typically observed reduced numbers of introductions along international routes covered by travel restrictions; however, the overall effects of this on controlling national transmission depended on the extent to which lineages were already locally well established. During the global expansion of SARS-CoV-2, international exportations were driven initially by dispersal from China; however, the number of exports declined rapidly following the cessation of China’s international flights in January 2020 (ref.11). Endemic transmission began in Italy during mid-February 2020, with establishment in other European countries soon thereafter12. The shift in global dissemination towards greater intercontinental exportation from Europe was associated with the expansion of a lineage bearing the D614G spike mutation13. Virus lineage migrations from Europe to North America increased until the declaration by WHO of a pandemic on 11 March 2020, suggesting that air travel restrictions subsequently slowed international spread14. In South Africa, international introductions plummeted after travel restrictions began on 26 March 2020 (ref.15). Similar observations were made in other nationally focused studies, including those from Italy16, New Zealand, Australia, Iceland, Taiwan17 and the UK18.

The impact of international travel restrictions depended on the level of domestic transmission control and whether restrictions were implemented before full establishment of local transmission. A study of 427 genomes from Brazil applied a discrete asymmetric phylogeographic model and estimated at least 104 international introductions during March and April 2020; these fell into three monophyletic clusters (Box 1) of apparently European origin, and a molecular clock approach indicated that they arrived in late February 2020. Domestic transmission in Brazil was already well established by early March, suggesting that international restrictions implemented thereafter may have had little impact19. In the USA, an early study investigated the efficacy of international travel restrictions in Connecticut20. Seven of nine Connecticut genomes fell into a clade of mostly Washington State genomes, whereas two clustered with genomes from China and Europe. As the Connecticut genomes were derived from people with no history of recent travel, their phylogenetic placement in a cluster of genetically similar genomes indicated community transmission of recently imported lineages; again, flight restrictions may have been more effective in reducing cases if they had been implemented earlier20. Similar patterns were observed in other countries, including Italy12 and the UK21.

Many countries strengthened travel restrictions later in 2020, aiming to slow the spread of variants associated with changes in transmissibility (see the section Tracking lineages of interest). In Brazil, a phylogeny of SARS-CoV-2 from cases detected in São Paulo in late December 2020 indicated two independent international introductions of lineage B.1.1.7 (the alpha VOC) from London, UK22. These introductions occurred despite the suspension of flights to and from the UK. Similarly, phylodynamics suggested multiple international introductions to the USA and hidden transmission of B.1.1.7 since November 2020, and that lineage B.1.1.7 expanded to 33 US states by January 2021 with a doubling time of 9.8 days23. Investigations have also considered the factors that drove the resurgence of transmission in Europe in late summer 2020. A recent study using a Bayesian time-scaled phylogeographic model (Box 1) found that by mid-August a large fraction of the lineages then circulating in European countries had been introduced after 15 June, the date when many countries in the Schengen area opened their borders24. The study also found that newly introduced lineages tended to expand more quickly when entering a region of low incidence, and that for most countries resurgence was driven by new introductions rather than persistence of lineages from the spring24.

Phylogeographic inference of SARS-CoV-2 migration patterns has typically used either a discrete trait analysis (DTA; for example, analyses with travel and/or mobility data)18,19,25 or structured event-based birth–death (BD) models12,17,26. The approaches differ in that DTA assigns discrete states (locations) to nodes on a phylogeny, whereas structured models explicitly model migration events and rates at a population level. The advantage of DTA is its relatively low computational demand26 and its ability to incorporate discrete metadata, such as travel histories (for example, ref.25), in a straightforward manner. However, DTA does not accommodate the interdependency of tree shape and migration rate or population size, and it is more difficult to interpret DTA model parameters26. Structured BD approaches are more computationally costly, but can model variable sampling between regions (DTA is less robust to sampling patterns27), and they infer parameters that can be more readily compared with those obtained from epidemiological or mobility data sets17.

Local transmission and interventions

Non-pharmaceutical interventions (NPIs) include travel restrictions, person-to-person distancing and mandatory mask wearing. Two phylogenetic approaches were typically adopted to investigate the impact of NPIs. First, the frequencies of lineage movement between regions within a country were assessed using phylogeographic analyses (as discussed above for international dissemination). Second, estimates of virus population size, epidemic doubling time and general reproduction number (Rt) were calculated from virus genome sequences using phylodynamic approaches.

Molecular clock dating of SARS-CoV-2 lineages indicated multiple introductions from Wuhan to Guangdong in early January 2020, with a fall in lineage diversity thereafter, suggesting that within-country travel restrictions combined with comprehensive tracing and isolation in Guangdong were effective in controlling transmission28. A phylogenetic study of transmission in Boston, USA, also reported a drop in importations to Boston from other domestic locations after national restrictions began14, and phylodynamic methods estimated a reduction in Rt in Israel of at least two-thirds, coincident with the imposition of quarantine measures29. By contrast, a study of NPIs in Italy16 suggested that domestic travel restrictions failed to prevent community transmission, although there was evidence that transmission was inhibited, and the relatively low genome sampling density means that NPIs could have greatly restricted the transmission of lower-frequency lineages. One global study of 29,000 SARS-CoV-2 genome sequences used a compartmental structured coalescent model to estimate the time of epidemic seeding in 57 different locations30. The authors found that locations with early implementation of strong NPIs experienced less severe morbidity and mortality during the study30 and that stringent interventions 2 weeks earlier would have approximately halved cumulative deaths in the immediate post-intervention period.

Phylodynamic BD models can estimate growth rates and other population parameters from genetic sequence data. A BD-skyline approach estimated a national fall in Rt from 1.63 to 0.48 in Australia after the introduction of travel restrictions and social distancing on 27 March 2020 (ref.31). Similar approaches were used to show that Rt for a New Zealand transmission cluster fell from 7.0 to 0.2 during March 2020, demonstrating the impact of NPIs targeting this cluster32. A multitype BD model applied to data from Taiwan showed a decrease in Rt throughout the early pandemic even in the absence of substantially decreased local human mobility or stay-at-home orders17, suggesting that interventions such as effective contact tracing and widespread face mask use may be sufficient for adequate outbreak control. Phylodynamic studies have provided other parameter estimates that are useful for understanding virus biology and transmission, or for use as statistical priors (Box 1) in further Bayesian modelling (Table 1).

Table 1 SARS-CoV-2 epidemiological parameter estimates using phylodynamic approaches

Phylodynamic analyses have repeatedly demonstrated hidden circulation of SARS-CoV-2 for days to months before first-case detection. Such results are important in determining whether existing surveillance adequately captures ongoing community transmission33. A US study of 346 genomes, covering January to mid-March 2020, examined the establishment of community transmission in Washington State. A phylogeny consistent with community transmission was reported, with most genomes clustered in a clade containing WA1 (USA-WA1-2020, the genome of the first detected US case). The estimated date of origin for the major clade was 18 January to 9 February 2020. This date was used to parameterize a stochastic epidemiological model that suggested 1,600 active infections in Washington State by mid-March34. Similarly, a molecular clock analysis of genomes from Scotland estimated transmission began around 19 February 2020, predating first-case detection by almost 2 weeks21.

Phylogenetics and phylodynamics have also contributed near real-time insights that are suitable for guiding the responses of public health authorities. Many investigations of hospital or event-associated outbreaks during the pandemic employed phylogenetic methods and rapidly provided actionable information. For example, phylogenetics supported public health examinations of numerous outbreaks in New Zealand35, and in the Netherlands excluded a church service, initially implicated, as a source of an outbreak in a care home36. Other studies influenced policy changes: a phylogeographic study of the impact of travel restrictions on lineage imports and transmission37 was used to support the re-introduction of restrictions in Wales in October 2020. A phylogenetic investigation of the June 2020 re-emergence of epidemic transmission in Australia implicated the national mandatory hotel quarantine system38, and the findings led to reform of the quarantine programme. The study38 also used phylodynamics to show the initial growth rate of the second wave to be similar to that of B.1.1.7 emerging in the UK. Furthermore, phylodynamic detection of the increased transmissibility of B.1.1.7 in England39 contributed to the evidence base that informed responses to the first VOC.

Outbreak phylogenetics

Evolutionary approaches can help to refute or confirm suspected transmission routes, supplementing our understanding from contact tracing of cases. Phylogenetic insights can reveal factors associated with transmission, help to establish the polarity of transmission between individuals and estimate outbreak parameters. Genetic analyses have reconstructed events in travel-associated outbreaks and can be used to cross-validate epidemiological records, helping to rule out spurious connections between cases.

Nosocomial transmission

Studies of health-care settings are used to determine whether personal protective equipment (PPE) guidelines are sufficient to prevent nosocomial transmission. In Durban, South Africa, routine surveillance identified a clade of cases in co-workers at a city hospital, suggesting a nosocomial outbreak. The observation of community cases with additional mutations implied community transmission beyond the hospital40. A lack of phylogenetic clustering of cases by ward among health-care workers in a hospital in the Netherlands showed community transmission to be more likely than nosocomial transmission41. Furthermore, a study in Australia ruled out associations between 54 cases across four health services, where shared health-care workers had been initially implicated in dissemination. Phylogenetics revealed that the cases instead actually clustered according to a common social event31.

At a UK renal unit, virus genomes were used to assign responsibility for an outbreak to a shared bus service used to transport outpatients, rather than to transmission from in-patients. Rapid and extensive sequencing resulted in timely revision of the hospital’s infection control procedures42. In a second UK study, phylogenetic analysis of infections from 31 care home staff and 61 residents indicated transmission within, and possibly between, care homes, as well as from staff to staff — the study supported the case against the use of locum staff in such settings43. Policy change was also called for in a Boston hospital study: virus genomes with shared substitutions suggested at least two patient-to-staff transmission events, despite an apparent lack of aerosol-generating procedures and the staff wearing masks and face shields44. In South Korea, the observation of eight near-identical B.2.1 lineage genomes across two Seoul hospitals suggested that the outbreak in one hospital was seeded by a patient transferred from the other45. Multiple introductions were inferred for an outbreak at a San Francisco nursing facility, with one worker, who had also worked in Washington State, implicated in the introduction of WA1-related virus genomes46. Other applications of phylogenetics in investigations of outbreaks in medical or care settings are found in reports from Chile47, France48, Minnesota49 and the Netherlands36. Nevertheless, although phylogenetics has supported the confirmation of nosocomial transmission in some cases, it has also helped reveal the contributions of wider social contact, outside of hospitals and care homes, in the maintenance of transmission networks that span nosocomial settings.

Public gatherings and super-spreading

Epidemiological studies of SARS-CoV-2 have indicated a relatively high attack rate50,51, and phylogenetics has corroborated this finding. For example, an outbreak that affected 11 workers in a large open-plan office in Sweden was supported by a phylogenetic clade of virus genomes from eight workers (six genomes were identical and two near-identical)52. In some cases, local bursts of transmission seem to precede national-scale transmission. Phylogenetic analysis of the early epidemic in Boston identified 28 cases from an international business conference that formed a cluster. All cases shared a novel C2416T substitution, and by November 2020 genomes containing this substitution seemed to underlie 35% of Boston’s cases and 1.9% of US genomes14. This finding showed that individual mass-infection events can facilitate transmission and virus dissemination.

The role of large celebrations in triggering super-spreading has also been explored. Discrete-state phylogeography was used to suggest that a Mardi Gras-associated super-spreading event led to outward (inter-State) dissemination in the southern USA and the acceleration of the early epidemic there53. Resurgence of an early outbreak in Japan was hypothesized initially to be linked to increased travel to cherry blossom sites during the national holiday of 20–22 March 2020. Clarification through sequencing later showed that the late March cases were not directly related to cases from the first epidemic ‘wave’54. In Germany, three events at a Berlin nightclub in early March 2020 led to a series of outbreaks. Phylogenetics confirmed the club as a potential focus of super-spreading; Germany decided to prohibit such events from 16 March55. In the USA, phylodynamics linked the establishment of B.1.1.7 to the Thanksgiving holiday travel surge in November 2020 (ref.23).

Travel and transport

The contribution of transport settings to SARS-CoV-2 transmission has been keenly debated. Virus genomes supported the case for in-flight transmission on a Massachusetts to Hong Kong flight; two flight attendants and two related passengers were detected with B.1 lineage infections, despite B.1 being unknown in Hong Kong at that time56. A similar indication of in-flight transmission was reported for a flight between Dubai (United Arab Emirates) and Auckland (New Zealand)57.

The predominance of one major clade in the February 2020 Diamond Princess cruise ship outbreak suggested that most passengers became infected while attending on-board events, with a single introduction before quarantine measures58. Similarly, a phylogenetic study involving samples from northern California and outbreaks on two consecutive cruises of the Grand Princess ship, with a common crew, found that infected passengers carried three substitutions characteristic of WA1. WA1 at that time was dominant in Washington State, and all cases sampled from the Grand Princess also shared two substitutions that were common in WA1 viruses then circulating in Washington and California. This finding suggested that the source or sources of infection on the cruise were more likely local (that is, California) than either of the cruise destinations. The second cruise, immediately following the first outbreak, shared a subset of passengers with the first cruise. The outbreak phylogeny indicated that one of the first-cruise genomes was ancestral to the second-cruise genomes and also to Californian WA1 genomes in general. This suggested that the shared cohort of passengers seeded the outbreak on the second cruise59. The patterns of shared, derived, mutations in the Grand Princess outbreaks imply large numbers of infections from probably a single infected passenger or crew member (or related transmission cluster) and that the source of infection was local, for example, a crew member, rather than from any station of disembarkation. The implications are that revision of infection management procedures and practices was essential to protect passengers early during the pandemic.

Genomic analyses aided the tracing of transmission during a Chinese–German business meeting in greater Munich (19–22 January 2020), which began an outbreak in Bavaria and involved 16 cases (detected from 27 January to 11 February 2020). Genomes indicated that transmission may have occurred in the pre-symptomatic phase of infection between two individuals who sat briefly back-to-back in a canteen. Sequencing helped to refine estimates of incubation periods and attack rate, and revealed the order of transmissions in a subsequent household cluster60.

The ability of virus genomes to distinguish prolonged infection from cases of re-infection clarifies the reconstruction of transmission chains, and is crucial to understanding why some people repeatedly test virus positive. Similarly, co-infection with more than one virus phylogenetic lineage in a host at the same time could mask an international lineage introduction. Sequencing supported re-infection of an air traveller to Hong Kong (from Spain, via the UK) who had a high viral load and a B.1.79 lineage infection in August 2020; the same passenger had a B.2 lineage infection in March 2020 and was reverse transcription–PCR (RT–PCR) negative in mid-April 2020 (ref.61) (see also a similar case from the USA62).

Tracking lineages of interest

VOCs are genetic variants of SARS-CoV-2 that carry mutations that are known or suspected to affect key virus phenotypes such as increased transmissibility or immune escape. Phylogenetic analysis has revealed the independent emergence of VOCs, some of which share identical mutations (evolutionary convergence), and has reconstructed the accumulation of substitutions in time and space, shedding light on virus evolutionary or adaptive strategies.

The end of 2020 saw the discovery of the first VOCs, with multiple instances of convergent molecular evolution among them (Fig. 2; see the next section). For example, lineage B.1.1.7 (first labelled VOC-202012/01 and now termed VOC alpha) (Table 2) was determined by Public Health England to be a VOC on 21 December 2020 because its increase in frequency seemed to be related to the presence of particular genetic changes in the virus’s spike protein that had already been implicated in greater transmissibility (for example, N501Y and P681H) and antibody escape (for example, deletion Δ69/70)63. Lineage B.1.1.7 became dominant in the UK just a few months after its emergence, and phylodynamic studies showed it to have an estimated growth rate 40–70% higher than previous lineages25. In the global SARS-CoV-2 phylogeny, B.1.1.7 descends from the B.1.1 parental lineage via a long branch, suggesting that either the immediate ancestors of B.1.1.7 were unsampled or that the variant arose through a discrete evolutionary event during which multiple mutations were acquired, possibly during protracted infection of a single patient64. Slightly before the emergence of B.1.1.7, the N501Y spike mutation was detected in an independent lineage in South Africa. This lineage, B.1.351 (VOC-501Y.V2, now named VOC beta) also carried mutation E484K in the receptor-binding domain (RBD) of its spike protein15.

Fig. 2: The emergence of E484-bearing lineages from late 2020 to March 2021.
figure 2

Spike amino acid mutations and deletions are shown as symbols on the pins marking the approximate locations of first detection. The symbols include only those mutations that were implicated in possible immune escape or as suspected drivers of lineage growth and that were shared by two or more lineages. The locality of first detection may not be that of the lineage’s origin; however, the intercontinental spread of first detections is consistent with multiple independent origins. The B.1.1.7 lineage coloured in red differs from the other B.1.1.7 viruses in that it bears S494P rather than a substitution at E484. Lineage B.1.617 bears E484Q rather than E484K. Some lineages (B.1.1.7 and A.23.1) also have members that lack E484K, and some virus genotypes may have arisen multiple times (for example, B.1.1.7 with E484K). The near coincidental first detection of the same variants in genomes of phylogenetically distant lineages in countries worldwide, in early 2020, is a clear sign of convergent evolution and was a major factor leading to numerous studies aimed at detecting any selective advantage of the variants of concern (VOCs), including the search for vaccine escape phenotypes. Lineages and variants are based on the following publications: A.23.1 (ref.141); B.1.1.318, B.1.1.7 + E484K, B.1.1.7 + S494P, B.1.324.1 (ref.74); B.1.351 (refs15,71,74); B.1.525 (ref.74); B.1.617 (ref.142); P.1 (ref.73); P.2 (refs143,144); P.3 (ref.145). Note that B.1.324.1 was not designated as a sublineage of B.1.324, and reference here is to the variant described as B.1.324.1 in the Technical briefing Table 17 of ref.74. Pin heights indicate time relative to detection of the first lineage, that is, P.2 in Rio de Janeiro, 13 October 2020 (not to scale, but ranked in time, with days since detection of P.2 marked on each pin).

Table 2 Pango lineages of interest or concern during the first year of the pandemic

Phylogenetics can help to reveal the order in which variants accrue substitutions, which could provide clues to the functional advantages of convergent variants. For example, a phylogeny for the then emerging P.1 VOC (now named VOC gamma) indicated that the lineage’s characteristic mutations were gained in two phases, with a molecular clock analysis suggesting an intervening gap of several months65. Similarly, the nascent lineage B.1.351 detected in samples taken in South Africa during October 2020, lacked L18F, R246I and K417N; the latter substitution is among the nine changes that define B.1.351 and appeared in samples from the lineage in November 2020 (ref.66). Nevertheless, it is sometimes impossible to resolve the order of evolutionary events, because either genome sampling through time is insufficiently frequent or several mutations occurred very quickly. For example, ΔH69/V70 has arisen independently in several lineages (Fig. 3a) and is thought to compensate for decreased infectivity due to antibody escape substitutions such as N501Y; however, it is currently not clear whether or not the deletion preceded the RBD substitutions in B.1.1.7 (ref.67). The sudden appearance of lineages with constellations of 30 or so key substitutions relative to ancestral genomes is unlikely a priori given the low long-term substitution rate of SARS-CoV-2. The recent emergence of BA.1 (VOC omicron) has reignited interest in this phenomenon; evolution during a prolonged infection of an immunocompromised patient, or isolation within and then re-introduction from an unsampled human or animal population, are being considered as hypotheses for the origins of omicron68.

Fig. 3: Convergent evolution of SARS-CoV-2 spike protein.
figure 3

a | Phylogenies for the first year of the pandemic show the independent emergence of spike ΔH69/V70, indicated in red, in genomes of the B.1.1.7 and B.1.258 lineages respectively — note, the B.1.258 clade in red includes some branches without the deletion. Phylogeny from Nextstrain146,147 (which used data from the Europe ncov GISAID data set148), visualized in Figtree. Acknowledgements of authors responsible for the genetic sequence data generated, shared via the GISAID initiative and used to generate the Nextstrain tree, may be found in Supplementary Table 1. For clarity, not all Pango lineages are shown. b | By the start of 2020 several commonly occurring spike substitutions and deletions had been recognized as shared between lineages. The illustrated substitutions are found in the exposed (that is, outermost on the surface of the virion) subunit of spike, termed S1, or in the spike N-terminal domain (NTD), and are those shared by variants of interest or concern, excluding those shared sporadically or in minor sublineages. B.1.351 and P.1 share K417T/N and (in some B.1.351 sublineages) L18F, as well as two other recurrent substitutions; this is indicated by the overlap of their extended shading. ‘Mink’ refers to the severe acute respiratory syndrome coronavirus 2 (SARS-CoV-2) mink–human sublineage, termed ‘cluster 5’, which exhibited ΔH69/V70 and N501T (and other spike substitutions)91; the second B.1.1.7 lineage (VOC-202102/02, the grey ellipse with broken-line border) is a cluster of B.1.1.7 that also bears E484K71. N501T is a homoplasy that emerged in mink and may have transferred to humans; it is relatively uncommon, as it was found in only five mink in the original mink farm epidemic in Denmark. Nevertheless, N501T seemed to have emerged independently four times and has been detected in ten human cases149. L18F is an NTD substitution found in some B.1.351 and several of its sublineages, and it is increasing in frequency in B.1.1.7 (ref.67). As in Fig. 2, we see that the same substitutions appear in multiple lineages, implying that they arose independently at different times and places. Here, we also see that not only are individual substitutions shared, but constellations of several changes also seem to co-occur in more than one lineage; this suggests epistatic interactions, with perhaps compensatory changes following immune escape variants.

The E484K mutation in B.1.351 has been associated with antibody escape and potential resistance to convalescent plasma therapies63,69. In vitro, B.1.351 exhibits improved ability to escape antibody responses targeted at VOCs that arose earlier in the pandemic, such as B.1.1.7 (an escape phenotype attributed mostly to E484K and K417N)70,71, and shows increased transmissibility70. Although the B.1.1.7 lineage did not carry E484K when it first emerged, by 1 February 2021 this mutation had appeared in 13 English and two Welsh B.1.1.7 genomes. The phylogenetic relationships between these genomes suggested at least two independent acquisitions of E484K in the UK. Lentiviral and vesicular stomatitis virus (VSV) pseudotyping experiments indicate that the E484K mutation on the B.1.1.7 lineage backbone results in a reduction of neutralizing activity by vaccine sera71,72. The P.1 lineage was first reported in international travellers from Brazil entering Japan73 and showed 11 amino acid substitutions relative to its ancestral lineage B.1.1.28. Three of these substitutions fall within the RBD (K417T, E484K and N501Y), and all three sites are also modified in B.1.351 and some B.1.1.7 lineages63. P.1 seems to have originated in Brazil73,74 and also shows signs of increased transmissibility relative to its parental lineage B.1.1.28 (ref.75).

Although the phenotypic effect of mutations carried by VOCs can be investigated in vitro (see refs76,77,78,79 for examples), their epidemiological significance is harder to evaluate. Changes in mutation frequency during an emerging epidemic may not always directly reflect transmission potential or selective advantage, because they can also be influenced by founder effects, genetic linkage to other mutations, ascertainment bias and uneven sampling across regions33. Studies with a phylogenetic or phylodynamic basis have the potential to ameliorate some of these issues. The first amino acid replacement substitution to show a marked change in prevalence was D614G. Globally, SARS-CoV-2 with glycine (G) at spike position 614 rose from 10% prevalence before 1 March 2020, to overall global predominance by April 2020 (ref.77). Relative growth rates for D614G and other substitutions were estimated by phylogenetic diversification; this suggested that most variants were weakly deleterious and not more transmissible80. Sequence data from repeated international introductions of SARS-CoV-2 to the UK were leveraged to provide replicate observations of the growth of 614D and 614G lineages81. Modelling and phylodynamic analyses of 307 independent introductions between 29 January and 16 June 2020 suggested a genuine (that is, not a sampling effect) replacement of D by G in the UK, with a growth effect of around 20% and phylogenetic estimates of the basic reproduction number (R0) of 2.7–3.5 for 614D and 3.1–4.8 for 614G; however, indications of positive selection for 614G were not significant in all analyses. A separate analysis suggested that founder effects were responsible for the apparent selective advantage of 614G82, noting that the expansion of 614G coincided with a shift in the nexus of global dispersal from Asia to Europe.

Deep sequencing and phylogenetics have also been used to track virus evolution during co-infections83 (Fig. 4a) and prolonged infections84 (Fig. 4b) and to distinguish chronic infections from re-infections. In addition, phylogenetic assumptions regarding the distribution and independence of mutations can be violated by virus mutational patterns related to host antiviral defences85 (Fig. 4c). Within-host variation also has implications for SARS-CoV-2 phylogenetics, as co-infections may complicate tracing of transmission networks83 (Fig. 4d).

Fig. 4: Effects of within-host evolution and dynamics on epidemiological observations.
figure 4

Phylogenetic and phylodynamic approaches help detect and understand complex infections, measure within-patient lineage turnover and explore how host-induced mutation affects outbreak investigations. a | Co-infections may confound transmissibility and aetiological studies, but they can be detected using phylogenetics. Specifically, co-infections are identified when viral genomes sequenced from multiple isolates from the same patient are not monophyletic. b | Lineage turnover can occur if within-host lineages share a recent common ancestor and arise from evolution within the host itself. Lineage turnover may complicate patient treatment, as a lineage with lesser susceptibility to host immune responses may give way to a more transmissible lineage after apparently successful completion of a course of therapy. Nevertheless, phylogenetic features, such as longitudinal samples falling into different sister lineages and relative branch lengths, can help detect and account for lineage turnover. c | The antiviral activities of host APOBEC cytidine deaminases, which promote C → U hypermutation, adenosine deaminases that act on RNA (ADARs) and similar host systems, can lead to biases such as C → U homoplasies (convergent evolution) in the case of APOBECs, and changes in virus genome CpG content as a response150,151,152. Phylogenetics can highlight such convergent changes, which will be seen arising in lineages that are not closely related, and phylogenetic and phylodynamic approaches can be adjusted to account for the elevated rate of particular transitions. d | Co-infections and superinfections can complicate attempts to trace transmission chains, through either lineage turnover or sampling bias (for example, differential PCR amplification or through effects of organotropy). The result can be failure to connect two related transmission chains. A superinfected individual could also cryptically contribute to more than one heterochronous outbreak. The schema shows potential transmission events within households, or similar units (for example, workplaces), in a simplified transmission scenario. The dashed lines indicate transmission events between households. Circles represent individuals, with empty circles indicating infection chains involving lineage 1 and filled circles those involving lineage 2. The red asterisk indicates a co-infected individual who carries both lineages. The phylogeny shows that the true relationship between individuals X and Y may be unclear if lineage 1 dominates the co-infection at the time of sampling.

Homoplasy and recombination

Lineages bearing N501Y and E484K appeared independently in Brazil, South Africa, Canada and the UK in late 2020. Evolutionary convergence was observed, with the same changes being acquired independently on several branches scattered across the virus phylogeny (homoplasy) (Fig. 3a), and several lineages may share one or more substitutions (Fig. 3b). For example, both B.1.351 and P.1 (VOCs beta and gamma) showed escape-associated RBD substitutions at sites 417, 484 and 501 (Fig. 3b), as well as at positions 614 and 701 in the spike protein, but these two lineages do not share immediate common ancestry. The concurrent emergence and spread of the same mutations in different places and on different genomic backgrounds suggests that there were shared selective pressures acting on the virus86, such as the need to increase intrinsic transmissibility, extend the duration of infection or evade host immune responses (whether elicited by natural infection or by vaccination)87. The parallel emergence of constellations of functionally relevant mutations88 further suggests the existence of fitness interactions (epistasis) between them. Some mutations may only grow to a detectable population frequency if preceded, or closely followed by, a second permissive or compensatory mutation — several such mutations have been suggested in SARS-CoV-2 (refs67,89,90,91).

The epidemiological and phylogenetic context of these convergent changes indicates that they arose through independent, parallel mutation. However, it is known that such changes (homoplasies) can arise also through recombination, and evolutionary analyses suggest that recombination could now be relevant to SARS-CoV-2 evolution92. The level, scale and consequences of recombination during the pandemic are unclear; one earlier study of phylogenetic inconsistency found no clear signals of recombination93, whereas a more recent analysis of UK sequence data discovered at least four groups of natural recombinants of B.1.1.7 and other parental lineages94. The increasing co-circulation in 2021 of genetically diverse viruses increases the likelihood that further SARS-CoV-2 recombinants will be detected95.

Tackling sampling bias in genomic epidemiology

Uneven sampling of genomes has emerged as an issue for SARS-CoV-2 phylogenetics. Sampling was effectively absent during the first days and weeks, becoming more extensive as the pandemic progressed96, and often concentrated towards particularly large outbreaks or the radiation of VOCs. Some countries sequence routinely, others only for outbreak investigation and some not at all. The UK and Danish virus sequencing programmes are examples of large-scale, sustained sampling intensity, and they have allowed quantitative assessments of virus properties and public health interventions.

Although these and other countries have generated and openly shared many virus genomes97, even large-scale genomic programmes sometimes achieve only moderate sampling densities. For example, models of total UK infections to 4 May 2020 suggested that 3.4 million people were infected; of these, around 0.3% had their viral genomes sequenced98. The impact of low and uneven sampling intensity on SARS-CoV-2 phylogeography has been recognized25,99, and is commonly addressed through downsampling or a bespoke sampling regimen99. Some phylodynamic models allow for explicit modelling of sampling bias, and this has enhanced several SARS-CoV-2 studies12,25,100. In some cases, a single uniquely shared variant may be sufficient to determine the origin of a transmission chain60,101, and in studies in which one or several regions or periods were poorly sampled, unsampled individuals have been modelled and added to analyses102, or unobserved ancestral locations jointly inferred using phylogenies103. Delays between case detection104 and genome sequence availability, as well as insufficiently dense sampling, can hinder outbreak analysis. Further study of sampling effects is required and there is a need for standard definitions of sampling schemes to minimize bias in large-scale analyses while maintaining practicality99.

Ascertainment bias towards symptomatic cases potentially complicated attempts to determine any greater transmissibility conferred by spike 614G over 614D81, and undersampling in a region of high incidence has been suggested as a cause of overestimation of size and duration of SARS-CoV-2 transmission chains28. A generalized approach to assessing the impact of incomplete sampling on outbreak analyses has been proposed, based on the probability of correctly detecting linked transmission pairs given a certain level of sampling effort105. Genetic detection of viruses in wastewater has the potential to augment genome sampling before first-case detection and to estimate the impact of NPIs on local and regional prevalence106. Recently, SARS-CoV-2 sequences from wastewater samples have indicated genomes with novel constellations of mutations not seen in any known circulating lineage or VOC, including ‘cryptic lineages’ with multiple immune escape mutations107, with some studies using phylogenetics to place cryptic wastewater lineages in relation to known lineages108.

Uneven sampling through time can be addressed by adding an explicit sampling model to phylodynamic inference. One current solution uses structured (epoch-based) models to condition on the rate of genomic sampling relative to all PCR-confirmed SARS-CoV-2 cases, and reportedly improves molecular clock accuracy100. Methods that can accommodate changing rate of sequencing through time have been developed, for example, the coalescent-based Bayesian Epoch Skyline Plot (ESP)109 (see also Box 1), an approach analogous to the classic BD-skyline110. An alternative is to model sampling while linking sample location to regional variations in sampling effort; this has improved estimation of population size history for at least some data sets111, and the BD-skyline with variable sampling rate has also been applied to SARS-CoV-2 (refs100,112). The relationship between genetic variation and transmission patterns is one of interdependence, and therefore combining phylodynamic estimation with epidemiological data should generate stronger inferences. There has been notable progress on such integrated approaches. One recent method113 allowed the incorporation of non-genomic incidence data and epidemic dynamics models with a novel phylodynamic approach that represented both original and downstream members of transmission chains (that is, phylogenies with extant internal nodes). This joint epidemiological and phylodynamic approach is reportedly less susceptible to bias arising through undiagnosed cases, imported cases and changes in sampling levels, and so produces more reliable estimates of transmission rates than epidemiological data alone113. Analytical methods for a priori estimation of appropriate sampling intensity, sizes and strategy for analyses in virus genomic epidemiology are urgently required but not well developed. Such methods would reduce sequencing costs for longer-term initiatives as well as help ameliorate sampling bias105. Guidance is being developed to ensure that project objectives are considered and addressed using economically efficient genome sampling and sequencing approaches114.

Undetected SARS-CoV-2 transmission has been incorporated into standard Bayesian phylogenetic inference using an epidemiological model that includes data on confirmed but unsequenced cases and combines molecular sequence, case count data and temporal information. The model can infer undetected transmission or track changes following interventions and can also incorporate changes in sampling strategy, such as a decision to begin testing asymptomatic individuals112. The approach has been successfully used to infer the R0 and cumulative case count trajectories for the Diamond Princess cruise ship outbreak, a closed system for which reliable epidemiological (non-genetic) data are available to validate corresponding phylodynamic estimates112,115. These, and similar efforts to unite phylogenetics and epidemiology, are promising tools for the study of virus epidemics, including SARS-CoV-2 (ref.116).

Some phylodynamics approaches, especially those based on Bayesian inference, require substantial computation time, and population genetics could provide simpler evaluations of estimated population size117. One such method applied to SARS-CoV-2 used the relationship between number of haplotypes and number of sequences deviating from a reference to estimate effective population size (Ne) from sequence data118. The size trajectories estimated were similar to those inferred using phylodynamic methods. Such methods may be useful for rapid assessment, as they can be computed within high-throughput pipelines or in simulations. Subsampling of large data sets in virus genomic epidemiology has been a popular solution to reduce computational cost and to ameliorate differences in sampling intensity between countries. Although a few studies have explored how to optimize downsampling for parameter estimation (for example, see ref.119), there is a current lack of formal methods for downsampling model selection and implementation, and alternative strategies have been discussed99.

Further innovation is nevertheless required, particularly concerning the estimation of large virus phylogenies. Several studies have noted that the comparatively low substitution rate of SARS-CoV-2 reduces phylogenetic signal, potentially hampering studies of events early in the pandemic or of local epidemics120, hindering estimation of substitution rates121 and lowering phylogenetic resolution122,123. The effect was seen in the accumulation of sequence data for Washington State that overturned the initial suggestion of prolonged cryptic transmission in the state, by linking the first and second outbreaks there34, to favour independent international introductions as founding the later outbreaks102. Solutions to these problems include the use of sets of plausible trees, rather than a single tree120, testing of alternative root-placements123 and randomization tests for phylogenetic signal120. Consideration is also being given to rapid generation of maximally stable topologies, from multiple studies based on different data, methods and assumptions, implemented using entropy-weighted tree distances to highlight the least stable clades124. Some methods in virus genomic epidemiology have inherent assumptions — such as negligible variation within patients and absence of superinfection — that may not hold for SARS-CoV-2. The analyses being applied to the virus matured throughout the first year of the pandemic, and solutions arose from across diverse biological science disciplines, often in a highly collaborative manner. For example, approaches to quasispecies deconvolution were adopted from practices in oncology125,126. In addition, there is a vast literature devoted to the reduction of technical error and improvement of genome sequence quality127,128.

Conclusions and the way forwards

The contributions of evolutionary analyses to the global pandemic response are substantial and varied. The first year of the SARS-CoV-2 pandemic highlighted the progress that has been made over the past decade in virus genomics and phylodynamic analyses, while revealing technical and social challenges that remain to be addressed. The rapid, open sharing of protocols and data has been critically important, and more extensive for SARS-CoV-2 than ever before, yet hesitancy to share sequencing data before publication remains104 because of concerns that data may be used elsewhere without appropriate credit being given to producers129. Greater insights into SARS-CoV-2 transmission could be gained through the incorporation of more and varied data (for example, mobility data); however, this must be balanced with privacy and anonymization concerns. Flexible and robust methods for incorporation of diverse metadata into phylodynamic analyses are also required, as are standards for their collection and availability99.

In addition, the nomenclature of lineages and variants was initially inconsistent; this complicated scientific discussion, and encouraged the media to adopt simple but inappropriate naming of lineages based on the location of their first detection (for example, ‘South Africa variant’)130. The problem of toponymic naming in the popular literature has been partly overcome by the adoption of Greek letter designation for VOCs and VOIs by the WHO, with the Pango nomenclature adopted by researchers requiring a systematic nomenclature or for epidemiologically relevant lineages. Nevertheless, some confusion can still arise between the possible naming of recurring constellations of variants by the WHO, and their phylogenetic context as indicated by a Pango designation131.

In many countries, current research recruitment, evaluation and funding frameworks disincentivize the long-term participation of researchers with phylodynamic analysis skills in public health surveillance and control, because such participation diverts from those activities that are used to evaluate career progress (for example, research publications and grants)132. Consequently, new career pathways or evaluation systems are required to encourage greater embedding of evolutionary genomic approaches in public health. Investment in the training and retention of those with bioinformatic and phylogenetic expertise is required in many low and middle income countries, where the capacity for computational analysis sometimes lags behind that for genetic sequencing114. Further investigation into these ethical and technical challenges is needed to prepare for future pandemics and to sustain our tracking of SARS-CoV-2, transmission, new VOCs, new recombinants and cross-species transmission events.

Phylodynamics has demonstrated the impact of interventions and highlighted cases where they could have been applied more effectively or their use better timed. Phylogenetics has distinguished local onward transmission from new introductions and thereby informed infection control and planning. The history of pandemic transmission is recorded in virus genomes, allowing a global overview of virus epidemiology to be obtained even with samples taken in limited geographical areas or unevenly through time. Accordingly, phylogenetic concepts are likely to continue to play an important part in efforts to combat SARS-CoV-2 and in the prediction of the virus’s next move.