Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Frameshift variants may interfere with typing in type_vcf.py, even with low ALT_FREQ #92

Open
dfornika opened this issue Jan 16, 2021 · 7 comments
Open

Frameshift variants may interfere with typing in type_vcf.py, even with low ALT_FREQ #92

dfornika opened this issue Jan 16, 2021 · 7 comments

Comments

@dfornika
Copy link
Contributor

dfornika commented Jan 16, 2021 •
edited
Loading

Running type_vcf.py -dp 10 -af 0.75 on the following input:

REGION	POS	REF	ALT	REF_DP	REF_RV	REF_QUAL	ALT_DP	ALT_RV	ALT_QUAL	ALT_FREQ	TOTAL_DP	PVAL	PASS	GFF_FEATURE	REF_CODON	REF_AA	ALT_CODON	ALT_AA
MN908947.3	241	C	T	0	0	0	653	118	44	1	653	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	913	C	T	1	0	39	590	289	42	0.998308	591	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	3037	C	T	0	0	0	247	162	43	1	247	5.21659e-158	TRUE	NA	NA	NA	NA	NA
MN908947.3	3267	C	T	0	0	0	341	2	42	1	341	2.44084e-216	TRUE	NA	NA	NA	NA	NA
MN908947.3	5388	C	A	0	0	0	153	2	40	0.993506	154	5.24955e-95	TRUE	NA	NA	NA	NA	NA
MN908947.3	5986	C	T	0	0	0	155	70	40	1	155	6.34799e-96	TRUE	NA	NA	NA	NA	NA
MN908947.3	6954	T	C	0	0	0	250	167	40	1	250	3.29933e-154	TRUE	NA	NA	NA	NA	NA
MN908947.3	11287	G	-TCTGGTTTT	618	522	41	655	0	20	0.967504	677	1.10773e-145	TRUE	NA	NA	NA	NA	NA
MN908947.3	14408	C	T	0	0	0	260	163	43	1	260	2.14099e-166	TRUE	NA	NA	NA	NA	NA
MN908947.3	14676	C	T	0	0	0	110	1	42	1	110	7.17302e-69	TRUE	NA	NA	NA	NA	NA
MN908947.3	15279	C	T	0	0	0	127	79	41	1	127	2.25986e-79	TRUE	NA	NA	NA	NA	NA
MN908947.3	16176	T	C	0	0	0	217	110	40	1	217	6.3057e-135	TRUE	NA	NA	NA	NA	NA
MN908947.3	17066	A	T	329	163	43	195	99	43	0.372137	524	7.79486e-76	TRUE	NA	NA	NA	NA	NA
MN908947.3	17609	A	-T	346	292	44	135	0	20	0.335821	402	1.56485e-34	TRUE	NA	NA	NA	NA	NA
MN908947.3	17615	A	G	0	0	0	337	283	44	1	337	2.76713e-218	TRUE	NA	NA	NA	NA	NA
MN908947.3	21701	G	+T	21	1	51	10	0	20	0.357143	28	0.000849992	TRUE	NA	NA	NA	NA	NA
MN908947.3	21764	A	-TACATG	18	6	47	22	0	20	1	22	5.12135e-06	TRUE	NA	NA	NA	NA	NA
MN908947.3	21834	A	-T	23	8	39	7	0	20	0.291667	24	0.0132072	FALSE	NA	NA	NA	NA	NA
MN908947.3	21990	T	-TTA	19	8	41	20	0	20	0.909091	22	1.64422e-05	TRUE	NA	NA	NA	NA	NA
MN908947.3	22000	C	+A	21	10	39	9	0	20	0.409091	22	0.00470272	TRUE	NA	NA	NA	NA	NA
MN908947.3	23063	A	T	0	0	0	235	104	39	1	235	9.76466e-144	TRUE	NA	NA	NA	NA	NA
MN908947.3	23271	C	A	0	0	0	209	132	41	1	209	2.13804e-130	TRUE	NA	NA	NA	NA	NA
MN908947.3	23403	A	G	0	0	0	199	137	44	1	199	1.41954e-129	TRUE	NA	NA	NA	NA	NA
MN908947.3	23604	C	A	1	0	77	92	26	39	0.989247	93	2.18401e-54	TRUE	NA	NA	NA	NA	NA
MN908947.3	23709	C	T	0	0	0	102	25	44	1	102	5.7399e-66	TRUE	NA	NA	NA	NA	NA
MN908947.3	24506	T	G	0	0	0	120	86	41	1	120	5.37115e-76	TRUE	NA	NA	NA	NA	NA
MN908947.3	24914	G	C	0	0	0	273	95	39	1	273	9.11315e-170	TRUE	NA	NA	NA	NA	NA
MN908947.3	27972	C	T	0	0	0	513	96	44	1	513	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	28048	G	T	0	0	0	744	188	43	1	744	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	28111	A	G	0	0	0	839	212	42	1	839	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	28270	T	-A	731	344	41	686	0	20	0.873885	785	2.59338e-153	TRUE	NA	NA	NA	NA	NA
MN908947.3	28280	G	C	0	0	0	618	293	40	1	618	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	28281	A	T	0	0	0	620	297	41	1	620	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	28282	T	A	0	0	0	634	309	41	1	634	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	28881	G	A	1	1	33	739	446	42	0.998649	740	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	28882	G	A	1	1	33	743	452	41	0.998656	744	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	28883	G	C	1	1	33	747	454	42	0.998663	748	0	TRUE	NA	NA	NA	NA	NA
MN908947.3	28977	C	T	0	0	0	501	351	40	1	501	9.95898e-309	TRUE	NA	NA	NA	NA	NA

