Host exogenous nucleic acid identification method based on high-throughput sequencing technology

Abstract

The invention discloses a host exogenous nucleic acid identification method based on a high-throughput sequencing technology. The method comprises the following steps: acquiring sequencing data of a to-be-detected sample; filtering the sequencing data twice by using a host reference genome database to remove a host sequence; species classification is performed on the sequencing data after host filtering, and interference of homologous species sequences is reduced; performing sequence comparison by utilizing an identified species reference genome to obtain an accurate species sequence; low-complexity sequence removal is carried out on the compared sequences, so that the influence of sequence preference is reduced; and carrying out secondary comparison on the species reference genome by using BLAST or comparing the species reference genome to an NT library. By applying the method, the exogenous nucleic acid sequences of different host samples can be accurately and specifically identified, and effective analysis of biomass samples of high-host background and low-exogenous species is realized.

Description

Technical Field

[0001] The present invention relates to the field of bioinformatics, and specifically to a method for identifying host exogenous nucleic acids based on high-throughput sequencing technology. Background Art

[0002] With the development of high-throughput sequencing technology and the reduction of costs, it has become possible to screen clinical diseases, identify sequences of different species, and explore the source of host exogenous nucleic acid sequences through sequencing. Sequencing different host samples, such as lung cell lavage fluid, tissue samples, urine, etc., to identify exogenous nucleic acid sequences, can effectively assist in the diagnosis of clinical infections, including infections of bacteria, viruses, fungi and parasites (echinococcosis and liver fluke disease, etc.), and can detect more pathogens than traditional PCR methods. However, due to the large background value of host sequences in these samples from the host, for example, in plasma samples, human sequences account for more than 98% of the sequences, which will be affected by host sequences in the diagnosis of blood infection diseases, significantly reducing the accuracy of detection. In addition to microbial sequences, plants and animals, which belong to eukaryotes, are more susceptible to homologous sequences, low-complexity sequences, and sequences of other closely related species. The current analysis method still has residues when removing the host's own sequences and homologous sequences. In the case of high host background and low abundance of exogenous species sequences, it is easy to produce inaccurate results. In addition, current analysis methods, when dealing with mixed sequences of multiple species (mixed samples of humans and parasites, humans and microorganisms, and humans and plants), especially samples with animal and plant sequences, find it difficult to distinguish the accurate animal and plant sequences of each species. Therefore, the bioinformatics analysis method for identifying the source of exogenous nucleic acids in the host based on high-throughput sequencing technology needs to be further improved.

[0003] At present, when using high-throughput sequencing technology to detect the exogenous nucleic acid sequences of the host, the main research is on the exogenous nucleic acid sequences of human body infection bacteria and viruses that are stored in different human samples (blood, sputum and body fluids, etc.). Such as methods and devices for detecting microorganisms in blood (CN105525033A), methods and devices for microbial analysis of host samples (CN111009286A), methods for determining the presence of predetermined species in environmental samples by metagenomic sequencing (CN114574565A), pathogenic microorganism analysis and identification systems and applications (CN111462821A), a metagenomic data analysis method and system (CN108334750B), and pathogen analysis methods, analysis devices, equipment and storage media based on metagenomics (CN111951895A). Their methods mainly perform simple quality control on the offline data of high-throughput sequencing technology, and after filtering the host sequence in a single step, directly compare or blast to the constructed database of pathogenic bacteria or viruses, and count the nucleic acid sequence information of related bacteria or virus species. Although this method can detect as many sequences of infected species as possible, it may cause false positives or fail to detect when the content of exogenous nucleic acid sequences is low and the background noise of host nucleic acid sequences is large in the early stage of infection. In addition, host humans and microorganisms are far apart in biological classification and are relatively easy to distinguish, but when there is a parasitic infection in the host human, it is more difficult to distinguish the nucleic acid sequences of the two. For example, during the inventor's research and development process, it was found that the blood samples of patients suspected of echinococcosis infection had a low content of nucleic acid sequences released in the early stage of parasitic infection, so it was difficult to capture specific signals using current conventional analysis methods. In addition, few people are currently focusing on studying whether plant nucleic acid sequences may exist in different hosts. Studies have shown that plant RNA sequences exist in mammals. Therefore, in the face of nucleic acids with stronger homology between animal and animal sequences, and animal and plant sequences, accurate nucleic acid sequence identification methods are more needed to obtain accurate and reliable results. Summary of the invention

[0004] The present invention aims to solve one of the technical problems in the related art to a certain extent. To this end, the purpose of the present invention is to propose a bioinformatics analysis method based on high-throughput sequencing technology, which can quickly and accurately identify the species of microorganisms, animals, and plants in samples with high host backgrounds, and improve the credibility of detecting non-host species sequences. The present invention effectively improves the credibility and accuracy of classified species by combining multiple sequence comparison strategy algorithms, and can assist in the extensive identification of different types of host samples, such as pathogenic microorganisms and other exogenous species sequences.

[0005] The purpose of the present invention is achieved by at least one of the following technical solutions.

