A method for analyzing segmented RNA viruses based on macrotranscriptome sequencing and its application

Abstract

The present invention discloses a method for analyzing segmented RNA viruses based on macrotranscriptome sequencing, wherein the method annotates the macrotranscriptome sequencing data of the sample to be tested with the NR database and the RdRP database, and performs Spearman correlation analysis on each sequence in the macrotranscriptome sequencing data, and screens based on the annotation results, the correlation analysis results, and the sequence information in the macrotranscriptome sequencing data to obtain segmented RNA virus data. The method can discover segmented RNA viruses and obtain complete fragments of segmented RNA virus genomes, and can not only realize accurate detection of existing segmented RNA viruses, but also discover new segmented RNA viruses in biological samples or environmental samples with unknown segmented RNA virus information.

Description

Technical Field

[0001] The present invention relates to the field of biological detection technology, and in particular to a method for analyzing segmented RNA viruses based on macrotranscriptome sequencing and an application thereof. Background Art

[0002] Segmented RNA viruses refer to viruses whose genetic material is not composed of a continuous nucleic acid molecule, but is composed of two or more independent nucleic acid fragments (usually RNA). Each nucleic acid fragment of these viruses can encode one or more viral proteins. During the life cycle of the virus, these fragments must exist in the same host cell at the same time in order to assemble into complete virus particles for replication and transmission.

[0003] Before the emergence of high-throughput sequencing technology, the search for segmented RNA viruses was mainly through experimental expression of corresponding proteins, followed by sequencing and homology analysis for identification and classification. With the popularization and development of high-throughput sequencing technology, researchers have gradually adopted two more efficient methods to search for segmented viruses: (1) similar load method and (2) co-occurrence method. The similar load method utilizes the load similarity of different genomes of segmented RNA viruses. By comparing the loads of different fragments, the similarities between different fragments are determined, and then the segments belonging to the same segmented RNA virus are found (10.1073 / pnas.1324194111); the co-occurrence method is based on the fact that the various fragments of segmented RNA viruses in different samples appear together to discover segmented viruses (10.1101 / 2020.02.10.942854).

[0004] These two methods combine high-throughput sequencing technology and bioinformatics analysis. Although they can detect segmented RNA viruses and their segments, they can only find the RdRP fragments of segmented RNA viruses, and cannot find the complete genome of segmented RNA viruses. At the same time, these two methods are complex to operate, have a long operation cycle, and cannot repeatedly search for segmented RNA viruses in batches.

[0005] Therefore, there is an urgent need for a new method for analyzing segmented RNA viruses, which can shorten the analysis time while achieving accurate detection of segmented RNA viruses and their segments, and find out the complete genome of the segmented RNA virus as much as possible. Summary of the invention

[0006] The purpose of the present invention is to overcome the above-mentioned deficiencies of the prior art and to provide a method for analyzing segmented RNA viruses based on macrotranscriptome sequencing and its application.

[0007] The first object of the present invention is to provide a method for analyzing segmented RNA viruses based on macrotranscriptome sequencing.

[0008] The second object of the present invention is to provide application of the above method in constructing a model for analyzing segmented RNA viruses.

[0009] A third object of the present invention is to provide a model for analyzing segmented RNA viruses.

[0010] A fourth object of the present invention is to provide application of the above method and / or the above model in analyzing segmented RNA viruses.

[0011] In order to achieve the above object, the present invention is implemented by the following scheme:

[0012] A method for analyzing segmented RNA viruses based on macrotranscriptome sequencing comprises the following steps:

[0013] S1. Take no less than 3 different samples as test samples, perform macrotranscriptome sequencing on the test samples, collect the macrotranscriptome sequencing data in FastQ format from the offline data for quality control, remove the ribosomal RNA data from the macrotranscriptome sequencing data after quality control, and obtain clean fq data;

