set.seed(1337)
n <- 100
x1 <- runif(n)
x2 <- runif(n)
x3 <- runif(n)
y <- rnorm(n, 0 + 1 * x1 + 1 * x3, exp(-3 + 1 * x2 + 1 * x3))
m <- lmls(y ~ x1 + x3, ~ x2 + x3, light = FALSE)
# score -----------------------------------------------------------------------
beta <- gamma <- c(0, 0, 0)
m <- set_coef(m, "location", beta)
m <- set_coef(m, "scale", gamma)
f <- function(x, predictor) {
loglik(set_coef(m, predictor, x))
}
test_that("score of beta is correct", {
num_score <- numDeriv::grad(f, beta, predictor = "location")
expect_exactly(score(m, "location"), num_score)
})
test_that("score of gamma is correct", {
num_score <- numDeriv::grad(f, gamma, predictor = "scale")
expect_exactly(score(m, "scale"), num_score)
})
# fisher info -----------------------------------------------------------------
beta <- c(0, 1, 1)
gamma <- c(-3, 1, 1)
m <- set_coef(m, "location", beta)
m <- set_coef(m, "scale", gamma)
nsim <- 5000
reps <- replicate(nsim, {
y <- rnorm(n, 0 + 1 * x1 + 1 * x3, exp(-3 + 1 * x2 + 1 * x3))
m <- lmls(y ~ x1 + x3, ~ x2 + x3, light = FALSE)
m <- set_coef(m, "location", beta)
m <- set_coef(m, "scale", gamma)
c(score_beta(m), score_gamma(m))
})
num_info <- cov(t(reps))
test_that("fisher info of beta is correct", {
expect_roughly(info_beta(m), num_info[1:3, 1:3])
})
test_that("fisher info of gamma is correct", {
expect_roughly(info_gamma(m), num_info[4:6, 4:6])
})
test_that("chol_info() works", {
expect_equal(chol_info(m, "location"), chol(info_beta(m)))
expect_equal(chol_info(m, "scale"), chol(info_gamma(m)))
})
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.