单细胞分析(Signac): PBMC scATAC-seq 整合

2024-05-24 12:36

本文主要是介绍单细胞分析(Signac): PBMC scATAC-seq 整合,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

引言

在本教学指南中,我们将探讨由10x Genomics公司提供的人类外周血单核细胞(PBMCs)的单细胞ATAC-seq数据集。

加载包

首先加载 Signac、Seurat 和我们将用于分析人类数据的其他一些包。

if (!requireNamespace("EnsDb.Hsapiens.v75", quietly = TRUE))
    BiocManager::install("EnsDb.Hsapiens.v75")

library(Signac)
library(Seurat)
library(EnsDb.Hsapiens.v75)

library(ggplot2)
library(patchwork)

与 scRNA-seq 数据整合

为了更好地理解 scATAC-seq 测序数据,我们可以通过与人类外周血单核细胞(PBMC)相关的单细胞 RNA 测序(scRNA-seq)实验来对细胞进行分类。我们采用了一种跨数据类型整合和标签迁移的方法,这种方法的详细教程可以在这里找到。我们的目标是找出基因活动矩阵与 scRNA-seq 数据集之间的共同相关性模式,以此来确定两种技术中对应的生物学状态。通过这种方法,我们可以为每个细胞在 scRNA-seq 定义的每个聚类标签下得到一个分类评分。

在本例中,我们加载了一个由 10x Genomics 提供的经过预处理的人类 PBMC 的 scRNA-seq 数据集。你可以从 10x Genomics 的官方网站下载这项实验的原始数据,也可以在 GitHub 上查看构建该数据对象所使用的代码。此外,你还可以直接下载已经预处理过的 Seurat 数据对象。

# Load the pre-processed scRNA-seq data for PBMCs
pbmc_rna <- readRDS("../vignette_data/pbmc_10k_v3.rds")
pbmc_rna <- UpdateSeuratObject(pbmc_rna)

transfer.anchors <- FindTransferAnchors(
  reference = pbmc_rna,
  query = pbmc,
  reduction = 'cca'
)

## Running CCA
## Merging objects
## Finding neighborhoods
## Finding anchors
##  Found 17505 anchors

predicted.labels <- TransferData(
  anchorset = transfer.anchors,
  refdata = pbmc_rna$celltype,
  weight.reduction = pbmc[['lsi']],
  dims = 2:30
)

## Finding integration vectors
## Finding integration vector weights
## Predicting cell labels

pbmc <- AddMetaData(object = pbmc, metadata = predicted.labels)

plot1 <- DimPlot(
  object = pbmc_rna,
  group.by = 'celltype',
  label = TRUE,
  repel = TRUE) + NoLegend() + ggtitle('scRNA-seq')

plot2 <- DimPlot(
  object = pbmc,
  group.by = 'predicted.id',
  label = TRUE,
  repel = TRUE) + NoLegend() + ggtitle('scATAC-seq')

plot1 + plot2
alt

观察到基于单细胞 RNA 测序(scRNA)的分类结果与利用单细胞染色质可及性测序(scATAC-seq)数据生成的 UMAP(均匀流形近似和投影)图呈现出一致性。我们注意到,在 scATAC-seq 数据集中并未预测出血小板细胞,这是符合预期的,因为血小板是无核细胞,无法通过 scATAC-seq 技术被检测到。

为了整合我们的 scATAC-seq 聚类分析和标签转移结果,我们可以选择将每个聚类的名称重新指定为该聚类中最常见的预测标签。另外,我们也可以选择直接使用每个细胞的预测标签来进行后续分析。

# replace each label with its most likely prediction
for(i in levels(pbmc)) {
  cells_to_reid <- WhichCells(pbmc, idents = i)
  newid <- names(which.max(table(pbmc$predicted.id[cells_to_reid])))
  Idents(pbmc, cells = cells_to_reid) <- newid
}

查找细胞类型之间差异可及的峰

