yzhlinscau的个人博客分享 http://blog.sciencenet.cn/u/yzhlinscau

博文

双性状分析之MCMCglmm篇

已有 8823 次阅读 2013-10-25 14:12 |个人分类:R|系统分类:科研笔记

对于多性状分析,R中的程序包MCMCglmm也可进行asreml-r的类似分析,只是代码更为复杂。

 MCMCglmm程序包代码如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
############ 代码清单 ##########
library( MCMCglmm ) 
phen.var <- matrix( c( var( df.2$dj, na.rm = T) , 0, 0,
                  var( df.2$h5, na.rm = T )), 2, 2)
prior<-list(G=list( G1 = list( V = phen.var / 2, n = 2 )), # G结构初始值
              R = list( V = phen.var / 2, n = 2 ))     # R结构初始值
 
m2.glmm <- MCMCglmm( cbind(dj, h5) ~ trait -1 + trait : Rep,
                     random = ~ us(trait):Fam,
                     rcov = ~ us(trait) : units,
                     data = df.2, family = c( "gaussian", "gaussian" ),
                     nitt = 130000, thin = 100, burnin = 30000,
                     prior = prior, verbose = F ) 
# 结果提取
summary(m2.glmm)
posterior.mode(m2.glmm$VCV)
 
# 计算遗传力
A.h2.glmm <- 4 * m2.glmm$VCV[,'dj:dj.Fam'] / ( m2.glmm$VCV[,'dj:dj.Fam']  +  m2.glmm$VCV[,'dj:dj.units'])
posterior.mode(A.h2.glmm)
 
B.h2.glmm <- 4 * m2.glmm$VCV[,'h5:h5.Fam'] / (m2.glmm$VCV[, 'h5:h5.Fam'] +  m2.glmm$VCV[,'h5:h5.units'])
 
# 计算遗传相关
gCorr.glmm <-m2.glmm$VCV[,'h5:dj.Fam'] / sqrt( m2.glmm$VCV[, 'dj:dj.Fam']  * m2.glmm$VCV[, 'h5:h5.Fam'])
 
# 计算表型相关

pCorr.glmm <- (m2.glmm$VCV[,'h5:dj.Fam'] + m2.glmm$VCV[, 'h5:dj.units'] ) / sqrt(( m2.glmm$VCV[, 'dj:dj.Fam'] + m2.glmm$VCV[, 'dj:dj.units'] ) * ( m2.glmm$VCV[, 'h5:h5.Fam'] + m2.glmm$VCV[,'h5:h5.units']))

posterior.mode(B.h2.glmm);posterior.mode(gCorr.glmm);posterior.mode(pCorr.glmm)


运行结果如下:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
> posterior.mode(m2.glmm$VCV)   # 方差分量
dj:dj.Fam     h5:dj.Fam     dj:h5.Fam     h5:h5.Fam   dj:dj.units
   7.723708e-05  4.709730e-02  4.709730e-02  6.967567e+02  4.945213e-04
h5:dj.units   dj:h5.units   h5:h5.units
-1.887868e-01 -1.887868e-01  5.417873e+03
 
> posterior.mode(A.h2.glmm)   # dj遗传力
    var1
0.4764733
 
> posterior.mode(B.h2.glmm)   # h5遗传力
    var1
0.5959831
 
> posterior.mode(gCorr.glmm)  # 遗传相关
    var1
0.2416088
 
> posterior.mode(pCorr.glmm)  # 表型相关
     var1
-0.0980515

   比较asreml-rMCMCglmm分析的结果可以发现,性状dj的方差分量、遗传力结果都比较接近,而对于h5MCMCglmm所得的方差分量偏大。两性状之间的遗传协方差和误差协方差、表型相关结果基本一致,但遗传相关,MCMCglmm所得结果偏小。




https://wap.sciencenet.cn/blog-1114360-736009.html

上一篇:G×E分析之ASReml-R篇
下一篇:主成分分析之R篇
收藏 IP: 58.254.92.*| 热度|

0

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

数据加载中...
扫一扫,分享此博文

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

GMT+8, 2024-6-2 13:05

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部