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

# ── Empagliflozin PopPK Model ──
# Baron et al., Diabetes Ther (2016) 7:455-471
# 2-compartment model with lagged first-order absorption
# Key covariate: eGFR on CL/F

model_code <- '
$PARAM @annotated
TVCL   : 10.6   : Typical CL/F (L/h) - reference patient
TVV2   : 3.14   : Typical central volume V2/F (L)
TVQ    : 6.34   : Typical intercompartmental CL Q/F (L/h)
TVV3   : 70.6   : Typical peripheral volume V3/F (L)
TVKA   : 0.196  : Typical absorption rate constant (1/h)
APTS   : 0.5    : Absorption lag time (h) - fixed
EGFR   : 100    : eGFR (mL/min/1.73m2)
REF_EGFR : 100  : Reference eGFR
THETA_EGFR_CL : 0.449 : eGFR effect on CL/F (power)
SEX    : 0      : Sex (0=male, 1=female)
THETA_SEX_CL  : 0.066 : Female effect on CL/F (proportional)
AGE    : 50     : Age (years)
BMI    : 25     : BMI (kg/m2)

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

$MAIN
double CLi = TVCL * pow(EGFR / REF_EGFR, THETA_EGFR_CL);
if (SEX == 1) CLi = CLi * (1.0 + THETA_SEX_CL);
double V2i = TVV2;
double V3i = TVV3;
double Qi  = TVQ;
double KAi = TVKA;
ALAG_DEPOT = APTS;

$ODE
dxdt_DEPOT  = -KAi * DEPOT;
dxdt_CENT   = KAi * DEPOT - (CLi/V2i + Qi/V2i) * CENT + (Qi/V3i) * PERIPH;
dxdt_PERIPH = (Qi/V2i) * CENT - (Qi/V3i) * PERIPH;

$TABLE
double CP = CENT / V2i;
double CP_nmol = CP * 1000.0 / 450.91;

$CAPTURE @annotated
CP       : Plasma concentration (mg/L)
CP_nmol  : Plasma concentration (nmol/L)
CLi      : Individual CL/F (L/h)
'

mod <- mcode("empagliflozin", 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; }
")