为了识别不同细胞聚类间的差异开放区域,我们可以进行差异可及性(DA)分析。一种简便的方法是应用 Wilcoxon 秩和检验,Presto 软件包已经内嵌了一个非常高效的 Wilcoxon 检验程序,适用于 Seurat 数据对象的分析。

我们还可以使用逻辑回归进行差异性分析,这是 Ntranos 等人在 2018 年针对单细胞 RNA 测序(scRNA-seq)数据提出的建议,同时将片段总数作为一个隐变量来考虑,以减少测序深度差异对分析结果的影响。本研究将重点比较未激活的 CD4+ T 细胞和 CD14+ 单核细胞,但这些方法也适用于任何细胞群体的比较。此外,我们还可以使用 Seurat 提供的小提琴图、特征图、散点图、热图或其他任何可视化工具来展示这些特征峰。

# change back to working with peaks instead of gene activities
DefaultAssay(pbmc) <- 'peaks'

da_peaks <- FindMarkers(
  object = pbmc,
  ident.1 = "CD4 Naive",
  ident.2 = "CD14+ Monocytes",
  test.use = 'LR',
  latent.vars = 'nCount_peaks'
)
head(da_peaks)

##                                  p_val avg_log2FC pct.1 pct.2     p_val_adj
## chr14-99721608-99741934  1.414340e-280   5.571161 0.868 0.022 1.238410e-275
## chr14-99695477-99720910  5.529507e-222   5.100310 0.797 0.021 4.841691e-217
## chr17-80084198-80086094  8.962660e-221   7.032939 0.668 0.005 7.847795e-216
## chr7-142501666-142511108 6.506936e-212   4.757277 0.754 0.029 5.697538e-207
## chr2-113581628-113594911 5.746054e-188  -5.051770 0.035 0.663 5.031302e-183
## chr6-44025105-44028184   2.328105e-179  -4.394199 0.046 0.616 2.038512e-174

plot1 <- VlnPlot(
  object = pbmc,
  features = rownames(da_peaks)[1],
  pt.size = 0.1,
  idents = c("CD4 Naive","CD14+ Monocytes")
)
plot2 <- FeaturePlot(
  object = pbmc,
  features = rownames(da_peaks)[1],
  pt.size = 0.1
)

plot1 | plot2
alt

要识别两个细胞群体间的差异可访问区域,我们还可以通过比较这两组细胞的可访问性变化倍数来进行。这种方法比执行更复杂的差异可访问性测试要快捷,但它无法考虑到如细胞间总测序深度差异等潜在因素,并且不进行任何统计检验。尽管如此,这种方法仍然是一个有用的数据快速探索手段,可以通过 Seurat 软件包中的 FoldChange() 函数来实现。

fc <- FoldChange(pbmc, ident.1 = "CD4 Naive", ident.2 = "CD14+ Monocytes")
# order by fold change
fc <- fc[order(fc$avg_log2FC, decreasing = TRUE), ]
head(fc)

##                          avg_log2FC pct.1 pct.2
## chr6-28416849-28417227     11.45411 0.067     0
## chr7-110665002-110665493   11.41080 0.073     0
## chr8-19317420-19317942     11.22146 0.061     0
## chr1-172836553-172836955   11.20312 0.058     0
## chr2-191380525-191380926   11.15811 0.056     0
## chr8-90255778-90256179     11.03921 0.050     0

仅仅依靠峰值坐标可能难以进行分析解读。我们可以通过 ClosestFeature() 函数来识别每个峰值所对应的最近基因。在对基因列表进行深入分析时,可以发现,在未激活的 T 细胞(Naive T cells)中开放的峰值通常与如 BCL11B 和 GATA3(T 细胞分化的关键调控基因)等基因相邻,而在单核细胞(monocytes)中开放的峰值则通常与如 CEBPB(单核细胞分化的关键调控基因)等基因相邻。为了进一步深入分析,我们可以对 ClosestFeature() 函数返回的基因集进行基因本体(Gene Ontology, GO)富集分析,而目前有许多 R 语言包能够执行此类分析(例如,可以参考 GOstats 包)。

