(C) PLOS One This story was originally published by PLOS One and is unaltered. . . . . . . . . . . COVID-19 infection and transmission includes complex sequence diversity [1] ['Ernest R. Chan', 'Institute For Computational Biology', 'Case Western Reserve University', 'Cleveland', 'Ohio', 'United States Of America', 'Department Of Population', 'Quantitative Health Sciences', 'Lucas D. Jones', 'Department Of Molecular Biology'] Date: 2022-11 SARS-CoV-2 whole genome sequencing has played an important role in documenting the emergence of polymorphisms in the viral genome and its continuing evolution during the COVID-19 pandemic. Here we present data from over 360 patients to characterize the complex sequence diversity of individual infections identified during multiple variant surges (e.g., Alpha and Delta). Across our survey, we observed significantly increasing SARS-CoV-2 sequence diversity during the pandemic and frequent occurrence of multiple biallelic sequence polymorphisms in all infections. This sequence polymorphism shows that SARS-CoV-2 infections are heterogeneous mixtures. Convention for reporting microbial pathogens guides investigators to report a majority consensus sequence. In our study, we found that this approach would under-report sequence variation in all samples tested. As we find that this sequence heterogeneity is efficiently transmitted from donors to recipients, our findings illustrate that infection complexity must be monitored and reported more completely to understand SARS-CoV-2 infection and transmission dynamics. Many of the nucleotide changes that would not be reported in a majority consensus sequence have now been observed as lineage defining SNPs in Omicron BA.1 and/or BA.2 variants. This suggests that minority alleles in earlier SARS-CoV-2 infections may play an important role in the continuing evolution of new variants of concern. Evolution of the virus causing COVID-19 (SARS-CoV-2) has been associated with significant transmission surges. With evolution of SARS-CoV-2, evidence has accumulated regarding increased transmissibility of lineages, varying severity of illness, evasion of vaccines and diagnostic tests. Continuous tracking of SARS-CoV-2 lineage evolution distills very large and complex viral sequence data sets down to consensus sequences that report the majority nucleotide at each of over 29,000 positions in the SARS-CoV-2 genome. We observe that this eliminates considerable sequence variation and leads to a significant underestimation of SARS-CoV-2 infection diversity and transmission complexity. Additionally, concentration on the majority consensus sequence diverts attention from genetic variation that may contribute significantly to the continuing evolution of the COVID-19 pandemic. Funding: This work was supported by a Merit Review grant (CX001848) from the US Department of Veterans Affairs to CJD and by a grant from the US Department of Veterans Affairs Office of Research and Development as part of funding for VASeqCURE (grant number N/A) which in turn received funding from the American Rescue Plan Act funds. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. This is an open access article, free of all copyright, and may be freely reproduced, distributed, transmitted, modified, built upon, or otherwise used by anyone for any lawful purpose. The work is made available under the Creative Commons CC0 public domain dedication. Here we evaluate SARS-CoV-2 whole genome data from a sample of patients, administrative and medical staff experiencing COVID-19 in the VA Northeast Ohio Healthcare System. We have applied stringent inclusion criteria based on SARS-CoV-2 genome coverage and read-depth to normalize comparisons across sequences derived independently from individual patient isolates as well as patients who have had direct prolonged contact with other infected individuals. We examine multiple perspectives, including details of individual infection diversity and transmission outcomes which are critical for the successful management of the COVID-19 pandemic. The findings reveal details about SARS-CoV-2 infection and transmission complexity that are significantly under-reported. Data analyzed for lineage tracking relies on extensive curation and validation of SARS-CoV-2 genomic data. This includes minimization and resolution of ambiguous sequence data. In fact, most sequence data used to monitor the COVID-19 pandemic distills the whole genome sequence generated from an infection sample, to a single consensus sequence (the majority nucleotide present at each genomic position) [ 21 – 24 ] to national and international databases [ 25 , 26 ]. Out of the wide range of SARS-CoV-2 genomics studies, sequence variation within infected individuals has been described in relatively few studies as intra-host single nucleotide variants (iSNVs) [ 27 – 35 ]. Outcomes of these studies appear to be mixed in their assessments of individual infection diversity and transmission outcomes but inevitably reveal and acknowledge that these variations within a single infection exist. From earlier studies of the SARS-CoV-1 virus, it was suggested that the error-correcting capacity of the RNA-dependent RNA polymerase (RdRP) conferred by nonstructural protein 14 (NSP14) [ 12 , 13 ] would contribute to replication error rates significantly lower than other RNA viruses (e.g. HIV-1 and influenza virus) [ 14 , 15 ]. Studies comparing SARS-CoV-1 and SARS-CoV-2 estimate mutation rates that range from 1–3 × 10 −3 substitutions/site/year [ 7 , 16 – 18 ], respectively (extrapolating to 30–90 mutations/genome/year) are consistent with this prediction. Data provided through the Global Initiative on Sharing All Influenza Data ( GISAID ), and NextStrain provides a continuously updating report on SARS-CoV-2 lineage classification and illustrates that an apparently low mutation rate is no guarantee that the virus will exhibit limited capacity for variation, particularly once it has become so widely dispersed through millions of daily infections [ 5 , 19 , 20 ]. A further potential source of sequence complexity includes involvement of more than one variant in an individual infection. Early investigation into the outbreak of the novel pneumonia of unknown cause, from the city of Wuhan in late December 2019, quickly identified the genomic sequence [ 1 ] associated with a coronavirus, closely related to the SARS-CoV [ 2 , 3 ]. The virus, SARS-CoV-2, is now understood to be the causative agent of COVID-19. Sequence variation across the ~30kb positive-strand RNA SARS-CoV-2 genome was reported in the earliest studies, including single nucleotide polymorphisms (SNPs), insertions and deletions (indels) [ 1 , 4 , 5 ]. Significant efforts have since focused on characterizing sequence variation and documenting the evolution of SARS-CoV-2 strains during the pandemic through over 1.8 million sequences submitted from 172 countries as of April 2021 [ 6 ]. Extensive effort [ 7 – 10 ] has been and continues to be performed to define and track emergence of lineages across the ongoing evolution of SARS-CoV-2 populations around the world. This effort has captured and presented evidence of the continuous emergence of the major SARS-CoV-2 lineages sharing sequence variations across the global landscape during the pandemic. Many notable lineages including, B.1.1.7 (Alpha), B.1.351 (Beta), P.1 (Gamma), B.1.617.2 (Delta) have been associated with important surges in COVID-19 transmission and on-going phases of the pandemic [ 8 , 11 ]. With the surge of each new variant, concern naturally arises (and evidence accumulates) about strain-specific transmission and virulence, as well as the effectiveness of COVID-19 diagnostics, therapeutics, and vaccines. Finally, given published reports of extensive SARS-CoV-2 sequence diversity in association with infected immunocompromised (IM-) individuals [ 41 – 44 ], we have evaluated the complexity of SARS-CoV-2 sequence in the context of the immunological conditions of the infected individuals in our study. Because the infections studied included asymptomatic individuals and those who receive primary care outside the VA Northeast Ohio Healthcare System, it was not possible to conduct a complete association study on this focus and represents a limitation of our study. Additionally, we were not able to observe emergence of increasing numbers of SARS-CoV-2 mutations that have been described in longitudinal studies. Nevertheless, SARS-CoV-2 sequence meeting our analysis inclusion threshold was possible for 16 individuals with a range of clinical characteristics consistent with immunodeficiencies (including those classified as mild to severe by the CDC–see Methods , Patient Clinical Characteristics and Outcomes). Given the observations shown in Fig 1F , that genome-wide SNPs have increased during the pandemic, we compared SNPs detected in viruses from IM- with immunocompetent (IM+) individuals in relation to the time when their samples were collected. Results summarized in S2 Fig showed that there were no obvious differences in the number of variations harbored by the virus from IM- compared to IM+ patients. Interestingly, the individuals identified to carry a B.1.1.7/B.1.617.2 mixed infection were discordant with respect to their immune status (one IM- the other IM+). a Summary of bottleneck evaluations for 11 donor-recipient pairs. Dates of transmission include June/July 2020 –MT1-RN1, MT2-RN2, NA-MT2, NA-RN2; November 2020/January 2021 – P11-P13, P11-E30, P13-E30; December 3, 2020 –D-P1, D1-P2, D1-P3; January 23, 2021 – D2-P4. Individual results provide the mean predicted number of virions for each transmission event and 95% CI. b-d Transmission graphs showing AAFs between the 11 infected donor and exposed recipient pairs. Infection complexity is demonstrated by an array of AAFs from 0.05 to 1.0 (y-axes). Transmission settings include a VA Clinical Ward b ; a VA Emergency Department c ; two independent commuter transport vans d . Inset date for each transmission graph show the maximum change in allele frequency (Max AFD), the median allele frequency change (Med AFD), and the average allele frequency change (Avg AFD). Allele sharing assessments were also performed for two randomly chosen pairs who did not work together or have known contact—CLCCB1 & EDCB6: Max: AFD = 99.95%, Med AFD = 0.081%, Avg AFD = 0.221%; SURVCB26 & EDCB5: Max AFD = 100.00%, Med AFD = 0.077%, Avg AFD = 0.249%. Examination of the total infection sequence composition has also been considered in our study of infection transmission complexity from multiple exposure events. We performed bottleneck assessment using the exact beta-binomial method of Sobel-Leonard et al. [ 38 ], that recognizes AAF and RAF proportions instead of treating nucleotide sequence variation on a presence/absence basis taken by alternative models [ 39 , 40 ]. These analyses focused on 11 donor-recipient pairs, some that have been studied previously [ 33 – 35 ]. While donors were known from transmission events that occurred during commuter van transport events, we performed reciprocal calculations (A to B and B to A) for all transmission pairs. Results from our overall analyses show the maximum likelihood estimate (MLE) of the bottlenecks among our transmission pairs was between 6 to 123 virions ( Fig 5A ). These results include transmission events that occurred toward the end of pandemic Year-1 (December 2020 and January 2021) with a viral population characterized by greater genetic diversity compared to the first quarter of 2020 ( Fig 1F ). Paired SNP proportions in the “transmission plots” ( Fig 5B–5D ) show the AAF for each position in the SARS-CoV-2 genome between the respective infected donor (left) and subsequent recipient (right). In these graphs, most of the nucleotide-specific data are at or below the 5% alternate allele threshold and therefore would be designated the reference allele (>29,000 nucleotide positions). Results further indicate that among the majority of biallelic positions, the allele frequencies were transmitted in very similar proportions (AAF:RAF). As such, the average allele frequency differences (AFΔ) between donor and recipient pairs were minimal in all transmission scenarios, regardless of exposure details (0.09 to 0.17%; Fig 5B to 5D ). Counts of samples observed to show variation (among the 140 meeting analytical threshold) across the SARS-CoV-2 genome. Black histogram components include those positions where biallelic variation was observed at or below 50% and may not have been included in reported consensus sequences. Red histogram components identify positions and occurrences where the alternate allele was observed at proportions from 51–100% and would have been included in the reported consensus sequence. The multicolored bar below the graph identifies genomic segments encoding different SARS-CoV-2 genes designated in the legend (ORF–Open Reading Frame; S–Spike; E–Envelope; M = Membrane; N–Nucleocapsid). Quantitative details summarized in Fig 4 are available in S2A to S2B Table . In addition to the accumulation of increased polymorphism across the viral genome over time, closer examination of the genomic variations in the summary zebra plot (so-called due to the gray-scale patterns; Fig 3A ) revealed significant pattern shifts in the sequences that characterize the major variants in our sample population that correspond to epidemiological transition periods. Following the sample collection time series ( Fig 3A , April 2020 (top) to August 2021 (bottom)), four approximate variant signature time frames are observed—Pangolin lineages B.1 (no shading; April–November 2020), B.1.2 (yellow shading; November 2020 –March 2021), B.1.1.7 (Alpha; blue shading; March–June 2021) and B.1.617.2 (Delta; red shading; May–August 2021); lineage-specific references sequences compiled from frequencies of mutations deposited in GISAID are identified by lineage designations (e.g., B.1.2, B.1.1.7, B.1.617.2). Sequences within green and purple shades identify transition periods, B.1.2 to B.1.1.7 and B.1.1.7 to B.1.617.2, respectively. Of particular interest, the two exceptional samples identified in Fig 1D–1F (samples 106 –May 14, 2021; 110 –May 21, 2021) are found in the B.1.1.7 to B.1.617.2 transition period. In evaluation of their genomic sequence data, variants are observed with AAFs > 5% but not reaching major consensus, that are consistent with B.1.617.2-specific SNPs ( Fig 3B ). Following conventional practices, the major variants (with frequencies above 50% allele frequency cutoff) identified for these infections would have been reported as B.1.1.7 and not a mixed strain infection of B.1.1.7/B.1.617.2. Observation of these minor alleles above the 5% cutoff demonstrates early detection of the next major variant of concern ( Fig 1F ). A final summary of the distribution of SNVs observed in our study population is presented in Fig 4 . The data show that variation is concentrated among the genes encoding structural (mean = 4.2%; range– 1.3% to 7.7%) and accessory proteins (mean = 5.3% range– 3.1% to 7.4%), compared to non-structural proteins (mean = 3.1% range– 2.9% to 3.3%). This is particularly true for positions characterized by biallelic variation where AAF may not rise above 50% (Black components of stacked histogram bars); of which many are observed in the SARS-CoV-2 spike protein. Our observations of biallelic variant positions are consistent with those from other studies describing subconsensus SNPs or iSNVs [ 27 – 35 ]. To determine how often biallelic variations were included in SARS-CoV-2 sequences submitted to international repositories, we selected and analyzed 110 additional COVID-19 infections ( S4 Table ) reported from 18 different laboratories or consortia with raw sequence data available in the sequence read archive (SRA; total whole genome SARS-CoV-2 sequences available = 41,525) [ 36 ]; the same ≥80% coverage, ≥100X read-depth inclusion threshold was applied. Results summarized in S1 Fig confirmed that all SRA files from the 110 external infections showed evidence of multiple biallelic iSNVs that would reflect the complexity of the SARS-CoV-2 infection (nucleotide positions characterized by biallelic iSNVs in our VA study population and the samples queried in the SRA are designated by a green dot in Figs 2 and S1 ). Similar to observations from the VA Northeast Ohio Healthcare System patient population, 46 alternate biallelic iSNVs from a total of 87 variable positions did not clear the 50% consensus threshold. Among the 41 SNPs where AAFs did rise above the 50% threshold, the alternate allele would be reported in 668 occurrences and not reported in 312. Following the consensus sequence reporting convention across the 87 positions, 36% of the SARS-CoV-2 variations with AAFs between 5–50% (369 iSNVs) would not have been reported across the infections studied. Given that all SRA data queried showed this complexity, we asked whether attempts to report biallelic sequence variation were being made through the use of standard ambiguity codes [ 37 ] used to identify heterozygosity in diploid genetic systems. We first observed that the international data repositories understandably discourage the use of ambiguity codes given the substantial effort in attempting to curate sequence and determine lineage relationships as the pandemic continues to emerge. When we queried the database of 1,801,208 SARS-Cov-2 sequences accessible through Nextstrain/GISAID (available up to October 3, 2021), while the repository reports extensive nucleotide substitutions, results showed that 85.9% reported no ambiguity sequence codes suggesting no biallelic polymorphism in their SARS-CoV-2 sequences; 0.46% reported >10 suggesting potential for a complex multi-strain infection ( S6 Table ). Given our presentation that most, if not all, SARS-CoV-2 full genome sequences are characterized by the occurrence of multiple biallelic positions, we have integrated appropriate ambiguity codes at positions where the AAF is >5% and observed that sequences are characterized by a minimum of 12 and maximum of 132 biallelic ambiguities. In head-to-head comparisons, nucleotide-specific substitution was higher in Nextstrain/GISAID compared to our study (average substitution rate = 1.69, st. dev. = 0.52), while biallelic ambiguity code usage is 100-fold lower (average substitution rate = 0.0097, st. dev. = 0.018); overall identification of mutations and ambiguity usage is summarized in S5C to S5F Table . Due to the increasing levels of SNVs observed over the time of sample collection, we report the average, minimum and maximum numbers of biallelic positions in quartiles (April 1 to August 1, 2020; September 28 to December 9, 2020; December 10, 2020 to April 13, 2021; April 15 to August 18, 2021), corresponding to breaks in sample availability ( Table 1 ). a Zebra-plot alignment of 140 sequences from the Cleveland VA study population, collected from April 2020 to August 2021 at 200 lineage defining signature positions. The grayscale bar (upper right) indicates the proportion of the alternate allele at each biallelic variant position. Integrated into the alignment of individual patient infection sequences are the Wuhan reference sequence (NC_045512.2), and aggregate plots for B.1.2 (51,069 sequences), B.1.1.7 (479,770 sequences), B.1.617.2 (258,910 sequences), AY2 (1,938 sequences), AY7 (3,443 sequences), AY11 (785 sequences), AY19 (106 sequences), AY22 (70 sequences), AY3.1 (5,363 sequences). Lineage defining signature positions are defined as a position in which >50% of samples from that lineage exhibits the alternate allele. b Zebra-plot and AAF proportion assessment of a SARS-CoV-2 infections 106 (May 14, 2021) and 110 (May 21, 2021) showing evidence of a B.1.1.7(Alpha)/B.1.617.2(Delta) infection. Box and whisker plots illustrate the spread of RAF:AAF proportions for each mutation (median (bar); 25 th and 75 th percentiles (bottom and top of box, respectively); 10 th and 90 th percentiles (bottom and top whiskers, respectively)). The red dashed line at 50% identifies the demarcation at which the alternate allele is recognized and reported as the conventional consensus nucleotide that characterizes the infecting virus. These 75 mutations were chosen for further analysis because they all showed evidence of biallelic states at greater than the 5% background threshold and were shared in more than 10 infections (blue dashed line). Among these varying positions, Illumina sequence reads for 30 mutations did not reach an AAF:RAF >50% for any of the 140 infections (shaded gray in aligned table); 45 positions exceeded the RAF:AAF 50% threshold for a portion of the 140 infections and would therefore have contributed to characterization of the infecting viral strain. The accompanying table enumerates the number of samples with AAF between 5 and 95%, <50% and >50%. Overall, among the 5,390 alternate alleles detected at a frequency >5%, just 975 (18%) were observed at >50% frequency necessary to be reported in consensus sequence. Nucleotide positions characterized by biallelic iSNVs in our VA study population and the samples queried in the sequence read archive (SRA) ( S1 Fig ) are designated by a green dot. The variant positions across the SARS-CoV-2 genomic sequences detected with allele frequencies >5% reveal extensive biallelic polymorphism throughout our results ( Figs 2 and 3 ). We first examine these data by assessing alternate vs reference allele frequency proportions (AAF:RAF) for the 75 SNV positions that were shared most commonly (in ≥10 samples) ( Fig 2 ). The red dashed line in the figure at 50% allele frequency is the commonly accepted cutoff for variant detection indicating major allele. Only those alternate allelic SNPs with frequencies above 50% would be reported as part of the consensus sequence characterizing the virus present in the infection sample. Following this convention, 30 alternate allelic positions observed in our surveyed population would not have been reported in any of the infection-specific consensus sequences, totaling 2,865 SNP events. Among the 45 SNPs where AAFs did rise above the 50% threshold, the alternate allele would be reported in 1,143 occurrences (AAF >50%) and not reported in 1,553 (AAF 5–50%). Following the consensus sequence reporting convention across the 75 positions monitored, most of the SARS-CoV-2 variations (4,418 SNP out of 5,561, 79%) with alternate allele frequencies between 5–50% would not have been reported across the infections studied. These findings shown previously in Fig 1D–1F , predicting lower levels of polymorphism at a 50% vs 5% threshold, are consistent with this analysis. Enumeration of SNPs across the SARS-CoV-2 genes and a full listing of all SNPs observed in our study are presented in S1 Table (SNPs observed in ≥10 infections) and S2 Table (all SNPs), respectively. Alignment of the 140 SARS-CoV-2 consensus sequences ( S3 Table ) that met the coverage thresholds, revealed a total of 1,096 SNV positions. Among these polymorphisms, 406 SNPs were observed in multiple infection samples. a RT-PCR cycle thresholds predicting SARS-CoV-2 sequence read quality. Comparison of RT-PCR cycle thresholds for samples according to technical breadth and depth thresholds ((≥80% of SARS-CoV-2 genome covered at ≥100X read-depth; met thresholds, n = 140; did not meet thresholds n = 114). b SARS-CoV-2 full genome sequence depth comparisons. Among 140 sequences meeting our technical thresholds, the average coverage was 1,081 X (range 223 to 2,235 X); for the 114 sequences not meeting our technical thresholds, the average coverage was 71 X (range 0 to 1,281 X). Among all samples, average uniquely mapped read count was 292,430 reads (min = 21, max = 2,680,239). Among samples that pass >80% 100X coverage, the average uniquely mapped read count was 279,408 reads (min = 63,296, max = 2,680,239). Among samples that did not pass our cutoff, the average uniquely mapped read count was 175,895 reads (min = 21, max = 1,283,414). c Correlation between RT-PCR Ct value and average viral genome coverage; black circles identify samples that met the coverage and read-depth thresholds, and gray circles identify samples that did not meet these thresholds. d Correlation between detection of alternate alleles and RT-PCR Ct value where alternate allele frequency (AAF) was >5% (black circles >80%100x and open gray circles <80%100x), or where AAF was >50% (blue diamonds >80% 100x, open blue diamonds <80% 100x). e Correlation between detection of alternate alleles and average SARS-CoV-2 genome coverage where alternate allele frequency (AAF) was >5% (black circles >80%100x and open gray circles <80%100x), or where AAF was >50% (blue diamonds >80% 100x, open blue diamonds <80% 100x). f Correlation between the number of alternate alleles identified in the 140 >80% 100x samples and their collection dates—AAF >5% (black circles), or AAF >50% (blue diamonds). Significant correlation was observed between an increase in AAFs during the sample collection time-period (AAF >5%, Pearson R = 0.76112; AAF >50%, Pearson R = 0.89253; both P-values <0.0001). Red arrows in d-e identify samples that carried significantly more sequence variation compared to expectation (Sample 106–110 SNPs; Sample 110–132 SNPs). Green arrows in f link the position of the two samples with additional sequence variation at the AAF >50% and AAF >5% thresholds. Each sample was initially evaluated by RT-PCR to approximate levels of infection that would enable high-resolution evaluation of individual infection complexities. The Ct values for samples that resulted in SARS-CoV-2 full genome coverage (≥80% of the genome covered at ≥100X read-depth) were significantly lower (Average Ct value of 27.96 [n = 138]), implying higher viral RNA copy number than samples that did not meet this coverage threshold (Average Ct value of 35.85 [n = 103]) ( Fig 1A ; Mann-Whitney U = 1,971.5, P<0.0001). Coverage across the SARS-CoV-2 genome for samples with sequence generated above (black line) and below (gray line) the coverage threshold is summarized in Fig 1B and direct correlation between Ct and average coverage for each individual sample is provided in Fig 1C . Among the 140 (two samples did not have RT-PCR data) sequences meeting our technical thresholds, the average coverage was 1,081X (range 223 to 2,235X) with an average of 279,408 uniquely mapped reads per sample (min = 63,296, max = 2,680,239). Further examination of our data showed that the ability to detect SNPs was not influenced by Ct or average coverage if a minimum coverage was reached ( Fig 1D and 1E ). Samples that passed our inclusion criteria, ≥80% of the genome covered at ≥100X read depth, showed comparable number of SNPs defined by alternate allele frequency (AAF) > 5% (filled black circles) or SNPs defined by AAF > 50% (filled diamonds), regardless of Ct value. As expected, fewer SNPs were detected in samples that did not meet our inclusion (open circles and diamonds). With read depth of ≥100X, we observed consistency in detecting expected sequence variations within our sampled population and across genomic sequences reported into the international databases. From these analyses we established confidence in detecting minor alleles at frequencies ≥5% (whether reference or alternate allele). For the 114 samples not meeting our criteria for inclusion, the average coverage was 71X (range 0 to 1,281X) and these samples were removed from subsequent analyses. A comparison of the genome sequence data across the total sampling period of our survey showed a significant increase in the number SNPs detected over time ( Fig 1F and Table 1 ) with average coverage remaining consistent through this period; approximately doubling (30 to 65 AAF ≥5%) or tripling (10 to 35 AAF ≥50%) from April 2020 to August 2021. The two samples demonstrating substantially higher SNP counts ( Fig 1D–1F red and green arrows ) will be covered in their epidemiological context below. Full genome sequences from 364 SARS-CoV-2 infections (254 from present study; 110 from international data repositories) were analyzed in comparison to the Wuhan reference genome (GenBank NC_045512 ). A total of 254 samples from patients (N = 179) and healthcare personnel (N = 75) with COVID-19 were collected to (a) perform high-resolution contact tracing of COVID-19 infections and (b) gain insight into SARS-CoV-2 transmission dynamics between donor and recipient individuals [ 33 – 35 ]. Because of the infection control and outbreak monitoring nature of this work, the individuals represented a wide demographic range, including gender, age, race, and health status. Of 175 patients with available clinical information, the mean age was 61.4 years (range, 24 to 98), 31 (18%) were asymptomatic, 58 (33%) were hospitalized for COVID-19 or non-COVID-19 conditions, 32 (18%) had severe COVID-19, and 11 (6%) died due to complications of COVID-19. Twenty-four (14%) patients with COVID-19 were moderately or severely immunocompromised, including 17 with active treatment for malignancy with chemotherapy or severely immunosuppressive medications, 2 with advanced or untreated HIV infection, 1 with an organ transplant requiring immunosuppressive therapy, and 4 with other chronic medical conditions requiring immunosuppressive therapy. All patient samples have been analyzed as part of an effort to assess the extent of sequence variability within SARS-CoV-2 infections across different surges associated with variants of interest and concern from June 2020 to October 2021. An additional 110 samples from the Sequence Read Archive (SRA) have been analyzed for assessment of diversity outside the local region. Discussion Our study has revealed significant within host infection diversity of SARS-CoV-2 through full genome sequence analysis. We greatly appreciate the work of international data repositories such as GISAID and NextStrain for curating the millions of reported infections and believe that these repositories are a tremendous and invaluable asset for the monitoring of SARS-CoV-2. The shear amount of data that needs to be collected and annotated makes it understandable and perhaps necessary to represent these infections by a single consensus sequence. While the practice of reporting a single unambiguous sequence assists with tracking SARS-CoV-2 in regard to evolution and global distribution of lineages, emphasis on reporting the majority nucleotide signature significantly reduces our ability to detect emerging new variants that could ignite a new surge with public health impact. Additionally, the majority consensus approach limits our ability to understand the complexity of infections, dynamic changes in the viral population over time within infected individuals, and between donor and recipient individuals across transmission events. In this study, we addressed the majority consensus problem as it relates to SARS-CoV-2 infection through the use of standard ambiguity codes [37]. All consensus sequences reported to GenBank from this study have included ambiguity codes. Taking this approach, our study accurately identified two individuals who harbored mixed B.1.1.7/B.1.617.2 infections in mid-May 2021. The collection of these samples preceded by one month the first reports from the State of Ohio that B.1.617.2 was present in the Cleveland area and occurred at a time when the City lifted its Civil Emergency Proclamation (May 28, 2021) and the State of Ohio rescinded its masking mandate and social distancing orders (June 2, 2021) [45]. Communicating knowledge of B.1.617.2 infection and transmission in Cleveland ahead of relaxing mitigation measures may have contributed to a more cautious public health strategy. Additionally, nucleotides with biallelic mixtures observed repeatedly in multiple samples have included many of the functionally significant polymorphisms that have been reported in variants of concern. While widely distributed across the SARS-CoV-2 genome, biallelic mixtures are observed at amino acid positions within the SARS-CoV-2 spike protein’s receptor binding domain (e.g., L452R and T478K; B.1.617.2) that are prominent in lineage determination and associated with increased transmission efficiency (e.g., D614G). Moreover, there are numerous examples of alternate alleles observed in this study that did not clear the 50% threshold for inclusion in conventional consensus sequences that have subsequently become included as lineage defining SNPs in Omicron BA.1 (ORF1a:T3255I; S:T95I; S:K417N; S:H655Y; M:D3G; N:P13L; N:R203K; N:G204R) and/or BA.2 (ORF1a:T3255I; ORF6:D61L; S:T19I; S:G142D; S:K417N; S:H655Y; N:P13L; N:R203K; N:G204R) variants. These latter observations suggest that minority alleles in earlier SARS-CoV-2 infections may contribute to the continuing evolution of new variants of concern [46]. Confirming this potential will require consideration of minor alternate alleles as suggested here in more extensive whole genome analysis of SARS-CoV-2. In our assessment of transmission from infected donors to recipients we were also interested to see that the overall sequence diversity was largely retained. This appeared to be true regardless of the AAF:RAF mixtures (Fig 5B to 5D). These transmission plots provide both quantitative and qualitative perspectives on SARS-CoV-2 complexity. Many of the AAF:RAF mixtures were observed to be present in donor and recipient sequences at nearly identical ratios. For the events captured in commuter vans (two independent events) where the direction of transmission was known to be from infected drivers to uninfected passengers, this suggests that the majority of SARS-CoV-2 strain diversity is included during transmission and subsequent infection complexity [35]. This appeared to be true even though transmission occurred within the context of circulating air. Maximum likelihood estimates of the bottlenecks among our transmission pairs (van driver to passengers as well as others Fig 5A) were similar or higher than those observed by Braun [32] and Lythgoe [31], but did not reach the highest values reported by Popa [29]. Finally, given reports that have suggested evolution of complexed populations of SARS-CoV-2 in immunocompromised people [41–44], we were interested to compare infections between IM+ and IM- individuals. While normalizing for overall complexity by comparing infections that occurred within time-similar cohorts, we did not observe increased numbers of iSNVs in IM- patients. This appears to differ from the observations of Al Khatib et al. who reported increased numbers of biallelic subconsensus variants (indicating AAF:RAF mixtures) in individuals classified as severe cases compared to individuals with mild cases [27] with comorbidities similar to the IM- patients studied here. As we indicated previously, we acknowledge limitations in this part of our study. We were not able to study single patients over extended time series known to have extended for up to one year by other investigators [43]. Therefore, we did not have the opportunity to compare potential evolution of SARS-CoV-2 variants within individual IM+ and IM- patients over time. Our conclusions are limited by the fact that we studied a small sample size and were unable to control for potential confounding factors such as duration of infection, genetic differences, and COVID-19 treatment. The continuing dynamics of the COVID-19 pandemic have become increasingly complex and unpredictable. Underlying the challenge of COVID-19 transmission among humans, Sender et al. estimate that every infected person carries 1 billion to 100 billion virions during peak infection [47]. As a result, despite proof-reading activity within the SARS-CoV-2 RdRp that would limit mutation, efficient virus replication and transmission increasingly favors dispersal of mutations across abundant globally distributed infections (>460,000,000 to-date [19,20]). As levels of immunity stimulated by infection and vaccination fluctuate in the human population ((a) no infection or vaccine exposure; (b) previous infection(s) and no vaccination; (c) no infections and vaccination(s); (d) previous infection(s) and vaccination(s)), the global SARS-CoV-2 population will encounter significant heterogeneity in acquired and natural human host defense mechanisms. Given the global disparities influencing access to, and uptake of vaccines, as well as availability and practice of non-pharmaceutical interventions, the virus has a free range for making thousands of deleterious as well as advantageous mutations. Therefore, significant opportunity exists for random chance to determine what future chapters of the pandemic will present. This emphasizes the need to study genetic complexity of SARS-CoV-2 infections more intensively and requires increased examination to achieve a greater understanding of the capacity of this virus to evade our best efforts at mitigation. While there are efforts to sequence the virus and document how it is evolving, most of the scientific community recognizes that this effort needs to be increased [48–53]. Coupled with this effort, approaches to report characteristics of multi-strain complex infections must be considered to enable closer monitoring and full understanding of the virus’s evolutionary capacity. Lest it be thought that under-reporting of the genetic diversity of infections occurs only with SARS-CoV-2 genome sequence data, consensus sequence reporting is the submission practice for most (if not all) infectious disease pathogens to national and international data repositories (https://VEuPathDB.org) from the US National Institute of Allergy and Infectious Diseases (NIAID) and the Wellcome Trust (example pathogens and vectors include biosecurity pathogens (e.g. Yersinia pestis (plague), Marburg virus, Ebola virus; pathogens of global concern [e.g. Plasmodium, Mycobacterium tuberculosis). Adequately addressing the under-reporting of infection complexity represented in data repositories has significant potential to alter vitally important characteristics of infectious diseases, including assessment of drug treatment efficacy and resistance and vaccine escape. [END] --- [1] Url: https://journals.plos.org/plosgenetics/article?id=10.1371/journal.pgen.1010200 Published and (C) by PLOS One Content appears here under this condition or license: Creative Commons - Attribution BY 4.0. via Magical.Fish Gopher News Feeds: gopher://magical.fish/1/feeds/news/plosone/