A method for pathogen metagenomic analysis based on nanopore sequencing data

Abstract

The present application relates to the technical field of genome sequencing, and particularly relates to a pathogen metagenome analysis method based on nanopore sequencing data. In the first round of alignment, four kinds of rapid alignment software and algorithms, namely Kraken2, Centrifuge, Pandora and Minimap2, are used to realize preliminary identification of species for typical data (1Gb nanopore sequencing data). In the second and third round of verification analysis, the read of each species identified in the last round is cross-verified with the reference sequence of the species and other pathogen reference sequences identified in the sample by using the accurate alignment algorithm BLAST, so that the calculation amount of accurate alignment is effectively reduced. The Minimap2 is aligned with the database of drug-resistant genes / point mutations, and the drug-resistant genes / point mutations are identified by using the blast accurate alignment algorithm. Through the construction of the method and system, the accurate analysis and report generation of typical nanopore sequencing pathogen metagenome data within 30 minutes are realized, and the second-generation sequencing data analysis is effectively compatible.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of genomic sequencing technology, in particular to a pathogen metagenomic analysis method based on nanopore sequencing data. BACKGROUND

[0002] Clinical samples will inevitably have colonizing bacteria during the collection process, and reagent background bacteria will be introduced during the experiment process. These bacteria may introduce drug-resistant genes, causing false positives and causing difficulties for clinical precision medicine. Based on short read data of second-generation sequencing for drug resistance analysis, it is difficult to distinguish which drug-resistant genes are carried by colonizing bacteria / reagent background bacteria and which are carried by pathogenic bacteria, making it difficult to establish the relationship between pathogenic microorganisms and drug-resistant genes. Although various algorithm models are used to analyze the attribution of pathogenic microorganisms and drug-resistant genes, they are not as accurate as direct comparison of long read data of nanopore sequencing.

[0003] Pathogen metagenomic detection generates a large amount of sequencing data, which contains not only pathogen nucleic acid information but also a large amount of human sequence information, and is inevitably mixed with background (such as environmental, engineering reagent bacteria) contaminated sequences. The existing second-generation pathogen metagenomic sequencing data analysis method is prone to multiple alignment and false alignment, with the risk of false positives, making it difficult to ensure the accuracy of the identification results.

[0004] The construction of a database is an important aspect of ensuring the accuracy of the identification results. At the present stage, a plurality of genomic sequence information databases have been established worldwide, such as GenBank, SRA, ENA, etc. These databases contain a large amount of pathogenic microorganism genomic DNA sequence and annotation information, but part of the data in these databases is derived from user-submitted sequence information, so the database needs to be sorted out, verified, classified, extracted and integrated. In addition, a database of background contamination sequences can also be self-built. The degree of database construction will also limit the accuracy of the identification results. A single analysis process is difficult to achieve the requirement of greatly speeding up the analysis speed on the basis of ensuring the accuracy of the identification results. Multiple analysis processes run in parallel, and different databases are used for gene sequence alignment to ensure the accuracy and timeliness of the identification results.

[0005] The bwa algorithm based on the BWT conversion algorithm mainly analyzes short read data with a read length < 300 bp, and cannot analyze long read (> 2 kb) data of nanopore sequencing. The BLAST algorithm is slow when aligning a large amount of data. Therefore, exploring the analysis method of long read data of nanopore sequencing is a crucial link in making good use of long read data at the present time.

[0006] In view of the above defects, the present inventors have finally obtained the present application after a long period of research and practice. SUMMARY

[0007] The present application aims to solve the problem that the existing sequence alignment software used for nanopore sequencing can reduce the timeliness, accuracy and sensitivity of analysis, and provides a pathogen metagenome analysis method based on nanopore sequencing data.

[0008] In order to achieve the above-mentioned purpose, the present application discloses a pathogen metagenome analysis method based on nanopore sequencing data, comprising the following steps:

[0009] S1, reading off-line data, using the internationally accepted sequencing data standard FASTQ format, supporting gzip method for data compression to reduce storage occupation;

[0010] S2, according to the data type, using Nanofilt software and fastQC software and corresponding parameters to filter low-quality data;

[0011] S3, using samtools software to extract high-quality sequencing data passing data quality control for subsequent analysis;