...results in the following .csq.vcf:

##fileformat=VCFv4.2
##FILTER=<ID=PASS,Description="All filters passed">
##source=iVar
##contig=<ID=MN908947.3>
##INFO=<ID=DP,Number=1,Type=Integer,Description="Total Depth">
##FILTER=<ID=FAIL,Description="Result of p-value > 0.05">
##FORMAT=<ID=GT,Number=1,Type=String,Description="Genotype">
##FORMAT=<ID=REF_DP,Number=1,Type=Integer,Description="Depth of reference base">
##FORMAT=<ID=REF_RV,Number=1,Type=Integer,Description="Depth of reference base on reverse reads">
##FORMAT=<ID=REF_QUAL,Number=1,Type=Integer,Description="Mean quality of reference base">
##FORMAT=<ID=ALT_DP,Number=1,Type=Integer,Description="Depth of alternate base">
##FORMAT=<ID=ALT_RV,Number=1,Type=Integer,Description="Depth of alternate base on reverse reads">
##FORMAT=<ID=ALT_QUAL,Number=1,Type=String,Description="Mean quality of alternate base">
##FORMAT=<ID=ALT_FREQ,Number=1,Type=String,Description="Frequency of alternate base">
##bcftools/csqVersion=1.11+htslib-1.11
##bcftools/csqCommand=csq -v 0 -f MN908947.3.fasta -g MN908947.3.gff; Date=Fri Jan 15 15:08:49 2021
##INFO=<ID=BCSQ,Number=.,Type=String,Description="Haplotype-aware consequence annotation from BCFtools/csq, see http://samtools.github.io/bcftools/howtos/csq-calling.html for details. Format: Consequence|gene|transcript|biotype|strand|amino_acid_change|dna_change">
##FORMAT=<ID=BCSQ,Number=.,Type=Integer,Description="Bitmask of indexes to INFO/BCSQ, with interleaved first/second haplotype. Use \"bcftools query -f'[%CHROM\t%POS\t%SAMPLE\t%TBCSQ\n]'\" to translate.">
#CHROM	POS	ID	REF	ALT	QUAL	FILTER	INFO	FORMAT	SAMPLE
MN908947.3	241	.	C	T	.	PASS	DP=653	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ	1:0:0:0:653:118:44:1
MN908947.3	913	.	C	T	.	PASS	DP=591;BCSQ=synonymous|ORF1ab|ENSSAST00005000003|protein_coding|+|216S|913C>T,synonymous|ORF1ab|ENSSAST00005000002|protein_coding|+|216S|913C>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:1:0:39:590:289:42:0.998308:5
MN908947.3	3037	.	C	T	.	PASS	DP=247;BCSQ=synonymous|ORF1ab|ENSSAST00005000003|protein_coding|+|924F|3037C>T,synonymous|ORF1ab|ENSSAST00005000002|protein_coding|+|924F|3037C>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:247:162:43:1:5
MN908947.3	3267	.	C	T	.	PASS	DP=341;BCSQ=missense|ORF1ab|ENSSAST00005000003|protein_coding|+|1001T>1001I|3267C>T,missense|ORF1ab|ENSSAST00005000002|protein_coding|+|1001T>1001I|3267C>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:341:2:42:1:5
MN908947.3	5388	.	C	A	.	PASS	DP=154;BCSQ=missense|ORF1ab|ENSSAST00005000003|protein_coding|+|1708A>1708D|5388C>A,missense|ORF1ab|ENSSAST00005000002|protein_coding|+|1708A>1708D|5388C>A	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:153:2:40:0.993506:5
MN908947.3	5986	.	C	T	.	PASS	DP=155;BCSQ=synonymous|ORF1ab|ENSSAST00005000003|protein_coding|+|1907F|5986C>T,synonymous|ORF1ab|ENSSAST00005000002|protein_coding|+|1907F|5986C>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:155:70:40:1:5
MN908947.3	6954	.	T	C	.	PASS	DP=250;BCSQ=missense|ORF1ab|ENSSAST00005000003|protein_coding|+|2230I>2230T|6954T>C,missense|ORF1ab|ENSSAST00005000002|protein_coding|+|2230I>2230T|6954T>C	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:250:167:40:1:5
MN908947.3	11287	.	GTCTGGTTTT	G	.	PASS	DP=677;BCSQ=inframe_deletion|ORF1ab|ENSSAST00005000003|protein_coding|+|3674LSGF>3674L|11287GTCTGGTTTT>G,inframe_deletion|ORF1ab|ENSSAST00005000002|protein_coding|+|3674LSGF>3674L|11287GTCTGGTTTT>G	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:618:522:41:655:0:20:0.967504:5
MN908947.3	14408	.	C	T	.	PASS	DP=260;BCSQ=synonymous|ORF1ab|ENSSAST00005000002|protein_coding|+|4715L|14408C>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:260:163:43:1:1
MN908947.3	14676	.	C	T	.	PASS	DP=110;BCSQ=missense|ORF1ab|ENSSAST00005000002|protein_coding|+|4804P>4801L|14676C>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:110:1:42:1:1
MN908947.3	15279	.	C	T	.	PASS	DP=127;BCSQ=missense|ORF1ab|ENSSAST00005000002|protein_coding|+|5005T>5002I|15279C>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:127:79:41:1:1
MN908947.3	16176	.	T	C	.	PASS	DP=217;BCSQ=missense|ORF1ab|ENSSAST00005000002|protein_coding|+|5304L>5301P|16176T>C	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:217:110:40:1:1
MN908947.3	17066	.	A	T	.	PASS	DP=524;BCSQ=missense|ORF1ab|ENSSAST00005000002|protein_coding|+|5601I>5598F|17066A>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:329:163:43:195:99:43:0.372137:1
MN908947.3	17609	.	AT	A	.	PASS	DP=402;BCSQ=frameshift|ORF1ab|ENSSAST00005000002|protein_coding|+|5782IISLKHIKTNQLNALKCFIRVLSRMMFHLQLTGHK*>5779K*|17609AT>A+17615A>G	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:346:292:44:135:0:20:0.335821:1
MN908947.3	17615	.	A	G	.	PASS	DP=337;BCSQ=@17609	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:337:283:44:1:1
MN908947.3	21701	.	G	GT	.	PASS	DP=28;BCSQ=inframe_deletion|S|ENSSAST00005000004|protein_coding|+|47VLHSTQDLFLPFFSNVTWFHAIHVSGTNGTKRFDNPVLPFNDGVY>47VFTFNSGLVLTFLFQCYLVPCYLWDQWY*|21701G>GT+21764ATACATG>A+21834AT>A	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:21:1:51:10:0:20:0.357143:1
MN908947.3	21764	.	ATACATG	A	.	PASS	DP=22;BCSQ=@21701	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:18:6:47:22:0:20:1:1
MN908947.3	21834	.	AT	A	.	FAIL	DP=24;BCSQ=@21701	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:23:8:39:7:0:20:0.291667:1
MN908947.3	21990	.	TTTA	T	.	PASS	DP=22;BCSQ=*inframe_deletion|S|ENSSAST00005000004|protein_coding|+|143VY>141V|21990TTTA>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:19:8:41:20:0:20:0.909091:1
MN908947.3	22000	.	C	CA	.	PASS	DP=22;BCSQ=*frameshift|S|ENSSAST00005000004|protein_coding|+|146HKNNKSWMESEFRVYSSANNCTFEYVSQPFLMDLEGKQGNFKNLREFVFKNIDGYFKIYSKHTPINLVRDLPQGFSALEPLVDLPIGINITRFQTLLALHRSYLTPGDSSSGWTAGAAAYYVGYLQPRTFLLKYNENGTITDAVDCALDPLSETKCTLKSFTVEKGIYQTSNFRVQPTESIVRFPNITNLCPFGEVFNATRFASVYAWNRKRISNCVADYSVLYNSASFSTFKCYGVSPTKLNDLCFTNVYADSFVIRGDEVRQIAPGQTGKIADYNYKLPDDFTGCVIAWNSNNLDSKVGGNYNYLYRLFRKSNLKPFERDISTEIYQAGSTPCNGVEGFNCYFPLQSYGFQPTNGVGYQPYRVVVLSFELLHAPATVCGPKKSTNLVKNKCVNFNFNGLTGTGVLTESNKKFLPFQQFGRDIADTTDAVRDPQTLEILDITPCSFGGVSVITPGTNTSNQVAVLYQDVNCTEVPVAIHADQLTPTWRVYSTGSNVFQTRAGCLIGAEHVNNSYECDIPIGAGICASYQTQTNSPRRARSVASQSIIAYTMSLGAENSVAYSNNSIAIPTNFTISVTTEILPVSMTKTSVDCTMYICGDSTECSNLLLQYGSFCTQLNRALTGIAVEQDKNTQEVFAQVKQIYKTPPIKDFGGFNFSQILPDPSKPSKRSFIEDLLFNKVTLADAGFIKQYGDCLGDIAARDLICAQKFNGLTVLPPLLTDEMIAQYTSALLAGTITSGWTFGAGAALQIPFAMQMAYRFNGIGVTQNVLYENQKLIANQFNSAIGKIQDSLSSTASALGKLQDVVNQNAQALNTLVKQLSSNFGAISSVLNDILSRLDKVEAEVQIDRLITGRLQSLQTYVTQQLIRAAEIRASANLAATKMSECVLGQSKRVDFCGKGYHLMSFPQSAPHGVVFLHVTYVPAQEKNFTTAPAICHDGKAHFPREGVFVSNGTHWFVTQRNFYEPQIITTDNTFVSGNCDVVIGIVNNTVYDPLQPELDSFKEELDKYFKNHTSPDVDLGDISGINASVVNIQKEIDRLNEVAKNLNESLIDLQELGKYEQYIKWPWYIWLGFIAGLIAIVMVTIMLCCMTSCCSCLKGCCSCGSCCKFDEDDSEPVLKGVKLHYT*>143HKKQQKLDGK*|22000C>CA+23063A>T+23271C>A+23403A>G+23604C>A+23709C>T+24506T>G+24914G>C	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:21:10:39:9:0:20:0.409091:1
MN908947.3	23063	.	A	T	.	PASS	DP=235;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:235:104:39:1:1
MN908947.3	23271	.	C	A	.	PASS	DP=209;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:209:132:41:1:1
MN908947.3	23403	.	A	G	.	PASS	DP=199;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:199:137:44:1:1
MN908947.3	23604	.	C	A	.	PASS	DP=93;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:1:0:77:92:26:39:0.989247:1
MN908947.3	23709	.	C	T	.	PASS	DP=102;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:102:25:44:1:1
MN908947.3	24506	.	T	G	.	PASS	DP=120;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:120:86:41:1:1
MN908947.3	24914	.	G	C	.	PASS	DP=273;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:273:95:39:1:1
MN908947.3	27972	.	C	T	.	PASS	DP=513;BCSQ=stop_gained|ORF8|ENSSAST00005000008|protein_coding|+|27Q>27*|27972C>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:513:96:44:1:1
MN908947.3	28048	.	G	T	.	PASS	DP=744;BCSQ=*missense|ORF8|ENSSAST00005000008|protein_coding|+|52R>52I|28048G>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:744:188:43:1:1
MN908947.3	28111	.	A	G	.	PASS	DP=839;BCSQ=*missense|ORF8|ENSSAST00005000008|protein_coding|+|73Y>73C|28111A>G	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:839:212:42:1:1
MN908947.3	28270	.	TA	T	.	PASS	DP=785	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ	1:731:344:41:686:0:20:0.873885
MN908947.3	28280	.	G	C	.	PASS	DP=618;BCSQ=missense|N|ENSSAST00005000005|protein_coding|+|3D>3L|28280G>C+28281A>T+28282T>A	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:618:293:40:1:1
MN908947.3	28281	.	A	T	.	PASS	DP=620;BCSQ=@28280	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:620:297:41:1:1
MN908947.3	28282	.	T	A	.	PASS	DP=634;BCSQ=@28280	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:634:309:41:1:1
MN908947.3	28881	.	G	A	.	PASS	DP=740;BCSQ=missense|N|ENSSAST00005000005|protein_coding|+|203R>203K|28881G>A+28882G>A	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:1:1:33:739:446:42:0.998649:1
MN908947.3	28882	.	G	A	.	PASS	DP=744;BCSQ=@28881	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:1:1:33:743:452:41:0.998656:1
MN908947.3	28883	.	G	C	.	PASS	DP=748;BCSQ=missense|N|ENSSAST00005000005|protein_coding|+|204G>204R|28883G>C	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:1:1:33:747:454:42:0.998663:1
MN908947.3	28977	.	C	T	.	PASS	DP=501;BCSQ=missense|N|ENSSAST00005000005|protein_coding|+|235S>235F|28977C>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:501:351:40:1:1

