test/cmdlinetest

   1 #!/bin/bash
   2
   3 # integrative subcmds testing
   4
   5 if [ $# -lt 1 ];then
   6     echo "Run all tests for subcommands of MACS2. Need 1 parameter for a tag name! A unique string combining date, time and MACS2 version is recommended. ./test.sh <TAG>"
   7     exit
   8 fi
   9
  10 # test all sub-commands
  11 TAG=$1
  12
  13 CHIP=CTCF_SE_ChIP_chr22_50k.bed.gz
  14 CTRL=CTCF_SE_CTRL_chr22_50k.bed.gz
  15
  16 CHIPPE=CTCF_PE_ChIP_chr22_50k.bam
  17 CTRLPE=CTCF_PE_CTRL_chr22_50k.bam
  18
  19 CHIPBEDPE=CTCF_PE_ChIP_chr22_50k.bedpe.gz
  20 CTRLBEDPE=CTCF_PE_CTRL_chr22_50k.bedpe.gz
  21
  22 CHIPCONTIGS50K=contigs50k.bed.gz
  23
  24 CHIPBIGSPEEDTEST=CTCF_12878_5M.bed.gz
  25 CTRLBIGSPEEDTEST=input_12878_5M.bed.gz
  26
  27 # callpeak
  28 echo "1. callpeak"
  29 echo "1.1 callpeak narrow"
  30
  31 mkdir ${TAG}_run_callpeak_narrow
  32
  33 macs2 callpeak -g 52000000 -t $CHIP -c $CTRL -n run_callpeak_narrow0 -B --outdir ${TAG}_run_callpeak_narrow --cutoff-analysis &> ${TAG}_run_callpeak_narrow/run_callpeak_narrow0.log
  34 macs2 callpeak -g 52000000 -t $CHIP -c $CTRL -n run_callpeak_narrow1 -B --outdir ${TAG}_run_callpeak_narrow --d-min 15 --call-summits &> ${TAG}_run_callpeak_narrow/run_callpeak_narrow1.log
  35 macs2 callpeak -g 52000000 -t $CHIP -c $CTRL -n run_callpeak_narrow2 -B --outdir ${TAG}_run_callpeak_narrow --nomodel --extsize 100 &> ${TAG}_run_callpeak_narrow/run_callpeak_narrow2.log
  36 macs2 callpeak -g 52000000 -t $CHIP -c $CTRL -n run_callpeak_narrow3 -B --outdir ${TAG}_run_callpeak_narrow --nomodel --extsize 100 --shift -50 &> ${TAG}_run_callpeak_narrow/run_callpeak_narrow3.log
  37 macs2 callpeak -g 52000000 -t $CHIP -c $CTRL -n run_callpeak_narrow4 -B --outdir ${TAG}_run_callpeak_narrow --nomodel --nolambda --extsize 100 --shift -50 &> ${TAG}_run_callpeak_narrow/run_callpeak_narrow4.log
  38 macs2 callpeak -g 52000000 -t $CHIP -c $CTRL -n run_callpeak_narrow5 -B --outdir ${TAG}_run_callpeak_narrow --scale-to large &> ${TAG}_run_callpeak_narrow/run_callpeak_narrow5.log
  39
  40 echo "1.2 callpeak broad"
  41 mkdir ${TAG}_run_callpeak_broad
  42
  43 macs2 callpeak -g 52000000 -t $CHIP -c $CTRL -n run_callpeak_broad -B --outdir ${TAG}_run_callpeak_broad --broad &> ${TAG}_run_callpeak_broad/run_callpeak_broad.log
  44
  45 echo "1.3 callpeak on PE narrow/broad"
  46
  47 mkdir  ${TAG}_run_callpeak_pe_narrow
  48
  49 macs2 callpeak -g 52000000 -f BAMPE -t $CHIPPE -c $CTRLPE -n run_callpeak_bampe_narrow -B --outdir ${TAG}_run_callpeak_pe_narrow --call-summits &> ${TAG}_run_callpeak_pe_narrow/run_callpeak_bampe_narrow.log
  50 macs2 callpeak -g 52000000 -f BEDPE -t $CHIPBEDPE -c $CTRLBEDPE -n run_callpeak_bedpe_narrow -B --outdir ${TAG}_run_callpeak_pe_narrow --call-summits &> ${TAG}_run_callpeak_pe_narrow/run_callpeak_bedpe_narrow.log
  51 macs2 callpeak -g 52000000 -f BEDPE -t $CHIPBEDPE -n run_callpeak_pe_narrow_onlychip -B --outdir ${TAG}_run_callpeak_pe_narrow &> ${TAG}_run_callpeak_pe_narrow/run_callpeak_pe_narrow_onlychip.log
  52
  53 mkdir  ${TAG}_run_callpeak_pe_broad
  54
  55 macs2 callpeak -g 52000000 -f BAMPE -t $CHIPPE -c $CTRLPE -n run_callpeak_bampe_broad -B --outdir ${TAG}_run_callpeak_pe_broad --broad &> ${TAG}_run_callpeak_pe_broad/run_callpeak_bampe_broad.log
  56 macs2 callpeak -g 52000000 -f BEDPE -t $CHIPBEDPE -c $CTRLBEDPE -n run_callpeak_bedpe_broad -B --outdir ${TAG}_run_callpeak_pe_broad --broad &> ${TAG}_run_callpeak_pe_broad/run_callpeak_bedpe_broad.log
  57
  58 # pileup
  59 echo "2. pileup"
  60
  61 mkdir ${TAG}_run_pileup
  62
  63 macs2 pileup -f BED -i $CHIP --extsize 200 --outdir ${TAG}_run_pileup -o run_pileup_ChIP.bed.bdg &> ${TAG}_run_pileup/run_pileup_ChIP.bed.log
  64 macs2 pileup -f BED -i $CTRL --extsize 200 --outdir ${TAG}_run_pileup -o run_pileup_CTRL.bed.bdg &> ${TAG}_run_pileup/run_pileup_CTRL.bed.log
  65
  66 macs2 pileup -f BAMPE -i $CHIPPE --outdir ${TAG}_run_pileup -o run_pileup_ChIPPE.bampe.bdg &> ${TAG}_run_pileup/run_pileup_ChIPPE.bampe.log
  67 macs2 pileup -f BAMPE -i $CTRLPE --outdir ${TAG}_run_pileup -o run_pileup_CTRLPE.bampe.bdg &> ${TAG}_run_pileup/run_pileup_CTRLPE.bampe.log
  68
  69 macs2 pileup -f BEDPE -i $CHIPBEDPE --outdir ${TAG}_run_pileup -o run_pileup_ChIPPE.bedpe.bdg &> ${TAG}_run_pileup/run_pileup_ChIPPE.bedpe.log
  70 macs2 pileup -f BEDPE -i $CTRLBEDPE --outdir ${TAG}_run_pileup -o run_pileup_CTRLPE.bedpe.bdg &> ${TAG}_run_pileup/run_pileup_CTRLPE.bedpe.log
  71
  72 # filterdup
  73 echo "3. filterdup"
  74
  75 mkdir ${TAG}_run_filterdup
  76
  77 macs2 filterdup -g 52000000 -i $CHIP --outdir ${TAG}_run_filterdup -o run_filterdup_result.bed --dry-run &> ${TAG}_run_filterdup/run_filterdup_d.log
  78 macs2 filterdup -g 52000000 -i $CHIP --outdir ${TAG}_run_filterdup -o run_filterdup_result.bed &> ${TAG}_run_filterdup/run_filterdup.log
  79
  80 macs2 filterdup -g 52000000 -i $CHIPPE -f BAMPE --outdir ${TAG}_run_filterdup -o run_filterdup_result_pe.bedpe --dry-run &> ${TAG}_run_filterdup/run_filterdup_pe_d.log
  81 macs2 filterdup -g 52000000 -i $CHIPPE -f BAMPE --outdir ${TAG}_run_filterdup -o run_filterdup_result_pe.bedpe &> ${TAG}_run_filterdup/run_filterdup_pe.log
  82
  83 # predictd
  84 echo "4. predictd"
  85
  86 mkdir ${TAG}_run_predictd
  87
  88 macs2 predictd -g 52000000 -i $CHIP --d-min 10 --outdir ${TAG}_run_predictd --rfile run_predictd.R &> ${TAG}_run_predictd/run_predictd.log
  89
  90 # randsample
  91 echo "5. randsample"
  92
  93 mkdir ${TAG}_run_randsample
  94
  95 macs2 randsample -i $CHIP -n 10000 --seed 31415926 --outdir ${TAG}_run_randsample -o run_randsample.bed &> ${TAG}_run_randsample/run_randsample.log
  96
  97 macs2 randsample -f BAMPE -i $CHIPPE -n 10000 --seed 31415926 --outdir ${TAG}_run_randsample -o run_randsample_bampe.bedpe &> ${TAG}_run_randsample/run_randsample_bampe.log
  98
  99 macs2 randsample -f BEDPE -i $CHIPBEDPE -n 10000 --seed 31415926 --outdir ${TAG}_run_randsample -o run_randsample_bedpe.bedpe &> ${TAG}_run_randsample/run_randsample_bedpe.log
 100
 101 # refinepeak
 102 echo "6. refinepeak"
 103
 104 mkdir ${TAG}_run_refinepeak
 105
 106 macs2 refinepeak -b ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_peaks.narrowPeak -i $CHIP --outdir ${TAG}_run_refinepeak --o-prefix run_refinepeak_w_prefix &> ${TAG}_run_refinepeak/run_refinepeak_w_prefix.log
 107
 108 macs2 refinepeak -b ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_peaks.narrowPeak -i $CHIP --outdir ${TAG}_run_refinepeak -o run_refinepeak_w_ofile.bed &> ${TAG}_run_refinepeak/run_refinepeak_w_ofile.log
 109
 110 # bdgcmp
 111 echo "7. bdgcmp"
 112
 113 mkdir ${TAG}_run_bdgcmp
 114
 115 macs2 bdgcmp -t ${TAG}_run_pileup/run_pileup_ChIP.bed.bdg -c ${TAG}_run_pileup/run_pileup_CTRL.bed.bdg  -m ppois FE -p 1 --outdir ${TAG}_run_bdgcmp --o-prefix run_bdgcmp &> ${TAG}_run_bdgcmp/run_bdgcmp.log
 116
 117 # bdgpeakcall
 118 echo "8. bdgpeakcall"
 119
 120 mkdir ${TAG}_run_bdgpeakcall
 121
 122 macs2 bdgpeakcall -i ${TAG}_run_bdgcmp/run_bdgcmp_FE.bdg -c 2 --outdir ${TAG}_run_bdgpeakcall --o-prefix run_bdgpeakcall_w_prefix &> ${TAG}_run_bdgpeakcall/run_bdgpeakcall_w_prefix.log
 123
 124 macs2 bdgpeakcall -i ${TAG}_run_bdgcmp/run_bdgcmp_FE.bdg -c 2 --outdir ${TAG}_run_bdgpeakcall -o run_bdgpeakcall_cutoff.txt --cutoff-analysis &> ${TAG}_run_bdgpeakcall/run_bdgpeakcall_cutoff.log
 125
 126 # bdgbroadcall
 127 echo "9. bdgbroadcall"
 128
 129 mkdir ${TAG}_run_bdgbroadcall
 130
 131 macs2 bdgbroadcall -i ${TAG}_run_bdgcmp/run_bdgcmp_FE.bdg -c 2 -C 1.5 --outdir ${TAG}_run_bdgbroadcall --o-prefix run_bdgbroadcall_w_prefix &> ${TAG}_run_bdgbroadcall/run_bdgbroadcall_w_prefix.log
 132
 133 # bdgdiff
 134 echo "10. bdgdiff"
 135
 136 mkdir ${TAG}_run_callpeak_narrow_revert
 137 mkdir ${TAG}_run_bdgdiff
 138
 139 macs2 callpeak -g 10000000 --nomodel --extsize 250 -c $CHIP -t $CTRL -n run_callpeak_narrow_revert -B --outdir ${TAG}_run_callpeak_narrow_revert &> ${TAG}_run_callpeak_narrow_revert/run_callpeak_narrow_revert.log
 140
 141 macs2 bdgdiff --t1 ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_treat_pileup.bdg --c1 ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_control_lambda.bdg --t2 ${TAG}_run_callpeak_narrow_revert/run_callpeak_narrow_revert_treat_pileup.bdg --c2 ${TAG}_run_callpeak_narrow_revert/run_callpeak_narrow_revert_control_lambda.bdg --o-prefix run_bdgdiff_prefix --outdir ${TAG}_run_bdgdiff &> ${TAG}_run_bdgdiff/run_bdgdiff_w_prefix.log
 142
 143 macs2 bdgdiff --t1 ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_treat_pileup.bdg --c1 ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_control_lambda.bdg --t2 ${TAG}_run_callpeak_narrow_revert/run_callpeak_narrow_revert_treat_pileup.bdg --c2 ${TAG}_run_callpeak_narrow_revert/run_callpeak_narrow_revert_control_lambda.bdg -o cond1.bed cond2.bed common.bed --outdir ${TAG}_run_bdgdiff &> ${TAG}_run_bdgdiff/run_bdgdiff_w_o_file.log
 144
 145 # cmbreps
 146 echo "11. cmbreps"
 147
 148 mkdir ${TAG}_run_cmbreps
 149
 150 macs2 cmbreps -i ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_treat_pileup.bdg ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_control_lambda.bdg ${TAG}_run_bdgcmp/run_bdgcmp_ppois.bdg -m max -o run_cmbreps_max.bdg --outdir ${TAG}_run_cmbreps &> ${TAG}_run_cmbreps/run_cmbreps_max.log
 151 macs2 cmbreps -i ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_treat_pileup.bdg ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_control_lambda.bdg ${TAG}_run_bdgcmp/run_bdgcmp_ppois.bdg -m mean -o run_cmbreps_mean.bdg --outdir ${TAG}_run_cmbreps &> ${TAG}_run_cmbreps/run_cmbreps_mean.log
 152 macs2 cmbreps -i ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_treat_pileup.bdg ${TAG}_run_callpeak_narrow/run_callpeak_narrow0_control_lambda.bdg ${TAG}_run_bdgcmp/run_bdgcmp_ppois.bdg -m fisher -o run_cmbreps_fisher.bdg --outdir ${TAG}_run_cmbreps &> ${TAG}_run_cmbreps/run_cmbreps_fisher.log
 153
 154 # test large amount of contigs
 155 echo "12. 50k contigs with buffersize"
 156
 157 mkdir ${TAG}_run_50kcontigs
 158
 159 echo "12.1 callpeak"
 160 macs2 callpeak -g 10000000 -t $CHIPCONTIGS50K -n run_callpeak_50kcontigs --outdir ${TAG}_run_50kcontigs --buffer-size 1000 --nomodel --extsize 200 &> ${TAG}_run_50kcontigs/run_callpeak_50kcontigs.log
 161
 162 echo "12.2 filterdup"
 163 macs2 filterdup -g 10000000 -i $CHIPCONTIGS50K --outdir ${TAG}_run_50kcontigs -o run_filterdup_result.bed --buffer-size 1000 &> ${TAG}_run_50kcontigs/run_filterdup_50kcontigs.log
 164
 165 echo "12.3 pileup"
 166 macs2 pileup -f BED -i $CHIPCONTIGS50K --extsize 200 --outdir ${TAG}_run_50kcontigs -o run_pileup_ChIP.bed.bdg --buffer-size 1000 &> ${TAG}_run_50kcontigs/run_pileup_ChIP.bed.log
 167
 168 echo "12.4 randsample"
 169 macs2 randsample -i $CHIPCONTIGS50K -n 100000 --seed 31415926 --outdir ${TAG}_run_50kcontigs -o run_randsample.bed --buffer-size 1000 &> ${TAG}_run_50kcontigs/run_randsample.log
 170
 171 echo "13. search for errors or warnings in log files"
 172 flag=0
 173 for i in `ls ${TAG}_run_*/*.log`;do
 174     echo " checking $i..."
 175     egrep -A 2 -B 2 -i "error|warning" $i;
 176     if [[ $? == 0 ]];then
 177         echo " ... error/warning found!"
 178         let "flag=1";
 179     else
 180         echo " ... clear!"
 181     fi
 182 done;
 183
 184 if [[ $flag == 1 ]]; then
 185     echo " Error or Warning can be found! Quit the test!"
 186     exit 1;
 187 fi
 188
 189 echo "14. compare with standard outputs"
 190
 191 flag=0
 192 # we do not check the consistency of random subsampling results
 193 subfolders=(50kcontigs bdgbroadcall bdgcmp bdgdiff bdgpeakcall callpeak_broad callpeak_narrow callpeak_narrow_revert callpeak_pe_broad callpeak_pe_narrow cmbreps filterdup pileup refinepeak)
 194
 195 m=0
 196 n=0
 197 for i in ${subfolders[@]};do
 198     let "m+=1"
 199     let "n=0"
 200     for j in `ls standard_results_${i}/*`; do
 201         k=`basename $j`
 202         let "n+=1"
 203         echo "14.${m}.${n} checking $i $k ..."
 204         fq=${TAG}_run_${i}/${k}
 205         fs=standard_results_${i}/${k}
 206         sort $fq > tmp_fq.txt
 207         sort $fs > tmp_fs.txt
 208         d=`diff tmp_fq.txt tmp_fs.txt`
 209         if [ -z "$d" ]; then
 210             echo " ... success!"
 211         else
 212             echo " ... failed! Difference:"
 213             diff tmp_fq.txt tmp_fs.txt | head -10
 214             let "flag=1";
 215         fi
 216     done
 217 done
 218 rm -f tmp_fq.txt tmp_fs.txt
 219
 220 # exit with 1 if test fails
 221 if [ $flag -eq 1 ];then
 222     exit $flag;
 223 fi
 224
 225 echo "15. brief speed test"
 226
 227 macs2 callpeak -t CTCF_12878_5M.bed.gz -c Input_12878_5M.bed.gz -n speedtest.py37 -B &> speedtest.py37.log &
 228 # monitor the macs2 callpeak run and output CPU and mem usage
 229 ./prockreport 1 $!
 230 wait;
 231
 232 # END of test