MACS3—探索基因组调控的钥匙

生信菜鸟团

发布于 2024-05-30 16:26:00

89700

代码可运行

文章被收录于专栏：生信菜鸟团生信菜鸟团

运行总次数：0

代码可运行

工欲善其事必先利其器

1MACS3

MACS3 是由Liu Tao 主导开发的一款ChIP-Seq 数据分析工具，作为 MACS 系列软件的第三代版本，MACS3 继承并扩展了前两代软件的功能，成为 ChIP-Seq 数据分析中不可或缺的工具，因其具有以下特性，被广泛应用于基因组学、转录调控和表观遗传学研究等领域。

高精度峰值检测
- 模型驱动方法：MACS3 采用动态泊松分布模型，能够更精确地识别 ChIP-Seq 数据中的显著峰值。这种方法有效地区分了真实信号和背景噪音，显著提高了峰值检测的准确性。
- 自动片段大小估计：MACS3 能自动估计 DNA 片段的大小和位移大小，这对于优化峰值检测至关重要，有助于更准确地定位蛋白质-DNA 结合区域。
支持多种峰值类型
- 窄峰和宽峰检测：MACS3 支持识别不同类型的峰值，包括窄峰（如转录因子的结合位点）和宽峰（如组蛋白修饰区域），以满足不同实验的需求。
- 宽峰模式（Broadpeak）：通过 --broad 选项，用户可以检测宽峰，这在组蛋白修饰研究中尤其有用。
高效的数据处理
- 利用对照样本：MACS3 能够利用对照样本，通过比较实验样本和对照样本，减少背景噪音，提高峰值检测的特异性。
- 批量处理能力：MACS3 能处理大规模数据集，适用于高通量测序数据的分析需求。
灵活的参数设置
- 自定义参数：MACS3 提供了丰富的参数选项，用户可以根据具体需求调整参数设置。例如，--keep-dup 选项用于处理重复读取，--cutoff-analysis 选项用于分析不同的统计显著性水平。
- 多样的输出格式：MACS3 生成多种格式的输出文件，包括标准峰值文件（如 BED 格式）、统计报告和用于可视化的 bedGraph 或 bigWig 文件，方便后续的分析和展示。

2发表文章

题目：Model-based Analysis of ChIP-Seq (MACS) 期刊：Genome Biology 日期：2008年9月17日 作者&单位：Liu Tao & Dana-Farber Cancer Institute and Harvard School of Public Health DOI：https://doi.org/10.1186/gb-2008-9-9-r137 GitHub: https://github.com/macs3-project/MACS 文档：https://macs3-project.github.io/MACS/ 编程语言：核心主要使用 Python 语言编写，高性能计算的代码部分利用了 Cython 来提高性能 支持平台：仅支持 Linux 和 Mac OS系统

3如何安装

推荐使用conda 安装，方便快捷。

conda actiavte chipseq
conda install bioconda::macs3

4功能简述

MACS3功能

5最小化使用

callpeak

MACS3 中的主要函数，它用于处理各种格式的对齐文件，进而识别基因组中的显著富集区域，也就是所谓的“峰值”。这些峰值区域通常对应于蛋白质（如转录因子）与DNA的结合位置，对于理解基因表达调控机制极为重要。

#### 常规转录因子 ChIP-seq 的峰值识别
macs3 callpeak -t ChIP.bam -c Control.bam -f BAM -g hs -n test -B -q 0.01 --outdir outdir_name
#### 组蛋白标记 ChIP-seq 的宽峰识别
macs3 callpeak -t ChIP.bam -c Control.bam --broad -g hs --broad-cutoff 0.1
#### ATAC-seq 的峰值识别（成对末端模式）
macs3 callpeak -f BAMPE -t ATAC.bam -g hs -n test -B -q 0.01
#### ATAC-seq 的峰值识别（关注插入位点，使用单端模式）
macs3 callpeak -f BAM -t ATAC.bam -g hs -n test -B -q 0.01 --shift -50 --extension 100

-t #处理组文件，通常是ChIP-seq 实验数据的文件。可以同时指定多个文件，如 `-t A B C`，MACS3 会将这些文件中的数据合并处理。
-c #对照组文件，通常是背景或对照样本的数据文件。也可以同时指定多个文件，如 `-c A B C`，这些文件将被合并以估计 ChIP-seq 的背景噪声。
-f #指定输入文件的格式。支持格式：{AUTO,BAM,SAM,BED,ELAND,ELANDMULTi,ELANDEXPORT,BOWTIE,BAMPE,BEDPE} 。默认选项 `AUTO` 将由 MACS 自动决定文件格式。（注： `BAMPE` 和 `BEDPE` 需要明确设置）
-g #指定有效的基因组大小。可以是具体数字，如 `1.0e+9` 或 `1000000000`，或者使用简写，如 'hs' 代表人类（2,913,022,398），'mm' 代表小鼠（2,652,783,500），'ce' 代表秀丽隐杆线虫（100,286,401）和 'dm' 代表果蝇（142,573,017）。数字来源于 Deeptools ：https://deeptools.readthedocs.io/en/develop/content/feature/effectiveGenomeSize.html
-n #指定输出文件前缀
-q #设定富集区域检测的统计显著性阈值。阈值越低，找到的峰值质量越高，但可能会错过一些真实的低信号峰。
-B #设定是否分别保存extended fragment pileup 和 local lambda tracks[局部背景信号估计] 到bedGraph文件，这些文件对于在基因组浏览器中可视化峰值和背景信号非常有用。默认关闭。
--outdir #指定输出文件目录。默认是当前工作目录
--broad #启用检测宽峰模式
--shift -50 #定义一个任意的位移量（以碱基对为单位），用于在没有使用模型的情况下调整读取的切割端（5' 端）。默认值0
--extension 100 #定义一个任意的延伸大小（以碱基对为单位），当没有使用模型（`nomodel`）时，MACS 将使用此值作为片段大小将每个读取向 3' 端延伸，然后将它们堆叠起来。默认值200

