From 6395e8fc1fc6b33c3b305d74d9be3642dfca8417 Mon Sep 17 00:00:00 2001 From: gmungoc Date: Wed, 6 Jul 2016 13:08:28 +0100 Subject: [PATCH] JAL-2110 ensuring mapped dbrefs between protein and cds --- src/jalview/analysis/AlignmentUtils.java | 308 +++++++++++++++++++----- src/jalview/analysis/CrossRef.java | 24 +- src/jalview/datamodel/AlignedCodonFrame.java | 16 ++ src/jalview/datamodel/Alignment.java | 5 +- src/jalview/ext/ensembl/EnsemblSeqProxy.java | 12 +- test/jalview/analysis/AlignmentUtilsTests.java | 82 +++++-- 6 files changed, 351 insertions(+), 96 deletions(-) diff --git a/src/jalview/analysis/AlignmentUtils.java b/src/jalview/analysis/AlignmentUtils.java index aa7cb18..35aa91d 100644 --- a/src/jalview/analysis/AlignmentUtils.java +++ b/src/jalview/analysis/AlignmentUtils.java @@ -22,6 +22,7 @@ package jalview.analysis; import static jalview.io.gff.GffConstants.CLINICAL_SIGNIFICANCE; +import jalview.api.DBRefEntryI; import jalview.datamodel.AlignedCodon; import jalview.datamodel.AlignedCodonFrame; import jalview.datamodel.AlignedCodonFrame.SequenceToSequenceMapping; @@ -39,6 +40,7 @@ import jalview.io.gff.SequenceOntologyFactory; import jalview.io.gff.SequenceOntologyI; import jalview.schemes.ResidueProperties; import jalview.util.Comparison; +import jalview.util.DBRefUtils; import jalview.util.MapList; import jalview.util.MappingUtils; import jalview.util.StringUtils; @@ -850,6 +852,11 @@ public class AlignmentUtils */ public static int alignProteinAsDna(AlignmentI protein, AlignmentI dna) { + if (protein.isNucleotide() || !dna.isNucleotide()) + { + System.err.println("Wrong alignment type in alignProteinAsDna"); + return 0; + } List unmappedProtein = new ArrayList(); Map> alignedCodons = buildCodonColumnsMap( protein, dna, unmappedProtein); @@ -857,6 +864,68 @@ public class AlignmentUtils } /** + * Realigns the given dna to match the alignment of the protein, using codon + * mappings to translate aligned peptide positions to codons. + * + * @param dna + * the alignment whose sequences are realigned by this method + * @param protein + * the protein alignment whose alignment we are 'copying' + * @return the number of sequences that were realigned + */ + public static int alignCdsAsProtein(AlignmentI dna, AlignmentI protein) + { + if (protein.isNucleotide() || !dna.isNucleotide()) + { + System.err.println("Wrong alignment type in alignProteinAsDna"); + return 0; + } + // todo: implement this + List mappings = protein.getCodonFrames(); + int alignedCount = 0; + for (SequenceI dnaSeq : dna.getSequences()) + { + if (alignCdsSequenceAsProtein(dnaSeq, protein, mappings)) + { + alignedCount++; + } + } + return alignedCount; + } + + /** + * Helper method to align (if possible) the dna sequence (cds only) to match + * the alignment of a mapped protein sequence + * + * @param dnaSeq + * @param protein + * @param mappings + * @return + */ + static boolean alignCdsSequenceAsProtein(SequenceI dnaSeq, + AlignmentI protein, List mappings) + { + List dnaMappings = MappingUtils + .findMappingsForSequence(dnaSeq, mappings); + for (AlignedCodonFrame mapping : dnaMappings) + { + SequenceI peptide = mapping.findAlignedSequence(dnaSeq, protein); + if (peptide != null) + { + Mapping map = mapping.getMappingBetween(dnaSeq, peptide); + if (map != null) + { + /* + * traverse peptide adding gaps or codons to new cds sequence + */ + MapList mapList = map.getMap(); + } + } + } + return false; + } + + /** * Builds a map whose key is an aligned codon position (3 alignment column * numbers base 0), and whose value is a map from protein sequence to each * protein's peptide residue for that codon. The map generates an ordering of @@ -1404,9 +1473,9 @@ public class AlignmentUtils * added to the alignment dataset. * * @param dna - * aligned dna sequences + * aligned nucleotide (dna or cds) sequences * @param dataset - * - throws error if not given a dataset + * the alignment dataset the sequences belong to * @param products * (optional) to restrict results to CDS that map to specified * protein products @@ -1414,20 +1483,21 @@ public class AlignmentUtils * sequences (or null if no mappings are found) */ public static AlignmentI makeCdsAlignment(SequenceI[] dna, - AlignmentI dataset, AlignmentI products) + AlignmentI dataset, SequenceI[] products) { - if (dataset.getDataset() != null) + if (dataset == null || dataset.getDataset() != null) { - throw new Error( + throw new IllegalArgumentException( "IMPLEMENTATION ERROR: dataset.getDataset() must be null!"); } + List foundSeqs = new ArrayList(); List cdsSeqs = new ArrayList(); List mappings = dataset.getCodonFrames(); HashSet productSeqs = null; if (products != null) { productSeqs = new HashSet(); - for (SequenceI seq : products.getSequences()) + for (SequenceI seq : products) { productSeqs.add(seq.getDatasetSequence() == null ? seq : seq .getDatasetSequence()); @@ -1435,24 +1505,32 @@ public class AlignmentUtils } /* - * construct CDS sequences from the (cds-to-protein) mappings made earlier; - * this makes it possible to model multiple products from dna (e.g. EMBL); - * however it does mean we don't have the EMBL protein_id (a property on - * the CDS features) in order to make the CDS sequence name :-( + * Construct CDS sequences from mappings on the alignment dataset. + * The logic is: + * - find the protein product(s) mapped to from each dna sequence + * - if the mapping covers the whole dna sequence (give or take start/stop + * codon), take the dna as the CDS sequence + * - else search dataset mappings for a suitable dna sequence, i.e. one + * whose whole sequence is mapped to the protein + * - if no sequence found, construct one from the dna sequence and mapping + * (and add it to dataset so it is found if this is repeated) */ - for (SequenceI seq : dna) + for (SequenceI dnaSeq : dna) { - SequenceI seqDss = seq.getDatasetSequence() == null ? seq : seq - .getDatasetSequence(); + SequenceI dnaDss = dnaSeq.getDatasetSequence() == null ? dnaSeq + : dnaSeq.getDatasetSequence(); + List seqMappings = MappingUtils - .findMappingsForSequence(seq, mappings); + .findMappingsForSequence(dnaSeq, mappings); for (AlignedCodonFrame mapping : seqMappings) { - List mappingsFromSequence = mapping.getMappingsFromSequence(seq); + List mappingsFromSequence = mapping + .getMappingsFromSequence(dnaSeq); for (Mapping aMapping : mappingsFromSequence) { - if (aMapping.getMap().getFromRatio() == 1) + MapList mapList = aMapping.getMap(); + if (mapList.getFromRatio() == 1) { /* * not a dna-to-protein mapping (likely dna-to-cds) @@ -1461,55 +1539,33 @@ public class AlignmentUtils } /* - * check for an existing CDS sequence i.e. a 3:1 mapping to - * the dna mapping's product + * skip if mapping is not to one of the target set of proteins */ - SequenceI cdsSeq = null; - - // TODO better mappings collection data model so we can do - // a direct lookup instead of double loops to find mappings - SequenceI proteinProduct = aMapping.getTo(); - - /* - * skip if not mapped to one of a specified set of proteins - */ if (productSeqs != null && !productSeqs.contains(proteinProduct)) { continue; } - for (AlignedCodonFrame acf : MappingUtils - .findMappingsForSequence(proteinProduct, mappings)) + /* + * try to locate the CDS from the dataset mappings; + * guard against duplicate results (for the case that protein has + * dbrefs to both dna and cds sequences) + */ + SequenceI cdsSeq = findCdsForProtein(mappings, dnaSeq, + seqMappings, aMapping); + if (cdsSeq != null) { - for (SequenceToSequenceMapping map : acf.getMappings()) + if (!foundSeqs.contains(cdsSeq)) { - if (map.getMapping().getMap().getFromRatio() == 3 - && proteinProduct == map.getMapping().getTo() - && seqDss != map.getFromSeq()) + foundSeqs.add(cdsSeq); + SequenceI derivedSequence = cdsSeq.deriveSequence(); + cdsSeqs.add(derivedSequence); + if (!dataset.getSequences().contains(cdsSeq)) { - /* - * found a 3:1 mapping to the protein product which is not - * from the dna sequence...assume it is from the CDS sequence - * TODO mappings data model that brings together related - * dna-cds-protein mappings in one object - */ - cdsSeq = map.getFromSeq(); + dataset.addSequence(cdsSeq); } } - } - if (cdsSeq != null) - { - /* - * mappings are always to dataset sequences so create an aligned - * sequence to own it; add the dataset sequence to the dataset - */ - SequenceI derivedSequence = cdsSeq.deriveSequence(); - cdsSeqs.add(derivedSequence); - if (!dataset.getSequences().contains(cdsSeq)) - { - dataset.addSequence(cdsSeq); - } continue; } @@ -1517,7 +1573,7 @@ public class AlignmentUtils * didn't find mapped CDS sequence - construct it and add * its dataset sequence to the dataset */ - cdsSeq = makeCdsSequence(seq.getDatasetSequence(), aMapping); + cdsSeq = makeCdsSequence(dnaSeq.getDatasetSequence(), aMapping); SequenceI cdsSeqDss = cdsSeq.createDatasetSequence(); cdsSeqs.add(cdsSeq); if (!dataset.getSequences().contains(cdsSeqDss)) @@ -1530,11 +1586,10 @@ public class AlignmentUtils */ List cdsRange = Collections.singletonList(new int[] { 1, cdsSeq.getLength() }); - MapList map = new MapList(cdsRange, aMapping.getMap() - .getToRanges(), aMapping.getMap().getFromRatio(), - aMapping.getMap().getToRatio()); + MapList cdsToProteinMap = new MapList(cdsRange, mapList.getToRanges(), + mapList.getFromRatio(), mapList.getToRatio()); AlignedCodonFrame cdsToProteinMapping = new AlignedCodonFrame(); - cdsToProteinMapping.addMap(cdsSeq, proteinProduct, map); + cdsToProteinMapping.addMap(cdsSeq, proteinProduct, cdsToProteinMap); /* * guard against duplicating the mapping if repeating this action @@ -1545,22 +1600,59 @@ public class AlignmentUtils } /* + * copy protein's dbrefs to CDS sequence + * this enables Get Cross-References from CDS alignment + */ + DBRefEntry[] proteinRefs = DBRefUtils.selectDbRefs(false, + proteinProduct.getDBRefs()); + if (proteinRefs != null) + { + for (DBRefEntry ref : proteinRefs) + { + DBRefEntry cdsToProteinRef = new DBRefEntry(ref); + cdsToProteinRef.setMap(new Mapping(proteinProduct, + cdsToProteinMap)); + cdsSeqDss.addDBRef(cdsToProteinRef); + } + } + + /* * add another mapping from original 'from' range to CDS */ - AlignedCodonFrame dnaToProteinMapping = new AlignedCodonFrame(); - map = new MapList(aMapping.getMap().getFromRanges(), cdsRange, 1, + AlignedCodonFrame dnaToCdsMapping = new AlignedCodonFrame(); + MapList dnaToCdsMap = new MapList(mapList.getFromRanges(), + cdsRange, 1, 1); - dnaToProteinMapping.addMap(seq.getDatasetSequence(), cdsSeq, map); - if (!mappings.contains(dnaToProteinMapping)) + dnaToCdsMapping.addMap(dnaSeq.getDatasetSequence(), cdsSeq, + dnaToCdsMap); + if (!mappings.contains(dnaToCdsMapping)) { - mappings.add(dnaToProteinMapping); + mappings.add(dnaToCdsMapping); } + /* + * add DBRef with mapping from protein to CDS + * (this enables Get Cross-References from protein alignment) + * This is tricky because we can't have two DBRefs with the + * same source and accession, so need a different accession for + * the CDS from the dna sequence + */ + DBRefEntryI dnaRef = dnaDss.getSourceDBRef(); + if (dnaRef != null) + { + // assuming cds version same as dna ?!? + DBRefEntry proteinToCdsRef = new DBRefEntry(dnaRef.getSource(), + dnaRef.getVersion(), cdsSeq.getName()); + proteinToCdsRef.setMap(new Mapping(cdsSeqDss, cdsToProteinMap + .getInverse())); + proteinProduct.addDBRef(proteinToCdsRef); + } /* * transfer any features on dna that overlap the CDS */ - transferFeatures(seq, cdsSeq, map, null, SequenceOntologyI.CDS); + transferFeatures(dnaSeq, cdsSeq, cdsToProteinMap, null, + SequenceOntologyI.CDS); } } } @@ -1573,6 +1665,83 @@ public class AlignmentUtils } /** + * A helper method that finds a CDS sequence in the alignment dataset that is + * mapped to the given protein sequence, and either is, or has a mapping from, + * the given dna sequence. + * + * @param mappings + * set of all mappings on the dataset + * @param dnaSeq + * a dna (or cds) sequence we are searching from + * @param seqMappings + * the set of mappings involving dnaSeq + * @param aMapping + * an initial candidate from seqMappings + * @return + */ + static SequenceI findCdsForProtein(List mappings, + SequenceI dnaSeq, List seqMappings, + Mapping aMapping) + { + /* + * TODO a better dna-cds-protein mapping data representation to allow easy + * navigation; until then this clunky looping around lists of mappings + */ + SequenceI seqDss = dnaSeq.getDatasetSequence() == null ? dnaSeq + : dnaSeq.getDatasetSequence(); + SequenceI proteinProduct = aMapping.getTo(); + + /* + * is this mapping from the whole dna sequence (i.e. CDS)? + * allowing for possible stop codon on dna but not peptide + */ + int mappedFromLength = MappingUtils.getLength(aMapping.getMap() + .getFromRanges()); + int dnaLength = seqDss.getLength(); + if (mappedFromLength == dnaLength || mappedFromLength == dnaLength - 3) + { + return seqDss; + } + + /* + * looks like we found the dna-to-protein mapping; search for the + * corresponding cds-to-protein mapping + */ + List mappingsToPeptide = MappingUtils + .findMappingsForSequence(proteinProduct, mappings); + for (AlignedCodonFrame acf : mappingsToPeptide) + { + for (SequenceToSequenceMapping map : acf.getMappings()) + { + Mapping mapping = map.getMapping(); + if (mapping != aMapping && mapping.getMap().getFromRatio() == 3 + && proteinProduct == mapping.getTo() + && seqDss != map.getFromSeq()) + { + mappedFromLength = MappingUtils.getLength(mapping.getMap() + .getFromRanges()); + if (mappedFromLength == map.getFromSeq().getLength()) + { + /* + * found a 3:1 mapping to the protein product which covers + * the whole dna sequence i.e. is from CDS; finally check it + * is from the dna start sequence + */ + SequenceI cdsSeq = map.getFromSeq(); + List dnaToCdsMaps = MappingUtils + .findMappingsForSequence(cdsSeq, seqMappings); + if (!dnaToCdsMaps.isEmpty()) + { + return cdsSeq; + } + } + } + } + } + return null; + } + + /** * Helper method that makes a CDS sequence as defined by the mappings from the * given sequence i.e. extracts the 'mapped from' ranges (which may be on * forward or reverse strand). @@ -1609,8 +1778,15 @@ public class AlignmentUtils } } - SequenceI newSeq = new Sequence(seq.getName() + "|" - + mapping.getTo().getName(), newSeqChars, 1, newPos); + /* + * assign 'from id' held in the mapping if set (e.g. EMBL protein_id), + * else generate a sequence name + */ + String mapFromId = mapping.getMappedFromId(); + String seqId = "CDS|" + (mapFromId != null ? mapFromId : seq.getName()); + SequenceI newSeq = new Sequence(seqId, newSeqChars, 1, newPos); + // newSeq.setDescription(mapFromId); + return newSeq; } diff --git a/src/jalview/analysis/CrossRef.java b/src/jalview/analysis/CrossRef.java index 2b5a0e2..df5ca10 100644 --- a/src/jalview/analysis/CrossRef.java +++ b/src/jalview/analysis/CrossRef.java @@ -295,12 +295,14 @@ public class CrossRef if (fromDna) { // map is from dna seq to a protein product - cf.addMap(dss, rsq, xref.getMap().getMap()); + cf.addMap(dss, rsq, xref.getMap().getMap(), xref.getMap() + .getMappedFromId()); } else { // map should be from protein seq to its coding dna - cf.addMap(rsq, dss, xref.getMap().getMap().getInverse()); + cf.addMap(rsq, dss, xref.getMap().getMap().getInverse(), + xref.getMap().getMappedFromId()); } } } @@ -621,14 +623,14 @@ public class CrossRef void updateDbrefMappings(SequenceI mapFrom, DBRefEntry[] xrefs, SequenceI[] retrieved, AlignedCodonFrame acf, boolean fromDna) { - SequenceIdMatcher matcher = new SequenceIdMatcher(retrieved); + SequenceIdMatcher idMatcher = new SequenceIdMatcher(retrieved); for (DBRefEntry xref : xrefs) { if (!xref.hasMap()) { String targetSeqName = xref.getSource() + "|" + xref.getAccessionId(); - SequenceI[] matches = matcher.findAllIdMatches(targetSeqName); + SequenceI[] matches = idMatcher.findAllIdMatches(targetSeqName); if (matches == null) { return; @@ -705,6 +707,20 @@ public class CrossRef return false; } xref.setMap(new Mapping(mapTo, mapping)); + + /* + * and add a reverse DbRef with the inverse mapping + */ + if (mapFrom.getDatasetSequence() != null + && mapFrom.getDatasetSequence().getSourceDBRef() != null) + { + DBRefEntry dbref = new DBRefEntry(mapFrom.getDatasetSequence() + .getSourceDBRef()); + dbref.setMap(new Mapping(mapFrom.getDatasetSequence(), mapping + .getInverse())); + mapTo.addDBRef(dbref); + } + if (fromDna) { AlignmentUtils.computeProteinFeatures(mapFrom, mapTo, mapping); diff --git a/src/jalview/datamodel/AlignedCodonFrame.java b/src/jalview/datamodel/AlignedCodonFrame.java index bb705b6..326cc4e 100644 --- a/src/jalview/datamodel/AlignedCodonFrame.java +++ b/src/jalview/datamodel/AlignedCodonFrame.java @@ -128,6 +128,21 @@ public class AlignedCodonFrame */ public void addMap(SequenceI dnaseq, SequenceI aaseq, MapList map) { + addMap(dnaseq, aaseq, map, null); + } + + /** + * Adds a mapping between the dataset sequences for the associated dna and + * protein sequence objects + * + * @param dnaseq + * @param aaseq + * @param map + * @param mapFromId + */ + public void addMap(SequenceI dnaseq, SequenceI aaseq, MapList map, + String mapFromId) + { // JBPNote DEBUG! THIS ! // dnaseq.transferAnnotation(aaseq, mp); // aaseq.transferAnnotation(dnaseq, new Mapping(map.getInverse())); @@ -155,6 +170,7 @@ public class AlignedCodonFrame * otherwise, add a new sequence mapping */ Mapping mp = new Mapping(toSeq, map); + mp.setMappedFromId(mapFromId); mappings.add(new SequenceToSequenceMapping(fromSeq, mp)); } diff --git a/src/jalview/datamodel/Alignment.java b/src/jalview/datamodel/Alignment.java index 2b194cc..f5665db 100755 --- a/src/jalview/datamodel/Alignment.java +++ b/src/jalview/datamodel/Alignment.java @@ -30,7 +30,6 @@ import java.util.Collections; import java.util.Enumeration; import java.util.HashSet; import java.util.Hashtable; -import java.util.Iterator; import java.util.List; import java.util.Map; import java.util.Set; @@ -1724,6 +1723,10 @@ public class Alignment implements AlignmentI { return AlignmentUtils.alignProteinAsDna(this, al); } + else if (thatIsProtein && thisIsNucleotide) + { + return AlignmentUtils.alignCdsAsProtein(this, al); + } return AlignmentUtils.alignAs(this, al); } diff --git a/src/jalview/ext/ensembl/EnsemblSeqProxy.java b/src/jalview/ext/ensembl/EnsemblSeqProxy.java index c86469f..31552af 100644 --- a/src/jalview/ext/ensembl/EnsemblSeqProxy.java +++ b/src/jalview/ext/ensembl/EnsemblSeqProxy.java @@ -5,7 +5,6 @@ import jalview.analysis.Dna; import jalview.datamodel.Alignment; import jalview.datamodel.AlignmentI; import jalview.datamodel.DBRefEntry; -import jalview.datamodel.DBRefSource; import jalview.datamodel.Mapping; import jalview.datamodel.SequenceFeature; import jalview.datamodel.SequenceI; @@ -315,13 +314,6 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient for (DBRefEntry xref : xrefs) { seq.addDBRef(xref); - /* - * Save any Uniprot xref to be the reference for SIFTS mapping - */ - if (DBRefSource.UNIPROT.equals(xref.getSource())) - { - seq.setSourceDBRef(xref); - } } /* @@ -330,6 +322,7 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient DBRefEntry self = new DBRefEntry(getDbSource(), getEnsemblDataVersion(), seq.getName()); seq.addDBRef(self); + seq.setSourceDBRef(self); } /** @@ -387,8 +380,9 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient if (ids.contains(name) || ids.contains(name.replace("ENSP", "ENST"))) { - DBRefUtils.parseToDbRef(sq, getDbSource(), + DBRefEntry dbref = DBRefUtils.parseToDbRef(sq, getDbSource(), getEnsemblDataVersion(), name); + sq.setSourceDBRef(dbref); } } if (alignment == null) diff --git a/test/jalview/analysis/AlignmentUtilsTests.java b/test/jalview/analysis/AlignmentUtilsTests.java index d704ec6..727cf22 100644 --- a/test/jalview/analysis/AlignmentUtilsTests.java +++ b/test/jalview/analysis/AlignmentUtilsTests.java @@ -22,6 +22,7 @@ package jalview.analysis; import static org.testng.AssertJUnit.assertEquals; import static org.testng.AssertJUnit.assertFalse; +import static org.testng.AssertJUnit.assertNotNull; import static org.testng.AssertJUnit.assertNull; import static org.testng.AssertJUnit.assertSame; import static org.testng.AssertJUnit.assertTrue; @@ -973,10 +974,17 @@ public class AlignmentUtilsTests @Test(groups = { "Functional" }) public void testMakeCdsAlignment() { + /* + * scenario: + * dna1 --> [4, 6] [10,12] --> pep1 + * dna2 --> [1, 3] [7, 9] [13,15] --> pep1 + */ SequenceI dna1 = new Sequence("dna1", "aaaGGGcccTTTaaa"); SequenceI dna2 = new Sequence("dna2", "GGGcccTTTaaaCCC"); SequenceI pep1 = new Sequence("pep1", "GF"); SequenceI pep2 = new Sequence("pep2", "GFP"); + pep1.addDBRef(new DBRefEntry("UNIPROT", "0", "pep1")); + pep2.addDBRef(new DBRefEntry("UNIPROT", "0", "pep2")); dna1.createDatasetSequence(); dna2.createDatasetSequence(); pep1.createDatasetSequence(); @@ -984,6 +992,19 @@ public class AlignmentUtilsTests AlignmentI dna = new Alignment(new SequenceI[] { dna1, dna2 }); dna.setDataset(null); + /* + * need a sourceDbRef if we are to construct dbrefs to the CDS + * sequence + */ + DBRefEntry dbref = new DBRefEntry("ENSEMBL", "0", "dna1"); + dna1.getDatasetSequence().setSourceDBRef(dbref); + dbref = new DBRefEntry("ENSEMBL", "0", "dna2"); + dna2.getDatasetSequence().setSourceDBRef(dbref); + + /* + * CDS sequences are 'discovered' from dna-to-protein mappings on the alignment + * dataset (e.g. added from dbrefs by CrossRef.findXrefSequences) + */ MapList map = new MapList(new int[] { 4, 6, 10, 12 }, new int[] { 1, 2 }, 3, 1); AlignedCodonFrame acf = new AlignedCodonFrame(); @@ -1001,6 +1022,9 @@ public class AlignmentUtilsTests AlignmentI cds = AlignmentUtils.makeCdsAlignment(new SequenceI[] { dna1, dna2 }, dna.getDataset(), null); + /* + * verify cds sequences + */ assertEquals(2, cds.getSequences().size()); assertEquals("GGGTTT", cds.getSequenceAt(0).getSequenceAsString()); assertEquals("GGGTTTCCC", cds.getSequenceAt(1).getSequenceAsString()); @@ -1009,19 +1033,45 @@ public class AlignmentUtilsTests * verify shared, extended alignment dataset */ assertSame(dna.getDataset(), cds.getDataset()); - assertTrue(dna.getDataset().getSequences() - .contains(cds.getSequenceAt(0).getDatasetSequence())); - assertTrue(dna.getDataset().getSequences() - .contains(cds.getSequenceAt(1).getDatasetSequence())); + SequenceI cds1Dss = cds.getSequenceAt(0).getDatasetSequence(); + SequenceI cds2Dss = cds.getSequenceAt(1).getDatasetSequence(); + assertTrue(dna.getDataset().getSequences().contains(cds1Dss)); + assertTrue(dna.getDataset().getSequences().contains(cds2Dss)); + + /* + * verify CDS has a dbref with mapping to peptide + */ + assertNotNull(cds1Dss.getDBRefs()); + assertEquals(1, cds1Dss.getDBRefs().length); + dbref = cds1Dss.getDBRefs()[0]; + assertEquals("UNIPROT", dbref.getSource()); + assertEquals("0", dbref.getVersion()); + assertEquals("pep1", dbref.getAccessionId()); + assertNotNull(dbref.getMap()); + assertSame(pep1.getDatasetSequence(), dbref.getMap().getTo()); + MapList cdsMapping = new MapList(new int[] { 1, 6 }, + new int[] { 1, 2 }, 3, 1); + assertEquals(cdsMapping, dbref.getMap().getMap()); /* - * Verify mappings from CDS to peptide, cDNA to CDS, and cDNA to peptide - * the mappings are on the shared alignment dataset + * verify peptide has added a dbref with reverse mapping to CDS */ - List cdsMappings = cds.getDataset().getCodonFrames(); + assertNotNull(pep1.getDBRefs()); + assertEquals(2, pep1.getDBRefs().length); + dbref = pep1.getDBRefs()[1]; + assertEquals("ENSEMBL", dbref.getSource()); + assertEquals("0", dbref.getVersion()); + assertEquals("CDS|dna1", dbref.getAccessionId()); + assertNotNull(dbref.getMap()); + assertSame(cds1Dss, dbref.getMap().getTo()); + assertEquals(cdsMapping.getInverse(), dbref.getMap().getMap()); + /* + * Verify mappings from CDS to peptide, cDNA to CDS, and cDNA to peptide + * the mappings are on the shared alignment dataset * 6 mappings, 2*(DNA->CDS), 2*(DNA->Pep), 2*(CDS->Pep) */ + List cdsMappings = cds.getDataset().getCodonFrames(); assertEquals(6, cdsMappings.size()); /* @@ -1048,13 +1098,13 @@ public class AlignmentUtilsTests SearchResults sr = MappingUtils.buildSearchResults(pep1, 1, mappings); assertEquals(1, sr.getResults().size()); Match m = sr.getResults().get(0); - assertSame(cds.getSequenceAt(0).getDatasetSequence(), m.getSequence()); + assertSame(cds1Dss, m.getSequence()); assertEquals(1, m.getStart()); assertEquals(3, m.getEnd()); // map F to TTT sr = MappingUtils.buildSearchResults(pep1, 2, mappings); m = sr.getResults().get(0); - assertSame(cds.getSequenceAt(0).getDatasetSequence(), m.getSequence()); + assertSame(cds1Dss, m.getSequence()); assertEquals(4, m.getStart()); assertEquals(6, m.getEnd()); @@ -1072,19 +1122,19 @@ public class AlignmentUtilsTests sr = MappingUtils.buildSearchResults(pep2, 1, mappings); assertEquals(1, sr.getResults().size()); m = sr.getResults().get(0); - assertSame(cds.getSequenceAt(1).getDatasetSequence(), m.getSequence()); + assertSame(cds2Dss, m.getSequence()); assertEquals(1, m.getStart()); assertEquals(3, m.getEnd()); // map F to TTT sr = MappingUtils.buildSearchResults(pep2, 2, mappings); m = sr.getResults().get(0); - assertSame(cds.getSequenceAt(1).getDatasetSequence(), m.getSequence()); + assertSame(cds2Dss, m.getSequence()); assertEquals(4, m.getStart()); assertEquals(6, m.getEnd()); // map P to CCC sr = MappingUtils.buildSearchResults(pep2, 3, mappings); m = sr.getResults().get(0); - assertSame(cds.getSequenceAt(1).getDatasetSequence(), m.getSequence()); + assertSame(cds2Dss, m.getSequence()); assertEquals(7, m.getStart()); assertEquals(9, m.getEnd()); } @@ -1169,7 +1219,7 @@ public class AlignmentUtilsTests SequenceI cdsSeq = cds.get(0); assertEquals("GGGTTT", cdsSeq.getSequenceAsString()); // assertEquals("dna1|A12345", cdsSeq.getName()); - assertEquals("dna1|pep1", cdsSeq.getName()); + assertEquals("CDS|dna1", cdsSeq.getName()); // assertEquals(1, cdsSeq.getDBRefs().length); // DBRefEntry cdsRef = cdsSeq.getDBRefs()[0]; // assertEquals("EMBLCDS", cdsRef.getSource()); @@ -1179,7 +1229,7 @@ public class AlignmentUtilsTests cdsSeq = cds.get(1); assertEquals("aaaccc", cdsSeq.getSequenceAsString()); // assertEquals("dna1|A12346", cdsSeq.getName()); - assertEquals("dna1|pep2", cdsSeq.getName()); + assertEquals("CDS|dna1", cdsSeq.getName()); // assertEquals(1, cdsSeq.getDBRefs().length); // cdsRef = cdsSeq.getDBRefs()[0]; // assertEquals("EMBLCDS", cdsRef.getSource()); @@ -1189,7 +1239,7 @@ public class AlignmentUtilsTests cdsSeq = cds.get(2); assertEquals("aaaTTT", cdsSeq.getSequenceAsString()); // assertEquals("dna1|A12347", cdsSeq.getName()); - assertEquals("dna1|pep3", cdsSeq.getName()); + assertEquals("CDS|dna1", cdsSeq.getName()); // assertEquals(1, cdsSeq.getDBRefs().length); // cdsRef = cdsSeq.getDBRefs()[0]; // assertEquals("EMBLCDS", cdsRef.getSource()); @@ -2255,7 +2305,7 @@ public class AlignmentUtilsTests * execute method under test to find CDS for EMBL peptides only */ AlignmentI cds = AlignmentUtils.makeCdsAlignment(new SequenceI[] { - dna1, dna2 }, dna.getDataset(), emblPeptides); + dna1, dna2 }, dna.getDataset(), emblPeptides.getSequencesArray()); assertEquals(2, cds.getSequences().size()); assertEquals("GGGTTT", cds.getSequenceAt(0).getSequenceAsString()); -- 1.7.10.2