Advertisement
jhangyu

VCF to Circos

Jun 22nd, 2018
117
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
Bash 1.44 KB | None | 0 0
  1. And here is my script to process vcf file into Circos Plot image:
  2.  
  3. #Turn VCF Spots into number counts in region
  4. for ID in $(ls | grep -P -o '.*(?=(.vcf))')
  5. do
  6.     if [ ! -f ${ID}.dat.txt ] ;then
  7.         awk '/^#/ {next} {printf("%s\t%d\n",$1,$2-$2%10000);}' ${ID}.vcf | sort | uniq -c | awk '{printf("%s\t%s\t%d\t%s\n",$2,$3,$3+10000,$1);}' > ${ID}.dat.txt
  8.     fi
  9. done
  10.  
  11. for ID in $(ls | grep -P -o '.*(?=(.dat.txt))')
  12. do
  13.     if [ ! -f ${ID}.chrd.txt ] ;then
  14.         echo process ${ID} changing chromosome id to name
  15.         cp ${ID}.dat.txt ${ID}.chrd.txt
  16.         pids="$pids $!"
  17.         wait $pids
  18.         while read ref_seq chr
  19.         do
  20.             #Convert ref_seq ID into chromosome number
  21.             sed -i s/${ref_seq}/${chr}/g ${ID}.chrd.txt
  22.             pids="$pids $!"
  23.             wait $pids
  24.         done < ref_chr.list
  25.     fi
  26.     if [ ! -f ${ID}.png ] ;then
  27.         echo Create ${ID} config file for circos
  28.         #Input the plot file and the max point
  29.         maxpoint=$(sort -rnk4 ${ID}.chrd.txt | awk '!x[$2]++' | head -1 | cut -f 4)
  30.         sed s/plotfileinput.txt/${ID}.chrd.txt/g vcf_density.conf > ${ID}_vcf_density.conf
  31.         sed -i s/max_point_vale/${maxpoint}/g ${ID}_vcf_density.conf
  32.         pids="$pids $!"
  33.         wait $pids
  34.         echo Drawing ${ID} circos plot
  35.         circos -outputdir ./ -outputfile ${ID} -conf ${ID}_vcf_density.conf
  36.         pids="$pids $!"
  37.         wait $pids
  38.     fi
  39.     rm ${ID}.dat.txt
  40. done
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement