library(shiny)
library(bslib)
library(mrgsolve)
library(dplyr)
library(ggplot2)

# Sacituzumab govitecan (SG) PopPK model — three analytes
# Source: Sathe et al., Clinical Pharmacokinetics 2024 (doi:10.1007/s40262-024-01366-3)
# SG: 2-compartment IV infusion
# Free SN-38: sequential 2-compartment model driven by first-order SG release
# tAB: 2-compartment model with time-dependent clearance

model_code <- '
$PARAM @annotated
CLSG_ref  : 0.133   : SG clearance (L/h)
QSG_ref   : 0.006   : SG inter-compartmental clearance (L/h)
V1SG_ref  : 2.770   : SG central volume (L)
V2SG_ref  : 0.908   : SG peripheral volume (L)
WTEXP_CL_SG : 0.508 : Weight exponent on SG CL and Q
WTEXP_V_SG  : 0.532 : Weight exponent on SG volumes
ALB_EXP_SG  : -0.355 : Albumin exponent on SG clearance

KREL       : 0.0961 : First-order SG release rate constant (1/h)
CLSN38_ref : 409.0  : SN-38 apparent clearance (L/h)
QSN38_ref  : 247.0  : SN-38 apparent inter-compartmental clearance (L/h)
V1SN38     : 49.0   : SN-38 central volume (L, fixed)
V2SN38     : 2177.0 : SN-38 peripheral volume (L, fixed)
WTEXP_CL_SN38 : 0.500 : Weight exponent on SN-38 CL and Q
F_CONV_cal : 87.0   : Calibrated SG-to-SN-38 conversion (ng/ug)

CLtAB_ref  : 0.016   : tAB clearance (L/h)
QtAB_ref   : 0.010   : tAB inter-compartmental clearance (L/h)
V1tAB_ref  : 3.06    : tAB central volume (L)
V2tAB_ref  : 1.20    : tAB peripheral volume (L)
WTEXP_CL_tAB : 0.372 : Weight exponent on tAB CL and Q
WTEXP_V_tAB  : 0.446 : Weight exponent on tAB volumes
ALB_EXP_tAB  : -0.735 : Albumin exponent on tAB clearance
MAXRED     : 0.167   : Max relative reduction in tAB clearance
KRED       : 6.08e-4 : Rate constant for time-dependent clearance reduction (1/h)

WT         : 70.0    : Body weight (kg)
MALE       : 0.0     : Sex (0=female, 1=male)
ALB        : 38.0    : Baseline albumin (g/L)

$CMT @annotated
CENT_SG    : SG central amount (ug)
PERI_SG    : SG peripheral amount (ug)
CENT_SN38  : Free SN-38 central amount (ng)
PERI_SN38  : Free SN-38 peripheral amount (ng)
CENT_tAB   : tAB central amount (ug)
PERI_tAB   : tAB peripheral amount (ug)

$MAIN
double CLSG = CLSG_ref * pow(WT / 70.0, WTEXP_CL_SG) * pow(ALB / 38.0, ALB_EXP_SG);
double QSG  = QSG_ref  * pow(WT / 70.0, WTEXP_CL_SG);
double V1SG = V1SG_ref * pow(WT / 70.0, WTEXP_V_SG);
double V2SG = V2SG_ref * pow(WT / 70.0, WTEXP_V_SG);

double CLSN38 = CLSN38_ref * pow(WT / 70.0, WTEXP_CL_SN38);
double QSN38  = QSN38_ref  * pow(WT / 70.0, WTEXP_CL_SN38);

double CLtAB_base = CLtAB_ref * pow(WT / 70.0, WTEXP_CL_tAB) * pow(ALB / 38.0, ALB_EXP_tAB);
double QtAB = QtAB_ref * pow(WT / 70.0, WTEXP_CL_tAB);
double V1tAB = V1tAB_ref * pow(WT / 70.0, WTEXP_V_tAB) * (MALE > 0.5 ? 1.121 : 1.0);
double V2tAB = V2tAB_ref * pow(WT / 70.0, WTEXP_V_tAB);

