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