Merge branch 'develop' into features/JAL-653_JAL-1766_htslib_refseqsupport
[jalview.git] / test / jalview / analysis / AlignmentUtilsTests.java
index 09bd64e..8bdd740 100644 (file)
@@ -48,6 +48,7 @@ import java.util.ArrayList;
 import java.util.Arrays;
 import java.util.HashSet;
 import java.util.Iterator;
+import java.util.LinkedHashMap;
 import java.util.List;
 import java.util.Map;
 import java.util.Set;
@@ -571,9 +572,15 @@ public class AlignmentUtilsTests
   @Test(groups = { "Functional" })
   public void testTranslatesAs()
   {
+    // null arguments check
+    assertFalse(AlignmentUtils.translatesAs(null, 0, null));
+    assertFalse(AlignmentUtils.translatesAs(new char[] { 't' }, 0, null));
+    assertFalse(AlignmentUtils.translatesAs(null, 0, new char[] { 'a' }));
+
+    // straight translation
     assertTrue(AlignmentUtils.translatesAs("tttcccaaaggg".toCharArray(), 0,
             "FPKG".toCharArray()));
-    // with start codon (not in protein)
+    // with extra start codon (not in protein)
     assertTrue(AlignmentUtils.translatesAs("atgtttcccaaaggg".toCharArray(),
             3, "FPKG".toCharArray()));
     // with stop codon1 (not in protein)
@@ -601,7 +608,7 @@ public class AlignmentUtilsTests
     assertTrue(AlignmentUtils.translatesAs(
             "atgtttcccaaagggtga".toCharArray(), 3, "FPKG".toCharArray()));
 
-    // with embedded stop codon
+    // with embedded stop codons
     assertTrue(AlignmentUtils.translatesAs(
             "atgtttTAGcccaaaTAAgggtga".toCharArray(), 3,
             "F*PK*G".toCharArray()));
@@ -609,6 +616,26 @@ public class AlignmentUtilsTests
     // wrong protein
     assertFalse(AlignmentUtils.translatesAs("tttcccaaaggg".toCharArray(),
             0, "FPMG".toCharArray()));
+
+    // truncated dna
+    assertFalse(AlignmentUtils.translatesAs("tttcccaaagg".toCharArray(), 0,
+            "FPKG".toCharArray()));
+
+    // truncated protein
+    assertFalse(AlignmentUtils.translatesAs("tttcccaaaggg".toCharArray(),
+            0, "FPK".toCharArray()));
+
+    // overlong dna (doesn't end in stop codon)
+    assertFalse(AlignmentUtils.translatesAs(
+            "tttcccaaagggttt".toCharArray(), 0, "FPKG".toCharArray()));
+
+    // dna + stop codon + more
+    assertFalse(AlignmentUtils.translatesAs(
+            "tttcccaaagggttaga".toCharArray(), 0, "FPKG".toCharArray()));
+
+    // overlong protein
+    assertFalse(AlignmentUtils.translatesAs("tttcccaaaggg".toCharArray(),
+            0, "FPKGQ".toCharArray()));
   }
 
   /**
@@ -1014,6 +1041,18 @@ public class AlignmentUtilsTests
     dna2.createDatasetSequence();
     pep1.createDatasetSequence();
     pep2.createDatasetSequence();
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds1", 4, 6, 0f,
+            null));
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds2", 10, 12, 0f,
+            null));
+    dna2.addSequenceFeature(new SequenceFeature("CDS", "cds3", 1, 3, 0f,
+            null));
+    dna2.addSequenceFeature(new SequenceFeature("CDS", "cds4", 7, 9, 0f,
+            null));
+    dna2.addSequenceFeature(new SequenceFeature("CDS", "cds5", 13, 15, 0f,
+            null));
+    AlignmentI dna = new Alignment(new SequenceI[] { dna1, dna2 });
+    dna.setDataset(null);
 
     List<AlignedCodonFrame> mappings = new ArrayList<AlignedCodonFrame>();
     MapList map = new MapList(new int[] { 4, 6, 10, 12 },
@@ -1028,10 +1067,21 @@ public class AlignmentUtilsTests
     mappings.add(acf);
 
     AlignmentI cds = AlignmentUtils.makeCdsAlignment(new SequenceI[] {
-        dna1, dna2 }, mappings);
+        dna1, dna2 }, mappings, dna);
     assertEquals(2, cds.getSequences().size());
-    assertEquals("GGGTTT", cds.getSequenceAt(0).getSequenceAsString());
-    assertEquals("GGGTTTCCC", cds.getSequenceAt(1).getSequenceAsString());
+    assertEquals("---GGG---TTT---", cds.getSequenceAt(0)
+            .getSequenceAsString());
+    assertEquals("GGG---TTT---CCC", cds.getSequenceAt(1)
+            .getSequenceAsString());
+
+    /*
+     * verify shared, extended alignment dataset
+     */
+    assertSame(dna.getDataset(), cds.getDataset());
+    assertTrue(dna.getDataset().getSequences()
+            .contains(cds.getSequenceAt(0).getDatasetSequence()));
+    assertTrue(dna.getDataset().getSequences()
+            .contains(cds.getSequenceAt(1).getDatasetSequence()));
 
     /*
      * Verify updated mappings
@@ -1048,14 +1098,14 @@ public class AlignmentUtilsTests
     SearchResults sr = MappingUtils.buildSearchResults(pep1, 1, mappings);
     assertEquals(1, sr.getResults().size());
     Match m = sr.getResults().get(0);
-    assertEquals(cds.getSequenceAt(0).getDatasetSequence(),
+    assertSame(cds.getSequenceAt(0).getDatasetSequence(),
             m.getSequence());
     assertEquals(1, m.getStart());
     assertEquals(3, m.getEnd());
     // map F to TTT
     sr = MappingUtils.buildSearchResults(pep1, 2, mappings);
     m = sr.getResults().get(0);
-    assertEquals(cds.getSequenceAt(0).getDatasetSequence(),
+    assertSame(cds.getSequenceAt(0).getDatasetSequence(),
             m.getSequence());
     assertEquals(4, m.getStart());
     assertEquals(6, m.getEnd());
@@ -1070,21 +1120,21 @@ public class AlignmentUtilsTests
     sr = MappingUtils.buildSearchResults(pep2, 1, mappings);
     assertEquals(1, sr.getResults().size());
     m = sr.getResults().get(0);
-    assertEquals(cds.getSequenceAt(1).getDatasetSequence(),
+    assertSame(cds.getSequenceAt(1).getDatasetSequence(),
             m.getSequence());
     assertEquals(1, m.getStart());
     assertEquals(3, m.getEnd());
     // map F to TTT
     sr = MappingUtils.buildSearchResults(pep2, 2, mappings);
     m = sr.getResults().get(0);
-    assertEquals(cds.getSequenceAt(1).getDatasetSequence(),
+    assertSame(cds.getSequenceAt(1).getDatasetSequence(),
             m.getSequence());
     assertEquals(4, m.getStart());
     assertEquals(6, m.getEnd());
     // map P to CCC
     sr = MappingUtils.buildSearchResults(pep2, 3, mappings);
     m = sr.getResults().get(0);
-    assertEquals(cds.getSequenceAt(1).getDatasetSequence(),
+    assertSame(cds.getSequenceAt(1).getDatasetSequence(),
             m.getSequence());
     assertEquals(7, m.getStart());
     assertEquals(9, m.getEnd());
@@ -1118,8 +1168,12 @@ public class AlignmentUtilsTests
     mappings.add(acf);
 
     AlignedCodonFrame newMapping = new AlignedCodonFrame();
+    List<int[]> ungappedColumns = new ArrayList<int[]>();
+    ungappedColumns.add(new int[] { 4, 6 });
+    ungappedColumns.add(new int[] { 10, 12 });
     List<SequenceI> cdsSeqs = AlignmentUtils.makeCdsSequences(dna1, acf,
-            newMapping);
+            ungappedColumns,
+            newMapping, '-');
     assertEquals(1, cdsSeqs.size());
     SequenceI cdsSeq = cdsSeqs.get(0);
 
@@ -1148,6 +1202,18 @@ public class AlignmentUtilsTests
     pep1.createDatasetSequence();
     pep2.createDatasetSequence();
     pep3.createDatasetSequence();
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds1", 4, 6, 0f,
+            null));
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds2", 10, 12, 0f,
+            null));
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds3", 1, 3, 0f,
+            null));
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds4", 7, 9, 0f,
+            null));
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds5", 1, 3, 0f,
+            null));
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds6", 10, 12, 0f,
+            null));
     pep1.getDatasetSequence().addDBRef(
             new DBRefEntry("EMBLCDS", "2", "A12345"));
     pep2.getDatasetSequence().addDBRef(
@@ -1156,9 +1222,7 @@ public class AlignmentUtilsTests
             new DBRefEntry("EMBLCDS", "4", "A12347"));
 
     /*
-     * Make the mappings from dna to protein. Using LinkedHashset is a
-     * convenience so results are in the input order. There is no assertion that
-     * the generated exon sequences are in any particular order.
+     * Make the mappings from dna to protein
      */
     List<AlignedCodonFrame> mappings = new ArrayList<AlignedCodonFrame>();
     // map ...GGG...TTT to GF
