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

# ── Rucaparib Clinical PK Model (2-compartment oral) ───────────────────────
pk_code <- '
$PARAM @annotated
CL   : 37.6  : Apparent oral clearance CL/F (L/h)
V1   : 310   : Central volume V1/F (L)
V2   : 230   : Peripheral volume V2/F (L)
Q    : 16    : Inter-compartmental clearance Q/F (L/h)
KA   : 1.15  : Absorption rate constant (1/h)
WT   : 70    : Body weight (kg)
HEPF : 1.0   : Hepatic function (fraction of normal)

$CMT @annotated
DEPOT  : Oral depot (mg)
CENT   : Central compartment (mg)
PERI   : Peripheral compartment (mg)

$MAIN
double CLi  = CL * pow(WT/70.0, 0.75) * HEPF;
double V1i  = V1 * (WT/70.0);
double V2i  = V2 * (WT/70.0);
double Qi   = Q  * pow(WT/70.0, 0.75);

$ODE
dxdt_DEPOT = -KA * DEPOT;
dxdt_CENT  = KA * DEPOT - (CLi/V1i)*CENT - (Qi/V1i)*CENT + (Qi/V2i)*PERI;
dxdt_PERI  = (Qi/V1i)*CENT - (Qi/V2i)*PERI;

$TABLE
double CP = CENT / V1i;     // mg/L
double CPng = CP * 1000.0;  // ng/mL

$CAPTURE CP CPng
'

# ── DDR-based Tumor Growth Inhibition Model (simplified from paper) ────────
# Uses Gompertz growth + PK-driven kill with HR deficiency modulating sensitivity
# Based on the DDR framework from Villette et al. (2025) Br J Cancer
tgi_code <- '
$PARAM @annotated
// Mouse PK (1-comp oral, rucaparib in mouse)
CLm   : 3.0   : Mouse apparent clearance (L/h/kg)
Vm    : 2.5   : Mouse apparent volume (L/kg)
KAm   : 2.0   : Mouse absorption rate (1/h)

// Drug effect params (calibrated to paper behavior)
kmax    : 0.015  : Maximum kill rate constant (1/h)
KC50    : 5.0    : Concentration for half-max kill (mg/L)
sens    : 4.0    : HR deficiency sensitivity amplifier (fold)
ip_fork : 0.5    : PARPi fork trapping potency (/(mg/L)/h)
clpf    : 0.45   : Clearance of fork trapping effect (1/h)

// Tumor biology
lambda  : 0.004  : Gompertz growth rate (1/h)
TVmax   : 3.0    : Carrying capacity (cm3)
def4    : 0.2    : HR repair deficiency pathway 4 (0-1)
TV0     : 0.1    : Initial tumor volume (cm3)

$CMT @annotated
GDEPOT   : Drug depot (mg/kg)
GCENT    : Drug central (mg/kg)
EFF_FORK : PARPi fork trapping effect
TV       : Tumor volume (cm3)

$MAIN
if(NEWIND <= 1) {
  TV_0 = TV0;
}

$ODE
// Mouse PK (1-comp oral)
dxdt_GDEPOT = -KAm * GDEPOT;
dxdt_GCENT  = KAm * GDEPOT - (CLm/Vm) * GCENT;
double Cplasma = GCENT / Vm;
if(Cplasma < 0) Cplasma = 0;

// Fork trapping effect (indirect response, Eq 8 from paper)
dxdt_EFF_FORK = ip_fork * Cplasma - clpf * EFF_FORK;

// PK-driven kill rate with HR deficiency sensitization
// Synthetic lethality: def4 amplifies drug kill (BRCA-mut tumors are more sensitive)
// Fork trapping adds to effective kill
double drug_kill = kmax * Cplasma / (KC50 + Cplasma) * (1.0 + sens * def4) * (1.0 + 0.5 * EFF_FORK);

// Gompertz tumor growth - drug kill
double growth = 0.0;
if(TV > 1e-6 && TV < TVmax) {
  growth = lambda * TV * log(TVmax / TV);
} else if(TV >= TVmax) {
  growth = 0.0;
}
dxdt_TV = growth - drug_kill * TV;

// Floor at tiny value
if(TV < 1e-6 && dxdt_TV < 0) dxdt_TV = 0;

$TABLE
double TVout = TV;
double Cpl = GCENT / Vm;

$CAPTURE TVout Cpl
'

