XCRs and tumour clones

library(ithi.utils)
load_base_libs()

library(ithi.meta)
library(ithi.xcr)
library(ithi.clones)

Colour palettes

pal_patient <- select_palette("patient")

Parameters

db_path <- snakemake@params$db

xcr_table_path <- snakemake@input$xcr_table
prevalence_option <- snakemake@params$prevalence_option
distance_method <- snakemake@params$xcr_distance_method

clone_tree_file <- snakemake@input$clone_tree_file
clone_prevalence_file <- snakemake@input$clone_prevalence_file
clone_branch_length_file <- snakemake@input$clone_branch_length_file

Metadata

db <- src_sqlite(db_path, create = FALSE)

duplicates <- collect(tbl(db, "duplicates"))

Analysis

XCR repertoire similarity

xcr_table <- read_clonotypes(xcr_table_path, duplicates = FALSE, db_path, verbose = 1)

Read 19.7% of 304822 rows
Read 59.1% of 304822 rows
Read 91.9% of 304822 rows
Read 304822 rows and 18 (of 18) columns from 0.070 GB file in 00:00:05

tcr_segment_type <- "TRB"
bcr_segment_type <- "IGH"

id_type <- "condensed_id"

patients <- unique(xcr_table$patient_id)

dists <- lapply(patients, function(patient) {
    tcr_clonotypes <- subset(xcr_table, type == tcr_segment_type & patient_id == 
        patient)
    bcr_clonotypes <- subset(xcr_table, type == bcr_segment_type & patient_id == 
        patient)
    tcr_cross_table <- cross_tabulate(tcr_clonotypes, id_type = id_type)
    bcr_cross_table <- cross_tabulate(bcr_clonotypes, id_type = id_type)
    
    cross_tables <- list(tcr = tcr_cross_table, bcr = bcr_cross_table)
    
    distance_matrices <- lapply(cross_tables, function(cross_table) {
        distmat <- compute_immune_distance_matrix(cross_table, method = distance_method)
        mat <- as.matrix(distmat)
        return(mat)
    })
    return(distance_matrices)
})
tcr_dists <- lapply(dists, function(x) x$tcr)
bcr_dists <- lapply(dists, function(x) x$bcr)

xcr_dists <- tibble(patient_id = patients, tcr = tcr_dists, bcr = bcr_dists)
xcr_dists$patient_id <- factor(xcr_dists$patient_id)
xcr_dists <- xcr_dists[order(xcr_dists$patient_id), ]

The distance metric being used is horn.

Tumour clones

tree_branch_data <- read_clone_tree_data(clone_tree_file, clone_branch_length_file, 
    clone_prevalence_file, db_path)
trees <- lapply(tree_branch_data, function(x) x$tree)
prevalences <- rbind.fill(lapply(tree_branch_data, function(x) x$prevalence_dat))
branch_lengths <- rbind.fill(lapply(tree_branch_data, function(x) x$branch_dat))

ccfs <- compute_ccf(prevalences, trees, id_type = id_type)
ccfs_labeled <- merge(ccfs, branch_lengths, by = c("label", "node", "patient_id"))

## Need non-normalized distances
clone_dists <- clone_distances(ccfs_labeled, normalize = FALSE, id_type = "condensed_id")

Combination

clone_dists$patient_id <- factor(clone_dists$patient_id)

total_dists <- inner_join(clone_dists, xcr_dists)

compute_overlap_similarities <- function(vals) {
    common_samples <- lapply(vals, function(x) {
        colnames(x)
    })
    common_samples <- Reduce(intersect, common_samples)
    
    sim_dfs <- lapply(names(vals), function(category) {
        mat <- vals[[category]]
        mat <- mat[common_samples, common_samples]
        
        res <- setNames(melt(as.matrix(mat)), c("sample1", "sample2", category))
        return(res)
    })
    similarities <- Reduce(merge, sim_dfs)
    result <- subset(similarities, as.numeric(sample1) < as.numeric(sample2))
    return(result)
}

pairwise_similarities <- rbind.fill(lapply(1:nrow(total_dists), function(i) {
    clone_distmat <- total_dists[i, ]$dist_clones_weighted[[1]]
    tcr_distmat <- total_dists[i, ]$tcr[[1]]
    bcr_distmat <- total_dists[i, ]$bcr[[1]]
    
    inputs <- list(clone_distmat, tcr_distmat, bcr_distmat)
    names(inputs) <- c("clones", "tcr", "bcr")
    
    patient_id <- total_dists[i, ]$patient_id
    
    pairwise_sims <- cbind(patient_id = patient_id, compute_overlap_similarities(inputs))
    return(pairwise_sims)
}))
pairwise_similarities$patient_id <- factor(pairwise_similarities$patient_id)

