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

# ── VRC07-523LS PopPK Model ──
# Two-compartment model with zero-order SC absorption
# Source: Huynh et al. J Antimicrob Chemother 2026; dkaf449
# Population: HIV-exposed infants + healthy adults
# NONMEM ADVAN4 TRANS4 parameterization

code <- '
$PARAM @annotated
WT    : 3.0  : Body weight (kg)
DOSE  : 80   : Dose (mg)
D1HR  : 36   : Zero-order absorption duration (hours) - infant
TVCL  : 47.76 : Typical CL (mL/day/70 kg) - infant base
TVVC  : 1.48  : Typical Vc (L/70 kg)
TVVP  : 2.28  : Typical Vp (L/70 kg)
TVQ   : 0.5832 : Typical Q (L/day/70 kg) [0.0243 L/h * 24]
TVF   : 0.30  : Bioavailability (SC)
ISADULT : 0   : Adult flag (0=infant, 1=adult)
ISREPEAT : 0  : Repeat dose flag (0=first, 1=subsequent)
ADULCL : 1.69 : Adult factor on CL
ADULTD : 2.79 : Adult factor on D1
REPVSS : 1.49 : Repeat dose factor on Vss

$CMT @annotated
DEPOT  : SC depot (zero-order input)
CENT   : Central compartment
PERIPH : Peripheral compartment

$MAIN
double CLi = TVCL * pow(WT/70.0, 0.85);
if(ISADULT > 0.5) CLi = CLi * ADULCL;
double CLF = CLi / TVF;  // CL/F in mL/day

double VCi = TVVC * pow(WT/70.0, 1.0);
double VPi = TVVP * pow(WT/70.0, 1.0);
if(ISREPEAT > 0.5) {
  VCi = VCi * REPVSS;
  VPi = VPi * REPVSS;
}
double VCF = VCi / TVF;  // Vc/F in L
double VPF = VPi / TVF;  // Vp/F in L

double Qi = TVQ * pow(WT/70.0, 0.85);
double QF = Qi / TVF;    // Q/F in L/day

double D1i = D1HR / 24.0;  // Convert to days
if(ISADULT > 0.5) D1i = D1i * ADULTD;

D_CENT = D1i;  // Zero-order duration in days
F_CENT = 1.0;  // F already accounted for in CL/F, V/F

$ODE
double ke = CLF / 1000.0 / VCF;   // CLF mL/day -> L/day
double k12 = QF / VCF;             // QF already in L/day
double k21 = QF / VPF;

dxdt_DEPOT = 0;
dxdt_CENT = -ke * CENT - k12 * CENT + k21 * PERIPH;
dxdt_PERIPH = k12 * CENT - k21 * PERIPH;

$TABLE
double CP = CENT / VCF;  // mg / L = mcg/mL

$CAPTURE @annotated
CP : Serum concentration (mcg/mL)
'

mod <- mcode("VRC07_523LS", code, compile = TRUE)

# ── WHO Growth Chart (boys, 50th percentile, birth to 6 months) ──
# Reference: CDC/WHO growth standards (ref 16 in paper)
who_days <- c(0, 7, 14, 28, 42, 56, 70, 84, 98, 112, 126, 140, 154, 168)
who_wts_50th <- c(3.3, 3.6, 3.9, 4.5, 5.1, 5.6, 6.0, 6.4, 6.7, 7.0, 7.3, 7.5, 7.7, 7.9)

# ── Helper: Simulate ──
run_sim <- function(mod, population, dose1, dose2, dose2_week,
                    weight, second_dose_enabled) {

  is_adult <- ifelse(population == "Adult", 1, 0)

  if (population == "Infant") {
    # Scale WHO growth chart to actual birth weight
    scale_factor <- weight / 3.3
    who_wts <- who_wts_50th * scale_factor

    # Interpolation function for weight at any time
    wt_at <- function(t) approx(who_days, who_wts, xout = t, rule = 2)$y

    dose2_time <- if (second_dose_enabled) dose2_week * 7 else 999

    # Build dataset with weight updates every week + dose events
    wt_times <- seq(7, 168, by = 7)
    dat <- data.frame(
      ID = 1, time = wt_times, amt = 0, cmt = 0, evid = 0, rate = 0,
      WT = sapply(wt_times, wt_at),
      ISREPEAT = ifelse(wt_times >= dose2_time, 1, 0),
      ISADULT = 0
    )
    # Dose 1
    dat <- rbind(
      data.frame(ID = 1, time = 0, amt = dose1, cmt = 2, evid = 1, rate = -2,
                 WT = weight, ISREPEAT = 0, ISADULT = 0),
      dat
    )
    # Dose 2
    if (second_dose_enabled) {
      dat <- rbind(dat,
        data.frame(ID = 1, time = dose2_time, amt = dose2, cmt = 2, evid = 1, rate = -2,
                   WT = wt_at(dose2_time), ISREPEAT = 1, ISADULT = 0)
      )
    }
    dat <- dat %>% arrange(time, desc(evid))

    out <- mod %>%
      data_set(dat) %>%
      mrgsim(end = 168, delta = 0.25) %>%
      as.data.frame()
  } else {
    # Adults — constant weight, single simulation
    events <- data.frame(
      ID = 1, time = 0, amt = dose1, cmt = 2, evid = 1, rate = -2
    )
    if (second_dose_enabled) {
      dose2_time <- dose2_week * 7
      events <- rbind(events, data.frame(
        ID = 1, time = dose2_time, amt = dose2, cmt = 2, evid = 1, rate = -2
      ))
    }
    out <- mod %>%
      param(WT = weight, ISADULT = is_adult, ISREPEAT = 0) %>%
      ev(as.ev(events)) %>%
      mrgsim(end = 168, delta = 0.25) %>%
      as.data.frame()
  }

  out$week <- out$time / 7
  out
}

