software

用 Signac 合并两份 ATAC 数据,为什么同一个峰不能直接把矩阵拼起来?

作者:发布于 2026-11-14约 7 分钟阅读
两个样本原来是五和六,重数呢。教学算例|比较区域从各自原峰改为统一的并集或交集。原峰,甲为五,乙为六,样本各自使用的边界不同。并集,甲为六,乙为六,同一区域内片段总数相等。交集,甲为二,乙为二,还需检查并集左右两侧信号是否重新分布。

甲样本的峰标成一百到二百,乙样本在附近也有一个峰,只是边界向右移了五十个碱基。分析人员想把两行合成一行,没测到的位置填零,再比较两组开放程度。这样做藏着一个难以从图上发现的错误。两个数字对应的不是同一段基因组,即使给两行贴上相同名称,也不会让过去的计数自动覆盖新的边界。Signac 可以从各自的片段文件重新定量,前提是坐标体系、峰集合和细胞的身份都对得上。

先画两个半重叠的峰,再决定比较什么

为避免边界写法混淆,下面的教学算例使用左闭右开的坐标约定,所有片段都放在区间内部。甲样本原先的峰是第一百到第二百,乙样本的峰是第一百五十到第二百五十,两峰有一段重叠。若科学问题是这一大片候选区域是否整体更开放,可以选择它们的并集作为共同区域。若问题是重叠核心是否改变,就应使用交集。选哪一个直接决定分母所对应的空间范围,应在计算前说明。

设甲来源的片段在左侧独有区、共同区、右侧独有区分别有三、二、一个。乙来源分别有零、二、四个。下表特意假定每条片段完全落在一个小区间,便于读者手算。甲的原峰计五个,乙的原峰计六个。按共同并集重新定量,两者其实都是六个,按交集重新定量,则都是两个。

在哪一段比较 甲样本 乙样本 它们是否指同一段坐标
各自原来的峰 3 + 2 = 5 2 + 4 = 6 左右边界不同
两峰并集 3 + 2 + 1 = 6 0 + 2 + 4 = 6 同一并集
两峰交集 2 2 同一交集
两个样本原来是五和六,重数呢。教学算例|比较区域从各自原峰改为统一的并集或交集。原峰,甲为五,乙为六,样本各自使用的边界不同。并集,甲为六,乙为六,同一区域内片段总数相等。交集,甲为二,乙为二,还需检查并集左右两侧信号是否重新分布。
教学算例|比较区域从各自原峰改为统一的并集或交集。

这个例子还保留了一个真实的生物学可能。尽管两个样本的并集总量相同,甲的左侧更强,乙的右侧更强。若只汇报合并后的总数,就会错过局部位置变化。对于特别宽、包含多个局部高点的区域,可以查看分供体的覆盖曲线,按有依据的窄区间另作比较。多个互相搭边的峰还可能被串成一个很长的并集,合并后先查宽度分布、异常宽区域及其附近的高点。Signac 作者的多样本合并教程

从片段文件重数,而不是从两个成品矩阵猜

用共同的区域集在每个样本的片段文件上计数,可获得定义一致的行。只拥有原峰矩阵时,甲峰的五个计数究竟有几个落在交集、几个落在左侧,从一个总数无法反推。物理上未观察到的区段也不能一律解释成零。尤其当两个数据来自不同参考版本、不同染色体命名,或有一方缺某条染色体的片段时,强行补零会制造生物学差异。先确认二者使用同一参考坐标,检查主要染色体和目标区间能否真实查询。

片段文件应采用可按区间查询的格式,并有同目录对应的索引。模板里的 CreateFragmentObject 还可核对指定细胞是否在文件中存在。如果为合并对象给细胞加上样本前缀,构造映射时,向量的名字是对象里的细胞名,向量的值才是文件第四列存放的原始条形码。把两者倒置,即使所有字符串都眼熟也可能查不到对应细胞。两个文库恰好有相同原始条形码时,更需要各自带样本前缀保留身份。片段对象及条形码映射

