Objective 1 - Step 7: Preliminary Clinical Analysis (Figure 1)

Overview

This script builds Figure 1 (“Preliminary clinical and molecular characteristics of the PDAC cohort”), using the demographic, clinical, and pathology variables available in this dataset:

Panel A - boxplots comparing PDAC and BBP on the four laboratory variables that most strongly differentiated the two groups in Table 1 (total bilirubin, ALP, albumin, CRP).
Panel B - sex distribution, BBP vs. PDAC, reflecting the significant demographic difference between the groups identified in Table 1.
Panel C - PDAC disease category distribution (RPC/LAPC/MPC), the resectability-based staging established in Table 2.
Panel D - a preliminary Kaplan-Meier survival comparison across disease categories, restricted to the subset of patients with complete follow-up data.

clinical_clean <- readRDS("clinical_clean.rds")
nms <- names(clinical_clean)

cat("CA19-9 column present:", any(grepl("ca19|ca_19", nms, ignore.case = TRUE)), "\n")

## CA19-9 column present: FALSE

cat("Stage/grade/resection/margin/TNM column present:",
    any(grepl("stage|grade|resect|margin|tnm", nms, ignore.case = TRUE)), "\n")

## Stage/grade/resection/margin/TNM column present: FALSE

library(dplyr)
library(ggplot2)
library(patchwork)
library(survival)

theme_set(theme_minimal(base_size = 10))

cc1 <- clinical_clean %>%
  mutate(group2 = case_when(
    pathology == "BBP" ~ "BBP",
    pathology %in% c("RPC", "LAPC", "MPC") ~ "PDAC",
    TRUE ~ NA_character_
  )) %>%
  filter(!is.na(group2))

table(cc1$group2)

## 
##  BBP PDAC 
##   83  102

Panel A: PDAC vs. BBP — key differentiating clinical variables

Bilirubin and ALP are shown on a log10 scale given their wide range; CRP has a minimum of 0 in this dataset, so it is left on a linear scale rather than log-transformed (log of zero is undefined).

make_box <- function(varlab, varname, logscale = FALSE) {
  d <- cc1 %>% filter(!is.na(.data[[varname]]))
  p <- ggplot(d, aes(x = group2, y = .data[[varname]], fill = group2)) +
    geom_boxplot(outlier.size = 0.8, width = 0.6) +
    scale_fill_manual(values = c(BBP = "#4C72B0", PDAC = "#DD8452"), guide = "none") +
    labs(x = NULL, y = NULL, title = varlab) +
    theme(plot.title = element_text(size = 8, face = "bold"))
  if (logscale) p <- p + scale_y_log10()
  p
}

pA1 <- make_box("Total bilirubin (umol/L)", "total_bilirubin_tbil_umol_l", logscale = TRUE)
pA2 <- make_box("ALP (U/L)", "alkaline_phosphatase_alp_u_l", logscale = TRUE)
pA3 <- make_box("Albumin (g/L)", "albumin_g_l")
pA4 <- make_box("CRP (mg/L)", "c_reactive_protein_mg_l")

rowA <- (pA1 | pA2 | pA3 | pA4)
rowA

Panel B: Sex distribution (BBP vs. PDAC)

sex_d <- cc1 %>%
  filter(!is.na(gender)) %>%
  count(group2, gender) %>%
  group_by(group2) %>%
  mutate(pct = 100 * n / sum(n))

pB <- ggplot(sex_d, aes(x = group2, y = pct, fill = gender)) +
  geom_col(width = 0.6) +
  scale_fill_manual(values = c(female = "#C44E52", male = "#55A868")) +
  labs(x = NULL, y = "Percent", title = "Sex distribution (BBP vs PDAC)", fill = NULL) +
  theme(plot.title = element_text(size = 8, face = "bold"), legend.position = "bottom")

pB

Panel C: PDAC disease category distribution

cat_d <- cc1 %>%
  filter(group2 == "PDAC") %>%
  count(pathology) %>%
  mutate(pathology = factor(pathology, levels = c("RPC", "LAPC", "MPC")))

pC <- ggplot(cat_d, aes(x = pathology, y = n, fill = pathology)) +
  geom_col(width = 0.6) +
  geom_text(aes(label = n), vjust = -0.4, size = 3) +
  scale_fill_manual(values = c(RPC = "#4C72B0", LAPC = "#DD8452", MPC = "#C44E52"), guide = "none") +
  labs(x = NULL, y = "n patients", title = "PDAC disease category distribution") +
  theme(plot.title = element_text(size = 8, face = "bold")) +
  ylim(0, max(cat_d$n) * 1.15)

pC

Panel D: Preliminary survival by disease category

Caveat, stated plainly before the result: only 33 of 102 PDAC patients (21 RPC, 8 LAPC, 4 MPC) have both a recorded follow-up time and vital status. This is a small, non-random subset (patients with complete follow-up capture are unlikely to be representative of the full PDAC cohort), and the LAPC/MPC arms in particular are too small for a reliable survival estimate. This panel is preliminary and hypothesis-generating only, not a definitive survival analysis.

dead_alive is coded 0 = Dead, 1 = Alive in the source data (verified against date_of_death in 00_Clinical_Data_Import_QC.Rmd); this is converted to the standard event = 1 (death observed) / event = 0 (censored) coding expected by the survival package.

pdac_surv <- clinical_clean %>%
  filter(pathology %in% c("RPC", "LAPC", "MPC")) %>%
  filter(!is.na(followup), !is.na(dead_alive)) %>%
  mutate(
    event = ifelse(dead_alive == 0, 1, 0),
    pathology = factor(pathology, levels = c("RPC", "LAPC", "MPC"))
  )

