作者，Evil Genius

今天分享的內(nèi)容來(lái)自于文獻(xiàn)Single-cell and spatial transcriptome analysis reveals the cellular heterogeneity of liver metastatic colorectal cancer, 2023年6月發(fā)表于SCIENCE ADVANCES，IF 9.79分，可見(jiàn)隨著單細(xì)胞空間聯(lián)合分析的越來(lái)越深入，發(fā)文的難度也在提升。

首先我們需要確立一個(gè)觀(guān)點(diǎn)，那就是之前我們都在強(qiáng)調(diào)單細(xì)胞空間的聯(lián)合分析注釋空間轉(zhuǎn)錄組，但很多時(shí)候過(guò)于精細(xì)往往會(huì)忽略結(jié)構(gòu)的整體作用，所以空間聚類(lèi)的作用凸顯了出來(lái)，尤其在組織形態(tài)的加持下，當(dāng)然了，細(xì)胞類(lèi)型的精確匹配也是非常重要，這也是為什么都在提時(shí)空組學(xué)的原因。

華大的時(shí)空組學(xué)是目前國(guó)內(nèi)最好的平臺(tái)，而且匹配了華大研究院專(zhuān)門(mén)分析研究時(shí)空組學(xué)數(shù)據(jù)，在這里雖然我不是華大的人，但還是很希望有朝一日投入到國(guó)產(chǎn)的時(shí)空組學(xué)研究中去，同時(shí)也能離家近一點(diǎn)，希望有機(jī)會(huì)吧，天津華大研究院，不過(guò)目前還是主要10X時(shí)空組項(xiàng)目居多。

好了，開(kāi)始我們今天的結(jié)直腸癌的時(shí)空研究，其中文獻(xiàn)的思路是最主要的。

結(jié)直腸癌(CRC)是一種惡性腫瘤，其中部分CRC存在轉(zhuǎn)移性，尤其是肝轉(zhuǎn)移，這也是預(yù)后較差的一個(gè)指標(biāo)。轉(zhuǎn)移過(guò)程涉及多個(gè)步驟，包括癌細(xì)胞從原發(fā)腫瘤部位逃逸、在血流中存活、在遠(yuǎn)處部位播散，最終生長(zhǎng)為轉(zhuǎn)移性腫瘤。癌細(xì)胞和基質(zhì)細(xì)胞之間的交流在促進(jìn)轉(zhuǎn)移擴(kuò)散中起著至關(guān)重要的作用。這種交流通過(guò)分泌細(xì)胞因子、生長(zhǎng)因子和蛋白酶來(lái)重塑腫瘤微環(huán)境(TME)。TME是一種復(fù)雜的細(xì)胞組成，包括不同群體的成纖維細(xì)胞和免疫細(xì)胞，所有這些細(xì)胞在癌癥逃避、轉(zhuǎn)移和對(duì)治療的反應(yīng)中發(fā)揮重要作用。所以時(shí)空組學(xué)的研究就在于發(fā)現(xiàn)兩種腫瘤之間的區(qū)別。

首先第一步，仍然是CRC單細(xì)胞全圖譜（BD平臺(tái)），包括原發(fā)性CRC和肝轉(zhuǎn)移瘤。包括原發(fā)性CC、鄰近正常結(jié)直腸黏膜(CN)、肝轉(zhuǎn)移(LM)、鄰近正常肝組織(LN)和外周血(PB)。

我們這里注意下細(xì)胞注釋的細(xì)節(jié)，根據(jù)EPCAM和SOX9的表達(dá)鑒定腫瘤細(xì)胞，根據(jù)COL1A1和COL1A2的表達(dá)鑒定成纖維細(xì)胞，根據(jù)PECAM1和CD34的表達(dá)鑒定內(nèi)皮細(xì)胞。免疫細(xì)胞鑒定為T(mén)細(xì)胞(CD3)、自然殺傷(NK)細(xì)胞(CD56)、B細(xì)胞/漿細(xì)胞(CD19)、單核細(xì)胞/巨噬細(xì)胞(CD14)、樹(shù)突狀細(xì)胞(dc) (HLA.DRA)和肥大細(xì)胞(TPSAB1)。當(dāng)然這里也要注意細(xì)胞分布的樣本特異性。