mod_pk  <- mcode("rucaparib_pk", pk_code)
mod_tgi <- mcode("rucaparib_tgi", tgi_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; }
  .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: #f0f4ff; border-left: 4px solid #8b5cf6; padding: 12px 16px; border-radius: 4px; margin-top: 10px; font-size: 13px; }
  .ref-box a { color: #8b5cf6; }
")

# ── UI ─────────────────────────────────────────────────────────────────────
ui <- page_sidebar(
  title = "Rucaparib PK/PD Simulator — PARP Inhibitor",
  theme = app_theme,
  sidebar = sidebar(
    title = "Simulation Settings", width = 360,
    conditionalPanel(
      condition = "input.tabs == 'Clinical PK'",
      h6("Dosing Regimen"),
      sliderInput("dose_pk", "Dose (mg)", min = 200, max = 800, value = 600, step = 100),
      radioButtons("regimen", "Regimen", choices = c("BID (twice daily)" = 12, "QD (once daily)" = 24), selected = 12, inline = TRUE),
      numericInput("n_days_pk", "Duration (days)", value = 14, min = 1, max = 28),
      hr(), h6("Patient Characteristics"),
      sliderInput("wt", "Weight (kg)", min = 40, max = 120, value = 70, step = 5),
      sliderInput("hepf", "Hepatic Function", min = 0.3, max = 1.0, value = 1.0, step = 0.1),
      hr(),
      checkboxInput("log_scale", "Log Scale (Y-axis)", value = FALSE),
      checkboxInput("show_target", "Show Target Ctrough (1500 ng/mL)", value = TRUE)
    ),
    conditionalPanel(
      condition = "input.tabs == 'Tumor Growth Inhibition'",
      h6("Rucaparib Dosing (Preclinical)"),
      sliderInput("dose_tgi", "Rucaparib Dose (mg/kg)", min = 0, max = 150, value = 50, step = 10),
      radioButtons("schedule_tgi", "Schedule",
        choices = c("QD (daily)" = 24, "BID (twice daily)" = 12), selected = 24, inline = TRUE),
      numericInput("n_days_tgi", "Treatment Duration (days)", value = 35, min = 7, max = 56),
      hr(), h6("Tumor Characteristics"),
      sliderInput("def4_tgi", "HR Deficiency (def4)", min = 0, max = 1.0, value = 0.2, step = 0.05,
        post = ""),
      helpText(style = "font-size:11px; margin-top:-8px;",
        "0% = HR proficient (HRD−) | 20-30% = HRD+ | 80-100% = BRCA mutant"),
      sliderInput("tv0_tgi", "Initial Tumor Volume (cm³)", min = 0.05, max = 0.3, value = 0.1, step = 0.01),
      hr(),
      checkboxInput("show_vehicle", "Show Vehicle (Untreated)", value = TRUE)
    )
  ),
  navset_card_tab(
    id = "tabs",
    title = "Rucaparib PK/PD",
    full_screen = TRUE,

    # ── Tab 1: Clinical PK ──────────────────────────────────────────────
    nav_panel("Clinical PK",
      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,ss (ng/mL)")),
        div(class = "metric-card metric-warning",
          div(class = "metric-value", textOutput("ctrough")),
          div(class = "metric-label", "Ctrough,ss (ng/mL)")),
        div(class = "metric-card metric-primary",
          div(class = "metric-value", textOutput("auc")),
          div(class = "metric-label", HTML("AUC<sub>0-&tau;</sub> (&mu;g&middot;h/mL)"))),
        div(class = "metric-card metric-info",
          div(class = "metric-value", textOutput("thalf")),
          div(class = "metric-label", "t½ (h)"))
      ),
      plotOutput("pkPlot", height = "500px")
    ),

    # ── Tab 2: Tumor Growth Inhibition ──────────────────────────────────
    nav_panel("Tumor Growth Inhibition",
      layout_column_wrap(
        width = 1/4,
        fill = FALSE,
        div(class = "metric-card metric-success",
          div(class = "metric-value", textOutput("tgi_percent")),
          div(class = "metric-label", "TGI (%)")),
        div(class = "metric-card metric-warning",
          div(class = "metric-value", textOutput("tv_final")),
          div(class = "metric-label", "Tumor Vol at EOT (cm³)")),
        div(class = "metric-card metric-primary",
          div(class = "metric-value", textOutput("tv_vehicle")),
          div(class = "metric-label", "Vehicle Vol at EOT (cm³)")),
        div(class = "metric-card metric-info",
          div(class = "metric-value", textOutput("hrd_status")),
          div(class = "metric-label", "HRD Status"))
      ),
      plotOutput("tgiPlot", height = "500px")
    ),

    # ── Tab 3: Model Information ────────────────────────────────────────
    nav_panel("Model Information", markdown("
## Rucaparib — Semi-Mechanistic PARP Inhibitor PK/PD Model

### Clinical Pharmacokinetics
**Structure:** Two-compartment model with first-order oral absorption and first-order elimination  
**CL/F:** 37.6 L/h (apparent oral clearance, allometric scaling on weight)  
**V1/F:** 310 L (central volume)  
**V2/F:** 230 L (peripheral volume)  
**Q/F:** 16 L/h (inter-compartmental clearance)  
**Ka:** 1.15 h⁻¹  
**Bioavailability:** 36% (absolute, range 30–45%)  
**t½:** ~17 h (terminal half-life)  
**Tmax:** ~1.9 h (median)  
**Steady state:** ~7 days with BID dosing  

### Mechanism of Action: DNA Damage Response (DDR) Model
Rucaparib inhibits PARP enzymes, which play a key role in:
- **SSB repair** — blocking single-strand break repair pathway 1
- **DSB repair** — impairing double-strand break repair pathway 3
- **PARP trapping** — accelerating conversion of SSBs to DSBs at replication forks

#### Synthetic Lethality
In tumors with **HR deficiency** (e.g., BRCA mutation), the alternative DSB repair
pathway (pathway 4) is already impaired. When PARP inhibitor blocks pathway 3,
both DSB repair pathways are compromised → DSB accumulation → cell death.

#### Key DDR Model Equations
- SSB dynamics: dSSB/dt = r_ssb − R_conv·SSB − R_ssb·SSB  
- DSB dynamics: dDSB/dt = r_dsb + R_conv·SSB − R_dsb·DSB  
- Drug effect: dX/dt = u·C_drug − Cl_x·X (indirect response)  
- Tumor growth: dTV/dt = k_g·TV − k_kill·DSB·TV  

### Therapeutic Context
- **Approved dose:** 600 mg BID (oral tablets)  
- **Indications:** Recurrent ovarian cancer (maintenance & treatment), mCRPC with BRCA mutation  
- **Target Ctrough:** ≥1500 ng/mL associated with clinical efficacy  
- **Metabolism:** CYP2D6 (primary), CYP1A2 and CYP3A4 (minor)  
- **Key toxicities:** Fatigue, nausea, anemia, thrombocytopenia  
    ")),

    # ── Tab 4: References ───────────────────────────────────────────────
    nav_panel("References", div(class = "ref-box",
      tags$h5(HTML("&#128218; Key References")),
      tags$ol(
        tags$li(HTML("Villette CC et al. (2025) <em>Semi-mechanistic efficacy model for PARP + ATR inhibitors — application to rucaparib and talazoparib in combination with gartisertib in breast cancer PDXs.</em> Br J Cancer 132:481–491. ",
          tags$a(href = "https://doi.org/10.1038/s41416-024-02935-w", target = "_blank", "doi:10.1038/s41416-024-02935-w"))),
        tags$li(HTML("Liao M et al. (2022) <em>Clinical Pharmacokinetics and Pharmacodynamics of Rucaparib.</em> Clin Pharmacokinet 61:1477–1493. ",
          tags$a(href = "https://doi.org/10.1007/s40262-022-01157-8", target = "_blank", "doi:10.1007/s40262-022-01157-8"))),
        tags$li(HTML("Xiao JJ et al. (2022) <em>Population pharmacokinetics of rucaparib in patients with advanced ovarian cancer or other solid tumors.</em> Cancer Chemother Pharmacol 90:55–65. ",
          tags$a(href = "https://doi.org/10.1007/s00280-022-04413-7", target = "_blank", "doi:10.1007/s00280-022-04413-7"))),
        tags$li(HTML("Shapiro GI et al. (2019) <em>Pharmacokinetic study of rucaparib in patients with advanced solid tumors.</em> Clin Pharmacol Drug Dev 8(1):107–118. ",
          tags$a(href = "https://doi.org/10.1002/cpdd.575", target = "_blank", "doi:10.1002/cpdd.575")))
      ),
      tags$h5(HTML("&#128138; Therapeutic Context")),
      tags$ul(
        tags$li(tags$strong("Class:"), " PARP1/2/3 inhibitor (DDR-targeted agent)"),
        tags$li(tags$strong("Brand name:"), " Rubraca"),
        tags$li(tags$strong("Approved dose:"), " 600 mg BID, continuous dosing"),
        tags$li(tags$strong("Route:"), " Oral (tablet)"),
        tags$li(tags$strong("Indications:"), " Recurrent ovarian cancer (maintenance/treatment), mCRPC (BRCA-mutant)"),
        tags$li(tags$strong("Biomarkers:"), " BRCA1/2 mutation, HRD status, PAR levels"),
        tags$li(tags$strong("Metabolism:"), " CYP2D6 (primary), with CYP1A2/3A4 (minor); inhibits CYP1A2, MATE1/2-K, OCT1/2")
      )
    ))
  ),

  # ── Footer ─────────────────────────────────────────────────────────────
  div(style = "text-align: center; padding: 20px; margin-top: 30px; 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's 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) {

  # ── Clinical PK reactive ───────────────────────────────────────────────
  pk_data <- reactive({
    tau <- as.numeric(input$regimen)
    n_doses <- ceiling(input$n_days_pk * 24 / tau)
    ev1 <- ev(amt = input$dose_pk, cmt = 1, ii = tau, addl = n_doses - 1)
    mod_pk %>%
      param(WT = input$wt, HEPF = input$hepf) %>%
      ev(ev1) %>%
      mrgsim(end = input$n_days_pk * 24, delta = 0.25) %>%
      as.data.frame() %>%
      mutate(time_h = time, time_d = time / 24)
  })

  # ── PK metrics ─────────────────────────────────────────────────────────
  output$cmax <- renderText({
    d <- pk_data()
    tau <- as.numeric(input$regimen)
    last_start <- (input$n_days_pk - 1) * 24
    last <- d %>% filter(time >= last_start)
    sprintf("%.0f", max(last$CPng, na.rm = TRUE))
  })

  output$ctrough <- renderText({
    d <- pk_data()
    tau <- as.numeric(input$regimen)
    last_start <- (input$n_days_pk - 1) * 24
    last <- d %>% filter(time >= last_start)
    sprintf("%.0f", min(last$CPng[last$CPng > 0], na.rm = TRUE))
  })

  output$auc <- renderText({
    d <- pk_data()
    tau <- as.numeric(input$regimen)
    last_start <- (input$n_days_pk - 1) * 24
    last <- d %>% filter(time >= last_start, time <= last_start + tau)
    if (nrow(last) < 2) return("—")
    auc <- sum(diff(last$time) * (head(last$CP, -1) + tail(last$CP, -1)) / 2)
    sprintf("%.1f", auc)
  })

  output$thalf <- renderText({
    CLi <- 37.6 * (input$wt / 70)^0.75 * input$hepf
    V1i <- 130 * (input$wt / 70)
    V2i <- 250 * (input$wt / 70)
    Qi  <- 15  * (input$wt / 70)^0.75
    ab_sum  <- CLi/V1i + Qi/V1i + Qi/V2i
    ab_prod <- CLi * Qi / (V1i * V2i)
    beta_val <- (ab_sum - sqrt(ab_sum^2 - 4 * ab_prod)) / 2
    sprintf("%.1f", log(2) / beta_val)
  })

  # ── PK Plot ────────────────────────────────────────────────────────────
  output$pkPlot <- renderPlot({
    d <- pk_data()
    tau <- as.numeric(input$regimen)
    reg_label <- if (tau == 12) "BID" else "QD"

    p <- ggplot(d, aes(x = time_d, y = CPng)) +
      geom_line(color = "#8b5cf6", linewidth = 0.9)

    if (input$show_target) {
      p <- p +
        geom_hline(yintercept = 1500, linetype = "dashed", color = "#10b981", alpha = 0.7) +
        annotate("text", x = max(d$time_d) * 0.02, y = 1600,
          label = "Target Ctrough: 1500 ng/mL", hjust = 0, size = 3.5, color = "#10b981")
    }

    p <- p +
      labs(
        x = "Time (days)", y = "Plasma Concentration (ng/mL)",
        title = paste0("Rucaparib ", input$dose_pk, " mg ", reg_label,
          " | Wt: ", input$wt, " kg | Hepatic fn: ", input$hepf * 100, "%")
      ) +
      theme_minimal(base_size = 14) +
      theme(plot.title = element_text(face = "bold", size = 13))

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

  # ── TGI reactive ──────────────────────────────────────────────────────
  tgi_data <- reactive({
    tau <- as.numeric(input$schedule_tgi)
    n_doses <- ceiling(input$n_days_tgi * 24 / tau)
    total_h <- (input$n_days_tgi + 14) * 24  # simulate 2 weeks beyond treatment

    # Treated
    ev_treat <- ev(amt = input$dose_tgi, cmt = 1, ii = tau, addl = n_doses - 1)
    sim_treat <- mod_tgi %>%
      param(def4 = input$def4_tgi, TV0 = input$tv0_tgi) %>%
      ev(ev_treat) %>%
      mrgsim(end = total_h, delta = 2) %>%
      as.data.frame() %>%
      mutate(time_d = time / 24, group = "Rucaparib")

    # Vehicle (no treatment)
    ev_veh <- ev(amt = 0, cmt = 1)
    sim_veh <- mod_tgi %>%
      param(def4 = input$def4_tgi, TV0 = input$tv0_tgi) %>%
      ev(ev_veh) %>%
      mrgsim(end = total_h, delta = 2) %>%
      as.data.frame() %>%
      mutate(time_d = time / 24, group = "Vehicle")

    list(treat = sim_treat, vehicle = sim_veh)
  })

  # ── TGI metrics ────────────────────────────────────────────────────────
  output$tgi_percent <- renderText({
    d <- tgi_data()
    eot <- input$n_days_tgi
    tv_treat <- d$treat %>% filter(abs(time_d - eot) < 0.5) %>% pull(TVout) %>% tail(1)
    tv_veh   <- d$vehicle %>% filter(abs(time_d - eot) < 0.5) %>% pull(TVout) %>% tail(1)
    if (length(tv_veh) == 0 || tv_veh <= 0) return("\u2014")
    tgi <- (1 - tv_treat / tv_veh) * 100
    if (tgi < 0) tgi <- 0
    sprintf("%.0f%%", tgi)
  })

  output$tv_final <- renderText({
    d <- tgi_data()
    eot <- input$n_days_tgi
    tv <- d$treat %>% filter(abs(time_d - eot) < 0.5) %>% pull(TVout) %>% tail(1)
    if (length(tv) == 0) return("\u2014")
    sprintf("%.2f", tv)
  })

  output$tv_vehicle <- renderText({
    d <- tgi_data()
    eot <- input$n_days_tgi
    tv <- d$vehicle %>% filter(abs(time_d - eot) < 0.5) %>% pull(TVout) %>% tail(1)
    if (length(tv) == 0) return("\u2014")
    sprintf("%.2f", tv)
  })

  output$hrd_status <- renderText({
    d4 <- input$def4_tgi
    if (d4 >= 0.7) return("BRCA-mutant")
    if (d4 >= 0.2) return("HRD+")
    return("HRD\u2212")
  })

  # ── TGI Plot ───────────────────────────────────────────────────────────
  output$tgiPlot <- renderPlot({
    d <- tgi_data()
    eot <- input$n_days_tgi
    tau <- as.numeric(input$schedule_tgi)
    sched_label <- if (tau == 12) "BID" else "QD"
    def_label <- paste0("def4=", round(input$def4_tgi * 100), "%")

    plot_data <- d$treat
    if (input$show_vehicle) {
      plot_data <- bind_rows(d$treat, d$vehicle)
    }

    colors <- c("Rucaparib" = "#8b5cf6", "Vehicle" = "#95a5a6")

    p <- ggplot(plot_data, aes(x = time_d, y = TVout, color = group)) +
      geom_vline(xintercept = eot, linetype = "dotted", color = "gray50") +
      annotate("text", x = eot + 0.5, y = max(plot_data$TVout, na.rm = TRUE) * 0.95,
        label = "End of\ntreatment", size = 3, hjust = 0, color = "gray50") +
      geom_line(linewidth = 1) +
      scale_color_manual(values = colors) +
      labs(
        x = "Time (days)", y = expression(Tumor~Volume~(cm^3)),
        title = paste0("Rucaparib ", input$dose_tgi, " mg/kg ", sched_label,
          " | HRD: ", def_label, " | Tdoub: ", input$tdoub_tgi, "h"),
        color = NULL
      ) +
      theme_minimal(base_size = 14) +
      theme(
        plot.title = element_text(face = "bold", size = 13),
        legend.position = "top"
      )
    p
  })
}

shinyApp(ui = ui, server = server)