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

# ── mrgsolve PK/PD model ─────────────────────────────────────────────────────
model_code <- '
$PARAM @annotated
Vmax    : 207.0 : Max DPD-mediated elimination rate (mg/h, at BSA 1.73 m2)
Km      : 0.70  : Michaelis-Menten half-saturation constant (mg/L)
V1      : 10.0  : Central volume of distribution (L)
V2      : 15.0  : Peripheral volume of distribution (L)
Q       : 5.0   : Inter-compartmental clearance (L/h)
BSA     : 1.73  : Body surface area (m2)
DPD_ACT : 1.0   : DPD activity fraction (1.0 = normal; 0.5 = deficiency)
Circ0   : 5000  : Baseline ANC (cells/uL)
Slope   : 0.009 : Linear drug-effect slope on proliferating cells (L/mg)
MTT     : 120   : Mean transit time of neutrophil maturation (h)
Gamma   : 0.27  : Feedback power exponent

$CMT @annotated
CENT  : Central compartment (mg)
PERI  : Peripheral compartment (mg)
PROL  : Proliferating cell precursors
TRAN1 : Maturation transit compartment 1
TRAN2 : Maturation transit compartment 2
TRAN3 : Maturation transit compartment 3
CIRC  : Circulating neutrophils (cells/uL)

$MAIN
double Vmaxi = Vmax * DPD_ACT * (BSA / 1.73);
double ktr   = 4.0 / MTT;

PROL_0  = Circ0;
TRAN1_0 = Circ0;
TRAN2_0 = Circ0;
TRAN3_0 = Circ0;
CIRC_0  = Circ0;

$ODE
double CPo      = CENT / V1;
double CPp      = PERI / V2;
double CIRCsafe = (CIRC > 10.0) ? CIRC : 10.0;

dxdt_CENT = -(Vmaxi * CPo) / (Km + CPo) - Q * (CPo - CPp);
dxdt_PERI = Q * (CPo - CPp);

double Edrug    = Slope * CPo;
if (Edrug > 0.95) Edrug = 0.95;
double feedback = pow(Circ0 / CIRCsafe, Gamma);

dxdt_PROL  = ktr * PROL  * (1.0 - Edrug) * feedback - ktr * PROL;
dxdt_TRAN1 = ktr * PROL  - ktr * TRAN1;
dxdt_TRAN2 = ktr * TRAN1 - ktr * TRAN2;
dxdt_TRAN3 = ktr * TRAN2 - ktr * TRAN3;
dxdt_CIRC  = ktr * TRAN3 - ktr * CIRC;

$TABLE
double CP  = CENT / V1;
double ANC = CIRC;

$CAPTURE @annotated
CP  : Plasma 5-FU concentration (mg/L)
ANC : Absolute neutrophil count (cells/uL)
'

mod <- mcode("fu5_pkpd", model_code)

