main.nf 18.7 KB
Newer Older
Beibei Chen's avatar
Beibei Chen committed
1
#!/usr/bin/env nextflow
Venkat Malladi's avatar
Venkat Malladi committed
2

Venkat Malladi's avatar
Venkat Malladi committed
3
4
5
6
7
8
9
10
11
12
/*

BICF ChIP-seq Analysis Workflow
#### Homepage / Documentation
https://git.biohpc.swmed.edu/BICF/Astrocyte/chipseq_analysis/
Licensed under MIT (https://git.biohpc.swmed.edu/BICF/Astrocyte/chipseq_analysis/LICENSE.md)


*/

Venkat Malladi's avatar
Venkat Malladi committed
13
14
15
16
17
// Path to an input file, or a pattern for multiple inputs
// Note - $baseDir is the location of this workflow file main.nf

// Define Input variables
params.reads = "$baseDir/../test_data/*.fastq.gz"
18
params.pairedEnd = false
Venkat Malladi's avatar
Venkat Malladi committed
19
params.designFile = "$baseDir/../test_data/design_ENCSR238SGC_SE.txt"
20
params.genome = 'GRCm38'
Venkat Malladi's avatar
Venkat Malladi committed
21
params.cutoffRatio = 1.2
Venkat Malladi's avatar
Venkat Malladi committed
22
params.outDir= "$baseDir/output"
23
params.extendReadsLen = 100
24
params.topPeakCount = 600
25
params.astrocyte = false
26
27
params.skipDiff = false
params.skipMotif = false
28
params.skipPlotProfile = false
29
params.references = "$baseDir/../docs/references.md"
Venkat Malladi's avatar
Venkat Malladi committed
30
params.multiqc =  "$baseDir/conf/multiqc_config.yaml"
31

32
// Assign variables if astrocyte
33
34
if (params.astrocyte) {
  print("Running under astrocyte")
35
  referenceLocation = "/project/shared/bicf_workflow_ref"
36
37
38
39
40
41
  if (params.genome == 'GRCh37') {
    params.bwaIndex = "$referenceLocation/human/$params.genome"
    params.chromSizes = "$referenceLocation/human/$params.genome/genomefile.txt"
    params.fasta = "$referenceLocation/human/$params.genome/genome.fa"
    params.gtf = "$referenceLocation/human/$params.genome/gencode.v19.chr_patch_hapl_scaff.annotation.gtf"
    params.geneNames = "$referenceLocation/human/$params.genome/genenames.txt"
42
    params.genomeSize = 'hs'
Venkat Malladi's avatar
Venkat Malladi committed
43
  } else if (params.genome == 'GRCm38') {
44
45
46
47
48
    params.bwaIndex = "$referenceLocation/mouse/$params.genome"
    params.chromSizes = "$referenceLocation/mouse/$params.genome/genomefile.txt"
    params.fasta = "$referenceLocation/mouse/$params.genome/genome.fa"
    params.gtf = "$referenceLocation/mouse/$params.genome/gencode.vM20.annotation.gtf"
    params.geneNames = "$referenceLocation/mouse/$params.genome/genenames.txt"
49
    params.genomeSize = 'mm'
50
51
52
53
54
55
56
  } else if (params.genome == 'GRCh38') {
    params.bwaIndex = "$referenceLocation/human/$params.genome"
    params.chromSizes = "$referenceLocation/human/$params.genome/genomefile.txt"
    params.fasta = "$referenceLocation/human/$params.genome/genome.fa"
    params.gtf = "$referenceLocation/human/$params.genome/gencode.v25.chr_patch_hapl_scaff.annotation.gtf"
    params.geneNames = "$referenceLocation/human/$params.genome/genenames.txt"
    params.genomeSize = 'hs'
57
  }
58
59
60
61
62
} else {
    params.bwaIndex = params.genome ? params.genomes[ params.genome ].bwa ?: false : false
    params.genomeSize = params.genome ? params.genomes[ params.genome ].genomesize ?: false : false
    params.chromSizes = params.genome ? params.genomes[ params.genome ].chromsizes ?: false : false
    params.fasta = params.genome ? params.genomes[ params.genome ].fasta ?: false : false
63
    params.gtf = params.genome ? params.genomes[ params.genome ].gtf ?: false : false
64
    params.geneNames = params.genome ? params.genomes[ params.genome ].geneNames ?: false : false
65
66
67
}


