test/jalview/analysis/AlignmentUtilsTests.java

   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.analysis;
  22
  23 import static org.testng.Assert.assertNotEquals;
  24 import static org.testng.AssertJUnit.assertEquals;
  25 import static org.testng.AssertJUnit.assertFalse;
  26 import static org.testng.AssertJUnit.assertNotNull;
  27 import static org.testng.AssertJUnit.assertNull;
  28 import static org.testng.AssertJUnit.assertSame;
  29 import static org.testng.AssertJUnit.assertTrue;
  30
  31 import java.awt.Color;
  32 import java.io.IOException;
  33 import java.util.ArrayList;
  34 import java.util.Arrays;
  35 import java.util.HashMap;
  36 import java.util.LinkedHashMap;
  37 import java.util.List;
  38 import java.util.Map;
  39 import java.util.SortedMap;
  40 import java.util.TreeMap;
  41
  42 import org.testng.Assert;
  43 import org.testng.annotations.BeforeClass;
  44 import org.testng.annotations.DataProvider;
  45 import org.testng.annotations.Test;
  46
  47 import jalview.datamodel.AlignedCodonFrame;
  48 import jalview.datamodel.Alignment;
  49 import jalview.datamodel.AlignmentAnnotation;
  50 import jalview.datamodel.AlignmentI;
  51 import jalview.datamodel.Annotation;
  52 import jalview.datamodel.ContactListI;
  53 import jalview.datamodel.ContactMatrixI;
  54 import jalview.datamodel.DBRefEntry;
  55 import jalview.datamodel.GeneLociI;
  56 import jalview.datamodel.Mapping;
  57 import jalview.datamodel.SearchResultMatchI;
  58 import jalview.datamodel.SearchResultsI;
  59 import jalview.datamodel.SeqDistanceContactMatrix;
  60 import jalview.datamodel.Sequence;
  61 import jalview.datamodel.SequenceFeature;
  62 import jalview.datamodel.SequenceGroup;
  63 import jalview.datamodel.SequenceI;
  64 import jalview.gui.JvOptionPane;
  65 import jalview.io.AppletFormatAdapter;
  66 import jalview.io.DataSourceType;
  67 import jalview.io.FileFormat;
  68 import jalview.io.FileFormatI;
  69 import jalview.io.FormatAdapter;
  70 import jalview.io.gff.SequenceOntologyI;
  71 import jalview.util.Comparison;
  72 import jalview.util.MapList;
  73 import jalview.util.MappingUtils;
  74
  75 public class AlignmentUtilsTests
  76 {
  77   private static Sequence ts = new Sequence("short",
  78           "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklm");
  79
  80   @BeforeClass(alwaysRun = true)
  81   public void setUpJvOptionPane()
  82   {
  83     JvOptionPane.setInteractiveMode(false);
  84     JvOptionPane.setMockResponse(JvOptionPane.CANCEL_OPTION);
  85
  86     AlignmentAnnotation ann1 = new AlignmentAnnotation("Secondary Structure", "Secondary Structure",
  87             new Annotation[] {});
  88     AlignmentAnnotation ann2 = new AlignmentAnnotation("jnetpred", "jnetpred",
  89             new Annotation[] {});
  90     AlignmentAnnotation ann3 = new AlignmentAnnotation("Temp", "Temp",
  91             new Annotation[] {});
  92     AlignmentAnnotation ann4 = new AlignmentAnnotation("Temp", "Temp",
  93             new Annotation[] {});
  94
  95     AlignmentAnnotation[] anns1 = new AlignmentAnnotation[] {ann1, ann3, ann4};
  96
  97     AlignmentAnnotation[] anns2 = new AlignmentAnnotation[] {ann2, ann3, ann4};
  98
  99     AlignmentAnnotation[] anns3 = new AlignmentAnnotation[] {ann3, ann4};
 100
 101     AlignmentAnnotation[] anns4 = new AlignmentAnnotation[0];
 102
 103     AlignmentAnnotation[] anns5 = new AlignmentAnnotation[] {ann1, ann2, ann3, ann4};
 104   }
 105
 106   @Test(groups = { "Functional" })
 107   public void testExpandContext()
 108   {
 109     AlignmentI al = new Alignment(new Sequence[] {});
 110     for (int i = 4; i < 14; i += 2)
 111     {
 112       SequenceI s1 = ts.deriveSequence().getSubSequence(i, i + 7);
 113       al.addSequence(s1);
 114     }
 115     System.out.println(new AppletFormatAdapter()
 116             .formatSequences(FileFormat.Clustal, al, true));
 117     for (int flnk = -1; flnk < 25; flnk++)
 118     {
 119       AlignmentI exp = AlignmentUtils.expandContext(al, flnk);
 120       System.out.println("\nFlank size: " + flnk);
 121       System.out.println(new AppletFormatAdapter()
 122               .formatSequences(FileFormat.Clustal, exp, true));
 123       if (flnk == -1)
 124       {
 125         /*
 126          * Full expansion to complete sequences
 127          */
 128         for (SequenceI sq : exp.getSequences())
 129         {
 130           String ung = sq.getSequenceAsString().replaceAll("-+", "");
 131           final String errorMsg = "Flanking sequence not the same as original dataset sequence.\n"
 132                   + ung + "\n"
 133                   + sq.getDatasetSequence().getSequenceAsString();
 134           assertTrue(errorMsg, ung.equalsIgnoreCase(
 135                   sq.getDatasetSequence().getSequenceAsString()));
 136         }
 137       }
 138       else if (flnk == 24)
 139       {
 140         /*
 141          * Last sequence is fully expanded, others have leading gaps to match
 142          */
 143         assertTrue(exp.getSequenceAt(4).getSequenceAsString()
 144                 .startsWith("abc"));
 145         assertTrue(exp.getSequenceAt(3).getSequenceAsString()
 146                 .startsWith("--abc"));
 147         assertTrue(exp.getSequenceAt(2).getSequenceAsString()
 148                 .startsWith("----abc"));
 149         assertTrue(exp.getSequenceAt(1).getSequenceAsString()
 150                 .startsWith("------abc"));
 151         assertTrue(exp.getSequenceAt(0).getSequenceAsString()
 152                 .startsWith("--------abc"));
 153       }
 154     }
 155   }
 156
 157   /**
 158    * Test that annotations are correctly adjusted by expandContext
 159    */
 160   @Test(groups = { "Functional" })
 161   public void testExpandContext_annotation()
 162   {
 163     AlignmentI al = new Alignment(new Sequence[] {});
 164     SequenceI ds = new Sequence("Seq1", "ABCDEFGHI");
 165     // subsequence DEF:
 166     SequenceI seq1 = ds.deriveSequence().getSubSequence(3, 6);
 167     al.addSequence(seq1);
 168
 169     /*
 170      * Annotate DEF with 4/5/6 respectively
 171      */
 172     Annotation[] anns = new Annotation[] { new Annotation(4),
 173         new Annotation(5), new Annotation(6) };
 174     AlignmentAnnotation ann = new AlignmentAnnotation("SS",
 175             "secondary structure", anns);
 176     seq1.addAlignmentAnnotation(ann);
 177
 178     /*
 179      * The annotations array should match aligned positions
 180      */
 181     assertEquals(3, ann.annotations.length);
 182     assertEquals(4, ann.annotations[0].value, 0.001);
 183     assertEquals(5, ann.annotations[1].value, 0.001);
 184     assertEquals(6, ann.annotations[2].value, 0.001);
 185
 186     /*
 187      * Check annotation to sequence position mappings before expanding the
 188      * sequence; these are set up in Sequence.addAlignmentAnnotation ->
 189      * Annotation.setSequenceRef -> createSequenceMappings
 190      */
 191     assertNull(ann.getAnnotationForPosition(1));
 192     assertNull(ann.getAnnotationForPosition(2));
 193     assertNull(ann.getAnnotationForPosition(3));
 194     assertEquals(4, ann.getAnnotationForPosition(4).value, 0.001);
 195     assertEquals(5, ann.getAnnotationForPosition(5).value, 0.001);
 196     assertEquals(6, ann.getAnnotationForPosition(6).value, 0.001);
 197     assertNull(ann.getAnnotationForPosition(7));
 198     assertNull(ann.getAnnotationForPosition(8));
 199     assertNull(ann.getAnnotationForPosition(9));
 200
 201     /*
 202      * Expand the subsequence to the full sequence abcDEFghi
 203      */
 204     AlignmentI expanded = AlignmentUtils.expandContext(al, -1);
 205     assertEquals("abcDEFghi",
 206             expanded.getSequenceAt(0).getSequenceAsString());
 207
 208     /*
 209      * Confirm the alignment and sequence have the same SS annotation,
 210      * referencing the expanded sequence
 211      */
 212     ann = expanded.getSequenceAt(0).getAnnotation()[0];
 213     assertSame(ann, expanded.getAlignmentAnnotation()[0]);
 214     assertSame(expanded.getSequenceAt(0), ann.sequenceRef);
 215
 216     /*
 217      * The annotations array should have null values except for annotated
 218      * positions
 219      */
 220     assertNull(ann.annotations[0]);
 221     assertNull(ann.annotations[1]);
 222     assertNull(ann.annotations[2]);
 223     assertEquals(4, ann.annotations[3].value, 0.001);
 224     assertEquals(5, ann.annotations[4].value, 0.001);
 225     assertEquals(6, ann.annotations[5].value, 0.001);
 226     assertNull(ann.annotations[6]);
 227     assertNull(ann.annotations[7]);
 228     assertNull(ann.annotations[8]);
 229
 230     /*
 231      * sequence position mappings should be unchanged
 232      */
 233     assertNull(ann.getAnnotationForPosition(1));
 234     assertNull(ann.getAnnotationForPosition(2));
 235     assertNull(ann.getAnnotationForPosition(3));
 236     assertEquals(4, ann.getAnnotationForPosition(4).value, 0.001);
 237     assertEquals(5, ann.getAnnotationForPosition(5).value, 0.001);
 238     assertEquals(6, ann.getAnnotationForPosition(6).value, 0.001);
 239     assertNull(ann.getAnnotationForPosition(7));
 240     assertNull(ann.getAnnotationForPosition(8));
 241     assertNull(ann.getAnnotationForPosition(9));
 242   }
 243
 244   /**
 245    * Test method that returns a map of lists of sequences by sequence name.
 246    *
 247    * @throws IOException
 248    */
 249   @Test(groups = { "Functional" })
 250   public void testGetSequencesByName() throws IOException
 251   {
 252     final String data = ">Seq1Name\nKQYL\n" + ">Seq2Name\nRFPW\n"
 253             + ">Seq1Name\nABCD\n";
 254     AlignmentI al = loadAlignment(data, FileFormat.Fasta);
 255     Map<String, List<SequenceI>> map = AlignmentUtils
 256             .getSequencesByName(al);
 257     assertEquals(2, map.keySet().size());
 258     assertEquals(2, map.get("Seq1Name").size());
 259     assertEquals("KQYL", map.get("Seq1Name").get(0).getSequenceAsString());
 260     assertEquals("ABCD", map.get("Seq1Name").get(1).getSequenceAsString());
 261     assertEquals(1, map.get("Seq2Name").size());
 262     assertEquals("RFPW", map.get("Seq2Name").get(0).getSequenceAsString());
 263   }
 264
 265   /**
 266    * Helper method to load an alignment and ensure dataset sequences are set up.
 267    *
 268    * @param data
 269    * @param format
 270    *          TODO
 271    * @return
 272    * @throws IOException
 273    */
 274   protected AlignmentI loadAlignment(final String data, FileFormatI format)
 275           throws IOException
 276   {
 277     AlignmentI a = new FormatAdapter().readFile(data, DataSourceType.PASTE,
 278             format);
 279     a.setDataset(null);
 280     return a;
 281   }
 282
 283   /**
 284    * Test mapping of protein to cDNA, for the case where we have no sequence
 285    * cross-references, so mappings are made first-served 1-1 where sequences
 286    * translate.
 287    *
 288    * @throws IOException
 289    */
 290   @Test(groups = { "Functional" })
 291   public void testMapProteinAlignmentToCdna_noXrefs() throws IOException
 292   {
 293     List<SequenceI> protseqs = new ArrayList<>();
 294     protseqs.add(new Sequence("UNIPROT|V12345", "EIQ"));
 295     protseqs.add(new Sequence("UNIPROT|V12346", "EIQ"));
 296     protseqs.add(new Sequence("UNIPROT|V12347", "SAR"));
 297     AlignmentI protein = new Alignment(protseqs.toArray(new SequenceI[3]));
 298     protein.setDataset(null);
 299
 300     List<SequenceI> dnaseqs = new ArrayList<>();
 301     dnaseqs.add(new Sequence("EMBL|A11111", "TCAGCACGC")); // = SAR
 302     dnaseqs.add(new Sequence("EMBL|A22222", "GAGATACAA")); // = EIQ
 303     dnaseqs.add(new Sequence("EMBL|A33333", "GAAATCCAG")); // = EIQ
 304     dnaseqs.add(new Sequence("EMBL|A44444", "GAAATTCAG")); // = EIQ
 305     AlignmentI cdna = new Alignment(dnaseqs.toArray(new SequenceI[4]));
 306     cdna.setDataset(null);
 307
 308     assertTrue(AlignmentUtils.mapProteinAlignmentToCdna(protein, cdna));
 309
 310     // 3 mappings made, each from 1 to 1 sequence
 311     assertEquals(3, protein.getCodonFrames().size());
 312     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(0)).size());
 313     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(1)).size());
 314     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(2)).size());
 315
 316     // V12345 mapped to A22222
 317     AlignedCodonFrame acf = protein.getCodonFrame(protein.getSequenceAt(0))
 318             .get(0);
 319     assertEquals(1, acf.getdnaSeqs().length);
 320     assertEquals(cdna.getSequenceAt(1).getDatasetSequence(),
 321             acf.getdnaSeqs()[0]);
 322     Mapping[] protMappings = acf.getProtMappings();
 323     assertEquals(1, protMappings.length);
 324     MapList mapList = protMappings[0].getMap();
 325     assertEquals(3, mapList.getFromRatio());
 326     assertEquals(1, mapList.getToRatio());
 327     assertTrue(
 328             Arrays.equals(new int[]
 329             { 1, 9 }, mapList.getFromRanges().get(0)));
 330     assertEquals(1, mapList.getFromRanges().size());
 331     assertTrue(
 332             Arrays.equals(new int[]
 333             { 1, 3 }, mapList.getToRanges().get(0)));
 334     assertEquals(1, mapList.getToRanges().size());
 335
 336     // V12346 mapped to A33333
 337     acf = protein.getCodonFrame(protein.getSequenceAt(1)).get(0);
 338     assertEquals(1, acf.getdnaSeqs().length);
 339     assertEquals(cdna.getSequenceAt(2).getDatasetSequence(),
 340             acf.getdnaSeqs()[0]);
 341
 342     // V12347 mapped to A11111
 343     acf = protein.getCodonFrame(protein.getSequenceAt(2)).get(0);
 344     assertEquals(1, acf.getdnaSeqs().length);
 345     assertEquals(cdna.getSequenceAt(0).getDatasetSequence(),
 346             acf.getdnaSeqs()[0]);
 347
 348     // no mapping involving the 'extra' A44444
 349     assertTrue(protein.getCodonFrame(cdna.getSequenceAt(3)).isEmpty());
 350   }
 351
 352   /**
 353    * Test for the alignSequenceAs method that takes two sequences and a mapping.
 354    */
 355   @Test(groups = { "Functional" })
 356   public void testAlignSequenceAs_withMapping_noIntrons()
 357   {
 358     MapList map = new MapList(new int[] { 1, 6 }, new int[] { 1, 2 }, 3, 1);
 359
 360     /*
 361      * No existing gaps in dna:
 362      */
 363     checkAlignSequenceAs("GGGAAA", "-A-L-", false, false, map,
 364             "---GGG---AAA");
 365
 366     /*
 367      * Now introduce gaps in dna but ignore them when realigning.
 368      */
 369     checkAlignSequenceAs("-G-G-G-A-A-A-", "-A-L-", false, false, map,
 370             "---GGG---AAA");
 371
 372     /*
 373      * Now include gaps in dna when realigning. First retaining 'mapped' gaps
 374      * only, i.e. those within the exon region.
 375      */
 376     checkAlignSequenceAs("-G-G--G-A--A-A-", "-A-L-", true, false, map,
 377             "---G-G--G---A--A-A");
 378
 379     /*
 380      * Include all gaps in dna when realigning (within and without the exon
 381      * region). The leading gap, and the gaps between codons, are subsumed by
 382      * the protein alignment gap.
 383      */
 384     checkAlignSequenceAs("-G-GG--AA-A---", "-A-L-", true, true, map,
 385             "---G-GG---AA-A---");
 386
 387     /*
 388      * Include only unmapped gaps in dna when realigning (outside the exon
 389      * region). The leading gap, and the gaps between codons, are subsumed by
 390      * the protein alignment gap.
 391      */
 392     checkAlignSequenceAs("-G-GG--AA-A-", "-A-L-", false, true, map,
 393             "---GGG---AAA---");
 394   }
 395
 396   /**
 397    * Test for the alignSequenceAs method that takes two sequences and a mapping.
 398    */
 399   @Test(groups = { "Functional" })
 400   public void testAlignSequenceAs_withMapping_withIntrons()
 401   {
 402     /*
 403      * Exons at codon 2 (AAA) and 4 (TTT)
 404      */
 405     MapList map = new MapList(new int[] { 4, 6, 10, 12 },
 406             new int[]
 407             { 1, 2 }, 3, 1);
 408
 409     /*
 410      * Simple case: no gaps in dna
 411      */
 412     checkAlignSequenceAs("GGGAAACCCTTTGGG", "--A-L-", false, false, map,
 413             "GGG---AAACCCTTTGGG");
 414
 415     /*
 416      * Add gaps to dna - but ignore when realigning.
 417      */
 418     checkAlignSequenceAs("-G-G-G--A--A---AC-CC-T-TT-GG-G-", "--A-L-", false,
 419             false, map, "GGG---AAACCCTTTGGG");
 420
 421     /*
 422      * Add gaps to dna - include within exons only when realigning.
 423      */
 424     checkAlignSequenceAs("-G-G-G--A--A---A-C-CC-T-TT-GG-G-", "--A-L-", true,
 425             false, map, "GGG---A--A---ACCCT-TTGGG");
 426
 427     /*
 428      * Include gaps outside exons only when realigning.
 429      */
 430     checkAlignSequenceAs("-G-G-G--A--A---A-C-CC-T-TT-GG-G-", "--A-L-",
 431             false, true, map, "-G-G-GAAAC-CCTTT-GG-G-");
 432
 433     /*
 434      * Include gaps following first intron if we are 'preserving mapped gaps'
 435      */
 436     checkAlignSequenceAs("-G-G-G--A--A---A-C-CC-T-TT-GG-G-", "--A-L-", true,
 437             true, map, "-G-G-G--A--A---A-C-CC-T-TT-GG-G-");
 438
 439     /*
 440      * Include all gaps in dna when realigning.
 441      */
 442     checkAlignSequenceAs("-G-G-G--A--A---A-C-CC-T-TT-GG-G-", "--A-L-", true,
 443             true, map, "-G-G-G--A--A---A-C-CC-T-TT-GG-G-");
 444   }
 445
 446   /**
 447    * Test for the case where not all of the protein sequence is mapped to cDNA.
 448    */
 449   @Test(groups = { "Functional" })
 450   public void testAlignSequenceAs_withMapping_withUnmappedProtein()
 451   {
 452     /*
 453      * Exons at codon 2 (AAA) and 4 (TTT) mapped to A and P
 454      */
 455     final MapList map = new MapList(new int[] { 4, 6, 10, 12 },
 456             new int[]
 457             { 1, 1, 3, 3 }, 3, 1);
 458
 459     /*
 460      * -L- 'aligns' ccc------
 461      */
 462     checkAlignSequenceAs("gggAAAcccTTTggg", "-A-L-P-", false, false, map,
 463             "gggAAAccc------TTTggg");
 464   }
 465
 466   /**
 467    * Helper method that performs and verifies the method under test.
 468    *
 469    * @param alignee
 470    *          the sequence to be realigned
 471    * @param alignModel
 472    *          the sequence whose alignment is to be copied
 473    * @param preserveMappedGaps
 474    * @param preserveUnmappedGaps
 475    * @param map
 476    * @param expected
 477    */
 478   protected void checkAlignSequenceAs(final String alignee,
 479           final String alignModel, final boolean preserveMappedGaps,
 480           final boolean preserveUnmappedGaps, MapList map,
 481           final String expected)
 482   {
 483     SequenceI alignMe = new Sequence("Seq1", alignee);
 484     alignMe.createDatasetSequence();
 485     SequenceI alignFrom = new Sequence("Seq2", alignModel);
 486     alignFrom.createDatasetSequence();
 487     AlignedCodonFrame acf = new AlignedCodonFrame();
 488     acf.addMap(alignMe.getDatasetSequence(), alignFrom.getDatasetSequence(),
 489             map);
 490
 491     AlignmentUtils.alignSequenceAs(alignMe, alignFrom, acf, "---", '-',
 492             preserveMappedGaps, preserveUnmappedGaps);
 493     assertEquals(expected, alignMe.getSequenceAsString());
 494   }
 495
 496   /**
 497    * Test for the alignSequenceAs method where we preserve gaps in introns only.
 498    */
 499   @Test(groups = { "Functional" })
 500   public void testAlignSequenceAs_keepIntronGapsOnly()
 501   {
 502
 503     /*
 504      * Intron GGGAAA followed by exon CCCTTT
 505      */
 506     MapList map = new MapList(new int[] { 7, 12 }, new int[] { 1, 2 }, 3,
 507             1);
 508
 509     checkAlignSequenceAs("GG-G-AA-A-C-CC-T-TT", "AL", false, true, map,
 510             "GG-G-AA-ACCCTTT");
 511   }
 512
 513   /**
 514    * Test the method that realigns protein to match mapped codon alignment.
 515    */
 516   @Test(groups = { "Functional" })
 517   public void testAlignProteinAsDna()
 518   {
 519     // seq1 codons are [1,2,3] [4,5,6] [7,8,9] [10,11,12]
 520     SequenceI dna1 = new Sequence("Seq1", "TGCCATTACCAG-");
 521     // seq2 codons are [1,3,4] [5,6,7] [8,9,10] [11,12,13]
 522     SequenceI dna2 = new Sequence("Seq2", "T-GCCATTACCAG");
 523     // seq3 codons are [1,2,3] [4,5,7] [8,9,10] [11,12,13]
 524     SequenceI dna3 = new Sequence("Seq3", "TGCCA-TTACCAG");
 525     AlignmentI dna = new Alignment(new SequenceI[] { dna1, dna2, dna3 });
 526     dna.setDataset(null);
 527
 528     // protein alignment will be realigned like dna
 529     SequenceI prot1 = new Sequence("Seq1", "CHYQ");
 530     SequenceI prot2 = new Sequence("Seq2", "CHYQ");
 531     SequenceI prot3 = new Sequence("Seq3", "CHYQ");
 532     SequenceI prot4 = new Sequence("Seq4", "R-QSV"); // unmapped, unchanged
 533     AlignmentI protein = new Alignment(
 534             new SequenceI[]
 535             { prot1, prot2, prot3, prot4 });
 536     protein.setDataset(null);
 537
 538     MapList map = new MapList(new int[] { 1, 12 }, new int[] { 1, 4 }, 3,
 539             1);
 540     AlignedCodonFrame acf = new AlignedCodonFrame();
 541     acf.addMap(dna1.getDatasetSequence(), prot1.getDatasetSequence(), map);
 542     acf.addMap(dna2.getDatasetSequence(), prot2.getDatasetSequence(), map);
 543     acf.addMap(dna3.getDatasetSequence(), prot3.getDatasetSequence(), map);
 544     ArrayList<AlignedCodonFrame> acfs = new ArrayList<>();
 545     acfs.add(acf);
 546     protein.setCodonFrames(acfs);
 547
 548     /*
 549      * Translated codon order is [1,2,3] [1,3,4] [4,5,6] [4,5,7] [5,6,7] [7,8,9]
 550      * [8,9,10] [10,11,12] [11,12,13]
 551      */
 552     AlignmentUtils.alignProteinAsDna(protein, dna);
 553     assertEquals("C-H--Y-Q-", prot1.getSequenceAsString());
 554     assertEquals("-C--H-Y-Q", prot2.getSequenceAsString());
 555     assertEquals("C--H--Y-Q", prot3.getSequenceAsString());
 556     assertEquals("R-QSV", prot4.getSequenceAsString());
 557   }
 558
 559   /**
 560    * Test the method that tests whether a CDNA sequence translates to a protein
 561    * sequence
 562    */
 563   @Test(groups = { "Functional" })
 564   public void testTranslatesAs()
 565   {
 566     // null arguments check
 567     assertFalse(AlignmentUtils.translatesAs(null, 0, null));
 568     assertFalse(AlignmentUtils.translatesAs(new char[] { 't' }, 0, null));
 569     assertFalse(AlignmentUtils.translatesAs(null, 0, new char[] { 'a' }));
 570
 571     // straight translation
 572     assertTrue(AlignmentUtils.translatesAs("tttcccaaaggg".toCharArray(), 0,
 573             "FPKG".toCharArray()));
 574     // with extra start codon (not in protein)
 575     assertTrue(AlignmentUtils.translatesAs("atgtttcccaaaggg".toCharArray(),
 576             3, "FPKG".toCharArray()));
 577     // with stop codon1 (not in protein)
 578     assertTrue(AlignmentUtils.translatesAs("tttcccaaagggtaa".toCharArray(),
 579             0, "FPKG".toCharArray()));
 580     // with stop codon1 (in protein as *)
 581     assertTrue(AlignmentUtils.translatesAs("tttcccaaagggtaa".toCharArray(),
 582             0, "FPKG*".toCharArray()));
 583     // with stop codon2 (not in protein)
 584     assertTrue(AlignmentUtils.translatesAs("tttcccaaagggtag".toCharArray(),
 585             0, "FPKG".toCharArray()));
 586     // with stop codon3 (not in protein)
 587     assertTrue(AlignmentUtils.translatesAs("tttcccaaagggtga".toCharArray(),
 588             0, "FPKG".toCharArray()));
 589     // with start and stop codon1
 590     assertTrue(AlignmentUtils.translatesAs(
 591             "atgtttcccaaagggtaa".toCharArray(), 3, "FPKG".toCharArray()));
 592     // with start and stop codon1 (in protein as *)
 593     assertTrue(AlignmentUtils.translatesAs(
 594             "atgtttcccaaagggtaa".toCharArray(), 3, "FPKG*".toCharArray()));
 595     // with start and stop codon2
 596     assertTrue(AlignmentUtils.translatesAs(
 597             "atgtttcccaaagggtag".toCharArray(), 3, "FPKG".toCharArray()));
 598     // with start and stop codon3
 599     assertTrue(AlignmentUtils.translatesAs(
 600             "atgtttcccaaagggtga".toCharArray(), 3, "FPKG".toCharArray()));
 601
 602     // with embedded stop codons
 603     assertTrue(AlignmentUtils.translatesAs(
 604             "atgtttTAGcccaaaTAAgggtga".toCharArray(), 3,
 605             "F*PK*G".toCharArray()));
 606
 607     // wrong protein
 608     assertFalse(AlignmentUtils.translatesAs("tttcccaaaggg".toCharArray(), 0,
 609             "FPMG".toCharArray()));
 610
 611     // truncated dna
 612     assertFalse(AlignmentUtils.translatesAs("tttcccaaagg".toCharArray(), 0,
 613             "FPKG".toCharArray()));
 614
 615     // truncated protein
 616     assertFalse(AlignmentUtils.translatesAs("tttcccaaaggg".toCharArray(), 0,
 617             "FPK".toCharArray()));
 618
 619     // overlong dna (doesn't end in stop codon)
 620     assertFalse(AlignmentUtils.translatesAs("tttcccaaagggttt".toCharArray(),
 621             0, "FPKG".toCharArray()));
 622
 623     // dna + stop codon + more
 624     assertFalse(AlignmentUtils.translatesAs(
 625             "tttcccaaagggttaga".toCharArray(), 0, "FPKG".toCharArray()));
 626
 627     // overlong protein
 628     assertFalse(AlignmentUtils.translatesAs("tttcccaaaggg".toCharArray(), 0,
 629             "FPKGQ".toCharArray()));
 630   }
 631
 632   /**
 633    * Test mapping of protein to cDNA, for cases where the cDNA has start and/or
 634    * stop codons in addition to the protein coding sequence.
 635    *
 636    * @throws IOException
 637    */
 638   @Test(groups = { "Functional" })
 639   public void testMapProteinAlignmentToCdna_withStartAndStopCodons()
 640           throws IOException
 641   {
 642     List<SequenceI> protseqs = new ArrayList<>();
 643     protseqs.add(new Sequence("UNIPROT|V12345", "EIQ"));
 644     protseqs.add(new Sequence("UNIPROT|V12346", "EIQ"));
 645     protseqs.add(new Sequence("UNIPROT|V12347", "SAR"));
 646     AlignmentI protein = new Alignment(protseqs.toArray(new SequenceI[3]));
 647     protein.setDataset(null);
 648
 649     List<SequenceI> dnaseqs = new ArrayList<>();
 650     // start + SAR:
 651     dnaseqs.add(new Sequence("EMBL|A11111", "ATGTCAGCACGC"));
 652     // = EIQ + stop
 653     dnaseqs.add(new Sequence("EMBL|A22222", "GAGATACAATAA"));
 654     // = start +EIQ + stop
 655     dnaseqs.add(new Sequence("EMBL|A33333", "ATGGAAATCCAGTAG"));
 656     dnaseqs.add(new Sequence("EMBL|A44444", "GAAATTCAG"));
 657     AlignmentI cdna = new Alignment(dnaseqs.toArray(new SequenceI[4]));
 658     cdna.setDataset(null);
 659
 660     assertTrue(AlignmentUtils.mapProteinAlignmentToCdna(protein, cdna));
 661
 662     // 3 mappings made, each from 1 to 1 sequence
 663     assertEquals(3, protein.getCodonFrames().size());
 664     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(0)).size());
 665     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(1)).size());
 666     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(2)).size());
 667
 668     // V12345 mapped from A22222
 669     AlignedCodonFrame acf = protein.getCodonFrame(protein.getSequenceAt(0))
 670             .get(0);
 671     assertEquals(1, acf.getdnaSeqs().length);
 672     assertEquals(cdna.getSequenceAt(1).getDatasetSequence(),
 673             acf.getdnaSeqs()[0]);
 674     Mapping[] protMappings = acf.getProtMappings();
 675     assertEquals(1, protMappings.length);
 676     MapList mapList = protMappings[0].getMap();
 677     assertEquals(3, mapList.getFromRatio());
 678     assertEquals(1, mapList.getToRatio());
 679     assertTrue(
 680             Arrays.equals(new int[]
 681             { 1, 9 }, mapList.getFromRanges().get(0)));
 682     assertEquals(1, mapList.getFromRanges().size());
 683     assertTrue(
 684             Arrays.equals(new int[]
 685             { 1, 3 }, mapList.getToRanges().get(0)));
 686     assertEquals(1, mapList.getToRanges().size());
 687
 688     // V12346 mapped from A33333 starting position 4
 689     acf = protein.getCodonFrame(protein.getSequenceAt(1)).get(0);
 690     assertEquals(1, acf.getdnaSeqs().length);
 691     assertEquals(cdna.getSequenceAt(2).getDatasetSequence(),
 692             acf.getdnaSeqs()[0]);
 693     protMappings = acf.getProtMappings();
 694     assertEquals(1, protMappings.length);
 695     mapList = protMappings[0].getMap();
 696     assertEquals(3, mapList.getFromRatio());
 697     assertEquals(1, mapList.getToRatio());
 698     assertTrue(
 699             Arrays.equals(new int[]
 700             { 4, 12 }, mapList.getFromRanges().get(0)));
 701     assertEquals(1, mapList.getFromRanges().size());
 702     assertTrue(
 703             Arrays.equals(new int[]
 704             { 1, 3 }, mapList.getToRanges().get(0)));
 705     assertEquals(1, mapList.getToRanges().size());
 706
 707     // V12347 mapped to A11111 starting position 4
 708     acf = protein.getCodonFrame(protein.getSequenceAt(2)).get(0);
 709     assertEquals(1, acf.getdnaSeqs().length);
 710     assertEquals(cdna.getSequenceAt(0).getDatasetSequence(),
 711             acf.getdnaSeqs()[0]);
 712     protMappings = acf.getProtMappings();
 713     assertEquals(1, protMappings.length);
 714     mapList = protMappings[0].getMap();
 715     assertEquals(3, mapList.getFromRatio());
 716     assertEquals(1, mapList.getToRatio());
 717     assertTrue(
 718             Arrays.equals(new int[]
 719             { 4, 12 }, mapList.getFromRanges().get(0)));
 720     assertEquals(1, mapList.getFromRanges().size());
 721     assertTrue(
 722             Arrays.equals(new int[]
 723             { 1, 3 }, mapList.getToRanges().get(0)));
 724     assertEquals(1, mapList.getToRanges().size());
 725
 726     // no mapping involving the 'extra' A44444
 727     assertTrue(protein.getCodonFrame(cdna.getSequenceAt(3)).isEmpty());
 728   }
 729
 730   /**
 731    * Test mapping of protein to cDNA, for the case where we have some sequence
 732    * cross-references. Verify that 1-to-many mappings are made where
 733    * cross-references exist and sequences are mappable.
 734    *
 735    * @throws IOException
 736    */
 737   @Test(groups = { "Functional" })
 738   public void testMapProteinAlignmentToCdna_withXrefs() throws IOException
 739   {
 740     List<SequenceI> protseqs = new ArrayList<>();
 741     protseqs.add(new Sequence("UNIPROT|V12345", "EIQ"));
 742     protseqs.add(new Sequence("UNIPROT|V12346", "EIQ"));
 743     protseqs.add(new Sequence("UNIPROT|V12347", "SAR"));
 744     AlignmentI protein = new Alignment(protseqs.toArray(new SequenceI[3]));
 745     protein.setDataset(null);
 746
 747     List<SequenceI> dnaseqs = new ArrayList<>();
 748     dnaseqs.add(new Sequence("EMBL|A11111", "TCAGCACGC")); // = SAR
 749     dnaseqs.add(new Sequence("EMBL|A22222", "ATGGAGATACAA")); // = start + EIQ
 750     dnaseqs.add(new Sequence("EMBL|A33333", "GAAATCCAG")); // = EIQ
 751     dnaseqs.add(new Sequence("EMBL|A44444", "GAAATTCAG")); // = EIQ
 752     dnaseqs.add(new Sequence("EMBL|A55555", "GAGATTCAG")); // = EIQ
 753     AlignmentI cdna = new Alignment(dnaseqs.toArray(new SequenceI[5]));
 754     cdna.setDataset(null);
 755
 756     // Xref A22222 to V12345 (should get mapped)
 757     dnaseqs.get(1).addDBRef(new DBRefEntry("UNIPROT", "1", "V12345"));
 758     // Xref V12345 to A44444 (should get mapped)
 759     protseqs.get(0).addDBRef(new DBRefEntry("EMBL", "1", "A44444"));
 760     // Xref A33333 to V12347 (sequence mismatch - should not get mapped)
 761     dnaseqs.get(2).addDBRef(new DBRefEntry("UNIPROT", "1", "V12347"));
 762     // as V12345 is mapped to A22222 and A44444, this leaves V12346 unmapped.
 763     // it should get paired up with the unmapped A33333
 764     // A11111 should be mapped to V12347
 765     // A55555 is spare and has no xref so is not mapped
 766
 767     assertTrue(AlignmentUtils.mapProteinAlignmentToCdna(protein, cdna));
 768
 769     // 4 protein mappings made for 3 proteins, 2 to V12345, 1 each to V12346/7
 770     assertEquals(3, protein.getCodonFrames().size());
 771     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(0)).size());
 772     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(1)).size());
 773     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(2)).size());
 774
 775     // one mapping for each of the first 4 cDNA sequences
 776     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(0)).size());
 777     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(1)).size());
 778     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(2)).size());
 779     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(3)).size());
 780
 781     // V12345 mapped to A22222 and A44444
 782     AlignedCodonFrame acf = protein.getCodonFrame(protein.getSequenceAt(0))
 783             .get(0);
 784     assertEquals(2, acf.getdnaSeqs().length);
 785     assertEquals(cdna.getSequenceAt(1).getDatasetSequence(),
 786             acf.getdnaSeqs()[0]);
 787     assertEquals(cdna.getSequenceAt(3).getDatasetSequence(),
 788             acf.getdnaSeqs()[1]);
 789
 790     // V12346 mapped to A33333
 791     acf = protein.getCodonFrame(protein.getSequenceAt(1)).get(0);
 792     assertEquals(1, acf.getdnaSeqs().length);
 793     assertEquals(cdna.getSequenceAt(2).getDatasetSequence(),
 794             acf.getdnaSeqs()[0]);
 795
 796     // V12347 mapped to A11111
 797     acf = protein.getCodonFrame(protein.getSequenceAt(2)).get(0);
 798     assertEquals(1, acf.getdnaSeqs().length);
 799     assertEquals(cdna.getSequenceAt(0).getDatasetSequence(),
 800             acf.getdnaSeqs()[0]);
 801
 802     // no mapping involving the 'extra' A55555
 803     assertTrue(protein.getCodonFrame(cdna.getSequenceAt(4)).isEmpty());
 804   }
 805
 806   /**
 807    * Test mapping of protein to cDNA, for the case where we have some sequence
 808    * cross-references. Verify that once we have made an xref mapping we don't
 809    * also map un-xrefd sequeces.
 810    *
 811    * @throws IOException
 812    */
 813   @Test(groups = { "Functional" })
 814   public void testMapProteinAlignmentToCdna_prioritiseXrefs()
 815           throws IOException
 816   {
 817     List<SequenceI> protseqs = new ArrayList<>();
 818     protseqs.add(new Sequence("UNIPROT|V12345", "EIQ"));
 819     protseqs.add(new Sequence("UNIPROT|V12346", "EIQ"));
 820     AlignmentI protein = new Alignment(
 821             protseqs.toArray(new SequenceI[protseqs.size()]));
 822     protein.setDataset(null);
 823
 824     List<SequenceI> dnaseqs = new ArrayList<>();
 825     dnaseqs.add(new Sequence("EMBL|A11111", "GAAATCCAG")); // = EIQ
 826     dnaseqs.add(new Sequence("EMBL|A22222", "GAAATTCAG")); // = EIQ
 827     AlignmentI cdna = new Alignment(
 828             dnaseqs.toArray(new SequenceI[dnaseqs.size()]));
 829     cdna.setDataset(null);
 830
 831     // Xref A22222 to V12345 (should get mapped)
 832     // A11111 should then be mapped to the unmapped V12346
 833     dnaseqs.get(1).addDBRef(new DBRefEntry("UNIPROT", "1", "V12345"));
 834
 835     assertTrue(AlignmentUtils.mapProteinAlignmentToCdna(protein, cdna));
 836
 837     // 2 protein mappings made
 838     assertEquals(2, protein.getCodonFrames().size());
 839     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(0)).size());
 840     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(1)).size());
 841
 842     // one mapping for each of the cDNA sequences
 843     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(0)).size());
 844     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(1)).size());
 845
 846     // V12345 mapped to A22222
 847     AlignedCodonFrame acf = protein.getCodonFrame(protein.getSequenceAt(0))
 848             .get(0);
 849     assertEquals(1, acf.getdnaSeqs().length);
 850     assertEquals(cdna.getSequenceAt(1).getDatasetSequence(),
 851             acf.getdnaSeqs()[0]);
 852
 853     // V12346 mapped to A11111
 854     acf = protein.getCodonFrame(protein.getSequenceAt(1)).get(0);
 855     assertEquals(1, acf.getdnaSeqs().length);
 856     assertEquals(cdna.getSequenceAt(0).getDatasetSequence(),
 857             acf.getdnaSeqs()[0]);
 858   }
 859
 860   /**
 861    * Test the method that shows or hides sequence annotations by type(s) and
 862    * selection group.
 863    */
 864   @Test(groups = { "Functional" })
 865   public void testShowOrHideSequenceAnnotations()
 866   {
 867     SequenceI seq1 = new Sequence("Seq1", "AAA");
 868     SequenceI seq2 = new Sequence("Seq2", "BBB");
 869     SequenceI seq3 = new Sequence("Seq3", "CCC");
 870     Annotation[] anns = new Annotation[] { new Annotation(2f) };
 871     AlignmentAnnotation ann1 = new AlignmentAnnotation("Structure", "ann1",
 872             anns);
 873     ann1.setSequenceRef(seq1);
 874     AlignmentAnnotation ann2 = new AlignmentAnnotation("Structure", "ann2",
 875             anns);
 876     ann2.setSequenceRef(seq2);
 877     AlignmentAnnotation ann3 = new AlignmentAnnotation("Structure", "ann3",
 878             anns);
 879     AlignmentAnnotation ann4 = new AlignmentAnnotation("Temp", "ann4",
 880             anns);
 881     ann4.setSequenceRef(seq1);
 882     AlignmentAnnotation ann5 = new AlignmentAnnotation("Temp", "ann5",
 883             anns);
 884     ann5.setSequenceRef(seq2);
 885     AlignmentAnnotation ann6 = new AlignmentAnnotation("Temp", "ann6",
 886             anns);
 887     AlignmentI al = new Alignment(new SequenceI[] { seq1, seq2, seq3 });
 888     al.addAnnotation(ann1); // Structure for Seq1
 889     al.addAnnotation(ann2); // Structure for Seq2
 890     al.addAnnotation(ann3); // Structure for no sequence
 891     al.addAnnotation(ann4); // Temp for seq1
 892     al.addAnnotation(ann5); // Temp for seq2
 893     al.addAnnotation(ann6); // Temp for no sequence
 894     List<String> types = new ArrayList<>();
 895     List<SequenceI> scope = new ArrayList<>();
 896
 897     /*
 898      * Set all sequence related Structure to hidden (ann1, ann2)
 899      */
 900     types.add("Structure");
 901     AlignmentUtils.showOrHideSequenceAnnotations(al, types, null, false,
 902             false);
 903     assertFalse(ann1.visible);
 904     assertFalse(ann2.visible);
 905     assertTrue(ann3.visible); // not sequence-related, not affected
 906     assertTrue(ann4.visible); // not Structure, not affected
 907     assertTrue(ann5.visible); // "
 908     assertTrue(ann6.visible); // not sequence-related, not affected
 909
 910     /*
 911      * Set Temp in {seq1, seq3} to hidden
 912      */
 913     types.clear();
 914     types.add("Temp");
 915     scope.add(seq1);
 916     scope.add(seq3);
 917     AlignmentUtils.showOrHideSequenceAnnotations(al, types, scope, false,
 918             false);
 919     assertFalse(ann1.visible); // unchanged
 920     assertFalse(ann2.visible); // unchanged
 921     assertTrue(ann3.visible); // not sequence-related, not affected
 922     assertFalse(ann4.visible); // Temp for seq1 hidden
 923     assertTrue(ann5.visible); // not in scope, not affected
 924     assertTrue(ann6.visible); // not sequence-related, not affected
 925
 926     /*
 927      * Set Temp in all sequences to hidden
 928      */
 929     types.clear();
 930     types.add("Temp");
 931     scope.add(seq1);
 932     scope.add(seq3);
 933     AlignmentUtils.showOrHideSequenceAnnotations(al, types, null, false,
 934             false);
 935     assertFalse(ann1.visible); // unchanged
 936     assertFalse(ann2.visible); // unchanged
 937     assertTrue(ann3.visible); // not sequence-related, not affected
 938     assertFalse(ann4.visible); // Temp for seq1 hidden
 939     assertFalse(ann5.visible); // Temp for seq2 hidden
 940     assertTrue(ann6.visible); // not sequence-related, not affected
 941
 942     /*
 943      * Set all types in {seq1, seq3} to visible
 944      */
 945     types.clear();
 946     scope.clear();
 947     scope.add(seq1);
 948     scope.add(seq3);
 949     AlignmentUtils.showOrHideSequenceAnnotations(al, types, scope, true,
 950             true);
 951     assertTrue(ann1.visible); // Structure for seq1 set visible
 952     assertFalse(ann2.visible); // not in scope, unchanged
 953     assertTrue(ann3.visible); // not sequence-related, not affected
 954     assertTrue(ann4.visible); // Temp for seq1 set visible
 955     assertFalse(ann5.visible); // not in scope, unchanged
 956     assertTrue(ann6.visible); // not sequence-related, not affected
 957
 958     /*
 959      * Set all types in all scope to hidden
 960      */
 961     AlignmentUtils.showOrHideSequenceAnnotations(al, types, null, true,
 962             false);
 963     assertFalse(ann1.visible);
 964     assertFalse(ann2.visible);
 965     assertTrue(ann3.visible); // not sequence-related, not affected
 966     assertFalse(ann4.visible);
 967     assertFalse(ann5.visible);
 968     assertTrue(ann6.visible); // not sequence-related, not affected
 969   }
 970
 971   /**
 972    * Tests for the method that checks if one sequence cross-references another
 973    */
 974   @Test(groups = { "Functional" })
 975   public void testHasCrossRef()
 976   {
 977     assertFalse(AlignmentUtils.hasCrossRef(null, null));
 978     SequenceI seq1 = new Sequence("EMBL|A12345", "ABCDEF");
 979     assertFalse(AlignmentUtils.hasCrossRef(seq1, null));
 980     assertFalse(AlignmentUtils.hasCrossRef(null, seq1));
 981     SequenceI seq2 = new Sequence("UNIPROT|V20192", "ABCDEF");
 982     assertFalse(AlignmentUtils.hasCrossRef(seq1, seq2));
 983
 984     // different ref
 985     seq1.addDBRef(new DBRefEntry("UNIPROT", "1", "v20193"));
 986     assertFalse(AlignmentUtils.hasCrossRef(seq1, seq2));
 987
 988     // case-insensitive; version number is ignored
 989     seq1.addDBRef(new DBRefEntry("UNIPROT", "1", "v20192"));
 990     assertTrue(AlignmentUtils.hasCrossRef(seq1, seq2));
 991
 992     // right case!
 993     seq1.addDBRef(new DBRefEntry("UNIPROT", "1", "V20192"));
 994     assertTrue(AlignmentUtils.hasCrossRef(seq1, seq2));
 995     // test is one-way only
 996     assertFalse(AlignmentUtils.hasCrossRef(seq2, seq1));
 997   }
 998
 999   /**
1000    * Tests for the method that checks if either sequence cross-references the
1001    * other
1002    */
1003   @Test(groups = { "Functional" })
1004   public void testHaveCrossRef()
1005   {
1006     assertFalse(AlignmentUtils.hasCrossRef(null, null));
1007     SequenceI seq1 = new Sequence("EMBL|A12345", "ABCDEF");
1008     assertFalse(AlignmentUtils.haveCrossRef(seq1, null));
1009     assertFalse(AlignmentUtils.haveCrossRef(null, seq1));
1010     SequenceI seq2 = new Sequence("UNIPROT|V20192", "ABCDEF");
1011     assertFalse(AlignmentUtils.haveCrossRef(seq1, seq2));
1012
1013     seq1.addDBRef(new DBRefEntry("UNIPROT", "1", "V20192"));
1014     assertTrue(AlignmentUtils.haveCrossRef(seq1, seq2));
1015     // next is true for haveCrossRef, false for hasCrossRef
1016     assertTrue(AlignmentUtils.haveCrossRef(seq2, seq1));
1017
1018     // now the other way round
1019     seq1.setDBRefs(null);
1020     seq2.addDBRef(new DBRefEntry("EMBL", "1", "A12345"));
1021     assertTrue(AlignmentUtils.haveCrossRef(seq1, seq2));
1022     assertTrue(AlignmentUtils.haveCrossRef(seq2, seq1));
1023
1024     // now both ways
1025     seq1.addDBRef(new DBRefEntry("UNIPROT", "1", "V20192"));
1026     assertTrue(AlignmentUtils.haveCrossRef(seq1, seq2));
1027     assertTrue(AlignmentUtils.haveCrossRef(seq2, seq1));
1028   }
1029
1030   /**
1031    * Test the method that extracts the cds-only part of a dna alignment.
1032    */
1033   @Test(groups = { "Functional" })
1034   public void testMakeCdsAlignment()
1035   {
1036     /*
1037      * scenario:
1038      *     dna1 --> [4, 6] [10,12]        --> pep1
1039      *     dna2 --> [1, 3] [7, 9] [13,15] --> pep2
1040      */
1041     SequenceI dna1 = new Sequence("dna1", "aaaGGGcccTTTaaa");
1042     SequenceI dna2 = new Sequence("dna2", "GGGcccTTTaaaCCC");
1043     SequenceI pep1 = new Sequence("pep1", "GF");
1044     SequenceI pep2 = new Sequence("pep2", "GFP");
1045     pep1.addDBRef(new DBRefEntry("UNIPROT", "0", "pep1"));
1046     pep2.addDBRef(new DBRefEntry("UNIPROT", "0", "pep2"));
1047     dna1.createDatasetSequence();
1048     dna2.createDatasetSequence();
1049     pep1.createDatasetSequence();
1050     pep2.createDatasetSequence();
1051     AlignmentI dna = new Alignment(new SequenceI[] { dna1, dna2 });
1052     dna.setDataset(null);
1053
1054     /*
1055      * put a variant feature on dna2 base 8
1056      * - should transfer to cds2 base 5
1057      */
1058     dna2.addSequenceFeature(
1059             new SequenceFeature("variant", "hgmd", 8, 8, 0f, null));
1060
1061     /*
1062      * need a sourceDbRef if we are to construct dbrefs to the CDS
1063      * sequence from the dna contig sequences
1064      */
1065     DBRefEntry dbref = new DBRefEntry("ENSEMBL", "0", "dna1");
1066     dna1.getDatasetSequence().addDBRef(dbref);
1067     org.testng.Assert.assertEquals(dbref, dna1.getPrimaryDBRefs().get(0));
1068     dbref = new DBRefEntry("ENSEMBL", "0", "dna2");
1069     dna2.getDatasetSequence().addDBRef(dbref);
1070     org.testng.Assert.assertEquals(dbref, dna2.getPrimaryDBRefs().get(0));
1071
1072     /*
1073      * CDS sequences are 'discovered' from dna-to-protein mappings on the alignment
1074      * dataset (e.g. added from dbrefs by CrossRef.findXrefSequences)
1075      */
1076     MapList mapfordna1 = new MapList(new int[] { 4, 6, 10, 12 },
1077             new int[]
1078             { 1, 2 }, 3, 1);
1079     AlignedCodonFrame acf = new AlignedCodonFrame();
1080     acf.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(),
1081             mapfordna1);
1082     dna.addCodonFrame(acf);
1083     MapList mapfordna2 = new MapList(new int[] { 1, 3, 7, 9, 13, 15 },
1084             new int[]
1085             { 1, 3 }, 3, 1);
1086     acf = new AlignedCodonFrame();
1087     acf.addMap(dna2.getDatasetSequence(), pep2.getDatasetSequence(),
1088             mapfordna2);
1089     dna.addCodonFrame(acf);
1090
1091     /*
1092      * In this case, mappings originally came from matching Uniprot accessions
1093      * - so need an xref on dna involving those regions.
1094      * These are normally constructed from CDS annotation
1095      */
1096     DBRefEntry dna1xref = new DBRefEntry("UNIPROT", "ENSEMBL", "pep1",
1097             new Mapping(mapfordna1));
1098     dna1.addDBRef(dna1xref);
1099     assertEquals(2, dna1.getDBRefs().size()); // to self and to pep1
1100     DBRefEntry dna2xref = new DBRefEntry("UNIPROT", "ENSEMBL", "pep2",
1101             new Mapping(mapfordna2));
1102     dna2.addDBRef(dna2xref);
1103     assertEquals(2, dna2.getDBRefs().size()); // to self and to pep2
1104
1105     /*
1106      * execute method under test:
1107      */
1108     AlignmentI cds = AlignmentUtils
1109             .makeCdsAlignment(new SequenceI[]
1110             { dna1, dna2 }, dna.getDataset(), null);
1111
1112     /*
1113      * verify cds sequences
1114      */
1115     assertEquals(2, cds.getSequences().size());
1116     assertEquals("GGGTTT", cds.getSequenceAt(0).getSequenceAsString());
1117     assertEquals("GGGTTTCCC", cds.getSequenceAt(1).getSequenceAsString());
1118
1119     /*
1120      * verify shared, extended alignment dataset
1121      */
1122     assertSame(dna.getDataset(), cds.getDataset());
1123     SequenceI cds1Dss = cds.getSequenceAt(0).getDatasetSequence();
1124     SequenceI cds2Dss = cds.getSequenceAt(1).getDatasetSequence();
1125     assertTrue(dna.getDataset().getSequences().contains(cds1Dss));
1126     assertTrue(dna.getDataset().getSequences().contains(cds2Dss));
1127
1128     /*
1129      * verify CDS has a dbref with mapping to peptide
1130      */
1131     assertNotNull(cds1Dss.getDBRefs());
1132     assertEquals(2, cds1Dss.getDBRefs().size());
1133     dbref = cds1Dss.getDBRefs().get(0);
1134     assertEquals(dna1xref.getSource(), dbref.getSource());
1135     // version is via ensembl's primary ref
1136     assertEquals(dna1xref.getVersion(), dbref.getVersion());
1137     assertEquals(dna1xref.getAccessionId(), dbref.getAccessionId());
1138     assertNotNull(dbref.getMap());
1139     assertSame(pep1.getDatasetSequence(), dbref.getMap().getTo());
1140     MapList cdsMapping = new MapList(new int[] { 1, 6 }, new int[] { 1, 2 },
1141             3, 1);
1142     assertEquals(cdsMapping, dbref.getMap().getMap());
1143
1144     /*
1145      * verify peptide has added a dbref with reverse mapping to CDS
1146      */
1147     assertNotNull(pep1.getDBRefs());
1148     // FIXME pep1.getDBRefs() is 1 - is that the correct behaviour ?
1149     assertEquals(2, pep1.getDBRefs().size());
1150     dbref = pep1.getDBRefs().get(1);
1151     assertEquals("ENSEMBL", dbref.getSource());
1152     assertEquals("0", dbref.getVersion());
1153     assertEquals("CDS|dna1", dbref.getAccessionId());
1154     assertNotNull(dbref.getMap());
1155     assertSame(cds1Dss, dbref.getMap().getTo());
1156     assertEquals(cdsMapping.getInverse(), dbref.getMap().getMap());
1157
1158     /*
1159      * verify cDNA has added a dbref with mapping to CDS
1160      */
1161     assertEquals(3, dna1.getDBRefs().size());
1162     DBRefEntry dbRefEntry = dna1.getDBRefs().get(2);
1163     assertSame(cds1Dss, dbRefEntry.getMap().getTo());
1164     MapList dnaToCdsMapping = new MapList(new int[] { 4, 6, 10, 12 },
1165             new int[]
1166             { 1, 6 }, 1, 1);
1167     assertEquals(dnaToCdsMapping, dbRefEntry.getMap().getMap());
1168     assertEquals(3, dna2.getDBRefs().size());
1169     dbRefEntry = dna2.getDBRefs().get(2);
1170     assertSame(cds2Dss, dbRefEntry.getMap().getTo());
1171     dnaToCdsMapping = new MapList(new int[] { 1, 3, 7, 9, 13, 15 },
1172             new int[]
1173             { 1, 9 }, 1, 1);
1174     assertEquals(dnaToCdsMapping, dbRefEntry.getMap().getMap());
1175
1176     /*
1177      * verify CDS has added a dbref with mapping to cDNA
1178      */
1179     assertEquals(2, cds1Dss.getDBRefs().size());
1180     dbRefEntry = cds1Dss.getDBRefs().get(1);
1181     assertSame(dna1.getDatasetSequence(), dbRefEntry.getMap().getTo());
1182     MapList cdsToDnaMapping = new MapList(new int[] { 1, 6 },
1183             new int[]
1184             { 4, 6, 10, 12 }, 1, 1);
1185     assertEquals(cdsToDnaMapping, dbRefEntry.getMap().getMap());
1186     assertEquals(2, cds2Dss.getDBRefs().size());
1187     dbRefEntry = cds2Dss.getDBRefs().get(1);
1188     assertSame(dna2.getDatasetSequence(), dbRefEntry.getMap().getTo());
1189     cdsToDnaMapping = new MapList(new int[] { 1, 9 },
1190             new int[]
1191             { 1, 3, 7, 9, 13, 15 }, 1, 1);
1192     assertEquals(cdsToDnaMapping, dbRefEntry.getMap().getMap());
1193
1194     /*
1195      * Verify mappings from CDS to peptide, cDNA to CDS, and cDNA to peptide
1196      * the mappings are on the shared alignment dataset
1197      * 6 mappings, 2*(DNA->CDS), 2*(DNA->Pep), 2*(CDS->Pep)
1198      */
1199     List<AlignedCodonFrame> cdsMappings = cds.getDataset().getCodonFrames();
1200     assertEquals(6, cdsMappings.size());
1201
1202     /*
1203      * verify that mapping sets for dna and cds alignments are different
1204      * [not current behaviour - all mappings are on the alignment dataset]
1205      */
1206     // select -> subselect type to test.
1207     // Assert.assertNotSame(dna.getCodonFrames(), cds.getCodonFrames());
1208     // assertEquals(4, dna.getCodonFrames().size());
1209     // assertEquals(4, cds.getCodonFrames().size());
1210
1211     /*
1212      * Two mappings involve pep1 (dna to pep1, cds to pep1)
1213      * Mapping from pep1 to GGGTTT in first new exon sequence
1214      */
1215     List<AlignedCodonFrame> pep1Mappings = MappingUtils
1216             .findMappingsForSequence(pep1, cdsMappings);
1217     assertEquals(2, pep1Mappings.size());
1218     List<AlignedCodonFrame> mappings = MappingUtils
1219             .findMappingsForSequence(cds.getSequenceAt(0), pep1Mappings);
1220     assertEquals(1, mappings.size());
1221
1222     // map G to GGG
1223     SearchResultsI sr = MappingUtils.buildSearchResults(pep1, 1, mappings);
1224     assertEquals(1, sr.getResults().size());
1225     SearchResultMatchI m = sr.getResults().get(0);
1226     assertSame(cds1Dss, m.getSequence());
1227     assertEquals(1, m.getStart());
1228     assertEquals(3, m.getEnd());
1229     // map F to TTT
1230     sr = MappingUtils.buildSearchResults(pep1, 2, mappings);
1231     m = sr.getResults().get(0);
1232     assertSame(cds1Dss, m.getSequence());
1233     assertEquals(4, m.getStart());
1234     assertEquals(6, m.getEnd());
1235
1236     /*
1237      * Two mappings involve pep2 (dna to pep2, cds to pep2)
1238      * Verify mapping from pep2 to GGGTTTCCC in second new exon sequence
1239      */
1240     List<AlignedCodonFrame> pep2Mappings = MappingUtils
1241             .findMappingsForSequence(pep2, cdsMappings);
1242     assertEquals(2, pep2Mappings.size());
1243     mappings = MappingUtils.findMappingsForSequence(cds.getSequenceAt(1),
1244             pep2Mappings);
1245     assertEquals(1, mappings.size());
1246     // map G to GGG
1247     sr = MappingUtils.buildSearchResults(pep2, 1, mappings);
1248     assertEquals(1, sr.getResults().size());
1249     m = sr.getResults().get(0);
1250     assertSame(cds2Dss, m.getSequence());
1251     assertEquals(1, m.getStart());
1252     assertEquals(3, m.getEnd());
1253     // map F to TTT
1254     sr = MappingUtils.buildSearchResults(pep2, 2, mappings);
1255     m = sr.getResults().get(0);
1256     assertSame(cds2Dss, m.getSequence());
1257     assertEquals(4, m.getStart());
1258     assertEquals(6, m.getEnd());
1259     // map P to CCC
1260     sr = MappingUtils.buildSearchResults(pep2, 3, mappings);
1261     m = sr.getResults().get(0);
1262     assertSame(cds2Dss, m.getSequence());
1263     assertEquals(7, m.getStart());
1264     assertEquals(9, m.getEnd());
1265
1266     /*
1267      * check cds2 acquired a variant feature in position 5
1268      */
1269     List<SequenceFeature> sfs = cds2Dss.getSequenceFeatures();
1270     assertNotNull(sfs);
1271     assertEquals(1, sfs.size());
1272     assertEquals("variant", sfs.get(0).type);
1273     assertEquals(5, sfs.get(0).begin);
1274     assertEquals(5, sfs.get(0).end);
1275   }
1276
1277   /**
1278    * Test the method that makes a cds-only alignment from a DNA sequence and its
1279    * product mappings, for the case where there are multiple exon mappings to
1280    * different protein products.
1281    */
1282   @Test(groups = { "Functional" })
1283   public void testMakeCdsAlignment_multipleProteins()
1284   {
1285     SequenceI dna1 = new Sequence("dna1", "aaaGGGcccTTTaaa");
1286     SequenceI pep1 = new Sequence("pep1", "GF"); // GGGTTT
1287     SequenceI pep2 = new Sequence("pep2", "KP"); // aaaccc
1288     SequenceI pep3 = new Sequence("pep3", "KF"); // aaaTTT
1289     dna1.createDatasetSequence();
1290     pep1.createDatasetSequence();
1291     pep2.createDatasetSequence();
1292     pep3.createDatasetSequence();
1293     pep1.getDatasetSequence()
1294             .addDBRef(new DBRefEntry("EMBLCDS", "2", "A12345"));
1295     pep2.getDatasetSequence()
1296             .addDBRef(new DBRefEntry("EMBLCDS", "3", "A12346"));
1297     pep3.getDatasetSequence()
1298             .addDBRef(new DBRefEntry("EMBLCDS", "4", "A12347"));
1299
1300     /*
1301      * Create the CDS alignment
1302      */
1303     AlignmentI dna = new Alignment(new SequenceI[] { dna1 });
1304     dna.setDataset(null);
1305
1306     /*
1307      * Make the mappings from dna to protein
1308      */
1309     // map ...GGG...TTT to GF
1310     MapList map = new MapList(new int[] { 4, 6, 10, 12 },
1311             new int[]
1312             { 1, 2 }, 3, 1);
1313     AlignedCodonFrame acf = new AlignedCodonFrame();
1314     acf.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(), map);
1315     dna.addCodonFrame(acf);
1316
1317     // map aaa...ccc to KP
1318     map = new MapList(new int[] { 1, 3, 7, 9 }, new int[] { 1, 2 }, 3, 1);
1319     acf = new AlignedCodonFrame();
1320     acf.addMap(dna1.getDatasetSequence(), pep2.getDatasetSequence(), map);
1321     dna.addCodonFrame(acf);
1322
1323     // map aaa......TTT to KF
1324     map = new MapList(new int[] { 1, 3, 10, 12 }, new int[] { 1, 2 }, 3, 1);
1325     acf = new AlignedCodonFrame();
1326     acf.addMap(dna1.getDatasetSequence(), pep3.getDatasetSequence(), map);
1327     dna.addCodonFrame(acf);
1328
1329     /*
1330      * execute method under test
1331      */
1332     AlignmentI cdsal = AlignmentUtils
1333             .makeCdsAlignment(new SequenceI[]
1334             { dna1 }, dna.getDataset(), null);
1335
1336     /*
1337      * Verify we have 3 cds sequences, mapped to pep1/2/3 respectively
1338      */
1339     List<SequenceI> cds = cdsal.getSequences();
1340     assertEquals(3, cds.size());
1341
1342     /*
1343      * verify shared, extended alignment dataset
1344      */
1345     assertSame(cdsal.getDataset(), dna.getDataset());
1346     assertTrue(dna.getDataset().getSequences()
1347             .contains(cds.get(0).getDatasetSequence()));
1348     assertTrue(dna.getDataset().getSequences()
1349             .contains(cds.get(1).getDatasetSequence()));
1350     assertTrue(dna.getDataset().getSequences()
1351             .contains(cds.get(2).getDatasetSequence()));
1352
1353     /*
1354      * verify aligned cds sequences and their xrefs
1355      */
1356     SequenceI cdsSeq = cds.get(0);
1357     assertEquals("GGGTTT", cdsSeq.getSequenceAsString());
1358     // assertEquals("dna1|A12345", cdsSeq.getName());
1359     assertEquals("CDS|dna1", cdsSeq.getName());
1360     // assertEquals(1, cdsSeq.getDBRefs().length);
1361     // DBRefEntry cdsRef = cdsSeq.getDBRefs()[0];
1362     // assertEquals("EMBLCDS", cdsRef.getSource());
1363     // assertEquals("2", cdsRef.getVersion());
1364     // assertEquals("A12345", cdsRef.getAccessionId());
1365
1366     cdsSeq = cds.get(1);
1367     assertEquals("aaaccc", cdsSeq.getSequenceAsString());
1368     // assertEquals("dna1|A12346", cdsSeq.getName());
1369     assertEquals("CDS|dna1", cdsSeq.getName());
1370     // assertEquals(1, cdsSeq.getDBRefs().length);
1371     // cdsRef = cdsSeq.getDBRefs()[0];
1372     // assertEquals("EMBLCDS", cdsRef.getSource());
1373     // assertEquals("3", cdsRef.getVersion());
1374     // assertEquals("A12346", cdsRef.getAccessionId());
1375
1376     cdsSeq = cds.get(2);
1377     assertEquals("aaaTTT", cdsSeq.getSequenceAsString());
1378     // assertEquals("dna1|A12347", cdsSeq.getName());
1379     assertEquals("CDS|dna1", cdsSeq.getName());
1380     // assertEquals(1, cdsSeq.getDBRefs().length);
1381     // cdsRef = cdsSeq.getDBRefs()[0];
1382     // assertEquals("EMBLCDS", cdsRef.getSource());
1383     // assertEquals("4", cdsRef.getVersion());
1384     // assertEquals("A12347", cdsRef.getAccessionId());
1385
1386     /*
1387      * Verify there are mappings from each cds sequence to its protein product
1388      * and also to its dna source
1389      */
1390     List<AlignedCodonFrame> newMappings = cdsal.getCodonFrames();
1391
1392     /*
1393      * 6 mappings involve dna1 (to pep1/2/3, cds1/2/3)
1394      */
1395     List<AlignedCodonFrame> dnaMappings = MappingUtils
1396             .findMappingsForSequence(dna1, newMappings);
1397     assertEquals(6, dnaMappings.size());
1398
1399     /*
1400      * dna1 to pep1
1401      */
1402     List<AlignedCodonFrame> mappings = MappingUtils
1403             .findMappingsForSequence(pep1, dnaMappings);
1404     assertEquals(1, mappings.size());
1405     assertEquals(1, mappings.get(0).getMappings().size());
1406     assertSame(pep1.getDatasetSequence(),
1407             mappings.get(0).getMappings().get(0).getMapping().getTo());
1408
1409     /*
1410      * dna1 to cds1
1411      */
1412     List<AlignedCodonFrame> dnaToCds1Mappings = MappingUtils
1413             .findMappingsForSequence(cds.get(0), dnaMappings);
1414     Mapping mapping = dnaToCds1Mappings.get(0).getMappings().get(0)
1415             .getMapping();
1416     assertSame(cds.get(0).getDatasetSequence(), mapping.getTo());
1417     assertEquals("G(1) in CDS should map to G(4) in DNA", 4,
1418             mapping.getMap().getToPosition(1));
1419
1420     /*
1421      * dna1 to pep2
1422      */
1423     mappings = MappingUtils.findMappingsForSequence(pep2, dnaMappings);
1424     assertEquals(1, mappings.size());
1425     assertEquals(1, mappings.get(0).getMappings().size());
1426     assertSame(pep2.getDatasetSequence(),
1427             mappings.get(0).getMappings().get(0).getMapping().getTo());
1428
1429     /*
1430      * dna1 to cds2
1431      */
1432     List<AlignedCodonFrame> dnaToCds2Mappings = MappingUtils
1433             .findMappingsForSequence(cds.get(1), dnaMappings);
1434     mapping = dnaToCds2Mappings.get(0).getMappings().get(0).getMapping();
1435     assertSame(cds.get(1).getDatasetSequence(), mapping.getTo());
1436     assertEquals("c(4) in CDS should map to c(7) in DNA", 7,
1437             mapping.getMap().getToPosition(4));
1438
1439     /*
1440      * dna1 to pep3
1441      */
1442     mappings = MappingUtils.findMappingsForSequence(pep3, dnaMappings);
1443     assertEquals(1, mappings.size());
1444     assertEquals(1, mappings.get(0).getMappings().size());
1445     assertSame(pep3.getDatasetSequence(),
1446             mappings.get(0).getMappings().get(0).getMapping().getTo());
1447
1448     /*
1449      * dna1 to cds3
1450      */
1451     List<AlignedCodonFrame> dnaToCds3Mappings = MappingUtils
1452             .findMappingsForSequence(cds.get(2), dnaMappings);
1453     mapping = dnaToCds3Mappings.get(0).getMappings().get(0).getMapping();
1454     assertSame(cds.get(2).getDatasetSequence(), mapping.getTo());
1455     assertEquals("T(4) in CDS should map to T(10) in DNA", 10,
1456             mapping.getMap().getToPosition(4));
1457   }
1458
1459   @Test(groups = { "Functional" })
1460   public void testIsMappable()
1461   {
1462     SequenceI dna1 = new Sequence("dna1", "cgCAGtgGT");
1463     SequenceI aa1 = new Sequence("aa1", "RSG");
1464     AlignmentI al1 = new Alignment(new SequenceI[] { dna1 });
1465     AlignmentI al2 = new Alignment(new SequenceI[] { aa1 });
1466
1467     assertFalse(AlignmentUtils.isMappable(null, null));
1468     assertFalse(AlignmentUtils.isMappable(al1, null));
1469     assertFalse(AlignmentUtils.isMappable(null, al1));
1470     assertFalse(AlignmentUtils.isMappable(al1, al1));
1471     assertFalse(AlignmentUtils.isMappable(al2, al2));
1472
1473     assertTrue(AlignmentUtils.isMappable(al1, al2));
1474     assertTrue(AlignmentUtils.isMappable(al2, al1));
1475   }
1476
1477   /**
1478    * Test creating a mapping when the sequences involved do not start at residue
1479    * 1
1480    *
1481    * @throws IOException
1482    */
1483   @Test(groups = { "Functional" })
1484   public void testMapCdnaToProtein_forSubsequence() throws IOException
1485   {
1486     SequenceI prot = new Sequence("UNIPROT|V12345", "E-I--Q", 10, 12);
1487     prot.createDatasetSequence();
1488
1489     SequenceI dna = new Sequence("EMBL|A33333", "GAA--AT-C-CAG", 40, 48);
1490     dna.createDatasetSequence();
1491
1492     MapList map = AlignmentUtils.mapCdnaToProtein(prot, dna);
1493     assertEquals(10, map.getToLowest());
1494     assertEquals(12, map.getToHighest());
1495     assertEquals(40, map.getFromLowest());
1496     assertEquals(48, map.getFromHighest());
1497   }
1498
1499   /**
1500    * Test for the alignSequenceAs method where we have protein mapped to protein
1501    */
1502   @Test(groups = { "Functional" })
1503   public void testAlignSequenceAs_mappedProteinProtein()
1504   {
1505
1506     SequenceI alignMe = new Sequence("Match", "MGAASEV");
1507     alignMe.createDatasetSequence();
1508     SequenceI alignFrom = new Sequence("Query", "LQTGYMGAASEVMFSPTRR");
1509     alignFrom.createDatasetSequence();
1510
1511     AlignedCodonFrame acf = new AlignedCodonFrame();
1512     // this is like a domain or motif match of part of a peptide sequence
1513     MapList map = new MapList(new int[] { 6, 12 }, new int[] { 1, 7 }, 1,
1514             1);
1515     acf.addMap(alignFrom.getDatasetSequence(), alignMe.getDatasetSequence(),
1516             map);
1517
1518     AlignmentUtils.alignSequenceAs(alignMe, alignFrom, acf, "-", '-', true,
1519             true);
1520     assertEquals("-----MGAASEV-------", alignMe.getSequenceAsString());
1521   }
1522
1523   /**
1524    * Test for the alignSequenceAs method where there are trailing unmapped
1525    * residues in the model sequence
1526    */
1527   @Test(groups = { "Functional" })
1528   public void testAlignSequenceAs_withTrailingPeptide()
1529   {
1530     // map first 3 codons to KPF; G is a trailing unmapped residue
1531     MapList map = new MapList(new int[] { 1, 9 }, new int[] { 1, 3 }, 3, 1);
1532
1533     checkAlignSequenceAs("AAACCCTTT", "K-PFG", true, true, map,
1534             "AAA---CCCTTT---");
1535   }
1536
1537   /**
1538    * Tests for transferring features between mapped sequences
1539    */
1540   @Test(groups = { "Functional" })
1541   public void testTransferFeatures()
1542   {
1543     SequenceI dna = new Sequence("dna/20-34", "acgTAGcaaGCCcgt");
1544     SequenceI cds = new Sequence("cds/10-15", "TAGGCC");
1545
1546     // no overlap
1547     dna.addSequenceFeature(
1548             new SequenceFeature("type1", "desc1", 1, 2, 1f, null));
1549     // partial overlap - to [1, 1]
1550     dna.addSequenceFeature(
1551             new SequenceFeature("type2", "desc2", 3, 4, 2f, null));
1552     // exact overlap - to [1, 3]
1553     dna.addSequenceFeature(
1554             new SequenceFeature("type3", "desc3", 4, 6, 3f, null));
1555     // spanning overlap - to [2, 5]
1556     dna.addSequenceFeature(
1557             new SequenceFeature("type4", "desc4", 5, 11, 4f, null));
1558     // exactly overlaps whole mapped range [1, 6]
1559     dna.addSequenceFeature(
1560             new SequenceFeature("type5", "desc5", 4, 12, 5f, null));
1561     // no overlap (internal)
1562     dna.addSequenceFeature(
1563             new SequenceFeature("type6", "desc6", 7, 9, 6f, null));
1564     // no overlap (3' end)
1565     dna.addSequenceFeature(
1566             new SequenceFeature("type7", "desc7", 13, 15, 7f, null));
1567     // overlap (3' end) - to [6, 6]
1568     dna.addSequenceFeature(
1569             new SequenceFeature("type8", "desc8", 12, 12, 8f, null));
1570     // extended overlap - to [6, +]
1571     dna.addSequenceFeature(
1572             new SequenceFeature("type9", "desc9", 12, 13, 9f, null));
1573
1574     MapList map = new MapList(new int[] { 4, 6, 10, 12 },
1575             new int[]
1576             { 1, 6 }, 1, 1);
1577
1578     /*
1579      * transferFeatures() will build 'partial overlap' for regions
1580      * that partially overlap 5' or 3' (start or end) of target sequence
1581      */
1582     AlignmentUtils.transferFeatures(dna, cds, map, null);
1583     List<SequenceFeature> sfs = cds.getSequenceFeatures();
1584     assertEquals(6, sfs.size());
1585
1586     SequenceFeature sf = sfs.get(0);
1587     assertEquals("type2", sf.getType());
1588     assertEquals("desc2", sf.getDescription());
1589     assertEquals(2f, sf.getScore());
1590     assertEquals(1, sf.getBegin());
1591     assertEquals(1, sf.getEnd());
1592
1593     sf = sfs.get(1);
1594     assertEquals("type3", sf.getType());
1595     assertEquals("desc3", sf.getDescription());
1596     assertEquals(3f, sf.getScore());
1597     assertEquals(1, sf.getBegin());
1598     assertEquals(3, sf.getEnd());
1599
1600     sf = sfs.get(2);
1601     assertEquals("type4", sf.getType());
1602     assertEquals(2, sf.getBegin());
1603     assertEquals(5, sf.getEnd());
1604
1605     sf = sfs.get(3);
1606     assertEquals("type5", sf.getType());
1607     assertEquals(1, sf.getBegin());
1608     assertEquals(6, sf.getEnd());
1609
1610     sf = sfs.get(4);
1611     assertEquals("type8", sf.getType());
1612     assertEquals(6, sf.getBegin());
1613     assertEquals(6, sf.getEnd());
1614
1615     sf = sfs.get(5);
1616     assertEquals("type9", sf.getType());
1617     assertEquals(6, sf.getBegin());
1618     assertEquals(6, sf.getEnd());
1619   }
1620
1621   /**
1622    * Tests for transferring features between mapped sequences
1623    */
1624   @Test(groups = { "Functional" })
1625   public void testTransferFeatures_withOmit()
1626   {
1627     SequenceI dna = new Sequence("dna/20-34", "acgTAGcaaGCCcgt");
1628     SequenceI cds = new Sequence("cds/10-15", "TAGGCC");
1629
1630     MapList map = new MapList(new int[] { 4, 6, 10, 12 },
1631             new int[]
1632             { 1, 6 }, 1, 1);
1633
1634     // [5, 11] maps to [2, 5]
1635     dna.addSequenceFeature(
1636             new SequenceFeature("type4", "desc4", 5, 11, 4f, null));
1637     // [4, 12] maps to [1, 6]
1638     dna.addSequenceFeature(
1639             new SequenceFeature("type5", "desc5", 4, 12, 5f, null));
1640     // [12, 12] maps to [6, 6]
1641     dna.addSequenceFeature(
1642             new SequenceFeature("type8", "desc8", 12, 12, 8f, null));
1643
1644     // desc4 and desc8 are the 'omit these' varargs
1645     AlignmentUtils.transferFeatures(dna, cds, map, null, "type4", "type8");
1646     List<SequenceFeature> sfs = cds.getSequenceFeatures();
1647     assertEquals(1, sfs.size());
1648
1649     SequenceFeature sf = sfs.get(0);
1650     assertEquals("type5", sf.getType());
1651     assertEquals(1, sf.getBegin());
1652     assertEquals(6, sf.getEnd());
1653   }
1654
1655   /**
1656    * Tests for transferring features between mapped sequences
1657    */
1658   @Test(groups = { "Functional" })
1659   public void testTransferFeatures_withSelect()
1660   {
1661     SequenceI dna = new Sequence("dna/20-34", "acgTAGcaaGCCcgt");
1662     SequenceI cds = new Sequence("cds/10-15", "TAGGCC");
1663
1664     MapList map = new MapList(new int[] { 4, 6, 10, 12 },
1665             new int[]
1666             { 1, 6 }, 1, 1);
1667
1668     // [5, 11] maps to [2, 5]
1669     dna.addSequenceFeature(
1670             new SequenceFeature("type4", "desc4", 5, 11, 4f, null));
1671     // [4, 12] maps to [1, 6]
1672     dna.addSequenceFeature(
1673             new SequenceFeature("type5", "desc5", 4, 12, 5f, null));
1674     // [12, 12] maps to [6, 6]
1675     dna.addSequenceFeature(
1676             new SequenceFeature("type8", "desc8", 12, 12, 8f, null));
1677
1678     // "type5" is the 'select this type' argument
1679     AlignmentUtils.transferFeatures(dna, cds, map, "type5");
1680     List<SequenceFeature> sfs = cds.getSequenceFeatures();
1681     assertEquals(1, sfs.size());
1682
1683     SequenceFeature sf = sfs.get(0);
1684     assertEquals("type5", sf.getType());
1685     assertEquals(1, sf.getBegin());
1686     assertEquals(6, sf.getEnd());
1687   }
1688
1689   /**
1690    * Test the method that extracts the cds-only part of a dna alignment, for the
1691    * case where the cds should be aligned to match its nucleotide sequence.
1692    */
1693   @Test(groups = { "Functional" })
1694   public void testMakeCdsAlignment_alternativeTranscripts()
1695   {
1696     SequenceI dna1 = new Sequence("dna1", "aaaGGGCC-----CTTTaaaGGG");
1697     // alternative transcript of same dna skips CCC codon
1698     SequenceI dna2 = new Sequence("dna2", "aaaGGGCC-----cttTaaaGGG");
1699     // dna3 has no mapping (protein product) so should be ignored here
1700     SequenceI dna3 = new Sequence("dna3", "aaaGGGCCCCCGGGcttTaaaGGG");
1701     SequenceI pep1 = new Sequence("pep1", "GPFG");
1702     SequenceI pep2 = new Sequence("pep2", "GPG");
1703     dna1.createDatasetSequence();
1704     dna2.createDatasetSequence();
1705     dna3.createDatasetSequence();
1706     pep1.createDatasetSequence();
1707     pep2.createDatasetSequence();
1708
1709     AlignmentI dna = new Alignment(new SequenceI[] { dna1, dna2, dna3 });
1710     dna.setDataset(null);
1711
1712     MapList map = new MapList(new int[] { 4, 12, 16, 18 },
1713             new int[]
1714             { 1, 4 }, 3, 1);
1715     AlignedCodonFrame acf = new AlignedCodonFrame();
1716     acf.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(), map);
1717     dna.addCodonFrame(acf);
1718     map = new MapList(new int[] { 4, 8, 12, 12, 16, 18 },
1719             new int[]
1720             { 1, 3 }, 3, 1);
1721     acf = new AlignedCodonFrame();
1722     acf.addMap(dna2.getDatasetSequence(), pep2.getDatasetSequence(), map);
1723     dna.addCodonFrame(acf);
1724
1725     AlignmentI cds = AlignmentUtils
1726             .makeCdsAlignment(new SequenceI[]
1727             { dna1, dna2, dna3 }, dna.getDataset(), null);
1728     List<SequenceI> cdsSeqs = cds.getSequences();
1729     assertEquals(2, cdsSeqs.size());
1730     assertEquals("GGGCCCTTTGGG", cdsSeqs.get(0).getSequenceAsString());
1731     assertEquals("GGGCCTGGG", cdsSeqs.get(1).getSequenceAsString());
1732
1733     /*
1734      * verify shared, extended alignment dataset
1735      */
1736     assertSame(dna.getDataset(), cds.getDataset());
1737     assertTrue(dna.getDataset().getSequences()
1738             .contains(cdsSeqs.get(0).getDatasetSequence()));
1739     assertTrue(dna.getDataset().getSequences()
1740             .contains(cdsSeqs.get(1).getDatasetSequence()));
1741
1742     /*
1743      * Verify 6 mappings: dna1 to cds1, cds1 to pep1, dna1 to pep1
1744      * and the same for dna2/cds2/pep2
1745      */
1746     List<AlignedCodonFrame> mappings = cds.getCodonFrames();
1747     assertEquals(6, mappings.size());
1748
1749     /*
1750      * 2 mappings involve pep1
1751      */
1752     List<AlignedCodonFrame> pep1Mappings = MappingUtils
1753             .findMappingsForSequence(pep1, mappings);
1754     assertEquals(2, pep1Mappings.size());
1755
1756     /*
1757      * Get mapping of pep1 to cds1 and verify it
1758      * maps GPFG to 1-3,4-6,7-9,10-12
1759      */
1760     List<AlignedCodonFrame> pep1CdsMappings = MappingUtils
1761             .findMappingsForSequence(cds.getSequenceAt(0), pep1Mappings);
1762     assertEquals(1, pep1CdsMappings.size());
1763     SearchResultsI sr = MappingUtils.buildSearchResults(pep1, 1,
1764             pep1CdsMappings);
1765     assertEquals(1, sr.getResults().size());
1766     SearchResultMatchI m = sr.getResults().get(0);
1767     assertEquals(cds.getSequenceAt(0).getDatasetSequence(),
1768             m.getSequence());
1769     assertEquals(1, m.getStart());
1770     assertEquals(3, m.getEnd());
1771     sr = MappingUtils.buildSearchResults(pep1, 2, pep1CdsMappings);
1772     m = sr.getResults().get(0);
1773     assertEquals(4, m.getStart());
1774     assertEquals(6, m.getEnd());
1775     sr = MappingUtils.buildSearchResults(pep1, 3, pep1CdsMappings);
1776     m = sr.getResults().get(0);
1777     assertEquals(7, m.getStart());
1778     assertEquals(9, m.getEnd());
1779     sr = MappingUtils.buildSearchResults(pep1, 4, pep1CdsMappings);
1780     m = sr.getResults().get(0);
1781     assertEquals(10, m.getStart());
1782     assertEquals(12, m.getEnd());
1783
1784     /*
1785      * Get mapping of pep2 to cds2 and verify it
1786      * maps GPG in pep2 to 1-3,4-6,7-9 in second CDS sequence
1787      */
1788     List<AlignedCodonFrame> pep2Mappings = MappingUtils
1789             .findMappingsForSequence(pep2, mappings);
1790     assertEquals(2, pep2Mappings.size());
1791     List<AlignedCodonFrame> pep2CdsMappings = MappingUtils
1792             .findMappingsForSequence(cds.getSequenceAt(1), pep2Mappings);
1793     assertEquals(1, pep2CdsMappings.size());
1794     sr = MappingUtils.buildSearchResults(pep2, 1, pep2CdsMappings);
1795     assertEquals(1, sr.getResults().size());
1796     m = sr.getResults().get(0);
1797     assertEquals(cds.getSequenceAt(1).getDatasetSequence(),
1798             m.getSequence());
1799     assertEquals(1, m.getStart());
1800     assertEquals(3, m.getEnd());
1801     sr = MappingUtils.buildSearchResults(pep2, 2, pep2CdsMappings);
1802     m = sr.getResults().get(0);
1803     assertEquals(4, m.getStart());
1804     assertEquals(6, m.getEnd());
1805     sr = MappingUtils.buildSearchResults(pep2, 3, pep2CdsMappings);
1806     m = sr.getResults().get(0);
1807     assertEquals(7, m.getStart());
1808     assertEquals(9, m.getEnd());
1809   }
1810
1811   /**
1812    * Test the method that realigns protein to match mapped codon alignment.
1813    */
1814   @Test(groups = { "Functional" })
1815   public void testAlignProteinAsDna_incompleteStartCodon()
1816   {
1817     // seq1: incomplete start codon (not mapped), then [3, 11]
1818     SequenceI dna1 = new Sequence("Seq1", "ccAAA-TTT-GGG-");
1819     // seq2 codons are [4, 5], [8, 11]
1820     SequenceI dna2 = new Sequence("Seq2", "ccaAA-ttT-GGG-");
1821     // seq3 incomplete start codon at 'tt'
1822     SequenceI dna3 = new Sequence("Seq3", "ccaaa-ttt-GGG-");
1823     AlignmentI dna = new Alignment(new SequenceI[] { dna1, dna2, dna3 });
1824     dna.setDataset(null);
1825
1826     // prot1 has 'X' for incomplete start codon (not mapped)
1827     SequenceI prot1 = new Sequence("Seq1", "XKFG"); // X for incomplete start
1828     SequenceI prot2 = new Sequence("Seq2", "NG");
1829     SequenceI prot3 = new Sequence("Seq3", "XG"); // X for incomplete start
1830     AlignmentI protein = new Alignment(
1831             new SequenceI[]
1832             { prot1, prot2, prot3 });
1833     protein.setDataset(null);
1834
1835     // map dna1 [3, 11] to prot1 [2, 4] KFG
1836     MapList map = new MapList(new int[] { 3, 11 }, new int[] { 2, 4 }, 3,
1837             1);
1838     AlignedCodonFrame acf = new AlignedCodonFrame();
1839     acf.addMap(dna1.getDatasetSequence(), prot1.getDatasetSequence(), map);
1840
1841     // map dna2 [4, 5] [8, 11] to prot2 [1, 2] NG
1842     map = new MapList(new int[] { 4, 5, 8, 11 }, new int[] { 1, 2 }, 3, 1);
1843     acf.addMap(dna2.getDatasetSequence(), prot2.getDatasetSequence(), map);
1844
1845     // map dna3 [9, 11] to prot3 [2, 2] G
1846     map = new MapList(new int[] { 9, 11 }, new int[] { 2, 2 }, 3, 1);
1847     acf.addMap(dna3.getDatasetSequence(), prot3.getDatasetSequence(), map);
1848
1849     ArrayList<AlignedCodonFrame> acfs = new ArrayList<>();
1850     acfs.add(acf);
1851     protein.setCodonFrames(acfs);
1852
1853     /*
1854      * verify X is included in the aligned proteins, and placed just
1855      * before the first mapped residue
1856      * CCT is between CCC and TTT
1857      */
1858     AlignmentUtils.alignProteinAsDna(protein, dna);
1859     assertEquals("XK-FG", prot1.getSequenceAsString());
1860     assertEquals("--N-G", prot2.getSequenceAsString());
1861     assertEquals("---XG", prot3.getSequenceAsString());
1862   }
1863
1864   /**
1865    * Tests for the method that maps the subset of a dna sequence that has CDS
1866    * (or subtype) feature - case where the start codon is incomplete.
1867    */
1868   @Test(groups = "Functional")
1869   public void testFindCdsPositions_fivePrimeIncomplete()
1870   {
1871     SequenceI dnaSeq = new Sequence("dna", "aaagGGCCCaaaTTTttt");
1872     dnaSeq.createDatasetSequence();
1873     SequenceI ds = dnaSeq.getDatasetSequence();
1874
1875     // CDS for dna 5-6 (incomplete codon), 7-9
1876     SequenceFeature sf = new SequenceFeature("CDS", "", 5, 9, 0f, null);
1877     sf.setPhase("2"); // skip 2 bases to start of next codon
1878     ds.addSequenceFeature(sf);
1879     // CDS for dna 13-15
1880     sf = new SequenceFeature("CDS_predicted", "", 13, 15, 0f, null);
1881     ds.addSequenceFeature(sf);
1882
1883     List<int[]> ranges = AlignmentUtils.findCdsPositions(dnaSeq);
1884
1885     /*
1886      * check the mapping starts with the first complete codon
1887      */
1888     assertEquals(6, MappingUtils.getLength(ranges));
1889     assertEquals(2, ranges.size());
1890     assertEquals(7, ranges.get(0)[0]);
1891     assertEquals(9, ranges.get(0)[1]);
1892     assertEquals(13, ranges.get(1)[0]);
1893     assertEquals(15, ranges.get(1)[1]);
1894   }
1895
1896   /**
1897    * Tests for the method that maps the subset of a dna sequence that has CDS
1898    * (or subtype) feature.
1899    */
1900   @Test(groups = "Functional")
1901   public void testFindCdsPositions()
1902   {
1903     SequenceI dnaSeq = new Sequence("dna", "aaaGGGcccAAATTTttt");
1904     dnaSeq.createDatasetSequence();
1905     SequenceI ds = dnaSeq.getDatasetSequence();
1906
1907     // CDS for dna 10-12
1908     SequenceFeature sf = new SequenceFeature("CDS_predicted", "", 10, 12,
1909             0f, null);
1910     sf.setStrand("+");
1911     ds.addSequenceFeature(sf);
1912     // CDS for dna 4-6
1913     sf = new SequenceFeature("CDS", "", 4, 6, 0f, null);
1914     sf.setStrand("+");
1915     ds.addSequenceFeature(sf);
1916     // exon feature should be ignored here
1917     sf = new SequenceFeature("exon", "", 7, 9, 0f, null);
1918     ds.addSequenceFeature(sf);
1919
1920     List<int[]> ranges = AlignmentUtils.findCdsPositions(dnaSeq);
1921     /*
1922      * verify ranges { [4-6], [12-10] }
1923      * note CDS ranges are ordered ascending even if the CDS
1924      * features are not
1925      */
1926     assertEquals(6, MappingUtils.getLength(ranges));
1927     assertEquals(2, ranges.size());
1928     assertEquals(4, ranges.get(0)[0]);
1929     assertEquals(6, ranges.get(0)[1]);
1930     assertEquals(10, ranges.get(1)[0]);
1931     assertEquals(12, ranges.get(1)[1]);
1932   }
1933
1934   /**
1935    * Tests for the method that maps the subset of a dna sequence that has CDS
1936    * (or subtype) feature, with CDS strand = '-' (reverse)
1937    */
1938   // test turned off as currently findCdsPositions is not strand-dependent
1939   // left in case it comes around again...
1940   @Test(groups = "Functional", enabled = false)
1941   public void testFindCdsPositions_reverseStrand()
1942   {
1943     SequenceI dnaSeq = new Sequence("dna", "aaaGGGcccAAATTTttt");
1944     dnaSeq.createDatasetSequence();
1945     SequenceI ds = dnaSeq.getDatasetSequence();
1946
1947     // CDS for dna 4-6
1948     SequenceFeature sf = new SequenceFeature("CDS", "", 4, 6, 0f, null);
1949     sf.setStrand("-");
1950     ds.addSequenceFeature(sf);
1951     // exon feature should be ignored here
1952     sf = new SequenceFeature("exon", "", 7, 9, 0f, null);
1953     ds.addSequenceFeature(sf);
1954     // CDS for dna 10-12
1955     sf = new SequenceFeature("CDS_predicted", "", 10, 12, 0f, null);
1956     sf.setStrand("-");
1957     ds.addSequenceFeature(sf);
1958
1959     List<int[]> ranges = AlignmentUtils.findCdsPositions(dnaSeq);
1960     /*
1961      * verify ranges { [12-10], [6-4] }
1962      */
1963     assertEquals(6, MappingUtils.getLength(ranges));
1964     assertEquals(2, ranges.size());
1965     assertEquals(12, ranges.get(0)[0]);
1966     assertEquals(10, ranges.get(0)[1]);
1967     assertEquals(6, ranges.get(1)[0]);
1968     assertEquals(4, ranges.get(1)[1]);
1969   }
1970
1971   /**
1972    * Tests for the method that maps the subset of a dna sequence that has CDS
1973    * (or subtype) feature - reverse strand case where the start codon is
1974    * incomplete.
1975    */
1976   @Test(groups = "Functional", enabled = false)
1977   // test turned off as currently findCdsPositions is not strand-dependent
1978   // left in case it comes around again...
1979   public void testFindCdsPositions_reverseStrandThreePrimeIncomplete()
1980   {
1981     SequenceI dnaSeq = new Sequence("dna", "aaagGGCCCaaaTTTttt");
1982     dnaSeq.createDatasetSequence();
1983     SequenceI ds = dnaSeq.getDatasetSequence();
1984
1985     // CDS for dna 5-9
1986     SequenceFeature sf = new SequenceFeature("CDS", "", 5, 9, 0f, null);
1987     sf.setStrand("-");
1988     ds.addSequenceFeature(sf);
1989     // CDS for dna 13-15
1990     sf = new SequenceFeature("CDS_predicted", "", 13, 15, 0f, null);
1991     sf.setStrand("-");
1992     sf.setPhase("2"); // skip 2 bases to start of next codon
1993     ds.addSequenceFeature(sf);
1994
1995     List<int[]> ranges = AlignmentUtils.findCdsPositions(dnaSeq);
1996
1997     /*
1998      * check the mapping starts with the first complete codon
1999      * expect ranges [13, 13], [9, 5]
2000      */
2001     assertEquals(6, MappingUtils.getLength(ranges));
2002     assertEquals(2, ranges.size());
2003     assertEquals(13, ranges.get(0)[0]);
2004     assertEquals(13, ranges.get(0)[1]);
2005     assertEquals(9, ranges.get(1)[0]);
2006     assertEquals(5, ranges.get(1)[1]);
2007   }
2008
2009   @Test(groups = "Functional")
2010   public void testAlignAs_alternateTranscriptsUngapped()
2011   {
2012     SequenceI dna1 = new Sequence("dna1", "cccGGGTTTaaa");
2013     SequenceI dna2 = new Sequence("dna2", "CCCgggtttAAA");
2014     AlignmentI dna = new Alignment(new SequenceI[] { dna1, dna2 });
2015     ((Alignment) dna).createDatasetAlignment();
2016     SequenceI cds1 = new Sequence("cds1", "GGGTTT");
2017     SequenceI cds2 = new Sequence("cds2", "CCCAAA");
2018     AlignmentI cds = new Alignment(new SequenceI[] { cds1, cds2 });
2019     ((Alignment) cds).createDatasetAlignment();
2020
2021     AlignedCodonFrame acf = new AlignedCodonFrame();
2022     MapList map = new MapList(new int[] { 4, 9 }, new int[] { 1, 6 }, 1, 1);
2023     acf.addMap(dna1.getDatasetSequence(), cds1.getDatasetSequence(), map);
2024     map = new MapList(new int[] { 1, 3, 10, 12 }, new int[] { 1, 6 }, 1, 1);
2025     acf.addMap(dna2.getDatasetSequence(), cds2.getDatasetSequence(), map);
2026
2027     /*
2028      * verify CDS alignment is as:
2029      *   cccGGGTTTaaa (cdna)
2030      *   CCCgggtttAAA (cdna)
2031      *
2032      *   ---GGGTTT--- (cds)
2033      *   CCC------AAA (cds)
2034      */
2035     dna.addCodonFrame(acf);
2036     AlignmentUtils.alignAs(cds, dna);
2037     assertEquals("---GGGTTT", cds.getSequenceAt(0).getSequenceAsString());
2038     assertEquals("CCC------AAA",
2039             cds.getSequenceAt(1).getSequenceAsString());
2040   }
2041
2042   @Test(groups = { "Functional" })
2043   public void testAddMappedPositions()
2044   {
2045     SequenceI from = new Sequence("dna", "ggAA-ATcc-TT-g");
2046     SequenceI seq1 = new Sequence("cds", "AAATTT");
2047     from.createDatasetSequence();
2048     seq1.createDatasetSequence();
2049     Mapping mapping = new Mapping(seq1,
2050             new MapList(new int[]
2051             { 3, 6, 9, 10 }, new int[] { 1, 6 }, 1, 1));
2052     Map<Integer, Map<SequenceI, Character>> map = new TreeMap<>();
2053     AlignmentUtils.addMappedPositions(seq1, from, mapping, map);
2054
2055     /*
2056      * verify map has seq1 residues in columns 3,4,6,7,11,12
2057      */
2058     assertEquals(6, map.size());
2059     assertEquals('A', map.get(3).get(seq1).charValue());
2060     assertEquals('A', map.get(4).get(seq1).charValue());
2061     assertEquals('A', map.get(6).get(seq1).charValue());
2062     assertEquals('T', map.get(7).get(seq1).charValue());
2063     assertEquals('T', map.get(11).get(seq1).charValue());
2064     assertEquals('T', map.get(12).get(seq1).charValue());
2065
2066     /*
2067      *
2068      */
2069   }
2070
2071   /**
2072    * Test case where the mapping 'from' range includes a stop codon which is
2073    * absent in the 'to' range
2074    */
2075   @Test(groups = { "Functional" })
2076   public void testAddMappedPositions_withStopCodon()
2077   {
2078     SequenceI from = new Sequence("dna", "ggAA-ATcc-TT-g");
2079     SequenceI seq1 = new Sequence("cds", "AAATTT");
2080     from.createDatasetSequence();
2081     seq1.createDatasetSequence();
2082     Mapping mapping = new Mapping(seq1,
2083             new MapList(new int[]
2084             { 3, 6, 9, 10 }, new int[] { 1, 6 }, 1, 1));
2085     Map<Integer, Map<SequenceI, Character>> map = new TreeMap<>();
2086     AlignmentUtils.addMappedPositions(seq1, from, mapping, map);
2087
2088     /*
2089      * verify map has seq1 residues in columns 3,4,6,7,11,12
2090      */
2091     assertEquals(6, map.size());
2092     assertEquals('A', map.get(3).get(seq1).charValue());
2093     assertEquals('A', map.get(4).get(seq1).charValue());
2094     assertEquals('A', map.get(6).get(seq1).charValue());
2095     assertEquals('T', map.get(7).get(seq1).charValue());
2096     assertEquals('T', map.get(11).get(seq1).charValue());
2097     assertEquals('T', map.get(12).get(seq1).charValue());
2098   }
2099
2100   /**
2101    * Test for the case where the products for which we want CDS are specified.
2102    * This is to represent the case where EMBL has CDS mappings to both Uniprot
2103    * and EMBLCDSPROTEIN. makeCdsAlignment() should only return the mappings for
2104    * the protein sequences specified.
2105    */
2106   @Test(groups = { "Functional" })
2107   public void testMakeCdsAlignment_filterProducts()
2108   {
2109     SequenceI dna1 = new Sequence("dna1", "aaaGGGcccTTTaaa");
2110     SequenceI dna2 = new Sequence("dna2", "GGGcccTTTaaaCCC");
2111     SequenceI pep1 = new Sequence("Uniprot|pep1", "GF");
2112     SequenceI pep2 = new Sequence("Uniprot|pep2", "GFP");
2113     SequenceI pep3 = new Sequence("EMBL|pep3", "GF");
2114     SequenceI pep4 = new Sequence("EMBL|pep4", "GFP");
2115     dna1.createDatasetSequence();
2116     dna2.createDatasetSequence();
2117     pep1.createDatasetSequence();
2118     pep2.createDatasetSequence();
2119     pep3.createDatasetSequence();
2120     pep4.createDatasetSequence();
2121     AlignmentI dna = new Alignment(new SequenceI[] { dna1, dna2 });
2122     dna.setDataset(null);
2123     AlignmentI emblPeptides = new Alignment(new SequenceI[] { pep3, pep4 });
2124     emblPeptides.setDataset(null);
2125
2126     AlignedCodonFrame acf = new AlignedCodonFrame();
2127     MapList map = new MapList(new int[] { 4, 6, 10, 12 },
2128             new int[]
2129             { 1, 2 }, 3, 1);
2130     acf.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(), map);
2131     acf.addMap(dna1.getDatasetSequence(), pep3.getDatasetSequence(), map);
2132     dna.addCodonFrame(acf);
2133
2134     acf = new AlignedCodonFrame();
2135     map = new MapList(new int[] { 1, 3, 7, 9, 13, 15 }, new int[] { 1, 3 },
2136             3, 1);
2137     acf.addMap(dna2.getDatasetSequence(), pep2.getDatasetSequence(), map);
2138     acf.addMap(dna2.getDatasetSequence(), pep4.getDatasetSequence(), map);
2139     dna.addCodonFrame(acf);
2140
2141     /*
2142      * execute method under test to find CDS for EMBL peptides only
2143      */
2144     AlignmentI cds = AlignmentUtils
2145             .makeCdsAlignment(new SequenceI[]
2146             { dna1, dna2 }, dna.getDataset(),
2147                     emblPeptides.getSequencesArray());
2148
2149     assertEquals(2, cds.getSequences().size());
2150     assertEquals("GGGTTT", cds.getSequenceAt(0).getSequenceAsString());
2151     assertEquals("GGGTTTCCC", cds.getSequenceAt(1).getSequenceAsString());
2152
2153     /*
2154      * verify shared, extended alignment dataset
2155      */
2156     assertSame(dna.getDataset(), cds.getDataset());
2157     assertTrue(dna.getDataset().getSequences()
2158             .contains(cds.getSequenceAt(0).getDatasetSequence()));
2159     assertTrue(dna.getDataset().getSequences()
2160             .contains(cds.getSequenceAt(1).getDatasetSequence()));
2161
2162     /*
2163      * Verify mappings from CDS to peptide, cDNA to CDS, and cDNA to peptide
2164      * the mappings are on the shared alignment dataset
2165      */
2166     List<AlignedCodonFrame> cdsMappings = cds.getDataset().getCodonFrames();
2167     /*
2168      * 6 mappings, 2*(DNA->CDS), 2*(DNA->Pep), 2*(CDS->Pep)
2169      */
2170     assertEquals(6, cdsMappings.size());
2171
2172     /*
2173      * verify that mapping sets for dna and cds alignments are different
2174      * [not current behaviour - all mappings are on the alignment dataset]
2175      */
2176     // select -> subselect type to test.
2177     // Assert.assertNotSame(dna.getCodonFrames(), cds.getCodonFrames());
2178     // assertEquals(4, dna.getCodonFrames().size());
2179     // assertEquals(4, cds.getCodonFrames().size());
2180
2181     /*
2182      * Two mappings involve pep3 (dna to pep3, cds to pep3)
2183      * Mapping from pep3 to GGGTTT in first new exon sequence
2184      */
2185     List<AlignedCodonFrame> pep3Mappings = MappingUtils
2186             .findMappingsForSequence(pep3, cdsMappings);
2187     assertEquals(2, pep3Mappings.size());
2188     List<AlignedCodonFrame> mappings = MappingUtils
2189             .findMappingsForSequence(cds.getSequenceAt(0), pep3Mappings);
2190     assertEquals(1, mappings.size());
2191
2192     // map G to GGG
2193     SearchResultsI sr = MappingUtils.buildSearchResults(pep3, 1, mappings);
2194     assertEquals(1, sr.getResults().size());
2195     SearchResultMatchI m = sr.getResults().get(0);
2196     assertSame(cds.getSequenceAt(0).getDatasetSequence(), m.getSequence());
2197     assertEquals(1, m.getStart());
2198     assertEquals(3, m.getEnd());
2199     // map F to TTT
2200     sr = MappingUtils.buildSearchResults(pep3, 2, mappings);
2201     m = sr.getResults().get(0);
2202     assertSame(cds.getSequenceAt(0).getDatasetSequence(), m.getSequence());
2203     assertEquals(4, m.getStart());
2204     assertEquals(6, m.getEnd());
2205
2206     /*
2207      * Two mappings involve pep4 (dna to pep4, cds to pep4)
2208      * Verify mapping from pep4 to GGGTTTCCC in second new exon sequence
2209      */
2210     List<AlignedCodonFrame> pep4Mappings = MappingUtils
2211             .findMappingsForSequence(pep4, cdsMappings);
2212     assertEquals(2, pep4Mappings.size());
2213     mappings = MappingUtils.findMappingsForSequence(cds.getSequenceAt(1),
2214             pep4Mappings);
2215     assertEquals(1, mappings.size());
2216     // map G to GGG
2217     sr = MappingUtils.buildSearchResults(pep4, 1, mappings);
2218     assertEquals(1, sr.getResults().size());
2219     m = sr.getResults().get(0);
2220     assertSame(cds.getSequenceAt(1).getDatasetSequence(), m.getSequence());
2221     assertEquals(1, m.getStart());
2222     assertEquals(3, m.getEnd());
2223     // map F to TTT
2224     sr = MappingUtils.buildSearchResults(pep4, 2, mappings);
2225     m = sr.getResults().get(0);
2226     assertSame(cds.getSequenceAt(1).getDatasetSequence(), m.getSequence());
2227     assertEquals(4, m.getStart());
2228     assertEquals(6, m.getEnd());
2229     // map P to CCC
2230     sr = MappingUtils.buildSearchResults(pep4, 3, mappings);
2231     m = sr.getResults().get(0);
2232     assertSame(cds.getSequenceAt(1).getDatasetSequence(), m.getSequence());
2233     assertEquals(7, m.getStart());
2234     assertEquals(9, m.getEnd());
2235   }
2236
2237   /**
2238    * Test the method that just copies aligned sequences, provided all sequences
2239    * to be aligned share the aligned sequence's dataset
2240    */
2241   @Test(groups = "Functional")
2242   public void testAlignAsSameSequences()
2243   {
2244     SequenceI dna1 = new Sequence("dna1", "cccGGGTTTaaa");
2245     SequenceI dna2 = new Sequence("dna2", "CCCgggtttAAA");
2246     AlignmentI al1 = new Alignment(new SequenceI[] { dna1, dna2 });
2247     ((Alignment) al1).createDatasetAlignment();
2248
2249     SequenceI dna3 = new Sequence(dna1);
2250     SequenceI dna4 = new Sequence(dna2);
2251     assertSame(dna3.getDatasetSequence(), dna1.getDatasetSequence());
2252     assertSame(dna4.getDatasetSequence(), dna2.getDatasetSequence());
2253     String seq1 = "-cc-GG-GT-TT--aaa";
2254     dna3.setSequence(seq1);
2255     String seq2 = "C--C-Cgg--gtt-tAA-A-";
2256     dna4.setSequence(seq2);
2257     AlignmentI al2 = new Alignment(new SequenceI[] { dna3, dna4 });
2258     ((Alignment) al2).createDatasetAlignment();
2259
2260     /*
2261      * alignment removes gapped columns (two internal, two trailing)
2262      */
2263     assertTrue(AlignmentUtils.alignAsSameSequences(al1, al2));
2264     String aligned1 = "-cc-GG-GTTT-aaa";
2265     assertEquals(aligned1, al1.getSequenceAt(0).getSequenceAsString());
2266     String aligned2 = "C--C-Cgg-gtttAAA";
2267     assertEquals(aligned2, al1.getSequenceAt(1).getSequenceAsString());
2268
2269     /*
2270      * add another sequence to 'aligned' - should still succeed, since
2271      * unaligned sequences still share a dataset with aligned sequences
2272      */
2273     SequenceI dna5 = new Sequence("dna5", "CCCgggtttAAA");
2274     dna5.createDatasetSequence();
2275     al2.addSequence(dna5);
2276     assertTrue(AlignmentUtils.alignAsSameSequences(al1, al2));
2277     assertEquals(aligned1, al1.getSequenceAt(0).getSequenceAsString());
2278     assertEquals(aligned2, al1.getSequenceAt(1).getSequenceAsString());
2279
2280     /*
2281      * add another sequence to 'unaligned' - should fail, since now not
2282      * all unaligned sequences share a dataset with aligned sequences
2283      */
2284     SequenceI dna6 = new Sequence("dna6", "CCCgggtttAAA");
2285     dna6.createDatasetSequence();
2286     al1.addSequence(dna6);
2287     // JAL-2110 JBP Comment: what's the use case for this behaviour ?
2288     assertFalse(AlignmentUtils.alignAsSameSequences(al1, al2));
2289   }
2290
2291   @Test(groups = "Functional")
2292   public void testAlignAsSameSequencesMultipleSubSeq()
2293   {
2294     SequenceI dna1 = new Sequence("dna1", "cccGGGTTTaaa");
2295     SequenceI dna2 = new Sequence("dna2", "CCCgggtttAAA");
2296     SequenceI as1 = dna1.deriveSequence(); // cccGGGTTTaaa/1-12
2297     SequenceI as2 = dna1.deriveSequence().getSubSequence(3, 7); // GGGT/4-7
2298     SequenceI as3 = dna2.deriveSequence(); // CCCgggtttAAA/1-12
2299     as1.insertCharAt(6, 5, '-');
2300     assertEquals("cccGGG-----TTTaaa", as1.getSequenceAsString());
2301     as2.insertCharAt(6, 5, '-');
2302     assertEquals("GGGT-----", as2.getSequenceAsString());
2303     as3.insertCharAt(3, 5, '-');
2304     assertEquals("CCC-----gggtttAAA", as3.getSequenceAsString());
2305     AlignmentI aligned = new Alignment(new SequenceI[] { as1, as2, as3 });
2306
2307     // why do we need to cast this still ?
2308     ((Alignment) aligned).createDatasetAlignment();
2309     SequenceI uas1 = dna1.deriveSequence();
2310     SequenceI uas2 = dna1.deriveSequence().getSubSequence(3, 7);
2311     SequenceI uas3 = dna2.deriveSequence();
2312     AlignmentI tobealigned = new Alignment(
2313             new SequenceI[]
2314             { uas1, uas2, uas3 });
2315     ((Alignment) tobealigned).createDatasetAlignment();
2316
2317     /*
2318      * alignAs lines up dataset sequences and removes empty columns (two)
2319      */
2320     assertTrue(AlignmentUtils.alignAsSameSequences(tobealigned, aligned));
2321     assertEquals("cccGGG---TTTaaa", uas1.getSequenceAsString());
2322     assertEquals("GGGT", uas2.getSequenceAsString());
2323     assertEquals("CCC---gggtttAAA", uas3.getSequenceAsString());
2324   }
2325
2326   @Test(groups = { "Functional" })
2327   public void testTransferGeneLoci()
2328   {
2329     SequenceI from = new Sequence("transcript",
2330             "aaacccgggTTTAAACCCGGGtttaaacccgggttt");
2331     SequenceI to = new Sequence("CDS", "TTTAAACCCGGG");
2332     MapList map = new MapList(new int[] { 1, 12 }, new int[] { 10, 21 }, 1,
2333             1);
2334
2335     /*
2336      * first with nothing to transfer
2337      */
2338     AlignmentUtils.transferGeneLoci(from, map, to);
2339     assertNull(to.getGeneLoci());
2340
2341     /*
2342      * next with gene loci set on 'from' sequence
2343      */
2344     int[] exons = new int[] { 100, 105, 155, 164, 210, 229 };
2345     MapList geneMap = new MapList(new int[] { 1, 36 }, exons, 1, 1);
2346     from.setGeneLoci("human", "GRCh38", "7", geneMap);
2347     AlignmentUtils.transferGeneLoci(from, map, to);
2348
2349     GeneLociI toLoci = to.getGeneLoci();
2350     assertNotNull(toLoci);
2351     // DBRefEntry constructor upper-cases 'source'
2352     assertEquals("HUMAN", toLoci.getSpeciesId());
2353     assertEquals("GRCh38", toLoci.getAssemblyId());
2354     assertEquals("7", toLoci.getChromosomeId());
2355
2356     /*
2357      * transcript 'exons' are 1-6, 7-16, 17-36
2358      * CDS 1:12 is transcript 10-21
2359      * transcript 'CDS' is 10-16, 17-21
2360      * which is 'gene' 158-164, 210-214
2361      */
2362     MapList toMap = toLoci.getMapping();
2363     assertEquals(1, toMap.getFromRanges().size());
2364     assertEquals(2, toMap.getFromRanges().get(0).length);
2365     assertEquals(1, toMap.getFromRanges().get(0)[0]);
2366     assertEquals(12, toMap.getFromRanges().get(0)[1]);
2367     assertEquals(2, toMap.getToRanges().size());
2368     assertEquals(2, toMap.getToRanges().get(0).length);
2369     assertEquals(158, toMap.getToRanges().get(0)[0]);
2370     assertEquals(164, toMap.getToRanges().get(0)[1]);
2371     assertEquals(210, toMap.getToRanges().get(1)[0]);
2372     assertEquals(214, toMap.getToRanges().get(1)[1]);
2373     // or summarised as (but toString might change in future):
2374     assertEquals("[ [1, 12] ] 1:1 to [ [158, 164] [210, 214] ]",
2375             toMap.toString());
2376
2377     /*
2378      * an existing value is not overridden
2379      */
2380     geneMap = new MapList(new int[] { 1, 36 }, new int[] { 36, 1 }, 1, 1);
2381     from.setGeneLoci("inhuman", "GRCh37", "6", geneMap);
2382     AlignmentUtils.transferGeneLoci(from, map, to);
2383     assertEquals("GRCh38", toLoci.getAssemblyId());
2384     assertEquals("7", toLoci.getChromosomeId());
2385     toMap = toLoci.getMapping();
2386     assertEquals("[ [1, 12] ] 1:1 to [ [158, 164] [210, 214] ]",
2387             toMap.toString());
2388   }
2389
2390   /**
2391    * Tests for the method that maps nucleotide to protein based on CDS features
2392    */
2393   @Test(groups = "Functional")
2394   public void testMapCdsToProtein()
2395   {
2396     SequenceI peptide = new Sequence("pep", "KLQ");
2397
2398     /*
2399      * Case 1: CDS 3 times length of peptide
2400      * NB method only checks lengths match, not translation
2401      */
2402     SequenceI dna = new Sequence("dna", "AACGacgtCTCCT");
2403     dna.createDatasetSequence();
2404     dna.addSequenceFeature(new SequenceFeature("CDS", "", 1, 4, null));
2405     dna.addSequenceFeature(new SequenceFeature("CDS", "", 9, 13, null));
2406     MapList ml = AlignmentUtils.mapCdsToProtein(dna, peptide);
2407     assertEquals(3, ml.getFromRatio());
2408     assertEquals(1, ml.getToRatio());
2409     assertEquals("[[1, 3]]",
2410             Arrays.deepToString(ml.getToRanges().toArray()));
2411     assertEquals("[[1, 4], [9, 13]]",
2412             Arrays.deepToString(ml.getFromRanges().toArray()));
2413
2414     /*
2415      * Case 2: CDS 3 times length of peptide + stop codon
2416      * (note code does not currently check trailing codon is a stop codon)
2417      */
2418     dna = new Sequence("dna", "AACGacgtCTCCTCCC");
2419     dna.createDatasetSequence();
2420     dna.addSequenceFeature(new SequenceFeature("CDS", "", 1, 4, null));
2421     dna.addSequenceFeature(new SequenceFeature("CDS", "", 9, 16, null));
2422     ml = AlignmentUtils.mapCdsToProtein(dna, peptide);
2423     assertEquals(3, ml.getFromRatio());
2424     assertEquals(1, ml.getToRatio());
2425     assertEquals("[[1, 3]]",
2426             Arrays.deepToString(ml.getToRanges().toArray()));
2427     assertEquals("[[1, 4], [9, 13]]",
2428             Arrays.deepToString(ml.getFromRanges().toArray()));
2429
2430     /*
2431      * Case 3: CDS longer than 3 * peptide + stop codon - no mapping is made
2432      */
2433     dna = new Sequence("dna", "AACGacgtCTCCTTGATCA");
2434     dna.createDatasetSequence();
2435     dna.addSequenceFeature(new SequenceFeature("CDS", "", 1, 4, null));
2436     dna.addSequenceFeature(new SequenceFeature("CDS", "", 9, 19, null));
2437     ml = AlignmentUtils.mapCdsToProtein(dna, peptide);
2438     assertNull(ml);
2439
2440     /*
2441      * Case 4: CDS shorter than 3 * peptide - no mapping is made
2442      */
2443     dna = new Sequence("dna", "AACGacgtCTCC");
2444     dna.createDatasetSequence();
2445     dna.addSequenceFeature(new SequenceFeature("CDS", "", 1, 4, null));
2446     dna.addSequenceFeature(new SequenceFeature("CDS", "", 9, 12, null));
2447     ml = AlignmentUtils.mapCdsToProtein(dna, peptide);
2448     assertNull(ml);
2449
2450     /*
2451      * Case 5: CDS 3 times length of peptide + part codon - mapping is truncated
2452      */
2453     dna = new Sequence("dna", "AACGacgtCTCCTTG");
2454     dna.createDatasetSequence();
2455     dna.addSequenceFeature(new SequenceFeature("CDS", "", 1, 4, null));
2456     dna.addSequenceFeature(new SequenceFeature("CDS", "", 9, 15, null));
2457     ml = AlignmentUtils.mapCdsToProtein(dna, peptide);
2458     assertEquals(3, ml.getFromRatio());
2459     assertEquals(1, ml.getToRatio());
2460     assertEquals("[[1, 3]]",
2461             Arrays.deepToString(ml.getToRanges().toArray()));
2462     assertEquals("[[1, 4], [9, 13]]",
2463             Arrays.deepToString(ml.getFromRanges().toArray()));
2464
2465     /*
2466      * Case 6: incomplete start codon corresponding to X in peptide
2467      */
2468     dna = new Sequence("dna", "ACGacgtCTCCTTGG");
2469     dna.createDatasetSequence();
2470     SequenceFeature sf = new SequenceFeature("CDS", "", 1, 3, null);
2471     sf.setPhase("2"); // skip 2 positions (AC) to start of next codon (GCT)
2472     dna.addSequenceFeature(sf);
2473     dna.addSequenceFeature(new SequenceFeature("CDS", "", 8, 15, null));
2474     peptide = new Sequence("pep", "XLQ");
2475     ml = AlignmentUtils.mapCdsToProtein(dna, peptide);
2476     assertEquals("[[2, 3]]",
2477             Arrays.deepToString(ml.getToRanges().toArray()));
2478     assertEquals("[[3, 3], [8, 12]]",
2479             Arrays.deepToString(ml.getFromRanges().toArray()));
2480   }
2481
2482   /**
2483    * Tests for the method that locates the CDS sequence that has a mapping to
2484    * the given protein. That is, given a transcript-to-peptide mapping, find the
2485    * cds-to-peptide mapping that relates to both, and return the CDS sequence.
2486    */
2487   @Test(groups = "Functional")
2488   public void testFindCdsForProtein()
2489   {
2490     List<AlignedCodonFrame> mappings = new ArrayList<>();
2491     AlignedCodonFrame acf1 = new AlignedCodonFrame();
2492     mappings.add(acf1);
2493
2494     SequenceI dna1 = new Sequence("dna1", "cgatATcgGCTATCTATGacg");
2495     dna1.createDatasetSequence();
2496
2497     // NB we currently exclude STOP codon from CDS sequences
2498     // the test would need to change if this changes in future
2499     SequenceI cds1 = new Sequence("cds1", "ATGCTATCT");
2500     cds1.createDatasetSequence();
2501
2502     SequenceI pep1 = new Sequence("pep1", "MLS");
2503     pep1.createDatasetSequence();
2504     List<AlignedCodonFrame> seqMappings = new ArrayList<>();
2505     MapList mapList = new MapList(new int[] { 5, 6, 9, 15 },
2506             new int[]
2507             { 1, 3 }, 3, 1);
2508     Mapping dnaToPeptide = new Mapping(pep1.getDatasetSequence(), mapList);
2509
2510     // add dna to peptide mapping
2511     seqMappings.add(acf1);
2512     acf1.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(),
2513             mapList);
2514
2515     /*
2516      * first case - no dna-to-CDS mapping exists - search fails
2517      */
2518     SequenceI seq = AlignmentUtils.findCdsForProtein(mappings, dna1,
2519             seqMappings, dnaToPeptide);
2520     assertNull(seq);
2521
2522     /*
2523      * second case - CDS-to-peptide mapping exists but no dna-to-CDS
2524      * - search fails
2525      */
2526     // todo this test fails if the mapping is added to acf1, not acf2
2527     // need to tidy up use of lists of mappings in AlignedCodonFrame
2528     AlignedCodonFrame acf2 = new AlignedCodonFrame();
2529     mappings.add(acf2);
2530     MapList cdsToPeptideMapping = new MapList(new int[] { 1, 9 },
2531             new int[]
2532             { 1, 3 }, 3, 1);
2533     acf2.addMap(cds1.getDatasetSequence(), pep1.getDatasetSequence(),
2534             cdsToPeptideMapping);
2535     assertNull(AlignmentUtils.findCdsForProtein(mappings, dna1, seqMappings,
2536             dnaToPeptide));
2537
2538     /*
2539      * third case - add dna-to-CDS mapping - CDS is now found!
2540      */
2541     MapList dnaToCdsMapping = new MapList(new int[] { 5, 6, 9, 15 },
2542             new int[]
2543             { 1, 9 }, 1, 1);
2544     acf1.addMap(dna1.getDatasetSequence(), cds1.getDatasetSequence(),
2545             dnaToCdsMapping);
2546     seq = AlignmentUtils.findCdsForProtein(mappings, dna1, seqMappings,
2547             dnaToPeptide);
2548     assertSame(seq, cds1.getDatasetSequence());
2549   }
2550
2551   /**
2552    * Tests for the method that locates the CDS sequence that has a mapping to
2553    * the given protein. That is, given a transcript-to-peptide mapping, find the
2554    * cds-to-peptide mapping that relates to both, and return the CDS sequence.
2555    * This test is for the case where transcript and CDS are the same length.
2556    */
2557   @Test(groups = "Functional")
2558   public void testFindCdsForProtein_noUTR()
2559   {
2560     List<AlignedCodonFrame> mappings = new ArrayList<>();
2561     AlignedCodonFrame acf1 = new AlignedCodonFrame();
2562     mappings.add(acf1);
2563
2564     SequenceI dna1 = new Sequence("dna1", "ATGCTATCTTAA");
2565     dna1.createDatasetSequence();
2566
2567     // NB we currently exclude STOP codon from CDS sequences
2568     // the test would need to change if this changes in future
2569     SequenceI cds1 = new Sequence("cds1", "ATGCTATCT");
2570     cds1.createDatasetSequence();
2571
2572     SequenceI pep1 = new Sequence("pep1", "MLS");
2573     pep1.createDatasetSequence();
2574     List<AlignedCodonFrame> seqMappings = new ArrayList<>();
2575     MapList mapList = new MapList(new int[] { 1, 9 }, new int[] { 1, 3 }, 3,
2576             1);
2577     Mapping dnaToPeptide = new Mapping(pep1.getDatasetSequence(), mapList);
2578
2579     // add dna to peptide mapping
2580     seqMappings.add(acf1);
2581     acf1.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(),
2582             mapList);
2583
2584     /*
2585      * first case - transcript lacks CDS features - it appears to be
2586      * the CDS sequence and is returned
2587      */
2588     SequenceI seq = AlignmentUtils.findCdsForProtein(mappings, dna1,
2589             seqMappings, dnaToPeptide);
2590     assertSame(seq, dna1.getDatasetSequence());
2591
2592     /*
2593      * second case - transcript has CDS feature - this means it is
2594      * not returned as a match for CDS (CDS sequences don't have CDS features)
2595      */
2596     dna1.addSequenceFeature(
2597             new SequenceFeature(SequenceOntologyI.CDS, "cds", 1, 12, null));
2598     seq = AlignmentUtils.findCdsForProtein(mappings, dna1, seqMappings,
2599             dnaToPeptide);
2600     assertNull(seq);
2601
2602     /*
2603      * third case - CDS-to-peptide mapping exists but no dna-to-CDS
2604      * - search fails
2605      */
2606     // todo this test fails if the mapping is added to acf1, not acf2
2607     // need to tidy up use of lists of mappings in AlignedCodonFrame
2608     AlignedCodonFrame acf2 = new AlignedCodonFrame();
2609     mappings.add(acf2);
2610     MapList cdsToPeptideMapping = new MapList(new int[] { 1, 9 },
2611             new int[]
2612             { 1, 3 }, 3, 1);
2613     acf2.addMap(cds1.getDatasetSequence(), pep1.getDatasetSequence(),
2614             cdsToPeptideMapping);
2615     assertNull(AlignmentUtils.findCdsForProtein(mappings, dna1, seqMappings,
2616             dnaToPeptide));
2617
2618     /*
2619      * fourth case - add dna-to-CDS mapping - CDS is now found!
2620      */
2621     MapList dnaToCdsMapping = new MapList(new int[] { 1, 9 },
2622             new int[]
2623             { 1, 9 }, 1, 1);
2624     acf1.addMap(dna1.getDatasetSequence(), cds1.getDatasetSequence(),
2625             dnaToCdsMapping);
2626     seq = AlignmentUtils.findCdsForProtein(mappings, dna1, seqMappings,
2627             dnaToPeptide);
2628     assertSame(seq, cds1.getDatasetSequence());
2629   }
2630
2631   @Test(groups = "Functional")
2632   public void testAddReferenceAnnotations()
2633   {
2634     SequenceI longseq = new Sequence("longA", "ASDASDASDASDAASDASDASDASDA");
2635     Annotation[] aa = new Annotation[longseq.getLength()];
2636
2637     for (int p = 0; p < aa.length; p++)
2638     {
2639       aa[p] = new Annotation("P", "pos " + (p + 1), (char) 0,
2640               (float) p + 1);
2641     }
2642     AlignmentAnnotation refAnnot = new AlignmentAnnotation("LongSeqAnnot",
2643             "Annotations", aa);
2644     refAnnot.setCalcId("Test");
2645     longseq.addAlignmentAnnotation(refAnnot);
2646     verifyExpectedSequenceAnnotation(refAnnot);
2647
2648     Alignment ourAl = new Alignment(
2649             new SequenceI[]
2650             { longseq.getSubSequence(5, 10),
2651                 longseq.getSubSequence(7, 12) });
2652     ourAl.createDatasetAlignment();
2653
2654     // transfer annotation
2655     SortedMap<String, String> tipEntries = new TreeMap<>();
2656     Map<SequenceI, List<AlignmentAnnotation>> candidates = new LinkedHashMap<>();
2657
2658     AlignmentUtils.findAddableReferenceAnnotations(ourAl.getSequences(),
2659             tipEntries, candidates, ourAl);
2660     AlignmentUtils.addReferenceAnnotations(candidates, ourAl, null);
2661
2662     assertNotNull(ourAl.getAlignmentAnnotation());
2663     assertEquals(ourAl.getAlignmentAnnotation().length, 2);
2664
2665     for (AlignmentAnnotation alan : ourAl.getAlignmentAnnotation())
2666     {
2667       verifyExpectedSequenceAnnotation(alan);
2668     }
2669     // Everything above works for 2.11.3 and 2.11.2.x.
2670     // now simulate copy/paste to new alignment
2671     SequenceI[] newSeqAl = new SequenceI[2];
2672     // copy sequences but no annotation
2673     newSeqAl[0] = new Sequence(ourAl.getSequenceAt(0),
2674             ourAl.getSequenceAt(0).getAnnotation());
2675     newSeqAl[1] = new Sequence(ourAl.getSequenceAt(1),
2676             ourAl.getSequenceAt(1).getAnnotation());
2677
2678     Alignment newAl = new Alignment(newSeqAl);
2679     // delete annotation
2680     for (SequenceI sq : newAl.getSequences())
2681     {
2682       sq.setAlignmentAnnotation(new AlignmentAnnotation[0]);
2683     }
2684     // JAL-4182 scenario test
2685     SequenceGroup sg = new SequenceGroup(Arrays.asList(newSeqAl));
2686     sg.setStartRes(0);
2687     sg.setEndRes(newAl.getWidth());
2688     AlignmentUtils.addReferenceAnnotationTo(newAl, newSeqAl[0],
2689             newSeqAl[0].getDatasetSequence().getAnnotation()[0], sg);
2690     AlignmentUtils.addReferenceAnnotationTo(newAl, newSeqAl[1],
2691             newSeqAl[1].getDatasetSequence().getAnnotation()[0], sg);
2692     for (AlignmentAnnotation alan : newAl.getAlignmentAnnotation())
2693     {
2694       verifyExpectedSequenceAnnotation(alan);
2695     }
2696   }
2697
2698   /**
2699    * helper - tests annotation is mapped to position it was originally created
2700    * for
2701    *
2702    * @param alan
2703    */
2704   private void verifyExpectedSequenceAnnotation(AlignmentAnnotation alan)
2705   {
2706     for (int c = 0; c < alan.annotations.length; c++)
2707     {
2708       Annotation a = alan.annotations[c];
2709       if (a != null)
2710       {
2711         assertEquals("Misaligned annotation at " + c,
2712                 (float) alan.sequenceRef.findPosition(c), a.value);
2713       }
2714       else
2715       {
2716         assertTrue("Unexpected Null at position " + c,
2717                 c >= alan.sequenceRef.getLength()
2718                         || Comparison.isGap(alan.sequenceRef.getCharAt(c)));
2719       }
2720     }
2721   }
2722
2723   @Test(groups = "Functional")
2724   public void testAddReferenceContactMap()
2725   {
2726     SequenceI sq = new Sequence("a", "SSSQ");
2727     ContactMatrixI cm = new SeqDistanceContactMatrix(4);
2728     AlignmentAnnotation cm_aan = sq.addContactList(cm);
2729     cm_aan.description = cm_aan.description + " cm1";
2730     SequenceI dssq = sq.createDatasetSequence();
2731
2732     // remove annotation on our non-dataset sequence
2733     sq.removeAlignmentAnnotation(sq.getAnnotation()[0]);
2734     // test transfer
2735     Alignment al = new Alignment(new SequenceI[] { sq });
2736     SortedMap<String, String> tipEntries = new TreeMap<>();
2737     Map<SequenceI, List<AlignmentAnnotation>> candidates = new LinkedHashMap<>();
2738
2739     AlignmentUtils.findAddableReferenceAnnotations(al.getSequences(),
2740             tipEntries, candidates, al);
2741     AlignmentUtils.addReferenceAnnotations(candidates, al, null);
2742     assertTrue("No contact map annotation transferred",
2743             al.getAlignmentAnnotation() != null
2744                     && al.getAlignmentAnnotation().length == 1);
2745     AlignmentAnnotation alan = al.findAnnotations(sq, null, cm_aan.label)
2746             .iterator().next();
2747     ContactMatrixI t_cm = al.getContactMatrixFor(alan);
2748     assertNotNull("No contact map for the transferred annotation row.",
2749             t_cm);
2750     assertTrue(t_cm instanceof SeqDistanceContactMatrix);
2751     assertTrue(((SeqDistanceContactMatrix) t_cm).hasReferenceSeq());
2752
2753     ContactListI cl = al.getContactListFor(alan, 1);
2754     assertNotNull(
2755             "No contact matrix recovered after reference annotation transfer",
2756             cl);
2757     // semantics of sequence associated contact list is slightly tricky - column
2758     // 3 in alignment should have data
2759     cl = al.getContactListFor(alan, 3);
2760     assertNotNull(
2761             "Contact matrix should have data for last position in sequence",
2762             cl);
2763
2764     ContactMatrixI cm2 = new SeqDistanceContactMatrix(4);
2765     dssq.addContactList(cm2);
2766     tipEntries = new TreeMap<>();
2767     candidates = new LinkedHashMap<>();
2768
2769     AlignmentUtils.findAddableReferenceAnnotations(al.getSequences(),
2770             tipEntries, candidates, al);
2771     AlignmentUtils.addReferenceAnnotations(candidates, al, null);
2772     assertTrue("Expected two contact map annotation transferred",
2773             al.getAlignmentAnnotation() != null
2774                     && al.getAlignmentAnnotation().length == 2);
2775
2776   }
2777
2778     @Test(groups = "Functional", dataProvider = "SecondaryStructureAnnotations")
2779     public void testSecondaryStructurePresentAndSources(AlignmentAnnotation[] annotations, boolean expectedSSPresent, ArrayList<String> expectedSSSources) {
2780         Assert.assertEquals(expectedSSPresent, AlignmentUtils.isSecondaryStructurePresent(annotations));
2781         Assert.assertEquals(expectedSSSources, AlignmentUtils.getSecondaryStructureSources(annotations));
2782     }
2783
2784     @DataProvider(name = "SecondaryStructureAnnotations")
2785     public static Object[][] provideSecondaryStructureAnnotations() {
2786         AlignmentAnnotation ann1 = new AlignmentAnnotation("Secondary Structure", "Secondary Structure", new Annotation[]{});
2787         AlignmentAnnotation ann2 = new AlignmentAnnotation("jnetpred", "jnetpred", new Annotation[]{});
2788         AlignmentAnnotation ann3 = new AlignmentAnnotation("Temp", "Temp", new Annotation[]{});
2789         AlignmentAnnotation ann4 = new AlignmentAnnotation("Temp", "Temp", new Annotation[]{});
2790
2791         List<String> ssSources1 = new ArrayList<>(Arrays.asList("3D Structures"));
2792         List<String> ssSources2 = new ArrayList<>(Arrays.asList("JPred"));
2793         List<String> ssSources3 = new ArrayList<>(Arrays.asList("3D Structures", "JPred"));
2794         List<String> ssSources4 = new ArrayList<>();
2795
2796         return new Object[][]{
2797             {new AlignmentAnnotation[]{ann1, ann3, ann4}, true, ssSources1},
2798             {new AlignmentAnnotation[]{ann2, ann3, ann4}, true, ssSources2},
2799             {new AlignmentAnnotation[]{ann3, ann4}, false, ssSources4},
2800             {new AlignmentAnnotation[]{}, false, ssSources4},
2801             {new AlignmentAnnotation[]{ann1, ann2, ann3, ann4}, true, ssSources3}
2802         };
2803     }
2804
2805     @Test(dataProvider = "SecondaryStructureAnnotationColours")
2806     public void testSecondaryStructureAnnotationColour(char symbol, Color expectedColor) {
2807         Color actualColor = AlignmentUtils.getSecondaryStructureAnnotationColour(symbol);
2808         Assert.assertEquals(actualColor, expectedColor);
2809     }
2810
2811     @DataProvider(name = "SecondaryStructureAnnotationColours")
2812     public static Object[][] provideSecondaryStructureAnnotationColours() {
2813         return new Object[][]{
2814             {'C', Color.gray},
2815             {'E', Color.green},
2816             {'H', Color.red},
2817             {'-', Color.gray}
2818         };
2819     }
2820
2821     @Test(dataProvider = "SSAnnotationPresence")
2822     public void testIsSSAnnotationPresent(Map<SequenceI, List<AlignmentAnnotation>> annotations, boolean expectedPresence) {
2823         boolean actualPresence = AlignmentUtils.isSSAnnotationPresent(annotations);
2824         Assert.assertEquals(actualPresence, expectedPresence);
2825     }
2826
2827     @DataProvider(name = "SSAnnotationPresence")
2828     public static Object[][] provideSSAnnotationPresence() {
2829         Map<SequenceI, List<AlignmentAnnotation>> annotations1 = new HashMap<>();
2830         SequenceI seq1 = new Sequence("Seq1", "ASD---ASD---ASD", 37, 45);
2831         List<AlignmentAnnotation> annotationsList1 = new ArrayList<>();
2832         annotationsList1.add(new AlignmentAnnotation("Secondary Structure", "Secondary Structure", new Annotation[]{}));
2833         annotations1.put(seq1, annotationsList1); // Annotation present secondary structure for seq1
2834
2835         Map<SequenceI, List<AlignmentAnnotation>> annotations2 = new HashMap<>();
2836         SequenceI seq2 = new Sequence("Seq2", "ASD---ASD------", 37, 42);
2837         List<AlignmentAnnotation> annotationsList2 = new ArrayList<>();
2838         annotationsList2.add(new AlignmentAnnotation("Other Annotation", "Other Annotation", new Annotation[]{}));
2839         annotations2.put(seq2, annotationsList2); // Annotation not related to any of secondary structure for seq2
2840
2841         Map<SequenceI, List<AlignmentAnnotation>> annotations3 = new HashMap<>();
2842         // Empty annotation map
2843
2844         Map<SequenceI, List<AlignmentAnnotation>> annotations4 = new HashMap<>();
2845         SequenceI seq4 = new Sequence("Seq4", "ASD---ASD---AS-", 37, 44);
2846         List<AlignmentAnnotation> annotationsList4 = new ArrayList<>();
2847         annotationsList4.add(new AlignmentAnnotation("jnetpred", "jnetpred", new Annotation[]{}));
2848         annotations4.put(seq4, annotationsList4); // Annotation present from JPred for seq4
2849
2850
2851         return new Object[][]{
2852             {annotations1, true}, // Annotations present secondary structure present
2853             {annotations2, false}, // No annotations related to any of the secondary structure present
2854             {annotations3, false},  // Empty annotation map
2855             {annotations4, true}, // Annotations present from JPred secondary structure present
2856         };
2857     }
2858 //
2859 //    @Test(dataProvider = "SSSourceFromAnnotationDescription")
2860 //    public void testGetSSSourceFromAnnotationDescription(Map<SequenceI, List<AlignmentAnnotation>> annotations, String expectedSSSource) {
2861 //        String actualSSSource = AlignmentUtils.getSSSourceFromAnnotationDescription(annotations);
2862 //        Assert.assertEquals(actualSSSource, expectedSSSource);
2863 //    }
2864
2865     @DataProvider(name = "SSSourceFromAnnotationDescription")
2866     public static Object[][] provideSSSourceFromAnnotationDescription() {
2867         Map<SequenceI, List<AlignmentAnnotation>> annotations1 = new HashMap<>();
2868         SequenceI seq1 = new Sequence("Seq1", "ASD---ASD---ASD", 37, 45);
2869         List<AlignmentAnnotation> annotationsList1 = new ArrayList<>();
2870         annotationsList1.add(new AlignmentAnnotation("jnetpred", "JPred Output", new Annotation[]{}));
2871         annotations1.put(seq1, annotationsList1); // Annotation present from JPred for seq1
2872
2873         Map<SequenceI, List<AlignmentAnnotation>> annotations2 = new HashMap<>();
2874         SequenceI seq2 = new Sequence("Seq2", "ASD---ASD------", 37, 42);
2875         List<AlignmentAnnotation> annotationsList2 = new ArrayList<>();
2876         annotationsList2.add(new AlignmentAnnotation("Secondary Structure",
2877                 "Secondary Structure for af-q43517-f1A", new Annotation[]{}));
2878         annotations2.put(seq2, annotationsList2); // Annotation present secondary structure from Alphafold for seq2
2879
2880         Map<SequenceI, List<AlignmentAnnotation>> annotations3 = new HashMap<>();
2881         // Empty annotation map
2882
2883         Map<SequenceI, List<AlignmentAnnotation>> annotations4 = new HashMap<>();
2884         SequenceI seq4 = new Sequence("Seq4", "ASD---ASD---AS-", 37, 44);
2885         List<AlignmentAnnotation> annotationsList4 = new ArrayList<>();
2886         annotationsList4.add(new AlignmentAnnotation("Secondary Structure",
2887                 "Secondary Structure for 4zhpA", new Annotation[]{}));
2888         annotations4.put(seq4, annotationsList4); // Annotation present secondary structure from pdb for seq4
2889
2890         Map<SequenceI, List<AlignmentAnnotation>> annotations5 = new HashMap<>();
2891         SequenceI seq5 = new Sequence("Seq5", "ASD---ASD---AS-", 37, 44);
2892         List<AlignmentAnnotation> annotationsList5 = new ArrayList<>();
2893         annotationsList5.add(new AlignmentAnnotation("Secondary Structure",
2894                 "Secondary Structure for p09911_54-147__3a7wzn.1.p3502557454997462030P",
2895                 new Annotation[]{}));
2896         annotations5.put(seq5, annotationsList5); // Annotation present secondary structure from Swiss model for seq5
2897
2898
2899         //JPred Output - JPred
2900         //Secondary Structure for af-q43517-f1A - Alphafold
2901         //Secondary Structure for 4zhpA - Experimental
2902         //Secondary Structure for p09911_54-147__3a7wzn.1.p3502557454997462030P - Swiss Model
2903
2904         return new Object[][]{
2905             {annotations1, "JPred"},
2906             {annotations2, "Alphafold"},
2907             {annotations3, null},
2908             {annotations4, "PDB"},
2909             {annotations5, "Swiss Model"}
2910         };
2911     }
2912
2913 }