dat <- denis.missing
# View missingness structure
require(reshape2)
acast(dat, env~gen, value.var='yield')
# Stability variance (Table 3 in Piepho)
library("nlme")
m1 <- lme(yield ~ -1 + gen, data=dat, random= ~ 1|env,
weights = varIdent(form= ~ 1|gen),
na.action=na.omit)
svar <- m1$sigma^2 * c(1, coef(m1$modelStruct$varStruct, unc = FALSE))^2
round(svar, 2)
## G5 G3 G1 G2
## 39.25 22.95 54.36 12.17 23.77Run the code above in your browser using DataLab