68

Venkat Malladi's avatar
Venkat Malladi committed
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
// Check inputs
if( params.bwaIndex ){
  bwaIndex = Channel
    .fromPath(params.bwaIndex)
    .ifEmpty { exit 1, "BWA index not found: ${params.bwaIndex}" }
} else {
  exit 1, "No reference genome specified."
}

// Define List of Files
readsList = Channel
  .fromPath( params.reads )
  .flatten()
  .map { file -> [ file.getFileName().toString(), file.toString() ].join("\t")}
  .collectFile( name: 'fileList.tsv', newLine: true )

Venkat Malladi's avatar
Venkat Malladi committed
85
// Define regular variables
Venkat Malladi's avatar
Venkat Malladi committed
86
pairedEnd = params.pairedEnd
Venkat Malladi's avatar
Venkat Malladi committed
87
designFile = Channel.fromPath(params.designFile)
Venkat Malladi's avatar
Venkat Malladi committed
88
genomeSize = params.genomeSize
89
genome = params.genome
Venkat Malladi's avatar
Venkat Malladi committed
90
chromSizes = params.chromSizes
Venkat Malladi's avatar
Venkat Malladi committed
91
fasta = params.fasta
Venkat Malladi's avatar
Venkat Malladi committed
92
cutoffRatio = params.cutoffRatio
Venkat Malladi's avatar
Venkat Malladi committed
93
outDir = params.outDir
94
extendReadsLen = params.extendReadsLen
Venkat Malladi's avatar
Venkat Malladi committed
95
topPeakCount = params.topPeakCount
96
97
skipDiff = params.skipDiff
skipMotif = params.skipMotif
98
skipPlotProfile = params.skipPlotProfile
Venkat Malladi's avatar
Venkat Malladi committed
99
references = params.references
Venkat Malladi's avatar
Venkat Malladi committed
100
multiqc = params.multiqc
Venkat Malladi's avatar
Venkat Malladi committed
101
102
gtfFile = params.gtf
geneNames = params.geneNames
Venkat Malladi's avatar
Venkat Malladi committed
103

104
// Check design file for errors
Venkat Malladi's avatar
Venkat Malladi committed
105
106
process checkDesignFile {

Venkat Malladi's avatar
Venkat Malladi committed
107
  publishDir "$outDir/design", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
108
109
110

  input:

Venkat Malladi's avatar
Venkat Malladi committed
111
112
  file designFile
  file readsList
Venkat Malladi's avatar
Venkat Malladi committed
113
114
115
116
117
118
119
120
121

  output:

  file("design.tsv") into designFilePaths

  script:

  if (pairedEnd) {
    """
122
    module load python/3.6.1-2-anaconda
Venkat Malladi's avatar
Venkat Malladi committed
123
124
125
126
127
    python3 $baseDir/scripts/check_design.py -d $designFile -f $readsList -p
    """
  }
  else {
    """
128
    module load python/3.6.1-2-anaconda
Venkat Malladi's avatar
Venkat Malladi committed
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
    python $baseDir/scripts/check_design.py -d $designFile -f $readsList
    """
  }

}

// Define channel for raw reads
if (pairedEnd) {
  rawReads = designFilePaths
    .splitCsv(sep: '\t', header: true)
    .map { row -> [ row.sample_id, [row.fastq_read1, row.fastq_read2], row.experiment_id, row.biosample, row.factor, row.treatment, row.replicate, row.control_id ] }
} else {
rawReads = designFilePaths
  .splitCsv(sep: '\t', header: true)
  .map { row -> [ row.sample_id, [row.fastq_read1], row.experiment_id, row.biosample, row.factor, row.treatment, row.replicate, row.control_id ] }
}

