which one is better for paired-ends peak calling? MACS2 OR MACS14?

yueli · 03-20-2019, 03:28 PM

Hello,

I'am working on paired-ends chip-seq analysis.

Which one is better for peak calling? macs2 OR macs14?

I ran both of them, macs2 has less peaks than macs14.

Thanks in advance for any help!

With Best!

Yue Li

KB* · 06-21-2019, 08:48 AM

Hi, @Yue Li

Did you find an answer for your question?

I am analysing ChIP-Seq data and I used MACS2. How may peaks did you get? I had to use q=0.005 (minimum FDR (q-value) cutoff for peak detection) to reduce number of called peaks. And I have got ~30k. I worry it is too many.
I will appreciate your feedback.

yueli · 06-22-2019, 10:31 AM

Hello, @KB

Flollwing is my command, qvalue cutoff = 5.00e-02. I got peaks around 400.

If using paired end reads use “--format BAMPE” to let MACS2 pileup the whole
fragments in general. If you want to focus on looking for where the 'cutting sites' are, then
“--nomodel --shift -100 --extsize 200” should work.

Since the DNA wrapped on a nucleosome is about 147bp, for single nucleosome detection

use “ --nomodel --shift -37 --extsize 73

INFO @ Thu, 21 Mar 2019 12:30:48: #1 fragment size = 2729.7

Please free to contract me if you have any questions.

Best,

Yue

administrator@ACB-HuangLab-Ubuntu:~/MACS2-2.1.2$ macs2 callpeak -t myc1_marked_duplicates.bam -c control_marked_duplicates.bam -f BAMPE -g hs -n myc1 -B --nomodel --shift -100 --extsize 200 --slocal 3000

/usr/lib/python2.7/dist-packages/pkg_resources.py:1031: UserWarning: /home/administrator/.python-eggs is writable by group/others and vulnerable to attack when used with get_resource_filename. Consider a more secure location (set with .set_extraction_path or the PYTHON_EGG_CACHE environment variable).

warnings.warn(msg, UserWarning)

INFO @ Fri, 22 Mar 2019 12:34:52:

# Command line: callpeak -t myc1_marked_duplicates.bam -c control_marked_duplicates.bam -f BAMPE -g hs -n myc1 -B --nomodel --shift -100 --extsize 200 --slocal 3000

# ARGUMENTS LIST:

# name = myc1

# format = BAMPE

# ChIP-seq file = ['myc1_marked_duplicates.bam']

# control file = ['control_marked_duplicates.bam']

# effective genome size = 2.70e+09

# band width = 300

# model fold = [5, 50]

# qvalue cutoff = 5.00e-02

# The maximum gap between significant sites is assigned as the read length/tag size.

# The minimum length of peaks is assigned as the predicted fragment length "d".

# Larger dataset will be scaled towards smaller dataset.

# Range for calculating regional lambda is: 3000 bps and 10000 bps

# Broad region calling is off

# Paired-End mode is on

INFO @ Fri, 22 Mar 2019 12:34:52: #1 read fragment files...

INFO @ Fri, 22 Mar 2019 12:34:52: #1 read treatment fragments...

INFO @ Fri, 22 Mar 2019 12:34:58: 1000000

INFO @ Fri, 22 Mar 2019 12:35:10: 2000000

INFO @ Fri, 22 Mar 2019 12:35:17: 3000000

INFO @ Fri, 22 Mar 2019 12:35:24: 4000000

INFO @ Fri, 22 Mar 2019 12:35:31: 5000000

INFO @ Fri, 22 Mar 2019 12:35:37: 6000000

INFO @ Fri, 22 Mar 2019 12:35:44: 7000000

INFO @ Fri, 22 Mar 2019 12:35:51: 8000000

INFO @ Fri, 22 Mar 2019 12:35:57: 9000000

INFO @ Fri, 22 Mar 2019 12:36:04: 10000000

INFO @ Fri, 22 Mar 2019 12:36:12: 11000000