$ODE
double CSG1 = CENT_SG / V1SG;
double CSG2 = PERI_SG / V2SG;
double release_sn38 = KREL * CENT_SG * F_CONV_cal;
double CSN1 = CENT_SN38 / V1SN38;
double CSN2 = PERI_SN38 / V2SN38;
double CLtAB = CLtAB_base * (1.0 - MAXRED * (1.0 - exp(-KRED * SOLVERTIME)));
double CtAB1 = CENT_tAB / V1tAB;
double CtAB2 = PERI_tAB / V2tAB;

dxdt_CENT_SG   = -CLSG * CSG1 - QSG * (CSG1 - CSG2);
dxdt_PERI_SG   =                  QSG * (CSG1 - CSG2);
dxdt_CENT_SN38 = release_sn38 - CLSN38 * CSN1 - QSN38 * (CSN1 - CSN2);
dxdt_PERI_SN38 =                                QSN38 * (CSN1 - CSN2);
dxdt_CENT_tAB  = -CLtAB * CtAB1 - QtAB * (CtAB1 - CtAB2);
dxdt_PERI_tAB  =                    QtAB * (CtAB1 - CtAB2);

$TABLE
double CP_SG   = fmax((CENT_SG / V1SG) / 1000.0, 0.0);
double CP_SN38 = fmax((CENT_SN38 / V1SN38) / 1000.0, 0.0);
double CP_tAB  = fmax((CENT_tAB / V1tAB) / 1000.0, 0.0);

$CAPTURE CP_SG CP_SN38 CP_tAB
'

mod <- mcode("sg_popk_project_364", model_code, quiet = TRUE)

