Changes requested by CRAN

DESCRIPTION

@@ -10,7 +10,7 @@ Authors@R:
 	   person(given="Kofi Placid", family="Adragni", role="ctb", comment="Author of code taken from the ldr package"),
 	   person(given="Andrew", family="Raim", role="ctb",comment="Author of code taken from the ldr package")
 	   )
-Description: Methods to calculate the expected value of information from a decision-analytic model.  This includes the expected value of perfect information (EVPI), partial perfect information (EVPPI) and sample information (EVSI), and the expected net benefit of sampling (ENBS).  A range of alternative computational methods are provided under the same user interface.  The package accompanies a forthcoming book: "Value of Information for Healthcare Decision-Making" (CRC Press, eds. Heath, Jackson, Kunst).
+Description: Methods to calculate the expected value of information from a decision-analytic model.  This includes the expected value of perfect information (EVPI), partial perfect information (EVPPI) and sample information (EVSI), and the expected net benefit of sampling (ENBS).  A range of alternative computational methods are provided under the same user interface.  See Jackson et al. (2022) <doi:10.1146/annurev-statistics-040120-010730>.
 License: GPL-3
 Encoding: UTF-8
 LazyData: true

R/check_regression.R

@@ -113,6 +113,8 @@ check_plot_default <- function(mod){
   if (!is.numeric(res))
     warning("residuals() does not work on regression model object, so can't produce diagnostic plots")
   dat <- data.frame(fit = fit, res = res)
+  oldpar <- graphics::par(no.readonly=TRUE)
+  on.exit(par(oldpar))
   graphics::par(mfrow=c(2,1))
   bw <- 2 * IQR(dat$res) / length(dat$res)^(1/3)
   p1 <- ggplot2::ggplot(dat, aes(x=res)) +

R/evppi_earth.R

@@ -22,6 +22,8 @@ fitted_earth <- function(y, inputs, pars, verbose=FALSE, ...){
 }
 
 check_plot_earth <- function(mod){
+  oldpar <- graphics::par(no.readonly=TRUE)
+  on.exit(par(oldpar))
   graphics::par(mfrow=c(2,2))
   plot(mod)
 }

R/evppi_gam.R

@@ -25,8 +25,10 @@ fitted_rep_gam <- function(model, B) {
 }
 
 check_plot_gam <- function(mod){
-    graphics::par(mfrow=c(2,2))
-    mgcv::gam.check(mod)
+  oldpar <- graphics::par(no.readonly=TRUE)
+  on.exit(par(oldpar))
+  graphics::par(mfrow=c(2,2))
+  mgcv::gam.check(mod)
 }
 
 check_stats_gam <- function(mod){

R/evppi_gp.R

@@ -243,8 +243,10 @@ gp.check <- function(mod){
 }
 
 check_plot_gp <- function(mod){
-    graphics::par(mfrow=c(2,2))
-    gp.check(mod)
+  oldpar <- graphics::par(no.readonly=TRUE)
+  on.exit(par(oldpar))
+  graphics::par(mfrow=c(2,2))
+  gp.check(mod)
 }
 
 check_stats_gp <- function(mod){

R/evppi_inla.R

@@ -156,8 +156,10 @@ fit.inla <- function(parameter, inputs, x, mesh,
 }
 
 check_plot_inla <- function(mod){
-    graphics::par(mfrow=c(2,2))
-    gp.check(mod)
+  oldpar <- graphics::par(no.readonly=TRUE)
+  on.exit(par(oldpar))
+  graphics::par(mfrow=c(2,2))
+  gp.check(mod)
 }
 
 check_stats_inla <- function(mod){

R/evppivar.R

@@ -6,7 +6,8 @@
 ##'
 ##' Typically this will come from a Monte Carlo sample, where we first sample from the uncertainty distributions of the parameters, and then compute the quantity of interest as a function of the parameters.  It might also be produced by a Markov Chain Monte Carlo sample from the joint distribution of parameters and outputs. 
 ##'
-##'
+##' @return A data frame with a column \code{pars}, indicating the parameter(s), and a column \code{evppi}, giving the corresponding EVPPI.
+##' 
 ##' @inheritParams evppi
 ##' 
 ##' @references

inst/book_misc/he_graphs.R

@@ -0,0 +1,29 @@
+library(tidyverse)
+library(viridis)
+
+bet <- SHELF::fitdist(vals=c(0.2, 0.4, 0.6),
+                      probs=c(0.025, 0.5, 0.975),
+                      lower=0, upper=1)$Beta
+
+## https://en.wikipedia.org/wiki/Logit-normal_distribution
+dlogitnorm <- function(x, mu, sd){
+  dnorm(qlogis(x), mu, sd) / (x*(1-x))
+}
+
+pdat <- data.frame(x = seq(0, 1, by=0.01)) %>%
+  mutate(dbet = dbeta(x, bet[["shape1"]], bet[["shape2"]]),
+         dln = plogitnorm(x, qlogis(0.4), 
+                     (qlogis(0.6) - qlogis(0.2))/(2*qnorm(0.975))))
+
+green <- viridis(5)[4]
+purple <- viridis(5)[2]
+pdf("~/work/voibook/voibook/Figures/02-healthecon/elic_compare.pdf", width=6, height=4)
+ggplot(pdat, aes(x=x, y=dbet)) + 
+  geom_line(lwd=2, col=green) + 
+  geom_line(aes(y=dln), col=purple, lwd=2) + 
+  geom_text(aes(x=0.55, y=3), label="Beta(8.9, 13.1)", col=green, hjust=0) + 
+  geom_text(aes(x=0.55, y=2.6), label="Logit-normal(-0.4,0.46)", col=purple, hjust=0) + 
+  theme_bw() + 
+  xlab(expression(italic(p))) + ylab("Probability density") + 
+  scale_x_continuous(breaks = c(0, 0.2, 0.4, 0.6, 0.8, 1))
+dev.off()