|
0
|
1 #!/bin/bash
|
|
|
2 #$Id: bigwig2wig 23 2014-01-28 12:09:22Z jens $
|
|
|
3
|
|
|
4 #SCRIPT CONVERTS BIGWIG FILE TO FIXED-STEP WIGGLE FORMAT FILE
|
|
|
5 #RESOLUTION IS CONTROLLED THROUGH THE BIN SIZE
|
|
|
6
|
|
|
7 #default bin_size
|
|
|
8 bin_size=500
|
|
|
9 mylab="wiggle file"
|
|
|
10
|
|
|
11 #parse input
|
|
|
12 while getopts hf:b:l: myarg
|
|
|
13 do case "$myarg" in
|
|
|
14 h) echo "Usage: bigwig_correlation -f <bigwig_file> -b <bin_size>"
|
|
|
15 echo "Ex: bigwig_correlation -f <MYFILE.bw> -b 600"
|
|
|
16 exit ;;
|
|
|
17 f) bigwig_file="$OPTARG" ;; #required
|
|
|
18 l) mylab="$OPTARG" ;; #optional
|
|
|
19 b) bin_size="$OPTARG" ;; #optional
|
|
|
20 [?]) echo "Usage: bigwig_correlation -f <MYFILE.bw> -b <bin_size>"
|
|
|
21 exit 1 ;;
|
|
|
22 esac
|
|
|
23 done
|
|
|
24
|
|
|
25 ###################################################
|
|
|
26 ###VALIDATE INPUT
|
|
|
27 ###################################################
|
|
|
28
|
|
|
29 #make tmp-filename to hold chromosome info
|
|
|
30 org_assembly_file=`mktemp -u`
|
|
|
31 bigWigInfo -chroms $bigwig_file | perl -ne "/^\tchr/ && print" | perl -pe "s/ +/\t/g" | cut -f2,4 > $org_assembly_file
|
|
|
32
|
|
|
33 #check bin_size & define step_size
|
|
|
34 bin_size_mod=`expr $bin_size % 2` #determine modulus
|
|
|
35 if [ $bin_size_mod -ne 0 ]; then
|
|
|
36 echo "Chosen bin_size must be an even positive number, added +1 to bin_size"
|
|
|
37 bin_size=$(($bin_size + 1))
|
|
|
38 fi
|
|
|
39
|
|
|
40 if [ $bin_size -lt 100 ]; then
|
|
|
41 echo "ERROR: Chosen bin_size must be a positive number >=100"
|
|
|
42 exit 1
|
|
|
43 fi
|
|
|
44 #set stetp size equal to bin size i.e. non-overlapping intervals
|
|
|
45 step_size=$bin_size
|
|
|
46
|
|
|
47 ###################################################
|
|
|
48 ###EXTRACT DENSITIES FROM NORMALIZED BIGWIG FILES
|
|
|
49 ###################################################
|
|
|
50
|
|
|
51
|
|
|
52 #make track definition line
|
|
|
53 echo "track type=wiggle_0 name=$mylab description=\"fixedStep format\""
|
|
|
54
|
|
|
55 #for each chromsome
|
|
|
56 cat $org_assembly_file | while read line; do
|
|
|
57
|
|
|
58 cur_chr=`echo $line | cut --delimiter=" " -f1`
|
|
|
59 cur_length=`echo $line | cut --delimiter=" " -f2`
|
|
|
60
|
|
|
61 n_bins=`echo "scale=0; (${cur_length}-${step_size})/${bin_size}" | bc`
|
|
|
62
|
|
|
63 start=1
|
|
|
64 stop=`echo "$n_bins * $bin_size" | bc`
|
|
|
65
|
|
|
66 #write header line for each chromosome
|
|
|
67 echo "fixedStep chrom=$cur_chr start=$start step=$step_size span=$step_size"
|
|
|
68
|
|
|
69 #get densities along chr in n_bins with chosen bin_size and step_size (giving overlap in bins)
|
|
|
70 nice bigWigSummary $bigwig_file $cur_chr $start $stop $n_bins | perl -pe 's/\t/\n/g' | perl -pe "s/n\/a/0/"
|
|
|
71 #gives warning if no data in/for current chromosome
|
|
|
72
|
|
|
73 done
|
|
|
74
|
|
|
75 #rm tmp
|
|
|
76 rm $org_assembly_file
|
|
|
77
|
