单细胞转录组的细分亚群的降维聚类分群加上gsea或者gsva以及转录因子和拟时序流程（仅需8000元）

生信技能树jimmy

发布于 2021-10-20 13:23:01

2.3K0

发布于 2021-10-20 13:23:01

文章被收录于专栏：单细胞天地单细胞天地

网罗了一大波生物信息学数据分析方面的工程师，是时候官宣咱们的单细胞数据处理业务列表：

现在介绍10x的fastq文件的数据分析流程

我们的单细胞产品线还是蛮丰富的：

明码标价之10X技术单细胞（2.5万每个样品）（标准100G测序数据）
明码标价之10X转录组原始测序数据的cellranger流程（单个项目收费800起）
明码标价之单细胞转录组的质控降维聚类分群和生物学注释（单个项目收费800起）
扫描下面二维码即可添加微信咨询！（添加好友务必备注高校或者工作单位+姓名，方便后续认识）

只需8000元，我们会提供如下所示的全部分析流程和图表，包括代码！

细胞图谱时代早就过去了，不再是随便选取一个物种，一种组织或者一种疾病，挑选几个样品做单细胞，简单的说明清楚其单细胞亚群组成比例和生物学意义就足够的时机了。但是图谱既然是被做烂了的思路，就应该是有系统性梳理的价值，比如肿瘤样品的单细胞就应该是首先是按照如下所示的标记基因进行第一次分群：

immune (CD45+,PTPRC),
epithelial/cancer (EpCAM+,EPCAM),
stromal (CD10+,MME,fibo or CD31+,PECAM1,endo)

然后每个亚群进行第二层次细分亚群，甚至第三层次，第四次分群，结构清晰明了。我们以上皮细胞亚群的细分来举例说明每个分析点的工作量：

现在我们要介绍的是发表于2020的文章，标题是：《Single-cell transcriptomic architecture and intercellular crosstalk of human intrahepatic cholangiocarcinoma》，数据集在；https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE138709

是肝癌里面的 intrahepatic cholangiocarcinoma 这个细分疾病，本研究共5.6万个细胞，取样策略如下，

单细胞实验取样策略

可以看到研究者的第一层次降维聚类分群并不是我们前面的提到的上皮，免疫和间质细胞，而是直接进行了第二层次分群，一般来说背诵下面的基因即可：

# T Cells (CD3D, CD3E, CD8A), 
# B cells (CD19, CD79A, MS4A1 [CD20]), 
# Plasma cells (IGHG1, MZB1, SDC1, CD79A), 
# Monocytes and macrophages (CD68, CD163, CD14),
# NK Cells (FGFBP2, FCG3RA, CX3CR1),  
# Photoreceptor cells (RCVRN), 
# Fibroblasts (FGF7, MME), 
# Endothelial cells (PECAM1, VWF). 
# epi or tumor (EPCAM, KRT19, PROM1, ALDH1A1, CD24).
#   immune (CD45+,PTPRC), epithelial/cancer (EpCAM+,EPCAM), 
# stromal (CD10+,MME,fibo or CD31+,PECAM1,endo) 

library(ggplot2) 
genes_to_check = c('PTPRC', 'CD3D', 'CD3E', 'CD4','CD8A',
                   'CD19', 'CD79A', 'MS4A1' ,
                   'IGHG1', 'MZB1', 'SDC1',
                   'CD68', 'CD163', 'CD14', 
                   'TPSAB1' , 'TPSB2',  # mast cells,
                   'RCVRN','FPR1' , 'ITGAM' ,
                   'C1QA',  'C1QB',  # mac
                   'S100A9', 'S100A8', 'MMP19',# monocyte
                   'LAMP3', 'IDO1','IDO2',## DC3 
                   'CD1E','CD1C', # DC2
                   'KLRB1','NCR1', # NK 
                   'FGF7','MME', 'ACTA2', ## fibo 
                   'DCN', 'LUM',  'GSN' , ## mouse PDAC fibo 
                   'Amy1' , 'Amy2a2', # Acinar_cells
                   'PECAM1', 'VWF',  ## endo 
                   'EPCAM' , 'KRT19', 'PROM1', 'ALDH1A1' )

单细胞转录组数据分析的标准降维聚类分群，并且进行生物学注释后的结果。可以参考前面的例子：人人都能学会的单细胞聚类分群注释 ，我们演示了第一层次的分群。

