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

model_code <- '
$PARAM @annotated
CLREF    : 0.243 : Typical clearance at reference covariates (L/day)
VCREF    : 4.86  : Typical central volume at reference covariates (L)
Q        : 0.656 : Intercompartmental clearance (L/day)
VPREF    : 4.25  : Typical peripheral volume at reference covariates (L)
KA       : 0.229 : First-order absorption rate constant (day^-1)
F_SC_P3  : 0.890 : Absolute SC bioavailability, phase III formulation
F_SC_P12 : 0.710 : Absolute SC bioavailability, phase I/II formulation
WT       : 90    : Body weight (kg)
ALB      : 44    : Serum albumin (g/L)
SCR      : 76    : Serum creatinine (umol/L)
HSCRP    : 2.8   : High-sensitivity CRP (mg/L)
ADA128   : 0     : ADA titer >=128 effect flag (0/1)
PHASE3   : 1     : Phase III SC formulation flag (1=yes, 0=no)

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

$MAIN
double wt_ratio = pow(WT / 70.0, 0.933);
double alb_ratio = pow(ALB / 44.0, -0.715);
double scr_ratio = pow(SCR / 76.0, -0.253);
double hscrp_ratio = pow(HSCRP / 2.8, 0.044);
double ada_ratio = ADA128 > 0.5 ? 1.428 : 1.0;

double CLi = CLREF * wt_ratio * alb_ratio * scr_ratio * hscrp_ratio * ada_ratio;
double V1i = VCREF * pow(WT / 70.0, 1.17);
double V2i = VPREF * pow(WT / 70.0, 0.377);
double FSC = PHASE3 > 0.5 ? F_SC_P3 : F_SC_P12;

F_DEPOT = FSC;

$ODE
dxdt_DEPOT = -KA * DEPOT;
dxdt_CENT = KA * DEPOT - (CLi / V1i) * CENT - (Q / V1i) * CENT + (Q / V2i) * PERIPH;
dxdt_PERIPH = (Q / V1i) * CENT - (Q / V2i) * PERIPH;

$TABLE
double wt_ratio_tab = pow(WT / 70.0, 0.933);
double alb_ratio_tab = pow(ALB / 44.0, -0.715);
double scr_ratio_tab = pow(SCR / 76.0, -0.253);
double hscrp_ratio_tab = pow(HSCRP / 2.8, 0.044);
double ada_ratio_tab = ADA128 > 0.5 ? 1.428 : 1.0;

double CL_OUT = CLREF * wt_ratio_tab * alb_ratio_tab * scr_ratio_tab * hscrp_ratio_tab * ada_ratio_tab;
double V1_OUT = VCREF * pow(WT / 70.0, 1.17);
double V2_OUT = VPREF * pow(WT / 70.0, 0.377);
double CP = CENT / V1_OUT;
double K10 = CL_OUT / V1_OUT;
double K12 = Q / V1_OUT;
double K21 = Q / V2_OUT;
double disc = pow(K10 + K12 + K21, 2.0) - 4.0 * K10 * K21;
if (disc < 0.0) disc = 0.0;
double beta = 0.5 * ((K10 + K12 + K21) - sqrt(disc));
double THALF = beta > 0.0 ? log(2.0) / beta : 0.0;
double F_OUT = PHASE3 > 0.5 ? F_SC_P3 : F_SC_P12;

$CAPTURE @annotated
CP     : Plasma concentration (mg/L)
CL_OUT : Individualized clearance (L/day)
V1_OUT : Individualized central volume (L)
V2_OUT : Individualized peripheral volume (L)
F_OUT  : SC bioavailability fraction
THALF  : Terminal half-life estimate (days)
'

