Nothing
asr41.lib <- "D:\\Analyses\\R ASReml4.1"
cat("#### Test estimateV specials for Cane with asreml41\n")
test_that("HEB25_estimateV_asreml41", {
skip_if_not_installed("asreml")
skip_on_cran()
library(dae)
library(asreml, lib.loc = asr41.lib)
library(asremlPlus)
### Parana Local 851
### 200 test lines, 40 check plots
### 10 Columns x 24 rows
data(local851)
test.specials <- c( "ar1", "ar2", "ar3", "sar","sar2",
"ma1", "ma2", "arma", "exp", "gau",
"cor")
vpar.vals <- c(-1.274424e-01, 0.0669068, -2.265849e-01, -1.238088e-01, -1.57972e-01,
1.200424e-01, -4.0854110e-02, 0, 3.116814e-08, 2.721824e-08, 4.134244e-02)
# vpar.vals <- c(-1.278067e-01, 6.676302e-02, -2.265831e-01, -1.242148e-01, -1.580255e-01,
# 1.186700e-01, -4.063130e-02, 0, 1.937006e-07, 1.549264e-07, 4.164659e-02 )
names(vpar.vals) <- test.specials
V.el <- c(-14.6976, -14.82033, -12.26416, -13.75973, -13.78844,
-12.3181, -11.08808, 0, 10.94163, 10.94113, 10.41591)
# V.el <- c(-14.72583, -14.82721, -12.17283, -13.78792, -13.80899,
# -12.23504, -11.12781, 0, 10.91292, 10.91283, 10.48848)
names(V.el) <- test.specials
### Model with genetic variance only
for (func in test.specials[test.specials != "arma"])
{
ranform <- as.formula(paste("~ ar1(Col):", func, "(Row)", sep = ""))
asreml.obj <- asreml(tch ~ Control/Check,
random = ~ Col + Row + New,
residual = ranform,
data=site2,
na.action=na.method(y="include", x="include"))
print(asreml.obj$vparameters[length(asreml.obj$vparameters)])
cat("\n",func,": ", asreml.obj$vparameters[length(asreml.obj$vparameters)], " and ",vpar.vals[func],"\n\n")
testthat::expect_true(abs(asreml.obj$vparameters[length(asreml.obj$vparameters)] -
vpar.vals[func]) < 1e-03)
V <- estimateV(asreml.obj)
cat("\n",func,": ", V[2, 1], " and ",V.el[func],"\n\n")
testthat::expect_true(abs(V[2, 1] - V.el[func]) < 1e-03)
}
### This had a bug, but seems to be working - need to see if can get to converge
testthat::expect_error(asreml.obj <- asreml(tch ~ Control/Check,
random = ~ Col + Row + New,
residual = ~ ar1(Col):corb(Row, b = 3),
data=site2,
na.action=na.method(y="include", x="include")))
testthat::expect_true(asreml.obj$converge)
testthat::expect_equal(nrow(summary(asreml.obj)$varcomp), 6)
### Model with genetic competition without pedigree
asreml.options(design = TRUE)
asreml.obj <- asreml(tch ~ Control/Check,
random=~ Col + Row + str(~New+ grp(neighbour),~us(2):id(201)),
residual =~ ar1(Col):ar1(Row),
group=list(neighbour=c(11:211)),data=site2,
na.action=na.method(y="include", x="include"))
summary.asreml(asreml.obj)$varcomp
asreml.obj$vparameters
V.g <- with(site2, fac.vcmat(Col, asreml.obj$vparameters["Col"]) +
fac.vcmat(Row, asreml.obj$vparameters["Row"]))
G.g <- kronecker(matrix(asreml.obj$vparameters[c(3,4,4,5)], nrow = 2, ncol = 2), mat.I(201))
cols <- c(grep("New", colnames(asreml.obj$design), fixed = TRUE),
grep("grp", colnames(asreml.obj$design), fixed = TRUE))
V.g <- V.g + (asreml.obj$design[, cols] %*% G.g %*% t(as.matrix(asreml.obj$design[, cols])))
ar1C <- fac.ar1mat(site2$Col, asreml.obj$vparameters["Col:Row!Col!cor"])
ar1R <- fac.ar1mat(site2$Row, asreml.obj$vparameters["Col:Row!Row!cor"])
V.g <- asreml.obj$sigma2 * (V.g + ar1C * ar1R)
V <- estimateV(asreml.obj)
testthat::expect_true(all(abs(V - V.g) < 1e-06))
#Test for grp function
asreml.obj <- asreml(tch ~ Control/Check,
random=~ Col + Row + New+ grp(neighbour),
residual =~ ar1(Col):ar1(Row),
group=list(neighbour=c(11:211)),data=site2,
na.action=na.method(y="include", x="include"))
summary.asreml(asreml.obj)$varcomp
asreml.obj$vparameters
V.g <- with(site2, fac.vcmat(Col, asreml.obj$vparameters["Col"]) +
fac.vcmat(Row, asreml.obj$vparameters["Row"]) +
fac.vcmat(New, asreml.obj$vparameters["New"]))
G.g <- asreml.obj$vparameters["grp(neighbour)"] * mat.I(201)
cols <- grep("grp", colnames(asreml.obj$design), fixed = TRUE)
V.g <- V.g + (asreml.obj$design[, cols] %*% G.g %*% t(as.matrix(asreml.obj$design[, cols])))
ar1C <- fac.ar1mat(site2$Col, asreml.obj$vparameters["Col:Row!Col!cor"])
ar1R <- fac.ar1mat(site2$Row, asreml.obj$vparameters["Col:Row!Row!cor"])
V.g <- asreml.obj$sigma2 * (V.g + ar1C * ar1R)
V <- estimateV(asreml.obj)
testthat::expect_true(all(abs(V - V.g) < 1e-03))
asreml.options(design = FALSE)
})
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.