Classes

Class	Line #	Total Statements	Complexity	TOTAL Coverage	Actions
VCFLoaderTest	54	379	11	1.0100%	Show methods

Class VCFLoaderTest

Class VCFLoaderTest	Line # 54	Total Statements 379	Complexity 11	TOTAL Coverage 1.0100%
setUp() : void   setUp() : void	100100	5.05	1.01	1.0 1.0100%
setUpBeforeTest() : void   setUpBeforeTest() : void	113113	1.01	1.01	1.0 1.0100%
testDoLoad() : void   testDoLoad() : void	122122	82.082	2.02	1.0 1.0100%
makeVcfFile() : File   makeVcfFile() : File	250250	7.07	1.01	1.0 1.0100%
buildAlignment() : AlignmentI   buildAlignment() : AlignmentI	268268	46.046	1.01	1.0 1.0100%
testDoLoad_reverseStrand() : void   testDoLoad_reverseStrand() : void	381381	65.065	2.02	1.0 1.0100%
testDoLoad_vepCsq() : void   testDoLoad_vepCsq() : void	519519	134.0134	2.02	1.0 1.0100%
testLoadVCFContig() : void   testLoadVCFContig() : void	734734	39.039	1.01	1.0 1.0100%

 

Contributing tests

Contributing tests

Select all

	Test contribution	Test	Result	Actions
	0.46819338	jalview.io.vcf.VCFLoaderTest.testDoLoad_vepCsqjalview.io.vcf.VCFLoaderTest.testDoLoad_vepCsq	1PASS
	0.3562341	jalview.io.vcf.VCFLoaderTest.testDoLoadjalview.io.vcf.VCFLoaderTest.testDoLoad	1PASS
	0.3129771	jalview.io.vcf.VCFLoaderTest.testDoLoad_reverseStrandjalview.io.vcf.VCFLoaderTest.testDoLoad_reverseStrand	1PASS
	0.10178117	jalview.io.vcf.VCFLoaderTest.testLoadVCFContigjalview.io.vcf.VCFLoaderTest.testLoadVCFContig	1PASS

 

Source view

1		/*
2		* Jalview - A Sequence Alignment Editor and Viewer ($$Version-Rel$$)
3		* Copyright (C) $$Year-Rel$$ The Jalview Authors
4		*
5		* This file is part of Jalview.
6		*
7		* Jalview is free software: you can redistribute it and/or
8		* modify it under the terms of the GNU General Public License
9		* as published by the Free Software Foundation, either version 3
10		* of the License, or (at your option) any later version.
11		*
12		* Jalview is distributed in the hope that it will be useful, but
13		* WITHOUT ANY WARRANTY; without even the implied warranty
14		* of MERCHANTABILITY or FITNESS FOR A PARTICULAR
15		* PURPOSE. See the GNU General Public License for more details.
16		*
17		* You should have received a copy of the GNU General Public License
18		* along with Jalview. If not, see <http://www.gnu.org/licenses/>.
19		* The Jalview Authors are detailed in the 'AUTHORS' file.
20		*/
21		package jalview.io.vcf;
22
23		import static jalview.io.gff.SequenceOntologyI.SEQUENCE_VARIANT;
24		import static org.testng.Assert.assertEquals;
25		import static org.testng.Assert.assertNull;
26		import static org.testng.Assert.assertTrue;
27
28		import jalview.bin.Cache;
29		import jalview.bin.Console;
30		import jalview.datamodel.AlignmentI;
31		import jalview.datamodel.DBRefEntry;
32		import jalview.datamodel.Mapping;
33		import jalview.datamodel.Sequence;
34		import jalview.datamodel.SequenceFeature;
35		import jalview.datamodel.SequenceI;
36		import jalview.datamodel.features.FeatureAttributes;
37		import jalview.datamodel.features.SequenceFeatures;
38		import jalview.gui.AlignFrame;
39		import jalview.io.DataSourceType;
40		import jalview.io.FileLoader;
41		import jalview.io.gff.Gff3Helper;
42		import jalview.util.MapList;
43
44		import java.io.File;
45		import java.io.IOException;
46		import java.io.PrintWriter;
47		import java.util.List;
48		import java.util.Map;
49
50		import org.testng.annotations.BeforeClass;
51		import org.testng.annotations.BeforeTest;
52		import org.testng.annotations.Test;
53
		100% Uncovered Elements: 0 (393) Complexity: 11 Complexity Density: 0.03