// Trim raw reads using trimgalore
process trimReads {

  tag "$sampleId-$replicate"
Venkat Malladi's avatar
Venkat Malladi committed
150
  publishDir "$outDir/${task.process}/${sampleId}", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
151
152
153
154
155
156
157
158

  input:

  set sampleId, reads, experimentId, biosample, factor, treatment, replicate, controlId from rawReads

  output:

  set sampleId, file('*.fq.gz'), experimentId, biosample, factor, treatment, replicate, controlId into trimmedReads
159
160
  file('*trimming_report.txt') into trimgaloreResults
  file('version_*.txt') into trimReadsVersions
Venkat Malladi's avatar
Venkat Malladi committed
161
162
163
164
165

  script:

  if (pairedEnd) {
    """
166
167
    module load python/3.6.1-2-anaconda
    module load trimgalore/0.4.1
168
    python3 $baseDir/scripts/trim_reads.py -f ${reads[0]} ${reads[1]} -s $sampleId -p
Venkat Malladi's avatar
Venkat Malladi committed
169
170
171
172
    """
  }
  else {
    """
173
174
    module load python/3.6.1-2-anaconda
    module load trimgalore/0.4.1
175
    python3 $baseDir/scripts/trim_reads.py -f ${reads[0]} -s $sampleId
Venkat Malladi's avatar
Venkat Malladi committed
176
177
178
179
180
181
182
183
    """
  }

}

// Align trimmed reads using bwa
process alignReads {

184
  queue '128GB,256GB,256GBv1'
Venkat Malladi's avatar
Venkat Malladi committed
185
  tag "$sampleId-$replicate"
Venkat Malladi's avatar
Venkat Malladi committed
186
  publishDir "$outDir/${task.process}/${sampleId}", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
187
188
189
190
191
192
193
194
195

  input:

  set sampleId, reads, experimentId, biosample, factor, treatment, replicate, controlId from trimmedReads
  file index from bwaIndex.first()

  output:

  set sampleId, file('*.bam'), experimentId, biosample, factor, treatment, replicate, controlId into mappedReads
Venkat Malladi's avatar
Venkat Malladi committed
196
  file '*.flagstat.qc' into mappedReadsStats
197
  file('version_*.txt') into alignReadsVersions
Venkat Malladi's avatar
Venkat Malladi committed
198
199
200
201
202

  script:

  if (pairedEnd) {
    """
203
204
205
    module load python/3.6.1-2-anaconda
    module load bwa/intel/0.7.12
    module load samtools/1.6
206
    python3 $baseDir/scripts/map_reads.py -f ${reads[0]} ${reads[1]} -r ${index}/genome.fa -s $sampleId -p
Venkat Malladi's avatar
Venkat Malladi committed
207
208
209
210
    """
  }
  else {
    """
211
212
213
    module load python/3.6.1-2-anaconda
    module load bwa/intel/0.7.12
    module load samtools/1.6
214
    python3 $baseDir/scripts/map_reads.py -f $reads -r ${index}/genome.fa -s $sampleId
Venkat Malladi's avatar
Venkat Malladi committed
215
216
217
218
219
220
221
222
    """
  }

}

