Oxford Nanopore MinION Sequencing Enables Rapid Whole Genome Assembly of Rickettsia typhi in a Resource-Limited Setting

INTRODUCTION

Until recently, whole genome sequencing (WGS) has been the preserve of high-income settings. Although the costs of WGS have dramatically decreased over the past decade,¹ the need for initial investment in sequencing platforms and associated equipment and for the technical expertise to run and maintain them have prevented their introduction into lower income settings. Downstream processing of sequencing data is frequently hampered by poor Internet connectivity, a lack of local expertise, and the perceived requirement for substantial computational infrastructure.

Oxford Nanopore Technologies’ (ONT) MinION sequencer is a portable device for DNA and RNA sequencing that generates data for local analysis, in real time. The MinION weighs under 100 g and plugs directly into a laptop via a USB port, with no additional computing infrastructure required.² Only basic laboratory facilities are needed to extract DNA and prepare sequencing libraries. The MinION has been used to sequence viral,³ bacterial,^4,5 and eukaryotic^6–9 genomes. Despite its high raw error rate, nanopore data can in many cases produce highly contiguous assemblies. Illumina short read data can be added to improve consensus accuracy, in many cases to that of a high-quality draft assembly, without extra finishing steps.

The Lao People’s Democratic Republic (Laos) is a land-linked country of ∼7 million people in Southeast Asia. Lao People’s Democratic Republic is ranked 135th in the world on the Human Development Index, a composite statistic of life expectancy, education, and per capita income.¹⁰ There are few functioning molecular biology laboratories in Laos, no previous published WGS projects, and relatively poor Internet connectivity.¹¹

Rickettsia typhi is an obligate intracellular, Gram-negative bacterium in the family Rickettsiaceae that causes the disease murine typhus and has a worldwide tropical and subtropical distribution.^12,13 The pathogen is transmitted to humans primarily by flea fecal contamination of the bite of the oriental rat flea Xenopsylla cheopis.¹⁴ The disease is an important and grossly under-recognized global cause of febrile illness.¹⁵ Those affected suffer from high fever, headache, myalgia, arthralgia, nausea, and vomiting and may have a macular rash. Complications are infrequent but include myocarditis, meningoencephalitis, seizures, and renal failure.^13,16 When diagnosed, recovery is typically rapid after treatment with doxycycline and mortality is low (around 1–2%) with antibiotic treatment.^12,17

Despite its importance as a pathogen, to date just three R. typhi whole genome sequences, with wide geographic and temporal distribution, have been published. The type strain Wilmington was isolated from a patient in North Carolina in 1928, and the sequence was published in 2004.^18,19 Two further sequences originating from a patient in Northern Thailand collected in 1965 and the other from a bandicoot rat, Bandicota sp., collected in Burma (now Myanmar) in the 1970s^20–22 were published online in 2014 (Table 1).

We report the successful WGS of a human isolate of R. typhi using only the MinION platform in Laos. We assess the potential for use of ONT data alone to perform comparative analyses and discuss the challenges of undertaking WGS in resource-limited settings.

MATERIALS AND METHODS

RESULTS

The MinION flow cell and reagents were shipped from Oxford, United Kingdom, to Laos at +4°C. Flow cell quality control before dispatch recorded 1,242 active pores and on receipt after 72hours travel, 1,103 active pores were recorded as available. Sequencing using the MinKNOW platform was performed for 48 hours and generated approximately 250,000 reads, with a total fast5 file size of 20 GB. In total, 222,848 reads passed quality filters and were used in further analysis.

We assessed the species composition of the reads that passed quality filters. Because the R. typhi sample was grown in L929 fibroblast cells, we expected to see some contaminating reads from the mouse genome. In all, 134,802 (60%) of the reads were classified as belonging to the genus Rickettsia, 126,447 reads were classified as R. typhi at the species level, and 41,428 (19%) were classified as belonging to the genus Mus (Figure 1).

View Full Size

Pavian visualization showing the proportions of Oxford Nanopore Technologies′ reads assigned up to the genus level³².

Citation: The American Journal of Tropical Medicine and Hygiene 102, 2; 10.4269/ajtmh.19-0383

• Download Figure
• Download figure as PowerPoint slide

The FASTQ files produced by Albacore were assembled using Canu, an assembler designed for noisy single-molecule sequencing such as that produced by the MinION. This produced a single-contig assembly of 1,078,916 bp. Nanopolish, a tool to use the signal-level data from Oxford Nanopore sequencing, was used to improve the consensus sequence of the Canu assembly, producing a new assembly of 1,099,322 bp. To further correct errors in the assembly produced from long reads alone, we polished the Canu + Nanopolish genome assembly using Pilon, which uses short read sequencing data to correct errors. We repeated rounds of mapping and correction with Pilon until no further errors were corrected. After four rounds of correction, Pilon did not correct any further errors, with the exception of two short indels which were removed and then reinserted by successive rounds of polishing. In total, 94.2% of errors were corrected by the first round of polishing, with a further 4.8% corrected in the second round (Supplemental Table 1).

