|
0
|
1 #!/bin/bash
|
|
|
2
|
|
|
3 vcffile=$1
|
|
|
4 outputfile=$2
|
|
|
5
|
|
|
6 # vcf columns: CHROM-POS-ID-REF-ALT
|
|
|
7 # LV cloumns: variantId-chromosome-start-end-reference-alleleSeq-xRef
|
|
|
8
|
|
|
9
|
|
|
10 # add chr prefix if not present
|
|
|
11 # determine varType (snp, ins, del, sub)
|
|
|
12 # convert coordinates to 0-based halfopen
|
|
|
13 # calculate end coordinate from position and length
|
|
|
14 # remove leading reference base from the non-SNP variants, update position
|
|
|
15
|
|
|
16 awk 'BEGIN{
|
|
|
17 FS="\t";
|
|
|
18 OFS="\t";
|
|
|
19 count=0;
|
|
|
20
|
|
|
21 #output new header
|
|
4
|
22 print "variantId", "chromosome", "begin", "end", "varType", "reference", "alleleSeq", "xRef" > "headerline.txt"
|
|
0
|
23 }{
|
|
|
24
|
|
|
25 if(substr($0,1,1)!="#" && $5 != "."){ #skip header or nonvariant entries (period in ALT column)
|
|
|
26
|
|
|
27 # detect multivariants
|
|
|
28 chrom=$1
|
|
|
29 pos=$2
|
|
|
30 ref=$4
|
|
|
31 #alt=$5
|
|
|
32 reflen=length($4)
|
|
|
33
|
|
4
|
34 # excel adds quotes sometimes :s
|
|
|
35 gsub(/"/,"",ref)
|
|
|
36 gsub(/"/,"",alt)
|
|
|
37
|
|
0
|
38 # add chr prefix if needed
|
|
|
39 if(substr($1,1,3)!="chr")
|
|
|
40 chromosome="chr"$1
|
|
|
41 else
|
|
|
42 chromosome=chrom
|
|
|
43
|
|
|
44 # split ALT column in case of multiple variant alleles
|
|
|
45 split($5,alleles,",");
|
|
|
46
|
|
|
47 for (i in alleles) {
|
|
|
48 alt=alleles[i]
|
|
|
49
|
|
|
50
|
|
|
51 # determine varType
|
|
|
52 if(length(ref) == 1 && length(alt) == 1)
|
|
|
53 varType="snp"
|
|
4
|
54 else if (length(ref) == 1 && substr(ref,1,1)==substr(alt,1,1) )
|
|
0
|
55 varType="ins"
|
|
4
|
56 else if (length(alt) == 1 && substr(ref,1,1)==substr(alt,1,1) )
|
|
0
|
57 varType="del"
|
|
|
58 else
|
|
|
59 varType="sub"
|
|
|
60
|
|
|
61 # determine start and end coordinates in 0-based half-open coordinate system
|
|
|
62
|
|
|
63 if (varType=="snp"){
|
|
|
64 start=pos-1
|
|
|
65 end=pos
|
|
|
66 }
|
|
|
67 else if (varType=="ins"){
|
|
|
68 start=pos
|
|
4
|
69 end=pos
|
|
0
|
70 }
|
|
|
71 else if (varType=="del"){
|
|
|
72 start=pos
|
|
|
73 end=pos+(reflen-1)
|
|
|
74 }
|
|
|
75 else if (varType=="sub"){
|
|
4
|
76 start=pos-1
|
|
0
|
77 end=pos+(reflen-1)
|
|
|
78 }
|
|
|
79
|
|
|
80 # remove leading reference base
|
|
4
|
81 if ( varType!="snp" && substr(ref,1,1)==substr(alt,1,1) ){ #subs not mandatory leading reference base :s
|
|
0
|
82 reference=substr(ref,2)
|
|
|
83 alleleSeq=substr(alt,2)
|
|
4
|
84 if (varType =="sub"){
|
|
|
85 start+=1
|
|
|
86 }
|
|
0
|
87 }
|
|
|
88 else{
|
|
|
89 reference=ref
|
|
|
90 alleleSeq=alt
|
|
|
91 }
|
|
|
92
|
|
|
93 #print output variant(s)
|
|
|
94
|
|
4
|
95 if(chromosome == "chr1" || chromosome == "chr2" || chromosome == "chr3" || chromosome == "chr4" || chromosome == "chr5" || chromosome == "chr6" || chromosome == "chr7" || chromosome == "chr8" || chromosome == "chr9" || chromosome == "chr10" || chromosome == "chr11" || chromosome == "chr12" || chromosome == "chr13" || chromosome == "chr14" || chromosome == "chr15" || chromosome == "chr16" || chromosome == "chr17" ||chromosome == "chr18" ||chromosome == "chr19" ||chromosome == "chr20" ||chromosome == "chr21" ||chromosome == "chr22" ||chromosome == "chrX" ||chromosome == "chrY" )
|
|
|
96 print count, chromosome, start, end, varType, reference, alleleSeq, ""
|
|
0
|
97
|
|
|
98 count+=1
|
|
|
99 }
|
|
|
100 }
|
|
4
|
101 }END{}' $vcffile > $outputfile.almost
|
|
0
|
102
|
|
|
103
|
|
4
|
104
|
|
|
105 # due to overlapping variants that we reduce to more canonical forms, variants may have become out of order, so resort to be sure
|
|
|
106 sort -k2,2V -k3,3g $outputfile.almost > $outputfile.almost2
|
|
|
107
|
|
|
108 cat headerline.txt $outputfile.almost2 > $outputfile
|
|
|
109
|
|
|
110
|
|
|
111
|
|
0
|
112
|
|
|
113 #from 100Genomes site:
|
|
|
114
|
|
|
115 #CHROM chromosome: an identifier from the reference genome. All entries for a specific CHROM should form a contiguous block within the VCF file.(Alphanumeric String, Required)
|
|
|
116 #POS position: The reference position, with the 1st base having position 1. Positions are sorted numerically, in increasing order, within each reference sequence CHROM. (Integer, Required)
|
|
|
117 #ID semi-colon separated list of unique identifiers where available. If this is a dbSNP variant it is encouraged to use the rs number(s). No identifier should be present in more than one data record. If there is no identifier available, then the missing value should be used. (Alphanumeric String)
|
|
|
118 #REF reference base(s): Each base must be one of A,C,G,T,N. Bases should be in uppercase. Multiple bases are permitted. The value in the POS field refers to the position of the first base in the String. For InDels, the reference String must include the base before the event (which must be reflected in the POS field). (String, Required).
|
|
4
|
119 #ALT comma separated list of alternate non-reference alleles called on at least one of the samples. Options are base Strings made up of the bases A,C,G,T,N, or an angle-bracketed ID String (”<ID>”). If there are no alternative alleles, then the missing value should be used. Bases should be in uppercase. (Alphanumeric String; no whitespace, commas, or angle-brackets are permitted in the ID String itself) |