Note the variant at position 22000:

REGION	POS	REF	ALT	REF_DP	REF_RV	REF_QUAL	ALT_DP	ALT_RV	ALT_QUAL	ALT_FREQ	TOTAL_DP	PVAL	PASS	GFF_FEATURE	REF_CODON	REF_AA	ALT_CODON	ALT_AA
MN908947.3	22000	C	+A	21	10	39	9	0	20	0.409091	22	0.00470272	TRUE	NA	NA	NA	NA	NA

...is included in the .csq.vcf file despite having ALT_FREQ of 0.409091. This introduces a frameshift in the S gene that interferes with consequence assignment downstream:

MN908947.3	22000	.	C	CA	.	PASS	DP=22;BCSQ=*frameshift|S|ENSSAST00005000004|protein_coding|+|146HKNNKSWMESEFRVYSSANNCTFEYVSQPFLMDLEGKQGNFKNLREFVFKNIDGYFKIYSKHTPINLVRDLPQGFSALEPLVDLPIGINITRFQTLLALHRSYLTPGDSSSGWTAGAAAYYVGYLQPRTFLLKYNENGTITDAVDCALDPLSETKCTLKSFTVEKGIYQTSNFRVQPTESIVRFPNITNLCPFGEVFNATRFASVYAWNRKRISNCVADYSVLYNSASFSTFKCYGVSPTKLNDLCFTNVYADSFVIRGDEVRQIAPGQTGKIADYNYKLPDDFTGCVIAWNSNNLDSKVGGNYNYLYRLFRKSNLKPFERDISTEIYQAGSTPCNGVEGFNCYFPLQSYGFQPTNGVGYQPYRVVVLSFELLHAPATVCGPKKSTNLVKNKCVNFNFNGLTGTGVLTESNKKFLPFQQFGRDIADTTDAVRDPQTLEILDITPCSFGGVSVITPGTNTSNQVAVLYQDVNCTEVPVAIHADQLTPTWRVYSTGSNVFQTRAGCLIGAEHVNNSYECDIPIGAGICASYQTQTNSPRRARSVASQSIIAYTMSLGAENSVAYSNNSIAIPTNFTISVTTEILPVSMTKTSVDCTMYICGDSTECSNLLLQYGSFCTQLNRALTGIAVEQDKNTQEVFAQVKQIYKTPPIKDFGGFNFSQILPDPSKPSKRSFIEDLLFNKVTLADAGFIKQYGDCLGDIAARDLICAQKFNGLTVLPPLLTDEMIAQYTSALLAGTITSGWTFGAGAALQIPFAMQMAYRFNGIGVTQNVLYENQKLIANQFNSAIGKIQDSLSSTASALGKLQDVVNQNAQALNTLVKQLSSNFGAISSVLNDILSRLDKVEAEVQIDRLITGRLQSLQTYVTQQLIRAAEIRASANLAATKMSECVLGQSKRVDFCGKGYHLMSFPQSAPHGVVFLHVTYVPAQEKNFTTAPAICHDGKAHFPREGVFVSNGTHWFVTQRNFYEPQIITTDNTFVSGNCDVVIGIVNNTVYDPLQPELDSFKEELDKYFKNHTSPDVDLGDISGINASVVNIQKEIDRLNEVAKNLNESLIDLQELGKYEQYIKWPWYIWLGFIAGLIAIVMVTIMLCCMTSCCSCLKGCCSCGSCCKFDEDDSEPVLKGVKLHYT*>143HKKQQKLDGK*|22000C>CA+23063A>T+23271C>A+23403A>G+23604C>A+23709C>T+24506T>G+24914G>C	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:21:10:39:9:0:20:0.409091:1
MN908947.3	23063	.	A	T	.	PASS	DP=235;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:235:104:39:1:1
MN908947.3	23271	.	C	A	.	PASS	DP=209;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:209:132:41:1:1
MN908947.3	23403	.	A	G	.	PASS	DP=199;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:199:137:44:1:1
MN908947.3	23604	.	C	A	.	PASS	DP=93;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:1:0:77:92:26:39:0.989247:1
MN908947.3	23709	.	C	T	.	PASS	DP=102;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:102:25:44:1:1
MN908947.3	24506	.	T	G	.	PASS	DP=120;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:120:86:41:1:1
MN908947.3	24914	.	G	C	.	PASS	DP=273;BCSQ=@22000	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:273:95:39:1:1
MN908947.3	27972	.	C	T	.	PASS	DP=513;BCSQ=stop_gained|ORF8|ENSSAST00005000008|protein_coding|+|27Q>27*|27972C>T	GT:REF_DP:REF_RV:REF_QUAL:ALT_DP:ALT_RV:ALT_QUAL:ALT_FREQ:BCSQ	1:0:0:0:513:96:44:1:1
MN908947.3	28048	.	G	T	.	PASS	DP=744;

