# Colonoscopy follow-up measure example

#### Kenneth Nieser

This vignette includes an example of hospital profiling based on a measure of whether patients receive appropriate recommendations for a follow-up colonoscopy.

    library(QualityMeasure)
    #> Loading required package: ggplot2
    #> Loading required package: doParallel
    #> Loading required package: foreach
    #> Loading required package: iterators
    #> Loading required package: parallel
    #> Loading required package: lme4
    #> Loading required package: Matrix
    #> Package 'QualityMeasure' version 2.0.2
    #>
    #> If you have issues or feedback, please email me at nieser@stanford.edu, so I can make this package better!
    #>
    #> For more info about this package, see https://github.com/knieser/quality_measure_reliability.

First, we’ll load the dataset included with the `QualityMeasure` package.

    df <- colonoscopy
    knitr::kable(head(df), 'simple')

entity | p | n | x
--- | --- | --- | ---
010001 | 0.72 | 29 | 21
010005 | 1.00 | 210 | 210
010006 | 0.86 | 85 | 73
010007 | 0.65 | 52 | 34
010008 | 0.67 | 12 | 8
010011 | 0.83 | 42 | 35

Next, we will calculate reliability using the Beta-Binomial method for aggregated data.

    BB.results <- calcBetaBin(df = df, df.aggregate = T, n = 'n', x = 'x')
    #>  Currently, Beta-Binomial reliability estimates do not account for risk-adjustment (even if you specified a model). Updates to this function to account for risk-adjustment are in progress.
    #> Note that aggregated data are being used, so Beta-Binomial reliability estimates with random effects predictions cannot be calculated.

Beta-Binomial parameter estimates are: alpha = 4.516 and beta = 0.43. The between-entity variance in rates is 0.013.

Below is a summary of the distribution of entity-level reliability estimates.

    summary(BB.results$est.BB)
    #>    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
    #>  0.6898  0.9010  0.9381  0.9171  0.9546  0.9975

We can also plot entity-level reliability results by sample size.

    plot.df <- data.frame(
      n = df$n,
      rel = BB.results$est.BB
    )

    fig <- ggplot(data = plot.df, aes(x = n, y = rel)) +
      geom_point(size = 3) +
      geom_hline(yintercept = median(BB.results$est.BB), linetype = 'dashed', col = 'red', linewidth = 2) +
      annotate('text', x = 1700, y = 0.92, label = 'Median reliability', size = 6, col = 'red') +
      xlab('Entity sample size') +
      ylab('Reliability') +
      theme_classic() +
      theme(
        panel.grid.major = element_line(),
        panel.grid.minor = element_line(),
        axis.text = element_text(size = 16),
        axis.ticks.length = unit(.25, 'cm'),
        axis.title = element_text(size = 18, face = 'bold')
      )
    fig