54		public class VCFLoaderTest
55		{
56		private static final float DELTA = 0.00001f;
57
58		// columns 9717- of gene P30419 from Ensembl (much modified)
59		private static final String FASTA = "" +
60		/*
61		* forward strand 'gene' and 'transcript' with two exons
62		*/
63		">gene1/1-25 chromosome:GRCh38:17:45051610:45051634:1\n"
64		+ "CAAGCTGGCGGACGAGAGTGTGACA\n"
65		+ ">transcript1/1-18\n--AGCTGGCG----AGAGTGTGAC-\n"
66
67		/*
68		* reverse strand gene and transcript (reverse complement alleles!)
69		*/
70		+ ">gene2/1-25 chromosome:GRCh38:17:45051610:45051634:-1\n"
71		+ "TGTCACACTCTCGTCCGCCAGCTTG\n" + ">transcript2/1-18\n"
72		+ "-GTCACACTCT----CGCCAGCT--\n"
73
74		/*
75		* 'gene' on chromosome 5 with two transcripts
76		*/
77		+ ">gene3/1-25 chromosome:GRCh38:5:45051610:45051634:1\n"
78		+ "CAAGCTGGCGGACGAGAGTGTGACA\n"
79		+ ">transcript3/1-18\n--AGCTGGCG----AGAGTGTGAC-\n"
80		+ ">transcript4/1-18\n-----TGG-GGACGAGAGTGTGA-A\n";
81
82		private static final String[] VCF = { "##fileformat=VCFv4.2",
83		// fields other than AF are ignored when parsing as they have no INFO
84		// definition
85		"##INFO=<ID=AF,Number=A,Type=Float,Description=\"Allele Frequency, for each ALT allele, in the same order as listed\">",
86		"##INFO=<ID=AC_Female,Number=A,Type=Integer,Description=\"Allele count in Female genotypes\"",
87		"##INFO=<ID=AF_AFR,Number=A,Type=Float,Description=\"Allele Frequency among African/African American genotypes\"",
88		"##reference=Homo_sapiens/GRCh38",
89		"#CHROM\tPOS\tID\tREF\tALT\tQUAL\tFILTER\tINFO",
90		// A/T,C variants in position 2 of gene sequence (precedes transcript)
91		// should create 2 variant features with respective AF values
92		// malformed values for AC_Female and AF_AFR should be ignored
93		"17\t45051611\trs384765\tA\tT,C\t1666.64\tRF;XYZ\tAC=15;AF=5.0e-03,4.0e-03;AC_Female=12,3d;AF_AFR=low,2.3e-4",
94		// SNP G/C in position 4 of gene sequence, position 2 of transcript
95		// insertion G/GA is transferred to nucleotide but not to peptide
96		"17\t45051613\t.\tG\tGA,C\t1666.65\t.\tAC=15;AF=3.0e-03,2.0e-03",
97		// '.' in INFO field should be ignored
98		"17\t45051615\t.\tG\tC\t1666.66\tRF\tAC=16;AF=." };
99
		100% Uncovered Elements: 0 (5) Complexity: 1 Complexity Density: 0.2
100	1	@BeforeClass(alwaysRun = true)...
101		public void setUp()
102		{
103		/*
104		* configure to capture all available VCF and VEP (CSQ) fields
105		*/
106	1	Cache.loadProperties("test/jalview/io/testProps.jvprops");
107	1	Cache.setProperty("VCF_FIELDS", ".*");
108	1	Cache.setProperty("VEP_FIELDS", ".*");
109	1	Cache.setProperty("VCF_ASSEMBLY", "GRCh38=GRCh38");
110	1	Console.initLogger();
111		}
112
		100% Uncovered Elements: 0 (1) Complexity: 1 Complexity Density: 1
113	4	@BeforeTest(alwaysRun = true)...
114		public void setUpBeforeTest()
115		{
116		/*
117		* clear down feature attributes metadata
118		*/
119	4	FeatureAttributes.getInstance().clear();
120		}
121
		100% Uncovered Elements: 0 (84) Complexity: 2 Complexity Density: 0.02 1PASS
122	1	@Test(groups = "Functional")...
123		public void testDoLoad() throws IOException
124		{
125	1	AlignmentI al = buildAlignment();
126
127	1	File f = makeVcfFile();
128	1	VCFLoader loader = new VCFLoader(f.getPath());
129
130	1	loader.doLoad(al.getSequencesArray(), null);
131
132		/*
133		* verify variant feature(s) added to gene
134		* NB alleles at a locus may not be processed, and features added,
135		* in the order in which they appear in the VCF record as method
136		* VariantContext.getAlternateAlleles() does not guarantee order
137		* - order of assertions here matches what we find (is not important)
138		*/
139	1	List<SequenceFeature> geneFeatures = al.getSequenceAt(0)
140		.getSequenceFeatures();
141	1	SequenceFeatures.sortFeatures(geneFeatures, true);
142	1	assertEquals(geneFeatures.size(), 5);
143	1	SequenceFeature sf = geneFeatures.get(0);
144	1	assertEquals(sf.getFeatureGroup(), "VCF");
145	1	assertEquals(sf.getBegin(), 2);
146	1	assertEquals(sf.getEnd(), 2);
147	1	assertEquals(sf.getScore(), 0f);
148	1	assertEquals(sf.getValue("AF"), "4.0e-03");
149	1	assertEquals(sf.getValue("AF_AFR"), "2.3e-4");
150	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "A,C");
151	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
152	1	assertEquals(sf.getValue("POS"), "45051611");
153	1	assertEquals(sf.getValue("ID"), "rs384765");
154	1	assertEquals(sf.getValue("QUAL"), "1666.64");
155	1	assertEquals(sf.getValue("FILTER"), "RF;XYZ");
156		// malformed integer for AC_Female is ignored (JAL-3375)
157	1	assertNull(sf.getValue("AC_Female"));
158
159	1	sf = geneFeatures.get(1);
160	1	assertEquals(sf.getFeatureGroup(), "VCF");
161	1	assertEquals(sf.getBegin(), 2);
162	1	assertEquals(sf.getEnd(), 2);
163	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
164	1	assertEquals(sf.getScore(), 0f);
165	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 5.0e-03,
166		DELTA);
167	1	assertEquals(sf.getValue("AC_Female"), "12");
168		// malformed float for AF_AFR is ignored (JAL-3375)
169	1	assertNull(sf.getValue("AC_AFR"));
170	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "A,T");
171
172	1	sf = geneFeatures.get(2);
173	1	assertEquals(sf.getFeatureGroup(), "VCF");
174	1	assertEquals(sf.getBegin(), 4);
175	1	assertEquals(sf.getEnd(), 4);
176	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
177	1	assertEquals(sf.getScore(), 0f);
178	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 2.0e-03,
179		DELTA);
180	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,C");
181
182	1	sf = geneFeatures.get(3);
183	1	assertEquals(sf.getFeatureGroup(), "VCF");
184	1	assertEquals(sf.getBegin(), 4);
185	1	assertEquals(sf.getEnd(), 4);
186	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
187	1	assertEquals(sf.getScore(), 0f);
188	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 3.0e-03,
189		DELTA);
190	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,GA");
191	1	assertNull(sf.getValue("ID")); // '.' is ignored
192	1	assertNull(sf.getValue("FILTER")); // '.' is ignored
193
194	1	sf = geneFeatures.get(4);
195	1	assertEquals(sf.getFeatureGroup(), "VCF");
196	1	assertEquals(sf.getBegin(), 6);
197	1	assertEquals(sf.getEnd(), 6);
198	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
199	1	assertEquals(sf.getScore(), 0f);
200		// AF=. should not have been captured
201	1	assertNull(sf.getValue("AF"));
202	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,C");
203
204		/*
205		* verify variant feature(s) added to transcript
206		*/
207	1	List<SequenceFeature> transcriptFeatures = al.getSequenceAt(1)
208		.getSequenceFeatures();
209	1	assertEquals(transcriptFeatures.size(), 3);
210	1	sf = transcriptFeatures.get(0);
211	1	assertEquals(sf.getFeatureGroup(), "VCF");
212	1	assertEquals(sf.getBegin(), 2);
213	1	assertEquals(sf.getEnd(), 2);
214	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
215	1	assertEquals(sf.getScore(), 0f);
216	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 2.0e-03,
217		DELTA);
218	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,C");
219	1	sf = transcriptFeatures.get(1);
220	1	assertEquals(sf.getFeatureGroup(), "VCF");
221	1	assertEquals(sf.getBegin(), 2);
222	1	assertEquals(sf.getEnd(), 2);
223	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
224	1	assertEquals(sf.getScore(), 0f);
225	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 3.0e-03,
226		DELTA);
227	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,GA");
228
229		/*
230		* verify SNP variant feature(s) computed and added to protein
231		* first codon AGC varies to ACC giving S/T
232		*/
233	1	List<DBRefEntry> dbRefs = al.getSequenceAt(1).getDBRefs();
234	1	SequenceI peptide = null;
235	1	for (DBRefEntry dbref : dbRefs)
236		{
237	2	if (dbref.getMap().getMap().getFromRatio() == 3)
238		{
239	1	peptide = dbref.getMap().getTo();
240		}
241		}
242	1	List<SequenceFeature> proteinFeatures = peptide.getSequenceFeatures();
243
244		/*
245		* JAL-3187 don't precompute protein features, do dynamically instead
246		*/
247	1	assertTrue(proteinFeatures.isEmpty());
248		}
249
		100% Uncovered Elements: 0 (7) Complexity: 1 Complexity Density: 0.14