// Dedup reads using sambamba
process filterReads {

223
  queue '128GB,256GB,256GBv1'
Venkat Malladi's avatar
Venkat Malladi committed
224
  tag "$sampleId-$replicate"
Venkat Malladi's avatar
Venkat Malladi committed
225
  publishDir "$outDir/${task.process}/${sampleId}", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
226
227
228
229
230
231
232
233
234

  input:

  set sampleId, mapped, experimentId, biosample, factor, treatment, replicate, controlId from mappedReads

  output:

  set sampleId, file('*.bam'), file('*.bai'), experimentId, biosample, factor, treatment, replicate, controlId into dedupReads
  set sampleId, file('*.bam'), experimentId, biosample, factor, treatment, replicate, controlId into convertReads
Venkat Malladi's avatar
Venkat Malladi committed
235
236
237
  file '*.flagstat.qc' into dedupReadsStats
  file '*.pbc.qc' into dedupReadsComplexity
  file '*.dedup.qc' into dupReads
238
  file('version_*.txt') into filterReadsVersions
Venkat Malladi's avatar
Venkat Malladi committed
239
240
241
242
243

  script:

  if (pairedEnd) {
    """
244
245
246
247
    module load python/3.6.1-2-anaconda
    module load samtools/1.6
    module load sambamba/0.6.6
    module load bedtools/2.26.0
Venkat Malladi's avatar
Venkat Malladi committed
248
249
250
251
252
    python3 $baseDir/scripts/map_qc.py -b $mapped -p
    """
  }
  else {
    """
253
254
255
256
    module load python/3.6.1-2-anaconda
    module load samtools/1.6
    module load sambamba/0.6.6
    module load bedtools/2.26.0
Venkat Malladi's avatar
Venkat Malladi committed
257
258
259
260
261
262
263
264
265
266
    python3 $baseDir/scripts/map_qc.py -b $mapped
    """
  }

}

// Define channel collecting dedup reads into new design file
dedupReads
.map{ sampleId, bam, bai, experimentId, biosample, factor, treatment, replicate, controlId ->
"$sampleId\t$bam\t$bai\t$experimentId\t$biosample\t$factor\t$treatment\t$replicate\t$controlId\n"}
Venkat Malladi's avatar
Venkat Malladi committed
267
.collectFile(name:'design_dedup.tsv', seed:"sample_id\tbam_reads\tbam_index\texperiment_id\tbiosample\tfactor\ttreatment\treplicate\tcontrol_id\n", storeDir:"$outDir/design")
Venkat Malladi's avatar
Venkat Malladi committed
268
269
270
271
272
.into { dedupDesign; preDiffDesign }

// Quality Metrics using deeptools
process experimentQC {

273
  queue '128GB,256GB,256GBv1'
Venkat Malladi's avatar
Venkat Malladi committed
274
  publishDir "$outDir/${task.process}", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
275
276
277
278
279
280
281

  input:

  file dedupDesign

  output:

Venkat Malladi's avatar
Venkat Malladi committed
282
  file '*.{pdf,npz}' into experimentQCStats
283
  file('version_*.txt') into experimentQCVersions
Venkat Malladi's avatar
Venkat Malladi committed
284
285
286
287

  script:

  """
288
289
  module load python/3.6.1-2-anaconda
  module load deeptools/2.5.0.1
290
  python3 $baseDir/scripts/experiment_qc.py -d $dedupDesign -e $extendReadsLen
Venkat Malladi's avatar
Venkat Malladi committed
291
292
  """

293
}
Venkat Malladi's avatar
Venkat Malladi committed
294

Venkat Malladi's avatar
Venkat Malladi committed
295
296
297
// Convert reads to bam
process convertReads {

298
  queue '128GB,256GB,256GBv1'
Venkat Malladi's avatar
Venkat Malladi committed
299
  tag "$sampleId-$replicate"
Venkat Malladi's avatar
Venkat Malladi committed
300
  publishDir "$outDir/${task.process}/${sampleId}", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
301
302
303
304
305
306
307
308

  input:

  set sampleId, deduped, experimentId, biosample, factor, treatment, replicate, controlId from convertReads

  output:

  set sampleId, file('*.tagAlign.gz'), file('*.bed{pe,se}.gz'), experimentId, biosample, factor, treatment, replicate, controlId into tagReads
309
  file('version_*.txt') into convertReadsVersions
Venkat Malladi's avatar
Venkat Malladi committed
310
311
312
313
314

  script:

  if (pairedEnd) {
    """
315
316
317
    module load python/3.6.1-2-anaconda
    module load samtools/1.6
    module load bedtools/2.26.0
Venkat Malladi's avatar
Venkat Malladi committed
318
319
320
321
322
    python3 $baseDir/scripts/convert_reads.py -b $deduped -p
    """
  }
  else {
    """
323
324
325
    module load python/3.6.1-2-anaconda
    module load samtools/1.6
    module load bedtools/2.26.0
Venkat Malladi's avatar
Venkat Malladi committed
326
327
328
329
330
331
332
333
334
335
    python3 $baseDir/scripts/convert_reads.py -b $deduped
    """
  }

}