下面示意一个样本的计数。假设已有按同一组装版本建立的基因组区间对象,并且原条形码向量来自该捕获通过质控的细胞。对第二个样本用自己的片段文件和细胞标识执行同样过程,比较返回的行名和顺序,检查一致后再拼列。单独记下空区间和缺少的染色体,不要仅为绕过维度不同的报错而补齐数字。FeatureMatrix 的输入和返回值

library(Signac)
map <- setNames(raw_barcodes, paste0("A_", raw_barcodes))
fr <- CreateFragmentObject(path = "A/fragments.tsv.gz", cells = map)
mat_a <- FeatureMatrix(fragments = fr, features = common_ranges,
                       cells = names(map))
stopifnot(!anyDuplicated(rownames(mat_a)), !anyDuplicated(colnames(mat_a)))

最好挑两枚已知细胞,用带前缀和原始的两套名字分别人工追溯。若这个小测试不通过,先修正映射再做全量计数。若行列名都有、总矩阵却在某一条染色体全为零,先查片段文件和区域集合中的染色体命名是否一方带前缀、一方不带。

实跑经验|共同峰总数相等,右侧局部计数却是四比一

我们安装 Signac 1.17.1,把上面的教学例子写成两份真实的片段格式文件,各含六条片段,全部按坐标排序,再做分块压缩并建立 Tabix 索引。两个样本都故意使用相同的原始条形码 bar01,分别映射到带样本前缀的对象列名。用同一个 FeatureMatrix 调用分别对两份文件定量。按各自最初的峰观察,甲为五、乙为六。改用共同并集,两者都是六,共同核心两者都是二,而只看最右的独有区,甲为一、乙为四。因此“整段总开放相同”和“局部片段分布不同”可以从同一份原始数据同时得到。

另做一项容易忽略的对照,保持甲样本所有片段的坐标和条形码不变,只把第一条片段第五列的支持读段数由一改为四。相同的五个区间重新定量,计数依次仍为五、三、六、二、一。对本例所用入口,不能把片段表第五列直接当作按需累加的区间权重。若将甲的六条片段、乙的六条片段和含额外测序支持的甲副本同时制表,应明确后者是同一数据的输入对照,不能当成第三个独立样本。片段文件各列的官方定义

用一个小区域排除坐标偏移

拿一个覆盖清楚、能在基因组浏览器定位的短区域做核查。在原始片段文件中查询区间,确认返回的染色体名称、坐标边界、样本名和条形码符合预期。文件采用的起点约定与 R 的区间对象可能不同,读入时使用合适的格式转换,避免多加一次一或漏掉转换。坐标偏移一个碱基在宽峰中不一定明显,却可能改变窄区域边界附近的计数。把一个恰好贴近边缘的片段单独做对照比只看总峰数更容易识别这个问题。

同样检查总计数的分母。如果一个样本每核可用片段普遍较多,在相同区域看到更多片段并不直接等价于处理激活了这个调控元件。可以先按样本内相同细胞类型聚合,再保留每位供体的有效片段水平,按实验设计比较。两个条件各只有一个供体时,即使得到几千枚细胞,处理和个体来源也难以区分。共同峰定义解决的是量尺一致的问题,深度与生物学重复需要另外处理。

另做一次区域选择的敏感性检查。固定片段和细胞集合,在一个有实际争议的位点分别用共同并集、较窄的重叠核心和左右两侧区间计数。如果结论只能在某侧出现,回看局部覆盖、比对可信度和该处是否有两个相邻高点,再决定说“区域总开放增加”还是“信号在局部重新分布”。这两个描述对后续验证位置的选择有实际影响。

向合作者交付时附上共同峰集、参考版本、样本条形码映射和定量所用片段文件版本。日后要换一套增强子边界时,保留原片段就还有办法重新作答。需要在团队中共享多样本片段并反复比较不同区域定义时,可了解云生信的生物信息分析服务器。

把知识用于当前任务

下一步可以直接开始分析

按研究目标查找在线工具、课程、计算资源或完整分析服务。

按目标找方案仍有疑问,联系客服 ↗