随着测序技术的进步，染色质免疫沉淀技术被广泛用于研究全基因组蛋白-DNA互作。macs 基于一种新的模型可以很好的识别转录因子结合位点。macs 可以直接应用于ChIP-Seq 数据，也可以将ChIP-Seq数据与control结合起来提高特异性。
安装
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
-
control:
-
-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读入文件格式,"ELAND", "BED", "ELANDMULTI", "ELANDEXPORT", "ELANDMULTIPET" (for pair-end tags), "SAM", "BAM", "BOWTIE", "BAMPE" or "BEDPE"；默认自动检测输入文件格式，因此可以使用不同格式的文件。 -
-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)
- fold enrichment for this peak summit against random Poisson distribution with local lambda
- -log10(qvalue) at peak summit
- peak name
2.NAME_peaks.narrowPeak
BED6+4格式，包含peak位置信息，peak summit, pvalue and qvalue，可以使用UCSC genome browser查看。其中几列信息如下：
- 1th: 染色体名
- 2th: peak 起始位置
- 3th: peak 终止位置
- 4th: peak name
- 5th: integer score for display,
int(-10*log10(pvalue)) - 7th: fold-change
- 8th: -log10(pvalue)
- 9th: -log10qvalue
- 10th: 峰位与peak起点的距离
3.NAME_summits.bed
BED格式，包含peak summits(peak最高点)位置；如果想寻找结合位点的motifs ，建议使用此文件。
- 5th: -log10pvalue
4.NAME_peaks.broadPeak
ED6+3格式，与narrowPeak类似，除了没有第10列peak summit的注释信息。
5.NAME_peaks.gappedPeak
BED12+3格式，存放broad region 和 narrow peaks，可以使用UCSC genome browser查看。
6.NAME_model.r
R程序，运行后生成基于输入数据产生的模型图片
$ Rscript NAME_model.r
7. .bdg files
bedGraph 文件，可以导入UCSC genome browser查看，或转格式为bigWig 文件；
- treat_pileup ：实验组bedGraph 文件
- control_lambda ：对照组bedGraph 文件