[0006] A method for identifying host exogenous nucleic acid sequences based on high-throughput sequencing technology comprises the following steps:

[0007] S1. Obtain whole genome sequencing data, targeted sequencing data or metagenomic sequencing data of host samples, and perform routine quality control on the sequencing data;

[0008] S2. Remove the host sequence in the sequencing data: obtain the latest host reference genome, filter out the nucleic acid sequences in the sequencing data that are identical to the host reference genome nucleic acid sequences, and complete the first filtering;

[0009] S3, filtering the nucleic acid sequences identical to the host reference genome nucleic acid sequences in the sequencing data after the first filtration, completing the second filtration;

[0010] S4. Specify the species to be analyzed, including parasites, microorganisms and other species that do not come from the host, obtain the reference genome sequence of the specified analysis species, and build a species classification database; use the metagenomic analysis software to compare the sequencing data after the second filtration to the species classification database to obtain the accurate nucleic acid sequence of the specified analysis species;

[0011] S5. Re-align the accurate taxonomic nucleic acid sequence of the designated analysis species to the reference genome sequence of the designated analysis species, retain the double-end sequences that are successfully aligned to the reference genome of the designated analysis species, filter out the results that cannot be aligned and the results of single-end alignment, and obtain the accurate species nucleic acid sequence;

[0012] S6, removing low-complexity nucleic acid sequences from accurate species nucleic acid sequences;

[0013] S7. Using the BLAST alignment algorithm, align the accurate species nucleic acid sequence after removing the low-complexity nucleic acid sequence to the reference genome sequence of the designated analyzed species, retain the double-end sequence alignment results with a confidence level greater than 90 and a coverage level greater than 90% in the BLAST alignment, and obtain a nucleic acid sequence with high credibility;

[0014] S8. Align the highly reliable nucleic acid sequences to the NT library of NCBI, retain the double-end sequence alignment results with a reliability greater than 90 and a coverage greater than 90% in the BLAST alignment, further filter the nucleic acid sequences of distantly related species, and finally obtain the nucleic acid sequences that are truly from the specified analysis species;

[0015] S9. If the actual nucleic acid sequence from the designated analysis species obtained in step S8 is the same as the accurate sequence of the designated analysis species, it indicates that the host sample contains the nucleic acid sequence of the designated analysis species; otherwise, it indicates that the host sample does not contain the nucleic acid sequence of the designated analysis species.

[0016] Furthermore, the host sample includes a body fluid sample or a tissue sample.

[0017] This method can accurately identify the presence of exogenous nucleic acid sequences in the host, such as some parasites, bacteria, viruses and potential non-host nucleic acid sequences in the human body. This method can effectively distinguish mixed nucleic acid sequences of combinations of animals and plants, animals and microorganisms, plants and microorganisms, etc., and subdivide them into designated species. Here, the host of this method generally refers to the species with the largest proportion of sequences in the mixed sequence, and the rest can be called exogenous nucleic acid sequences. The specific delineation can be based on the needs of actual sample analysis.

[0018] Furthermore, the nucleic acid sequence includes a DNA sequence and an RNA sequence.

[0019] Furthermore, in step S1, the conventional quality control refers to removing sequencing residual adapter sequences and low-quality sequencing reads; the low-quality sequencing reads include sequencing reads with an N ratio greater than 5% and sequencing reads with a length less than 30.

[0020] Furthermore, during the second filtering, a different alignment algorithm than that used during the first filtering needs to be used for nucleic acid sequence filtering.

[0021] Furthermore, in step S4, if there is a species that may have homology with the sequence of the designated analysis species during the analysis process and affect the analysis results, the reference genome sequence of the species affecting the analysis is obtained, and the genomes of the designated analysis species and the species affecting the analysis are mixed to construct a species classification database.

[0022] Furthermore, in step S6, the low-complexity nucleic acid sequence refers to a nucleic acid sequence in which the sequencing reads have an interval in which the length of two or more consecutive identical bases exceeds 5%.

[0023] Furthermore, in step S6, different software is needed from that used in step S5 when aligning to the species reference genome to remove low-complexity nucleic acid sequences from the accurate species nucleic acid sequence; thereby once again verifying the accuracy of the sequence and effectively reducing the systematic errors and biases caused by different software.

[0024] Further, in step S8, the distant species refers to species that are different at the family level when the double-ended sequence is annotated using the NT library for species annotation. These species may have similar sequences in evolution, but their taxa are far apart, and such sequences should not be used as the true sequence of the specified species, such as the result of cross-border between eukaryotes and bacteria, eukaryotes and viruses. After filtering by distant species, the specific sequence of the species can be obtained, which greatly improves the credibility of the sequence information.

[0025] Furthermore, if the reference genome sequence of the specified analysis species is not included in the NT library of NCBI, no valid results will be produced. In this case, the comparison with the NT library can be ignored, that is, step S8 is not executed, and the highly credible sequence obtained in step S7 is directly used as the accurate sequence of the specified analysis species, and then step S9 is executed. At this time, the interference of sequences of distant species is not removed.

[0026] Compared with the prior art, the advantages of the present invention are:

[0027] Through multiple comparison and filtering, the present invention can almost completely remove host sequences from sequencing data containing a high proportion (99%) of host sequences, capture weak exogenous nucleic acid sequence signals, accurately identify nucleic acid sequences from exogenous sources, effectively improve the accuracy of analysis results, and reduce false positives. BRIEF DESCRIPTION OF THE DRAWINGS

[0028] Figure 1 This is a schematic diagram of an information analysis process of a method for identifying host exogenous nucleic acid sequences based on high-throughput sequencing technology in Example 1 of the present invention;

[0029] Figure 2 This is a schematic diagram of an information analysis process of a method for identifying host exogenous nucleic acid sequences based on high-throughput sequencing technology in Example 2 of the present invention;

[0030] Figure 3 This is a schematic diagram of the genome coverage of the exogenous sequences identified according to the results of Example 3 of the present invention. DETAILED DESCRIPTION

[0031] Embodiments of the present invention will be described in detail below with reference to the accompanying drawings, which are intended to explain the present invention but should not be construed as limiting the present invention.

[0032] Example 1 - Identification of Echinococcus Sequences in Patients Infected with Echinococcosis:

[0033] Human echinococcosis is a parasitic disease caused by tapeworms of the genus Echinococcus, also known as hydatid disease. When humans are infected with Echinococcus, hydatid cysts often form in the liver, and the incubation period of the disease is as long as several years. Currently, the diagnosis of echinococcosis mainly relies on imaging diagnosis, but it is not easy to detect in the early latent stage, which greatly increases the cost of treatment. The free nucleic acid sequence of Echinococcus in the blood can effectively detect early echinococcosis, but due to the excessive background noise of its host sequence, the sequence from Echinococcus is not easy to distinguish, which can easily lead to false positive results. A more accurate analysis method is needed to distinguish between parasite sequences and host sequences.

[0034] In order to evaluate the effectiveness of the inventive method in distinguishing parasite sequences from host sequences, in this example, 50 simulated data samples were generated using Wgsim (https: / / github.com / lh3 / wgsim) with default parameters using humans as hosts. In addition to 10 million paired-end sequenced host sequences, 1000 pairs of paired-end sequenced sequences of the two main infectious species, Echinococcus granulosus (GCA_021556725.1) and Echinococcus multilocularis (GCA_000469725.3), were also included.

[0035] A method for identifying host exogenous nucleic acid sequences based on high-throughput sequencing technology, such as Figure 1 As shown, the following steps are included:

[0036] S1. Obtain whole genome sequencing data, targeted sequencing data or metagenomic sequencing data of host samples, and perform routine quality control on the sequencing data;

[0037] The host sample includes a body fluid sample or a tissue sample.

[0038] In this embodiment, the comparison databases of BOWTIE2 and Kraken2 were first established based on the hg38 reference genome downloaded from NCBI. Subsequently, the reference genomes of Echinococcus granulosus and Echinococcus multilocularis were downloaded from NCBI to establish the comparison databases of Krakenuniq, BOWTIE2 and BLAST. The above databases will be used for subsequent bioinformatics analysis.

[0039] The conventional quality control refers to the removal of sequencing residual adapter sequences and low-quality sequencing reads; the low-quality sequencing reads include sequencing reads with an N ratio greater than 5% and sequencing reads with a length less than 30.

[0040] S2. Removing the host sequence in the sequencing data: Obtain the latest host reference genome. In this embodiment, the nucleic acid sequence in the sequencing data that is identical to the host reference genome nucleic acid sequence is filtered using a local alignment algorithm software to complete the first filtering;

[0041] The nucleic acid sequence includes a DNA sequence and an RNA sequence.

[0042] The preferred software is BOWTIE2, and the parameter is set to "—very-senstive-local" to use a local alignment method to remove sequences from the host as completely as possible.

[0043] S3. In this embodiment, the k-mer alignment algorithm software is used to filter out the nucleic acid sequences identical to the host reference genome nucleic acid sequences in the sequencing data after the first filtration, and the second filtration is completed;

[0044] The preferred software is Kraken2. Use the Kraken2-build function to build a comparison database based on the host reference genome. Set the parameter confidence to 0, which will classify as many sequences as possible to the host and remove host sequences as much as possible.

[0045] In this embodiment, the initial analysis data is aligned to the hg38 reference genome using BOWTIE2, the sequences on the aligned genome are filtered out, and the paired-end sequences that cannot be aligned are retained. Preferably, BOWTIE2 uses "-very-sensitive-local" as a parameter. The retained data is then aligned to the hg38 reference genome using Kraken2, and the paired-end sequence results that cannot be aligned are retained. Preferably, Kraken2 uses the parameter "-paired-confidence=0".

[0046] S4. Specify the species to be analyzed, including parasites, microorganisms and other species that do not come from the host, obtain the reference genome sequence of the specified analysis species, and build a species classification database; use the metagenomic analysis software to compare the sequencing data after the second filtration to the species classification database to obtain the accurate nucleic acid sequence of the specified analysis species;

