Code for Shinyapp for PIRO Model: Predicting AKI in ACLF

server.R

library(shiny)

shinyServer(function(input, output, session) {

  observe({
    input$calc
    isolate({
      creat_cat <- input$creat_cat
      urea <- as.numeric(input$urea)
      bil_cat <- input$bil
      k_cat <- input$pot
      nephrotoxic <- input$neph
      sirs <- input$sirs
      circFailure <- input$c_fail
      ci <- as.numeric(input$ci)

      nd <- data.frame(creat_cat = creat_cat, urea = urea, bil_cat = bil_cat, k_cat = k_cat, nephrotoxic = nephrotoxic, sirs = sirs, circFailure = circFailure, stringsAsFactors = FALSE)

      predy <- rms:::predict.lrm(mod_reduced, nd, se.fit = TRUE)

      z <- qnorm(1 - ((1-ci)/2))

      prob <- plogis(predy[["linear.predictors"]])
      lcl <- plogis(predy[["linear.predictors"]] - z*predy[["se.fit"]])
      ucl <- plogis(predy[["linear.predictors"]] + z*predy[["se.fit"]])

      output$result <- renderUI({
        list(
          h2(strong("Predicted Probability of AKI within 7 days")),
          h2(strong(textOutput("prob_kf"))),
          h2(paste0(100*ci, "% Confidence Interval")),
          h2(textOutput("ci_prob_kf")),
          h4("Linear Predictor"),
          textOutput("lp")
        )
      })

      output$prob_kf <- renderText(paste0(as.character(round(100 * prob, 3)), "%"))
      output$ci_prob_kf <- renderText(paste0(round(100 * lcl, 3), " - ", round(100 * ucl, 3), "%"))
      output$lp <- renderText(as.character(round(predy[["linear.predictors"]], 3)))
    })
  })

  output$formula <- renderPrint(body(Hmisc::Function(mod_reduced)))
  output$cov <- renderTable(mod_reduced$var, digits = 4)
  output$colnames <- renderPrint(c(mod_reduced$Design$name, "kidFailure"))
  output$colvalues <- renderPrint(list(mod_reduced$Design$values, kidFailure = c("N", "Y")))

  # observe({
  #
  #   inFile <- input$val_data
  #   if (is.null(inFile)) {
  #     return(NULL)
  #   }
  #   read.csv(inFile$datapath, stringsAsFactors = FALSE)
  # }, label = "valdata")

  valdata <- reactive({
    inFile <- input$val_data
    if (is.null(inFile)) {
      return(NULL)
    } else {
      read.csv(inFile$datapath)
    }
  })

  output$val_data_details <- renderText({

    if (!is.null(valdata())) {
      nrow <- nrow(valdata())
      cnames <- colnames(valdata())
      cnames_mod <- c(mod_reduced$Design$name, "kidFailure")

      if (length(setdiff(cnames, cnames_mod)) == 0L) {
        same_names <- "column names are same."
      } else {
        same_names <- "column names are different than what is required."
      }

      num_row <- paste0("number of rows are ", nrow, ".")

      paste(same_names, num_row)
    } else {
      return(NULL)
    }

  })

  observe({
    if (!is.null(valdata())) {
      val_val <- val.prob(logit = predict(mod_reduced, valdata()),
                          y = valdata()[["kidFailure"]] == "Y", pl = FALSE, g = 10, logistic.cal = FALSE)

      output$plot_val_raw <- renderPlot({
        val.prob(logit = predict(mod_reduced, valdata()),
                 y = valdata()[["kidFailure"]] == "Y", pl = TRUE, g = 10, logistic.cal = FALSE, xlab = "Predicted Probability Validation Dataset (Unshrunken)")
        #rms:::plot.val.prob(val_val, xlab = "Predicted Probability Validation Dataset (Unshrunken)")
      })

      sd_c_val <- sd_cidx(mod_reduced, data = valdata(), response_var = "kidFailure",
                          true_level = "Y")
      ret_df <- data_frame(`C index: unshrunken` = val_val["C (ROC)"],
                           `95% LCL: unshrunken` = `C index: unshrunken` - 1.96*sd_c_val,
                           `95% UCL: unshrunken` = `C index: unshrunken` + 1.96*sd_c_val)

      val_val_shrnk <- val.prob(logit = shrunked_lp(val_prob_obj = val_val,
                                                    model = mod_reduced, data = valdata()),
                                y = valdata()[["kidFailure"]] == "Y", pl = FALSE,
                                g = 10, logistic.cal = FALSE)

      output$plot_val_shrnk <- renderPlot({
        val.prob(logit = shrunked_lp(val_prob_obj = val_val,
                                     model = mod_reduced, data = valdata()),
                 y = valdata()[["kidFailure"]] == "Y", pl = TRUE,
                 xlab = "Predicted Probability Validation Dataset (Shrunken)",
                 g = 10, logistic.cal = FALSE)
      })

      ret_df <- bind_cols(ret_df, data_frame(`C Index: after shrinkage` = val_val_shrnk["C (ROC)"],
                                             `95% LCL: after shrinkage` = `C Index: after shrinkage` - 1.96*sd_c_val,
                                             `95% UCL: after shrinkage` = `C Index: after shrinkage` + 1.96*sd_c_val))

      output$perf_val_data <- renderTable({
        ret_df
      }, digits = 4)

      output$shrinkage_val <- renderPrint({
        list(Intercept = val_val["Intercept"], Slope = val_val["Slope"])
      })
    }

  })
})

