src/jalview/datamodel/AlignedCodonFrame.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.datamodel;
  22
  23 import jalview.util.MapList;
  24 import jalview.util.MappingUtils;
  25
  26 import java.util.AbstractList;
  27 import java.util.ArrayList;
  28 import java.util.List;
  29
  30 /**
  31  * Stores mapping between the columns of a protein alignment and a DNA alignment
  32  * and a list of individual codon to amino acid mappings between sequences.
  33  */
  34 public class AlignedCodonFrame
  35 {
  36
  37   /*
  38    * Data bean to hold mappings from one sequence to another
  39    */
  40   public class SequenceToSequenceMapping
  41   {
  42     private SequenceI fromSeq;
  43
  44     private Mapping mapping;
  45
  46     SequenceToSequenceMapping(SequenceI from, Mapping map)
  47     {
  48       this.fromSeq = from;
  49       this.mapping = map;
  50     }
  51
  52     /**
  53      * Readable representation for debugging only, not guaranteed not to change
  54      */
  55     @Override
  56     public String toString()
  57     {
  58       return String.format("From %s %s", fromSeq.getName(),
  59               mapping.toString());
  60     }
  61
  62     /**
  63      * Returns a hashCode derived from the hashcodes of the mappings and fromSeq
  64      *
  65      * @see SequenceToSequenceMapping#hashCode()
  66      */
  67     @Override
  68     public int hashCode()
  69     {
  70       return (fromSeq == null ? 0 : fromSeq.hashCode() * 31)
  71               + mapping.hashCode();
  72     }
  73
  74     /**
  75      * Answers true if the objects hold the same mapping between the same two
  76      * sequences
  77      *
  78      * @see Mapping#equals
  79      */
  80     @Override
  81     public boolean equals(Object obj)
  82     {
  83       if (!(obj instanceof SequenceToSequenceMapping))
  84       {
  85         return false;
  86       }
  87       SequenceToSequenceMapping that = (SequenceToSequenceMapping) obj;
  88       if (this.mapping == null)
  89       {
  90         return that.mapping == null;
  91       }
  92       // TODO: can simplify by asserting fromSeq is a dataset sequence
  93       return (this.fromSeq == that.fromSeq || (this.fromSeq != null
  94               && that.fromSeq != null
  95               && this.fromSeq.getDatasetSequence() != null && this.fromSeq
  96               .getDatasetSequence() == that.fromSeq
  97               .getDatasetSequence())) && this.mapping.equals(that.mapping);
  98     }
  99
 100     public SequenceI getFromSeq()
 101     {
 102       return fromSeq;
 103     }
 104
 105     public Mapping getMapping()
 106     {
 107       return mapping;
 108     }
 109   }
 110
 111   private List<SequenceToSequenceMapping> mappings;
 112
 113   /**
 114    * Constructor
 115    */
 116   public AlignedCodonFrame()
 117   {
 118     mappings = new ArrayList<SequenceToSequenceMapping>();
 119   }
 120
 121   /**
 122    * Adds a mapping between the dataset sequences for the associated dna and
 123    * protein sequence objects
 124    *
 125    * @param dnaseq
 126    * @param aaseq
 127    * @param map
 128    */
 129   public void addMap(SequenceI dnaseq, SequenceI aaseq, MapList map)
 130   {
 131     // JBPNote DEBUG! THIS !
 132     // dnaseq.transferAnnotation(aaseq, mp);
 133     // aaseq.transferAnnotation(dnaseq, new Mapping(map.getInverse()));
 134
 135     SequenceI fromSeq = (dnaseq.getDatasetSequence() == null) ? dnaseq
 136             : dnaseq.getDatasetSequence();
 137     SequenceI toSeq = (aaseq.getDatasetSequence() == null) ? aaseq : aaseq
 138             .getDatasetSequence();
 139
 140     /*
 141      * if we already hold a mapping between these sequences, just add to it
 142      * note that 'adding' a duplicate map does nothing; this protects against
 143      * creating duplicate mappings in AlignedCodonFrame
 144      */
 145     for (SequenceToSequenceMapping ssm : mappings)
 146     {
 147       if (ssm.fromSeq == fromSeq && ssm.mapping.to == toSeq)
 148       {
 149         ssm.mapping.map.addMapList(map);
 150         return;
 151       }
 152     }
 153
 154     /*
 155      * otherwise, add a new sequence mapping
 156      */
 157     Mapping mp = new Mapping(toSeq, map);
 158     mappings.add(new SequenceToSequenceMapping(fromSeq, mp));
 159   }
 160
 161   public SequenceI[] getdnaSeqs()
 162   {
 163     // TODO return a list instead?
 164     // return dnaSeqs;
 165     List<SequenceI> seqs = new ArrayList<SequenceI>();
 166     for (SequenceToSequenceMapping ssm : mappings)
 167     {
 168       seqs.add(ssm.fromSeq);
 169     }
 170     return seqs.toArray(new SequenceI[seqs.size()]);
 171   }
 172
 173   public SequenceI[] getAaSeqs()
 174   {
 175     // TODO not used - remove?
 176     List<SequenceI> seqs = new ArrayList<SequenceI>();
 177     for (SequenceToSequenceMapping ssm : mappings)
 178     {
 179       seqs.add(ssm.mapping.to);
 180     }
 181     return seqs.toArray(new SequenceI[seqs.size()]);
 182   }
 183
 184   public MapList[] getdnaToProt()
 185   {
 186     List<MapList> maps = new ArrayList<MapList>();
 187     for (SequenceToSequenceMapping ssm : mappings)
 188     {
 189       maps.add(ssm.mapping.map);
 190     }
 191     return maps.toArray(new MapList[maps.size()]);
 192   }
 193
 194   public Mapping[] getProtMappings()
 195   {
 196     List<Mapping> maps = new ArrayList<Mapping>();
 197     for (SequenceToSequenceMapping ssm : mappings)
 198     {
 199       maps.add(ssm.mapping);
 200     }
 201     return maps.toArray(new Mapping[maps.size()]);
 202   }
 203
 204   /**
 205    * Returns the first mapping found which is to or from the given sequence, or
 206    * null.
 207    *
 208    * @param seq
 209    * @return
 210    */
 211   public Mapping getMappingForSequence(SequenceI seq)
 212   {
 213     SequenceI seqDs = seq.getDatasetSequence();
 214     seqDs = seqDs != null ? seqDs : seq;
 215
 216     for (SequenceToSequenceMapping ssm : mappings)
 217     {
 218       if (ssm.fromSeq == seqDs || ssm.mapping.to == seqDs)
 219       {
 220         return ssm.mapping;
 221       }
 222     }
 223     return null;
 224   }
 225
 226   /**
 227    * Return the corresponding aligned or dataset aa sequence for given dna
 228    * sequence, null if not found.
 229    *
 230    * @param sequenceRef
 231    * @return
 232    */
 233   public SequenceI getAaForDnaSeq(SequenceI dnaSeqRef)
 234   {
 235     SequenceI dnads = dnaSeqRef.getDatasetSequence();
 236     for (SequenceToSequenceMapping ssm : mappings)
 237     {
 238       if (ssm.fromSeq == dnaSeqRef || ssm.fromSeq == dnads)
 239       {
 240         return ssm.mapping.to;
 241       }
 242     }
 243     return null;
 244   }
 245
 246   /**
 247    *
 248    * @param sequenceRef
 249    * @return null or corresponding aaSeq entry for dnaSeq entry
 250    */
 251   public SequenceI getDnaForAaSeq(SequenceI aaSeqRef)
 252   {
 253     SequenceI aads = aaSeqRef.getDatasetSequence();
 254     for (SequenceToSequenceMapping ssm : mappings)
 255     {
 256       if (ssm.mapping.to == aaSeqRef || ssm.mapping.to == aads)
 257       {
 258         return ssm.fromSeq;
 259       }
 260     }
 261     return null;
 262   }
 263
 264   /**
 265    * test to see if codon frame involves seq in any way
 266    *
 267    * @param seq
 268    *          a nucleotide or protein sequence
 269    * @return true if a mapping exists to or from this sequence to any translated
 270    *         sequence
 271    */
 272   public boolean involvesSequence(SequenceI seq)
 273   {
 274     return getAaForDnaSeq(seq) != null || getDnaForAaSeq(seq) != null;
 275   }
 276
 277   /**
 278    * Add search results for regions in other sequences that translate or are
 279    * translated from a particular position in seq
 280    *
 281    * @param seq
 282    * @param index
 283    *          position in seq
 284    * @param results
 285    *          where highlighted regions go
 286    */
 287   public void markMappedRegion(SequenceI seq, int index,
 288           SearchResults results)
 289   {
 290     int[] codon;
 291     SequenceI ds = seq.getDatasetSequence();
 292     for (SequenceToSequenceMapping ssm : mappings)
 293     {
 294       if (ssm.fromSeq == seq || ssm.fromSeq == ds)
 295       {
 296         codon = ssm.mapping.map.locateInTo(index, index);
 297         if (codon != null)
 298         {
 299           for (int i = 0; i < codon.length; i += 2)
 300           {
 301             results.addResult(ssm.mapping.to, codon[i], codon[i + 1]);
 302           }
 303         }
 304       }
 305       else if (ssm.mapping.to == seq || ssm.mapping.to == ds)
 306       {
 307         {
 308           codon = ssm.mapping.map.locateInFrom(index, index);
 309           if (codon != null)
 310           {
 311             for (int i = 0; i < codon.length; i += 2)
 312             {
 313               results.addResult(ssm.fromSeq, codon[i], codon[i + 1]);
 314             }
 315           }
 316         }
 317       }
 318     }
 319   }
 320
 321   /**
 322    * Returns the DNA codon positions (base 1) for the given position (base 1) in
 323    * a mapped protein sequence, or null if no mapping is found.
 324    *
 325    * Intended for use in aligning cDNA to match aligned protein. Only the first
 326    * mapping found is returned, so not suitable for use if multiple protein
 327    * sequences are mapped to the same cDNA (but aligning cDNA as protein is
 328    * ill-defined for this case anyway).
 329    *
 330    * @param seq
 331    *          the DNA dataset sequence
 332    * @param aaPos
 333    *          residue position (base 1) in a protein sequence
 334    * @return
 335    */
 336   public int[] getDnaPosition(SequenceI seq, int aaPos)
 337   {
 338     /*
 339      * Adapted from markMappedRegion().
 340      */
 341     MapList ml = null;
 342     int i = 0;
 343     for (SequenceToSequenceMapping ssm : mappings)
 344     {
 345       if (ssm.fromSeq == seq)
 346       {
 347         ml = getdnaToProt()[i];
 348         break;
 349       }
 350       i++;
 351     }
 352     return ml == null ? null : ml.locateInFrom(aaPos, aaPos);
 353   }
 354
 355   /**
 356    * Convenience method to return the first aligned sequence in the given
 357    * alignment whose dataset has a mapping with the given (aligned or dataset)
 358    * sequence.
 359    *
 360    * @param seq
 361    *
 362    * @param al
 363    * @return
 364    */
 365   public SequenceI findAlignedSequence(SequenceI seq, AlignmentI al)
 366   {
 367     /*
 368      * Search mapped protein ('to') sequences first.
 369      */
 370     for (SequenceToSequenceMapping ssm : mappings)
 371     {
 372       if (ssm.fromSeq == seq || ssm.fromSeq == seq.getDatasetSequence())
 373       {
 374         for (SequenceI sourceAligned : al.getSequences())
 375         {
 376           if (ssm.mapping.to == sourceAligned.getDatasetSequence()
 377                   || ssm.mapping.to == sourceAligned)
 378           {
 379             return sourceAligned;
 380           }
 381         }
 382       }
 383     }
 384
 385     /*
 386      * Then try mapped dna sequences.
 387      */
 388     for (SequenceToSequenceMapping ssm : mappings)
 389     {
 390       if (ssm.mapping.to == seq
 391               || ssm.mapping.to == seq.getDatasetSequence())
 392       {
 393         for (SequenceI sourceAligned : al.getSequences())
 394         {
 395           if (ssm.fromSeq == sourceAligned.getDatasetSequence())
 396           {
 397             return sourceAligned;
 398           }
 399         }
 400       }
 401     }
 402
 403     return null;
 404   }
 405
 406   /**
 407    * Returns the region in the target sequence's dataset that is mapped to the
 408    * given position (base 1) in the query sequence's dataset. The region is a
 409    * set of start/end position pairs.
 410    *
 411    * @param target
 412    * @param query
 413    * @param queryPos
 414    * @return
 415    */
 416   public int[] getMappedRegion(SequenceI target, SequenceI query,
 417           int queryPos)
 418   {
 419     SequenceI targetDs = target.getDatasetSequence() == null ? target
 420             : target.getDatasetSequence();
 421     SequenceI queryDs = query.getDatasetSequence() == null ? query : query
 422             .getDatasetSequence();
 423     if (targetDs == null || queryDs == null /*|| dnaToProt == null*/)
 424     {
 425       return null;
 426     }
 427     for (SequenceToSequenceMapping ssm : mappings)
 428     {
 429       /*
 430        * try mapping from target to query
 431        */
 432       if (ssm.fromSeq == targetDs && ssm.mapping.to == queryDs)
 433       {
 434         int[] codon = ssm.mapping.map.locateInFrom(queryPos, queryPos);
 435         if (codon != null)
 436         {
 437           return codon;
 438         }
 439       }
 440       /*
 441        * else try mapping from query to target
 442        */
 443       else if (ssm.fromSeq == queryDs && ssm.mapping.to == targetDs)
 444       {
 445         int[] codon = ssm.mapping.map.locateInTo(queryPos, queryPos);
 446         if (codon != null)
 447         {
 448           return codon;
 449         }
 450       }
 451     }
 452     return null;
 453   }
 454
 455   /**
 456    * Returns the mapped DNA codons for the given position in a protein sequence,
 457    * or null if no mapping is found. Returns a list of (e.g.) ['g', 'c', 't']
 458    * codons. There may be more than one codon mapped to the protein if (for
 459    * example), there are mappings to cDNA variants.
 460    *
 461    * @param protein
 462    *          the peptide dataset sequence
 463    * @param aaPos
 464    *          residue position (base 1) in the peptide sequence
 465    * @return
 466    */
 467   public List<char[]> getMappedCodons(SequenceI protein, int aaPos)
 468   {
 469     MapList ml = null;
 470     SequenceI dnaSeq = null;
 471     List<char[]> result = new ArrayList<char[]>();
 472
 473     for (SequenceToSequenceMapping ssm : mappings)
 474     {
 475       if (ssm.mapping.to == protein)
 476       {
 477         ml = ssm.mapping.map;
 478         dnaSeq = ssm.fromSeq;
 479
 480         int[] codonPos = ml.locateInFrom(aaPos, aaPos);
 481         if (codonPos == null)
 482         {
 483           return null;
 484         }
 485
 486         /*
 487          * Read off the mapped nucleotides (converting to position base 0)
 488          */
 489         codonPos = MappingUtils.flattenRanges(codonPos);
 490         char[] dna = dnaSeq.getSequence();
 491         int start = dnaSeq.getStart();
 492         result.add(new char[] { dna[codonPos[0] - start],
 493             dna[codonPos[1] - start], dna[codonPos[2] - start] });
 494       }
 495     }
 496     return result.isEmpty() ? null : result;
 497   }
 498
 499   /**
 500    * Returns any mappings found which are from the given sequence, and to
 501    * distinct sequences.
 502    *
 503    * @param seq
 504    * @return
 505    */
 506   public List<Mapping> getMappingsFromSequence(SequenceI seq)
 507   {
 508     List<Mapping> result = new ArrayList<Mapping>();
 509     List<SequenceI> related = new ArrayList<SequenceI>();
 510     SequenceI seqDs = seq.getDatasetSequence();
 511     seqDs = seqDs != null ? seqDs : seq;
 512
 513     for (SequenceToSequenceMapping ssm : mappings)
 514     {
 515       final Mapping mapping = ssm.mapping;
 516       if (ssm.fromSeq == seqDs)
 517       {
 518         if (!related.contains(mapping.to))
 519         {
 520           result.add(mapping);
 521           related.add(mapping.to);
 522         }
 523       }
 524     }
 525     return result;
 526   }
 527
 528   /**
 529    * Test whether the given sequence is substitutable for one or more dummy
 530    * sequences in this mapping
 531    *
 532    * @param map
 533    * @param seq
 534    * @return
 535    */
 536   public boolean isRealisableWith(SequenceI seq)
 537   {
 538     return realiseWith(seq, false) > 0;
 539   }
 540
 541   /**
 542    * Replace any matchable mapped dummy sequences with the given real one.
 543    * Returns the count of sequence mappings instantiated.
 544    *
 545    * @param seq
 546    * @return
 547    */
 548   public int realiseWith(SequenceI seq)
 549   {
 550     return realiseWith(seq, true);
 551   }
 552
 553   /**
 554    * Returns the number of mapped dummy sequences that could be replaced with
 555    * the given real sequence.
 556    *
 557    * @param seq
 558    *          a dataset sequence
 559    * @param doUpdate
 560    *          if true, performs replacements, else only counts
 561    * @return
 562    */
 563   protected int realiseWith(SequenceI seq, boolean doUpdate)
 564   {
 565     SequenceI ds = seq.getDatasetSequence() != null ? seq
 566             .getDatasetSequence() : seq;
 567     int count = 0;
 568
 569     /*
 570      * check for replaceable DNA ('map from') sequences
 571      */
 572     for (SequenceToSequenceMapping ssm : mappings)
 573     {
 574       SequenceI dna = ssm.fromSeq;
 575       if (dna instanceof SequenceDummy
 576               && dna.getName().equals(ds.getName()))
 577       {
 578         Mapping mapping = ssm.mapping;
 579         int mapStart = mapping.getMap().getFromLowest();
 580         int mapEnd = mapping.getMap().getFromHighest();
 581         boolean mappable = couldRealiseSequence(dna, ds, mapStart, mapEnd);
 582         if (mappable)
 583         {
 584           count++;
 585           if (doUpdate)
 586           {
 587             // TODO: new method ? ds.realise(dna);
 588             // might want to copy database refs as well
 589             ds.setSequenceFeatures(dna.getSequenceFeatures());
 590             // dnaSeqs[i] = ds;
 591             ssm.fromSeq = ds;
 592             System.out.println("Realised mapped sequence " + ds.getName());
 593           }
 594         }
 595       }
 596
 597       /*
 598        * check for replaceable protein ('map to') sequences
 599        */
 600       Mapping mapping = ssm.mapping;
 601       SequenceI prot = mapping.getTo();
 602       int mapStart = mapping.getMap().getToLowest();
 603       int mapEnd = mapping.getMap().getToHighest();
 604       boolean mappable = couldRealiseSequence(prot, ds, mapStart, mapEnd);
 605       if (mappable)
 606       {
 607         count++;
 608         if (doUpdate)
 609         {
 610           // TODO: new method ? ds.realise(dna);
 611           // might want to copy database refs as well
 612           ds.setSequenceFeatures(dna.getSequenceFeatures());
 613           ssm.mapping.setTo(ds);
 614         }
 615       }
 616     }
 617     return count;
 618   }
 619
 620   /**
 621    * Helper method to test whether a 'real' sequence could replace a 'dummy'
 622    * sequence in the map. The criteria are that they have the same name, and
 623    * that the mapped region overlaps the candidate sequence.
 624    *
 625    * @param existing
 626    * @param replacement
 627    * @param mapStart
 628    * @param mapEnd
 629    * @return
 630    */
 631   protected static boolean couldRealiseSequence(SequenceI existing,
 632           SequenceI replacement, int mapStart, int mapEnd)
 633   {
 634     if (existing instanceof SequenceDummy
 635             && !(replacement instanceof SequenceDummy)
 636             && existing.getName().equals(replacement.getName()))
 637     {
 638       int start = replacement.getStart();
 639       int end = replacement.getEnd();
 640       boolean mappingOverlapsSequence = (mapStart >= start && mapStart <= end)
 641               || (mapEnd >= start && mapEnd <= end);
 642       if (mappingOverlapsSequence)
 643       {
 644         return true;
 645       }
 646     }
 647     return false;
 648   }
 649
 650   /**
 651    * Change any mapping to the given sequence to be to its dataset sequence
 652    * instead. For use when mappings are created before their referenced
 653    * sequences are instantiated, for example when parsing GFF data.
 654    *
 655    * @param seq
 656    */
 657   public void updateToDataset(SequenceI seq)
 658   {
 659     if (seq == null || seq.getDatasetSequence() == null)
 660     {
 661       return;
 662     }
 663     SequenceI ds = seq.getDatasetSequence();
 664
 665     for (SequenceToSequenceMapping ssm : mappings)
 666     /*
 667      * 'from' sequences
 668      */
 669     {
 670       if (ssm.fromSeq == seq)
 671       {
 672         ssm.fromSeq = ds;
 673       }
 674
 675       /*
 676        * 'to' sequences
 677        */
 678       if (ssm.mapping.to == seq)
 679       {
 680         ssm.mapping.to = ds;
 681       }
 682     }
 683   }
 684
 685   /**
 686    * Answers true if this object contains no mappings
 687    *
 688    * @return
 689    */
 690   public boolean isEmpty()
 691   {
 692     return mappings.isEmpty();
 693   }
 694
 695   /**
 696    * Method for debug / inspection purposes only, may change in future
 697    */
 698   @Override
 699   public String toString()
 700   {
 701     return mappings == null ? "null" : mappings.toString();
 702   }
 703
 704   /**
 705    * Returns the first mapping found that is between 'fromSeq' and 'toSeq', or null
 706    * if none found
 707    *
 708    * @param fromSeq
 709    *          aligned or dataset sequence
 710    * @param toSeq
 711    *          aligned or dataset sequence
 712    * @return
 713    */
 714   public Mapping getMappingBetween(SequenceI fromSeq, SequenceI toSeq)
 715   {
 716     SequenceI dssFrom = fromSeq.getDatasetSequence() == null ? fromSeq
 717             : fromSeq.getDatasetSequence();
 718     SequenceI dssTo = toSeq.getDatasetSequence() == null ? toSeq : toSeq
 719             .getDatasetSequence();
 720
 721     for (SequenceToSequenceMapping mapping : mappings)
 722     {
 723       SequenceI from = mapping.fromSeq;
 724       SequenceI to = mapping.mapping.to;
 725       if ((from == dssFrom && to == dssTo)
 726               || (from == dssTo && to == dssFrom))
 727       {
 728         return mapping.mapping;
 729       }
 730     }
 731     return null;
 732   }
 733
 734   /**
 735    * Returns a hashcode derived from the list of sequence mappings
 736    *
 737    * @see SequenceToSequenceMapping#hashCode()
 738    * @see AbstractList#hashCode()
 739    */
 740   @Override
 741   public int hashCode()
 742   {
 743     return this.mappings.hashCode();
 744   }
 745
 746   /**
 747    * Two AlignedCodonFrame objects are equal if they hold the same ordered list
 748    * of mappings
 749    *
 750    * @see SequenceToSequenceMapping#
 751    */
 752   @Override
 753   public boolean equals(Object obj)
 754   {
 755     if (!(obj instanceof AlignedCodonFrame))
 756     {
 757       return false;
 758     }
 759     return this.mappings.equals(((AlignedCodonFrame) obj).mappings);
 760   }
 761
 762   public List<SequenceToSequenceMapping> getMappings()
 763   {
 764     return mappings;
 765   }
 766 }