250	2	private File makeVcfFile() throws IOException...
251		{
252	2	File f = File.createTempFile("Test", ".vcf");
253	2	f.deleteOnExit();
254	2	PrintWriter pw = new PrintWriter(f);
255	2	for (String vcfLine : VCF)
256		{
257	18	pw.println(vcfLine);
258		}
259	2	pw.close();
260	2	return f;
261		}
262
263		/**
264		* Make a simple alignment with one 'gene' and one 'transcript'
265		*
266		* @return
267		*/
		100% Uncovered Elements: 0 (46) Complexity: 1 Complexity Density: 0.02
268	3	private AlignmentI buildAlignment()...
269		{
270	3	AlignFrame af = new FileLoader().LoadFileWaitTillLoaded(FASTA,
271		DataSourceType.PASTE);
272
273		/*
274		* map gene1 sequence to chromosome (normally done when the sequence is fetched
275		* from Ensembl and transcripts computed)
276		*/
277	3	AlignmentI alignment = af.getViewport().getAlignment();
278	3	SequenceI gene1 = alignment.findName("gene1");
279	3	int[] to = new int[] { 45051610, 45051634 };
280	3	int[] from = new int[] { gene1.getStart(), gene1.getEnd() };
281	3	gene1.setGeneLoci("homo_sapiens", "GRCh38", "17",
282		new MapList(from, to, 1, 1));
283
284		/*
285		* map 'transcript1' to chromosome via 'gene1'
286		* transcript1/1-18 is gene1/3-10,15-24
287		* which is chromosome 45051612-45051619,45051624-45051633
288		*/
289	3	to = new int[] { 45051612, 45051619, 45051624, 45051633 };
290	3	SequenceI transcript1 = alignment.findName("transcript1");
291	3	from = new int[] { transcript1.getStart(), transcript1.getEnd() };
292	3	transcript1.setGeneLoci("homo_sapiens", "GRCh38", "17",
293		new MapList(from, to, 1, 1));
294
295		/*
296		* map gene2 to chromosome reverse strand
297		*/
298	3	SequenceI gene2 = alignment.findName("gene2");
299	3	to = new int[] { 45051634, 45051610 };
300	3	from = new int[] { gene2.getStart(), gene2.getEnd() };
301	3	gene2.setGeneLoci("homo_sapiens", "GRCh38", "17",
302		new MapList(from, to, 1, 1));
303
304		/*
305		* map 'transcript2' to chromosome via 'gene2'
306		* transcript2/1-18 is gene2/2-11,16-23
307		* which is chromosome 45051633-45051624,45051619-45051612
308		*/
309	3	to = new int[] { 45051633, 45051624, 45051619, 45051612 };
310	3	SequenceI transcript2 = alignment.findName("transcript2");
311	3	from = new int[] { transcript2.getStart(), transcript2.getEnd() };
312	3	transcript2.setGeneLoci("homo_sapiens", "GRCh38", "17",
313		new MapList(from, to, 1, 1));
314
315		/*
316		* add a protein product as a DBRef on transcript1
317		*/
318	3	SequenceI peptide1 = new Sequence("ENSP001", "SWRECD");
319	3	MapList mapList = new MapList(new int[] { 1, 18 }, new int[] { 1, 6 },
320		3, 1);
321	3	Mapping map = new Mapping(peptide1, mapList);
322	3	DBRefEntry product = new DBRefEntry("", "", "ENSP001", map);
323	3	transcript1.addDBRef(product);
324
325		/*
326		* add a protein product as a DBRef on transcript2
327		*/
328	3	SequenceI peptide2 = new Sequence("ENSP002", "VTLSPA");
329	3	mapList = new MapList(new int[] { 1, 18 }, new int[] { 1, 6 }, 3, 1);
330	3	map = new Mapping(peptide2, mapList);
331	3	product = new DBRefEntry("", "", "ENSP002", map);
332	3	transcript2.addDBRef(product);
333
334		/*
335		* map gene3 to chromosome
336		*/
337	3	SequenceI gene3 = alignment.findName("gene3");
338	3	to = new int[] { 45051610, 45051634 };
339	3	from = new int[] { gene3.getStart(), gene3.getEnd() };
340	3	gene3.setGeneLoci("homo_sapiens", "GRCh38", "5",
341		new MapList(from, to, 1, 1));
342
343		/*
344		* map 'transcript3' to chromosome
345		*/
346	3	SequenceI transcript3 = alignment.findName("transcript3");
347	3	to = new int[] { 45051612, 45051619, 45051624, 45051633 };
348	3	from = new int[] { transcript3.getStart(), transcript3.getEnd() };
349	3	transcript3.setGeneLoci("homo_sapiens", "GRCh38", "5",
350		new MapList(from, to, 1, 1));
351
352		/*
353		* map 'transcript4' to chromosome
354		*/
355	3	SequenceI transcript4 = alignment.findName("transcript4");
356	3	to = new int[] { 45051615, 45051617, 45051619, 45051632, 45051634,
357		45051634 };
358	3	from = new int[] { transcript4.getStart(), transcript4.getEnd() };
359	3	transcript4.setGeneLoci("homo_sapiens", "GRCh38", "5",
360		new MapList(from, to, 1, 1));
361
362		/*
363		* add a protein product as a DBRef on transcript3
364		*/
365	3	SequenceI peptide3 = new Sequence("ENSP003", "SWRECD");
366	3	mapList = new MapList(new int[] { 1, 18 }, new int[] { 1, 6 }, 3, 1);
367	3	map = new Mapping(peptide3, mapList);
368	3	product = new DBRefEntry("", "", "ENSP003", map);
369	3	transcript3.addDBRef(product);
370
371	3	return alignment;
372		}
373
374		/**
375		* Test with 'gene' and 'transcript' mapped to the reverse strand of the
376		* chromosome. The VCF variant positions (in forward coordinates) should get
377		* correctly located on sequence positions.
378		*
379		* @throws IOException
380		*/
		100% Uncovered Elements: 0 (67) Complexity: 2 Complexity Density: 0.03 1PASS