resulting SAMPLE.variants.csv file is:

sampleID,gene,aa_var,dna_var
SAMPLE,ORF1ab,Syn.216S,C913T
SAMPLE,ORF1ab,Syn.924F,C3037T
SAMPLE,ORF1ab,T1001I,C3267T
SAMPLE,ORF1ab,A1708D,C5388A
SAMPLE,ORF1ab,Syn.1907F,C5986T
SAMPLE,ORF1ab,I2230T,T6954C
SAMPLE,ORF1ab,LSGF3674L,GTCTGGTTTT11287G
SAMPLE,ORF1ab,Syn.4715L,C14408T
SAMPLE,ORF1ab,P4804L,C14676T
SAMPLE,ORF1ab,T5005I,C15279T
SAMPLE,ORF1ab,L5304P,T16176C
SAMPLE,S,VY143V,TTTA21990T
SAMPLE,ORF8,Q27*,C27972T
SAMPLE,ORF8,R52I,G28048T
SAMPLE,ORF8,Y73C,A28111G
SAMPLE,N,D3L,G28280C
SAMPLE,N,R203K,G28881A
SAMPLE,N,G204R,G28883C
SAMPLE,N,S235F,C28977T

...as a result, the VOC-202012/01 variant is not called.

Should variants be excluded from the .csq.vcf file if their ALT_FREQ is below --allele_freq_threshold? I've implemented that in another repo where I've extracted out the typing portion of this pipeline and it does recover the correct variant calling for this example.