@@ -1184,8 +1248,10 @@ public class AlignmentUtilsTests
      * Create the Exon alignment; also replaces the dna-to-protein mappings with
      * exon-to-protein and exon-to-dna mappings
      */
+    AlignmentI dna = new Alignment(new SequenceI[] { dna1 });
+    dna.setDataset(null);
     AlignmentI exal = AlignmentUtils.makeCdsAlignment(
-            new SequenceI[] { dna1 }, mappings);
+            new SequenceI[] { dna1 }, mappings, dna);
 
     /*
      * Verify we have 3 cds sequences, mapped to pep1/2/3 respectively
@@ -1193,8 +1259,22 @@ public class AlignmentUtilsTests
     List<SequenceI> cds = exal.getSequences();
     assertEquals(3, cds.size());
 
+    /*
+     * verify shared, extended alignment dataset
+     */
+    assertSame(exal.getDataset(), dna.getDataset());
+    assertTrue(dna.getDataset().getSequences()
+            .contains(cds.get(0).getDatasetSequence()));
+    assertTrue(dna.getDataset().getSequences()
+            .contains(cds.get(1).getDatasetSequence()));
+    assertTrue(dna.getDataset().getSequences()
+            .contains(cds.get(2).getDatasetSequence()));
+
+    /*
+     * verify aligned cds sequences and their xrefs
+     */
     SequenceI cdsSeq = cds.get(0);
