Merge branch 'release/v0.6.1'

CITATION.cff

@@ -8,7 +8,7 @@ message: 'To cite package "multinma" in publications use:'
 type: software
 license: GPL-3.0-only
 title: 'multinma: Bayesian Network Meta-Analysis of Individual and Aggregate Data'
-version: 0.6.0
+version: 0.6.1
 doi: 10.5281/zenodo.3904454
 abstract: Network meta-analysis and network meta-regression models for aggregate data,
   individual patient data, and mixtures of both individual and aggregate data using
@@ -28,7 +28,7 @@ preferred-citation:
     email: david.phillippo@bristol.ac.uk
     orcid: https://orcid.org/0000-0003-2672-7841
   year: '2023'
-  notes: R package version 0.6.0
+  notes: R package version 0.6.1
   url: https://dmphillippo.github.io/multinma/
   doi: 10.5281/zenodo.3904454
 repository: https://CRAN.R-project.org/package=multinma

DESCRIPTION

@@ -1,6 +1,6 @@
 Package: multinma
 Title: Bayesian Network Meta-Analysis of Individual and Aggregate Data
-Version: 0.6.0
+Version: 0.6.1
 Authors@R: 
     person(given = c("David", "M."),
            family = "Phillippo",

NEWS.md

@@ -1,3 +1,8 @@
+# multinma 0.6.1
+
+* Fix: Piecewise exponential hazard models no longer give errors during set-up. 
+Calculation of RW1 prior weights needed to be handled as a special case.
+
 # multinma 0.6.0
 
 ## Feature: Survival/time-to-event models are now supported

R/integration.R

@@ -515,14 +515,19 @@ distr <- function(qfun, ...) {
 
 #' The Bernoulli Distribution
 #'
-#' The quantile function `qbern` for a Bernoulli distribution, with success
-#' probability `prob`. This is equivalent to `qbinom(p, 1, prob)`.
+#' The density function `dbern()`, distribution function `pbern()`, and quantile
+#' function `qbern()` for a Bernoulli distribution, with success probability
+#' `prob`. These are equivalent to `dbinom(p, 1, prob)`, `pbinom(p, 1, prob)`
+#' and `qbinom(p, 1, prob)`.
 #'
 #' @param x,q vector of quantiles
 #' @param p vector of probabilities
 #' @param prob probability of success
 #' @param lower.tail,log.p,log see [stats::Binomial]
 #'
+#' @return Numeric vector of length equal to the maximum of the lengths of the
+#'   input arguments.
+#'
 #' @export
 #' @rdname Bernoulli
 #' @aliases qbern
@@ -556,6 +561,9 @@ dbern <- function(x, prob, log = FALSE) {
 #' @param mean,sd mean and standard deviation, overriding `shape` and
 #'   `rate` or `scale` if specified
 #'
+#' @return Numeric vector of length equal to the maximum of the lengths of the
+#'   input arguments.
+#'
 #' @export
 #' @rdname GammaDist
 #' @aliases qgamma
@@ -606,6 +614,9 @@ pgamma <- function(q, shape, rate = 1, scale = 1/rate, lower.tail = TRUE,
 #' @param mean,sd mean and standard deviation, overriding `mu` and `sigma` if
 #'   specified
 #'
+#' @return Numeric vector of length equal to the maximum of the lengths of the
+#'   input arguments.
+#'
 #' @export
 #' @rdname logitNormal
 #' @aliases qlogitnorm

R/mspline.R

@@ -268,8 +268,12 @@ rw1_prior_weights <- function(basis) {
   iknots <- attr(basis, "knots")
   bknots <- attr(basis, "Boundary.knots")
   knots <- c(rep(bknots[1], ord), iknots, rep(bknots[2], ord))
-  wts <- 1/(ord - 1) * (knots[(ord + 1):(nscoef + ord - 1)] - knots[2:nscoef])
-  return(sqrt(wts / sum(wts)))
+  if (ord == 1) {
+    wts <- (knots[2:nscoef] - knots[1:(nscoef - 1)]) / (knots[nscoef] - bknots[1])
+  } else {
+    wts <- (knots[(ord + 1):(nscoef + ord - 1)] - knots[2:nscoef]) / ((ord - 1) * (bknots[2] - bknots[1]))
+  }
+  return(sqrt(wts))
 }
 
 #' softmax transform

R/nma_data-class.R

@@ -66,6 +66,8 @@ NULL
 #' @param ... other options (not used)
 #' @param n number of studies of each type to print
 #'
+#' @return `x` is returned invisibly.
+#'
 #' @export
 print.nma_data <- function(x, ..., n = 10) {
   cwidth <- getOption("width")

R/nma_nodesplit-class.R

@@ -33,6 +33,8 @@ NULL
 #' @seealso The summary method [summary.nma_nodesplit_df()] summarises the
 #'   node-splitting results.
 #'
+#' @return `x` is returned invisibly.
+#'
 #' @export
 print.nma_nodesplit_df <- function(x, ...) {
   n_ns <- nrow(dplyr::filter(x, !is.na(.data$trt1) & !is.na(.data$trt2)))

R/predict.R

@@ -2089,12 +2089,11 @@ make_surv_predict <- function(eta, aux, times, likelihood,
   if (type %in% c("survival", "hazard", "cumhaz") && n_eta == 1) { # Multiple times, single linear predictor
     if (!is.null(aux)) aux <- matrix(aux, ncol = dim(aux)[3], dimnames = list(NULL, dimnames(aux)[[3]]))
     for (i in 1:length(times)) {
-      out[, , i] <- do.call(surv_predfuns[[likelihood]][[type]],
-                            args = list(times = times[i],
-                                        eta = as.vector(eta),
-                                        aux = aux,
-                                        quantiles = quantiles,
-                                        basis = basis))
+      out[, , i] <- surv_predfuns[[likelihood]][[type]](times = times[i],
+                                                        eta = as.vector(eta),
+                                                        aux = aux,
+                                                        quantiles = quantiles,
+                                                        basis = basis)
     }
   } else if (n_eta > 1) { # Single/multiple times, multiple linear predictors
     if (type == "quantile") {
@@ -2118,22 +2117,20 @@ make_surv_predict <- function(eta, aux, times, likelihood,
         }
       }
       out[ , , ((i-1)*iinc+1):(i*iinc)] <-
-        do.call(surv_predfuns[[likelihood]][[type]],
-                args = list(times = ti,
-                            eta = as.vector(eta[ , , i]),
-                            aux = auxi,
-                            quantiles = quantiles,
-                            basis = basis))
+        surv_predfuns[[likelihood]][[type]](times = ti,
+                                            eta = as.vector(eta[ , , i]),
+                                            aux = auxi,
+                                            quantiles = quantiles,
+                                            basis = basis)
     }
   } else { # Single time, single linear predictor
     if (!is.null(aux)) aux <- matrix(aux, ncol = dim(aux)[3], dimnames = list(NULL, dimnames(aux)[[3]]))
     if (type == "quantile") quantiles <- rep(quantiles, each = length(eta))
-    out <- array(do.call(surv_predfuns[[likelihood]][[type]],
-                         args = list(times = times,
-                                     eta = as.vector(eta),
-                                     aux = aux,
-                                     quantiles = quantiles,
-                                     basis = basis)),
+    out <- array(surv_predfuns[[likelihood]][[type]](times = times,
+                                                     eta = as.vector(eta),
+                                                     aux = aux,
+                                                     quantiles = quantiles,
+                                                     basis = basis),
                  dim = d_out, dimnames = dn_out)
   }
 
@@ -2284,11 +2281,10 @@ rmst_Sbar <- function(times, eta, weights, aux, likelihood, basis, start = 0) {
 pSbar <- function(times, eta, weights, aux, likelihood, basis) {
   Sbar <- numeric(length(times))
   for (i in 1:length(times)) {
-    S <- do.call(surv_predfuns[[likelihood]][["survival"]],
-                       args = list(times = times[i],
-                                   eta = eta,
-                                   aux = aux,
-                                   basis = basis))
+    S <- surv_predfuns[[likelihood]][["survival"]](times = times[i],
+                                                   eta = eta,
+                                                   aux = aux,
+                                                   basis = basis)
 
     Sbar[i] <- weighted.mean(S, weights)
   }
@@ -2346,11 +2342,10 @@ qSbar <- function(times, p, eta, weights, aux, likelihood, basis) {
         auxi <- as.vector(aux[ , , i])
       }
     }
-    S[ , , i] <- do.call(surv_predfuns[[likelihood]][["survival"]],
-                         args = list(times = times,
-                                     eta = as.vector(eta[ , , i]),
-                                     aux = auxi,
-                                     basis = basis))
+    S[ , , i] <- surv_predfuns[[likelihood]][["survival"]](times = times,
+                                                           eta = as.vector(eta[ , , i]),
+                                                           aux = auxi,
+                                                           basis = basis)
   }
 
   Sbar <- apply(S, MARGIN = 1:2, FUN = weighted.mean, w = weights)

R/stan_nma-class.R

@@ -40,6 +40,8 @@ NULL
 #' @param x A [stan_nma] object
 #' @param ... Further arguments passed to [print.stanfit()]
 #'
+#' @return `x` is returned invisibly.
+#'
 #' @export
 print.stan_nma <- function(x, ...) {
   if (inherits(x$network, "mlnmr_data")) type <- "ML-NMR"

README.md

@@ -78,7 +78,7 @@ papers:
 The `multinma` package can be cited as follows:
 
 > Phillippo, D. M. (2024). *multinma: Bayesian Network Meta-Analysis of
-> Individual and Aggregate Data*. R package version 0.6.0, doi:
+> Individual and Aggregate Data*. R package version 0.6.1, doi:
 > [10.5281/zenodo.3904454](https://doi.org/10.5281/zenodo.3904454).
 
 When fitting ML-NMR models, please cite the methods paper:

cran-comments.md

@@ -1,12 +1,35 @@
-This release ensures compatibility with the forthcoming update to rstan and
-StanHeaders, as well as providing extensive new features, improvements, and
-fixes.
+This release addresses an "additional checks" UBSAN failure, caused by a memory
+allocation bug in Stan (https://github.com/stan-dev/rstan/issues/1111).
+StanHeaders has now been patched which resolves the issue; no changes are
+necessary in multinma. Verified passing UBSAN checks with rocker/r-devel-san.
+
+In response to CRAN request, I have also:
+
+* Added the few missing \value tags to documentation
+
+Regarding the other issues raised:
+
+* \dontrun{} is used in man/pairs.stan_nma.Rd, as the example is designed to
+  raise warning messages. These are divergent transitions in Stan, which we then
+  demonstrate for the user how to examine using the pairs() plot.
+
+* The global environment is never modified. rm(list = ls()) in some large test
+  files is never run on CRAN or by the user, these lines only run when testing 
+  on GitHub Actions to avoid running out of memory. These tests are skipped on 
+  CRAN to keep check time manageable.
+
+* Similarly, .GlobalEnv is never modified. In the .Rd-files starting with 
+  aa_example_, the assign() to .GlobalEnv only runs and is only needed when 
+  building the package website using pkgdown, i.e. when pkgdown::in_pkgdown() is
+  TRUE. This code is never run otherwise.
+
 
 ## Test environments
-* local R installation (Windows 10), R 4.3.2
+* local rocker/r-devel-san image with sanitizers
+* local R installation (Windows 10), R 4.3.3
 * win-builder (release, devel)
-* Ubuntu 22.04.3 on GitHub Actions (release, devel, oldrel)
-* Mac OS 12.7.2 on GitHub Actions (release)
+* Ubuntu 22.04.4 on GitHub Actions (release, devel, oldrel)
+* Mac OS 12.7.3 on GitHub Actions (release)
 * Windows Server 2022 10.0.20348 on GitHub Actions (release)
 
 ## R CMD check results