patient_counts <- pairwise_similarities %>% group_by(patient_id) %>% summarise(n = n())
singletons <- subset(patient_counts, n == 1)$patient_id
pairwise_similarities <- subset(pairwise_similarities, !patient_id %in% singletons)

TCR-clones

stats <- setNames(ddply(pairwise_similarities, .(patient_id), function(x) {
    res <- with(x, cor.test(clones, tcr, method = "pearson"))
    pvalue <- res$p.value
    eq <- substitute(italic(P) == p, list(p = format(pvalue, digits = 3)))
    return(as.character(as.expression(eq)))
}), c("patient_id", "p.value"))

summary(lmer(tcr ~ clones + (1 | patient_id), data = pairwise_similarities))

Linear mixed model fit by REML t-tests use Satterthwaite approximations
  to degrees of freedom [lmerMod]
Formula: tcr ~ clones + (1 | patient_id)
   Data: pairwise_similarities

REML criterion at convergence: -66.5

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-4.2035 -0.3406  0.1428  0.4982  1.7980 

Random effects:
 Groups     Name        Variance Std.Dev.
 patient_id (Intercept) 0.03909  0.1977  
 Residual               0.01950  0.1396  
Number of obs: 116, groups:  patient_id, 13

Fixed effects:
             Estimate Std. Error        df t value Pr(>|t|)    
(Intercept) 7.045e-01  5.957e-02 1.417e+01   11.83 9.86e-09 ***
clones      2.708e-05  7.941e-06 1.136e+02    3.41 0.000901 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
       (Intr)
clones -0.313
fit warnings:
Some predictor variables are on very different scales: consider rescaling

ggplot(pairwise_similarities, aes(x = clones, y = tcr)) + geom_point(aes(colour = patient_id)) + 
    theme_bw() + theme_Publication() + scale_color_manual(values = pal_patient) + 
    xlab("Clone distance") + ylab("TCR distance") + facet_wrap(~patient_id, 
    scales = "free") + geom_text(data = stats, aes(x = Inf, y = Inf, label = p.value), 
    hjust = 1.1, vjust = 1.5, size = 3, parse = TRUE)

BCR-clones

stats <- setNames(ddply(pairwise_similarities, .(patient_id), function(x) {
    res <- with(x, cor.test(clones, bcr, method = "pearson"))
    pvalue <- res$p.value
    eq <- substitute(italic(P) == p, list(p = format(pvalue, digits = 3)))
    return(as.character(as.expression(eq)))
}), c("patient_id", "p.value"))

summary(lmer(bcr ~ clones + (1 | patient_id), data = pairwise_similarities))

Linear mixed model fit by REML t-tests use Satterthwaite approximations
  to degrees of freedom [lmerMod]
Formula: bcr ~ clones + (1 | patient_id)
   Data: pairwise_similarities

REML criterion at convergence: -88.1

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-5.2307 -0.3481  0.1803  0.4483  1.6530 

Random effects:
 Groups     Name        Variance Std.Dev.
 patient_id (Intercept) 0.01768  0.1330  
 Residual               0.01720  0.1312  
Number of obs: 116, groups:  patient_id, 13

Fixed effects:
             Estimate Std. Error        df t value Pr(>|t|)    
(Intercept) 7.949e-01  4.254e-02 1.532e+01  18.688 5.86e-12 ***
clones      2.192e-05  7.205e-06 1.121e+02   3.043  0.00292 ** 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Correlation of Fixed Effects:
       (Intr)
clones -0.394
fit warnings:
Some predictor variables are on very different scales: consider rescaling

ggplot(pairwise_similarities, aes(x = clones, y = bcr)) + geom_point(aes(colour = patient_id)) + 
    theme_bw() + theme_Publication() + scale_color_manual(values = pal_patient) + 
    xlab("Clone distance") + ylab("BCR distance") + facet_wrap(~patient_id, 
    scales = "free") + geom_text(data = stats, aes(x = Inf, y = Inf, label = p.value), 
    hjust = 1.1, vjust = 1.5, size = 3, parse = TRUE)

TCR-clone correlations are good for most patients, while BCR-clone correlations are poor for most patients. Patient 1 has good correlations for both – not exactly surprising as patient 1 is very immunoreactive (highest out of the ITH2 cohort). Patient 4 is a clear opposite to this trend. Is it possible that patient 4 harbors a B-cell driven response, rather than a T-cell driven one?

Do we have corroborating evidence for this?

Supports:

• Patient 4 has the highest intrapatient BCR similarity across the entire cohort (immune responses highly similar from site-to-site could imply that there’s a widespread antigen, and tends to be associated with higher CD20+ TIL)
• Patient 4 has relatively high CD20+ TIL density (see the Immune Variability page) – tied for 3rd/21

Does not support:

• Patient 3 also has high CD20+ TIL density (tied for 3rd/21), but does not exhibit this pattern