第二步就是CRC和肝轉(zhuǎn)移瘤的空間圖譜

1、通過(guò)蘇木素-伊紅(he)染色和基因表達(dá)特征鑒定每個(gè)樣本的腫瘤區(qū)域(T)和癌旁區(qū)域(PT)。
2、癌旁組織和腫瘤組織的差異分析，并對(duì)空間數(shù)據(jù)聚類(lèi)識(shí)別腫瘤的空間異質(zhì)性。
3、整合單細(xì)胞空間的數(shù)據(jù)，識(shí)別空間每個(gè)區(qū)域的細(xì)胞類(lèi)型組成（這個(gè)單細(xì)胞空間聯(lián)合的方法是SPOTlight），ST組織中鑒定出B細(xì)胞、T細(xì)胞、NK細(xì)胞、漿細(xì)胞、髓樣細(xì)胞、腫瘤細(xì)胞、成纖維細(xì)胞和內(nèi)皮細(xì)胞。這些細(xì)胞簇的比例和分?jǐn)?shù)在每個(gè)區(qū)域都不同。細(xì)胞類(lèi)型在癌區(qū)和癌旁有不同的分布模式。

第三步，研究原發(fā)性和肝轉(zhuǎn)移腫瘤細(xì)胞的異質(zhì)性

注意這里是研究同一種細(xì)胞類(lèi)型---腫瘤細(xì)胞，不同條件下的異質(zhì)性，之前的內(nèi)容講到過(guò)，轉(zhuǎn)錄因子（TF）分析和GSVA分析主要針對(duì)這種情況，以此來(lái)鑒別腫瘤細(xì)胞前后發(fā)生的調(diào)控變化。

微環(huán)境重塑髓系細(xì)胞組成

在腫瘤微環(huán)境(TME)中，髓系細(xì)胞是由單核/巨噬細(xì)胞、傳統(tǒng)樹(shù)突狀細(xì)胞(cdc)和肥大細(xì)胞組成的異質(zhì)性亞群。髓系細(xì)胞在不同部位的分布模式不同，具有器官特異性。

髓系細(xì)胞免疫圖譜

CXCL13+ T cells are enriched in liver metastatic tumors

在適應(yīng)性免疫應(yīng)答中，T細(xì)胞尤其是CD8+ T細(xì)胞和CD4+ T細(xì)胞發(fā)揮主導(dǎo)作用。分析中發(fā)現(xiàn)CD8+ T細(xì)胞和CD4+ T細(xì)胞都包含一個(gè)表達(dá)CXCL13的cluster，CXCL13是已知與CXCR5相互作用的趨化因子。與正常肝組織(LN)相比，CD8_CXCL13細(xì)胞在肝轉(zhuǎn)移(LM)中顯著升高，而在正常結(jié)直腸黏膜(CN)、LN和外周血(PB)中很少檢測(cè)到。進(jìn)一步分析CXCL13+ T細(xì)胞亞群的特征。基因集變異分析(Gene set variation analysis, GSVA)顯示CD8_CXCL13亞群富集T細(xì)胞增殖通路，而CD4_CXCL13亞群富集G2-M檢查點(diǎn)通路，表明這兩個(gè)亞群的增殖特性。對(duì)每個(gè)cluster內(nèi)的基因表達(dá)進(jìn)行的分析表明，CD8_CXCL13細(xì)胞高水平表達(dá)ITGAE，這是組織駐留記憶T (TRM)細(xì)胞的標(biāo)志物。

CD8_CXCL13和CD4_CXCL13細(xì)胞在結(jié)直腸癌LM中表達(dá)上調(diào)，其機(jī)制可能與高增殖能力有關(guān)。

CXCL13+ T細(xì)胞與結(jié)直腸癌患者的良好預(yù)后相關(guān)

