Variable Life Adjusted Display and Other Risk-Adjusted Quality Control Charts

Contains functions to set up risk-adjusted quality control charts in health care. For the variable life adjusted display (VLAD) proposed by Lovegrove et al. (1997) signaling rules derived in Wittenberg et al. (2018) are implemented. Additionally, for the risk-adjusted cumulative sum chart based on log-likelihood ratio statistic introduced by Steiner et al. (2000) average run length and control limits can be computed.


An R-package which contains functions to set up risk-adjusted quality control charts in health care.

Main features

  • Risk-adjusted CUSUM chart based on E-O by Wittenberg et al. (2018)
  • Risk-adjusted CUSUM chart based on log-likelihood ratio statistic by Steiner et al. (2000)
  • Algorithms are implemented using Rcpp and RcppArmadillo
  • High performance with parallel computation

Installation

You can install the released version of vlad from CRAN with:

install.packages("vlad")

And the development version from GitHub with:

# install.packages("devtools")
devtools::install_github("wittenberg/vlad")

Example

Load libraries:

library("vlad")
library("dplyr")
library("tidyr")
library("ggplot2")

Subset the dataset cardiacsurgery into Phase I (first two years) and Phase II (five years) and estimate a risk model based on phaseI.

data("cardiacsurgery", package = "spcadjust")
cardiacsurgery <- cardiacsurgery %>% rename(s = Parsonnet) %>%
  mutate(y = ifelse(status == 1 & time <= 30, 1, 0),
        phase = factor(ifelse(date < 2*365, "I", "II")))
head(cardiacsurgery)
#>   date time status  s surgeon y phase
#> 1    1   90      0 15       7 0     I
#> 2    1   90      0  9       3 0     I
#> 3    2   90      0  2       5 0     I
#> 4    3   90      0  8       7 0     I
#> 5    3   90      0  7       1 0     I
#> 6    3   90      0 40       1 0     I
phaseI <- filter(cardiacsurgery, phase == "I") %>% select(s, y)
coeff <- round(coef(glm(y ~ s, data = phaseI, family = "binomial")), 3)
print(coeff)
#> (Intercept)           s 
#>       -3.79        0.08

Create VLADs for seven surgeons

By using the estimated risk model coefficients coeff, for each pair of Parsonnet score s and operation outcome values y, the difference between expected and observed outcome is calculated with the function calceo(). Thereafter, differences are cummulated to create the VLAD. This is done for all seven surgeons of the cardiacsurgery dataset. Results are saved to the object vlads7.

vlads7 <- lapply(1:7, function(j){
  Si <- filter(cardiacsurgery, surgeon == j)
  EO <- sapply(seq_along(Si$s), function(i) calceo(df = Si[i, c("s", "y")], coeff = coeff))
  select(Si, surgeon, phase) %>%  mutate(n = 1:length(EO), cEO = cumsum(EO))
}) 

Create Variable life-adjusted Displays for each surgeon from the object vlads7.

vlads7 %>% 
  bind_rows() %>%  
  gather(key = "Surgeon", value = value, c(-n, -surgeon, -phase)) %>%
  ggplot(aes(x = n, y = value, colour = phase, group = Surgeon)) +
    geom_hline(yintercept = 0, colour = "darkgreen", linetype = "dashed") +
    geom_line(size = 1.1) + facet_wrap( ~ surgeon, ncol = 2, scales = "free") +
    labs(x="Patient number n", y="CUSUM E-O") + theme_classic() +
    scale_y_continuous(sec.axis = dup_axis(name = NULL, labels = NULL)) +
    scale_x_continuous(sec.axis = dup_axis(name = NULL, labels = NULL))

Create a VLAD for surgeon 2

S2 <- filter(cardiacsurgery, surgeon == 2) %>% select(phase, s, y)
S2I <- subset(S2, c(phase == "I"))
S2II <- subset(S2, c(phase == "II"))
coeff <- coef(glm(y ~ s, data = S2I, family = "binomial"))
EO <- sapply(1:nrow(S2), function(i) calceo(df = S2[i, c("s", "y")], coeff = coeff))
 
df1 <- select(S2, phase) %>% mutate(n = row_number(), cEO = cumsum(EO))
df2 <- gather(df1, variable, value, c(-n, -phase))
 
