关于PAGA轨迹分析的一些理解

PAGA(Partition-based graph abstraction)轨迹分析功能包含在Scanpy软件中,是一种基于python的软件,对于没有python使用经验的人来说用起来还是比较吃力的。2019年发表的文献对其工作原理进行了详细阐述,下面是我在做轨迹分析时的一些理解。

在聚类分析中我们寻找的是细胞间的离散属性,也就是细胞之间的不同,而估计分析中我们寻找的是连续属性,也就是细胞之间的相同之处。PAGA是一种基于图形的分析方法,先通过 Louvain algorithm算法对细胞进行降维,生成低纬度的聚类图,基于聚类图进一步分析不同细胞类群之间的关系。下边是我用自己的数据生成的肌细胞谱系的PAGA轨迹分析结果,其实与其说是PAGA轨迹分析图,不如说是PAGA轨迹关系图

1. 从seurat到PAGA

个人比较习惯R studio的操作环境,所以细胞质控、数据预处理、注释等全部是用Seurat软件分析的(当然scanpy也有同样的功能),将处理后的数据先转为example.h5Seurat,在转为example.5had文件,PAGA可直接分析.h5ad文件.

remotes::install_github("mojaveazure/seurat-disk")  ##安装seuratdisk软件包
library(SeuratDisk)
SaveH5Seurat(example.cell, "example.cell.h5Seurat") ##example.cell直接转为h5Seurat格式保存
Convert("example.cell.h5Seurat", dest = "h5ad")  ##h5Seurat转为h5ad

同样的
#####读取h5ad文件
##h5ad转为h5seurat
Convert('example.cell.h5ad', "h5seurat",overwrite = TRUE,assay = "RNA")
example.cell <- LoadH5Seurat("example.cell.h5seurat")

2. PAGA图

左图是PAGA软件对全部细胞的分群结果,右图是轨迹关系推断(n_neighbors=10;threshold=0.7),可以看到在两张图片中不同细胞类群的位置关系还是比较一致的。轨迹关系结果仅仅只是简单的点-线图,那么如何理解这张图?首先点代表一个细胞类群,两点之间的连线代表两个细胞类群之间有关系,线的长短反映两个细胞类群在聚类图上的位置关系(并没有太重要的作用),线的粗细表示得到的轨迹关系的置信度(可以理解为显著性),线越粗,置信度越高。

 3. 不同参数

做任何分析都会有不同的参数,PAGA轨迹分析也不例外,主要的参数有

resolution, n_pcs,n_neighbors,threshold

resulution, n_pcs两个参数与用seurat分析设置的参数感觉相似,可以设置为相同的值,应该不会对结果造成太大的影响。n_neighbors是计算参数,threshold是结果展示参数(个人理解,可能不太准确)。参数的选择是具有强烈主观性的,分析结果的好坏要看是否符合预期(主观性显著)。下面展示了两个参数不同的分析结果之间的比较。可以看到n_neighbors从5增加到15并不会对结果产生太大的影响。threshold分别设置成0.01,0.03,0.07和0.1,轨迹关系结构发生了显著的变化。从0.01到0.1,轨迹关系判别结果愈发严格,当我们对复杂组织、多种细胞类型的数据进行轨迹分析,可以调高threshold,而简单细胞群体的轨迹判别,可以设置成低的threshol。在这里选择n_neighbors=10,threshold=0.03的结果(主观)。

4. 几个有用的链接/教程

scanpy官网Trajectory inference for hematopoiesis in mouse — Scanpy documentation

文献官网PAGA: graph abstraction reconciles clustering with trajectory inference through a topology preserving map of single cells | Genome Biology | Full Text

两个教程

单细胞转录组之Scanpy - 轨迹推断/拟时序分析 - 简书

轨迹分析_Scanpy进行单细胞分析及发育轨迹推断_颜不良文丑的博客-CSDN博客

