笔记|R语言画图 | 如何看已知基因list的细胞类型特异性表达（）开发语言|r语言

自己学东西，觉得举步维艰。网上的资源太多了，反而让人很容易陷入一种焦躁和慌乱中。不知道从哪里理出头绪来。是碎片化的，不成体系的。自己需要对信息进行过滤，整合成自己的方法。

目标：加载作者的处理过的有细胞类型标签的数据，对数据进行差异表达分析，找到特定的组织中差异表达的基因。
注明以下使用到的数据来源：http://development.psychencode.org/files/processed_data/scRNA-seq/Sestan.fetalHuman.Psychencode.Rdata

该数据为这篇文章中作者用到的，胎儿的大脑数据。
在Rdata中存放有三个数据矩阵。

count2：原始count值

cpm2：归一化处理后的数据矩阵

meta2：关于细胞与样本的详细注释信息，里面有作者已经分类后的每一个细胞对应的数据类型
整体上这个数据有60155行，762列。

我们需要在上述已有数据的基础上，实现我们的目标。

load("F:/Rworkplace/celltype_specifc/Sestan.fetalHuman.Psychencode.Rdata") library(dplyr) library(Seurat) library(patchwork) dim(count2) [1] 60155 762 data<-CreateSeuratObject(counts=count2, project = "fetal", meta.data = https://www.it610.com/article/meta2, min.cells=0, min.features=0) #data@meta.data<-meta2 data<-NormalizeData(data,normalization.method="RC",scale.factor = 1000000) #这里就是计算了cpm的值 #这里实际上进行了cpm的转换#筛选标准：该基因至少在300个细胞中表达，每一个细胞中至少有300个基因表达，对这个matrix进行行列层面的筛选。

将原先meta.data中的celltype的变量赋值到Idents中。得到细胞的类型标签。

Idents(object=data)<-data@meta.data$ctype #这一步是成功的第一步，已有的细胞注释的结果 levels(x=data)

到这一步的话，似乎就是利用原始的数据来做的。

我之所以举步维艰，是因为：
（1）我不知道每一步的含义，是否可以视特殊情况而去除；
（2）我不知道更多的函数应该怎样使用；

接着进行差异表达基因分析。

markers<-FindAllMarkers(data,only.pos = TRUE ,logfc.threshold = 2,min.pct = 0.25,test.use="roc") write.csv(markers,"DEG_all2.csv") VlnPlot(object=data,feature=c("ENSG00000189238-LINC00943"),log = T)

到这里虽然还不知道如何使用cpm后的数据进行差异表达分析。
但是比较欣慰的是，可以确保找到的差异表达基因与作者的基因是一样的了。即，得到的DEG_all.csv这个list是有意义的。
使用roc进行差异表达分析，似乎效果更好。【首选：有时间要总结一下差异表达分析常用的方法】
得到这个list之后，我们与我们自己的list进行比照，看看我们的list中，有多少是在这个细胞类型特异性的差异表达基因中的。

DEG<-read.csv("DEG_all2.csv") library(stringr) gene_split<-as.data.frame(str_split_fixed(DEG$gene,"-",n=2)) head(gene_split) head(DEG) DEG_m<-cbind(DEG,gene_split) interest<-read.table("result_coding_gene_id_2.txt") DEG_interest<-merge(DEG_m,interest) #112/4205 write.csv(DEG_interest$gene,"DEG_interest_celltype.csv")

得到了一些基因在特定的细胞类型中表达。

table(DEG_interest_2$cluster) 112/4205

Astro Endo ExN1 ExN2 ExN3 InN1 IPC1
6 26 2 7 9 2 4
IPC2 Microglia NasN1 NasN5 NasN6 NEPRGC1 NEPRGC2
1 20 2 1 1 1 2
NEPRGC3 Oligo OPC1 Pericyte
2 4 2 20

根据这些基因，绘制heatmap。

DEG_interest_2<-DEG_interest[order(DEG_interest$cluster,decreasing =FALSE ),] DoHeatmap(data,features = DEG_interest_2$gene,draw.lines = TRUE,label = T,size = 2)+NoLegend()

可能这个工具还没有用熟悉，所以画出来的图是真的很丑（后面再进行润色和调整）。

文章图片
保存这次的处理结果。

DEG_interest_2<-DEG_interest[order(DEG_interest$cluster,decreasing =FALSE ),] a<-subset(DEG_interest_2,select=c("gene","cluster")) tissue<-rep("fetal",dim(DEG_interest_2)[1]) d<-cbind(a,tissue) colnames(d)<-c("gene","cell_type","label") write.csv(d,"DEG_interest_celltype2.csv")

最后测试一下分类标签是否正确。

VlnPlot(object=data,feature=c("ENSG00000122547-EEPD1"),log = T)

文章图片
完成目标。

现在存在的疑惑就是，是否需要吧ExN1、ExN2、ExN3这些都属于ExN细胞的细胞类型合并起来进行差异表达分析。
现在先这样，可能后期，还有再加上这些大的class（如，ExN/NasN/IPC）所共有的相对于其他class而言的差异基因。

接下来想处理adult的数据，没想到遇到了不能避免的数据太大的问题。

load("F://Rworkplace//celltype_specifc//phyEncode_adult//Sestan.adultHumanNuclei.Psychencode.Rdata") library(dplyr) library(Seurat) library(patchwork)

报错：