-    assertEquals("GGGTTT", cdsSeq.getSequenceAsString());
+    assertEquals("---GGG---TTT", cdsSeq.getSequenceAsString());
     assertEquals("dna1|A12345", cdsSeq.getName());
     assertEquals(1, cdsSeq.getDBRefs().length);
     DBRefEntry cdsRef = cdsSeq.getDBRefs()[0];
@@ -1203,7 +1283,7 @@ public class AlignmentUtilsTests
     assertEquals("A12345", cdsRef.getAccessionId());
 
     cdsSeq = cds.get(1);
-    assertEquals("aaaccc", cdsSeq.getSequenceAsString());
+    assertEquals("aaa---ccc---", cdsSeq.getSequenceAsString());
     assertEquals("dna1|A12346", cdsSeq.getName());
     assertEquals(1, cdsSeq.getDBRefs().length);
     cdsRef = cdsSeq.getDBRefs()[0];
@@ -1212,7 +1292,7 @@ public class AlignmentUtilsTests
     assertEquals("A12346", cdsRef.getAccessionId());
 
     cdsSeq = cds.get(2);
-    assertEquals("aaaTTT", cdsSeq.getSequenceAsString());
+    assertEquals("aaa------TTT", cdsSeq.getSequenceAsString());
     assertEquals("dna1|A12347", cdsSeq.getName());
     assertEquals(1, cdsSeq.getDBRefs().length);
     cdsRef = cdsSeq.getDBRefs()[0];
@@ -1289,7 +1369,7 @@ public class AlignmentUtilsTests
    * @throws IOException
    */
   @Test(groups = { "Functional" })
-  public void testMapProteinSequenceToCdna_forSubsequence()
+  public void testMapCdnaToProtein_forSubsequence()
           throws IOException
   {
     SequenceI prot = new Sequence("UNIPROT|V12345", "E-I--Q", 10, 12);
@@ -1298,7 +1378,7 @@ public class AlignmentUtilsTests
     SequenceI dna = new Sequence("EMBL|A33333", "GAA--AT-C-CAG", 40, 48);
     dna.createDatasetSequence();
 
-    MapList map = AlignmentUtils.mapProteinSequenceToCdna(prot, dna);
+    MapList map = AlignmentUtils.mapCdnaToProtein(prot, dna);
     assertEquals(10, map.getToLowest());
     assertEquals(12, map.getToHighest());
     assertEquals(40, map.getFromLowest());
@@ -1342,6 +1422,9 @@ public class AlignmentUtilsTests
             "AAA---CCCTTT---");
   }
 
+  /**
+   * Tests for transferring features between mapped sequences
+   */
   @Test(groups = { "Functional" })
   public void testTransferFeatures()
   {
@@ -1380,34 +1463,547 @@ public class AlignmentUtilsTests
             new int[] { 1, 6 }, 1, 1);
 
     /*
-     * behaviour of transferFeatures depends on MapList.locateInTo()
-     * if start and end positions are mapped, returns the mapped region
-     * if either is not mapped, does _not_ search for overlapped region 
+     * transferFeatures() will build 'partial overlap' for regions
+     * that partially overlap 5' or 3' (start or end) of target sequence
      */
-    AlignmentUtils.transferFeatures(dna, cds, map);
+    AlignmentUtils.transferFeatures(dna, cds, map, null);
     SequenceFeature[] sfs = cds.getSequenceFeatures();
-    assertEquals(4, sfs.length);
+    assertEquals(6, sfs.length);
 
     SequenceFeature sf = sfs[0];
+    assertEquals("type2", sf.getType());
+    assertEquals("desc2", sf.getDescription());
+    assertEquals(2f, sf.getScore());
+    assertEquals(1, sf.getBegin());
+    assertEquals(1, sf.getEnd());
+
+    sf = sfs[1];
     assertEquals("type3", sf.getType());
     assertEquals("desc3", sf.getDescription());
     assertEquals(3f, sf.getScore());
     assertEquals(1, sf.getBegin());
     assertEquals(3, sf.getEnd());
 
-    sf = sfs[1];
+    sf = sfs[2];
     assertEquals("type4", sf.getType());
     assertEquals(2, sf.getBegin());
     assertEquals(5, sf.getEnd());
 
-    sf = sfs[2];
+    sf = sfs[3];
     assertEquals("type5", sf.getType());
     assertEquals(1, sf.getBegin());
     assertEquals(6, sf.getEnd());
 