<think>好的,用户需要寻找用于评估斑马鱼T细胞状态的R代码示例。首先,我需要考虑斑马鱼T细胞分析常用的方法。单细胞RNA测序(scRNA-seq)是常用的技术,因此应该涉及相关R包,如Seurat或Scanpy。不过用户明确提到R语言,所以重点放在Seurat上。 接下来,用户提到的“状态评估”可能包括细胞聚类、标记基因识别和功能注释。根据引用[3],细胞注释通常需要参考数据库或文献中的标记基因。斑马鱼作为模式生物,可能有特定的标记基因,比如cd4、cd8等。需要确保代码中包含这些基因的检测步骤。 引用[4]中的PAGA图用于发育轨迹分析,这可能与细胞状态转变有关。虽然用户没有明确提到轨迹分析,但状态评估可能包括发育阶段或激活状态的确定。因此可以考虑添加轨迹推断的步骤,比如使用monocle3或类似工具。 另外,数据预处理步骤也很重要,包括标准化、降维和聚类。需要检查代码示例是否涵盖这些基础分析步骤。同时,可视化部分如UMAP图和热图展示标记基因,能帮助用户直观评估细胞状态。 考虑到用户可能不熟悉具体流程,代码示例应尽量详细并附有注释。同时,提醒用户注意数据质量控制,如线粒体基因比例过滤,确保分析可靠性。最后,引用相关文献和工具,如引用[3]中的参考数据库和标记基因比较方法,以及引用[4]中的轨迹分析,以增强答案的可信度。 需要确认提供的代码是否符合R语言的最新包版本,避免过时的函数。此外,斑马鱼可能使用Ensembl基因ID,需提醒用户注意基因ID转换问题。最后,生成相关问题,帮助用户深入理解后续步骤,如数据整合、轨迹分析和多组学方法,引用相关文献支持。</think>以下是一个基于单细胞转录组数据分析斑马鱼T细胞状态的R代码示例,整合了细胞注释、标记基因分析和可视化流程: ```r # 加载必要的包 library(Seurat) library(clustree) library(celldex) library(SingleR) library(scater) # 1. 数据预处理 tcell.data <- Read10X(data.dir = "path/to/your/filtered_feature_bc_matrix") tcell <- CreateSeuratObject(counts = tcell.data, project = "Zebrafish_Tcells", min.cells = 3, min.features = 200) # 质量控制 tcell[["percent.mt"]] <- PercentageFeatureSet(tcell, pattern = "^mt-") tcell <- subset(tcell, subset = nFeature_RNA > 200 & percent.mt < 20) # 标准化和数据缩放 tcell <- NormalizeData(tcell) tcell <- FindVariableFeatures(tcell, selection.method = "vst", nfeatures = 2000) tcell <- ScaleData(tcell) # 2. 降维和聚类 tcell <- RunPCA(tcell, features = VariableFeatures(object = tcell)) tcell <- FindNeighbors(tcell, dims = 1:15) tcell <- FindClusters(tcell, resolution = c(0.2, 0.4, 0.6)) tcell <- RunUMAP(tcell, dims = 1:15) # 3. 细胞类型注释(参考引用[3]方法) # 加载斑马鱼参考数据集 ref <- BlueprintEncodeData() pred <- SingleR(test = as.SingleCellExperiment(tcell), ref = ref, labels = ref$label.main) # 添加注释结果 tcell$celltype <- pred$labels # 4. T细胞亚群分析 tcell.markers <- FindAllMarkers(tcell, only.pos = TRUE, min.pct = 0.25) top10 <- tcell.markers %>% group_by(cluster) %>% top_n(n = 10, wt = avg_log2FC) # 可视化标记基因 VlnPlot(tcell, features = c("cd4", "cd8a", "foxp3a", "il2ra")) FeaturePlot(tcell, features = c("lck", "rag1", "stat3")) # 5. 状态评估(激活/抑制状态) tcell <- AddModuleScore(tcell, features = list(Activation = c("il2", "ifng", "tnfa"), Exhaustion = c("pdcd1", "lag3", "ctla4")), name = "State_") # 可视化功能状态 FeaturePlot(tcell, features = c("State_Activation1", "State_Exhaustion1")) ``` ### 关键分析步骤说明: 1. **数据预处理**:包含质量控制步骤,通过线粒体基因比例过滤低质量细胞[^3] 2. **细胞注释**:采用SingleR方法整合参考数据集进行自动注释[^3] 3. **标记基因分析**:识别T细胞亚群特征基因(如cd4、cd8等) 4. **功能状态评估**:通过基因模块评分量化细胞激活/耗竭状态 ### 注意事项: - 斑马鱼基因命名可能需要转换为同源人类基因名(使用biomaRt包) - T细胞发育轨迹分析可结合PAGA方法(引用[4]) - 建议样本量>1000个细胞以保证分析可靠性[^1]
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值