[0014] S2. Splicing the clean fq data obtained in step S1 to obtain clean spliced ​​data, annotating the clean spliced ​​data using the NR database and the RdRP database, and removing the sequences in the clean spliced ​​data that are aligned to the non-viral NT database to obtain the data to be analyzed and the annotation information of the data to be analyzed;

[0015] The non-viral NT database is an NT database after the virus sequences in the NT database are removed;

[0016] The RdRP database is a database that summarizes the RdRP data in Redefining the invertebrate RNA virosphere and The evolutionary history of vertebrate RNA viruses and annotates them using the pfamA database, extracts and merges data annotated with RdRP, RNA_pol, VP1 and / or NS5 characters in the RdRP data;

[0017] S3. De-redundancy processing is performed on the data to be analyzed obtained in step S2 to obtain de-redundancy data, and the TPM load of each gene in the de-redundancy data is calculated in the clean spliced ​​data obtained in step S2. The TPM load of the gene annotated with RdRP information in the de-redundancy data is used as a benchmark, and Spearman correlation analysis is performed based on the TPM load of each gene to obtain the correlation between each gene in the de-redundancy data and the gene annotated with RdRP information in the de-redundancy data. Genes with a correlation of >0.79 in the de-redundancy data are selected as segmented virus initial screening data;

[0018] S4. Based on the annotation information of the data to be analyzed obtained in step S2 and the TPM load of each gene in the de-redundant data obtained in step S3, the gene clusters that meet the removal conditions in the segmented RNA virus preliminary screening data obtained in step S3 are removed, that is, the segmented RNA virus data are scored;

[0019] The removal conditions are: gene clusters containing Virus, virus, Viruses, viruses, Phage, phage and / or Riboviria in the annotation information and with a consistency of >30% with the NR database, gene clusters consisting entirely of RdRP sequences, gene clusters containing multiple different RdRP sequences, gene clusters with TPM <200 and gene coverage <20, and / or gene clusters with a virus occurrence frequency <3.

[0020] Preferably, the three different samples in step S1 are samples from the same source and / or different sources.

[0021] More preferably, the three different samples in step S1 are samples from the same source.

[0022] Further preferably, the sample is a biological sample and / or a non-biological sample.

[0023] More preferably, the biological sample is an animal sample and / or a plant sample, and the non-biological sample is an environmental sample.

[0024] Preferably, the quality control in step S1 is performed using fastp software.

[0025] More preferably, quality control is performed using fastp software with version number 0.22.0 based on default parameters.

[0026] Preferably, in step S1, RiboDetector software is used to remove ribosomal RNA data from the quality-controlled macrotranscriptome sequencing data.

[0027] More preferably, RiboDetector software with version number 0.2.7 is used to remove ribosomal RNA data from the quality-controlled macrotranscriptome sequencing data based on default parameters.

[0028] Preferably, the step S2 of splicing the clean fq data obtained in step S1 is performed by using megahit software and / or spades software.

[0029] Preferably, the method for establishing the non-viral NT database in step S2 is: deleting sequence data whose sequence IDs contain virus, viruses, riboviria, phage, viridae and / or proviral in the NT database to obtain the remaining NT data; using prot.accession2taxid as a reference data set, annotating the remaining NT data with taxid, and using taxonkit software to extract sequence data annotated with taxid in the remaining NT data as original sequences; using the viral reference sequence database in the NCBI database as a reference database, deleting viral sequences in the sequence data not annotated with taxid in the remaining NT data, and obtaining sequences not annotated with taxid after deleting the viral sequences; merging the original sequences and the sequences not annotated with taxid after deleting the viral sequences to obtain the non-viral NT database.

[0030] Preferably, the DOI of Redefining the invertebrate RNA virosphere in step S2 is 10.1038 / nature20167, and the DOI of The evolutionary history of vertebrate RNA viruses is https: / / doi.org / 10.1038 / s41586-018-0012-7.

[0031] Preferably, the annotation information of the data to be analyzed in step S2 includes RdRP annotation information and NR annotation information of the data to be analyzed.

