Figure 5 (A)
# Load the library
library(eulerr)
library(grid)
# Define functions
compute_venn <- function(nv_genes, sf_genes, sample1 = "Lung NV", sample2 = "Lung SF", digits = 4) {
set1_only <- setdiff(nv_genes, sf_genes)
set2_only <- setdiff(sf_genes, nv_genes)
shared <- intersect(nv_genes, sf_genes)
counts <- c(length(set1_only), length(set2_only), length(shared))
names(counts) <- c(sample1, sample2, paste(sample1, sample2, sep = "&"))
total <- sum(counts)
labels <- paste0(counts, "\n(", signif(counts / total * 100, digits), "%)")
names(labels) <- names(counts)
list(counts = counts, labels = labels)
}
plot_marker_euler <- function(venn_res, title = "Sample", fills = c("#F9DF91", "#A8CFE8", "#E6D9AE")) {
venn_obj <- euler(venn_res$counts)
grid.newpage()
p <- plot(
venn_obj,
fills = list(fill = fills, alpha = 0.7),
labels = FALSE,
edges = NULL,
quantities = FALSE,
main = list(label = title, cex = 2.5, font = 2),
output = "grob"
)
grid.draw(p)
grid.text(venn_res$labels[1], x = 0.14, y = 0.5, gp = gpar(fontsize = 20, fontface = "bold"))
grid.text(venn_res$labels[2], x = 0.86, y = 0.5, gp = gpar(fontsize = 20, fontface = "bold"))
grid.text(venn_res$labels[3], x = 0.5, y = 0.5, gp = gpar(fontsize = 25, fontface = "bold"))
}
# Load data and show the figure
load("markers_lnvsf.RData")
celltype <- "Mesenchymal_Alveolar"
nv_genes <- markers_l_nv$gene[markers_l_nv$cluster == celltype]
sf_genes <- markers_l_sf$gene[markers_l_sf$cluster == celltype]
venn_res <- compute_venn(nv_genes, sf_genes)
plot_marker_euler(venn_res, title = celltype)
Figure 5 (B)
# Load the library
library(Seurat)
library(ggplot2)
library(org.Mm.eg.db)
library(dplyr)
library(clusterProfiler)
# Define functions
perform_KEGG <- function(marker_data, p_cutoff = 0.05, q_cutoff = 0.2) {
markers <- marker_data %>%
filter(p_val_adj < 0.001) %>%
group_by(cluster) %>%
ungroup()
gid <- bitr(unique(markers$gene), fromType = "SYMBOL", toType = "ENTREZID", OrgDb = org.Mm.eg.db)
markers <- full_join(markers, gid, by = c("gene" = "SYMBOL")) %>%
filter(!is.na(ENTREZID))
compareCluster(
ENTREZID ~ cluster,
data = markers,
fun = "enrichKEGG",
organism = "mmu",
keyType = "ncbi-geneid",
pvalueCutoff = p_cutoff,
qvalueCutoff = q_cutoff
)
}
clean_go_description <- function(df) {
df$Description <- gsub(" - Mus musculus \\(house mouse\\)", "", df$Description)
df
}
combine_results <- function(go_obj_list, sample_names) {
combined <- lapply(seq_along(go_obj_list), function(i) {
df <- clean_go_description(go_obj_list[[i]]@compareClusterResult)
df$sample <- sample_names[i]
df
}) %>% bind_rows()
combined_go <- go_obj_list[[1]]
combined_go@compareClusterResult <- combined
combined_go
}
# Load the data and combine
lnv_KEGG <- perform_KEGG(markers_l_nv, p_cutoff = 0.05)
lsf_KEGG <- perform_KEGG(markers_l_sf, p_cutoff = 0.05)
l_combined_KEGG <- combine_results(
go_obj_list = list(lnv_KEGG, lsf_KEGG),
sample_names = c("Lung NV", "Lung SF")
)
l_combined_KEGG@compareClusterResult$cluster <- factor(
l_combined_KEGG@compareClusterResult$cluster,
levels = unique(l_combined_KEGG@compareClusterResult$cluster)
)
# Show the figure
p <- dotplot(
l_combined_KEGG,
x = "cluster",
color = "p.adjust",
font.size = 18,
showCategory = 5,
label_format = 40
) +
theme(
axis.text.x = element_text(angle = 30, hjust = 1, size = 16, face = "bold"),
axis.text.y = element_text(size = 16, face = "bold"),
axis.title = element_blank(),
plot.title = element_text(size = 32, face = "bold"),
legend.title = element_text(size = 20, face = "bold"),
legend.text = element_text(size = 20, face = "bold"),
strip.text = element_text(size = 22, face = "bold")
) +
facet_grid(. ~ sample)
p
Figure 5 (D)
import warnings
warnings.filterwarnings('ignore')
import stereo as st
import pickle
import matplotlib.pyplot as plt
# Define functions
class SafeUnpickler(pickle.Unpickler):
def find_class(self, module, name):
if module == "stereo.core.result" and name == "Result":
class DummyResult(dict):
def __setitem__(self, key, value):
super().__setitem__(key, value)
return DummyResult
return super().find_class(module, name)
def load_stereo_data(file_path):
with open(file_path, "rb") as f:
return SafeUnpickler(f).load()
def subset_cells(data, platform):
cell_names = data.cells.cell_name
platforms = [n.split("_")[0].split(" ")[1] for n in cell_names]
cells = [c for c, p in zip(cell_names, platforms) if p == platform]
sub_data = data.tl.filter_cells(cell_list=cells, inplace=False)
sub_module = data.tl.result["spatial_hotspot"].module_scores.loc[cells, :]
sub_data.tl.result["spatial_hotspot"].module_scores = sub_module
return sub_data
def plot_module(data, platform, module_id, vmin=-2, vmax=3.5, cbar_ticks=None, title=None):
sub_data = subset_cells(data, platform)
scores = sub_data.tl.result['spatial_hotspot'].module_scores[module_id]
coords = sub_data.cells_matrix['spatial']
if cbar_ticks is None:
cbar_ticks = [vmin, 0, vmax]
plt.figure(figsize=(6, 5))
sc = plt.scatter(coords[:, 0], coords[:, 1],
c=scores, cmap="RdYlBu_r", s=5,
vmin=vmin, vmax=vmax)
ax = plt.gca()
ax.invert_yaxis()
ax.set_aspect('equal')
ax.set_xticks([]); ax.set_yticks([])
for spine in ax.spines.values():
spine.set_visible(False)
cbar = plt.colorbar(sc, shrink=0.4, aspect=17)
cbar.set_ticks(cbar_ticks)
plt.axis("off")
plt.show()
# Load the data and plot figure
data = load_stereo_data("Lung_NVSF.pkl")
plot_module(data, platform="NV", module_id=1,
vmin=-2, vmax=3.5, cbar_ticks=[-1, 0, 1, 2, 3])
# Show the figure
plot_module(data, platform="SF", module_id=1,
vmin=-2, vmax=3.5, cbar_ticks=[-1, 0, 1, 2, 3])
