如何进行MACS2 peak calling的实战

发布时间：2021-12-31 11:37:08 作者：柒染
来源：亿速云阅读：320

# 如何进行MACS2 peak calling的实战

## 1. 什么是MACS2及其应用场景

MACS2（Model-based Analysis of ChIP-Seq 2）是用于识别染色质免疫沉淀测序（ChIP-seq）数据中富集区域（peak calling）的主流工具。它通过统计建模解决以下关键问题：

- **噪声过滤**：区分真实信号与背景噪声
- **峰识别**：准确定位转录因子结合位点或组蛋白修饰区域
- **双峰检测**：特别适用于识别转录因子结合位点周围的核小体缺失区域

典型应用场景包括：
- 转录因子结合位点分析
- 组蛋白修饰研究（H3K27ac, H3K4me3等）
- ATAC-seq数据分析

## 2. 安装与数据准备

### 2.1 安装MACS2
推荐通过conda安装：
```bash
conda install -c bioconda macs2

2.2 输入数据要求

需要两种输入文件： 1. 处理后的BAM文件： - 已完成比对（推荐使用bowtie2/hisat2） - 已去除重复 reads（建议使用picard） - 示例文件：treatment.bam

对照样本（可选但强烈建议）：
- Input DNA或IgG对照
- 示例文件：control.bam

3. 基础分析流程

3.1 基础peak calling命令

macs2 callpeak -t treatment.bam -c control.bam \
    -f BAM -g hs -n experiment_name \
    --outdir peaks/ --bdg --qvalue 0.05

参数解析： - -t：处理样本 - -c：对照样本 - -f：输入格式（BAM/BED等） - -g：基因组大小（hs=human, mm=mouse） - --bdg：输出bedGraph格式信号文件 - --qvalue：显著性阈值（默认0.05）

3.2 输出文件说明

关键输出文件： - _peaks.narrowPeak：主结果文件（BED6+4格式） - _peaks.xls：包含详细统计信息的表格 - _summits.bed：每个peak的精确最高点 - _treat_pileup.bdg：处理样本信号轨迹 - _control_lambda.bdg：背景信号估计

4. 进阶参数调优

4.1 双端测序数据处理

macs2 callpeak -t treatment.bam -f BAMPE \
    --keep-dup all -g hs --broad

特殊参数： - BAMPE：自动计算片段长度 - --keep-dup：重复reads处理策略 - --broad：用于组蛋白修饰的宽峰

4.2 灵敏度调节

macs2 callpeak -t treatment.bam -c control.bam \
    --broad-cutoff 0.1 --min-length 200 \
    --max-gap 100

4.3 批次处理多个样本

for sample in TF1 TF2 TF3; do
    macs2 callpeak -t ${sample}.bam -c input.bam \
        -f BAM -g hs -n ${sample}
done

5. 结果解读与可视化

5.1 peak文件格式解析

.narrowPeak文件示例：

chr1 1000 2000 peak_1 500 . 8.32 45.67 32.14 150

各列含义： 1. 染色体 2. start 3. end 4. peak名称 5. 得分 6. 链信息 7. fold-change 8. -log10(pvalue) 9. -log10(qvalue) 10. summit位置

5.2 使用IGV可视化

加载BAM文件和bedGraph
添加.narrowPeak注释轨道
使用_summits.bed精确定位

5.3 质控指标解读

关键质控参数： - FRiP（Fraction of Reads in Peaks）：>1%为合格 - NSC（Normalized Strand Cross-correlation）：>1.05 - RSC（Relative Strand Cross-correlation）：>0.8

6. 常见问题解决方案

6.1 内存不足问题

macs2 callpeak -t large.bam --tempdir /tmp \
    --buffer-size 1000000

6.2 处理混合样本

macs2 callpeak -t rep1.bam rep2.bam -c input.bam \
    --nomodel --extsize 200

6.3 非模式生物分析

macs2 callpeak -t non_model.bam -g 1.2e8 \
    --keep-dup auto

7. 下游分析建议

差异peak分析：

使用DiffBind包（R语言）

library(DiffBind)
dba <- dba(sampleSheet="samples.csv")
dba <- dba.count(dba)
dba <- dba.contrast(dba)
dba <- dba.analyze(dba)

motif分析：

使用MEME Suite或HOMER

findMotifsGenome.pl peaks.bed genome.fa output_dir

功能注释：

使用ChIPseeker（R包）

library(ChIPseeker)
peaks <- readPeakFile("peaks.narrowPeak")
annotatePeak(peaks, TxDb=txdb)

8. 完整示例代码

# 质控
fastqc *.fastq
multiqc .

# 比对
bowtie2 -x hg19 -U input.fastq -S input.sam
samtools sort -o input.bam input.sam

# peak calling
macs2 callpeak -t ChIP.bam -c input.bam \
    -f BAM -g hs -n MY_TF \
    --outdir peaks/ --bdg --qvalue 0.01

# 可视化准备
bedSort MY_TF_treat_pileup.bdg MY_TF_sorted.bdg
bedGraphToBigWig MY_TF_sorted.bdg hg19.chrom.sizes MY_TF.bw

9. 参考文献

Zhang Y, et al. (2008) “Model-based Analysis of ChIP-Seq (MACS)”. Genome Biology
最新文档：https://github.com/macs3-project/MACS

”`