app_theme <- bs_theme(
  version = 5,
  bootswatch = "flatly",
  primary = "#8b5cf6"
) |> bs_add_rules("
  .metric-card {
    background: #f8f9fa;
    border-radius: 8px;
    padding: 15px;
    margin: 5px;
    text-align: center;
    border: 1px solid #dee2e6;
  }
  .metric-value {
    font-size: 24px;
    font-weight: bold;
    color: #2c3e50;
  }
  .metric-label {
    font-size: 12px;
    color: #7f8c8d;
  }
  .metric-success .metric-value { color: #10b981; }
  .metric-warning .metric-value { color: #f59e0b; }
  .metric-primary .metric-value { color: #8b5cf6; }
  .metric-info .metric-value { color: #0dcaf0; }
  .ref-box {
    background: #f8f9fa;
    border-left: 4px solid #8b5cf6;
    padding: 12px 16px;
    border-radius: 4px;
    font-size: 13px;
  }
")

ui <- page_sidebar(
  title = "Sacituzumab Govitecan (Trodelvy) PK Simulator",
  theme = app_theme,
  sidebar = sidebar(
    title = "Simulation Settings",
    width = 320,
    h6("Dosing", class = "text-muted mb-1"),
    sliderInput("dose_mgkg", "Dose (mg/kg)", min = 6, max = 14, value = 10, step = 0.5),
    sliderInput("n_cycles", "Cycles to Simulate", min = 1, max = 6, value = 3, step = 1),
    p(class = "small text-muted mt-0", "3-hour IV infusion on Days 1 and 8 of each 21-day cycle."),
    hr(class = "my-2"),
    h6("Patient Characteristics", class = "text-muted mb-1"),
    sliderInput("wt", "Body Weight (kg)", min = 37, max = 140, value = 70, step = 1),
    radioButtons("sex", "Sex", choices = c("Female" = 0, "Male" = 1), selected = 0, inline = TRUE),
    sliderInput("alb", "Baseline Albumin (g/L)", min = 19, max = 50, value = 38, step = 1),
    hr(class = "my-2"),
    selectInput(
      "analyte",
      "Displayed Analyte",
      choices = c(
        "Sacituzumab Govitecan (SG)" = "CP_SG",
        "Free SN-38" = "CP_SN38",
        "Total Antibody (tAB)" = "CP_tAB"
      ),
      selected = "CP_SG"
    ),
    checkboxInput("log_scale", "Log scale", value = FALSE)
  ),
  layout_column_wrap(
    width = 1 / 4,
    fill = FALSE,
    div(
      class = "metric-card metric-success",
      div(class = "metric-value", textOutput("cmax")),
      div(class = "metric-label", textOutput("cmax_label"))
    ),
    div(
      class = "metric-card metric-warning",
      div(class = "metric-value", textOutput("ctrough")),
      div(class = "metric-label", textOutput("ctrough_label"))
    ),
    div(
      class = "metric-card metric-primary",
      div(class = "metric-value", textOutput("auc")),
      div(class = "metric-label", textOutput("auc_label"))
    ),
    div(
      class = "metric-card metric-info",
      div(class = "metric-value", textOutput("thalf")),
      div(class = "metric-label", textOutput("thalf_label"))
    )
  ),
  card(
    full_screen = TRUE,
    height = "500px",
    card_header(textOutput("plot_title")),
    card_body(plotOutput("pk_plot", height = "430px"))
  ),
  navset_card_underline(
    title = "Study Details",
    nav_panel(
      "Model Information",
      div(
        class = "p-3",
        tags$h5("Population PK Structure"),
        tags$ul(
          tags$li("SG: two-compartment IV infusion model with linear elimination."),
          tags$li("Free SN-38: sequential two-compartment model generated from SG via first-order release (KREL = 0.0961 h^-1)."),
          tags$li("Total antibody: two-compartment model with time-dependent clearance reduction over repeated dosing."),
          tags$li("Body weight scales clearance and distribution for all three analytes; albumin affects SG and tAB clearance.")
        ),
        tags$h5("Key Parameter Estimates"),
        tags$table(
          class = "table table-sm table-striped",
          tags$thead(tags$tr(tags$th("Parameter"), tags$th("Estimate"), tags$th("Unit"))),
          tags$tbody(
            tags$tr(tags$td("CLSG"), tags$td("0.133"), tags$td("L/h")),
            tags$tr(tags$td("QSG"), tags$td("0.006"), tags$td("L/h")),
            tags$tr(tags$td("V1SG"), tags$td("2.77"), tags$td("L")),
            tags$tr(tags$td("V2SG"), tags$td("0.908"), tags$td("L")),
            tags$tr(tags$td("KREL"), tags$td("0.0961"), tags$td("1/h")),
            tags$tr(tags$td("CLSN38/F"), tags$td("409"), tags$td("L/h")),
            tags$tr(tags$td("QSN38/F"), tags$td("247"), tags$td("L/h")),
            tags$tr(tags$td("V1SN38/F"), tags$td("49"), tags$td("L (fixed)")),
            tags$tr(tags$td("V2SN38/F"), tags$td("2177"), tags$td("L (fixed)")),
            tags$tr(tags$td("CLtAB"), tags$td("0.016"), tags$td("L/h")),
            tags$tr(tags$td("QtAB"), tags$td("0.010"), tags$td("L/h")),
            tags$tr(tags$td("V1tAB"), tags$td("3.06"), tags$td("L")),
            tags$tr(tags$td("V2tAB"), tags$td("1.20"), tags$td("L"))
          )
        ),
        div(
          class = "ref-box",
          strong("Study population: "),
          "529 patients from IMMU-132-01 and ASCENT. Approved regimen: 10 mg/kg IV on Days 1 and 8 every 21 days."
        )
      )
    ),
    nav_panel(
      "References",
      div(
        class = "p-3",
        tags$ol(
          tags$li(
            "Sathe AG, Singh I, Singh P, et al. ",
            tags$em("Population Pharmacokinetics of Sacituzumab Govitecan in Patients with Metastatic Triple-Negative Breast Cancer and Other Solid Tumors."),
            " Clin Pharmacokinet. 2024;63:669-681. ",
            tags$a("DOI:10.1007/s40262-024-01366-3", href = "https://doi.org/10.1007/s40262-024-01366-3", target = "_blank")
          ),
          tags$li(
            "Trodelvy prescribing information. ",
            tags$a("FDA label", href = "https://www.accessdata.fda.gov/drugsatfda_docs/label/2025/761115s059lbl.pdf", target = "_blank")
          ),
          tags$li(
            "Klein CE, Gupta E, Reid JM, et al. Population pharmacokinetic model for irinotecan and SN-38 disposition used as literature support for fixed SN-38 distribution volumes."
          )
        ),
        div(
          class = "ref-box",
          tags$ul(
            tags$li(tags$strong("Drug class:"), " TROP2-directed antibody-drug conjugate"),
            tags$li(tags$strong("Therapeutic area:"), " Oncology"),
            tags$li(tags$strong("Route:"), " Intravenous infusion"),
            tags$li(tags$strong("Indications:"), " Metastatic TNBC, HR+/HER2-negative metastatic breast cancer, and metastatic urothelial cancer.")
          )
        )
      )
    )
  ),
  div(
    class = "text-center text-muted small py-3 mt-2 border-top",
    HTML('Built with <a href="https://pkpdbuilder.com" target="_blank" style="color:#8b5cf6;text-decoration:none;">PKPDBuilder</a> &bull; Parameters from Sathe et al. 2024 &bull; For research and educational use only')
  )
)

analyte_meta <- function(code) {
  switch(
    code,
    CP_SG = list(name = "Sacituzumab Govitecan (SG)", unit = "ug/mL", color = "#8b5cf6"),
    CP_SN38 = list(name = "Free SN-38", unit = "ng/mL", color = "#f97316"),
    CP_tAB = list(name = "Total Antibody (tAB)", unit = "ug/mL", color = "#10b981")
  )
}

trap_auc <- function(time, conc) {
  keep <- is.finite(time) & is.finite(conc)
  time <- time[keep]
  conc <- conc[keep]
  if (length(time) < 2) {
    return(NA_real_)
  }
  sum(diff(time) * (head(conc, -1) + tail(conc, -1)) / 2)
}

terminal_half_life <- function(analyte, wt, alb, sex_code) {
  if (identical(analyte, "CP_SG")) {
    cl <- 0.133 * (wt / 70)^0.508 * (alb / 38)^(-0.355)
    q <- 0.006 * (wt / 70)^0.508
    v1 <- 2.77 * (wt / 70)^0.532
    v2 <- 0.908 * (wt / 70)^0.532
  } else if (identical(analyte, "CP_SN38")) {
    cl <- 409 * (wt / 70)^0.5
    q <- 247 * (wt / 70)^0.5
    v1 <- 49
    v2 <- 2177
  } else {
    cl <- 0.016 * (wt / 70)^0.372 * (alb / 38)^(-0.735)
    q <- 0.010 * (wt / 70)^0.372
    v1 <- 3.06 * (wt / 70)^0.446 * ifelse(sex_code > 0.5, 1.121, 1)
    v2 <- 1.20 * (wt / 70)^0.446
  }

  k10 <- cl / v1
  k12 <- q / v1
  k21 <- q / v2
  disc <- (k10 + k12 + k21)^2 - 4 * k10 * k21
  beta <- 0.5 * ((k10 + k12 + k21) - sqrt(max(disc, 0)))
  if (!is.finite(beta) || beta <= 0) {
    return(NA_real_)
  }
  log(2) / beta
}

server <- function(input, output, session) {
  sim_data <- reactive({
    shiny::req(input$dose_mgkg, input$n_cycles, input$wt, input$sex, input$alb)

    ii_h <- 504
    tinf_h <- 3
    n_cycles <- as.integer(input$n_cycles)
    sim_end <- n_cycles * ii_h
    dose_ug <- input$dose_mgkg * input$wt * 1000

    dose_times <- sort(c(
      seq(0, by = ii_h, length.out = n_cycles),
      seq(168, by = ii_h, length.out = n_cycles)
    ))
    dose_times <- dose_times[dose_times <= sim_end]

    ev_list <- lapply(
      dose_times,
      function(t) {
        ev(
          amt = dose_ug,
          cmt = c("CENT_SG", "CENT_tAB"),
          time = t,
          rate = dose_ug / tinf_h
        )
      }
    )
    events <- do.call(c, ev_list)

    mod %>%
      param(WT = input$wt, MALE = as.numeric(input$sex), ALB = input$alb) %>%
      mrgsim(events = events, end = sim_end, delta = 1, carry_out = "evid") %>%
      as.data.frame() |>
      dplyr::filter(.data$evid == 0 | .data$time == 0) |>
      mutate(
        time_days = time / 24,
        cycle = floor(time / ii_h) + 1
      )
  })

  selected_meta <- reactive({
    analyte_meta(input$analyte)
  })

  metrics <- reactive({
    df <- sim_data()
    yvar <- input$analyte
    cycle1 <- df |>
      dplyr::filter(time <= 504)
    day8_trough <- cycle1 |>
      dplyr::filter(time < 168) |>
      dplyr::slice_max(order_by = time, n = 1)

    list(
      cmax = max(cycle1[[yvar]], na.rm = TRUE),
      ctrough = if (nrow(day8_trough) > 0) day8_trough[[yvar]][1] else NA_real_,
      auc = trap_auc(cycle1$time, cycle1[[yvar]]),
      thalf = terminal_half_life(input$analyte, input$wt, input$alb, as.numeric(input$sex))
    )
  })

  output$cmax <- renderText({
    sprintf("%.2f", metrics()$cmax)
  })
  output$ctrough <- renderText({
    sprintf("%.2f", metrics()$ctrough)
  })
  output$auc <- renderText({
    sprintf("%.1f", metrics()$auc)
  })
  output$thalf <- renderText({
    th <- metrics()$thalf
    if (is.na(th)) "NA" else sprintf("%.1f", th / 24)
  })

  output$cmax_label <- renderText({
    sprintf("Cmax (%s, Cycle 1)", selected_meta()$unit)
  })
  output$ctrough_label <- renderText({
    sprintf("Ctrough at Day 8 (%s)", selected_meta()$unit)
  })
  output$auc_label <- renderText({
    sprintf("AUC0-21d (%s·h/mL)", selected_meta()$unit)
  })
  output$thalf_label <- renderText({
    "Approx terminal t1/2 (days)"
  })

  output$plot_title <- renderText({
    sprintf("%s Concentration-Time Profile", selected_meta()$name)
  })

  output$pk_plot <- renderPlot({
    df <- sim_data()
    yvar <- input$analyte
    meta <- selected_meta()
    vlines <- data.frame(x = sort(unique(c(seq(7, by = 21, length.out = input$n_cycles), seq(21, by = 21, length.out = input$n_cycles - 1)))))

    plot_df <- df
    if (isTRUE(input$log_scale)) {
      plot_df <- plot_df |>
        dplyr::filter(.data[[yvar]] > 0)
    }

    p <- ggplot(plot_df, aes(x = .data$time_days, y = .data[[yvar]])) +
      geom_line(color = meta$color, linewidth = 1) +
      geom_vline(data = vlines, aes(xintercept = .data$x), linetype = "dashed", color = "#d1d5db", linewidth = 0.4) +
      labs(
        x = "Time (days)",
        y = sprintf("%s (%s)", meta$name, meta$unit)
      ) +
      theme_bw(base_size = 13) +
      theme(
        panel.grid.minor = element_blank(),
        panel.border = element_rect(color = "#dee2e6"),
        axis.title = element_text(size = 11),
        axis.text = element_text(size = 10)
      )

    if (isTRUE(input$log_scale)) {
      p <- p + scale_y_log10()
    }

    p
  })
}

shinyApp(ui, server)