在模型假定中,数据需要符合正态分布,在动物模型分析中,如果对数据进行转化,对育种值和方差组分变化的情况,进行测试。
1library(asreml)
2library(learnasreml)
3
4data("animalmodel.dat")
5data("animalmodel.ped")
6
7dat = animalmodel.dat
8ped = animalmodel.ped
9dat[dat==0] = NA
10
11# 原始数据分布
12hist(dat$BWT)
方差组分:
1# 原始数据,计算动物模型方差组分
2head(dat)
3ainv = ainverse(ped)
4m1 = asreml(BWT ~ BYEAR + SEX, random = ~ vm(ANIMAL,ainv),data=dat)
5summary(m1)$varcomp
育种值blup:
1# 计算原始数据育种值blup
2blup = tiqu_blup(coef(m1)$random)
3head(blup)
结论:
1# 原始数据减去一个常数:10
2# 方差组分不变,育种值不变
3dat$a = dat$BWT - 10
4hist(dat$a)
5m2 = asreml(a ~ BYEAR + SEX, random = ~ vm(ANIMAL,ainv),data=dat)
6summary(m2)$varcomp
7blup1 = tiqu_blup(coef(m2)$random)
8head(blup1)
9cor(blup$effect,blup1$effect)
表型数据分布:
方差组分:
结论:
1# 原始数据,log10(x)转化
2# 方差组分变化,育种值与原来育种值相关性系数:0.9
3dat$a = log10(dat$BWT)
4hist(dat$a) # 偏离正态分布
5m2 = asreml(a ~ BYEAR + SEX, random = ~ vm(ANIMAL,ainv),data=dat)
6summary(m2)$varcomp
7blup1 = tiqu_blup(coef(m2)$random)
8cor(blup$effect,blup1$effect)
9plot(blup$effect,blup1$effect)
表型数据分布:
方差组分:
和原育种值相关性分析:
结论:
1# 原始数据,log10(x +10)转化
2# 方差组分变化,育种值与原来育种值相关性系数:0.99
3dat$a = log10(dat$BWT+10)
4hist(dat$a) # 符合正态分布
5m2 = asreml(a ~ BYEAR + SEX, random = ~ vm(ANIMAL,ainv),data=dat)
6summary(m2)$varcomp
7blup1 = tiqu_blup(coef(m2)$random)
8cor(blup$effect,blup1$effect)
9plot(blup$effect,blup1$effect)
表型数据分布:
方差组分:
育种值相关性分析:
方差组分变化:
遗传力变化:
所以,不要轻易做数据转化,多多考虑加大数据量,或者提高数据质量吧。统计方法不能拯救失败的数据,只能告诉你数据分析为什么会失败,这是我说的。
下面是Fisher说的:
科学试验做完后再找统计学家分析数据,如同病人死了找医生进行尸体解剖,医生会告诉你病人死的原因。同样,统计学家会告诉你试验失败的原因。 ------Ronald Fisher
All models are wrong, but some models are useful | 模型皆有误,或尤建奇功。 George E. P. Box 在终极的分析中,一切知识都是历史;在抽象的意义下,一切科学都是数学;在理性的基础上,所有的判断都是统计学。 C.R.Rao 统计学是对令人困惑费解的问题做出数字设想的艺术。 David Freedma 统计具有非凡的能力处理各种复杂的问题,它需要非常精细的方法和小心翼翼的解释。当人类科学探索者在问题的丛林中遇到难以逾越的障碍时,唯有统计工具可为其开辟一条前进的通道。 Francis Galton 对统计学的一知半解常常会造成一些不必要的上当受骗,对统计学的一概排斥往往会造成某些不必要的愚昧无知。 C.R.Rao