[0047] If there are species that may have homology with the sequences of the designated analysis species and affect the analysis results during the analysis process, the reference genome sequence of the species affecting the analysis is obtained, and the genomes of the designated analysis species and the species affecting the analysis are mixed to construct a species classification database.

[0048] The preferred software is MetaPhlAn4 or Krakenuniq. In the present invention, Krakenuniq is used as an example to establish a species classification database based on species sequences and annotate and classify the sequences by species. This step can utilize unique k-mer to better classify most homologous sequences and reduce interference between model species.

[0049] In this example, the data retained after filtering the host sequences were aligned to the constructed species database using Krakenuniq using the parameter "-paired", and the sequences annotated to Echinococcus granulosus and Echinococcus multilocularis were extracted for subsequent analysis. This step effectively separated the sequences of the two homologous Echinococcus species.

[0050] S5. In this embodiment, the end-to-end alignment algorithm is used to re-align the accurate classification nucleic acid sequence of the specified analysis species to the reference genome sequence of the specified analysis species, retain the double-end sequence successfully aligned to the reference genome of the specified analysis species, filter out the results that cannot be aligned and the single-end alignment, and obtain the accurate species nucleic acid sequence;

[0051] In this example, the obtained sequences of Echinococcus granulosus and Echinococcus multilocularis were aligned to their respective genomes using BOWTIE2. The results of the double-end sequence alignment were then extracted using SAMTOOLS. Preferably, BOWTIE2 uses an "end-to-end" alignment mode.

[0052] S6. In this embodiment, the DUST algorithm is used to remove low-complexity nucleic acid sequences from the accurate species nucleic acid sequences;

[0053] The low-complexity nucleic acid sequence refers to a nucleic acid sequence in which the sequencing reads have more than two consecutive identical bases with a base length exceeding 5%, such as a protein sequence "PPCDPPPPPKDKKKKDDGPP" and a nucleic acid sequence in the form of "AATAAAAAAAATAAAAAAT". Since the sequence has a strong preference and a single base, it will correspond to multiple species when performing a sequence similarity search, and it is impossible to effectively determine the species. Therefore, it is necessary to remove these sequences with preferences to improve the accuracy of identification.

[0054] The preferred software is PRINSEQ. The parameters of the software are as follows: the low complexity sequence masking algorithm is DUST, and the range of calculated bases is 7.

[0055] In this example, the sequences aligned to Echinococcus are used to remove low-complexity sequences using PRINSEQ. Preferably, PRINSEQ uses the parameter "-lc_method dust-lc_threshold 7".

[0056] S7. Using the BLAST alignment algorithm, align the accurate species nucleic acid sequence after removing the low-complexity nucleic acid sequence to the reference genome sequence of the designated analyzed species, retain the double-end sequence alignment results with a confidence level greater than 90 and a coverage level greater than 90% in the BLAST alignment, and obtain a nucleic acid sequence with high credibility;

[0057] S8. Align the highly reliable nucleic acid sequences to the NT library of NCBI, retain the double-end sequence alignment results with a reliability greater than 90 and a coverage greater than 90% in the BLAST alignment, further filter the nucleic acid sequences of distantly related species, and finally obtain the nucleic acid sequences that are truly from the specified analysis species;

[0058] Use BLAST to compare the species reference genome or NT library, and give priority to the latest database on NCBI;

[0059] If the reference genome sequence of the specified analysis species is not included in the NT library of NCBI, no valid results will be produced. In this case, the comparison with the NT library can be ignored, that is, step S8 is not executed, and the highly reliable sequence obtained in step S7 is directly used as the accurate sequence of the specified analysis species, and then step S9 is executed. At this time, the interference of sequences of distant species is not removed.

[0060] S9. If the actual nucleic acid sequence from the designated analysis species obtained in step S8 is the same as the accurate sequence of the designated analysis species, it indicates that the host sample contains the nucleic acid sequence of the designated analysis species; otherwise, it indicates that the host sample does not contain the nucleic acid sequence of the designated analysis species.

[0061] In this embodiment, BLAST uses the parameters "-evalue 1e-5-perc_identity 90-qcov_hsp_perc90".

[0062] Through the present invention, the information analysis results of Example 1 are shown in Table 1. The results show that the present invention can completely remove the sequences from the host, and the final sequence recovery rates from Echinococcus granulosus and Echinococcus multilocularis are 95.76% and 95.89%, respectively. For plasma samples of patients with early Echinococcus infection, in the case of containing a low biomass of Echinococcus sequences, the present invention can also effectively reduce the interference of the host and retain a trace amount of Echinococcus sequence signals.

[0063] Table 1 Information analysis results of Example 1

[0064] Step Echinococcus_granulosus Echinococcus_multilocularis Homo_sapiens Raw 991 998 10000000 Bowtie2 991±0 998±0 4325852±1887.69 Kraken2 984.35±0.80 995.47±0.89 5013.80±1212.93 Krakenuniq 976.84±0.75 986.27±0.45 0 Bowite2 971.39±4.40 974.39±4.01 0 PRINSEQ 970.65±3.96 973.65±4.46 0 BLAST 949.04±1.76 957.27±0.45 0 BLAST recovery rate (%) 95.76％ 95.89％ 0

