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

# ── mrgsolve Model ──
# Delattre et al. (2010) - Amikacin PopPK in Critically Ill Septic Patients
# 2-compartment model with first-order elimination, IV infusion (30 min)
# Covariate: Creatinine clearance (Cockcroft-Gault) on CL
model_code <- '
$PARAM @annotated
V1     : 19.2  : Central volume of distribution (L)
V2     : 9.38  : Peripheral volume of distribution (L)
Q      : 4.38  : Intercompartmental clearance (L/h)
TVCL   : 0.77  : Base clearance intercept (L/h)
THETA2 : 1.42  : CLCR effect on CL (fractional)
CLCR   : 55.5  : Creatinine clearance - Cockcroft-Gault (mL/min)
REF_CLCR : 55.5 : Reference CLCR median (mL/min)

$CMT @annotated
CENT   : Central compartment (mg)
PERIPH : Peripheral compartment (mg)

$MAIN
double CLi = TVCL + THETA2 * (CLCR / REF_CLCR);
if(CLi < 0.1) CLi = 0.1;

$ODE
dxdt_CENT   = -(CLi/V1 + Q/V1) * CENT + (Q/V2) * PERIPH;
dxdt_PERIPH = (Q/V1) * CENT - (Q/V2) * PERIPH;

$TABLE
double CP = CENT / V1;
CP = CP > 0 ? CP : 0;

$CAPTURE @annotated
CP   : Plasma concentration (mg/L)
CLi  : Individual clearance (L/h)
'

mod <- mcode("amikacin_delattre2010", 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;
  }
  .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; }
  .target-line { font-size: 11px; margin-top: 2px; }
")