mod <- mcode("risankizumab_pop_pk_339", model_code)

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 <- page_sidebar(
  title = "Risankizumab PopPK Simulator",
  theme = app_theme,
  sidebar = sidebar(
    title = "Simulation Settings",
    width = 340,
    accordion(
      accordion_panel(
        "Dosing",
        radioButtons(
          "route",
          "Administration route",
          choices = c("Subcutaneous" = "sc", "Intravenous" = "iv"),
          selected = "sc",
          inline = TRUE
        ),
        numericInput("dose_mg", "Dose (mg)", value = 150, min = 1, max = 1200, step = 1),
        selectInput(
          "regimen",
          "Regimen",
          choices = c(
            "Phase III psoriasis: weeks 0, 4, then q12w" = "phase3",
            "Fixed q12w" = "q12w",
            "Fixed q8w" = "q8w",
            "Fixed q4w" = "q4w"
          ),
          selected = "phase3"
        ),
        numericInput("n_doses", "Total doses", value = 6, min = 1, max = 12, step = 1),
        conditionalPanel(
          "input.route === 'sc'",
          selectInput(
            "sc_formulation",
            "SC formulation",
            choices = c(
              "Phase III clinical supply (F = 0.89)" = "phase3",
              "Phase I/II clinical supply (F = 0.71)" = "phase12"
            ),
            selected = "phase3"
          )
        )
      ),
      accordion_panel(
        "Covariates",
        sliderInput("wt", "Body weight (kg)", min = 40, max = 200, value = 90, step = 1),
        sliderInput("albumin", "Serum albumin (g/L)", min = 34, max = 58, value = 44, step = 0.5),
        sliderInput("creatinine", "Serum creatinine (umol/L)", min = 35.4, max = 203, value = 76, step = 0.1),
        sliderInput("hscrp", "hs-CRP (mg/L)", min = 0.1, max = 131, value = 2.8, step = 0.1),
        checkboxInput("ada128", "ADA titer >= 128", value = FALSE)
      ),
      accordion_panel(
        "Display",
        checkboxInput("log_scale", "Log scale concentration axis", value = FALSE)
      ),
      open = c("Dosing", "Covariates")
    )
  ),
  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, last interval (ug/mL)")
    ),
    div(
      class = "metric-card metric-warning",
      div(class = "metric-value", textOutput("ctrough")),
      div(class = "metric-label", "Ctrough, end interval (ug/mL)")
    ),
    div(
      class = "metric-card metric-primary",
      div(class = "metric-value", textOutput("auc")),
      div(class = "metric-label", "AUC0-tau, last interval (ug*day/mL)")
    ),
    div(
      class = "metric-card metric-info",
      div(class = "metric-value", textOutput("thalf")),
      div(class = "metric-label", "Terminal t1/2 (days)")
    )
  ),
  card(
    card_header("Serum Risankizumab Concentration-Time Profile"),
    full_screen = TRUE,
    plotOutput("pk_plot", height = "500px")
  ),
  navset_card_underline(
    title = "Drug Information",
    nav_panel(
      "Model Information",
      div(
        style = "padding: 10px 4px;",
        tags$h4("Final population PK model"),
        tags$p("Two-compartment model with first-order absorption and elimination fit to 13,123 plasma concentrations from 1,899 subjects across phase I to III healthy volunteer and plaque psoriasis studies."),
        tags$ul(
          tags$li(tags$strong("Typical CL:"), " 0.243 L/day at reference covariates"),
          tags$li(tags$strong("Typical Vc:"), " 4.86 L"),
          tags$li(tags$strong("Typical Q:"), " 0.656 L/day"),
          tags$li(tags$strong("Typical Vp:"), " 4.25 L"),
          tags$li(tags$strong("Ka:"), " 0.229 day^-1"),
          tags$li(tags$strong("SC bioavailability:"), " 0.89 in phase III formulation, 0.71 in phase I/II formulation")
        ),
        tags$h5("Published covariates"),
        tags$p("The manuscript retained body weight, serum albumin, serum creatinine, hs-CRP, and ADA titer >=128 on clearance; body weight also affected Vc and Vp. The reference values used here are the published medians: weight 70 kg for scaling, albumin 44 g/L, creatinine 76 umol/L, and hs-CRP 2.8 mg/L."),
        tags$p("For a typical 90 kg psoriasis subject, the paper reported CL ~0.31 L/day, Vss ~11.2 L, terminal half-life ~28 days, and phase III SC bioavailability 89%."),
        tags$h5("Clinical context"),
        tags$p("The approved psoriasis regimen is 150 mg SC at weeks 0 and 4, then every 12 weeks thereafter. Model-predicted mean Ctrough values reported in the paper were 2.48 ug/mL at week 16 and 1.91 ug/mL at week 52.")
      )
    ),
    nav_panel(
      "References",
      div(
        class = "ref-box",
        tags$h5("Primary source"),
        tags$ol(
          tags$li("Suleiman AA, Minocha M, Khatri A, Pang Y, Othman AA. Population Pharmacokinetics of Risankizumab in Healthy Volunteers and Subjects with Moderate to Severe Plaque Psoriasis: Integrated Analyses of Phase I-III Clinical Trials. Clin Pharmacokinet. 2019.")
        ),
        tags$h5("Drug summary"),
        tags$ul(
          tags$li(tags$strong("Brand:"), " Skyrizi"),
          tags$li(tags$strong("Class:"), " Anti-IL-23 p19 monoclonal antibody"),
          tags$li(tags$strong("Therapeutic area:"), " Immunology / dermatology"),
          tags$li(tags$strong("Indication in source paper:"), " Moderate to severe plaque psoriasis"),
          tags$li(tags$strong("Routes evaluated:"), " Intravenous and subcutaneous")
        )
      )
    )
  ),
  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"
    ),
    " | Built by Sunny",
    br(),
    tags$span(
      style = "font-size: 10px;",
      "For research and educational purposes only. Not for clinical decision-making."
    )
  )
)