# ── UI ──
ui <- page_sidebar(
  title = "VRC07-523LS PK Simulator",
  theme = bs_theme(
    version = 5,
    bootswatch = "flatly",
    primary = "#8b5cf6"
  ),

  sidebar = sidebar(
    width = 320,
    title = "Parameters",

    selectInput("population", "Population",
                choices = c("Infant", "Adult"), selected = "Infant"),

    conditionalPanel(
      condition = "input.population == 'Infant'",
      sliderInput("birth_weight", "Birth Weight (kg)",
                  min = 2.0, max = 5.0, value = 2.8, step = 0.1)
    ),
    conditionalPanel(
      condition = "input.population == 'Adult'",
      sliderInput("adult_weight", "Body Weight (kg)",
                  min = 40, max = 120, value = 70, step = 1)
    ),

    hr(),
    h6("Dose 1"),
    sliderInput("dose1", "Dose 1 (mg)",
                min = 20, max = 200, value = 80, step = 10),

    hr(),
    checkboxInput("second_dose", "Give Second Dose", value = TRUE),
    conditionalPanel(
      condition = "input.second_dose",
      sliderInput("dose2", "Dose 2 (mg)",
                  min = 20, max = 200, value = 100, step = 10),
      sliderInput("dose2_week", "Second Dose Timing (weeks)",
                  min = 4, max = 20, value = 12, step = 1)
    ),

    hr(),
    h6("Display Options"),
    checkboxInput("show_target", "Show 10 µg/mL Target", value = TRUE),
    checkboxInput("show_log", "Log Scale", value = FALSE),

    hr(),
    tags$small(
      tags$b("Route:"), "Subcutaneous (SC)", tags$br(),
      tags$b("Model:"), "2-CMT + zero-order absorption", tags$br(),
      tags$b("Software:"), "NONMEM → mrgsolve"
    )
  ),

  # ── Main Panel ──
  layout_column_wrap(
    width = 1/3,
    fill = FALSE,
    value_box(
      title = "Cmax",
      value = textOutput("cmax_val"),
      showcase = icon("arrow-up"),
      theme = "primary"
    ),
    value_box(
      title = "C(12 wk)",
      value = textOutput("c12wk_val"),
      showcase = icon("clock"),
      theme = "info"
    ),
    value_box(
      title = "t½",
      value = textOutput("thalf_val"),
      showcase = icon("hourglass-half"),
      theme = "success"
    )
  ),

  navset_card_tab(
    title = "VRC07-523LS Pharmacokinetics",
    nav_panel(
      "Simulation",
      plotOutput("pk_plot", height = "500px")
    ),
    nav_panel(
      "PT80 Analysis",
      plotOutput("pt80_plot", height = "500px"),
      tags$p(class = "text-muted mt-2",
        "PT80 = Serum [VRC07-523LS] / IC80. ",
        "Median IC80 = 0.238 µg/mL (IQR: 0.103–0.545). ",
        "PT80 > 200 correlates with ~90% prevention efficacy; PT80 > 32 with ~50% efficacy (Gilbert et al., Nat Med 2022)."
      )
    ),
    nav_panel(
      "Model Information",
      tags$div(
        class = "p-3",
        tags$h5("Population Pharmacokinetic Model"),
        tags$p("Two-compartment model with zero-order subcutaneous absorption (NONMEM ADVAN4 TRANS4)."),
        tags$h6("Structural Parameters (Table 2)"),
        tags$table(
          class = "table table-sm table-bordered",
          tags$thead(tags$tr(
            tags$th("Parameter"), tags$th("Estimate"), tags$th("Description")
          )),
          tags$tbody(
            tags$tr(tags$td("Vc"), tags$td("1.48 L/70 kg"), tags$td("Central volume")),
            tags$tr(tags$td("Vp"), tags$td("2.28 L/70 kg"), tags$td("Peripheral volume")),
            tags$tr(tags$td("CL"), tags$td("47.76 mL/day/70 kg"), tags$td("Clearance (infant base)")),
            tags$tr(tags$td("Q"), tags$td("0.0243 L/h/70 kg"), tags$td("Intercompartmental clearance")),
            tags$tr(tags$td("F"), tags$td("0.30"), tags$td("SC bioavailability")),
            tags$tr(tags$td("D1"), tags$td("36 hours"), tags$td("Zero-order absorption duration (infant)")),
            tags$tr(tags$td("Adult CL factor"), tags$td("1.69"), tags$td("CL 69% higher in adults")),
            tags$tr(tags$td("Adult D1 factor"), tags$td("2.79"), tags$td("Absorption 2.79× slower in adults")),
            tags$tr(tags$td("Repeat dose Vss"), tags$td("1.49"), tags$td("Vss 49% higher with repeat dosing"))
          )
        ),
        tags$h6("Allometric Scaling"),
        tags$ul(
          tags$li("CL, Q ~ (WT/70)^0.85"),
          tags$li("Vc, Vp ~ (WT/70)^1.0")
        ),
        tags$h6("Between-Subject Variability"),
        tags$ul(
          tags$li("Vss: 35.7% CV"),
          tags$li("CL: 29.3% CV"),
          tags$li("D1: 10.0% CV")
        ),
        tags$h6("Apparent Parameters (Infant, 70 kg normalized)"),
        tags$ul(
          tags$li("CL/F = 159 mL/day/70 kg"),
          tags$li("Vss/F = 12.53 L/70 kg"),
          tags$li("t½ ≈ 39 days")
        ),
        tags$h6("Apparent Parameters (Adult, 70 kg normalized)"),
        tags$ul(
          tags$li("CL/F = 269 mL/day/70 kg"),
          tags$li("t½ ≈ 31 days")
        )
      )
    ),
    nav_panel(
      "References",
      tags$div(
        class = "p-3",
        tags$h5("Primary Source"),
        tags$p(
          "Huynh D, Nikanjam M, Cunningham CK, et al. ",
          tags$em("Model-based assessment of VRC07-523LS dosing in infants through population pharmacokinetic–pharmacodynamic modelling in adults and infants."),
          " J Antimicrob Chemother. 2026;81(2):dkaf449. ",
          tags$a(href = "https://doi.org/10.1093/jac/dkaf449", target = "_blank", "doi:10.1093/jac/dkaf449")
        ),
        tags$h5("Key Supporting References"),
        tags$ul(
          tags$li("Gaudinski MR, et al. Safety and pharmacokinetics of VRC07-523LS in healthy adults. Lancet HIV. 2019;6:e667–79."),
          tags$li("Gilbert PB, et al. Neutralization titer biomarker for antibody-mediated prevention of HIV-1 acquisition. Nat Med. 2022;28:1924–32."),
          tags$li("Cunningham CK, et al. Safety, tolerability, and PK of VRC07-523LS in newborn infants. J Pediatric Infect Dis Soc. 2025;14:piaf002."),
          tags$li("Chawana TD, et al. PK interaction assessment of VRC07-523LS: cross-protocol analysis. BMC Immunol. 2025;26:8."),
          tags$li("Walsh SR, et al. Safety and PK of VRC07-523LS via different routes (HVTN 127/HPTN 087). PLoS Med. 2024;21:e1004329.")
        ),
        tags$h5("Clinical Trials"),
        tags$ul(
          tags$li("IMPAACT P1112 (NCT02256631) — Infants exposed to HIV"),
          tags$li("VRC605 (NCT03015181) — Healthy adults, dose-escalation")
        ),
        tags$h5("About VRC07-523LS"),
        tags$p(
          "VRC07-523LS is a broadly neutralizing antibody (bNAb) targeting the CD4 binding site of HIV-1 gp120. ",
          "The LS mutation (M428L/N434S) in the Fc region enhances FcRn binding, reducing clearance and prolonging ",
          "the serum half-life compared to VRC01. It neutralizes >95% of Clade B and C HIV-1 isolates. ",
          "Currently under development for HIV prophylaxis in infants during breastfeeding and in adults."
        )
      )
    )
  ),

  # ── 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({
    wt <- ifelse(input$population == "Infant", input$birth_weight, input$adult_weight)
    run_sim(
      mod,
      population = input$population,
      dose1 = input$dose1,
      dose2 = ifelse(input$second_dose, input$dose2, 0),
      dose2_week = input$dose2_week,
      weight = wt,
      second_dose_enabled = input$second_dose
    )
  })

  output$cmax_val <- renderText({
    d <- sim_data()
    cmax <- max(d$CP, na.rm = TRUE)
    paste0(round(cmax, 1), " µg/mL")
  })

  output$c12wk_val <- renderText({
    d <- sim_data()
    idx <- which.min(abs(d$time - 84))  # 12 weeks = 84 days
    c12 <- d$CP[idx]
    paste0(round(c12, 1), " µg/mL")
  })

  output$thalf_val <- renderText({
    d <- sim_data()
    cmax <- max(d$CP, na.rm = TRUE)
    cmax_time <- d$time[which.max(d$CP)]
    # Find time to half of Cmax after Cmax (terminal phase)
    post_cmax <- d %>% filter(time > cmax_time)
    half_target <- cmax / 2
    idx <- which.min(abs(post_cmax$CP - half_target))
    thalf <- post_cmax$time[idx] - cmax_time
    paste0(round(thalf, 0), " days")
  })

  output$pk_plot <- renderPlot({
    d <- sim_data()
    p <- ggplot(d, aes(x = week, y = CP)) +
      geom_line(color = "#8b5cf6", linewidth = 1.2) +
      labs(
        x = "Time (weeks)",
        y = "VRC07-523LS Concentration (µg/mL)",
        title = paste0("VRC07-523LS PK Profile — ",
                       input$population,
                       ifelse(input$population == "Infant",
                              paste0(" (", input$birth_weight, " kg)"),
                              paste0(" (", input$adult_weight, " kg)")))
      ) +
      theme_minimal(base_size = 14) +
      theme(
        plot.title = element_text(face = "bold"),
        panel.grid.minor = element_blank()
      )

    if (input$show_target) {
      p <- p +
        geom_hline(yintercept = 10, linetype = "dashed", color = "#e74c3c", linewidth = 0.8) +
        annotate("text", x = max(d$week) * 0.95, y = 12,
                 label = "Target: 10 µg/mL", color = "#e74c3c", hjust = 1, size = 3.5)
    }

    if (input$second_dose) {
      p <- p +
        geom_vline(xintercept = input$dose2_week, linetype = "dotted",
                   color = "#3498db", linewidth = 0.6) +
        annotate("text", x = input$dose2_week + 0.3, y = max(d$CP) * 0.9,
                 label = paste0("Dose 2\n(", input$dose2, " mg)"),
                 color = "#3498db", hjust = 0, size = 3.5)
    }

    # Annotate Dose 1
    p <- p +
      annotate("text", x = 0.3, y = max(d$CP) * 0.95,
               label = paste0("Dose 1: ", input$dose1, " mg SC"),
               color = "#8b5cf6", hjust = 0, size = 3.5)

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

    p
  })

  output$pt80_plot <- renderPlot({
    d <- sim_data()
    # Calculate PT80 using median IC80 = 0.238 µg/mL
    ic80_median <- 0.238
    ic80_low <- 0.103  # IQR lower
    ic80_high <- 0.545  # IQR upper

    d$pt80_median <- d$CP / ic80_median
    d$pt80_low <- d$CP / ic80_high   # high IC80 = low PT80
    d$pt80_high <- d$CP / ic80_low   # low IC80 = high PT80

    ggplot(d, aes(x = week)) +
      geom_ribbon(aes(ymin = pt80_low, ymax = pt80_high),
                  fill = "#8b5cf6", alpha = 0.2) +
      geom_line(aes(y = pt80_median), color = "#8b5cf6", linewidth = 1.2) +
      geom_hline(yintercept = 200, linetype = "dashed", color = "#e74c3c", linewidth = 0.8) +
      geom_hline(yintercept = 32, linetype = "dashed", color = "#f39c12", linewidth = 0.8) +
      annotate("text", x = max(d$week) * 0.95, y = 220,
               label = "PT80 = 200 (~90% efficacy)", color = "#e74c3c", hjust = 1, size = 3.5) +
      annotate("text", x = max(d$week) * 0.95, y = 40,
               label = "PT80 = 32 (~50% efficacy)", color = "#f39c12", hjust = 1, size = 3.5) +
      scale_y_log10() +
      labs(
        x = "Time (weeks)",
        y = "PT80 (Concentration / IC80)",
        title = "Predicted PT80 — Serum Neutralization Titre Ratio",
        subtitle = "Shaded: IQR of IC80 distribution (0.103–0.545 µg/mL)"
      ) +
      theme_minimal(base_size = 14) +
      theme(
        plot.title = element_text(face = "bold"),
        panel.grid.minor = element_blank()
      )
  })
}

shinyApp(ui, server)