Datasets
Models
Downloads: 1
2140 lines (1836 with data), 90.1 kB
---
title: "Supplementary Figures"
author: Tobias Roider
date: "Last compiled on `r format(Sys.time(), '%d %B, %Y, %X')`"
output:
rmdformats::readthedown:
editor_options:
chunk_output_type: console
---
```{r options, include=FALSE, warning = FALSE}
library(knitr)
opts_chunk$set(echo=TRUE, tidy=FALSE, include=TRUE, message=FALSE, cache.lazy = FALSE,
dpi = 100, cache = FALSE, warning = FALSE)
opts_knit$set(root.dir = "../")
options(bitmapType = "cairo")
```
# Read data, functions and packages
```{r read data}
source("R/ReadPackages.R")
source("R/Functions.R")
source("R/ReadData.R")
source("R/ThemesColors.R")
source("R/Helpers.R")
```
# Supplementary Figure 1
## Patient characteristics
```{r patient characteristics SF1 part 1}
p1 <-
df_meta %>%
ggplot(aes(x=Order, y=Tcells, fill=Entity))+
geom_bar(stat = "identity", color="white", width=0.5, size=0.25)+
geom_text(aes(x=6, y=84, label=paste0("n = ", nrow(df_meta))), check_overlap = T, size=2.75)+
scale_y_continuous(name="% T-cells", limits=c(0, 90), expand = c(0.025, 0.025))+
scale_x_discrete(limits=as.character(1:101))+
scale_fill_manual(values = colors_characteristics)+
mytheme_1+
coord_cartesian(clip = "off")+
theme(axis.text.x = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_text(size=7),
axis.ticks.x = element_blank(),
plot.tag = element_text(margin = unit(c(0,0,0,-0.75), "cm")),
plot.margin = unit(c(0.1,0.1,0,1), "lines"))
cex_ <- 0.6
pos_ <- -0.25
size_ <- 2.4
p_ann <-
ggplot()+
geom_tile(data=df_meta, aes(y=0, x=Order, fill=Entity))+
geom_tile(data=df_meta, aes(y=-1.75, x=Order, fill=Status))+
geom_tile(data=df_meta, aes(y=-3.5, x=Order, fill=Sex))+
geom_text(data=df_meta, aes(y=-5.25, x=Order, label=`CITEseq`), size=size_)+
geom_text(data=df_meta, aes(y=-7, x=Order, label=`TCRseq`), size=size_)+
geom_text(data=df_meta, aes(y=-8.75, x=Order, label=`MultiIF`), size=size_)+
geom_text(data=df_meta, aes(y=-10.5, x=Order, label=`MultiFlow`), size=size_)+
annotation_custom(grob=textGrob(label = "Entity", gp = gpar(cex=cex_), just="right"), xmin = pos_, xmax = pos_, ymin = 0, ymax = 0)+
annotation_custom(grob=textGrob(label = "Collection", gp = gpar(cex=cex_), just="right"), xmin = pos_, xmax = pos_, ymin = -1.75, ymax = -1.75)+
annotation_custom(grob=textGrob(label = "Sex", gp = gpar(cex=cex_), just="right"), xmin = pos_, xmax = pos_, ymin = -3.5, ymax = -3.5)+
annotation_custom(grob=textGrob(label = "CITE-seq", gp = gpar(cex=cex_), just="right"), xmin = pos_, xmax = pos_, ymin = -5.25, ymax = -5.25)+
annotation_custom(grob=textGrob(label = "TCR-seq", gp = gpar(cex=cex_), just="right"), xmin = pos_, xmax = pos_, ymin = -7, ymax = -7)+
annotation_custom(grob=textGrob(label = "Multi-IF", gp = gpar(cex=cex_), just="right"), xmin = pos_, xmax = pos_, ymin = -8.75, ymax = -8.75)+
annotation_custom(grob=textGrob(label = "Multi-flow", gp = gpar(cex=cex_), just="right"), xmin = pos_, xmax = pos_, ymin = -10.5, ymax = -10.5)+
scale_x_discrete(limits=as.character(1:101))+
scale_y_continuous(expand = c(0,0), name=NULL)+
geom_vline(xintercept = seq(0.5, nrow(df_meta), 1), size=1, color="white")+
scale_fill_manual(values = colors_characteristics)+
coord_cartesian(clip = "off")+
theme_void()+
theme(legend.position = "none",
plot.margin = unit(c(0,0.1,0,0.1), "lines"))
plot_legend1 <-
ggplot()+
geom_tile(data=df_meta, aes(y=4, x=Order, fill=Entity))+
scale_x_discrete(expand = c(0,0), name=NULL)+
scale_y_discrete(expand = c(0,0), name=NULL)+
geom_vline(xintercept = seq(0.5, 51.5, 1), color="white")+
scale_fill_manual(values = colors_characteristics, limits=df_meta$Entity %>% unique())+
guides(fill=guide_legend(direction = "horizontal", keywidth = 0.3, keyheight = 0.4))+
theme_void()+
mytheme_1+
theme(legend.position = "bottom",
legend.title = element_text(size=7, face = "bold"),
legend.text = element_text(size=7))
legend1 <- get_legend(plot_legend1)
plot_legend2 <-
ggplot()+
geom_tile(data=df_meta, aes(y=4, x=Order, fill=Sex))+
scale_x_discrete(expand = c(0,0), name=NULL)+
scale_y_discrete(expand = c(0,0), name=NULL)+
geom_vline(xintercept = seq(0.5, 51.5, 1), color="white")+
scale_fill_manual(values = colors_characteristics, limits=df_meta$Sex %>% unique())+
guides(fill=guide_legend(direction = "horizontal", keywidth = 0.3, keyheight = 0.4))+
theme_void()+
mytheme_1+
theme(legend.position = "bottom",
legend.title = element_text(size=7, face = "bold"),
legend.text = element_text(size=7))
legend2 <- get_legend(plot_legend2)
plot_legend3 <-
ggplot()+
geom_tile(data=df_meta, aes(y=4, x=Order, fill=Status))+
scale_x_discrete(expand = c(0,0), name=NULL)+
scale_y_discrete(expand = c(0,0), name=NULL)+
geom_vline(xintercept = seq(0.5, 51.5, 1), color="white")+
scale_fill_manual(values = colors_characteristics, name="Collection",
limits=filter(df_meta, Status!="NA")$Status %>% unique())+
guides(fill=guide_legend(direction = "horizontal", keywidth = 0.3, keyheight = 0.4))+
theme_void()+
mytheme_1+
theme(legend.position = "bottom",
legend.title = element_text(size=7, face = "bold"),
legend.text = element_text(size=7))
legend3 <- get_legend(plot_legend3)
```
## Associations with overall T-cell frequencies
```{r patient characteristics SF1 part 2}
# Sex
p2 <- df_meta %>%
filter(!Entity %in% c("CLL", "rLN")) %>%
mutate(Tcells=ifelse(is.na(Tcells), TcellsDx, Tcells)) %>%
filter(Tcells>1) %>%
ggplot(aes(x=Sex, y=Tcells))+
geom_boxplot(width=0.35, size=0.25)+
#ggbeeswarm::geom_beeswarm(size=0.8, shape=21, stroke=0.25, cex = 2.75, aes(fill=Entity))+
stat_compare_means(comparisons = list(c("F", "M")), label.y = 78, size=2.5, bracket.size = 0.25)+
scale_fill_manual(values = hue_pal()(5))+
scale_x_discrete(labels=c("Female", "Male"))+
ylim(0,90)+
ylab("T-cells in %")+
mytheme_1+
theme(axis.title.x = element_blank(),
plot.tag = element_text(margin = unit(c(0,0,0,-0.75), "cm")))+
labs(tag = "B")
# Status
p3 <- df_meta %>%
filter(!Entity %in% c("CLL", "rLN")) %>%
mutate(Tcells=ifelse(is.na(Tcells), TcellsDx, Tcells)) %>%
filter(Tcells>1) %>%
ggplot(aes(x=Status, y=Tcells))+
geom_boxplot(width=0.35, size=0.25)+
#ggbeeswarm::geom_beeswarm(size=0.8, shape=21, stroke=0.25, cex = 2.75, aes(fill=Entity))+
stat_compare_means(comparisons = list(c("Initial diagnosis", "Relapse")), label.y = 78, bracket.size = 0.25, size=2.5)+
scale_x_discrete(labels=c("Initial \ndiagnosis", "Relapse"))+
scale_fill_manual(values = hue_pal()(5))+
ylim(0,90)+
ylab("T-cells in %")+
mytheme_1+
theme(axis.title.x = element_blank())+
labs(tag = "C")
# Age
p5 <- df_meta %>%
filter(!Entity %in% c("CLL", "rLN")) %>%
mutate(Tcells=ifelse(is.na(Tcells), TcellsDx, Tcells)) %>%
filter(Tcells>1) %>%
ggplot(aes(x=Age, y=Tcells))+
geom_point(size=1.25, shape=21, stroke=0.25, color="white", aes(fill=Entity))+
stat_cor(aes(label=..r.label..), size=2.5, label.y = 82)+
scale_fill_manual(values = colors_characteristics, name=NULL, limits=c("DLBCL", "FL", "MCL", "MZL"))+
ylim(0,90)+
ylab("T-cells in %")+
mytheme_1+
theme(legend.position = "right",
axis.title.x = element_text(vjust = 6.5),
legend.box.margin = unit(c(0,0,0,-0.38), "cm"),
legend.key.width = unit("cm", x = 0.15),
legend.key.height = unit("cm", x = 0.35))+
labs(tag = "D")
# Entity
df_tmp <- df_meta %>%
filter(!Entity %in% c("CLL", "rLN")) %>%
mutate(Tcells=ifelse(is.na(Tcells), TcellsDx, Tcells)) %>%
filter(Tcells>1)
p6 <- ggplot(df_meta, aes(x=Entity, y=Tcells))+
geom_boxplot(width=0.4, size=0.25)+
ggbeeswarm::geom_beeswarm(size=0.75, shape=21, stroke=0.25, cex = 2.25, aes(fill=Entity))+
stat_compare_means(data=df_tmp %>% filter(!Entity %in% c("rLN")),
label.y = 83, label.x = 2.5, size=2.5, hjust=0.5, bracket.size = 0.25)+
scale_fill_manual(values = colors_characteristics, name=NULL,
limits=c("DLBCL", "FL", "MCL", "MZL", "rLN"))+
geom_segment(data = NULL, aes(x=1, xend=4, y=80, yend=80), size=0.1)+
ylim(0,90)+
ylab("T-cells in %")+
mytheme_1+
theme(axis.title.x = element_blank())+
labs(tag = "E")
```
## Assemble plot
```{r assemble SF1, fig.height=5.5}
# Plot
wrap_plots(p1+labs(tag = "A")+p_ann+plot_layout(nrow = 2, heights = c(1, 0.5)))/
plot_spacer()/
wrap_plots(p2+p3+p5+p6+plot_layout(nrow = 1, widths = c(0.55,0.55,1,1)))+
plot_layout(heights = c(1.25,0.075,0.7))
#ggsave(width = 18.5, height = 12, units = "cm", filename = "SF1.pdf")
```
## Legend
```{r legend SF1, fig.height=0.5}
# Legend
as_ggplot(legend1)+as_ggplot(legend2)+as_ggplot(legend3)+plot_layout(nrow = 1, widths = c(1, 0.5, 0.75))
ggsave(width = 13, height = 1, units = "cm", filename = "S1_p1_legend.pdf")
```
# Supplementary Figure 2
## Confusion Matrix
```{r, confusion SF2, fig.height=4.5}
confmat <- confusionMatrix(GBclasses_surface$test, GBclasses_surface$predict)
conf_freq <- confmat$table %>%
data.frame() %>%
group_by(Reference) %>%
mutate(Prop=Freq/sum(Freq))
plot_surface <- ggplot(conf_freq, aes(x=Reference, y=Prediction, fill=Prop, label=round(Prop, 2)))+
geom_tile()+
geom_text(data=conf_freq %>% filter(Prop>0.4), color="white", size=1.75)+
scale_fill_gradientn(colours = colorRampPalette(RColorBrewer::brewer.pal(9, "BuPu"))(100), name="Sensitivity",
limits=c(0,0.9), breaks=seq(0,0.8,0.2))+
geom_hline(yintercept = seq(1.5, 13.5), color="black",linetype="solid", size=0.25, alpha=0.25)+
geom_vline(xintercept = seq(1.5, 13.5), color="black",linetype="solid", size=0.25, alpha=0.25)+
scale_y_discrete(expand = c(0,0), limits=factor(cluster_order), labels=unlist(labels_cl))+
scale_x_discrete(expand = c(0,0), limits=factor(cluster_order), labels=unlist(labels_cl))+
coord_fixed()+
ggtitle("Surface marker only")+
mytheme_1+
theme(legend.position = "none",
axis.text.x = element_text(angle=45, hjust = 1))+
labs(tag = "A")
confmat_plus <- confusionMatrix(GBclasses_surfaceplus$test, GBclasses_surfaceplus$predict)
conf_freq <- confmat_plus$table %>%
data.frame() %>%
group_by(Reference) %>%
mutate(Prop=Freq/sum(Freq))
plot_surfaceplus <- ggplot(conf_freq, aes(x=Reference, y=Prediction, fill=Prop, label=round(Prop, 2)))+
geom_tile()+
geom_text(data=conf_freq %>% filter(Prop>0.4), color="white", size=1.75)+
scale_fill_gradientn(colours = colorRampPalette(RColorBrewer::brewer.pal(9, "BuPu"))(100), name="Sensitivity",
limits=c(0,0.9), breaks=seq(0,0.8,0.2))+
geom_hline(yintercept = seq(1.5, 13.5), color="black",linetype="solid", size=0.25, alpha=0.25)+
geom_vline(xintercept = seq(1.5, 13.5), color="black",linetype="solid", size=0.25, alpha=0.25)+
scale_y_discrete(expand = c(0,0), limits=factor(cluster_order), labels=unlist(labels_cl))+
scale_x_discrete(expand = c(0,0), limits=factor(cluster_order), labels=unlist(labels_cl))+
coord_fixed()+
ggtitle("Surface plus FoxP3, IKZF3, and Ki67")+
mytheme_1+
theme(legend.position = "right",
legend.key.height = unit(0.35, "cm"),
legend.key.width = unit(0.28, "cm"),
legend.margin = margin(c(0,0.2,0,0), unit = "cm"),
axis.text.x = element_text(angle=45, hjust = 1))+
labs(tag = "B")
plot_surface+plot_surfaceplus+plot_layout(guides = "collect", widths = c(1,1))
#ggsave(width = 19, height = 9, units = "cm", filename = "SF2.pdf")
```
# Supplementary Figure 3
## Frequencies overview
```{r freq overview SF3}
df_pop
mat_complete <- rbind(
df_facs %>%
left_join(., df_meta %>% select(PatientID, `CITEseq`)) %>%
filter(`CITEseq`=="-") %>%
select(PatientID, Population, Prop=FACS),
df_freq %>%
mutate(RNA=ifelse(is.nan(RNA), 0, RNA)) %>%
select(PatientID, Population, Prop=RNA)) %>%
filter(Population %in% df_pop$Population) %>%
left_join(., df_pop) %>%
select(-Population) %>%
left_join(., df_meta %>% select(PatientID, Tcells)) %>%
mutate(Prop=Tcells*Prop/100) %>%
select(-Tcells) %>%
pivot_wider(names_from = "IdentI", values_from = "Prop", values_fill = 0) %>%
column_to_rownames("PatientID")
cl_order <- c(6,1,2,9,8,13,15,11,12,16,3,5,14,19)
names(labels_cl_parsed) <- as.character(cluster_order)
df_freqPlot <-
mat_complete %>%
rownames_to_column("PatientID") %>%
pivot_longer(cols=2:ncol(.), names_to = "IdentI", values_to = "Prop") %>%
add_entity() %>%
mutate(Entity=factor(Entity, levels = c("rLN", "DLBCL", "MCL", "FL", "MZL"))) %>%
mutate(IdentI=factor(IdentI, levels=cl_order)) %>%
mutate(Label=factor(IdentI, levels=cl_order, labels = labels_cl_parsed[as.character(cl_order)])) %>%
group_by(Entity, IdentI) %>%
mutate(outlier = (Prop > quantile(Prop, 0.75) + IQR(Prop) * 1.5) | (Prop < quantile(Prop, 0.25) - IQR(Prop) * 1.5))
df_medianLines <- df_freqPlot %>%
filter(Entity=="rLN") %>%
group_by(IdentI, Label) %>%
summarise(MedianProp=median(Prop))
df_freqPlot_pvalues <- df_freqPlot %>%
group_by(IdentI) %>%
wilcox_test(data=., formula = Prop ~ Entity, detailed = T, ref.group = "rLN") %>%
adjust_pvalue(method = "BH") %>%
select(IdentI, Entity=group2, p.adj, estimate) %>%
mutate(Entity=factor(Entity, levels = c("rLN", "DLBCL", "MCL", "FL", "MZL"))) %>%
mutate(p.adj=ifelse(p.adj>0.05, NA, p.adj)) %>%
mutate(p.adj_s=format(p.adj, scientific = TRUE, digits=1)) %>%
mutate(p.adj_f=case_when(p.adj > 0.05 ~ "NA",
p.adj==0.05 ~ "0.05",
p.adj < 0.05 & p.adj > 0.001 ~ as.character(round(p.adj, 3)),
p.adj==0.001 ~ "0.001",
p.adj < 0.001 ~ p.adj_s)) %>%
filter(!is.na(p.adj)) %>%
left_join(., df_freqPlot %>% select(IdentI, Label) %>% distinct) %>%
left_join(., data.frame(IdentI=factor(cl_order), height=c(28, 28, 21, 21, 10, 10, 12, 12, 14, 14, 50, 50, 12, 12)))
```
## Plot
```{r freq overview, fig.height=6, fig.width=5.2}
p <- list()
for(i in c(1:7)){
y <- list(c(1,6),c(2,9),c(8,13),c(15,11),c(12,16),c(3,5),c(14,19))[[i]]
ylim <- c(38,30,13,15.5,20,65,15)
p[[i]] <-
ggplot(data=df_freqPlot %>% filter(IdentI %in% y) %>%
mutate(Label=factor(Label, levels = labels_cl_parsed[as.character(y)])),
aes(y=Prop, x=Entity, fill=IdentI))+
geom_hline(data=df_medianLines %>%filter(IdentI %in% y), aes(yintercept=MedianProp),
size=0.25, linetype="dashed", color="grey60")+
geom_boxplot(width=0.4, outlier.shape = 21, outlier.size = 1, outlier.color = "white",
outlier.alpha = 0, show.legend = F, size=0.25)+
ggbeeswarm::geom_beeswarm(data = function(x) dplyr::filter_(x, ~ outlier), cex = 3, stroke=0.25,
groupOnX = TRUE, shape = 21, size = 1, color = "white", alpha = 1)+
geom_text(data=df_freqPlot_pvalues %>% filter(IdentI %in% y),
inherit.aes = F, aes(y=height, x=Entity, label=p.adj_f), hjust=0.1, size=2.3, angle=45)+
scale_fill_manual(values = colors_umap_cl)+
scale_y_continuous(name="% of total cells", limits=c(0,ylim[i]))+
scale_x_discrete(expand = c(0.17,0.17))+
facet_wrap(~Label, ncol = 2, labeller = label_parsed)+
mytheme_1+
theme(axis.title.x = element_blank(),
strip.background = element_rect(color=NA),
plot.margin = unit(c(0.1,0.25,0,0), units = "cm"),
panel.border = element_rect(size = 0.5),
axis.text.x = element_text(angle=45, hjust = 1))
if(i!=7){
p[[i]] <- p[[i]]+
theme(axis.text.x = element_blank(),
axis.ticks.x = element_blank())
}
if(!i %in% c(4,5)){
p[[i]] <- p[[i]]+
theme(axis.title.y = element_blank())
}
}
plot_freq <- wrap_plots(p, ncol = 2)
plot_freq
ggsave(width = 18.5, height = 15, units = "cm", filename = "SF3.pdf")
```
# Supplementary Figure 4
## Lasso model (beta coefficients)
```{r lasso model SF4}
mat_complete <- rbind(
df_facs %>%
left_join(., df_meta %>% select(PatientID, `CITEseq`)) %>%
filter(`CITEseq`=="-") %>%
select(PatientID, Population, Prop=FACS),
df_freq %>%
mutate(RNA=ifelse(is.nan(RNA), 0, RNA)) %>%
select(PatientID, Population, Prop=RNA)) %>%
filter(Population %in% df_pop$Population) %>%
left_join(., df_pop) %>%
select(-Population) %>%
pivot_wider(names_from = "IdentI", values_from = "Prop", values_fill = 0) %>%
column_to_rownames("PatientID")
total <- mat_complete %>%
rownames_to_column("PatientID") %>%
left_join(., df_meta %>% select(PatientID, Tcells)) %>%
add_entity() %>%
mutate_if(is.numeric, .funs = ~./100)
cell_types <- total %>% select(-Entity, -PatientID) %>% colnames()
gt <- my_glmnet(total)
limits_coefs=c(cluster_order, "Tcells", "(Intercept)")
plot_coefs <- gt$coefs %>%
mutate(Entity=factor(Entity, levels = c("rLN", "DLBCL", "MCL", "FL", "MZL"))) %>%
ggplot(aes(x=beta, y=cell_type, fill=cell_type, color=cell_type))+
geom_hline(yintercept = c(1.5, 5.5, 9.5, 10.5, 13.5)+2, linetype="solid", size=0.25, alpha=0.1)+
geom_bar(stat = "identity", width = 0.5)+
scale_y_discrete(limits=rev(limits_coefs), labels=c("Intercept", "T cells", rev(unlist(unname(labels_cl)))))+
scale_color_manual(limits=rev(limits_coefs), values = c("black", "black", rev(unname(colors_umap_cl[as.character(cluster_order)]))))+
scale_fill_manual(limits=rev(limits_coefs), values = c("black", "black", rev(unname(colors_umap_cl[as.character(cluster_order)]))))+
facet_wrap(~Entity, ncol=5)+
mytheme_1+
theme(axis.title.y = element_blank(),
panel.border = element_rect(size = 0.5))
```
## Multivariate Model
```{r multivariate model SF4, fig.height=7}
models <- list()
entities <- c("Overall", "FL", "MZL", "MCL")
for(i in colnames(mat_complete)){
for(e in entities) {
ent <- e
if(e=="Overall"){ent <- c("FL", "MZL", "MCL", "DLBCL")}
models[[paste0(i, "_", e)]] <- mat_complete %>%
rownames_to_column("PatientID") %>%
pivot_longer(cols=2:ncol(.), names_to = "IdentI", values_to = "Prop") %>%
left_join(., df_meta %>% select(PatientID, Status, Pretreatment, Sex, Entity, Department, Age, Subtype, Tcells, `CITEseq`) %>% distinct,
by="PatientID") %>%
filter(Entity %in% ent) %>%
filter(IdentI==i) %>%
glm(formula = Prop ~ Status + Age + Sex,
family = "gaussian") %>% summary()
models[[paste0(i, "_", e)]] <- models[[paste0(i, "_", e)]]$coefficients %>%
`colnames<-`(c("Estimate", "Std.error", "t.value", "p.value")) %>%
data.frame() %>%
dplyr::mutate(IdentI=as.character(i)) %>%
rownames_to_column("Parameter")
}
}
for(i in colnames(mat_complete)){
ent <- "DLBCL"
e <- "DLBCL"
models[[paste0(i, "_", e)]] <- mat_complete %>%
rownames_to_column("PatientID") %>%
pivot_longer(cols=2:ncol(.), names_to = "IdentI", values_to = "Prop") %>%
left_join(., df_meta %>% select(PatientID, Status, Pretreatment, Sex, Entity, Department, Age, Subtype, Tcells, `CITEseq`) %>% distinct,
by="PatientID") %>%
filter(Entity %in% ent) %>%
filter(IdentI==i) %>%
glm(formula = Prop ~ Status + Age + Sex + Subtype,
family = "gaussian") %>% summary()
models[[paste0(i, "_", e)]] <- models[[paste0(i, "_", e)]]$coefficients %>%
`colnames<-`(c("Estimate", "Std.error", "t.value", "p.value")) %>%
data.frame() %>%
dplyr::mutate(IdentI=as.character(i)) %>%
rownames_to_column("Parameter")
}
p.value <- bind_rows(models, .id = "model") %>%
filter(Parameter!="(Intercept)") %>%
mutate(Entity=strsplit(model, split = "_") %>% sapply(., "[[", 2)) %>%
group_by(IdentI, Entity) %>%
mutate(p.adj=p.adjust(p.value, method = "BH")) %>%
mutate(sign=ifelse(Estimate<0 & p.adj<0.05, "-", NA)) %>%
mutate(sign=ifelse(Estimate>0 & p.adj<0.05, "+", sign)) %>%
mutate(Entity=factor(Entity, levels = c("Overall", "DLBCL", "MCL", "FL", "MZL")))
plot_multi <- ggplot()+
geom_point(data=p.value %>% filter(Entity=="Overall"), aes(x=IdentI, y=Parameter, size=-log10(p.adj), fill=t.value), shape=21, stroke=0.25,
position = position_nudge(y=0.16))+
geom_point(data=p.value %>% filter(Entity!="Overall"), aes(x=IdentI, y=Parameter, size=-log10(p.adj), fill=t.value, group=Entity), shape=21, stroke=0.25,
position = ggpp::position_dodgenudge(y = -0.14, width = 0.78))+
geom_text(data=p.value %>% filter(Entity=="Overall", p.adj<0.05), aes(x=IdentI, y=Parameter, label=round(p.adj, 3)), hjust=0.5, nudge_y = 0.4, size=2.5)+
geom_text(data=p.value %>% filter(Entity=="Overall", p.adj<0.05), aes(x=IdentI, y=Parameter, label=sign), hjust=0.5, nudge_y = 0.18, size=2.5)+
geom_text(data=p.value %>% filter(Entity!="Overall"), aes(x=IdentI, y=Parameter, label=sign, group=Entity), size=2.5,
position = ggpp::position_dodgenudge(y = -0.11, width = 0.78))+
scale_x_discrete(limits=factor(cluster_order), labels=unlist(labels_cl))+
scale_size_continuous(range=c(0.25,3.5), limits=c(0,3), name=expression('-log'[10]~'p'))+
scale_fill_gradientn(colors=brewer.pal(9, name = "RdBu"), limits=c(-5, 5), breaks=c(-3,0,3), name="Statistic")+
guides(fill=guide_colorbar(barwidth = 2.5))+
guides(size=guide_legend(keywidth = 0.3))+
geom_hline(yintercept = seq(1,6,1), color="grey90", size=0.25)+
geom_vline(xintercept = seq(1.5,13.5,1), color="grey90", size=0.25)+
scale_y_discrete(limits=rev(c("StatusRelapse", "Age", "SexM", "Subtypenon-GCB")),
labels=rev(c("<b style='color:#B2182B'>Initial diagnosis </b><br>vs. <b style='color:#2166AC'>Relapse</b>",
"<b style='color:#B2182B'>Younger </b>vs.<br><b style='color:#2166AC'>Older</b>",
"<b style='color:#B2182B'>Female </b>vs.<br><b style='color:#2166AC'>Male</b>",
"<b style='color:#B2182B'>GCB</b> vs.<br><b style='color:#2166AC'>non-GCB</b><br>[only DLBCL]")))+
mytheme_1+
theme(axis.title = element_blank(),
legend.position = "top",
plot.title = element_text(margin = unit(c(0,0,0,0), units = "cm")),
axis.text.x = element_text(angle = 45, hjust = 1),
panel.border = element_rect(size=0.5),
legend.key.height = unit(x = 0.3, units = "cm"),
legend.box.margin = unit(c(0,-8.5,-0.35,0), "cm"),
plot.margin = unit(c(0,0.25,0,0.25), "cm"),
legend.key.width = unit(x = 0.35, units = "cm"),
axis.text.y = element_markdown(hjust = 1, lineheight = 1.2))+
labs(tag = "B")
plot_coefs+labs(tag = "A")+plot_multi+plot_layout(heights = c(1.2,1))
ggsave(width = 18, height = 15, units = "cm", filename = "SF4.pdf")
```
# Supplementary Figure 5
## Proportions of T-cell subsets in 5' scRNA versus CITE-seq
### Data handling
```{r data handling SF5}
df_freq_5prime <- DFtotal_5prime %>%
select(Barcode_full, IdentI, PatientID) %>%
distinct() %>%
add_prop(vars = c("PatientID", "IdentI"), group.vars = 1) %>%
mutate(Class="5prime") %>%
filter(PatientID %in% unique(df_comb$PatientID))
df_freq_3prime <- df_comb %>%
add_prop(vars = c( "PatientID", "IdentI"), group.vars = 1) %>%
mutate(Class="3prime") %>%
filter(PatientID %in% DFtotal_5prime$PatientID)
df_freq_prime <- rbind(df_freq_3prime, df_freq_5prime) %>%
mutate(Prop=100*Prop) %>%
pivot_wider(names_from = "Class", values_from = "Prop", values_fill = 0)
```
### Correlation plots
```{r cor plots SF5}
this_theme <- theme(plot.margin = unit(c(0,0.25,0.1,0), units = "cm"),
plot.title = element_text(vjust = -1),
axis.title = element_blank(),
axis.text = element_text(size=7, color="black"))
cor_plots_prime <- list()
cor_plots_prime[["TFH"]] <-
df_freq_prime %>%
filter(IdentI==6) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["6"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["6"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson", size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,60), breaks = c(0, 20, 40, 60))+
scale_y_continuous(limits = c(0,60), breaks = c(0, 20, 40, 60))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["6"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
this_theme
cor_plots_prime[["TPR"]] <- df_freq_prime %>%
filter(IdentI==14) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["14"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["14"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,12), breaks = c(0, 4, 8, 12))+
scale_y_continuous(limits = c(0,12), breaks = c(0, 4, 8, 12))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["14"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
this_theme+
theme(axis.title.y = element_text(size=7))
cor_plots_prime[["TDN"]] <- df_freq_prime %>%
filter(IdentI==19) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["19"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["19"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,8), breaks = c(0, 2, 4, 6, 8))+
scale_y_continuous(limits = c(0,8), breaks = c(0, 2, 4, 6, 8))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["19"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
theme_axis_sub3
cor_plots_prime[["THCM2"]] <- df_freq_prime %>%
filter(IdentI==9) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["9"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["9"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,15), breaks = c(0, 5, 10, 15))+
scale_y_continuous(limits = c(0,15), breaks = c(0, 5, 10, 15))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["9"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
this_theme
cor_plots_prime[["THCM1"]] <- df_freq_prime %>%
filter(IdentI==2) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["2"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["2"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,20), breaks = c(0, 6, 12, 18))+
scale_y_continuous(limits = c(0,20), breaks = c(0, 6, 12, 18))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["2"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
this_theme
cor_plots_prime[["THNaive"]] <- df_freq_prime %>%
filter(IdentI==1) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["1"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["1"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,45), breaks = c(0, 15, 30, 45))+
scale_y_continuous(limits = c(0,45), breaks = c(0, 15, 30, 45))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["1"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
this_theme
cor_plots_prime[["TREGCM1"]] <- df_freq_prime %>%
filter(IdentI==8) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["8"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["8"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["8"]])+
scale_x_continuous(limits = c(0,25), breaks = c(0, 8, 16, 24))+
scale_y_continuous(limits = c(0,25), breaks = c(0, 8, 16, 24))+
theme_bw()+
coord_fixed()+
mytheme_1+
this_theme
cor_plots_prime[["TREGEM2"]] <- df_freq_prime %>%
filter(IdentI==11) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["11"]],size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["11"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,16), breaks = c(0, 5, 10, 15))+
scale_y_continuous(limits = c(0,16), breaks = c(0, 5, 10, 15))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["11"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
theme_axis_sub3
cor_plots_prime[["TREGEM1"]] <- df_freq_prime %>%
filter(IdentI==15) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["15"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["15"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,20), breaks = c(0, 6, 12, 18))+
scale_y_continuous(limits = c(0,20), breaks = c(0, 6, 12, 18))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["15"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
this_theme+
theme(axis.title.y = element_text(size=7),
axis.title.x = element_text(size=7))
cor_plots_prime[["TREGCM2"]] <- df_freq_prime %>%
filter(IdentI==13) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["13"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["13"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,20), breaks = c(0, 6, 12, 18))+
scale_y_continuous(limits = c(0,20), breaks = c(0, 6, 12, 18))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["13"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
this_theme
cor_plots_prime[["TTOXNaive"]] <-df_freq_prime %>%
filter(IdentI==12) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["12"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["12"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,27), breaks = c(0, 8, 16, 24))+
scale_y_continuous(limits = c(0,27), breaks = c(0, 8, 16, 24))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["12"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
theme_axis_sub3
cor_plots_prime[["TTOXEM3"]] <- df_freq_prime %>%
filter(IdentI==5) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color="#b50923", size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color="#b50923", se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,80), breaks = c(0, 25, 50, 75))+
scale_y_continuous(limits = c(0,80), breaks = c(0, 25, 50, 75))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["5"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
theme_axis_sub3
cor_plots_prime[["TTOXEM1"]] <- df_freq_prime %>%
filter(IdentI==3) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color="#fea044", size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x,na.rm = T,
color="#fea044", se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,50), breaks = c(0, 15, 30, 45))+
scale_y_continuous(limits = c(0,50), breaks = c(0, 15, 30, 45))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["3"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
theme_axis_sub3
cor_plots_prime[["TTOXEM2"]] <- df_freq_prime %>%
filter(IdentI==16) %>%
ggplot(aes(y=`5prime`, x=`3prime`))+
geom_point(color=colors_umap_cl[["16"]], size=0.65, alpha=0.75)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_umap_cl[["16"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson",size=2.5, label.x.npc = c(0.45), label.y.npc = c(0.9))+
scale_x_continuous(limits = c(0,15), breaks = c(0, 4, 8, 12))+
scale_y_continuous(limits = c(0,15), breaks = c(0, 4, 8, 12))+
labs(x="3' scRNA - %", y="5' scRNA - %")+
ggtitle(labels_cl[["16"]])+
theme_bw()+
coord_fixed()+
mytheme_1+
theme_axis_sub3
cor_plots_prime_all <- cor_plots_prime$TPR+labs(tag ="A")+
cor_plots_prime$THNaive+cor_plots_prime$THCM1+cor_plots_prime$THCM2+cor_plots_prime$TFH+
cor_plots_prime$TREGCM1+cor_plots_prime$TREGCM2+cor_plots_prime$TREGEM1+cor_plots_prime$TREGEM2+
cor_plots_prime$TTOXNaive+cor_plots_prime$TTOXEM1+cor_plots_prime$TTOXEM2+cor_plots_prime$TTOXEM3+
cor_plots_prime$TDN+
plot_layout(nrow = 2)
df_freq_prime %>%
group_by(IdentI) %>%
summarise(R=cor.test(`3prime`, `5prime`)$estimate) %>% pull(R) %>% median()
```
## CD4 and CD8 expression in clonal T-cells
```{r CD4 and CD8 expression in SF5}
lim_CD4 <- 0.3
lim_CD8 <- 0.6
df_clone_expr <- lapply(sobjs_T_5prime, function(sobj){
FetchData(sobj, vars=c("CD4", "CD8A")) %>%
data.frame() %>%
rownames_to_column("Barcode_full")
}) %>% bind_rows() %>% remove_rownames()
df_clone_CD4 <- DFtotal_5prime %>%
select(Barcode_full, PatientID, refUMAP_1, refUMAP_2, Entity, IdentI, raw_clonotype_id) %>%
distinct() %>%
add_count(IdentI, PatientID, raw_clonotype_id) %>%
left_join(., df_clone_expr) %>%
mutate(IdentI_new=factor(IdentI, levels = cluster_order, labels = labels_cl_parsed)) %>%
filter(IdentI %in% c(3,5)) %>%
mutate(Expanded=n>2) %>%
mutate(CD4CD8=case_when(CD4 > lim_CD4 & CD8A < lim_CD8 ~ "CD4pos",
CD4 < lim_CD4 & CD8A > lim_CD8 ~ "CD8pos",
CD4 > lim_CD4 & CD8A > lim_CD8 ~ "CD4CD8pos",
CD4 < lim_CD4 & CD8A < lim_CD8 ~ "CD4CD8neg"))
df_clone_CD4_freq <- df_clone_CD4 %>%
drop_na() %>%
add_prop(vars = c("IdentI", "Expanded", "CD4CD8"), group.vars = c(1,2)) %>%
fill_zeros(names_from = "IdentI", values_from = "Prop") %>%
mutate(Prop=round(Prop, 2)) %>%
left_join(., data.frame(CD4CD8=c("CD4CD8neg", "CD4pos", "CD4CD8pos", "CD8pos"),
CD8=c(-0.35, -0.35, 3.7, 3.7),
CD4=c(-0.35, 1.9, 1.9, -0.35))) %>%
mutate(IdentI_new=factor(IdentI, levels = cluster_order, labels = labels_cl_parsed))
plot_expr1 <- ggplot()+
geom_point(data=df_clone_CD4 %>% filter(n<3), inherit.aes = F, aes(x=CD4, y=CD8A, color=IdentI),
position = position_jitter(width = 0.1, height = 0.1), na.rm = T,
size=0.25, alpha=0.25, stroke=0)+
scale_color_manual(values = colors_umap_cl, guide="none")+
scale_fill_manual(values = colors_umap_cl, guide="none")+
geom_text(data=df_clone_CD4_freq %>% filter(Expanded==F), inherit.aes = F, aes(x=CD4, y=CD8, label=Prop), size=2.5, alpha=1)+
geom_hline(yintercept = 0.6, linetype="dashed", size=0.25)+
geom_vline(xintercept = 0.3, linetype="dashed", size=0.25)+
scale_x_continuous(breaks = c(0,0.5,1,1.5,2), limits=c(-0.5, 2.4))+
scale_y_continuous(breaks = c(0,1,2,3), limits=c(-0.5, 4))+
labs(x="<i>CD4</i> - RNA expression",
y="<i>CD8</i> - RNA expression",
title="Clone size < 3",
tag = "B")+
facet_wrap(~IdentI_new, nrow = 1, labeller = label_parsed)+
mytheme_1+
theme(axis.title.x = element_textbox(size=7, halign = 0.5, margin = unit(units = "cm", c(0.1,0,0,0))),
axis.title.y = element_textbox(size=7, orientation = "left-rotated", margin = unit(units = "cm", c(0,0,0.1,0))),
panel.border = element_rect(size=0.4))
plot_expr2 <- ggplot()+
geom_point(data=df_clone_CD4 %>% filter(n>2), inherit.aes = F, aes(x=CD4, y=CD8A, color=IdentI),
position = position_jitter(width = 0.1, height = 0.1), na.rm = T,
size=0.25, alpha=0.25, stroke=0)+
scale_color_manual(values = colors_umap_cl, guide="none")+
scale_fill_manual(values = colors_umap_cl, guide="none")+
scale_size_continuous(range=c(1, 5), limits=c(3, 50), breaks=c(3, 20, 35, 50),
labels=c("3", "20", "35", "> 50"), name = "Clonotype size")+
geom_hline(yintercept = 0.6, linetype="dashed", size=0.25)+
geom_vline(xintercept = 0.3, linetype="dashed", size=0.25)+
geom_text(data=df_clone_CD4_freq %>% filter(Expanded==T), inherit.aes = F, aes(x=CD4, y=CD8, label=Prop), size=2.5, alpha=1)+
facet_wrap(~IdentI_new, nrow = 1, labeller = label_parsed)+
scale_x_continuous(breaks = c(0,0.5,1,1.5,2), limits=c(-0.5, 2.4))+
scale_y_continuous(breaks = c(0,1,2,3), limits=c(-0.5, 4))+
labs(x="<i>CD4</i> - RNA expression",
y="<i>CD8</i> - RNA expression",
title="Clone size > 2",
tag = "C")+
facet_wrap(~IdentI_new, nrow = 1, labeller = label_parsed)+
mytheme_1+
theme(axis.title.x = element_textbox(size=7, halign = 0.5, margin = unit(units = "cm", c(0.1,0,0,0))),
axis.title.y = element_textbox(size=7, orientation = "left-rotated", margin = unit(units = "cm", c(0,0,0.1,0))),
panel.border = element_rect(size=0.4))
```
## Assemble plot
```{r assemble SF5}
cor_plots_prime_all/wrap_plots(plot_expr1+plot_spacer()+plot_expr2+plot_layout(widths = c(1,0.1,1)))+
plot_layout(heights = c(1.6,1))
#ggsave(width = 18.5, height = 12.5, units = "cm", filename = "SF5.pdf")
```
## TCR diversity
### Read
```{r TCR diversity read}
# Single cell T-cell receptor data read by immunarch package
# RNA–seq, epitope and TCR raw and processed data have been deposited in the Gene Expression Omnibus (GEO) under accession codes GSE252608 and GSE252455.
DF_immunarchTCR <- repLoad(list.files(path = "countMatrices", pattern = "TCRrep", full.names = T))
DF_immunarchTCR$meta$Sample <- strsplit(DF_immunarchTCR$meta$Sample, split = "_") %>% sapply("[[", 1)
names(DF_immunarchTCR$data) <- DF_immunarchTCR$meta$Sample
```
### Plot
```{r TCR diversity plot, fig.height=3.5}
plots <- list()
for(i in unique(DFtotal_5prime$PatientID)){
set.seed(substr(i, 4,7) %>% as.numeric()+5)
plots[[i]] <-
DFtotal_5prime %>% filter(PatientID==i) %>%
select(Barcode_full, PatientID, raw_clonotype_id) %>%
distinct() %>%
dplyr::count(raw_clonotype_id) %>%
drop_na() %>%
mutate(Prop=n/sum(n)) %>%
dplyr::arrange(-n) %>%
mutate(Cumsum=cumsum(n)) %>%
mutate(Max=sum(n)) %>%
filter(Cumsum < 0.1*Max) %>%
mutate(new_id=as.character(1:nrow(.))) %>%
mutate(PatientID=i) %>%
add_entity() %>%
mutate(PatientID_new=paste0(PatientID, " (", Entity, ")")) %>%
ggplot(aes(x=new_id, y=n))+
geom_segment(aes(x=new_id, xend=new_id, y=0, yend=n), size=0.2)+
facet_wrap(~PatientID_new)+
xlab("Clonotype ID")+
scale_y_continuous(name = "Clonotype size")+
scale_x_discrete(limits=sample(as.character(1:195), 195), name="Unique clonotype ID")+
mytheme_1+
theme(legend.position = "none",
panel.border = element_rect(size=0.4),
strip.background = element_rect(colour = NA),
axis.title.x = element_text(vjust = 5),
axis.text.x = element_blank(),
axis.ticks.x = element_blank())
}
imm_raref <- repDiversity(DF_immunarchTCR$data, "raref", .verbose = F) %>%
rename(PatientID=Sample) %>%
add_entity() %>%
filter(!PatientID %in% c("LN0256", "LN0367"))
tops <- imm_raref %>%
group_by(PatientID) %>%
top_n(n=1, Size) %>%
#mutate(Mean=ifelse(PatientID=="LN0302", 5.0, Mean)) %>%
mutate(Size=ifelse(PatientID=="LN0132", 88, Size)) %>%
mutate(Size=ifelse(PatientID=="LN0110", 72, Size)) %>%
mutate(Mean=ifelse(PatientID=="LN0110", 19, Mean)) %>%
mutate(Size=ifelse(PatientID=="LN0259", 90, Size)) %>%
mutate(Mean=ifelse(PatientID=="LN0259", 20.75, Mean)) %>%
mutate(Size=ifelse(PatientID=="LN0264", 55, Size)) %>%
mutate(Mean=ifelse(PatientID=="LN0264", 13.5, Mean)) %>%
mutate(Size=ifelse(PatientID=="LN0217", Size-10, Size)) %>%
mutate(Size=ifelse(PatientID=="LN0144", 145, Size)) %>%
mutate(Mean=ifelse(PatientID=="LN0193", Mean+1, Mean)) %>%
mutate(Size=ifelse(PatientID=="LN0198", 115, Size)) %>%
mutate(Mean=ifelse(PatientID=="LN0198", Mean+0.2, Mean)) %>%
mutate(Mean=ifelse(PatientID=="LN0278", Mean+1, Mean)) %>%
mutate(Mean=ifelse(PatientID=="LN0144", Mean+0.75, Mean)) %>%
mutate(Mean=ifelse(PatientID=="LN0302", Mean-1.5, Mean)) %>%
mutate(Size=ifelse(PatientID=="LN0302", 50, Size)) %>%
filter(!PatientID %in% c("LN0417", "LN0104", "LN0249")) # Manuel labelling in Powerpoint
plot_rare <-
imm_raref %>%
ggplot(aes(x=Size, y=Mean, group=PatientID, color=Entity))+
geom_line(linetype="solid", size=0.25)+
geom_label(data=tops, aes(x=Size, y=Mean, label=PatientID), size=2.25,
show.legend = F, fill="white", color="white",
label.padding = unit(units = "cm", 0.02), label.size = 0)+
geom_text(data=tops, aes(x=Size, y=Mean, label=PatientID), size=2, show.legend = F)+
guides(color=guide_legend(override.aes = list(size=0.35, linetype="solid")))+
scale_color_brewer(palette = "Paired", limits=c("DLBCL", "MCL", "FL", "MZL", "rLN"))+
ylab("Estimated diversity")+
xlab("Clone size")+
#xlim(0, 70)+
mytheme_1+
theme(#legend.position = "top",
legend.title = element_blank(),
panel.grid = element_blank(),
legend.spacing.x = unit("cm", x = 0.05),
legend.box.margin = unit(c(0,0,-0.95,0), "cm"),
panel.border = element_rect(size=0.4),
legend.key.height = unit("cm", x = 0.36),
legend.key.width = unit("cm", x = 0.5))+
labs(tag = "D")
plot_rare+wrap_plots((plots$LN0132+theme(axis.title.x = element_blank())+labs(tag = "E"))/plots$LN0217+labs(tag = "F"))+
plot_layout(widths = c(1,1.7))
#ggsave(width = 18, height = 8, units = "cm", filename = "SF5.pdf")
```
# Supplementary Figure 6
## T-cell exhaustion UMAP
### Calculate exhaustion module
```{r calculate exhaustion module}
exhausted_cells <- ttox@meta.data %>%
mutate(Exhausted=ifelse( Pseudotime>=24, "yes", "no")) %>%
filter(Exhausted=="yes") %>% rownames()
Combined_T@meta.data$Exhausted <- Combined_T@meta.data %>%
mutate(Exhausted=Barcode_full %in% exhausted_cells) %>%
pull(Exhausted)
Idents(Combined_T) <- "Exhausted"
module_exhausted <- FindMarkers(Combined_T, ident.1 = "TRUE", assay = "integratedRNA", test.use = "roc") %>%
rownames_to_column("Gene") %>%
mutate(Module=paste0(Gene, ifelse(myAUC>0.5, "+", "-")),
Assay="RNA")
module_exhausted_prot <- FindMarkers(Combined_T, ident.1 = "TRUE", assay = "integratedADT", test.use = "roc") %>%
rownames_to_column("Gene") %>%
mutate(Module=paste0(Gene, ifelse(myAUC>0.5, "+", "-")),
Assay="Protein")
#module_exhausted --> Supplementary Table 5
#WriteXLS::WriteXLS(rbind(module_exhausted, module_exhausted_prot), ExcelFileName = "SuppTable5.xlsx")
module_exhausted <- list(module_exhausted$Module)
names(module_exhausted) <- "exhausted"
Combined_T <- UCell::AddModuleScore_UCell(Combined_T, features = module_exhausted)
```
### Plot exhaustion module
```{r plot exhaustion module}
set.seed(1)
plot_exh <- FetchData(Combined_T, vars = c("wnnUMAP_1", "wnnUMAP_2", "exhausted_UCell")) %>%
sample_frac(0.2) %>%
ggplot(aes(x=wnnUMAP_1, y=wnnUMAP_2, fill= exhausted_UCell))+
ggrastr::geom_point_rast(size=0.25, stroke=0, shape=21, raster.dpi = 600, alpha=0.75)+
scale_fill_gradientn(colours = brewer.pal(n = 9, name = "YlOrRd")[2:9],
name="Score")+
xlab("wnnUMAP-1")+
ylab("wnnUMAP-2")+
ggtitle("Exhaustion signature")+
mytheme_1+
theme(panel.border = element_rect(size = 0.2),
axis.title.x = element_text(margin = unit(units = "cm", c(-0.75,0,0,0))),
#legend.margin = margin(c(0,0,0,-0.35), unit = "cm"),
legend.box.margin = unit(c(0,0,0,-0.35), units = "cm"),
legend.title = element_text(size=6),
legend.text = element_text(size=6),
legend.position = "right",
plot.title = element_text(face = "plain", vjust = -0.5),
panel.background = element_rect(fill=NA),
legend.key.height = unit(units="cm", 0.2),
legend.key.width = unit(units="cm", 0.15),
legend.box.background = element_rect(fill=NA, color=NA),
legend.background = element_rect(fill=NA, color=NA)
)+
labs(tag="A")
```
## Association with cell-of-origin in DLBCL
```{r SF6 part 1}
### Schmitz
p1 <- left_join(df_surv_schmitz, df_ttoxcompl_schmitz) %>%
drop_na() %>%
ggplot(aes(x=Subtype, y=Exhausted/Absolute, group=Subtype))+
geom_boxplot(outlier.alpha = 0, width=0.4, size=0.25)+
stat_compare_means(size=2.25, vjust = 1, aes(label=paste0("p = ", ..p.format..)))+
stat_compare_means(comparisons = list(c("ABC", "GCB")), size=2.25)+
scale_y_continuous(limits=c(0, 0.65), name = "Exhausted T-cells")+
#ggtitle("Schmitz et al. 2018")+
xlab("Cell-of-origin")+
mytheme_1+
theme(axis.text.x = element_text(angle = 45, hjust = 1),
#axis.title.x = element_text(margin = unit(units = "cm", c(-0.1,0,0,0))),
plot.title = element_text(face = "plain", vjust = -0.5))+
labs(tag = "B")
### Chapuy
p2 <- left_join(df_surv_chapuy, df_ttoxcompl_chapuy) %>%
drop_na() %>%
ggplot(aes(x=Subtype, y=Exhausted/Absolute, group=Subtype))+
geom_boxplot(outlier.alpha = 0, width=0.4, size=0.25)+
stat_compare_means(size=2.25, vjust = 1, aes(label=paste0("p = ", ..p.format..)))+
scale_y_continuous(limits=c(0, 0.65), name="Exhausted T-cells")+
#ggtitle("Chapuy et al. 2018")+
xlab("Cell-of-origin")+
mytheme_1+
theme(axis.text.x = element_text(angle = 45, hjust = 1),
#axis.title.x = element_text(margin = unit(units = "cm", c(-0.1,0,0,0))),
plot.title = element_text(face = "plain", vjust = -0.5))+
labs(tag = "C")
```
## Association with genetic Subtype
```{r SF6 part 2}
### Schmitz
p3 <- left_join(df_surv_schmitz, df_ttoxcompl_schmitz) %>%
drop_na() %>%
ggplot(aes(x=GenSubtype, y=Exhausted/Absolute, group=GenSubtype))+
geom_boxplot(outlier.alpha = 0, width=0.4, size=0.25)+
stat_compare_means(size=2.25, vjust = 1, aes(label=paste0("p = ", ..p.format..)))+
scale_y_continuous(limits = c(0, 0.65), name="Exhausted T-cells")+
#ggtitle("Schmitz et al. 2018")+
xlab("Cluster")+
mytheme_1+
theme(axis.text.x = element_text(angle = 45, hjust = 1),
plot.title = element_text(face = "plain", vjust = -0.5),
axis.title.y = element_blank(),
axis.text.y = element_blank(),
axis.ticks.y = element_blank(),
#axis.title.x = element_text(margin = unit(units = "cm", c(-0.1,0,0,0)))
)
### Chapuy
p4 <- left_join(df_surv_chapuy, df_ttoxcompl_chapuy) %>%
drop_na() %>%
mutate(Cluster=paste0("C", Cluster)) %>%
ggplot(aes(x=Cluster, y=Exhausted/Absolute, group=Cluster))+
geom_boxplot(outlier.alpha = 0, width=0.4, size=0.25)+
stat_compare_means(size=2.25, vjust = 1, aes(label=paste0("p = ", ..p.format..)))+
scale_y_continuous(limits=c(0, 0.65), name="Exhausted T-cells")+
#ggtitle("Chapuy et al. 2018")+
xlab("Cluster")+
mytheme_1+
theme(axis.text.x = element_text(angle = 45, hjust = 1),
plot.title = element_text(face = "plain", vjust = -0.5),
axis.title.y = element_blank(),
axis.text.y = element_blank(),
axis.ticks.y = element_blank(),
#axis.title.x = element_text(margin = unit(units = "cm", c(-0.65,0,0,0)))
)
```
## Assemble plot
### Part 1
```{r assemble plot SF6 part 1, fig.height=2.5}
plot_exh+p1+p3+p2+p4+plot_layout(nrow = 1, widths = c(1.4,0.75,1,0.75,1))
#ggsave(width = 19, height = 6, units = "cm", filename = "SF6_p1.pdf")
```
## Associations with genetic features
### Mutations
```{r mutations SF6}
pairs <- list(c(1:42), c(43:85))
df_mut <- df_snvs_chapuy %>%
filter(Description=="Mutation", value %in% c(0,2)) %>%
mutate(value=factor(value, levels = c("0", "2"), labels = c("wt", "mut"))) %>%
group_by(Name, value) %>%
dplyr::summarise(MeanExhausted=mean(Exhausted/Absolute)) %>%
mutate(MeanExhausted=ifelse(MeanExhausted>0.15, 0.15, MeanExhausted),
Group=case_when(Name %in% unique(.$Name)[pairs[[1]]] ~ "group1",
Name %in% unique(.$Name)[pairs[[2]]] ~ "group2"))
p.values <- df_snvs_chapuy %>%
filter(Description=="Mutation", value %in% c(0,2)) %>%
mutate(Exhausted=Exhausted/Absolute) %>%
compare_means(formula = Exhausted ~ value, group.by = "Name", p.adjust.method = "BH") %>%
filter(p<0.05) %>%
left_join(., df_mut %>% select(Name, Group), by="Name")
p5 <- df_mut %>%
ggplot(aes(x=value, y=Name, fill=MeanExhausted))+
geom_tile()+
scale_fill_gradientn(name="Expression", colours = brewer.pal(5, "PRGn"), limits=c(0,0.17))+
geom_vline(xintercept = 1.5, color="black", size=0.25)+
geom_segment(data=data.frame(y=seq(1.5, 42.5,1)), inherit.aes = F,
aes(y = y, yend = y, x=0.05, xend=2.8), color="white", size=0.25)+
geom_text(data = p.values, inherit.aes = F, aes(y=Name, x=2.8, label=round(p, 2)), size=2.5)+
facet_wrap(~Group, ncol=3, scales = "free_y")+
ggtitle("Mutations")+
coord_cartesian(clip = "off")+
mytheme_1+
theme(strip.background = element_blank(),
strip.text = element_blank(),
axis.text.x = element_text(angle=45, hjust=1, size=7),
axis.ticks = element_blank(),
axis.text.y = element_text(size=7, margin = unit(c(0,-0.4,0,0), units = "cm")),
plot.margin = unit(units = "cm", c(0,1.25,0,0)),
axis.title = element_blank(),
panel.border = element_blank())+
labs(tag = "E")
```
### Copy number gain
```{r copy number SF6}
df_gain <- df_snvs_chapuy %>%
filter(Description=="CN gain") %>%
mutate(value=factor(value, levels = c("0", "1", "2"), labels = c("wt", "gain", "gain"))) %>%
group_by(Name, value) %>%
dplyr::summarise(MeanExhausted=mean(Exhausted/Absolute))
p.values <- df_snvs_chapuy %>%
filter(Description=="CN gain") %>%
mutate(Exhausted=Exhausted/Absolute) %>%
mutate(value=factor(value, levels = c("0", "1", "2"), labels = c("wt", "gain", "gain"))) %>%
compare_means(formula = Exhausted ~ value, group.by = "Name", p.adjust.method = "BH") %>%
filter(p<0.05)
p6 <- df_gain %>%
ggplot(aes(x=value, y=Name, fill=MeanExhausted))+
geom_tile()+
scale_fill_gradientn(name="Expression", colours = brewer.pal(5, "PRGn"), limits=c(0,0.17))+
geom_vline(xintercept = 1.5, color="black", size=0.25)+
geom_segment(data=data.frame(y=seq(1.5, 31.5,1)), inherit.aes = F,
aes(y = y, yend = y, x=0.25, xend=2.8), color="white", size=0.25)+
geom_text(data = p.values, inherit.aes = F, aes(y=Name, x=2.8, label=round(p, 2)), size=2.5)+
coord_cartesian(clip = "off")+
ggtitle("CN gain")+
mytheme_1+
theme(strip.background = element_blank(),
strip.text = element_blank(),
axis.text.x = element_text(angle=45, hjust=1, size=7),
axis.ticks = element_blank(),
axis.text.y = element_text(size=7, margin = unit(c(0,-0.15,0,0), units = "cm")),
plot.margin = unit(units = "cm", c(0,1.25,0,0)),
axis.title = element_blank(),
panel.border = element_blank())+
labs(tag = "F")
```
### Copy number loss
```{r copy number loss SF6}
df_loss <- df_snvs_chapuy %>%
filter(Description=="CN loss") %>%
mutate(value=factor(value, levels = c("0", "1", "2"), labels = c("wt", "loss", "loss"))) %>%
group_by(Name, value) %>%
dplyr::summarise(MeanExhausted=mean(Exhausted/Absolute))
p.values <- df_snvs_chapuy %>%
filter(Description=="CN loss") %>%
mutate(Exhausted=Exhausted/Absolute) %>%
mutate(value=factor(value, levels = c("0", "1", "2"), labels = c("wt", "loss", "loss"))) %>%
compare_means(formula = Exhausted ~ value, group.by = "Name", p.adjust.method = "BH") %>%
filter(p<0.05) %>%
mutate(p=ifelse(p<0.005, 0.0051, p))
p7 <- df_loss %>%
ggplot(aes(x=value, y=Name, fill=MeanExhausted))+
geom_tile()+
scale_fill_gradientn(name="Expression", colours = brewer.pal(5, "PRGn"), limits=c(0,0.17))+
geom_vline(xintercept = 1.5, color="black", size=0.25)+
geom_segment(data=data.frame(y=seq(1.5, 32.5,1)), inherit.aes = F,
aes(y = y, yend = y, x=0.25, xend=2.8), color="white", size=0.25)+
geom_text(data = p.values, inherit.aes = F, aes(y=Name, x=2.8, label=round(p, 2)), size=2.5)+
ggtitle("CN loss")+
coord_cartesian(clip = "off")+
mytheme_1+
theme(strip.background = element_blank(),
strip.text = element_blank(),
axis.text.x = element_text(angle=45, hjust=1, size=7),
axis.ticks = element_blank(),
axis.text.y = element_text(size=7, margin = unit(c(0,-0.15,0,0), units = "cm")),
axis.title = element_blank(),
panel.border = element_blank())+
labs(tag = "G")
```
### Structural variants
```{r structural variants SF6}
df_struct <- df_snvs_chapuy %>%
filter(Description=="SV") %>%
mutate(value=factor(value, levels = c("0", "3"), labels = c("wt", "mut")))%>%
group_by(Name, value) %>%
dplyr::summarise(MeanExhausted=mean(Exhausted/Absolute))
p.values <- df_snvs_chapuy %>%
filter(Description=="SV") %>%
mutate(Exhausted=Exhausted/Absolute) %>%
mutate(value=factor(value, levels = c("0", "3"), labels = c("wt", "mut"))) %>%
compare_means(formula = Exhausted ~ value, group.by = "Name", p.adjust.method = "BH") %>%
filter(p<0.05)
p8 <- df_struct %>%
ggplot(aes(x=value, y=Name, fill=MeanExhausted))+
geom_tile()+
scale_fill_gradientn(name="Expression", colours = brewer.pal(5, "PRGn"), limits=c(0,0.17))+
geom_vline(xintercept = 1.5, color="black", size=0.25)+
geom_segment(data=data.frame(y=seq(1.5, 7.5,1)), inherit.aes = F,
aes(y = y, yend = y, x=0.25, xend=2.75), color="white", size=0.25)+
ggtitle("Structural variants")+
coord_cartesian(clip = "off")+
mytheme_1+
theme(strip.background = element_blank(),
strip.text = element_blank(),
axis.text.x = element_text(angle=45, hjust=1, size=7),
axis.ticks = element_blank(),
axis.text.y = element_text(size=7, margin = unit(c(0,-0.15,0,0), units = "cm")),
plot.margin = unit(units = "cm", c(0,1.25,0,0)),
axis.title = element_blank(),
panel.border = element_blank())+
labs(tag = "H")
```
### Part 2
```{r assemble plot SF6 part 2, fig.height=7}
p5+wrap_plots(p6/p8+plot_layout(heights = c(2,0.5)))+(p7/plot_spacer()+plot_layout(heights = c(2,0.5)))+
plot_layout(nrow = 1, widths = c(2.9,1,1))
#ggsave(width = 18, height = 14, units = "cm", filename = "SF6_p2.pdf")
```
### Part 3 (Legend)
```{r assemble plot SF6 part 3, fig.height=1}
as_ggplot(get_legend(p8+guides(fill=guide_colorbar(nrow = 2, title = "Exhausted\nT-cells"))+
theme(legend.position = "right",
legend.key.height = unit(units = "cm", 0.3),
legend.key.width = unit(units = "cm", 0.3))))
ggsave(width = 2, height = 2.4, units = "cm", filename = "SF6_legend.pdf")
```
# Supplementary Figure 7
Panel A was generated using FlowJo.
## Flow cytometry: IKZF3
```{r, fig.height=3}
med <- df_ikzf3 %>% filter(Entity=="rLN") %>% pull(`FoxP3+/IKZF3+`) %>% median()
pvalues <- df_ikzf3 %>% rename(IKZF3=`FoxP3+/IKZF3+`) %>%
data.frame() %>%
compare_means(data=., formula = IKZF3 ~ Entity, ref.group = "rLN") %>%
select(Entity=group2, p) %>%
filter(p<0.05) %>%
mutate(p=round(p,3))
nrow(df_ikzf3)
plot_aiolos <- df_ikzf3 %>%
ggplot(aes(x=Entity,y=`FoxP3+/IKZF3+`))+
geom_hline(yintercept = med, size=0.25, linetype="dashed", color="grey60")+
geom_boxplot(width=0.5, outlier.alpha = 0, size=0.25)+
ggbeeswarm::geom_beeswarm(size=0.75, shape=21, stroke=0.25, cex = 2.25, aes(fill=Entity))+
geom_text(inherit.aes = F, data = pvalues %>% mutate(Y=c(75, 75)),
aes(x=Entity, y=Y, label=p), size=2.5, check_overlap = T)+
scale_fill_brewer(palette = "Paired", limits=c("DLBCL", "MCL", "FL", "MZL", "rLN"))+
scale_y_continuous(limits = c(0,80), name=expression('% IKZF3'^'+'~'of FoxP3'^'+'))+
scale_x_discrete(limits=c("rLN", "DLBCL", "MCL", "FL", "MZL"))+
ggtitle("Flow cytometry")+
mytheme_1+
theme(legend.position = "none",
strip.background = element_rect(color=NA),
axis.title.x = element_blank(),
plot.title = element_text(face = "plain", size=7),
panel.border = element_rect(size=0.5),
axis.text.x = element_text(angle=45, hjust = 1, size=7),
axis.text.y = element_text(size=7),
axis.title.y = element_text(size=7),
panel.background = element_rect(fill=NA),
plot.margin = unit(c(0,0.25,0,0.25), "cm"))
plot_spacer()+plot_aiolos+plot_layout(widths = c(3,1))
ggsave(width = 19, height = 5.7, units = "cm", filename = "SF7.pdf")
```
## Treg clonotypes
```{r treg clonotypes SF7, fig.height=3}
df_clonotypes_shared <-
left_join(DFtotal_5prime %>% filter(!is.na(raw_clonotype_id)) %>%
select(Barcode_fulla=Barcode_full, PatientID, refUMAP_1a=refUMAP_1, refUMAP_2a=refUMAP_2, IdentIa=IdentI, raw_clonotype_id) %>% distinct(),
DFtotal_5prime %>% filter(!is.na(raw_clonotype_id)) %>%
select(Barcode_fullb=Barcode_full, PatientID, refUMAP_1b=refUMAP_1, refUMAP_2b=refUMAP_2, IdentIb=IdentI, raw_clonotype_id) %>% distinct()
) %>%
filter(Barcode_fulla!=Barcode_fullb) %>%
filter(IdentIa!=IdentIb)
treg_shared <- list()
for(i in c(8,13,15)){
df_subset <-
df_clonotypes_shared %>%
add_entity() %>%
filter(IdentIb==i)
treg_shared[[i]] <- ggplot()+
geom_point_rast(data=DFtotal_5prime,
aes(x=refUMAP_1, y=refUMAP_2, fill=IdentI), size=0.25,
alpha=ifelse(DFtotal_5prime$IdentI==i, 0.4, 0.04), stroke=0, shape=21)+
geom_curve(data= df_subset,
aes(x=refUMAP_1a, y=refUMAP_2a, xend=refUMAP_1b, yend=refUMAP_2b, color=IdentIa,
group=paste(raw_clonotype_id, PatientID)), curvature = -0.4, size=0.15, alpha=0.4)+
scale_color_manual(values = colors_umap_cl, limits=factor(cluster_order),
labels=unlist(labels_cl), guide="none")+
scale_fill_manual(values = colors_umap_cl, limits=factor(cluster_order), guide="none",
labels=unlist(labels_cl))+
guides(fill=guide_legend(nrow = 7, byrow = F, override.aes = list(size=1.75, stroke=0, shape=21, alpha=1, color="white")))+
coord_cartesian(clip = "off")+
xlab("refUMAP-1")+
ylab("refUMAP-2")+
mytheme_1+
theme(legend.position = "right",
panel.border = element_rect(size=0.25),
plot.title = element_textbox_simple(size=7, width = NULL, face = "plain",
padding = margin(1.25, 0, 1, 0),
lineheight = 1.25,
halign=0.5),
legend.text = element_text(size=7),
legend.spacing.x = unit("cm", x = 0.13),
axis.title.x = element_text(size=7),
axis.title.y = element_text(size=7),
axis.text = element_text(size=7),
legend.spacing.y = unit("cm", x = 0.001),
legend.key.width = unit("cm", x = 0.05),
legend.key.height = unit("cm", x = 0.5),
legend.box.margin = margin(unit = "cm",c(0,-0.35,0,-1)),
legend.title = element_blank())
if(i==8)
treg_shared[[i]] <- treg_shared[[i]]+
labs(title="Paired clonotypes of <span style='color:#C6DBEF'>T<sub>REG</sub> CM<sub>1</sub></span>",
tag = "C")
if(i==13)
treg_shared[[i]] <- treg_shared[[i]]+
labs(title="Paired clonotypes of <span style='color:#6BAED6'>T<sub>REG</sub> CM<sub>2</sub></span>")+
theme(axis.ticks.y = element_blank(),
axis.title.y = element_blank(),
axis.text.y = element_blank())
if(i==15)
treg_shared[[i]] <- treg_shared[[i]]+
labs(title="Paired clonotypes of <span style='color:#2171B5'>T<sub>REG</sub> EM<sub>1</sub></span>")+
theme(axis.ticks.y = element_blank(),
axis.title.y = element_blank(),
axis.text.y = element_blank())
}
plot_treg <- treg_shared[[8]]+treg_shared[[13]]+treg_shared[[15]]+plot_layout(guides = "collect")
#plot_treg
```
## Survival analysis
```{r survival plots}
kmplot_tfh <- readRDS("data/SurvPlot_Tfh_Gallium.rds")
kmplot_tfh$plot$plot_env$legend <- c(0.32, 0.2)
kmplot_tfh$plot$theme$legend.position <- c(0.32, 0.2)
kmplot_tfh$plot$theme$legend.background$fill <- NA
kmplot_tfh$plot$theme$legend.text$colour <- NA
kmplot_tfh$plot$theme$legend.text$size <- 7
kmplot_tfh$plot <- kmplot_tfh$plot+annotation_custom(grob = textGrob(label = expression('T'[FH]~'High'), gp = gpar(cex=0.5), x=0.36, y=0.215))+
annotation_custom(grob = textGrob(label = expression('T'[FH]~'Low'), gp = gpar(cex=0.5), x=0.36, y=0.125))+
labs(tag = "D")
kmplot_treg <- readRDS("data/SurvPlot_TregEff2_Gallium.rds")
kmplot_treg$plot$plot_env$legend <- c(0.32, 0.2)
kmplot_treg$plot$theme$legend.position <- c(0.32, 0.2)
kmplot_treg$plot$theme$legend.background$fill <- NA
kmplot_treg$plot$theme$legend.text$colour <- NA
kmplot_tfh$plot$theme$legend.text$size <- 7
kmplot_treg$plot <- kmplot_treg$plot+annotation_custom(grob = textGrob(label = expression('T'[REG]~'EM'[2]~'High'), gp = gpar(cex=0.5), x=0.36, y=0.215))+
annotation_custom(grob = textGrob(label = expression('T'[REG]~'EM'[2]~'Low'), gp = gpar(cex=0.5), x=0.36, y=0.125))+
labs(tag = "E")
```
## Assemble plot
```{r}
plot_treg/
wrap_plots(kmplot_tfh$plot+kmplot_treg$plot+plot_spacer())+
plot_layout(heights = c(1.2,1))
#ggsave(width = 18, height = 11, units = "cm", filename = "SF7.pdf")
```
# Supplementary Figure 8
## Dendrogram T-cells
```{r dendrogram}
# Create data frame
data <- data.frame(
level1="_Tcells",
level2=c("_'T'[Pr]",
rep("_'T'[H]",2),
"_'T'[FH]",
rep("_'T'[REG]",1),
rep("_'T'[TOX]",3)),
level3=c("_'T'[Pr]",
"TH_'CD4'^'+'*' Naive'",
"TH_'non-Naive (CM'[1]*' + CM'[2]*')'",
"_'T'[FH]",
"_'T'[REG]",
"TTOX_'CD8'^'+'*' Naive'",
"TTOX_'non-Naive (EM'[1]*' + EM'[2]*')'",
"TTOX_' Exhausted (EM'[3]*')'")
)
# Data handling
edges_level1_2 <- data %>% select(level1, level2) %>% unique %>% rename(from=level1, to=level2)
edges_level2_3 <- data %>% select(level2, level3) %>% unique %>% rename(from=level2, to=level3)
edge_list=rbind(edges_level1_2, edges_level2_3)
vert <- data.frame(
name=unique(c(data$level1, data$level2, data$level3))) %>%
mutate(cluster=as.character(c(NA, 14, 'TH', 6, 'TREG', "TTOX", 1, 2, 12, 3, 5))) %>%
mutate(label=strsplit(name, split = "_") %>% sapply(., "[[", 2))
# Make ggraph object
mygraph_codex <- graph_from_data_frame( edge_list ,vertices = vert)
# Small codex dendrogramm
ggraph(mygraph_codex, layout = 'tree', circular = FALSE) +
geom_edge_diagonal(strength = 1.4, edge_width=0.25)+
geom_node_label(aes(label=label),
parse = T, nudge_y=0.11, label.padding = unit(units = "cm", 0.12),
size=2.35, alpha=1,
fill=c(rep("white", 4), "black", rep("white", 6)), vjust=1, color=NA,
label.size = 0, label.r = unit(units = "cm", 0))+
geom_node_text(aes(label=label, color=cluster),
vjust=1, nudge_y=0.05,
parse = T,
alpha=c(0,rep(1,5),rep(0,5)),
size=2.35)+
scale_color_manual(values = colors_dendrogramm_codex)+
coord_cartesian(clip = "off")+
ggtitle("T-cell subsets \nidentified by mIF")+
theme_void()+
theme(legend.position = "none")+
theme(plot.margin = unit(c(-0.5,0.35,0.25,0.25), units = "cm"),
plot.title = element_text(hjust = 0.4, vjust=-7.5, size=7))
#ggsave(filename = "Figure8_p2.pdf", width = 12.6, height = 4.25, units = "cm")
```
## Handle data
```{r handle data SF8}
# Read CODEX expression data
# Available at BioStudies database (https://www.ebi.ac.uk/biostudies/) under accession number S-BIAD565
codex_expression <- data.table::fread("data/cells_expression.csv") %>% tibble() %>%
rename(unique_cell_id=V1)
proteins_selected <- c("PAX5", "CD20", "CD79a", "CD21", "PDPN", "CD38", "MCT", "GRZB", "CD56", "CD163", "CD206", "CD11c",
"CD15", "CD34", "CD31", "CD90", "Ki67", "PD1", "CXCR5", "ICOS", "CD69", "CD45RO", "TIM3", "LAG3",
"CD57", "CD8", "CD45RA", "FOXP3", "CD4", "CD3")
codex_meanExp <- codex_expression %>%
left_join(., codex_annotation %>% select(unique_cell_id, Merged_final)) %>%
filter(!is.na(Merged_final)) %>%
select(-unique_cell_id) %>%
group_by(Merged_final) %>%
summarise_all(mean) %>%
pivot_longer(cols = 2:ncol(.), names_to = "Protein", values_to = "Expression") %>%
group_by(Protein) %>%
mutate(Expression=(Expression-min(Expression))/(max(Expression)-min(Expression))) %>%
filter(Protein %in% proteins_selected)
```
## Makers and cell types
```{r plot markers SF8}
p2 <- codex_meanExp %>%
ggplot(aes(x=Protein, y=Merged_final, fill=Expression))+
geom_tile()+
scale_fill_gradientn(name="Expression", colours = brewer.pal(5, "GnBu"), limits=c(0,1), breaks=c(0,0.5,1))+
geom_hline(yintercept = seq(1.5, 17.5, 1), size=0.25, color="white")+
geom_vline(xintercept = seq(1.5, 32.5, 1), size=0.25, color="white")+
scale_x_discrete( expand = c(0,0), limits=c(proteins_selected))+
scale_y_discrete(expand = c(0,0), limits=c("B", "FDC", "PC", "MC", "NK", "NKT", "Macro", "DC", "Granulo", "Stromal cells",
"TPR", "TFH", "TTOX_exh", "TTOX", "TTOXNaive", "Treg", "CD4T", "CD4TNaive"),
labels=unlist(list("B-cells", "FDC", "Plasma cells", "Mast cells", "NK cells", "NK T-cells", "Macrophages", "Dendritic cells", "Granulocytes",
"Stromal cells", labels_codex$TPR, labels_codex$TFH, labels_codex$TTOX_exh, labels_codex$TTOX, labels_codex$TTOXNaive,
labels_codex$Treg, labels_codex$CD4T, labels_codex$CD4TNaive)))+
guides(fill=guide_colorbar(ticks.colour = "black"))+
theme_bw()+
mytheme_1+
theme(axis.title = element_blank(),
legend.position = "top",
axis.text.y = element_text(size=7),
legend.text = element_text(size = 7, color="black"),
legend.title = element_text(size = 7, color="black", vjust = 0.8, margin = unit(units = "cm", c(0,0.2,0,0))),
legend.key.height = unit(0.25, "cm"),
plot.margin = unit(c(0,0,0,0), units = "cm"),
legend.key.width = unit(0.2, "cm"),
legend.box.spacing = unit(0.1, "cm"),
legend.box.margin = unit(c(0,0,0,0), units = "cm"),
plot.title = element_text(face = "plain", vjust = -1),
plot.tag = element_text(margin = unit(c(0,-0.5,-0.25,0), units = "cm")),
axis.text.x = element_text(size=6.5, angle = 45, hjust = 1))+
labs(tag = "C")
```
## T-cell numbers in codex
```{r T-cell numbers SF8}
df_codex_no <-
codex_annotation %>%
filter(Merged_final %in% c("TTOXNaive", "TTOX_exh", "TTOX", "Treg", "TPR", "TFH", "CD4T", "CDT4Naive")) %>%
count(PatientID, Merged_final, unique_region) %>%
group_by(PatientID) %>%
mutate(Sum=sum(n)) %>%
ungroup() %>%
mutate(No=dense_rank(desc(Sum)))
regions_random <- df_codex_no %>%
select(PatientID, unique_region) %>%
distinct() %>%
group_by(PatientID) %>%
sample_n(1)
p3 <- ggplot()+
geom_hline(yintercept = 47280, size=0.25, linetype="dashed")+
geom_bar(data=df_codex_no %>% filter(unique_region %in% regions_random$unique_region),
aes(x=No-0.15, y=n, fill=Merged_final), color="white",
stat = "identity", width=0.2, size=0.25, alpha=0.7)+
geom_bar(data=df_codex_no %>% filter(!unique_region %in% regions_random$unique_region),
aes(x=No+0.15, y=n, fill=Merged_final), color="white",
stat = "identity", width=0.25, size=0.25, alpha=0.7)+
scale_y_continuous(name = "Absolute number of cells - mIF", limits = c(0, 100000))+
scale_fill_manual(values = colors_codex[c(2:9)],
limits=limits_codex[c(2:9)],
labels=labels_codex[c(2:9)],
name=NULL)+
guides(fill=guide_legend(nrow = 2, default.unit = "cm", override.aes = list(color="white"),
keywidth = 0.3, keyheight = 0.3, byrow = T))+
scale_x_continuous(name="Patients",
breaks=unique(df_codex_no$No),
labels=unique(df_codex_no$PatientID),
expand = c(0.02,0.02))+
mytheme_1+
theme(axis.text.x = element_text(angle = 45, hjust = 1),
plot.margin = unit(c(0,0,0,0), units = "cm"),
legend.box.margin = unit(c(0,0,-0.3,0), units = "cm"),
axis.title.y = element_text(size=7, vjust = -15),
axis.title.x = element_blank(),
legend.spacing.y = unit(0.15, units = "cm"),
plot.tag = element_text(margin = unit(c(0,-0.5,0,0), units = "cm")),
legend.position = "top")+
labs(tag = "D")
```
## Correlations
```{r correlations SF8}
freq_codex <-
codex_annotation %>%
filter(Merged_final %in% c("TTOXNaive", "TTOX_exh", "TTOX", "Treg", "TPR", "TFH", "CD4T", "CD4TNaive")) %>%
filter(!unique_region %in% c("191_1reg004", "191_4reg004", "191_4reg005", "191_1reg006")) %>%
add_prop(vars = c("Merged_final", "PatientID"), group.vars = 2) %>%
rename(Prop_codex=Prop, IdentI=Merged_final)
freq_citeseq <-
Combined_T@meta.data %>%
add_prop(vars = c("IdentI", "PatientID"), group.vars = 2) %>%
mutate(IdentII=case_when(IdentI==1 ~ "CD4TNaive",
IdentI %in% c(2,9) ~ "CD4T",
IdentI==14 ~ "TPR",
IdentI==6 ~ "TFH",
IdentI %in% c(8,11,13,15) ~ "Treg",
IdentI==12 ~ "TTOXNaive",
IdentI %in% c(3,16) ~ "TTOX",
IdentI %in% c(5) ~ "TTOX_exh")) %>%
group_by(IdentII, PatientID) %>%
summarise(Prop=sum(Prop)) %>%
rename(Prop_citeseq=Prop, IdentI=IdentII) %>%
fill_zeros(names_from = "IdentI", values_from = "Prop_citeseq")
freq_joined <- left_join(freq_codex, freq_citeseq) %>%
mutate(Prop_codex=100*Prop_codex, Prop_citeseq=100*Prop_citeseq)
this_theme <-
theme_bw()+
mytheme_1+
theme(plot.margin = unit(c(0,0.2,0.2,0.35), units = "cm"),
axis.title.y = element_blank(),
axis.title.x = element_blank(),
plot.title = element_textbox_simple(face = "plain", halign=0.5, width = 2, padding = margin(0, 0, 0, 0)),
plot.tag = element_text(margin = unit(c(0,-0.5,0,0), units = "cm")),
axis.text = element_text(size=7, color="black"))
cor_plots_codex <- list()
cor_plots_codex[["TFH"]] <-
freq_joined %>%
filter(IdentI=="TFH") %>%
ggplot(aes(x=Prop_codex, y=Prop_citeseq))+
geom_point(fill=colors_codex[["TFH"]], size=1, alpha=0.75, shape=21, stroke=0.1)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_codex[["TFH"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson", size=2.5, label.x.npc = c(0.5), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,40), breaks = c(0, 12, 24, 36), name = "mIF")+
scale_y_continuous(limits = c(0,40), breaks = c(0, 12, 24, 36), name = "CITE-seq")+
labs(title="T<sub>FH</sub>")+
coord_fixed()+
this_theme+
theme(plot.margin = unit(c(0,0,0.2,0), units = "cm"))
cor_plots_codex[["TREG"]] <-
freq_joined %>%
filter(IdentI=="Treg") %>%
ggplot(aes(x=Prop_codex, y=Prop_citeseq))+
geom_point(fill=colors_codex[["Treg"]], size=1, alpha=0.75, shape=21, stroke=0.1)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color="black", se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson", size=2.5, label.x.npc = c(0.5), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,56), breaks = c(0, 16, 32, 48), name = "mIF")+
scale_y_continuous(limits = c(0,56), breaks = c(0, 16, 32, 48), name = "CITE-seq")+
labs(title="T<sub>REG</sub>")+
coord_fixed()+
this_theme+
theme(axis.title.y = element_text(size=7, angle = 90, vjust = 2.5))
cor_plots_codex[["TTOXNaive"]] <-
freq_joined %>%
filter(IdentI=="TTOXNaive") %>%
ggplot(aes(x=Prop_codex, y=Prop_citeseq))+
geom_point(fill=colors_codex[["TTOXNaive"]], size=1, alpha=0.75, shape=21, stroke=0.1)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_codex[["TTOXNaive"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson", size=2.5, label.x.npc = c(0.5), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,13), breaks = c(0, 4, 8, 12), name = "mIF")+
scale_y_continuous(limits = c(0,13), breaks = c(0, 4, 8, 12), name = "CITE-seq")+
labs(title="CD8<sup>+</sup> Naive")+
coord_fixed()+
this_theme+
theme(plot.margin = unit(c(0,0,0.2,0), units = "cm"))
cor_plots_codex[["THNaive"]] <-
freq_joined %>%
filter(IdentI=="CD4TNaive") %>%
ggplot(aes(x=Prop_codex, y=Prop_citeseq))+
geom_point(fill=colors_codex[["CD4TNaive"]], size=1, alpha=0.75, shape=21, stroke=0.1)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_codex[["CD4TNaive"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson", size=2.5, label.x.npc = c(0.53), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,30), breaks = c(0, 10, 20, 30), name = "mIF")+
scale_y_continuous(limits = c(0,30), breaks = c(0, 10, 20, 30), name = "CITE-seq")+
labs(title="CD4<sup>+</sup> Naive")+
coord_fixed()+
this_theme+
theme(plot.title = element_textbox_simple(face = "plain", halign=0.5, margin = unit(units = "cm", c(0,0,-1.75,0)),
width = 2, padding = margin(0, 0, 0, 0)),
plot.margin = unit(c(0,0,0.2,0), units = "cm"))
cor_plots_codex[["TPR"]] <-
freq_joined %>%
filter(IdentI=="TPR") %>%
ggplot(aes(x=Prop_codex, y=Prop_citeseq))+
geom_point(fill=colors_codex[["TPR"]], size=1, alpha=0.75, shape=21, stroke=0.1)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_codex[["TPR"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson", size=2.5, label.x.npc = c(0.5), label.y.npc = c(0.93))+
scale_x_continuous(limits = c(0,13), breaks = c(0, 4, 8, 12), name = "mIF")+
scale_y_continuous(limits = c(0,13), breaks = c(0, 4, 8, 12), name = "CITE-seq")+
labs(title="T<sub>Pr</sub>")+
coord_fixed()+
this_theme+
theme(axis.title.y = element_text(size=7, angle = 90, vjust = 2.5),
plot.title = element_textbox_simple(face = "plain", halign=0.5, margin = unit(units = "cm", c(0,0,-1.75,0)),
width = 2, padding = margin(0, 0, 0, 0)),
plot.tag = element_text(margin = unit(c(0,0,-0.25,0), units = "cm")))
cor_plots_codex[["TH"]] <-
freq_joined %>%
filter(IdentI=="CD4T") %>%
ggplot(aes(x=Prop_codex, y=Prop_citeseq))+
geom_point(fill=colors_codex[["CD4T"]], size=1, alpha=0.75, shape=21, stroke=0.1)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_codex[["CD4T"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson", size=2.5, label.x.npc = c(0.5), label.y.npc = c(0.93))+
scale_x_continuous(limits = c(0,60), breaks = c(0, 20, 40, 60), name = "mIF")+
scale_y_continuous(limits = c(0,60), breaks = c(0, 20, 40, 60), name = "CITE-seq")+
labs(title="Memory T<sub>H</sub>")+
coord_fixed()+
this_theme+
theme(axis.title.y = element_text(size=7, angle = 90, vjust = 2.5))
cor_plots_codex[["TTOX"]] <-
freq_joined %>%
filter(IdentI=="TTOX") %>%
ggplot(aes(x=Prop_codex, y=Prop_citeseq))+
geom_point(fill=colors_codex[["TTOX"]], size=1, alpha=0.75, shape=21, stroke=0.1)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_codex[["TTOX"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson", size=2.5, label.x.npc = c(0.5), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,45), breaks = c(0, 15, 30, 45), name = "mIF")+
scale_y_continuous(limits = c(0,45), breaks = c(0, 15, 30, 45), name = "CITE-seq")+
labs(title="Memory T<sub>TOX</sub>")+
coord_fixed()+
this_theme+
theme(axis.title.y = element_text(size=7, angle = 90, vjust = 2.5),
axis.title.x = element_text(size=7))
cor_plots_codex[["TTOX_exh"]] <-
freq_joined %>%
filter(IdentI=="TTOX_exh") %>%
ggplot(aes(x=Prop_codex, y=Prop_citeseq))+
geom_point(fill=colors_codex[["TTOX_exh"]], size=1, alpha=0.75, shape=21, stroke=0.1)+
geom_smooth(method = "lm", linetype="dashed", size=0.25, formula = y ~ x, na.rm = T,
color=colors_codex[["TTOX_exh"]], se=F, fullrange=T)+
stat_cor(aes(label=..r.label..), method = "pearson", size=2.5, label.x.npc = c(0.5), label.y.npc = c(0.1))+
scale_x_continuous(limits = c(0,65), breaks = c(0, 20, 40, 60), name = "mIF")+
scale_y_continuous(limits = c(0,65), breaks = c(0, 20, 40, 60), name = "CITE-seq")+
labs(title="PD1<sup>+</sup> TIM3<sup>+</sup> T<sub>TOX</sub>")+
coord_fixed()+
this_theme+
theme(axis.title.x = element_text(size=7),
plot.margin = unit(c(0,0,0.2,0), units = "cm"))
freq_joined %>%
group_by(IdentI) %>%
summarise(R=cor.test(Prop_codex, Prop_citeseq)$estimate) %>% pull(R) %>% median()
```
## Assemble plot
```{r assemble SF9, fig.height=5.5}
p_full <- wrap_plots(p2/p3+plot_layout(heights = c(1.25,1)))+wrap_plots(cor_plots_codex$TPR+labs(tag = "E")+cor_plots_codex$THNaive+cor_plots_codex$TH+cor_plots_codex$TFH+
cor_plots_codex$TREG+cor_plots_codex$TTOXNaive+cor_plots_codex$TTOX+cor_plots_codex$TTOX_exh+
plot_layout(ncol = 2))+
plot_layout(widths = c(2.25,1.15))
p_full
#ggsave(p_full, width = 18, height = 14, units = "cm", filename = "SF8.pdf")
```
# Supplementary Figure 9
## Immunofluorescence images
```{r images SF9, fig.height=3}
plots_codex <- list()
for(r in c("191_3reg008", "191_4reg004", "191_2reg007", "191_5reg002", "191_1reg003", "empty")) {
df_tmp <- codex_annotation %>% filter(unique_region== r) %>%
mutate(Merged_all_simple=ifelse(Merged_final %in% c("Granulo", "Macro", "DC"), "Myeloid", Merged_final)) %>%
mutate(Merged_all_simple=ifelse(Merged_all_simple %in% c("MC", "NKT", "PC", "NK"), "Other", Merged_all_simple)) %>%
filter(((x-mean(.$x))^2+(y-mean(.$y))^2)<2500^2)
plots_codex[[r]] <- ggplot()+
geom_point_rast(data=df_tmp %>% filter(Merged_all_simple=="B"), aes(x=x,y=y), shape=21, size=0.25, stroke=0, alpha=1, raster.dpi =300,
color=colors_codex[["B"]], fill=colors_codex[["B"]])+
geom_point_rast(data=df_tmp %>% filter(Merged_all_simple!="B"), aes(x=x,y=y, fill=Merged_all_simple, color=Merged_all_simple),
shape=21, size=0.25, stroke=0, alpha=1, raster.dpi=300)+
scale_color_manual(values = colors_codex, limits=limits_codex, labels=labels_codex, name="Cell type")+
scale_fill_manual(values = colors_codex, limits=limits_codex, labels=labels_codex, name="Cell type")+
ggtitle(unique(df_tmp$Entity))+
coord_fixed()+
theme_void()+
theme(legend.position = "none",
plot.title = element_text(color="white", hjust=0.1, size=8,
margin = unit(units = "cm", c(0,0,-0.6,0)), face = "bold"),
plot.margin = unit(units = "cm", c(0.1, 0.1, 0.1, 0.1)),
panel.background = element_rect(fill = "black", color="black"))
}
plots_codex
#ggsave(wrap_plots(plots_codex), width = 19, height = 12.5, units = "cm", filename = "SFigure9_p2.pdf")
legend_plot_codex <- ggplot()+
geom_point_rast(data=df_tmp %>% filter(Merged_all_simple!="B"), aes(x=x,y=y, fill=Merged_all_simple, color=Merged_all_simple),
shape=21, size=0.25, stroke=0, alpha=0, raster.dpi=300)+
scale_color_manual(values = colors_codex, limits=limits_codex, labels=labels_codex, name="Cell type")+
scale_fill_manual(values = colors_codex, limits=limits_codex, labels=labels_codex, name="Cell type")+
guides(fill=guide_legend(ncol = 2, override.aes = list(size=1.75, stroke=0, shape=21, alpha=1, color=NA)))+
guides(color=guide_legend(ncol = 2))+
ggtitle("")+
coord_fixed()+
mytheme_codex+
theme(panel.background = element_rect(fill = "black", color="black"),
legend.position = "right",
legend.text = element_text(size=6, color="white"),
legend.box.background = element_rect(fill = "black", color="black"),
legend.spacing.x = unit("cm", x = 0.13),
legend.spacing.y = unit("cm", x = 0.001),
legend.key.width = unit("cm", x = 0.05),
legend.key.height = unit("cm", x = 0.5),
legend.title = element_blank())
as_ggplot(get_legend(legend_plot_codex))
#ggsave(width = 4, height = 4, units = "cm", filename = "SFigure9_p2_legend.pdf")
```
# Supplementary Figure 10
## Load analysis
```{r run analysis SF10}
# Read results from neighborhood analysis
# Please run file: analysis/NeighborhoodAnalysis.Rmd
load("output/Neighborhood_results.RData")
# Add codex annotation
codex_annotation <- left_join(codex_annotation, nn_classes, by="unique_cell_id")
codex_annotation
```
## Tissue cores
### Images
```{r plots SF10, fig.height=6}
regions <- codex_annotation %>% pull(unique_region) %>% unique()
plots <- list()
df <- list()
for(r in regions){
df[[r]] <- codex_annotation %>%
filter(!is.na(Region), unique_region %in% r) %>%
filter(x>500, x<7500) %>%
filter(y>500, y<7500) %>%
filter(((x-mean(.$x))^2+(y-mean(.$y))^2)<2500^2)
plots[[r]] <- df[[r]] %>%
ggplot()+
ggrastr::geom_point_rast(aes(x=x,y=y,color=Region, fill=Region), shape=21, size=0.25, stroke=0, alpha=1, raster.dpi =400)+
scale_color_manual(values = colors_nn)+
scale_fill_manual(values = colors_nn)+
guides(color=guide_legend(override.aes = list(size=3)))+
ggtitle(unique(df[[r]]$PatientID))+
coord_fixed(clip = "off")+
theme_void()+
theme(legend.position = "none",
plot.margin = unit(units = "cm", c(0.1,0.1,0.1,0.1)),
plot.subtitle = element_text(size=7, face = "bold", hjust=0.5, margin = unit(units = "cm", c(0,0,0,0))),
plot.title = element_text(size=6.5, face = "plain", margin = unit(units = "cm", c(-0.1,0,-0.6,-0.1))),
panel.background = element_rect(fill = NA, color = NA),
plot.background = element_rect(fill = NA, color = NA))
}
empty <- codex_annotation %>%
filter(!is.na(Region), unique_region %in% "empty") %>%
ggplot()+
ggrastr::geom_point_rast(aes(x=x,y=y,color=Region, fill=Region), alpha=0, raster.dpi =400)+
guides(color=guide_legend(override.aes = list(size=3, alpha=1)))+
ggtitle("")+
coord_fixed(clip = "off")+
theme_void()+
theme(legend.position = "right",
plot.margin = unit(units = "cm", c(0.1,0.1,0.1,0.1)),
plot.title = element_text(size=8, face = "plain", margin = unit(units = "cm", c(0,0,-0.75,0))),
panel.background = element_rect(fill = NA, color = NA),
plot.background = element_rect(fill = NA, color = NA))
p_full <-
wrap_plots(plots$`191_1reg006`+labs(tag = "A", subtitle = "rLN")+
plots$`191_3reg007`+
plots$`191_5reg005`+plot_layout(ncol = 1))+
wrap_plots(plots$`191_4reg004`+labs(tag = "B", subtitle = "DLBCL")+#
plots$`191_3reg001`+
plots$`191_4reg006`+plot_layout(ncol = 1))+
wrap_plots(plots$`191_2reg007`+labs(tag = "C", subtitle = "MCL")+
plots$`191_2reg002`+
plots$`191_3reg006`+plot_layout(ncol = 1))+
wrap_plots(plots$`191_5reg002`+labs(tag = "D", subtitle = "FL")+
plots$`191_3reg002`+
plots$`191_5reg001`+plot_layout(ncol = 1))+
wrap_plots(plots$`191_1reg003`+labs(tag = "E", subtitle = "MZL")+
plots$`191_1reg004`+
empty+plot_layout(ncol = 1))+
plot_layout(ncol = 5)
p_full
#ggsave(p_full, width = 18, height = 11.5, units = "cm", filename = "SF11.pdf")
```
### Legend
```{r legend SF10, fig.height=1}
labels_nn <- c(
"N1: B-cells / FDC" ,
'N2: B-cells / FDC / T'[FH]~'',
'N3: T'[Pr]~'/ T'[REG]~'',
'N4: Macrophages / B-cells / Exh. T'[TOX]~'',
"N5: B-cells",
"N6: T-cell area I" ,
"N7: T-cell area II" ,
'N8: PC / NK / Memory T'[TOX]~'' ,
"N9: T-cell area III" ,
"N10: Stromal cells / Macrophages"
)
p_legend <- codex_annotation %>%
filter(!is.na(Region), unique_region %in% r) %>%
ggplot()+
ggrastr::geom_point_rast(aes(x=x,y=y,color=Region, fill=Region), shape=21, size=0.25,
stroke=0, alpha=1, raster.dpi =400)+
scale_color_manual(values = colors_nn, limits=factor(1:10), labels=labels_nn)+
scale_fill_manual(values = colors_nn, limits=factor(1:10), labels=labels_nn)+
guides(color=guide_legend(nrow = 2, override.aes = list(size=2, color="black", stroke=0.25)))+
ggtitle(unique(df[[r]]$PatientID))+
coord_fixed(clip = "off")+
theme_void()+
theme(legend.position = "right",
plot.margin = unit(units = "cm", c(0.1,0.1,0.1,0.1)),
legend.text = element_text(size=6.5),
legend.title = element_text(size=7, face = "bold"),
plot.subtitle = element_text(size=7, face = "bold", hjust=0.5, margin = unit(units = "cm", c(0,0,0,0))),
plot.title = element_text(size=6.5, face = "plain", margin = unit(units = "cm", c(-0.1,0,-0.6,-0.1))),
panel.background = element_rect(fill = NA, color = NA),
plot.background = element_rect(fill = NA, color = NA),
legend.box.margin = unit(c(0,0,0,-0.38), "cm"),
legend.key.width = unit("cm", x = 0.1),
legend.spacing.x = unit("cm", x = 0.1),
legend.key.height = unit("cm", x = 0.35))
as_ggplot(get_legend(p_legend))
#ggsave(width = 19, height = 1.5, units = "cm", filename = "SF11_legend.pdf")
```
## Composition of neighborhoods
```{r neighborhood composition, fig.height=3}
df_freq_nh <- codex_annotation %>%
add_prop(vars = c("Entity", "Region", "unique_region"), group.vars = 3) %>%
fill_zeros(names_from = "Region", values_from = "Prop") %>%
group_by(Entity, Region) %>%
mutate(Max=0.06+max(Prop),
Region_label=paste0("N", Region)) %>%
mutate(Region_label=factor(Region_label, levels = paste0("N", 1:10))) %>%
mutate(Max=ifelse(Region_label=="N8" & Entity=="MCL", 0.18, Max)) %>%
ungroup()
pvalues <- df_freq_nh %>%
compare_means(data=., formula = Prop ~ Entity, ref.group = "rLN",
group.by = "Region_label", p.adjust.method = "BH") %>%
filter(p.adj<0.05) %>%
mutate(p.adj_s=format(p.adj, scientific = TRUE, digits=1)) %>%
mutate(p.adj_f=case_when(p.adj > 0.01 ~ as.character(round(p.adj, 2)),
p.adj==0.01 ~ "0.01",
p.adj < 0.01 ~ p.adj_s),
Entity=group2) %>%
mutate(Entity=factor(Entity, levels = c("rLN", "DLBCL", "MCL", "FL", "MZL"))) %>%
left_join(., df_freq_nh %>% select(Region_label, Max, Entity) %>% distinct, by = c("Region_label", "Entity"))
df_medianLines <- df_freq_nh %>%
filter(Entity=="rLN") %>%
group_by(Region_label) %>%
dplyr::summarise(MedianProp=median(Prop))
df_freq_nh %>%
ggplot(aes(x=Entity, y=Prop)) +
geom_hline(data=df_medianLines, aes(yintercept=MedianProp),
size=0.25, linetype="dashed", color="grey60")+
geom_boxplot(width=0.5, outlier.alpha = 0, size=0.25)+
ggbeeswarm::geom_beeswarm(size=1, shape=21, stroke=0.1, cex = 2, aes(fill=Region))+
geom_text(data=pvalues, inherit.aes = F, aes(y=Max, x=Entity, label=p.adj_f), size=2.5)+
scale_fill_manual(values = colors_nn)+
scale_x_discrete(limits=c("rLN", "DLBCL", "MCL", "FL", "MZL"))+
facet_wrap(~Region_label, strip.position = "top", scales = "free_y", nrow = 2)+
scale_y_continuous(name = "% of total area", expand = c(0,0.075))+
mytheme_1+
theme(strip.text.y = element_text(angle = 0, size=6),
axis.text.x = element_text(angle=45, hjust=1),
axis.title.x = element_blank(),
strip.background = element_blank(),
plot.margin = unit(c(0,0.1,0,0.1), "cm"))+
labs(tag = "F")
ggsave(width = 18, height = 7, units = "cm", filename = "SF10.pdf")
```
# Session info
```{r session}
sessionInfo()
```