server <- function(input, output, session) {
  dose_schedule <- reactive({
    n_doses <- as.integer(input$n_doses)
    shiny::req(!is.na(n_doses), n_doses >= 1)

    if (identical(input$regimen, "phase3")) {
      if (n_doses == 1) {
        times <- 0
      } else if (n_doses == 2) {
        times <- c(0, 28)
      } else {
        times <- c(0, 28, 28 + seq(84, by = 84, length.out = n_doses - 2))
      }
    } else {
      interval_days <- switch(
        input$regimen,
        q12w = 84,
        q8w = 56,
        q4w = 28,
        84
      )
      times <- seq(0, by = interval_days, length.out = n_doses)
    }

    route_cmt <- if (identical(input$route, "iv")) 2 else 1
    data.frame(
      ID = rep(1L, length(times)),
      time = times,
      amt = rep(as.numeric(input$dose_mg), length(times)),
      cmt = rep(route_cmt, length(times)),
      evid = rep(1, length(times))
    )
  })

  interval_bounds <- reactive({
    sched <- dose_schedule()
    shiny::req(nrow(sched) >= 1)

    dose_times <- sched$time
    if (length(dose_times) == 1) {
      start_time <- 0
      end_time <- max(168, dose_times[1] + 168)
    } else {
      start_time <- dose_times[length(dose_times) - 1]
      end_time <- dose_times[length(dose_times)]
    }

    list(
      start = start_time,
      end = end_time,
      sim_end = max(end_time + 168, max(dose_times) + 168)
    )
  })

  sim_data <- reactive({
    sched <- dose_schedule()
    bounds <- interval_bounds()

    out <- mod %>%
      param(
        WT = input$wt,
        ALB = input$albumin,
        SCR = input$creatinine,
        HSCRP = input$hscrp,
        ADA128 = ifelse(isTRUE(input$ada128), 1, 0),
        PHASE3 = ifelse(identical(input$sc_formulation, "phase3"), 1, 0)
      ) %>%
      data_set(sched) %>%
      mrgsim(end = bounds$sim_end, delta = 0.25) %>%
      as.data.frame() |>
      mutate(
        time_days = time,
        conc_ugml = pmax(CP, 0),
        route_label = ifelse(input$route == "iv", "IV", "SC")
      )

    out
  })

  last_interval <- reactive({
    d <- sim_data()
    bounds <- interval_bounds()
    d |>
      dplyr::filter(time_days >= bounds$start, time_days <= bounds$end)
  })

  output$cmax <- renderText({
    d <- last_interval()
    shiny::req(nrow(d) > 0)
    sprintf("%.2f", max(d$conc_ugml, na.rm = TRUE))
  })

  output$ctrough <- renderText({
    d <- last_interval()
    shiny::req(nrow(d) > 0)
    sprintf("%.2f", d$conc_ugml[nrow(d)])
  })

  output$auc <- renderText({
    d <- last_interval()
    shiny::req(nrow(d) > 1)
    auc_val <- sum(diff(d$time_days) * (head(d$conc_ugml, -1) + tail(d$conc_ugml, -1)) / 2)
    sprintf("%.1f", auc_val)
  })

  output$thalf <- renderText({
    d <- sim_data()
    shiny::req(nrow(d) > 0)
    sprintf("%.1f", d$THALF[1])
  })

  output$pk_plot <- renderPlot({
    d <- sim_data()
    shiny::req(nrow(d) > 0)

    plot_data <- if (isTRUE(input$log_scale)) {
      d |>
        dplyr::filter(conc_ugml > 0)
    } else {
      d
    }

    p <- ggplot(plot_data, aes(x = time_days, y = conc_ugml)) +
      geom_line(color = "#8b5cf6", linewidth = 0.9) +
      labs(
        x = "Time (days)",
        y = "Risankizumab concentration (ug/mL)",
        title = paste0("Risankizumab ", input$dose_mg, " mg ", toupper(input$route))
      ) +
      theme_minimal(base_size = 14) +
      theme(
        panel.grid.minor = element_blank(),
        plot.title = element_text(face = "bold")
      )

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

    p
  })
}

shinyApp(ui = ui, server = server)