CN112750501B - 一种宏病毒组流程的优化分析方法 - Google Patents
一种宏病毒组流程的优化分析方法 Download PDFInfo
- Publication number
- CN112750501B CN112750501B CN202011596555.9A CN202011596555A CN112750501B CN 112750501 B CN112750501 B CN 112750501B CN 202011596555 A CN202011596555 A CN 202011596555A CN 112750501 B CN112750501 B CN 112750501B
- Authority
- CN
- China
- Prior art keywords
- sequence
- sequences
- virus
- contig
- quality
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
- 238000004458 analytical method Methods 0.000 title claims abstract description 15
- 241000700605 Viruses Species 0.000 title claims description 44
- 238000003908 quality control method Methods 0.000 claims abstract description 7
- 108090000623 proteins and genes Proteins 0.000 claims description 24
- 102000004169 proteins and genes Human genes 0.000 claims description 11
- 238000012216 screening Methods 0.000 claims description 11
- 238000000034 method Methods 0.000 claims description 9
- 241000894007 species Species 0.000 claims description 9
- 238000007621 cluster analysis Methods 0.000 claims description 5
- 238000001914 filtration Methods 0.000 claims description 4
- 230000003612 virological effect Effects 0.000 abstract description 7
- 239000000203 mixture Substances 0.000 description 4
- 102000004196 processed proteins & peptides Human genes 0.000 description 4
- 108090000765 processed proteins & peptides Proteins 0.000 description 4
- 244000000001 Virome Species 0.000 description 2
- 238000007405 data analysis Methods 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000012163 sequencing technique Methods 0.000 description 2
- 238000003556 assay Methods 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 238000012165 high-throughput sequencing Methods 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 238000012545 processing Methods 0.000 description 1
Classifications
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B40/00—ICT specially adapted for biostatistics; ICT specially adapted for bioinformatics-related machine learning or data mining, e.g. knowledge discovery or pattern finding
-
- G—PHYSICS
- G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
- G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
- G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
Landscapes
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Health & Medical Sciences (AREA)
- Engineering & Computer Science (AREA)
- Medical Informatics (AREA)
- Biophysics (AREA)
- Theoretical Computer Science (AREA)
- Spectroscopy & Molecular Physics (AREA)
- General Health & Medical Sciences (AREA)
- Evolutionary Biology (AREA)
- Biotechnology (AREA)
- Bioinformatics & Computational Biology (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Data Mining & Analysis (AREA)
- Public Health (AREA)
- Artificial Intelligence (AREA)
- Evolutionary Computation (AREA)
- Epidemiology (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Analytical Chemistry (AREA)
- Chemical & Material Sciences (AREA)
- Bioethics (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明公开了一种宏病毒组流程的优化分析方法，1)序列质控，2)序列拼接，3)病毒序列预测，4)病毒序列聚类处理，5)物种注释，本发明分析数据更为准确,能拼接出长度较长的contig，从而提高比对结果的准确性和注释结果的可靠性。
Description
技术领域
本发明涉及基因检测技术领域，尤其涉及一种宏病毒组流程的优化分析方法。
背景技术
宏病毒组测序是指对样本种所有病毒基因组进行高通量测序，通过分析测序序列的构成分析特定环境中病毒群体的构成情况或基因的组成以及功能。借助不同环境下病毒的构成差异分析我们可以分析病毒与宿主之间的关系，寻找标志性病毒或特定功能的基因。为了实现着一些目的。当前运用于宏病毒组分析的最常用方法，是基于reads分析方法。该方法的原理是根据reads与病毒数据库进行比对，将比对上的序列注释为该病毒。但是，由于reads比对的可靠性不高，准确性较差，所得到的结果往往会受到广大的质疑，因此，需要引入基于拼接的数据分析方法。
而现有的基于reads分析方法具有如下缺陷：
(1)Reads太短，往往只有150bp，由此去比对，往往得到的结果为假阳性；
(2)基于reads的分析方法，其结果往往包含了大量的病毒结果，而大部分结果都不可信。
发明内容
本发明的提供一种宏病毒组流程的优化分析方法。
本发明的方案是：
一种宏病毒组流程的优化分析方法，包括下列步骤：
1)序列质控，使用fastp软件，将原始序列进行质量过滤，获得高质量序列；
2)序列拼接，使用megahit软件包对上一步得到的高质量序列进行拼接，得到拼接好的contig文件；
3)病毒序列预测，将上一步得到的contig文件使用virsorter软件进行病毒序列的预测，得到预测为病毒的contig序列；
4)病毒序列聚类处理，将步骤3)中挑选预测为病毒的contig序列，使用vcontact2软件，对序列进行聚类分析，获得聚类的结果文件；
5)物种注释，将上一步得到的聚类结果，使用lca的算法，计算出病毒序列的注释结果。
作为优选的技术方案，所述步骤1)中的原始数据是采用Cutadapt识别3'端潜在的接头序列，并在识别的接头序列处截断；与接头序列的匹配长度≥3bp，允许≤20％的碱基错配率。
作为优选的技术方案，所述步骤1)中fast软件p采用滑动窗口法对序列进行质量筛查：窗口大小为5bp，从5'端第一个碱基位置开始移动，要求窗口中碱基平均质量≥Q20，从第一个平均质量值低于Q20的窗口的3'端碱基处截断序列；经上述质量筛查后，去除序列长度小于50bp的序列；经上述质量筛查后，去除序列中含有模糊碱基的序列。
作为优选的技术方案，所述步骤3)中所述contig长度大于1kb；Virsorter取出其中的category1，category2，category4，category5；余下大于5kb的contigs提取这部分contigs的蛋白序列做vHMM。
作为优选的技术方案，所述contigs的蛋白序列做vHMM如下：
一个contigs至少要求五个蛋白有HMM的hit，占10％以上，有KO注释的基因小于20％，PFAM注释的基因小于40％；
有HMM hit的基因数大于等于PFAM hit的基因数,；
HMM hit大于等于60％的总基因数。
作为优选的技术方案，所述步骤5)中使用vContact2软件对预测为病毒的contig序列进行物种注释。
本发明提供了一种宏病毒组流程的优化分析方法，1)序列质控，使用fastp软件，将原始序列进行质量过滤，获得高质量序列；2)序列拼接，使用megahit软件包对上一步得到的高质量序列进行拼接，得到拼接好的contig文件；3)病毒序列预测，将上一步得到的contig文件使用virsorter软件进行病毒序列的预测，得到预测为病毒的contig序列；4)病毒序列聚类处理，将步骤3)中挑选预测为病毒的contig序列，使用vcontact2软件，对序列进行聚类分析，获得聚类的结果文件；5)物种注释，将上一步得到的聚类结果，使用lca的算法，计算出病毒序列的注释结果。
本发明的优点：分析数据相较基于reads的方法更为准确,能拼接出长度较长的contig，从而提高比对结果的准确性和注释结果的可靠性；
以vcontact2软件将病毒contig进行聚类，基于lca算法对病毒contig进行注释，结果可靠，避免了假阳性问题，数据分析更精准。
附图说明
图1为本发明实施例的流程示意图。
具体实施方式
为了弥补以上不足，本发明提供了一种宏病毒组流程的优化分析方法以解决上述背景技术中的问题。
一种宏病毒组流程的优化分析方法，包括下列步骤：
1)序列质控，使用fastp软件，将原始序列进行质量过滤，获得高质量序列；
2)序列拼接，使用megahit软件包对上一步得到的高质量序列进行拼接，得到拼接好的contig文件；
3)病毒序列预测，将上一步得到的contig文件使用virsorter软件进行病毒序列的预测，得到预测为病毒的contig序列；
4)病毒序列聚类处理，将步骤3)中挑选预测为病毒的contig序列，使用vcontact2软件，对序列进行聚类分析，获得聚类的结果文件；
5)物种注释，将上一步得到的聚类结果，使用lca的算法，计算出病毒序列的注释结果。
所述步骤1)中的原始数据是采用Cutadapt识别3'端潜在的接头序列，并在识别的接头序列处截断；与接头序列的匹配长度≥3bp，允许≤20％的碱基错配率。
所述步骤1)中fast软件p采用滑动窗口法对序列进行质量筛查：窗口大小为5bp，从5'端第一个碱基位置开始移动，要求窗口中碱基平均质量≥Q20，从第一个平均质量值低于Q20的窗口的3'端碱基处截断序列；经上述质量筛查后，去除序列长度小于50bp的序列；经上述质量筛查后，去除序列中含有模糊碱基的序列。
所述步骤3)中所述contig长度大于1kb；Virsorter取出其中的category1，category2，category4，category5；余下大于5kb的contigs提取这部分contigs的蛋白序列做vHMM。
所述contigs的蛋白序列做vHMM如下：
一个contigs至少要求五个蛋白有HMM的hit，占10％以上，有KO注释的基因小于20％，PFAM注释的基因小于40％；
有HMM hit的基因数大于等于PFAM hit的基因数,；
HMM hit大于等于60％的总基因数。
所述步骤5)中使用vContact2软件对预测为病毒的contig序列进行物种注释。
为了使本发明实现的技术手段、创作特征、达成目的与功效易于明白了解，下面结合具体实施例，进一步阐述本发明。
实施例：
在步骤S101中，进行序列质控，将现有的原始序列进行质量控制处理，得到高质量序列。
在本发明实施例中，采用Cutadapt识别3'端潜在的接头序列(极少数测通的情况会出现)，并在识别的接头序列处截断。要求与接头序列(R1:AGATCGGAAGAGCACACGTCTGAACTCCAGTCA；R2:AGATCGGAAGAGCGTCGTGTAGGGAAAGAGTGT)的匹配长度至少达到3bp，且允许至多20％的碱基错配率。
去除3'端接头序列后，使用fastp采用滑动窗口法对序列进行质量筛查：窗口大小为5bp，从5'端第一个碱基位置开始移动，要求窗口中碱基平均质量大于等于Q20(即碱基平均测序准确率大于99％)，从第一个平均质量值低于Q20的窗口的3'端碱基处截断序列。
经上述质量筛查后，去除序列长度小于50bp的序列。
经上述质量筛查后，去除序列中含有模糊碱基的序列。
在步骤S102中，使用megahit对上一步得到的高质量蓄力进行拼接，得到拼接好的contig文件。
在本发明实施例中，首先以K-mer～[27,127]的参数设置，调用megahit(https://hku-bal.github.io/megabox/)，对每个样本的双端序列分别进行从头(De novo)组装拼接，通过De Bruijn图构建Contigs序列，并对生成的Contigs序列进行组装拼接效果评价。
在步骤S103中，将上一步得到的contig使用virsorter软件进行病毒序列的预测，得到预测为病毒的contig序列。
在本发明实施例中，先使用Virsorter预测contigs(只取长度大于1kb；Virsorter更新PFAM数据库为最新)，取出其中的category1，2，4，5(依据文献“Uncovering Earth’svirome”等)；余下的contigs(要求大于5kb)，提取这部分contigs的蛋白序列做vHMM(依据文献Uncovering Earth’s virome，下载HMM模型)，要求如下：
即：1.一个contigs至少要求五个蛋白有HMM的hit，且占10％以上，有KO注释的基因小于20％，PFAM注释的基因小于40％。
2.有HMM hit的基因数大于等于PFAM hit的基因数。(与以上并列)
3.HMM hit大于等于60％的总基因数。(与以上并列)。
在步骤S104中，挑选预测为病毒的contig序列，使用vcontact2软件，对序列进行聚类分析，获得聚类的结果文件。
在步骤S105中，将上一步得到的聚类结果，使用lca的算法，计算出病毒序列的注释结果。
在本发明实施例中，使用vContact2软件对上述预测为病毒的contig序列进行物种注释，由于每一条目标序列可能和多多条参考序列聚类在一起，而这些匹配的参考序列又分属不同的分类单元，为使分析严谨可靠，同时又不丢失生物学意义，我们采取“最低共同祖先(Lowest Common Ancestor，LCA)”算法，将参考序列分化为不同物种分枝前的最后一级共同分类，作为目标序列的物种分类注释信息。
以上显示和描述了本发明的基本原理、主要特征及本发明的优点。本行业的技术人员应该了解，本发明不受上述实施例的限制，上述实施例和说明书中描述的只是说明本发明的原理，在不脱离本发明精神和范围的前提下，本发明还会有各种变化和改进，这些变化和改进都落入要求保护的本发明范围内。本发明要求保护范围由所附的权利要求书及其等效物界定。
Claims (1)
1.一种宏病毒组流程的优化分析方法，其特征在于，包括下列步骤：
1）序列质控，使用fastp软件，将原始序列进行质量过滤，获得高质量序列；原始数据是采用 Cutadapt识别 3' 端潜在的接头序列 ，并在识别的接头序列处截断；与接头序列的匹配长度≥3 bp，允许≤ 20% 的碱基错配率；fastp软件采用滑动窗口法对序列进行质量筛查：窗口大小为 5 bp，从 5' 端第一个碱基位置开始移动，要求窗口中碱基平均质量≥Q20 ，从第一个平均质量值低于 Q20 的窗口的 3' 端碱基处截断序列；经上述质量筛查后，去除序列长度小于 50 bp 的序列；经上述质量筛查后，去除序列中含有模糊碱基的序列；
2）序列拼接，使用megahit软件包对上一步得到的高质量序列进行拼接，得到拼接好的contig文件；
3）病毒序列预测，将上一步得到的contig文件使用virsorter软件进行病毒序列的预测，得到预测为病毒的contig序列；所述contig长度大于1kb；Virsorter 取出其中的category1，category2，category4，category5；余下大于5kb的contigs提取这部分contigs的蛋白序列做vHMM；所述contigs的蛋白序列做vHMM如下：
一个contigs至少要求五个蛋白有HMM的hit，占10%以上，有KO注释的基因小于20%，PFAM注释的基因小于40%；
有HMM hit的基因数大于等于PFAM hit的基因数;
HMM hit大于等于60%的总基因数；
4）病毒序列聚类处理，将步骤3）中挑选预测为病毒的contig序列，使用vcontact2软件，对序列进行聚类分析，获得聚类的结果文件；
5）物种注释，将上一步得到的聚类结果，使用lca的算法，计算出病毒序列的注释结果。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011596555.9A CN112750501B (zh) | 2020-12-29 | 2020-12-29 | 一种宏病毒组流程的优化分析方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011596555.9A CN112750501B (zh) | 2020-12-29 | 2020-12-29 | 一种宏病毒组流程的优化分析方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN112750501A CN112750501A (zh) | 2021-05-04 |
| CN112750501B true CN112750501B (zh) | 2024-04-02 |
Family
ID=75646795
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202011596555.9A Active CN112750501B (zh) | 2020-12-29 | 2020-12-29 | 一种宏病毒组流程的优化分析方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN112750501B (zh) |
Families Citing this family (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN113539378A (zh) * | 2021-07-16 | 2021-10-22 | 明科生物技术(杭州)有限公司 | 一种病毒数据库的数据分析方法、系统、设备及存储介质 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2008000186A1 (fr) * | 2006-06-21 | 2008-01-03 | Beijing Bioway-Fortune Research Center For Gene Drugs Ltd. | Méthode d'identification d'un nouveau gène et nouveaux gènes résultants |
| CN107609347A (zh) * | 2017-08-21 | 2018-01-19 | 上海派森诺生物科技股份有限公司 | 一种基于高通量测序技术的宏转录组数据分析方法 |
| CN110033826A (zh) * | 2018-12-10 | 2019-07-19 | 上海派森诺生物科技股份有限公司 | 一种应用于宏病毒组高通量测序数据的分析方法 |
| WO2020250068A1 (en) * | 2019-06-14 | 2020-12-17 | University College Cork – National University Of Ireland, Cork | Materials and methods for assessing virome and microbiome matter |
-
2020
- 2020-12-29 CN CN202011596555.9A patent/CN112750501B/zh active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| WO2008000186A1 (fr) * | 2006-06-21 | 2008-01-03 | Beijing Bioway-Fortune Research Center For Gene Drugs Ltd. | Méthode d'identification d'un nouveau gène et nouveaux gènes résultants |
| CN107609347A (zh) * | 2017-08-21 | 2018-01-19 | 上海派森诺生物科技股份有限公司 | 一种基于高通量测序技术的宏转录组数据分析方法 |
| CN110033826A (zh) * | 2018-12-10 | 2019-07-19 | 上海派森诺生物科技股份有限公司 | 一种应用于宏病毒组高通量测序数据的分析方法 |
| WO2020250068A1 (en) * | 2019-06-14 | 2020-12-17 | University College Cork – National University Of Ireland, Cork | Materials and methods for assessing virome and microbiome matter |
Also Published As
| Publication number | Publication date |
|---|---|
| CN112750501A (zh) | 2021-05-04 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN111341383B (zh) | 一种检测拷贝数变异的方法、装置和存储介质 | |
| CN108197434B (zh) | 去除宏基因组测序数据中人源基因序列的方法 | |
| Pu et al. | Detection and analysis of ancient segmental duplications in mammalian genomes | |
| WO2022267867A1 (zh) | 基因测序分析方法、装置、存储介质和计算机设备 | |
| CN110021355B (zh) | 二倍体基因组测序片段的单倍体分型和变异检测方法和装置 | |
| Saheb Kashaf et al. | Recovering prokaryotic genomes from host-associated, short-read shotgun metagenomic sequencing data | |
| Lange et al. | Analysis pipelines for cancer genome sequencing in mice | |
| KR20140006846A (ko) | Dna 서열의 데이터 분석 | |
| CN112599198A (zh) | 一种用于宏基因组测序数据的微生物物种与功能组成分析方法 | |
| CN112786102A (zh) | 一种基于宏基因组学分析精准识别水体中未知微生物群落的方法 | |
| CN110621785A (zh) | 基于三代捕获测序对二倍体基因组单倍体分型的方法和装置 | |
| CN112750501B (zh) | 一种宏病毒组流程的优化分析方法 | |
| CN110648722B (zh) | 新生儿遗传病患病风险评估的装置 | |
| CN113539378A (zh) | 一种病毒数据库的数据分析方法、系统、设备及存储介质 | |
| CN116189763A (zh) | 一种基于二代测序的单样本拷贝数变异检测方法 | |
| CN113066532B (zh) | 基于高通量测序技术的宿主中病毒来源sRNA数据分析方法 | |
| CN107862177B (zh) | 一种区分鲤群体的单核苷酸多态性分子标记集的构建方法 | |
| KR20210040714A (ko) | 핵산 서열 분석에서 위양성 변이를 검출하는 방법 및 장치 | |
| CN115295084A (zh) | 一种肿瘤新抗原免疫组库数据可视化分析方法和系统 | |
| Lara et al. | A web tool to discover full-length sequences—Full-Lengther | |
| US20210363589A1 (en) | Immunotherapy using multi-omics data to extract microsatellite instability-based neoantigen | |
| JP2021503128A (ja) | アラインされていないシーケンシングデータの高速品質管理のためのｋ−ｍｅｒの使用 | |
| CN114496073B (zh) | 用于识别阳性重排的方法、计算设备和计算机存储介质 | |
| Deshpande et al. | RNA-seq data science: From raw data to effective interpretation | |
| CN115198036B (zh) | 一种基于纳米孔和高通量测序数据的噬菌体鉴定和宿主预测方法 |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant | ||
| EE01 | Entry into force of recordation of patent licensing contract |
Application publication date: 20210504 Assignee: Shanghai Boyio Gene Technology Co.,Ltd. Assignor: SHANGHAI PERSONAL BIOTECHNOLOGY Co.,Ltd. Contract record no.: X2024980004994 Denomination of invention: An Optimization Analysis Method for Macro Virus Group Process Granted publication date: 20240402 License type: Common License Record date: 20240428 |
|
| EE01 | Entry into force of recordation of patent licensing contract |