Last updated: 2022-02-22
Checks: 7 0
Knit directory: MelanomaIMC/
This reproducible R Markdown analysis was created with workflowr (version 1.7.0). The Checks tab describes the reproducibility checks that were applied when the results were created. The Past versions tab lists the development history.
Great! Since the R Markdown file has been committed to the Git repository, you know the exact version of the code that produced these results.
Great job! The global environment was empty. Objects defined in the global environment can affect the analysis in your R Markdown file in unknown ways. For reproduciblity it’s best to always run the code in an empty environment.
The command set.seed(20200728) was run prior to running the code in the R Markdown file. Setting a seed ensures that any results that rely on randomness, e.g. subsampling or permutations, are reproducible.
Great job! Recording the operating system, R version, and package versions is critical for reproducibility.
Nice! There were no cached chunks for this analysis, so you can be confident that you successfully produced the results during this run.
Great job! Using relative paths to the files within your workflowr project makes it easier to run your code on other machines.
Great! You are using Git for version control. Tracking code development and connecting the code version to the results is critical for reproducibility.
The results in this page were generated with repository version d246c15. See the Past versions tab to see a history of the changes made to the R Markdown and HTML files.
Note that you need to be careful to ensure that all relevant files for the analysis have been committed to Git prior to generating the results (you can use wflow_publish or wflow_git_commit). workflowr only checks the R Markdown file, but you know if there are other scripts or data files that it depends on. Below is the status of the Git repository when the results were generated:
Ignored files:
Ignored: .DS_Store
Ignored: .Rproj.user/
Ignored: Table_S4.csv
Ignored: code/.DS_Store
Ignored: code/._.DS_Store
Ignored: data/.DS_Store
Ignored: data/._.DS_Store
Ignored: data/data_for_analysis/
Ignored: data/full_data/
Unstaged changes:
Modified: .gitignore
Modified: analysis/Supp-Figure_10.rmd
Modified: analysis/_site.yml
Deleted: analysis/license.Rmd
Note that any generated files, e.g. HTML, png, CSS, etc., are not included in this status report because it is ok for generated content to have uncommitted changes.
These are the previous versions of the repository in which changes were made to the R Markdown (analysis/Figure_3.rmd) and HTML (docs/Figure_3.html) files. If you’ve configured a remote Git repository (see ?wflow_git_remote), click on the hyperlinks in the table below to view the files as they were in that past version.
| File | Version | Author | Date | Message |
|---|---|---|---|---|
| html | 73aa800 | toobiwankenobi | 2022-02-22 | add .html for static website |
| Rmd | a4bcb73 | toobiwankenobi | 2022-02-09 | clean repo |
| Rmd | dfe5f09 | toobiwankenobi | 2022-02-09 | change Figure order |
| Rmd | f9a3a83 | toobiwankenobi | 2022-02-08 | clean repo for release |
| Rmd | 588dbb1 | toobiwankenobi | 2022-02-06 | Figure Order |
| Rmd | 3da15db | toobiwankenobi | 2021-11-24 | changes for revision |
| Rmd | c4e2793 | toobiwankenobi | 2021-08-04 | rearrange figure order to match pre-print |
| html | 4109ff1 | toobiwankenobi | 2021-07-07 | delete html files and adapt gitignore |
| Rmd | 0f72ef1 | toobiwankenobi | 2021-05-11 | figure adaptations |
| html | 0f72ef1 | toobiwankenobi | 2021-05-11 | figure adaptations |
| html | 4affda4 | toobiwankenobi | 2021-04-14 | figure adaptations |
| Rmd | 3203891 | toobiwankenobi | 2021-02-19 | change celltype names |
| html | 3203891 | toobiwankenobi | 2021-02-19 | change celltype names |
| Rmd | ee1595d | toobiwankenobi | 2021-02-12 | clean repo and adapt files |
| html | ee1595d | toobiwankenobi | 2021-02-12 | clean repo and adapt files |
| html | 3f5af3f | toobiwankenobi | 2021-02-09 | add .html files |
| Rmd | afa7957 | toobiwankenobi | 2021-02-08 | minor changes on figures and figure order |
| Rmd | 20a1458 | toobiwankenobi | 2021-02-04 | adapt figure order |
| Rmd | f9bb33a | toobiwankenobi | 2021-02-04 | new Figure 5 and minor changes in figure order |
| Rmd | 2ac1833 | toobiwankenobi | 2021-01-08 | changes to Figures |
| Rmd | 545c207 | toobiwankenobi | 2020-12-22 | clean up branch |
| Rmd | 64d1f24 | toobiwankenobi | 2020-12-22 | start new branch with clean scripts |
| Rmd | 9442cb9 | toobiwankenobi | 2020-12-22 | add all new files |
| Rmd | 1af3353 | toobiwankenobi | 2020-10-16 | add stuff |
| Rmd | a6b51cd | toobiwankenobi | 2020-10-14 | clean scripts, add new subfigures |
| Rmd | d8819f2 | toobiwankenobi | 2020-10-08 | read new data (nuclei expansion) and adapt scripts |
| Rmd | 2c11d5c | toobiwankenobi | 2020-08-05 | add new scripts |
This script generates plots for Figure 3. Panel A and B were created manually.
knitr::opts_chunk$set(echo = TRUE, message= FALSE)
knitr::opts_knit$set(root.dir = rprojroot::find_rstudio_root_file())sapply(list.files("code/helper_functions", full.names = TRUE), source) code/helper_functions/calculateSummary.R
value ?
visible FALSE
code/helper_functions/censor_dat.R
value ?
visible FALSE
code/helper_functions/detect_mRNA_expression.R
value ?
visible FALSE
code/helper_functions/DistanceToClusterCenter.R
value ?
visible FALSE
code/helper_functions/findMilieu.R code/helper_functions/findPatch.R
value ? ?
visible FALSE FALSE
code/helper_functions/getInfoFromString.R
value ?
visible FALSE
code/helper_functions/getSpotnumber.R
value ?
visible FALSE
code/helper_functions/plotCellCounts.R
value ?
visible FALSE
code/helper_functions/plotCellFractions.R
value ?
visible FALSE
code/helper_functions/plotDist.R code/helper_functions/read_Data.R
value ? ?
visible FALSE FALSE
code/helper_functions/scatter_function.R
value ?
visible FALSE
code/helper_functions/sceChecks.R
value ?
visible FALSE
code/helper_functions/validityChecks.R
value ?
visible FALSE library(SingleCellExperiment)
library(reshape2)
library(tidyverse)
library(dplyr)
library(data.table)
library(ggplot2)
library(ComplexHeatmap)
library(colorRamps)
library(circlize)
library(RColorBrewer)
library(ggbeeswarm)
library(destiny)
library(scater)
library(dittoSeq)
library(gridExtra)
library(ggpmisc)
library(cowplot)
library(viridis)
library(ggpubr)
library(rstatix)
library(sf)
library(concaveman)
library(RANN)
library(pheatmap)
library(SummarizedExperiment)
library(imcRtools)sce_rna = readRDS(file = "data/data_for_analysis/sce_RNA.rds")
sce_prot = readRDS(file = "data/data_for_analysis/sce_protein.rds")
# meta data
dat_relation = fread(file = "data/data_for_analysis/rna/Object relationships.csv",stringsAsFactors = FALSE)
sce_rna <- sce_rna[,sce_rna$Location != "CTRL"]
sce_prot <- sce_prot[,sce_prot$Location != "CTRL"]example <- findPatch(sce_rna[,sce_rna$ImageNumber == 58], sce_rna[,sce_rna$CXCL10 == 1]$cellID,
'cellID',
'Center_X', 'Center_Y',
'ImageNumber',
distance = 20,
min_clust_size = 10,
output_colname = "example_cluster")Time difference of 1.163448 secs
[1] "patches successfully added to sce object"example <- findMilieu(example,
'cellID',
'Center_X', 'Center_Y',
'ImageNumber',
'example_cluster',
distance = 25,
output_colname = "chemokine_community_i",
plot = TRUE)
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
Time difference of 4.468902 secs
[1] "milieus successfully added to sce object"example <- findPatch(sce_rna[,sce_rna$ImageNumber == 58], sce_rna[,sce_rna$CXCL10 == 1]$cellID,
'cellID',
'Center_X', 'Center_Y',
'ImageNumber',
distance = 20,
min_clust_size = 10,
output_colname = "example_cluster")Time difference of 0.9644794 secs
[1] "patches successfully added to sce object"example <- findMilieu(example,
'cellID',
'Center_X', 'Center_Y',
'ImageNumber',
'example_cluster',
distance = 25,
output_colname = "chemokine_community_i",
plot = TRUE,
xlim = c(725,850),
ylim = c(500,675),
point_size = 14)Warning: Removed 533 rows containing missing values (geom_point).
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
Time difference of 3.071393 secs
[1] "milieus successfully added to sce object"# define fractions of chemokines present in community
cur_dt <- data.frame(colData(sce_rna[,sce_rna$Location != "CTRL"]))
plot_list <- list()
for(i in names(cur_dt[,grepl(glob2rx("*pure"),names(cur_dt))])) {
chemokine_name <- toupper(str_split(i, "_")[[1]][1])
# select all cells that are in a milieu
unique_comms <- unique(cur_dt[cur_dt[,i] > 0,i])
cur_dt_sub <- cur_dt[cur_dt[,i] %in% unique_comms,]
cur_dt_sub <- cbind(cur_dt_sub[,i],
cur_dt_sub[,grepl(glob2rx("C*L*"),names(cur_dt_sub))])
colnames(cur_dt_sub)[1] <- i
# add celltype and MM_location_simplified
cur_dt_sub$cellID <- rownames(cur_dt_sub)
cur_dt_sub <- left_join(cur_dt_sub, cur_dt[,c("cellID", "celltype", "MM_location_simplified")])
# melt the table
cur_dt_sub <- cur_dt_sub %>%
reshape2::melt(id.vars = c("cellID", "celltype", "MM_location_simplified", i), variable.name = "chemokine", value.name = "status")
# all cells that do not produce a chemokine
non_producer <- cur_dt_sub %>%
group_by(cellID) %>%
summarise(sum = sum(status)) %>%
filter(sum == 0) %>%
select(cellID)
# all cells that produce a chemokine - regardless of what chemokine
producer <- cur_dt_sub %>%
group_by(cellID) %>%
summarise(sum = sum(status)) %>%
filter(sum > 0) %>%
select(cellID)
# select non-producing cells and count
non_producer <- cur_dt_sub[cur_dt_sub$cellID %in% non_producer$cellID,] %>%
distinct(cellID, .keep_all = TRUE) %>%
group_by(celltype, MM_location_simplified, chemokine) %>%
summarise(n=n())
non_producer$chemokine <- "no chemokine"
# select producing cells and count chemokines
producer <- cur_dt_sub[cur_dt_sub$cellID %in% producer$cellID,] %>%
filter(status == 1) %>%
group_by(celltype, MM_location_simplified, chemokine) %>%
summarise(n=n())
summary <- rbind(producer, non_producer)
# celltypes numbers
summary_celltypes <- summary %>%
group_by(celltype) %>%
summarise(n=sum(n))
# chemokines per celltype numbers
summary_chemokines <- summary %>%
group_by(chemokine) %>%
summarise(n=sum(n))
# color_vector for cells and chemokines
col_vector_cells <- metadata(sce_rna)$colour_vector$celltype
col_vector_chemokines <- metadata(sce_rna)$colour_vectors$chemokine_single
col_vector <- c(col_vector_cells, col_vector_chemokines)
# add "no chemokine" to col_vector
col_vector <- c(col_vector, "white")
names(col_vector) <- c(names(col_vector[-length(col_vector)]), "no chemokine")
# create labels for middle of sunburst plot
# Number of detected Patches
numberOfPatches <- paste(length(unique_comms), ifelse(length(unique_comms)>1," Milieus", " Milieu"), sep = "")
# Median Number of Chemokine XY Producing Cells in a Patch
medianCells <- cur_dt[cur_dt[,i] > 0 & cur_dt[,chemokine_name] == 1,] %>%
group_by_at(i) %>%
summarise(n=n()) %>%
mutate(median = median(n))
medianCells <- paste(round(unique(medianCells$median)), " Cells", sep = "")
# Percentage of chemokines produced by milieu cells
percentageInPatches <- cur_dt[cur_dt[,chemokine_name] == 1,] %>%
mutate(in_patch = ifelse(.[,i] > 0, 1, 0)) %>%
group_by(in_patch) %>%
summarise(n=n()) %>%
mutate(percentage = n / sum(n) * 100) %>%
filter(in_patch == 1)
percentageInPatches <- paste(round(unique(percentageInPatches$percentage)), "% in Patch", sep = "")
# Number of Patients showing this type of patch
numberPatients <- cur_dt[cur_dt[,i] > 0,] %>%
distinct(PatientID) %>%
summarise(n=n()) %>%
mutate(fraction = n / length(unique(sce_prot[,sce_prot$Location!="CTRL"]$PatientID))*100)
numberPatients <- paste0(round(numberPatients$fraction,0), "% Patients")
CMratio <- cur_dt[cur_dt[,i] %in% unique_comms,] %>%
distinct_at(i, .keep_all = T) %>%
group_by(Location) %>%
summarise(n=n()) %>%
mutate(percentage = n / sum(n) * 100)
percentMargin <- paste(round(CMratio[CMratio$Location == "M",]$percentage), "%", sep = "")
# paste all numbers
label <- paste(paste(paste(numberOfPatches, medianCells, sep = "\n"), numberPatients , sep = "\n"), percentageInPatches, sep = "\n")
# sunburst plot
plt1 <- ggplot() +
geom_col(aes(x = 0, y = percentage, fill = Location),
data = CMratio, col="black") +
theme_void() +
scale_fill_manual(values = c("black", "grey"),
breaks = c("C", "M"),
labels = c("C", "M")) +
theme(axis.ticks=element_blank(),
plot.margin = unit(c(0,0,0,0), "cm"),
axis.text=element_blank(),
axis.title=element_blank(),
legend.position = "none",
text = element_text(size = 10),
plot.title = element_text(hjust = 0.5))
plt2 <- ggplot() +
geom_text(aes(x=0,y=0, label = label), size=5) +
geom_col(aes(x = 8, y = n, fill = celltype),
data = summary_celltypes,
color = "white",
width = 3,
lwd = 0.5) +
xlim(0, 11) + labs(x = NULL, y = NULL) +
scale_fill_manual(values = unname(col_vector),
breaks = names(col_vector),
labels = names(col_vector)) +
theme_void() +
theme(axis.ticks=element_blank(),
plot.margin = unit(c(0,0,0,0), "cm"),
axis.text=element_blank(),
axis.title=element_blank(),
legend.position = "none") +
coord_polar(theta = "y")
plt3 <- ggplot() +
geom_col(aes(x=0,y = n, fill = chemokine),
data = summary_chemokines, col="black") +
theme_void() +
scale_fill_manual(values = unname(col_vector),
breaks = names(col_vector),
labels = names(col_vector)) +
theme(axis.ticks=element_blank(),
plot.margin = unit(c(0,0,0,0), "cm"),
axis.text=element_blank(),
axis.title=element_blank(),
legend.position = "none",
text = element_text(size = 18),
plot.title = element_text(hjust = 0.5))
plt <- grid.arrange(plt1,plt2,plt3, nrow=1, widths=c(0.2,1,0.2),
padding = unit(0,"cm"),
top = textGrob(chemokine_name,gp=gpar(fontsize=15)))
# add to list
plot_list[[i]] <- plot_grid(plt)
}
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
# plot sunburst plots (without CCL4, CCL22, CCL8 - low abundance communities)
plot_grid(plot_list$cxcl8_pure, plot_list$ccl2_pure,
plot_list$cxcl10_pure, plot_list$cxcl9_pure,
plot_list$ccl18_pure, plot_list$ccl19_pure,
plot_list$cxcl12_pure, plot_list$cxcl13_pure,
plot_list$ccl4_pure, plot_list$ccl22_pure,
plot_list$ccl8_pure,
ncol = 2, aligh = "hv")Warning in as_grob.default(plot): Cannot convert object of class character into
a grob.
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
lgd1 = Legend(labels = c("C", "M"), title = "Core/Margin",
legend_gp = gpar(fill = c("black", "grey")))
# create legend for celltypes
lgd2 = Legend(labels = names(col_vector_cells), title = "Cell Type ", legend_gp = gpar(fill = unname(col_vector_cells)))
# create legend for chemokines
lgd3 = Legend(labels = names(col_vector_chemokines), title = "Chemokines", legend_gp = gpar(fill = unname(col_vector_chemokines)))
draw(packLegend(lgd1, lgd2, lgd3, direction = "horizontal"))
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
milieus <- data.frame(colData(sce_rna)) %>%
filter(celltype == "CD8+ T cell") %>%
select(cellID, contains("pure")) %>%
mutate_if(is.numeric, ~1 * (. > 0))
milieus$number_of_milieus <- rowSums(milieus[,-1])
# keep CD8+ T cells that are part of at least one milieu
milieus <- milieus %>%
filter(number_of_milieus > 0) %>%
select(-number_of_milieus) %>%
reshape2::melt(id.vars = "cellID", variable.name = "milieu", value.name = "is_part") %>%
filter(is_part > 0) %>%
select(cellID, milieu)marker_rna <- c("Lag3", "T8_CXCL13", "T5_CCL4")
# rna data
dat_rna <- data.frame(t(assay(sce_rna[marker_rna, sce_rna$celltype == "CD8+ T cell"], "asinh")))
dat_rna$cellID <- rownames(dat_rna)
dat_rna <- left_join(milieus, dat_rna)
# melt
dat_rna <- dat_rna %>%
reshape2::melt(id.vars = c("cellID", "milieu"), variable.name = "channel", value.name = "asinh")
# remove CCL4/CCL8/CXCL8 milieus due to too few data points
dat_rna <- dat_rna %>%
filter(!(milieu %in% c("ccl4_pure", "ccl8_pure", "cxcl8_pure")))
# rename milieus
dat_rna <- dat_rna %>%
mutate(milieu_short = toupper(str_split(milieu, "_", n = 2, simplify = TRUE)[,1]))
col_vector_chemokines <- metadata(sce_rna)$colour_vectors$chemokine_single
# add channel medium
dat_rna <- dat_rna %>%
group_by(channel) %>%
mutate(channel_median = median(asinh))
# one-sample t test
stat.test <- data.frame()
# loop through all channels (each has a different µ)
for(j in unique(dat_rna$channel)){
cur.mu <- unique(dat_rna[dat_rna$channel == j, ]$channel_median)
# calculate p-value for different milieus in one channel and adjust pvalue
cur.test <- dat_rna[dat_rna$channel == j, ] %>%
group_by(channel) %>%
wilcox_test(asinh ~ milieu_short, ref.group = ".all.") %>%
adjust_pvalue(method = "BH") %>%
add_x_position(x="milieu_short")
stat.test <- rbind(stat.test, cur.test)
}
# adjust again for testing across different channels
stat.test <- stat.test %>%
group_by(channel) %>%
adjust_pvalue(method = "BH") %>%
add_significance("p.adj",cutpoints = c(0, 1e-04, 0.001, 0.01, 0.1, 1))
# plot
plot_list <- list()
ylim_list <- list("Lag3" = c(0,0.9), "T8_CXCL13" = c(0,2.3), "T5_CCL4" = c(0,1))
for(i in unique(dat_rna$channel)){
cur.stat.test <- stat.test[stat.test$channel == i, ]
plot_list[[i]] <- ggplot(dat_rna[dat_rna$channel == i,], aes(x=milieu_short, y=asinh)) +
geom_boxplot(alpha=1, lwd=0.5, outlier.shape = NA, position = position_dodge(1.1), aes(fill=milieu_short)) +
theme_bw() +
theme(text = element_text(size=18),
axis.title.x=element_blank(),
axis.text.x=element_blank(),
axis.ticks.x=element_blank()) +
guides(fill=guide_legend("Milieu", override.aes = list(alpha=1)), col="none") +
stat_pvalue_manual(
cur.stat.test, x = "xmin", y.position = ylim_list[[i]][2]-0.05,
label = "p.adj.signif",
position = position_dodge(0.8),
size=2.5) +
ylab("Mean Expression (asinh)") +
xlab("") +
geom_hline(aes(yintercept = channel_median, group = channel), colour = 'black', linetype = 2, size=0.5) +
scale_fill_manual(values = unname(col_vector_chemokines),
breaks = names(col_vector_chemokines),
labels = names(col_vector_chemokines)) +
coord_cartesian(ylim=ylim_list[[i]]) +
facet_wrap(~channel)
}
leg_c <- cowplot::get_legend(plot_list[[1]])
grid.arrange(plot_list[[1]] + theme(legend.position = "none"),
plot_list[[2]] + theme(legend.position = "none") + ylab(""),
plot_list[[3]] + theme(legend.position = "none") + ylab(""),
ncol=2)
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
Here, we will perform an enrichment analysis checking if certain cell types are enriched or avoid in chemokine milieus. We will test the enrichment of each cell type within this community type using a Fisher’s exact test. We will also exclude cells that already express the cytokine.
cur_res_list <- list()
cur_df <- colData(sce_rna)
# Set thresholds
perc_cells <- 0
# We remove ccl8, ccl4 and ccl22 due to few images with milieus
for (i in c("ccl18_pure", "cxcl8_pure", "cxcl10_pure",
"cxcl12_pure", "cxcl13_pure", "ccl2_pure",
"cxcl9_pure", "ccl19_pure")) {
cur_perc <- cur_df %>%
as.data.frame() %>%
group_by(ImageNumber) %>%
dplyr::summarize(perc_cells = sum(!!sym(i) != 0)/n())
cur_chemo <- toupper(sub("_pure", "", i))
for (j in unique(cur_df$celltype)) {
cur_res <- cur_df %>%
as.data.frame() %>%
filter(ImageNumber %in% cur_perc$ImageNumber[cur_perc$perc_cells > perc_cells]) %>%
filter(!(!!sym(cur_chemo) == 1 & celltype == j)) %>%
group_by(ImageNumber) %>%
dplyr::summarize(celltype_inside = sum(celltype == j & !!sym(i) != 0),
other_inside = sum(celltype != j & !!sym(i) != 0),
celltype_outside = sum(celltype == j & !!sym(i) == 0),
other_outside = sum(celltype != j & !!sym(i) == 0))
if (nrow(cur_res) == 0) {
next
}
cur_tests <- as.data.frame(t(apply(as.matrix(cur_res), 1, function(x){
cur_mat <- matrix(x[2:5], ncol = 2, nrow = 2, byrow = FALSE)
rownames(cur_mat) <- c("celltype", "other")
colnames(cur_mat) <- c("inside", "outside")
cur_test <- fisher.test(cur_mat)
c(cur_test$p.value, cur_test$estimate)
})))
colnames(cur_tests)[1] <- "p.value"
cur_tests$adj.p <- p.adjust(cur_tests$p.value, method = "BH")
cur_tests$community <- i
cur_tests$celltype <- j
cur_tests$ImageNumber <- cur_res$ImageNumber
cur_res_list[[paste0(i, "_", j)]] <- cur_tests
}
}
out <- do.call("rbind", cur_res_list)
out$adj.p.all <- p.adjust(out$p.value, method = "BH")
final <- out %>% mutate(sigval = ifelse(`odds ratio` > 1 & adj.p.all < 0.1, 1,
ifelse(`odds ratio` < 1 & adj.p.all < 0.1, -1, 0))) %>%
group_by(community, celltype) %>%
dplyr::summarize(mean = mean(sigval))
# Number of tested images
out %>%
group_by(community, celltype) %>%
dplyr::summarize(count = n()) %>%
as.data.frame() community celltype count
1 ccl18_pure CD38 14
2 ccl18_pure CD8- T cell 14
3 ccl18_pure CD8+ T cell 14
4 ccl18_pure HLA-DR 14
5 ccl18_pure Macrophage 14
6 ccl18_pure Neutrophil 14
7 ccl18_pure Stroma 14
8 ccl18_pure Tumor 14
9 ccl18_pure unknown 14
10 ccl18_pure Vasculature 14
11 ccl19_pure CD38 19
12 ccl19_pure CD8- T cell 19
13 ccl19_pure CD8+ T cell 19
14 ccl19_pure HLA-DR 19
15 ccl19_pure Macrophage 19
16 ccl19_pure Neutrophil 19
17 ccl19_pure Stroma 19
18 ccl19_pure Tumor 19
19 ccl19_pure unknown 19
20 ccl19_pure Vasculature 19
21 ccl2_pure CD38 15
22 ccl2_pure CD8- T cell 15
23 ccl2_pure CD8+ T cell 15
24 ccl2_pure HLA-DR 15
25 ccl2_pure Macrophage 15
26 ccl2_pure Neutrophil 15
27 ccl2_pure Stroma 15
28 ccl2_pure Tumor 15
29 ccl2_pure unknown 15
30 ccl2_pure Vasculature 15
31 cxcl10_pure CD38 43
32 cxcl10_pure CD8- T cell 43
33 cxcl10_pure CD8+ T cell 43
34 cxcl10_pure HLA-DR 43
35 cxcl10_pure Macrophage 43
36 cxcl10_pure Neutrophil 43
37 cxcl10_pure Stroma 43
38 cxcl10_pure Tumor 43
39 cxcl10_pure unknown 43
40 cxcl10_pure Vasculature 43
41 cxcl12_pure CD38 10
42 cxcl12_pure CD8- T cell 10
43 cxcl12_pure CD8+ T cell 10
44 cxcl12_pure HLA-DR 10
45 cxcl12_pure Macrophage 10
46 cxcl12_pure Neutrophil 10
47 cxcl12_pure Stroma 10
48 cxcl12_pure Tumor 10
49 cxcl12_pure unknown 10
50 cxcl12_pure Vasculature 10
51 cxcl13_pure CD38 23
52 cxcl13_pure CD8- T cell 23
53 cxcl13_pure CD8+ T cell 23
54 cxcl13_pure HLA-DR 23
55 cxcl13_pure Macrophage 23
56 cxcl13_pure Neutrophil 23
57 cxcl13_pure Stroma 23
58 cxcl13_pure Tumor 23
59 cxcl13_pure unknown 23
60 cxcl13_pure Vasculature 23
61 cxcl8_pure CD38 8
62 cxcl8_pure CD8- T cell 8
63 cxcl8_pure CD8+ T cell 8
64 cxcl8_pure HLA-DR 8
65 cxcl8_pure Macrophage 8
66 cxcl8_pure Neutrophil 8
67 cxcl8_pure Stroma 8
68 cxcl8_pure Tumor 8
69 cxcl8_pure unknown 8
70 cxcl8_pure Vasculature 8
71 cxcl9_pure CD38 37
72 cxcl9_pure CD8- T cell 37
73 cxcl9_pure CD8+ T cell 37
74 cxcl9_pure HLA-DR 37
75 cxcl9_pure Macrophage 37
76 cxcl9_pure Neutrophil 37
77 cxcl9_pure Stroma 37
78 cxcl9_pure Tumor 37
79 cxcl9_pure unknown 37
80 cxcl9_pure Vasculature 37for_plot <- final %>%
pivot_wider(names_from = community, values_from = mean) %>%
as.data.frame()
rownames(for_plot) <- for_plot$celltype
for_plot <- for_plot[,-1]pheatmap(as.matrix(for_plot),
color = colorRampPalette(c("dark blue", "white", "dark red"))(100),
breaks = seq(-1, 1, length.out = 100))
| Version | Author | Date |
|---|---|---|
| 3697a9b | toobiwankenobi | 2022-02-22 |
sessionInfo()R version 4.1.2 (2021-11-01)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: Ubuntu 20.04.3 LTS
Matrix products: default
BLAS/LAPACK: /usr/lib/x86_64-linux-gnu/openblas-pthread/libopenblasp-r0.3.8.so
locale:
[1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
[3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
[5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
[7] LC_PAPER=en_US.UTF-8 LC_NAME=C
[9] LC_ADDRESS=C LC_TELEPHONE=C
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
attached base packages:
[1] grid stats4 stats graphics grDevices utils datasets
[8] methods base
other attached packages:
[1] imcRtools_1.0.2 SpatialExperiment_1.4.0
[3] pheatmap_1.0.12 RANN_2.6.1
[5] concaveman_1.1.0 sf_1.0-5
[7] rstatix_0.7.0 ggpubr_0.4.0
[9] viridis_0.6.2 viridisLite_0.4.0
[11] cowplot_1.1.1 ggpmisc_0.4.5
[13] ggpp_0.4.3 gridExtra_2.3
[15] dittoSeq_1.6.0 scater_1.22.0
[17] scuttle_1.4.0 destiny_3.8.1
[19] ggbeeswarm_0.6.0 RColorBrewer_1.1-2
[21] circlize_0.4.13 colorRamps_2.3
[23] ComplexHeatmap_2.10.0 data.table_1.14.2
[25] forcats_0.5.1 stringr_1.4.0
[27] purrr_0.3.4 readr_2.1.2
[29] tidyr_1.2.0 tibble_3.1.6
[31] ggplot2_3.3.5 tidyverse_1.3.1
[33] reshape2_1.4.4 SingleCellExperiment_1.16.0
[35] SummarizedExperiment_1.24.0 Biobase_2.54.0
[37] GenomicRanges_1.46.1 GenomeInfoDb_1.30.1
[39] IRanges_2.28.0 S4Vectors_0.32.3
[41] BiocGenerics_0.40.0 MatrixGenerics_1.6.0
[43] matrixStats_0.61.0 dplyr_1.0.7
[45] workflowr_1.7.0
loaded via a namespace (and not attached):
[1] SparseM_1.81 ggthemes_4.2.4
[3] R.methodsS3_1.8.1 bit64_4.0.5
[5] knitr_1.37 irlba_2.3.5
[7] DelayedArray_0.20.0 R.utils_2.11.0
[9] RCurl_1.98-1.5 doParallel_1.0.16
[11] generics_0.1.2 ScaledMatrix_1.2.0
[13] terra_1.5-17 callr_3.7.0
[15] proxy_0.4-26 bit_4.0.4
[17] tzdb_0.2.0 xml2_1.3.3
[19] lubridate_1.8.0 httpuv_1.6.5
[21] assertthat_0.2.1 xfun_0.29
[23] hms_1.1.1 jquerylib_0.1.4
[25] evaluate_0.14 promises_1.2.0.1
[27] DEoptimR_1.0-10 fansi_1.0.2
[29] dbplyr_2.1.1 readxl_1.3.1
[31] htmlwidgets_1.5.4 igraph_1.2.11
[33] DBI_1.1.2 ellipsis_0.3.2
[35] RSpectra_0.16-0 backports_1.4.1
[37] V8_4.0.0 svgPanZoom_0.3.4
[39] sparseMatrixStats_1.6.0 vctrs_0.3.8
[41] quantreg_5.87 TTR_0.24.3
[43] abind_1.4-5 RcppEigen_0.3.3.9.1
[45] withr_2.4.3 ggforce_0.3.3
[47] cytomapper_1.6.0 robustbase_0.93-9
[49] vroom_1.5.7 vcd_1.4-9
[51] xts_0.12.1 svglite_2.0.0
[53] cluster_2.1.2 laeken_0.5.2
[55] crayon_1.4.2 labeling_0.4.2
[57] edgeR_3.36.0 pkgconfig_2.0.3
[59] units_0.7-2 tweenr_1.0.2
[61] vipor_0.4.5 nnet_7.3-17
[63] rlang_1.0.0 lifecycle_1.0.1
[65] MatrixModels_0.5-0 modelr_0.1.8
[67] rsvd_1.0.5 polyclip_1.10-0
[69] cellranger_1.1.0 rprojroot_2.0.2
[71] RcppHNSW_0.3.0 lmtest_0.9-39
[73] tiff_0.1-11 Matrix_1.4-0
[75] raster_3.5-15 carData_3.0-5
[77] Rhdf5lib_1.16.0 boot_1.3-28
[79] zoo_1.8-9 RTriangle_1.6-0.10
[81] reprex_2.0.1 beeswarm_0.4.0
[83] whisker_0.4 ggridges_0.5.3
[85] GlobalOptions_0.1.2 processx_3.5.2
[87] png_0.1-7 rjson_0.2.21
[89] shinydashboard_0.7.2 bitops_1.0-7
[91] getPass_0.2-2 R.oo_1.24.0
[93] KernSmooth_2.23-20 rhdf5filters_1.6.0
[95] DelayedMatrixStats_1.16.0 shape_1.4.6
[97] classInt_0.4-3 jpeg_0.1-9
[99] ggsignif_0.6.3 beachmat_2.10.0
[101] scales_1.1.1 magrittr_2.0.2
[103] plyr_1.8.6 hexbin_1.28.2
[105] zlibbioc_1.40.0 compiler_4.1.2
[107] dqrng_0.3.0 pcaMethods_1.86.0
[109] clue_0.3-60 cli_3.1.1
[111] XVector_0.34.0 ps_1.6.0
[113] ggplot.multistats_1.0.0 MASS_7.3-55
[115] tidyselect_1.1.1 stringi_1.7.6
[117] highr_0.9 yaml_2.2.2
[119] BiocSingular_1.10.0 locfit_1.5-9.4
[121] ggrepel_0.9.1 sass_0.4.0
[123] EBImage_4.36.0 tools_4.1.2
[125] parallel_4.1.2 rstudioapi_0.13
[127] foreach_1.5.2 git2r_0.29.0
[129] smoother_1.1 farver_2.1.0
[131] scatterplot3d_0.3-41 ggraph_2.0.5
[133] DropletUtils_1.14.2 digest_0.6.29
[135] shiny_1.7.1 Rcpp_1.0.8
[137] car_3.0-12 broom_0.7.12
[139] later_1.3.0 httr_1.4.2
[141] colorspace_2.0-2 rvest_1.0.2
[143] fs_1.5.2 ranger_0.13.1
[145] graphlayouts_0.8.0 sp_1.4-6
[147] systemfonts_1.0.3 xtable_1.8-4
[149] jsonlite_1.7.3 tidygraph_1.2.0
[151] R6_2.5.1 mime_0.12
[153] pillar_1.7.0 htmltools_0.5.2
[155] DT_0.20 glue_1.6.1
[157] fastmap_1.1.0 VIM_6.1.1
[159] BiocParallel_1.28.3 fftwtools_0.9-11
[161] BiocNeighbors_1.12.0 class_7.3-20
[163] codetools_0.2-18 utf8_1.2.2
[165] lattice_0.20-45 bslib_0.3.1
[167] curl_4.3.2 magick_2.7.3
[169] limma_3.50.0 rmarkdown_2.11
[171] munsell_0.5.0 e1071_1.7-9
[173] GetoptLong_1.0.5 rhdf5_2.38.0
[175] GenomeInfoDbData_1.2.7 iterators_1.0.13
[177] HDF5Array_1.22.1 haven_2.4.3
[179] gtable_0.3.0