[0032] Preferably, the de-redundancy processing in step S3 is to perform nucleic acid de-redundancy and amino acid de-redundancy on the data to be analyzed obtained in step S2 in sequence; before performing amino acid de-redundancy, the data to be analyzed obtained in step S2 is translated into a corresponding amino acid sequence.

[0033] More preferably, the redundancy removal is performed using cd-hit software.

[0034] Further preferably, when the cd-hit software is used to deduplicate the sequences annotated with RdRP data information in the data to be analyzed obtained in step S2, the similarity threshold of the cd-hit software is 0.999; when the cd-hit software is used to deduplicate the sequences not annotated with RdRP data information in the data to be analyzed obtained in step S2, the similarity threshold of the cd-hit software is 0.8.

[0035] Preferably, the gene coverage rate in step S4 is the multi value in the output data when splicing is performed using megahit software in step S2 or the cov value in the output data when splicing is performed using spades software.

[0036] The present invention also claims to protect the use of any of the above methods in constructing a model for analyzing segmented RNA viruses.

[0037] The present invention also claims a model for analyzing segmented RNA viruses, comprising a data acquisition module, an analysis module, and a result output module;

[0038] The data acquisition module is used to obtain the macrotranscriptome sequencing data in FastQ format of the sample to be tested; the sample to be tested is no less than 3 different samples;

[0039] The analysis module uses the macrotranscriptome sequencing data acquired by the data acquisition module as input data, and executes any of the above methods to obtain segmented RNA virus data of the sample to be tested;

[0040] The result output module is used to output the segmented RNA virus data of the sample to be tested obtained by the analysis module.

[0041] The present invention also claims the use of any of the above methods and / or models in analyzing segmented RNA viruses.

[0042] Preferably, the analyzing of segmented RNA viruses is obtaining complete fragments of segmented RNA virus genomes and / or discovering segmented RNA viruses.

[0043] Compared with the prior art, the present invention has the following beneficial effects:

[0044] The present invention provides a method for analyzing segmented RNA viruses based on macrotranscriptome sequencing, wherein the method annotates the macrotranscriptome sequencing data of the sample to be tested with the NR database and the RdRP database, and performs Spearman correlation analysis on each sequence in the macrotranscriptome sequencing data, and screens based on the annotation results, the correlation analysis results, and the sequence information in the macrotranscriptome sequencing data to obtain segmented RNA virus data. The method can discover segmented RNA viruses and obtain complete fragments of segmented RNA virus genomes, and can not only realize accurate detection of existing segmented RNA viruses, but also discover new segmented RNA viruses in biological samples or environmental samples with unknown segmented RNA virus information. BRIEF DESCRIPTION OF THE DRAWINGS

[0045] Figure 1 A flowchart of a method for analyzing segmented RNA viruses based on metatranscriptome sequencing;

[0046] Figure 2 The flowchart of establishing the non-viral NT database in Example 1 is as follows;

[0047] Figure 3 This is a bar graph of the segmented RNA virus of Drosophila melanogaster obtained by analyzing according to the method shown in Example 1. DETAILED DESCRIPTION

[0048] The present invention is further described in detail below in conjunction with the accompanying drawings and specific examples of the specification. The examples are only used to explain the present invention and are not used to limit the scope of the present invention. The test methods used in the following examples are conventional methods unless otherwise specified; the materials and reagents used are reagents and materials that can be obtained from commercial channels unless otherwise specified.

[0049] Example 1 A method for analyzing segmented RNA viruses based on macrotranscriptome sequencing

[0050] A method for analyzing segmented RNA viruses based on macrotranscriptome sequencing, the flow chart is as follows Figure 1 As shown, at least 3 different animal, plant or environmental samples are used as test samples, and the test samples are sequenced by macrotranscriptome sequencing, and the macrotranscriptome sequencing data in FastQ format in the offline data is obtained, as follows:

[0051] S1. Use fastp software with version number 0.22.0 based on default parameters to perform quality control on the macrotranscriptome sequencing data in FastQ format to obtain the quality-controlled macrotranscriptome sequencing data, and then use RiboDetector software with version number 0.2.7 based on default parameters to remove ribosomal RNA data from the quality-controlled macrotranscriptome sequencing data to obtain clean fq data;

[0052] S2. The clean fq data obtained in step S1 are spliced ​​using megahit software to obtain clean spliced ​​data in fasta format, and then the clean spliced ​​data are annotated using the NR database and the RdRP database in the NCBI database, and the sequences in the clean spliced ​​data that are aligned to the non-viral NT database are removed to obtain the data to be analyzed and the annotation information of the data to be analyzed; the annotation information of the data to be analyzed includes the RdRP annotation information and the NR annotation information of the data to be analyzed;

[0053] The non-viral NT database is the NT database after the virus sequences in the NT database are removed. The flow chart of establishing the non-viral NT database is as follows: Figure 2 As shown, specifically: delete the sequence data with sequence IDs containing virus, viruses, riboviria, phage, viridae and / or proviral in the NT database to obtain the remaining NT data; use prot.accession2taxid as a reference data set, annotate the remaining NT data with taxid, and use taxonkit software to extract the sequence data annotated with taxid in the remaining NT data as the original sequence; use the virus reference sequence database in the NCBI database as a reference database, delete the virus sequence in the sequence data not annotated with taxid in the remaining NT data, and obtain the sequence not annotated with taxid after deleting the virus sequence; merge the original sequence and the sequence not annotated with taxid after deleting the virus sequence, that is, the non-viral NT database;

[0054] The RdRP database was constructed by summarizing the RdRP data in Redefining the invertebrate RNA virosphere (DOI: 10.1038 / nature20167) and The evolutionary history of vertebrate RNA viruses (DOI: https: / / doi.org / 10.1038 / s41586-018-0012-7) and annotating them using the pfamA database, extracting and merging the data annotated with RdRP, RNA_pol, VP1 or NS5 characters in the RdRP data to obtain the RdRP database.

[0055] S3. The data to be analyzed obtained in step S2 is sequentially subjected to nucleic acid de-redundancy and amino acid de-redundancy using cd-hit software to obtain de-redundant data, and the TPM loading of each gene in the de-redundant data is calculated from the clean spliced ​​data obtained in step S2;

[0056] Taking the TPM load of genes annotated with RdRP information in the de-redundant data as the benchmark, Spearman correlation analysis was used to analyze the TPM load of each gene to obtain the correlation between each gene in the de-redundant data and the genes annotated with RdRP information in the de-redundant data, and genes with correlation > 0.79 in the de-redundant data were selected as the initial screening data for segmented RNA viruses;

[0057] When performing amino acid de-redundancy on the data to be analyzed, the data to be analyzed are translated into the corresponding amino acid sequence; when performing nucleic acid de-redundancy and amino acid de-redundancy on the sequences annotated with RdRP data information in the data to be analyzed, the similarity threshold of the cd-hit software is 0.999; when performing de-redundancy on the sequences not annotated with RdRP data information in the data to be analyzed, the similarity threshold of the cd-hit software is 0.8;

[0058] The calculation formula of TPM loading of each gene is shown in Formula I,

[0059] Formula I: TPM=(Transcripts Per Kilobase of exon model per Million mapped reads) / (Transcript length in kilobases)×10 6 ;

[0060] Where (Transcripts Per Kilobase of exon model per Million mapped reads) is the mapping value of the clean spliced ​​data in fasta format in step S2, and (Transcript length in kilobases) is the length of the clean spliced ​​data in fasta format in step S2; the mapping value of the clean spliced ​​data and the length of the clean spliced ​​data are both output when splicing using the megahit software.

[0061] S4. Based on the annotation information of each gene cluster in the data to be analyzed obtained in step S2 and the TPM load of each gene in the de-redundant data obtained in step S3, the gene clusters that meet the removal conditions in the segmented RNA virus preliminary screening data obtained in step S3 are removed, that is, the segmented RNA virus data are scored;

