We sought to improve the coverage evenness of multiplex PCR panel targeting SARS-CoV-2 by increasing the size of amplicons. To find a favorable size of amplicons, we designed and synthesized three two-pool panels respectively consisted of ~1,000-bp, and ~2,000-bp, and ~3,000-bp amplicon (Additional file). One nasopharyngeal swab sample was used for preliminary validation of the three panels. The multiplex PCR products were sequenced with MinION. The result shows that both the panels of 1000-bp and 2000-bp amplicons could obtain full viral genome, but the amplification of the 3000-bp panel failed due to that most RNA/cDNA degraded to short fragments. Then, the coverage profiles of the 1000-bp and 2000-bp panel were compared. The 2000-bp panel had a much larger bias of coverage among amplicons than the 1000-bp panel (Figure 1A). For the 2000-bp panel, 30.6% of sequencing data were assigned to one amplicon. Therefore, we selected the panel of ~1000-bp amplicons for the following analyses, referred to as 1K-panel. The 1K-panel contains 36 primer pairs in two pools (18 pairs each), and the sizes of the amplicons ranged from 880-bp to 1027-bp with an averagely 112-bp overlap.
We next compared the evenness of viral coverage by using the 1K-panel and a widely used 98-plex primer panel provided by the ARTIC network (http://artic.network/ncov-2019, version V3). RNA was extracted from six nasopharyngeal swabs with varied SARS-CoV-2 titers (Ct values, 22.9 - 31.0), and aliquots were respectively amplified with the 1K-panel and the ARTIC panel, followed by sequencing with MiSeq. We defined the viral genome region with a >20% mean depth as the high-coverage region, and the proportion of the high-coverage regions was used to quantify coverage evenness.
We found that the 1K-panel generated a more even sequencing coverage than the ARTIC panel (Figure 1B). With the 1K-panel, the proportions of the high-coverage regions were averagely 93.0% (SD=1.7%) for the six samples; whereas it was 80.6% (SD=8.7%) for the ARTIC panel. The coverage evenness was dependent on the viral titer of samples, and we found that 1K-panel was less affected by low viral titers than the ARTIC panel (Figure 2A). As the Ct values increasing, the proportions of the high-coverage region slightly decreased from 95.3% to 90.5% for the 1K-panel, compared with 88.1% to 62.4% for the ARTIC panel. We also evaluated a higher cutoff (>30% mean depth) to define the high-coverage region, and the 1K-panel maintained the advantage (Figure 2A).
Next, we included another 49 nasopharyngeal swabs to further evaluate the performance of the 1K-panel to recover the viral genome. We observed a comparable efficiency of viral amplification and genome coverage evenness (Figure 2B). Among 29 of the 49 samples with a <30 Ct value, except one sample, all had a viral genome recovery with >90% high-coverage region (mean 93.8%, SD 8.1%). The evenness decreases for the samples with a >30 Ct value, but five of six samples with a 30-33 Ct value had >70% high-coverage regions.
The amplicon sizes of the 1K-panel exceeded the maximum read length of MPS, but are favorable for nanopore sequencing in which fragmentation of PCR products is not needed. Recent studies have shown that nanopore sequencing by using MinION sequencing could generate accurate consensus genome of SARS-CoV-2 (12, 16). Thus, we focused on assessing the level of intrinsic noise to identify viral single nucleotide variations (SNVs), especially the sub-clonal variations, via the 1K-panel amplification followed by MinION sequencing. We sequenced the PCR products of the 1K-panel from eight samples by using both MinION (ligation-based kit and R9.4.1 flowcell) and MiSeq (enzymatic fragmentation and pair-end 2×300 bp flowcell). A high depth coverage was obtained for both sequencing devices (Table 1) which is required to identify sub-clonal SNVs. Then we implemented a home-mode bioinformatics workflow based on piled sequencing bases (see Methods) to identify SNVs. We included the SNV with a mutated allele frequency (MuAF) larger than 0.2 for the assessments. With this cutoff of MuAF, the MiSeq-based identification of SNVs should be reliable as previous studies shown [17,22].
The SNVs respectively identified with MinION and MiSeq sequences are shown in Figure 3A. Based on the distributions of MuAFs of these SNVs, we found that, with a MuAF≥0.8 threshold to identify SNVs, 96.8% (30/31) SNVs based on MiSeq could be identified by using MinION. One SNV (G26526T, in sample S1) was missed due to a much reduced MuAF (MiSeq 0.95, MinION 0.51). No false positive SNV with a ≥0.8 MuAF was observed through MinION sequencing, which is consistent with previous studies [12,16]. However, MinION-based SNVs with a lower MuAF were highly error-prone. Benchmarked by the MiSeq-based SNVs with a 0.2-0.8 MuAF, the MinION-based SNVs had a recall rate of 100% (4 of 4) and a precision as low as 12.9% (4/31). We illustrated these artificial SNVs and the true positives along the viral genome (Figure 2B). Most artificial SNVs were shared among samples, which is consistent with the feature of systematic error of nanopore sequencing [23,24]. No obvious enrichment of artificial SNVs at certain viral genome region was observed. The result suggests that the sub-clonal SNVs could not be reliably identified based on MinION sequencing.
As the 1K-panel could generate an even sequencing coverage on the SARS-CoV-2 genome and no need to detect minor viral populations, a large amount of MinION sequencing data per sample (about 300 to 500 Mbp) seems over-abundant. Therefore, to test the limit of the 1K-panel to recover SARS-CoV-2 genome, we re-sequenced the eight samples with an ultra-low throughput by using MinION (Flongle R9.4.1). We obtained 1.21 to 6.14 Mbp data per sample, the majority of viral genomes were recovered with sequencing reads (Figure 4). With a cutoff of ≥5-fold sequencing depth, 29 (93.5%) of the 31 high-MuAF SNVs could be identified, and no false positive SNV was observed. The two SNVs (T11158C and T28144C) failed to be identified due to insufficient or no coverage in sample S7, which had the lowest 1.21 Mbp sequencing data. We furthermore examined the coverage evenness. In 7 of the 8 samples, >90% (90.5% - 96.6%) of the viral genome had high coverage, and for sample S8 the proportion was 81.6%. Thus, the coverage evenness by using the 1K-panel and MinION sequencing was not affected by the extremely low data input.