// Calculate Cross-correlation using phantompeaktools
process crossReads {

  tag "$sampleId-$replicate"
Venkat Malladi's avatar
Venkat Malladi committed
336
  publishDir "$outDir/${task.process}/${sampleId}", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
337
338
339
340
341
342
343
344

  input:

  set sampleId, seTagAlign, tagAlign, experimentId, biosample, factor, treatment, replicate, controlId from tagReads

  output:

  set sampleId, seTagAlign, tagAlign, file('*.cc.qc'), experimentId, biosample, factor, treatment, replicate, controlId into xcorReads
345
346
  set file('*.cc.qc'), file('*.cc.plot.pdf') into crossReadsStats
  file('version_*.txt') into crossReadsVersions
Venkat Malladi's avatar
Venkat Malladi committed
347
348
349
350
351

  script:

  if (pairedEnd) {
    """
352
353
    module load python/3.6.1-2-anaconda
    module load phantompeakqualtools/1.2
Venkat Malladi's avatar
Venkat Malladi committed
354
355
356
357
358
    python3 $baseDir/scripts/xcor.py -t $seTagAlign -p
    """
  }
  else {
    """
359
360
    module load python/3.6.1-2-anaconda
    module load phantompeakqualtools/1.2
Venkat Malladi's avatar
Venkat Malladi committed
361
362
363
364
365
366
367
368
369
370
    python3 $baseDir/scripts/xcor.py -t $seTagAlign
    """
  }

}

// Define channel collecting tagAlign and xcor into design file
xcorDesign = xcorReads
              .map{ sampleId, seTagAlign, tagAlign, xcor, experimentId, biosample, factor, treatment, replicate, controlId ->
              "$sampleId\t$seTagAlign\t$tagAlign\t$xcor\t$experimentId\t$biosample\t$factor\t$treatment\t$replicate\t$controlId\n"}
Venkat Malladi's avatar
Venkat Malladi committed
371
              .collectFile(name:'design_xcor.tsv', seed:"sample_id\tse_tag_align\ttag_align\txcor\texperiment_id\tbiosample\tfactor\ttreatment\treplicate\tcontrol_id\n", storeDir:"$outDir/design")
Venkat Malladi's avatar
Venkat Malladi committed
372
373
374
375

// Make Experiment design files to be read in for downstream analysis
process defineExpDesignFiles {

Venkat Malladi's avatar
Venkat Malladi committed
376
  publishDir "$outDir/design", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
377
378
379
380
381
382
383

  input:

  file xcorDesign

  output:

Venkat Malladi's avatar
Venkat Malladi committed
384
  file '*.tsv' into experimentObjs mode flatten
Venkat Malladi's avatar
Venkat Malladi committed
385
386
387
388

  script:

  """
389
  module load python/3.6.1-2-anaconda
Venkat Malladi's avatar
Venkat Malladi committed
390
391
392
393
394
395
396
397
398
399
400
  python3 $baseDir/scripts/experiment_design.py -d $xcorDesign
  """

}


// Make Experiment design files to be read in for downstream analysis
process poolAndPsuedoReads {


  tag "${experimentObjs.baseName}"
Venkat Malladi's avatar
Venkat Malladi committed
401
  publishDir "$outDir/design", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
402
403
404
405
406
407
408

  input:

  file experimentObjs

  output:

Venkat Malladi's avatar
Venkat Malladi committed
409
  file '*.tsv' into experimentPoolObjs
Venkat Malladi's avatar
Venkat Malladi committed
410
411
412
413
414

  script:

  if (pairedEnd) {
    """
415
    module load python/3.6.1-2-anaconda
Venkat Malladi's avatar
Venkat Malladi committed
416
417
418
419
420
    python3 $baseDir/scripts/pool_and_psuedoreplicate.py -d $experimentObjs -c $cutoffRatio -p
    """
  }
  else {
    """
421
    module load python/3.6.1-2-anaconda
Venkat Malladi's avatar
Venkat Malladi committed
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
    python3 $baseDir/scripts/pool_and_psuedoreplicate.py -d $experimentObjs -c $cutoffRatio
    """
  }

}

