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

# ── PBPK-PD Model: T-DXd Multi-Tissue + Tumor PD ─────────────────────────
# Based on: Wang Z et al. (2025) Eur J Pharm Sci 213:107234
# Full tissue PBPK: lung, tumor, liver, kidney, rest-of-body
# ILD risk endpoint: DXd lung Cmax vs IC thresholds
# Efficacy endpoint: Tumor growth inhibition by DXd
# ──────────────────────────────────────────────────────────────────────────
model_code <- '
$PARAM @annotated
// ── T-DXd ADC pharmacokinetics (Yin et al., 2021 PopPK) ──
CL_ADC   : 0.0175  : T-DXd clearance (L/h)
V1_ADC   : 2.68    : Central volume (L)
V2_ADC   : 2.06    : Peripheral volume (L)
Q_ADC    : 0.00958 : Intercompartmental CL (L/h)

// ── Lung endosomal processing (Wang et al. 2025 calibrated) ──
KIN_L    : 0.0058  : ADC uptake into lung endosomes (L/h)
KPROC    : 0.0085  : Endosomal processing rate (1/h)
FGEN     : 4.10    : DXd yield per ADC degraded (ug/mg)
KGEN_SYS : 0.0025  : Systemic DXd generation (ug/h per mg ADC)

// ── Tumor endosomal processing (HER2-mediated) ──
KIN_T    : 0.025   : ADC tumor uptake rate (L/h) [overridden by HER2 level]
KPROC_T  : 0.025   : Tumor endosomal processing (1/h)

// ── DXd elimination ──
CL_DXD   : 3.0     : DXd systemic clearance (L/h)

// ── Physiological parameters (65y female, Wang 2025 Table S6) ──
QC       : 289     : Cardiac output (L/h)
FQL      : 0.255   : Liver blood flow fraction
FQK      : 0.190   : Kidney blood flow fraction
FQT      : 0.005   : Tumor blood flow fraction

// ── Tissue volumes (L) ──
VL_T     : 0.95    : Lung tissue volume (L)
VLV      : 1.45    : Liver volume (L)
VK       : 0.28    : Kidney volume (L)
VR       : 38.0    : Rest-of-body volume (L)
VP_DXD   : 2.04    : DXd plasma volume (L)
VT       : 0.05    : Tumor volume (L) [overridden by user input]

// ── DXd tissue partition coefficients (Kp, tissue:plasma) ──
KP_L     : 2.50    : Lung Kp (calibrated; lipophilic payload accumulates in lung)
KP_LV    : 3.5     : Liver Kp
KP_K     : 3.2     : Kidney Kp
KP_R     : 0.8     : Rest-of-body Kp
KP_T     : 8.0     : Tumor Kp [overridden by HER2 level]

// ── Tumor PD (growth inhibition) ──
KG       : 0.005   : Tumor growth rate (1/h)
KMAX     : 0.08    : Max DXd kill rate (1/h)
IC50_T   : 15.0    : DXd tumor IC50 (ng/mL)
HILL_T   : 1.5     : Hill coefficient

// ── Patient factor ──
AGE_GRP  : 2       : Age group (1=50-60y 2=60-70y 3=70-80y)

$CMT @annotated
ADC_C  : T-DXd central compartment (mg)
ADC_P  : T-DXd peripheral compartment (mg)
ENDO_L : ADC in lung endosomes (mg)
ENDO_T : ADC in tumor endosomes (mg)
DXD_L  : DXd in lung tissue (ug)
DXD_T  : DXd in tumor tissue (ug)
DXD_LV : DXd in liver (ug)
DXD_K  : DXd in kidney (ug)
DXD_R  : DXd in rest-of-body (ug)
DXD_P  : DXd in systemic plasma (ug)
TUMOR  : Tumor size (normalized)

$MAIN
// Age-dependent lung physiological adjustments (Wang 2025 Table S6)
double gen_f    = 1.0;
double kproc_i  = KPROC;

if (AGE_GRP == 1) {
  // 50-60 years: higher cardiac output, larger lungs
  gen_f   = 0.871;
} else if (AGE_GRP == 2) {
  // 60-70 years: reference
  gen_f   = 1.0;
} else if (AGE_GRP == 3) {
  // 70-80 years: reduced CO, slower metabolism
  gen_f   = 0.953;
  kproc_i = KPROC * 0.93;
}

// Tumor initial condition
TUMOR_0 = 1.0;