[0065] In Example 1, under the mixed sequences of all 9 model species, the present invention can still distinguish the sequences of the 9 species well, and maintain consistency and specificity, with a good sequence recovery rate, which confirms the effectiveness of the analysis method.

[0066] Example 2 - Identification of contaminating sequences from human host sequences mixed with different model species:

[0067] During the experiment and sequencing process, there may be sequence contamination from some model species of animals, plants or microorganisms. If you need to obtain a pure sequence after removing the contamination, this method can be used to effectively separate the specified sequence of the corresponding species and recover the sequence.

[0068] like Figure 2 As shown, a method for identifying host exogenous nucleic acid sequences based on high-throughput sequencing technology comprises the following steps:

[0069] S1. Obtain whole genome sequencing data, targeted sequencing data or metagenomic sequencing data of host samples, and perform routine quality control on the sequencing data;

[0070] In this embodiment, in order to verify the accuracy and consistency of the analysis method, 50 simulated data samples were generated using Wgsim (https: / / github.com / lh3 / wgsim) with default parameters using humans as hosts. In addition to 10 million double-end sequenced host sequences, it also contains mixed sequences of 2 plant species sequences (potatoes and tomatoes), 2 animal species sequences (domestic pigs and chickens), and 2 microbial species sequences (new coronavirus and Escherichia coli). In the simulated data samples, the sequence of each species accounts for only 0.01% of the host sequence, all of which are 1000 double-end sequencing reads, and the specific reference genome information is shown in Table 2 below.

[0071] Table 2 Simulation data construction information

[0072] Species Total of reads (PE) Assembly_accession* Sample number Solanum tuberosum (potato) 1000 GCA_014189475.1 50 Solanum lycopersicum (tomato) 1000 GCF_000188115.5 50 Susscrofa (domestic pig) 1000 GCF_000003025.6 50 Gallusgallu (chicken) 1000 GCF_016699485.2 50 SARS-COV-2 (New Coronavirus) 1000 GCF_009858895.2 50 Escherichia coli 1000 GCF_000005845.2 50 Homosapiens 10000000 GCF_000001405.40 50

[0073] *Assembly accession: Assembly accession number, the sequence accession number assigned to the genome assembly format by NCBI database staff.

[0074] S2. Removing host sequences from sequencing data: Obtain the latest host reference genome, and use local alignment algorithm software to filter out nucleic acid sequences in the sequencing data that are identical to the host reference genome nucleic acid sequences, completing the first filtering;

[0075] In this embodiment, the generated simulated data (simulated data does not require quality control) is aligned to the human reference genome using alignment software, the aligned sequences are removed, and the sequences that cannot be aligned are retained. Preferably, the "-very-sensitive-local" parameter of BOWTIE2 is used for local alignment. In addition, alignment software such as BWA can also be used for sequence alignment, and then SAMTOOLS is used to extract sequences that may not be aligned. The human reference genome here uses the hg38 reference genome downloaded from NCBI, and the human T2T genome can also be selected. After this step of processing, a large number of host sequences can be reduced, reducing the interference caused to subsequent analysis.

[0076] S3, using k-mer alignment algorithm software to filter out the nucleic acid sequences in the sequencing data after the first filtration that are identical to the host reference genome nucleic acid sequences, and completing the second filtration;

[0077] In the present embodiment, the metagenomic analysis software Kraken2 is selected to establish a host reference database based on the human reference genome. Kraken2 uses the k-mer alignment algorithm, which can further filter the host sequence. Then the sequence preliminarily filtered in step 2 is aligned to the host database by Kraken2, and the parameter "-paired confidence 0" is used to filter the sequence classified to the host, retain the sequence that cannot be classified, and output the unclassified sequence as a fastq format file. Through the second filtering, the sequences similar to and homologous to the host can be further filtered out, and the remaining sequences almost do not have sequences belonging to the host. The human reference genome of the host database established in this step is hg38, and the T2T genome of the human can also be selected to construct the database.

[0078] S4. Specify the species to be analyzed, including parasites, microorganisms and other species that do not come from the host, obtain the reference genome sequence of the specified analysis species, and build a species classification database; use the metagenomic analysis software to compare the sequencing data after the second filtration to the species classification database to obtain the accurate nucleic acid sequence of the specified analysis species;

[0079] In the present embodiment, the metagenomic analysis software Krakenuniq is selected, and according to the above-mentioned 6 specified species reference genomes, the species classification database is established based on the above-mentioned species using the Krakenuniq-build function. In the process of establishing the classification database, it is necessary to download the taxdump.tar.gz classification tree information file from NCBI to assist in establishing a complete classification result. Subsequently, the sequence after step three filtering is subjected to detailed species classification, and the sequence of the specified species is extracted by the extract_kraken_reads.py script of Krakentools, and the output is a fastq format result. By establishing a homologous database of multiple species and performing sequence classification of species, the interference of homologous species sequences can be effectively reduced, and the accuracy and consistency of the analysis can be improved. The reference genomes for establishing the database in this step are all from NCBI, and the specific information is referred to Table 2.