open_cd4naive <- rownames(da_peaks[da_peaks$avg_log2FC > 3, ])
open_cd14mono <- rownames(da_peaks[da_peaks$avg_log2FC < -3, ])

closest_genes_cd4naive <- ClosestFeature(pbmc, regions = open_cd4naive)
closest_genes_cd14mono <- ClosestFeature(pbmc, regions = open_cd14mono)

head(closest_genes_cd4naive)
##                           tx_id gene_name         gene_id   gene_biotype type
## ENST00000443726 ENST00000443726    BCL11B ENSG00000127152 protein_coding  cds
## ENST00000357195 ENST00000357195    BCL11B ENSG00000127152 protein_coding  cds
## ENST00000583593 ENST00000583593    CCDC57 ENSG00000176155 protein_coding  cds
## ENSE00002456092 ENST00000463701     PRSS1 ENSG00000204983 protein_coding exon
## ENST00000546420 ENST00000546420    CCDC64 ENSG00000135127 protein_coding  cds
## ENST00000455990 ENST00000455990     HOOK1 ENSG00000134709 protein_coding  cds
##                            closest_region              query_region distance
## ENST00000443726   chr14-99737498-99737555   chr14-99721608-99741934        0
## ENST00000357195   chr14-99697682-99697894   chr14-99695477-99720910        0
## ENST00000583593   chr17-80085568-80085694   chr17-80084198-80086094        0
## ENSE00002456092  chr7-142460719-142460923  chr7-142501666-142511108    40742
## ENST00000546420 chr12-120427684-120428101 chr12-120426014-120428613        0
## ENST00000455990    chr1-60280790-60280852    chr1-60279767-60281364        0

head(closest_genes_cd14mono)
##                           tx_id     gene_name         gene_id   gene_biotype
## ENST00000432018 ENST00000432018          IL1B ENSG00000125538 protein_coding
## ENSE00001638912 ENST00000455005 RP5-1120P11.3 ENSG00000231881        lincRNA
## ENST00000445003 ENST00000445003 RP11-290F20.3 ENSG00000224397        lincRNA
## ENST00000568649 ENST00000568649         PPCDC ENSG00000138621 protein_coding
## ENST00000409245 ENST00000409245         TTC7A ENSG00000068724 protein_coding
## ENST00000484822 ENST00000484822          RXRA ENSG00000186350 protein_coding
##                 type           closest_region             query_region distance
## ENST00000432018  cds chr2-113593760-113593806 chr2-113581628-113594911        0
## ENSE00001638912 exon   chr6-44041650-44042535   chr6-44025105-44028184    13465
## ENST00000445003  gap  chr20-48884201-48894027  chr20-48889794-48893313        0
## ENST00000568649  cds  chr15-75335782-75335877  chr15-75334903-75336779        0
## ENST00000409245  cds   chr2-47300841-47301062   chr2-47297968-47301173        0
## ENST00000484822  gap chr9-137211331-137293477 chr9-137263243-137268534        0

本文由 mdnice 多平台发布

这篇关于单细胞分析(Signac): PBMC scATAC-seq 整合的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/998427

相关文章

性能分析之MySQL索引实战案例

文章目录 一、前言二、准备三、MySQL索引优化四、MySQL 索引知识回顾五、总结 一、前言 在上一讲性能工具之 JProfiler 简单登录案例分析实战中已经发现SQL没有建立索引问题,本文将一起从代码层去分析为什么没有建立索引? 开源ERP项目地址:https://gitee.com/jishenghua/JSH_ERP 二、准备 打开IDEA找到登录请求资源路径位置

SWAP作物生长模型安装教程、数据制备、敏感性分析、气候变化影响、R模型敏感性分析与贝叶斯优化、Fortran源代码分析、气候数据降尺度与变化影响分析

查看原文>>>全流程SWAP农业模型数据制备、敏感性分析及气候变化影响实践技术应用 SWAP模型是由荷兰瓦赫宁根大学开发的先进农作物模型,它综合考虑了土壤-水分-大气以及植被间的相互作用;是一种描述作物生长过程的一种机理性作物生长模型。它不但运用Richard方程,使其能够精确的模拟土壤中水分的运动,而且耦合了WOFOST作物模型使作物的生长描述更为科学。 本文让更多的科研人员和农业工作者

