CN110211635A - 用于畜禽基因组选择分析的方法及畜禽育种方法 - Google Patents

用于畜禽基因组选择分析的方法及畜禽育种方法 Download PDF

Info

Publication number
CN110211635A
CN110211635A CN201910505483.3A CN201910505483A CN110211635A CN 110211635 A CN110211635 A CN 110211635A CN 201910505483 A CN201910505483 A CN 201910505483A CN 110211635 A CN110211635 A CN 110211635A
Authority
CN
China
Prior art keywords
data
livestock
group
pedigree
poultry
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
Application number
CN201910505483.3A
Other languages
English (en)
Other versions
CN110211635B (zh
Inventor
刘继强
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Beijing Compson Agricultural Technology Co Ltd
Original Assignee
Beijing Compson Agricultural Technology Co Ltd
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Beijing Compson Agricultural Technology Co Ltd filed Critical Beijing Compson Agricultural Technology Co Ltd
Priority to CN201910505483.3A priority Critical patent/CN110211635B/zh
Publication of CN110211635A publication Critical patent/CN110211635A/zh
Application granted granted Critical
Publication of CN110211635B publication Critical patent/CN110211635B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

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
    • G16B25/00—ICT specially adapted for hybridisation; ICT specially adapted for gene or protein expression
    • G16B25/10—Gene or protein expression profiling; Expression-ratio estimation or normalisation
    • 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

Landscapes

  • Health & Medical Sciences (AREA)
  • Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Medical Informatics (AREA)
  • Biophysics (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Genetics & Genomics (AREA)
  • Theoretical Computer Science (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • General Health & Medical Sciences (AREA)
  • Evolutionary Biology (AREA)
  • Biotechnology (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Epidemiology (AREA)
  • Software Systems (AREA)
  • Public Health (AREA)
  • Artificial Intelligence (AREA)
  • Evolutionary Computation (AREA)
  • Bioethics (AREA)
  • Databases & Information Systems (AREA)
  • Data Mining & Analysis (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Molecular Biology (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)

Abstract

本发明公开一种用于畜禽基因组选择分析的方法及畜禽育种方法。其中，用于畜禽基因组选择分析的方法包括选择具有表型数据、基因型数据和系谱数据的多个畜禽个体组成参考群；对于来自参考群的个体的表型数据、基因型数据和系谱数据分别进行处理；基于处理后的数据构建参考群模型，得到最优模型；利用最优模型预测候选群；和结果处理及输出。本发明的方法模型是基于系谱信息的传统BLUP法和基于SNP标记信息的GBLUP法的合并，它在模型形式上与BLUP及GBLUP法并无明显区别，但是本发明的方法通过特定的数据处理和清洗方法和原则，并通过构建特定模型大大提高了基因组选择的准确性。另外，本发明的结果输出更符合企业用户的需求，分析周期更短。

Description

用于畜禽基因组选择分析的方法及畜禽育种方法
技术领域
本发明涉及遗传育种领域，具体地涉及用于畜禽基因组选择分析的方法及畜禽育种的方法。
背景技术
生物的遗传基础来自于细胞核内的染色体，由脱氧核糖核酸(DNA)构成，在DNA的双螺旋结构上，不同碱基(A/T/G/C)的核苷酸以不同的排列顺序在DNA的某一段落上构成某个功能片段，即基因。不同的个体出现变异主要源于基因片段上核苷酸组合的变化。由多态、缺失、插入及交换几种变化方式组成，称之为“SNP多态性”。
基因芯片(gene chip)(又称DNA芯片、生物芯片)的原型是80年代中期提出的。基因芯片的测序原理是杂交测序方法，即通过与一组已知序列的核酸探针杂交进行核酸序列测定的方法，在一块基片表面固定了序列已知的靶核苷酸的探针。当溶液中带有荧光标记的核酸序列与基因芯片上对应位置的核酸探针产生互补匹配时，通过确定荧光强度最强的探针位置，获得一组序列完全互补的探针序列。据此可重组出靶核酸的序列。
在2002年，基因组选择(Genomic selection，GS)被首次提出，在很大程度上实现了标记辅助选择的优势。这种方法是利用覆盖全基因组的高密度分子标记进行标记辅助选择，可以追溯到大量影响不同数量性状的基因，从而实现对数量性状进行更准确的评定。
在动植物基因组中存在大量遗传标记(Single Nucleotide Polymorphism，SNP)，影响性状的所有基因都至少与一个标记紧密连锁，通过对所有标记效应的估计，实现对全基因组所有基因效应的估计。利用估计的标记效应计算个体育种值，即基因组育种值(Genomic estimated breeding values，GEBV)，然后根据GEBV的大小进行选择。
现有技术中已公开了一些常用的育种技术，其中包括BLUP法和GBLUP两类方法。关于BLUP法，例如，吴丽丽等在云南畜牧兽医(2004年第2期)中公开了一种先进的动物育种值估计方法-动物模型BLUP法，该方法对育种值估计有实用价值。再例如，CN 104106493A公开了一种黔北麻羊提纯复壮快速繁育方法。该方法通过组建基础群体，以BLUP(Best LinearUnbiased Prediction)技术为核心现代育种方法对基础群体选育2～3个世代，选育出整齐度好、生产性能稳定的优秀个体组成黔北麻羊核心群体，以人工授精技术为核心的现代繁殖技术对核心群体生产出来优秀种用个体进行人工授精扩散，迅速扩大优秀群体的数量；最终实现黔北麻羊的提纯复壮及快繁。
关于GBLUP方法，例如，CN 105512510 A公开了一种通过基因组数据对遗传力进行评估的方法，对于某一数量性状，通过使用不同数量的参考群个体利用GBLUP算法进行全基因组的标记效应的估计，进而得到估计群的育种值，并计算出估计准确度；通过基因组估计准确度与参考群体大小进行曲线直线化拟合，拟合出的回归方程的截距的倒数为遗传力的估计值；该方法通过基因组的数据对数量性状的遗传力进行评估，所研究的成果可直接应用于动植物数量性状育种中，本发明的算法不对个体进行系谱记录而是对个体基因组进行测序，通过全基因组标记来预测性状的遗传力，遗传力估计结果主要用于将来的育种工作中，另外，测序可以捕获到孟德尔抽样误差，相对记录系谱数据能够获得更准确的系谱信息。再例如，CN 107338321 A公开了一种确定最佳SNP数量及通过筛选标记对大黄鱼生产性能进行基因组选择育种的方法。先对参考群个体进行生产性能的表型测定和基因组测序，获得SNP位点；筛选出合格的SNP位点，并将缺失的基因型补齐；将参考群分为训练集和验证集进行杂交验证；通过单标记分析筛选与性状最显著关联的SNP位点，然后只使用这些位点通过GBLUP方法计算验证集个体的GEBV；进一步得到各个筛选SNP数量下的育种值估计准确度；最后确定SNP筛选的最佳数量。再根据该最佳数量，通过GBLUP方法计算出GEBV，进一步得到育种值估计准确度，根据该值的高低进行基因组选择育种。该方法可显著节省对大黄鱼生产性能的基因组选择费用。
上述现有技术均推进了基因组选择在育种领域的应用，发挥了基因组选择在动物育种领域的优势。然而，上述基因组选择的准确性仍需进一步提高。另外，目前的基因组选择分析中各个步骤都是分开的，缺少完整的育种解决方案，分析周期长，结果不稳定，影响因素大。这些方法只适用于科研院所。由于缺乏专业的育种技术人员，是无法独立将全套流程进行实操的，而且对于专业性很强的结果也不能很好地理解，这对将全基因组选择技术真正地应用到企业生产实际中造成了很大的阻碍。
发明内容
为解决现有技术中的至少部分技术问题，本发明提供一种准确性更高，且便于企业用户使用的畜禽基因组选择分析及畜禽育种方案。具体地，本发明包括以下内容。
本发明的第一方面，提供一种用于畜禽基因组选择分析的方法，其包括：
(1)选择具有表型数据、基因型数据和系谱数据的多个畜禽个体组成参考群；
(2)对于来自所述参考群的个体的表型数据进行处理，删除其中错误值和异常值，并选择数据呈正态分布的表型数据；
(3)对于来自所述参考群的个体的基因型数据进行处理，去除没有定位在染色体上的SNP位点、性染色体上的SNP位点、缺失率大于10％的SNP位点、最小等位基因频率小于1％的SNP位点，并且SNP基因分型检出率大于90％；
(4)对于来自所述参考群的个体的系谱数据进行处理，以使系谱根据出生日期排序，且至少包括三个世代的系谱数据，并排除重复ID的数据，将缺失数据定义为0；
(5)构建参考群模型，其包括下述子步骤：
(5-1)利用步骤(3)和(4)得到的处理数据构建下述H矩阵：
其中：
A11：A矩阵中非测序个体与非测序个体构成的系谱矩阵
A12：A矩阵中非测序个体与测序个体构成的系谱矩阵
A21：A矩阵中测序个体与非测序个体构成的系谱矩阵
A22：A矩阵中测序个体与测序个体构成的系谱矩阵
G：测序个体与测序个体构成的G矩阵
(5-2)进一步得到H逆矩阵：
其中，A-1为全部系谱关系的逆矩阵，G-1为基因组关系逆矩阵，为测序的个体系谱关系逆矩阵，
(5-3)进一步利用下式调整H逆矩阵的参数：
其中，w是A22的权重，其默认为0.05，τ和ω的默认值各自分别为1，a是变异中不能由基因组信息解释的百分比,默认是5％,即0.05，b＝1-a,默认是95％,即0.95；
将G矩阵矫正至A矩阵的尺度，调整G矩阵和A22的相对权重，并设置参数值τ和ω；
(5-4)确定固定因子和随机因子，并确定协变量及权重变量。本发明中，固定因子根据混合线性模型中Wald显著性检验,如果不显著,则去掉该固定因子,只保留Wald检验达到显著性的固定因子，随机因子使用似然比(LRT)检验,如果不显著,则去掉该随机因子,只保留似然比(LRT)检验达到显著性的随机因子；
(5-5)将所述参考群中的数据随机的分为验证样本和训练样本，用以验证模型，从而确定最优模型；
(6)利用所述最优模型预测候选群；和
(7)结果处理及输出。
优选地，所述固定因子包括场年季，所述随机因子包括加性效应、母体效应和永久环境效应。
优选地，在所述步骤(5-5)中，包括将参考群中的数据分为n份，取其中的第一份去除表型数据后作为验证样本，其余的n-1份全部作为训练样本，利用训练样本预测验证样本的育种值，并计算所述育种值与TBV(真实育种值)的相关系数和准确性，并评估计算预测的可靠性，可靠性的计算方法可根据评估的育种值,与真实育种值做线性回归分析,所得到的回归系数即为可靠性；
取n份数据中的第二份去除表型数据后作为验证样本，其余的n-1份全部作为训练样本，重复进行n次验证。
优选地，所述步骤(6)包括对基因型数据的处理：去除没有定位在染色体上的SNP位点、性染色体上的SNP位点、缺失率大于10％的SNP位点、最小等位基因频率小于1％的SNP位点，同时保证SNP基因分型检出率大于90％，由此得到处理数据；和使用Beagle对所述处理数据进行填充，得到待分析的基因型数据。
优选地，所述步骤(6)还包括对系谱数据的处理：根据出生日期排序以使祖代在前，子代在后，系谱数据至少包括三个世代，并且去除ID重复和父母本间交叉重复的数据，同时将缺失数据定义0，由此得到待分析的系谱数据。
优选地，所述步骤(6)进一步包括运行下述混合线性方程组计算固定因子效应值和随机因子效应值，并根据方差组分(例如使用限制性极大似然法(REML)进行评估得到的方差组分)计算单性状遗传力：
其中，Y为待分析的性状，X为固定因子，Z为随机因子。
优选地，所述步骤(6)进一步包括运行下述混合线性方程组计算固定因子效应值和随机因子效应值，并根据方差组分计算多性状遗传力和遗传相关性：
其中，Y为待分析的性状，X为固定因子，Z为随机因子。
优选地，步骤(7)包括表型数据的统计及展示，单性状模型和多性状模型各自中固定因子和随时因子的设定，各性状交叉验证提高准确性和无偏性，和将各性状的单性状和多性状模型育种值排名。
优选地，在步骤(7)中以图表的形式展示表型数据的统计结果。
本发明的第二方面，提供一种畜禽育种的方法，其使用第一方面所述的方法进行畜禽基因组选择分析，或者使用所述畜禽基因组选择分析的结果进行选择育种。
本发明的方法模型是基于系谱信息的传统BLUP法和基于SNP标记信息的GBLUP法的合并，它在模型形式上与BLUP及GBLUP法并无明显区别，但是本发明的方法通过特定的数据处理和清洗方法和原则，并通过构建特定模型大大提高了基因组选择的准确性。
另外，本发明的方法针对于目前畜禽基因组选择分析技术存在各个实现手段都是分开的，没有一个完整的解决方案，针对于企业用户，现有方案不具有实践性，时间周期太长，无法有效解决企业用户关注的问题的缺点，而创建了一步法基因组选择的流程，结果输出更符合企业用户的需求，分析周期更短，解决了企业因为缺乏专业育种技术人员而无法开展育种工作的问题，帮助企业应用基因组选择技术加快育种进程。同时也可以让企业能够及时根据结果去筛选出优良的个体进行选育选配工作，从而达到及时有效解决企业用户关注的问题的目的。
附图说明
图1参考群构建的流程图。
图2交叉验证的示意图。
图3预测候选群的流程图。
图4以图像展示的表型数据。
具体实施方式
现详细说明本发明的多种示例性实施方式，该详细说明不应认为是对本发明的限制，而应理解为是对本发明的某些方面、特性和实施方案的更详细的描述。
应理解本发明中所述的术语仅仅是为描述特别的实施方式，并非用于限制本发明。另外，对于本发明中的数值范围，应理解为具体公开了该范围的上限和下限以及它们之间的每个中间值。在任何陈述值或陈述范围内的中间值以及任何其他陈述值或在所述范围内的中间值之间的每个较小的范围也包括在本发明内。这些较小范围的上限和下限可独立地包括或排除在范围内。
除非另有说明，否则本文使用的所有技术和科学术语具有本发明所述领域的常规技术人员通常理解的相同含义。虽然本发明仅描述了优选的方法和材料，但是在本发明的实施或测试中也可以使用与本文所述相似或等同的任何方法和材料。本说明书中提到的所有文献通过引用并入，用以公开和描述与所述文献相关的方法和/或材料。在与任何并入的文献冲突时，以本说明书的内容为准。除非另有说明，否则“％”为基于重量的百分数。
本发明的“畜禽”是指可供发展畜牧业的牲畜、家禽。其中牲畜的实例包括但不限于猪、牛、羊和兔等哺乳动物。家禽的实例包括但不限于鸡、鸭、鹅和鸽子等。本发明的畜禽优选为猪和牛。
本发明的基因组选择分析不仅包括利用覆盖全基因组的高密度遗传标记(SNP)计算个体的基因组估计育种值(GEBV)，而且还包括基于系谱数据或信息估计育种值(estimated breeding value,EBV)。与仅基于系谱数据或信息的方法以及仅基于SNP数据的方法相比，本发明的准确性更高，从而为早期选择提供了可能。
[用于畜禽基因组选择分析的方法]
本发明的第一方面，提供一种用于畜禽基因组选择分析的方法，其至少包括以下七个步骤。
(1)选择具有表型数据、基因型数据和系谱数据的多个畜禽个体组成参考群；
(2)对于来自所述参考群的个体的表型数据进行处理，删除其中错误值和异常值，并选择数据呈正态分布的表型数据；
(3)对于来自所述参考群的个体的基因型数据进行处理，去除没有定位在染色体上的SNP位点、性染色体上的SNP位点、缺失率大于10％的SNP位点、最小等位基因频率小于1％的SNP位点，并且SNP基因分型检出率大于90％；
(4)对于来自所述参考群的个体的系谱数据进行处理，以使系谱根据出生日期排序，且至少包括三个世代的系谱数据，并排除重复ID的数据，将缺失数据定义为0；
(5)构建参考群模型；
(6)利用所述最优模型预测候选群；和
(7)结果处理及输出。
下面详细说明各步骤。
步骤(1)
本发明的步骤(1)为选择具有表型数据、基因型数据和系谱数据的多个畜禽个体组成参考群的步骤。
本发明的参考群中每个畜禽个体均具有表型数据、基因型数据和系谱数据，不包括只具有表型数据、基因型数据或系谱数据，或者仅具有上述两者的个体。
本发明的参考群中畜禽个体的数量不特别限定，但由于畜禽个体的数量越多，则所得到的模型的可靠性和准确性越高。因此，本发明中畜禽个体的数量一般为500以上，优选1000以上，更优选2000以上。
步骤(2)
本发明的步骤(2)为对于来自参考群的个体的表型数据进行处理以提高分析准确性。表型数据的处理包括删除错误值和异常值，并选择呈正态分布的表型数据。表型数据是否属于错误值和异常值可由已知方法而容易地确认。本发明中将三倍标准差以外的数据作为错误值和异常值予以删除。表型数据是否符合正态分布可通过已知方法确认，例如可以使用条形图、箱线图和QQ图查看数据分布。对于不满足正态分布的数据还可包括数据转化步骤，例如通过log转化、Box-Cox转化以使这些数据符合正态分布。
步骤(3)
本发明的步骤(3)为对于来自参考群的个体的基因型数据进行处理，去除没有定位在染色体上的SNP位点、性染色体上的SNP位点、缺失率大于10％的SNP位点、最小等位基因频率小于1％的SNP位点，并且SNP基因分型检出率大于90％。经过上述处理，特别是去除性染色体上的SNP位点，发明人发现可以大大提高预测结果的可靠性。
在某些实施方案中，进一步包括对于来自参考群的个体的基因型数据进行处理后的处理数据进行填充，得到待分析的基因型数据。可使用已知的方法进行填充。例如使用Beagle进行填充，其示例性流程包括：将Plink筛选后的基因型数据转化为VCF文件；使用Beagle软件对VCF文件进行填充；将结果转化为Plink格式。
步骤(4)
本发明的步骤(4)为对于来自参考群的个体的系谱数据进行处理，以使系谱根据出生日期排序，并排除重复ID的数据，将缺失数据定义为0。系谱数据至少包括3个世代的数据，优选包括5代以上的数据，更优选包括7代以上的数据。
步骤(5)
本发明的步骤(5)为构建参考群模型的步骤。一般情况下，构建参考群模型具有多个子步骤，优选地子步骤如下：
(5-1)利用步骤(3)和(4)得到的处理数据构建H矩阵：
其中，A11：A矩阵中非测序个体与非测序个体构成的系谱矩阵
A12：A矩阵中非测序个体与测序个体构成的系谱矩阵
A21：A矩阵中测序个体与非测序个体构成的系谱矩阵
A22：A矩阵中测序个体与测序个体构成的系谱矩阵
G：测序个体与测序个体构成的G矩阵。
(5-2)由H矩阵进一步得到H逆矩阵：
其中，A-1为全部系谱关系的逆矩阵，G-1为基因组关系逆矩阵，为测序的个体系谱关系逆矩阵。可使用已知方法来构建A矩阵和G矩阵。对A矩阵和G矩阵求逆即得到A-1和G-1。
(5-3)进一步利用下式调整H逆矩阵的参数：
其中，w是A22的权重，其默认为0.05，τ和ω的默认值各自分别为1，
将G矩阵矫正至A矩阵的尺度，调整G矩阵和A22的相对权重，并设置参数值τ和ω；
(5-4)确定固定因子和随机因子，并确定协变量及权重变量。固定因子包括场年季；随机因子包括加性效应、母体效应和永久环境效应。
(5-5)将参考群中的数据随机地分为验证样本和训练样本，用以验证模型，从而确定最优模型。最优模型的确定一般包括将参考群中的数据分为n份，取其中的第一份去除表型数据后作为验证样本，其余的n-1份全部作为训练样本，利用训练样本预测验证样本的育种值，并计算育种值与TBV的相关系数和准确性，并评估计算预测的可靠性。取n份数据中的第二份去除表型数据后作为验证样本，其余的n-1份全部作为训练样本，重复进行n次验证。其中，n为5-100之间的自然数，优选10-50次。
在某些实施方案中，最优模型的确定包括：将原始数据分为五份；每一次取一份作为验证样本(去除表型值)，其余作为训练样本；使用训练群体预测测试群体的育种值(GEBV)；计算GEBV和TBV的相关系数、准确性，并计算可靠性；共重复五次以上。
步骤(6)
本发明的步骤(6)为利用最优模型预测候选群。在利用已知的基因型数据和系谱数据预测之前，优选包括对基因型数据和系谱数据进行处理。其中，对于基因型数据的处理可与步骤(3)中的处理类似，其可包括去除没有定位在染色体上的SNP位点、性染色体上的SNP位点、缺失率大于10％的SNP位点、最小等位基因频率小于1％的SNP位点，同时保证SNP基因分型检出率大于90％，由此得到处理数据。需要说明的是，本发明中对于基因型数据的处理必去除性染色体上的SNP位点，这对于本发明准确性的提高而言是重要的。另外，本发明中对于基因型数据的处理还包括对数据的填充。填充方法可使用已知的方法。例如使用Beagle对处理后数据进行填充，得到待分析的基因型数据。
本发明的步骤(6)还包括对系谱数据的处理。其中，对于系谱数据的处理可与步骤(4)相同。可包括根据出生日期排序以使祖代在前，子代在后，系谱数据至少包括三个世代，并且去除ID重复和父母本间交叉重复的数据，同时将缺失数据定义0，由此得到待分析的系谱数据。
本发明的步骤(6)中，利用最优模型预测候选群时包括运行单性状模型，即运行下述混合线性方程组计算固定因子效应值和随机因子效应值，并根据方差组分计算单性状遗传力：
其中，Y为待分析的性状，X为固定因子，Z为随机因子。
本发明的步骤(6)中，利用最优模型预测候选群时还包括运行多性状模型，即运行下述混合线性方程组计算固定因子效应值和随机因子效应值，并根据方差组分计算多性状遗传力和遗传相关性：
其中，Y为待分析的性状，X为固定因子，Z为随机因子。
本发明的步骤(6)中，可使用已知的方法进行单性状模型和多性状模型的运算。例如，使用ASREML软件。
步骤(7)
本发明的步骤(7)为结果处理及输出。步骤(7)包括表型数据的统计及展示。表型数据的统计包括统计各个性状的总数、平均值、偏差、SD和CV值等。优选地，以图表的形式展示表型数据的统计结果。步骤(7)的结果输出还可包括以报表的形式设定单性状模型和多性状模型各自中固定因子和随时因子。步骤(7)的结果输出还包括交叉验证各性状，以提高准确性和无偏性。步骤(7)的结果输出还包括将各性状的单性状和多性状模型育种值排名。包括SSBLUP、标准误、可靠性、性状育种值和遗传力，方差组分。
交叉验证过程:
将原始数据分为五份
每一次取一份作为验证样本(去除表型值)，其余作为训练样本
使用训练群体预测测试群体的育种值(GEBV)
计算GEBV和TBV的相关系数：准确性
可靠性的计算
共重复五次
[畜禽育种的方法]
本发明的第二方面，提供一种畜禽育种的方法，其包括使用第一方面所述的方法进行畜禽基因组选择分析，或者使用该畜禽基因组选择分析的结果进行选择育种。
实施例
本实施例用于示例性说明本发明的方法。具体地，本发明包括：
S1，构建参考群，确定模型；
S2，利用构建好的参考群进行候选群表型值的预测，然后根据表型值进行选择，达到提前选择提升育种效率的目的。
下面详细说明各步骤。
一、S1参考群构建流程：
参考群的构建流程参见图1。
S1.1芯片数据清洗
芯片数据下机后，是plink格式的数据，name.map和name.ped；
S1.1.1使用Plink进行基因组数据的清洗
清洗标准：
去除没有定位在染色体上的SNP位点和性染色体上的SNP位点；
去除缺失率大于10％的SNP位点；
SNP基因分型检出率需要大于90％(call rate>90％)；
去除最小等位基因频率(minor allele frequency,MAF)小于1％的SNP位点。
S1.1.2使用Beagle进行基因组数据的填充
使用Beagle填充流程:
将Plink筛选后的基因型数据转化为VCF文件；
使用Beagle软件对VCF文件进行填充；
将结果转化为Plink格式。
S1.2表型数据清洗
数据删除明显的错误值和异常值；
分析的性状符合正态分布，可以使用条形图、箱线图和QQ图查看数据分布；
对于不满足正态分布的数据，进行数据转化(log转化，Box-Cox转化等)，以保证数据要满足正态分布，或者删除三倍标准差以外的数据。
S1.3系谱数据清洗
系谱应当根据出生日期排序(祖代在前，子代在后)；
系谱数据至少包括三个世代；
数据清洗(ID不能重复，父母本间没有交叉重复)，缺失数据定义：“0”。
S1.4参考群构建模型
S1.4.1根据系谱和基因型数据构建H逆矩阵
由下式构建H矩阵：
构建H逆矩阵：
其中，A-1为全部系谱关系的逆矩阵，G-1为基因组关系逆矩阵，为测序的个体系谱关系逆矩阵。可使用已知方法来构建A矩阵和G矩阵。对A矩阵和G矩阵求逆即得到A-1和G-1。
对清洗后的系谱数据和芯片数据构建H逆矩阵，使用R语言编程。
S1.4.2调整H逆矩阵不同参数
H逆矩阵参数设置汇总
其中，w是A22的权重，其默认为0.05，τ和ω的默认值各自分别为1，a、b的定义如前所述。
将G矩阵矫正到A矩阵的尺度；
调整G矩阵和A22的相对权重；
设置参数值τ和ω。
S1.4.3考虑不同固定因子和随机因子：
确定固定因子(场年季等)；
确定随机因子(加性效应，母体效应，永久环境效应等)；
确定协变量及权重变量。
S1.4.4交叉验证，选择最优模型：如图2所示，进行交叉验证。其中，每行中浅色部分为验证样本，其余四个为训练样本。
将原始数据分为五份；
每一次取一份作为验证样本(去除表型值)，其余作为训练样本；
使用训练群体预测测试群体的育种值(GEBV)；
计算GEBV和TBV的相关系数：准确性；
可靠性的计算；
共重复五次。
二、S2预测候选群流程：
S2.1芯片数据清洗：
芯片数据下机后，是plink格式的数据，name.map和name.ped。
S2.1.1使用Plink进行基因组数据的清洗：
清洗标准：
去除没有定位在染色体上的SNP位点和性染色体上的SNP位点；
去除缺失率大于10％的SNP位点；
SNP基因分型检出率需要大于90％(call rate>90％)；
去除最小等位基因频率(minor allele frequency,MAF)小于1％的SNP位点。
S2.1.2使用Beagle进行基因组数据的填充
用Beagle填充流程:
将Plink筛选后的基因型数据转化为VCF文件；
使用Beagle软件对VCF文件进行填充；
将结果转化为Plink格式。
S2.2系谱数据清洗
系谱应当根据出生日期排序(祖代在前，子代在后)；
系谱数据至少包括三个世代；
数据清洗(ID不能重复，父母本间没有交叉重复)，缺失数据定义：“0”。
S2.3候选群基因组选择
S2.3.1运行单性状模型
单性状混合线性方程组如下：
单性状模型，软件使用ASREML，处理流程如下：
确定需要分析的性状，比如产仔数(Y)
确定固定因子，比如场年季(X)
确定随机因子，比如加性效应，窝别效应(Z)
使用Wald-F test计算固定因子的显著性
使用REML方法估计随机因子的方差组分及标准误
求解混合线性方程组，计算固定因子效应值BLUE
求解混合线性方程组，计算随机因子效应值BLUP
根据方差组分，计算性状遗传力
S2.3.2运行多性状模型
多性状混合线性方程组如下：
多性状模型，软件使用ASREML，处理流程如下：
确定需要分析的性状，比如产仔数，背膘厚(Y)
确定固定因子，比如场年季(X)
确定随机因子，比如加性效应，窝别效应(Z)
使用Wald-F test计算固定因子的显著性
使用REML方法估计随机因子的方差组分及标准误
求解混合线性方程组，计算固定因子效应值BLUE
求解混合线性方程组，计算随机因子效应值BLUP
根据方差组分，计算性状遗传力和遗传相关性。
S2.3.3输出结果报表
S2.3.3.1GS分析报表1-表型数据汇总统计
将表型数据汇总统计如下表1所示：
表1
ID Total_number Missing_number Mean Variance SD CV
y1 50 0 28.507047 76.894049 8.7689252 30.760553
y2 50 0 28.507047 76.894049 8.7689252 30.760553
y3 50 0 28.162378 112.04075 10.58493 37.585357
y4 50 0 31.761742 146.5482 12.105709 38.114123
y5 50 0 30.763959 129.05688 11.36032 36.927368
表型数据以图像展。如图4所示。
S2.3.3.2 GS分析报表2-模型设定
包括单性状模型和多性状模型的固定因子和随机因子设定：
单性状模型固定因子的设定：
固定因子：F1+F2+F3
随机因子：加性效应
其中，y1、y3、y4、y5、y6固定因子为三个，y2无固定因子。
多性状模型：
多性状组合：y1、y2、y3、y4、y5、y6
固定因子：F1+F2+F3
随机因子：加性效应
其中，y1、y2、y3、y4、y5、y68固定因子为三个。
模型中H矩阵的各个参数设置：
H矩阵构建时，可以进行不同参数的设置，进行模型的调整，可以设置的参数：
a和b是将G矩阵和A矩阵进行尺度的调整，它是根据两个方程进行的计算，默认参数。
d是G矩阵不能解释的方差变异，这部分由A22解释，默认为0.05；
τ和ω是可以调整的参数，其默认值各自分别为1。
S2.3.3.3 GS分析报表3-交叉验证：
包括各个性状交叉验证的准确性和无偏性，相对于系谱动物模型预测准确性的提升。
性状 类型 准确性 标准误 可靠性 标准误
y1 ablup 0.133709758 0.0080841 0.285908 0.0302332
y1 hblup 0.530340341 0.0100877 0.6931378 0.0225313
y2 ablup 0.209960298 0.0138328 0.5131189 0.0645642
y2 hblup 0.745133171 0.0095429 0.9384778 0.024603
y3 ablup 0.386504166 0.010035 0.8162168 0.0277909
y3 hblup 0.649089088 0.0097295 0.8746612 0.0297613
y4 ablup 0.255854127 0.0161493 0.5413297 0.0425037
y4 hblup 0.598049463 0.0073989 0.8121373 0.0339824
y5 ablup 0.336049113 0.0064751 0.714417 0.0261479
y5 hblup 0.631972658 0.0080424 0.8948797 0.0298291
y6 ablup 0.109601697 0.0100453 0.2704414 0.0256011
y6 hblup 0.388222768 0.0164567 0.611726 0.033442
S2.3.3.4 GS分析报表4-结果文件：
包括各个性状的单性状和多性状模型育种值排名，包括SSBLUP、标准误，可靠性，性状育种值和遗传力，方差组分。
ID HBLUP se Rank
ID669 13.7 0.7997 1
ID833 11.28 0.7571 2
ID823 10.77 0.7914 3
ID198 10.15 0.7799 4
ID179 10.02 0.7907 5
ID665 9.274 0.7935 6
ID42 9.207 2.616 7
ID667 8.884 0.798 8
ID858 8.846 0.792 9
ID668 8.84 0.7979 10
ID666 8.688 0.7988 11
ID664 8.673 0.7944 12
ID682 8.502 0.7917 13
ID678 8.433 0.7976 14
ID681 8.384 0.7895 15
在不背离本发明的范围或精神的情况下，可对本发明说明书的具体实施方式做多种改进和变化，这对本领域技术人员而言是显而易见的。由本发明的说明书得到的其他实施方式对技术人员而言是显而易见得的。本申请说明书和实施例仅是示例性的。

Claims (10)

1.一种用于畜禽基因组选择分析的方法，其特征在于，包括以下步骤：
(1)选择具有表型数据、基因型数据和系谱数据的多个畜禽个体组成参考群；
(2)对于来自所述参考群的个体的表型数据进行处理，删除其中错误值和异常值，并选择呈正态分布的表型数据；
(3)对于来自所述参考群的基因型数据进行处理，去除没有定位在染色体上的SNP位点、性染色体上的SNP位点、缺失率大于10％的SNP位点、最小等位基因频率小于1％的SNP位点，并且使SNP基因分型检出率大于90％；
(4)对于来自所述参考群的系谱数据进行处理，以使系谱根据出生日期排序，且至少包括三个世代的系谱数据，排除重复ID的数据，将缺失数据定义为0；
(5)构建参考群模型，其包括下述子步骤：
(5-1)利用步骤(3)和(4)得到的处理数据构建H矩阵：
(5-2)进一步得到H逆矩阵：
其中，A-1为全部系谱关系的逆矩阵，G-1为基因组关系逆矩阵，为测序的个体系谱关系逆矩阵，
(5-3)进一步利用下式调整H逆矩阵的参数：
其中，w是A22的权重，其默认为0.05，τ和ω的默认值各自分别为1，a是变异中不能由基因组信息解释的百分比,默认是5％,即0.05，b＝1-a,默认是95％,即0.95；
将G矩阵矫正至A矩阵的尺度，调整G矩阵和A22的相对权重，并设置参数值τ和ω；
(5-4)确定固定因子和随机因子，并确定协变量及权重变量；
(5-5)将所述参考群中的数据随机的分为验证样本和训练样本，用以验证模型，从而确定最优模型；
(6)利用所述最优模型预测候选群；
(7)结果处理及输出。
2.根据权利要求1所述的用于畜禽基因组选择分析的方法，其特征在于，所述固定因子包括场年季，所述随机因子包括加性效应、母体效应和永久环境效应。
3.根据权利要求2所述的用于畜禽基因组选择分析的方法，其特征在于，所述步骤(5-5)包括将参考群中的数据分为n份，取其中的第一份去除表型数据后作为验证样本，其余的n-1份全部作为训练样本，利用训练样本预测验证样本的育种值，并计算所述育种值与TBV的相关系数和准确性，并评估计算预测的可靠性；
取n份数据中的第二份去除表型数据后作为验证样本，其余的n-1份全部作为训练样本，重复进行n次验证。
4.根据权利要求3所述的用于畜禽基因组选择分析的方法，其特征在于，所述步骤(6)包括对基因型数据的处理：去除没有定位在染色体上的SNP位点、性染色体上的SNP位点、缺失率大于10％的SNP位点、最小等位基因频率小于1％的SNP位点，同时保证SNP基因分型检出率大于90％，由此得到处理数据；和使用Beagle对所述处理数据进行填充，得到待分析的基因型数据。
5.根据权利要求4所述的用于畜禽基因组选择分析的方法，其特征在于，所述步骤(6)还包括对系谱数据的处理：根据出生日期排序以使祖代在前，子代在后，系谱数据至少包括三个世代，并且去除ID重复和父母本间交叉重复的数据，同时将缺失数据定义0，由此得到待分析的系谱数据。
6.根据权利要求5所述的用于畜禽基因组选择分析的方法，其特征在于，所述步骤(6)进一步包括运行下述混合线性方程组计算固定因子效应值和随机因子效应值，并根据方差组分计算单性状遗传力：
其中，Y为待分析的性状，X为固定因子，Z为随机因子。
7.根据权利要求6所述的用于畜禽基因组选择分析的方法，其特征在于，所述步骤(6)进一步包括运行下述混合线性方程组计算固定因子效应值和随机因子效应值，并根据方差组分计算多性状遗传力和遗传相关性：
其中，Y为待分析的性状，X为固定因子，Z为随机因子。
8.根据权利要求1所述的用于畜禽基因组选择分析的方法，其特征在于，步骤(7)包括表型数据的统计及展示，单性状模型和多性状模型各自中固定因子和随时因子的设定，各性状交叉验证提高准确性和无偏性，和将各性状的单性状和多性状模型育种值排名。
9.根据权利要求8所述的用于畜禽基因组选择分析的方法，其特征在于，在步骤(7)中以图表的形式展示表型数据的统计结果。
10.一种畜禽育种的方法，其特征在于，使用根据权利要求1-9任一项所述的方法进行畜禽基因组选择分析，或者使用所述畜禽基因组选择分析的结果进行选择育种。
CN201910505483.3A 2019-06-12 2019-06-12 用于畜禽基因组选择分析的方法及畜禽育种方法 Active CN110211635B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910505483.3A CN110211635B (zh) 2019-06-12 2019-06-12 用于畜禽基因组选择分析的方法及畜禽育种方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910505483.3A CN110211635B (zh) 2019-06-12 2019-06-12 用于畜禽基因组选择分析的方法及畜禽育种方法

Publications (2)

Publication Number Publication Date
CN110211635A true CN110211635A (zh) 2019-09-06
CN110211635B CN110211635B (zh) 2020-07-21

Family

ID=67792213

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910505483.3A Active CN110211635B (zh) 2019-06-12 2019-06-12 用于畜禽基因组选择分析的方法及畜禽育种方法

Country Status (1)

Country Link
CN (1) CN110211635B (zh)

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110610744A (zh) * 2019-09-11 2019-12-24 华中农业大学 一种高效可并行运算且高准确性的全基因组选择方法
CN110853711A (zh) * 2019-11-20 2020-02-28 云南省烟草农业科学研究院 一种预测烟草果糖含量的全基因组选择模型及其应用
CN111797544A (zh) * 2020-07-24 2020-10-20 新疆七色花信息科技有限公司 基于blup算法的遗传育种信息处理系统及处理方法
CN112273291A (zh) * 2020-10-28 2021-01-29 厦门大学 基于全基因组选择的大黄鱼抗刺激隐核虫病的选育方法
CN113517027A (zh) * 2020-04-09 2021-10-19 杭州锘崴信息科技有限公司 基于隐私保护并实现全基因组关联分析的联盟学习系统及方法
CN113951169A (zh) * 2021-12-16 2022-01-21 山东新希望六和集团有限公司 生长性能测定模型的训练方法、测定方法及装置
CN116072226A (zh) * 2023-01-17 2023-05-05 中国农业大学 一种用于蛋鸡产蛋性状基因组选择的机器学习方法及系统
CN116064846A (zh) * 2023-01-30 2023-05-05 中国海洋大学三亚海洋研究院 一种评估花鲈生长和抗性性状综合育种值的方法及应用

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20150156996A1 (en) * 2013-08-27 2015-06-11 Recombinetics, Inc. Livestock with genetically modified prolactin receptor
CN107338321A (zh) * 2017-08-29 2017-11-10 集美大学 一种确定最佳snp数量及其通过筛选标记对大黄鱼生产性能进行基因组选择育种的方法
CN107563147A (zh) * 2017-08-02 2018-01-09 中国农业大学 一种估计基因组育种值的方法及装置

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20150156996A1 (en) * 2013-08-27 2015-06-11 Recombinetics, Inc. Livestock with genetically modified prolactin receptor
CN107563147A (zh) * 2017-08-02 2018-01-09 中国农业大学 一种估计基因组育种值的方法及装置
CN107338321A (zh) * 2017-08-29 2017-11-10 集美大学 一种确定最佳snp数量及其通过筛选标记对大黄鱼生产性能进行基因组选择育种的方法

Non-Patent Citations (10)

* Cited by examiner, † Cited by third party
Title
刘瑞: ""非正态和多元过程能力分析及其在流程和离散制造业中的应用"", 《中国优秀硕士学位论文全文数据库 工程科技Ⅱ辑》 *
张璟等: ""地方政府"掠夺之手"、金融发展对地区收入差距的影响——来自中国省际面板数据的经验证据"", 《软科学》 *
朱才业: ""不同尾型绵羊全基因组关联分析、拷贝数变异及选择信号检测"", 《中国博士学位论文全文数据库 农业科技辑》 *
朱波: ""一步法和多性状基因组选择在西门塔尔牛群体中的应用研究"", 《中国博士学位论文全文数据库 农业科技辑》 *
李景密: ""跟驰状态下城市主干道交通流特征研究"", 《中国优秀硕士学位论文全文数据库 工程科技Ⅱ辑》 *
杨泽明: ""猪育种信息管理与分析系统软件的研究与开发"", 《中国优秀硕博士学位论文全文数据库（博士） 农业科技辑》 *
白皓: ""鸡喙畸形性状全基因组关联分析及拷贝数变异研究"", 《中国博士学位论文全文数据库 农业科技辑》 *
简银巧: ""热带玉米全长泛转录组和基因组大小变异及应用"", 《中国博士学位论文全文数据库 农业科技辑》 *
赵德胜等: ""家畜全基因组关联分析研究进展"", 《中国畜牧兽医》 *
赵海燕: ""基于高密度SNP芯片的深县猪体型性状全基因组关联分析"", 《中国优秀硕士学位论文全文数据库 农业科技辑》 *

Cited By (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN110610744A (zh) * 2019-09-11 2019-12-24 华中农业大学 一种高效可并行运算且高准确性的全基因组选择方法
CN110853711B (zh) * 2019-11-20 2023-09-12 云南省烟草农业科学研究院 一种预测烟草果糖含量的全基因组选择模型及其应用
CN110853711A (zh) * 2019-11-20 2020-02-28 云南省烟草农业科学研究院 一种预测烟草果糖含量的全基因组选择模型及其应用
CN113517027A (zh) * 2020-04-09 2021-10-19 杭州锘崴信息科技有限公司 基于隐私保护并实现全基因组关联分析的联盟学习系统及方法
CN113517027B (zh) * 2020-04-09 2024-05-24 杭州锘崴信息科技有限公司 基于隐私保护并实现全基因组关联分析的联盟学习系统及方法
CN111797544A (zh) * 2020-07-24 2020-10-20 新疆七色花信息科技有限公司 基于blup算法的遗传育种信息处理系统及处理方法
CN111797544B (zh) * 2020-07-24 2024-05-31 新疆七色花信息科技有限公司 基于blup算法的遗传育种信息处理系统及处理方法
CN112273291A (zh) * 2020-10-28 2021-01-29 厦门大学 基于全基因组选择的大黄鱼抗刺激隐核虫病的选育方法
CN112273291B (zh) * 2020-10-28 2021-09-07 厦门大学 基于全基因组选择的大黄鱼抗刺激隐核虫病的选育方法
CN113951169A (zh) * 2021-12-16 2022-01-21 山东新希望六和集团有限公司 生长性能测定模型的训练方法、测定方法及装置
CN113951169B (zh) * 2021-12-16 2022-04-22 山东新希望六和集团有限公司 生长性能测定模型的训练方法、测定方法及装置
CN116072226A (zh) * 2023-01-17 2023-05-05 中国农业大学 一种用于蛋鸡产蛋性状基因组选择的机器学习方法及系统
CN116064846A (zh) * 2023-01-30 2023-05-05 中国海洋大学三亚海洋研究院 一种评估花鲈生长和抗性性状综合育种值的方法及应用

Also Published As

Publication number Publication date
CN110211635B (zh) 2020-07-21

Similar Documents

Publication Publication Date Title
CN110211635A (zh) 用于畜禽基因组选择分析的方法及畜禽育种方法
VanRaden et al. Selecting sequence variants to improve genomic predictions for dairy cattle
CN106779076A (zh) 基于生物信息的选育良种系统及其算法
Ma et al. Functional mapping of quantitative trait loci underlying the character process: a theoretical framework
García-Gámez et al. Linkage disequilibrium and inbreeding estimation in Spanish Churra sheep
CN107278877B (zh) 一种玉米出籽率的全基因组选择育种方法
AU2011261447B2 (en) Methods and compositions for predicting unobserved phenotypes (PUP)
Zhang et al. Advances in genomic selection in domestic animals
CN107338321A (zh) 一种确定最佳snp数量及其通过筛选标记对大黄鱼生产性能进行基因组选择育种的方法
CN104615912B (zh) 一种改进的基于通路的全基因组关联分析算法
Yin et al. Strategy for the simulation and analysis of longitudinal phenotypic and genomic data in the context of a temperature× humidity-dependent covariate
CN105010233A (zh) 一种应用snp辅助选择育种技术培育高繁殖性能种兔的方法
CN114921561B (zh) 杜洛克猪全基因组低密度snp芯片及其制备方法和应用
Guillaume et al. Estimation by simulation of the efficiency of the French marker-assisted selection program in dairy cattle (Open Access publication)
CN111370058B (zh) 一种基于全基因组snp信息追溯水牛血统来源以及进行基因组选配的方法
CN114410746B (zh) 一种东星斑分子溯源选择育种方法及其应用
CN110144414A (zh) 与公猪精子畸形率相关的分子遗传标记及其应用和获取方法
Bondoc Animal breeding: Principles and practice in the Philippine context
CN115305289A (zh) 一种整合snp点集先验信息的鸡腹脂率降低的基因组选择方法
Wei et al. Optimizing the construction and update strategies for the genomic selection of pig reference and Candidate populations in China
CN110195116A (zh) 一种与公猪精子活力相关的分子遗传标记及其应用和获取方法
De Koning et al. Designs for QTL detection in livestock and their implications for MAS
CN112514790B (zh) 水稻分子导航育种方法及应用
Alekya et al. Chapter-7 whole genome strategies for marker assisted selection in plant breeding
CN112837750B (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