INFO @ Fri, 22 Mar 2019 12:36:19: 12000000

INFO @ Fri, 22 Mar 2019 12:36:25: 13000000

INFO @ Fri, 22 Mar 2019 12:36:32: 14000000

INFO @ Fri, 22 Mar 2019 12:36:39: 15000000

INFO @ Fri, 22 Mar 2019 12:36:46: 16000000

INFO @ Fri, 22 Mar 2019 12:36:53: 17000000

INFO @ Fri, 22 Mar 2019 12:37:01: 18000000

INFO @ Fri, 22 Mar 2019 12:37:08: 19000000

INFO @ Fri, 22 Mar 2019 12:37:14: 20000000

INFO @ Fri, 22 Mar 2019 12:37:21: 21000000

INFO @ Fri, 22 Mar 2019 12:37:28: 22000000

INFO @ Fri, 22 Mar 2019 12:37:35: 23000000

INFO @ Fri, 22 Mar 2019 12:37:42: 24000000

INFO @ Fri, 22 Mar 2019 12:37:49: 25000000

INFO @ Fri, 22 Mar 2019 12:37:56: 26000000

INFO @ Fri, 22 Mar 2019 12:38:04: 27000000

INFO @ Fri, 22 Mar 2019 12:38:11: 28000000

INFO @ Fri, 22 Mar 2019 12:38:19: 29000000

INFO @ Fri, 22 Mar 2019 12:38:28: 30000000

INFO @ Fri, 22 Mar 2019 12:38:36: 31000000

INFO @ Fri, 22 Mar 2019 12:39:05: #1.2 read input fragments...

INFO @ Fri, 22 Mar 2019 12:39:14: 1000000

INFO @ Fri, 22 Mar 2019 12:39:21: 2000000

INFO @ Fri, 22 Mar 2019 12:39:27: 3000000

INFO @ Fri, 22 Mar 2019 12:39:34: 4000000

INFO @ Fri, 22 Mar 2019 12:39:41: 5000000

INFO @ Fri, 22 Mar 2019 12:39:48: 6000000

INFO @ Fri, 22 Mar 2019 12:39:55: 7000000

INFO @ Fri, 22 Mar 2019 12:40:02: 8000000

INFO @ Fri, 22 Mar 2019 12:40:09: 9000000

INFO @ Fri, 22 Mar 2019 12:40:16: 10000000

INFO @ Fri, 22 Mar 2019 12:40:23: 11000000

INFO @ Fri, 22 Mar 2019 12:40:30: 12000000

INFO @ Fri, 22 Mar 2019 12:40:37: 13000000

INFO @ Fri, 22 Mar 2019 12:40:46: 14000000

INFO @ Fri, 22 Mar 2019 12:41:02: #1 mean fragment size is determined as 2508.6 bp from treatment

INFO @ Fri, 22 Mar 2019 12:41:02: #1 note: mean fragment size in control is 2606.8 bp -- value ignored

INFO @ Fri, 22 Mar 2019 12:41:02: #1 fragment size = 2508.6

INFO @ Fri, 22 Mar 2019 12:41:02: #1 total fragments in treatment: 31689440

INFO @ Fri, 22 Mar 2019 12:41:02: #1 user defined the maximum fragments...

INFO @ Fri, 22 Mar 2019 12:41:02: #1 filter out redundant fragments by allowing at most 1 identical fragment(s)

INFO @ Fri, 22 Mar 2019 12:41:56: #1 fragments after filtering in treatment: 31667534

INFO @ Fri, 22 Mar 2019 12:41:56: #1 Redundant rate of treatment: 0.00

INFO @ Fri, 22 Mar 2019 12:41:56: #1 total fragments in control: 14684539

INFO @ Fri, 22 Mar 2019 12:41:56: #1 user defined the maximum fragments...

INFO @ Fri, 22 Mar 2019 12:41:56: #1 filter out redundant fragments by allowing at most 1 identical fragment(s)

INFO @ Fri, 22 Mar 2019 12:42:22: #1 fragments after filtering in control: 14676434

