CN111798922A - 基于重测序数据中多态性位点密度鉴定小麦育种的基因组选择利用区间的方法 - Google Patents
基于重测序数据中多态性位点密度鉴定小麦育种的基因组选择利用区间的方法 Download PDFInfo
- Publication number
- CN111798922A CN111798922A CN202010742489.5A CN202010742489A CN111798922A CN 111798922 A CN111798922 A CN 111798922A CN 202010742489 A CN202010742489 A CN 202010742489A CN 111798922 A CN111798922 A CN 111798922A
- Authority
- CN
- China
- Prior art keywords
- snp
- interval
- variation
- crop
- genome
- 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.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 71
- 241000209140 Triticum Species 0.000 title claims abstract description 37
- 235000021307 Triticum Nutrition 0.000 title claims abstract description 37
- 238000009395 breeding Methods 0.000 title claims abstract description 34
- 230000001488 breeding effect Effects 0.000 title claims abstract description 34
- 239000000463 material Substances 0.000 claims abstract description 124
- 238000009826 distribution Methods 0.000 claims abstract description 55
- 238000012163 sequencing technique Methods 0.000 claims abstract description 32
- 238000001914 filtration Methods 0.000 claims abstract description 13
- 108020004414 DNA Proteins 0.000 claims abstract description 11
- 230000002068 genetic effect Effects 0.000 claims abstract description 8
- 108091028043 Nucleic acid sequence Proteins 0.000 claims abstract description 6
- 241000196324 Embryophyta Species 0.000 claims description 16
- 244000098338 Triticum aestivum Species 0.000 claims description 8
- 230000008569 process Effects 0.000 claims description 8
- 241000209504 Poaceae Species 0.000 claims description 6
- 238000004422 calculation algorithm Methods 0.000 claims description 6
- 244000038559 crop plants Species 0.000 claims description 5
- 230000005540 biological transmission Effects 0.000 claims description 3
- 238000012937 correction Methods 0.000 claims description 3
- 238000005457 optimization Methods 0.000 claims description 3
- 241000894007 species Species 0.000 claims description 3
- 238000011161 development Methods 0.000 claims description 2
- 238000012545 processing Methods 0.000 claims description 2
- 239000002773 nucleotide Substances 0.000 abstract description 4
- 125000003729 nucleotide group Chemical group 0.000 abstract description 4
- 238000003908 quality control method Methods 0.000 abstract description 4
- 238000004458 analytical method Methods 0.000 description 11
- 230000003321 amplification Effects 0.000 description 7
- 238000003199 nucleic acid amplification method Methods 0.000 description 7
- 238000012549 training Methods 0.000 description 6
- 239000011159 matrix material Substances 0.000 description 5
- 238000012360 testing method Methods 0.000 description 4
- 238000012795 verification Methods 0.000 description 4
- 238000012800 visualization Methods 0.000 description 4
- 108091027544 Subgenomic mRNA Proteins 0.000 description 3
- 210000000349 chromosome Anatomy 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 3
- 230000006870 function Effects 0.000 description 3
- 238000009499 grossing Methods 0.000 description 3
- 239000000203 mixture Substances 0.000 description 3
- 238000011160 research Methods 0.000 description 3
- LZZYPRNAOMGNLH-UHFFFAOYSA-M Cetrimonium bromide Chemical compound [Br-].CCCCCCCCCCCCCCCC[N+](C)(C)C LZZYPRNAOMGNLH-UHFFFAOYSA-M 0.000 description 2
- 238000000137 annealing Methods 0.000 description 2
- 238000004364 calculation method Methods 0.000 description 2
- 238000006243 chemical reaction Methods 0.000 description 2
- 238000002474 experimental method Methods 0.000 description 2
- 239000012634 fragment Substances 0.000 description 2
- 238000005192 partition Methods 0.000 description 2
- 230000007704 transition Effects 0.000 description 2
- 108700028369 Alleles Proteins 0.000 description 1
- 238000007400 DNA extraction Methods 0.000 description 1
- 238000001712 DNA sequencing Methods 0.000 description 1
- 238000000729 Fisher's exact test Methods 0.000 description 1
- 238000012408 PCR amplification Methods 0.000 description 1
- 239000003242 anti bacterial agent Substances 0.000 description 1
- 229940088710 antibiotic agent Drugs 0.000 description 1
- 239000003153 chemical reaction reagent Substances 0.000 description 1
- 238000013461 design Methods 0.000 description 1
- 238000001514 detection method Methods 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000011156 evaluation Methods 0.000 description 1
- 238000000605 extraction Methods 0.000 description 1
- 235000013305 food Nutrition 0.000 description 1
- 238000007689 inspection Methods 0.000 description 1
- 238000012986 modification Methods 0.000 description 1
- 230000004048 modification Effects 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 230000003252 repetitive effect Effects 0.000 description 1
- 238000007480 sanger sequencing Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 239000000344 soap Substances 0.000 description 1
- 238000012358 sourcing Methods 0.000 description 1
- 230000000007 visual effect Effects 0.000 description 1
Images
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
- G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
- G16B20/20—Allele or variant detection, e.g. single nucleotide polymorphism [SNP] detection
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
- G06F18/214—Generating training patterns; Bootstrap methods, e.g. bagging or boosting
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/29—Graphical models, e.g. Bayesian networks
- G06F18/295—Markov models or related models, e.g. semi-Markov models; Markov random fields; Networks embedding Markov models
-
- 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
- G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
- G16B20/30—Detection of binding sites or motifs
-
- 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
- G16B30/10—Sequence alignment; Homology search
Landscapes
- Engineering & Computer Science (AREA)
- Life Sciences & Earth Sciences (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Bioinformatics & Cheminformatics (AREA)
- Bioinformatics & Computational Biology (AREA)
- Health & Medical Sciences (AREA)
- Evolutionary Biology (AREA)
- Data Mining & Analysis (AREA)
- Analytical Chemistry (AREA)
- Biophysics (AREA)
- Chemical & Material Sciences (AREA)
- Biotechnology (AREA)
- Spectroscopy & Molecular Physics (AREA)
- General Health & Medical Sciences (AREA)
- Proteomics, Peptides & Aminoacids (AREA)
- Medical Informatics (AREA)
- Molecular Biology (AREA)
- Genetics & Genomics (AREA)
- Artificial Intelligence (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Evolutionary Computation (AREA)
- General Engineering & Computer Science (AREA)
- General Physics & Mathematics (AREA)
- Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)
Abstract
本发明公开了基于重测序数据中多态性位点密度鉴定小麦育种的基因组选择利用区间的方法。该方法包括提取待鉴定作物材料基因组DNA并进行全基因组重测序；质控重测序数据并将其比对到作物参考基因组上得到DNA序列比对文件；分别过滤鉴定待鉴定作物材料得到待鉴定作物材料与作物参考基因组相比的排除CNV变异区间的SNP变异位点；计算所述得到的待鉴定作物材料SNP变异位点之间的差异SNP变异位点的密度分布；根据密度分布拟合得到混合正态分布后根据分布边界确定阈值，得到待鉴定作物材料之间的SNP多态性热点区间和SNP低频变异区间。本发明可用于确定不同品种间的遗传差异，从基因组水平对种质资源的相似性进行评估。
Description
技术领域
本发明涉及植物基因组学领域，具体涉及一种基于重测序数据中多态性位点密度鉴定小麦育种的基因组选择利用区间的方法。
背景技术
近年来，关于定位和鉴定作物育种热点选择区间的研究受到广泛关注，同时不同品种间的育种利用相似性鉴定也是品种保护和差异分析的重要评价指标。与之相关的鉴定方法一般是通过比对不同样本内的分子标记扩增结果或者芯片鉴定结果以确定特定区段在不同品种内的选择结果与传递关系。但是芯片技术中检测的位点数量有限，且基于芯片的实验具有较大的批次效应，芯片数据具有较高噪音，鉴定和分析能力有限。此外，芯片数据也不易判断基因组序列的丢失和重复的情况。随着普通小麦参考基因组的公布和基于第二代测序的测序成本的降低，利用大规模重测序数据对样本基因组间进行比较提供了高质量的数据基础，也大大提高了分析的精度，为从全基因组水平准确鉴定相同来源的育种利用区间提供了可能。
普通小麦是一种重要的粮食作物，属于异源六倍体，基因组庞大且重复序列高。研究并确定普通小麦品种间的育种利用异同区间在科研、育种等领域都有重要价值。因此开发一种鉴定小麦育种选择利用的基因组区间的方法对于研究小麦基因组资源的鉴定和利用具有重要实际意义和应用价值。
发明内容
本发明所要解决的技术问题是如何鉴定作物不同材料间的基因组异同区间或如何鉴定作物材料间育种的基因组选择区间或如何鉴定作物材料间的“SNP多态性热点区间”。
为了解决上述技术问题，本发明首先提供了一种方法：
一种鉴定作物育种材料间基因组选择利用区间和/或“SNP多态性热点区间”的方法，包括如下步骤：
A1)提取待鉴定作物材料基因组DNA并进行全基因组重测序获得待鉴定作物重测序数据；
A2)分别过滤名称为材料A和材料B的两个所述待鉴定作物材料重测序数据以对其进行质量控制，并将其分别比对到作物参考基因组上，得到名称为文件A的所述材料A基因组与所述作物参考基因组的DNA序列比对文件，和名称为文件B的所述材料B基因组与所述作物参考基因组的DNA序列比对文件；
A3)对所述文件A进行过滤低质量变异位点，得到所述材料A与所述作物参考基因组相比存在的SNP变异位点和CNV变异区间，然后将所述CNV变异区间中存在的SNP变异位点去除，得到所述材料A基因组与所述作物参考基因组相比存在的排除CNV变异区间的SNP变异位点，将该SNP变异位点称为A-参考SNP变异位点；
和对所述文件B进行过滤低质量变异位点，得到所述材料B与所述作物参考基因组相比存在的SNP变异位点和CNV变异区间，然后将所述CNV变异区间中存在的SNP变异位点去除，得到所述材料B基因组与所述作物参考基因组相比存在的排除CNV变异区间的SNP变异位点，将该SNP变异位点称为B-参考SNP变异位点；
A4)将所述A-参考SNP变异位点和B-参考SNP变异位点进行比较得到二者之间的差异SNP变异位点，计算所述差异SNP变异位点在所述待鉴定作物参考基因组上的密度分布；
A5)根据所述密度分布拟合得到混合正态分布后根据边界确定阈值，并使用所述阈值对混合正态分布进行分割，得到所述材料A和所述材料B之间的“SNP多态性热点区间”和所述材料A和所述材料B之间的“SNP低频变异区间”；所述材料A和所述材料B之间的“SNP多态性热点区间”即为所述材料A和所述材料B之间的基因组差异区间，所述材料A和所述材料B之间的“SNP低频变异区间”即为所述材料A和所述材料B之间的基因组相似区间或育种利用区间。
上述方法中，还包括如下步骤：
A6)对所述材料A和所述材料B之间的“SNP多态性热点区间”和所述材料A和所述材料B之间的“SNP低频变异区间”进行调整和优化，得到优化的所述材料A和所述材料B之间的“SNP多态性热点区间”和优化的所述材料A和所述材料B之间的“SNP低频变异区间”；优化的所述材料A和所述材料B之间的“SNP多态性热点区间”即为优化的所述材料A和所述材料B之间的基因组差异区间，优化的所述材料A和所述材料B之间的“SNP低频变异区间”即为优化的所述材料A和所述材料B之间的基因组相似区间或育种利用区间。
上述方法中，所述材料可为品种和/或品系。
上述方法A5)步骤中的所述阈值是以拟合得到的正态分布的交叉边界确定的，是对基因组区间进行初步划分的硬性标准；使用所述阈值对所述正态分布进行分割得到的密度较大的正态分布中的SNP位点所在区间为所述材料A和所述材料B之间的“SNP多态性热点区间”，密度较小的正态分布中的SNP位点所在区间为所述材料A和所述材料B之间的“SNP低频变异区间”。
上述方法A4)步骤中所述计算所述材料A和所述材料B之间的差异SNP变异位点的密度分布，包括以下步骤：
A4-1)在0.5-5Mbp范围内选择固定区间长度为单位将所述材料A和所述材料B之间全基因组进行划分小窗，
A4-2)对所述A-参考SNP变异位点和所述B-参考SNP变异位点在全基因组每个小窗中的SNP频率进行分析，
A4-3)统计小窗中所述差异SNP变异位点数，依据所述得到的差异SNP变异位点频率数与非丢失位点比例的比值得到所述差异SNP变异位点的密度分布。
上述A4-3)中所述密度分布具体是以10为底数取密度值+1的对数得到的分布；所述密度值为在所述固定区间长度中所述差异SNP位点的个数除以非丢失位点的比例的比值；所述非丢失位点比例为所述材料A和所述材料B的所述固定区间长度中的所有通过测序判定基因型的位点的数量和区间总长度的比例，包括因CNV变异或局部区间覆盖度低等原因产生。所述密度值+1是避免个别小窗在实际计算时获得0值而无法取对数的情形进行矫正。
上述方法A5)步骤中所述拟合得到的正态分布包括两个正态分布，所述拟合得到正态分布的过程，是通过利用“最大期望算法”(Expectation-Maximization algorithm,EM)对所述A4)步骤中得到的所述差异SNP变异位点的密度分布进行计算后拟合得到两个正态分布的均值和方差。
上文所述作物为下述任一种：
D1)双子叶植物，
D2)单子叶植物，
D3)禾本目，
D4)禾本科，
D5)小麦属，
D6)普通小麦。
上述方法A5)步骤中，所述作物为普通小麦，所述阈值为32/1Mbp区间。
上述方法A6)步骤中对两个材料间的“SNP多态性热点区间”和“SNP低频变异区间”进行调整和优化过程，是利用“隐马尔可夫模型”(Hidden Markov Model，HMM)对“多态性热点区间”和“低频率变异区间”进行平滑校正处理，完成对区间划分的微调，减少因硬阈值产生的随机误差，具体包括以下步骤：
C1)确定观测序列和初始参数，
C2)对初始参数进行参数训练得到模型，
C3)运用维特比算法使用观测序列推测状态序列(即相似区间与多态性区间的排列顺序)。
上述方法所述A3)步骤中对所述过滤低质量位点或噪音位点，可使用GATK软件进行，具体可为GATK软件中的HaplotypeCaller、GenotypeGVCF、SelectVariants和/或VariantFiltration功能进行过滤；也可使用其他软件进行，如Samtools、bcftools、Atlas-SNP2、SOAPsnp、SLIDERII、dindel、DeepVariant、SpeedSeq软件等，本发明不做限制。
上述方法A4)步骤中所述作物材料为小麦，所述差异SNP变异位点为纯合的差异SNP位点；杂合的差异SNP位点不被考虑在内。
上述方法A2)步骤中，所述比对过程可使用BWA软件中的BWA-MEM工具，也可使用BWA软件地其他工具或其他比对软件进行，如SOAP、Bowtie2、minimap、SAM/BAM flagexplain、lastz等，本发明不做限制。
上述方法A3)步骤中，所述CNV变异区间可通过如下步骤获得：分别以0.5-5Mbp范围内选择以固定长度为单位将所述材料A和所述材料B全基因组划分为小窗，通过计算所述材料A和所述材料B每个小窗的“平均相对覆盖深度”，选择值低于一定阈值的小窗内的DNA测序片段为相应材料相对于作物参考基因组的CNV变异区间；所述阈值可为30％-50％范围中的一个值(相对于全基因组的平局覆盖深度)；所述步骤也可通过软件过滤分析获得；所述软件可为freec、CNVnator、Genome STRiP软件等，本发明不做限制。
上述方法A6)步骤中所述进行调整和优化得到所述材料A和所述材料B间的“SNP多态性热点区间”和“SNP低频变异区间”，可通过使用R语言编写可视化脚本实现可视化，以进行更直观地观察。
为了解决上述技术问题，本发明还提供了上述方法的下述任一种应用：
P1、上述方法在定位和鉴定作物育种热点选择区间中的应用；
P2、上述方法在鉴定作物品种基因组间关系中的应用；
P3、上述方法在开发作物新种质资源中的应用；
P4、上述方法在作物育种中的应用；
P5、上述方法在分析作物品种育种来源中的应用；
P6、上述方法在鉴定不同作物品种的相同来源育种区间中的应用；
P7、上述方法在鉴定作物的不同品种间的遗传差异中的应用；
P8、上述方法在鉴定作物育种过程中世代间基因组遗传信息的传递情况中的应用；
P9、上述方法在鉴定作物不同品种间的遗传区间利用关系中的应用。
上文中所述作物为下述任一种：
D1)双子叶植物，
D2)单子叶植物，
D3)禾本目，
D4)禾本科，
D5)小麦属，
D6)普通小麦。
本发明的实施例中，通过使用本发明所述的方法得到两个待鉴定小麦材料济麦22和良星99间的“SNP多态性热点区间”和“SNP低频变异区间”，其中“SNP多态性热点区间”为两个样本间的基因组差异区间。通过对所得到的小麦济麦22与良星99间的“SNP多态性热点区间”选择区段设计引物进行PCR检测，发现在所选择的基因组差异区间中济麦22和良星99间存在核苷酸序列差异。本发明所提供的方法可有效地鉴定作物的不同材料间的“SNP低频变异区间”和“SNP多态性热点区间”。
附图说明
图1为本方法的流程示意图。
图2为全基因组小窗“平均相对覆盖深度”的密度分布。
图3中A为全基因组小窗中济麦22与良星99间存在差异的纯合SNP与非丢失位点的密度分布；B为使用EM算法计算A中密度分布的均值与方差后模拟得出的两个正态分布并以此得出区分两者的阈值；该图的分布均以小窗定长为1Mb的前提下计算得出的分布。
图4为使用隐马尔可夫模型对差异区间与相似区间平滑后的全基因组可视化结果，白色区段为CNV变异区间。
图5中为济麦22与良星99的“SNP多态性热点区间”的PCR产物测序验证结果。其中Chinese Spring代表参考基因组中国春，Jimai 22代表济麦22，Liangxing 99代表良星99，Consensus代表济麦22、良星99和中国春三者碱基相同的位点,“★”代表差异SNP位点的位置。
具体实施方式
下面结合具体实施方式对本发明进行进一步的详细描述，给出的实施例仅为了阐明本发明，而不是为了限制本发明的范围。以下提供的实施例可作为本技术领域普通技术人员进行进一步改进的指南，并不以任何方式构成对本发明的限制。
下述实施例中的实验方法，如无特殊说明，均为常规方法，按照本领域内的文献所描述的技术或条件或者按照产品说明书进行。下述实施例中所用的材料、试剂等，如无特殊说明，均可从商业途径得到。
供试材料：待测小麦为国审小麦品种良星99(国审麦2006016)与济麦22(国审麦2006018)品种，分别来源于山东良星种业有限公司和安徽皖垦种业股份有限公司。
实施例1、鉴定小麦育种中基因组选择利用区间和/或“SNP多态性热点区间”的方法的建立
1.待测小麦基因组DNA的提取
取材供试品种材料济麦22和良星99的幼根采用CTAB法提取基因组DNA，得到济麦22和良星99的DNA样品。
2.待测小麦全基因组重测序
将步骤1提取的两个品种的DNA样品，送由北京诺禾致源科技股份有限公司使用Illumina HiSeq2500平台进行二代测序技术全基因组双端重测序，测序深度不低于5X，得到两个品种的原始测序数据Rawdata文件。
3.待测小麦全基因组重测序数据数据分析
3.1重测序数据质量控制和基因组比对
3.1.1重测序数据质量控制
利用Trimmomatic v0.36软件，分别对两个待测小麦的原始测序数据RawData进行质量控制，软件来源为(http://www.usadellab.org/cms/index.php？page＝trimmomatic)，去除低质量的测序数据，具体参数为LEADING:3(reads开头碱基质量值低于该阈值则丢弃)；TRAILING:3(reads末尾碱基质量值低于该阈值则丢弃)；SLIDINGWINDOW:4:15(从reads的5'端开始的滑窗大小)；MINLEN:36(经过剪切后reads的长度该阈值则丢弃)，最后得到两个品种过滤后的数据Cleandata用于后续分析。
3.1.2重测序数据基因组比对
使用BWA软件中的BWA-MEM工具(https://sourceforge.net/projects/bio-bwa/files/)分别将过滤后的两个品种的Cleandata数据比对到中国春小麦参考基因组(IWGSCv1，https://urgi.versailles.inra.fr/download/iwgsc/IWGSC_RefSeq_Annotations/v1.1/)，选择保留存在“唯一最优匹配(unique best hit)”的读段对的数据。
然后利用bamtools(v2.5.1,https://github.com/pezmaster31/bamtools)与samtools(v1.4，http://www.htslib.org/)工具分别对两个品种比对到参考基因组的数据进行过滤。使用bamtools去除读段对间的距离大于10,000bp或为0bp的读段以及回帖质量小于1的读段。使用samtools去除PCR duplicate数据。通过比对和过滤，最后获得两个品种的BAM文件用于后续分析。
3.2过滤鉴定待测小麦的变异位点
3.2.1查找SNP变异位点
使用GATK软件(v3.8，https://gatk.broadinstitute.org/hc/en-us)的HaplotypeCaller、GenotypeGVCF、SelectVariants和VariantFiltration等功能分别对3.1中得到的两个品种的BAM文件进行过滤，以查找两个品种相对于小麦中国春的参考基因组的单核苷酸多态性信息(single nucleotide polymorphism，SNP)。其中，SNP位点过滤参数设定为“QD<2.0，FS>60.0，MQRankSum<-12.5，ReadPosRankSum<-8.0，SOR>3.0，MQ<40.0，DP>30||DP<3”。其中各参数含义如下：QD：QualByDepth由等位基因位点深度标准化后的质量分数；FS：FisherStrand链方向Fisher精确检验值；MQRankSum：MappingQualityRankSumTest回帖质量的秩和检验值；SOR：StrandOddsRatio链比值比；ReadPosRankSum：Read PositionRank Sum Test读段内位点所在位置秩和检验；MQ：RMSMappingQuality位点所有读段回帖质量的均方根；DP：Depth位点读段测序深度。最后分别得到两个品种相对于中国春参考基因组的SNP变异信息的名称为VCF的文件，即济麦22-中国春VCF文件和良星99-中国春VCF文件。
3.2.2鉴定CNV变异
分别将两个品种以1Mbp为单位将全基因组划分为小窗，利用bedtools(v2.26.0，https://bedtools.readthedocs.io/en/latest/)软件计算每个品种的重测序比对读段在每个窗口内的“平均覆盖深度”(Depbin)；并结合该材料的全基因组“平均读段覆盖深度”(Depave)将两者相除进行归一化，得到每个小窗的“平均相对覆盖深度”：Depbin/Depave。Depbin/Depave值低于0.5的小窗被视为待测小麦品种相对于中国春参考基因组的“CNV变异区间”(图2)，最终分别得到两个品种的含有CNV变异区间信息的BED文件，即济麦22-中国春BED文件和良星99-中国春BED文件。在济麦22-中国春BED文件和良星99-中国春BED文件中同时存在的CNV变异区间被称为两个品种的“共有CNV变异区间”。
将3.1.1中鉴定得到的SNP变异统计文件济麦22-中国春VCF文件和良星99-中国春VCF文件分别去除3.2.2中鉴定出的CNV变异区间中的SNP变异位点，最后分别得到的排除CNV变异区间的SNP变异位点的文件VCF文件，即济麦22-中国春VCF文件(去除CNV区间SNP)和良星99-中国春VCF文件(去除CNV区间SNP)用于后续分析。此方法仅对拷贝数正常的区间进行分析。同时将所得济麦22-中国春VCF文件(去除CNV区间SNP)和良星99-中国春VCF文件(去除CNV区间SNP)进行比较得到二者之间的差异SNP变异位点。
4.计算SNP多态性位点地密度分布
利用3中得到的差异SNP变异位点计算其在小麦参考基因组上的密度分布，以鉴定济麦22和良星99两个品种的基因组序列差异区间。
首先分别对两个品种以1Mbp为单位对全基因组进行划分小窗，以对每个小窗中SNP频率进行分析。由于小麦是自花授粉作物，其基因组中的杂合位点较少，此处仅考虑纯合的SNP变异位点，经测序和计算得到的杂合SNP位点不被考虑在内。统计每个小窗中纯合SNP位点数量以及非丢失位点(或可判断基因型位点)的数量。
然后统计小窗中两品种间存在差异的纯合SNP数与非丢失位点比例比值的密度，以10为底数取密度值+1的对数，使用R语言软件(https://cran.r-project.org/)将所得密度值+1的对数作图，得到密度分布，该密度分布通常可被视为两个类似正态分布的混合分布(图3中A所示)。其中非丢失位点比例可通过计算定长的小窗区间中可判定基因型的位点的数量和区间总长度的比值获得；该部分通常为个别小窗中的覆盖度较低或出现CNV丢失情况时对密度值进行适当校正。
5.确定多态性热点区间和低频变异区间
由于步骤4中所得的密度分布具有明显的高斯混合分布的特点，因此利用最大期望算法(Expectation-Maximization algorithm,EM)通过计算步骤4中所得密度分布拟合出混合的两个正态分布的均值与方差，该拟合过程采用R语言mixtools包中的normalmixEM功能按默认参数进行，按拟合结果，根据混合正态分布的交叉边界确定划分阈值，选择阈值32(/1Mbp区间)作为两个分布的分割点(图3中B所示)。对拟合的两个正态分布进行分割后，得到两个类基因组区间，其中，密度较大的正态分布中的SNP位点被认为是两个品种间“SNP多态性热点区间”的SNP位点，由其组成的SNP位点区间为“SNP多态性热点区间”；密度较小的正态分布中的SNP位点被认为是两个品种间“SNP低频变异区间”的SNP位点，由其组成的SNP位点区间为“SNP低频变异区间”。
6.多态性热点区间和低频变异区间进行调整和优化
对步骤5中确定的两个品种间的相似区间与多态性热点区间利用基于python中hmmlearn包的hmm.MultinomialHMM功能利用隐马尔可夫模型(Hidden Markov Model，HMM)进行平滑，得到优化的两个品种的“SNP多态性热点区间”和优化的“SNP低频变异区间”。
首先将步骤5中的相似区间(SNP低频变异区间)与多态性区间(SNP多态性热点区间)的排列顺序作为观测序列，统计初始参数即相似区间之间、相似区间与多态性热点区间之间、多态性热点区间之间的转换频率。
对初始状态概率向量π按0.33333进行平均分配，并基于初始参数对得到的观测概率矩阵(发射矩阵)与状态转移矩阵进行参数训练，A亚基因组训练16次，B亚基因组训练19次，D亚基因组训练32次，得到最终的状态转移矩阵，发射矩阵，初始状态向量得到模型。
利用上述三个训练结果得到隐马尔可夫模型，之后运用维特比算法使用观测序列推测状态序列，即按照当前的模型与观测序列，确定每个区间概率最高的状态。状态序列就是利用隐马尔科夫模型平滑后的相似区间与多态性区间在小麦基因组染色体上的排列顺序。
7.分析结果可视化
使用R语言软件(https://cran.r-project.org/)对步骤7的结果进行可视化，得到更直观的分布于各个染色体上的两个品种的“SNP低频变异区间”与“SNP多态性热点区间”的位置信息(图4中所示)。图4中，灰色区段为优化的“SNP低频变异区间”，即所述两个品种在作物育种中基因组选择利用区间；黑色区段为优化的两个品种的“SNP多态性热点区间”。优化的“SNP低频变异区间”为两个样本间的遗传相似区段(基因组相似区间)，即为育种基因组选择利用区间(简称育种利用区间)，也可作为育种热点区间；优化的“SNP多态性差异热点区间”为两个品种的基因组差异区间，可用于分析小麦品种的不同育种来源。
上述鉴定小麦育种中基因组选择利用区间和/或“SNP多态性热点区间”的方法的流程图如图1所示。
实施例2、鉴定小麦育种中基因组选择利用区间和/或“SNP多态性热点区间”的方法的验证
针对实施例1中得到的“SNP多态性热点区间”，使用软件设计的基因组特异引物进行扩增和扩增产物测序验证。
1.DNA提取
使用CTAB法分别提取供试小麦济麦22、良星99和中国春的DNA，其中中国春作为未使用实施例1中方法进行分析的对照。
2.PCR扩增及产物测序
选取实施例1中对小麦济麦22与良星99分析得到的位于基因组D组2号染色体(chr2D)上的优化两个品种间的“SNP多态性热点区间”(图4所示黑色区域)，利用基因组特异引物对济麦22、良星99以及中国春的DNA进行片段扩增，将所得PCR产物送由天一辉远生物科技有限公司进行Sanger测序验证所选DNA区段内核苷酸序列的差异(图5)。
2.1引物设计
扩增引物由在线软件Primer 3.0设计，其中引物长度范围为18-24bp，GC含量范围在40％-60％之间，退火温度范围为54-60℃，扩增引物扩增产物大小100-500bp。(引物序列：F：CAGGTCGAGACAGAGAACAA；R：ATCGAGCCCCTCAATTTCAT)
2.2PCR操作流程
PCR反应体系为20μL，包括10μL的2X M5 HiPer plus Taq HiFi PCR mix，正、反向引物(10μmol L–1)各1μL，150ngμL–1模板DNA 2μL，用ddH2O补至20μL。PCR扩增程序为95℃3min；95℃30s，56-57.4℃30-60s(依引物退火温度以及目标序列而定)，72℃2min，35个循环；72℃5min。
2.3PCR产物测序结果分析
PCR产物测序结果如图5所示，在济麦22(图5中Jimai 22所示)与良星99(图5中Liangxing 99所示)的“SNP多态性热点区间”中，两个品种在所选区段确实存在核苷酸序列的变异，表明实施例1中鉴定“SNP多态性热点区间”的方法是有效的。
上述验证结果表明实施例1中的方法可有效鉴定和区分小麦两个品种间的“SNP低频变异区间”(即基因组选择利用区间)和/或“SNP多态性热点区间”(即基因组差异区间)。
Claims (10)
1.一种鉴定作物育种材料间基因组选择利用区间和/或SNP多态性热点区间的方法，包括如下步骤：
A1)提取待鉴定作物材料基因组DNA并进行全基因组重测序获得待鉴定作物重测序数据；
A2)分别过滤名称为材料A和材料B的两个所述待鉴定作物材料重测序数据以对其进行质量控制，并将其分别比对到作物参考基因组上，得到名称为文件A的所述材料A基因组与所述作物参考基因组的DNA序列比对文件，和名称为文件B的所述材料B基因组与所述作物参考基因组的DNA序列比对文件；
A3)对所述文件A进行过滤，得到所述材料A与所述作物参考基因组相比存在的SNP变异位点和CNV变异区间，然后将所述CNV变异区间中存在的SNP变异位点去除，得到所述材料A基因组与所述作物参考基因组相比存在的排除CNV变异区间的SNP变异位点，将该SNP变异位点称为A-参考SNP变异位点；
和对所述文件B进行过滤，得到所述材料B与所述作物参考基因组相比存在的SNP变异位点和CNV变异区间，然后将所述CNV变异区间中存在的SNP变异位点去除，得到所述材料B基因组与所述作物参考基因组相比存在的排除CNV变异区间的SNP变异位点，将该SNP变异位点称为B-参考SNP变异位点；
A4)将所述A-参考SNP变异位点和B-参考SNP变异位点进行比较得到二者之间的差异SNP变异位点，计算所述差异SNP变异位点在所述待鉴定作物参考基因组上的密度分布；
A5)根据所述密度分布拟合得到正态分布后根据边界确定阈值，并使用所述阈值对正态分布进行分割，得到所述材料A和所述材料B之间的SNP多态性热点区间和SNP低频变异区间；所述材料A和所述材料B之间的SNP多态性热点区间即为所述材料A和所述材料B之间的基因组差异区间，所述材料A和所述材料B之间的SNP低频变异区间即为所述材料A和所述材料B之间的基因组相似区间或育种利用区间。
2.根据权利要求1所述的方法，其特征在于：所述方法还包括：
A6)对所述材料A和所述材料B之间的SNP多态性热点区间和SNP低频变异区间进行调整和优化，得到优化的所述材料A和所述材料B之间的SNP多态性热点区间和优化的所述材料A和所述材料B之间的SNP低频变异区间；优化的所述材料A和所述材料B之间的SNP多态性热点区间即为优化的所述材料A和所述材料B之间的基因组差异区间，优化的所述材料A和所述材料B之间的SNP低频变异区间即为优化的所述材料A和所述材料B之间的基因组相似区间或育种利用区间。
3.根据权利要求1所述的方法，其特征在于：所述A5)步骤中的所述阈值是以拟合得到的混合正态分布的交叉边界作为硬阈值而确定的；使用所述阈值对所述混合正态分布进行分割得到的密度较大的正态分布中的SNP位点所在区间为所述材料A和所述材料B之间的SNP多态性热点区间，密度较小的正态分布中的SNP位点所在区间为所述材料A和所述材料B之间的SNP低频变异区间。
4.根据权利要求1或2或3所述的方法，其特征在于：A4)步骤中计算所述差异SNP变异位点的密度分布，包括以下步骤：
A4-1)在0.5-5Mbp范围内选择以固定长度为单位将所述材料A和所述材料B全基因组进行划分小窗，
A4-2)对所述A-参考SNP变异位点和所述B-参考SNP变异位点在全基因组每个小窗中的SNP频率进行分析，
A4-3)统计小窗中所述差异SNP变异位点数，依据所述得到的差异SNP变异位点数与非丢失位点比例的比值得到所述差异SNP变异位点的密度分布。
5.根据权利要求1-4任一权利要求所述的方法，其特征在于：所述A5)步骤中的正态分布，是通过利用EM算法对所述A4)步骤中得到的差异SNP变异位点的密度分布进行拟合得到的。
6.根据权利要求1-5任一权利要求所述的方法，其特征在于：所述作物为普通小麦，所述阈值为32/1Mbp区间。
7.根据权利要求2-6任一权利要求所述的方法，其特征在于：所述A6)中对所述材料A和所述材料B间的SNP多态性热点区间和SNP低频变异区间进行调整和优化过程，是利用隐马尔可夫模型对所述材料A和所述材料B间的SNP多态性热点区间和SNP低频变异区间进行平滑校正处理。
8.根据权利要求1-7任一权利要求所述的方法，其特征在于：所述A3)步骤中的过滤，是使用GATK软件进行的。
9.权利要求1-8任一权利要求所述的方法的下述任意一种应用：
P1、权利要求1-8任一权利要求所述的方法在定位和鉴定作物育种热点选择区间中的应用；
P2、权利要求1-8任一权利要求所述的方法在鉴定作物品种基因组间关系中的应用；
P3、权利要求1-8任一权利要求所述的方法在开发作物新种质资源中的应用；
P4、权利要求1-8任一权利要求所述的方法在作物育种中的应用；
P5、权利要求1-8任一权利要求所述的方法在分析作物品种育种来源中的应用；
P6、权利要求1-8任一权利要求所述的方法在鉴定不同作物品种的相同来源育种区间中的应用；
P7、权利要求1-8任一权利要求所述的方法在鉴定作物的不同品种间的遗传差异中的应用；
P8、权利要求1-8任一权利要求所述的方法在鉴定作物育种过程中世代间基因组遗传信息的传递情况中的应用；
P9、权利要求1-8任一权利要求所述的方法在鉴定作物不同品种间的遗传区间利用关系中的应用。
10.根据权利要求1-8任一权利要求所述的方法和/或权利要求9所述的应用，其特征在于：所述作物为下述任一种：
D1)双子叶植物，
D2)单子叶植物，
D3)禾本目，
D4)禾本科，
D5)小麦属，
D6)普通小麦。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010742489.5A CN111798922B (zh) | 2020-07-29 | 2020-07-29 | 基于重测序数据中多态性位点密度鉴定小麦育种的基因组选择利用区间的方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202010742489.5A CN111798922B (zh) | 2020-07-29 | 2020-07-29 | 基于重测序数据中多态性位点密度鉴定小麦育种的基因组选择利用区间的方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN111798922A true CN111798922A (zh) | 2020-10-20 |
| CN111798922B CN111798922B (zh) | 2024-04-02 |
Family
ID=72828053
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202010742489.5A Active CN111798922B (zh) | 2020-07-29 | 2020-07-29 | 基于重测序数据中多态性位点密度鉴定小麦育种的基因组选择利用区间的方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN111798922B (zh) |
Cited By (1)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115831223A (zh) * | 2023-02-20 | 2023-03-21 | 吉林工商学院 | 一种挖掘近源物种间染色体结构变异的分析方法及系统 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20130110407A1 (en) * | 2011-09-16 | 2013-05-02 | Complete Genomics, Inc. | Determining variants in genome of a heterogeneous sample |
| CN106202995A (zh) * | 2016-07-13 | 2016-12-07 | 北京麦美瑞生物科技有限公司 | 小麦BSR‑Seq基因定位的方法 |
| CN110016497A (zh) * | 2018-01-09 | 2019-07-16 | 北京大学 | 一种检测肿瘤单细胞基因组拷贝数变异的方法 |
| CN111264858A (zh) * | 2019-10-30 | 2020-06-12 | 成都图径生物科技有限公司 | 海藻燕窝及其制备方法和用途 |
-
2020
- 2020-07-29 CN CN202010742489.5A patent/CN111798922B/zh active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20130110407A1 (en) * | 2011-09-16 | 2013-05-02 | Complete Genomics, Inc. | Determining variants in genome of a heterogeneous sample |
| CN106202995A (zh) * | 2016-07-13 | 2016-12-07 | 北京麦美瑞生物科技有限公司 | 小麦BSR‑Seq基因定位的方法 |
| CN110016497A (zh) * | 2018-01-09 | 2019-07-16 | 北京大学 | 一种检测肿瘤单细胞基因组拷贝数变异的方法 |
| CN111264858A (zh) * | 2019-10-30 | 2020-06-12 | 成都图径生物科技有限公司 | 海藻燕窝及其制备方法和用途 |
Cited By (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN115831223A (zh) * | 2023-02-20 | 2023-03-21 | 吉林工商学院 | 一种挖掘近源物种间染色体结构变异的分析方法及系统 |
| CN115831223B (zh) * | 2023-02-20 | 2023-06-13 | 吉林工商学院 | 一种挖掘近源物种间染色体结构变异的分析方法及系统 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN111798922B (zh) | 2024-04-02 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN109346130B (zh) | 一种直接从全基因组重测序数据中得到微单体型及其分型的方法 | |
| CN109196123B (zh) | 用于水稻基因分型的snp分子标记组合及其应用 | |
| CN109545278B (zh) | 一种鉴定植物lncRNA与基因互作的方法 | |
| CN111304303B (zh) | 微卫星不稳定的预测方法及其应用 | |
| US20210285063A1 (en) | Genome-wide maize snp array and use thereof | |
| CN108304694B (zh) | 基于二代测序数据分析基因突变的方法 | |
| US20190287646A1 (en) | Identifying copy number aberrations | |
| CN105046105B (zh) | 染色体跨度的单体型图及其构建方法 | |
| CN115678979A (zh) | 菠萝液相芯片及其应用 | |
| CN111485026A (zh) | 一种与绵羊出生重相关的snp位点、应用、分子标记和引物 | |
| CN111798922A (zh) | 基于重测序数据中多态性位点密度鉴定小麦育种的基因组选择利用区间的方法 | |
| Roy et al. | NGS-μsat: bioinformatics framework supporting high throughput microsatellite genotyping from next generation sequencing platforms | |
| CN102154452B (zh) | 一种鉴定顺式和反式调控作用的方法和系统 | |
| CN114420214A (zh) | 核酸测序数据的质量评估方法和筛选方法 | |
| CN111411167B (zh) | 烟草品种的dna指纹图谱库及其应用 | |
| JP2020517304A (ja) | Ｄｎａ分析のためのオフターゲット配列の使用 | |
| CN114300045A (zh) | 基于对照组的半监督snp分型方法、装置和电子设备 | |
| JP7166638B2 (ja) | 多型検出法 | |
| KR20220050296A (ko) | 배추 계통 구분을 위한 단일 염기 다형성 기반 마커 및 이의 용도 | |
| CN118109605B (zh) | 与尼里-拉菲水牛生长性状相关的snp分子标记组合及应用 | |
| CN111304358B (zh) | 基于冬瓜转录组序列开发的est-ssr引物及其应用 | |
| CN116590435B (zh) | 一种与猪背膘厚相关的因果候选基因及其鉴定方法和应用 | |
| CN117637020B (zh) | 一种基于深度学习的四倍体牡蛎全基因组snp分型方法 | |
| CN116855596A (zh) | 水稻品种同质性评价方法 | |
| CN114959068A (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 |