1. 为何要计算BLUE值?
一年多点或者多年多点的植物数据中,一个基因型(品种)往往有多个表型数据,但只有一个基因型,在GWAS关联分析中,就需要一个基因型对应一个表型数据。
之所以有多个表型数据的原因:
- 或者是多个重复
- 或者是多个地点的数据
- 或者是多个年份的数据
问题:如何计算得到一个表型数据呢?
解答:可以使用多个表型值的平均值,作为品种的表型值,现在有更好的方法:BLUE值。
2. 为何使用BLUE值?
一般,有两个选择,BLUE值或者BLUP值,在GWAS中大都使用的BLUE值。
BLUE和BLUP的区别:
- BLUE值是混合线性模型中固定因子的估计效应值
- BLUP值是混合线性模型中随机因子的估计效应值
BLUE和BLUP的代表:
- BLUE值着重在于评估品种现在的表现
- BLUP值着重在于预测品种将来的表现
BLUE和BLUP的方差变化
- BLUE只是对表型值根据地点,年份进行矫正,得到的数据和原来数据尺度一样
- BLUP值会对表型数据进行压缩
3. 示例数据
数据为learnasreml
中的MET数据集。数据包括2年,5个地点,每个地点4个重复,共有10品种,观测值为产量(yield)
代码
library(learnasreml)
library(magrittr)
data("MET")
head(MET)
数据
4. lme4包如何分析
模型:
- 固定因子:Cul
- 随机因子:Year + Location + Location:Rep
代码:
dat = MET
m1 = lmer(Yield ~ Cul + (1|Location) + (1|Location:Rep) + (1|Year), data=dat)
as.data.frame(fixef(m1))
注意: 植物中,一般的BLUE值需要加上截距(Intercept)。因为BLUE值中,第一个水平会当做0,其它为相对值,可以手动进行相加,也可以使用lsmeans
包中的lsmeans
。
library(lsmeans)
re = lsmeans(m1,"Cul")
re
数据中的lsmeans即为品种的BLUE值,可以作为GWAS或者GS的表型值进行后续的计算。
5. asreml对比结果
总所周知,asreml
是一个非常强大的商业软件,如果用asreml进行结果对比,可以判断lme4计算是否正确。
代码
library(asreml)
m2 = asreml(Yield ~ Cul , random = ~Location/Rep + Year, data=dat)
pre = as.data.frame(predict(m2,classify = "Cul")$pvals)
pre
结果
结果中的predicted.value即为品种的BLUE值。
两者结果对比
library(tidyverse)
ls_value = as.data.frame(re)
as_value = pre
head(ls_value)
head(as_value)
merge(ls_value,as_value,by = "Cul") %>% select(lsmean ,predicted.value) %>% cor
相关性分析结果:
散点图:
可以看出,两者结果完全一致。