381	1	@Test(groups = "Functional")...
382		public void testDoLoad_reverseStrand() throws IOException
383		{
384	1	AlignmentI al = buildAlignment();
385
386	1	File f = makeVcfFile();
387
388	1	VCFLoader loader = new VCFLoader(f.getPath());
389
390	1	loader.doLoad(al.getSequencesArray(), null);
391
392		/*
393		* verify variant feature(s) added to gene2
394		* gene2/1-25 maps to chromosome 45051634- reverse strand
395		*/
396	1	List<SequenceFeature> geneFeatures = al.getSequenceAt(2)
397		.getSequenceFeatures();
398	1	SequenceFeatures.sortFeatures(geneFeatures, true);
399	1	assertEquals(geneFeatures.size(), 5);
400	1	SequenceFeature sf;
401
402		/*
403		* insertion G/GA at 45051613 maps to an insertion at
404		* the preceding position (21) on reverse strand gene
405		* reference: CAAGC -> GCTTG/21-25
406		* genomic variant: CAAGAC (G/GA)
407		* gene variant: GTCTTG (G/GT at 21)
408		*/
409	1	sf = geneFeatures.get(1);
410	1	assertEquals(sf.getFeatureGroup(), "VCF");
411	1	assertEquals(sf.getBegin(), 21);
412	1	assertEquals(sf.getEnd(), 21);
413	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
414	1	assertEquals(sf.getScore(), 0f);
415	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,GT");
416	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 3.0e-03,
417		DELTA);
418
419		/*
420		* variant G/C at 45051613 maps to C/G at gene position 22
421		*/
422	1	sf = geneFeatures.get(2);
423	1	assertEquals(sf.getFeatureGroup(), "VCF");
424	1	assertEquals(sf.getBegin(), 22);
425	1	assertEquals(sf.getEnd(), 22);
426	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
427	1	assertEquals(sf.getScore(), 0f);
428	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "C,G");
429	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 2.0e-03,
430		DELTA);
431
432		/*
433		* variant A/C at 45051611 maps to T/G at gene position 24
434		*/
435	1	sf = geneFeatures.get(3);
436	1	assertEquals(sf.getFeatureGroup(), "VCF");
437	1	assertEquals(sf.getBegin(), 24);
438	1	assertEquals(sf.getEnd(), 24);
439	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
440	1	assertEquals(sf.getScore(), 0f);
441	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "T,G");
442	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 4.0e-03,
443		DELTA);
444
445		/*
446		* variant A/T at 45051611 maps to T/A at gene position 24
447		*/
448	1	sf = geneFeatures.get(4);
449	1	assertEquals(sf.getFeatureGroup(), "VCF");
450	1	assertEquals(sf.getBegin(), 24);
451	1	assertEquals(sf.getEnd(), 24);
452	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
453	1	assertEquals(sf.getScore(), 0f);
454	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "T,A");
455	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 5.0e-03,
456		DELTA);
457
458		/*
459		* verify 3 variant features added to transcript2
460		*/
461	1	List<SequenceFeature> transcriptFeatures = al.getSequenceAt(3)
462		.getSequenceFeatures();
463	1	assertEquals(transcriptFeatures.size(), 3);
464
465		/*
466		* insertion G/GT at position 21 of gene maps to position 16 of transcript
467		*/
468	1	sf = transcriptFeatures.get(1);
469	1	assertEquals(sf.getFeatureGroup(), "VCF");
470	1	assertEquals(sf.getBegin(), 16);
471	1	assertEquals(sf.getEnd(), 16);
472	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
473	1	assertEquals(sf.getScore(), 0f);
474	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,GT");
475	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 3.0e-03,
476		DELTA);
477
478		/*
479		* SNP C/G at position 22 of gene maps to position 17 of transcript
480		*/
481	1	sf = transcriptFeatures.get(2);
482	1	assertEquals(sf.getFeatureGroup(), "VCF");
483	1	assertEquals(sf.getBegin(), 17);
484	1	assertEquals(sf.getEnd(), 17);
485	1	assertEquals(sf.getType(), SEQUENCE_VARIANT);
486	1	assertEquals(sf.getScore(), 0f);
487	1	assertEquals(sf.getValue(Gff3Helper.ALLELES), "C,G");
488	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 2.0e-03,
489		DELTA);
490
491		/*
492		* verify variant feature(s) computed and added to protein
493		* last codon GCT varies to GGT giving A/G in the last peptide position
494		*/
495	1	List<DBRefEntry> dbRefs = al.getSequenceAt(3).getDBRefs();
496	1	SequenceI peptide = null;
497	1	for (DBRefEntry dbref : dbRefs)
498		{
499	2	if (dbref.getMap().getMap().getFromRatio() == 3)
500		{
501	1	peptide = dbref.getMap().getTo();
502		}
503		}
504	1	List<SequenceFeature> proteinFeatures = peptide.getSequenceFeatures();
505
506		/*
507		* JAL-3187 don't precompute protein features, do dynamically instead
508		*/
509	1	assertTrue(proteinFeatures.isEmpty());
510		}
511
512		/**
513		* Tests that if VEP consequence (CSQ) data is present in the VCF data, then
514		* it is added to the variant feature, but restricted where possible to the
515		* consequences for a specific transcript
516		*
517		* @throws IOException
518		*/
		100% Uncovered Elements: 0 (136) Complexity: 2 Complexity Density: 0.01 1PASS