MOLE 2.5 分析分子通道和孔隙

软件介绍 生物大分子通道和孔隙在生物学中发挥着重要作用,例如在分子识别和酶底物特异性方面。 我们介绍了一种名为 MOLE 2.5 的高级软件工具,该工具旨在分析分子通道和孔隙。 与其他可用软件工具的基准测试表明,MOLE 2.5 相比更快、更强大、功能更丰富。作为一项新功能,MOLE 2.5 可以估算已识别通道的物理化学性质。 软件下载 https://pan.quark.cn/s/57

衡石分析平台使用手册-单机安装及启动

单机安装及启动​ 本文讲述如何在单机环境下进行 HENGSHI SENSE 安装的操作过程。 在安装前请确认网络环境,如果是隔离环境,无法连接互联网时,请先按照 离线环境安装依赖的指导进行依赖包的安装,然后按照本文的指导继续操作。如果网络环境可以连接互联网,请直接按照本文的指导进行安装。 准备工作​ 请参考安装环境文档准备安装环境。 配置用户与安装目录。 在操作前请检查您是否有 sud

线性因子模型 - 独立分量分析(ICA)篇

序言 线性因子模型是数据分析与机器学习中的一类重要模型,它们通过引入潜变量( latent variables \text{latent variables} latent variables)来更好地表征数据。其中,独立分量分析( ICA \text{ICA} ICA)作为线性因子模型的一种,以其独特的视角和广泛的应用领域而备受关注。 ICA \text{ICA} ICA旨在将观察到的复杂信号

【软考】希尔排序算法分析

目录 1. c代码2. 运行截图3. 运行解析 1. c代码 #include <stdio.h>#include <stdlib.h> void shellSort(int data[], int n){// 划分的数组,例如8个数则为[4, 2, 1]int *delta;int k;// i控制delta的轮次int i;// 临时变量,换值int temp;in

三相直流无刷电机(BLDC)控制算法实现:BLDC有感启动算法思路分析

一枚从事路径规划算法、运动控制算法、BLDC/FOC电机控制算法、工控、物联网工程师,爱吃土豆。如有需要技术交流或者需要方案帮助、需求:以下为联系方式—V 方案1:通过霍尔传感器IO中断触发换相 1.1 整体执行思路 霍尔传感器U、V、W三相通过IO+EXIT中断的方式进行霍尔传感器数据的读取。将IO口配置为上升沿+下降沿中断触发的方式。当霍尔传感器信号发生发生信号的变化就会触发中断在中断

kubelet组件的启动流程源码分析

概述 摘要: 本文将总结kubelet的作用以及原理,在有一定基础认识的前提下,通过阅读kubelet源码,对kubelet组件的启动流程进行分析。 正文 kubelet的作用 这里对kubelet的作用做一个简单总结。 节点管理 节点的注册 节点状态更新 容器管理(pod生命周期管理) 监听apiserver的容器事件 容器的创建、删除(CRI) 容器的网络的创建与删除

PostgreSQL核心功能特性与使用领域及场景分析

PostgreSQL有什么优点? 开源和免费 PostgreSQL是一个开源的数据库管理系统,可以免费使用和修改。这降低了企业的成本,并为开发者提供了一个活跃的社区和丰富的资源。 高度兼容 PostgreSQL支持多种操作系统(如Linux、Windows、macOS等)和编程语言(如C、C++、Java、Python、Ruby等),并提供了多种接口(如JDBC、ODBC、ADO.NET等

RabbitMQ使用及与spring boot整合

1.MQ   消息队列(Message Queue,简称MQ)——应用程序和应用程序之间的通信方法   应用:不同进程Process/线程Thread之间通信   比较流行的中间件:     ActiveMQ     RabbitMQ(非常重量级,更适合于企业级的开发)     Kafka(高吞吐量的分布式发布订阅消息系统)     RocketMQ   在高并发、可靠性、成熟度等