# ── Theme ────────────────────────────────────────────────────────────────────
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; margin-top: 4px; }
  .metric-success .metric-value { color: #10b981; }
  .metric-warning .metric-value { color: #f59e0b; }
  .metric-primary .metric-value { color: #8b5cf6; }
  .metric-info    .metric-value { color: #0dcaf0; }
  .auc-badge { display: inline-block; padding: 2px 9px; border-radius: 20px;
               font-size: 10px; font-weight: 700; margin-left: 4px; }
  .auc-ok   { background: #d1fae5; color: #065f46; }
  .auc-low  { background: #fef3c7; color: #92400e; }
  .auc-high { background: #fee2e2; color: #991b1b; }
  .ref-box { background: #f0f4ff; border-left: 4px solid #8b5cf6;
             padding: 14px 18px; border-radius: 4px; margin-top: 10px; font-size: 13px; }
  .footer-bar { text-align: center; padding: 16px; margin-top: 8px;
                border-top: 1px solid #e9ecef; color: #6c757d; font-size: 12px; }
")

# ── Regimen presets ──────────────────────────────────────────────────────────
regimen_presets <- list(
  "FOLFOX6 (q2w)"      = list(dose = 2400, dur = 46, ii = 14),
  "FOLFIRI (q2w)"      = list(dose = 2400, dur = 46, ii = 14),
  "LV5FU2 (q2w)"       = list(dose = 600,  dur = 22, ii = 14),
  "Continuous (5-day)" = list(dose = 1000, dur = 120, ii = 28),
  "Custom"             = list(dose = 2000, dur = 46,  ii = 14)
)

# ── UI ───────────────────────────────────────────────────────────────────────
ui <- page_sidebar(
  title = "5-Fluorouracil (5-FU) PK/PD Simulator",
  theme = app_theme,

  sidebar = sidebar(
    title = "Simulation Settings", width = 340,

    h6("Chemotherapy Regimen"),
    selectInput("regimen", NULL,
                choices = names(regimen_presets),
                selected = "FOLFOX6 (q2w)"),
    hr(),

    h6("Dosing Parameters"),
    sliderInput("bsa",         "Body Surface Area (m\u00b2)", min = 1.2, max = 2.4,  value = 1.73, step = 0.01),
    sliderInput("dose_per_m2", "5-FU Dose (mg/m\u00b2)",     min = 400, max = 4000, value = 2400, step = 100),
    sliderInput("inf_dur",     "Infusion Duration (h)",       min = 2,   max = 120,  value = 46,   step = 1),
    sliderInput("n_cycles",    "Number of Cycles",            min = 1,   max = 6,    value = 3,    step = 1),
    sliderInput("cycle_ii",    "Cycle Interval (days)",       min = 7,   max = 28,   value = 14,   step = 7),
    hr(),

    h6("Pharmacogenomics"),
    sliderInput("dpd_act", "DPD Activity (fraction of normal)",
                min = 0.3, max = 1.5, value = 1.0, step = 0.05),
    helpText("\u26a0\ufe0f DPD deficiency (<0.7) \u2192 severe toxicity risk at standard doses"),
    hr(),
    checkboxInput("log_scale", "Log scale (concentration axis)", 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", "Cmax (mg/L)")),
    div(class = "metric-card metric-primary",
        div(class = "metric-value", textOutput("auc_c1")),
        div(class = "metric-label", htmlOutput("auc_label"))),
    div(class = "metric-card metric-warning",
        div(class = "metric-value", textOutput("nadir_anc")),
        div(class = "metric-label", "Nadir ANC (cells/\u03bcL)")),
    div(class = "metric-card metric-info",
        div(class = "metric-value", textOutput("t_half")),
        div(class = "metric-label", "Effective t\u00bd (h)"))
  ),

  navset_card_underline(
    full_screen = TRUE,

    nav_panel("Pharmacokinetics",
      plotOutput("pkPlot", height = "470px")),

    nav_panel("Myelosuppression (Friberg)",
      plotOutput("ancPlot", height = "470px")),

    nav_panel("Model Information", div(style = "padding:12px;",
      h4("5-Fluorouracil (5-FU) \u2014 2-CMT Michaelis-Menten + Friberg PK/PD Model"),
      tags$hr(),
      h5("Pharmacokinetic Model"),
      tags$p("Structure: 2-compartment with Michaelis-Menten (saturable) elimination.
             5-FU is primarily catabolized by hepatic DPD (\u226580% of elimination).
             Capacity-limited kinetics mean small dose changes cause disproportionate AUC changes
             \u2014 the pharmacological basis for TDM."),
      tags$table(class="table table-sm table-bordered",
        tags$thead(tags$tr(tags$th("Parameter"), tags$th("Typical Value"), tags$th("Description"))),
        tags$tbody(
          tags$tr(tags$td("Vmax"), tags$td("207 mg/h per 1.73 m\u00b2"), tags$td("Max DPD-mediated elimination")),
          tags$tr(tags$td("Km"),   tags$td("0.70 mg/L"),                  tags$td("Half-saturation constant")),
          tags$tr(tags$td("V1"),   tags$td("10 L"),                       tags$td("Central volume")),
          tags$tr(tags$td("V2"),   tags$td("15 L"),                       tags$td("Peripheral volume")),
          tags$tr(tags$td("Q"),    tags$td("5 L/h"),                      tags$td("Inter-compartmental CL"))
        )
      ),
      tags$hr(),
      h5("Myelosuppression Model (Friberg Semi-Physiological)"),
      tags$p("Proliferating cells \u2192 3 maturation transit compartments \u2192 circulating neutrophils.
             Drug linearly inhibits proliferation. Feedback from circulating cell count
             drives compensatory response."),
      tags$table(class="table table-sm table-bordered",
        tags$thead(tags$tr(tags$th("Parameter"), tags$th("Value"), tags$th("Description"))),
        tags$tbody(
          tags$tr(tags$td("Circ\u2080"), tags$td("5,000 cells/\u03bcL"), tags$td("Baseline ANC")),
          tags$tr(tags$td("MTT"),        tags$td("120 h"),                tags$td("Mean transit time")),
          tags$tr(tags$td("Slope"),      tags$td("0.009 L/mg"),           tags$td("Linear drug effect")),
          tags$tr(tags$td("\u03b3"),     tags$td("0.27"),                 tags$td("Feedback exponent"))
        )
      ),
      tags$hr(),
      h5("Therapeutic Drug Monitoring (TDM)"),
      tags$p(HTML("<strong>Target AUC: 20\u201330 mg\u00b7h/L</strong> per infusion cycle.<br>
             IATDMCT (2018) strongly recommends TDM of 5-FU in clinical practice.")),
      tags$ul(
        tags$li(tags$strong("< 20 mg\u00b7h/L:"), " Underdosing \u2014 dose escalation recommended"),
        tags$li(tags$strong("20\u201330 mg\u00b7h/L:"), " Therapeutic \u2014 optimal exposure-efficacy"),
        tags$li(tags$strong("> 30 mg\u00b7h/L:"), " Overdosing \u2014 excess toxicity, dose reduction required")
      )
    )),

    nav_panel("References", div(class = "ref-box",
      tags$h5("\U0001f4da Primary References"),
      tags$ol(
        tags$li(tags$strong("Kobuchi S & Ito Y. (2020)"),
          " Pharmacometrics approaches for personalized 5-FU-based chemotherapy (Review). ",
          tags$em("Anticancer Research"), " 40:6585\u20136597."),
        tags$li(tags$strong("Terret C et al. (2000)"),
          " Dose and time dependencies of 5-fluorouracil pharmacokinetics. ",
          tags$em("Clin Pharmacol Ther"), " 68:270\u2013279."),
        tags$li(tags$strong("Friberg LE et al. (2002)"),
          " Model of chemotherapy-induced myelosuppression. ",
          tags$em("J Clin Oncol"), " 20:4713\u20134721."),
        tags$li(tags$strong("Gamelin E et al. (1998)"),
          " Individual fluorouracil dose adjustment based on pharmacokinetic follow-up. ",
          tags$em("J Clin Oncol"), " 16:2099\u20132110."),
        tags$li(tags$strong("Kobuchi S et al. (2016)"),
          " Population PK-PD model for myelosuppression in 5-FU-treated rats. ",
          tags$em("Eur J Drug Metab Pharmacokinet"), " 42:703\u2013714.")
      ),
      tags$h5("\U0001f489 Clinical Context"),
      tags$ul(
        tags$li(tags$strong("Drug class:"), " Fluoropyrimidine antimetabolite"),
        tags$li(tags$strong("Indications:"), " Colorectal, gastric, breast, pancreatic cancer"),
        tags$li(tags$strong("Key toxicities:"), " Myelosuppression, mucositis, hand-foot syndrome, diarrhea"),
        tags$li(tags$strong("Pharmacogenomics:"), " DPYD variants predict DPD deficiency \u2014 screening recommended before treatment"),
        tags$li(tags$strong("TDM:"), " IATDMCT strongly recommends TDM in clinical practice (2018)")
      )
    ))
  ),

  div(class = "footer-bar",
    "Powered by ",
    tags$a(href = "https://www.pkpdbuilder.com", target = "_blank",
           style = "color: #8b5cf6; font-weight: 600;", "PKPDBuilder.com"),
    " \u2022 Based on: Kobuchi & Ito (2020), Anticancer Research 40:6585\u20136597", br(),
    tags$span(style = "font-size: 10px;",
      "For research and educational use only. Not for clinical decision-making.")
  )
)

# ── Server ────────────────────────────────────────────────────────────────────
server <- function(input, output, session) {

  # Regimen presets
  observeEvent(input$regimen, {
    p <- regimen_presets[[input$regimen]]
    updateSliderInput(session, "dose_per_m2", value = p$dose)
    updateSliderInput(session, "inf_dur",     value = p$dur)
    updateSliderInput(session, "cycle_ii",    value = p$ii)
  }, ignoreInit = TRUE)

  sim_data <- reactive({
    shiny::req(input$dose_per_m2, input$bsa, input$inf_dur, input$n_cycles, input$cycle_ii, input$dpd_act)

    dose_total <- input$dose_per_m2 * input$bsa
    cycle_h    <- input$cycle_ii * 24
    n_cyc      <- as.integer(input$n_cycles)
    dur_h      <- input$inf_dur
    rate_val   <- dose_total / dur_h

    # Build event table: one infusion per cycle
    ev_list <- lapply(seq_len(n_cyc) - 1, function(i) {
      ev(amt = dose_total, rate = rate_val, cmt = 1, time = i * cycle_h)
    })
    ev_combined <- Reduce("+", ev_list)

    end_time <- (n_cyc - 1) * cycle_h + max(dur_h + 72, cycle_h)

    mod |>
      param(BSA = input$bsa, DPD_ACT = input$dpd_act) |>
      ev(ev_combined) |>
      mrgsim(end = end_time, delta = 0.25) |>
      as.data.frame()
  })

  # PK metrics ---------------------------------------------------------------
  output$cmax <- renderText({
    d <- sim_data()
    sprintf("%.2f", max(d$CP, na.rm = TRUE))
  })

  output$auc_c1 <- renderText({
    d <- sim_data()
    dur_h <- input$inf_dur
    d1 <- dplyr::filter(d, time >= 0, time <= dur_h + 48)
    auc <- sum(diff(d1$time) * (head(d1$CP, -1) + tail(d1$CP, -1)) / 2)
    sprintf("%.1f", auc)
  })

  output$auc_label <- renderUI({
    d <- sim_data()
    dur_h <- input$inf_dur
    d1 <- dplyr::filter(d, time >= 0, time <= dur_h + 48)
    auc <- sum(diff(d1$time) * (head(d1$CP, -1) + tail(d1$CP, -1)) / 2)
    badge_class <- if (auc < 20) "auc-badge auc-low" else if (auc > 30) "auc-badge auc-high" else "auc-badge auc-ok"
    badge_label <- if (auc < 20) "Underdosed" else if (auc > 30) "Overdosed" else "Therapeutic"
    tagList(
      span("AUC Cycle 1 (mg\u00b7h/L)"),
      span(class = badge_class, badge_label)
    )
  })

  output$nadir_anc <- renderText({
    d <- sim_data()
    sprintf("%.0f", min(d$ANC, na.rm = TRUE))
  })

  output$t_half <- renderText({
    d <- sim_data()
    cp <- d$CP[d$CP > 0.01]
    if (length(cp) < 2) return("—")
    # Estimate from terminal slope during last cycle decay
    dur_h    <- input$inf_dur
    cycle_h  <- input$cycle_ii * 24
    last_cyc <- (input$n_cycles - 1) * cycle_h
    decay    <- dplyr::filter(d, time > last_cyc + dur_h, time < last_cyc + dur_h + 48, CP > 0.001)
    if (nrow(decay) < 4) return("—")
    fit <- tryCatch(
      lm(log(CP) ~ time, data = decay),
      error = function(e) NULL
    )
    if (is.null(fit)) return("—")
    k <- -coef(fit)[["time"]]
    if (k <= 0) return("—")
    sprintf("%.1f", log(2) / k)
  })

  # PK Plot ------------------------------------------------------------------
  output$pkPlot <- renderPlot({
    d <- sim_data()
    target_lo <- 20; target_hi <- 30

    p <- ggplot(d, aes(x = time, y = CP)) +
      annotate("rect", xmin = -Inf, xmax = Inf,
               ymin = target_lo, ymax = target_hi,
               fill = "#10b981", alpha = 0.12) +
      geom_hline(yintercept = c(target_lo, target_hi),
                 linetype = "dashed", color = "#10b981", alpha = 0.6) +
      geom_line(color = "#8b5cf6", linewidth = 0.8) +
      annotate("text", x = max(d$time) * 0.02, y = target_hi * 1.05,
               label = "TDM target AUC zone: 20\u201330 mg\u00b7h/L",
               hjust = 0, size = 3.5, color = "#10b981") +
      labs(x = "Time (h)",
           y = "5-FU Concentration (mg/L)",
           title = sprintf("5-FU %.0f mg/m\u00b2 over %gh \u00d7 %d cycles (BSA %.2f m\u00b2, DPD %.0f%%)",
                           input$dose_per_m2, input$inf_dur,
                           input$n_cycles, input$bsa, input$dpd_act * 100)) +
      theme_minimal(base_size = 14) +
      theme(plot.title = element_text(size = 12))

    if (input$log_scale) {
      p <- p + scale_y_log10(
        limits = c(1e-3, max(d$CP[d$CP > 0], na.rm = TRUE) * 2)
      )
    }
    p
  })

  # ANC Plot -----------------------------------------------------------------
  output$ancPlot <- renderPlot({
    d <- sim_data()
    nadir_t <- d$time[which.min(d$ANC)]
    nadir_v <- min(d$ANC, na.rm = TRUE)

    sev_grade <- dplyr::case_when(
      nadir_v < 500   ~ "Grade 4 (Severe)",
      nadir_v < 1000  ~ "Grade 3",
      nadir_v < 1500  ~ "Grade 2",
      nadir_v < 2000  ~ "Grade 1",
      TRUE            ~ "No toxicity"
    )
    grade_col <- dplyr::case_when(
      nadir_v < 500  ~ "#dc2626",
      nadir_v < 1000 ~ "#f59e0b",
      nadir_v < 1500 ~ "#facc15",
      TRUE           ~ "#10b981"
    )

    ggplot(d, aes(x = time, y = ANC)) +
      annotate("rect", xmin = -Inf, xmax = Inf, ymin = 0,    ymax = 500,
               fill = "#fee2e2", alpha = 0.5) +
      annotate("rect", xmin = -Inf, xmax = Inf, ymin = 500,  ymax = 1000,
               fill = "#fef3c7", alpha = 0.5) +
      annotate("rect", xmin = -Inf, xmax = Inf, ymin = 1000, ymax = 1500,
               fill = "#fefce8", alpha = 0.4) +
      geom_hline(yintercept = c(500, 1000, 1500),
                 linetype = "dotted", color = "grey50", alpha = 0.7) +
      annotate("text", x = max(d$time) * 0.98, y = 250,  label = "Grade 4 ANC <500",  hjust = 1, size = 3, color = "#dc2626") +
      annotate("text", x = max(d$time) * 0.98, y = 750,  label = "Grade 3 ANC <1000", hjust = 1, size = 3, color = "#d97706") +
      annotate("text", x = max(d$time) * 0.98, y = 1250, label = "Grade 2 ANC <1500", hjust = 1, size = 3, color = "#ca8a04") +
      geom_line(color = "#3b82f6", linewidth = 0.8) +
      geom_point(aes(x = nadir_t, y = nadir_v), color = grade_col, size = 4) +
      annotate("label", x = nadir_t, y = nadir_v + max(d$ANC) * 0.08,
               label = paste0("Nadir: ", round(nadir_v), " cells/\u03bcL\n", sev_grade),
               size = 3.3, fill = "white", color = grade_col) +
      labs(x = "Time (h)",
           y = "ANC (cells/\u03bcL)",
           title = "Friberg Semi-Physiological Myelosuppression Model \u2014 Neutrophil Dynamics") +
      theme_minimal(base_size = 14) +
      theme(plot.title = element_text(size = 12))
  })
}

shinyApp(ui = ui, server = server)