pathways.Rmd

---
title: "Figure 1C: Pathway point mutation frequency comparison"
output: html_notebook
---

# Setup
```{r}
library(tidyverse)
library(janitor)
library(ggpubr)
library(rstatix)
```

```{r}
patient_sheet_paired <- read_delim('data/paired_participant_cohort_sheet.txt', show_col_types = FALSE) %>% 
  mutate(Cohort = 'Paired')

patient_sheet_unpaired <- read_delim('data/unpaired_participant_cohort_sheet.txt', show_col_types = FALSE)

patient_sheet <- bind_rows(patient_sheet_paired, patient_sheet_unpaired) %>% 
  mutate(harmonized_tki_duration = ifelse(str_equal(Cohort, 'Paired'), sampled_tki_duration, earliest_tki_duration),
         harmonized_os_tki_dfd = ifelse(str_equal(Cohort, 'Paired'), os_sampled_tki_dfd, os_earliest_tki_dfd))

maf_paired <- read_delim('data/paired.maf',
                         show_col_types = FALSE,
                         col_types = cols(level_of_evidence = col_character())) %>% 
  filter((!str_equal(trajectory, 'low_ccf')) & (!is.na(gene_set)) & (!is.na(level_of_evidence))) %>% 
  distinct(Patient_ID, unique_mut_id, .keep_all = TRUE)

maf_unpaired <- read_delim('data/unpaired.maf',
                           show_col_types = FALSE,
                           col_types = cols(level_of_evidence = col_character())) %>% 
  filter(postDP_ccf_mode >= 0.1 & (!is.na(gene_set)) & (!is.na(level_of_evidence))) %>% 
  distinct(Patient_ID, unique_mut_id, .keep_all = TRUE)

maf_combined <- bind_rows(maf_paired, maf_unpaired) %>% 
  mutate(gene_set = ifelse(str_equal(gene_set, 'ESR1'), 'ESR1 + regulators', gene_set)) %>% 
  left_join(patient_sheet %>% select(participant_id, Cohort, resistance_type), by = c('Patient_ID' = 'participant_id'))
```

# Calculate pathway mutation frequencies
```{r}
to_plot_muts <- maf_combined %>%
  mutate(Cohort = factor(case_when(
    str_equal(Cohort, 'Paired') ~ str_c('Paired ', resistance_type),
    str_equal(Cohort, 'Exposed') ~ 'Post-TKI',
    str_equal(Cohort, 'Non-exposed') ~ 'Pre-TKI'
  ), levels = c('Post-TKI', 'Pre-TKI', 'Paired Acquired', 'Paired Intrinsic'))) %>%
  group_by(gene_set, Cohort) %>% 
  summarize(n_pts_mut = n_distinct(Patient_ID), .groups = 'drop') %>%
  complete(gene_set, nesting(Cohort),
           fill = list(n_pts_mut = 0)) %>% 
  mutate(pct_pts_mut = 100 * case_when(
    str_equal(Cohort, 'Paired Acquired') ~ n_pts_mut / 8,
    str_equal(Cohort, 'Paired Intrinsic') ~ n_pts_mut / 18,
    str_equal(Cohort, 'Post-TKI') ~ n_pts_mut / 30,
    str_equal(Cohort, 'Pre-TKI') ~ n_pts_mut / 55
  )) %>% 
  group_by(gene_set) %>% 
  mutate(order_sum = sum(pct_pts_mut)) %>% 
  ungroup %>% 
  arrange(order_sum) %>% 
  mutate(gene_set = replace(gene_set, str_equal(gene_set, 'Cell Cycle'), 'Cell cycle')) %>% 
  mutate(label = paste0(round(pct_pts_mut, 1), '% (n = ', n_pts_mut, ')'))

to_plot_muts
```

# Plot pathway mutation frequencies
```{r}
plot_mut_freq <- ggbarplot(to_plot_muts, 'gene_set', 'pct_pts_mut',
          color = 'Cohort', fill = 'Cohort', label = to_plot_muts %>% pull(label),
          orientation = 'horiz',
          position = position_dodge(0.75), lab.vjust = 0.5, lab.hjust = -0.1,
          lab.col = 'black',
          ylab = 'Frequency of patients with mutation (%)',
          xlab = 'Pathway',
          lab.size = 3) +
  theme_classic(base_size = 16) +
  theme(legend.position = c(0.8, 0.22)) +
  scale_y_continuous(limits = c(0, 100)) +
  coord_flip(clip = 'off') +
  scale_color_manual(values = c('#0072B2','#F0E442','#009E73','#D55E00')) +
  scale_fill_manual(values = c('#0072B2','#F0E442','#009E73','#D55E00')) +
  guides(color = guide_legend(reverse = TRUE),
         fill = guide_legend(reverse = TRUE))

plot_mut_freq

ggsave(plot_mut_freq, filename = '../results/2024-07-02_figure-1/mut_freq_pathways.2024-07-02.pdf',
       width = 6, height = 8, units = 'in')
```

# Compare frequencies
(intrinsic vs. acquired and pre-TKI vs. post-TKI, with Fisher's Exact test)
```{r}
to_fisher <- to_plot_muts %>% 
  select(-pct_pts_mut, -order_sum) %>% 
  mutate(n_pts_not_mut = case_when(
    str_equal(Cohort, 'Paired Acquired') ~ 8 - n_pts_mut,
    str_equal(Cohort, 'Paired Intrinsic') ~ 18 - n_pts_mut,
    str_equal(Cohort, 'Post-TKI') ~ 30 - n_pts_mut,
    str_equal(Cohort, 'Pre-TKI') ~ 55 - n_pts_mut
  ))

fisher_results <- bind_rows(
  to_fisher %>% 
    filter(Cohort %in% c('Paired Intrinsic', 'Paired Acquired')) %>% 
    group_by(gene_set) %>% 
    summarize(data = list(row_wise_fisher_test(
      as.matrix(select(cur_data(), starts_with('n_'))), p.adjust.method = 'BH'))) %>%
    unnest_wider(data) %>% 
    unnest(c(group:p.adj.signif)) %>% 
    mutate(comparison = 'Intrinsic vs. Acquired') %>% 
    distinct(gene_set, p, comparison),
  to_fisher %>% 
    filter(Cohort %in% c('Post-TKI', 'Pre-TKI')) %>% 
    group_by(gene_set) %>% 
    summarize(data = list(row_wise_fisher_test(
      as.matrix(select(cur_data(), starts_with('n_'))), p.adjust.method = 'BH'))) %>%
    unnest_wider(data) %>% 
    unnest(c(group:p.adj.signif)) %>% 
    mutate(comparison = 'Pre-TKI vs. Post-TKI') %>% 
    distinct(gene_set, p, comparison)
)

fisher_results %>% arrange(p)
```

```{r}
sessionInfo()
```