[0012] S4, first round of rapid alignment, using four kinds of rapid alignment software and algorithms of Kraken2, Centrifuge, Pandora and Minimap2 to preliminarily identify species;

[0013] S5, based on the results of step S3 and step S4, using bcftools software to extract suspected pathogen reads and split according to the preliminarily identified species;

[0014] S6, BLAST aligning the suspected pathogen reads split in step S5 with each species reference sequence, extracting each read optimal alignment, and screening according to the following conditions: alignment length greater than 80% of read length, sequence similarity higher than 90%, and alignment expectation value less than 1E-5;

[0015] S7, performing species statistical analysis at each classification level according to the pathogen type and genome type for the pathogen reads verified in step S6, to determine the pathogen type, genome type, genus level and species level read quantity and constituent ratio in the sample;

[0016] S8, rapidly aligning the data generated in step S3 with the drug resistance gene and mutation database, and preliminarily identifying species by using Minimap2 rapid alignment algorithm;

[0017] S9, based on the results of step S3 and S8, using samtools software to extract suspected pathogen drug resistance gene reads and generating fasta sequence according to the preliminary identification result, and then using Bcftools to obtain consistent sequence;

[0018] S10, the suspected pathogen resistance gene reads obtained in step S9 are subjected to BLAST alignment with each species reference sequence, the optimal alignment of each read is extracted, and the alignment length is greater than 80% of the read length, the sequence similarity is higher than 90%, and the alignment expectation value is less than 1E-5, for screening;

[0019] S11, the pathogen resistance gene reads verified in step S10 are associated with the read taxid number generated in step S7, and the abundance statistical analysis is performed according to the species, so as to determine the types of drug-resistant bacteria, drug-resistant gene and point mutation in the sample, the number of reads and the constituent ratio;

[0020] S12, according to the pathogen statistical information obtained in steps S7 and S11, and the drug-resistant bacteria and drug-resistant gene information, the detection report is generated by supplementing the pathogen annotation information, sample information and patient information.

[0021] In the step S1, the cat or zcat command is used to combine multiple FASTQ data files or compressed FASTQ.gz files from the same sample.

[0022] In the step S2, the Nanofilt software is used for quality control detection and filtering of the FASTQ sequencing data file of the third-generation sequencing, and the sequencing reads with average sequencing data quality greater than Q7 are screened.

[0023] In the step S2, the fastQC software is used for the second-generation data, and the sequencing reads with average sequencing data quality greater than Q30 are screened.

[0024] In the step S4, Kraken2 is an algorithm for identifying bacteria, fungi, viruses and parasites based on species-specific k-mer sequences, Centrifuge is a fast metagenomic sequence classification method for identifying bacteria, fungi, viruses and parasites, Pandora is a fast sequence identification method, Minimap2 is a supplementary alignment for pathogenic organisms with large genomic variations such as viruses, based on minimizer hash table search, chaining algorithm and dynamic programming algorithm.

[0025] In the step S7, the pathogen types include bacteria, fungi, viruses and parasites, and the genome types include DNA and RNA.

[0026] Compared with existing technologies, the advantages of this invention are as follows: For long-read data generated by nanopore sequencing, this invention integrates multiple databases and constructs multiple parallel analysis algorithms for pathogen identification based on sequence alignment, ensuring the accuracy and timeliness of pathogen identification and drug resistance gene analysis. Firstly, in the first round of alignment, four parallel fast alignment algorithms are used to achieve preliminary species identification for typical data (1Gb nanopore sequencing data). In the second and third rounds of verification analysis, the precise alignment algorithm BLAST is used to cross-align the reads of each species identified in the previous round with the reference sequence of that species and other pathogen reference sequences identified in this sample, effectively reducing the computational load of precise alignment. Alignment with drug resistance gene / point mutation databases is performed using Minimap2, and the precise alignment algorithm BLAST is used to identify drug resistance genes / point mutations. By constructing this method and system, accurate analysis and report generation of typical nanopore sequencing pathogen metagenomic data are achieved within 30 minutes, and it is effectively compatible with next-generation sequencing data analysis. Attached Figure Description

[0027] Figure 1 This is a flowchart of the present invention;

[0028] Figure 2 The computing environment and time consumption statistics for the tests of this invention;

[0029] Figure 3 This is a comparison of the time consumption of this invention with common second-generation analysis software;

[0030] Figure 4 This is the quality control chart for Example 1;

