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.junit.Assert.assertEquals;
  24 import static org.junit.Assert.assertFalse;
  25 import static org.junit.Assert.assertSame;
  26 import static org.junit.Assert.assertTrue;
  27
  28 import java.io.IOException;
  29 import java.util.ArrayList;
  30 import java.util.Arrays;
  31 import java.util.Collections;
  32 import java.util.HashSet;
  33 import java.util.List;
  34 import java.util.Map;
  35 import java.util.Set;
  36
  37 import org.junit.Test;
  38
  39 import jalview.datamodel.AlignedCodonFrame;
  40 import jalview.datamodel.Alignment;
  41 import jalview.datamodel.AlignmentAnnotation;
  42 import jalview.datamodel.AlignmentI;
  43 import jalview.datamodel.Annotation;
  44 import jalview.datamodel.DBRefEntry;
  45 import jalview.datamodel.Mapping;
  46 import jalview.datamodel.SearchResults;
  47 import jalview.datamodel.SearchResults.Match;
  48 import jalview.datamodel.Sequence;
  49 import jalview.datamodel.SequenceI;
  50 import jalview.io.AppletFormatAdapter;
  51 import jalview.io.FormatAdapter;
  52 import jalview.util.MapList;
  53 import jalview.util.MappingUtils;
  54
  55 public class AlignmentUtilsTests
  56 {
  57   // @formatter:off
  58   private static final String TEST_DATA =
  59           "# STOCKHOLM 1.0\n" +
  60           "#=GS D.melanogaster.1 AC AY119185.1/838-902\n" +
  61           "#=GS D.melanogaster.2 AC AC092237.1/57223-57161\n" +
  62           "#=GS D.melanogaster.3 AC AY060611.1/560-627\n" +
  63           "D.melanogaster.1          G.AGCC.CU...AUGAUCGA\n" +
  64           "#=GR D.melanogaster.1 SS  ................((((\n" +
  65           "D.melanogaster.2          C.AUUCAACU.UAUGAGGAU\n" +
  66           "#=GR D.melanogaster.2 SS  ................((((\n" +
  67           "D.melanogaster.3          G.UGGCGCU..UAUGACGCA\n" +
  68           "#=GR D.melanogaster.3 SS  (.(((...(....(((((((\n" +
  69           "//";
  70
  71   private static final String AA_SEQS_1 =
  72           ">Seq1Name\n" +
  73           "K-QY--L\n" +
  74           ">Seq2Name\n" +
  75           "-R-FP-W-\n";
  76
  77   private static final String CDNA_SEQS_1 =
  78           ">Seq1Name\n" +
  79           "AC-GG--CUC-CAA-CT\n" +
  80           ">Seq2Name\n" +
  81           "-CG-TTA--ACG---AAGT\n";
  82
  83   private static final String CDNA_SEQS_2 =
  84           ">Seq1Name\n" +
  85           "GCTCGUCGTACT\n" +
  86           ">Seq2Name\n" +
  87           "GGGTCAGGCAGT\n";
  88   // @formatter:on
  89
  90   public static Sequence ts=new Sequence("short","ASDASDASDASDASDASDASDASDASDASDASDASDASD");
  91
  92   @Test
  93   public void testExpandFlanks()
  94   {
  95     AlignmentI al = new Alignment(new Sequence[] {});
  96     for (int i=4;i<14;i+=3)
  97     {
  98       SequenceI s1=ts.deriveSequence().getSubSequence(i, i+7);
  99       al.addSequence(s1);
 100     }
 101     System.out.println(new AppletFormatAdapter().formatSequences("Clustal", al, true));
 102     for (int flnk=-1;flnk<25; flnk++)
 103     {
 104       AlignmentI exp;
 105       System.out.println("\nFlank size: "+flnk);
 106       System.out.println(new AppletFormatAdapter().formatSequences("Clustal", exp=AlignmentUtils.expandContext(al, flnk), true));
 107       if (flnk==-1) {
 108         for (SequenceI sq:exp.getSequences())
 109       {
 110           String ung = sq.getSequenceAsString().replaceAll("-+", "");
 111           assertTrue("Flanking sequence not the same as original dataset sequence.\n"+ung+"\n"+sq.getDatasetSequence().getSequenceAsString(),ung.equalsIgnoreCase(sq.getDatasetSequence().getSequenceAsString()));
 112       }
 113       }
 114     }
 115     }
 116
 117   /**
 118    * Test method that returns a map of lists of sequences by sequence name.
 119    *
 120    * @throws IOException
 121    */
 122   @Test
 123   public void testGetSequencesByName() throws IOException
 124   {
 125     final String data = ">Seq1Name\nKQYL\n" + ">Seq2Name\nRFPW\n"
 126             + ">Seq1Name\nABCD\n";
 127     AlignmentI al = loadAlignment(data, "FASTA");
 128     Map<String, List<SequenceI>> map = AlignmentUtils
 129             .getSequencesByName(al);
 130     assertEquals(2, map.keySet().size());
 131     assertEquals(2, map.get("Seq1Name").size());
 132     assertEquals("KQYL", map.get("Seq1Name").get(0).getSequenceAsString());
 133     assertEquals("ABCD", map.get("Seq1Name").get(1).getSequenceAsString());
 134     assertEquals(1, map.get("Seq2Name").size());
 135     assertEquals("RFPW", map.get("Seq2Name").get(0).getSequenceAsString());
 136   }
 137   /**
 138    * Helper method to load an alignment and ensure dataset sequences are set up.
 139    *
 140    * @param data
 141    * @param format TODO
 142    * @return
 143    * @throws IOException
 144    */
 145   protected AlignmentI loadAlignment(final String data, String format) throws IOException
 146   {
 147     Alignment a = new FormatAdapter().readFile(data,
 148             AppletFormatAdapter.PASTE, format);
 149     a.setDataset(null);
 150     return a;
 151   }
 152
 153   /**
 154    * Test mapping of protein to cDNA, for the case where we have no sequence
 155    * cross-references, so mappings are made first-served 1-1 where sequences
 156    * translate.
 157    *
 158    * @throws IOException
 159    */
 160   @Test
 161   public void testMapProteinToCdna_noXrefs() throws IOException
 162   {
 163     List<SequenceI> protseqs = new ArrayList<SequenceI>();
 164     protseqs.add(new Sequence("UNIPROT|V12345", "EIQ"));
 165     protseqs.add(new Sequence("UNIPROT|V12346", "EIQ"));
 166     protseqs.add(new Sequence("UNIPROT|V12347", "SAR"));
 167     AlignmentI protein = new Alignment(protseqs.toArray(new SequenceI[3]));
 168     protein.setDataset(null);
 169
 170     List<SequenceI> dnaseqs = new ArrayList<SequenceI>();
 171     dnaseqs.add(new Sequence("EMBL|A11111", "TCAGCACGC")); // = SAR
 172     dnaseqs.add(new Sequence("EMBL|A22222", "GAGATACAA")); // = EIQ
 173     dnaseqs.add(new Sequence("EMBL|A33333", "GAAATCCAG")); // = EIQ
 174     dnaseqs.add(new Sequence("EMBL|A44444", "GAAATTCAG")); // = EIQ
 175     AlignmentI cdna = new Alignment(dnaseqs.toArray(new SequenceI[4]));
 176     cdna.setDataset(null);
 177
 178     assertTrue(AlignmentUtils.mapProteinToCdna(protein, cdna));
 179
 180     // 3 mappings made, each from 1 to 1 sequence
 181     assertEquals(3, protein.getCodonFrames().size());
 182     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(0)).size());
 183     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(1)).size());
 184     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(2)).size());
 185
 186     // V12345 mapped to A22222
 187     AlignedCodonFrame acf = protein.getCodonFrame(
 188             protein.getSequenceAt(0)).get(0);
 189     assertEquals(1, acf.getdnaSeqs().length);
 190     assertEquals(cdna.getSequenceAt(1).getDatasetSequence(),
 191             acf.getdnaSeqs()[0]);
 192     Mapping[] protMappings = acf.getProtMappings();
 193     assertEquals(1, protMappings.length);
 194     MapList mapList = protMappings[0].getMap();
 195     assertEquals(3, mapList.getFromRatio());
 196     assertEquals(1, mapList.getToRatio());
 197     assertTrue(Arrays.equals(new int[]
 198     { 1, 9 }, mapList.getFromRanges().get(0)));
 199     assertEquals(1, mapList.getFromRanges().size());
 200     assertTrue(Arrays.equals(new int[]
 201     { 1, 3 }, mapList.getToRanges().get(0)));
 202     assertEquals(1, mapList.getToRanges().size());
 203
 204     // V12346 mapped to A33333
 205     acf = protein.getCodonFrame(protein.getSequenceAt(1)).get(0);
 206     assertEquals(1, acf.getdnaSeqs().length);
 207     assertEquals(cdna.getSequenceAt(2).getDatasetSequence(),
 208             acf.getdnaSeqs()[0]);
 209
 210     // V12347 mapped to A11111
 211     acf = protein.getCodonFrame(protein.getSequenceAt(2)).get(0);
 212     assertEquals(1, acf.getdnaSeqs().length);
 213     assertEquals(cdna.getSequenceAt(0).getDatasetSequence(),
 214             acf.getdnaSeqs()[0]);
 215
 216     // no mapping involving the 'extra' A44444
 217     assertTrue(protein.getCodonFrame(cdna.getSequenceAt(3)).isEmpty());
 218   }
 219
 220   /**
 221    * Test for the alignSequenceAs method that takes two sequences and a mapping.
 222    */
 223   @Test
 224   public void testAlignSequenceAs_withMapping_noIntrons()
 225   {
 226     MapList map = new MapList(new int[]
 227     { 1, 6 }, new int[]
 228     { 1, 2 }, 3, 1);
 229
 230     /*
 231      * No existing gaps in dna:
 232      */
 233     checkAlignSequenceAs("GGGAAA", "-A-L-", false, false, map,
 234             "---GGG---AAA");
 235
 236     /*
 237      * Now introduce gaps in dna but ignore them when realigning.
 238      */
 239     checkAlignSequenceAs("-G-G-G-A-A-A-", "-A-L-", false, false, map,
 240             "---GGG---AAA");
 241
 242     /*
 243      * Now include gaps in dna when realigning. First retaining 'mapped' gaps
 244      * only, i.e. those within the exon region.
 245      */
 246     checkAlignSequenceAs("-G-G--G-A--A-A-", "-A-L-", true, false, map,
 247             "---G-G--G---A--A-A");
 248
 249     /*
 250      * Include all gaps in dna when realigning (within and without the exon
 251      * region). The leading gap, and the gaps between codons, are subsumed by
 252      * the protein alignment gap.
 253      */
 254     checkAlignSequenceAs("-G-GG--AA-A-", "-A-L-", true, true, map,
 255             "---G-GG---AA-A-");
 256
 257     /*
 258      * Include only unmapped gaps in dna when realigning (outside the exon
 259      * region). The leading gap, and the gaps between codons, are subsumed by
 260      * the protein alignment gap.
 261      */
 262     checkAlignSequenceAs("-G-GG--AA-A-", "-A-L-", false, true, map,
 263             "---GGG---AAA-");
 264   }
 265
 266   /**
 267    * Test for the alignSequenceAs method that takes two sequences and a mapping.
 268    */
 269   @Test
 270   public void testAlignSequenceAs_withMapping_withIntrons()
 271   {
 272     /*
 273      * Exons at codon 2 (AAA) and 4 (TTT)
 274      */
 275     MapList map = new MapList(new int[]
 276     { 4, 6, 10, 12 }, new int[]
 277     { 1, 2 }, 3, 1);
 278
 279     /*
 280      * Simple case: no gaps in dna
 281      */
 282     checkAlignSequenceAs("GGGAAACCCTTTGGG", "--A-L-", false, false, map,
 283             "GGG---AAACCCTTTGGG");
 284
 285     /*
 286      * Add gaps to dna - but ignore when realigning.
 287      */
 288     checkAlignSequenceAs("-G-G-G--A--A---AC-CC-T-TT-GG-G-", "--A-L-",
 289             false, false, map, "GGG---AAACCCTTTGGG");
 290
 291     /*
 292      * Add gaps to dna - include within exons only when realigning.
 293      */
 294     checkAlignSequenceAs("-G-G-G--A--A---A-C-CC-T-TT-GG-G-", "--A-L-",
 295             true, false, map, "GGG---A--A---ACCCT-TTGGG");
 296
 297     /*
 298      * Include gaps outside exons only when realigning.
 299      */
 300     checkAlignSequenceAs("-G-G-G--A--A---A-C-CC-T-TT-GG-G-", "--A-L-",
 301             false, true, map, "-G-G-GAAAC-CCTTT-GG-G-");
 302
 303     /*
 304      * Include gaps following first intron if we are 'preserving mapped gaps'
 305      */
 306     checkAlignSequenceAs("-G-G-G--A--A---A-C-CC-T-TT-GG-G-", "--A-L-",
 307             true, true, map, "-G-G-G--A--A---A-C-CC-T-TT-GG-G-");
 308
 309     /*
 310      * Include all gaps in dna when realigning.
 311      */
 312     checkAlignSequenceAs("-G-G-G--A--A---A-C-CC-T-TT-GG-G-", "--A-L-",
 313             true, true, map, "-G-G-G--A--A---A-C-CC-T-TT-GG-G-");
 314   }
 315
 316   /**
 317    * Test for the case where not all of the protein sequence is mapped to cDNA.
 318    */
 319   @Test
 320   public void testAlignSequenceAs_withMapping_withUnmappedProtein()
 321   {
 322
 323     /*
 324      * Exons at codon 2 (AAA) and 4 (TTT) mapped to A and P
 325      */
 326     final MapList map = new MapList(new int[]
 327     { 4, 6, 10, 12 }, new int[]
 328     { 1, 1, 3, 3 }, 3, 1);
 329
 330
 331     /*
 332      * Expect alignment does nothing (aborts realignment). Change this test
 333      * first if different behaviour wanted.
 334      */
 335     checkAlignSequenceAs("GGGAAACCCTTTGGG", "-A-L-P-", false,
 336             false, map, "GGGAAACCCTTTGGG");
 337   }
 338
 339   /**
 340    * Helper method that performs and verifies the method under test.
 341    *
 342    * @param dnaSeq
 343    * @param proteinSeq
 344    * @param preserveMappedGaps
 345    * @param preserveUnmappedGaps
 346    * @param map
 347    * @param expected
 348    */
 349   protected void checkAlignSequenceAs(final String dnaSeq,
 350           final String proteinSeq, final boolean preserveMappedGaps,
 351           final boolean preserveUnmappedGaps, MapList map,
 352           final String expected)
 353   {
 354     SequenceI dna = new Sequence("Seq1", dnaSeq);
 355     dna.createDatasetSequence();
 356     SequenceI protein = new Sequence("Seq1", proteinSeq);
 357     protein.createDatasetSequence();
 358     AlignedCodonFrame acf = new AlignedCodonFrame();
 359     acf.addMap(dna.getDatasetSequence(), protein.getDatasetSequence(), map);
 360
 361     AlignmentUtils.alignSequenceAs(dna, protein, acf, "---", '-',
 362             preserveMappedGaps, preserveUnmappedGaps);
 363     assertEquals(expected, dna.getSequenceAsString());
 364   }
 365
 366   /**
 367    * Test for the alignSequenceAs method where we preserve gaps in introns only.
 368    */
 369   @Test
 370   public void testAlignSequenceAs_keepIntronGapsOnly()
 371   {
 372
 373     /*
 374      * Intron GGGAAA followed by exon CCCTTT
 375      */
 376     MapList map = new MapList(new int[]
 377     { 7, 12 }, new int[]
 378     { 1, 2 }, 3, 1);
 379
 380     checkAlignSequenceAs("GG-G-AA-A-C-CC-T-TT", "AL",
 381             false, true, map, "GG-G-AA-ACCCTTT");
 382   }
 383
 384   /**
 385    * Test for the method that generates an aligned translated sequence from one
 386    * mapping.
 387    */
 388   @Test
 389   public void testGetAlignedTranslation_dnaLikeProtein()
 390   {
 391     // dna alignment will be replaced
 392     SequenceI dna = new Sequence("Seq1", "T-G-CC-A--T-TAC-CAG-");
 393     dna.createDatasetSequence();
 394     // protein alignment will be 'applied' to dna
 395     SequenceI protein = new Sequence("Seq1", "-CH-Y--Q-");
 396     protein.createDatasetSequence();
 397     MapList map = new MapList(new int[]
 398     { 1, 12 }, new int[]
 399     { 1, 4 }, 3, 1);
 400     AlignedCodonFrame acf = new AlignedCodonFrame();
 401     acf.addMap(dna.getDatasetSequence(), protein.getDatasetSequence(), map);
 402
 403     final SequenceI aligned = AlignmentUtils
 404                 .getAlignedTranslation(protein, '-', acf);
 405     assertEquals("---TGCCAT---TAC------CAG---", aligned.getSequenceAsString());
 406     assertSame(aligned.getDatasetSequence(), dna.getDatasetSequence());
 407   }
 408
 409   /**
 410    * Test the method that realigns protein to match mapped codon alignment.
 411    */
 412   @Test
 413   public void testAlignProteinAsDna()
 414   {
 415     // seq1 codons are [1,2,3] [4,5,6] [7,8,9] [10,11,12]
 416     SequenceI dna1 = new Sequence("Seq1", "TGCCATTACCAG-");
 417     // seq2 codons are [1,3,4] [5,6,7] [8,9,10] [11,12,13]
 418     SequenceI dna2 = new Sequence("Seq2", "T-GCCATTACCAG");
 419     // seq3 codons are [1,2,3] [4,5,7] [8,9,10] [11,12,13]
 420     SequenceI dna3 = new Sequence("Seq3", "TGCCA-TTACCAG");
 421     AlignmentI dna = new Alignment(new SequenceI[]
 422     { dna1, dna2, dna3 });
 423     dna.setDataset(null);
 424
 425     // protein alignment will be realigned like dna
 426     SequenceI prot1 = new Sequence("Seq1", "CHYQ");
 427     SequenceI prot2 = new Sequence("Seq2", "CHYQ");
 428     SequenceI prot3 = new Sequence("Seq3", "CHYQ");
 429     AlignmentI protein = new Alignment(new SequenceI[]
 430     { prot1, prot2, prot3 });
 431     protein.setDataset(null);
 432
 433     MapList map = new MapList(new int[]
 434     { 1, 12 }, new int[]
 435     { 1, 4 }, 3, 1);
 436     AlignedCodonFrame acf = new AlignedCodonFrame();
 437     acf.addMap(dna1.getDatasetSequence(), prot1.getDatasetSequence(), map);
 438     acf.addMap(dna2.getDatasetSequence(), prot2.getDatasetSequence(), map);
 439     acf.addMap(dna3.getDatasetSequence(), prot3.getDatasetSequence(), map);
 440     protein.setCodonFrames(Collections.singleton(acf));
 441
 442     /*
 443      * Translated codon order is [1,2,3] [1,3,4] [4,5,6] [4,5,7] [5,6,7] [7,8,9]
 444      * [8,9,10] [10,11,12] [11,12,13]
 445      */
 446     AlignmentUtils.alignProteinAsDna(protein, dna);
 447     assertEquals("C-H--Y-Q-", prot1.getSequenceAsString());
 448     assertEquals("-C--H-Y-Q", prot2.getSequenceAsString());
 449     assertEquals("C--H--Y-Q", prot3.getSequenceAsString());
 450   }
 451
 452   /**
 453    * Test the method that tests whether a CDNA sequence translates to a protein
 454    * sequence
 455    */
 456   @Test
 457   public void testTranslatesAs()
 458   {
 459     assertTrue(AlignmentUtils.translatesAs("tttcccaaaggg".toCharArray(), 0,
 460             "FPKG".toCharArray()));
 461     // with start codon
 462     assertTrue(AlignmentUtils.translatesAs("atgtttcccaaaggg".toCharArray(),
 463             3, "FPKG".toCharArray()));
 464     // with stop codon1
 465     assertTrue(AlignmentUtils.translatesAs("tttcccaaagggtaa".toCharArray(),
 466             0, "FPKG".toCharArray()));
 467     // with stop codon2
 468     assertTrue(AlignmentUtils.translatesAs("tttcccaaagggtag".toCharArray(),
 469             0, "FPKG".toCharArray()));
 470     // with stop codon3
 471     assertTrue(AlignmentUtils.translatesAs("tttcccaaagggtga".toCharArray(),
 472             0, "FPKG".toCharArray()));
 473     // with start and stop codon1
 474     assertTrue(AlignmentUtils.translatesAs(
 475             "atgtttcccaaaggtaa".toCharArray(), 3, "FPKG".toCharArray()));
 476     // with start and stop codon2
 477     assertTrue(AlignmentUtils.translatesAs(
 478             "atgtttcccaaaggtag".toCharArray(), 3, "FPKG".toCharArray()));
 479     // with start and stop codon3
 480     assertTrue(AlignmentUtils.translatesAs(
 481             "atgtttcccaaaggtga".toCharArray(), 3, "FPKG".toCharArray()));
 482
 483     // wrong protein
 484     assertFalse(AlignmentUtils.translatesAs("tttcccaaaggg".toCharArray(),
 485             0,
 486             "FPMG".toCharArray()));
 487   }
 488
 489   /**
 490    * Test mapping of protein to cDNA, for cases where the cDNA has start and/or
 491    * stop codons in addition to the protein coding sequence.
 492    *
 493    * @throws IOException
 494    */
 495   @Test
 496   public void testMapProteinToCdna_withStartAndStopCodons()
 497           throws IOException
 498   {
 499     List<SequenceI> protseqs = new ArrayList<SequenceI>();
 500     protseqs.add(new Sequence("UNIPROT|V12345", "EIQ"));
 501     protseqs.add(new Sequence("UNIPROT|V12346", "EIQ"));
 502     protseqs.add(new Sequence("UNIPROT|V12347", "SAR"));
 503     AlignmentI protein = new Alignment(protseqs.toArray(new SequenceI[3]));
 504     protein.setDataset(null);
 505
 506     List<SequenceI> dnaseqs = new ArrayList<SequenceI>();
 507     // start + SAR:
 508     dnaseqs.add(new Sequence("EMBL|A11111", "ATGTCAGCACGC"));
 509     // = EIQ + stop
 510     dnaseqs.add(new Sequence("EMBL|A22222", "GAGATACAATAA"));
 511     // = start +EIQ + stop
 512     dnaseqs.add(new Sequence("EMBL|A33333", "ATGGAAATCCAGTAG"));
 513     dnaseqs.add(new Sequence("EMBL|A44444", "GAAATTCAG"));
 514     AlignmentI cdna = new Alignment(dnaseqs.toArray(new SequenceI[4]));
 515     cdna.setDataset(null);
 516
 517     assertTrue(AlignmentUtils.mapProteinToCdna(protein, cdna));
 518
 519     // 3 mappings made, each from 1 to 1 sequence
 520     assertEquals(3, protein.getCodonFrames().size());
 521     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(0)).size());
 522     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(1)).size());
 523     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(2)).size());
 524
 525     // V12345 mapped from A22222
 526     AlignedCodonFrame acf = protein.getCodonFrame(
 527             protein.getSequenceAt(0)).get(0);
 528     assertEquals(1, acf.getdnaSeqs().length);
 529     assertEquals(cdna.getSequenceAt(1).getDatasetSequence(),
 530             acf.getdnaSeqs()[0]);
 531     Mapping[] protMappings = acf.getProtMappings();
 532     assertEquals(1, protMappings.length);
 533     MapList mapList = protMappings[0].getMap();
 534     assertEquals(3, mapList.getFromRatio());
 535     assertEquals(1, mapList.getToRatio());
 536     assertTrue(Arrays.equals(new int[]
 537     { 1, 9 }, mapList.getFromRanges().get(0)));
 538     assertEquals(1, mapList.getFromRanges().size());
 539     assertTrue(Arrays.equals(new int[]
 540     { 1, 3 }, mapList.getToRanges().get(0)));
 541     assertEquals(1, mapList.getToRanges().size());
 542
 543     // V12346 mapped from A33333 starting position 4
 544     acf = protein.getCodonFrame(protein.getSequenceAt(1)).get(0);
 545     assertEquals(1, acf.getdnaSeqs().length);
 546     assertEquals(cdna.getSequenceAt(2).getDatasetSequence(),
 547             acf.getdnaSeqs()[0]);
 548     protMappings = acf.getProtMappings();
 549     assertEquals(1, protMappings.length);
 550     mapList = protMappings[0].getMap();
 551     assertEquals(3, mapList.getFromRatio());
 552     assertEquals(1, mapList.getToRatio());
 553     assertTrue(Arrays.equals(new int[]
 554     { 4, 12 }, mapList.getFromRanges().get(0)));
 555     assertEquals(1, mapList.getFromRanges().size());
 556     assertTrue(Arrays.equals(new int[]
 557     { 1, 3 }, mapList.getToRanges().get(0)));
 558     assertEquals(1, mapList.getToRanges().size());
 559
 560     // V12347 mapped to A11111 starting position 4
 561     acf = protein.getCodonFrame(protein.getSequenceAt(2)).get(0);
 562     assertEquals(1, acf.getdnaSeqs().length);
 563     assertEquals(cdna.getSequenceAt(0).getDatasetSequence(),
 564             acf.getdnaSeqs()[0]);
 565     protMappings = acf.getProtMappings();
 566     assertEquals(1, protMappings.length);
 567     mapList = protMappings[0].getMap();
 568     assertEquals(3, mapList.getFromRatio());
 569     assertEquals(1, mapList.getToRatio());
 570     assertTrue(Arrays.equals(new int[]
 571     { 4, 12 }, mapList.getFromRanges().get(0)));
 572     assertEquals(1, mapList.getFromRanges().size());
 573     assertTrue(Arrays.equals(new int[]
 574     { 1, 3 }, mapList.getToRanges().get(0)));
 575     assertEquals(1, mapList.getToRanges().size());
 576
 577     // no mapping involving the 'extra' A44444
 578     assertTrue(protein.getCodonFrame(cdna.getSequenceAt(3)).isEmpty());
 579   }
 580
 581   /**
 582    * Test mapping of protein to cDNA, for the case where we have some sequence
 583    * cross-references. Verify that 1-to-many mappings are made where
 584    * cross-references exist and sequences are mappable.
 585    *
 586    * @throws IOException
 587    */
 588   @Test
 589   public void testMapProteinToCdna_withXrefs() throws IOException
 590   {
 591     List<SequenceI> protseqs = new ArrayList<SequenceI>();
 592     protseqs.add(new Sequence("UNIPROT|V12345", "EIQ"));
 593     protseqs.add(new Sequence("UNIPROT|V12346", "EIQ"));
 594     protseqs.add(new Sequence("UNIPROT|V12347", "SAR"));
 595     AlignmentI protein = new Alignment(protseqs.toArray(new SequenceI[3]));
 596     protein.setDataset(null);
 597
 598     List<SequenceI> dnaseqs = new ArrayList<SequenceI>();
 599     dnaseqs.add(new Sequence("EMBL|A11111", "TCAGCACGC")); // = SAR
 600     dnaseqs.add(new Sequence("EMBL|A22222", "ATGGAGATACAA")); // = start + EIQ
 601     dnaseqs.add(new Sequence("EMBL|A33333", "GAAATCCAG")); // = EIQ
 602     dnaseqs.add(new Sequence("EMBL|A44444", "GAAATTCAG")); // = EIQ
 603     dnaseqs.add(new Sequence("EMBL|A55555", "GAGATTCAG")); // = EIQ
 604     AlignmentI cdna = new Alignment(dnaseqs.toArray(new SequenceI[5]));
 605     cdna.setDataset(null);
 606
 607     // Xref A22222 to V12345 (should get mapped)
 608     dnaseqs.get(1).addDBRef(new DBRefEntry("UNIPROT", "1", "V12345"));
 609     // Xref V12345 to A44444 (should get mapped)
 610     protseqs.get(0).addDBRef(new DBRefEntry("EMBL", "1", "A44444"));
 611     // Xref A33333 to V12347 (sequence mismatch - should not get mapped)
 612     dnaseqs.get(2).addDBRef(new DBRefEntry("UNIPROT", "1", "V12347"));
 613     // as V12345 is mapped to A22222 and A44444, this leaves V12346 unmapped.
 614     // it should get paired up with the unmapped A33333
 615     // A11111 should be mapped to V12347
 616     // A55555 is spare and has no xref so is not mapped
 617
 618     assertTrue(AlignmentUtils.mapProteinToCdna(protein, cdna));
 619
 620     // 4 protein mappings made for 3 proteins, 2 to V12345, 1 each to V12346/7
 621     assertEquals(3, protein.getCodonFrames().size());
 622     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(0)).size());
 623     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(1)).size());
 624     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(2)).size());
 625
 626     // one mapping for each of the first 4 cDNA sequences
 627     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(0)).size());
 628     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(1)).size());
 629     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(2)).size());
 630     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(3)).size());
 631
 632     // V12345 mapped to A22222 and A44444
 633     AlignedCodonFrame acf = protein.getCodonFrame(
 634             protein.getSequenceAt(0)).get(0);
 635     assertEquals(2, acf.getdnaSeqs().length);
 636     assertEquals(cdna.getSequenceAt(1).getDatasetSequence(),
 637             acf.getdnaSeqs()[0]);
 638     assertEquals(cdna.getSequenceAt(3).getDatasetSequence(),
 639             acf.getdnaSeqs()[1]);
 640
 641     // V12346 mapped to A33333
 642     acf = protein.getCodonFrame(protein.getSequenceAt(1)).get(0);
 643     assertEquals(1, acf.getdnaSeqs().length);
 644     assertEquals(cdna.getSequenceAt(2).getDatasetSequence(),
 645             acf.getdnaSeqs()[0]);
 646
 647     // V12347 mapped to A11111
 648     acf = protein.getCodonFrame(protein.getSequenceAt(2)).get(0);
 649     assertEquals(1, acf.getdnaSeqs().length);
 650     assertEquals(cdna.getSequenceAt(0).getDatasetSequence(),
 651             acf.getdnaSeqs()[0]);
 652
 653     // no mapping involving the 'extra' A55555
 654     assertTrue(protein.getCodonFrame(cdna.getSequenceAt(4)).isEmpty());
 655   }
 656
 657   /**
 658    * Test mapping of protein to cDNA, for the case where we have some sequence
 659    * cross-references. Verify that once we have made an xref mapping we don't
 660    * also map un-xrefd sequeces.
 661    *
 662    * @throws IOException
 663    */
 664   @Test
 665   public void testMapProteinToCdna_prioritiseXrefs() throws IOException
 666   {
 667     List<SequenceI> protseqs = new ArrayList<SequenceI>();
 668     protseqs.add(new Sequence("UNIPROT|V12345", "EIQ"));
 669     protseqs.add(new Sequence("UNIPROT|V12346", "EIQ"));
 670     AlignmentI protein = new Alignment(
 671             protseqs.toArray(new SequenceI[protseqs.size()]));
 672     protein.setDataset(null);
 673
 674     List<SequenceI> dnaseqs = new ArrayList<SequenceI>();
 675     dnaseqs.add(new Sequence("EMBL|A11111", "GAAATCCAG")); // = EIQ
 676     dnaseqs.add(new Sequence("EMBL|A22222", "GAAATTCAG")); // = EIQ
 677     AlignmentI cdna = new Alignment(dnaseqs.toArray(new SequenceI[dnaseqs
 678             .size()]));
 679     cdna.setDataset(null);
 680
 681     // Xref A22222 to V12345 (should get mapped)
 682     // A11111 should then be mapped to the unmapped V12346
 683     dnaseqs.get(1).addDBRef(new DBRefEntry("UNIPROT", "1", "V12345"));
 684
 685     assertTrue(AlignmentUtils.mapProteinToCdna(protein, cdna));
 686
 687     // 2 protein mappings made
 688     assertEquals(2, protein.getCodonFrames().size());
 689     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(0)).size());
 690     assertEquals(1, protein.getCodonFrame(protein.getSequenceAt(1)).size());
 691
 692     // one mapping for each of the cDNA sequences
 693     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(0)).size());
 694     assertEquals(1, protein.getCodonFrame(cdna.getSequenceAt(1)).size());
 695
 696     // V12345 mapped to A22222
 697     AlignedCodonFrame acf = protein.getCodonFrame(protein.getSequenceAt(0))
 698             .get(0);
 699     assertEquals(1, acf.getdnaSeqs().length);
 700     assertEquals(cdna.getSequenceAt(1).getDatasetSequence(),
 701             acf.getdnaSeqs()[0]);
 702
 703     // V12346 mapped to A11111
 704     acf = protein.getCodonFrame(protein.getSequenceAt(1)).get(0);
 705     assertEquals(1, acf.getdnaSeqs().length);
 706     assertEquals(cdna.getSequenceAt(0).getDatasetSequence(),
 707             acf.getdnaSeqs()[0]);
 708   }
 709
 710   /**
 711    * Test the method that shows or hides sequence annotations by type(s) and
 712    * selection group.
 713    */
 714   @Test
 715   public void testShowOrHideSequenceAnnotations()
 716   {
 717     SequenceI seq1 = new Sequence("Seq1", "AAA");
 718     SequenceI seq2 = new Sequence("Seq2", "BBB");
 719     SequenceI seq3 = new Sequence("Seq3", "CCC");
 720     Annotation[] anns = new Annotation[]
 721     { new Annotation(2f) };
 722     AlignmentAnnotation ann1 = new AlignmentAnnotation("Structure", "ann1",
 723             anns);
 724     ann1.setSequenceRef(seq1);
 725     AlignmentAnnotation ann2 = new AlignmentAnnotation("Structure", "ann2",
 726             anns);
 727     ann2.setSequenceRef(seq2);
 728     AlignmentAnnotation ann3 = new AlignmentAnnotation("Structure", "ann3",
 729             anns);
 730     AlignmentAnnotation ann4 = new AlignmentAnnotation("Temp", "ann4", anns);
 731     ann4.setSequenceRef(seq1);
 732     AlignmentAnnotation ann5 = new AlignmentAnnotation("Temp", "ann5", anns);
 733     ann5.setSequenceRef(seq2);
 734     AlignmentAnnotation ann6 = new AlignmentAnnotation("Temp", "ann6", anns);
 735     AlignmentI al = new Alignment(new SequenceI[] {seq1, seq2, seq3});
 736     al.addAnnotation(ann1); // Structure for Seq1
 737     al.addAnnotation(ann2); // Structure for Seq2
 738     al.addAnnotation(ann3); // Structure for no sequence
 739     al.addAnnotation(ann4); // Temp for seq1
 740     al.addAnnotation(ann5); // Temp for seq2
 741     al.addAnnotation(ann6); // Temp for no sequence
 742     List<String> types = new ArrayList<String>();
 743     List<SequenceI> scope = new ArrayList<SequenceI>();
 744
 745     /*
 746      * Set all sequence related Structure to hidden (ann1, ann2)
 747      */
 748     types.add("Structure");
 749     AlignmentUtils.showOrHideSequenceAnnotations(al, types, null, false,
 750             false);
 751     assertFalse(ann1.visible);
 752     assertFalse(ann2.visible);
 753     assertTrue(ann3.visible); // not sequence-related, not affected
 754     assertTrue(ann4.visible); // not Structure, not affected
 755     assertTrue(ann5.visible); // "
 756     assertTrue(ann6.visible); // not sequence-related, not affected
 757
 758     /*
 759      * Set Temp in {seq1, seq3} to hidden
 760      */
 761     types.clear();
 762     types.add("Temp");
 763     scope.add(seq1);
 764     scope.add(seq3);
 765     AlignmentUtils.showOrHideSequenceAnnotations(al, types, scope, false,
 766             false);
 767     assertFalse(ann1.visible); // unchanged
 768     assertFalse(ann2.visible); // unchanged
 769     assertTrue(ann3.visible); // not sequence-related, not affected
 770     assertFalse(ann4.visible); // Temp for seq1 hidden
 771     assertTrue(ann5.visible); // not in scope, not affected
 772     assertTrue(ann6.visible); // not sequence-related, not affected
 773
 774     /*
 775      * Set Temp in all sequences to hidden
 776      */
 777     types.clear();
 778     types.add("Temp");
 779     scope.add(seq1);
 780     scope.add(seq3);
 781     AlignmentUtils.showOrHideSequenceAnnotations(al, types, null, false,
 782             false);
 783     assertFalse(ann1.visible); // unchanged
 784     assertFalse(ann2.visible); // unchanged
 785     assertTrue(ann3.visible); // not sequence-related, not affected
 786     assertFalse(ann4.visible); // Temp for seq1 hidden
 787     assertFalse(ann5.visible); // Temp for seq2 hidden
 788     assertTrue(ann6.visible); // not sequence-related, not affected
 789
 790     /*
 791      * Set all types in {seq1, seq3} to visible
 792      */
 793     types.clear();
 794     scope.clear();
 795     scope.add(seq1);
 796     scope.add(seq3);
 797     AlignmentUtils.showOrHideSequenceAnnotations(al, types, scope, true,
 798             true);
 799     assertTrue(ann1.visible); // Structure for seq1 set visible
 800     assertFalse(ann2.visible); // not in scope, unchanged
 801     assertTrue(ann3.visible); // not sequence-related, not affected
 802     assertTrue(ann4.visible); // Temp for seq1 set visible
 803     assertFalse(ann5.visible); // not in scope, unchanged
 804     assertTrue(ann6.visible); // not sequence-related, not affected
 805
 806     /*
 807      * Set all types in all scope to hidden
 808      */
 809     AlignmentUtils.showOrHideSequenceAnnotations(al, types, null, true,
 810             false);
 811     assertFalse(ann1.visible);
 812     assertFalse(ann2.visible);
 813     assertTrue(ann3.visible); // not sequence-related, not affected
 814     assertFalse(ann4.visible);
 815     assertFalse(ann5.visible);
 816     assertTrue(ann6.visible); // not sequence-related, not affected
 817   }
 818
 819   /**
 820    * Tests for the method that checks if one sequence cross-references another
 821    */
 822   @Test
 823   public void testHasCrossRef()
 824   {
 825     assertFalse(AlignmentUtils.hasCrossRef(null, null));
 826     SequenceI seq1 = new Sequence("EMBL|A12345", "ABCDEF");
 827     assertFalse(AlignmentUtils.hasCrossRef(seq1, null));
 828     assertFalse(AlignmentUtils.hasCrossRef(null, seq1));
 829     SequenceI seq2 = new Sequence("UNIPROT|V20192", "ABCDEF");
 830     assertFalse(AlignmentUtils.hasCrossRef(seq1, seq2));
 831
 832     // different ref
 833     seq1.addDBRef(new DBRefEntry("UNIPROT", "1", "v20193"));
 834     assertFalse(AlignmentUtils.hasCrossRef(seq1, seq2));
 835
 836     // case-insensitive; version number is ignored
 837     seq1.addDBRef(new DBRefEntry("UNIPROT", "1", "v20192"));
 838     assertTrue(AlignmentUtils.hasCrossRef(seq1, seq2));
 839
 840     // right case!
 841     seq1.addDBRef(new DBRefEntry("UNIPROT", "1", "V20192"));
 842     assertTrue(AlignmentUtils.hasCrossRef(seq1, seq2));
 843     // test is one-way only
 844     assertFalse(AlignmentUtils.hasCrossRef(seq2, seq1));
 845   }
 846
 847   /**
 848    * Tests for the method that checks if either sequence cross-references the
 849    * other
 850    */
 851   @Test
 852   public void testHaveCrossRef()
 853   {
 854     assertFalse(AlignmentUtils.hasCrossRef(null, null));
 855     SequenceI seq1 = new Sequence("EMBL|A12345", "ABCDEF");
 856     assertFalse(AlignmentUtils.haveCrossRef(seq1, null));
 857     assertFalse(AlignmentUtils.haveCrossRef(null, seq1));
 858     SequenceI seq2 = new Sequence("UNIPROT|V20192", "ABCDEF");
 859     assertFalse(AlignmentUtils.haveCrossRef(seq1, seq2));
 860
 861     seq1.addDBRef(new DBRefEntry("UNIPROT", "1", "V20192"));
 862     assertTrue(AlignmentUtils.haveCrossRef(seq1, seq2));
 863     // next is true for haveCrossRef, false for hasCrossRef
 864     assertTrue(AlignmentUtils.haveCrossRef(seq2, seq1));
 865
 866     // now the other way round
 867     seq1.setDBRef(null);
 868     seq2.addDBRef(new DBRefEntry("EMBL", "1", "A12345"));
 869     assertTrue(AlignmentUtils.haveCrossRef(seq1, seq2));
 870     assertTrue(AlignmentUtils.haveCrossRef(seq2, seq1));
 871
 872     // now both ways
 873     seq1.addDBRef(new DBRefEntry("UNIPROT", "1", "V20192"));
 874     assertTrue(AlignmentUtils.haveCrossRef(seq1, seq2));
 875     assertTrue(AlignmentUtils.haveCrossRef(seq2, seq1));
 876   }
 877
 878   /**
 879    * Test the method that extracts the exon-only part of a dna alignment.
 880    */
 881   @Test
 882   public void testMakeExonAlignment()
 883   {
 884     SequenceI dna1 = new Sequence("dna1", "aaaGGGcccTTTaaa");
 885     SequenceI dna2 = new Sequence("dna2", "GGGcccTTTaaaCCC");
 886     SequenceI pep1 = new Sequence("pep1", "GF");
 887     SequenceI pep2 = new Sequence("pep2", "GFP");
 888     dna1.createDatasetSequence();
 889     dna2.createDatasetSequence();
 890     pep1.createDatasetSequence();
 891     pep2.createDatasetSequence();
 892
 893     Set<AlignedCodonFrame> mappings = new HashSet<AlignedCodonFrame>();
 894     MapList map = new MapList(new int[]
 895     { 4, 6, 10, 12 }, new int[]
 896     { 1, 2 }, 3, 1);
 897     AlignedCodonFrame acf = new AlignedCodonFrame();
 898     acf.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(), map);
 899     mappings.add(acf);
 900     map = new MapList(new int[]
 901     { 1, 3, 7, 9, 13, 15 }, new int[]
 902     { 1, 3 }, 3, 1);
 903     acf = new AlignedCodonFrame();
 904     acf.addMap(dna2.getDatasetSequence(), pep2.getDatasetSequence(), map);
 905     mappings.add(acf);
 906
 907     AlignmentI exons = AlignmentUtils.makeExonAlignment(new SequenceI[]
 908     { dna1, dna2 }, mappings);
 909     assertEquals(2, exons.getSequences().size());
 910     assertEquals("GGGTTT", exons.getSequenceAt(0).getSequenceAsString());
 911     assertEquals("GGGTTTCCC", exons.getSequenceAt(1).getSequenceAsString());
 912
 913     /*
 914      * Verify updated mappings
 915      */
 916     assertEquals(2, mappings.size());
 917
 918     /*
 919      * Mapping from pep1 to GGGTTT in first new exon sequence
 920      */
 921     List<AlignedCodonFrame> pep1Mapping = MappingUtils
 922             .findMappingsForSequence(pep1, mappings);
 923     assertEquals(1, pep1Mapping.size());
 924     // map G to GGG
 925     SearchResults sr = MappingUtils.buildSearchResults(pep1, 1, mappings);
 926     assertEquals(1, sr.getResults().size());
 927     Match m = sr.getResults().get(0);
 928     assertEquals(exons.getSequenceAt(0).getDatasetSequence(),
 929             m.getSequence());
 930     assertEquals(1, m.getStart());
 931     assertEquals(3, m.getEnd());
 932     // map F to TTT
 933     sr = MappingUtils.buildSearchResults(pep1, 2, mappings);
 934     m = sr.getResults().get(0);
 935     assertEquals(exons.getSequenceAt(0).getDatasetSequence(),
 936             m.getSequence());
 937     assertEquals(4, m.getStart());
 938     assertEquals(6, m.getEnd());
 939
 940     /*
 941      * Mapping from pep2 to GGGTTTCCC in second new exon sequence
 942      */
 943     List<AlignedCodonFrame> pep2Mapping = MappingUtils
 944             .findMappingsForSequence(pep2, mappings);
 945     assertEquals(1, pep2Mapping.size());
 946     // map G to GGG
 947     sr = MappingUtils.buildSearchResults(pep2, 1, mappings);
 948     assertEquals(1, sr.getResults().size());
 949     m = sr.getResults().get(0);
 950     assertEquals(exons.getSequenceAt(1).getDatasetSequence(),
 951             m.getSequence());
 952     assertEquals(1, m.getStart());
 953     assertEquals(3, m.getEnd());
 954     // map F to TTT
 955     sr = MappingUtils.buildSearchResults(pep2, 2, mappings);
 956     m = sr.getResults().get(0);
 957     assertEquals(exons.getSequenceAt(1).getDatasetSequence(),
 958             m.getSequence());
 959     assertEquals(4, m.getStart());
 960     assertEquals(6, m.getEnd());
 961     // map P to CCC
 962     sr = MappingUtils.buildSearchResults(pep2, 3, mappings);
 963     m = sr.getResults().get(0);
 964     assertEquals(exons.getSequenceAt(1).getDatasetSequence(),
 965             m.getSequence());
 966     assertEquals(7, m.getStart());
 967     assertEquals(9, m.getEnd());
 968   }
 969
 970   /**
 971    * Test the method that makes an exon-only sequence from a DNA sequence and
 972    * its product mapping. Test includes the expected case that the DNA sequence
 973    * already has a protein product (Uniprot translation) which in turn has an
 974    * x-ref to the EMBLCDS record.
 975    */
 976   @Test
 977   public void testMakeExonSequence()
 978   {
 979     SequenceI dna1 = new Sequence("dna1", "aaaGGGcccTTTaaa");
 980     SequenceI pep1 = new Sequence("pep1", "GF");
 981     dna1.createDatasetSequence();
 982     pep1.createDatasetSequence();
 983     pep1.getDatasetSequence().addDBRef(
 984             new DBRefEntry("EMBLCDS", "2", "A12345"));
 985
 986     /*
 987      * Make the mapping from dna to protein. The protein sequence has a DBRef to
 988      * EMBLCDS|A12345.
 989      */
 990     Set<AlignedCodonFrame> mappings = new HashSet<AlignedCodonFrame>();
 991     MapList map = new MapList(new int[]
 992     { 4, 6, 10, 12 }, new int[]
 993     { 1, 2 }, 3, 1);
 994     AlignedCodonFrame acf = new AlignedCodonFrame();
 995     acf.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(), map);
 996     mappings.add(acf);
 997
 998     AlignedCodonFrame newMapping = new AlignedCodonFrame();
 999     SequenceI exon = AlignmentUtils.makeExonSequence(dna1, acf, newMapping);
1000
1001     assertEquals("GGGTTT", exon.getSequenceAsString());
1002     assertEquals("dna1|A12345", exon.getName());
1003     assertEquals(1, exon.getDBRef().length);
1004     DBRefEntry cdsRef = exon.getDBRef()[0];
1005     assertEquals("EMBLCDS", cdsRef.getSource());
1006     assertEquals("2", cdsRef.getVersion());
1007     assertEquals("A12345", cdsRef.getAccessionId());
1008
1009   }
1010 }