{"id":27993,"date":"2024-02-04T17:23:06","date_gmt":"2024-02-04T09:23:06","guid":{"rendered":"http:\/\/www.biocloudservice.com\/wordpress\/?p=27993"},"modified":"2024-02-04T17:23:08","modified_gmt":"2024-02-04T09:23:08","slug":"%e6%83%b3%e5%b1%95%e7%a4%ba%e7%8b%ac%e7%89%b9%e7%9a%84%e5%af%8c%e9%9b%86%e7%bb%93%e6%9e%9c%ef%bc%9f%e8%bf%99%e4%bb%bdgo%e5%88%86%e6%9e%90%e5%86%b2%e7%a7%af%e5%9b%be%e7%9a%84%e4%bd%bf%e7%94%a8","status":"publish","type":"post","link":"http:\/\/www.biocloudservice.com\/wordpress\/?p=27993","title":{"rendered":"\u60f3\u5c55\u793a\u72ec\u7279\u7684\u5bcc\u96c6\u7ed3\u679c\uff1f\u8fd9\u4efdGO\u5206\u6790+\u51b2\u79ef\u56fe\u7684\u4f7f\u7528\u6559\u7a0b\u8bf7\u6536\u597d\uff01"},"content":{"rendered":"

\u6211\u4eec\u90fd\u77e5\u9053\u5bcc\u96c6\u5206\u6790\u662f\u4e00\u79cd\u5e7f\u6cdb\u5e94\u7528\u4e8e\u751f\u7269\u4fe1\u606f\u5b66\u7814\u7a76\u7684\u7edf\u8ba1\u65b9\u6cd5\uff0c\u4e3b\u8981\u7528\u4e8e\u68c0\u9a8c\u4e00\u4e2a\u57fa\u56e0\u96c6\u5408\u4e2d\u67d0\u4e9b\u529f\u80fd\u6216\u7279\u5f81\u7684\u5bcc\u96c6\u7a0b\u5ea6\uff0c\u8fdb\u800c\u5e2e\u52a9\u6211\u4eec\u89e3\u8bfb\u4e00\u7ec4\u57fa\u56e0\u80cc\u540e\u6240\u4ee3\u8868\u7684\u751f\u7269\u5b66\u77e5\u8bc6\uff0c\u5e76\u63ed\u793a\u5176\u5728\u7ec6\u80de\u5185\u6216\u7ec6\u80de\u5916\u626e\u6f14\u4e86\u4ec0\u4e48\u6837\u7684\u89d2\u8272\u3002\u5f53\u6211\u4eec\u5c55\u793a\u5bcc\u96c6\u5206\u6790\u7684\u7ed3\u679c\u7684\u65f6\u5019\uff0c\u6709\u591a\u79cd\u7684\u8868\u793a\u65b9\u6cd5\uff0c\u5c0f\u679c\u5076\u5c14\u4e5f\u60f3\u7ed8\u5236\u4e00\u4e9b\u65b0\u5947\u7cbe\u5f69\u7684\u56fe\u7247\uff0c\u6240\u4ee5\u7279\u5730\u5b66\u4e60\u4e86<\/a>ggalluvial\u8fd9\u4e2a\u5305\uff0c\u4e0d\u4f46\u53ef\u4ee5\u5c55\u793a\u5bcc\u96c6\u5f97\u5230\u7684\u57fa\u56e0\uff0c\u8fd8\u53ef\u4ee5\u5c06\u914d\u4f53\u4e0e\u53d7\u4f53\u7684\u5173\u7cfb\u4e5f\u4e00\u5e76\u7ed8\u5236\u5728\u56fe\u4e0a\uff01\u4eca\u5929\u7684\u6570\u636e\u5904\u7406\u91cf\u6bd4\u8f83\u5927\uff0c\u5982\u679c\u5c0f\u4f19\u4f34\u4eec\u81ea\u5df1\u7684\u670d\u52a1\u5668\u8fd0\u884c\u7a7a\u95f4\u6709\u9650\uff0c\u6b22\u8fce\u8ba2\u9605\u6211\u4eec\u7684\u670d\u52a1\u5668\u54e6\uff5e<\/p>\n

\u5728\u6559\u5927\u5bb6\u5982\u4f55\u4f7f\u7528ggalluvial\u7ed8\u5236\u51b2\u79ef\u56fe\u4e4b\u524d\uff0c\u5c0f\u679c\u9700\u8981\u5148\u7ed9\u5c0f\u4f19\u4f34\u4eec\u4ecb\u7ecd\u4e00\u4e0b\uff0c\u4ec0\u4e48\u662fGO\u5bcc\u96c6\u5206\u6790\uff08Gene Ontology Enrichment Analysis\uff09\u3002\u5f53\u6211\u4eec\u83b7\u5f97\u4e00\u6279\u8f6c\u5f55\u7ec4\u7684\u6570\u636e\u540e\uff0c\u5fc5\u7136\u8981\u505a\u7684\u662f\u4e0d\u540c\u6837\u672c\u4e4b\u95f4\u7684\u5dee\u5f02\u57fa\u56e0\u5206\u6790\uff0c\u5728\u5f97\u5230\u5dee\u5f02\u7684\u57fa\u56e0\u6216\u86cb\u767d\u5217\u8868\u540e\uff0c\u8fd8\u4f1a\u9047\u5230\u4e00\u4e2a\u95ee\u9898\uff0c\u90a3\u5c31\u662f\u5dee\u5f02\u57fa\u56e0\u662f\u975e\u5e38\u591a\u7684\uff0c\u9700\u8981\u6709\u4e00\u79cd\u65b9\u6cd5\u53ef\u4ee5\u6279\u91cf\u5904\u7406\u3002<\/p>\n

\u6b64\u65f6\u6211\u4eec\u5c31\u9700\u8981\u628a\u8fd9\u4e9b\u5dee\u5f02\u57fa\u56e0\u8fdb\u884c\u6ce8\u91ca\uff0c\u628a\u8fd9\u4e9b\u57fa\u56e0\u6216\u86cb\u767d\u5206\u6210\u51e0\u5927\u7c7b\uff08\u4e00\u4e2a\u7c7b\u522b\u5c31\u76f8\u5f53\u4e8e\u4e00\u4e2aGO term\uff09\uff0c\u6b64\u65f6\u53ea\u8981\u5224\u65ad\u8fd9\u51e0\u4e2aGO term\u7684\u533a\u522b\uff0c\u5c31\u53ef\u4ee5\u66f4\u5bb9\u6613\u5730\u5bf9\u5927\u6279\u91cf\u7684\u6570\u636e\u8fdb\u884c\u5206\u6790\u4e86\uff0c\u800c\u8fd9\u4e2a\u8fc7\u7a0b\u5c31\u662f\u5bcc\u96c6\u5206\u6790\u3002\u5bcc\u96c6\u5206\u6790\u5c5e\u4e8e\u5dee\u5f02\u57fa\u56e0\u7684\u4e0b\u6e38\u5206\u6790\uff0c\u5e38\u5e38\u6d89\u53ca\u5230\u4e24\u4e2a\u6982\u5ff5\uff0c\u524d\u666f\u57fa\u56e0\u548c\u80cc\u666f\u57fa\u56e0\uff0c\u524d\u666f\u57fa\u56e0\u5c31\u662f\u6211\u4eec\u5173\u6ce8\u7684\u91cd\u70b9\u7814\u7a76\u7684\u57fa\u56e0\u96c6\uff0c\u80cc\u666f\u57fa\u56e0\u5c31\u662f\u6240\u6709\u7684\u57fa\u56e0\u96c6\u3002<\/p>\n

\u5728\u4e86\u89e3\u4e86\u5bcc\u96c6\u5206\u6790\u539f\u7406\u7684\u57fa\u7840\u4e0a\uff0c\u6211\u4eec\u5c31\u53ef\u4ee5\u8fdb\u884cGO\u5bcc\u96c6\u5206\u6790\uff0c\u5e76\u4f7f\u7528\u51b2\u79ef\u56fe\u5c55\u793a\u5176\u72ec\u7279\u7684\u7ed3\u679c\u5566\uff01<\/p>\n

\u7b2c\u4e00\u6b65\uff1a\u5bf9\u6240\u6709\u80cc\u666f\u57fa\u56e0\u7684DEG\uff08\u5dee\u5f02\u8868\u8fbe\u57fa\u56e0\uff09\u8fdb\u884cGO\u5206\u6790<\/p>\n

\u5728\u5df2\u7ecf\u8fa8\u660e\u80bf\u7624\u7c07\u4e0e\u57fa\u8d28\u7c07\u7684\u57fa\u7840\u4e0a\uff0c\u4f9d\u636e\u5c0f\u679c\u63d0\u4f9b\u7684\u4ee3\u7801\uff0c\u5bfc\u5165dplyr\uff0cstats\uff0cggplot2\uff0cgplots\uff0cRColorBrewer\u5305\u3002\u9996\u5148\u8bbe\u5b9aGO\u57fa\u56e0\u96c6\u4e3a\u9ad8\u7b49\u7ea7\u7684GO terms\uff0c\u63a5\u7740\u8fc7\u6ee4\u5c11\u4e8e10\u4e2a\u57fa\u56e0\u7684\u57fa\u56e0\u96c6\uff0c\u7b5b\u9009\u5dee\u5f02\u8868\u8fbe\u57fa\u56e0\uff0c\u5e76\u5bf9\u4e00\u4e2a\u57fa\u56e0\u96c6\u5408\u8fdb\u884cGO\u5bcc\u96c6\u5206\u6790\uff0c\u4f7f\u7528\u8d85\u51e0\u4f55\u5206\u5e03\u68c0\u9a8c\u6765\u8ba1\u7b97\u6bcf\u4e2aGO\u672f\u8bed\u5728\u6bcf\u4e2a\u7ec6\u80de\u7fa4\u4e2d\u7684\u5bcc\u96c6\u7a0b\u5ea6\uff0c\u7136\u540e\u4f7f\u7528FDR\uff08false discovery rate\uff09\u65b9\u6cd5\u6765\u6821\u6b63p\u503c\uff0c\u5e76\u8fc7\u6ee4\u6389\u6240\u6709FDR<0.01\u7684\u57fa\u56e0\uff0c\u6700\u540e\u627e\u5230\u6bcf\u7c07\u7684\u8868\u8fbe\u6700\u9ad8\u768420\u4e2a\u57fa\u56e0\u7ec4\uff0c\u7ed8\u5236\u70ed\u56fe\uff0c\u5c06\u4ece\u6570\u636e\u4e2d\u9009\u62e9\u7279\u5b9a\u7684\u751f\u7269\u8fc7\u7a0b\u8fd0\u7528\u5b8c\u5168\u94fe\u63a5\u8fdb\u884c\u805a\u7c7b\u3002<\/p>\n

<\/p>\n

\u5982\u56fe\u662f\u6211\u4eec\u7684\u8f93\u5165\u6570\u636e\uff0c\u5305\u542b\u4e86GO\u7684\u5206\u7c7b\u96c6\u3002<\/p>\n

# \u5bfc\u5165\u5305<\/p>\n

library(dplyr)<\/p>\n

library(stats)<\/p>\n

library(ggplot2)<\/p>\n

library(gplots)<\/p>\n

library(RColorBrewer)<\/p>\n

data_dir = ‘\/data’<\/p>\n

save_dir = ‘\/results’<\/p>\n

fnm = paste0(save_dir,<\/p>\n

‘\/seurat_object_TM_SM.RData’)<\/p>\n

load(fnm)<\/p>\n

## \u8bbe\u5b9aGO\u57fa\u56e0\u96c6\u4e3a\u9ad8\u7b49\u7ea7\u7684GO terms<\/p>\n

read_and_merge <- function(data_dir, file_names) {<\/p>\n

df <- tibble()<\/p>\n

for (file_name in file_names) {<\/p>\n

file_path <- glue(“{data_dir}\/DAVIDKnowledgebase\/{file_name}”)<\/p>\n

df_temp <- read_tsv(file_path)<\/p>\n

df <- rbind(df, df_temp)<\/p>\n

}<\/p>\n

return(df)<\/p>\n

}<\/p>\n

file_names <- c(“OFFICIAL_GENE_SYMBOL2GOTERM_BP_2.txt”,<\/p>\n

“OFFICIAL_GENE_SYMBOL2GOTERM_BP_3.txt”,<\/p>\n

“OFFICIAL_GENE_SYMBOL2GOTERM_BP_4.txt”)<\/p>\n

df <- read_and_merge(data_dir, file_names)<\/p>\n

GO <- unique(df$GO_term)<\/p>\n

GO_set <- df %>%<\/p>\n

group_by(GO_term) %>%<\/p>\n

summarise(gene = toString(unique(gene)))<\/p>\n

## \u8fc7\u6ee4\u5c11\u4e8e10\u4e2a\u57fa\u56e0\u7684\u57fa\u56e0\u96c6<\/p>\n

indx = c()<\/p>\n

for (i in seq(1,dim(GO_set)[1])){<\/p>\n

temp = GO_set[i,2][[1]]<\/p>\n

x = strsplit(temp,’, ‘)<\/p>\n

geneset = x[[1]]<\/p>\n

if (length(geneset)>10){<\/p>\n

indx = c(indx,i)<\/p>\n

}<\/p>\n

}<\/p>\n

GO_set = GO_set[indx,]<\/p>\n

GO_set = GO_set[!duplicated(GO_set), ] # remove duplicate<\/p>\n

GO_set = GO_set[!is.na(GO_set$GO_term),] # remove NA<\/p>\n

## \u7b5b\u9009\u5dee\u5f02\u8868\u8fbe\u57fa\u56e0<\/p>\n

dir1 = paste0(save_dir,’\/ST differentially expresssed genes up and down’)<\/p>\n

files = list.files(dir1)<\/p>\n

clusters = list()<\/p>\n

clu_nm = c()<\/p>\n

for (i in seq(1,length(files))){<\/p>\n

fnm = files[i]<\/p>\n

clu = unlist(strsplit(fnm,’.x’))[1]<\/p>\n

tmp1 = read.table(paste0(dir1,’\/’,fnm),header = TRUE, stringsAsFactors = F, check.names=F)<\/p>\n

tmp1 = tmp1[tmp1$avg_logFC>0,]<\/p>\n

tmp1 = tmp1[tmp1$p_val<0.01,]<\/p>\n

clusters[[clu]] = list(tmp1$gene)<\/p>\n

clu_nm = c(clu_nm,clu)<\/p>\n

}<\/p>\n

## \u5bf9\u4e00\u4e2a\u57fa\u56e0\u96c6\u5408\u8fdb\u884cGO\u5bcc\u96c6\u5206\u6790\uff0c\u4f7f\u7528\u8d85\u51e0\u4f55\u5206\u5e03\u68c0\u9a8c\u6765\u8ba1\u7b97\u6bcf\u4e2aGO\u672f\u8bed\u5728\u6bcf\u4e2a\u7ec6\u80de\u7fa4\u4e2d\u7684\u5bcc\u96c6\u7a0b\u5ea6\uff0c\u7136\u540e\u4f7f\u7528FDR\uff08false discovery rate\uff09\u65b9\u6cd5\u6765\u6821\u6b63p\u503c<\/p>\n

enrich_matrix = data.frame(matrix(1, ncol = length(clusters), nrow = dim(GO_set)[1]))<\/p>\n

colnames(enrich_matrix) = clu_nm<\/p>\n

temp = as.character(GO_set$GO_term)<\/p>\n

GO_names = c()<\/p>\n

for (i in seq(1,dim(GO_set)[1])){<\/p>\n

temp = GO_set[i,2][[1]]<\/p>\n

x = strsplit(temp,’, ‘)<\/p>\n

geneset = x[[1]]<\/p>\n

temp = as.character(GO_set[i,1][[1]])<\/p>\n

y = strsplit(temp,’~’)<\/p>\n

go_name = y[[1]][2]<\/p>\n

GO_names = c(GO_names,go_name)<\/p>\n

for (j in seq(1,length(clusters))){<\/p>\n

cluster_DE_genes = unlist(clusters[[j]])<\/p>\n

N_ST = length(cluster_DE_genes)<\/p>\n

N_geneset = length(geneset)<\/p>\n

N_intersect = length(intersect(cluster_DE_genes, geneset))<\/p>\n

N_tot = 20000<\/p>\n

# N_tot = 16522<\/p>\n

P = phyper(N_intersect, N_geneset, N_tot-N_geneset, N_ST, lower.tail = FALSE, log.p = FALSE)<\/p>\n

enrich_matrix[i,j] = P<\/p>\n

}<\/p>\n

}<\/p>\n

row.names(enrich_matrix) = GO_names<\/p>\n

enrich_matrix_FDR = data.frame(matrix(1, ncol = length(clusters), nrow = dim(GO_set)[1]))<\/p>\n

colnames(enrich_matrix_FDR) = clu_nm<\/p>\n

row.names(enrich_matrix_FDR) = GO_names<\/p>\n

for (i in seq(1,length(clusters))){<\/p>\n

temp = enrich_matrix[,i]<\/p>\n

enrich_matrix_FDR[,i] = p.adjust(temp,method=”fdr”)<\/p>\n

}<\/p>\n

## \u8fc7\u6ee4\u6389\u6240\u6709FDR<0.01\u7684\u57fa\u56e0\uff08FDR-\u9519\u8bef\u53d1\u73b0\u7387\uff09<\/p>\n

indx = c()<\/p>\n

for (i in seq(1,dim(enrich_matrix_FDR)[1])){<\/p>\n

if (min(enrich_matrix_FDR[i,])<0.01){<\/p>\n

indx = c(indx,i)<\/p>\n

}<\/p>\n

}<\/p>\n

enrich_matrix_f = enrich_matrix_FDR[indx,]<\/p>\n

## \u627e\u5230\u6bcf\u7c07\u7684\u8868\u8fbe\u6700\u9ad8\u768420\u4e2a\u57fa\u56e0\u7ec4<\/p>\n

all_genesets = row.names(enrich_matrix_f)<\/p>\n

top_genesets = c()<\/p>\n

for (i in seq(1,length(clusters))){<\/p>\n

temp = enrich_matrix_f[,i]<\/p>\n

a = sort(temp,index.return = TRUE)<\/p>\n

top = a$ix[1:10]<\/p>\n

top_genesets = c(top_genesets,all_genesets[top])<\/p>\n

}<\/p>\n

top_genesets = unique(top_genesets)<\/p>\n

enrich_matrix_f = enrich_matrix_f[top_genesets,]<\/p>\n

enrich_matrix_log = -log2(enrich_matrix_f+1e-100)<\/p>\n

## \u7ed8\u5236\u70ed\u56fe<\/p>\n

my_palette <- colorRampPalette(brewer.pal(10, “RdBu”))(256)<\/p>\n

data = data.matrix(enrich_matrix_log)<\/p>\n

<\/p>\n

GO-BP-enrichment-top20<\/p>\n

\u5982\u56fe\uff0c\u662f\u6211\u4eec\u5bcc\u96c6\u4e86\u8868\u8fbe\u6c34\u5e73\u6700\u9ad8\u768420\u4e2a\u57fa\u56e0\u7684\u70ed\u56fe\uff0c\u5176\u4e2d\u7ea2\u8272\u4ee3\u8868\u9ad8\u8868\u8fbe\u6c34\u5e73\uff0c\u84dd\u8272\u4ee3\u8868\u4f4e\u8868\u8fbe\u6c34\u5e73\uff0c\u8fd8\u5217\u51fa\u4e86\u4e0d\u540c\u7684\u751f\u7269\u5b66\u8fc7\u7a0b\u548c\u76f8\u5173\u57fa\u56e0\uff0c\u6bd4\u5982\u8bf4\u5728\u7ec6\u80de\u514d\u75ab\u53cd\u5e94\u8fc7\u7a0b\u5728A4\uff0cA5\u7ec6\u80de\u4e2d\u7684\u8868\u8fbe\u91cf\u660e\u663e\u6bd4A10\u548cA12\u9ad8\u3002<\/p>\n

## \u8fd0\u7528\u5b8c\u5168\u94fe\u63a5\u8fdb\u884c\u805a\u7c7b<\/p>\n

hr <- hclust(as.dist(1-cor(t(data), method=”spearman”)), method=”complete”)<\/p>\n

hc <- hclust(as.dist(1-cor(data, method=”spearman”)), method=”complete”)<\/p>\n

data[data>30] = 30<\/p>\n

filename = paste0(save_dir,’\/heatmap\/GO-BP-enrichment-top20.pdf’)<\/p>\n

pdf(filename,12,12)<\/p>\n

heatmap.2(data,<\/p>\n

col = rev(my_palette),<\/p>\n

Rowv=as.dendrogram(hr), Colv=as.dendrogram(hc),<\/p>\n

trace=”none”,<\/p>\n

scale = “none”,<\/p>\n

density.info=”none”,<\/p>\n

margins=c(5,25),<\/p>\n

lhei=c(1,8), lwid=c(1,4),<\/p>\n

keysize=0.2, key.par = list(cex=0.5),<\/p>\n

cexRow=0.9,cexCol=1,srtCol=90, # rotate column label<\/p>\n

key.title = ‘-log(FDR)’,<\/p>\n

key.xlab = ‘-log(FDR)’<\/p>\n

)<\/p>\n

dev.off()<\/p>\n

## \u4ece\u6570\u636e\u4e2d\u9009\u62e9\u7279\u5b9a\u7684\u751f\u7269\u8fc7\u7a0b<\/p>\n

selected = c(‘extracellular matrix organization’, ‘cell adhesion’, ‘cell migration’,<\/p>\n

‘cell cycle process’,’mitotic cell cycle process’,’protein transport’,<\/p>\n

‘immune response’,’RNA metabolic process’,’peptide metabolic process’,<\/p>\n

‘vasculature development’,<\/p>\n

‘extracellular matrix disassembly’,’cellular metabolic process’,<\/p>\n

‘cellular respiration’,’respiratory electron transport chain’,<\/p>\n

‘response to stress’,’cell motility’)<\/p>\n

temp1 = row.names(data)<\/p>\n

selected = intersect(selected,temp1)<\/p>\n

data1 = data[selected,]<\/p>\n

\u5728\u83b7\u5f97\u7684\u70ed\u56fe\u4e2d\uff0c\u4f9d\u636e\u4ee3\u7801\uff0c\u6211\u4eec\u4ece\u521d\u6b65\u5bcc\u96c6\u5206\u6790\u7684\u7ed3\u679c\u4e2d\uff0c\u518d\u8fdb\u884c\u5bcc\u96c6\u5206\u6790\uff0c\u5e76\u4e14\u4ece\u4e2d\u9009\u51fa\u8868\u8fbe\u91cf\u9ad8\u7684\u751f\u7269\u5b66\u8fc7\u7a0b\u8fdb\u884c\u805a\u7c7b\uff0c\u4f7f\u5f97\u6574\u4e2a\u5bcc\u96c6\u7684\u7ed3\u679c\u66f4\u52a0\u7b80\u6d01\u3002<\/p>\n

<\/p>\n

GO-BP-enrichment-selected<\/p>\n

\u7b2c\u4e8c\u6b65\uff1a\u5bf9\u6240\u6709\u7684\u524d\u666f\u57fa\u56e0\uff08\u764c\u75c7\u76f8\u5173\u6210\u7ea4\u7ef4\u57fa\u56e0\uff09\u8fdb\u884cGO\u5206\u6790<\/p>\n

## \u8bbe\u5b9aGO\u57fa\u56e0\u96c6\u4e3a\u9ad8\u7b49\u7ea7\u7684GO terms<\/p>\n

database = paste0(data_dir,’\/DAVIDKnowledgebase\/OFFICIAL_GENE_SYMBOL2GOTERM_BP_2.txt’)<\/p>\n

df1 = read.delim2(database, header = FALSE, sep = “\\t”)<\/p>\n

colnames(df1) = c(‘gene’,’GO_term’)<\/p>\n

database = paste0(data_dir,’\/DAVIDKnowledgebase\/OFFICIAL_GENE_SYMBOL2GOTERM_BP_3.txt’)<\/p>\n

df2 = read.delim2(database, header = FALSE, sep = “\\t”)<\/p>\n

colnames(df2) = c(‘gene’,’GO_term’)<\/p>\n

df = rbind(df1,df2)<\/p>\n

database = paste0(data_dir,’\/DAVIDKnowledgebase\/OFFICIAL_GENE_SYMBOL2GOTERM_BP_4.txt’)<\/p>\n

df2 = read.delim2(database, header = FALSE, sep = “\\t”)<\/p>\n

colnames(df2) = c(‘gene’,’GO_term’)<\/p>\n

df = rbind(df,df2)<\/p>\n

GO = unique(df$`GO_term`)<\/p>\n

GO_set = df %>%<\/p>\n

group_by(GO_term) %>%<\/p>\n

summarise(gene = toString(unique(gene)))<\/p>\n

## \u8fc7\u6ee4\u591a\u4e8e10\u4e2a\u57fa\u56e0\u7684\u57fa\u56e0\u96c6<\/p>\n

indx = c()<\/p>\n

for (i in seq(1,dim(GO_set)[1])){<\/p>\n

temp = GO_set[i,2][[1]]<\/p>\n

x = strsplit(temp,’, ‘)<\/p>\n

geneset = x[[1]]<\/p>\n

if (length(geneset)>10){<\/p>\n

indx = c(indx,i)<\/p>\n

}<\/p>\n

}<\/p>\n

GO_set = GO_set[indx,]<\/p>\n

## \u7b5b\u9009\u524d\u666f\u57fa\u56e0\u7684\u5dee\u5f02\u8868\u8fbe\u57fa\u56e0<\/p>\n

dir1 = paste0(save_dir,’\/ST differentially expresssed genes up and down’)<\/p>\n

files = list.files(dir1)<\/p>\n

clusters = list()<\/p>\n

clu_nm = c()<\/p>\n

clusters_select = c(“A4_c0”, “A4_c7”, “A5_c0”, “A5_c1”, “A10_c0”, “A10_c1”, “A10_c4”, “A12_c2”, “A12_c3”, “A12_c4”)<\/p>\n

for (i in seq(1,length(files))){<\/p>\n

fnm = files[i]<\/p>\n

clu = unlist(strsplit(fnm,’.x’))[1]<\/p>\n

if (clu %in% clusters_select){<\/p>\n

tmp1 = read.table(paste0(dir1,’\/’,fnm),header = TRUE, stringsAsFactors = F, check.names=F)<\/p>\n

tmp1 = tmp1[tmp1$avg_logFC>0.15,]<\/p>\n

tmp1 = tmp1[tmp1$p_val<0.01,]<\/p>\n

clusters[[clu]] = list(tmp1$gene)<\/p>\n

clu_nm = c(clu_nm,clu)<\/p>\n

}<\/p>\n

}<\/p>\n

## \u8d85\u51e0\u4f55\u5206\u5e03\u68c0\u9a8c<\/p>\n

calc_p_fdr <- function(cluster, geneset, N_tot = 20000) {<\/p>\n

cluster_DE_genes <- unlist(cluster)<\/p>\n

N_ST <- length(cluster_DE_genes)<\/p>\n

N_geneset <- length(geneset)<\/p>\n

N_intersect <- length(intersect(cluster_DE_genes, geneset))<\/p>\n

P <- phyper(N_intersect, N_geneset, N_tot – N_geneset, N_ST, lower.tail = FALSE, log.p = FALSE)<\/p>\n

return(P)<\/p>\n

}<\/p>\n

# \u5b9a\u4e49\u4e00\u4e2a\u51fd\u6570\uff0c\u7528\u4e8e\u63d0\u53d6 GO \u540d\u79f0\u548c\u57fa\u56e0\u96c6<\/p>\n

extract_go_geneset <- function(go_term) {<\/p>\n

temp <- as.character(go_term)<\/p>\n

x <- strsplit(temp, “, “)<\/p>\n

geneset <- x[[1]]<\/p>\n

y <- strsplit(temp, “~”)<\/p>\n

go_name <- y[[1]][2]<\/p>\n

return(list(go_name = go_name, geneset = geneset))<\/p>\n

}<\/p>\n

GO_list <- map(GO_set, extract_go_geneset)<\/p>\n

# \u4f7f\u7528 map \u51fd\u6570\u5bf9 GO_list \u7684\u6bcf\u4e2a\u5143\u7d20\u8fdb\u884c\u8ba1\u7b97\u64cd\u4f5c\uff0c\u5f97\u5230 P \u503c\u7684\u77e9\u9635<\/p>\n

enrich_matrix <- map_dfc(clusters, ~ {<\/p>\n

map_dbl(GO_list, ~ {<\/p>\n

calc_p_fdr(.x, .y$geneset)<\/p>\n

})<\/p>\n

}) %>%<\/p>\n

as.data.frame()<\/p>\n

colnames(enrich_matrix) <- clu_nm<\/p>\n

GO_names <- map_chr(GO_list, “go_name”)<\/p>\n

row.names(enrich_matrix) <- GO_names<\/p>\n

enrich_matrix_FDR <- enrich_matrix %>%<\/p>\n

mutate(across(everything(), ~<\/p>\n

p.adjust(., method = “fdr”)<\/p>\n

}))<\/p>\n

row.names(enrich_matrix_FDR) <- GO_names<\/p>\n

## \u8fc7\u6ee4\u6389\u6240\u6709FDR<0.01\u7684\u57fa\u56e0\uff08FDR-\u9519\u8bef\u53d1\u73b0\u7387\uff09<\/p>\n

indx = c()<\/p>\n

for (i in seq(1,dim(enrich_matrix_FDR)[1])){<\/p>\n

if (min(enrich_matrix_FDR[i,])<0.01){<\/p>\n

indx = c(indx,i)<\/p>\n

}<\/p>\n

}<\/p>\n

enrich_matrix_f = enrich_matrix_FDR[indx,]<\/p>\n

## \u627e\u5230\u6bcf\u7c07\u7684\u8868\u8fbe\u6700\u9ad8\u768420\u4e2a\u57fa\u56e0\u7ec4<\/p>\n

all_genesets = row.names(enrich_matrix_f)<\/p>\n

top_genesets = c()<\/p>\n

for (i in seq(1,length(clusters))){<\/p>\n

temp = enrich_matrix_f[,i]<\/p>\n

a = sort(temp,index.return = TRUE)<\/p>\n

top = a$ix[1:20]<\/p>\n

top_genesets = c(top_genesets,all_genesets[top])<\/p>\n

}<\/p>\n

top_genesets = unique(top_genesets)<\/p>\n

enrich_matrix_f = enrich_matrix_f[top_genesets,]<\/p>\n

enrich_matrix_log = -log2(enrich_matrix_f+1e-100)<\/p>\n

## \u7ed8\u5236\u70ed\u56fe<\/p>\n

my_palette <- colorRampPalette(brewer.pal(10, “RdBu”))(256)<\/p>\n

data = data.matrix(enrich_matrix_log)<\/p>\n

\u8fd9\u4e00\u6b65\u7684\u5206\u6790\u4ee3\u7801\u4e0e\u4e0a\u4e00\u6b65\u57fa\u672c\u7c7b\u4f3c\uff0c\u533a\u522b\u5728\u4e8e\uff0c\u4e0a\u4e00\u6b65\u5206\u6790\u7684\u662f\u6240\u6709\u6837\u672c\uff0c\u800c\u8fd9\u4e00\u6b65\u53ea\u5206\u6790\u4e86\u524d\u666f\u57fa\u56e0\u8fd9\u90e8\u5206\u6837\u672c\u3002<\/p>\n

<\/p>\n

GO-BP-enrichment-CAF-detail<\/p>\n

\u7b2c\u4e09\u6b65\uff0c\u52a0\u8f7d\u914d\u4f53\u53d7\u4f53\u53c2\u8003\u6570\u636e\u5e93<\/p>\n

\u9996\u5148\u8bbe\u5b9a\u914d\u4f53\u53d7\u4f53\u6570\u636e\u5e93\uff0c\u5982\u56fe\u6240\u793a\uff0c\u9700\u8981\u5bfc\u5165\u76f8\u5173\u7684\u53d7\u4f53\u4e0e\u914d\u4f53\u7684\u6570\u636e\uff0c\u7136\u4e4e\u83b7\u53d6\u80bf\u7624\u548c\u57fa\u8d28\u7a7a\u95f4\u4e0a\u7684\u76f8\u90bb\u7684\u914d\u4f53\u53d7\u4f53\u5bf9\uff0c\u5e76\u63d0\u53d6\u4f4d\u70b9\u5750\u6807\uff0c\u8fdb\u884c\u53ef\u89c6\u5316\u3002<\/p>\n

## \u8bbe\u5b9a\u53d7\u4f53\u914d\u4f53\u6570\u636e\u5e93<\/p>\n

dir1 = paste0(data_dir,’\/LR_manual_revised.txt’)<\/p>\n

file = read.delim(dir1,sep = ‘\\t’,stringsAsFactors = FALSE)<\/p>\n

folder = paste0(save_dir,’\/ligand_receptor’)<\/p>\n

if (!file.exists(folder)){<\/p>\n

dir.create(folder)<\/p>\n

}<\/p>\n

dir2 = paste0(folder,’\/nearest neighbour ligand receptor interactions of two clusters’)<\/p>\n

if (!file.exists(dir2)){<\/p>\n

dir.create(dir2)<\/p>\n

}<\/p>\n

for (sampleid in c(‘A4′,’A5′,’A10′,’A12′)){<\/p>\n

fnm = paste0(dir2,’\/’,sampleid)<\/p>\n

if (!file.exists(fnm)){<\/p>\n

dir.create(fnm)<\/p>\n

}<\/p>\n

}<\/p>\n

dir3 = paste0(folder,’\/nearest neighbours of two clusters images’)<\/p>\n

if (!file.exists(dir3)){<\/p>\n

dir.create(dir3)<\/p>\n

}<\/p>\n

## \u83b7\u53d6\u4e24\u4e2a\u7ec6\u80de\u7fa4\uff08\u80bf\u7624\uff0c\u57fa\u8d28\uff09\u4e2d\u7a7a\u95f4\u4e0a\u76f8\u90bb\u6591\u70b9\u7684\u914d\u4f53\u53d7\u4f53\u5bf9<\/p>\n

# \u5b9a\u4e49\u4e00\u4e2a\u51fd\u6570\uff0c\u7528\u4e8e\u8ba1\u7b97\u57fa\u56e0\u8868\u8fbe\u7684\u5747\u503c\u548c\u6807\u51c6\u5dee<\/p>\n

calc_expression <- function(SOE, gene, tumor_loc, stroma_loc, stroma_ids) {<\/p>\n

# \u8c03\u7528 min_distance \u51fd\u6570\uff0c\u5f97\u5230\u57fa\u8d28\u5230\u80bf\u7624\u7684\u6700\u5c0f\u8ddd\u79bb<\/p>\n

dist_to_tumor <- min_distance(tumor_loc, stroma_loc)<\/p>\n

distance <- sort(unique(dist_to_tumor))<\/p>\n

data_mat <- SOE@assays$SCT@scale.data[, stroma_ids]<\/p>\n

expression_dist_mean <- map_dbl(distance, ~ {<\/p>\n

indx <- which(dist_to_tumor == .x)<\/p>\n

mean_exp <- mean(data_mat[gene, indx])<\/p>\n

return(mean_exp)<\/p>\n

})<\/p>\n

expression_dist_std <- map_dbl(distance, ~ {<\/p>\n

indx <- which(dist_to_tumor == .x)<\/p>\n

std_exp <- sd(data_mat[gene, indx])<\/p>\n

return(std_exp)<\/p>\n

}<\/p>\n

return(list(expression_dist_mean, expression_dist_std))<\/p>\n

}<\/p>\n

# \u5b9a\u4e49\u4e00\u4e2a\u51fd\u6570\uff0c\u7528\u4e8e\u8ba1\u7b97\u6700\u5c0f\u8ddd\u79bb<\/p>\n

min_distance <- function(tumor_loc, stroma_loc) {<\/p>\n

mini_dist <- map_dbl(stroma_loc, ~ {<\/p>\n

xF <- .x[1]<\/p>\n

yF <- .x[2]<\/p>\n

min_distance <- 100000<\/p>\n

dist <- map_dbl(tumor_loc, ~ {<\/p>\n

# \u63d0\u53d6\u5f53\u524d\u70b9\u7684 x \u548c y \u5750\u6807<\/p>\n

xT <- .x[1]<\/p>\n

yT <- .x[2]<\/p>\n

sqrt((xF – xT)^2 + (yF – yT)^2)<\/p>\n

})<\/p>\n

min_distance <- min(dist)<\/p>\n

return(min_distance)<\/p>\n

})<\/p>\n

return(mini_dist)<\/p>\n

}<\/p>\n

## \u63d0\u53d6\u6240\u6709\u4f4d\u70b9\u7684\u5750\u6807<\/p>\n

extract_loc <- function(SOE,SID){<\/p>\n

tmp = as.matrix(GetAssayData(object = SOE, assay = “SCT”, slot = “data”))<\/p>\n

tmp1 = colnames(tmp)<\/p>\n

if (nchar(SID)==1){<\/p>\n

tmp1 = substr(tmp1,4,100)<\/p>\n

} else{<\/p>\n

tmp1 = substr(tmp1,5,100)<\/p>\n

}<\/p>\n

x = as.numeric(unlist(lapply(tmp1, function(x){strsplit(x,”x”)[[1]][1]})))<\/p>\n

y = as.numeric(unlist(lapply(tmp1, function(x){strsplit(x,”x”)[[1]][2]})))<\/p>\n

d = data.frame(x,y)<\/p>\n

return(d)<\/p>\n

}<\/p>\n

## \u627e\u5230\u4e24\u4e2a\u8ddd\u79bb\u6700\u8fd1\u7684\u7c07<\/p>\n

find_nearest_spots = function(SOE,loc,clu1,clu2){<\/p>\n

Location = loc<\/p>\n

temp = SOE@meta.data$seurat_clusters<\/p>\n

clu1ID = which(temp==clu1)<\/p>\n

clu2ID = which(temp==clu2)<\/p>\n

Nearest_clu_ID = data.frame(c1 = integer(0),c2 = integer(0))<\/p>\n

count = 1<\/p>\n

for (i in seq(1,length(clu1ID))){<\/p>\n

xF = Location[clu1ID[i],1]<\/p>\n

yF = Location[clu1ID[i],2]<\/p>\n

for (j in seq(1,length(clu2ID))){<\/p>\n

xT = Location[clu2ID[j],1]<\/p>\n

yT = Location[clu2ID[j],2]<\/p>\n

distance = (xF-xT)^2+(yF-yT)^2<\/p>\n

if (distance<=1){<\/p>\n

Nearest_clu_ID[count,1] = clu1ID[i]<\/p>\n

Nearest_clu_ID[count,2] = clu2ID[j]<\/p>\n

count = count+1<\/p>\n

}<\/p>\n

}<\/p>\n

}<\/p>\n

return(Nearest_clu_ID)<\/p>\n

}<\/p>\n

## \u53ef\u89c6\u5316<\/p>\n

# \u5b9a\u4e49\u4e00\u4e2a\u51fd\u6570\uff0c\u7528\u4e8e\u8ba1\u7b97\u70b9\u7684\u4f4d\u7f6e<\/p>\n

get_points <- function(d, topleft, topright, bottomleft, bottomright) {<\/p>\n

delta_x = mean(abs(topright[1] – topleft[1]), abs(bottomright[1] – bottomleft[1])) \/ 32<\/p>\n

delta_y = mean(abs(topright[2] – bottomright[2]), abs(topleft[2] – bottomleft[2])) \/ 34<\/p>\n

adjust_x = function(step) {return((step – 1) * (bottomleft[1] – topleft[1]) \/ 34)}<\/p>\n

adjust_y = function(step) {return((step – 1) * (topright[2] – topleft[2]) \/ 32)}<\/p>\n

x = topleft[1] + (d[, 1] – 1) * delta_x + apply(d, 1, function(x) {return(adjust_x(x[2]))})<\/p>\n

y = topleft[2] – (d[, 2] – 1) * delta_y + apply(d, 1, function(x) {return(adjust_y(x[1]))})<\/p>\n

return(data.frame(x, y))<\/p>\n

}<\/p>\n

# \u5b9a\u4e49\u4e00\u4e2a\u51fd\u6570\uff0c\u7528\u4e8e\u7ed8\u5236\u56fe\u7247<\/p>\n

draw_image <- function(im, im1, points, mycols) {<\/p>\n

png(width = 512, height = 512)<\/p>\n

par(mar = c(0, 0, 0, 0)) # \u8bbe\u7f6e\u96f6\u8fb9\u8ddd<\/p>\n

plot(x = NULL, y = NULL, xlim = c(0, dim(im)[1]), ylim = c(0, dim(im)[2]), pch = ”,<\/p>\n

xaxt = ‘n’, yaxt = ‘n’, xlab = ”, ylab = ”, xaxs = ‘i’, yaxs = ‘i’,<\/p>\n

bty = ‘n’) # \u7ed8\u5236\u7a7a\u767d\u56fe<\/p>\n

rasterImage(im1, xleft = 0, ybottom = 0, xright = dim(im1)[1], ytop = dim(im1)[2]) # \u6dfb\u52a0\u80bf\u7624\u8f6e\u5ed3<\/p>\n

points(points$x, points$y, col = mycols, pch = 20, cex = 3) # \u6dfb\u52a0\u805a\u7c7b\u70b9<\/p>\n

dev.off()<\/p>\n

}<\/p>\n

# \u5b9a\u4e49\u4e00\u4e2a\u51fd\u6570\uff0c\u7528\u4e8e\u9009\u62e9\u805a\u7c7b<\/p>\n

select_cluster <- function(clu1, clu2, SID, SOE) {<\/p>\n

d = extract_loc(SOE, SID) # \u63d0\u53d6\u4f4d\u7f6e\u4fe1\u606f<\/p>\n

temp = SOE@meta.data$seurat_cluster # \u63d0\u53d6\u805a\u7c7b\u4fe1\u606f<\/p>\n

nspots = length(temp)<\/p>\n

mycols <- rep(‘snow2’, nspots)<\/p>\n

mycolpalette = c(‘orange’, ‘paleturquoise4′)<\/p>\n

tmp = find_nearest_spots(SOE, d, clu1, clu2)<\/p>\n

clu1ID = tmp[, 1] # \u63d0\u53d6\u7b2c\u4e00\u4e2a\u805a\u7c7b\u7684ID<\/p>\n

clu2ID = tmp[, 2]<\/p>\n

mycols[clu1ID] = mycolpalette[1] # \u8bbe\u7f6e\u7b2c\u4e00\u4e2a\u805a\u7c7b\u7684\u989c\u8272<\/p>\n

mycols[clu2ID] = mycolpalette[2]<\/p>\n

return(mycols)<\/p>\n

}<\/p>\n

# \u5b9a\u4e49\u4e00\u4e2a\u51fd\u6570\uff0c\u7528\u4e8e\u53ef\u89c6\u5316<\/p>\n

visualize <- function(SID, clu1, clu2) {<\/p>\n

if (SID == ’10’) {<\/p>\n

fnm = paste0(image_dir, ‘\/Sample_’, SID, ‘_clear.jpg’) # downsample to 0.25x<\/p>\n

im = load.image(fnm)<\/p>\n

xd = dim(im)[1]<\/p>\n

yd = dim(im)[2]<\/p>\n

topleft = c(18, yd – 22)<\/p>\n

topright = c(1046, yd – 18)<\/p>\n

bottomleft = c(18, yd – 1116)<\/p>\n

bottomright = c(1043, yd – 1112)<\/p>\n

} else if (SID == ’12’) {<\/p>\n

fnm = paste0(image_dir, ‘\/Sample_’, SID, ‘_clear.jpg’) # downsample to 0.25x<\/p>\n

im = load.image(fnm)<\/p>\n

xd = dim(im)[1]<\/p>\n

yd = dim(im)[2]<\/p>\n

topleft = c(32, yd – 30)<\/p>\n

topright = c(984, yd – 22)<\/p>\n

bottomleft = c(40, yd – 1114)<\/p>\n

bottomright = c(992, yd – 1106)<\/p>\n

} else if (SID == ‘4’) {<\/p>\n

fnm = paste0(image_dir, ‘\/Sample_’, SID, ‘_ds.jpg’) # downsample to 0.25x<\/p>\n

im = load.image(fnm)<\/p>\n

xd = dim(im)[1]<\/p>\n

yd = dim(im)[2]<\/p>\n

topleft = c(58, yd – 50)<\/p>\n

topright = c(1412, yd – 64)<\/p>\n

bottomleft = c(48, yd – 1512)<\/p>\n

bottomright = c(1400, yd – 1524)<\/p>\n

} else if (SID == ‘5’) {<\/p>\n

fnm = paste0(image_dir, ‘\/Sample_’, SID, ‘_ds.jpg’) # downsample to 0.25x<\/p>\n

im = load.image(fnm)<\/p>\n

xd = dim(im)[1]<\/p>\n

yd = dim(im)[2]<\/p>\n

topleft = c(30, yd – 12)<\/p>\n

topright = c(1382, yd – 20)<\/p>\n

bottomleft = c(22, yd – 1458)<\/p>\n

bottomright = c(1368, yd – 1464)<\/p>\n

}<\/p>\n

# \u52a0\u8f7d\u80bf\u7624\u8f6e\u5ed3\u56fe\u7247<\/p>\n

fnm = paste0(save_dir, ‘\/map_tumor_outline’, ‘\/A’, SID, ‘_tumor_outline.png’) # downsample to 0.25x<\/p>\n

im1 = load.image(fnm)<\/p>\n

# \u9009\u62e9\u805a\u7c7b<\/p>\n

mycols = select_cluster(clu1, clu2, SID, SOE)<\/p>\n

# \u8ba1\u7b97\u70b9\u7684\u4f4d\u7f6e<\/p>\n

points = get_points(d, topleft, topright, bottomleft, bottomright)<\/p>\n

# \u7ed8\u5236\u56fe\u7247<\/p>\n

draw_image(im, im1, points, mycols)<\/p>\n

}<\/p>\n

<\/p>\n

A4_c2_c0<\/p>\n

<\/p>\n

A5_c5_c0<\/p>\n

<\/p>\n

A10_c2_c0<\/p>\n

<\/p>\n

A12_c0_c2<\/p>\n

\u7ed3\u679c\u5982\u56fe\u6240\u793a\uff0c\u662f\u53d7\u4f53\u4e0e\u914d\u4f53\u7684\u76f8\u5173\u914d\u5bf9\u60c5\u51b5\uff0c\u8fd9\u4e2a\u7ed3\u679c\u5bf9\u4e8e\u6211\u4eec\u540e\u7eed\u7ed8\u5236\u51b2\u79ef\u56fe\u662f\u5fc5\u987b\u7684\u3002<\/p>\n

## \u8bbe\u5b9a\u4e3b\u51fd\u6570\uff0c\u5faa\u73af<\/p>\n

SIDs = c(rep(‘4’,6),rep(‘5′,4),rep(’10’,5),rep(’12’,9)) # sample ID<\/p>\n

clu1s = c(rep(2,5),rep(5,5),rep(2,2),rep(5,1),rep(6,2),rep(0,3),rep(1,3),rep(5,3)) # tumor clusters<\/p>\n

clu2s = c(0,1,3,6,7,7,0,1,2,6,0,4,0,0,1,2,3,4,2,3,4,2,3,4) # stroma clusters<\/p>\n

for (i in seq(1,length(clu1s))){<\/p>\n

if (SIDs[i]==’4′){<\/p>\n

SOE = sample4<\/p>\n

}<\/p>\n

if (SIDs[i]==’5′){<\/p>\n

SOE = sample5<\/p>\n

}<\/p>\n

if (SIDs[i]==’10’){<\/p>\n

SOE = sample10<\/p>\n

}<\/p>\n

if (SIDs[i]==’12’){<\/p>\n

SOE = sample12<\/p>\n

}<\/p>\n

SID = SIDs[i]<\/p>\n

clu1 = clu1s[i]<\/p>\n

clu2 = clu2s[i]<\/p>\n

<\/p>\n

A4_ligand_c0_receptor_c2<\/p>\n

\u9664\u4e86\u53ef\u89c6\u5316\u7684\u7ed3\u679c\u4e4b\u5916\uff0c\u6211\u4eec\u8fd8\u53ef\u4ee5\u4ece\u6587\u4ef6\u4e2d\u89c2\u5bdf\u5230\u56db\u4e2a\u6837\u672c\u7684\u53d7\u4f53\u914d\u4f53\u914d\u5bf9\u60c5\u51b5\u3002<\/p>\n

\u7b2c\u56db\u6b65\uff0c\u7ed8\u5236\u6bcf\u4e2a\u6837\u54c1\u7684\u914d\u4f53-\u53d7\u4f53\u76f8\u4e92\u4f5c\u7528\u7684 ggalluvial \u56fe<\/p>\n

\u4f9d\u636e\u5c0f\u679c\u63d0\u4f9b\u7684\u4ee3\u7801\uff0c\u5206\u522b\u4e3a\u4e3aLTS\uff0cSTS\u5206\u522b\u7ed8\u5236 ggalluvial \u56fe\uff0c\u5982\u56fe\u6240\u793a\uff0c\u4f9d\u636eggalluvial \u56fe\uff0c\u53ef\u4ee5\u975e\u5e38\u6e05\u6670\u7684\u770b\u5230\u4e0d\u540c\u7684\u57fa\u56e0\u8868\u8fbe\u4ea7\u7269\u7684\u53d7\u4f53\u914d\u4f53\u7684\u914d\u5bf9\u5173\u7cfb\uff0c\u6839\u636e\uff21\uff14\uff0c\uff21\uff15\uff0c\uff21\uff11\uff10\uff0c\uff21\uff11\uff12\uff0c\u6c47\u96c6\u4e3aLTS\uff08\u957f\u751f\u5b58\u5468\u671f\u60a3\u8005\uff09\uff0cSTS\uff08\u77ed\u751f\u5b58\u5468\u671f\u60a3\u8005\uff09\uff0c\u5373\u53ef\u8fa8\u660e\u5bfc\u81f4\u75c5\u4eba\u751f\u5b58\u5468\u671f\u53d8\u5316\u7684\u5173\u952e\u57fa\u56e0\uff0c\u4ece\u800c\u4e3a\u764c\u7684\u6cbb\u7406\u63d0\u4f9b\u6307\u5bfc\u3002<\/p>\n

## \u7ed8\u5236 ggalluvial \u56fe<\/p>\n

ggall_plot <- function(df){<\/p>\n

ggplot(data = df,<\/p>\n

aes(axis1 = ligand, axis2 = receptor, y = sqrt(connection))) +<\/p>\n

geom_alluvium(aes(fill = ligand)) +<\/p>\n

geom_stratum(alpha = .5, width = 1\/4) +<\/p>\n

geom_text(stat = “stratum”,size = 1.9,<\/p>\n

aes(label = after_stat(stratum))) +<\/p>\n

scale_x_discrete(limits = c(“ligand”, “receptor”),<\/p>\n

expand = c(0.1, 0)) +<\/p>\n

theme_void()<\/p>\n

}<\/p>\n

## \u4e3a\u4e2a\u4f53\u6837\u672c\u7ed8\u5236 ggalluvial \u56fe<\/p>\n

for (j in seq(1,length(dirs))){<\/p>\n

files = list.files(dirs[j])<\/p>\n

df = data.frame()<\/p>\n

for (i in seq(1,length(files))){<\/p>\n

file = paste0(dirs[j],’\/’,files[i])<\/p>\n

tmp = read.delim2(file,sep = ‘,’)<\/p>\n

tmp$avg_logFC.of.ligand = as.numeric(as.character(tmp$avg_logFC.of.ligand))<\/p>\n

tmp$avg_logFC.of.receptor = as.numeric(as.character(tmp$avg_logFC.of.receptor))<\/p>\n

tmp$fraction.of.spots.the.ligand.is.detected = as.numeric(as.character(tmp$fraction.of.spots.the.ligand.is.detected))<\/p>\n

tmp$fraction.of.spots.the.receptor.is.detected = as.numeric(as.character(tmp$fraction.of.spots.the.receptor.is.detected))<\/p>\n

tmp$connection = sqrt((tmp$avg_logFC.of.ligand)*(tmp$avg_logFC.of.receptor)*(tmp$fraction.of.spots.the.ligand.is.detected)*(tmp$fraction.of.spots.the.receptor.is.detected))<\/p>\n

df = rbind(df,tmp)<\/p>\n

}<\/p>\n

df = df[,c(‘ligand’,’receptor’,’connection’)]<\/p>\n

df1 = aggregate(df$connection,by = list(df$ligand,df$receptor),FUN = mean)<\/p>\n

colnames(df1) = c(‘ligand’,’receptor’,’connection’)<\/p>\n

temp = str_split(dirs[j],’\/’)<\/p>\n

fnm = temp[[1]][length(temp[[1]])]<\/p>\n

print(fnm)<\/p>\n

#plot.new()<\/p>\n

pdf(paste0(folder,’\/’,fnm,’.pdf’),8,8)<\/p>\n

fig = ggall_plot(df1)<\/p>\n

print(fig)<\/p>\n

dev.off()<\/p>\n

# chordDiagram(df1)<\/p>\n

}<\/p>\n

## \u4e3aLTS\/STS\u7ed8\u5236 ggalluvial \u56fe<\/p>\n

## STS<\/p>\n

df = data.frame()<\/p>\n

for (j in seq(1,2)){<\/p>\n

files = list.files(dirs[j])<\/p>\n

for (i in seq(1,length(files))){<\/p>\n

file = paste0(dirs[j],’\/’,files[i])<\/p>\n

tmp = read.delim2(file,sep = ‘,’)<\/p>\n

tmp$avg_logFC.of.ligand = as.numeric(as.character(tmp$avg_logFC.of.ligand))<\/p>\n

tmp$avg_logFC.of.receptor = as.numeric(as.character(tmp$avg_logFC.of.receptor))<\/p>\n

tmp$fraction.of.spots.the.ligand.is.detected = as.numeric(as.character(tmp$fraction.of.spots.the.ligand.is.detected))<\/p>\n

tmp$fraction.of.spots.the.receptor.is.detected = as.numeric(as.character(tmp$fraction.of.spots.the.receptor.is.detected))<\/p>\n

tmp$connection = sqrt((tmp$avg_logFC.of.ligand)*(tmp$avg_logFC.of.receptor)*(tmp$fraction.of.spots.the.ligand.is.detected)*(tmp$fraction.of.spots.the.receptor.is.detected))<\/p>\n

df = rbind(df,tmp)<\/p>\n

}<\/p>\n

}<\/p>\n

df = df[,c(‘ligand’,’receptor’,’connection’)]<\/p>\n

df1 = aggregate(df$connection,by = list(df$ligand,df$receptor),FUN = mean)<\/p>\n

colnames(df1) = c(‘ligand’,’receptor’,’connection’)<\/p>\n

pdf(paste0(folder,’\/STS.pdf’),8,10)<\/p>\n

fig = ggall_plot(df1)<\/p>\n

print(fig)<\/p>\n

dev.off()<\/p>\n

## LTS<\/p>\n

df = data.frame()<\/p>\n

for (j in seq(3,4)){<\/p>\n

files = list.files(dirs[j])<\/p>\n

for (i in seq(1,length(files))){<\/p>\n

file = paste0(dirs[j],’\/’,files[i])<\/p>\n

tmp = read.delim2(file,sep = ‘,’)<\/p>\n

tmp$avg_logFC.of.ligand = as.numeric(as.character(tmp$avg_logFC.of.ligand))<\/p>\n

tmp$avg_logFC.of.receptor = as.numeric(as.character(tmp$avg_logFC.of.receptor))<\/p>\n

tmp$fraction.of.spots.the.ligand.is.detected = as.numeric(as.character(tmp$fraction.of.spots.the.ligand.is.detected))<\/p>\n

tmp$fraction.of.spots.the.receptor.is.detected = as.numeric(as.character(tmp$fraction.of.spots.the.receptor.is.detected))<\/p>\n

tmp$connection = sqrt((tmp$avg_logFC.of.ligand)*(tmp$avg_logFC.of.receptor)*(tmp$fraction.of.spots.the.ligand.is.detected)*(tmp$fraction.of.spots.the.receptor.is.detected))<\/p>\n

df = rbind(df,tmp)<\/p>\n

}<\/p>\n

}<\/p>\n

df = df[,c(‘ligand’,’receptor’,’connection’)]<\/p>\n

df1 = aggregate(df$connection,by = list(df$ligand,df$receptor),FUN = mean)<\/p>\n

colnames(df1) = c(‘ligand’,’receptor’,’connection’)<\/p>\n

pdf(paste0(folder,’\/LTS.pdf’),8,10)<\/p>\n

fig = ggall_plot(df1)<\/p>\n

print(fig)<\/p>\n

dev.off()<\/p>\n

<\/p>\n

LTS<\/p>\n

STS<\/p>\n

\u5982\u56fe\u662f\u6211\u4eec\u7ed8\u5236\u7684\u51b2\u79ef\u56fe\uff0c\u53ef\u4ee5\u975e\u5e38\u76f4\u89c2\u7684\u89c2\u5bdf\u5230\u53d7\u4f53\u4e0e\u914d\u4f53\u7684\u5206\u5e03\u60c5\u51b5\u4ee5\u53ca\u53cd\u5e94\u60c5\u51b5\uff0c\u540c\u65f6\uff0c\u4e5f\u53ef\u4ee5\u6839\u636e\u6807\u8bb0\u7684\u53d7\u4f53\u4e0e\u914d\u4f53\u540d\u79f0\u5bf9\u5dee\u5f02\u57fa\u56e0\u8fdb\u884c\u533a\u5206\uff0c\u53ef\u4ee5\u8bf4\u662f\u4e00\u4e3e\u4e24\u5f97\u3002<\/p>\n

\u4ee5\u4e0a\u5c31\u662f\u5c0f\u679c\u4eca\u5929\u7ed9\u5927\u5bb6\u5e26\u6765\u7684\u4f7f\u7528\u5bcc\u96c6\u5206\u6790\u7684\u7ed3\u679c\u4e0eggalluvial\u7ed8\u5236\u51b2\u79ef\u56fe\u7684\u6559\u7a0b\u5566\uff0c\u5982\u679c\u5927\u5bb6\u5bf9\u4eca\u5929\u7684\u5206\u4eab\u6709\u4efb\u4f55\u7591\u95ee\uff0c\u6b22\u8fce\u5173\u6ce8\u516c\u4f17\u53f7\u5411\u5c0f\u679c\u63d0\u95ee\uff0c\u6211\u4eec\u6709\u4e13\u95e8\u7684\u751f\u4fe1\u5206\u6790\u670d\u52a1\u54e6\u3002\u540c\u65f6\uff0c\u5982\u679c\u5927\u5bb6\u6ca1\u6709\u505a\u51fa\u6765\u56fe\u4e5f\u4e0d\u8981\u7070\u5fc3\uff0c\u53ef\u4ee5\u8bd5\u4e00\u8bd5\u6211\u4eec\u7684\u4e91\u751f\u4fe1\u5c0f\u5de5\u5177\u54e6\uff0c\u53ea\u8981\u8f93\u5165\u5408\u9002\u7684\u6570\u636e\u4ee5\u53ca\u6307\u4ee4\u5c31\u53ef\u4ee5\u76f4\u63a5\u7ed8\u5236\u60f3\u8981\u7684\u56fe\u5462\uff0c\u94fe\u63a5\uff1ahttp:\/\/www.biocloudservice.com\/home.html\u3002<\/p>\n","protected":false},"excerpt":{"rendered":"

\u6211\u4eec\u90fd\u77e5\u9053\u5bcc\u96c6\u5206\u6790\u662f\u4e00\u79cd\u5e7f\u6cdb\u5e94\u7528\u4e8e\u751f\u7269\u4fe1\u606f\u5b66\u7814\u7a76\u7684\u7edf\u8ba1\u65b9\u6cd5\uff0c\u4e3b\u8981\u7528\u4e8e\u68c0\u9a8c\u4e00\u4e2a\u57fa\u56e0\u96c6\u5408\u4e2d\u67d0\u4e9b\u529f\u80fd\u6216\u7279\u5f81\u7684\u5bcc\u96c6\u7a0b\u5ea6 […]<\/p>\n","protected":false},"author":1,"featured_media":0,"comment_status":"closed","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":{"jetpack_post_was_ever_published":false,"_jetpack_newsletter_access":""},"categories":[1],"tags":[],"jetpack_featured_media_url":"","_links":{"self":[{"href":"http:\/\/www.biocloudservice.com\/wordpress\/index.php?rest_route=\/wp\/v2\/posts\/27993"}],"collection":[{"href":"http:\/\/www.biocloudservice.com\/wordpress\/index.php?rest_route=\/wp\/v2\/posts"}],"about":[{"href":"http:\/\/www.biocloudservice.com\/wordpress\/index.php?rest_route=\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"http:\/\/www.biocloudservice.com\/wordpress\/index.php?rest_route=\/wp\/v2\/users\/1"}],"replies":[{"embeddable":true,"href":"http:\/\/www.biocloudservice.com\/wordpress\/index.php?rest_route=%2Fwp%2Fv2%2Fcomments&post=27993"}],"version-history":[{"count":1,"href":"http:\/\/www.biocloudservice.com\/wordpress\/index.php?rest_route=\/wp\/v2\/posts\/27993\/revisions"}],"predecessor-version":[{"id":28005,"href":"http:\/\/www.biocloudservice.com\/wordpress\/index.php?rest_route=\/wp\/v2\/posts\/27993\/revisions\/28005"}],"wp:attachment":[{"href":"http:\/\/www.biocloudservice.com\/wordpress\/index.php?rest_route=%2Fwp%2Fv2%2Fmedia&parent=27993"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"http:\/\/www.biocloudservice.com\/wordpress\/index.php?rest_route=%2Fwp%2Fv2%2Fcategories&post=27993"},{"taxonomy":"post_tag","embeddable":true,"href":"http:\/\/www.biocloudservice.com\/wordpress\/index.php?rest_route=%2Fwp%2Fv2%2Ftags&post=27993"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}