[0031] Figure 5 This is the detection quality control chart for Example 2. Detailed Implementation

[0032] The above-mentioned and other technical features and resource consumption of the present invention will be described in more detail below with reference to the accompanying drawings.

[0033] Data detection process:

[0034] 1. Nucleic acid extraction

[0035] Transfer the sample to a 50 mL sterile centrifuge tube and centrifuge at 3500 rpm for 10 min in a benchtop low-speed centrifuge. After centrifugation, remove the supernatant and mix the precipitate obtained by centrifugation with a sterile 1 mL plastic dropper. Transfer the precipitate to a new 1.5 mL centrifuge tube and extract nucleic acid using the kit.

[0036] 2. Determine the nucleic acid concentration using Qubit 4.0. The specific procedure is as follows:

[0037] 1) Prepare working solution with buffer and nucleic acid dye in the ratio of 199:1;

[0038] 2) Prepare standard solution 1 and 2 with working solution and standard in the ratio of 190:10;

[0039] 3) Prepare sample solution with working solution and sample in the ratio of 199:1;

[0040] 4) Draw standard curve with standard solution 1 and 2;

[0041] 5) Determine the concentration of sample in sample solution under the standard curve of step 4).

[0042] 3. Construction of metagenomic sequencing library, specifically:

[0043] 1) Library preparation, in a PCR tube, add 3ul template DNA, 1ul Fragmentation mixture, 20ul nuclease-free water, 1ul RLB(01-12A) and 25ul 2x master mix, put the PCR tube into the PCR instrument, the program is as follows: 98℃ 3min, (98℃ 15s---54℃ 20s---75℃ 3min) 25 cycle---75℃ 7min;

[0044] 2) Magnetic bead purification: transfer the PCR product into a new 1.5ml centrifuge tube, add 40ul magnetic beads and mix well, adsorb at room temperature for 5min, after instant centrifugation, place the centrifuge tube on the magnetic stand, stand at room temperature for 10min, after the solution is clear, remove the supernatant with a pipette, rinse with 75% ethanol and remove the supernatant, air dry at room temperature for 10min, add DNAse / RNAse-free H2O to dissolve, take 1ul for qubit quantification;

[0045] 3) Sequencing: mix barcode with sample, after centrifugation, take 200ul mixed sample, slowly add to the chip, cover the cover, select the program parameters, set the file name, sample name and sample type, select real-time base calling for real-time sequencing.

[0046] 4. Analysis by nanopore pathogen metagenome analysis system

[0047] As shown in Figure 1 and Figure 2 , the pathogen metagenome analysis method based on nanopore sequencing data and the time consumption under fixed computing resources proposed by the application are as follows:

[0048] (1) Read the off-machine data, use the international standard FASTQ format for sequencing data, support gzip method for data compression to reduce storage occupation; use cat or zcat command to merge multiple FASTQ data files or compressed FASTQ.gz files from the same sample; adopt the way of data batch uploading and analysis, this step takes about 2 minutes;

[0049] (2) According to the data type, use different data quality control software and parameters for low-quality data filtering:

[0050] ① For the FASTQ sequencing data file of the third-generation sequencing, use Nanofilt software for quality control detection and filtering, and select the sequencing reads with average sequencing data quality greater than Q7;

[0051] ② For the second-generation data, use fastQC software to select sequencing reads with average sequencing data quality greater than Q30;

[0052] (3) Extract high-quality sequencing data passing data quality control for subsequent analysis;

[0053] (4) First round of rapid alignment, use the following four rapid alignment algorithms for preliminary species identification:

[0054] ① Kraken2, based on species-specific k-mer sequence for bacteria, fungi, viruses and parasites identification;

[0055] ② Centrifuge, a fast metagenomic sequence classification method, for bacteria, fungi, viruses and parasites identification;

[0056] ③ Pandora, a fast and sequence identification method;

[0057] ④ Minimap2, supplementary alignment for pathogens with large genomic variation such as viruses, based on minimizer (a seed of the smallest hash value within a sequence) hash table search, chaining algorithm and dynamic programming algorithm;

[0058] ⑤ Run the first round of four rapid alignment in parallel, this step takes about 10 minutes;

[0059] (5) Based on the results of steps (3), (4), use bcftools to extract suspected pathogenic reads and split them according to the preliminary identified species, this extraction target sequence step takes about 5 minutes;