table(pdac_surv$pathology)

## 
##  RPC LAPC  MPC 
##   21    8    4

fit <- survfit(Surv(followup, event) ~ pathology, data = pdac_surv)
fit

## Call: survfit(formula = Surv(followup, event) ~ pathology, data = pdac_surv)
## 
##                 n events median 0.95LCL 0.95UCL
## pathology=RPC  21     12    356     126      NA
## pathology=LAPC  8      5    159      85      NA
## pathology=MPC   4      4    100      14      NA

km_df <- data.frame(time = fit$time, surv = fit$surv,
                     strata = rep(names(fit$strata), fit$strata))
km_df$strata <- gsub("pathology=", "", km_df$strata)

logrank <- survdiff(Surv(followup, event) ~ pathology, data = pdac_surv)
p_logrank <- 1 - pchisq(logrank$chisq, length(logrank$n) - 1)
cat("Log-rank p-value:", round(p_logrank, 4), "\n")

## Log-rank p-value: 0.1495

pD <- ggplot(km_df, aes(x = time, y = surv, color = strata)) +
  geom_step(linewidth = 0.8) +
  scale_color_manual(values = c(RPC = "#4C72B0", LAPC = "#DD8452", MPC = "#C44E52")) +
  labs(x = "Follow-up (days)", y = "Survival probability",
       title = "Preliminary survival by disease category",
       subtitle = paste0("n=33 complete cases; log-rank p=", round(p_logrank, 3)),
       color = NULL) +
  theme(plot.title = element_text(size = 8, face = "bold"),
        plot.subtitle = element_text(size = 6.5, color = "grey40"),
        legend.position = "bottom") +
  ylim(0, 1)

pD

The log-rank test does not reach significance (p = 0.15), consistent with the very small LAPC and MPC sample sizes in the complete-case subset; the curves separate in the expected direction (MPC declining fastest, RPC most favourable), but this should be read as directionally suggestive only, not confirmatory, given the sample size.

Combined Figure 1

rowBCD <- (pB | pC | pD)

figure1 <- rowA / rowBCD +
  plot_annotation(
    tag_levels = list(c("A", "", "", "", "B", "C", "D")),
    title = "Figure 1. Preliminary clinical and molecular characteristics of the PDAC cohort"
  ) &
  theme(plot.tag = element_text(face = "bold"))

figure1

ggsave("Figure1_Preliminary_Clinical.png", figure1, width = 12, height = 7, dpi = 200)

Results text (draft)

PDAC and BBP were first compared with respect to demographic and clinical variables (Table 1; Figure 1A-B). PDAC patients were significantly older and more frequently male than BBP patients, and showed significantly elevated total bilirubin, ALP, and CRP alongside significantly reduced albumin, consistent with biliary obstruction, hepatobiliary compromise, and a systemic inflammatory response, respectively (Figure 1A). Resectability-based disease category (RPC/LAPC/MPC) is reported in Table 2 and Figure 1C, with the cohort predominantly resectable at presentation (79/102, 77.5%).

A preliminary, hypothesis-generating survival comparison across disease categories (Figure 1D) did not reach statistical significance (log-rank p = 0.15), reflecting the small subset of patients (n = 33/102) with complete follow-up data rather than necessarily the absence of a true survival difference; this is revisited with the full available survival data in Objective 2.

Session info

sessionInfo()

## R version 4.6.1 (2026-06-24 ucrt)
## Platform: x86_64-w64-mingw32/x64
## Running under: Windows 11 x64 (build 26200)
## 
## Matrix products: default
##   LAPACK version 3.12.1
## 
## locale:
## [1] LC_COLLATE=English_United States.utf8 
## [2] LC_CTYPE=English_United States.utf8   
## [3] LC_MONETARY=English_United States.utf8
## [4] LC_NUMERIC=C                          
## [5] LC_TIME=English_United States.utf8    
## 
## time zone: Africa/Johannesburg
## tzcode source: internal
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
## [1] survival_3.8-6  patchwork_1.3.2 ggplot2_4.0.3   dplyr_1.2.1    
## 
## loaded via a namespace (and not attached):
##  [1] Matrix_1.7-5       gtable_0.3.6       jsonlite_2.0.0     compiler_4.6.1    
##  [5] tidyselect_1.2.1   jquerylib_0.1.4    textshaping_1.0.5  systemfonts_1.3.2 
##  [9] splines_4.6.1      scales_1.4.0       yaml_2.3.12        fastmap_1.2.0     
## [13] lattice_0.22-9     R6_2.6.1           labeling_0.4.3     generics_0.1.4    
## [17] knitr_1.51         tibble_3.3.1       bslib_0.11.0       pillar_1.11.1     
## [21] RColorBrewer_1.1-3 rlang_1.2.0        cachem_1.1.0       xfun_0.59         
## [25] sass_0.4.10        S7_0.2.2           otel_0.2.0         cli_3.6.6         
## [29] withr_3.0.3        magrittr_2.0.5     digest_0.6.39      grid_4.6.1        
## [33] rstudioapi_0.19.0  lifecycle_1.0.5    vctrs_0.7.3        evaluate_1.0.5    
## [37] glue_1.8.1         farver_2.1.2       ragg_1.5.2         rmarkdown_2.31    
## [41] tools_4.6.1        pkgconfig_2.0.3    htmltools_0.5.9