X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;f=test%2Fjalview%2Fanalysis%2FAlignmentUtilsTests.java;h=cf6ef136c5e475fe613928d49210a8eb070744ff;hb=44d634a7732851685a37a502e48fbf2b967c6baf;hp=d229a3917a423d34c1152d6366efb9b842792c6c;hpb=e464da0ea20b0322846ce623598a16f2a1cba7ae;p=jalview.git diff --git a/test/jalview/analysis/AlignmentUtilsTests.java b/test/jalview/analysis/AlignmentUtilsTests.java index d229a39..cf6ef13 100644 --- a/test/jalview/analysis/AlignmentUtilsTests.java +++ b/test/jalview/analysis/AlignmentUtilsTests.java @@ -20,6 +20,7 @@ */ package jalview.analysis; +import static org.junit.Assert.assertNotEquals; import static org.testng.AssertJUnit.assertEquals; import static org.testng.AssertJUnit.assertFalse; import static org.testng.AssertJUnit.assertNotNull; @@ -41,6 +42,7 @@ import jalview.datamodel.SearchResultsI; import jalview.datamodel.Sequence; import jalview.datamodel.SequenceFeature; import jalview.datamodel.SequenceI; +import jalview.datamodel.AlignedCodonFrame.SequenceToSequenceMapping; import jalview.datamodel.features.SequenceFeatures; import jalview.gui.JvOptionPane; import jalview.io.AppletFormatAdapter; @@ -48,12 +50,14 @@ import jalview.io.DataSourceType; import jalview.io.FileFormat; import jalview.io.FileFormatI; import jalview.io.FormatAdapter; +import jalview.io.gff.SequenceOntologyI; import jalview.util.MapList; import jalview.util.MappingUtils; import java.io.IOException; import java.util.ArrayList; import java.util.Arrays; +import java.util.Iterator; import java.util.LinkedHashMap; import java.util.List; import java.util.Map; @@ -64,6 +68,8 @@ import org.testng.annotations.Test; public class AlignmentUtilsTests { + private static Sequence ts = new Sequence("short", + "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklm"); @BeforeClass(alwaysRun = true) public void setUpJvOptionPane() @@ -72,9 +78,6 @@ public class AlignmentUtilsTests JvOptionPane.setMockResponse(JvOptionPane.CANCEL_OPTION); } - private static Sequence ts = new Sequence("short", - "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklm"); - @Test(groups = { "Functional" }) public void testExpandContext() { @@ -264,14 +267,14 @@ public class AlignmentUtilsTests @Test(groups = { "Functional" }) public void testMapProteinAlignmentToCdna_noXrefs() throws IOException { - List protseqs = new ArrayList(); + List protseqs = new ArrayList<>(); protseqs.add(new Sequence("UNIPROT|V12345", "EIQ")); protseqs.add(new Sequence("UNIPROT|V12346", "EIQ")); protseqs.add(new Sequence("UNIPROT|V12347", "SAR")); AlignmentI protein = new Alignment(protseqs.toArray(new SequenceI[3])); protein.setDataset(null); - List dnaseqs = new ArrayList(); + List dnaseqs = new ArrayList<>(); dnaseqs.add(new Sequence("EMBL|A11111", "TCAGCACGC")); // = SAR dnaseqs.add(new Sequence("EMBL|A22222", "GAGATACAA")); // = EIQ dnaseqs.add(new Sequence("EMBL|A33333", "GAAATCCAG")); // = EIQ @@ -508,7 +511,7 @@ public class AlignmentUtilsTests acf.addMap(dna1.getDatasetSequence(), prot1.getDatasetSequence(), map); acf.addMap(dna2.getDatasetSequence(), prot2.getDatasetSequence(), map); acf.addMap(dna3.getDatasetSequence(), prot3.getDatasetSequence(), map); - ArrayList acfs = new ArrayList(); + ArrayList acfs = new ArrayList<>(); acfs.add(acf); protein.setCodonFrames(acfs); @@ -606,14 +609,14 @@ public class AlignmentUtilsTests public void testMapProteinAlignmentToCdna_withStartAndStopCodons() throws IOException { - List protseqs = new ArrayList(); + List protseqs = new ArrayList<>(); protseqs.add(new Sequence("UNIPROT|V12345", "EIQ")); protseqs.add(new Sequence("UNIPROT|V12346", "EIQ")); protseqs.add(new Sequence("UNIPROT|V12347", "SAR")); AlignmentI protein = new Alignment(protseqs.toArray(new SequenceI[3])); protein.setDataset(null); - List dnaseqs = new ArrayList(); + List dnaseqs = new ArrayList<>(); // start + SAR: dnaseqs.add(new Sequence("EMBL|A11111", "ATGTCAGCACGC")); // = EIQ + stop @@ -698,14 +701,14 @@ public class AlignmentUtilsTests @Test(groups = { "Functional" }) public void testMapProteinAlignmentToCdna_withXrefs() throws IOException { - List protseqs = new ArrayList(); + List protseqs = new ArrayList<>(); protseqs.add(new Sequence("UNIPROT|V12345", "EIQ")); protseqs.add(new Sequence("UNIPROT|V12346", "EIQ")); protseqs.add(new Sequence("UNIPROT|V12347", "SAR")); AlignmentI protein = new Alignment(protseqs.toArray(new SequenceI[3])); protein.setDataset(null); - List dnaseqs = new ArrayList(); + List dnaseqs = new ArrayList<>(); dnaseqs.add(new Sequence("EMBL|A11111", "TCAGCACGC")); // = SAR dnaseqs.add(new Sequence("EMBL|A22222", "ATGGAGATACAA")); // = start + EIQ dnaseqs.add(new Sequence("EMBL|A33333", "GAAATCCAG")); // = EIQ @@ -775,14 +778,14 @@ public class AlignmentUtilsTests public void testMapProteinAlignmentToCdna_prioritiseXrefs() throws IOException { - List protseqs = new ArrayList(); + List protseqs = new ArrayList<>(); protseqs.add(new Sequence("UNIPROT|V12345", "EIQ")); protseqs.add(new Sequence("UNIPROT|V12346", "EIQ")); AlignmentI protein = new Alignment( protseqs.toArray(new SequenceI[protseqs.size()])); protein.setDataset(null); - List dnaseqs = new ArrayList(); + List dnaseqs = new ArrayList<>(); dnaseqs.add(new Sequence("EMBL|A11111", "GAAATCCAG")); // = EIQ dnaseqs.add(new Sequence("EMBL|A22222", "GAAATTCAG")); // = EIQ AlignmentI cdna = new Alignment(dnaseqs.toArray(new SequenceI[dnaseqs @@ -849,8 +852,8 @@ public class AlignmentUtilsTests al.addAnnotation(ann4); // Temp for seq1 al.addAnnotation(ann5); // Temp for seq2 al.addAnnotation(ann6); // Temp for no sequence - List types = new ArrayList(); - List scope = new ArrayList(); + List types = new ArrayList<>(); + List scope = new ArrayList<>(); /* * Set all sequence related Structure to hidden (ann1, ann2) @@ -1784,10 +1787,13 @@ public class AlignmentUtilsTests map = new MapList(new int[] { 9, 11 }, new int[] { 2, 2 }, 3, 1); acf.addMap(dna3.getDatasetSequence(), prot3.getDatasetSequence(), map); - ArrayList acfs = new ArrayList(); + ArrayList acfs = new ArrayList<>(); acfs.add(acf); protein.setCodonFrames(acfs); - + Iterator protseq = protein.getSequences().iterator(); + for (SequenceI dnaseq:dna.getSequences()) { + assertCanResolveProteinCDS(dnaseq,protseq.next(),protein); + } /* * verify X is included in the aligned proteins, and placed just * before the first mapped residue @@ -1800,6 +1806,43 @@ public class AlignmentUtilsTests } /** + * assert that we can resolve the protein product in the given alignment given a DNA sequence with CDS mapping + * @param dnaseq + * @param protein + */ + private void assertCanResolveProteinCDS(SequenceI dnaseq, SequenceI expProtein, AlignmentI protein) + { + // try a few different methods to check all work + SequenceI aprot=null; + for (AlignedCodonFrame cf:protein.getCodonFrame(dnaseq)) + { + aprot=cf.getAaForDnaSeq(dnaseq); + if (aprot!=null) + { + assertTrue("getAaForDnaSeq didn't return expected protein sequence",aprot!=expProtein); + break; + } + } + assertNotNull("Didn't locate any proteins via AlignmentI.getCodonFrame .. AlignCodonFrame.getAaForDnaSeq", aprot); + // try mapping utils - + List mu_mappings=MappingUtils.findMappingsForSequence(dnaseq, protein.getCodonFrames()); + assertNotNull("No mappings found for dnaseq in protein alignment via MappingUtils.findMappingsForSequence",mu_mappings); + assertNotEquals("No mappings found for dnaseq in protein alignment via MappingUtils.findMappingsForSequence",0,mu_mappings.size()); + SequenceI mu_alignedprot=null; + List foundMap=null; + for (AlignedCodonFrame cf:mu_mappings) + { + foundMap=new ArrayList<>(); + mu_alignedprot = cf.findAlignedSequence(dnaseq, protein,foundMap); + if (mu_alignedprot!=null) { + break; + } + } + assertNotNull("Didn't locate proteins via MappingUtils.findMappingsForSequence",mu_alignedprot); + assertTrue("findAlignedSequence didn't return expected protein sequence",mu_alignedprot==expProtein); + } + + /** * Tests for the method that maps the subset of a dna sequence that has CDS * (or subtype) feature - case where the start codon is incomplete. */ @@ -1870,330 +1913,6 @@ public class AlignmentUtilsTests } /** - * Test the method that computes a map of codon variants for each protein - * position from "sequence_variant" features on dna - */ - @Test(groups = "Functional") - public void testBuildDnaVariantsMap() - { - SequenceI dna = new Sequence("dna", "atgAAATTTGGGCCCtag"); - MapList map = new MapList(new int[] { 1, 18 }, new int[] { 1, 5 }, 3, 1); - - /* - * first with no variants on dna - */ - LinkedHashMap[]> variantsMap = AlignmentUtils - .buildDnaVariantsMap(dna, map); - assertTrue(variantsMap.isEmpty()); - - /* - * single allele codon 1, on base 1 - */ - SequenceFeature sf1 = new SequenceFeature("sequence_variant", "", 1, 1, - 0f, null); - sf1.setValue("alleles", "T"); - sf1.setValue("ID", "sequence_variant:rs758803211"); - dna.addSequenceFeature(sf1); - - /* - * two alleles codon 2, on bases 2 and 3 (distinct variants) - */ - SequenceFeature sf2 = new SequenceFeature("sequence_variant", "", 5, 5, - 0f, null); - sf2.setValue("alleles", "T"); - sf2.setValue("ID", "sequence_variant:rs758803212"); - dna.addSequenceFeature(sf2); - SequenceFeature sf3 = new SequenceFeature("sequence_variant", "", 6, 6, - 0f, null); - sf3.setValue("alleles", "G"); - sf3.setValue("ID", "sequence_variant:rs758803213"); - dna.addSequenceFeature(sf3); - - /* - * two alleles codon 3, both on base 2 (one variant) - */ - SequenceFeature sf4 = new SequenceFeature("sequence_variant", "", 8, 8, - 0f, null); - sf4.setValue("alleles", "C, G"); - sf4.setValue("ID", "sequence_variant:rs758803214"); - dna.addSequenceFeature(sf4); - - // no alleles on codon 4 - - /* - * alleles on codon 5 on all 3 bases (distinct variants) - */ - SequenceFeature sf5 = new SequenceFeature("sequence_variant", "", 13, - 13, 0f, null); - sf5.setValue("alleles", "C, G"); // (C duplicates given base value) - sf5.setValue("ID", "sequence_variant:rs758803215"); - dna.addSequenceFeature(sf5); - SequenceFeature sf6 = new SequenceFeature("sequence_variant", "", 14, - 14, 0f, null); - sf6.setValue("alleles", "g, a"); // should force to upper-case - sf6.setValue("ID", "sequence_variant:rs758803216"); - dna.addSequenceFeature(sf6); - SequenceFeature sf7 = new SequenceFeature("sequence_variant", "", 15, - 15, 0f, null); - sf7.setValue("alleles", "A, T"); - sf7.setValue("ID", "sequence_variant:rs758803217"); - dna.addSequenceFeature(sf7); - - /* - * build map - expect variants on positions 1, 2, 3, 5 - */ - variantsMap = AlignmentUtils.buildDnaVariantsMap(dna, map); - assertEquals(4, variantsMap.size()); - - /* - * protein residue 1: variant on codon (ATG) base 1, not on 2 or 3 - */ - List[] pep1Variants = variantsMap.get(1); - assertEquals(3, pep1Variants.length); - assertEquals(1, pep1Variants[0].size()); - assertEquals("A", pep1Variants[0].get(0).base); // codon[1] base - assertSame(sf1, pep1Variants[0].get(0).variant); // codon[1] variant - assertEquals(1, pep1Variants[1].size()); - assertEquals("T", pep1Variants[1].get(0).base); // codon[2] base - assertNull(pep1Variants[1].get(0).variant); // no variant here - assertEquals(1, pep1Variants[2].size()); - assertEquals("G", pep1Variants[2].get(0).base); // codon[3] base - assertNull(pep1Variants[2].get(0).variant); // no variant here - - /* - * protein residue 2: variants on codon (AAA) bases 2 and 3 - */ - List[] pep2Variants = variantsMap.get(2); - assertEquals(3, pep2Variants.length); - assertEquals(1, pep2Variants[0].size()); - // codon[1] base recorded while processing variant on codon[2] - assertEquals("A", pep2Variants[0].get(0).base); - assertNull(pep2Variants[0].get(0).variant); // no variant here - // codon[2] base and variant: - assertEquals(1, pep2Variants[1].size()); - assertEquals("A", pep2Variants[1].get(0).base); - assertSame(sf2, pep2Variants[1].get(0).variant); - // codon[3] base was recorded when processing codon[2] variant - // and then the variant for codon[3] added to it - assertEquals(1, pep2Variants[2].size()); - assertEquals("A", pep2Variants[2].get(0).base); - assertSame(sf3, pep2Variants[2].get(0).variant); - - /* - * protein residue 3: variants on codon (TTT) base 2 only - */ - List[] pep3Variants = variantsMap.get(3); - assertEquals(3, pep3Variants.length); - assertEquals(1, pep3Variants[0].size()); - assertEquals("T", pep3Variants[0].get(0).base); // codon[1] base - assertNull(pep3Variants[0].get(0).variant); // no variant here - assertEquals(1, pep3Variants[1].size()); - assertEquals("T", pep3Variants[1].get(0).base); // codon[2] base - assertSame(sf4, pep3Variants[1].get(0).variant); // codon[2] variant - assertEquals(1, pep3Variants[2].size()); - assertEquals("T", pep3Variants[2].get(0).base); // codon[3] base - assertNull(pep3Variants[2].get(0).variant); // no variant here - - /* - * three variants on protein position 5 - */ - List[] pep5Variants = variantsMap.get(5); - assertEquals(3, pep5Variants.length); - assertEquals(1, pep5Variants[0].size()); - assertEquals("C", pep5Variants[0].get(0).base); // codon[1] base - assertSame(sf5, pep5Variants[0].get(0).variant); // codon[1] variant - assertEquals(1, pep5Variants[1].size()); - assertEquals("C", pep5Variants[1].get(0).base); // codon[2] base - assertSame(sf6, pep5Variants[1].get(0).variant); // codon[2] variant - assertEquals(1, pep5Variants[2].size()); - assertEquals("C", pep5Variants[2].get(0).base); // codon[3] base - assertSame(sf7, pep5Variants[2].get(0).variant); // codon[3] variant - } - - /** - * Tests for the method that computes all peptide variants given codon - * variants - */ - @Test(groups = "Functional") - public void testComputePeptideVariants() - { - /* - * scenario: AAATTTCCC codes for KFP - * variants: - * GAA -> E source: Ensembl - * CAA -> Q source: dbSNP - * AAG synonymous source: COSMIC - * AAT -> N source: Ensembl - * ...TTC synonymous source: dbSNP - * ......CAC,CGC -> H,R source: COSMIC - * (one variant with two alleles) - */ - SequenceI peptide = new Sequence("pep/10-12", "KFP"); - - /* - * two distinct variants for codon 1 position 1 - * second one has clinical significance - */ - String ensembl = "Ensembl"; - String dbSnp = "dbSNP"; - String cosmic = "COSMIC"; - SequenceFeature sf1 = new SequenceFeature("sequence_variant", "", 1, 1, - 0f, ensembl); - sf1.setValue("alleles", "A,G"); // GAA -> E - sf1.setValue("ID", "var1.125A>G"); - SequenceFeature sf2 = new SequenceFeature("sequence_variant", "", 1, 1, - 0f, dbSnp); - sf2.setValue("alleles", "A,C"); // CAA -> Q - sf2.setValue("ID", "var2"); - sf2.setValue("clinical_significance", "Dodgy"); - SequenceFeature sf3 = new SequenceFeature("sequence_variant", "", 3, 3, - 0f, cosmic); - sf3.setValue("alleles", "A,G"); // synonymous - sf3.setValue("ID", "var3"); - sf3.setValue("clinical_significance", "None"); - SequenceFeature sf4 = new SequenceFeature("sequence_variant", "", 3, 3, - 0f, ensembl); - sf4.setValue("alleles", "A,T"); // AAT -> N - sf4.setValue("ID", "sequence_variant:var4"); // prefix gets stripped off - sf4.setValue("clinical_significance", "Benign"); - SequenceFeature sf5 = new SequenceFeature("sequence_variant", "", 6, 6, - 0f, dbSnp); - sf5.setValue("alleles", "T,C"); // synonymous - sf5.setValue("ID", "var5"); - sf5.setValue("clinical_significance", "Bad"); - SequenceFeature sf6 = new SequenceFeature("sequence_variant", "", 8, 8, - 0f, cosmic); - sf6.setValue("alleles", "C,A,G"); // CAC,CGC -> H,R - sf6.setValue("ID", "var6"); - sf6.setValue("clinical_significance", "Good"); - - List codon1Variants = new ArrayList(); - List codon2Variants = new ArrayList(); - List codon3Variants = new ArrayList(); - List codonVariants[] = new ArrayList[3]; - codonVariants[0] = codon1Variants; - codonVariants[1] = codon2Variants; - codonVariants[2] = codon3Variants; - - /* - * compute variants for protein position 1 - */ - codon1Variants.add(new DnaVariant("A", sf1)); - codon1Variants.add(new DnaVariant("A", sf2)); - codon2Variants.add(new DnaVariant("A")); - codon2Variants.add(new DnaVariant("A")); - codon3Variants.add(new DnaVariant("A", sf3)); - codon3Variants.add(new DnaVariant("A", sf4)); - AlignmentUtils.computePeptideVariants(peptide, 1, codonVariants); - - /* - * compute variants for protein position 2 - */ - codon1Variants.clear(); - codon2Variants.clear(); - codon3Variants.clear(); - codon1Variants.add(new DnaVariant("T")); - codon2Variants.add(new DnaVariant("T")); - codon3Variants.add(new DnaVariant("T", sf5)); - AlignmentUtils.computePeptideVariants(peptide, 2, codonVariants); - - /* - * compute variants for protein position 3 - */ - codon1Variants.clear(); - codon2Variants.clear(); - codon3Variants.clear(); - codon1Variants.add(new DnaVariant("C")); - codon2Variants.add(new DnaVariant("C", sf6)); - codon3Variants.add(new DnaVariant("C")); - AlignmentUtils.computePeptideVariants(peptide, 3, codonVariants); - - /* - * verify added sequence features for - * var1 K -> E Ensembl - * var2 K -> Q dbSNP - * var4 K -> N Ensembl - * var6 P -> H COSMIC - * var6 P -> R COSMIC - */ - List sfs = peptide.getSequenceFeatures(); - SequenceFeatures.sortFeatures(sfs, true); - assertEquals(5, sfs.size()); - - /* - * features are sorted by start position ascending, but in no - * particular order where start positions match; asserts here - * simply match the data returned (the order is not important) - */ - SequenceFeature sf = sfs.get(0); - assertEquals(1, sf.getBegin()); - assertEquals(1, sf.getEnd()); - assertEquals("p.Lys1Asn", sf.getDescription()); - assertEquals("var4", sf.getValue("ID")); - assertEquals("Benign", sf.getValue("clinical_significance")); - assertEquals("ID=var4;clinical_significance=Benign", sf.getAttributes()); - assertEquals(1, sf.links.size()); - assertEquals( - "p.Lys1Asn var4|http://www.ensembl.org/Homo_sapiens/Variation/Summary?v=var4", - sf.links.get(0)); - assertEquals(ensembl, sf.getFeatureGroup()); - - sf = sfs.get(1); - assertEquals(1, sf.getBegin()); - assertEquals(1, sf.getEnd()); - assertEquals("p.Lys1Gln", sf.getDescription()); - assertEquals("var2", sf.getValue("ID")); - assertEquals("Dodgy", sf.getValue("clinical_significance")); - assertEquals("ID=var2;clinical_significance=Dodgy", sf.getAttributes()); - assertEquals(1, sf.links.size()); - assertEquals( - "p.Lys1Gln var2|http://www.ensembl.org/Homo_sapiens/Variation/Summary?v=var2", - sf.links.get(0)); - assertEquals(dbSnp, sf.getFeatureGroup()); - - sf = sfs.get(2); - assertEquals(1, sf.getBegin()); - assertEquals(1, sf.getEnd()); - assertEquals("p.Lys1Glu", sf.getDescription()); - assertEquals("var1.125A>G", sf.getValue("ID")); - assertNull(sf.getValue("clinical_significance")); - assertEquals("ID=var1.125A>G", sf.getAttributes()); - assertEquals(1, sf.links.size()); - // link to variation is urlencoded - assertEquals( - "p.Lys1Glu var1.125A>G|http://www.ensembl.org/Homo_sapiens/Variation/Summary?v=var1.125A%3EG", - sf.links.get(0)); - assertEquals(ensembl, sf.getFeatureGroup()); - - sf = sfs.get(3); - assertEquals(3, sf.getBegin()); - assertEquals(3, sf.getEnd()); - assertEquals("p.Pro3Arg", sf.getDescription()); - assertEquals("var6", sf.getValue("ID")); - assertEquals("Good", sf.getValue("clinical_significance")); - assertEquals("ID=var6;clinical_significance=Good", sf.getAttributes()); - assertEquals(1, sf.links.size()); - assertEquals( - "p.Pro3Arg var6|http://www.ensembl.org/Homo_sapiens/Variation/Summary?v=var6", - sf.links.get(0)); - assertEquals(cosmic, sf.getFeatureGroup()); - - // var5 generates two distinct protein variant features - sf = sfs.get(4); - assertEquals(3, sf.getBegin()); - assertEquals(3, sf.getEnd()); - assertEquals("p.Pro3His", sf.getDescription()); - assertEquals("var6", sf.getValue("ID")); - assertEquals("Good", sf.getValue("clinical_significance")); - assertEquals("ID=var6;clinical_significance=Good", sf.getAttributes()); - assertEquals(1, sf.links.size()); - assertEquals( - "p.Pro3His var6|http://www.ensembl.org/Homo_sapiens/Variation/Summary?v=var6", - sf.links.get(0)); - assertEquals(cosmic, sf.getFeatureGroup()); - } - - /** * Tests for the method that maps the subset of a dna sequence that has CDS * (or subtype) feature, with CDS strand = '-' (reverse) */ @@ -2309,7 +2028,7 @@ public class AlignmentUtilsTests seq1.createDatasetSequence(); Mapping mapping = new Mapping(seq1, new MapList( new int[] { 3, 6, 9, 10 }, new int[] { 1, 6 }, 1, 1)); - Map> map = new TreeMap>(); + Map> map = new TreeMap<>(); AlignmentUtils.addMappedPositions(seq1, from, mapping, map); /* @@ -2341,7 +2060,7 @@ public class AlignmentUtilsTests seq1.createDatasetSequence(); Mapping mapping = new Mapping(seq1, new MapList( new int[] { 3, 6, 9, 10 }, new int[] { 1, 6 }, 1, 1)); - Map> map = new TreeMap>(); + Map> map = new TreeMap<>(); AlignmentUtils.addMappedPositions(seq1, from, mapping, map); /* @@ -2513,9 +2232,16 @@ public class AlignmentUtilsTests AlignmentI al2 = new Alignment(new SequenceI[] { dna3, dna4 }); ((Alignment) al2).createDatasetAlignment(); + /* + * alignment removes gapped columns (two internal, two trailing) + */ assertTrue(AlignmentUtils.alignAsSameSequences(al1, al2)); - assertEquals(seq1, al1.getSequenceAt(0).getSequenceAsString()); - assertEquals(seq2, al1.getSequenceAt(1).getSequenceAsString()); + String aligned1 = "-cc-GG-GTTT-aaa"; + assertEquals(aligned1, + al1.getSequenceAt(0).getSequenceAsString()); + String aligned2 = "C--C-Cgg-gtttAAA"; + assertEquals(aligned2, + al1.getSequenceAt(1).getSequenceAsString()); /* * add another sequence to 'aligned' - should still succeed, since @@ -2525,8 +2251,8 @@ public class AlignmentUtilsTests dna5.createDatasetSequence(); al2.addSequence(dna5); assertTrue(AlignmentUtils.alignAsSameSequences(al1, al2)); - assertEquals(seq1, al1.getSequenceAt(0).getSequenceAsString()); - assertEquals(seq2, al1.getSequenceAt(1).getSequenceAsString()); + assertEquals(aligned1, al1.getSequenceAt(0).getSequenceAsString()); + assertEquals(aligned2, al1.getSequenceAt(1).getSequenceAsString()); /* * add another sequence to 'unaligned' - should fail, since now not @@ -2544,15 +2270,15 @@ public class AlignmentUtilsTests { SequenceI dna1 = new Sequence("dna1", "cccGGGTTTaaa"); SequenceI dna2 = new Sequence("dna2", "CCCgggtttAAA"); - SequenceI as1 = dna1.deriveSequence(); - SequenceI as2 = dna1.deriveSequence().getSubSequence(3, 7); - SequenceI as3 = dna2.deriveSequence(); + SequenceI as1 = dna1.deriveSequence(); // cccGGGTTTaaa/1-12 + SequenceI as2 = dna1.deriveSequence().getSubSequence(3, 7); // GGGT/4-7 + SequenceI as3 = dna2.deriveSequence(); // CCCgggtttAAA/1-12 as1.insertCharAt(6, 5, '-'); - String s_as1 = as1.getSequenceAsString(); + assertEquals("cccGGG-----TTTaaa", as1.getSequenceAsString()); as2.insertCharAt(6, 5, '-'); - String s_as2 = as2.getSequenceAsString(); - as3.insertCharAt(6, 5, '-'); - String s_as3 = as3.getSequenceAsString(); + assertEquals("GGGT-----", as2.getSequenceAsString()); + as3.insertCharAt(3, 5, '-'); + assertEquals("CCC-----gggtttAAA", as3.getSequenceAsString()); AlignmentI aligned = new Alignment(new SequenceI[] { as1, as2, as3 }); // why do we need to cast this still ? @@ -2564,10 +2290,13 @@ public class AlignmentUtilsTests uas3 }); ((Alignment) tobealigned).createDatasetAlignment(); + /* + * alignAs lines up dataset sequences and removes empty columns (two) + */ assertTrue(AlignmentUtils.alignAsSameSequences(tobealigned, aligned)); - assertEquals(s_as1, uas1.getSequenceAsString()); - assertEquals(s_as2, uas2.getSequenceAsString()); - assertEquals(s_as3, uas3.getSequenceAsString()); + assertEquals("cccGGG---TTTaaa", uas1.getSequenceAsString()); + assertEquals("GGGT", uas2.getSequenceAsString()); + assertEquals("CCC---gggtttAAA", uas3.getSequenceAsString()); } @Test(groups = { "Functional" }) @@ -2606,19 +2335,19 @@ public class AlignmentUtilsTests * transcript 'CDS' is 10-16, 17-21 * which is 'gene' 158-164, 210-214 */ - MapList toMap = toLoci.getMap(); + MapList toMap = toLoci.getMapping(); assertEquals(1, toMap.getFromRanges().size()); assertEquals(2, toMap.getFromRanges().get(0).length); assertEquals(1, toMap.getFromRanges().get(0)[0]); assertEquals(12, toMap.getFromRanges().get(0)[1]); - assertEquals(1, toMap.getToRanges().size()); - assertEquals(4, toMap.getToRanges().get(0).length); + assertEquals(2, toMap.getToRanges().size()); + assertEquals(2, toMap.getToRanges().get(0).length); assertEquals(158, toMap.getToRanges().get(0)[0]); assertEquals(164, toMap.getToRanges().get(0)[1]); - assertEquals(210, toMap.getToRanges().get(0)[2]); - assertEquals(214, toMap.getToRanges().get(0)[3]); + assertEquals(210, toMap.getToRanges().get(1)[0]); + assertEquals(214, toMap.getToRanges().get(1)[1]); // or summarised as (but toString might change in future): - assertEquals("[ [1, 12] ] 1:1 to [ [158, 164, 210, 214] ]", + assertEquals("[ [1, 12] ] 1:1 to [ [158, 164] [210, 214] ]", toMap.toString()); /* @@ -2629,8 +2358,248 @@ public class AlignmentUtilsTests AlignmentUtils.transferGeneLoci(from, map, to); assertEquals("GRCh38", toLoci.getAssemblyId()); assertEquals("7", toLoci.getChromosomeId()); - toMap = toLoci.getMap(); - assertEquals("[ [1, 12] ] 1:1 to [ [158, 164, 210, 214] ]", + toMap = toLoci.getMapping(); + assertEquals("[ [1, 12] ] 1:1 to [ [158, 164] [210, 214] ]", toMap.toString()); } + + /** + * Tests for the method that maps nucleotide to protein based on CDS features + */ + @Test(groups = "Functional") + public void testMapCdsToProtein() + { + SequenceI peptide = new Sequence("pep", "KLQ"); + + /* + * Case 1: CDS 3 times length of peptide + * NB method only checks lengths match, not translation + */ + SequenceI dna = new Sequence("dna", "AACGacgtCTCCT"); + dna.createDatasetSequence(); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 1, 4, null)); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 9, 13, null)); + MapList ml = AlignmentUtils.mapCdsToProtein(dna, peptide); + assertEquals(3, ml.getFromRatio()); + assertEquals(1, ml.getToRatio()); + assertEquals("[[1, 3]]", + Arrays.deepToString(ml.getToRanges().toArray())); + assertEquals("[[1, 4], [9, 13]]", + Arrays.deepToString(ml.getFromRanges().toArray())); + + /* + * Case 2: CDS 3 times length of peptide + stop codon + * (note code does not currently check trailing codon is a stop codon) + */ + dna = new Sequence("dna", "AACGacgtCTCCTCCC"); + dna.createDatasetSequence(); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 1, 4, null)); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 9, 16, null)); + ml = AlignmentUtils.mapCdsToProtein(dna, peptide); + assertEquals(3, ml.getFromRatio()); + assertEquals(1, ml.getToRatio()); + assertEquals("[[1, 3]]", + Arrays.deepToString(ml.getToRanges().toArray())); + assertEquals("[[1, 4], [9, 13]]", + Arrays.deepToString(ml.getFromRanges().toArray())); + + /* + * Case 3: CDS longer than 3 * peptide + stop codon - no mapping is made + */ + dna = new Sequence("dna", "AACGacgtCTCCTTGATCA"); + dna.createDatasetSequence(); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 1, 4, null)); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 9, 19, null)); + ml = AlignmentUtils.mapCdsToProtein(dna, peptide); + assertNull(ml); + + /* + * Case 4: CDS shorter than 3 * peptide - no mapping is made + */ + dna = new Sequence("dna", "AACGacgtCTCC"); + dna.createDatasetSequence(); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 1, 4, null)); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 9, 12, null)); + ml = AlignmentUtils.mapCdsToProtein(dna, peptide); + assertNull(ml); + + /* + * Case 5: CDS 3 times length of peptide + part codon - mapping is truncated + */ + dna = new Sequence("dna", "AACGacgtCTCCTTG"); + dna.createDatasetSequence(); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 1, 4, null)); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 9, 15, null)); + ml = AlignmentUtils.mapCdsToProtein(dna, peptide); + assertEquals(3, ml.getFromRatio()); + assertEquals(1, ml.getToRatio()); + assertEquals("[[1, 3]]", + Arrays.deepToString(ml.getToRanges().toArray())); + assertEquals("[[1, 4], [9, 13]]", + Arrays.deepToString(ml.getFromRanges().toArray())); + + /* + * Case 6: incomplete start codon corresponding to X in peptide + */ + dna = new Sequence("dna", "ACGacgtCTCCTTGG"); + dna.createDatasetSequence(); + SequenceFeature sf = new SequenceFeature("CDS", "", 1, 3, null); + sf.setPhase("2"); // skip 2 positions (AC) to start of next codon (GCT) + dna.addSequenceFeature(sf); + dna.addSequenceFeature(new SequenceFeature("CDS", "", 8, 15, null)); + peptide = new Sequence("pep", "XLQ"); + ml = AlignmentUtils.mapCdsToProtein(dna, peptide); + assertEquals("[[2, 3]]", + Arrays.deepToString(ml.getToRanges().toArray())); + assertEquals("[[3, 3], [8, 12]]", + Arrays.deepToString(ml.getFromRanges().toArray())); + } + + /** + * Tests for the method that locates the CDS sequence that has a mapping to + * the given protein. That is, given a transcript-to-peptide mapping, find the + * cds-to-peptide mapping that relates to both, and return the CDS sequence. + */ + @Test + public void testFindCdsForProtein() + { + List mappings = new ArrayList<>(); + AlignedCodonFrame acf1 = new AlignedCodonFrame(); + mappings.add(acf1); + + SequenceI dna1 = new Sequence("dna1", "cgatATcgGCTATCTATGacg"); + dna1.createDatasetSequence(); + + // NB we currently exclude STOP codon from CDS sequences + // the test would need to change if this changes in future + SequenceI cds1 = new Sequence("cds1", "ATGCTATCT"); + cds1.createDatasetSequence(); + + SequenceI pep1 = new Sequence("pep1", "MLS"); + pep1.createDatasetSequence(); + List seqMappings = new ArrayList<>(); + MapList mapList = new MapList( + new int[] + { 5, 6, 9, 15 }, new int[] { 1, 3 }, 3, 1); + Mapping dnaToPeptide = new Mapping(pep1.getDatasetSequence(), mapList); + + // add dna to peptide mapping + seqMappings.add(acf1); + acf1.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(), + mapList); + + /* + * first case - no dna-to-CDS mapping exists - search fails + */ + SequenceI seq = AlignmentUtils.findCdsForProtein(mappings, dna1, + seqMappings, dnaToPeptide); + assertNull(seq); + + /* + * second case - CDS-to-peptide mapping exists but no dna-to-CDS + * - search fails + */ + // todo this test fails if the mapping is added to acf1, not acf2 + // need to tidy up use of lists of mappings in AlignedCodonFrame + AlignedCodonFrame acf2 = new AlignedCodonFrame(); + mappings.add(acf2); + MapList cdsToPeptideMapping = new MapList(new int[] + { 1, 9 }, new int[] { 1, 3 }, 3, 1); + acf2.addMap(cds1.getDatasetSequence(), pep1.getDatasetSequence(), + cdsToPeptideMapping); + assertNull(AlignmentUtils.findCdsForProtein(mappings, dna1, seqMappings, + dnaToPeptide)); + + /* + * third case - add dna-to-CDS mapping - CDS is now found! + */ + MapList dnaToCdsMapping = new MapList(new int[] { 5, 6, 9, 15 }, + new int[] + { 1, 9 }, 1, 1); + acf1.addMap(dna1.getDatasetSequence(), cds1.getDatasetSequence(), + dnaToCdsMapping); + seq = AlignmentUtils.findCdsForProtein(mappings, dna1, seqMappings, + dnaToPeptide); + assertSame(seq, cds1.getDatasetSequence()); + } + + /** + * Tests for the method that locates the CDS sequence that has a mapping to + * the given protein. That is, given a transcript-to-peptide mapping, find the + * cds-to-peptide mapping that relates to both, and return the CDS sequence. + * This test is for the case where transcript and CDS are the same length. + */ + @Test + public void testFindCdsForProtein_noUTR() + { + List mappings = new ArrayList<>(); + AlignedCodonFrame acf1 = new AlignedCodonFrame(); + mappings.add(acf1); + + SequenceI dna1 = new Sequence("dna1", "ATGCTATCTTAA"); + dna1.createDatasetSequence(); + + // NB we currently exclude STOP codon from CDS sequences + // the test would need to change if this changes in future + SequenceI cds1 = new Sequence("cds1", "ATGCTATCT"); + cds1.createDatasetSequence(); + + SequenceI pep1 = new Sequence("pep1", "MLS"); + pep1.createDatasetSequence(); + List seqMappings = new ArrayList<>(); + MapList mapList = new MapList( + new int[] + { 1, 9 }, new int[] { 1, 3 }, 3, 1); + Mapping dnaToPeptide = new Mapping(pep1.getDatasetSequence(), mapList); + + // add dna to peptide mapping + seqMappings.add(acf1); + acf1.addMap(dna1.getDatasetSequence(), pep1.getDatasetSequence(), + mapList); + + /* + * first case - transcript lacks CDS features - it appears to be + * the CDS sequence and is returned + */ + SequenceI seq = AlignmentUtils.findCdsForProtein(mappings, dna1, + seqMappings, dnaToPeptide); + assertSame(seq, dna1.getDatasetSequence()); + + /* + * second case - transcript has CDS feature - this means it is + * not returned as a match for CDS (CDS sequences don't have CDS features) + */ + dna1.addSequenceFeature( + new SequenceFeature(SequenceOntologyI.CDS, "cds", 1, 12, null)); + seq = AlignmentUtils.findCdsForProtein(mappings, dna1, seqMappings, + dnaToPeptide); + assertNull(seq); + + /* + * third case - CDS-to-peptide mapping exists but no dna-to-CDS + * - search fails + */ + // todo this test fails if the mapping is added to acf1, not acf2 + // need to tidy up use of lists of mappings in AlignedCodonFrame + AlignedCodonFrame acf2 = new AlignedCodonFrame(); + mappings.add(acf2); + MapList cdsToPeptideMapping = new MapList(new int[] + { 1, 9 }, new int[] { 1, 3 }, 3, 1); + acf2.addMap(cds1.getDatasetSequence(), pep1.getDatasetSequence(), + cdsToPeptideMapping); + assertNull(AlignmentUtils.findCdsForProtein(mappings, dna1, seqMappings, + dnaToPeptide)); + + /* + * fourth case - add dna-to-CDS mapping - CDS is now found! + */ + MapList dnaToCdsMapping = new MapList(new int[] { 1, 9 }, + new int[] + { 1, 9 }, 1, 1); + acf1.addMap(dna1.getDatasetSequence(), cds1.getDatasetSequence(), + dnaToCdsMapping); + seq = AlignmentUtils.findCdsForProtein(mappings, dna1, seqMappings, + dnaToPeptide); + assertSame(seq, cds1.getDatasetSequence()); + } }