6实例测试

/usr/bin/time -v macs3 callpeak -t SRR19436492_H3K27ac-Comb_sort.bam SRR19436493_H3K27ac-Palb_sort.bam SRR19436494_H3K27ac-OXA_sort.bam SRR19436495_H3K27ac-CTRL_sort.bam \
-c SRR19436496_Input-Comb_sort.bam SRR19436497_Input-Palb_sort.bam SRR19436498_Input-OXA_sort.bam SRR19436499_Input-CTRL_sort.bam \
-f BAM -g hs -n test -B -q 0.01 \
--outdir /home/data/t020559/chip_seq/GSE205035_PRJNA843319/f_peaks 1>macs3_test.log 2>&1 

##单样本
macs3 callpeak \
-t SRR19436492_H3K27ac-Comb_sort.bam \
-c SRR19436496_Input-Comb_sort.bam \
-f BAM -g hs -n H3K27ac-Comb -B -q 0.01 \
--outdir /home/data/t020559/chip_seq/macs3_out 1>test_H3K27ac.log 2>&1

结果文件

7结果文件

NAME_peaks.xls

至少包含9列信息，分别为：

染色体名称 (Chromosome Name)
Peak起始位置 (Start Position of Peak)
Peak结束位置 (End Position of Peak)
Peak区间长度 (Length of Peak Region)
绝对峰值位置 (Absolute Peak Summit Position)：Peak区域内信号最强的点，通常代表结合位点或功能区域的精确位置
峰顶的堆叠高度 (Pileup Height at Peak Summit)：反映该位置的测序深度或信号强度
峰顶的 p 值的负对数值 (-log10(pvalue) for the Peak Summit)：p 值越小，-log10(pvalue) 值越大，峰的显著性越高
fold_enrichment ：折叠富集值，表示峰顶信号强度相对于局部背景的富集程度。这是通过比较峰顶的信号与基于局部 lambda（局部平均信号水平）计算的预期信号强度得出的比值
峰顶的 q 值的负对数值(-log10(qvalue) at Peak Summit)：值越小，-log10(qvalue) 值越大，表示峰值的显著性越高且假阳性风险越低

NAME_peaks.xls

pileup 如何理解？

Pileup 是指在基因组的某一特定位置或区域，多个读取叠加在一起的现象，描述了在某个特定的基因组位置，有多少测序读取覆盖了这个位置。MACS3 会计算每个基因组位置上的 pileup 值，并使用统计模型来评估这些值是否显著高于背景水平（即对照组的测序数据或基于局部序列复杂度的预期水平）。

pileup

NAME_peaks.narrowPeak 和 NAME_peaks.broadPeak

NAME_peaks.narrowPeak 文件是BED6+4格式，可直接加载到UCSC基因组浏览器中查看，每列内容为：

chrom - 染色体名称（或者是contig, scaffold等）。
chromStart - 特征在染色体或基座上的起始位置。染色体的第一个碱基编号为0。
chromEnd - 特征在染色体或基座上的结束位置。在显示特征时不包括chromEnd所指的碱基。例如，染色体的前100个碱基定义为chromStart=0, chromEnd=100，覆盖编号为0-99的碱基。
name - 分配给区域的名称（最好是唯一的）。如果没有指定名称，则使用"."。
score - 指示该峰值在浏览器中显示的深浅（0-1000）。
strand - 用+/-表示链或方向。如果没有指定方向，则使用"."。
signalValue - 折叠富集值，表示峰顶信号强度相对于局部背景的富集程度。（同NAME_peaks.xls第8列）
pValue - 统计显著性的度量，p 值的负对数值。如果没有指定p值，则使用-1。
qValue - q值的负对数值false discovery rate，-log10）度量的统计显著性。如果没有指定q值，则使用-1。
peak - 表示峰值最高点相对于起始位置的偏移，用于标定峰顶的精确位置。如果一个峰的起始位置（chromStart）是 100，并且峰顶位于 150，那么峰顶相对于起始位置的偏移量就是 50。那么在narrowPeak文件中这个偏移量会被记录为 50

NAME_peaks.narrowPeak

NAME_peaks.broadPeak 是 BED6+3格式，与 NAME_peaks.narrowPeak相比少了第10列的值

NAME_peaks.broadPeak

NAME_summits.bed

NAME_summits.bed

NAME_treat_pileup.bdg 和 NAME_control_lambda.bdg

pileup.bdg 和 lambda.bdg 文件都采用 bedGraph 格式。通常包含四列

Chromosome（染色体）
Start Position（起始位置），bedGraph 格式中，这个位置是从零开始计数的。
End Position（结束位置），实际的区域是从起始位置到结束位置前一个碱基对。
Data Value（数据值）

对于 NAME_treat_pileup.bdg 文件，这一列通常表示处理样本中每个基因组位置的 pileup 信号强度，这个信号通常会根据 --scale-to 选项进行归一化处理。
对于 NAME_control_lambda.bdg 文件，这一列表示从对照样本（或在缺少对照样本的情况下从处理样本）估计出的每个基因组位置的局部偏差（local biases）。