本文小编为大家详细介绍“R语言混合线性模型中BLUE值和BLUP值实例分析”,内容详细,步骤清晰,细节处理妥当,希望这篇“R语言混合线性模型中BLUE值和BLUP值实例分析”文章能帮助大家解决疑惑,下面跟着小编的思路慢慢深入,一起来学习新知识吧。
最佳线性无偏预测(best linear unbiased prediction, 简称BLUP),又音译为“布拉普”[1],是统计学上用于线性混合模型对随机效应进行预测的一种方法。最佳线性无偏预测由C.R. Henderson提出。随机效应的最佳线性无偏预测(BLUP)等同于固定效应的最佳线性无偏估计(best linear unbiased estimates, BLUE)(参见高斯-马尔可夫定理)。因为对固定效应使用估计一词,而对随机效应使用预测,这两个术语基本是等同的。BLUP被大量使用于动物育种。——wiki
BLUP值,相当于是对混合线性模型中随机因子的预测;
BLUE值,相当于是对混合线性模型中固定因子的估算
predict means:预测均值,固定因子和随机因子都可以预测均值,它的尺度和表型值尺度一致
将处理作为固定因子
将处理作为固定因子
setwd("D:\\02 ASReml\\blue VS blup")
library(asreml)
library(tidyverse)
dat <- read.csv("MaizeRILs.csv",head=T)
for (i in 1:4) dat[,i] <- as.factor(dat[,i])
as1 <- asreml(height ~ location/rep + location*RIL,data=dat)
ASReml: Tue May 08 11:07:55 2018
LogLik S2 DF wall cpu
-723.8797 64.8862 244 11:07:55 0.1
-723.8797 64.8862 244 11:07:55 0.0
Finished on: Tue May 08 11:07:55 2018
LogLikelihood Converged
#计算品种的BLUE值
ablue <- coef(as1)$fixed
blue1 <- ablue[grep("^RIL_RIL*",rownames(ablue)),] %>% as.data.frame()
head(blue1)
#计算品种的预测均值(predict means)
pv1 <- predict(as1,"RIL")$predictions$pvals
ASReml: Tue May 08 11:13:33 2018
LogLik S2 DF wall cpu
-723.8797 64.8862 244 11:13:33 0.0
-723.8797 64.8862 244 11:13:33 0.0
Finished on: Tue May 08 11:13:33 2018
LogLikelihood Converged
head(pv1) # 类似SAS中的lsmeans
#运行模型:因素作为随机因子
as2 <- asreml(height ~ 1,random = ~location/rep + location*RIL,data=dat)
ASReml: Tue May 08 11:13:34 2018
LogLik S2 DF wall cpu
-1646.5302 233.5135 495 11:13:34 0.0
-1569.0397 137.9186 495 11:13:34 0.0
-1507.3257 94.6888 495 11:13:34 0.0
-1471.3354 74.5149 495 11:13:34 0.0
-1462.9209 67.6142 495 11:13:34 0.0
-1461.7649 65.3553 495 11:13:34 0.0
-1461.7228 64.9069 495 11:13:34 0.0
-1461.7228 64.8863 495 11:13:34 0.0
-1461.7228 64.8862 495 11:13:34 0.0
Finished on: Tue May 08 11:13:34 2018
LogLikelihood Converged
blup <- coef(as2)$random
blup2 <- blup[grep("^RIL_RIL-*",rownames(blup)),] %>% as.data.frame()
head(blup2)
#预测均值
pv2 <- predict(as2,"RIL")$predictions$pvals
ASReml: Tue May 08 11:13:34 2018
LogLik S2 DF wall cpu
-1461.7228 64.8862 495 11:13:34 0.0
-1461.7228 64.8862 495 11:13:34 0.0
-1461.7228 64.8862 495 11:13:34 0.0
-1461.7228 64.8862 495 11:13:34 0.0
Finished on: Tue May 08 11:13:34 2018
LogLikelihood Converged
head(pv2)
#计算遗传力
summary(as2)$varcomp
str(dat)
'data.frame': 496 obs. of 9 variables:
$ location: Factor w/ 4 levels "ARC","CLY","PPAC",..: 1 1 2 2 3 3 4 4 1 1 ...
$ rep : Factor w/ 2 levels "1","2": 1 2 1 2 1 2 1 2 1 2 ...
$ block : Factor w/ 8 levels "1","2","3","4",..: 4 6 5 4 8 5 1 4 1 2 ...
$ plot : Factor w/ 122 levels "1","2","3","4",..: 28 47 36 92 64 40 7 27 6 9 ...
$ RIL : Factor w/ 62 levels "RIL-1","RIL-11",..: 1 1 1 1 1 1 1 1 2 2 ...
$ pollen : int 73 74 71 73 97 95 72 72 69 69 ...
$ silking : int 77 79 74 77 101 100 78 78 71 72 ...
$ ASI : int 4 5 3 4 4 5 6 6 2 3 ...
$ height : num 182 169 213 203 156 ...
VSNR::pin(as2,h3 ~ V3/(V3 + V4/4 + V5/(2*4)))
将数据保存到excel中
library(openxlsx)
write.xlsx(blue1,"blue.xlsx")
write.xlsx(blup2,"blup.xlsx")
write.xlsx(pv1,"pm1.xlsx")
write.xlsx(pv2,"pm2.xlsx")
结果解析
RIL是基因型
pm2-random是RIL作为随机因子的预测均值
pm1-fixed是RIL作为固定因子时的预测均值
blue是RIL作为固定因子的BLUE值
blup是RIL作为随机因子的BLUP值
pm2-blup 是随机因子的预测均值 减去 随机因子的BLUP值,可以看到得到的是一个常数(均值)
pm1-mu-random 是固定因子的预测均值 减去 固定依着你的BLUE值, 可以看到不是一个常数
blup/blue_effect=heritibility 是BLUP值 除以 BLUE效应值,得到的是遗传力常数
备注:blue_effect是用固定因子的预测均值 减去 整体均值
读到这里,这篇“R语言混合线性模型中BLUE值和BLUP值实例分析”文章已经介绍完毕,想要掌握这篇文章的知识点还需要大家自己动手实践使用过才能领会,如果想了解更多相关内容的文章,欢迎关注亿速云行业资讯频道。
亿速云「云服务器」,即开即用、新一代英特尔至强铂金CPU、三副本存储NVMe SSD云盘,价格低至29元/月。点击查看>>
免责声明:本站发布的内容(图片、视频和文字)以原创、转载和分享为主,文章观点不代表本网站立场,如果涉及侵权请联系站长邮箱:is@yisu.com进行举报,并提供相关证据,一经查实,将立刻删除涉嫌侵权内容。