519	1	@Test(groups = "Functional")...
520		public void testDoLoad_vepCsq() throws IOException
521		{
522	1	AlignmentI al = buildAlignment();
523
524	1	VCFLoader loader = new VCFLoader("test/jalview/io/vcf/testVcf.vcf");
525
526		/*
527		* VCF data file with variants at gene3 positions
528		* 1 C/A
529		* 5 C/T
530		* 9 CGT/C (deletion)
531		* 13 C/G, C/T
532		* 17 A/AC (insertion), A/G
533		*/
534	1	loader.doLoad(al.getSequencesArray(), null);
535
536		/*
537		* verify variant feature(s) added to gene3
538		*/
539	1	List<SequenceFeature> geneFeatures = al.findName("gene3")
540		.getSequenceFeatures();
541	1	SequenceFeatures.sortFeatures(geneFeatures, true);
542	1	assertEquals(geneFeatures.size(), 7);
543	1	SequenceFeature sf = geneFeatures.get(0);
544	1	assertEquals(sf.getBegin(), 1);
545	1	assertEquals(sf.getEnd(), 1);
546	1	assertEquals(sf.getScore(), 0f);
547	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.1f, DELTA);
548	1	assertEquals(sf.getValue("alleles"), "C,A");
549		// gene features include Consequence for all transcripts
550	1	Map map = (Map) sf.getValue("CSQ");
551	1	assertEquals(map.size(), 9);
552	1	assertEquals(map.get("PolyPhen"), "Bad");
553
554	1	sf = geneFeatures.get(1);
555	1	assertEquals(sf.getBegin(), 5);
556	1	assertEquals(sf.getEnd(), 5);
557	1	assertEquals(sf.getScore(), 0f);
558	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.2f, DELTA);
559	1	assertEquals(sf.getValue("alleles"), "C,T");
560	1	map = (Map) sf.getValue("CSQ");
561	1	assertEquals(map.size(), 9);
562	1	assertEquals(map.get("PolyPhen"), "Bad;;"); // %3B%3B decoded
563
564	1	sf = geneFeatures.get(2);
565	1	assertEquals(sf.getBegin(), 9);
566	1	assertEquals(sf.getEnd(), 11); // deletion over 3 positions
567	1	assertEquals(sf.getScore(), 0f);
568	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.3f, DELTA);
569	1	assertEquals(sf.getValue("alleles"), "CGG,C");
570	1	map = (Map) sf.getValue("CSQ");
571	1	assertEquals(map.size(), 9);
572
573	1	sf = geneFeatures.get(3);
574	1	assertEquals(sf.getBegin(), 13);
575	1	assertEquals(sf.getEnd(), 13);
576	1	assertEquals(sf.getScore(), 0f);
577	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.5f, DELTA);
578	1	assertEquals(sf.getValue("alleles"), "C,T");
579	1	map = (Map) sf.getValue("CSQ");
580	1	assertEquals(map.size(), 9);
581
582	1	sf = geneFeatures.get(4);
583	1	assertEquals(sf.getBegin(), 13);
584	1	assertEquals(sf.getEnd(), 13);
585	1	assertEquals(sf.getScore(), 0f);
586	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.4f, DELTA);
587	1	assertEquals(sf.getValue("alleles"), "C,G");
588	1	map = (Map) sf.getValue("CSQ");
589	1	assertEquals(map.size(), 9);
590
591	1	sf = geneFeatures.get(5);
592	1	assertEquals(sf.getBegin(), 17);
593	1	assertEquals(sf.getEnd(), 17);
594	1	assertEquals(sf.getScore(), 0f);
595	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.7f, DELTA);
596	1	assertEquals(sf.getValue("alleles"), "A,G");
597	1	map = (Map) sf.getValue("CSQ");
598	1	assertEquals(map.size(), 9);
599
600	1	sf = geneFeatures.get(6);
601	1	assertEquals(sf.getBegin(), 17);
602	1	assertEquals(sf.getEnd(), 17); // insertion
603	1	assertEquals(sf.getScore(), 0f);
604	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.6f, DELTA);
605	1	assertEquals(sf.getValue("alleles"), "A,AC");
606	1	map = (Map) sf.getValue("CSQ");
607	1	assertEquals(map.size(), 9);
608
609		/*
610		* verify variant feature(s) added to transcript3
611		* at columns 5 (1), 17 (2), positions 3, 11
612		* note the deletion at columns 9-11 is not transferred since col 11
613		* has no mapping to transcript 3
614		*/
615	1	List<SequenceFeature> transcriptFeatures = al.findName("transcript3")
616		.getSequenceFeatures();
617	1	SequenceFeatures.sortFeatures(transcriptFeatures, true);
618	1	assertEquals(transcriptFeatures.size(), 3);
619	1	sf = transcriptFeatures.get(0);
620	1	assertEquals(sf.getBegin(), 3);
621	1	assertEquals(sf.getEnd(), 3);
622	1	assertEquals(sf.getScore(), 0f);
623	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.2f, DELTA);
624	1	assertEquals(sf.getValue("alleles"), "C,T");
625		// transcript features only have Consequence for that transcripts
626	1	map = (Map) sf.getValue("CSQ");
627	1	assertEquals(map.size(), 9);
628	1	assertEquals(sf.getValueAsString("CSQ", "Feature"), "transcript3");
629
630	1	sf = transcriptFeatures.get(1);
631	1	assertEquals(sf.getBegin(), 11);
632	1	assertEquals(sf.getEnd(), 11);
633	1	assertEquals(sf.getScore(), 0f);
634	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.7f, DELTA);
635	1	assertEquals(sf.getValue("alleles"), "A,G");
636	1	assertEquals(map.size(), 9);
637	1	assertEquals(sf.getValueAsString("CSQ", "Feature"), "transcript3");
638
639	1	sf = transcriptFeatures.get(2);
640	1	assertEquals(sf.getBegin(), 11);
641	1	assertEquals(sf.getEnd(), 11);
642	1	assertEquals(sf.getScore(), 0f);
643	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.6f, DELTA);
644	1	assertEquals(sf.getValue("alleles"), "A,AC");
645	1	assertEquals(map.size(), 9);
646	1	assertEquals(sf.getValueAsString("CSQ", "Feature"), "transcript3");
647
648		/*
649		* verify variants computed on protein product for transcript3
650		* peptide is SWRECD
651		* codon variants are AGC/AGT position 1 which is synonymous
652		* and GAG/GGG which is E/G in position 4
653		* the insertion variant is not transferred to the peptide
654		*/
655	1	List<DBRefEntry> dbRefs = al.findName("transcript3").getDBRefs();
656	1	SequenceI peptide = null;
657	1	for (DBRefEntry dbref : dbRefs)
658		{
659	2	if (dbref.getMap().getMap().getFromRatio() == 3)
660		{
661	1	peptide = dbref.getMap().getTo();
662		}
663		}
664	1	List<SequenceFeature> proteinFeatures = peptide.getSequenceFeatures();
665		/*
666		* JAL-3187 don't precompute protein features, do dynamically instead
667		*/
668	1	assertTrue(proteinFeatures.isEmpty());
669		// SequenceFeatures.sortFeatures(proteinFeatures, true);
670		// assertEquals(proteinFeatures.size(), 2);
671		// sf = proteinFeatures.get(0);
672		// assertEquals(sf.getFeatureGroup(), "VCF");
673		// assertEquals(sf.getBegin(), 1);
674		// assertEquals(sf.getEnd(), 1);
675		// assertEquals(sf.getType(), SequenceOntologyI.SYNONYMOUS_VARIANT);
676		// assertEquals(sf.getDescription(), "agC/agT");
677		// sf = proteinFeatures.get(1);
678		// assertEquals(sf.getFeatureGroup(), "VCF");
679		// assertEquals(sf.getBegin(), 4);
680		// assertEquals(sf.getEnd(), 4);
681		// assertEquals(sf.getType(), SequenceOntologyI.NONSYNONYMOUS_VARIANT);
682		// assertEquals(sf.getDescription(), "p.Glu4Gly");
683
684		/*
685		* verify variant feature(s) added to transcript4
686		* at columns 13 (2) and 17 (2), positions 7 and 11
687		*/
688	1	transcriptFeatures = al.findName("transcript4").getSequenceFeatures();
689	1	SequenceFeatures.sortFeatures(transcriptFeatures, true);
690	1	assertEquals(transcriptFeatures.size(), 4);
691	1	sf = transcriptFeatures.get(0);
692	1	assertEquals(sf.getBegin(), 7);
693	1	assertEquals(sf.getEnd(), 7);
694	1	assertEquals(sf.getScore(), 0f);
695	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.5f, DELTA);
696	1	assertEquals(sf.getValue("alleles"), "C,T");
697	1	assertEquals(map.size(), 9);
698	1	assertEquals(sf.getValueAsString("CSQ", "Feature"), "transcript4");
699
700	1	sf = transcriptFeatures.get(1);
701	1	assertEquals(sf.getBegin(), 7);
702	1	assertEquals(sf.getEnd(), 7);
703	1	assertEquals(sf.getScore(), 0f);
704	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.4f, DELTA);
705	1	assertEquals(sf.getValue("alleles"), "C,G");
706	1	assertEquals(map.size(), 9);
707	1	assertEquals(sf.getValueAsString("CSQ", "Feature"), "transcript4");
708
709	1	sf = transcriptFeatures.get(2);
710	1	assertEquals(sf.getBegin(), 11);
711	1	assertEquals(sf.getEnd(), 11);
712	1	assertEquals(sf.getScore(), 0f);
713	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.7f, DELTA);
714	1	assertEquals(sf.getValue("alleles"), "A,G");
715	1	assertEquals(map.size(), 9);
716	1	assertEquals(sf.getValueAsString("CSQ", "Feature"), "transcript4");
717
718	1	sf = transcriptFeatures.get(3);
719	1	assertEquals(sf.getBegin(), 11);
720	1	assertEquals(sf.getEnd(), 11);
721	1	assertEquals(sf.getScore(), 0f);
722	1	assertEquals(Float.parseFloat((String) sf.getValue("AF")), 0.6f, DELTA);
723	1	assertEquals(sf.getValue("alleles"), "A,AC");
724	1	assertEquals(map.size(), 9);
725	1	assertEquals(sf.getValueAsString("CSQ", "Feature"), "transcript4");
726		}
727
728		/**
729		* A test that demonstrates loading a contig sequence from an indexed sequence
730		* database which is the reference for a VCF file
731		*
732		* @throws IOException
733		*/
		100% Uncovered Elements: 0 (39) Complexity: 1 Complexity Density: 0.03 1PASS
