Skip to content

Commit

Permalink
Added 'collapse' argument to cumul() functions.
Browse files Browse the repository at this point in the history
  • Loading branch information
@wviechtb
wviechtb committed Jul 9, 2024
1 parent fa074d5 commit 7ec4567
Show file tree
Hide file tree
Showing 312 changed files with 5,071 additions and 3,807 deletions.
5 changes: 3 additions & 2 deletions DESCRIPTION
View file
Original file line number Diff line number Diff line change
@@ -1,6 +1,6 @@
Package: metafor
Version: 4.7-23
Date: 2024-07-02
Version: 4.7-24
Date: 2024-07-09
Title: Meta-Analysis Package for R
Authors@R: person(given = "Wolfgang", family = "Viechtbauer", role = c("aut","cre"), email = "wvb@metafor-project.org", comment = c(ORCID = "0000-0003-3463-4063"))
Depends: R (>= 4.0.0), methods, Matrix, metadat, numDeriv
Expand All @@ -12,5 +12,6 @@ ByteCompile: TRUE
Encoding: UTF-8
RdMacros: mathjaxr
VignetteBuilder: R.rsp
BuildManual: TRUE
URL: https://www.metafor-project.org https://github.com/wviechtb/metafor https://wviechtb.github.io/metafor/ https://www.wvbauer.com
BugReports: https://github.com/wviechtb/metafor/issues
4 changes: 3 additions & 1 deletion NEWS.md
View file
Original file line number Diff line number Diff line change
@@ -1,4 +1,6 @@
# metafor 4.7-23 (2024-07-02)
# metafor 4.7-24 (2024-07-09)

- added `collapse` argument to the various `cumul()` functions

- the `predict.rma()` and `predict.rma.ls()` functions now also accept a matrix as input that includes a column for the intercept term (in which case the `intercept` argument is ignored)

Expand Down
171 changes: 104 additions & 67 deletions R/cumul.rma.mh.r
View file
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
cumul.rma.mh <- function(x, order, digits, transf, targs, progbar=FALSE, ...) {
cumul.rma.mh <- function(x, order, digits, transf, targs, collapse=FALSE, progbar=FALSE, ...) {

mstyle <- .get.mstyle()

Expand All @@ -9,6 +9,8 @@ cumul.rma.mh <- function(x, order, digits, transf, targs, progbar=FALSE, ...) {
if (!is.element(na.act, c("na.omit", "na.exclude", "na.fail", "na.pass")))
stop(mstyle$stop("Unknown 'na.action' specified under options()."))

if (na.act == "na.fail" && any(!x$not.na))
stop(mstyle$stop("Missing values in data."))

if (missing(digits)) {
digits <- .get.digits(xdigits=x$digits, dmiss=TRUE)
Expand All @@ -34,85 +36,122 @@ cumul.rma.mh <- function(x, order, digits, transf, targs, progbar=FALSE, ...) {
#########################################################################

if (grepl("^order\\(", deparse1(substitute(order))))
warning(mstyle$warning("Use of order() in 'order' argument is probably erroneous."), call.=FALSE)
warning(mstyle$warning("Use of order() in the 'order' argument is probably erroneous."), call.=FALSE)

if (missing(order)) {
order <- seq_len(x$k.all)

orvar <- seq_len(x$k.all)
collapse <- FALSE

} else {

mf <- match.call()
order <- .getx("order", mf=mf, data=x$data)
}
orvar <- .getx("order", mf=mf, data=x$data)

if (length(orvar) != x$k.all)
stop(mstyle$stop(paste0("Length of the 'order' argument (", length(orvar), ") does not correspond to the size of the original dataset (", x$k.all, ").")))

if (length(order) != x$k.all)
stop(mstyle$stop(paste0("Length of the 'order' argument (", length(order), ") does not correspond to the size of the original dataset (", x$k.all, ").")))
}

### note: order variable must be of the same length as the original dataset
### so we have to apply the same subsetting (if necessary)
### as was done during model fitting

order <- .getsubset(order, x$subset)

order <- order(order, decreasing=decreasing)

ai.f <- x$outdat.f$ai[order]
bi.f <- x$outdat.f$bi[order]
ci.f <- x$outdat.f$ci[order]
di.f <- x$outdat.f$di[order]
x1i.f <- x$outdat.f$x1i[order]
x2i.f <- x$outdat.f$x2i[order]
t1i.f <- x$outdat.f$t1i[order]
t2i.f <- x$outdat.f$t2i[order]
yi.f <- x$yi.f[order]
vi.f <- x$vi.f[order]
### so apply the same subsetting as was done during the model fitting

orvar <- .getsubset(orvar, x$subset)

### order data by the order variable (NAs in order variable are dropped)

order <- base::order(orvar, decreasing=decreasing, na.last=NA)

ai <- x$outdat.f$ai[order]
bi <- x$outdat.f$bi[order]
ci <- x$outdat.f$ci[order]
di <- x$outdat.f$di[order]
x1i <- x$outdat.f$x1i[order]
x2i <- x$outdat.f$x2i[order]
t1i <- x$outdat.f$t1i[order]
t2i <- x$outdat.f$t2i[order]
yi <- x$yi.f[order]
vi <- x$vi.f[order]
not.na <- x$not.na[order]
slab <- x$slab[order]
ids <- x$ids[order]
orvar <- orvar[order]

if (inherits(x$data, "environment")) {
data <- NULL
} else {
data <- x$data[order,]
}

beta <- rep(NA_real_, x$k.f)
se <- rep(NA_real_, x$k.f)
zval <- rep(NA_real_, x$k.f)
pval <- rep(NA_real_, x$k.f)
ci.lb <- rep(NA_real_, x$k.f)
ci.ub <- rep(NA_real_, x$k.f)
QE <- rep(NA_real_, x$k.f)
QEp <- rep(NA_real_, x$k.f)
I2 <- rep(NA_real_, x$k.f)
H2 <- rep(NA_real_, x$k.f)
if (collapse) {
uorvar <- unique(orvar)
} else {
uorvar <- orvar
}

### elements that need to be returned
k.o <- length(uorvar)

k <- rep(NA_integer_, k.o)
beta <- rep(NA_real_, k.o)
se <- rep(NA_real_, k.o)
zval <- rep(NA_real_, k.o)
pval <- rep(NA_real_, k.o)
ci.lb <- rep(NA_real_, k.o)
ci.ub <- rep(NA_real_, k.o)
QE <- rep(NA_real_, k.o)
QEp <- rep(NA_real_, k.o)
I2 <- rep(NA_real_, k.o)
H2 <- rep(NA_real_, k.o)
show <- rep(TRUE, k.o)

outlist <- "beta=beta, se=se, zval=zval, pval=pval, ci.lb=ci.lb, ci.ub=ci.ub, QE=QE, QEp=QEp, tau2=tau2, I2=I2, H2=H2"
### elements that need to be returned

### note: skipping NA cases
outlist <- "k=k, beta=beta, se=se, zval=zval, pval=pval, ci.lb=ci.lb, ci.ub=ci.ub, QE=QE, QEp=QEp, I2=I2, H2=H2"

if (progbar)
pbar <- pbapply::startpb(min=0, max=x$k.f)
pbar <- pbapply::startpb(min=0, max=k.o)

for (i in seq_len(x$k.f)) {
for (i in seq_len(k.o)) {

if (progbar)
pbapply::setpb(pbar, i)

if (!not.na[i])
next
if (collapse) {

if (all(!not.na[is.element(orvar, uorvar[i])])) {
if (na.act == "na.omit")
show[i] <- FALSE # if all studies to be added are !not.na, don't show (but a fit failure is still shown)
next
}

incl <- is.element(orvar, uorvar[1:i])

} else {

if (!not.na[i]) {
if (na.act == "na.omit")
show[i] <- FALSE # if study to be added is !not.na, don't show (but a fit failure is still shown)
next
}

incl <- 1:i

}

if (is.element(x$measure, c("RR","OR","RD"))) {
args <- list(ai=ai.f, bi=bi.f, ci=ci.f, di=di.f, measure=x$measure, add=x$add, to=x$to, drop00=x$drop00,
correct=x$correct, level=x$level, subset=seq_len(i), outlist=outlist)
args <- list(ai=ai, bi=bi, ci=ci, di=di, measure=x$measure, add=x$add, to=x$to, drop00=x$drop00,
correct=x$correct, level=x$level, subset=incl, outlist=outlist)
} else {
args <- list(x1i=x1i.f, x2i=x2i.f, t1i=t1i.f, t2i=t2i.f, measure=x$measure, add=x$add, to=x$to, drop00=x$drop00,
correct=x$correct, level=x$level, subset=seq_len(i), outlist=outlist)
args <- list(x1i=x1i, x2i=x2i, t1i=t1i, t2i=t2i, measure=x$measure, add=x$add, to=x$to, drop00=x$drop00,
correct=x$correct, level=x$level, subset=incl, outlist=outlist)
}

res <- try(suppressWarnings(.do.call(rma.mh, args)), silent=TRUE)

if (inherits(res, "try-error"))
next

k[i] <- res$k
beta[i] <- res$beta
se[i] <- res$se
zval[i] <- res$zval
Expand All @@ -139,12 +178,12 @@ cumul.rma.mh <- function(x, order, digits, transf, targs, progbar=FALSE, ...) {
if (is.function(transf)) {
if (is.null(targs)) {
beta <- sapply(beta, transf)
se <- rep(NA_real_, x$k.f)
se <- rep(NA_real_, k.o)
ci.lb <- sapply(ci.lb, transf)
ci.ub <- sapply(ci.ub, transf)
} else {
beta <- sapply(beta, transf, targs)
se <- rep(NA_real_, x$k.f)
se <- rep(NA_real_, k.o)
ci.lb <- sapply(ci.lb, transf, targs)
ci.ub <- sapply(ci.ub, transf, targs)
}
Expand All @@ -159,31 +198,29 @@ cumul.rma.mh <- function(x, order, digits, transf, targs, progbar=FALSE, ...) {

#########################################################################

if (na.act == "na.omit") {
out <- list(estimate=beta[not.na], se=se[not.na], zval=zval[not.na], pval=pval[not.na], ci.lb=ci.lb[not.na], ci.ub=ci.ub[not.na], Q=QE[not.na], Qp=QEp[not.na], I2=I2[not.na], H2=H2[not.na])
out$slab <- slab[not.na]
out$ids <- ids[not.na]
out$data <- data[not.na,]
}
out <- list(k=k[show], estimate=beta[show], se=se[show], zval=zval[show], pval=pval[show], ci.lb=ci.lb[show], ci.ub=ci.ub[show], Q=QE[show], Qp=QEp[show], I2=I2[show], H2=H2[show])

if (na.act == "na.exclude" || na.act == "na.pass") {
out <- list(estimate=beta, se=se, zval=zval, pval=pval, ci.lb=ci.lb, ci.ub=ci.ub, Q=QE, Qp=QEp, I2=I2, H2=H2)
out$slab <- slab
out$ids <- ids
out$data <- data
if (collapse) {
out$slab <- uorvar[show]
out$slab.null <- FALSE
} else {
out$slab <- slab[show]
out$ids <- ids[show]
out$data <- data[show,,drop=FALSE]
out$slab.null <- x$slab.null
}

if (na.act == "na.fail" && any(!x$not.na))
stop(mstyle$stop("Missing values in results."))
out$order <- uorvar[show]

out$digits <- digits
out$transf <- transf
out$slab.null <- x$slab.null
out$level <- x$level
out$measure <- x$measure
out$test <- x$test
out$digits <- digits
out$transf <- transf
out$level <- x$level
out$test <- x$test

attr(out$estimate, "measure") <- x$measure
if (!transf) {
out$measure <- x$measure
attr(out$estimate, "measure") <- x$measure
}

if (.isTRUE(ddd$time)) {
time.end <- proc.time()
Expand Down
Loading

0 comments on commit 7ec4567

Please sign in to comment.