Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

ael wip #7

Merged
merged 5 commits into from
Apr 12, 2024
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension


Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion .github/workflows/test-coverage.yaml
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
# Need help debugging build failures? Start at https://github.com/r-lib/actions#where-to-find-help
on:
push:
branches: [main, master]
branches: [main, master, dev]
pull_request:
branches: [main, master]

Expand Down
2 changes: 1 addition & 1 deletion DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,7 +1,7 @@
Type: Package
Package: melt
Title: Multiple Empirical Likelihood Tests
Version: 1.11.2
Version: 1.11.3
Authors@R: c(
person("Eunseop", "Kim", , "markean@pm.me", role = c("aut", "cph", "cre")),
person("Steven", "MacEachern", role = c("ctb", "ths")),
Expand Down
1 change: 1 addition & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -29,6 +29,7 @@ exportMethods(sigTests)
exportMethods(summary)
exportMethods(weights)
importFrom(Rcpp,sourceCpp)
importFrom(checkmate,assert_choice)
importFrom(checkmate,assert_class)
importFrom(checkmate,assert_int)
importFrom(checkmate,assert_logical)
Expand Down
8 changes: 8 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
@@ -1,3 +1,11 @@
# melt 1.11.3
## NEW FEATURES
* New `"ael"` option has been added in the `calibrate` argument of `elt()` for adjusted empirical likelihood calibration.

## MINOR IMPROVEMENTS
* The package vignette has been updated.


# melt 1.11.2
## MINOR IMPROVEMENTS
* The package vignette has been updated.
Expand Down
4 changes: 3 additions & 1 deletion R/AllClasses.R
Original file line number Diff line number Diff line change
Expand Up @@ -24,6 +24,8 @@
#' @slot nthreads A single integer for the number of threads for parallel
#' computation via OpenMP (if available).
#' @slot seed A single integer for the seed for random number generation.
#' @slot an A single numeric representing the scaling factor for adjusted
#' empirical likelihood calibration.
#' @slot b A single integer for the number of bootstrap replicates.
#' @slot m A single integer for the number of Monte Carlo samples.
#' @aliases ControlEL
Expand All @@ -33,7 +35,7 @@ setClass("ControlEL",
slots = c(
maxit = "integer", maxit_l = "integer", tol = "numeric", tol_l = "numeric",
step = "ANY", th = "ANY", verbose = "logical", keep_data = "logical",
nthreads = "integer", seed = "ANY", b = "integer", m = "integer"
nthreads = "integer", seed = "ANY", an = "ANY", b = "integer", m = "integer"
)
)

Expand Down
34 changes: 23 additions & 11 deletions R/AllGenerics.R
Original file line number Diff line number Diff line change
Expand Up @@ -283,7 +283,7 @@ setGeneric("elmt", function(object,

#' Empirical likelihood test
#'
#' Tests a linear hypothesis.
#' Tests a linear hypothesis with various calibration options.
#'
#' @param object An object that inherits from \linkS4class{EL}.
#' @param rhs A numeric vector or a column matrix for the right-hand side of
Expand All @@ -295,9 +295,9 @@ setGeneric("elmt", function(object,
#' of parameters. Or a character vector with a symbolic description of the
#' hypothesis is allowed. Defaults to `NULL`. See ‘Details’.
#' @param alpha A single numeric for the significance level. Defaults to `0.05`.
#' @param calibrate A single character for the calibration method. It is
#' case-insensitive and must be one of `"chisq"`, `"boot"`, or `"f"`.
#' Defaults to `"chisq"`. See ‘Details’.
#' @param calibrate A single character representing the calibration method. It
#' is case-insensitive and must be one of `"ael"`, `"boot"`, `"chisq"`, or
#' `"f"`. Defaults to `"chisq"`. See ‘Details’.
#' @param control An object of class \linkS4class{ControlEL} constructed by
#' [el_control()]. Defaults to `NULL` and inherits the `control` slot in
#' `object`.
Expand All @@ -318,12 +318,13 @@ setGeneric("elmt", function(object,
#' 1. If `lhs` is `NULL`, \eqn{L} is set to the identity matrix and the
#' problem reduces to evaluating at \eqn{r} as \eqn{l(r)}.
#'
#' `calibrate` specifies the calibration method used. Three methods are
#' available: `"chisq"` (chi-square calibration), `"boot"` (bootstrap
#' calibration), and `"f"` (\eqn{F} calibration). `"boot"` is applicable only
#' when `lhs` is `NULL`. The `nthreads`, `seed`, and `B` slots in `control`
#' apply to the bootstrap procedure. `"f"` is applicable only to the mean
#' parameter when `lhs` is `NULL`.
#' `calibrate` specifies the calibration method used. Four methods are
#' available: `"ael"` (adjusted empirical likelihood calibration), `"boot"`
#' (bootstrap calibration), `"chisq"` (chi-square calibration), and `"f"`
#' (\eqn{F} calibration). When `lhs` is not `NULL`, only `"chisq"` is
#' available. `"f"` is applicable only to the mean parameter. The `an` slot in
#' `control` applies specifically to `"ael"`, while the `nthreads`, `seed`,
#' and `B` slots apply to `"boot"`.
#' @return An object of class of \linkS4class{ELT}. If `lhs` is non-`NULL`, the
#' `optim` slot corresponds to that of \linkS4class{CEL}. Otherwise, it
#' corresponds to that of \linkS4class{EL}.
Expand All @@ -333,6 +334,11 @@ setGeneric("elmt", function(object,
#' \emph{Statistical Methods & Applications}, **19**(4), 463--476.
#' \doi{10.1007/s10260-010-0137-9}.
#' @references
#' Chen J, Variyath AM, Abraham B (2008).
#' ``Adjusted Empirical Likelihood and Its Properties.''
#' \emph{Journal of Computational and Graphical Statistics}, **17**(2),
#' 426--443.
#' @references
#' Kim E, MacEachern SN, Peruggia M (2024).
#' ``melt: Multiple Empirical Likelihood Tests in R.''
#' \emph{Journal of Statistical Software}, **108**(5), 1--33.
Expand All @@ -345,9 +351,15 @@ setGeneric("elmt", function(object,
#' @seealso \linkS4class{EL}, \linkS4class{ELT}, [elmt()], [el_control()]
#' @usage NULL
#' @examples
#' ## F calibration for the mean
#' ## Adjusted empirical likelihood calibration
#' data("precip")
#' fit <- el_mean(precip, 32)
#' elt(fit, rhs = 100, calibrate = "ael")
#'
#' ## Bootstrap calibration
#' elt(fit, rhs = 32, calibrate = "boot")
#'
#' ## F calibration
#' elt(fit, rhs = 32, calibrate = "f")
#'
#' ## Test of no treatment effect
Expand Down
8 changes: 7 additions & 1 deletion R/calibrate.R
Original file line number Diff line number Diff line change
Expand Up @@ -14,7 +14,7 @@
#' @noRd
calibrate <- function(calibrate, alpha, statistic, p, par, object, control) {
switch(calibrate,
"chisq" = {
"ael" = {
c(
cv = qchisq(1 - alpha, df = p),
pval = pchisq(statistic, df = p, lower.tail = FALSE)
Expand All @@ -32,6 +32,12 @@ calibrate <- function(calibrate, alpha, statistic, p, par, object, control) {
control@maxit_l, control@tol_l, control@th, getWeights(object)
)
},
"chisq" = {
c(
cv = qchisq(1 - alpha, df = p),
pval = pchisq(statistic, df = p, lower.tail = FALSE)
)
},
"f" = {
n <- nrow(getData(object))
c(
Expand Down
9 changes: 8 additions & 1 deletion R/el_control.R
Original file line number Diff line number Diff line change
Expand Up @@ -41,6 +41,9 @@
#' non-overlapping but reproducible sequence of random numbers. The
#' Xoshiro256+ pseudo-random number generator is used internally to work with
#' OpenMP.
#' @param an A single numeric representing the scaling factor for adjusted
#' empirical likelihood calibration. It only applies to [elt()] when
#' `calibrate` is set to `"ael"`. Defaults to `NULL`.
#' @param b A single integer for the number of bootstrap replicates. It only
#' applies to [elt()] when `calibrate` is set to `"boot"`. Defaults to
#' `10000`.
Expand All @@ -61,6 +64,7 @@ el_control <- function(maxit = 200L,
keep_data = TRUE,
nthreads,
seed = NULL,
an = NULL,
b = 10000L,
m = 1e+06L) {
maxit <- assert_int(maxit, lower = 1L, coerce = TRUE)
Expand Down Expand Up @@ -92,11 +96,14 @@ el_control <- function(maxit = 200L,
if (isFALSE(is.null(seed))) {
seed <- assert_int(seed, coerce = TRUE)
}
if (isFALSE(is.null(an))) {
an <- assert_number(an, lower = .Machine$double.eps, finite = TRUE)
}
b <- assert_int(b, lower = 1L, coerce = TRUE)
m <- assert_int(m, lower = 1L, coerce = TRUE)
new("ControlEL",
maxit = maxit, maxit_l = maxit_l, tol = tol, tol_l = tol_l, step = step,
th = th, verbose = verbose, keep_data = keep_data, nthreads = nthreads,
seed = seed, b = b, m = m
seed = seed, an = an, b = b, m = m
)
}
66 changes: 54 additions & 12 deletions R/elt-methods.R
Original file line number Diff line number Diff line change
Expand Up @@ -18,15 +18,16 @@ setMethod("elt", "EL", function(object,
onames <- names(logProb(object))
pnames <- names(getOptim(object)$par)
h <- validate_hypothesis(rhs, lhs, getNumPar(object), pnames)
alpha <- assert_number(alpha, lower = 0, upper = 1, finite = TRUE)
calibrate <- validate_calibrate(calibrate)
assert_number(alpha, lower = 0, upper = 1, finite = TRUE)
calibrate <- assert_choice(tolower(calibrate), c("ael", "boot", "chisq", "f"))
method <- getMethodEL(object)
maxit <- control@maxit
maxit_l <- control@maxit_l
tol <- control@tol
tol_l <- control@tol_l
step <- control@step
th <- control@th
an <- control@an
w <- getWeights(object)
if (is.null(lhs)) {
stopifnot(
Expand All @@ -35,6 +36,17 @@ setMethod("elt", "EL", function(object,
)
par <- h$r
out <- compute_EL(method, par, getData(object), maxit_l, tol_l, th, w)
if (isTRUE(calibrate == "ael")) {
if (is.null(an)) {
n <- nobs(object)
an <- log(n) / 2
}
g <- out$g
pg <- -an * colMeans(g)
g <- rbind(g, pg)
out <- compute_generic_EL(g, maxit_l, tol_l, th, w)
out$optim <- c(par = list(par), out$optim)
}
optim <- validate_optim(out$optim)
names(optim$par) <- pnames
optim$cstr <- FALSE
Expand All @@ -50,9 +62,11 @@ setMethod("elt", "EL", function(object,
}
# Proceed with chi-square calibration for non-NULL `lhs`
stopifnot(
"Bootstrap calibration is applicable only when `lhs` is NULL." =
"AEL calibration is applicable only when `lhs` is `NULL`." =
(calibrate != "ael"),
"Bootstrap calibration is applicable only when `lhs` is `NULL`." =
(calibrate != "boot"),
"F calibration is applicable only when `lhs` is NULL." =
"F calibration is applicable only when `lhs` is `NULL`." =
(calibrate != "f")
)
out <- test_hypothesis(
Expand Down Expand Up @@ -95,15 +109,16 @@ setMethod("elt", "QGLM", function(object,
nm <- names(getOptim(object)$par)
pnames <- nm[-getNumPar(object)]
h <- validate_hypothesis(rhs, lhs, getNumPar(object) - 1L, pnames)
alpha <- assert_number(alpha, lower = 0, upper = 1, finite = TRUE)
calibrate <- validate_calibrate(calibrate)
assert_number(alpha, lower = 0, upper = 1, finite = TRUE)
assert_choice(tolower(calibrate), c("ael", "boot", "chisq", "f"))
method <- getMethodEL(object)
maxit <- control@maxit
maxit_l <- control@maxit_l
tol <- control@tol
tol_l <- control@tol_l
step <- control@step
th <- control@th
an <- control@an
w <- getWeights(object)
if (is.null(lhs)) {
stopifnot(
Expand All @@ -115,6 +130,17 @@ setMethod("elt", "QGLM", function(object,
par <- c(h$r, object@dispersion)
}
out <- compute_EL(method, par, getData(object), maxit_l, tol_l, th, w)
if (isTRUE(calibrate == "ael")) {
if (is.null(an)) {
n <- nobs(object)
an <- log(n) / 2
}
g <- out$g
pg <- -an * colMeans(g)
g <- rbind(g, pg)
out <- compute_generic_EL(g, maxit_l, tol_l, th, w)
out$optim <- c(par = list(par), out$optim)
}
optim <- validate_optim(out$optim)
names(optim$par) <- nm
optim$cstr <- TRUE
Expand All @@ -124,12 +150,14 @@ setMethod("elt", "QGLM", function(object,
return(new("ELT",
optim = optim, logp = setNames(out$logp, onames), logl = out$logl,
loglr = out$loglr, statistic = out$statistic, df = length(par),
pval = unname(cal["pval"]), cv = unname(cal["cv"]), rhs = h$r, lhs = h$l,
pval = unname(cal["pval"]), cv = unname(cal["cv"]), rhs = par, lhs = h$l,
alpha = alpha, calibrate = calibrate, control = control
))
}
stopifnot(
"Bootstrap calibration is applicable only when `lhs` is NULL." =
"AEL calibration is applicable only when `lhs` is `NULL`." =
(calibrate != "ael"),
"Bootstrap calibration is applicable only when `lhs` is `NULL`." =
(calibrate != "boot"),
"F calibration is applicable only to the mean." = (calibrate != "f")
)
Expand Down Expand Up @@ -172,11 +200,12 @@ setMethod("elt", "SD", function(object,
onames <- names(logProb(object))
pnames <- names(getOptim(object)$par)
h <- validate_hypothesis(rhs, lhs, getNumPar(object), pnames)
alpha <- assert_number(alpha, lower = 0, upper = 1, finite = TRUE)
calibrate <- validate_calibrate(calibrate)
assert_number(alpha, lower = 0, upper = 1, finite = TRUE)
assert_choice(tolower(calibrate), c("ael", "boot", "chisq", "f"))
maxit_l <- control@maxit_l
tol_l <- control@tol_l
th <- control@th
an <- control@an
w <- getWeights(object)
if (is.null(lhs)) {
par <- h$r
Expand All @@ -188,6 +217,17 @@ setMethod("elt", "SD", function(object,
(par >= .Machine$double.eps)
)
out <- compute_EL("sd", par, getData(object), maxit_l, tol_l, th, w)
if (isTRUE(calibrate == "ael")) {
if (is.null(an)) {
n <- nobs(object)
an <- log(n) / 2
}
g <- out$g
pg <- -an * colMeans(g)
g <- rbind(g, pg)
out <- compute_generic_EL(g, maxit_l, tol_l, th, w)
out$optim <- c(par = list(par), out$optim)
}
optim <- validate_optim(out$optim)
names(optim$par) <- pnames
optim$cstr <- FALSE
Expand All @@ -200,6 +240,8 @@ setMethod("elt", "SD", function(object,
))
}
stopifnot(
"AEL calibration is applicable only when `lhs` is `NULL`." =
(calibrate != "ael"),
"Bootstrap calibration is applicable only when `lhs` is `NULL`." =
(calibrate != "boot"),
"F calibration is applicable only when `lhs` is `NULL`." =
Expand Down Expand Up @@ -233,8 +275,8 @@ setMethod("elt", "missing", function(object,
alpha = 0.05,
calibrate = "chisq",
control = NULL) {
alpha <- assert_number(alpha, lower = 0, upper = 1, finite = TRUE)
calibrate <- validate_calibrate(calibrate)
assert_number(alpha, lower = 0, upper = 1, finite = TRUE)
assert_choice(tolower(calibrate), c("ael", "boot", "chisq", "f"))
if (is.null(control)) {
control <- el_control()
} else {
Expand Down
1 change: 1 addition & 0 deletions R/melt-package.R
Original file line number Diff line number Diff line change
@@ -1,4 +1,5 @@
## usethis namespace: start
#' @importFrom checkmate assert_choice
#' @importFrom checkmate assert_class
#' @importFrom checkmate assert_int
#' @importFrom checkmate assert_logical
Expand Down
10 changes: 6 additions & 4 deletions R/print-methods.R
Original file line number Diff line number Diff line change
Expand Up @@ -57,8 +57,9 @@ setMethod("print", "ELT", function(x,
x@rhs, x@lhs, names(getOptim(x)$par)[seq_len(ncol(x@lhs))], digits
), sep = "\n")
method <- switch(x@calibrate,
"chisq" = "Chi-square",
"ael" = "Adjusted EL",
"boot" = "Bootstrap",
"chisq" = "Chi-square",
"f" = "F"
)
cat(paste0(
Expand Down Expand Up @@ -209,8 +210,9 @@ setMethod(
x@rhs, x@lhs, names(getOptim(x)$par)[seq_len(ncol(x@lhs))], digits
), sep = "\n")
method <- switch(x@calibrate,
"chisq" = "Chi-square",
"ael" = "Adjusted EL",
"boot" = "Bootstrap",
"chisq" = "Chi-square",
"f" = "F"
)
cat(paste0(
Expand Down Expand Up @@ -265,7 +267,7 @@ setMethod(
cat("\nLagrange multipliers:\n")
print.default(getOptim(x)$lambda, digits = digits, ...)
cat("\nMaximum EL estimates:\n")
print.default(coef(x)[,1L], digits = digits, ...)
print.default(coef(x)[, 1L], digits = digits, ...)
cat(paste(
"\nlogL:", format.default(logL(x), digits = digits, ...),
", logLR:", format.default(logLR(x), digits = digits, ...)
Expand Down Expand Up @@ -337,7 +339,7 @@ setMethod(
cat("\nLagrange multipliers:\n")
print.default(getOptim(x)$lambda, digits = digits, ...)
cat("\nMaximum EL estimates:\n")
print.default(coef(x)[,1L], digits = digits, ...)
print.default(coef(x)[, 1L], digits = digits, ...)
cat(paste(
"\nlogL:", format.default(logL(x), digits = digits, ...),
", logLR:", format.default(logLR(x), digits = digits, ...)
Expand Down
Loading