[0062] The removal conditions are: gene clusters with annotation information containing Virus, virus, Viruses, viruses, Phage, phage, or Riboviria and consistency with the NR database >30%, gene clusters with all RdRP sequences, gene clusters containing multiple different RdRP sequences, gene clusters with TPM load <200 and gene coverage <20, or gene clusters with virus occurrence frequency <3;

[0063] The gene coverage is the multi value in the output data when the megahit software is used for splicing in step S2.

[0064] Example 2 A method for analyzing segmented RNA viruses of Drosophila melanogaster

[0065] 1. Experimental methods

[0066] The macrotranscriptome data of Drosophila melanogaster (PRJNA454554) were downloaded from the NCBI database, and the transcriptome data of Drosophila melanogaster were analyzed according to the method shown in Example 1 to obtain the segmented RNA viruses of Drosophila melanogaster; and a bar graph was drawn with the name of the segmented RNA virus as the horizontal axis and the number of the segmented RNA virus as the vertical axis.

[0067] The segmented RNA virus of Drosophila melanogaster obtained according to the method shown in Example 1 was verified by using the publicly available segmented RNA virus of Drosophila melanogaster (based on the prior art "No detectable effect of Wolbachia wMel on the prevalence and abundance of the RNA virome of Drosophila melanogaster").

[0068] 2. Experimental results

[0069] The segmented RNA viruses of Drosophila melanogaster obtained by analysis according to the method shown in Example 1 are shown in Table 1.

[0070] Table 1 Information on segmented RNA viruses of Drosophila melanogaster

[0071]

[0072]

[0073] The bar graph of segmented RNA viruses of Drosophila melanogaster obtained by analysis according to the method shown in Example 1 is as follows Figure 3 shown.

[0074] The results showed that: by analyzing the macrotranscriptome data of Drosophila melanogaster (PRJNA454554) according to the method shown in Example 1, a total of 15 segmented RNA virus data were obtained, including segmented RNA virus data of 2 Craigies Hill viruses, 2 Dansomanviruses, 3 Galbut viruses, 2 Motts Mill viruses and 6 Torrey Pines viruses; compared with the publicly available segmented RNA viruses of Drosophila melanogaster, the segmented RNA viruses of Drosophila melanogaster can be accurately discovered using the method shown in Example 1, and no false discoveries occur, indicating that the method shown in Example 1 can be used to accurately discover segmented RNA viruses.

[0075] Example 3 A model for analyzing segmented RNA viruses

[0076] A model for analyzing segmented RNA viruses, comprising a data acquisition module, an analysis module, and a result output module;

[0077] The data acquisition module is used to obtain the macrotranscriptome sequencing data in FastQ format of the sample to be tested; the sample to be tested is no less than 3 different animal, plant or environmental samples;

[0078] The analysis module uses the macrotranscriptome sequencing data obtained by the data acquisition module as input data, executes the method shown in Example 1, and obtains segmented RNA virus data of the sample to be tested;

[0079] The result output module is used to output the segmented RNA virus data of the sample to be tested obtained by the analysis module.

[0080] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention rather than to limit the protection scope of the present invention. For ordinary technicians in this field, other different forms of changes or modifications can be made based on the above descriptions and ideas. It is not necessary and impossible to list all the implementation methods here. Any modifications, equivalent substitutions and improvements made within the spirit and principles of the present invention should be included in the protection scope of the claims of the present invention.

Claims