[0060] (6) The suspected pathogenic reads split in step (5) are subjected to BLAST alignment with each species reference sequence, the optimal alignment of each read is extracted, and the filtering is performed according to the alignment length greater than 80% of the read length, the sequence similarity higher than 90%, and the alignment expectation value (Evalue) less than 1E-5;

[0061] (7) The pathogenic reads verified in step (6) are subjected to statistical analysis according to the species, so as to determine the pathogenic type (bacteria, fungi, viruses, and parasites), the genomic type (DNA and RNA), the read quantity and constituent ratio at the genus level and the species level in the sample; the BLAST alignment with the target sequence in step (6) and step (7) takes about 10 minutes;

[0062] (8) The data generated in step (3) are subjected to rapid alignment with the drug resistance gene and mutation database, and the species are preliminarily identified by using the Minimap2 rapid alignment algorithm; this step and the species alignment in step (4) are run in parallel, and the time consumption is about 10 minutes;

[0063] (9) Based on the results of steps (3) and (8), the read of the suspected pathogen drug resistance gene is extracted, and the fasta sequence is generated according to the preliminary identification result, and the consistency sequence is obtained by using bcftools; this step and step (5) are run in parallel, and the time consumption is 5 minutes;

[0064] (10) The suspected pathogen drug resistance gene reads obtained in step (9) are subjected to BLAST alignment with each species reference sequence, the optimal alignment of each read is extracted, and the filtering is performed according to the alignment length greater than 80% of the read length, the sequence similarity higher than 90%, and the alignment expectation value (Evalue) less than 1E-5;

[0065] (11) The pathogenic drug resistance gene reads verified in step (10) are associated with the read taxid number generated in step S7, and subjected to statistical analysis according to the species, so as to determine the type of the drug-resistant bacteria, the drug resistance gene, and the point mutation, the read quantity and constituent ratio in the sample; steps (9) and (10) are run in parallel with steps (6) and (7), and the time consumption is about 10 minutes

[0066] (12) According to the statistical information of the pathogens in the sample obtained in steps (7) and (11), and the drug-resistant bacteria and drug resistance gene information, the detection report is generated by supplementing and adding the pathogen annotation information, the sample information, and the patient information, etc.; the statistical and automatic report generation process takes about 3 minutes.

[0067] Figure 3The time-consuming comparison of using the application and the commonly used software MetaPhlan3 to analyze the same batch of data is shown. In a Dell T7920 workstation, we configured 512G memory, Intel 5220R CPU 2 blocks, 2T SSD storage 10 blocks, Nvidia RTX A2000(12G) graphics card one block, under this computing resource, the application takes 25-30 minutes to analyze 1G of nanopore data, which is significantly better than the 45-50 minutes of MetaPhlan3 software.

[0068] Embodiment 1

[0069] On July 23, 2021, the lung lavage fluid of a patient with pulmonary infection collected from a hospital was detected by using this pathogen metagenomic analysis method based on nanopore sequencing data, and the detection results are as follows:

[0070] 1. List of detected viruses

[0071] (1) List of DNA viruses

[0072]

[0073] (2) List of RNA viruses

[0074]

[0075] (3) List of detected Mycobacterium tuberculosis / non-tuberculosis

[0076]

[0077] The quality control chart of the detection is shown in Figure 4 , it can be seen from Figure 4 that the sequencing data quality of this detection is qualified and the results are reliable.

[0078] Embodiment 2

[0079] On September 1, 2021, the lung lavage fluid of a patient with pulmonary infection collected from a hospital was detected by using this pathogen metagenomic analysis method based on nanopore sequencing data, and the detection results are as follows:

[0080] 1. List of detected bacteria

[0081]

[0082]

[0083] 2. List of detected fungi

[0084]

[0085] 3. List of detected viruses

[0086] (1) List of DNA viruses

[0087]

[0088] (2) List of RNA viruses

[0089]

[0090] 4. List of detected tuberculous / non-tuberculous mycobacteria

[0091]

[0092] 5. List of detected atypical pathogens

[0093]

[0094] 6. List of detected parasites

[0095]

[0096]

[0097] 7. List of oral and respiratory flora in humans

[0098]

[0099] 8. List of drug resistance genes detected

[0100]