734	1	@Test(groups = "Functional")...
735		public void testLoadVCFContig() throws IOException
736		{
737	1	VCFLoader loader = new VCFLoader("test/jalview/io/vcf/testVcf2.vcf");
738
739	1	SequenceI seq = loader.loadVCFContig("contig123");
740	1	assertEquals(seq.getLength(), 15);
741	1	assertEquals(seq.getSequenceAsString(), "AAAAACCCCCGGGGG");
742	1	List<SequenceFeature> features = seq.getSequenceFeatures();
743	1	SequenceFeatures.sortFeatures(features, true);
744	1	assertEquals(features.size(), 2);
745	1	SequenceFeature sf = features.get(0);
746	1	assertEquals(sf.getBegin(), 8);
747	1	assertEquals(sf.getEnd(), 8);
748	1	assertEquals(sf.getDescription(), "C,A");
749	1	sf = features.get(1);
750	1	assertEquals(sf.getBegin(), 12);
751	1	assertEquals(sf.getEnd(), 12);
752	1	assertEquals(sf.getDescription(), "G,T");
753
754	1	seq = loader.loadVCFContig("contig789");
755	1	assertEquals(seq.getLength(), 25);
756	1	assertEquals(seq.getSequenceAsString(), "GGGGGTTTTTAAAAACCCCCGGGGG");
757	1	features = seq.getSequenceFeatures();
758	1	SequenceFeatures.sortFeatures(features, true);
759	1	assertEquals(features.size(), 2);
760	1	sf = features.get(0);
761	1	assertEquals(sf.getBegin(), 2);
762	1	assertEquals(sf.getEnd(), 2);
763	1	assertEquals(sf.getDescription(), "G,T");
764	1	sf = features.get(1);
765	1	assertEquals(sf.getBegin(), 21);
766	1	assertEquals(sf.getEnd(), 21);
767	1	assertEquals(sf.getDescription(), "G,A");
768
769	1	seq = loader.loadVCFContig("contig456");
770	1	assertEquals(seq.getLength(), 20);
771	1	assertEquals(seq.getSequenceAsString(), "CCCCCGGGGGTTTTTAAAAA");
772	1	features = seq.getSequenceFeatures();
773	1	SequenceFeatures.sortFeatures(features, true);
774	1	assertEquals(features.size(), 1);
775	1	sf = features.get(0);
776	1	assertEquals(sf.getBegin(), 15);
777	1	assertEquals(sf.getEnd(), 15);
778	1	assertEquals(sf.getDescription(), "T,C");
779		}
780		}