1. A method for analyzing segmented RNA viruses based on macrotranscriptome sequencing, characterized in that: The following steps are involved: S1. Take no less than 3 different samples as test samples, perform macrotranscriptome sequencing on the test samples, collect the macrotranscriptome sequencing data in FastQ format from the offline data for quality control, remove the ribosomal RNA data from the macrotranscriptome sequencing data after quality control, and obtain clean fq data; S2. Use megahit software and / or spades software to splice the clean fq data obtained in step S1 to obtain clean spliced ​​data, use the NR database and the RdRP database to annotate the clean spliced ​​data, and remove the sequences in the clean spliced ​​data that are aligned to the non-viral NT database to obtain the data to be analyzed and the annotation information of the data to be analyzed; The non-viral NT database is an NT database after the virus sequences in the NT database are removed; The RdRP database is a database that summarizes and redefines the RdRP data in the invertebrate RNA virus circle and the evolutionary history of vertebrate RNA viruses and annotates them using the pfamA database, extracts and merges data annotated with RdRP, RNA_pol, VP1 and / or NS5 characters in the RdRP data; S3. De-redundancy processing is performed on the data to be analyzed obtained in step S2 to obtain de-redundancy data, and the TPM load of each gene in the de-redundancy data is calculated in the clean spliced ​​data obtained in step S2. The TPM load of the gene annotated with RdRP information in the de-redundancy data is used as a benchmark, and Spearman correlation analysis is used to analyze based on the TPM load of each gene to obtain the correlation between each gene in the de-redundancy data and the gene annotated with RdRP information in the de-redundancy data. Genes with a correlation of >0.79 in the de-redundancy data are selected as segmented virus initial screening data; S4. Based on the annotation information of the data to be analyzed obtained in step S2 and the TPM load of each gene in the de-redundant data obtained in step S3, the gene clusters that meet the removal conditions in the segmented RNA virus preliminary screening data obtained in step S3 are removed, that is, the segmented RNA virus data are scored; The removal conditions are: gene clusters containing Virus, virus, Viruses, viruses, Phage, phage and / or Riboviria in the annotation information and having a consistency with the NR database of >30%, gene clusters consisting entirely of RdRP sequences, gene clusters containing multiple different RdRP sequences, gene clusters with TPM <200 and gene coverage <20, and / or gene clusters with a virus occurrence frequency <3; The gene coverage rate is the multi value in the output data when the megahit software is used for splicing in step S2 or the cov value in the output data when the spades software is used for splicing.

2. The method according to claim 1, characterized in that The quality control described in step S1 is performed using fastp software.

3. The method according to claim 1, characterized in that In step S1, RiboDetector software is used to remove ribosomal RNA data from the metatranscriptome sequencing data after quality control.

4. The method according to claim 1, characterized in that: The de-redundancy processing in step S3 is to sequentially perform nucleic acid de-redundancy and amino acid de-redundancy on the data to be analyzed obtained in step S2; Before amino acid redundancy removal, the data to be analyzed obtained in step S2 is translated into a corresponding amino acid sequence.

5. The method according to claim 4, characterized in that The redundancy removal is performed by using cd-hit software.

6. The method according to claim 5, characterized in that When the cd-hit software is used to deduplicate the sequences annotated with RdRP data information in the data to be analyzed obtained in step S2, the similarity threshold of the cd-hit software is 0.999; when the cd-hit software is used to deduplicate the sequences not annotated with RdRP data information in the data to be analyzed obtained in step S2, the similarity threshold of the cd-hit software is 0.

8.

7. Use of the method according to any one of claims 1 to 6 in constructing a model for analyzing segmented RNA viruses.

8. A model for analyzing segmented RNA viruses, characterized in that: It includes data acquisition module, analysis module and result output module; The data acquisition module is used to obtain the macrotranscriptome sequencing data in FastQ format of the sample to be tested; the sample to be tested is no less than 3 different samples; The analysis module uses the macrotranscriptome sequencing data acquired by the data acquisition module as input data, and executes the method according to any one of claims 1 to 6 to obtain segmented RNA virus data of the sample to be tested; The result output module is used to output the segmented RNA virus data of the sample to be tested obtained by the analysis module.

9. Use of the method according to any one of claims 1 to 6 and / or the model according to claim 8 in analyzing segmented RNA viruses.

Images