$ODE
// ── T-DXd ADC: 2-compartment IV ──
dxdt_ADC_C = -(CL_ADC/V1_ADC + Q_ADC/V1_ADC)*ADC_C + (Q_ADC/V2_ADC)*ADC_P;
dxdt_ADC_P =  (Q_ADC/V1_ADC)*ADC_C - (Q_ADC/V2_ADC)*ADC_P;

double Cp_ADC = ADC_C / V1_ADC;          // mg/L (ug/mL)
double Cp_DXd = DXD_P / VP_DXD;          // ug/L (= ng/mL)

// ── Lung endosomal processing ──
dxdt_ENDO_L = KIN_L * gen_f * Cp_ADC - kproc_i * ENDO_L;

// ── Tumor endosomal processing (HER2-mediated) ──
dxdt_ENDO_T = KIN_T * Cp_ADC - KPROC_T * ENDO_T;

// ── Blood flows ──
// Systemic flows (sum = QC = one cardiac output)
double QLV  = FQL * QC;
double QK   = FQK * QC;
double QT_f = FQT * QC;
double QR   = QC - QLV - QK - QT_f;
// Lung receives ALL venous return (= QC) — same circuit, not additive

// ── DXd generation ──
double gen_lung  = FGEN * kproc_i * ENDO_L;   // ug/h
double gen_tumor = FGEN * KPROC_T * ENDO_T;   // ug/h

// ── Correct venous-arterial PBPK structure ──
// Systemic tissues receive ARTERIAL blood (Cp_DXd), drain to VENOUS
// Lung receives mixed VENOUS, outputs to ARTERIAL (= plasma pool)
// This correctly conserves mass: d(total DXd)/dt = generation - clearance

// Mixed venous concentration from systemic tissues
double C_ven = (QLV * DXD_LV/(VLV*KP_LV)
              + QK  * DXD_K /(VK *KP_K )
              + QT_f* DXD_T /(VT *KP_T )
              + QR  * DXD_R /(VR *KP_R )) / QC;  // ug/L = ng/mL

// Lung: perfused by venous blood, generates DXd from endosomal processing
dxdt_DXD_L  = QC  * (C_ven  - DXD_L /(VL_T * KP_L)) + gen_lung;

// Systemic tissues: receive arterial (Cp_DXd), drain to venous
dxdt_DXD_T  = QT_f* (Cp_DXd - DXD_T /(VT   * KP_T )) + gen_tumor;
dxdt_DXD_LV = QLV * (Cp_DXd - DXD_LV/(VLV  * KP_LV));
dxdt_DXD_K  = QK  * (Cp_DXd - DXD_K /(VK   * KP_K ));
dxdt_DXD_R  = QR  * (Cp_DXd - DXD_R /(VR   * KP_R ));

// Arterial plasma: receives lung output, distributes to systemic, eliminated
dxdt_DXD_P  = QC  * DXD_L/(VL_T*KP_L)   // lung → arterial
             - QC  * Cp_DXd              // arterial → systemic
             - CL_DXD * Cp_DXd           // elimination
             + KGEN_SYS * ADC_C;         // systemic ADC deconjugation (ug/h)

// ── Tumor growth inhibition ──
double CT_ngml = DXD_T / VT;   // ug/L = ng/mL
double kill_r  = KMAX * pow(CT_ngml, HILL_T) /
                (pow(IC50_T, HILL_T) + pow(CT_ngml, HILL_T));
dxdt_TUMOR = KG * TUMOR - kill_r * TUMOR;

$TABLE
double CP_TDXd   = ADC_C / V1_ADC;       // T-DXd plasma (ug/mL)
double CP_DXd_pl = DXD_P / VP_DXD;       // DXd plasma (ng/mL)
double CT_lung   = DXD_L / VL_T;         // DXd lung (ng/mL)
double CT_tumor  = DXD_T / VT;           // DXd tumor (ng/mL)
double CT_liver  = DXD_LV / VLV;         // DXd liver (ng/mL)
double CT_kidney = DXD_K / VK;           // DXd kidney (ng/mL)
double TUMOR_SZ  = TUMOR;                // normalized (1 = baseline)

$CAPTURE @annotated
CP_TDXd   : T-DXd plasma (ug/mL)
CP_DXd_pl : DXd plasma (ng/mL)
CT_lung   : DXd lung tissue (ng/mL)
CT_tumor  : DXd tumor (ng/mL)
CT_liver  : DXd liver (ng/mL)
CT_kidney : DXd kidney (ng/mL)
TUMOR_SZ  : Tumor size (normalized to baseline)
'

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

