anchor.diff:  ctl_deeptools_bedgraph = deeptools( sam: ctl_sam );
   
   % extract fasta
  -peaks_fa     = bedtools-getfasta( fa: fa, bed: peaks );
  +peaks_fa     = bedtools-getfasta( fa: fa, fai: fai, bed: peaks );
   
  -% restrict peak regions to +/- 100 ps around the peak and extract fasta
  -peaks_100    = restrict-peaks( bed: peaks );
  -peaks_100_fa = bedtools-getfasta( fa: fa, bed: peaks_100 );
  +% extend regions to 100 bps around the peak region and extract fasta
  +peaks_100    = extend-peaks( bed: peaks );
  +peaks_100_fa = bedtools-getfasta( fa: fa, fai: fai, bed: peaks_100 );
   
   
   %% =========================================================
  
positive.code: % quality control
  qc = fastqc( fq: tag_fq ctl_fq );
  
  % read mapping
  idx     = bowtie-build( fa: fa );
  tag_sam = bowtie-align( idx: idx, fq: tag_fq );
  
negative.code:   bedgraph = bamcoverage( bam: dedup_sorted_bam, bai: dedup_sorted_bai );
  }
  
  deftask bedtools-getfasta( bed_fa( File ) : fa( File ) fai( File ) bed( File ) )in bash *{
    bed_fa=$bed.fa
    ln -sf $fai $fa.fai
    bedtools getfasta -fi $fa -bed $bed -fo $bed_fa
  }*
  
  deftask extend-peaks( restricted_bed( File ) : bed( File ) )in bash *{
    restricted_bed=$bed.100.bed
    perl -lane '$start=$F[1] > 100 ? $F[1]-100 : 0; $end = $F[2]+100 ; print "$F[0]\t$start\t$end"' \
      $bed > $restricted_bed
  }*
  
  
  %% =========================================================
  %% INPUT DATA
  