diff --git a/CHANGELOG.md b/CHANGELOG.md index ff5abf0c..778a6d75 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -15,12 +15,15 @@ and this project adheres to [Semantic Versioning](http://semver.org/spec/v2.0.0. * [nf-core/atacseq#76](https://github.com/nf-core/atacseq/issues/76) - featureCounts coupled to DESeq2 * [nf-core/atacseq#79](https://github.com/nf-core/atacseq/issues/79) - Parallelize DESeq2 * [nf-core/atacseq#97](https://github.com/nf-core/atacseq/issues/97) - PBC1, PBC2 from pipeline? +* [nf-core/atacseq#107](https://github.com/nf-core/atacseq/issues/107) - Add options to change MACS2 parameters * [nf-core/atacseq#109](https://github.com/nf-core/atacseq/issues/109) - Specify custom gtf but gene bed is not generated from that gtf? * Regenerated screenshots and added collapsible sections for output files in `docs/output.md` * Update template to tools `1.9` * Replace `set` with `tuple` and `file()` with `path()` in all processes * Capitalise process names * Parameters: + * `--macs_fdr` to provide FDR threshold for MACS2 peak calling + * `--macs_pvalue` to provide p-value threshold for MACS2 peak calling * `--skip_peak_qc` to skip MACS2 peak QC plot generation * `--skip_peak_annotation` to skip annotation of MACS2 and consensus peaks with HOMER * `--skip_consensus_peaks` to skip consensus peak generation diff --git a/docs/usage.md b/docs/usage.md index b632e2cd..0c956cdd 100644 --- a/docs/usage.md +++ b/docs/usage.md @@ -35,6 +35,8 @@ * [Peaks](#peaks) * [`--narrow_peak`](#--narrow_peak) * [`--broad_cutoff`](#--broad_cutoff) + * [`--macs_fdr`](#--macs_fdr) + * [`--macs_pvalue`](#--macs_pvalue) * [`--min_reps_consensus`](#--min_reps_consensus) * [`--save_macs_pileup`](#--save_macs_pileup) * [`--skip_peak_qc`](#--skip_peak_qc) @@ -389,6 +391,14 @@ MACS2 is run by default with the [`--broad`](https://github.com/taoliu/MACS#--br Specifies broad cut-off value for MACS2. Only used when `--narrow_peak` isnt specified (Default: `0.1`). +### `--macs_fdr` + +Minimum FDR (q-value) cutoff for peak detection, `--macs_fdr` and `--macs_pvalue` are mutually exclusive (Default: false). + +### `--macs_pvalue` + +p-value cutoff for peak detection, `--macs_fdr` and `--macs_pvalue` are mutually exclusive (Default: false). If `--macs_pvalue` cutoff is set, q-value will not be calculated and reported as -1 in the final .xls file. + ### `--min_reps_consensus` Number of biological replicates required from a given condition for a peak to contribute to a consensus peak . If you are confident you have good reproducibility amongst your replicates then you can increase the value of this parameter to create a "reproducible" set of consensus of peaks. For example, a value of 2 will mean peaks that have been called in at least 2 replicates will contribute to the consensus set of peaks, and as such peaks that are unique to a given replicate will be discarded. diff --git a/main.nf b/main.nf index c5ebbf41..0514bf50 100755 --- a/main.nf +++ b/main.nf @@ -56,6 +56,8 @@ def helpMessage() { Peaks --narrow_peak [bool] Run MACS2 in narrowPeak mode (Default: false) --broad_cutoff [float] Specifies broad cutoff value for MACS2. Only used when --narrow_peak isnt specified (Default: 0.1) + --macs_fdr [float] Minimum FDR (q-value) cutoff for peak detection, --macs_fdr and --macs_pvalue are mutually exclusive (Default: false) + --macs_pvalue [float] p-value cutoff for peak detection, --macs_fdr and --macs_pvalue are mutually exclusive (Default: false) --min_reps_consensus [int] Number of biological replicates required from a given condition for a peak to contribute to a consensus peak (Default: 1) --save_macs_pileup [bool] Instruct MACS2 to create bedGraph files normalised to signal per million reads (Default: false) --skip_peak_qc [bool] Skip MACS2 peak QC plot generation (Default: false) @@ -230,6 +232,8 @@ summary['MACS2 Genome Size'] = params.macs_gsize ?: 'Not supplied' summary['Min Consensus Reps'] = params.min_reps_consensus if (params.macs_gsize) summary['MACS2 Narrow Peaks'] = params.narrow_peak ? 'Yes' : 'No' if (!params.narrow_peak) summary['MACS2 Broad Cutoff'] = params.broad_cutoff +if (params.macs_fdr) summary['MACS2 FDR'] = params.macs_fdr +if (params.macs_pvalue) summary['MACS2 P-value'] = params.macs_pvalue if (params.skip_trimming) { summary['Trimming Step'] = 'Skipped' } else { @@ -1121,6 +1125,8 @@ process MACS2 { broad = params.narrow_peak ? '' : "--broad --broad-cutoff ${params.broad_cutoff}" format = params.single_end ? 'BAM' : 'BAMPE' pileup = params.save_macs_pileup ? '-B --SPMR' : '' + fdr = params.macs_fdr ? "--qvalue ${params.macs_fdr}" : '' + pvalue = params.macs_pvalue ? "--pvalue ${params.macs_pvalue}" : '' """ macs2 callpeak \\ -t ${ipbam[0]} \\ @@ -1130,6 +1136,8 @@ process MACS2 { -g $params.macs_gsize \\ -n $ip \\ $pileup \\ + $fdr \\ + $pvalue \\ --keep-dup all cat ${ip}_peaks.${PEAK_TYPE} | wc -l | awk -v OFS='\t' '{ print "${ip}", \$1 }' | cat $peak_count_header - > ${ip}_peaks.count_mqc.tsv diff --git a/nextflow.config b/nextflow.config index 23aa20c5..287363ed 100644 --- a/nextflow.config +++ b/nextflow.config @@ -36,6 +36,8 @@ params { // Options: Peaks narrow_peak = false broad_cutoff = 0.1 + macs_fdr = false + macs_pvalue = false min_reps_consensus = 1 save_macs_pileup = false skip_peak_qc = false