育种数据分析之放飞自我分享 http://blog.sciencenet.cn/u/yijiaobai 关注:生物统计,数量遗传,混合线性模型,生物信息,R,Perl,Python,GWAS,GS相关方法,文章及代码

博文

全基因组选择介绍及实践-1

已有 5432 次阅读 2019-3-2 11:57 |个人分类:数量遗传学|系统分类:论文交流

1, 序言

这篇是基因组选择的理论加实践,因为我看到一句话,Talk is cheap. Show me the code,很有感触,有感而写。使用的包是R的sommer和asreml,其实强健的还是成熟的软件,比如DMU,BLUPF90,PIBLUP,ASreml等,但sommer作为基本功能的演示,非常合适。

2, 定义

基因组选择(Genomic Selection, GS), 利用覆盖全基因组的高密度分子遗传标记进行的标记辅助选择.

在这里插入图片描述

  • 表型选择, 这对于遗传力高的性状选择有效, 即表型可以遗传比例较大的性状, 比如植物里面千粒重, 比如动物里面的初生重等.
  • EBV selection, 应用系谱的动物模型
  • Marker assisted selection(MAS), 分子标记辅助选择, 对于有主效QTL或者主效基因的性状, 有优势.
  • Genomic selection, 全基因组选择

选择进展的定义
$$ 选择进展 = \frac{选择强度 \times 选择准确性\times 遗传标准差}{世代间隔}$$

  • 选择强度: Intensity of selection
  • 选择准确性: Accuracy of selection
  • 遗传标准差: Genetic standard deviation
  • 世代间隔: Generation interval

3, 不同选择方法的比较

基于系谱的动物模型
局限:

  • 对低遗传力性状效果较差
  • 不容易度量的性状(胴体性状, 肉质性状)效果较差
  • 不能早期度量的性状, 效果较差

分子标记辅助育种(MAS)
局限:

  • 需要先对主效基因或者QTL进行检测
  • 不同群体变化较大
  • 标记可解释的遗传变异百分比较低
  • 在动物育种中的应用非常有限

全基因组选择
优点:

  • 无需进行主效基因或者QTL的检测
  • 不依赖于表型信息(候选群)
  • 能够捕获基因组中的全部变异
  • 对于低遗传力, 难以度量的性状提升效果明显

4, 基因组选择流程

在这里插入图片描述

  • 1, 建立参考群体
    • 获得每个个体的性状表型值
    • 测定每个个体的SNP基因型(使用SNP芯片)
    • 估计SNP效应值
  • 2, 在候选群体中进行基因组选择
    • 测定候选个体的SNP基因型
    • 计算个体的GEBV
    • 依据GEBV进行选择

动物模型, GBLUP方法和Single-setp方法

在这里插入图片描述
三种方法的区别

  • 动物模型是利用的系谱构建的A矩阵
  • GBLUP是利用基因组信息构建的G矩阵
  • 一步法(single-setp)是利用系谱和基因组信息构建的H矩阵
    在这里插入图片描述

    5,其它方法

    除了GBLUP和Single-step, 还有其它方法用于基因组选择, 不过应用最广最强建最有优势的还是GBLUP和Single-step.
  • RR-BLUP(Rirdge Regression BLUP)
  • Bayesian method
    • Bayes A
    • Bayes B
    • Bayes C\pi
    • Bayesian Lasso
  • Principal component regression
  • Semiparametric procedures
  • Machine learning methods

6,Talk is cheap, show me the code

别说话, 上代码

数据描述

  • 数据使用QMSim软件模拟而成
  • 共有2000个体
  • 1200个有基因型
  • 共有5个世代, 3,4为参考群, 5为候选群, 每个世代400个体
  • 系谱数据为50k
  • 遗传力为0.3

软件使用

  • QMSim模拟数据
  • R包: learnasreml
  • R包: sommer
  • R包:asreml
# 载入软件包
if (!requireNamespace("devtools")) install.packages("devtools")
library(devtools)
if (!requireNamespace("learnasreml"))  install_github("dengfei2013/learnasreml")
library(learnasreml)
if (!requireNamespace("sommer")) install.packages("sommer")
library(sommer)

# 载入数据
data("gs_geno")
data("gs_phe")

# 动物模型animal model
gphe = gs_phe[gs_phe$G %in%c(3:5), ] # 这里提取测序个体
gphe$Progeny = as.factor(as.character(gphe$Progeny)) # 将Progeny转化为因子
ainv = asreml.Ainverse(gphe[,1:3])$ginv # 计算A逆矩阵
moda_as = asreml(Phen ~ Sex + G, random = ~ ped(Progeny),ginverse= list(Progeny=ainv),data = gphe) # 构建模型
summary(moda_as)$varcomp #结果

