JAL-3375 ignore '.' values for VCF data

diff --git a/src/jalview/datamodel/features/FeatureAttributes.java b/src/jalview/datamodel/features/FeatureAttributes.java
index 10249f3..26e6f8b 100644 (file)
--- a/src/jalview/datamodel/features/FeatureAttributes.java
+++ b/src/jalview/datamodel/features/FeatureAttributes.java
@@ -371,4 +371,27 @@ public class FeatureAttributes
     }
     return null;
   }
+
+  /**
+   * Resets all attribute metadata
+   */
+  public void clear()
+  {
+    attributes.clear();
+  }
+
+  /**
+   * Resets attribute metadata for one feature type
+   * 
+   * @param featureType
+   */
+  public void clear(String featureType)
+  {
+    Map<String[], AttributeData> map = attributes.get(featureType);
+    if (map != null)
+    {
+      map.clear();
+    }
+
+  }
 }

diff --git a/src/jalview/io/vcf/VCFLoader.java b/src/jalview/io/vcf/VCFLoader.java
index 7bf7791..053b52f 100644 (file)
--- a/src/jalview/io/vcf/VCFLoader.java
+++ b/src/jalview/io/vcf/VCFLoader.java
@@ -51,6 +51,8 @@ import htsjdk.variant.vcf.VCFInfoHeaderLine;
  */
 public class VCFLoader
 {
+  private static final String NO_VALUE = ".";
+
   private static final String DEFAULT_SPECIES = "homo_sapiens";
 
   /**
@@ -877,7 +879,7 @@ public class VCFLoader
 
     if (att instanceof String)
     {
-      return (String) att;
+      return NO_VALUE.equals(att) ? null : (String) att;
     }
     else if (att instanceof ArrayList)
     {

diff --git a/test/jalview/io/vcf/VCFLoaderTest.java b/test/jalview/io/vcf/VCFLoaderTest.java
index a87c160..6ba6fbe 100644 (file)
--- a/test/jalview/io/vcf/VCFLoaderTest.java
+++ b/test/jalview/io/vcf/VCFLoaderTest.java
@@ -1,6 +1,9 @@
 package jalview.io.vcf;
 
+import static jalview.io.gff.SequenceOntologyI.SEQUENCE_VARIANT;
 import static org.testng.Assert.assertEquals;
+import static org.testng.Assert.assertNull;
+import static org.testng.Assert.assertSame;
 
 import jalview.bin.Cache;
 import jalview.datamodel.AlignmentI;
@@ -9,6 +12,8 @@ import jalview.datamodel.Mapping;
 import jalview.datamodel.Sequence;
 import jalview.datamodel.SequenceFeature;
 import jalview.datamodel.SequenceI;
+import jalview.datamodel.features.FeatureAttributes;
+import jalview.datamodel.features.FeatureAttributes.Datatype;
 import jalview.datamodel.features.SequenceFeatures;
 import jalview.gui.AlignFrame;
 import jalview.io.DataSourceType;
@@ -24,6 +29,7 @@ import java.util.List;
 import java.util.Map;
 
 import org.testng.annotations.BeforeClass;
+import org.testng.annotations.BeforeTest;
 import org.testng.annotations.Test;
 
 public class VCFLoaderTest
@@ -56,15 +62,18 @@ public class VCFLoaderTest
           + ">transcript4/1-18\n-----TGG-GGACGAGAGTGTGA-A\n";
 
   private static final String[] VCF = { "##fileformat=VCFv4.2",
+      // fields other than AF are ignored when parsing as they have no INFO definition
       "##INFO=<ID=AF,Number=A,Type=Float,Description=\"Allele Frequency, for each ALT allele, in the same order as listed\">",
       "##reference=Homo_sapiens/GRCh38",
       "#CHROM\tPOS\tID\tREF\tALT\tQUAL\tFILTER\tINFO",
       // A/T,C variants in position 2 of gene sequence (precedes transcript)
-      // should create 2 variant features with respective scores
+      // should create 2 variant features with respective AF values
       "17\t45051611\t.\tA\tT,C\t1666.64\tRF\tAC=15;AF=5.0e-03,4.0e-03",
       // SNP G/C in position 4 of gene sequence, position 2 of transcript
       // insertion G/GA is transferred to nucleotide but not to peptide
-      "17\t45051613\t.\tG\tGA,C\t1666.64\tRF\tAC=15;AF=3.0e-03,2.0e-03" };
+      "17\t45051613\t.\tG\tGA,C\t1666.65\tRF\tAC=15;AF=3.0e-03,2.0e-03",
+      // '.' in INFO field should be ignored
+      "17\t45051615\t.\tG\tC\t1666.66\tRF\tAC=16;AF=." };
 
   @BeforeClass(alwaysRun = true)
   public void setUp()
@@ -79,12 +88,21 @@ public class VCFLoaderTest
     Cache.initLogger();
   }
 
+  @BeforeTest(alwaysRun = true)
+  public void setUpBeforeTest()
+  {
+    /*
+     * clear down feature attributes metadata
+     */
+    FeatureAttributes.getInstance().clear();
+  }
+
   @Test(groups = "Functional")
   public void testDoLoad() throws IOException
   {
     AlignmentI al = buildAlignment();
 
-    File f = makeVcf();
+    File f = makeVcfFile();
     VCFLoader loader = new VCFLoader(f.getPath());
 
     loader.doLoad(al.getSequencesArray(), null);
@@ -99,7 +117,7 @@ public class VCFLoaderTest
     List<SequenceFeature> geneFeatures = al.getSequenceAt(0)
             .getSequenceFeatures();
     SequenceFeatures.sortFeatures(geneFeatures, true);
-    assertEquals(geneFeatures.size(), 4);
+    assertEquals(geneFeatures.size(), 5);
     SequenceFeature sf = geneFeatures.get(0);
     assertEquals(sf.getFeatureGroup(), "VCF");
     assertEquals(sf.getBegin(), 2);
@@ -108,12 +126,12 @@ public class VCFLoaderTest
     assertEquals(Float.parseFloat((String) sf.getValue("AF")), 4.0e-03,
             DELTA);
     assertEquals(sf.getValue(Gff3Helper.ALLELES), "A,C");
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     sf = geneFeatures.get(1);
     assertEquals(sf.getFeatureGroup(), "VCF");
     assertEquals(sf.getBegin(), 2);
     assertEquals(sf.getEnd(), 2);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
     assertEquals(Float.parseFloat((String) sf.getValue("AF")), 5.0e-03,
             DELTA);
@@ -123,7 +141,7 @@ public class VCFLoaderTest
     assertEquals(sf.getFeatureGroup(), "VCF");
     assertEquals(sf.getBegin(), 4);
     assertEquals(sf.getEnd(), 4);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
     assertEquals(Float.parseFloat((String) sf.getValue("AF")), 2.0e-03,
             DELTA);
@@ -133,23 +151,33 @@ public class VCFLoaderTest
     assertEquals(sf.getFeatureGroup(), "VCF");
     assertEquals(sf.getBegin(), 4);
     assertEquals(sf.getEnd(), 4);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
     assertEquals(Float.parseFloat((String) sf.getValue("AF")), 3.0e-03,
             DELTA);
     assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,GA");
 
+    sf = geneFeatures.get(4);
+    assertEquals(sf.getFeatureGroup(), "VCF");
+    assertEquals(sf.getBegin(), 6);
+    assertEquals(sf.getEnd(), 6);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
+    assertEquals(sf.getScore(), 0f);
+    // AF=. should not have been captured
+    assertNull(sf.getValue("AF"));
+    assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,C");
+
     /*
      * verify variant feature(s) added to transcript
      */
     List<SequenceFeature> transcriptFeatures = al.getSequenceAt(1)
             .getSequenceFeatures();
-    assertEquals(transcriptFeatures.size(), 2);
+    assertEquals(transcriptFeatures.size(), 3);
     sf = transcriptFeatures.get(0);
     assertEquals(sf.getFeatureGroup(), "VCF");
     assertEquals(sf.getBegin(), 2);
     assertEquals(sf.getEnd(), 2);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
     assertEquals(Float.parseFloat((String) sf.getValue("AF")), 2.0e-03,
             DELTA);
@@ -158,7 +186,7 @@ public class VCFLoaderTest
     assertEquals(sf.getFeatureGroup(), "VCF");
     assertEquals(sf.getBegin(), 2);
     assertEquals(sf.getEnd(), 2);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
     assertEquals(Float.parseFloat((String) sf.getValue("AF")), 3.0e-03,
             DELTA);
@@ -178,16 +206,27 @@ public class VCFLoaderTest
       }
     }
     List<SequenceFeature> proteinFeatures = peptide.getSequenceFeatures();