// Collect list of experiment design files into a single channel
experimentRows = experimentPoolObjs
                .splitCsv(sep:'\t', header:true)
                .map { row -> [ row.sample_id, row.tag_align, row.xcor, row.experiment_id, row.biosample, row.factor, row.treatment, row.replicate, row.control_id, row.control_tag_align] }

// Call Peaks using MACS
process callPeaksMACS {

  tag "$sampleId-$replicate"
Venkat Malladi's avatar
Venkat Malladi committed
437
  publishDir "$outDir/${task.process}/${experimentId}/${replicate}", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
438
439
440
441
442
443
444

  input:
  set sampleId, tagAlign, xcor, experimentId, biosample, factor, treatment, replicate, controlId, controlTagAlign from experimentRows

  output:

  set sampleId, file('*.narrowPeak'), file('*.fc_signal.bw'), file('*.pvalue_signal.bw'), experimentId, biosample, factor, treatment, replicate, controlId into experimentPeaks
Venkat Malladi's avatar
Venkat Malladi committed
445
  file '*.xls' into callPeaksMACSsummit
446
  file('version_*.txt') into callPeaksMACSVersions
447
  file("*.fc_signal.bw") into bigwigs
Venkat Malladi's avatar
Venkat Malladi committed
448
449
450
451
452

  script:

  if (pairedEnd) {
    """
453
454
455
456
457
    module load python/3.6.1-2-anaconda
    module load macs/2.1.0-20151222
    module load UCSC_userApps/v317
    module load bedtools/2.26.0
    module load phantompeakqualtools/1.2
Venkat Malladi's avatar
Venkat Malladi committed
458
459
460
461
462
    python3 $baseDir/scripts/call_peaks_macs.py -t $tagAlign -x $xcor -c $controlTagAlign -s $sampleId -g $genomeSize -z $chromSizes -p
    """
  }
  else {
    """
463
464
465
466
467
    module load python/3.6.1-2-anaconda
    module load macs/2.1.0-20151222
    module load UCSC_userApps/v317
    module load bedtools/2.26.0
    module load phantompeakqualtools/1.2
Venkat Malladi's avatar
Venkat Malladi committed
468
469
470
471
472
473
474
475
476
477
    python3 $baseDir/scripts/call_peaks_macs.py -t $tagAlign -x $xcor -c $controlTagAlign -s $sampleId -g $genomeSize -z $chromSizes
    """
  }

}

// Define channel collecting peaks into design file
peaksDesign = experimentPeaks
              .map{ sampleId, peak, fcSignal, pvalueSignal, experimentId, biosample, factor, treatment, replicate, controlId ->
              "$sampleId\t$peak\t$fcSignal\t$pvalueSignal\t$experimentId\t$biosample\t$factor\t$treatment\t$replicate\t$controlId\n"}
Venkat Malladi's avatar
Venkat Malladi committed
478
              .collectFile(name:'design_peak.tsv', seed:"sample_id\tpeaks\tfc_signal\tpvalue_signal\texperiment_id\tbiosample\tfactor\ttreatment\treplicate\tcontrol_id\n", storeDir:"$outDir/design")
Venkat Malladi's avatar
Venkat Malladi committed
479

480
481
//plotProfile
process plotProfile {
482
  publishDir "$outDir/experimentQC", mode: 'copy'
483
484
485

  input:

Venkat Malladi's avatar
Venkat Malladi committed
486
  file bigWigList from bigwigs.collect()
487

Jeremy Mathews's avatar
Jeremy Mathews committed
488
489
  output:

Jeremy Mathews's avatar
Jeremy Mathews committed
490
  file '*.{png,gz}' into plotProfile
Jeremy Mathews's avatar
Jeremy Mathews committed
491

492
493
494
495
496
497
498
  when:

  !skipPlotProfile

  script:
  """
  module load deeptools/2.5.0.1
Venkat Malladi's avatar
Venkat Malladi committed
499
  bash $baseDir/scripts/plot_profile.sh -g $gtfFile
500
501
502
  """
}