如果你对单细胞数据分析还没有基础认知，可以看基础10讲：

本次我们介绍的重点是对上皮细胞进行降维聚类分群，以及单细胞高级分析，文章里面的主要的图表我们分开介绍：

A) tSNE plots and CNV box plots for 6 distinct epithelial cell subclusters.
B) Heatmap showing the top 10 DEGs (Wilcoxon test) in 6 epithelial cell subclusters.
C) Differences in pathway activity (scored per cell by GSVA) in 4 malignant cell sub- clusters.
D) Heatmap of the t-value for the area under the curve score of expression regulation by transcription factors, as estimated using SCENIC.
E) Violin plots showing the expression of marker genes in distinct malignant subclusters.
F) IHC staining showing the expression of CK19, TM4SF4, PROM1, CDH6, and SPINK1 in recurrent sample (ICC-20T) and adjacent liver tissue.
G) Kaplan–Meier survival curve of TM4SF4, SPINK1, and KRT19 expression using the optimal group cut-off point in patients with ICC (from TCGA).

标准降维聚类分群是6个细胞亚群（上皮细胞的细分）

其中0，1，2，3是有拷贝数变异的癌细胞，而4和5是正常水平细胞，分别是cholangiocytes和hepatocytes，如下所示：

第二层次降维聚类分群

可以很明显的看到标号为4和5的亚群在癌旁部分，而且CNV程度超级低，所以是正常细胞：

Subclusters 4 (cholangiocytes, top DEGs: FXYD2 and RHOB)
Subclusters 5 (hepatocytes, top DEGs: ASGR1 and APOE)

而有拷贝数变异的癌细胞就异质性很大了，下面是有拷贝数变异的癌细胞的4个亚群的各自介绍：

ubcluster 0 malignant cells were char- acterized by high expression levels of mesenchymal markers such as COL1A1, fibronectin, and IGFBP7, indicating epithelial- mesenchymal transition (EMT) characteristics.
Subcluster 1 displayed high expression levels of the malignancy-promoting factors S100P and FABP5;
subcluster 2 exhibited high expres- sion levels of the immune-associated genes CD74 and HLA-DRA;
subcluster 3 malignant cells, from the recurrent patient, displayed high expression levels of SPINK1

各个上皮细胞亚群的特异性基因

如下所示，各自单细胞亚群的特异性基因分别是：

展示标记基因

当然了，这样的每个亚群挑选自己生物学背景范围内的基因的策略，对纯粹的生物信息学工程师来说，非常的不友好，因为大家能背诵的基因屈指可数。

gsva对各个细胞亚群进行生物学功能数据库注释

不过，我们生信工程师有大杀器，就是生物学功能数据库注释，包括go和kegg的，通常是使用 clusterProfiler 包进行：

# 首先前面的降维聚类分群找到了 sce.markers 各个亚群的标记基因 ： 
library(clusterProfiler)
library(org.Hs.eg.db)
library(ggplot2)
ids=bitr(sce.markers$gene,'SYMBOL','ENTREZID','org.Hs.eg.db')
sce.markers=merge(sce.markers,ids,by.x='gene',by.y='SYMBOL')
gcSample=split(sce.markers$ENTREZID, sce.markers$cluster)
gcSample # entrez id , compareCluster 
xx <- compareCluster(gcSample, fun="enrichKEGG",
                     organism="hsa", 
                     pvalueCutoff= 1,
                     qvalueCutoff=1)
p=dotplot(xx) 
p+ theme(axis.text.x = element_text(angle = 45, 
                                    vjust = 0.5, hjust=0.5))
p

ggsave('compareCluster-sce.markers.pdf')

这里作者使用的gsva和gsea这样的策略，GSEA分析，我也多次讲解：

但实际上，绝大部分读者并没有去细看这个统计学原理，也不需要知道gsea分析的nes值如何计算，需要知道的是nes值的生物学意义。

而单细胞转录组数据的批量GSVA代码我在：借鉴escape包的一些可视化GSVA或者ssGSEA结果矩阵的方法 剧透过，其实跟 对单细胞表达矩阵做gsea分析的代码实现环节大同小异。无论是gsva还是gsea，本质上都是看某个生物学基因集的得分。