第3篇：用MACS2软件call peaks - 云+社区 - 腾讯云 (tencent.com)

MACS2软件call m6A甲基化peaks - 简书 (jianshu.com)

macs2 - 简书 (jianshu.com)

https://mp.weixin.qq.com/s/61VRtPkuRaHSxWeMxvispA  

MeRIP-seq 数据分析之 callpeak 及 peak 可视化

#######下载（新集群py37）

conda activate py37

pip install MACS2

#####################

conda install -c bioconda macs2

conda install -c bioconda/label/cf201901 macs2

了解相关参数：

1.输入文件参数：

-t: 实验组，IP的数据文件

c: 对照组

f：指定输入文件的格式，默认是自动检测输入数据是什么格式，支持bam,sam,bed等

g:有效基因组大小，由于基因组序列的重复性，基因组实际可以mapping的大小小于原始的基因组。这个参数要根据实际物种计算基因组的有效大小。软件里也给出了几个默认的-g 值：hs -- 2.7e9表示人类的基因组有效大小(UCSC human hg18 assembly).

hs: 2.7e9；mm: 1.87e9；ce: 9e7；dm: 1.2e8

2.输出文件参数：

--outdir:输出结果的存储路径 -n:输出文件名的前缀

-B/--bdg:输出bedgraph格式的文件，输出文件以NAME+'_treat_pileup.bdg' for treatment data, NAME+'_control_lambda.bdg' for local lambda values from control显示。

3.peak calling 参数：

-q/--qvalue和-p/--pvalue      

q value默认值是0.05，与pvalue不能同时使用。 

q值与峰宽有一定的联系。理想情况下，如果放宽阈值，您将简单地获得更多的峰值，但是使用MACS2放松阈值也会导致更宽的峰值。

--broad                peak有narrow peak和broad peak, 设置时可以call broad peak 的结果文件。

--broad-cutoff     和pvalue、以及qvalue相似

--nolambda:     不要考虑在峰值候选区域的局部偏差/λ

4.Shift 模型参数：

--nomodel     这个参数和extsize、shift是配套使用的，有这个参数才可以设置extsize和shift。

--extsize          当设置了nomodel时，MACS会用--extsize这个参数从5'->3'方向扩展reads修复fragments。比如说你的转录因子结合范围200bp，就设置这个参数是200。

--shift               当设置了--nomodel，MACS用这个参数从5' 端移动剪切，然后用--extsize延伸，如果--shift是负值表示从3'端方向移动。建议ChIP-seq数据集这个值保持默认值为0，对于检测富集剪切位点如DNAsel数据集设置为EXTSIZE的一半。

用法

#####################

awk '{print "macs2 callpeak -c 5.STAR_allignment/"$5"_Input/"$5"_InputAligned.sortedByCoord.out.bam -t 5.STAR_allignment/"$5"_IP/"$5"_IPAligned.sortedByCoord.out.bam -q 0.05 -f BAM -g hs -n 7.Macs2CallPeak/"$5"/"$5"_IP_Input"}' sampleinfo.input | xargs -iCMD -P30 bash -c CMD

################

常规峰值调用的示例：macs2 callpeak -t ChIP.bam -c Control.bam -f BAM -g hs -n test -B -q 0.01

宽峰调用示例：macs2 callpeak -t ChIP.bam -c Control.bam --broad -g hs --broad-cutoff 0.1

#######################

macs2 callpeak -t 116si_p.ip_siWTAP.bam \

           -c 116si_p.input_siWTAP.bam \

           -g hs -q 0.05 -f BAM

#有control数据

fori in *.bamdomacs2 callpeak -f BAM -c control.bam -t$i-n$i-g hs --outdir ../macs2/ --bdg -q 0.05done

基本选项

-t/--treatment FILENAME，IP 样本

这是MACS唯一的REQUIRED参数。文件可以是任何文件--format选项指定的受支持格式。检查 - 格式化详情。如果您有多个对齐文件，则可以指定

他们是-t A B C。 MACS会将所有这些文件汇集在一起。

-c/--control，Input 样本

控件或模拟数据文件。请与-t对应的control

-n/--name 【结果】文件的前缀

实验的名称字符串。 MACS将使用此字符串NAME创建输出文件，如NAME_peaks.xls，NAME_negative_peaks.xls，NAME_peaks.bed，NAME_summits.bed，NAME_model.r等等。所以请避免这些文件名和您的文件名之间的任何冲突现有文件。

--outdir  【输出文件夹】

MACS2会将所有输出文件保存到speficied文件夹中选项。

-f/--format FORMAT

标签文件的格式，可以是“ELAND”，“BED”，“ELANDMULTI”，“ELANDEXPORT”，“ELANDMULTIPET”（用于对端标签），“SAM”，“BAM”，“BOWTIE”，“BAMPE”或“BEDPE”。默认为“AUTO”，这将允许MACS自动决定格式。当您使用“AUTO”时也会使用

结合不同格式的文件。请注意，MACS无法检测到“BAMPE”或“BEDPE”格式带有“AUTO”，你必须隐含指定“BAMPE”和“BEDPE”的格式。

如今，最常见的格式是BED或BAM / SAM。

-g / --gsize,有效基因组大小，可以是数字，也可以是：hs，mm，ce，dm。

这是可映射的基因组大小或有效的基因组大小。定义为可以测序的基因组大小。因为在chromsomes上的重复特征，实际可映射的基因组大小将小于原始大小，约为基因组的90％或70％尺寸。对于UCSC人类hg18，建议使用默认的hs  -  2.7e9。以下是有效基因组的所有预编译参数。

size：

hs：2.7e9

mm：1.87e9

ce：9e7

dm：1.2e8

-q/--qvalue

调用重要区域的qvalue（最小FDR）默认是0.05。对于广泛的标记，您可以尝试0.05作为cutoff。 Q值是使用Benjamini-Hochberg程序从p值计算。

-p/--pvalue

如果指定了-p，则MACS2将使用pvalue代替qvalue。

--broad

设置该参数，MACS将尝试合并广泛的区域， The maximum length of broad region length is 4 times of d from MACS. DEFAULT:False。

-B/--bdg,保留the fragment pileup,control lambda到 bedGraph 文件。

bedGraph files will be stored in current directory named

--fe-cutoff：fold-enrichment。

--keep-dup：保留重复的 tag。

--outdir：文件输出路径。

-n：输出文件名称前缀。

-B：保留the fragment pileup,control lambda到 bedGraph 文件。

--SPMR：normalized pileup tracks across many datasets，可以在样本间进行比较 peaks。

--nomodel：不使用 macs2 建立 shifting model。

--extsize：与--nomodel 一起使用，从 5'->3'延伸 fragment 到指定长度。

--fe-cutoff：fold-enrichment。