INFO @ Fri, 22 Mar 2019 12:42:22: #1 Redundant rate of control: 0.00

INFO @ Fri, 22 Mar 2019 12:42:22: #1 finished!

INFO @ Fri, 22 Mar 2019 12:42:22: #2 Build Peak Model...

INFO @ Fri, 22 Mar 2019 12:42:22: #2 Skipped...

INFO @ Fri, 22 Mar 2019 12:42:22: #3 Call peaks...

INFO @ Fri, 22 Mar 2019 12:42:22: #3 Pre-compute pvalue-qvalue table...

INFO @ Fri, 22 Mar 2019 12:48:42: #3 In the peak calling step, the following will be performed simultaneously:

INFO @ Fri, 22 Mar 2019 12:48:42: #3 Write bedGraph files for treatment pileup (after scaling if necessary)... myc1_treat_pileup.bdg

INFO @ Fri, 22 Mar 2019 12:48:42: #3 Write bedGraph files for control lambda (after scaling if necessary)... myc1_control_lambda.bdg

INFO @ Fri, 22 Mar 2019 12:48:42: #3 Pileup will be based on sequencing depth in control.

INFO @ Fri, 22 Mar 2019 12:48:42: #3 Call peaks for each chromosome...

INFO @ Fri, 22 Mar 2019 12:53:34: #4 Write output xls file... myc1_peaks.xls

INFO @ Fri, 22 Mar 2019 12:53:34: #4 Write peak in narrowPeak format file... myc1_peaks.narrowPeak

INFO @ Fri, 22 Mar 2019 12:53:35: #4 Write summits bed file... myc1_summits.bed

INFO @ Fri, 22 Mar 2019 12:53:35: Done!

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
MACS2 broad peak calling	deep639	Bioinformatics	0	03-24-2016 04:05 PM
cut 5' and 3' ends of paired-ends reads	sebl	Bioinformatics	4	12-08-2015 02:48 AM
peak identification with macs2 for small RNA data	skingan	Bioinformatics	0	05-11-2015 09:12 AM
macs2 -- peak calling error on non-model genome	grgmncr	Bioinformatics	0	07-29-2014 12:08 PM

03-20-2019, 03:28 PM	#1
yueli Member Location: china Join Date: May 2013 Posts: 71	which one is better for paired-ends peak calling? MACS2 OR MACS14? Hello, I'am working on paired-ends chip-seq analysis. Which one is better for peak calling? macs2 OR macs14? I ran both of them, macs2 has less peaks than macs14. Thanks in advance for any help! With Best! Yue Li

06-21-2019, 08:48 AM	#2
KB* Member Location: UK Join Date: May 2018 Posts: 20	Hi, @Yue Li Did you find an answer for your question? I am analysing ChIP-Seq data and I used MACS2. How may peaks did you get? I had to use q=0.005 (minimum FDR (q-value) cutoff for peak detection) to reduce number of called peaks. And I have got ~30k. I worry it is too many. I will appreciate your feedback.