Venkat Malladi's avatar
Venkat Malladi committed
503
504
505
// Calculate Consensus Peaks
process consensusPeaks {

Venkat Malladi's avatar
Venkat Malladi committed
506
  publishDir "$outDir/${task.process}", mode: 'copy'
507
  publishDir "$outDir/design", mode: 'copy',  pattern: '*.{csv|tsv}'
Venkat Malladi's avatar
Venkat Malladi committed
508
509
510
511
512
513
514
515

  input:

  file peaksDesign
  file preDiffDesign

  output:

Venkat Malladi's avatar
Venkat Malladi committed
516
517
518
519
520
  file '*.replicated.*' into consensusPeaks
  file '*.rejected.*' into rejectedPeaks
  file 'design_diffPeaks.csv'  into designDiffPeaks
  file 'design_annotatePeaks.tsv'  into designAnnotatePeaks, designMotifSearch
  file 'unique_experiments.csv' into uniqueExperiments
521
  file('version_*.txt') into consensusPeaksVersions
Venkat Malladi's avatar
Venkat Malladi committed
522
523
524
525

  script:

  """
526
527
  module load python/3.6.1-2-anaconda
  module load bedtools/2.26.0
Venkat Malladi's avatar
Venkat Malladi committed
528
529
530
531
532
533
534
535
  python3 $baseDir/scripts/overlap_peaks.py -d $peaksDesign -f $preDiffDesign
  """

}

// Annotate Peaks
process peakAnnotation {

Venkat Malladi's avatar
Venkat Malladi committed
536
  publishDir "$outDir/${task.process}", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
537
538
539
540
541
542
543

  input:

  file designAnnotatePeaks

  output:

Venkat Malladi's avatar
Venkat Malladi committed
544
  file "*chipseeker*" into peakAnnotation
545
  file('version_*.txt') into peakAnnotationVersions
Venkat Malladi's avatar
Venkat Malladi committed
546
547
548
549

  script:

  """
550
  module load R/3.3.2-gccmkl
Venkat Malladi's avatar
Venkat Malladi committed
551
  Rscript $baseDir/scripts/annotate_peaks.R $designAnnotatePeaks $gtfFile $geneNames
Venkat Malladi's avatar
Venkat Malladi committed
552
553
554
555
  """

}

556
// Motif Search Peaks
557
process motifSearch {
Venkat Malladi's avatar
Venkat Malladi committed
558

Venkat Malladi's avatar
Venkat Malladi committed
559
  publishDir "$outDir/${task.process}", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
560
561
562

  input:

563
  file designMotifSearch
Venkat Malladi's avatar
Venkat Malladi committed
564
565
566

  output:

Venkat Malladi's avatar
Venkat Malladi committed
567
568
  file "*memechip" into motifSearch
  file "*narrowPeak" into filteredPeaks
569
  file('version_*.txt') into motifSearchVersions
Venkat Malladi's avatar
Venkat Malladi committed
570

571
  when:
572

573
  !skipMotif
Venkat Malladi's avatar
Venkat Malladi committed
574
575
576

  script:

577
  """
578
  module load python/3.6.1-2-anaconda
579
580
  module load meme/4.11.1-gcc-openmpi
  module load bedtools/2.26.0
581
582
  python3 $baseDir/scripts/motif_search.py -d $designMotifSearch -g $fasta -p $topPeakCount
  """
Venkat Malladi's avatar
Venkat Malladi committed
583
584
}

585
// Define channel to find number of unique experiments
586
uniqueExperimentsList = uniqueExperiments
587
                      .splitCsv(sep: '\t', header: true)
588
589
590