[0080] S5. Use an end-to-end alignment algorithm to re-align the accurate taxonomic nucleic acid sequence of the specified analysis species to the reference genome sequence of the specified analysis species, retain the double-end sequences that are successfully aligned to the reference genome of the specified analysis species, filter out the results that cannot be aligned and the results of single-end alignment, and obtain the accurate species nucleic acid sequence;

[0081] In the present embodiment, the sequence after species classification is compared to the reference genome of the specified species using the alignment software BOWTIE2 using default parameters. Preferably, using SAMTOOLS, select the double-ended sequence successfully compared to the species genome, filter the result (flag value includes 4 and 8 results) that cannot be compared and the single-ended comparison. Preferably, the double-ended comparison result compared to the species genome is retained, i.e., the specified species sequence. By comparing back to the genome of the species, the sequence unrelated to the species can be filtered out, and the pollution caused by the host and the remaining species can also be further reduced. The alignment software can also use the software of different precisions such as BWA or SNAP-aligner to improve the reliability of the alignment result. Selecting the reference species genome with higher assembly quality and integrity at the same time can also improve the accuracy of the result.

[0082] S6. Use the DUST algorithm to remove low-complexity nucleic acid sequences from accurate species nucleic acid sequences;

[0083] In this embodiment, low-complexity sequences exist in multiple species, have preferences, and cannot effectively distinguish the source. It is necessary to filter the low-complexity sequences to reduce the impact on subsequent analysis. Preferably, the software PRINSEQ is used, and the parameters "-lc_meathod dust-lc_threshold 7" are set to filter the low-complexity sequences, and the sequences with higher specificity are retained as the results of subsequent analysis.

[0084] S7. Using the BLAST alignment algorithm, align the accurate species nucleic acid sequence after removing the low-complexity nucleic acid sequence to the reference genome sequence of the designated analyzed species, retain the double-end sequence alignment results with a confidence level greater than 90 and a coverage level greater than 90% in the BLAST alignment, and obtain a nucleic acid sequence with high credibility;

[0085] In the present embodiment, in order to further improve the accuracy of the sequence screened out in step 4, a BLAST database of the reference genome of the specified species is constructed. Preferably, the sequence of step 4 is compared twice using BLAST comparison software, and the double-end sequence comparison results included with the credibility of BLAST comparison greater than 90 and the coverage greater than 90% are selected, and the output result is in fasta format. The results that cannot be compared and the single-end comparison are filtered to obtain a clean sequence of the specified species. The accuracy of the sequence is further verified by different comparison algorithms. Selecting multiple comparison algorithms can further reduce the bias brought by the software itself, improve the accuracy and uniqueness of the species sequence, and improve the accuracy for subsequent analysis.

[0086] S8. Align the highly reliable nucleic acid sequences to the NT library of NCBI, retain the double-end sequence alignment results with a reliability greater than 90 and a coverage greater than 90% in the BLAST alignment, further filter the nucleic acid sequences of distantly related species, and finally obtain the nucleic acid sequences that are truly from the specified analysis species;

[0087] S9. If the actual nucleic acid sequence from the designated analysis species obtained in step S8 is the same as the accurate sequence of the designated analysis species, it indicates that the host sample contains the nucleic acid sequence of the designated analysis species; otherwise, it indicates that the host sample does not contain the nucleic acid sequence of the designated analysis species.

[0088] In this embodiment, in order to further exclude the interference of similar sequences of homologous species. The clean sequence of the specified species is preferably aligned to the NT (May 6, 2023) library downloaded from NCBI using BLAST software. Remove the results with a credibility of less than 90 and a coverage of less than 90% in the alignment results. Retain the double-end results that can be aligned to the classification of the specified species. Finally, the results with high credibility and accuracy of the specified species are obtained, which are sufficient to support the subsequent further analysis and downstream application of the sequence of the species. Through this step, the interference of distant species sequences is removed, and the most accurate results are obtained compared with previous methods. Table 3 is the analysis results of each step of the simulated 2 plant sequences. In the results, through the two-step host sequence removal, the host sequence can be removed by 99.999713%, which can be almost completely removed. In the third step of the analysis step, the host sequence is no longer included in the sequence of the specified species. After the final secondary alignment, the sequence recovery rates of BLAST to the specified species were 89.2% and 91.29%, respectively, and the sequence consistency was high. At the same time, even after strict screening conditions of BLAST alignment to NT, the sequence recovery rates reached 71.1% and 91.39% respectively, excluding the interference of most distant species and homologous sequences, and still had a high recovery rate. This confirms that this method has relatively accurate analysis results in mixed animal and plant sequences.

[0089] Table 3 Sequence analysis results of two plants in Example 2

[0090]

[0091]