https://github.com/BCCDC-PHL/type-variants-nf/pull/4/files
@m-bull
Copy link
Contributor

m-bull commented Jan 18, 2021

Good catch - I guess the problem is that I run bcftools csq on the unfiltered vcf (which already breaks calls downstream of a low-AF frameshift), then filter the broken calls for AF later - by this point the damage is already done.

I think I'd like to output the un-AF-filtered vcf somewhere, but as in https://github.com/BCCDC-PHL/type-variants-nf/pull/4/files, i'll add the AF filter to the vcf before calling bcftools csq - does that seem like an OK solution?

@dfornika
Copy link
Contributor Author

That sounds good to me.

@ahmedmagds
Copy link

@dfornika It is my first time to use the pipeline and I tried to annotate using the -gff option but I got an error. I simply downloaded the reference.gff file from NCBI so not sure if this is the right file or there is something else that I need to use. It would be great if you can share the gff file that you use or guide me where I can find it. Thanks

@dfornika
Copy link
Contributor Author

@ahmedmagds look in the typing directory of this repository

@ahmedmagds
Copy link

@dfornika Thank you!

@pmenzel
Copy link

pmenzel commented Mar 4, 2021

Having the same issue with a low-frequency frameshift upstream of N501Y etc, which in turn are not called.

@pmenzel
Copy link

pmenzel commented Mar 4, 2021

Good catch - I guess the problem is that I run bcftools csq on the unfiltered vcf (which already breaks calls downstream of a low-AF frameshift), then filter the broken calls for AF later - by this point the damage is already done.

I think I'd like to output the un-AF-filtered vcf somewhere, but as in https://github.com/BCCDC-PHL/type-variants-nf/pull/4/files, i'll add the AF filter to the vcf before calling bcftools csq - does that seem like an OK solution?

Yeah, that would be great to have the original vcf file with all detected variants kept somewhere.

pmenzel added a commit to LaborBerlin/LB-ncov2019-artic-nf that referenced this issue Mar 18, 2021
@pmenzel
ignore frameshifts below AF treshold
b364ffd 
this patch fixes issue connor-lab#92
by applying the modifications from BCCDC-PHL/type-variants-nf#4

this commit also includes a patch for issue connor-lab#90
using the modification in connor-lab#91
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Assignees
No one assigned
Labels
None yet
Projects
None yet
Milestone
No milestone
Development

No branches or pull requests
4 participants
@dfornika @pmenzel @m-bull @ahmedmagds