-    sf = sfs[3];
+    sf = sfs[4];
     assertEquals("type8", sf.getType());
     assertEquals(6, sf.getBegin());
     assertEquals(6, sf.getEnd());
+
+    sf = sfs[5];
+    assertEquals("type9", sf.getType());
+    assertEquals(6, sf.getBegin());
+    assertEquals(6, sf.getEnd());
+  }
+
+  /**
+   * Tests for transferring features between mapped sequences
+   */
+  @Test(groups = { "Functional" })
+  public void testTransferFeatures_withOmit()
+  {
+    SequenceI dna = new Sequence("dna/20-34", "acgTAGcaaGCCcgt");
+    SequenceI cds = new Sequence("cds/10-15", "TAGGCC");
+
+    MapList map = new MapList(new int[] { 4, 6, 10, 12 },
+            new int[] { 1, 6 }, 1, 1);
+  
+    // [5, 11] maps to [2, 5]
+    dna.addSequenceFeature(new SequenceFeature("type4", "desc4", 5, 11, 4f,
+            null));
+    // [4, 12] maps to [1, 6]
+    dna.addSequenceFeature(new SequenceFeature("type5", "desc5", 4, 12, 5f,
+            null));
+    // [12, 12] maps to [6, 6]
+    dna.addSequenceFeature(new SequenceFeature("type8", "desc8", 12, 12,
+            8f, null));
+  
+    // desc4 and desc8 are the 'omit these' varargs
+    AlignmentUtils.transferFeatures(dna, cds, map, null, "type4", "type8");
+    SequenceFeature[] sfs = cds.getSequenceFeatures();
+    assertEquals(1, sfs.length);
+  
+    SequenceFeature sf = sfs[0];
+    assertEquals("type5", sf.getType());
+    assertEquals(1, sf.getBegin());
+    assertEquals(6, sf.getEnd());
+  }
+
+  /**
+   * Tests for transferring features between mapped sequences
+   */
+  @Test(groups = { "Functional" })
+  public void testTransferFeatures_withSelect()
+  {
+    SequenceI dna = new Sequence("dna/20-34", "acgTAGcaaGCCcgt");
+    SequenceI cds = new Sequence("cds/10-15", "TAGGCC");
+  
+    MapList map = new MapList(new int[] { 4, 6, 10, 12 },
+            new int[] { 1, 6 }, 1, 1);
+  
+    // [5, 11] maps to [2, 5]
+    dna.addSequenceFeature(new SequenceFeature("type4", "desc4", 5, 11, 4f,
+            null));
+    // [4, 12] maps to [1, 6]
+    dna.addSequenceFeature(new SequenceFeature("type5", "desc5", 4, 12, 5f,
+            null));
+    // [12, 12] maps to [6, 6]
+    dna.addSequenceFeature(new SequenceFeature("type8", "desc8", 12, 12,
+            8f, null));
+  
+    // "type5" is the 'select this type' argument
+    AlignmentUtils.transferFeatures(dna, cds, map, "type5");
+    SequenceFeature[] sfs = cds.getSequenceFeatures();
+    assertEquals(1, sfs.length);
+  
+    SequenceFeature sf = sfs[0];
+    assertEquals("type5", sf.getType());
+    assertEquals(1, sf.getBegin());
+    assertEquals(6, sf.getEnd());
+  }
+
+  /**
+   * Test the method that extracts the cds-only part of a dna alignment, for the
+   * case where the cds should be aligned to match its nucleotide sequence.
+   */
+  @Test(groups = { "Functional" })
+  public void testMakeCdsAlignment_alternativeTranscripts()
+  {
+    SequenceI dna1 = new Sequence("dna1", "aaaGGGCC-----CTTTaaaGGG");
+    // alternative transcript of same dna skips CCC codon
+    SequenceI dna2 = new Sequence("dna2", "aaaGGGCC-----cttTaaaGGG");
+    // dna3 has no mapping (protein product) so should be ignored here
+    SequenceI dna3 = new Sequence("dna3", "aaaGGGCCCCCGGGcttTaaaGGG");
+    SequenceI pep1 = new Sequence("pep1", "GPFG");
+    SequenceI pep2 = new Sequence("pep2", "GPG");
+    dna1.createDatasetSequence();
+    dna2.createDatasetSequence();
+    dna3.createDatasetSequence();
+    pep1.createDatasetSequence();
+    pep2.createDatasetSequence();
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds1", 4, 8, 0f,
+            null));
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds2", 9, 12, 0f,
+            null));
+    dna1.addSequenceFeature(new SequenceFeature("CDS", "cds3", 16, 18, 0f,
+            null));
+    dna2.addSequenceFeature(new SequenceFeature("CDS", "cds", 4, 8, 0f,
+            null));
+    dna2.addSequenceFeature(new SequenceFeature("CDS", "cds", 12, 12, 0f,
+            null));
+    dna2.addSequenceFeature(new SequenceFeature("CDS", "cds", 16, 18, 0f,
+            null));
+  
+    List<AlignedCodonFrame> mappings = new ArrayList<AlignedCodonFrame>();
+    MapList map = new MapList(new int[] { 4, 12, 16, 18 },
+            new int[] { 1, 4 }, 3, 1);
+    AlignedCodonFrame acf = new AlignedCodonFrame();
+    acf.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(), map);
+    mappings.add(acf);
+    map = new MapList(new int[] { 4, 8, 12, 12, 16, 18 },
+            new int[] { 1, 3 },
+            3, 1);
+    acf = new AlignedCodonFrame();
+    acf.addMap(dna2.getDatasetSequence(), pep2.getDatasetSequence(), map);
+    mappings.add(acf);
+  
+    AlignmentI dna = new Alignment(new SequenceI[] { dna1, dna2, dna3 });
+    dna.setDataset(null);
+    AlignmentI cds = AlignmentUtils.makeCdsAlignment(new SequenceI[] {
+        dna1, dna2, dna3 }, mappings, dna);
+    List<SequenceI> cdsSeqs = cds.getSequences();
+    assertEquals(2, cdsSeqs.size());
+    assertEquals("GGGCCCTTTGGG", cdsSeqs.get(0).getSequenceAsString());
+    assertEquals("GGGCC---TGGG", cdsSeqs.get(1).getSequenceAsString());
+  
+    /*
+     * verify shared, extended alignment dataset
+     */
+    assertSame(dna.getDataset(), cds.getDataset());
+    assertTrue(dna.getDataset().getSequences()
+            .contains(cdsSeqs.get(0).getDatasetSequence()));
+    assertTrue(dna.getDataset().getSequences()
+            .contains(cdsSeqs.get(1).getDatasetSequence()));
+
+    /*
+     * Verify updated mappings
+     */
+    assertEquals(2, mappings.size());
+  
+    /*
+     * Mapping from pep1 to GGGTTT in first new CDS sequence
+     */
+    List<AlignedCodonFrame> pep1Mapping = MappingUtils
+            .findMappingsForSequence(pep1, mappings);
+    assertEquals(1, pep1Mapping.size());
+    /*
+     * maps GPFG to 1-3,4-6,7-9,10-12
+     */
+    SearchResults sr = MappingUtils.buildSearchResults(pep1, 1, mappings);
+    assertEquals(1, sr.getResults().size());
+    Match m = sr.getResults().get(0);
+    assertEquals(cds.getSequenceAt(0).getDatasetSequence(),
+            m.getSequence());
+    assertEquals(1, m.getStart());
+    assertEquals(3, m.getEnd());
+    sr = MappingUtils.buildSearchResults(pep1, 2, mappings);
+    m = sr.getResults().get(0);
+    assertEquals(4, m.getStart());
+    assertEquals(6, m.getEnd());
+    sr = MappingUtils.buildSearchResults(pep1, 3, mappings);
+    m = sr.getResults().get(0);
+    assertEquals(7, m.getStart());
+    assertEquals(9, m.getEnd());
+    sr = MappingUtils.buildSearchResults(pep1, 4, mappings);
+    m = sr.getResults().get(0);
+    assertEquals(10, m.getStart());
+    assertEquals(12, m.getEnd());
+  
+    /*
+     * GPG in pep2 map to 1-3,4-6,7-9 in second CDS sequence
+     */
+    List<AlignedCodonFrame> pep2Mapping = MappingUtils
+            .findMappingsForSequence(pep2, mappings);
+    assertEquals(1, pep2Mapping.size());
+    sr = MappingUtils.buildSearchResults(pep2, 1, mappings);
+    assertEquals(1, sr.getResults().size());
+    m = sr.getResults().get(0);
+    assertEquals(cds.getSequenceAt(1).getDatasetSequence(),
+            m.getSequence());
+    assertEquals(1, m.getStart());
+    assertEquals(3, m.getEnd());
+    sr = MappingUtils.buildSearchResults(pep2, 2, mappings);
+    m = sr.getResults().get(0);
+    assertEquals(4, m.getStart());
+    assertEquals(6, m.getEnd());
+    sr = MappingUtils.buildSearchResults(pep2, 3, mappings);
+    m = sr.getResults().get(0);
+    assertEquals(7, m.getStart());
+    assertEquals(9, m.getEnd());
+  }
+
+  /**
+   * Tests for gapped column in sequences
+   */
+  @Test(groups = { "Functional" })
+  public void testIsGappedColumn()
+  {
+    SequenceI seq1 = new Sequence("Seq1", "a--c.tc-a-g");
+    SequenceI seq2 = new Sequence("Seq2", "aa---t--a-g");
+    SequenceI seq3 = new Sequence("Seq3", "ag-c t-g-");
+    List<SequenceI> seqs = Arrays
+            .asList(new SequenceI[] { seq1, seq2, seq3 });
+    // the column number is base 1
+    assertFalse(AlignmentUtils.isGappedColumn(seqs, 1));
+    assertFalse(AlignmentUtils.isGappedColumn(seqs, 2));
+    assertTrue(AlignmentUtils.isGappedColumn(seqs, 3));
+    assertFalse(AlignmentUtils.isGappedColumn(seqs, 4));
+    assertTrue(AlignmentUtils.isGappedColumn(seqs, 5));
+    assertFalse(AlignmentUtils.isGappedColumn(seqs, 6));
+    assertFalse(AlignmentUtils.isGappedColumn(seqs, 7));
+    assertFalse(AlignmentUtils.isGappedColumn(seqs, 8));
+    assertFalse(AlignmentUtils.isGappedColumn(seqs, 9));
+    assertTrue(AlignmentUtils.isGappedColumn(seqs, 10));
+    assertFalse(AlignmentUtils.isGappedColumn(seqs, 11));
+    // out of bounds:
+    assertTrue(AlignmentUtils.isGappedColumn(seqs, 0));
+    assertTrue(AlignmentUtils.isGappedColumn(seqs, 100));
+    assertTrue(AlignmentUtils.isGappedColumn(seqs, -100));
+    assertTrue(AlignmentUtils.isGappedColumn(null, 0));
+  }
+
+  @Test(groups = { "Functional" })
+  public void testFindCdsColumns()
+  {
+    // TODO target method belongs in a general-purpose alignment
+    // analysis method to find columns for feature
+
+    /*
+     * NB this method assumes CDS ranges are contiguous (no introns)
+     */
+    SequenceI gene = new Sequence("gene", "aaacccgggtttaaacccgggttt");
+    SequenceI seq1 = new Sequence("Seq1", "--ac-cgGG-GGaaACC--GGtt-");
+    SequenceI seq2 = new Sequence("Seq2", "AA--CCGG--g-AAA--cG-GTTt");
+    seq1.createDatasetSequence();
+    seq2.createDatasetSequence();
+    seq1.addSequenceFeature(new SequenceFeature("CDS", "cds", 5, 6, 0f,
+            null));
+    seq1.addSequenceFeature(new SequenceFeature("CDS", "cds", 7, 8, 0f,
+            null));
+    seq1.addSequenceFeature(new SequenceFeature("CDS", "cds", 11, 13, 0f,
+            null));
+    seq1.addSequenceFeature(new SequenceFeature("CDS", "cds", 14, 15, 0f,
+            null));
+    seq2.addSequenceFeature(new SequenceFeature("CDS", "cds", 1, 2, 0f,
+            null));
+    seq2.addSequenceFeature(new SequenceFeature("CDS", "cds", 3, 6, 0f,
+            null));
+    seq2.addSequenceFeature(new SequenceFeature("CDS", "cds", 8, 10, 0f,
+            null));
+    seq2.addSequenceFeature(new SequenceFeature("CDS", "cds", 12, 12, 0f,
+            null));
+    seq2.addSequenceFeature(new SequenceFeature("CDS", "cds", 13, 15, 0f,
+            null));
+
+    List<int[]> cdsColumns = AlignmentUtils.findCdsColumns(new SequenceI[] {
+        seq1, seq2 });
+    assertEquals(4, cdsColumns.size());
+    assertEquals("[1, 2]", Arrays.toString(cdsColumns.get(0)));
+    assertEquals("[5, 9]", Arrays.toString(cdsColumns.get(1)));
+    assertEquals("[11, 17]", Arrays.toString(cdsColumns.get(2)));
+    assertEquals("[19, 23]", Arrays.toString(cdsColumns.get(3)));
+  }
+
+  /**
+   * Test the method that realigns protein to match mapped codon alignment.
+   */
+  @Test(groups = { "Functional" })
+  public void testAlignProteinAsDna_incompleteStartCodon()
+  {
+    // seq1: incomplete start codon (not mapped), then [3, 11]
+    SequenceI dna1 = new Sequence("Seq1", "ccAAA-TTT-GGG-");
+    // seq2 codons are [4, 5], [8, 11]
+    SequenceI dna2 = new Sequence("Seq2", "ccaAA-ttT-GGG-");
+    // seq3 incomplete start codon at 'tt'
+    SequenceI dna3 = new Sequence("Seq3", "ccaaa-ttt-GGG-");
+    AlignmentI dna = new Alignment(new SequenceI[] { dna1, dna2, dna3 });
+    dna.setDataset(null);
+  
+    // prot1 has 'X' for incomplete start codon (not mapped)
+    SequenceI prot1 = new Sequence("Seq1", "XKFG"); // X for incomplete start
+    SequenceI prot2 = new Sequence("Seq2", "NG");
+    SequenceI prot3 = new Sequence("Seq3", "XG"); // X for incomplete start
+    AlignmentI protein = new Alignment(new SequenceI[] { prot1, prot2,
+        prot3 });
+    protein.setDataset(null);
+  
+    // map dna1 [3, 11] to prot1 [2, 4] KFG
+    MapList map = new MapList(new int[] { 3, 11 }, new int[] { 2, 4 }, 3, 1);
+    AlignedCodonFrame acf = new AlignedCodonFrame();
+    acf.addMap(dna1.getDatasetSequence(), prot1.getDatasetSequence(), map);
+
+    // map dna2 [4, 5] [8, 11] to prot2 [1, 2] NG
+    map = new MapList(new int[] { 4, 5, 8, 11 }, new int[] { 1, 2 }, 3, 1);
+    acf.addMap(dna2.getDatasetSequence(), prot2.getDatasetSequence(), map);
+
+    // map dna3 [9, 11] to prot3 [2, 2] G
+    map = new MapList(new int[] { 9, 11 }, new int[] { 2, 2 }, 3, 1);
+    acf.addMap(dna3.getDatasetSequence(), prot3.getDatasetSequence(), map);
+
+    ArrayList<AlignedCodonFrame> acfs = new ArrayList<AlignedCodonFrame>();
+    acfs.add(acf);
+    protein.setCodonFrames(acfs);
+
+    /*
+     * verify X is included in the aligned proteins, and placed just
+     * before the first mapped residue 
+     * CCT is between CCC and TTT
+     */
+    AlignmentUtils.alignProteinAsDna(protein, dna);
+    assertEquals("XK-FG", prot1.getSequenceAsString());
+    assertEquals("--N-G", prot2.getSequenceAsString());
+    assertEquals("---XG", prot3.getSequenceAsString());
+  }
+
+  /**
+   * Tests for the method that maps the subset of a dna sequence that has CDS
+   * (or subtype) feature - case where the start codon is incomplete.
+   */
+  @Test(groups = "Functional")
+  public void testGetCdsRanges_fivePrimeIncomplete()
+  {
+    SequenceI dnaSeq = new Sequence("dna", "aaagGGCCCaaaTTTttt");
+    dnaSeq.createDatasetSequence();
+    SequenceI ds = dnaSeq.getDatasetSequence();
+  
+    // CDS for dna 5-6 (incomplete codon), 7-9
+    SequenceFeature sf = new SequenceFeature("CDS", "", 5, 9, 0f, null);
+    sf.setPhase("2"); // skip 2 bases to start of next codon
+    ds.addSequenceFeature(sf);
+    // CDS for dna 13-15
+    sf = new SequenceFeature("CDS_predicted", "", 13, 15, 0f, null);
+    ds.addSequenceFeature(sf);
+  
+    List<int[]> ranges = AlignmentUtils.findCdsPositions(dnaSeq);
+  
+    /*
+     * check the mapping starts with the first complete codon
+     */
+    assertEquals(6, MappingUtils.getLength(ranges));
+    assertEquals(2, ranges.size());
+    assertEquals(7, ranges.get(0)[0]);
+    assertEquals(9, ranges.get(0)[1]);
+    assertEquals(13, ranges.get(1)[0]);
+    assertEquals(15, ranges.get(1)[1]);
+  }
+
+  /**
+   * Tests for the method that maps the subset of a dna sequence that has CDS
+   * (or subtype) feature.
+   */
+  @Test(groups = "Functional")
+  public void testGetCdsRanges()
+  {
+    SequenceI dnaSeq = new Sequence("dna", "aaaGGGcccAAATTTttt");
+    dnaSeq.createDatasetSequence();
+    SequenceI ds = dnaSeq.getDatasetSequence();
+  
+    // CDS for dna 3-6
+    SequenceFeature sf = new SequenceFeature("CDS", "", 4, 6, 0f, null);
+    ds.addSequenceFeature(sf);
+    // exon feature should be ignored here
+    sf = new SequenceFeature("exon", "", 7, 9, 0f, null);
+    ds.addSequenceFeature(sf);
+    // CDS for dna 10-12
+    sf = new SequenceFeature("CDS_predicted", "", 10, 12, 0f, null);
+    ds.addSequenceFeature(sf);
+  
+    List<int[]> ranges = AlignmentUtils.findCdsPositions(dnaSeq);
+    assertEquals(6, MappingUtils.getLength(ranges));
+    assertEquals(2, ranges.size());
+    assertEquals(4, ranges.get(0)[0]);
+    assertEquals(6, ranges.get(0)[1]);
+    assertEquals(10, ranges.get(1)[0]);
+    assertEquals(12, ranges.get(1)[1]);
+  }
+
+  /**
+   * Test the method that computes a map of codon variants for each protein
+   * position from "sequence_variant" features on dna
+   */
+  @Test(groups = "Functional")
+  public void testBuildDnaVariantsMap()
+  {
+    SequenceI dna = new Sequence("dna", "atgAAATTTGGGCCCtag");
+    MapList map = new MapList(new int[] { 1, 18 }, new int[] { 1, 5 }, 3, 1);
+
+    /*
+     * first with no variants on dna
+     */
+    LinkedHashMap<Integer, String[][]> variantsMap = AlignmentUtils
+            .buildDnaVariantsMap(dna, map);
+    assertTrue(variantsMap.isEmpty());
+
+    // single allele codon 1, on base 1
+    SequenceFeature sf = new SequenceFeature("sequence_variant", "", 1, 1,
+            0f, null);
+    sf.setValue("alleles", "T");
+    dna.addSequenceFeature(sf);
+
+    // two alleles codon 2, on bases 2 and 3
+    sf = new SequenceFeature("sequence_variant", "", 5, 5, 0f, null);
+    sf.setValue("alleles", "T");
+    dna.addSequenceFeature(sf);
+    sf = new SequenceFeature("sequence_variant", "", 6, 6, 0f, null);
+    sf.setValue("alleles", "G");
+    dna.addSequenceFeature(sf);
+
+    // two alleles codon 3, both on base 2
+    sf = new SequenceFeature("sequence_variant", "", 8, 8, 0f, null);
+    sf.setValue("alleles", "C, G");
+    dna.addSequenceFeature(sf);
+
+    // no alleles on codon 4
+    // alleles on codon 5 on all 3 bases
+    sf = new SequenceFeature("sequence_variant", "", 13, 13, 0f, null);
+    sf.setValue("alleles", "C, G"); // (C duplicates given base value)
+    dna.addSequenceFeature(sf);
+    sf = new SequenceFeature("sequence_variant", "", 14, 14, 0f, null);
+    sf.setValue("alleles", "g, a"); // should force to upper-case
+    dna.addSequenceFeature(sf);
+    sf = new SequenceFeature("sequence_variant", "", 15, 15, 0f, null);
+    sf.setValue("alleles", "A, T");
+    dna.addSequenceFeature(sf);
+
+    variantsMap = AlignmentUtils.buildDnaVariantsMap(dna, map);
+    assertEquals(4, variantsMap.size());
+    assertTrue(Arrays.deepEquals(new String[][] { { "A", "T" }, { "T" },
+        { "G" } }, variantsMap.get(1)));
+    assertTrue(Arrays.deepEquals(new String[][] { { "A" }, { "A", "T" },
+        { "A", "G" } }, variantsMap.get(2)));
+    assertTrue(Arrays.deepEquals(new String[][] { { "T" },
+        { "T", "C", "G" }, { "T" } }, variantsMap.get(3)));
+    // duplicated bases are not removed here, handled in computePeptideVariants
+    assertTrue(Arrays.deepEquals(new String[][] { { "C", "C", "G" },
+        { "C", "G", "A" }, { "C", "A", "T" } }, variantsMap.get(5)));
+  }
+
+  /**
+   * Tests for the method that computes all peptide variants given codon
+   * variants
+   */
+  @Test(groups = "Functional")
+  public void testComputePeptideVariants()
+  {
+    String[][] codonVariants = new String[][] { { "A" }, { "G" }, { "T" } };
+  
+    /*
+     * AGT codes for S - this is not included in the variants returned
+     */
+    List<String> variants = AlignmentUtils.computePeptideVariants(codonVariants, "S");
+    assertEquals("[]", variants.toString());
+  
+    // S is reported if it differs from the current value (A):
+    variants = AlignmentUtils.computePeptideVariants(codonVariants, "A");
+    assertEquals("[S]", variants.toString());
+  
+    /*
+     * synonymous variant is not reported
+     */
+    codonVariants = new String[][] { { "A" }, { "G" }, { "C", "T" } };
+    // AGC and AGT both code for S
+    variants = AlignmentUtils.computePeptideVariants(codonVariants, "s");
+    assertEquals("[]", variants.toString());
+  
+    /*
+     * equivalent variants are only reported once
+     */
+    codonVariants = new String[][] { { "C" }, { "T" },
+        { "A", "C", "G", "T" } };
+    // CTA CTC CTG CTT all code for L
+    variants = AlignmentUtils.computePeptideVariants(codonVariants, "S");
+    assertEquals("[L]", variants.toString());
+  
+    /*
+     * vary codons 1 and 2; variant products are sorted and non-redundant
+     */
+    codonVariants = new String[][] { { "a", "C" }, { "g", "T" }, { "A" } };
+    // aga ata cga cta code for R, I, R, L
+    variants = AlignmentUtils.computePeptideVariants(codonVariants, "S");
+    assertEquals("[I, L, R]", variants.toString());
+  
+    /*
+     * vary codons 2 and 3
+     */
+    codonVariants = new String[][] { { "a" }, { "g", "T" }, { "A", "c" } };
+    // aga agc ata atc code for R, S, I, I
+    variants = AlignmentUtils.computePeptideVariants(codonVariants, "S");
+    assertEquals("[I, R]", variants.toString());
+  
+    /*
+     * vary codons 1 and 3
+     */
+    codonVariants = new String[][] { { "a", "t" }, { "a" }, { "t", "g" } };
+    // aat aag tat tag code for N, K, Y, STOP - STOP sorted to end
+    variants = AlignmentUtils.computePeptideVariants(codonVariants, "S");
+    assertEquals("[K, N, Y, STOP]", variants.toString());
+  
+    /*
+     * vary codons 1, 2 and 3
+     */
+    codonVariants = new String[][] { { "a", "t" }, { "G", "C" },
+        { "t", "g" } };
+    // agt agg act acg tgt tgg tct tcg code for S, R, T, T, C, W, S, S
+    variants = AlignmentUtils.computePeptideVariants(codonVariants, "S");
+    assertEquals("[C, R, T, W]", variants.toString());
   }
 }