[0092] Table 4 shows the analysis results of each step of the simulated 2 animal sequences. In step S3, the sequence of the host is no longer included in the sequence of the specified species. After the final secondary alignment, the sequence recovery rates of BLAST to the specified species were 91.4% and 79.18%, respectively, and the sequence consistency was high. At the same time, even after the strict screening conditions of BLAST alignment to NT, the sequence recovery rates reached 57% and 24.12%, respectively, excluding the interference of most distant species and homologous sequences. Although the recovery rate is not very high, there are a large number of indistinguishable homologous and contaminated sequences in the NT library of higher animals. Therefore, through filtering and quality control of this method, relatively accurate results can still be obtained after secondary alignment. It is confirmed that this method can find species-specific sequences in animal and animal mixed sequences, thereby improving the accuracy of the identification results.

[0093] Table 4 Sequence analysis results of two animals in Example 2

[0094] Step Gallus_gallus(PE) Sus_scrofa(PE) Homo sapiens (PE) Raw 1000 999 1000000 Bowtie2 997.4±0.73 974.4±4.19 4320973.26±2174.37 Kraken2 968.62±4.63 852.2±11.59 287.04±12.82 Krakenuniq 968.54±4.61 851.7±11.59 0 Bowtie2 956.54±4.93 837.74±12.72 0 PRINSEQ 950.48±3.91 831.56±12.43 0 BLAST 914.6±8.46 791.56±14.37 0 BLAST recovery rate (%) 91.40％ 79.18％ 0 BLAST_NT 570.24±6.99 241.76±10.76 0 BLAST_NT recovery rate (%) 57.00％ 24.12％ 0

[0095] Table 5 shows the analysis results of each step of the simulated two microbial sequences. Microorganisms and hosts belong to different kingdoms in biological taxonomy and have a far species gap, so a higher recovery rate can be obtained. Despite this, pathogen detection of microorganisms still requires more accurate identification methods to reduce false positive results in clinical diagnosis. After secondary alignment, the sequence recovery rates of BLAST to the specified species were 97% and 97.1%, respectively, and the sequence consistency was high. At the same time, even after strict screening conditions of BLAST alignment to NT, the sequence recovery rates reached 97.3% and 98.7%, respectively, excluding the interference of most distant species and homologous sequences, and there was still a high recovery rate. It is confirmed that this method has relatively accurate analysis results in mixed sequences of animals and microorganisms.

[0096] Table 5 Results of sequence analysis of two microorganisms in Example 2

[0097] Step Escherichia_coli(PE) SRAS-COV-2(PE) Homo sapiens (PE) Raw 1000 1000 1000000 Bowtie2 1000±0 1000±0 4320973.26±2174.37 Kraken2 999.96±0.2 999.88±0.33 287.04±12.82 Krakenuniq 999.44±0.86 998.38±0.53 0 Bowtie2 987.3±3.54 987.98±1.53 0 PRINSEQ 987.3±3.54 987.98±1.53 0 BLAST 970.7±4.62 971.04±5.98 0 BLAST recovery rate (%) 97.00％ 97.10％ 0 BLAST_NT 973.52±4.68 987.84±1.77 0 BLAST_NT recovery rate (%) 97.30％ 98.70％ 0

[0098] Example 3 - Identification of Endobacterium Cardiomycetum infection in patients with endocarditis:

[0099] Identification of pathogens is the basis for the diagnosis of endocarditis (IE) and can guide clinical treatment and the choice of antibiotics. Currently, 31% of IE patients have negative microbial culture results in blood samples or heart valve samples. Metagenomic analysis using high-throughput sequencing technology can effectively identify the pathogens of infection and improve the effectiveness of clinical microbial diagnosis. By downloading sequencing data of a patient with endocarditis from the European Nucleotide Archive (ENA), project number PRJEB25228, including 2 negative samples (NEC1 and NEC2) and 2 positive samples (Sample1 and Sample2), we identified whether there were exogenous nucleic acid sequences of Endocardial Bacillus in the patient.

[0100] Through the analysis of the present invention, as shown in Table 6, the table shows the number of pair endreads remaining after each analysis step S1-S9. In the negative sample, after analysis, all host human sequences were finally removed by 100%, and a clean control result without Endobacter cardiogensis was obtained. In the positive sample, after analysis by the method of the invention, Sample1 and Sample2 obtained 265834 and 2884 pairs of double-end sequencing sequences of Endobacter cardiogensis, respectively, accounting for 19.03% and 0.19% of the original sequence, respectively. The identified sample1 and sample2 were aligned to the Endobacter cardiogensis reference genome (NCBI RefSeqassembly: GCF_900637305.1), and the alignment rates were 99.99% and 100%, respectively. Then, we calculated the coverage and average sequencing depth of the Endobacter cardiogensis identification sequences of the two samples on their genomes, and found that the coverage of sample1 was 85.49% and the average sequencing depth was 49.84; the coverage of sample2 was 30.09% and the average sequencing depth was 0.54. We plotted the distribution of Endobacterium cardiogensis sequences on the genome based on the reads identified by sample1, as shown in Figure 3 The above results confirm that the present invention can effectively identify exogenous nucleic acids in the human body and assist in clinical diagnosis.

[0101] Table 6 Information analysis results of Example 3

[0102]

[0103]