# ── UI ──
ui <- page_sidebar(
  title = "Amikacin PK Simulator — Critically Ill Septic Patients",
  theme = app_theme,

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

    h6("Patient Characteristics"),
    sliderInput("weight", "Body Weight (kg)", min = 38, max = 130, value = 70, step = 1),
    sliderInput("clcr", "Creatinine Clearance (mL/min)", min = 10, max = 200, value = 55.5, step = 0.5),
    helpText("Cockcroft-Gault estimate"),

    hr(),
    h6("Dosing"),
    numericInput("dose_mg_kg", "Dose (mg/kg)", value = 25, min = 5, max = 50, step = 1),
    numericInput("tinf", "Infusion Duration (min)", value = 30, min = 10, max = 120, step = 5),
    numericInput("n_doses", "Number of Doses", value = 1, min = 1, max = 10, step = 1),
    numericInput("interval", "Dosing Interval (h)", value = 24, min = 8, max = 72, step = 4),

    hr(),
    checkboxInput("log_scale", "Log Scale (Y-axis)", value = FALSE),
    checkboxInput("show_targets", "Show Therapeutic Targets", value = TRUE),

    hr(),
    card(
      card_header("Population PK Parameters", class = "bg-light"),
      p(strong("V1:"), "19.2 L (CV 39.1%)"),
      p(strong("V2:"), "9.38 L (CV 43.6%)"),
      p(strong("Q:"), "4.38 L/h"),
      p(strong("CL:"), "0.77 + 1.42 × (CLCR/55.5) L/h"),
      p(strong("CL IIV:"), "59.2%"),
      p(em("Delattre et al., Ther Drug Monit 2010"))
    )
  ),

  # ── Metrics ──
  layout_column_wrap(
        width = 1/4,
        fill = FALSE,
    div(class = "metric-card metric-success",
      div(class = "metric-value", textOutput("cmax_val")),
      div(class = "metric-label", "Cmax (mg/L)"),
      div(class = "target-line", "Target: 20–40 mg/L")
    ),
    div(class = "metric-card metric-warning",
      div(class = "metric-value", textOutput("ctrough_val")),
      div(class = "metric-label", "Ctrough (mg/L)"),
      div(class = "target-line", "Target: <5 mg/L")
    ),
    div(class = "metric-card metric-primary",
      div(class = "metric-value", textOutput("auc_val")),
      div(class = "metric-label", "AUC₀₋₂₄ (mg·h/L)")
    ),
    div(class = "metric-card metric-info",
      div(class = "metric-value", textOutput("thalf_val")),
      div(class = "metric-label", "Terminal t½ (h)")
    )
  ),

  # ── PK Plot ──
  card(
    card_header("Concentration-Time Profile"),
    full_screen = TRUE,
    plotOutput("pkPlot", height = "450px")
  ),

  # ── Tabs ──
  navset_card_tab(
    title = "Analysis",
    nav_panel("Dose Comparison",
      layout_columns(
        col_widths = c(6, 6),
        div(
          numericInput("compare_dose1", "Dose 1 (mg/kg)", value = 20, min = 5, max = 50),
          numericInput("compare_dose2", "Dose 2 (mg/kg)", value = 25, min = 5, max = 50),
          numericInput("compare_dose3", "Dose 3 (mg/kg)", value = 30, min = 5, max = 50)
        ),
        plotOutput("comparePlot", height = "400px")
      )
    ),
    nav_panel("CLCR Impact",
      plotOutput("clcrPlot", height = "400px")
    ),
    nav_panel("Model Information", 
      markdown(paste0(
        "## Amikacin Population PK Model\n\n",
        "**Source:** Delattre IK, Musuamba FT, Nyberg J, et al. ",
        "*Population Pharmacokinetic Modeling and Optimal Sampling Strategy ",
        "for Bayesian Estimation of Amikacin Exposure in Critically Ill Septic Patients.* ",
        "Ther Drug Monit. 2010;32(6):749-756.\n\n",
        "### Study Population\n",
        "- 88 critically ill septic patients (severe sepsis/septic shock)\n",
        "- Age: 22-89 years (median 65)\n",
        "- Weight: 38-125 kg (median 70)\n",
        "- CLCR: 12.3-408.3 mL/min (median 55.5)\n",
        "- 71% septic shock, 52% mechanically ventilated\n",
        "- 507 blood samples during first 24h of treatment\n\n",
        "### Model Structure\n",
        "- Two-compartment model with first-order elimination\n",
        "- IV infusion over 30 minutes\n",
        "- Covariate: Creatinine clearance (Cockcroft-Gault) on CL\n",
        "- CL = 0.77 + 1.42 × (CLCR / 55.5) L/h\n\n",
        "### Therapeutic Targets\n",
        "- **Peak (Cmax):** 20-40 mg/L (Cmax/MIC ≥ 8-10 for efficacy)\n",
        "- **Trough:** <5 mg/L (minimize nephro/ototoxicity)\n",
        "- Once-daily dosing (extended interval) preferred in sepsis\n\n",
        "### References\n",
        "- [Amikacin FDA Label (DailyMed)](https://dailymed.nlm.nih.gov/dailymed/drugInfo.cfm?setid=c0f57839-1c9b-49e5-8c7a-708e2d16495d)\n",
        "- [Amikacin Prescribing Information (Drugs.com)](https://www.drugs.com/pro/amikacin.html)\n",
        "- Delattre et al., Ther Drug Monit 2010;32:749-756\n"
      ))
    )
  ),

  # ── PKPDBuilder Branding 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)",
    br(),
    tags$span(style = "font-size: 10px;", 
      "For research and educational purposes only. Not for clinical decision-making.")
  )
)

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

  sim_data <- reactive({
    shiny::req(input$weight, input$clcr, input$dose_mg_kg, input$tinf, input$n_doses)

    dose_mg <- input$dose_mg_kg * input$weight
    tinf_h <- input$tinf / 60
    rate <- dose_mg / tinf_h

    ev <- ev(amt = dose_mg, rate = rate, cmt = 1, ii = input$interval, addl = input$n_doses - 1)
    sim_hours <- input$interval * input$n_doses + 12

    out <- mod %>%
      param(CLCR = input$clcr) %>%
      ev(ev) %>%
      mrgsim(end = sim_hours, delta = 0.1) %>%
      as.data.frame()

    out
  })

  # ── Metrics ──
  output$cmax_val <- renderText({
    d <- sim_data()
    sprintf("%.1f", max(d$CP, na.rm = TRUE))
  })

  output$ctrough_val <- renderText({
    d <- sim_data()
    # Trough = concentration at end of first dosing interval
    t_trough <- input$interval
    idx <- which.min(abs(d$time - t_trough))
    sprintf("%.2f", d$CP[idx])
  })

  output$auc_val <- renderText({
    d <- sim_data()
    d24 <- d %>% dplyr::filter(time <= 24) %>% dplyr::arrange(time)
    # Trapezoidal rule without pracma
    n <- nrow(d24)
    if(n < 2) return("--")
    auc <- sum(diff(d24$time) * (d24$CP[-n] + d24$CP[-1]) / 2)
    sprintf("%.0f", auc)
  })

  output$thalf_val <- renderText({
    # Terminal half-life from eigenvalues of 2-cpt model
    cl <- 0.77 + 1.42 * (input$clcr / 55.5)
    if(cl < 0.1) cl <- 0.1
    v1 <- 19.2; v2 <- 9.38; q <- 4.38
    k10 <- cl / v1; k12 <- q / v1; k21 <- q / v2
    beta_sum <- k10 + k12 + k21
    beta_prod <- k10 * k21
    disc <- beta_sum^2 - 4 * beta_prod
    disc <- max(disc, 0)
    beta <- (beta_sum - sqrt(disc)) / 2
    thalf <- ifelse(beta > 0, log(2) / beta, NA)
    sprintf("%.1f", thalf)
  })

  # ── PK Plot ──
  output$pkPlot <- renderPlot({
    d <- sim_data()

    if (input$log_scale) {
      d <- d %>% dplyr::filter(CP > 0.001)
    }

    p <- ggplot(d, aes(x = time, y = CP)) +
      geom_line(color = "#8b5cf6", linewidth = 1.2) +
      labs(x = "Time (hours)", y = "Amikacin Concentration (mg/L)",
           title = paste0("Amikacin ", input$dose_mg_kg, " mg/kg IV (",
                          input$tinf, "-min infusion) | CLCR = ",
                          input$clcr, " mL/min | WT = ", input$weight, " kg")) +
      theme_minimal(base_size = 14) +
      theme(plot.title = element_text(size = 13))

    if (input$show_targets) {
      p <- p +
        geom_hline(yintercept = 20, linetype = "dashed", color = "#10b981", alpha = 0.7) +
        geom_hline(yintercept = 40, linetype = "dashed", color = "#10b981", alpha = 0.7) +
        geom_hline(yintercept = 5, linetype = "dashed", color = "#ef4444", alpha = 0.7) +
        annotate("text", x = max(d$time) * 0.85, y = 21, label = "Peak target (20 mg/L)",
                 color = "#10b981", size = 3.5) +
        annotate("text", x = max(d$time) * 0.85, y = 41, label = "Peak target (40 mg/L)",
                 color = "#10b981", size = 3.5) +
        annotate("text", x = max(d$time) * 0.85, y = 6, label = "Trough limit (5 mg/L)",
                 color = "#ef4444", size = 3.5)
    }

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

    p
  })

  # ── Dose Comparison ──
  output$comparePlot <- renderPlot({
    shiny::req(input$compare_dose1, input$compare_dose2, input$compare_dose3)

    doses <- c(input$compare_dose1, input$compare_dose2, input$compare_dose3)
    tinf_h <- input$tinf / 60
    sim_hours <- input$interval + 12

    all_data <- lapply(doses, function(d_mg_kg) {
      dose_mg <- d_mg_kg * input$weight
      rate <- dose_mg / tinf_h
      ev_obj <- ev(amt = dose_mg, rate = rate, cmt = 1)
      out <- mod %>%
        param(CLCR = input$clcr) %>%
        ev(ev_obj) %>%
        mrgsim(end = sim_hours, delta = 0.1) %>%
        as.data.frame() %>%
        mutate(Dose = paste0(d_mg_kg, " mg/kg"))
      out
    })
    combined <- bind_rows(all_data)

    p <- ggplot(combined, aes(x = time, y = CP, color = Dose)) +
      geom_line(linewidth = 1.1) +
      geom_hline(yintercept = 20, linetype = "dashed", color = "#10b981", alpha = 0.5) +
      geom_hline(yintercept = 5, linetype = "dashed", color = "#ef4444", alpha = 0.5) +
      labs(x = "Time (hours)", y = "Concentration (mg/L)",
           title = "Dose Comparison") +
      scale_color_manual(values = c("#8b5cf6", "#f59e0b", "#0dcaf0")) +
      theme_minimal(base_size = 14)

    p
  })

  # ── CLCR Impact Plot ──
  output$clcrPlot <- renderPlot({
    clcr_vals <- c(20, 40, 60, 80, 120)
    dose_mg <- input$dose_mg_kg * input$weight
    tinf_h <- input$tinf / 60
    rate <- dose_mg / tinf_h
    sim_hours <- 36

    all_data <- lapply(clcr_vals, function(cr) {
      ev_obj <- ev(amt = dose_mg, rate = rate, cmt = 1)
      out <- mod %>%
        param(CLCR = cr) %>%
        ev(ev_obj) %>%
        mrgsim(end = sim_hours, delta = 0.1) %>%
        as.data.frame() %>%
        mutate(CLCR_label = paste0("CLCR = ", cr, " mL/min"))
      out
    })
    combined <- bind_rows(all_data)

    ggplot(combined, aes(x = time, y = CP, color = CLCR_label)) +
      geom_line(linewidth = 1) +
      geom_hline(yintercept = 5, linetype = "dashed", color = "#ef4444", alpha = 0.5) +
      labs(x = "Time (hours)", y = "Concentration (mg/L)",
           title = "Impact of Creatinine Clearance on Amikacin PK",
           color = "CLCR") +
      theme_minimal(base_size = 14)
  })
}

shinyApp(ui = ui, server = server)