After this process, a total of 6,082 errors were corrected, with over 95% (5,856) being small insertions and deletions, with the result that 12 kb was added to the genome during polishing for a final length of 1,111,939 bp. In comparison, when Pilon was run on the Canu draft genome before Nanopolish polishing, 19,970 errors were corrected, confirming that Nanopolish improves the draft genome.

We assessed the ONT read coverage and proximity to a homopolymer run (five bases or longer) for each of the errors (both single-nucleotide polymorphism [SNP] and indel) we corrected with Pilon. We observed no difference in the coverage distribution at sites with errors compared with the overall coverage distribution (Supplemental Figure 2). However, the positions around homopolymer runs were enriched for errors (two-tailed Kolmogorov–Smirnov test statistic=0.019, P = 0.035) (Figure 2). Whereas 15% of the final polished genome is within 5 bp of a homopolymer run (165 kb), 43% of errors (2,641) are within 5 bp of a homopolymer run.

View Full Size

The observed distance from a homopolymer run for the errors corrected by Pilon in the draft genome, compared with the expected distribution if errors were randomly distributed across the genome.

Citation: The American Journal of Tropical Medicine and Hygiene 102, 2; 10.4269/ajtmh.19-0383

• Download Figure
• Download figure as PowerPoint slide

The Illumina data alone were assembled using Unicycler to compare with assemblies using ONT only and the combined ONT and Illumina data. The assembly generated by Unicycler had 29 contigs, with a total length of 1,161,002 bp and an n50 of 166 kb. However, this included 20 short contigs under 6 kb, which are closely related to the mouse genome by BLAST search. These contigs are identified by their low coverage depth relative to R. typhi contigs. After removal of contaminating mouse contigs, the Unicycler assembly has eight contigs with a total length of 1,111,455 bp.

Furthermore, BUSCO was used to quantitatively measure the genome completeness of our assemblies based on a set of 48 conserved single-copy genes from the OrthoDB v9 database³³ (Table 2). None of the genes were found in the initial Canu assembly. Inspection of the TBLASTN results shows that while there were hits to many of the conserved genes, they fall below the threshold where BUSCO will call them as present. After polishing with Nanopolish using ONT data only, 14.2% of the genes could be found as a complete copy, with a further 32.4% found as fragmented genes. Two rounds of Pilon polishing improve this to 88.5% of genes found as complete copies. Running BUSCO on the Wilmington reference strain gives the same number of complete genes, suggesting this is the complete set of the orthologous reference genes that exists in R. typhi. The Unicycler assembly also contains this complete set of genes.

The final polished genome assembly is 1,111,939 bp, extremely similar to the three available genomes, which range from 1,111,496 to 1,112,957 bp in length (Table 1). Alignment of the four whole genomes in Mauve³⁴ shows that the genomes are colinear and that no rearrangements have taken place between strains (Supplemental Figure 3).

The TM2540 genome contained 866 predicted genes, very close to the 865 genes predicted in all other strains. Roary was used to define the core genome (genes present in all strains at 90% identity) of the four R. typhi strains, which included 863 genes. Of the five accessory genes not found in all strains, two were found in three of the four strains and three were unique to a single strain. Two of the genes not found in all samples were the same gene (RT_RS01750 in strain Wilmington), but clustered separately because of a 150-bp insertion in the B991CWPP strain that caused the sequence to fall below the 90% identity threshold to cluster with the other three strains. In two of the five accessory genes, a single base insertion in one strain has created a premature stop codon, splitting a gene into two sequences in that strain (RT_RS02175 in strain Wilmington, and TM2540_RS03925 in strain TM2540). In the final case, the gene is present in three strains next to a duplicate pseudogene which has lost the start and end of the gene, whereas in the fourth strain, a deletion has removed the end of the pseudogene and the start of the gene, leaving two incomplete and presumably inactive copies of the gene (RT_RS03370 in strain Wilmington). All five variably present genes are annotated as hypothetical proteins in all strains.

SNPs and indels in the four genomes were called from the alignment of the four reference genomes. We found 32 SNPs and 30 small indels (size range 1–782 bp). A full listing of the variants and their location in the Wilmington strain genome is given in Supplemental Table 2. Table 3 shows the number of differences between strains (counting each indel as a single event).

One of the small indels is in the RARP2 gene, which was previously noted to have a variable number of ankyrin repeats. Our new strain TM2540 has the same variant previously found in B9991CWPP and TH1527 and differs from the reference Wilmington strain. We also note a deletion in penicillin-binding protein 4 (pbpE) in TM2540 which leads to a truncated protein product. Variants are also present in the surface cell antigen (sca) gene families, with SNPs in sca1 and sca3 and indels in sca2 and sca4.

We compared the results of the final polished assembly to the best assemblies generated using Illumina data alone (filtered Unicycler assembly) or ONT data alone (Canu + Nanopolish assembly). The Illumina assembly has 866 predicted genes, and the core genome using this assembly has 862 core genes and seven accessory genes. This includes one gene unique to the Illumina assembly, which is part of a gene split into two pieces at a contig boundary. The ONT assembly has 1,745 predicted genes, but this number is inflated by a number of genes which have been split into multiple pieces likely because of sequence errors introducing premature stop codons in the sequence, and the core genome using this assembly has only 774 genes.