[0101]

[0102] The quality control chart is as follows Figure 5 As shown, by Figure 5 It can be concluded that the sequencing data from this test is of acceptable quality and the results are reliable.

[0103] The above description is merely a preferred embodiment of the present invention and is illustrative rather than restrictive. Those skilled in the art will understand that many changes, modifications, and even equivalents can be made within the spirit and scope defined by the claims of the present invention, all of which will fall within the protection scope of the present invention.

Claims

1. A method for pathogen metagenomic analysis based on nanopore sequencing data, characterized in that, Includes the following steps: S1 reads the data after sequencing, using the internationally recognized FASTQ sequencing data standard, and supports gzip compression to reduce storage space. S2, based on the data type, use Nanofilt software and fastQC software with corresponding parameters to filter low-quality data; S3: Use samtools software to extract high-quality sequencing data that has passed data quality control for subsequent analysis; S4, the first round of rapid comparison, uses four rapid comparison software and algorithms in parallel—Kraken2, Centrifuge, Pandora, and Minimap2—to make preliminary species identification; S5. Based on the results of steps S3 and S4, use samtools software to extract the readings of suspected pathogens and split them according to the preliminarily identified species. S6. The suspected pathogen reads split in step S5 are compared with the reference sequences of each species using BLAST. The best alignment of each read is extracted and then filtered according to the following criteria: alignment length greater than 80% of read length, sequence similarity greater than 90%, and alignment expectation value less than 1E-5. S7. Perform species statistical analysis at each taxonomic level on the pathogen reads verified in step S6 according to pathogen type and genome type to clarify the number and composition ratio of pathogen types, genome types, genus level and species level reads in the sample; S8. The data generated in step S3 is quickly compared with the drug resistance gene and mutation database, and the species is initially identified by using the Minimap2 fast comparison algorithm. S9. Based on the results of steps S3 and S8, use samtools software to extract reads of suspected pathogen drug resistance genes and generate FASTA sequences according to the preliminary identification results. Then use Bcftools to obtain the consistent sequence. S10, the suspected pathogen drug resistance gene reads obtained in step S9 are compared with the reference sequences of each species using BLAST, the best alignment of each read is extracted, and the reads are screened according to the following criteria: alignment length is greater than 80% of the read length, sequence similarity is greater than 90%, and alignment expectation value is less than 1E-5. S11, the pathogen drug resistance gene reads verified in step S10 are associated with the read taxid numbers generated in step S7, and abundance statistical analysis is performed by species to identify the types of drug-resistant bacteria, drug-resistant genes and their point mutations, as well as the number and composition ratio of reads in the sample; S12: Based on the pathogen statistical information, drug-resistant bacteria and drug-resistant gene information obtained in steps S7 and S11, supplement the pathogen annotation information, sample information and patient information, and generate a test report.

2. The pathogen metagenomic analysis method based on nanopore sequencing data as described in claim 1, characterized in that, In step S1, the cat or zcat command is used to merge multiple FASTQ data files or compressed FASTQ.gz files from the same sample.

3. The pathogen metagenomic analysis method based on nanopore sequencing data as described in claim 1, characterized in that, In step S2, quality control and filtering are performed on the FASTQ sequencing data files of the third-generation sequencing. The Nanofilt software is used to screen sequencing reads with an average sequencing data quality greater than Q7.

4. The pathogen metagenomic analysis method based on nanopore sequencing data as described in claim 1, characterized in that, In step S2, fastQC software is used to screen sequencing reads with an average sequencing data quality greater than Q30 for second-generation data.

5. The pathogen metagenomic analysis method based on nanopore sequencing data as described in claim 1, characterized in that, In step S4, Kraken2 is an algorithm for identifying bacteria, fungi, viruses, and parasites based on species-specific k-mer sequences; Centrifuge is a fast metagenomic sequence classification method for identifying bacteria, fungi, viruses, and parasites; Pandora is a fast alignment sequence identification method; and Minimap2 is a supplementary alignment method for pathogens with large viral genome variations, based on minimizer hash table search, chaining algorithm, and dynamic programming algorithm.

6. The pathogen metagenomic analysis method based on nanopore sequencing data as described in claim 1, characterized in that, In step S7, the pathogen types include bacteria, fungi, viruses, and parasites, and the genome types include DNA and RNA.

Images