在TME中富集的CXCL13+ T細(xì)胞是一個(gè)腫瘤反應(yīng)性的亞群，有助于改善CRC患者的預(yù)后。

Distinct subsets of fibroblasts exist in the primary and liver metastatic tumors of CRC

成纖維細(xì)胞是腫瘤微環(huán)境(TME)中基質(zhì)非免疫細(xì)胞的主要類(lèi)型，在多種腫瘤中表現(xiàn)出顯著的異質(zhì)性。在這項(xiàng)研究中，研究人員旨在表征原發(fā)性結(jié)直腸癌(CC)和肝轉(zhuǎn)移(LM)腫瘤中的成纖維細(xì)胞亞群。觀(guān)察到所有的成纖維細(xì)胞都高表達(dá)ACTA2，這是一種廣泛認(rèn)可的癌癥相關(guān)成纖維細(xì)胞(CAFs)的標(biāo)志物。

分析揭示了結(jié)直腸癌原發(fā)腫瘤和肝轉(zhuǎn)移腫瘤之間不同的表型特征和不同比例的成纖維細(xì)胞亞群。這些發(fā)現(xiàn)強(qiáng)調(diào)了TME在不同癌癥環(huán)境中的細(xì)胞異質(zhì)性。

The F3-expressing fibroblast subset enriched in the primary tumors secretes protumor factors and is associated with poor prognosis of CRC patients

在CC中富集的F4_F3成纖維細(xì)胞亞群與腫瘤因子的分泌、參與NRG1-ERBB3通路以及CRC患者較差的預(yù)后相關(guān)。這些發(fā)現(xiàn)闡明了成纖維細(xì)胞在TME形成中的作用及其對(duì)腫瘤進(jìn)展和患者結(jié)局的影響。

The MCAM-expressing fibroblast in TME of LM modulates the generation of CD8_CXCL13 cells through the Notch signaling pathway

研究探討了Notch信號(hào)通路在CD8_CXCL13和CD4_CXCL13亞群中對(duì)CXCL13表達(dá)的調(diào)節(jié)作用。我們發(fā)現(xiàn)Notch信號(hào)轉(zhuǎn)錄因子RBPJ在兩個(gè)亞群中均優(yōu)先表達(dá)。RBPJ在LM中的表達(dá)與CXCL13、ITGAE呈正相關(guān)，而在CC中未觀(guān)察到。CD8_CXCL13和CD4_CXCL13亞群主要表達(dá)NOTCH1受體。對(duì)Notch配體的分析表明，它們主要表達(dá)于成纖維細(xì)胞和內(nèi)皮細(xì)胞。具體而言，F(xiàn)2_MCAM子集富集了JAG1, F5_COCH子集富集了DLL1, E2_DLL4子集富集了DLL4、JAG1和JAG2。

交互作用分析顯示，E2_DLL4內(nèi)皮細(xì)胞亞群與CXCL13+ T細(xì)胞的交互作用最強(qiáng)，而F2_MCAM成纖維細(xì)胞亞群與CD8_CXCL13和CD4_CXCL13通過(guò)JAG1-NOTCH1相互作用。研究人員推測(cè)，F(xiàn)2_MCAM亞群由于其在TME中的分布模式，參與了CXCL13+ T細(xì)胞中Notch信號(hào)的激活?；颊叻治鲲@示，F(xiàn)2_MCAM-high組在LM中CD8_CXCL13亞群比例更高。同樣，F(xiàn)2_MCAM在LM中的浸潤(rùn)評(píng)分與CD8_CXCL13亞群的浸潤(rùn)評(píng)分相關(guān)，而在CC中則不相關(guān)。

空間轉(zhuǎn)錄組學(xué)(ST)分析證實(shí)了F2_MCAM和CD8_CXCL13在LM中的浸潤(rùn)相關(guān)。此外，與CC相比，Notch信號(hào)通路在LM中的相互作用強(qiáng)度更高，可能是由于在LM中F2_MCAM亞群的比例更高。此外，本研究還通過(guò)預(yù)測(cè)RBPJ在CXCL13啟動(dòng)子上的結(jié)合位點(diǎn)，探討了Notch信號(hào)通路對(duì)CXCL13表達(dá)的潛在調(diào)節(jié)作用，提示RBPJ可能是一個(gè)影響CXCL13表達(dá)的轉(zhuǎn)錄因子。

