单细胞转录组实战02：数据整理与之质控.docx

单细胞转录组实战02：数据整理与之质控 整理数据 从cellranger的输出目录中读取filtered_feature_bc_matrix.h5表达量矩阵，并把多个样本合并为1个anndata对象。 导库建立工作目录 from?pathlib?import?Path import?scanpy?as?sc OUTPUT_DIR=output/01.preprocess Path(OUTPUT_DIR).mkdir(parents=True,exist_ok=True) adata_dict?=?{} for?i?in?Path(quantify).glob(**/filtered_feature_bc_matrix.h5): adata?=?sc.read_10x_h5(i) adata.var_names_make_unique() adata.obs_names_make_unique() adata_dict[str(i.parent).split(/)[1]]?=?adata adata?=?sc.concat(adata_dict,label=Sample,axis=0) adata.obs_names_make_unique() Quality control 数据质控包括过滤细胞和过滤基因。 过滤细胞 评估一个细胞的质量：总的UMIs数量, 基因数量, 线粒体比例 过滤基因 一个基因至少在几个细胞中表达 简单的过滤一下细胞和基因，一个细胞至少表达300个基因，一个基因至少在10个细胞中表达。 sc.pp.filter_cells(adata,min_genes=300) sc.pp.filter_genes(adata,min_cells=10) 计算线粒体基因的比例，线粒体基因以MT开头 adata.var[Mito]?=?adata.var_names.str.startswith(rMT-,rmt-) sc.pp.calculate_qc_metrics(adata,?qc_vars=[Mito],?percent_top=None,?log1p=False,?inplace=True) 可视化count数量、基因数量、线粒体基因百分比 ks?=?(total_counts,n_genes_by_counts,pct_counts_Mito) _,?axes?=?plt.subplots(1,?3,?figsize=(6,?3)) axes?=?axes.flatten() for?a,k?in?enumerate(ks): sc.pl.violin(adata,?keys=k,jitter=False,show=False,ax=axes[a],ylabel=) plt.subplots_adjust(wspace?=?0.5); 每个样本的count数量、基因数量、线粒体基因百分比 _,?axes?=?plt.subplots(1,?2,?figsize=(6,?3)) for?a,k?in?enumerate(ks[1:]): sc.pl.scatter(adata,?x=total_counts,?y=k,color=pct_counts_Mito,ax=axes[a],legend_loc=right?margin,show=False) plt.subplots_adjust(wspace?=?0.3); _,?axes?=?plt.subplots(1,?3,?figsize=(6,?3)) axes?=?axes.flatten() for?a,k?in?enumerate(ks): sc.pl.violin(adata,keys=k,groupby=Sample,rotation=45,jitter=False,show=False,stripplot=False,ax=axes[a]) plt.subplots_adjust(wspace?=?0.5); 细胞中有较高的线粒体比例有两种可能，1是该细胞受到刺激或者是正在死亡的细胞，2是细胞活性比较强。 简单的过滤 线粒体基因比例小于20%，n_genes_by_counts小于整体的98%，total_counts_upper小于整体的98% n_genes_by_counts_upper_lim?=?np.quantile(adata.obs.n_genes_by_counts.values,?.98) total_counts_upper_lim?=?np.quantile(adata.obs.total_counts.values,?.98) afilter?=?{ pct_counts_Mito:x20, n_genes_by_