# ── ILD Risk Thresholds (IC values from Wang 2025 Table S3) ──────────────
IC5  <- 8.05
IC10 <- 17.00
IC20 <- 38.25
IC50 <- 153.0

get_ild_grade <- function(cmax_ngml) {
  if (is.na(cmax_ngml) || cmax_ngml < IC10) {
    list(grade = "Grade 1", label = "Low Risk", color = "#10b981",
         action = "Below IC\u2081\u2080 threshold. Monitor; temporary interruption if symptomatic.")
  } else if (cmax_ngml < IC20) {
    list(grade = "Grade 2", label = "Moderate Risk", color = "#f59e0b",
         action = "IC\u2081\u2080\u2013IC\u2082\u2080 range. Consider dose reduction; discontinue if Grade \u22652 symptoms.")
  } else if (cmax_ngml < IC50) {
    list(grade = "Grade 3", label = "High Risk", color = "#ef4444",
         action = "Above IC\u2082\u2080. Permanent discontinuation recommended; immediate intervention.")
  } else {
    list(grade = "Grade 4", label = "Severe Risk", color = "#dc2626",
         action = "Above IC\u2085\u2080. Permanent discontinuation; emergency intensive care.")
  }
}

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

# ── UI ───────────────────────────────────────────────────────────────────
ui <- page_sidebar(
  title = "T-DXd PBPK-PD: Multi-Tissue + ILD + Tumor Simulator",
  theme = app_theme,

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

    accordion(
      open = c("Dosing", "Patient"),

      accordion_panel(
        "Dosing",
        icon = bsicons::bs_icon("capsule"),
        selectInput("dose_mgkg", "T-DXd Dose (mg/kg)",
                    choices = c("3.2" = 3.2, "4.4" = 4.4,
                                "5.4 (FDA approved)" = 5.4, "6.4" = 6.4),
                    selected = 5.4),
        numericInput("n_cycles", "Number of Cycles (Q3W)", value = 4, min = 1, max = 8),
        helpText("Q3W = every 21 days (504 hours)")
      ),

      accordion_panel(
        "Patient & Tumor",
        icon = bsicons::bs_icon("person"),
        selectInput("age_grp", "Age Group",
                    choices = c("50-60 years" = 1, "60-70 years" = 2, "70-80 years" = 3),
                    selected = 2),
        sliderInput("wt", "Body Weight (kg)", min = 40, max = 120, value = 65, step = 5),
        selectInput("her2_lvl", "HER2 Expression",
                    choices = c("Low (IHC 1+)" = 1, "Medium (IHC 2+)" = 2, "High (IHC 3+)" = 3),
                    selected = 3),
        numericInput("tumor_vol_ml", "Tumor Volume (mL)", value = 50, min = 1, max = 500),
        helpText("Tumor volume affects DXd distribution and PD response")
      ),

      accordion_panel(
        "Display",
        icon = bsicons::bs_icon("gear"),
        checkboxInput("log_scale", "Log Scale (Y-axis)", value = FALSE),
        checkboxInput("show_ic", "Show IC Thresholds", value = TRUE)
      )
    )
  ),

  # ── Metric Cards ──
  layout_column_wrap(
    width = 1/5, fill = FALSE,

    div(class = "metric-card metric-primary",
        div(class = "metric-value", textOutput("cmax_adc")),
        div(class = "metric-label", "T-DXd Plasma Cmax"),
        div(class = "metric-unit", "\u00B5g/mL")),

    div(class = "metric-card metric-warning",
        div(class = "metric-value", textOutput("cmax_lung")),
        div(class = "metric-label", "DXd Lung Cmax"),
        div(class = "metric-unit", "ng/mL")),

    div(class = "metric-card metric-danger",
        div(class = "metric-value", textOutput("cmax_tumor")),
        div(class = "metric-label", "DXd Tumor Cmax"),
        div(class = "metric-unit", "ng/mL")),

    div(id = "ild_card", class = "metric-card",
        div(class = "metric-value", textOutput("ild_grade")),
        div(class = "metric-label", "ILD Risk Grade"),
        div(class = "metric-unit", htmlOutput("ild_label"))),

    div(id = "tumor_card", class = "metric-card",
        div(class = "metric-value", textOutput("tumor_reduction")),
        div(class = "metric-label", "Tumor Reduction"),
        div(class = "metric-unit", "% vs baseline at end"))
  ),

  # ── Main Tabs ──
  navset_card_underline(
    title = "T-DXd PBPK-PD Analysis",

    # 1. PK Profile
    nav_panel(
      "PK Profile",
      icon = bsicons::bs_icon("graph-up"),
      layout_columns(
        col_widths = c(6, 6),
        card(
          card_header("T-DXd Plasma Concentration"),
          full_screen = TRUE,
          plotOutput("pk_adc_plot", height = "380px")
        ),
        card(
          card_header("DXd Plasma Concentration"),
          full_screen = TRUE,
          plotOutput("pk_dxd_plot", height = "380px")
        )
      )
    ),

    # 2. Tissue Distribution
    nav_panel(
      "Tissue Distribution",
      icon = bsicons::bs_icon("diagram-3"),
      layout_columns(
        col_widths = c(8, 4),
        card(
          card_header("DXd Concentration by Tissue"),
          full_screen = TRUE,
          plotOutput("tissue_plot", height = "420px")
        ),
        card(
          card_header("Tissue Cmax Summary"),
          tableOutput("tissue_table"),
          hr(),
          markdown("
**Tissue ranking** reflects both perfusion and partitioning (Kp).
High Kp tissues (liver, tumor) accumulate more DXd relative to plasma.
Lung Cmax drives ILD risk; tumor Cmax drives anti-tumor efficacy.
          ")
        )
      )
    ),

    # 3. ILD Risk Assessment
    nav_panel(
      "ILD Risk Assessment",
      icon = bsicons::bs_icon("exclamation-triangle"),
      layout_columns(
        col_widths = c(7, 5),
        card(
          card_header("DXd Lung Cmax vs IC Thresholds"),
          full_screen = TRUE,
          plotOutput("ild_plot", height = "420px")
        ),
        card(
          card_header("Risk Stratification"),
          htmlOutput("risk_detail"),
          hr(),
          markdown("
**ILD Grading** (Wang et al., 2025):

| Grade | DXd Lung Conc. | Severity |
|-------|---------------|----------|
| 1 | < IC10 (17.0 ng/mL) | Low |
| 2 | IC10\u2013IC20 (17\u201338 ng/mL) | Moderate |
| 3 | IC20\u2013IC50 (38\u2013153 ng/mL) | High |
| 4 | > IC50 (>153 ng/mL) | Severe |
          ")
        )
      )
    ),

    # 4. Tumor PD
    nav_panel(
      "Tumor PD",
      icon = bsicons::bs_icon("graph-down-arrow"),
      layout_columns(
        col_widths = c(6, 6),
        card(
          card_header("Tumor Size Over Time"),
          full_screen = TRUE,
          plotOutput("tumor_pd_plot", height = "400px")
        ),
        card(
          card_header("DXd Tumor Concentration"),
          full_screen = TRUE,
          plotOutput("tumor_dxd_plot", height = "400px"),
          htmlOutput("tumor_pd_summary")
        )
      )
    ),

    # 5. Model Parameters
    nav_panel(
      "Model Parameters",
      icon = bsicons::bs_icon("table"),
      div(style = "padding: 16px;",
        markdown("
## PBPK-PD Model Structure

**Model Type:** Whole-body PBPK-PD (multi-tissue)
**Tissues:** Plasma, Lung, Tumor, Liver, Kidney, Rest-of-body
**ADC:** Trastuzumab deruxtecan (T-DXd, DS-8201a) | **Payload:** DXd | **DAR:** 8

### T-DXd ADC PK (Yin et al. 2021)

| Parameter | Value | Units |
|-----------|-------|-------|
| CL | 0.42 | L/day |
| V1 (central) | 2.68 | L |
| V2 (peripheral) | 2.06 | L |
| Q | 0.23 | L/day |

### Physiological Parameters (65y female, Wang 2025 Table S6)

| Tissue | Volume (L) | Blood Flow (% CO) | Kp |
|--------|-----------|-------------------|-----|
| Plasma | 2.04 | — | — |
| Lung | 0.95 | 100% (all CO) | 0.12 |
| Liver | 1.45 | 25.5% | 3.5 |
| Kidney | 0.28 | 19.0% | 3.2 |
| Tumor | user-defined | 0.5% | HER2-dep. |
| Rest | 38.0 | remainder | 0.8 |

### HER2-Dependent Tumor Parameters

| HER2 Level | KIN_T (L/h) | KP_T |
|-----------|------------|------|
| Low (IHC 1+) | 0.010 | 3.0 |
| Medium (IHC 2+) | 0.025 | 5.0 |
| High (IHC 3+) | 0.050 | 8.0 |

### Tumor PD Model

| Parameter | Value | Units |
|-----------|-------|-------|
| KG (growth rate) | 0.005 | 1/h |
| KMAX (max kill) | 0.080 | 1/h |
| IC50_T | 15.0 | ng/mL |
| Hill | 1.5 | — |
        ")
      )
    ),

    # 6. References
    nav_panel(
      "References",
      icon = bsicons::bs_icon("journal-text"),
      div(class = "ref-box", style = "margin: 16px;",
        tags$h5("\U0001F4DA Primary Reference"),
        tags$ol(
          tags$li(tags$strong("Wang Z et al."),
                  " (2025) PBPK-PD model for predicting pharmacokinetics, tumor growth inhibition, and toxicity risks of topoisomerase inhibitor ADCs. ",
                  tags$em("Eur J Pharm Sci"), " 213:107234. ",
                  tags$a(href = "https://doi.org/10.1016/j.ejps.2025.107234",
                         target = "_blank", "DOI"))
        ),
        tags$h5("\U0001F4DA Supporting References"),
        tags$ol(start = 2,
          tags$li("Yin O et al. (2021) Population PK of T-DXd. ", tags$em("Clin Pharmacol Ther"), " 109:1314."),
          tags$li("Doi T et al. (2017) Safety, PK, antitumour activity of T-DXd. ", tags$em("Lancet Oncol"), " 18:1512."),
          tags$li("Ogitani Y et al. (2016) DS-8201a: novel HER2-targeting ADC. ", tags$em("Clin Cancer Res"), " 22:5097."),
          tags$li("Scheuher B et al. (2023) Platform QSP model for ADCs. ", tags$em("J Pharmacokinet Pharmacodyn.")),
          tags$li("Powell CA et al. (2022) Drug-related ILD in T-DXd studies. ", tags$em("ESMO Open"), " 7:100554."),
          tags$li("Modi S et al. (2022) T-DXd in HER2-low breast cancer. ", tags$em("N Engl J Med"), " 387:9.")
        )
      )
    )
  ),

  # ── Footer ──
  div(style = "text-align: center; padding: 20px; margin-top: 20px;
               border-top: 1px solid #e9ecef; color: #6c757d; font-size: 12px;",
    "Powered by ",
    tags$a(href = "https://www.pkpdbuilder.com", target = "_blank",
           style = "color: #8b5cf6; font-weight: 500;", "PKPDBuilder.com"),
    " \u2022 Built by Sunny \u2600\uFE0F (Husain Attarwala\u2019s AI Assistant)",
    tags$br(),
    tags$span(style = "font-size: 10px;",
      "For research and educational purposes only. Not for clinical decision-making.")
  )
)

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

  # HER2-level lookup
  her2_params <- list(
    "1" = list(kin = 0.010, kp = 3.0),
    "2" = list(kin = 0.025, kp = 5.0),
    "3" = list(kin = 0.050, kp = 8.0)
  )

  # ── Reactive simulation ──
  sim_data <- reactive({
    shiny::req(input$dose_mgkg, input$n_cycles, input$age_grp,
               input$wt, input$her2_lvl, input$tumor_vol_ml)

    dose_mg  <- as.numeric(input$dose_mgkg) * input$wt
    ii_h     <- 504  # Q3W
    her2_key <- as.character(input$her2_lvl)
    kin_t    <- her2_params[[her2_key]]$kin
    kp_t     <- her2_params[[her2_key]]$kp
    vt       <- max(input$tumor_vol_ml, 1) / 1000  # mL → L

    ev1 <- ev(amt = dose_mg, cmt = 1, ii = ii_h,
              addl = as.integer(input$n_cycles) - 1L)
    end_time <- as.integer(input$n_cycles) * 504 + 200

    mod %>%
      param(AGE_GRP = as.numeric(input$age_grp),
            KIN_T   = kin_t,
            KP_T    = kp_t,
            VT      = vt) %>%
      init(TUMOR = 1.0) %>%
      ev(ev1) %>%
      mrgsim(end = end_time, delta = 1) %>%
      as.data.frame()
  })

  # ── Derived metrics ──
  metrics <- reactive({
    d <- sim_data()
    list(
      cmax_adc    = max(d$CP_TDXd,   na.rm = TRUE),
      cmax_lung   = max(d$CT_lung,   na.rm = TRUE),
      cmax_tumor  = max(d$CT_tumor,  na.rm = TRUE),
      cmax_liver  = max(d$CT_liver,  na.rm = TRUE),
      cmax_kidney = max(d$CT_kidney, na.rm = TRUE),
      cmax_dxd_pl = max(d$CP_DXd_pl, na.rm = TRUE),
      tmax_lung   = d$time[which.max(d$CT_lung)],
      tumor_end   = tail(d$TUMOR_SZ, 1)
    )
  })

  # ── Metric card outputs ──
  output$cmax_adc <- renderText(sprintf("%.1f", metrics()$cmax_adc))
  output$cmax_lung <- renderText(sprintf("%.1f", metrics()$cmax_lung))
  output$cmax_tumor <- renderText(sprintf("%.1f", metrics()$cmax_tumor))

  output$ild_grade <- renderText({
    get_ild_grade(metrics()$cmax_lung)$grade
  })
  output$ild_label <- renderUI({
    ild <- get_ild_grade(metrics()$cmax_lung)
    tags$span(style = paste0("color:", ild$color, "; font-weight:bold;"), ild$label)
  })

  output$tumor_reduction <- renderText({
    pct <- (1 - metrics()$tumor_end) * 100
    sprintf("%.0f%%", pct)
  })

  # ── T-DXd PK Plot ──
  output$pk_adc_plot <- renderPlot({
    d <- sim_data()
    d_plot <- if (input$log_scale) filter(d, CP_TDXd > 0.01) else d
    ggplot(d_plot, aes(x = time/24, y = CP_TDXd)) +
      geom_line(color = "#8b5cf6", linewidth = 0.9) +
      labs(x = "Time (days)", y = "T-DXd Plasma (\u00B5g/mL)",
           title = paste0("T-DXd ", input$dose_mgkg, " mg/kg Q3W"),
           subtitle = paste0("Cmax = ", sprintf("%.1f", metrics()$cmax_adc), " \u00B5g/mL")) +
      theme_minimal(base_size = 13) +
      theme(plot.title = element_text(face = "bold"),
            plot.subtitle = element_text(color = "#6c757d", size = 10),
            panel.grid.minor = element_blank()) +
      { if (input$log_scale) scale_y_log10(labels = scales::label_comma()) }
  })

  # ── DXd Plasma Plot ──
  output$pk_dxd_plot <- renderPlot({
    d <- sim_data()
    d_plot <- if (input$log_scale) filter(d, CP_DXd_pl > 0.001) else d
    ggplot(d_plot, aes(x = time/24, y = CP_DXd_pl)) +
      geom_line(color = "#0ea5e9", linewidth = 0.9) +
      labs(x = "Time (days)", y = "DXd Plasma (ng/mL)",
           title = "Free DXd Plasma Concentration",
           subtitle = paste0("Cmax = ", sprintf("%.1f", metrics()$cmax_dxd_pl), " ng/mL")) +
      theme_minimal(base_size = 13) +
      theme(plot.title = element_text(face = "bold"),
            plot.subtitle = element_text(color = "#6c757d", size = 10),
            panel.grid.minor = element_blank()) +
      { if (input$log_scale) scale_y_log10(labels = scales::label_comma()) }
  })

  # ── Tissue Distribution Plot ──
  output$tissue_plot <- renderPlot({
    d <- sim_data()

    tissue_colors <- c(
      "Lung"   = "#ef4444",
      "Tumor"  = "#f97316",
      "Liver"  = "#3b82f6",
      "Kidney" = "#10b981"
    )

    d_long <- d %>%
      select(time, CT_lung, CT_tumor, CT_liver, CT_kidney) %>%
      pivot_longer(-time, names_to = "tissue", values_to = "conc") %>%
      mutate(tissue = recode(tissue,
        CT_lung   = "Lung",
        CT_tumor  = "Tumor",
        CT_liver  = "Liver",
        CT_kidney = "Kidney"
      ))

    if (input$log_scale) d_long <- filter(d_long, conc > 0.001)

    p <- ggplot(d_long, aes(x = time/24, y = conc, color = tissue)) +
      geom_line(linewidth = 0.9) +
      scale_color_manual(values = tissue_colors) +
      labs(x = "Time (days)", y = "DXd Concentration (ng/mL)",
           title = "DXd Multi-Tissue Distribution",
           color = "Tissue") +
      theme_minimal(base_size = 13) +
      theme(plot.title = element_text(face = "bold"),
            panel.grid.minor = element_blank(),
            legend.position = "bottom")

    if (input$show_ic) {
      p <- p +
        geom_hline(yintercept = IC10, linetype = "dashed", color = "#f59e0b", alpha = 0.6) +
        geom_hline(yintercept = IC20, linetype = "dashed", color = "#ef4444", alpha = 0.6) +
        annotate("text", x = 0.3, y = IC10 + 1, label = "Lung IC10 (ILD threshold)",
                 hjust = 0, size = 3, color = "#f59e0b") +
        annotate("text", x = 0.3, y = IC20 + 1, label = "Lung IC20",
                 hjust = 0, size = 3, color = "#ef4444")
    }

    if (input$log_scale) p <- p + scale_y_log10(labels = scales::label_comma())
    p
  })

  # ── Tissue Cmax Table ──
  output$tissue_table <- renderTable({
    m <- metrics()
    her2_key <- as.character(input$her2_lvl)
    her2_label <- c("1" = "Low (IHC 1+)", "2" = "Medium (IHC 2+)", "3" = "High (IHC 3+)")[her2_key]
    data.frame(
      Tissue = c("Plasma", "Lung", "Tumor", "Liver", "Kidney"),
      `Cmax (ng/mL)` = sprintf("%.1f", c(
        m$cmax_dxd_pl, m$cmax_lung, m$cmax_tumor, m$cmax_liver, m$cmax_kidney
      )),
      Note = c("Systemic free DXd", "ILD endpoint", paste0("Efficacy (", her2_label, ")"),
               "High Kp tissue", "High Kp tissue"),
      check.names = FALSE
    )
  }, striped = TRUE, hover = TRUE, bordered = TRUE, spacing = "s")

  # ── ILD Risk Plot ──
  output$ild_plot <- renderPlot({
    cmax   <- metrics()$cmax_lung
    ild    <- get_ild_grade(cmax)
    ic_df  <- data.frame(
      threshold = c("IC5", "IC10", "IC20", "IC50"),
      value     = c(IC5, IC10, IC20, IC50),
      fill_color = c("#10b981", "#f59e0b", "#f97316", "#ef4444")
    )
    y_max <- max(IC50, cmax) * 1.15

    ggplot() +
      geom_col(data = ic_df, aes(x = threshold, y = value, fill = threshold),
               width = 0.5, alpha = 0.3) +
      scale_fill_manual(values = c("IC5" = "#10b981", "IC10" = "#f59e0b",
                                    "IC20" = "#f97316", "IC50" = "#ef4444")) +
      geom_hline(yintercept = cmax, color = ild$color, linewidth = 2) +
      annotate("label", x = 2.5, y = cmax,
               label = paste0("Lung DXd Cmax\n", sprintf("%.1f ng/mL", cmax)),
               fill = ild$color, color = "white", fontface = "bold", size = 4,
               label.padding = unit(0.5, "lines")) +
      geom_text(data = ic_df, aes(x = threshold, y = value + y_max*0.03,
                label = paste0(sprintf("%.1f", value), " ng/mL")), size = 3.5, color = "#495057") +
      scale_x_discrete(limits = c("IC5", "IC10", "IC20", "IC50")) +
      scale_y_continuous(limits = c(0, y_max), expand = c(0, 0)) +
      labs(title = paste0("ILD Risk: ", input$dose_mgkg, " mg/kg Q3W, ", input$wt, " kg"),
           subtitle = paste0("Predicted: ", ild$grade, " (", ild$label, ")"),
           x = "DXd Inhibitory Concentration Threshold",
           y = "DXd Lung Concentration (ng/mL)") +
      theme_minimal(base_size = 13) +
      theme(plot.title = element_text(face = "bold"),
            plot.subtitle = element_text(color = ild$color, face = "bold", size = 12),
            legend.position = "none",
            panel.grid.major.x = element_blank(),
            panel.grid.minor = element_blank())
  })

  # ── ILD Risk Detail ──
  output$risk_detail <- renderUI({
    cmax <- metrics()$cmax_lung
    ild  <- get_ild_grade(cmax)
    pct_ic20 <- cmax / IC20 * 100
    tagList(
      div(class = "risk-box", style = paste0("background:", ild$color, ";"),
          tags$h3(style = "margin:0;", ild$grade),
          tags$p(style = "margin:5px 0 0; font-size:14px;", ild$label)
      ),
      div(style = "padding: 10px;",
        tags$p(tags$strong("DXd Lung Cmax: "), sprintf("%.1f ng/mL", cmax)),
        tags$p(tags$strong("% of IC20: "),
               tags$span(style = paste0("color:", if (pct_ic20 > 80) "#ef4444" else "#495057",
                                        "; font-weight:bold;"),
                         sprintf("%.0f%%", pct_ic20))),
        tags$p(tags$strong("Action: "), ild$action),
        tags$p(style = "font-size: 11px; color: #6c757d; margin-top:8px;",
          tags$em("Powell et al. (2022): 15.4% ILD incidence with T-DXd, 2.2% fatal."))
      )
    )
  })

  # ── Tumor PD Plot ──
  output$tumor_pd_plot <- renderPlot({
    d <- sim_data()
    tumor_end_pct <- round((1 - tail(d$TUMOR_SZ, 1)) * 100, 0)
    her2_key  <- as.character(input$her2_lvl)
    her2_label <- c("1" = "HER2-Low", "2" = "HER2-Med", "3" = "HER2-High")[her2_key]

    ggplot(d, aes(x = time/24, y = TUMOR_SZ * 100)) +
      geom_line(color = "#8b5cf6", linewidth = 1.0) +
      geom_hline(yintercept = 50, linetype = "dashed", color = "#6c757d", alpha = 0.5) +
      annotate("text", x = 0.5, y = 52, label = "50% size", size = 3, color = "#6c757d") +
      labs(x = "Time (days)", y = "Tumor Size (% of baseline)",
           title = paste0("Tumor Growth Inhibition — ", her2_label),
           subtitle = paste0("Size at end: ", 100 - tumor_end_pct, "% reduction from baseline")) +
      theme_minimal(base_size = 13) +
      theme(plot.title = element_text(face = "bold"),
            plot.subtitle = element_text(color = "#8b5cf6", size = 10),
            panel.grid.minor = element_blank()) +
      scale_y_continuous(limits = c(0, 105), labels = function(x) paste0(x, "%"))
  })

  # ── Tumor DXd Concentration Plot ──
  output$tumor_dxd_plot <- renderPlot({
    d <- sim_data()
    d_plot <- if (input$log_scale) filter(d, CT_tumor > 0.001) else d

    p <- ggplot(d_plot, aes(x = time/24, y = CT_tumor)) +
      geom_line(color = "#f97316", linewidth = 0.9)

    if (input$show_ic) {
      p <- p +
        geom_hline(yintercept = 15.0, linetype = "dashed", color = "#8b5cf6", alpha = 0.8) +
        annotate("text", x = 0.5, y = 16.5, label = "IC50 tumor (15 ng/mL)",
                 hjust = 0, size = 3.2, color = "#8b5cf6")
    }

    p <- p +
      labs(x = "Time (days)", y = "DXd Tumor (ng/mL)",
           title = "DXd Tumor Tissue Concentration",
           subtitle = paste0("Cmax = ", sprintf("%.1f", metrics()$cmax_tumor), " ng/mL")) +
      theme_minimal(base_size = 13) +
      theme(plot.title = element_text(face = "bold"),
            plot.subtitle = element_text(color = "#f97316", size = 10),
            panel.grid.minor = element_blank())

    if (input$log_scale) p <- p + scale_y_log10(labels = scales::label_comma())
    p
  })

  # ── Tumor PD Summary ──
  output$tumor_pd_summary <- renderUI({
    m <- metrics()
    pct_red <- (1 - m$tumor_end) * 100
    cmax_t  <- m$cmax_tumor
    ic50_t  <- 15.0

    response_label <- if (pct_red >= 30) {
      tags$span(style = "color:#10b981; font-weight:bold;",
                if (pct_red >= 80) "Deep Response (CR/PR)" else "Partial Response")
    } else if (pct_red > 0) {
      tags$span(style = "color:#f59e0b; font-weight:bold;", "Stable Disease")
    } else {
      tags$span(style = "color:#ef4444; font-weight:bold;", "Progressive Disease")
    }

    div(style = "padding: 10px; font-size: 13px;",
      tags$p(tags$strong("Tumor Cmax: "), sprintf("%.1f ng/mL", cmax_t),
             " vs IC50 = ", sprintf("%.1f ng/mL", ic50_t)),
      tags$p(tags$strong("Cmax/IC50: "),
             tags$span(style = "font-weight:bold;", sprintf("%.1f\u00D7", cmax_t/ic50_t))),
      tags$p(tags$strong("Predicted response: "), response_label),
      tags$p(tags$strong("Tumor reduction: "), sprintf("%.0f%%", pct_red)),
      tags$p(style = "font-size: 11px; color: #6c757d; margin-top: 8px;",
        tags$em("Tumor PD based on Simeoni growth inhibition model. ",
                "Response correlates with HER2 expression and tumor volume."))
    )
  })
}

shinyApp(ui = ui, server = server)