# GBLUP sommer的操作代码
gdat = gs_geno #  重命名数据
gdat[1:10,1:10] # 预览数据
G = A.mat(gdat-1) # sommer计算G矩阵时, 支持的是-1,0,1, 所以需要转化
time_sommer = system.time({
  modg_mm = mmer2(Phen ~ Sex + G, random = ~ g(Progeny),G = list(Progeny=G),data = gphe) # 构建模型
})

# GBLUP asreml的操作代码
library(asreml)
diag(G) = diag(G)+0.01 # 矩阵奇异, 对角线+0.01防止奇异
ginv = write_relation_matrix(G,type="ginv") # 转化为三元组形式
attr(ginv,"rowNames") = rownames(G) # 对ginv进行rowNames的设置
time_asreml = system.time({
  modg_as = asreml(Phen ~ Sex + G, random = ~ giv(Progeny),ginverse= list(Progeny=ginv),data = gphe,workspace=1e9) # 构建模型
})

# 比较sommer和asreml结果 
summary(modg_mm)
summary(modg_as)$varcomp

# 比较运行时间
time_sommer
time_asreml

# 计算准确性:asreml
blup_as = as.data.frame(coef(modg_as)$random)
re1 = cbind(gphe,blup_as)
cor(re1$Polygene,re1$effect)

# 计算准确性:sommer
blup_mm = as.data.frame(randef(modg_mm))
head(blup_mm)
names(blup_mm) = "blup_mm"
re2 = cbind(gphe,blup_mm)
cor(re1$Polygene,re2$blup_mm)

方差组分结果比较

> # 比较sommer和asreml结果 
> summary(modg_mm)$var.comp.table
                       VarComp  VarCompSE    Zratio
g(Progeny).Phen-Phen 0.1947930 0.04421259  4.405827
units.Phen-Phen      0.8076997 0.04320582 18.694236

> summary(modg_as)$varcomp
                     gamma component  std.error   z.ratio constraint
giv(Progeny).giv 0.2417449 0.1947875 0.04518119  4.311252   Positive
R!variance       1.0000000 0.8057563 0.04384776 18.376227   Positive

可以看出, sommer和asreml方差组分结果基本一致, Vg=0.194, Ve=0.805

运行时间比较

> # 比较运行时间
> time_sommer
 用户  系统  流逝 
33.53  3.88  9.47 
> time_asreml
 用户  系统  流逝 
19.00  3.90 23.03

可以看到,sommer运行时间很快,sommer支持的是G矩阵类型,如果使用pedigree,有点麻烦。asreml-r版慢一点,个人体会,asreml-w更快一点。

准确性比较

> # 计算准确性:asreml
> blup_as = as.data.frame(coef(modg_as)$random)
> re1 = cbind(gphe,blup_as)
> cor(re1$Polygene,re1$effect)
[1] 0.6338289
> 
> # 计算准确性:sommer
> blup_mm = as.data.frame(randef(modg_mm))
Returning object of class 'list' where each element correspond to one random effect.
> names(blup_mm) = "blup_mm"
> re2 = cbind(gphe,blup_mm)
> cor(re1$Polygene,re2$blup_mm)
[1] 0.6327671

因为这是模拟数据,有真值(True breeding value),如果是真实数据,需要用交叉验证的手段评价准确性。

7,后续更新,欢迎持续关注微信公众号

  • 1,H矩阵构建及一步法计算
  • 2,H矩阵构建不同参数的设置
  • 3,多性状GBLUP及多性状一步法计算
  • 4,如何进行交叉验证
  • 5,如何筛选以及填充基因组数据
  • 6,基于贝叶斯的全基因组选择
  • 7,基于一步法的GWAS分析

公众号ID:R-breeding

在这里插入图片描述

8, 参考资料

张勤老师: 动物基因组选择-理论与应用

Meuwissen T H , Hayes B J , Goddard M E . Prediction of total genetic value using genome-wide dense marker maps[J]. Genetics, 2001, 157(4):1819-29.
Vanraden P M. Efficient methods to compute genomic predictions[J]. Journal of Dairy Science, 2008, 91(11):4414-4423.
Legarra A, Aguilar I, Misztal I. A relationship matrix including full pedigree and genomic information[J]. Journal of Dairy Science, 2009, 92(9):4656-4663.
Aguilar I, Misztal I, Johnson D L, et al. Hot topic: a unified approach to utilize phenotypic, full pedigree, and genomic information for genetic evaluation of Holstein final score.[J]. Journal of Dairy Science, 2010, 93(2):743-752.
Ishizaka A, Labib A. Review of the main developments in the analytic hierarchy process[J]. Expert Systems with Applications, 2011, 38(11):14336-14345.



http://blog.sciencenet.cn/blog-2577109-1165172.html

上一篇:R语言与独孤九剑
下一篇:全基因组选择介绍及实践-2:构建H矩阵

0

该博文允许注册用户评论 请点击登录 评论 (0 个评论)

数据加载中...

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2021-3-5 18:36

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部