// Calculate Differential Binding Activity
process diffPeaks {
Venkat Malladi's avatar
Venkat Malladi committed
591

Venkat Malladi's avatar
Venkat Malladi committed
592
  publishDir "$outDir/${task.process}", mode: 'copy'
Venkat Malladi's avatar
Venkat Malladi committed
593
594
595

  input:

596
  file designDiffPeaks
597
  val noUniqueExperiments from uniqueExperimentsList.count()
Venkat Malladi's avatar
Venkat Malladi committed
598
599
600

  output:

Venkat Malladi's avatar
Venkat Malladi committed
601
602
603
604
  file '*_diffbind.bed' into diffPeaks
  file '*_diffbind.csv' into diffPeaksCounts
  file '*.pdf' into diffPeaksStats
  file 'normcount_peaksets.txt' into normCountPeaks
605
  file('version_*.txt') into diffPeaksVersions
606
607

  when:
608

609
  noUniqueExperiments > 1 && !skipDiff
Venkat Malladi's avatar
Venkat Malladi committed
610

611
  script:
612

Venkat Malladi's avatar
Venkat Malladi committed
613
  """
614
  module load python/3.6.1-2-anaconda
615
  module load R/3.3.2-gccmkl
616
  Rscript $baseDir/scripts/diff_peaks.R $designDiffPeaks
Venkat Malladi's avatar
Venkat Malladi committed
617
618
  """
}
619

Venkat Malladi's avatar
Venkat Malladi committed
620
621
// Generate Multiqc Report, gerernate Software Versions and references
process multiqcReport {
Venkat Malladi's avatar
Venkat Malladi committed
622
623
  publishDir "$outDir/${task.process}", mode: 'copy'

624
625
  input:

Venkat Malladi's avatar
Venkat Malladi committed
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
  file ('trimReads_vf/*') from trimReadsVersions.first()
  file ('alignReads_vf/*') from alignReadsVersions.first()
  file ('filterReads_vf/*') from filterReadsVersions.first()
  file ('convertReads_vf/*') from convertReadsVersions.first()
  file ('crossReads_vf/*') from crossReadsVersions.first()
  file ('callPeaksMACS_vf/*') from callPeaksMACSVersions.first()
  file ('consensusPeaks_vf/*') from consensusPeaksVersions.first()
  file ('peakAnnotation_vf/*') from peakAnnotationVersions.first()
  file ('motifSearch_vf/*') from motifSearchVersions.first().ifEmpty()
  file ('diffPeaks_vf/*') from diffPeaksVersions.first().ifEmpty()
  file ('experimentQC_vf/*') from experimentQCVersions.first()
  file ('trimReads/*') from trimgaloreResults.collect()
  file ('alignReads/*') from mappedReadsStats.collect()
  file ('filterReads/*') from dedupReadsComplexity.collect()
  file ('crossReads/*') from crossReadsStats.collect()
Venkat Malladi's avatar
Venkat Malladi committed
641

642
  output:
Venkat Malladi's avatar
Venkat Malladi committed
643

Venkat Malladi's avatar
Venkat Malladi committed
644
645
  file('software_versions_mqc.yaml') into softwareVersions
  file('software_references_mqc.yaml') into softwareReferences
Venkat Malladi's avatar
Venkat Malladi committed
646
  file "multiqc_report.html" into multiqcReport
Venkat Malladi's avatar
Venkat Malladi committed
647
  file "*_data" into multiqcData
648
649

  script:
650

651
  """
Venkat Malladi's avatar
Venkat Malladi committed
652
  echo $workflow.nextflow.version > version_nextflow.txt
Venkat Malladi's avatar
Venkat Malladi committed
653
  singularity exec /project/shared/bicf_workflow_ref/singularity_images/bicf-multiqc-2.0.0.img multiqc --version > version_multiqc.txt
654
  python --version &> version_python.txt
Venkat Malladi's avatar
Venkat Malladi committed
655
  python3 $baseDir/scripts/generate_references.py -r $references -o software_references
656
  python3 $baseDir/scripts/generate_versions.py -o software_versions
Venkat Malladi's avatar
Venkat Malladi committed
657
  singularity exec /project/shared/bicf_workflow_ref/singularity_images/bicf-multiqc-2.0.0.img multiqc -c $multiqc .
658
659
  """
}