# ── UI ──
ui <- page_sidebar(
  title = "Empagliflozin PopPK Simulator",
  theme = app_theme,

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

    h6("Patient Characteristics"),
    sliderInput("egfr", "eGFR (mL/min/1.73m\u00B2)",
                min = 15, max = 150, value = 100, step = 5),
    selectInput("sex", "Sex", choices = c("Male" = 0, "Female" = 1)),
    sliderInput("age", "Age (years)", min = 18, max = 90, value = 50),

    hr(),
    h6("Dosing"),
    selectInput("dose", "Empagliflozin Dose",
                choices = c("10 mg" = 10, "25 mg" = 25), selected = 25),
    numericInput("n_days", "Treatment Duration (days)", value = 14, min = 1, max = 90),

    hr(),
    checkboxInput("log_scale", "Log Scale (Y-axis)", value = FALSE),
    checkboxInput("show_both", "Compare 10 mg vs 25 mg", value = FALSE),

    hr(),
    card(
      card_header("Model Parameters", class = "bg-light"),
      p(strong("CL/F:"), "10.6 L/h (ref)"),
      p(strong("V2/F:"), "3.14 L (central)"),
      p(strong("V3/F:"), "70.6 L (peripheral)"),
      p(strong("Q/F:"), "6.34 L/h"),
      p(strong("Ka:"), "0.196 h\u207B\u00B9"),
      p(strong("Tlag:"), "0.5 h"),
      p(em("Key covariate: eGFR on CL/F"))
    )
  ),

  # Metric cards
  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,ss (mg/L)")
    ),
    div(class = "metric-card metric-warning",
      div(class = "metric-value", textOutput("ctrough_val")),
      div(class = "metric-label", "Ctrough,ss (mg/L)")
    ),
    div(class = "metric-card metric-primary",
      div(class = "metric-value", textOutput("auc_val")),
      div(class = "metric-label", "AUCss,24h (mg\u00B7h/L)")
    ),
    div(class = "metric-card metric-info",
      div(class = "metric-value", textOutput("cl_val")),
      div(class = "metric-label", "CL/F (L/h)")
    )
  ),

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

  # Tabbed content
  navset_card_tab(
    title = "Analysis",
    nav_panel("eGFR Impact",
      plotOutput("egfrPlot", height = "400px")
    ),
    nav_panel("Dose Comparison",
      plotOutput("dosePlot", height = "400px")
    ),
    nav_panel("Model Information",
      markdown(paste0(
        "## Empagliflozin Population PK Model\n\n",
        "**Source:** Baron KT, Macha S, Broedl UC, Nock V, Retlich S, Riggs MT. ",
        "Population Pharmacokinetics and Exposure-Response (Efficacy and Safety/Tolerability) ",
        "of Empagliflozin in Patients with Type 2 Diabetes. *Diabetes Ther*. 2016;7:455-471.\n\n",
        "**DOI:** [10.1007/s13300-016-0174-y](https://doi.org/10.1007/s13300-016-0174-y)\n\n",
        "### Model Structure\n",
        "- Two-compartment disposition\n",
        "- Lagged first-order absorption (Tlag = 0.5 h)\n",
        "- First-order elimination\n\n",
        "### Key Findings\n",
        "- eGFR is the most prominent covariate on CL/F\n",
        "- AUCss increases ~88% with severe renal impairment (eGFR 15 mL/min/1.73m\u00B2)\n",
        "- No clinically relevant effects of age, sex, BMI, or race\n",
        "- Doses of 10 mg and 25 mg QD achieve near-maximal efficacy\n\n",
        "### Therapeutic Context\n",
        "- **Brand Name:** Jardiance\u00AE (Boehringer Ingelheim)\n",
        "- **Class:** SGLT2 Inhibitor\n",
        "- **Indications:** Type 2 Diabetes, Heart Failure, CKD\n",
        "- **Approved Doses:** 10 mg and 25 mg once daily\n",
        "- [FDA Label](https://www.accessdata.fda.gov/drugsatfda_docs/label/2025/204629s063lbl.pdf)\n\n",
        "### References\n",
        "1. Baron et al. (2016) Diabetes Ther 7:455-471 (this model)\n",
        "2. Macha et al. (2014) Clin Pharmacokinet 53:949-957 (prior PopPK)\n",
        "3. [FDA Label - Jardiance](https://www.accessdata.fda.gov/drugsatfda_docs/label/2025/204629s063lbl.pdf)\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({
    dose_mg <- as.numeric(input$dose)
    n_days <- input$n_days

    out <- mod %>%
      param(EGFR = input$egfr, SEX = as.numeric(input$sex), AGE = input$age) %>%
      ev(amt = dose_mg, cmt = 1, ii = 24, addl = n_days - 1) %>%
      mrgsim(end = n_days * 24, delta = 0.25) %>%
      as.data.frame()

    out$time_days <- out$time / 24
    out
  })

  compare_data <- reactive({
    shiny::req(input$show_both)
    n_days <- input$n_days

    bind_rows(
      mod %>%
        param(EGFR = input$egfr, SEX = as.numeric(input$sex)) %>%
        ev(amt = 10, cmt = 1, ii = 24, addl = n_days - 1) %>%
        mrgsim(end = n_days * 24, delta = 0.25) %>%
        as.data.frame() %>%
        mutate(Dose = "10 mg"),
      mod %>%
        param(EGFR = input$egfr, SEX = as.numeric(input$sex)) %>%
        ev(amt = 25, cmt = 1, ii = 24, addl = n_days - 1) %>%
        mrgsim(end = n_days * 24, delta = 0.25) %>%
        as.data.frame() %>%
        mutate(Dose = "25 mg")
    ) %>%
      mutate(time_days = time / 24)
  })

  # Metrics from last dosing interval (steady-state)
  ss_data <- reactive({
    d <- sim_data()
    n_days <- input$n_days
    last_start <- (n_days - 1) * 24
    d %>% dplyr::filter(time >= last_start)
  })

  output$cmax_val <- renderText({
    d <- ss_data()
    sprintf("%.2f", max(d$CP, na.rm = TRUE))
  })

  output$ctrough_val <- renderText({
    d <- ss_data()
    sprintf("%.3f", min(d$CP[d$CP > 0], na.rm = TRUE))
  })

  output$auc_val <- renderText({
    d <- ss_data()
    # Trapezoidal AUC over last 24h
    d <- d[order(d$time), ]
    auc <- sum(diff(d$time) * (head(d$CP, -1) + tail(d$CP, -1)) / 2)
    sprintf("%.1f", auc)
  })

  output$cl_val <- renderText({
    d <- ss_data()
    sprintf("%.1f", d$CLi[1])
  })

  output$pkPlot <- renderPlot({
    if (input$show_both) {
      d <- compare_data()
      d_plot <- if (input$log_scale) d %>% dplyr::filter(CP > 0.001) else d
      p <- ggplot(d_plot, aes(x = time_days, y = CP, color = Dose)) +
        geom_line(linewidth = 1) +
        scale_color_manual(values = c("10 mg" = "#8b5cf6", "25 mg" = "#10b981")) +
        labs(x = "Time (days)", y = "Concentration (mg/L)",
             title = "Empagliflozin PK: 10 mg vs 25 mg QD") +
        theme_minimal(base_size = 14) +
        theme(legend.position = "top")
    } else {
      d <- sim_data()
      d_plot <- if (input$log_scale) d %>% dplyr::filter(CP > 0.001) else d
      p <- ggplot(d_plot, aes(x = time_days, y = CP)) +
        geom_line(linewidth = 1, color = "#8b5cf6") +
        labs(x = "Time (days)", y = "Concentration (mg/L)",
             title = paste("Empagliflozin", input$dose, "mg QD PK Profile")) +
        theme_minimal(base_size = 14)
    }
    if (input$log_scale) p <- p + scale_y_log10()
    p
  })

  output$egfrPlot <- renderPlot({
    egfr_vals <- c(15, 30, 45, 60, 90, 120)
    dose_mg <- as.numeric(input$dose)

    all_sims <- bind_rows(lapply(egfr_vals, function(eg) {
      mod %>%
        param(EGFR = eg, SEX = as.numeric(input$sex)) %>%
        ev(amt = dose_mg, cmt = 1, ii = 24, addl = 13) %>%
        mrgsim(end = 14 * 24, delta = 0.5) %>%
        as.data.frame() %>%
        mutate(eGFR = paste0(eg, " mL/min/1.73m\u00B2"),
               time_days = time / 24)
    }))

    all_sims$eGFR <- factor(all_sims$eGFR,
      levels = paste0(egfr_vals, " mL/min/1.73m\u00B2"))

    d_plot <- if (input$log_scale) all_sims %>% dplyr::filter(CP > 0.001) else all_sims

    p <- ggplot(d_plot, aes(x = time_days, y = CP, color = eGFR)) +
      geom_line(linewidth = 0.8) +
      scale_color_viridis_d(option = "plasma") +
      labs(x = "Time (days)", y = "Concentration (mg/L)",
           title = paste("Impact of eGFR on Empagliflozin", dose_mg, "mg QD Exposure"),
           color = "eGFR") +
      theme_minimal(base_size = 14) +
      theme(legend.position = "right")
    if (input$log_scale) p <- p + scale_y_log10()
    p
  })

  output$dosePlot <- renderPlot({
    doses <- c(10, 25)
    all_sims <- bind_rows(lapply(doses, function(d_mg) {
      mod %>%
        param(EGFR = input$egfr, SEX = as.numeric(input$sex)) %>%
        ev(amt = d_mg, cmt = 1, ii = 24, addl = 13) %>%
        mrgsim(end = 14 * 24, delta = 0.5) %>%
        as.data.frame() %>%
        mutate(Dose = paste0(d_mg, " mg"), time_days = time / 24)
    }))

    # Get steady-state metrics
    ss_metrics <- all_sims %>%
      dplyr::filter(time >= 13 * 24) %>%
      group_by(Dose) %>%
      summarise(Cmax = max(CP), Ctrough = min(CP[CP > 0]),
                AUC24 = sum(diff(time[order(time)]) *
                  (head(CP[order(time)], -1) + tail(CP[order(time)], -1)) / 2),
                .groups = "drop")

    d_plot <- if (input$log_scale) all_sims %>% dplyr::filter(CP > 0.001) else all_sims

    p <- ggplot(d_plot, aes(x = time_days, y = CP, color = Dose)) +
      geom_line(linewidth = 1) +
      scale_color_manual(values = c("10 mg" = "#8b5cf6", "25 mg" = "#10b981")) +
      labs(x = "Time (days)", y = "Concentration (mg/L)",
           title = "Empagliflozin Dose Comparison at Steady State",
           subtitle = sprintf("eGFR = %d mL/min/1.73m\u00B2", input$egfr)) +
      theme_minimal(base_size = 14) +
      theme(legend.position = "top")
    if (input$log_scale) p <- p + scale_y_log10()
    p
  })
}

shinyApp(ui = ui, server = server)