
你的位置:pg电子麻辣火锅爆浆 > 新闻动态 >
发布日期:2025-03-07 13:02 点击次数:197
上次给大家简单整理了一下细胞鉴定曲线图理解,里面使用nCount_RNA或者nFeature_RNA在R语言里面绘制细胞鉴定曲线,找到一个合适的cutoff值,进行了一个初步的质控。
图片
结尾也提到了,很少会有下游是原始的rawcounts的数据,一般我们都是使用cellranger质控后的数据进行分析。不过对于处理后的数据集我们可以可视化一下nFeature_RNA和nCount_RNA来辅助进行质控
那首先我们基于Seurat官网的教程来了解回顾一下nFeature_RNA和nCount_RNA,并且可视化判断一下阈值,然后了解一下实际分析情况中的应用。
nFeature_RNA和nCount_RNA简介创建完seurat对象之后,在不进行任何操作时,seurat会为每个细胞创建一个元数据,保存在meta.data里面
#读取数据创建seurat对象pbmc.data <- Read10X(data.dir = "./filtered_gene_bc_matrices/hg19/")pbmc <- CreateSeuratObject(counts = pbmc.data, project = "pbmc3k", min.cells = 3) > dim(pbmc)[1] 13714 2700
图片
每一列的内容:
orig.ident:通常包含所知的样品名,默认为我们赋给project的值,如果不赋值那就是SeuratProjectnCount_RNA:每个细胞的UMI数目nFeature_RNA:每个细胞所检测到的基因数目可以看到nCount_RNA和nFeature_RNA还是有差异的,这就与它们的计算方法有关
图片
#nCount_RNA:总的UMI数即转录本数量colSums(sce@assays$RNA$counts)#nFeature_RNA:总的基因数目colSums(sce@assays$RNA$counts>0)可视化及阈值判断
可以使用小提琴图来简单可视化一下nFeature_RNA和nCount_RNA
VlnPlot(pbmc, features = c("nFeature_RNA", "nCount_RNA"))
图片
过滤前nFeature_RNA图:反映的是样品中每个细胞表达的基因数量,表达过高可能是双细胞或者多细胞,表达过低可能是空液滴或者包裹的是环境RNA
nCount_RNA图:反映的是每个细胞中包含的UMI数量也就是转录本的数量
在10X Genomics测序数据分析过程中,通过UMI对测序得到的reads进行简并之后,就可以看到一个细胞中被读到多少个基因。一般一个细胞可以得到40000-80000个有效的UMI,平均一个细胞的一个基因有10个左右的UMI。
所以我们在进行阈值判断的时候,可以直接基于nFeature_RNA值也就是基因的数量
图片
阈值判断We filter cells that have unique feature counts over 2,500 or less than 200
官网给的是大于200和小于2500,但可视化之后我们可以看到上限定在2000其实也可以。
不过pbmc是比较早期的数据了,测到的细胞数量比较少,上限设置的也比较低,如果是现在的单细胞数据还是要具体数据具体分析
#基于阈值过滤并且可视化pbmc <- subset(pbmc, subset = nFeature_RNA > 200 & nFeature_RNA < 2000)VlnPlot(pbmc, features = c("nFeature_RNA", "nCount_RNA"))> dim(pbmc)[1] 13714 2692
图片
过滤后过滤后,细胞从最开始的2700变为现在的2692,过滤掉了部分细胞。
以上是seurat官网pbmc3k_tutorial中QC的部分内容,接下来我们看看在实际数据中的应用。
实际分析中应用如果大家手里有技能树单细胞分析的标准分析代码,如果需要的话可以获取一下链接: https://pan.baidu.com/s/1bIBG9RciAzDhkTKKA7hEfQ?pwd=y4eh
那在我们的scRNA_scripts文件夹中有个qc.R的质控脚本文件,就是对读取进来的数据进行质控的。
图片
脚本函数首先是计算了线粒体、核糖体以及血红细胞的比例(下期给大家详细介绍),然后就可视化了细胞中这些参数的情况。我们还是先重点看看nFeature_RNA和nCount_RNA
#qc.R脚本中nFeature_RNA和nCount_RNA部分内容feats <- c("nFeature_RNA", "nCount_RNA")p1=VlnPlot(input_sce, group.by = "orig.ident", features = feats, pt.size = 0, ncol = 2) + NoLegend()p1 w=length(unique(input_sce$orig.ident))/3+5;wggsave(filename="Vlnplot1.pdf",plot=p1,width = w,height = 5)
图片
质控前一般走标准流程的时候,在创建seurat对象时候就会基于min.cells = 5和min.features = 300进行过滤,所以在qc脚本中是不进行这一步的过滤操作的。不过为了看一下过滤前后变化,我们可以基于可视化的结果进行一个简单的过滤操作。
图片
#简单过滤 if(T){ selected_c <- WhichCells(input_sce, expression = nFeature_RNA > 500 & nFeature_RNA < 2500) selected_f <- rownames(input_sce)[Matrix::rowSums(input_sce@assays$RNA$counts > 0 ) > 3] input_sce.filt <- subset(input_sce, features = selected_f, cells = selected_c) dim(input_sce) dim(input_sce.filt) } #可视化过滤后的情况 feats <- c("nFeature_RNA", "nCount_RNA") p1_filtered=VlnPlot(input_sce.filt, group.by = "orig.ident", features = feats, pt.size = 0, ncol = 2) + NoLegend() w=length(unique(input_sce.filt$orig.ident))/3+5;w ggsave(filename="Vlnplot1_filtered.pdf",plot=p1_filtered,width = w,height = 5)
图片
过滤后基本质控意义:可以去除掉每个样品中,一些表达量过高或者过低的基因。
除了在基本质控环节我们会可视化一下细胞中nFeature_RNA和nCount_RNA的情况,在进行降维聚类分群的时候,我们也会对nFeature_RNA和nCount_RNA进行可视化。
细胞降维聚类分群中应用在选择对应的阈值进行可视化的时候,我们会用到check-all-markers.R脚本,基于常见Marker基因进行一下可视化,以及绘制umap图
图片
在check-all-markers.R脚本,帮助我们查看确认每个细胞亚群中基因的表达情况,从而帮助我们判断是否是双细胞。
具体推文:如何排除双细胞
我们在进行亚群简单命名的时候,一般选择比较低的分辨率0.1,那在GSE208706数据的0.1分群中,我们可以很明显的看到第9群比较狭长,且包含了两个不同细胞亚群的Marker基因。
图片
为了判断是否是双细胞,我们就需要结合每个亚群的单个细胞的总的RNA数量进行判断
if("percent_mito" %in% colnames(sce.all.int@meta.data ) ){ #可视化细胞的上述比例情况 feats <- c("nFeature_RNA", "nCount_RNA", "percent_mito", "percent_ribo", "percent_hb") feats <- c("nFeature_RNA", "nCount_RNA") p1=VlnPlot(sce.all.int , features = feats, pt.size = 0, ncol = 2) + NoLegend() w=length(unique(sce.all.int$orig.ident))/3+5;w ggsave(filename=paste0(pro,"Vlnplot1.pdf"),plot=p1,width = w,height = 5) feats <- c("percent_mito", "percent_ribo", "percent_hb") p2=VlnPlot(sce.all.int, features = feats, pt.size = 0, ncol = 3, same.y.lims=T) + scale_y_continuous(breaks=seq(0, 100, 5)) + NoLegend() w=length(unique(sce.all.int$orig.ident))/2+5;w ggsave(filename=paste0(pro,"Vlnplot2.pdf"),plot=p2,width = w,height = 5) }
图片
nFeature_RNA可视化结果发现反而第8群表达量高,而第9群正常。基于Marker基因推断第8群是处于增殖期的细胞,所以表达量高是合理的。
并且提高分辨率之后,发现9群被细分为两个亚群,也不是双细胞。
图片
一般我们会根据中位线以及最高值来进行判断,再提高分辨率看亚群有没有分开,再确定是否是双细胞。
线粒体比例在官网以及我们的标准质控流程中,都会计算线粒体比例
图片
我们的qc.R脚本中还对核糖体以及血红细胞的比例进行了计算和可视化,那下期一起来了解一下这些内容吧!
图片
本站仅提供存储服务,所有内容均由用户发布,如发现有害或侵权内容,请点击举报。