-    assertEquals(proteinFeatures.size(), 1);
+    assertEquals(proteinFeatures.size(), 3);
     sf = proteinFeatures.get(0);
     assertEquals(sf.getFeatureGroup(), "VCF");
     assertEquals(sf.getBegin(), 1);
     assertEquals(sf.getEnd(), 1);
     assertEquals(sf.getType(), SequenceOntologyI.NONSYNONYMOUS_VARIANT);
     assertEquals(sf.getDescription(), "p.Ser1Thr");
+
+    /*
+     * check that sequence_variant attribute AF has been clocked as
+     * numeric with correct min and max values
+     * (i.e. invalid values have been ignored - JAL-3375)
+     */
+    FeatureAttributes fa = FeatureAttributes.getInstance();
+    assertSame(fa.getDatatype(SEQUENCE_VARIANT, "AF"), Datatype.Number);
+    float[] minmax = fa.getMinMax(SEQUENCE_VARIANT, "AF");
+    assertEquals(minmax[0], 0.002f);
+    assertEquals(minmax[1], 0.005f);
   }
 
-  private File makeVcf() throws IOException
+  private File makeVcfFile() throws IOException
   {
     File f = File.createTempFile("Test", ".vcf");
     f.deleteOnExit();
@@ -327,7 +366,7 @@ public class VCFLoaderTest
   {
     AlignmentI al = buildAlignment();
 
-    File f = makeVcf();
+    File f = makeVcfFile();
 
     VCFLoader loader = new VCFLoader(f.getPath());
 
@@ -340,96 +379,97 @@ public class VCFLoaderTest
     List<SequenceFeature> geneFeatures = al.getSequenceAt(2)
             .getSequenceFeatures();
     SequenceFeatures.sortFeatures(geneFeatures, true);
-    assertEquals(geneFeatures.size(), 4);
+    assertEquals(geneFeatures.size(), 5);
+    SequenceFeature sf;
 
     /*
-     * variant A/T at 45051611 maps to T/A at gene position 24
+     * insertion G/GA at 45051613 maps to an insertion at
+     * the preceding position (21) on reverse strand gene
+     * reference: CAAGC -> GCTTG/21-25
+     * genomic variant: CAAGAC (G/GA)
+     * gene variant: GTCTTG (G/GT at 21)
      */
-    SequenceFeature sf = geneFeatures.get(3);
+    sf = geneFeatures.get(1);
     assertEquals(sf.getFeatureGroup(), "VCF");
-    assertEquals(sf.getBegin(), 24);
-    assertEquals(sf.getEnd(), 24);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getBegin(), 21);
+    assertEquals(sf.getEnd(), 21);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
-    assertEquals(Float.parseFloat((String) sf.getValue("AF")), 5.0e-03,
+    assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,GT");
+    assertEquals(Float.parseFloat((String) sf.getValue("AF")), 3.0e-03,
             DELTA);
-    assertEquals(sf.getValue(Gff3Helper.ALLELES), "T,A");
 
     /*
-     * variant A/C at 45051611 maps to T/G at gene position 24
+     * variant G/C at 45051613 maps to C/G at gene position 22
      */
     sf = geneFeatures.get(2);
     assertEquals(sf.getFeatureGroup(), "VCF");
-    assertEquals(sf.getBegin(), 24);
-    assertEquals(sf.getEnd(), 24);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getBegin(), 22);
+    assertEquals(sf.getEnd(), 22);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
-    assertEquals(Float.parseFloat((String) sf.getValue("AF")), 4.0e-03,
+    assertEquals(sf.getValue(Gff3Helper.ALLELES), "C,G");
+    assertEquals(Float.parseFloat((String) sf.getValue("AF")), 2.0e-03,
             DELTA);
-    assertEquals(sf.getValue(Gff3Helper.ALLELES), "T,G");
 
     /*
-     * variant G/C at 45051613 maps to C/G at gene position 22
+     * variant A/C at 45051611 maps to T/G at gene position 24
      */
-    sf = geneFeatures.get(1);
+    sf = geneFeatures.get(3);
     assertEquals(sf.getFeatureGroup(), "VCF");
-    assertEquals(sf.getBegin(), 22);
-    assertEquals(sf.getEnd(), 22);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getBegin(), 24);
+    assertEquals(sf.getEnd(), 24);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
-    assertEquals(Float.parseFloat((String) sf.getValue("AF")), 2.0e-03,
+    assertEquals(sf.getValue(Gff3Helper.ALLELES), "T,G");
+    assertEquals(Float.parseFloat((String) sf.getValue("AF")), 4.0e-03,
             DELTA);
-    assertEquals(sf.getValue(Gff3Helper.ALLELES), "C,G");
 
     /*
-     * insertion G/GA at 45051613 maps to an insertion at
-     * the preceding position (21) on reverse strand gene
-     * reference: CAAGC -> GCTTG/21-25
-     * genomic variant: CAAGAC (G/GA)
-     * gene variant: GTCTTG (G/GT at 21)
+     * variant A/T at 45051611 maps to T/A at gene position 24
      */
-    sf = geneFeatures.get(0);
+    sf = geneFeatures.get(4);
     assertEquals(sf.getFeatureGroup(), "VCF");
-    assertEquals(sf.getBegin(), 21);
-    assertEquals(sf.getEnd(), 21);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getBegin(), 24);
+    assertEquals(sf.getEnd(), 24);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
-    assertEquals(Float.parseFloat((String) sf.getValue("AF")), 3.0e-03,
+    assertEquals(sf.getValue(Gff3Helper.ALLELES), "T,A");
+    assertEquals(Float.parseFloat((String) sf.getValue("AF")), 5.0e-03,
             DELTA);
-    assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,GT");
 
     /*
-     * verify 2 variant features added to transcript2
+     * verify 3 variant features added to transcript2
      */
     List<SequenceFeature> transcriptFeatures = al.getSequenceAt(3)
             .getSequenceFeatures();
-    assertEquals(transcriptFeatures.size(), 2);
+    assertEquals(transcriptFeatures.size(), 3);
 
     /*
      * insertion G/GT at position 21 of gene maps to position 16 of transcript
      */
-    sf = transcriptFeatures.get(0);
+    sf = transcriptFeatures.get(1);
     assertEquals(sf.getFeatureGroup(), "VCF");
     assertEquals(sf.getBegin(), 16);
     assertEquals(sf.getEnd(), 16);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
+    assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,GT");
     assertEquals(Float.parseFloat((String) sf.getValue("AF")), 3.0e-03,
             DELTA);
-    assertEquals(sf.getValue(Gff3Helper.ALLELES), "G,GT");
 
     /*
      * SNP C/G at position 22 of gene maps to position 17 of transcript
      */
-    sf = transcriptFeatures.get(1);
+    sf = transcriptFeatures.get(2);
     assertEquals(sf.getFeatureGroup(), "VCF");
     assertEquals(sf.getBegin(), 17);
     assertEquals(sf.getEnd(), 17);
-    assertEquals(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT);
+    assertEquals(sf.getType(), SEQUENCE_VARIANT);
     assertEquals(sf.getScore(), 0f);
+    assertEquals(sf.getValue(Gff3Helper.ALLELES), "C,G");
     assertEquals(Float.parseFloat((String) sf.getValue("AF")), 2.0e-03,
             DELTA);
-    assertEquals(sf.getValue(Gff3Helper.ALLELES), "C,G");
 
     /*
      * verify variant feature(s) computed and added to protein
@@ -445,7 +485,7 @@ public class VCFLoaderTest
       }
     }
     List<SequenceFeature> proteinFeatures = peptide.getSequenceFeatures();
-    assertEquals(proteinFeatures.size(), 1);
+    assertEquals(proteinFeatures.size(), 3);
     sf = proteinFeatures.get(0);
     assertEquals(sf.getFeatureGroup(), "VCF");
     assertEquals(sf.getBegin(), 6);