研究結(jié)果表明Notch信號(hào)通路在CD8_CXCL13和CD4_CXCL13亞群中調(diào)節(jié)CXCL13的表達(dá)，其中F2_MCAM成纖維細(xì)胞亞群在CXCL13+ T細(xì)胞中激活Notch信號(hào)通路中起關(guān)鍵作用。成纖維細(xì)胞、Notch信號(hào)和CXCL13表達(dá)之間的相互作用為T(mén)ME的復(fù)雜性及其對(duì)抗腫瘤免疫的影響提供了見(jiàn)解。

Intercellular interaction network in ST tissues（通訊分析居然放在了最后，分析方法采用了cellphoneDB）

原發(fā)腫瘤和肝轉(zhuǎn)移腫瘤間質(zhì)室內(nèi)的細(xì)胞相互作用。VEGFA-NRP1和VEGFA-NRP2配體-受體對(duì)在原發(fā)性和肝轉(zhuǎn)移性腫瘤中均有富集，表明它們參與了腫瘤的發(fā)生和進(jìn)展。此外，分析發(fā)現(xiàn)原發(fā)腫瘤(CC)和肝轉(zhuǎn)移腫瘤(LM)之間富集的配體-受體對(duì)不同。具體而言，ERBB3-NRG1對(duì)在原發(fā)腫瘤亞群C1和C3中富集，而在肝轉(zhuǎn)移腫瘤亞群L1和L2中缺失，表明其在原發(fā)腫瘤的發(fā)展和轉(zhuǎn)移中具有潛在作用。

方法需要注意的地方

SCENIC分析

最后，附上圖中繪制堆疊小提琴圖的代碼

source 腳本

library(Seurat)
library(Rmisc)
library(ggplot2)