[0104] The above embodiments are only preferred implementation modes of the present invention and are only used to explain the present invention rather than to limit the present invention. Any changes, substitutions, modifications, etc. made by those skilled in the art without departing from the spirit of the present invention should fall within the scope of protection of the present invention.

Claims

1. A method for identifying host exogenous nucleic acid sequences based on high-throughput sequencing technology, characterized in that: The steps include: S1. Obtain whole genome sequencing data, targeted sequencing data or metagenomic sequencing data of host samples, and perform routine quality control on the sequencing data; S2. Remove the host sequence in the sequencing data: obtain the latest host reference genome, filter out the nucleic acid sequences in the sequencing data that are identical to the host reference genome nucleic acid sequences, and complete the first filtering; S3, filtering the nucleic acid sequences identical to the host reference genome nucleic acid sequences in the sequencing data after the first filtration, completing the second filtration; S4. Specify the species to be analyzed, including parasites, microorganisms and other species that do not come from the host, obtain the reference genome sequence of the specified analysis species, and build a species classification database; use the metagenomic analysis software to compare the sequencing data after the second filtration to the species classification database to obtain the accurate nucleic acid sequence of the specified analysis species; S5. Re-align the accurate taxonomic nucleic acid sequence of the designated analysis species to the reference genome sequence of the designated analysis species, retain the double-end sequences that are successfully aligned to the reference genome of the designated analysis species, filter out the results that cannot be aligned and the results of single-end alignment, and obtain the accurate species nucleic acid sequence; S6, removing low-complexity nucleic acid sequences from accurate species nucleic acid sequences; S7. Using the BLAST alignment algorithm, align the accurate species nucleic acid sequence after removing the low-complexity nucleic acid sequence to the reference genome sequence of the designated analyzed species, retain the double-end sequence alignment results with a confidence level greater than 90 and a coverage level greater than 90% in the BLAST alignment, and obtain a nucleic acid sequence with high credibility; S8. Align the highly reliable nucleic acid sequences to the NT library of NCBI, retain the double-end sequence alignment results with a reliability greater than 90 and a coverage greater than 90% in the BLAST alignment, further filter the nucleic acid sequences of distantly related species, and finally obtain the nucleic acid sequences that are truly from the specified analysis species; S9. If the actual nucleic acid sequence from the designated analysis species obtained in step S8 is the same as the accurate sequence of the designated analysis species, it indicates that the host sample contains the nucleic acid sequence of the designated analysis species; otherwise, it indicates that the host sample does not contain the nucleic acid sequence of the designated analysis species.

2. The method for identifying host exogenous nucleic acid based on high-throughput sequencing technology according to claim 1, characterized in that: The host sample includes a body fluid sample or a tissue sample.

3. The method for identifying host exogenous nucleic acid sequences based on high-throughput sequencing technology according to claim 1, characterized in that: The nucleic acid sequence includes a DNA sequence and an RNA sequence.

4. The method for identifying host exogenous nucleic acid sequences based on high-throughput sequencing technology according to claim 1, characterized in that: In step S1, the conventional quality control refers to removing sequencing residual adapter sequences and low-quality sequencing reads; the low-quality sequencing reads include sequencing reads with an N ratio greater than 5% and sequencing reads with a length less than 30.

5. The method for identifying host exogenous nucleic acid sequences based on high-throughput sequencing technology according to claim 1, characterized in that: During the second filtering, a different alignment algorithm than that used during the first filtering needs to be used for nucleic acid sequence filtering.

6. The method for identifying host exogenous nucleic acid sequences based on high-throughput sequencing technology according to claim 1, characterized in that: In step S4, if there is a species that may have homology with the sequence of the designated analysis species and affect the analysis results during the analysis process, the reference genome sequence of the species affecting the analysis is obtained, and the genomes of the designated analysis species and the species affecting the analysis are mixed to construct a species classification database.

7. The method for identifying host exogenous nucleic acid based on high-throughput sequencing technology according to claim 1, characterized in that: In step S6, the low-complexity nucleic acid sequence refers to a nucleic acid sequence in which the sequencing reads have an interval with more than 2 consecutive identical bases exceeding 5% in length.

8. The method for identifying host exogenous nucleic acid based on high-throughput sequencing technology according to claim 1, characterized in that: In step S6, different software is needed from that used in step S5 when aligning to the species reference genome to remove low-complexity nucleic acid sequences from the accurate species nucleic acid sequence.

9. The method for identifying host exogenous nucleic acid based on high-throughput sequencing technology according to claim 1, characterized in that: In step S8, the distantly related species refers to species that are different at the family level when the paired-end sequences are annotated for species annotation using the NT library.

10. The method for identifying host exogenous nucleic acid based on high-throughput sequencing technology according to claim 1, characterized in that: If the reference genome sequence of the specified analysis species is not included in the NT library of NCBI, no valid results will be produced. At this time, the comparison with the NT library is not considered, that is, step S8 is not executed, and the highly reliable sequence obtained in step S7 is directly used as the accurate sequence of the specified analysis species, and then step S9 is executed. At this time, the interference of the sequences of distant species is not removed.

Images