Mercurial > repos > davidvanzessen > mutation_analysis

comparison wrapper.sh @ 62:4262e880472d draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author davidvanzessen
date Fri, 25 Mar 2016 04:39:18 -0400
parents 8bb4d6009e08
children a7381fd96dad
comparison
equal deleted inserted replaced
61:64e6a7803e07 62:4262e880472d
12 naive_output=$9 12 naive_output=$9
13 filter_unique=${10} 13 filter_unique=${10}
14 mkdir $outdir 14 mkdir $outdir
15 15
16 echo "---------------- read parameters ----------------" 16 echo "---------------- read parameters ----------------"
17 echo "---------------- read parameters ----------------" > $output
17 18
18 echo "unpacking IMGT file" 19 echo "unpacking IMGT file"
19 20
20 type="`file $input`" 21 type="`file $input`"
21 if [[ "$type" == *"Zip archive"* ]] ; then 22 if [[ "$type" == *"Zip archive"* ]] ; then
40 #BLASTN_DIR="/home/galaxy/tmp/blast/ncbi-blast-2.2.30+/bin" 41 #BLASTN_DIR="/home/galaxy/tmp/blast/ncbi-blast-2.2.30+/bin"
41 42
42 echo "${BLASTN_DIR}" 43 echo "${BLASTN_DIR}"
43 44
44 echo "identification ($method)" 45 echo "identification ($method)"
46 echo "identification ($method)" >> $output
45 47
46 echo "blast or custom" 48 echo "blast or custom"
47 49
48 if [[ "${method}" == "custom" ]] ; then 50 if [[ "${method}" == "custom" ]] ; then
49 echo "custom" 51 echo "custom"
62 echo -e "qseqid\tsseqid\tpident\tlength\tmismatch\tgapopen\tqstart\tqend\tsstart\tsend\tevalue\tbitscore" > $outdir/identified_genes.txt 64 echo -e "qseqid\tsseqid\tpident\tlength\tmismatch\tgapopen\tqstart\tqend\tsstart\tsend\tevalue\tbitscore" > $outdir/identified_genes.txt
63 ${BLASTN_DIR}/blastn -task blastn -db $dir/subclass_definition.db -query $PWD/sequences.fasta -outfmt 6 >> $outdir/identified_genes.txt 65 ${BLASTN_DIR}/blastn -task blastn -db $dir/subclass_definition.db -query $PWD/sequences.fasta -outfmt 6 >> $outdir/identified_genes.txt
64 fi 66 fi
65 67
66 echo "---------------- merge_and_filter.r ----------------" 68 echo "---------------- merge_and_filter.r ----------------"
69 echo "---------------- merge_and_filter.r ----------------" >> $output
67 70
68 Rscript $dir/merge_and_filter.r $PWD/summary.txt $PWD/sequences.txt $PWD/mutationanalysis.txt $PWD/mutationstats.txt $PWD/hotspots.txt $outdir/identified_genes.txt $outdir/merged.txt $outdir/unmatched.txt $method $functionality $unique ${filter_unique} 71 Rscript $dir/merge_and_filter.r $PWD/summary.txt $PWD/sequences.txt $PWD/mutationanalysis.txt $PWD/mutationstats.txt $PWD/hotspots.txt $outdir/identified_genes.txt $outdir/merged.txt $outdir/unmatched.txt $method $functionality $unique ${filter_unique}
69 72
70 echo "---------------- mutation_analysis.r ----------------" 73 echo "---------------- mutation_analysis.r ----------------"
74 echo "---------------- mutation_analysis.r ----------------" >> $output
71 75
72 genes="ca,ca1,ca2,cg,cg1,cg2,cg3,cg4,cm" 76 genes="ca,ca1,ca2,cg,cg1,cg2,cg3,cg4,cm"
73 echo "R mutation analysis" 77 echo "R mutation analysis"
74 Rscript $dir/mutation_analysis.r $outdir/merged.txt $genes $outdir ${include_fr1} 2>&1 78 Rscript $dir/mutation_analysis.r $outdir/merged.txt $genes $outdir ${include_fr1} 2>&1
75 79
77 #exit 0 81 #exit 0
78 82
79 83
80 84
81 echo "---------------- mutation_analysis.py ----------------" 85 echo "---------------- mutation_analysis.py ----------------"
86 echo "---------------- mutation_analysis.py ----------------" >> $output
82 87
83 python $dir/mutation_analysis.py --input $outdir/merged.txt --genes $genes --includefr1 "${include_fr1}" --output $outdir/hotspot_analysis.txt 88 python $dir/mutation_analysis.py --input $outdir/merged.txt --genes $genes --includefr1 "${include_fr1}" --output $outdir/hotspot_analysis.txt
84 echo "R AA histogram" 89 echo "R AA histogram"
85 90
86 echo "---------------- aa_histogram.r ----------------" 91 echo "---------------- aa_histogram.r ----------------"
92 echo "---------------- aa_histogram.r ----------------" >> $output
87 93
88 Rscript $dir/aa_histogram.r $outdir/aa_mutations.txt $outdir/aa_histogram.png 2>&1 94 Rscript $dir/aa_histogram.r $outdir/aa_mutations.txt $outdir/aa_histogram.png 2>&1
89 95
90 genes=(ca ca1 ca2 cg cg1 cg2 cg3 cg4 cm) 96 genes=(ca ca1 ca2 cg cg1 cg2 cg3 cg4 cm)
91 97
92 funcs=(sum mean median) 98 funcs=(sum mean median)
93 99
94 100
95 echo "<html><center><h1>$title</h1></center>" >> $output 101 echo "<html><center><h1>$title</h1></center>" > $output
102
103 #display the matched/unmatched for clearity
104
105 matched_count="`cat $outdir/merged.txt | tail -n +2 | wc -l`"
106 unmatched_count="`cat $outdir/unmatched.txt | tail -n +2 | wc -l`"
107 total_count=$((matched_count + unmatched_count))
108 perc_count=$((unmatched_count / total_count * 100))
109 perc_count=`bc -l <<< "scale=2; ${unmatched_count} / ${total_count} * 100"`
110 perc_count=`bc -l <<< "scale=2; (${unmatched_count} / ${total_count} * 100 ) / 1"`
111
112 echo "<center><h2>Total: ${total_count}</h2></center>" >> $output
113 echo "<center><h2>Matched: ${matched_count} Unmatched: ${unmatched_count}</h2></center>" >> $output
114 echo "<center><h2>Percentage unmatched: ${perc_count}</h2></center>" >> $output
115
96 echo "---------------- main tables ----------------" 116 echo "---------------- main tables ----------------"
97 for func in ${funcs[@]} 117 for func in ${funcs[@]}
98 do 118 do
99 119
100 echo "---------------- $func table ----------------" 120 echo "---------------- $func table ----------------"