####**** function part ****####
stacked_violin_plot=function(gene,seurat_object,cluster=NULL,limits.max=7,
         width=13,height=10.3,flip=T,filename="",text.size=10,Mean=T,
         col=colors){
  if(length(cluster)>0){
    seurat_object=subset(seurat_object,idents = cluster)
  }
  ave_expression=AverageExpression(seurat_object,assays = "RNA")$RNA
  ave_expression=log2(ave_expression+1)
  data_matrix=seurat_object@assays$RNA@data
  plot.list=list()
  #g=gene[5]
  no=1
  gene=c(gene[1],gene,gene[length(gene)])
  for (g in gene) {
      ave_gene_choose=ave_expression[which(rownames(ave_expression)==g),]
      data_matrix_choose=data_matrix[which(rownames(data_matrix)==g),]
      #data_matrix_choose=as.data.frame(data_matrix_choose)
      df=data.frame(expression=as.numeric(data_matrix_choose),cluster=as.character(seurat_object@active.ident))
      mean=vector()
      ave_gene_choose=as.data.frame(ave_gene_choose)
      if(Mean==T){
        for (i in df$cluster) {
          mean=c(mean,ave_gene_choose[i,])
        }
        df[,"mean"]=as.data.frame(mean)
        df=as.data.frame(df)
      }
    
  if(!flip){
    df[,"cluster"]=factor(df$cluster,levels=levels(seurat_object))
    if(Mean==F){
      p <- ggplot(df, aes(x=cluster, y=expression, fill= cluster, color=cluster))+
        geom_violin(scale="width") +
        labs(title=paste(g), y ="Expression", x="Cluster")+
        #theme_classic() +
        scale_fill_manual(values = col)+
        scale_color_manual(values = col)+
        theme(axis.title.y =  element_blank())+
        #theme(axis.ticks.y =  element_blank())+
        #theme(axis.line.y =   element_blank())+
        #theme(axis.text.y =   element_blank())+
        theme(axis.title.x = element_blank())+
        theme(legend.position="none" )
    }else{
      p <- ggplot(df, aes(x=cluster, y=expression, fill= mean, color=mean))+
        geom_violin(scale="width") +
        labs(title=paste(g), y ="Expression", x="Cluster")+
        #theme_classic() +
        scale_color_gradientn(colors = c("#FFFF00", "#FFD000","#FF0000","#360101"),
                              limits=c(0,limits.max))+
        scale_fill_gradientn(colors = c("#FFFF00", "#FFD000","#FF0000","#360101"),
                             limits=c(0,limits.max))+
        theme(axis.title.y =  element_blank())+
        #theme(axis.ticks.y =  element_blank())+
        #theme(axis.line.y =   element_blank())+
        #theme(axis.text.y =   element_blank())+
        theme(axis.title.x = element_blank())+
        theme(legend.position="none" )
    }
    
    if(no!=length(gene)){
      p<-p+
        theme( axis.line.x=element_blank(),
               axis.text.x = element_blank(),
               axis.ticks.x = element_blank())
    }else{
      p<-p+
        theme(axis.text.x = element_text(size = 10,vjust = 0.5,face = "bold",color = "black"))
    }
    #p<-p+theme(plot.title = element_text(size=text.size,face="bold",hjust = 0.5))
    p=p+theme(panel.border = element_rect(fill = "NA",size = 0.5,color = "black"),
              panel.background = element_blank(),
              panel.grid = element_blank())
    if(no==1){
      plot.margin=unit(c(0.2, 0.5, 0.2, 0.5), "cm")
      p=p+theme(legend.title = element_text(size = 9,face = "bold"), 
            legend.text  = element_text(size = 9,face = "bold"),
            legend.key.size = unit(1, "lines"))
      legend.position="none"
    }else if(no==length(gene)){
      plot.margin=unit(c(-0.3, 0.5, 0.2, 0.5), "cm")
      legend.position="none"
    }else{
      plot.margin=unit(c(-0.74, 0.5, 0, 0.5), "cm")
      legend.position="none"
    }
    p=p+xlab("") + ylab(g) + ggtitle("") +
      theme(legend.position = legend.position,
            #axis.text.x = element_blank(),
            axis.text.y = element_blank(),
            #axis.ticks.x = element_blank(),
            axis.ticks.y = element_blank(),
            axis.title.y = element_text(size=text.size,face="bold",hjust = 0.5),
            plot.margin = plot.margin )
    
    if(length(plot.list)==0){
      plot.list=list(p)
    }else{
      plot.list=c(plot.list,list(p))
    }
    no=no+1
  }else{
    ####******########
    df$cluster=factor(df$cluster,levels=rev(levels(seurat_object)))
    if(Mean==F){
      p <- ggplot(df, aes(x=cluster, y=expression, fill= cluster, color=cluster))+
        geom_violin(scale="width") +
        labs(title=paste(g), y ="Expression", x="Cluster")+
        #theme_classic() +
        scale_fill_manual(values = col)+
        scale_color_manual(values = col)+
        theme(axis.title.y =  element_blank())+
        #theme(axis.ticks.y =  element_blank())+
        #theme(axis.line.y =   element_blank())+
        #theme(axis.text.y =   element_blank())+
        theme(axis.title.x = element_blank())+
        theme(legend.position="right")
    }else{
      p <- ggplot(df, aes(x=cluster, y=expression, fill= mean, color=mean))+
        geom_violin(scale="width") +
        labs(title=paste(g), y ="Expression", x="Cluster")+
        #theme_classic() +
        scale_color_gradientn(colors = c("#FFFF00", "#FFD000","#FF0000","#360101"),
                              limits=c(0,limits.max))+
        scale_fill_gradientn(colors = c("#FFFF00", "#FFD000","#FF0000","#360101"),
                             limits=c(0,limits.max))+
        theme(axis.title.y =  element_blank())+
        #theme(axis.ticks.y =  element_blank())+
        #theme(axis.line.y =   element_blank())+
        #theme(axis.text.y =   element_blank())+
        theme(axis.title.x = element_blank())+
        theme(legend.position="right")
    }

    
    if(no!=length(gene)){
      p<-p+
        theme( axis.line.x=element_blank(),
               axis.text.x = element_blank(),
               axis.ticks.x = element_blank())
    }else{
      p<-p+
        theme(axis.text.x = element_text(size = 10,vjust = 0.2,face = "bold",color = "black"))
    }
    p<-p+theme(plot.title = element_text(size=4,face="bold",hjust = 0.5,color = "black"))
    p=p+theme(panel.border = element_rect(fill = "NA",size = 0.5,color = "black"),
              panel.background = element_blank(),
              panel.grid = element_blank())
    if(no==1){
      plot.margin=unit(c(0.5, 0.5, 0.5, 0.2), "cm")
      p=p+theme(legend.title = element_text(size = 9,face = "bold"), 
                legend.text  = element_text(size = 9,face = "bold"),
                legend.key.size = unit(0.5, "lines"))
      #legend.position="left"
    }else if(no==length(gene)){
      plot.margin=unit(c(0.5, 0.2, 0.5, 0.3), "cm")
      #legend.position="none"
    }else{
      plot.margin=unit(c(0.5, 0, 0.5, -0.11), "cm")
      #legend.position="none"
    }
    if(no==1){
      p=p+xlab("") + ylab("") +
        theme(legend.position = "none",
              axis.text.x = element_blank(),
              #axis.text.y = element_blank(),
              axis.ticks.x = element_blank(),
              #axis.ticks.y = element_blank(),
              #axis.title.x = element_text(size=10,face="bold",hjust = 0.5),
              plot.title = element_text(colour = "black", face = "bold", 
                                        size = text.size, vjust = 0.2),
              axis.text.y = element_text(size=10,face="bold",hjust = 1,color = "black"),
              plot.margin = plot.margin )+coord_flip()
    }else if(no==length(gene)){
      p=p+xlab("") + ylab("") +
        theme(legend.position = "right",
              axis.text.x = element_blank(),
              axis.text.y = element_blank(),
              axis.ticks.x = element_blank(),
              axis.ticks.y = element_blank(),
              #axis.title.x = element_text(size=10,face="bold",hjust = 0.5),
              plot.title = element_text(colour = "black", face = "bold", 
                                        size = text.size, vjust = 0.2),
              #axis.text.y = element_text(size=10,face="bold",hjust = 1,color = "black"),
              plot.margin = plot.margin )+coord_flip()
    }else{
      p=p+xlab("") + ylab("") +
        theme(legend.position = "none",
              axis.text.x = element_blank(),
              axis.text.y = element_blank(),
              axis.ticks.x = element_blank(),
              axis.ticks.y = element_blank(),
              #axis.title.x = element_text(size=10,face="bold",hjust = 0.5),
              plot.title = element_text(colour = "black", face = "bold", 
                                        size = text.size, vjust = 0.2),
              #axis.text.y = element_text(size=10,face="bold",hjust = 1,color = "black"),
              plot.margin = plot.margin )+coord_flip()
    }
    
    
    if(length(plot.list)==0){
      plot.list=list(p)
    }else{
      plot.list=c(plot.list,list(p))
    }
    no=no+1
  }
  }
  pdf(paste0(filename,"_vln_manual",".pdf"),width = width,height = height)
  if(flip){
    m=multiplot(plotlist = plot.list,cols = length(gene))
  }else{
    m=multiplot(plotlist = plot.list,cols=1)
  }
  dev.off()
}

用法

stacked_violin_plot(gene = c("CD40","CD80","CD86"),col = colors_myeloid,
                    filename = "Figure 3H.CD40",seurat_object = myeloid.combined,
                    Mean = T,flip = F,width = 10,height = 4,legend.position = "right",
                    limits.max = 3)

文章所有的代碼放在了這里，大家隨意查看下載

鏈接：https://pan.baidu.com/s/1jlVm6Ojx9Tx-Sh7Zyrv31g?pwd=ZyfO
提取碼：ZyfO

生活很好，有你更好