MACS2 参数详细学习
安装
pip install MACS2
MACS2功能：
macs2 callpeak 是macs2最主要的一个功能，能够利用bam文件寻找chip peak；
macs2 callpeak 使用：
# regular peak calling：
macs2 callpeak -t ChIP.bam -c Control.bam -f BAM -g hs -n test -B -q 0.01
# broad peak calling:
macs2 callpeak -t ChIP.bam -c Control.bam --broad -g hs --broad-cutoff 0.1
参数介绍
-T/–TREATMENT FILENAME：treat组
-C/–CONTROL：control 或 mock(非特异性抗体，如IgG)组
control:
input DNA，没有经过免疫共沉淀处理；
mock:
1）未使用抗体富集与蛋白结合的DNA片段
2）非特异性抗体，如IgG
-N/–NAME：为MACS2输出文件命名
‘NAME_peaks.xls’, ‘NAME_negative_peaks.xls’, ‘NAME_peaks.bed’ , ‘NAME_summits.bed’, ‘NAME_model.r’
–OUTDIR：MACS2结果文件保存路径
-F/–FORMAT FORMAT：MACS2读入文件格式(Tag文件),SAM，BAM,BED等，默认自动检测输入文件格式；
-G/–GSIZE：有效基因组大小(可比对基因组大小)；基因组中有大量重复序列测序测不到，实际上可比对的基因组大小只有原基因组90% 或 70%；人类默认值是– 2.7e9（UCSC human hg18 assembly）
| hs: | 2.7e9 |
|---|---|
| mm: | 1.87e9 |
| ce: | 9e7 |
| dm: | 1.2e8 |
-S/–TSIZE:测序读长；如果不设定，MACS 利用输入的前10个序列自动检测；
–BW：湿实验中，声波打断基因组的片段长度，用来建立模型；
-Q/–QVALUE：qvalue (minimum FDR)设定call significant regions的阈值；默认，0.01，对于 broad marks(组蛋白修饰的chipseq)，可以使用0.05；Q-values are calculated from p-values using Benjamini-Hochberg procedure.
-P/–PVALUE：设定p值时， qvalue不再起作用。
-M/–MFOLD：构建模型时，enrichment regions 选用标准（MFOLD range of high-confidence enrichment ratio against background to build model）;DEFAULT:5,50 means using all regions not too low (>5) and not too high (<50) to build paired-peaks model. MACS 无法找到超过100 regions 用来构建模型时，只有设定–fix-bimodal情况下，MACS 会调用参数–extsize。
–NOLAMBDA：不考虑peak 候选区域的偏差，使用背景λ作为 localλ。
–SLOCAL, –LLOCAL：设定两个水平检测peak 区域，从而计算最大λ作为local λ。默认，MACS 采用1000bp为small local region(–slocal)，10000bps为large local region(–llocal)计算开放染色体区域的偏差。区域设置的太小，尖峰会掩盖掉旁边显著性的峰。
–NOMODEL：MACS 不构建模型。
–EXTSIZE：设定–nomodel，MACS 会沿着 5’->3’方向延伸reads；如果转录因子结合区域长200bp，你也不想MACS建模，你就可以设定此参数为200.
–SHIFT：–shiftsize已经被 –extsize所替代；–nomodel设定之后，MACS 会用这个参数剪切reads5’，利用–extsize 延伸reads 3’端；如果设为负数，方向相反(3’->5’ );ChIP-Seq建议设置为0；当检测富集切割位点时，例如DNAseI-Seq datasets，此参数应该设为 -1 * half of EXTSIZE( EXTSIZE设为200，此参数为-100).
两个例子：
DNAse-Seq，想将平滑窗口设为200bps时，使用参数‘–nomodel –shift -100 –extsize 200’。
nucleosome-seq，使用核小体一半大小进行小波分析获得核小体中心的峰；当缠绕核小体DNA长度为147bps，可使用参数‘–nomodel –shift 37 –extsize 73’。
–KEEP-DUP：默认使用pvalue( 1e-5)基于二项式分布计算每个位置maximum tags；‘all’表示保留所有tags，如果设定了一个整数，那就是同一位置保留tags 的最大数。默认值为1，同一位置保留1 tag。
–BROAD：此参数会依据一个低的阈值(–broad-cutoff)将peaK附近富集区域归类到 broad region输出到BED12 格式文件。broad region最大长度是MACS计算的d的4倍。DEFAULT: False
–BROAD-CUTOFF：broad region阈值；pvalue 设定就是pvalue ,未设定就是qvalue；DEFAULT: 0.1。
–TO-LARGE：此参数设定后，线性放大小样本到大样本一样的深度；默认是缩小大样本到小样本深度。
注意：放大小样本可能产生更多的假阳性。
–DOWN-SAMPLE：设定此参数，使用随机抽样方法缩小大样本。随机抽样会使记过不稳定和不可重复。
-B/–BDG：保留the fragment pileup, control lambda, -log10pvalue 和 -log10qvalue scores到bedGraph 文件。
NAME+’_treat_pileup.bdg’：实验组数据
NAME+’_control_lambda.bdg’：对照组local lambda values
NAME+’_treat_pvalue.bdg’： Poisson pvalue scores (in -log10(pvalue) form)
NAME+’_treat_qvalue.bdg’ ： q-value scores from Benjamini–Hochberg–Yekutieli procedure
–CALL-SUMMITS：重新分析信号峰，从而获得主峰的临近峰；当要检测主峰周围的结合事件时，可使用此参数；结果中，同一主峰的临近峰有一样的范围 和不一样的分数，位置。
–VERBOSE：隐藏MACS运行过程信息，设置0；想了解各条染色体peak信息，设置为3或>3的数。
结果文件
1.NAME_peaks.xls
存放peak信息的文件
- 染色体名
- peak起始位置
- peak终止位置
- peak区域长度
- peak summit位置
- peak summit位置堆积信号, -log10(pvalue) for the peak summit (e.g. pvalue =1e-10, then this value should be 10)
- fold enrichment for this peak summit against random Poisson distribution with local lambda, -log10(qvalue) at peak summit
2.NAME_peaks.narrowPeak
BED6+4格式，包含peak位置信息，peak summit, pvalue and qvalue，可以使用UCSC genome browser查看。其中几列信息如下：
- 5th: integer score for display
- 7th: fold-change
- 8th: -log10pvalue
- 9th: -log10qvalue
- 10th: relative summit position to peak start
3.NAME_summits.bed
BED格式，包含peak summits(peak最高点)位置；如果想寻找结合位点的motifs ，建议使用此文件。
- 5th: -log10pvalue the same as NAME_peaks.bed
4.NAME_peaks.broadPeak
ED6+3格式，与narrowPeak类似，除了没有第10列peak summit的注释信息。
5.NAME_peaks.gappedPeak
BED12+3格式，存放broad region 和 narrow peaks，可以使用UCSC genome browser查看。
5th：10*-log10qvalue
7th：the start of the first narrow peak in the region
8th： the end of the first narrow peak in the region
9th： RGB color key, however, 0 for default color
10th：how many blocks including the starting 1bp and ending 1bp of broad regions.
11th：the length of each blocks
12th：the starts of each blocks.
13th: fold-change
14th: -log10pvalue
15th: -log10qvalue.
6.NAME_model.r
R程序，运行后生成基于输入数据产生的模型图片
$ Rscript NAME_model.r
7. .bdg files
bedGraph 文件，可以使用UCSC genome browser查看，或转格式为bigWig 文件； treat_pileup, and control_lambda
