所以我正在分析4种环境下水稻育种品系的多地点试验。简化后的数据如下:
https://drive.google.com/file/d/1jilVXX8JMkZCDVtIRmrwzB55kgR2GtYB/view?usp=sharing
我在sommer上运行一个简单的非结构化方差模型。我已经在lme4和lme3上做过了,但我想继续用summer。模型为:
m3 <- mmer(RDM ~ ENV ,
random=~ vsr(usr(ENV),GEN),
rcov=~ units,
data=d)
很简单,不是吗?然而,很快我就得到了错误:
System is singular (V). Stopping the job. Try a bigger number of tolparinv.
所以,好吧,我尝试一个更大的tolparinv数(因为我不能做一个更简单的模型)。但是使这个函数成立的最小的数是1000。所以,我的问题是:这意味着什么?
此外,我们假设这样运行模型是可以的。现在的情况是很多方差分量都是负的。这没有多大意义。
有人能解释一下吗?所以,具体的问题是:
- 为什么奇点来得这么快?
- 如果我增加这么多的parinv会发生什么? 这就是为什么我的方差是负的吗?
- 最重要的是:这个是可以修复的吗?如何?
谢谢! !
非结构化模型由于协方差参数较多,在REML过程中得到非全秩协方差矩阵的可能性较大,容易出现收敛性和奇异性问题。您可能想尝试降低秩的因子分析模型,以保证收敛并提高速度。普通的美国模型是这样的
library(sommer)
data(DT_example)
DT <- DT_example
A <- A_example
head(DT)
DT=DT[with(DT, order(Env)), ]
ans2us <- mmec(Yield~Env,
random= ~ vsc(usc(Env),isc(Name)),
rcov= ~ vsc(dsc(Env),isc(units)),
data=DT)
plot(ans2us$llik[1,])
summary(ans2us)$varcomp
pp <- predict(ans2us, D="Name")
pp$pvals
FA实现看起来像这样
# less accurate but is guaranteed to converge
M <- with(DT, rrc(Env, Name, Yield, nPC=2))
ans2fa <- mmec(Yield~Env,
random= ~ vsc(usc( M ),isc(Name)),
rcov= ~ vsc(dsc(Env),isc(units)),
nIters = 100,
data=DT)
summary(ans2fa)$varcomp
plot(ans2fa$llik[1,])
ppfa <- predict(ans2fa, D="Name")
Lam <- with(DT, rrc(Env, Name, Yield, nPC=2, returnLam = TRUE))$Lam
scores <- ans2fa$uList[[1]]
BLUP <- scores %*% t(Lam) # environment specific predictions
# highly correlated predictions
plot(pp$pvals$predicted.value, ppfa$pvals$predicted.value)
欢呼