diff --git a/R/anova.rma.r b/R/anova.rma.r index c683a7dcd..8e82232f4 100644 --- a/R/anova.rma.r +++ b/R/anova.rma.r @@ -215,8 +215,8 @@ anova.rma <- function(object, object2, btt, X, att, Z, rhs, digits, refit=FALSE, ### omnibus test of all hypotheses (only possible if 'Z' is of full rank) - QS <- NA_real_ ### need this in case QS cannot be calculated below - QSp <- NA_real_ ### need this in case QSp cannot be calculated below + QS <- NA_real_ # need this in case QS cannot be calculated below + QSp <- NA_real_ # need this in case QSp cannot be calculated below if (rankMatrix(Z) == m) { @@ -240,24 +240,24 @@ anova.rma <- function(object, object2, btt, X, att, Z, rhs, digits, refit=FALSE, hyp <- rep("", m) for (j in seq_len(m)) { - Zj <- round(Z[j,], digits[["est"]]) ### coefficients for the jth contrast - sel <- Zj != 0 ### TRUE if coefficient is != 0 - hyp[j] <- paste(paste(Zj[sel], rownames(x$alpha)[sel], sep="*"), collapse=" + ") ### coefficient*variable + coefficient*variable ... - hyp[j] <- gsub("1*", "", hyp[j], fixed=TRUE) ### turn '+1' into '+' and '-1' into '-' - hyp[j] <- gsub("+ -", "- ", hyp[j], fixed=TRUE) ### turn '+ -' into '-' + Zj <- round(Z[j,], digits[["est"]]) # coefficients for the jth contrast + sel <- Zj != 0 # TRUE if coefficient is != 0 + hyp[j] <- paste(paste(Zj[sel], rownames(x$alpha)[sel], sep="*"), collapse=" + ") # coefficient*variable + coefficient*variable ... + hyp[j] <- gsub("1*", "", hyp[j], fixed=TRUE) # turn '+1' into '+' and '-1' into '-' + hyp[j] <- gsub("+ -", "- ", hyp[j], fixed=TRUE) # turn '+ -' into '-' } if (identical(rhs, rep(0,m))) { - hyp <- paste0(hyp, " = 0") ### add '= 0' at the right + hyp <- paste0(hyp, " = 0") # add '= 0' at the right } else { if (length(unique(rhs)) == 1L) { - hyp <- paste0(hyp, " = ", round(rhs, digits=digits[["est"]])) ### add '= rhs' at the right + hyp <- paste0(hyp, " = ", round(rhs, digits=digits[["est"]])) # add '= rhs' at the right } else { - hyp <- paste0(hyp, " = ", fmtx(rhs, digits=digits[["est"]])) ### add '= rhs' at the right + hyp <- paste0(hyp, " = ", fmtx(rhs, digits=digits[["est"]])) # add '= rhs' at the right } } hyp <- data.frame(hyp, stringsAsFactors=FALSE) colnames(hyp) <- "" - rownames(hyp) <- paste0(seq_len(m), ":") ### add '1:', '2:', ... as row names + rownames(hyp) <- paste0(seq_len(m), ":") # add '1:', '2:', ... as row names res <- list(QS=QS, QSdf=QSdf, QSp=QSp, hyp=hyp, Za=Za, se=se, zval=zval, pval=pval, k=x$k, q=x$q, m=m, test=x$test, ddf=x$ddf.alpha, digits=digits, type="Wald.Za") @@ -400,24 +400,24 @@ anova.rma <- function(object, object2, btt, X, att, Z, rhs, digits, refit=FALSE, hyp <- rep("", m) for (j in seq_len(m)) { - Xj <- round(X[j,], digits[["est"]]) ### coefficients for the jth contrast - sel <- Xj != 0 ### TRUE if coefficient is != 0 - hyp[j] <- paste(paste(Xj[sel], rownames(x$beta)[sel], sep="*"), collapse=" + ") ### coefficient*variable + coefficient*variable ... - hyp[j] <- gsub("1*", "", hyp[j], fixed=TRUE) ### turn '+1' into '+' and '-1' into '-' - hyp[j] <- gsub("+ -", "- ", hyp[j], fixed=TRUE) ### turn '+ -' into '-' + Xj <- round(X[j,], digits[["est"]]) # coefficients for the jth contrast + sel <- Xj != 0 # TRUE if coefficient is != 0 + hyp[j] <- paste(paste(Xj[sel], rownames(x$beta)[sel], sep="*"), collapse=" + ") # coefficient*variable + coefficient*variable ... + hyp[j] <- gsub("1*", "", hyp[j], fixed=TRUE) # turn '+1' into '+' and '-1' into '-' + hyp[j] <- gsub("+ -", "- ", hyp[j], fixed=TRUE) # turn '+ -' into '-' } if (identical(rhs, rep(0,m))) { - hyp <- paste0(hyp, " = 0") ### add '= 0' at the right + hyp <- paste0(hyp, " = 0") # add '= 0' at the right } else { if (length(unique(rhs)) == 1L) { - hyp <- paste0(hyp, " = ", round(rhs, digits=digits[["est"]])) ### add '= rhs' at the right + hyp <- paste0(hyp, " = ", round(rhs, digits=digits[["est"]])) # add '= rhs' at the right } else { - hyp <- paste0(hyp, " = ", fmtx(rhs, digits=digits[["est"]])) ### add '= rhs' at the right + hyp <- paste0(hyp, " = ", fmtx(rhs, digits=digits[["est"]])) # add '= rhs' at the right } } hyp <- data.frame(hyp, stringsAsFactors=FALSE) colnames(hyp) <- "" - rownames(hyp) <- paste0(seq_len(m), ":") ### add '1:', '2:', ... as row names + rownames(hyp) <- paste0(seq_len(m), ":") # add '1:', '2:', ... as row names res <- list(QM=QM, QMdf=QMdf, QMp=QMp, hyp=hyp, Xb=Xb, se=se, zval=zval, pval=pval, k=x$k, p=x$p, m=m, test=x$test, ddf=ddf, digits=digits, type="Wald.Xb") diff --git a/R/bldiag.r b/R/bldiag.r index d62870286..7996ff0a3 100644 --- a/R/bldiag.r +++ b/R/bldiag.r @@ -22,7 +22,7 @@ bldiag <- function(..., order) { if (any(is00)) mlist <- mlist[!is00] - csdim <- rbind(c(0,0), apply(sapply(mlist,dim), 1, cumsum)) ### consider using rowCumsums() from matrixStats package + csdim <- rbind(c(0,0), apply(sapply(mlist,dim), 1, cumsum)) # consider using rowCumsums() from matrixStats package out <- array(0, dim=csdim[length(mlist) + 1,]) add1 <- matrix(rep(1:0, 2L), ncol=2) @@ -31,9 +31,9 @@ bldiag <- function(..., order) { indx <- apply(csdim[i:(i+1),] + add1, 2, function(x) x[1]:x[2]) - if (is.null(dim(indx))) { ### non-square matrix + if (is.null(dim(indx))) { # non-square matrix out[indx[[1]],indx[[2]]] <- mlist[[i]] - } else { ### square matrix + } else { # square matrix out[indx[,1],indx[,2]] <- mlist[[i]] } diff --git a/R/cooks.distance.rma.mv.r b/R/cooks.distance.rma.mv.r index a6e6ba491..8cc15aa8c 100644 --- a/R/cooks.distance.rma.mv.r +++ b/R/cooks.distance.rma.mv.r @@ -4,7 +4,7 @@ cooks.distance.rma.mv <- function(model, progbar=FALSE, cluster, reestimate=TRUE .chkclass(class(model), must="rma.mv") - #if (inherits(model, "robust.rma")) ### can compute Cook's distance also for 'robust.rma' objects + #if (inherits(model, "robust.rma")) # can compute Cook's distance also for 'robust.rma' objects # stop(mstyle$stop("Method not available for objects of class \"robust.rma\".")) na.act <- getOption("na.action") diff --git a/R/cumul.rma.mh.r b/R/cumul.rma.mh.r index c0fe1b029..1b58e6350 100644 --- a/R/cumul.rma.mh.r +++ b/R/cumul.rma.mh.r @@ -133,7 +133,7 @@ cumul.rma.mh <- function(x, order, digits, transf, targs, progbar=FALSE, ...) { ### if requested, apply transformation function - if (.isTRUE(transf) && is.element(x$measure, c("OR","RR","IRR"))) ### if transf=TRUE, apply exp transformation to ORs, RRs, and IRRs + if (.isTRUE(transf) && is.element(x$measure, c("OR","RR","IRR"))) # if transf=TRUE, apply exp transformation to ORs, RRs, and IRRs transf <- exp if (is.function(transf)) { diff --git a/R/cumul.rma.peto.r b/R/cumul.rma.peto.r index 9049429ec..11b1501da 100644 --- a/R/cumul.rma.peto.r +++ b/R/cumul.rma.peto.r @@ -123,7 +123,7 @@ cumul.rma.peto <- function(x, order, digits, transf, targs, progbar=FALSE, ...) ### if requested, apply transformation function - if (.isTRUE(transf)) ### if transf=TRUE, apply exp transformation to ORs + if (.isTRUE(transf)) # if transf=TRUE, apply exp transformation to ORs transf <- exp if (is.function(transf)) { diff --git a/R/forest.rma.r b/R/forest.rma.r index 08d0d1fa3..356403a1c 100644 --- a/R/forest.rma.r +++ b/R/forest.rma.r @@ -920,7 +920,7 @@ lty, fonts, cex, cex.lab, cex.axis, ...) { if (length(addpred) == 1L) addpred <- c(addpred, addpred) temp <- predict(x, level=level, tau2.levels=addpred[1], gamma2.levels=addpred[2], pi.type=pi.type) - addpred <- TRUE ### set addpred to TRUE, so if (!is.element(x$method, c("FE","EE","CE")) && addpred) further below works + addpred <- TRUE # set addpred to TRUE, so if (!is.element(x$method, c("FE","EE","CE")) && addpred) further below works } else { if (addpred) { ### here addpred=TRUE, but user has not specified the level, so throw an error diff --git a/R/hatvalues.rma.mv.r b/R/hatvalues.rma.mv.r index b8bfc44bc..c31599b87 100644 --- a/R/hatvalues.rma.mv.r +++ b/R/hatvalues.rma.mv.r @@ -19,7 +19,7 @@ hatvalues.rma.mv <- function(model, type="diagonal", ...) { W <- chol2inv(chol(x$M)) stXWX <- chol2inv(chol(as.matrix(t(x$X) %*% W %*% x$X))) H <- as.matrix(x$X %*% stXWX %*% crossprod(x$X,W)) - #H <- as.matrix(x$X %*% x$vb %*% crossprod(x$X,W)) ### x$vb may have been changed through robust() + #H <- as.matrix(x$X %*% x$vb %*% crossprod(x$X,W)) # x$vb may have been changed through robust() } else { A <- x$W stXAX <- chol2inv(chol(as.matrix(t(x$X) %*% A %*% x$X))) @@ -33,7 +33,7 @@ hatvalues.rma.mv <- function(model, type="diagonal", ...) { hii <- rep(NA_real_, x$k.f) hii[x$not.na] <- as.vector(diag(H)) - hii[hii > 1 - 10 * .Machine$double.eps] <- 1 ### as in lm.influence() + hii[hii > 1 - 10 * .Machine$double.eps] <- 1 # as in lm.influence() names(hii) <- x$slab if (na.act == "na.omit") diff --git a/R/hatvalues.rma.uni.r b/R/hatvalues.rma.uni.r index 3ef243a67..4206cdb7c 100644 --- a/R/hatvalues.rma.uni.r +++ b/R/hatvalues.rma.uni.r @@ -20,7 +20,7 @@ hatvalues.rma.uni <- function(model, type="diagonal", ...) { W <- diag(1/(x$vi + x$tau2), nrow=x$k, ncol=x$k) stXWX <- .invcalc(X=x$X, W=W, k=x$k) H <- x$X %*% stXWX %*% crossprod(x$X,W) - #H <- x$X %*% (x$vb / x$s2w) %*% crossprod(x$X,W) ### x$vb may be changed through robust() (and when test="knha") + #H <- x$X %*% (x$vb / x$s2w) %*% crossprod(x$X,W) # x$vb may be changed through robust() (and when test="knha") } else { A <- diag(x$weights, nrow=x$k, ncol=x$k) stXAX <- .invcalc(X=x$X, W=A, k=x$k) @@ -37,7 +37,7 @@ hatvalues.rma.uni <- function(model, type="diagonal", ...) { hii <- rep(NA_real_, x$k.f) hii[x$not.na] <- diag(H) - hii[hii > 1 - 10 * .Machine$double.eps] <- 1 ### as in lm.influence() + hii[hii > 1 - 10 * .Machine$double.eps] <- 1 # as in lm.influence() names(hii) <- x$slab if (na.act == "na.omit") diff --git a/R/hc.rma.uni.r b/R/hc.rma.uni.r index cd347cded..9b6c23e38 100644 --- a/R/hc.rma.uni.r +++ b/R/hc.rma.uni.r @@ -48,7 +48,7 @@ hc.rma.uni <- function(object, digits, transf, targs, control, ...) { ### original code by Henmi & Copas (2012), modified by Michael Dewey, small adjustments ### for consistency with other functions in the metafor package by Wolfgang Viechtbauer - wi <- 1/vi ### fixed effects weights + wi <- 1/vi # fixed effects weights W1 <- sum(wi) W2 <- sum(wi^2) / W1 diff --git a/R/labbe.rma.r b/R/labbe.rma.r index 10d533ed9..3c7a895d6 100644 --- a/R/labbe.rma.r +++ b/R/labbe.rma.r @@ -137,7 +137,7 @@ add=x$add, to=x$to, transf, targs, pch=21, psize, plim=c(0.5,3.5), col, bg, grid t1i[yi.is.na] <- NA_real_ t2i[yi.is.na] <- NA_real_ - options(na.action = "na.pass") ### to make sure dat.t and dat.c are of the same length + options(na.action = "na.pass") # to make sure dat.t and dat.c are of the same length measure <- switch(x$measure, "RR"="PLN", "OR"="PLO", "RD"="PR", "AS"="PAS", "IRR"="IRLN", "IRD"="IR", "IRSD"="IRS") diff --git a/R/leave1out.rma.mh.r b/R/leave1out.rma.mh.r index a3d479f36..f69fb05e4 100644 --- a/R/leave1out.rma.mh.r +++ b/R/leave1out.rma.mh.r @@ -100,7 +100,7 @@ leave1out.rma.mh <- function(x, digits, transf, targs, progbar=FALSE, ...) { ### if requested, apply transformation function - if (.isTRUE(transf) && is.element(x$measure, c("OR","RR","IRR"))) ### if transf=TRUE, apply exp transformation to ORs, RRs, and IRRs + if (.isTRUE(transf) && is.element(x$measure, c("OR","RR","IRR"))) # if transf=TRUE, apply exp transformation to ORs, RRs, and IRRs transf <- exp if (is.function(transf)) { diff --git a/R/leave1out.rma.peto.r b/R/leave1out.rma.peto.r index d1f3991f7..4616b6b1d 100644 --- a/R/leave1out.rma.peto.r +++ b/R/leave1out.rma.peto.r @@ -92,7 +92,7 @@ leave1out.rma.peto <- function(x, digits, transf, targs, progbar=FALSE, ...) { ### if requested, apply transformation function - if (.isTRUE(transf)) ### if transf=TRUE, apply exp transformation to ORs + if (.isTRUE(transf)) # if transf=TRUE, apply exp transformation to ORs transf <- exp if (is.function(transf)) { diff --git a/R/llplot.r b/R/llplot.r index 15990c0d3..13d9e080f 100644 --- a/R/llplot.r +++ b/R/llplot.r @@ -96,7 +96,7 @@ lty, lwd, col, level=99.99, refline=0, ...) { if (!.equal.length(yi, vi)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(yi) ### number of outcomes before subsetting + k <- length(yi) # number of outcomes before subsetting ### subsetting @@ -146,7 +146,7 @@ lty, lwd, col, level=99.99, refline=0, ...) { if (any(c(n1i < 0, n2i < 0), na.rm=TRUE)) stop(mstyle$stop("One or more group sizes are negative.")) - k <- length(ai) ### number of outcomes before subsetting + k <- length(ai) # number of outcomes before subsetting ### note studies that have at least one zero cell @@ -179,8 +179,8 @@ lty, lwd, col, level=99.99, refline=0, ...) { dat <- .do.call(escalc, measure="OR", ai=ai, bi=bi, ci=ci, di=di, drop00=drop00, onlyo1=onlyo1, addyi=addyi, addvi=addvi) - yi <- dat$yi ### one or more yi/vi pairs may be NA/NA - vi <- dat$vi ### one or more yi/vi pairs may be NA/NA + yi <- dat$yi # one or more yi/vi pairs may be NA/NA + vi <- dat$vi # one or more yi/vi pairs may be NA/NA } diff --git a/R/methods.list.rma.r b/R/methods.list.rma.r index fc528c247..6e911f813 100644 --- a/R/methods.list.rma.r +++ b/R/methods.list.rma.r @@ -19,7 +19,7 @@ if (length(out[[1]]) == 0L) return(NULL) - #out <- out[j] ### this causes all kinds of problems, so left out for now (TODO: check if this is really a problem) + #out <- out[j] # this causes all kinds of problems, so left out for now (TODO: check if this is really a problem) out$slab <- x$slab[i] diff --git a/R/misc.func.hidden.escalc.r b/R/misc.func.hidden.escalc.r index 255e4f5be..208cb46ac 100644 --- a/R/misc.func.hidden.escalc.r +++ b/R/misc.func.hidden.escalc.r @@ -256,8 +256,8 @@ } #integrate(function(x) .dcor(x, n=5, rho=.8), lower=-1, upper=1) -#integrate(function(x) x*.dcor(x, n=5, rho=.8), lower=-1, upper=1) ### should not be rho due to bias! -#integrate(function(x) x*.Fcalc(1/2, 1/2, (5-2)/2, 1-x^2)*.dcor(x, n=5, rho=.8), lower=-1, upper=1) ### should be ~rho +#integrate(function(x) x*.dcor(x, n=5, rho=.8), lower=-1, upper=1) # should not be rho due to bias! +#integrate(function(x) x*.Fcalc(1/2, 1/2, (5-2)/2, 1-x^2)*.dcor(x, n=5, rho=.8), lower=-1, upper=1) # should be ~rho ### pdf of ZCOR diff --git a/R/misc.func.hidden.glmm.r b/R/misc.func.hidden.glmm.r index d8eaf9d8c..5ee8daaed 100644 --- a/R/misc.func.hidden.glmm.r +++ b/R/misc.func.hidden.glmm.r @@ -4,7 +4,7 @@ ### Liao, J. G. & Rosen, O. (2001). Fast and stable algorithms for computing and sampling from the ### noncentral hypergeometric distribution. The American Statistician, 55, 366-369. -.dnoncenhypergeom <- function (x=NA_real_, n1, n2, m1, psi) { ### x=ai, n1=ai+bi, n2=ci+di, m1=ai+ci, psi=ORi +.dnoncenhypergeom <- function (x=NA_real_, n1, n2, m1, psi) { # x=ai, n1=ai+bi, n2=ci+di, m1=ai+ci, psi=ORi mstyle <- .get.mstyle() @@ -114,8 +114,8 @@ p <- ncol(X.fit) k <- length(ai) - beta <- parms[seq_len(p)] ### first p elemenets in parms are the model coefficients - tau2 <- ifelse(random, exp(parms[p+1]), 0) ### next value is tau^2 -- optimize over exp(tau^2) value or hold at 0 if random=FALSE + beta <- parms[seq_len(p)] # first p elemenets in parms are the model coefficients + tau2 <- ifelse(random, exp(parms[p+1]), 0) # next value is tau^2 -- optimize over exp(tau^2) value or hold at 0 if random=FALSE mu.i <- X.fit %*% cbind(beta) lli <- rep(NA_real_, k) diff --git a/R/misc.func.hidden.mv.r b/R/misc.func.hidden.mv.r index dc3b42da5..fb7cfbe29 100644 --- a/R/misc.func.hidden.mv.r +++ b/R/misc.func.hidden.mv.r @@ -74,7 +74,7 @@ ### get variables names in mf.g - g.names <- names(mf.g) ### names for inner and outer factors/variables + g.names <- names(mf.g) # names for inner and outer factors/variables ### check that inner variable is a factor (or character variable) for structures that require this (no longer required) @@ -107,7 +107,7 @@ ### get number of levels of each variable in mf.g (vector with two values, for the inner and outer factor) - #g.nlevels <- c(nlevels(mf.g[[1]]), nlevels(mf.g[[2]])) ### works only for factors + #g.nlevels <- c(nlevels(mf.g[[1]]), nlevels(mf.g[[2]])) # works only for factors if (is.element(struct, c("SPEXP","SPGAU","SPLIN","SPRAT","SPSPH","PHYBM","PHYPL","PHYPD","GEN","GDIAG"))) { g.nlevels <- c(length(unique(apply(mf.g[-nvars], 1, paste, collapse=" + "))), length(unique(mf.g[[nvars]]))) } else { @@ -116,7 +116,7 @@ ### get levels of each variable in mf.g - #g.levels <- list(levels(mf.g[[1]]), levels(mf.g[[2]])) ### works only for factors + #g.levels <- list(levels(mf.g[[1]]), levels(mf.g[[2]])) # works only for factors if (is.element(struct, c("SPEXP","SPGAU","SPLIN","SPRAT","SPSPH","PHYBM","PHYPL","PHYPD","GEN","GDIAG"))) { g.levels <- list(sort(unique(apply(mf.g[-nvars], 1, paste, collapse=" + "))), sort(unique((mf.g[[nvars]])))) } else { @@ -279,7 +279,7 @@ ### redo: get number of levels of each variable in mf.g (vector with two values, for the inner and outer factor) - #g.nlevels <- c(nlevels(mf.g[[1]]), nlevels(mf.g[[2]])) ### works only for factors + #g.nlevels <- c(nlevels(mf.g[[1]]), nlevels(mf.g[[2]])) # works only for factors if (is.element(struct, c("SPEXP","SPGAU","SPLIN","SPRAT","SPSPH","PHYBM","PHYPL","PHYPD","GEN","GDIAG"))) { g.nlevels <- c(length(unique(apply(mf.g[-nvars], 1, paste, collapse=" + "))), length(unique(mf.g[[nvars]]))) } else { @@ -288,7 +288,7 @@ ### redo: get levels of each variable in mf.g - #g.levels <- list(levels(mf.g[[1]]), levels(mf.g[[2]])) ### works only for factors + #g.levels <- list(levels(mf.g[[1]]), levels(mf.g[[2]])) # works only for factors if (is.element(struct, c("SPEXP","SPGAU","SPLIN","SPRAT","SPSPH","PHYBM","PHYPL","PHYPD","GEN","GDIAG"))) { g.levels <- list(sort(unique(apply(mf.g[-nvars], 1, paste, collapse=" + "))), sort(unique((mf.g[[nvars]])))) } else { @@ -613,7 +613,7 @@ if (any(g.levels.r) && struct == "UNR") { G[g.levels.r,] <- 0 G[,g.levels.r] <- 0 - diag(G) <- tau2 ### don't really need this + diag(G) <- tau2 # don't really need this rho <- G[lower.tri(G)] warning(mstyle$warning(paste0("Fixed ", ifelse(isG, 'rho', 'phi'), " value(s) to 0 for removed level(s).")), call.=FALSE) } diff --git a/R/misc.func.hidden.r b/R/misc.func.hidden.r index fcad079d5..ac704505c 100644 --- a/R/misc.func.hidden.r +++ b/R/misc.func.hidden.r @@ -8,14 +8,14 @@ if (missing(btt) || is.null(btt)) { - if (p > 1L) { ### if the model matrix has more than one column + if (p > 1L) { # if the model matrix has more than one column if (int.incl) { - btt <- seq.int(from=2, to=p) ### and the model has an intercept term, test all coefficients except the intercept + btt <- seq.int(from=2, to=p) # and the model has an intercept term, test all coefficients except the intercept } else { - btt <- seq_len(p) ### and the model does not have an intercept term, test all coefficients + btt <- seq_len(p) # and the model does not have an intercept term, test all coefficients } } else { - btt <- 1L ### if the model matrix has a single column, test that single coefficient + btt <- 1L # if the model matrix has a single column, test that single coefficient } } else { diff --git a/R/print.list.rma.r b/R/print.list.rma.r index aed9a4000..674e3d9fe 100644 --- a/R/print.list.rma.r +++ b/R/print.list.rma.r @@ -56,7 +56,7 @@ print.list.rma <- function(x, digits=x$digits, ...) { sav <- out[,seq_len(min.pos-1)] for (i in seq_len(min.pos-1)) { - if (inherits(out[,i], c("integer","logical","factor","character"))) { ### do not apply formating to these classes + if (inherits(out[,i], c("integer","logical","factor","character"))) { # do not apply formating to these classes out[,i] <- out[,i] } else { if (names(out)[i] %in% c("pred", "resid")) diff --git a/R/print.rma.mv.r b/R/print.rma.mv.r index 9a040d7ac..68eb3e30e 100644 --- a/R/print.rma.mv.r +++ b/R/print.rma.mv.r @@ -180,7 +180,7 @@ print.rma.mv <- function(x, digits, showfit=FALSE, signif.stars=getOption("show. diag(G.info) <- "-" vc <- cbind(G, "", G.info) - colnames(vc) <- c(paste0("rho.", abbreviate(x$g.levels.f[[1]])), "", abbreviate(x$g.levels.f[[1]])) ### FIXME: x$g.levels.f[[1]] may be numeric, in which case a wrapping 'header' is not recognized + colnames(vc) <- c(paste0("rho.", abbreviate(x$g.levels.f[[1]])), "", abbreviate(x$g.levels.f[[1]])) # FIXME: x$g.levels.f[[1]] may be numeric, in which case a wrapping 'header' is not recognized rownames(vc) <- x$g.levels.f[[1]] tmp <- capture.output(print(vc, quote=FALSE, right=right, print.gap=2)) .print.table(tmp, mstyle) @@ -311,7 +311,7 @@ print.rma.mv <- function(x, digits, showfit=FALSE, signif.stars=getOption("show. diag(H.info) <- "-" vc <- cbind(H, "", H.info) - colnames(vc) <- c(paste0("phi.", abbreviate(x$h.levels.f[[1]])), "", abbreviate(x$h.levels.f[[1]])) ### FIXME: x$h.levels.f[[1]] may be numeric, in which case a wrapping 'header' is not recognized + colnames(vc) <- c(paste0("phi.", abbreviate(x$h.levels.f[[1]])), "", abbreviate(x$h.levels.f[[1]])) # FIXME: x$h.levels.f[[1]] may be numeric, in which case a wrapping 'header' is not recognized rownames(vc) <- x$h.levels.f[[1]] tmp <- capture.output(print(vc, quote=FALSE, right=right, print.gap=2)) .print.table(tmp, mstyle) diff --git a/R/qqnorm.rma.uni.r b/R/qqnorm.rma.uni.r index 923b055c0..082a71fc2 100644 --- a/R/qqnorm.rma.uni.r +++ b/R/qqnorm.rma.uni.r @@ -35,7 +35,7 @@ label=FALSE, offset=0.3, pos=13, lty, ...) { bg <- .coladj(par("bg","fg"), dark=0.35, light=-0.35) if (missing(lty)) { - lty <- c("solid", "dotted") ### 1st value = diagonal line, 2nd value = pseudo confidence envelope + lty <- c("solid", "dotted") # 1st value = diagonal line, 2nd value = pseudo confidence envelope } else { if (length(lty) == 1L) lty <- c(lty, lty) @@ -94,11 +94,11 @@ label=FALSE, offset=0.3, pos=13, lty, ...) { ei <- ImH %*% dat ei <- apply(ei, 2, sort) if (bonferroni) { - lb <- apply(ei, 1, quantile, (level/2)/x$k) ### consider using rowQuantiles() from matrixStats package - ub <- apply(ei, 1, quantile, 1-(level/2)/x$k) ### consider using rowQuantiles() from matrixStats package + lb <- apply(ei, 1, quantile, (level/2)/x$k) # consider using rowQuantiles() from matrixStats package + ub <- apply(ei, 1, quantile, 1-(level/2)/x$k) # consider using rowQuantiles() from matrixStats package } else { - lb <- apply(ei, 1, quantile, (level/2)) ### consider using rowQuantiles() from matrixStats package - ub <- apply(ei, 1, quantile, 1-(level/2)) ### consider using rowQuantiles() from matrixStats package + lb <- apply(ei, 1, quantile, (level/2)) # consider using rowQuantiles() from matrixStats package + ub <- apply(ei, 1, quantile, 1-(level/2)) # consider using rowQuantiles() from matrixStats package } temp.lb <- qqnorm(lb, plot.it=FALSE) diff --git a/R/radial.rma.r b/R/radial.rma.r index e6db443b5..9c21f7286 100644 --- a/R/radial.rma.r +++ b/R/radial.rma.r @@ -43,8 +43,8 @@ transf, targs, pch=21, col, bg, back, arc.res=100, cex, cex.lab, cex.axis, ...) beta <- c(x$beta) ci.lb <- x$ci.lb ci.ub <- x$ci.ub - tau2 <- 1/mean(1/x$tau2) ### geometric mean of tau^2 values (hackish solution for models with multiple tau^2 values) - ### note: this works for 1/mean(1/0) = 0; TODO: consider something more sophisticated here + tau2 <- 1/mean(1/x$tau2) # geometric mean of tau^2 values (hackish solution for models with multiple tau^2 values) + # note: this works for 1/mean(1/0) = 0; TODO: consider something more sophisticated here if (is.null(aty)) { atyis <- range(yi) } else { @@ -82,18 +82,18 @@ transf, targs, pch=21, col, bg, back, arc.res=100, cex, cex.lab, cex.axis, ...) ### set x-axis limits if none are specified if (missing(xlim)) { - xlims <- c(0, (1.30*max(xi))) ### add 30% to upper bound + xlims <- c(0, (1.30*max(xi))) # add 30% to upper bound } else { xlims <- sort(xlim) } ### x-axis position of the confidence interval - ci.xpos <- xlims[2] + 0.12*(xlims[2]-xlims[1]) ### add 12% of range to upper bound + ci.xpos <- xlims[2] + 0.12*(xlims[2]-xlims[1]) # add 12% of range to upper bound ### x-axis position of the y-axis on the right - ya.xpos <- xlims[2] + 0.14*(xlims[2]-xlims[1]) ### add 14% of range to upper bound + ya.xpos <- xlims[2] + 0.14*(xlims[2]-xlims[1]) # add 14% of range to upper bound xaxismax <- xlims[2] diff --git a/R/ranef.rma.mv.r b/R/ranef.rma.mv.r index 74ce6cd7e..0d08ad86b 100644 --- a/R/ranef.rma.mv.r +++ b/R/ranef.rma.mv.r @@ -60,7 +60,7 @@ ranef.rma.mv <- function(object, level, digits, transf, targs, verbose=FALSE, .. ### compute residuals - ei <- c(x$yi - x$X %*% x$beta) ### use this instead of resid(), since this guarantees that the length is correct + ei <- c(x$yi - x$X %*% x$beta) # use this instead of resid(), since this guarantees that the length is correct ### create identity matrix diff --git a/R/ranef.rma.uni.r b/R/ranef.rma.uni.r index 45de4c14a..db16766d8 100644 --- a/R/ranef.rma.uni.r +++ b/R/ranef.rma.uni.r @@ -52,7 +52,7 @@ ranef.rma.uni <- function(object, level, digits, transf, targs, ...) { li <- ifelse(is.infinite(x$tau2.f), 1, x$tau2.f / (x$tau2.f + x$vi.f)) - for (i in seq_len(x$k.f)[x$not.na]) { ### note: skipping NA cases + for (i in seq_len(x$k.f)[x$not.na]) { # note: skipping NA cases Xi <- matrix(x$X.f[i,], nrow=1) if (is.element(x$method, c("FE","EE","CE"))) { pred[i] <- 0 diff --git a/R/rma.mh.r b/R/rma.mh.r index b86ef1f9e..ee7dfc2ea 100644 --- a/R/rma.mh.r +++ b/R/rma.mh.r @@ -1,6 +1,6 @@ rma.mh <- function(ai, bi, ci, di, n1i, n2i, x1i, x2i, t1i, t2i, measure="OR", data, slab, subset, -add=1/2, to="only0", drop00=TRUE, ### for add/to/drop00, 1st element for escalc(), 2nd for MH method +add=1/2, to="only0", drop00=TRUE, # for add/to/drop00, 1st element for escalc(), 2nd for MH method correct=TRUE, level=95, verbose=FALSE, digits, ...) { ######################################################################### @@ -117,7 +117,7 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { if (is.null(di)) di <- n2i - ci ni <- ai + bi + ci + di - k <- length(ai) ### number of outcomes before subsetting + k <- length(ai) # number of outcomes before subsetting k.all <- k if (length(ai)==0L || length(bi)==0L || length(ci)==0L || length(di)==0L) @@ -179,8 +179,8 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { ### calculate observed effect estimates and sampling variances dat <- .do.call(escalc, measure=measure, ai=ai, bi=bi, ci=ci, di=di, add=add[1], to=to[1], drop00=drop00[1], onlyo1=onlyo1, addyi=addyi, addvi=addvi) - yi <- dat$yi ### one or more yi/vi pairs may be NA/NA - vi <- dat$vi ### one or more yi/vi pairs may be NA/NA + yi <- dat$yi # one or more yi/vi pairs may be NA/NA + vi <- dat$vi # one or more yi/vi pairs may be NA/NA ### if drop00[2]=TRUE, set counts to NA for studies that have no events (or all events) in both arms @@ -201,7 +201,7 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { vi.f <- vi ni.f <- ni - k.f <- k ### total number of tables including all NAs + k.f <- k # total number of tables including all NAs ### check for NAs in table data and act accordingly @@ -249,8 +249,8 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { ni <- ni[not.na.yivi] warning(mstyle$warning("Some yi/vi values are NA."), call.=FALSE) - attr(yi, "measure") <- measure ### add measure attribute back - attr(yi, "ni") <- ni ### add ni attribute back + attr(yi, "measure") <- measure # add measure attribute back + attr(yi, "ni") <- ni # add ni attribute back } @@ -259,7 +259,7 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { } - k.yi <- length(yi) ### number of yi/vi pairs that are not NA (needed for QE df and fit.stats calculation) + k.yi <- length(yi) # number of yi/vi pairs that are not NA (needed for QE df and fit.stats calculation) ### add/to procedures for the 2x2 tables for the actual meta-analysis ### note: technically, nothing needs to be added, but Stata/RevMan add 1/2 by default for only0 studies (but drop studies with no/all events) @@ -326,7 +326,7 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { ni <- t1i + t2i - k <- length(x1i) ### number of outcomes before subsetting + k <- length(x1i) # number of outcomes before subsetting k.all <- k if (length(x1i)==0L || length(x2i)==0L || length(t1i)==0L || length(t2i)==0L) @@ -385,8 +385,8 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { ### calculate observed effect estimates and sampling variances dat <- .do.call(escalc, measure=measure, x1i=x1i, x2i=x2i, t1i=t1i, t2i=t2i, add=add[1], to=to[1], drop00=drop00[1], onlyo1=onlyo1, addyi=addyi, addvi=addvi) - yi <- dat$yi ### one or more yi/vi pairs may be NA/NA - vi <- dat$vi ### one or more yi/vi pairs may be NA/NA + yi <- dat$yi # one or more yi/vi pairs may be NA/NA + vi <- dat$vi # one or more yi/vi pairs may be NA/NA ### if drop00[2]=TRUE, set counts to NA for studies that have no events in both arms @@ -405,7 +405,7 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { vi.f <- vi ni.f <- ni - k.f <- k ### total number of tables including all NAs + k.f <- k # total number of tables including all NAs ### check for NAs in table data and act accordingly @@ -453,8 +453,8 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { ni <- ni[not.na.yivi] warning(mstyle$warning("Some yi/vi values are NA."), call.=FALSE) - attr(yi, "measure") <- measure ### add measure attribute back - attr(yi, "ni") <- ni ### add ni attribute back + attr(yi, "measure") <- measure # add measure attribute back + attr(yi, "ni") <- ni # add ni attribute back } @@ -463,7 +463,7 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { } - k.yi <- length(yi) ### number of yi/vi pairs that are not NA (needed for QE df and fitstats calculation) + k.yi <- length(yi) # number of yi/vi pairs that are not NA (needed for QE df and fitstats calculation) ### add/to procedures for the 2x2 tables for the actual meta-analysis ### note: technically, nothing needs to be added @@ -539,7 +539,7 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { } else { beta.exp <- R/S beta <- log(beta.exp) - se <- sqrt(1/2 * (sum(Pi*Ri)/R^2 + sum(Pi*Si + Qi*Ri)/(R*S) + sum(Qi*Si)/S^2)) ### based on Robins et al. (1986) + se <- sqrt(1/2 * (sum(Pi*Ri)/R^2 + sum(Pi*Si + Qi*Ri)/(R*S) + sum(Qi*Si)/S^2)) # based on Robins et al. (1986) zval <- beta / se pval <- 2*pnorm(abs(zval), lower.tail=FALSE) ci.lb <- beta - qnorm(level/2, lower.tail=FALSE) * se @@ -632,8 +632,8 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { if (measure == "RD") { beta <- sum(ai*(n2i/Ni) - ci*(n1i/Ni)) / sum(n1i*(n2i/Ni)) - se <- sqrt((beta * (sum(ci*(n1i/Ni)^2 - ai*(n2i/Ni)^2 + (n1i/Ni)*(n2i/Ni)*(n2i-n1i)/2)) + sum(ai*(n2i-ci)/Ni + ci*(n1i-ai)/Ni)/2) / sum(n1i*(n2i/Ni))^2) ### equation in: Sato, Greenland, & Robins (1989) - #se <- sqrt(sum(((ai/Ni^2)*bi*(n2i^2/n1i) + (ci/Ni^2)*di*(n1i^2/n2i))) / sum(n1i*(n2i/Ni))^2) ### equation in: Greenland & Robins (1985) + se <- sqrt((beta * (sum(ci*(n1i/Ni)^2 - ai*(n2i/Ni)^2 + (n1i/Ni)*(n2i/Ni)*(n2i-n1i)/2)) + sum(ai*(n2i-ci)/Ni + ci*(n1i-ai)/Ni)/2) / sum(n1i*(n2i/Ni))^2) # equation in: Sato, Greenland, & Robins (1989) + #se <- sqrt(sum(((ai/Ni^2)*bi*(n2i^2/n1i) + (ci/Ni^2)*di*(n1i^2/n2i))) / sum(n1i*(n2i/Ni))^2) # equation in: Greenland & Robins (1985) zval <- beta / se pval <- 2*pnorm(abs(zval), lower.tail=FALSE) ci.lb <- beta - qnorm(level/2, lower.tail=FALSE) * se @@ -686,7 +686,7 @@ correct=TRUE, level=95, verbose=FALSE, digits, ...) { if (measure == "IRD") { beta <- sum((x1i*t2i - x2i*t1i)/Ti) / sum((t1i/Ti)*t2i) - se <- sqrt(sum(((t1i/Ti)*t2i)^2*(x1i/t1i^2+x2i/t2i^2))) / sum((t1i/Ti)*t2i) ### from Rothland et al. (2008), chapter 15 + se <- sqrt(sum(((t1i/Ti)*t2i)^2*(x1i/t1i^2+x2i/t2i^2))) / sum((t1i/Ti)*t2i) # from Rothland et al. (2008), chapter 15 zval <- beta / se pval <- 2*pnorm(abs(zval), lower.tail=FALSE) ci.lb <- beta - qnorm(level/2, lower.tail=FALSE) * se diff --git a/R/rma.mv.r b/R/rma.mv.r index 12ca03e8f..8bb290f81 100644 --- a/R/rma.mv.r +++ b/R/rma.mv.r @@ -854,11 +854,11 @@ cvvc=FALSE, sparse=FALSE, verbose=FALSE, digits, control, ...) { W.f <- A ni.f <- ni mods.f <- mods - #mf.g.f <- mf.g ### copied further below - #mf.h.f <- mf.h ### copied further below - #mf.s.f <- mf.s ### copied further below + #mf.g.f <- mf.g # copied further below + #mf.h.f <- mf.h # copied further below + #mf.s.f <- mf.s # copied further below - k.f <- k ### total number of observed outcomes including all NAs + k.f <- k # total number of observed outcomes including all NAs ######################################################################### ######################################################################### @@ -1074,7 +1074,7 @@ cvvc=FALSE, sparse=FALSE, verbose=FALSE, digits, control, ...) { } - #mf.s.f <- mf.s ### not needed at the moment + #mf.s.f <- mf.s # not needed at the moment ### copy s.nlevels and s.levels (needed for ranef()) @@ -1122,7 +1122,7 @@ cvvc=FALSE, sparse=FALSE, verbose=FALSE, digits, control, ...) { } - mf.g.f <- mf.g ### needed for predict() + mf.g.f <- mf.g # needed for predict() ######################################################################### @@ -1165,7 +1165,7 @@ cvvc=FALSE, sparse=FALSE, verbose=FALSE, digits, control, ...) { } - mf.h.f <- mf.h ### needed for predict() + mf.h.f <- mf.h # needed for predict() # return(list(Z.G1=Z.G1, Z.G2=Z.G2, g.nlevels=g.nlevels, g.levels=g.levels, g.values=g.values, tau2=tau2, rho=rho, # Z.H1=Z.H1, Z.H2=Z.H2, h.nlevels=h.nlevels, h.levels=h.levels, h.values=h.values, gamma2=gamma2, phi=phi)) diff --git a/R/rma.peto.r b/R/rma.peto.r index 03b24a98a..5504858cb 100644 --- a/R/rma.peto.r +++ b/R/rma.peto.r @@ -1,6 +1,6 @@ rma.peto <- function(ai, bi, ci, di, n1i, n2i, data, slab, subset, -add=1/2, to="only0", drop00=TRUE, ### for add/to/drop00, 1st element for escalc(), 2nd for Peto's method +add=1/2, to="only0", drop00=TRUE, # for add/to/drop00, 1st element for escalc(), 2nd for Peto's method level=95, verbose=FALSE, digits, ...) { ######################################################################### @@ -48,7 +48,7 @@ level=95, verbose=FALSE, digits, ...) { .chkdots(ddd, c("outlist", "time")) - measure <- "PETO" ### set measure here so that it can be added below + measure <- "PETO" # set measure here so that it can be added below ### set defaults for digits @@ -104,7 +104,7 @@ level=95, verbose=FALSE, digits, ...) { if (is.null(di)) di <- n2i - ci ni <- ai + bi + ci + di - k <- length(ai) ### number of outcomes before subsetting + k <- length(ai) # number of outcomes before subsetting k.all <- k if (length(ai)==0L || length(bi)==0L || length(ci)==0L || length(di)==0L) @@ -166,8 +166,8 @@ level=95, verbose=FALSE, digits, ...) { ### calculate observed effect estimates and sampling variances dat <- .do.call(escalc, measure="PETO", ai=ai, bi=bi, ci=ci, di=di, add=add[1], to=to[1], drop00=drop00[1]) - yi <- dat$yi ### one or more yi/vi pairs may be NA/NA - vi <- dat$vi ### one or more yi/vi pairs may be NA/NA + yi <- dat$yi # one or more yi/vi pairs may be NA/NA + vi <- dat$vi # one or more yi/vi pairs may be NA/NA ### if drop00[2]=TRUE, set counts to NA for studies that have no events (or all events) in both arms @@ -188,7 +188,7 @@ level=95, verbose=FALSE, digits, ...) { vi.f <- vi ni.f <- ni - k.f <- k ### total number of tables including all NAs + k.f <- k # total number of tables including all NAs ### check for NAs in table data and act accordingly @@ -236,8 +236,8 @@ level=95, verbose=FALSE, digits, ...) { ni <- ni[not.na.yivi] warning(mstyle$warning("Some yi/vi values are NA."), call.=FALSE) - attr(yi, "measure") <- measure ### add measure attribute back - attr(yi, "ni") <- ni ### add ni attribute back + attr(yi, "measure") <- measure # add measure attribute back + attr(yi, "ni") <- ni # add ni attribute back } @@ -246,7 +246,7 @@ level=95, verbose=FALSE, digits, ...) { } - k.yi <- length(yi) ### number of yi/vi pairs that are not NA (needed for QE df and fit.stats calculation) + k.yi <- length(yi) # number of yi/vi pairs that are not NA (needed for QE df and fit.stats calculation) ### add/to procedures for the 2x2 tables for the actual meta-analysis ### note: technically, nothing needs to be added, but Stata/RevMan add 1/2 by default for only0 studies (but drop studies with no/all events) @@ -305,14 +305,14 @@ level=95, verbose=FALSE, digits, ...) { if (verbose) message(mstyle$message("Model fitting ...")) - xt <- ai + ci ### frequency of outcome1 in both groups combined - yt <- bi + di ### frequency of outcome2 in both groups combined + xt <- ai + ci # frequency of outcome1 in both groups combined + yt <- bi + di # frequency of outcome2 in both groups combined Ei <- xt * n1i / Ni - Vi <- xt * yt * (n1i/Ni) * (n2i/Ni) / (Ni - 1) ### 0 when xt = 0 or yt = 0 in a table + Vi <- xt * yt * (n1i/Ni) * (n2i/Ni) / (Ni - 1) # 0 when xt = 0 or yt = 0 in a table sumVi <- sum(Vi) - if (sumVi == 0L) ### sumVi = 0 when xt or yt = 0 in *all* tables + if (sumVi == 0L) # sumVi = 0 when xt or yt = 0 in *all* tables stop(mstyle$stop("One of the two outcomes never occurred in any of the tables. Peto's method cannot be used.")) beta <- sum(ai - Ei) / sumVi @@ -332,8 +332,8 @@ level=95, verbose=FALSE, digits, ...) { if (verbose) message(mstyle$message("Heterogeneity testing ...")) - k.pos <- sum(Vi > 0) ### number of tables with positive sampling variance - Vi[Vi == 0] <- NA_real_ ### set 0 sampling variances to NA + k.pos <- sum(Vi > 0) # number of tables with positive sampling variance + Vi[Vi == 0] <- NA_real_ # set 0 sampling variances to NA QE <- max(0, sum((ai - Ei)^2 / Vi, na.rm=TRUE) - sum(ai - Ei)^2 / sum(Vi, na.rm=TRUE)) if (k.pos > 1L) { diff --git a/R/rma.uni.r b/R/rma.uni.r index 8da9bb219..c5cdd950e 100644 --- a/R/rma.uni.r +++ b/R/rma.uni.r @@ -337,7 +337,7 @@ test="z", level=95, btt, att, tau2, verbose=FALSE, digits, control, ...) { ### allow easy setting of vi to a single value if (length(vi) == 1L) - vi <- rep(vi, k) ### note: k is number of outcomes before subsetting + vi <- rep(vi, k) # note: k is number of outcomes before subsetting ### check length of yi and vi @@ -900,7 +900,7 @@ test="z", level=95, btt, att, tau2, verbose=FALSE, digits, control, ...) { mods.f <- mods Z.f <- Z - k.f <- k ### total number of observed outcomes including all NAs + k.f <- k # total number of observed outcomes including all NAs ### check for NAs and act accordingly @@ -1017,7 +1017,7 @@ test="z", level=95, btt, att, tau2, verbose=FALSE, digits, control, ...) { int.incl <- FALSE } - p <- NCOL(X) ### number of columns in X (including the intercept if it is included) + p <- NCOL(X) # number of columns in X (including the intercept if it is included) ### make sure variable names in X and Z are unique @@ -1052,7 +1052,7 @@ test="z", level=95, btt, att, tau2, verbose=FALSE, digits, control, ...) { ### set/check 'btt' argument btt <- .set.btt(btt, p, int.incl, colnames(X)) - m <- length(btt) ### number of betas to test (m = p if all betas are tested) + m <- length(btt) # number of betas to test (m = p if all betas are tested) ######################################################################### @@ -2187,7 +2187,7 @@ test="z", level=95, btt, att, tau2, verbose=FALSE, digits, control, ...) { ### set/check 'att' argument att <- .set.btt(att, q, Z.int.incl, colnames(Z)) - m.alpha <- length(att) ### number of alphas to test (m = q if all alphas are tested) + m.alpha <- length(att) # number of alphas to test (m = q if all alphas are tested) ### ddf calculation @@ -2453,7 +2453,7 @@ test="z", level=95, btt, att, tau2, verbose=FALSE, digits, control, ...) { if (k > p) { wi <- 1/vi - W.FE <- diag(wi, nrow=k, ncol=k) ### note: ll.REML below involves W, so cannot overwrite W + W.FE <- diag(wi, nrow=k, ncol=k) # note: ll.REML below involves W, so cannot overwrite W stXWX <- .invcalc(X=X, W=W.FE, k=k) P <- W.FE - W.FE %*% X %*% stXWX %*% crossprod(X,W.FE) # need P below for calculation of I^2 QE <- max(0, c(crossprod(Ymc,P) %*% Ymc)) diff --git a/R/rstandard.rma.mh.r b/R/rstandard.rma.mh.r index e6f9f3534..ba5542c54 100644 --- a/R/rstandard.rma.mh.r +++ b/R/rstandard.rma.mh.r @@ -22,7 +22,7 @@ rstandard.rma.mh <- function(model, digits, ...) { resid <- c(x$yi.f - x$beta) resid[abs(resid) < 100 * .Machine$double.eps] <- 0 - #resid[abs(resid) < 100 * .Machine$double.eps * median(abs(resid), na.rm=TRUE)] <- 0 ### see lm.influence + #resid[abs(resid) < 100 * .Machine$double.eps * median(abs(resid), na.rm=TRUE)] <- 0 # see lm.influence ### note: these are like Pearson (or semi-standardized) residuals diff --git a/R/rstandard.rma.mv.r b/R/rstandard.rma.mv.r index b6ca93834..ea73a4284 100644 --- a/R/rstandard.rma.mv.r +++ b/R/rstandard.rma.mv.r @@ -65,7 +65,7 @@ rstandard.rma.mv <- function(model, digits, cluster, ...) { ei <- c(x$yi - x$X %*% x$beta) ei[abs(ei) < 100 * .Machine$double.eps] <- 0 - #ei[abs(ei) < 100 * .Machine$double.eps * median(abs(ei), na.rm=TRUE)] <- 0 ### see lm.influence + #ei[abs(ei) < 100 * .Machine$double.eps * median(abs(ei), na.rm=TRUE)] <- 0 # see lm.influence ### don't allow this; the SEs of the residuals cannot be estimated consistently for "robust.rma" objects #if (inherits(x, "robust.rma")) { diff --git a/R/rstandard.rma.peto.r b/R/rstandard.rma.peto.r index 8fcd665d4..484510328 100644 --- a/R/rstandard.rma.peto.r +++ b/R/rstandard.rma.peto.r @@ -22,7 +22,7 @@ rstandard.rma.peto <- function(model, digits, ...) { resid <- c(x$yi.f - x$beta) resid[abs(resid) < 100 * .Machine$double.eps] <- 0 - #resid[abs(resid) < 100 * .Machine$double.eps * median(abs(resid), na.rm=TRUE)] <- 0 ### see lm.influence + #resid[abs(resid) < 100 * .Machine$double.eps * median(abs(resid), na.rm=TRUE)] <- 0 # see lm.influence ### note: these are like Pearson (or semi-standardized) residuals diff --git a/R/rstandard.rma.uni.r b/R/rstandard.rma.uni.r index c6e942e89..06d3c3197 100644 --- a/R/rstandard.rma.uni.r +++ b/R/rstandard.rma.uni.r @@ -42,7 +42,7 @@ rstandard.rma.uni <- function(model, digits, type="marginal", ...) { ei <- c(x$yi - x$X %*% x$beta) ei[abs(ei) < 100 * .Machine$double.eps] <- 0 - #ei[abs(ei) < 100 * .Machine$double.eps * median(abs(ei), na.rm=TRUE)] <- 0 ### see lm.influence + #ei[abs(ei) < 100 * .Machine$double.eps * median(abs(ei), na.rm=TRUE)] <- 0 # see lm.influence ### don't allow this; the SEs of the residuals cannot be estimated consistently for "robust.rma" objects #if (inherits(x, "robust.rma")) { diff --git a/R/rstudent.rma.mh.r b/R/rstudent.rma.mh.r index 4f8618ba2..0d73bf71f 100644 --- a/R/rstudent.rma.mh.r +++ b/R/rstudent.rma.mh.r @@ -68,7 +68,7 @@ rstudent.rma.mh <- function(model, digits, progbar=FALSE, ...) { resid <- x$yi.f - delpred resid[abs(resid) < 100 * .Machine$double.eps] <- 0 - #resid[abs(resid) < 100 * .Machine$double.eps * median(abs(resid), na.rm=TRUE)] <- 0 ### see lm.influence + #resid[abs(resid) < 100 * .Machine$double.eps * median(abs(resid), na.rm=TRUE)] <- 0 # see lm.influence seresid <- sqrt(x$vi.f + vdelpred) stresid <- resid / seresid diff --git a/R/rstudent.rma.peto.r b/R/rstudent.rma.peto.r index 810dd993d..f41aef070 100644 --- a/R/rstudent.rma.peto.r +++ b/R/rstudent.rma.peto.r @@ -62,7 +62,7 @@ rstudent.rma.peto <- function(model, digits, progbar=FALSE, ...) { resid <- x$yi.f - delpred resid[abs(resid) < 100 * .Machine$double.eps] <- 0 - #resid[abs(resid) < 100 * .Machine$double.eps * median(abs(resid), na.rm=TRUE)] <- 0 ### see lm.influence + #resid[abs(resid) < 100 * .Machine$double.eps * median(abs(resid), na.rm=TRUE)] <- 0 # see lm.influence seresid <- sqrt(x$vi.f + vdelpred) stresid <- resid / seresid diff --git a/R/to.long.r b/R/to.long.r index eac809f27..aba5936d7 100644 --- a/R/to.long.r +++ b/R/to.long.r @@ -116,7 +116,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (!.all.specified(ai, bi, ci, di)) stop(mstyle$stop("Cannot compute outcomes. Check that all of the required information is specified\n via the appropriate arguments (i.e., ai, bi, ci, di or ai, n1i, ci, n2i).")) - k <- length(ai) ### number of outcomes before subsetting + k <- length(ai) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -138,7 +138,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (any(c(n1i < 0, n2i < 0), na.rm=TRUE)) stop(mstyle$stop("One or more group sizes are negative.")) - ni.u <- ai + bi + ci + di ### unadjusted total sample sizes + ni.u <- ai + bi + ci + di # unadjusted total sample sizes ### if drop00=TRUE, set counts to NA for studies that have no events (or all events) in both arms @@ -211,7 +211,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (!.equal.length(x1i, x2i, t1i, t2i)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(x1i) ### number of outcomes before subsetting + k <- length(x1i) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -227,7 +227,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (any(c(t1i, t2i) <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more person-times are <= 0.")) - ni.u <- t1i + t2i ### unadjusted total sample sizes + ni.u <- t1i + t2i # unadjusted total sample sizes ### if drop00=TRUE, set counts to NA for studies that have no events in both arms @@ -294,7 +294,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (!.equal.length(m1i, m2i, sd1i, sd2i, n1i, n2i)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(n1i) ### number of outcomes before subsetting + k <- length(n1i) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -312,7 +312,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (any(c(n1i, n2i) <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more group sizes are <= 0.")) - ni.u <- n1i + n2i ### unadjusted total sample sizes + ni.u <- n1i + n2i # unadjusted total sample sizes } @@ -332,7 +332,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (!.all.specified(ri, ni)) stop(mstyle$stop("Cannot compute outcomes. Check that all of the required information is specified\n via the appropriate arguments (i.e., ri, ni).")) - k <- length(ri) ### number of outcomes before subsetting + k <- length(ri) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k=k) @@ -346,7 +346,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (any(ni <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more sample sizes are <= 0.")) - ni.u <- ni ### unadjusted total sample sizes + ni.u <- ni # unadjusted total sample sizes } @@ -371,7 +371,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (!.all.specified(xi, mi)) stop(mstyle$stop("Cannot compute outcomes. Check that all of the required information is specified\n via the appropriate arguments (i.e., xi, mi or xi, ni).")) - k <- length(xi) ### number of outcomes before subsetting + k <- length(xi) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -390,7 +390,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (any(ni <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more group sizes are <= 0.")) - ni.u <- ni ### unadjusted total sample sizes + ni.u <- ni # unadjusted total sample sizes if (to == "all") { @@ -444,7 +444,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (!.equal.length(xi, ti)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(xi) ### number of outcomes before subsetting + k <- length(xi) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -458,7 +458,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (any(ti <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more person-times are <= 0.")) - ni.u <- ti ### unadjusted total sample sizes + ni.u <- ti # unadjusted total sample sizes if (to == "all") { @@ -510,7 +510,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (!.equal.length(mi, sdi, ni)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(ni) ### number of outcomes before subsetting + k <- length(ni) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -528,7 +528,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (is.element(measure, c("MNLN","CVLN")) && any(mi < 0, na.rm=TRUE)) stop(mstyle$stop("One or more means are negative.")) - ni.u <- ni ### unadjusted total sample sizes + ni.u <- ni # unadjusted total sample sizes } @@ -540,10 +540,10 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, m2i <- .getx("m2i", mf=mf, data=data, checknumeric=TRUE) sd1i <- .getx("sd1i", mf=mf, data=data, checknumeric=TRUE) sd2i <- .getx("sd2i", mf=mf, data=data, checknumeric=TRUE) - ri <- .getx("ri", mf=mf, data=data, checknumeric=TRUE) ### for SMCR, do not need to supply this + ri <- .getx("ri", mf=mf, data=data, checknumeric=TRUE) # for SMCR, do not need to supply this ni <- .getx("ni", mf=mf, data=data, checknumeric=TRUE) - k <- length(m1i) ### number of outcomes before subsetting + k <- length(m1i) # number of outcomes before subsetting if (is.element(measure, c("MC","SMCC","SMCRH","ROMC","CVRC"))) { @@ -587,7 +587,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (any(ni <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more sample sizes are <= 0.")) - ni.u <- ni ### unadjusted total sample sizes + ni.u <- ni # unadjusted total sample sizes } @@ -605,7 +605,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (!.equal.length(ai, mi, ni)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(ai) ### number of outcomes before subsetting + k <- length(ai) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -623,7 +623,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, vlong=FALSE, append=TRUE, if (any(ni <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more sample sizes are <= 0.")) - ni.u <- ni ### unadjusted total sample sizes + ni.u <- ni # unadjusted total sample sizes } diff --git a/R/to.table.r b/R/to.table.r index b565ffc13..f7e08afdd 100644 --- a/R/to.table.r +++ b/R/to.table.r @@ -88,7 +88,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (!.all.specified(ai, bi, ci, di)) stop(mstyle$stop("Cannot compute outcomes. Check that all of the required information is specified\n via the appropriate arguments (i.e., ai, bi, ci, di or ai, n1i, ci, n2i).")) - k <- length(ai) ### number of outcomes before subsetting + k <- length(ai) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -110,7 +110,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (any(c(n1i < 0, n2i < 0), na.rm=TRUE)) stop(mstyle$stop("One or more group sizes are negative.")) - ni.u <- ai + bi + ci + di ### unadjusted total sample sizes + ni.u <- ai + bi + ci + di # unadjusted total sample sizes ### if drop00=TRUE, set counts to NA for studies that have no events (or all events) in both arms @@ -183,7 +183,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (!.equal.length(x1i, x2i, t1i, t2i)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(x1i) ### number of outcomes before subsetting + k <- length(x1i) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -199,7 +199,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (any(c(t1i, t2i) <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more person-times are <= 0.")) - ni.u <- t1i + t2i ### unadjusted total sample sizes + ni.u <- t1i + t2i # unadjusted total sample sizes ### if drop00=TRUE, set counts to NA for studies that have no events in both arms @@ -266,7 +266,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (!.equal.length(m1i, m2i, sd1i, sd2i, n1i, n2i)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(n1i) ### number of outcomes before subsetting + k <- length(n1i) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -284,7 +284,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (any(c(n1i, n2i) <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more group sizes are <= 0.")) - ni.u <- n1i + n2i ### unadjusted total sample sizes + ni.u <- n1i + n2i # unadjusted total sample sizes } @@ -304,7 +304,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (!.all.specified(ri, ni)) stop(mstyle$stop("Cannot compute outcomes. Check that all of the required information is specified\n via the appropriate arguments (i.e., ri, ni).")) - k <- length(ri) ### number of outcomes before subsetting + k <- length(ri) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -318,7 +318,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (any(ni <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more sample sizes are <= 0.")) - ni.u <- ni ### unadjusted total sample sizes + ni.u <- ni # unadjusted total sample sizes } @@ -343,7 +343,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (!.all.specified(xi, mi)) stop(mstyle$stop("Cannot compute outcomes. Check that all of the required information is specified\n via the appropriate arguments (i.e., xi, mi or xi, ni).")) - k <- length(xi) ### number of outcomes before subsetting + k <- length(xi) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -362,7 +362,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (any(ni <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more group sizes are <= 0.")) - ni.u <- ni ### unadjusted total sample sizes + ni.u <- ni # unadjusted total sample sizes if (to == "all") { @@ -416,7 +416,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (!.equal.length(xi, ti)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(xi) ### number of outcomes before subsetting + k <- length(xi) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -430,7 +430,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (any(ti <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more person-times are <= 0.")) - ni.u <- ti ### unadjusted total sample sizes + ni.u <- ti # unadjusted total sample sizes if (to == "all") { @@ -482,7 +482,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (!.equal.length(mi, sdi, ni)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(ni) ### number of outcomes before subsetting + k <- length(ni) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -500,7 +500,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (is.element(measure, c("MNLN","CVLN")) && any(mi < 0, na.rm=TRUE)) stop(mstyle$stop("One or more means are negative.")) - ni.u <- ni ### unadjusted total sample sizes + ni.u <- ni # unadjusted total sample sizes } @@ -512,10 +512,10 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { m2i <- .getx("m2i", mf=mf, data=data, checknumeric=TRUE) sd1i <- .getx("sd1i", mf=mf, data=data, checknumeric=TRUE) sd2i <- .getx("sd2i", mf=mf, data=data, checknumeric=TRUE) - ri <- .getx("ri", mf=mf, data=data, checknumeric=TRUE) ### for SMCR, do not need to supply this + ri <- .getx("ri", mf=mf, data=data, checknumeric=TRUE) # for SMCR, do not need to supply this ni <- .getx("ni", mf=mf, data=data, checknumeric=TRUE) - k <- length(m1i) ### number of outcomes before subsetting + k <- length(m1i) # number of outcomes before subsetting if (is.element(measure, c("MC","SMCC","SMCRH","ROMC","CVRC"))) { @@ -559,7 +559,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (any(ni <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more sample sizes are <= 0.")) - ni.u <- ni ### unadjusted total sample sizes + ni.u <- ni # unadjusted total sample sizes } @@ -577,7 +577,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (!.equal.length(ai, mi, ni)) stop(mstyle$stop("Supplied data vectors are not all of the same length.")) - k <- length(ai) ### number of outcomes before subsetting + k <- length(ai) # number of outcomes before subsetting if (!is.null(subset)) { subset <- .chksubset(subset, k) @@ -595,7 +595,7 @@ data, slab, subset, add=1/2, to="none", drop00=FALSE, rows, cols) { if (any(ni <= 0, na.rm=TRUE)) stop(mstyle$stop("One or more sample sizes are <= 0.")) - ni.u <- ni ### unadjusted total sample sizes + ni.u <- ni # unadjusted total sample sizes } diff --git a/R/trimfill.rma.uni.r b/R/trimfill.rma.uni.r index 3aa93519b..f8dca6070 100644 --- a/R/trimfill.rma.uni.r +++ b/R/trimfill.rma.uni.r @@ -57,15 +57,15 @@ trimfill.rma.uni <- function(x, side, estimator="L0", maxiter=100, verbose=FALSE k <- length(yi) k0.sav <- -1 - k0 <- 0 ### estimated number of missing studies - iter <- 0 ### iteration counter + k0 <- 0 # estimated number of missing studies + iter <- 0 # iteration counter if (verbose) cat("\n") while (abs(k0 - k0.sav) > 0) { - k0.sav <- k0 ### save current value of k0 + k0.sav <- k0 # save current value of k0 iter <- iter + 1 @@ -85,9 +85,9 @@ trimfill.rma.uni <- function(x, side, estimator="L0", maxiter=100, verbose=FALSE beta <- c(res$beta) - yi.c <- yi - beta ### centered values - yi.c.r <- rank(abs(yi.c), ties.method="first") ### ranked absolute centered values - yi.c.r.s <- sign(yi.c) * yi.c.r ### signed ranked centered values + yi.c <- yi - beta # centered values + yi.c.r <- rank(abs(yi.c), ties.method="first") # ranked absolute centered values + yi.c.r.s <- sign(yi.c) * yi.c.r # signed ranked centered values ### estimate the number of missing studies with the R0 estimator diff --git a/R/zzz.r b/R/zzz.r index fa631dca9..b4de64050 100644 --- a/R/zzz.r +++ b/R/zzz.r @@ -1,6 +1,6 @@ .onAttach <- function(libname, pkgname) { - ver <- "4.5-13" + ver <- "4.6-0" loadmsg <- paste0("\nLoading the 'metafor' package (version ", ver, "). For an\nintroduction to the package please type: help(metafor)\n")