ui.R

library(shiny)

shinyUI(fluidPage(
  titlePanel("PIRO Model for predicting development of Acute Kidney Injury within 7 days in ACLF"),
  #img(src = "ilbs.png", align = "right"),
  tabsetPanel(
    # 1st tab: calculate piro score --------------------
    tabPanel("Probability of developing Kidney Failure", value = "calc_piro",
             sidebarPanel(sliderInput("urea", "Blood Urea (mg/dL) (P component)", 0, 350, 35, 0.1),
                          radioButtons("creat_cat", "Serum Creatinine (mg/dL) (P component)", c("< 2" = "< 2", ">= 2" = ">= 2")),
                          radioButtons("pot", "Serum Potassium (mmoL/L) (P component)", c("< 3" = "< 3", "3 - 5" = "3 - 5", ">= 5" = ">= 5")),
                          radioButtons("bil", "Serum Bilirubin (mg/dL) (P component)", c("< 12" = "< 12", "12 - 30" = "12 - 30", ">= 30" = ">= 30")),
                          hr(),
                          radioButtons("neph", "Nephrotoxicity (I component)", c("No" = "N", "Yes" = "Y")),
                          hr(),
                          radioButtons("sirs", "Systemic Inflammatory Response Syndrome (SIRS) (R component)", c("No" = "N", "Yes" = "Y")),
                          hr(),
                          radioButtons("c_fail", "Circulatory Failure (O component)", c("No" = "N", "Yes" = "Y")),
                          hr(),
                          sliderInput("ci", "Confidence Interval", 0.80, 0.99, 0.95, 0.01),
                          actionButton("calc", "Calculate ...")
             ),
             mainPanel(uiOutput("result"),
                       br(), br(), br(), br(),
                       h5(strong("Source:")),
                       p("Maiwall R, Sarin SK, Kumar S, et al. Development of predisposition, injury, response, organ failure model for predicting acute kidney injury in acute on chronic liver failure. Liver Int, 2017; 00:1-11", a("https://doi.org/10.1111/liv.13443", href = "https://doi.org/10.1111/liv.13443", target = "_blank"))
             )
    ),
    # tab 2: Variable definitions ------------
    tabPanel("Variable definitions", value = "defn",
             h3(strong("Definitions of variables")),
             p(strong("Nephrotoxity: "), "Use of nephrotoxic drugs as ACE inhibitors, ARBs, NSAIDs (aspirin, ibuprofen, etc), radiocontrast agents, antibiotics (aminoglycoside, vancomycin, amphotericin B, sulphonamides), acyclovir, lithium, methotrexate"), br(),
             p(strong("SIRS: "), "Presence of two or more of the following i.e.fever (temperature > 38.3°C) or hypothermia (temperature < 35.6°C); tachycardia (heart rate of > 90 beats/min); tachypnea (respiratory rate > 20 breaths/min or PaCO2 < 32 mm Hg) and alteration of the white cell count > 12,000 cells/mm3, < 4,000 cells/mm3 or > 10% immature neutrophils (bands)"), br(),
             p(strong("Circulatory Failure: "), "Mean arterial pressure (MAP) < 60 mm Hg or a reduction of 40 mmHg in systolic blood pressure from baseline despite adequate fluid resuscitation and cardiac output or need for vasopressors for at least one hour (infusion of dopamine ≥ 5 µg/kg/min or any dose of epinephrine, norepinephrine, dobutamine or terlipressin) to maintain a SBP ≥ 90 mmHg")
    ),
    # tab 3: model specification --------
    tabPanel("Model specification", value = "mod_sp",
             h1("Formula"),
             h4("Linear Predictor (LP) = "),
             p(textOutput("formula")),
             h4("Probability of developing Kidney Failure = "),
             p("1/(1 + exp(- LP))"),
             br(),
             h4("Covariance Matrix"),
             p(tableOutput("cov"))
    ),
    tabPanel("Download Prediction Nomogram", value = "nomo",
             h1("Prediction Nomogram"),
             h4("Introduction"),
             p("Using prediction nomogram is another way to apply the predictive
               model. You can download the nomogram as image file by right
               clicking the image and saving it and use it in your bedside
               by taking printouts of the same."),
             h4("How to use the nomogram"),
             p("For each variable, get points by drawing a vertical line up from
               the value in variable axis to the point axis. Add all the points
               and draw a vertical line from Total Points axis down to the
               Linear Predictor axis and Predicted probability of AKI after 7
               days axis to get the linear predictor and Predicted Probability."),
             img(src = "nomogram.png", align = "center", width = "70%")
             ),
    # tab 4: procedure of shrinkage of coefficients --------
    tabPanel("Validation dataset", value = "valdata",
             h1("Preparing validation dataset"),
             h3(strong("Introduction")),
             p("When a predictive model is made from the training dataset, it
               fits well to the underlying data, as the model is made out
               of it. When the model is tried on an external validation dataset,
               the performance of model usually suffers due to overfitting."),
             p("The performance of predictive model is classified into:"),
             shiny::tags$ol(strong("Discrimination: "), "The ability of model to rank the
                            patients correctly with respect to developing kidney
                            failure."),
             shiny::tags$ol(strong("Calibration: "), "The ability of the model to correctly
                            assign probability of developing kidney failure."),
             h3(strong("Sourcing validation dataset")),
             p("The validation dataset should be a csv file with the following format"),
             p(strong("Column Names")),
             p(verbatimTextOutput("colnames")),
             p(strong("Column Values (for strings)")),
             p(verbatimTextOutput("colvalues")),
             p("creat_cat: in mg/dl, k_cat: in mmol/l, bil_cat: mg/dl."),
             p(strong("Column Values (for numeric)")),
             p("urea: numbers in mg/dl"),
             p(strong("After making the csv file, upload it in the next tab and check for
                      predictive performance. NO MISSING VALUES ARE ACCEPTABLE.
                      Minimum of 200 patients are required for good analysis."))
             ),
    tabPanel("Model performance on validation dataset", value = "mod_perf",
             h2("Upload csv file"),
             fileInput("val_data", "Choose Validation Data", accept = c("text/csv", ".csv")),
             textOutput("val_data_details"),
             h2("Model performance on Validation Dataset"),
             h3("Calibration Plots"),
             h4("Unshrunken data"),
             plotOutput("plot_val_raw", width = "35%", height = "400px"),
             h4("Data after shrinkage of linear predictors"),
             plotOutput("plot_val_shrnk", width = "35%", height = "400px"),
             h3("C Index"),
             tableOutput("perf_val_data"),
             h2("Shrinkage parameters"),
             verbatimTextOutput("shrinkage_val"),
             h3("Shrinking linear predictor (LP)"),
             p("Shrunk LP = Intercept + (Slope * LP)"),
             p("Linear predictor is obtained from tab: Probability of developing Kidney Failure. For
             conversion of LP to probability of developing kidney failure, see
             tab on Model Specification.")
    )

             )))