p1 <- ggplot(df2, aes(x = n, y = value, colour = phase)) +
  geom_hline(yintercept = 0, linetype = "dashed") + geom_line() + geom_point() + 
  labs(x = "Patient number", y = "CUSUM E-O") + theme_classic() +
  scale_y_continuous(sec.axis = dup_axis(name = NULL, labels = NULL)) +
  scale_x_continuous(sec.axis = dup_axis(name = NULL, labels = NULL))
p1

Compute thresholds of a risk-adjusted CUSUM chart for surgeon 2

Upper and lower control limits of the risk-adjusted CUSUM chart based on log-likelihood ratio statistic can be computed with the function racusum_arl_h_sim(). The implemention uses parallel simulation and a multi-stage search procedure.

# set a random number generator for parallel computations
RNGkind("L'Ecuyer-CMRG")
# number of simulation runs
m <- 10^4
# assign cores
nc <- parallel::detectCores()
# verbose calculation 
UCL_sim <- racusum_crit_sim(L0 = 740, df = S2I[, c("s", "y")], coeff = coeff, m = m, RA = 2, nc = nc, verbose = TRUE)
#> h = 1    ARL = 75.006 
#> h = 2    ARL = 383.3554 
#> h = 3    ARL = 1312.8564 
#> h = 2.9  ARL = 1181.3768 
#> h = 2.8  ARL = 1052.3355 
#> h = 2.7  ARL = 928.6587 
#> h = 2.6  ARL = 822.3063 
#> h = 2.5  ARL = 730.0184 
#> h = 2.51     ARL = 739.8392 
#> h = 2.52     ARL = 747.3412 
#> h = 2.519    ARL = 746.4942 
#> h = 2.518    ARL = 745.9933 
#> h = 2.517    ARL = 745.6282 
#> h = 2.516    ARL = 744.6935 
#> h = 2.515    ARL = 744.0614 
#> h = 2.514    ARL = 743.2221 
#> h = 2.513    ARL = 742.2656 
#> h = 2.512    ARL = 741.7172 
#> h = 2.511    ARL = 741.3188 
#> h = 2.51     ARL = 739.8392 
#> h = 2.5101   ARL = 739.862 
#> h = 2.5102   ARL = 739.862 
#> h = 2.5103   ARL = 740.0315
# quite calculation
LCL_sim <- racusum_crit_sim(L0 = 740, df = S2I[, c("s", "y")], coeff = coeff, m = m, RA = 1/2, nc = nc, verbose = FALSE)
round(cbind(UCL_sim, LCL_sim), 3)
#>      UCL_sim LCL_sim
#> [1,]    2.51   2.281

References

Wittenberg et al. (2018). A simple signaling rule for variable life-adjusted display derived from an equivalent risk-adjusted CUSUM chart

Steiner et al. (2000). Monitoring surgical performance using risk-adjusted cumulative sum charts

Authors

Philipp Wittenberg and Sven Knoth

License

GPL (>= 2)

News

vlad 0.2.0

Changes

  • Renamed functions:

        original name | new name
    

:-----------------------|:------------- cusum_arl_h_sim() | bcusum_crit_sim() cusum_arl_sim() | bcusum_arl_sim() eocusum_arl_h_sim() | eocusum_crit_sim() eocusum_adoc_sim() | eocusum_ad_sim() racusum_arl_h_sim() | racusum_crit_sim() racusum_adoc_sim() | racusum_ad_sim()

Improvements

  • Include Markov Chain approximation for risk-adjusted CUSUM chart based on log-likelihood ratio statistic in function racusum_arl_mc() and racusum_crit_mc()
  • Add option of calculation of optimal k based on mean outcome in function optimal_k()
  • Add function for calculation of CUSUM scores for log-likelihood ratio CUSUM llr_cusum_scores() and E-O CUSUM eo_cusum_scores()
  • Extend unit tests for functions

vlad 0.1.0

Improvements

  • Preparation for CRAN release

Reference manual

It appears you don't have a PDF plugin for this browser. You can click here to download the reference manual.

install.packages("vlad")

0.2.0 by Philipp Wittenberg, 10 months ago


Report a bug at https://github.com/wittenberg/vlad/issues


Browse source code at https://github.com/cran/vlad


Authors: Philipp Wittenberg [aut, cre] , Sven Knoth [aut, ths]


Documentation:   PDF Manual  


GPL (>= 2) license


Imports Rcpp

Suggests dplyr, ggplot2, parallel, rmarkdown, spcadjust, testthat, tidyr

Linking to Rcpp, RcppArmadillo

System requirements: C++11


See at CRAN