06-22-2019, 10:31 AM	#3
yueli Member Location: china Join Date: May 2013 Posts: 71	Hello, @KB Flollwing is my command, qvalue cutoff = 5.00e-02. I got peaks around 400. If using paired end reads use “--format BAMPE” to let MACS2 pileup the whole fragments in general. If you want to focus on looking for where the 'cutting sites' are, then “--nomodel --shift -100 --extsize 200” should work. Since the DNA wrapped on a nucleosome is about 147bp, for single nucleosome detection use “ --nomodel --shift -37 --extsize 73 INFO @ Thu, 21 Mar 2019 12:30:48: #1 fragment size = 2729.7 Please free to contract me if you have any questions. Best, Yue administrator@ACB-HuangLab-Ubuntu:~/MACS2-2.1.2$ macs2 callpeak -t myc1_marked_duplicates.bam -c control_marked_duplicates.bam -f BAMPE -g hs -n myc1 -B --nomodel --shift -100 --extsize 200 --slocal 3000 /usr/lib/python2.7/dist-packages/pkg_resources.py:1031: UserWarning: /home/administrator/.python-eggs is writable by group/others and vulnerable to attack when used with get_resource_filename. Consider a more secure location (set with .set_extraction_path or the PYTHON_EGG_CACHE environment variable). warnings.warn(msg, UserWarning) INFO @ Fri, 22 Mar 2019 12:34:52: # Command line: callpeak -t myc1_marked_duplicates.bam -c control_marked_duplicates.bam -f BAMPE -g hs -n myc1 -B --nomodel --shift -100 --extsize 200 --slocal 3000 # ARGUMENTS LIST: # name = myc1 # format = BAMPE # ChIP-seq file = ['myc1_marked_duplicates.bam'] # control file = ['control_marked_duplicates.bam'] # effective genome size = 2.70e+09 # band width = 300 # model fold = [5, 50] # qvalue cutoff = 5.00e-02 # The maximum gap between significant sites is assigned as the read length/tag size. # The minimum length of peaks is assigned as the predicted fragment length "d". # Larger dataset will be scaled towards smaller dataset. # Range for calculating regional lambda is: 3000 bps and 10000 bps # Broad region calling is off # Paired-End mode is on INFO @ Fri, 22 Mar 2019 12:34:52: #1 read fragment files... INFO @ Fri, 22 Mar 2019 12:34:52: #1 read treatment fragments... INFO @ Fri, 22 Mar 2019 12:34:58: 1000000 INFO @ Fri, 22 Mar 2019 12:35:10: 2000000 INFO @ Fri, 22 Mar 2019 12:35:17: 3000000 INFO @ Fri, 22 Mar 2019 12:35:24: 4000000 INFO @ Fri, 22 Mar 2019 12:35:31: 5000000 INFO @ Fri, 22 Mar 2019 12:35:37: 6000000 INFO @ Fri, 22 Mar 2019 12:35:44: 7000000 INFO @ Fri, 22 Mar 2019 12:35:51: 8000000 INFO @ Fri, 22 Mar 2019 12:35:57: 9000000 INFO @ Fri, 22 Mar 2019 12:36:04: 10000000 INFO @ Fri, 22 Mar 2019 12:36:12: 11000000 INFO @ Fri, 22 Mar 2019 12:36:19: 12000000 INFO @ Fri, 22 Mar 2019 12:36:25: 13000000 INFO @ Fri, 22 Mar 2019 12:36:32: 14000000 INFO @ Fri, 22 Mar 2019 12:36:39: 15000000 INFO @ Fri, 22 Mar 2019 12:36:46: 16000000 INFO @ Fri, 22 Mar 2019 12:36:53: 17000000 INFO @ Fri, 22 Mar 2019 12:37:01: 18000000 INFO @ Fri, 22 Mar 2019 12:37:08: 19000000 INFO @ Fri, 22 Mar 2019 12:37:14: 20000000 INFO @ Fri, 22 Mar 2019 12:37:21: 21000000 INFO @ Fri, 22 Mar 2019 12:37:28: 22000000 INFO @ Fri, 22 Mar 2019 12:37:35: 23000000 INFO @ Fri, 22 Mar 2019 12:37:42: 24000000 INFO @ Fri, 22 Mar 2019 12:37:49: 25000000 INFO @ Fri, 22 Mar 2019 12:37:56: 26000000 INFO @ Fri, 22 Mar 2019 12:38:04: 27000000 INFO @ Fri, 22 Mar 2019 12:38:11: 28000000 INFO @ Fri, 22 Mar 2019 12:38:19: 29000000 INFO @ Fri, 22 Mar 2019 12:38:28: 30000000 INFO @ Fri, 22 Mar 2019 12:38:36: 31000000 INFO @ Fri, 22 Mar 2019 12:39:05: #1.2 read input fragments... INFO @ Fri, 22 Mar 2019 12:39:14: 1000000 INFO @ Fri, 22 Mar 2019 12:39:21: 2000000 INFO @ Fri, 22 Mar 2019 12:39:27: 3000000 INFO @ Fri, 22 Mar 2019 12:39:34: 4000000 INFO @ Fri, 22 Mar 2019 12:39:41: 5000000 INFO @ Fri, 22 Mar 2019 12:39:48: 6000000 INFO @ Fri, 22 Mar 2019 12:39:55: 7000000 INFO @ Fri, 22 Mar 2019 12:40:02: 8000000 INFO @ Fri, 22 Mar 2019 12:40:09: 9000000 INFO @ Fri, 22 Mar 2019 12:40:16: 10000000 INFO @ Fri, 22 Mar 2019 12:40:23: 11000000 INFO @ Fri, 22 Mar 2019 12:40:30: 12000000 INFO @ Fri, 22 Mar 2019 12:40:37: 13000000 INFO @ Fri, 22 Mar 2019 12:40:46: 14000000 INFO @ Fri, 22 Mar 2019 12:41:02: #1 mean fragment size is determined as 2508.6 bp from treatment INFO @ Fri, 22 Mar 2019 12:41:02: #1 note: mean fragment size in control is 2606.8 bp -- value ignored INFO @ Fri, 22 Mar 2019 12:41:02: #1 fragment size = 2508.6 INFO @ Fri, 22 Mar 2019 12:41:02: #1 total fragments in treatment: 31689440 INFO @ Fri, 22 Mar 2019 12:41:02: #1 user defined the maximum fragments... INFO @ Fri, 22 Mar 2019 12:41:02: #1 filter out redundant fragments by allowing at most 1 identical fragment(s) INFO @ Fri, 22 Mar 2019 12:41:56: #1 fragments after filtering in treatment: 31667534 INFO @ Fri, 22 Mar 2019 12:41:56: #1 Redundant rate of treatment: 0.00 INFO @ Fri, 22 Mar 2019 12:41:56: #1 total fragments in control: 14684539 INFO @ Fri, 22 Mar 2019 12:41:56: #1 user defined the maximum fragments... INFO @ Fri, 22 Mar 2019 12:41:56: #1 filter out redundant fragments by allowing at most 1 identical fragment(s) INFO @ Fri, 22 Mar 2019 12:42:22: #1 fragments after filtering in control: 14676434 INFO @ Fri, 22 Mar 2019 12:42:22: #1 Redundant rate of control: 0.00 INFO @ Fri, 22 Mar 2019 12:42:22: #1 finished! INFO @ Fri, 22 Mar 2019 12:42:22: #2 Build Peak Model... INFO @ Fri, 22 Mar 2019 12:42:22: #2 Skipped... INFO @ Fri, 22 Mar 2019 12:42:22: #3 Call peaks... INFO @ Fri, 22 Mar 2019 12:42:22: #3 Pre-compute pvalue-qvalue table... INFO @ Fri, 22 Mar 2019 12:48:42: #3 In the peak calling step, the following will be performed simultaneously: INFO @ Fri, 22 Mar 2019 12:48:42: #3 Write bedGraph files for treatment pileup (after scaling if necessary)... myc1_treat_pileup.bdg INFO @ Fri, 22 Mar 2019 12:48:42: #3 Write bedGraph files for control lambda (after scaling if necessary)... myc1_control_lambda.bdg INFO @ Fri, 22 Mar 2019 12:48:42: #3 Pileup will be based on sequencing depth in control. INFO @ Fri, 22 Mar 2019 12:48:42: #3 Call peaks for each chromosome... INFO @ Fri, 22 Mar 2019 12:53:34: #4 Write output xls file... myc1_peaks.xls INFO @ Fri, 22 Mar 2019 12:53:34: #4 Write peak in narrowPeak format file... myc1_peaks.narrowPeak INFO @ Fri, 22 Mar 2019 12:53:35: #4 Write summits bed file... myc1_summits.bed INFO @ Fri, 22 Mar 2019 12:53:35: Done!