**I need to add a context dependent custom tags in the FORMAT and SAMPLE column of a vcf file.**

**Here is my sample file:**

    #CHROM	POS	ID	REF	ALT	QUAL	FILTER	INFO	FORMAT	2ms01e	
    2	5187	.	G	C	2535.90	PASS		GT:AD:DP:GQ:PL:PG:PB:PI:PW:PC	0/1:47,17:64:99:965,0,1925:1|0:.,.,.,.,.,.,.,.:662:|:0.5	
    2	5207	.	T	A	2814.90	PASS		GT:AD:DP:GQ:PL:PG:PB:PI:PW:PC	0/1:41,30:71:99:955,0,1628:1|0:.,.,.,.,.,.,.,.:662:|:0.5		
    2	5238	.	G	C	8321.94	PASS		GT:AD:DP:GQ:PL:PG:PB:PI:PW:PC	0/1:48,59:107:99:1986,0,1336:0|1:.,.,.,.,.,.,.,.:662:|:0.5		
    2	5258	.	T	C	10415.90	PASS		GT:AD:DP:GQ:PL:PG:PB:PI:PW:PC	0/1:75,45:120:99:1369,0,2434:1|0:.,.,.,.,.,.,.,.:662:|:0.5	

**Step 01:** Look for phasing information at a particular variant site. For line #1, GT is 0/1 and PG is 1|0

**Step 02:** I need to take the information from the PG field (phased genotype) and update the vcf by creating a new field HP (haplotype phase) at the end of the FORMAT. The HP is the default output from GATK RBphasing but my other tool cannot output this field.


**The output file should be this:** 

    #CHROM	POS	ID	REF	ALT	QUAL	FILTER	INFO	FORMAT	2ms01e	
    2	5187	.	G	C	2535.90	PASS		GT:AD:DP:GQ:PL:PG:PB:PI:PW:PC:HP 0/1:47,17:64:99:965,0,1925:1|0:.,.,.,.,.,.,.,.:662:|:0.5:5187-2,5187-1	
    2	5207	.	T	A	2814.90	PASS		GT:AD:DP:GQ:PL:PG:PB:PI:PW:PC:HP 0/1:41,30:71:99:955,0,1628:1|0:.,.,.,.,.,.,.,.:662:|:0.5:5207-2,5207-1		
    2	5238	.	G	C	8321.94	PASS		GT:AD:DP:GQ:PL:PG:PB:PI:PW:PC:HP 0/1:48,59:107:99:1986,0,1336:0|1:.,.,.,.,.,.,.,.:662:|:0.5:5238-1,5238-2		
    2	5258	.	T	C	10415.90	PASS		GT:AD:DP:GQ:PL:PG:PB:PI:PW:PC:HP 0/1:75,45:120:99:1369,0,2434:1|0:.,.,.,.,.,.,.,.:662:|:0.5:5258-2,5258-1	

**Explanation of the process:**
1) We need to look at the PG field first and see if the particular locus is phased with any variants before-after it. If its phased the value for the PG field will have a pipe (|). In the first line the GT (genotype) is 0/1 (heterozygous) and the PG is 1|0.

2) Identify the positions (POS) of the heterozygous site. For the first line its 5187. So, we need to pick this value and connect to value with '-' to the value obtained in step 3.
- this steps store the value of the PG field as ':5187-'

3) So, with PG 1|0, the HP field should indicate the alternate allele (1) is reported first with value 2.
- this step store the value for the field as ':5178-2'

4) Now, complete the HP field to indicate the value of another allele, with reported phase value as 1.
- this step stores the value for the PG field as ':5178-2,5187-1'

So, `GT:AD:DP:GQ:PL:PG:PB:PI:PW:PC	0/1:47,17:64:99:965,0,1925:1|0:.,.,.,.,.,.,.,.:662:|:0.5`	
translates to 
`GT:AD:DP:GQ:PL:PG:PB:PI:PW:PC	:HP 0/1:47,17:64:99:965,0,1925:1|0:.,.,.,.,.,.,.,.:662:|:0.5:5187-2,5187-1`

	

**Note:**
- At some site the variant might be 1/2 (or say 1|2 in PG), so this needs to be considered too.
- I tried to looks for several variant annotation software including PESEQ, Annnovar, vcftools etc. but not finding any thing that can specifically do this. If there is any such tool that is even more great.
- It would be best if the program can be in python with description of what is happening at each step. I am trying to learn python so that would be of great help too.


Thanks much in advance !
- Bishwa K.



https://github.com/dshulgin2015/computing_python/blob/master/issue_1.py

Works well, if someone else wants to do the same.