Merge branch 'bug/JAL-3806_mappingCoversSequence_for2112' into develop bug/JAL-3806_without_shared_DS_for_2_11_2
authorJim Procter <j.procter@dundee.ac.uk>
Fri, 4 Feb 2022 14:34:20 +0000 (14:34 +0000)
committerJim Procter <j.procter@dundee.ac.uk>
Fri, 4 Feb 2022 14:34:20 +0000 (14:34 +0000)
Combines patches for JAL-3806 prior to redaction of shared CDS dataset sequences with patches for JAL-3673
 Conflicts:
src/jalview/datamodel/AlignedCodonFrame.java
src/jalview/util/MappingUtils.java
test/jalview/util/MappingUtilsTest.java

1  2 
src/jalview/datamodel/AlignedCodonFrame.java
src/jalview/util/MappingUtils.java
test/jalview/util/MappingUtilsTest.java

   */
  package jalview.datamodel;
  
 -import jalview.util.MapList;
 -import jalview.util.MappingUtils;
 -
  import java.util.AbstractList;
  import java.util.ArrayList;
  import java.util.List;
  
 +import jalview.util.MapList;
 +import jalview.util.MappingUtils;
 +
  /**
   * Stores mapping between the columns of a protein alignment and a DNA alignment
   * and a list of individual codon to amino acid mappings between sequences.
@@@ -107,7 -107,7 +107,7 @@@ public class AlignedCodonFram
      {
        return mapping;
      }
 -    
 +
      /**
       * Returns true if the mapping covers the full length of the given sequence.
       * This allows us to distinguish the CDS that codes for a protein from
       */
      public boolean covers(SequenceI seq)
      {
-       List<int[]> mappedRanges = null;
+       return covers(seq,false,false);
+     }
+     /**
+      * 
+      * @param seq
+      * @param localCover - when true - compare extent of seq's dataset sequence rather than the local extent
+      * @param either - when true coverage is required for either seq or the mapped sequence 
+      * @return true if mapping covers full length of given sequence (or the other if either==true)
+      */
+     public boolean covers(SequenceI seq, boolean localCover,boolean either)
+     {
+       List<int[]> mappedRanges = null,otherRanges=null;
        MapList mapList = mapping.getMap();
+       int mstart=seq.getStart(),mend=seq.getEnd(),ostart,oend;
+               ;
        if (fromSeq == seq || fromSeq == seq.getDatasetSequence())
        {
+         if (localCover && fromSeq !=seq)
+         {
+           mstart=fromSeq.getStart();
+           mend=fromSeq.getEnd();
+         }
          mappedRanges = mapList.getFromRanges();
+         otherRanges=mapList.getToRanges();
+         ostart=mapping.to.getStart();
+         oend=mapping.to.getEnd();
        }
        else if (mapping.to == seq || mapping.to == seq.getDatasetSequence())
        {
+         if (localCover && mapping.to !=seq)
+         {
+           mstart=mapping.to.getStart();
+           mend=mapping.to.getEnd();
+         }
          mappedRanges = mapList.getToRanges();
+         otherRanges=mapList.getFromRanges();
+         ostart=fromSeq.getStart();
+         oend=fromSeq.getEnd();
        }
        else
        {
         * (necessary for circular CDS - example EMBL:J03321:AAA91567)
         * and mapped length covers (at least) sequence length
         */
-       int length = 0;
+       int length = countRange(mappedRanges,mstart,mend);
+       if (length != -1)
+       {
+         // add 3 to mapped length to allow for a mapped stop codon
+         if (length + 3 >= (mend - mstart + 1))
+         {
+           return true;
+         }
+       }
+       if (either)
+       {
+         // also check coverage of the other range
+         length = countRange(otherRanges, ostart, oend);
+         if (length != -1)
+         {
+           if (length + 1 >= (oend - ostart + 1))
+           {
+             return true;
+           }
+         }
+       }
+       return false;
+     }
+     private int countRange(List<int[]> mappedRanges,int mstart,int mend)
+     {
+       int length=0;
        for (int[] range : mappedRanges)
        {
          int from = Math.min(range[0], range[1]);
          int to = Math.max(range[0], range[1]);
-         if (from < seq.getStart() || to > seq.getEnd())
+         if (from < mstart || to > mend)
          {
-           return false;
+           return -1;
          }
          length += (to - from + 1);
        }
-       // add 1 to mapped length to allow for a mapped stop codon
-       if (length + 1 < (seq.getEnd() - seq.getStart() + 1))
-       {
-         return false;
-       }
-       return true;
+       return length;
      }
  
      /**
       * Adds any regions mapped to or from position {@code pos} in sequence
       * {@code seq} to the given search results
--     * 
++     * Note: recommend first using the .covers(,true,true) to ensure mapping covers both sequences
       * @param seq
       * @param pos
       * @param sr
    }
  
    /**
 +   * Return the corresponding aligned or dataset dna sequence for given amino
 +   * acid sequence, or null if not found. returns the sequence from the first
 +   * mapping found that involves the protein sequence.
     * 
 -   * @param sequenceRef
 -   * @return null or corresponding aaSeq entry for dnaSeq entry
 +   * @param aaSeqRef
 +   * @return
     */
    public SequenceI getDnaForAaSeq(SequenceI aaSeqRef)
    {
  
    /**
     * Add search results for regions in other sequences that translate or are
 -   * translated from a particular position in seq
 +   * translated from a particular position in seq (which may be an aligned or
 +   * dataset sequence)
     * 
     * @param seq
     * @param index
    public void markMappedRegion(SequenceI seq, int index,
            SearchResultsI results)
    {
 -    int[] codon;
      SequenceI ds = seq.getDatasetSequence();
 -    for (SequenceToSequenceMapping ssm : mappings)
 +    if (ds == null)
      {
 -      if (ssm.covers(seq,true,true))
 -      {
 -      if ((ssm.fromSeq == seq || ssm.fromSeq == ds))
 -      {
 -        codon = ssm.mapping.map.locateInTo(index, index);
 -        if (codon != null)
 -        {
 -          for (int i = 0; i < codon.length; i += 2)
 -          {
 -            results.addResult(ssm.mapping.to, codon[i], codon[i + 1]);
 -          }
 -        }
 -      }
 -      else if ((ssm.mapping.to == seq || ssm.mapping.to == ds))
 -      {
 -        {
 -          codon = ssm.mapping.map.locateInFrom(index, index);
 -          if (codon != null)
 -          {
 -            for (int i = 0; i < codon.length; i += 2)
 -            {
 -              results.addResult(ssm.fromSeq, codon[i], codon[i + 1]);
 -            }
 -          }
 -        }
 -      }}
 +      ds = seq;
      }
 -  }
 -
 -  /**
 -   * Returns the DNA codon positions (base 1) for the given position (base 1) in
 -   * a mapped protein sequence, or null if no mapping is found.
 -   * 
 -   * Intended for use in aligning cDNA to match aligned protein. Only the first
 -   * mapping found is returned, so not suitable for use if multiple protein
 -   * sequences are mapped to the same cDNA (but aligning cDNA as protein is
 -   * ill-defined for this case anyway).
 -   * 
 -   * @param seq
 -   *          the DNA dataset sequence
 -   * @param aaPos
 -   *          residue position (base 1) in a protein sequence
 -   * @return
 -   */
 -  public int[] getDnaPosition(SequenceI seq, int aaPos)
 -  {
 -    /*
 -     * Adapted from markMappedRegion().
 -     */
 -    MapList ml = null;
 -    int i = 0;
      for (SequenceToSequenceMapping ssm : mappings)
      {
-       ssm.markMappedRegion(ds, index, results);
 -      if (ssm.fromSeq == seq)
 -      {
 -        ml = getdnaToProt()[i];
 -        break;
++      if (ssm.covers(seq,true,true)) {
++        ssm.markMappedRegion(ds, index, results);
+       }
 -      i++;
      }
 -    return ml == null ? null : ml.locateInFrom(aaPos, aaPos);
    }
  
    /**
     * Two AlignedCodonFrame objects are equal if they hold the same ordered list
     * of mappings
     * 
 -   * @see SequenceToSequenceMapping#
 +   * @see SequenceToSequenceMapping#equals
     */
    @Override
    public boolean equals(Object obj)
    {
      return mappings;
    }
 +
 +  /**
 +   * Returns the first mapping found which is between the two given sequences,
 +   * and covers the full extent of both.
 +   * 
 +   * @param seq1
 +   * @param seq2
 +   * @return
 +   */
 +  public SequenceToSequenceMapping getCoveringMapping(SequenceI seq1,
 +          SequenceI seq2)
 +  {
 +    for (SequenceToSequenceMapping mapping : mappings)
 +    {
 +      if (mapping.covers(seq2) && mapping.covers(seq1))
 +      {
 +        return mapping;
 +      }
 +    }
 +    return null;
 +  }
 +
 +  /**
 +   * Returns the first mapping found which is between the given dataset sequence
 +   * and another, is a triplet mapping (3:1 or 1:3), and covers the full extent
 +   * of both sequences involved
 +   * 
 +   * @param seq
 +   * @return
 +   */
 +  public SequenceToSequenceMapping getCoveringCodonMapping(SequenceI seq)
 +  {
 +    for (SequenceToSequenceMapping mapping : mappings)
 +    {
 +      if (mapping.getMapping().getMap().isTripletMap()
 +              && mapping.covers(seq))
 +      {
 +        if (mapping.fromSeq == seq
 +                && mapping.covers(mapping.getMapping().getTo()))
 +        {
 +          return mapping;
 +        }
 +        else if (mapping.getMapping().getTo() == seq
 +                && mapping.covers(mapping.fromSeq))
 +        {
 +          return mapping;
 +        }
 +      }
 +    }
 +    return null;
 +  }
  }
@@@ -41,6 -41,7 +41,7 @@@ import jalview.datamodel.AlignmentI
  import jalview.datamodel.AlignmentOrder;
  import jalview.datamodel.ColumnSelection;
  import jalview.datamodel.HiddenColumns;
+ import jalview.datamodel.Mapping;
  import jalview.datamodel.SearchResultMatchI;
  import jalview.datamodel.SearchResults;
  import jalview.datamodel.SearchResultsI;
@@@ -363,49 -364,53 +364,55 @@@ public final class MappingUtil
         */
        int startResiduePos = selected.findPosition(firstUngappedPos);
        int endResiduePos = selected.findPosition(lastUngappedPos);
-       for (AlignedCodonFrame acf : codonFrames)
+       for (SequenceI seq : mapTo.getAlignment().getSequences())
        {
-         for (SequenceI seq : mapTo.getAlignment().getSequences())
+         int mappedStartResidue = 0;
+         int mappedEndResidue = 0;
+         for (AlignedCodonFrame acf : codonFrames)
          {
-           SequenceI peptide = targetIsNucleotide ? selected : seq;
-           SequenceI cds = targetIsNucleotide ? seq : selected;
-           SequenceToSequenceMapping s2s = acf.getCoveringMapping(cds,
-                   peptide);
-           if (s2s == null)
-           {
-             continue;
-           }
-           int mappedStartResidue = 0;
-           int mappedEndResidue = 0;
-           List<AlignedCodonFrame> mapping = Arrays.asList(acf);
-           SearchResultsI sr = buildSearchResults(selected, startResiduePos,
-                   mapping);
-           for (SearchResultMatchI m : sr.getResults())
++          // rather than use acf.getCoveringMapping() we iterate through all
++          // mappings to make sure all CDS are selected for a protein
+           for (SequenceToSequenceMapping map: acf.getMappings())
            {
-             mappedStartResidue = m.getStart();
-             mappedEndResidue = m.getEnd();
-           }
-           sr = buildSearchResults(selected, endResiduePos, mapping);
-           for (SearchResultMatchI m : sr.getResults())
+           if (map.covers(selected) && map.covers(seq))
            {
-             mappedStartResidue = Math.min(mappedStartResidue, m.getStart());
-             mappedEndResidue = Math.max(mappedEndResidue, m.getEnd());
-           }
+             /*
+              * Found a sequence mapping. Locate the start/end mapped residues.
+              */
+             List<AlignedCodonFrame> mapping = Arrays
+                     .asList(new AlignedCodonFrame[]
+                     { acf });
+             // locate start 
+             SearchResultsI sr = buildSearchResults(selected,
+                     startResiduePos, mapping);
+             for (SearchResultMatchI m : sr.getResults())
+             {
+               mappedStartResidue = m.getStart();
+               mappedEndResidue = m.getEnd();
+             }
+             // locate end - allowing for adjustment of start range
+             sr = buildSearchResults(selected, endResiduePos, mapping);
+             for (SearchResultMatchI m : sr.getResults())
+             {
+               mappedStartResidue = Math.min(mappedStartResidue,
+                       m.getStart());
+               mappedEndResidue = Math.max(mappedEndResidue, m.getEnd());
+             }
  
-           /*
-            * Find the mapped aligned columns, save the range. Note findIndex
-            * returns a base 1 position, SequenceGroup uses base 0
-            */
-           int mappedStartCol = seq.findIndex(mappedStartResidue) - 1;
-           minStartCol = minStartCol == -1 ? mappedStartCol
-                   : Math.min(minStartCol, mappedStartCol);
-           int mappedEndCol = seq.findIndex(mappedEndResidue) - 1;
-           maxEndCol = maxEndCol == -1 ? mappedEndCol
-                   : Math.max(maxEndCol, mappedEndCol);
-           mappedGroup.addSequence(seq, false);
-           break;
-         }
+             /*
+              * Find the mapped aligned columns, save the range. Note findIndex
+              * returns a base 1 position, SequenceGroup uses base 0
+              */
+             int mappedStartCol = seq.findIndex(mappedStartResidue) - 1;
+             minStartCol = minStartCol == -1 ? mappedStartCol
+                     : Math.min(minStartCol, mappedStartCol);
+             int mappedEndCol = seq.findIndex(mappedEndResidue) - 1;
+             maxEndCol = maxEndCol == -1 ? mappedEndCol
+                     : Math.max(maxEndCol, mappedEndCol);
+             mappedGroup.addSequence(seq, false);
+             break;
+           }
+         }}
        }
      }
      mappedGroup.setStartRes(minStartCol < 0 ? 0 : minStartCol);
      {
        for (AlignedCodonFrame acf : mappings)
        {
 -        SequenceI mappedSeq = mappingToNucleotide ? acf.getDnaForAaSeq(seq)
 -                : acf.getAaForDnaSeq(seq);
 -        if (mappedSeq != null)
 -        {
            for (SequenceI seq2 : mapTo.getSequences())
            {
 -            if (seq2.getDatasetSequence() == mappedSeq)
 +            /*
 +             * the corresponding peptide / CDS is the one for which there is
 +             * a complete ('covering') mapping to 'seq'
 +             */
 +            SequenceI peptide = mappingToNucleotide ? seq2 : seq;
 +            SequenceI cds = mappingToNucleotide ? seq : seq2;
 +            SequenceToSequenceMapping s2s = acf.getCoveringMapping(cds,
 +                    peptide);
 +            if (s2s != null)
              {
                mappedOrder.add(seq2);
                j++;
                break;
              }
            }
 -        }
        }
      }
  
  
      if (colsel == null)
      {
 -      return; // mappedColumns;
 +      return; 
      }
  
      char fromGapChar = mapFrom.getAlignment().getGapCharacter();
        mapHiddenColumns(regions.next(), codonFrames, newHidden,
                fromSequences, toSequences, fromGapChar);
      }
 -    return; // mappedColumns;
 +    return; 
    }
  
    /**
           */
          for (SequenceI toSeq : toSequences)
          {
 -          if (toSeq.getDatasetSequence() == mappedSeq)
 +          if (toSeq.getDatasetSequence() == mappedSeq
 +                  && mappedStartResidue >= toSeq.getStart()
 +                  && mappedEndResidue <= toSeq.getEnd())
            {
              int mappedStartCol = toSeq.findIndex(mappedStartResidue);
              int mappedEndCol = toSeq.findIndex(mappedEndResidue);
@@@ -37,16 -37,6 +37,6 @@@ import java.util.List
  import org.testng.annotations.BeforeClass;
  import org.testng.annotations.Test;
  
- import java.awt.Color;
- import java.io.IOException;
- import java.util.ArrayList;
- import java.util.Arrays;
- import java.util.Iterator;
- import java.util.List;
- import org.testng.annotations.BeforeClass;
- import org.testng.annotations.Test;
  import jalview.api.AlignViewportI;
  import jalview.bin.Cache;
  import jalview.commands.EditCommand;
@@@ -254,7 -244,7 +244,7 @@@ public class MappingUtilsTes
      protein.setCodonFrames(acfList);
  
      /*
 -     * Select Seq1 and Seq3 in the protein (startRes=endRes=0)
 +     * Select Seq1 and Seq3 in the protein
       */
      SequenceGroup sg = new SequenceGroup();
      sg.setColourText(true);
      sg.setOutlineColour(Color.LIGHT_GRAY);
      sg.addSequence(protein.getSequenceAt(0), false);
      sg.addSequence(protein.getSequenceAt(2), false);
 +    sg.setEndRes(protein.getWidth() - 1);
  
      /*
       * Verify the mapped sequence group in dna
      assertSame(cdna.getSequenceAt(0), mappedGroup.getSequences().get(0));
      assertSame(cdna.getSequenceAt(2), mappedGroup.getSequences().get(1));
      assertEquals(0, mappedGroup.getStartRes());
 -    assertEquals(2, mappedGroup.getEndRes());
 +    assertEquals(2, mappedGroup.getEndRes()); // 3 columns (1 codon)
  
      /*
       * Verify mapping sequence group from dna to protein
      overlap = MappingUtils.findOverlap(ranges, 13, 15);
      assertNull(overlap);
    }
 +  
 +  /**
 +   * Test mapping a sequence group where sequences in and outside the group
 +   * share a dataset sequence (e.g. alternative CDS for the same gene)
 +   * <p>
 +   * This scenario doesn't arise after JAL-3763 changes, but test left as still valid
 +   * @throws IOException
 +   */
 +  @Test(groups = { "Functional" })
 +  public void testMapSequenceGroup_sharedDataset() throws IOException
 +  {
 +    /*
 +     * Set up dna and protein Seq1/2/3 with mappings (held on the protein
 +     * viewport). CDS sequences share the same 'gene' dataset sequence.
 +     */
 +    SequenceI dna = new Sequence("dna", "aaatttgggcccaaatttgggccc");
 +    SequenceI cds1 = new Sequence("cds1/1-6", "aaattt");
 +    SequenceI cds2 = new Sequence("cds1/4-9", "tttggg");
 +    SequenceI cds3 = new Sequence("cds1/19-24", "gggccc");
 +
 +    cds1.setDatasetSequence(dna);
 +    cds2.setDatasetSequence(dna);
 +    cds3.setDatasetSequence(dna);
 +
 +    SequenceI pep1 = new Sequence("pep1", "KF");
 +    SequenceI pep2 = new Sequence("pep2", "FG");
 +    SequenceI pep3 = new Sequence("pep3", "GP");
 +    pep1.createDatasetSequence();
 +    pep2.createDatasetSequence();
 +    pep3.createDatasetSequence();
 +
 +    /*
 +     * add mappings from coding positions of dna to respective peptides
 +     */
 +    AlignedCodonFrame acf = new AlignedCodonFrame();
 +    acf.addMap(dna, pep1,
 +            new MapList(new int[]
 +            { 1, 6 }, new int[] { 1, 2 }, 3, 1));
 +    acf.addMap(dna, pep2,
 +            new MapList(new int[]
 +            { 4, 9 }, new int[] { 1, 2 }, 3, 1));
 +    acf.addMap(dna, pep3,
 +            new MapList(new int[]
 +            { 19, 24 }, new int[] { 1, 2 }, 3, 1));
 +
 +    List<AlignedCodonFrame> acfList = Arrays
 +            .asList(new AlignedCodonFrame[]
 +            { acf });
 +
 +    AlignmentI cdna = new Alignment(new SequenceI[] { cds1, cds2, cds3 });
 +    AlignmentI protein = new Alignment(
 +            new SequenceI[]
 +            { pep1, pep2, pep3 });
 +    AlignViewportI cdnaView = new AlignViewport(cdna);
 +    AlignViewportI peptideView = new AlignViewport(protein);
 +    protein.setCodonFrames(acfList);
 +
 +    /*
 +     * Select pep1 and pep3 in the protein alignment
 +     */
 +    SequenceGroup sg = new SequenceGroup();
 +    sg.setColourText(true);
 +    sg.setIdColour(Color.GREEN);
 +    sg.setOutlineColour(Color.LIGHT_GRAY);
 +    sg.addSequence(pep1, false);
 +    sg.addSequence(pep3, false);
 +    sg.setEndRes(protein.getWidth() - 1);
 +
 +    /*
 +     * Verify the mapped sequence group in dna is cds1 and cds3
 +     */
 +    SequenceGroup mappedGroup = MappingUtils.mapSequenceGroup(sg,
 +            peptideView, cdnaView);
 +    assertTrue(mappedGroup.getColourText());
 +    assertSame(sg.getIdColour(), mappedGroup.getIdColour());
 +    assertSame(sg.getOutlineColour(), mappedGroup.getOutlineColour());
 +    assertEquals(2, mappedGroup.getSequences().size());
 +    assertSame(cds1, mappedGroup.getSequences().get(0));
 +    assertSame(cds3, mappedGroup.getSequences().get(1));
 +    // columns 1-6 selected (0-5 base zero)
 +    assertEquals(0, mappedGroup.getStartRes());
 +    assertEquals(5, mappedGroup.getEndRes());
 +
 +    /*
 +     * Select mapping sequence group from dna to protein
 +     */
 +    sg.clear();
 +    sg.addSequence(cds2, false);
 +    sg.addSequence(cds1, false);
 +    sg.setStartRes(0);
 +    sg.setEndRes(cdna.getWidth() - 1);
 +    mappedGroup = MappingUtils.mapSequenceGroup(sg, cdnaView, peptideView);
 +    assertTrue(mappedGroup.getColourText());
 +    assertSame(sg.getIdColour(), mappedGroup.getIdColour());
 +    assertSame(sg.getOutlineColour(), mappedGroup.getOutlineColour());
 +    assertEquals(2, mappedGroup.getSequences().size());
 +    assertSame(protein.getSequenceAt(1), mappedGroup.getSequences().get(0));
 +    assertSame(protein.getSequenceAt(0), mappedGroup.getSequences().get(1));
 +    assertEquals(0, mappedGroup.getStartRes());
 +    assertEquals(1, mappedGroup.getEndRes()); // two columns
 +  }
  }