X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;f=src%2Fjalview%2Fanalysis%2FCrossRef.java;h=7e77fc1e0903799660db67205966ddcf768a1897;hb=0ebcd487f08f9873987d9a32edce7df47aa27927;hp=68f6c93ed1f3eb760348f8267f57baba7e824564;hpb=ecc50d775a514f9840cdcc10a9d5b1e49c60582c;p=jalview.git diff --git a/src/jalview/analysis/CrossRef.java b/src/jalview/analysis/CrossRef.java index 68f6c93..7e77fc1 100644 --- a/src/jalview/analysis/CrossRef.java +++ b/src/jalview/analysis/CrossRef.java @@ -29,20 +29,13 @@ import jalview.datamodel.Mapping; import jalview.datamodel.Sequence; import jalview.datamodel.SequenceFeature; import jalview.datamodel.SequenceI; -import jalview.io.gff.SequenceOntology; -import jalview.schemes.ResidueProperties; import jalview.util.DBRefUtils; import jalview.util.MapList; -import jalview.util.MappingUtils; -import jalview.util.StringUtils; import jalview.ws.SequenceFetcher; import jalview.ws.seqfetcher.ASequenceFetcher; import java.util.ArrayList; -import java.util.Collections; -import java.util.LinkedHashMap; import java.util.List; -import java.util.Map.Entry; import java.util.Vector; /** @@ -54,6 +47,27 @@ import java.util.Vector; */ public class CrossRef { + /* + * A sub-class that ignores Parent attribute when comparing sequence + * features. This avoids 'duplicate' CDS features that only + * differ in their parent Transcript ids. + */ + class MySequenceFeature extends SequenceFeature + { + private SequenceFeature feat; + + MySequenceFeature(SequenceFeature sf) + { + this.feat = sf; + } + + @Override + public boolean equals(Object o) + { + return feat.equals(o, true); + } + } + /** * Select just the DNA or protein references for a protein or dna sequence * @@ -107,7 +121,7 @@ public class CrossRef { dss = dss.getDatasetSequence(); } - DBRefEntry[] rfs = findXDbRefs(dna, dss.getDBRef()); + DBRefEntry[] rfs = findXDbRefs(dna, dss.getDBRefs()); if (rfs != null) { for (DBRefEntry ref : rfs) @@ -121,13 +135,13 @@ public class CrossRef if (dataset != null) { // search for references to this sequence's direct references. - DBRefEntry[] lrfs = CrossRef.findXDbRefs(!dna, seq.getDBRef()); + DBRefEntry[] lrfs = CrossRef.findXDbRefs(!dna, seq.getDBRefs()); List rseqs = new ArrayList(); CrossRef.searchDatasetXrefs(seq, !dna, lrfs, dataset, rseqs, null); // don't need to specify codon frame for mapping here for (SequenceI rs : rseqs) { - DBRefEntry[] xrs = findXDbRefs(dna, rs.getDBRef()); + DBRefEntry[] xrs = findXDbRefs(dna, rs.getDBRefs()); if (xrs != null) { for (DBRefEntry ref : xrs) @@ -179,7 +193,7 @@ public class CrossRef Vector cseqs = new Vector(); for (int s = 0; s < seqs.length; s++) { - DBRefEntry[] cdna = findXDbRefs(true, seqs[s].getDBRef()); + DBRefEntry[] cdna = findXDbRefs(true, seqs[s].getDBRefs()); for (int c = 0; c < cdna.length; c++) { if (cdna[c].getSource().equals(DBRefSource.EMBLCDS)) @@ -206,30 +220,20 @@ public class CrossRef /** * - * @param dna - * @param seqs - * @return - */ - public static Alignment findXrefSequences(SequenceI[] seqs, boolean dna, - String source) - { - return findXrefSequences(seqs, dna, source, null); - } - - /** - * * @param seqs * sequences whose xrefs are being retrieved * @param dna * true if sequences are nucleotide * @param source - * @param dataset - * alignment to search for product sequences. + * @param al + * alignment to search for cross-referenced sequences (and possibly + * add to) * @return products (as dataset sequences) */ - public static Alignment findXrefSequences(SequenceI[] seqs, boolean dna, - String source, AlignmentI dataset) + public static Alignment findXrefSequences(SequenceI[] seqs, + final boolean dna, final String source, AlignmentI al) { + AlignmentI dataset = al.getDataset() == null ? al : al.getDataset(); List rseqs = new ArrayList(); AlignedCodonFrame cf = new AlignedCodonFrame(); for (SequenceI seq : seqs) @@ -240,12 +244,12 @@ public class CrossRef dss = dss.getDatasetSequence(); } boolean found = false; - DBRefEntry[] xrfs = CrossRef.findXDbRefs(dna, dss.getDBRef()); + DBRefEntry[] xrfs = CrossRef.findXDbRefs(dna, dss.getDBRefs()); if ((xrfs == null || xrfs.length == 0) && dataset != null) { System.out.println("Attempting to find ds Xrefs refs."); // FIXME should be dss not seq here? - DBRefEntry[] lrfs = CrossRef.findXDbRefs(!dna, seq.getDBRef()); + DBRefEntry[] lrfs = CrossRef.findXDbRefs(!dna, seq.getDBRefs()); // less ambiguous would be a 'find primary dbRefEntry' method. // filter for desired source xref here found = CrossRef.searchDatasetXrefs(dss, !dna, lrfs, dataset, @@ -278,12 +282,6 @@ public class CrossRef // map should be from protein seq to its coding dna cf.addMap(rsq, dss, xref.getMap().getMap().getInverse()); } - - /* - * compute peptide variants from dna variants - */ - rsq.createDatasetSequence(); - computeProteinVariants(seq, rsq, xref.getMap().getMap()); } found = true; } @@ -294,7 +292,8 @@ public class CrossRef // xrefs on this sequence. if (dataset != null) { - found |= searchDataset(dss, xref, dataset, rseqs, cf); // ,false,!dna); + found |= searchDataset(dss, xref, dataset, rseqs, cf, false, + !dna); if (found) { xrfs[r] = null; // we've recovered seqs for this one. @@ -328,8 +327,8 @@ public class CrossRef } if (l > 0) { - System.out - .println("Attempting to retrieve cross referenced sequences."); + // System.out + // .println("Attempting to retrieve cross referenced sequences."); DBRefEntry[] t = new DBRefEntry[l]; l = 0; for (int r = 0; r < xrfs.length; r++) @@ -352,12 +351,19 @@ public class CrossRef + seq.getName()); e.printStackTrace(); } + if (retrieved != null) { + updateDbrefMappings(dna, seq, xrfs, retrieved, cf); + + SequenceIdMatcher matcher = new SequenceIdMatcher( + dataset.getSequences()); + List copiedFeatures = new ArrayList(); + CrossRef me = new CrossRef(); for (int rs = 0; rs < retrieved.length; rs++) { // TODO: examine each sequence for 'redundancy' - DBRefEntry[] dbr = retrieved[rs].getDBRef(); + DBRefEntry[] dbr = retrieved[rs].getDBRefs(); if (dbr != null && dbr.length > 0) { for (int di = 0; di < dbr.length; di++) @@ -369,8 +375,25 @@ public class CrossRef { if (map.getTo() != null && map.getMap() != null) { - // should search the local dataset to find any existing - // candidates for To ! + SequenceI matched = matcher + .findIdMatch(map.getTo()); + if (matched != null) + { + /* + * already got an xref to this sequence; update this + * map to point to the same sequence, and add + * any new dbrefs to it + */ + for (DBRefEntry ref : map.getTo().getDBRefs()) + { + matched.addDBRef(ref); // add or update mapping + } + map.setTo(matched); + } + else + { + matcher.add(map.getTo()); + } try { // compare ms with dss and replace with dss in mapping @@ -379,17 +402,54 @@ public class CrossRef int sf = map.getMap().getToLowest(); int st = map.getMap().getToHighest(); SequenceI mappedrg = ms.getSubSequence(sf, st); - SequenceI loc = dss.getSubSequence(sf, st); + // SequenceI loc = dss.getSubSequence(sf, st); if (mappedrg.getLength() > 0 - && mappedrg.getSequenceAsString().equals( - loc.getSequenceAsString())) + && ms.getSequenceAsString().equals( + dss.getSequenceAsString())) + // && mappedrg.getSequenceAsString().equals( + // loc.getSequenceAsString())) { - System.err - .println("Mapping updated for retrieved crossreference"); + String msg = "Mapping updated from " + + ms.getName() + + " to retrieved crossreference " + + dss.getName(); + System.out.println(msg); // method to update all refs of existing To on // retrieved sequence with dss and merge any props // on To onto dss. map.setTo(dss); + /* + * copy sequence features as well, avoiding + * duplication (e.g. same variation from 2 + * transcripts) + */ + SequenceFeature[] sfs = ms + .getSequenceFeatures(); + if (sfs != null) + { + for (SequenceFeature feat : sfs) + { + /* + * we override SequenceFeature.equals here (but + * not elsewhere) to ignore Parent attribute + * TODO not quite working yet! + */ + if (!copiedFeatures + .contains(me.new MySequenceFeature( + feat))) + { + dss.addSequenceFeature(feat); + copiedFeatures.add(feat); + } + } + } + cf.addMap(retrieved[rs].getDatasetSequence(), + dss, map.getMap()); + } + else + { + cf.addMap(retrieved[rs].getDatasetSequence(), + map.getTo(), map.getMap()); } } catch (Exception e) { @@ -413,9 +473,7 @@ public class CrossRef Alignment ral = null; if (rseqs.size() > 0) { - SequenceI[] rsqs = new SequenceI[rseqs.size()]; - rseqs.toArray(rsqs); - ral = new Alignment(rsqs); + ral = new Alignment(rseqs.toArray(new SequenceI[rseqs.size()])); if (cf != null && !cf.isEmpty()) { ral.addCodonFrame(cf); @@ -425,6 +483,69 @@ public class CrossRef } /** + * Updates any empty mappings in the cross-references with one to a compatible + * retrieved sequence if found, and adds any new mappings to the + * AlignedCodonFrame + * + * @param dna + * @param mapFrom + * @param xrefs + * @param retrieved + * @param acf + */ + static void updateDbrefMappings(boolean dna, SequenceI mapFrom, + DBRefEntry[] xrefs, SequenceI[] retrieved, AlignedCodonFrame acf) + { + SequenceIdMatcher matcher = new SequenceIdMatcher(retrieved); + for (DBRefEntry xref : xrefs) + { + if (!xref.hasMap()) + { + String targetSeqName = xref.getSource() + "|" + + xref.getAccessionId(); + SequenceI[] matches = matcher.findAllIdMatches(targetSeqName); + if (matches == null) + { + return; + } + for (SequenceI seq : matches) + { + MapList mapping = null; + if (dna) + { + mapping = AlignmentUtils.mapCdnaToProtein(seq, mapFrom); + } + else + { + mapping = AlignmentUtils.mapCdnaToProtein(mapFrom, seq); + if (mapping != null) + { + mapping = mapping.getInverse(); + } + } + if (mapping != null) + { + xref.setMap(new Mapping(seq, mapping)); + if (dna) + { + AlignmentUtils.computeProteinFeatures(mapFrom, seq, mapping); + } + if (dna) + { + acf.addMap(mapFrom, seq, mapping); + } + else + { + acf.addMap(seq, mapFrom, mapping.getInverse()); + } + continue; + } + } + } + } + } + + /** * find references to lrfs in the cross-reference set of each sequence in * dataset (that is not equal to sequenceI) Identifies matching DBRefEntry * based on source and accession string only - Map and Version are nulled. @@ -529,7 +650,7 @@ public class CrossRef } // look for direct or indirect references in common - DBRefEntry[] poss = nxt.getDBRef(), cands = null; + DBRefEntry[] poss = nxt.getDBRefs(), cands = null; if (direct) { cands = jalview.util.DBRefUtils.searchRefs(poss, xrf); @@ -585,195 +706,6 @@ public class CrossRef } /** - * Computes variants in peptide product generated by variants in dna, and adds - * them as sequence_variant features on the protein sequence. Returns the - * number of variant features added. - * - * @param dnaSeq - * @param peptide - * @param dnaToProtein - */ - protected static int computeProteinVariants(SequenceI dnaSeq, - SequenceI peptide, MapList dnaToProtein) - { - /* - * map from peptide position to all variant features of the codon for it - * LinkedHashMap ensures we add the peptide features in sequence order - */ - LinkedHashMap variants = new LinkedHashMap(); - SequenceOntology so = SequenceOntology.getInstance(); - - SequenceFeature[] dnaFeatures = dnaSeq.getSequenceFeatures(); - if (dnaFeatures == null) - { - return 0; - } - - int[] lastCodon = null; - int lastPeptidePostion = 0; - - /* - * build a map of codon variations for peptides - */ - for (SequenceFeature sf : dnaFeatures) - { - int dnaCol = sf.getBegin(); - if (dnaCol != sf.getEnd()) - { - // not handling multi-locus variant features - continue; - } - if (so.isSequenceVariant(sf.getType())) - { - int[] mapsTo = dnaToProtein.locateInTo(dnaCol, dnaCol); - if (mapsTo == null) - { - // feature doesn't lie within coding region - continue; - } - int peptidePosition = mapsTo[0]; - String[][] codonVariants = variants.get(peptidePosition); - if (codonVariants == null) - { - codonVariants = new String[3][]; - variants.put(peptidePosition, codonVariants); - } - - /* - * extract dna variants to a string array - */ - String alls = (String) sf.getValue("alleles"); - if (alls == null) - { - continue; - } - String[] alleles = alls.split(","); - - /* - * get this peptides codon positions e.g. [3, 4, 5] or [4, 7, 10] - */ - int[] codon = peptidePosition == lastPeptidePostion ? lastCodon - : MappingUtils.flattenRanges(dnaToProtein.locateInFrom( - peptidePosition, peptidePosition)); - lastPeptidePostion = peptidePosition; - lastCodon = codon; - - /* - * save nucleotide (and this variant) for each codon position - */ - for (int codonPos = 0; codonPos < 3; codonPos++) - { - String nucleotide = String.valueOf(dnaSeq - .getCharAt(codon[codonPos] - 1)); - if (codon[codonPos] == dnaCol) - { - /* - * record current dna base and its alleles - */ - String[] dnaVariants = new String[alleles.length + 1]; - dnaVariants[0] = nucleotide; - System.arraycopy(alleles, 0, dnaVariants, 1, alleles.length); - codonVariants[codonPos] = dnaVariants; - } - else if (codonVariants[codonPos] == null) - { - /* - * record current dna base only - * (at least until we find any variation and overwrite it) - */ - codonVariants[codonPos] = new String[] { nucleotide }; - } - } - } - } - - /* - * scan codon variations, compute peptide variants and add to peptide sequence - */ - int count = 0; - for (Entry variant : variants.entrySet()) - { - int peptidePos = variant.getKey(); - String[][] codonVariants = variant.getValue(); - String residue = String.valueOf(peptide.getCharAt(peptidePos - 1)); // 0-based - List peptideVariants = computePeptideVariants(codonVariants, - residue); - if (!peptideVariants.isEmpty()) - { - Collections.sort(peptideVariants); - String desc = StringUtils.listToDelimitedString(peptideVariants, - ", "); - SequenceFeature sf = new SequenceFeature( - SequenceOntology.SEQUENCE_VARIANT, desc, peptidePos, - peptidePos, Float.NaN, null); - peptide.getDatasetSequence().addSequenceFeature(sf); - count++; - } - } - return count; - } - - /** - * Returns a non-redundant list of all peptide translations generated by the - * given dna variants, excluding the current residue value - * - * @param codonVariants - * an array of base values for codon positions 1, 2, 3 - * @param residue - * the current residue translation - * @return - */ - protected static List computePeptideVariants( - String[][] codonVariants, String residue) - { - List result = new ArrayList(); - for (String base1 : codonVariants[0]) - { - for (String base2 : codonVariants[1]) - { - for (String base3 : codonVariants[2]) - { - String codon = base1 + base2 + base3; - // TODO: report frameshift/insertion/deletion - // and multiple-base variants?! - String peptide = codon.contains("-") ? "-" : ResidueProperties - .codonTranslate(codon); - if (peptide != null && !result.contains(peptide) - && !peptide.equals(residue)) - { - result.add(peptide); - } - } - } - } - return result; - } - - /** - * Computes a list of all peptide variants given dna variants - * - * @param dnaSeq - * the coding dna sequence - * @param codonVariants - * variant features for each codon position (null if no variant) - * @param residue - * the canonical protein translation - * @return - */ - protected static List computePeptideVariants(SequenceI dnaSeq, - SequenceFeature[] codonVariants, String residue) - { - List result = new ArrayList(); - int[][] dnaVariants = new int[3][]; - for (int i = 0; i < 3; i++) - { - - } - // TODO Auto-generated method stub - return null; - } - - /** * precalculate different products that can be found for seqs in dataset and * return them. *