DISCUSSION

WGS of pathogens is an increasingly important tool in both research and clinical settings, whose use has increased over the last decade as its cost has decreased and availability improved. However, the initial investment in sequencing and computing equipment and the technical expertise needed to produce and analyze data remain a barrier to its introduction into resource-poor settings. Although our research laboratory can perform molecular diagnostics and cell culture, to date there has been no sequencing capacity, necessitating the slow and costly shipment of samples to other countries to undertake sequencing projects. The portable nature of the MinION sequencer has allowed the first bacterial WGS to be performed in Laos. We were able to run all steps of our bioinformatic analysis (except Centrifuge) on a laptop computer (MacBook Pro 2017) without continuous Internet access, demonstrating that sequencing and analysis can take place in relatively remote settings.

Murine typhus is an important and severely neglected tropical and subtropical disease of worldwide distribution. In Vientiane, an estimated 10% of non-malarial fevers in adult inpatients are caused by R. typhi.³⁵ Surprisingly, for a disease to which many millions of people are potentially exposed, there are only three published whole genome sequences, of which the last human isolate dates back more than 50 years.

The R. typhi genome is 1.1 Mb and contains almost no repeats, which makes genome assembly for this species relatively simple. Running Canu with the default ONT parameters was sufficient to assemble the genome into a single contig, and we polished the assemblies using both Nanopolish, to use the signal-level data from ONT sequencing, and Pilon, to use the extra information available from Illumina sequencing with its lower error rate. Although its synteny with other genomes suggests that there are no large misassemblies in our genome, the two methods of polishing corrected many small SNPs and indel errors in the ONT only draft genome. Thus, even for relatively simple bacterial genomes, a combination of ONT with data from another technology with a low error rate is currently necessary to produce an accurate sequence. In the case of R. typhi, the strains are so closely related that distinguishing between them requires a highly accurate sequence, and using ONT data alone, the sequencing errors are more numerous than the true differences between strains. Removing regions of the genome close to homopolymer runs would remove a large number of positions with errors by filtering out a comparatively small amount of the genome, but would still leave over 2,000 more errors than true differences. Oxford Nanopore Technologies data alone may suffice for applications which can work with sequences which still retain some base pair–level errors, such as determining large-scale genome rearrangements,³⁶ determining the species of an unknown sample,^37,38 or detecting antimicrobial resistance genes.³⁹ Increases in sequencing accuracy, combined with bioinformatics advances in assembly and polishing software, may allow for the future use of ONT data alone to give complete and accurate sequences, increasing the applications for this device.

Comparative analysis of the four R. typhi genomes shows very little variation between strains. Very few genes are variably present between strains, and the gene order is completely conserved. Although these four samples were collected from humans and a rodent over a long timeframe, we see only 23 SNPs between the two most divergent strains, no large-scale genome rearrangements, and very few indels. The low number of differences suggests that WGS may not be a useful tool to analyze the relationships between strains of R. typhi, as there may be few differences between strains and no recent shared ancestry. Recent phylogenomic studies of multiple rickettsial species have identified very little divergence within species.^40,41 Despite this small amount of variation, there are SNPs and indels which are private to each strain and could potentially form the basis of a molecular typing scheme to differentiate lineages of R. typhi. There are also variants in four of the five sca family genes present in R. typhi, which are important in adhesion and invasion of mammalian cells,⁴² and variation may be expected in these genes because of the selective pressure they are under as cell surface antigens; however, they are also large proteins and may be more likely to accumulate variants because of their size.

Some limitations remain on the use of ONT’s MinION platform. Currently, individual flow cells remain costly at USD500–USD900,⁴³ depending on the number purchased. However, the cost per sample can be reduced because multiplex sequencing now allows for multiple bacterial strains to be sequenced on the same flow cell.⁴⁴ Ongoing improvements in the technology continue to increase overall output and per-run capacity. Costs for setting up MinION sequencing facilities such as an appropriate laptop computer, shipping reagents, and flow cells are expected to be low. Shipping DNA or RNA at −80°C may cost more than USD1,000 per shipment and frequently requires supporting documentation and permissions. These costs remain a fraction of those required to set up and maintain non-portable sequencing platforms. Novel technologies are likely to improve access to sequencing when computational resources and support are limited and allow more cost-efficient sequencing of single samples.⁴⁵

MinION instruments continue to be used in an ever-widening set of challenging circumstances. The platform has been used to detect and monitor nosocomial tuberculosis in Zambia,⁴⁶ to sequence Ebola and Zika virus in Guinea and Brazil^3,47 at field sites, to sequence amphibian DNA in montane rainforest in Tanzania,⁴⁸ and for off-line use in the high Arctic and on the International Space Station.^49–51 Here, we demonstrate the current applicability of the MinION to resource-poor settings where some laboratory infrastructure exists, but WGS capacity is unavailable. With concerns in some countries about the export of biological samples for WGS in other countries, MinION systems could facilitate countries without current WGS facilities to undertake such work.

Supplemental tables and figures