X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;f=src%2Fjalview%2Fext%2Fensembl%2FEnsemblSeqProxy.java;h=869a7028f1d2c8f5da9462a36877ad202359d1de;hb=1e8c7a9ab9f5da589d0aa2482fd2e3361c320d57;hp=8698b78bdb75acb09837ece49adfdfa70aaf0847;hpb=a9f0472fe6fd4737b47d7955d198e76923e6aabc;p=jalview.git diff --git a/src/jalview/ext/ensembl/EnsemblSeqProxy.java b/src/jalview/ext/ensembl/EnsemblSeqProxy.java index 8698b78..869a702 100644 --- a/src/jalview/ext/ensembl/EnsemblSeqProxy.java +++ b/src/jalview/ext/ensembl/EnsemblSeqProxy.java @@ -11,7 +11,8 @@ import jalview.datamodel.SequenceI; import jalview.exceptions.JalviewException; import jalview.io.FastaFile; import jalview.io.FileParse; -import jalview.io.gff.SequenceOntology; +import jalview.io.gff.SequenceOntologyFactory; +import jalview.io.gff.SequenceOntologyI; import jalview.schemes.ResidueProperties; import jalview.util.DBRefUtils; import jalview.util.MapList; @@ -36,39 +37,39 @@ import java.util.Map.Entry; */ public abstract class EnsemblSeqProxy extends EnsemblRestClient { + private static final List CROSS_REFERENCES = Arrays + .asList(new String[] { "CCDS", "Uniprot/SWISSPROT" }); + protected static final String CONSEQUENCE_TYPE = "consequence_type"; protected static final String PARENT = "Parent"; protected static final String ID = "ID"; - /* - * this needs special handling, as it isA sequence_variant in the - * Sequence Ontology, but behaves in Ensembl as if it isA transcript - */ - protected static final String NMD_VARIANT = "NMD_transcript_variant"; - protected static final String NAME = "Name"; + /* + * enum for 'type' parameter to the /sequence REST service + */ public enum EnsemblSeqType { /** - * type=genomic for the full dna including introns + * type=genomic to fetch full dna including introns */ GENOMIC("genomic"), /** - * type=cdna for transcribed dna including UTRs + * type=cdna to fetch dna including UTRs */ CDNA("cdna"), /** - * type=cds for coding dna excluding UTRs + * type=cds to fetch coding dna excluding UTRs */ CDS("cds"), /** - * type=protein for the peptide product sequence + * type=protein to fetch peptide product sequence */ PROTEIN("protein"); @@ -124,7 +125,6 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient @Override public AlignmentI getSequenceRecords(String query) throws Exception { - long now = System.currentTimeMillis(); // TODO use a String... query vararg instead? // danger: accession separator used as a regex here, a string elsewhere @@ -153,17 +153,15 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient + " chunks. Unexpected problem (" + r.getLocalizedMessage() + ")"; System.err.println(msg); - if (alignment != null) - { - break; // return what we got - } - else - { - throw new JalviewException(msg, r); - } + break; } } + if (alignment == null) + { + return null; + } + /* * fetch and transfer genomic sequence features, * fetch protein product and add as cross-reference @@ -173,9 +171,11 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient addFeaturesAndProduct(accId, alignment); } - inProgress = false; - System.out.println(getClass().getName() + " took " - + (System.currentTimeMillis() - now) + "ms to fetch"); + for (SequenceI seq : alignment.getSequences()) + { + getCrossReferences(seq); + } + return alignment; } @@ -201,7 +201,7 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient * get 'dummy' genomic sequence with exon, cds and variation features */ SequenceI genomicSequence = null; - EnsemblOverlap gffFetcher = new EnsemblOverlap(); + EnsemblFeatures gffFetcher = new EnsemblFeatures(); EnsemblFeatureType[] features = getFeaturesToFetch(); AlignmentI geneFeatures = gffFetcher.getSequenceRecords(accId, features); @@ -265,12 +265,13 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient proteinSeq.createDatasetSequence(); querySeq.createDatasetSequence(); - getProteinCrossReferences(proteinSeq); - MapList mapList = mapCdsToProtein(querySeq, proteinSeq); if (mapList != null) { - Mapping map = new Mapping(proteinSeq.getDatasetSequence(), mapList); + // clunky: ensure Uniprot xref if we have one is on mapped sequence + SequenceI ds = proteinSeq.getDatasetSequence(); + ds.setSourceDBRef(proteinSeq.getSourceDBRef()); + Mapping map = new Mapping(ds, mapList); DBRefEntry dbr = new DBRefEntry(getDbSource(), getDbVersion(), accId, map); querySeq.getDatasetSequence().addDBRef(dbr); @@ -290,29 +291,44 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient } /** - * Get Uniprot and PDB xrefs from Ensembl, and attach them to the protein - * sequence + * Get database xrefs from Ensembl, and attach them to the sequence * - * @param proteinSeq + * @param seq */ - protected void getProteinCrossReferences(SequenceI proteinSeq) + protected void getCrossReferences(SequenceI seq) { - while (proteinSeq.getDatasetSequence() != null) + while (seq.getDatasetSequence() != null) { - proteinSeq = proteinSeq.getDatasetSequence(); + seq = seq.getDatasetSequence(); } EnsemblXref xrefFetcher = new EnsemblXref(); - List xrefs = xrefFetcher.getCrossReferences( - proteinSeq.getName(), "PDB", "Uniprot/SPTREMBL", - "Uniprot/SWISSPROT"); + List xrefs = xrefFetcher.getCrossReferences(seq.getName(), + getCrossReferenceDatabases()); for (DBRefEntry xref : xrefs) { - proteinSeq.addDBRef(xref); + seq.addDBRef(xref); + /* + * Save any Uniprot xref to be the reference for SIFTS mapping + */ + if (DBRefSource.UNIPROT.equals(xref.getSource())) + { + seq.setSourceDBRef(xref); + } } } /** + * Returns a list of database names to be used when fetching cross-references. + * + * @return + */ + protected List getCrossReferenceDatabases() + { + return CROSS_REFERENCES; + } + + /** * Returns a mapping from dna to protein by inspecting sequence features of * type "CDS" on the dna. * @@ -327,7 +343,7 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient int mappedDnaLength = getCdsRanges(dnaSeq, ranges); int proteinLength = proteinSeq.getLength(); - List proteinRange = new ArrayList(); + int proteinEnd = proteinLength; int proteinStart = 1; /* @@ -339,26 +355,32 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient proteinStart = 2; proteinLength--; } - proteinRange.add(new int[] { proteinStart, proteinLength }); + List proteinRange = new ArrayList(); /* * dna length should map to protein (or protein plus stop codon) */ int codesForResidues = mappedDnaLength / 3; - if (codesForResidues == proteinLength - || codesForResidues == (proteinLength + 1)) + if (codesForResidues == (proteinLength + 1)) + { + MappingUtils.unmapStopCodon(ranges, mappedDnaLength); + codesForResidues--; + } + if (codesForResidues == proteinLength) { + proteinRange.add(new int[] { proteinStart, proteinEnd }); return new MapList(ranges, proteinRange, 3, 1); } return null; } /** - * Adds CDS ranges to the ranges list, and returns the total length mapped. + * Adds CDS ranges to the ranges list, and returns the total length mapped + * from. * - * No need to worry about reverse strand dna here since the retrieved sequence - * is as transcribed (reverse complement for reverse strand), i.e in the same - * sense as the peptide. + * No need to worry about reverse strand dna, here since the retrieved + * sequence is as transcribed (reverse complement for reverse strand), i.e in + * the same sense as the peptide. * * @param dnaSeq * @param ranges @@ -371,13 +393,14 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient { return 0; } + SequenceOntologyI so = SequenceOntologyFactory.getInstance(); int mappedDnaLength = 0; for (SequenceFeature sf : sfs) { /* * process a CDS feature (or a sub-type of CDS) */ - if (SequenceOntology.getInstance().isA(sf.getType(), SequenceOntology.CDS)) + if (so.isA(sf.getType(), SequenceOntologyI.CDS)) { int phase = 0; try { @@ -392,7 +415,7 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient */ int begin = sf.getBegin(); int end = sf.getEnd(); - if (ranges.isEmpty() && phase > 0) + if (ranges.isEmpty()) { begin += phase; if (begin > end) @@ -565,7 +588,7 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient * the start position of the sequence we are mapping to * @return */ - protected MapList getGenomicRanges(SequenceI sourceSequence, + protected MapList getGenomicRangesFromFeatures(SequenceI sourceSequence, String accId, int start) { SequenceFeature[] sfs = sourceSequence.getSequenceFeatures(); @@ -591,11 +614,12 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient */ if (identifiesSequence(sf, accId)) { - int strand = sf.getStrand(); - - if (directionSet && strand != direction) - { - // abort - mix of forward and backward + int strand = sf.getStrand(); + strand = strand == 0 ? 1 : strand; // treat unknown as forward + + if (directionSet && strand != direction) + { + // abort - mix of forward and backward System.err.println("Error: forward and backward strand for " + accId); return null; @@ -640,8 +664,8 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient */ Collections.sort(regions, new RangeSorter(direction == 1)); - List to = new ArrayList(); - to.add(new int[] { start, start + mappedLength - 1 }); + List to = Arrays.asList(new int[] { start, + start + mappedLength - 1 }); return new MapList(regions, to, 1, 1); } @@ -696,17 +720,18 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient /* * for sequence_variant, make an additional feature with consequence */ - if (SequenceOntology.getInstance().isSequenceVariant(sf.getType())) - { - String consequence = (String) sf.getValue(CONSEQUENCE_TYPE); - if (consequence != null) - { - SequenceFeature sf2 = new SequenceFeature("consequence", - consequence, copy.getBegin(), copy.getEnd(), 0f, - null); - targetSequence.addSequenceFeature(sf2); - } - } + // if (SequenceOntologyFactory.getInstance().isA(sf.getType(), + // SequenceOntologyI.SEQUENCE_VARIANT)) + // { + // String consequence = (String) sf.getValue(CONSEQUENCE_TYPE); + // if (consequence != null) + // { + // SequenceFeature sf2 = new SequenceFeature("consequence", + // consequence, copy.getBegin(), copy.getEnd(), 0f, + // null); + // targetSequence.addSequenceFeature(sf2); + // } + // } } } @@ -726,15 +751,22 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient return false; } + // long start = System.currentTimeMillis(); SequenceFeature[] sfs = sourceSequence.getSequenceFeatures(); - MapList mapping = getGenomicRanges(sourceSequence, accessionId, + MapList mapping = getGenomicRangesFromFeatures(sourceSequence, accessionId, targetSequence.getStart()); if (mapping == null) { return false; } - return transferFeatures(sfs, targetSequence, mapping, accessionId); + boolean result = transferFeatures(sfs, targetSequence, mapping, + accessionId); + // System.out.println("transferFeatures (" + (sfs.length) + " --> " + // + targetSequence.getSequenceFeatures().length + ") to " + // + targetSequence.getName() + // + " took " + (System.currentTimeMillis() - start) + "ms"); + return result; } /** @@ -836,7 +868,7 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient SequenceFeature[] sfs = sequence.getSequenceFeatures(); if (sfs != null) { - SequenceOntology so = SequenceOntology.getInstance(); + SequenceOntologyI so = SequenceOntologyFactory.getInstance(); for (SequenceFeature sf :sfs) { if (so.isA(sf.getType(), type)) { @@ -874,7 +906,7 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient } AlignmentUtils.transferFeatures(dnaSeq, peptide, dnaToProtein, - SequenceOntology.EXON); + SequenceOntologyI.EXON); LinkedHashMap variants = buildDnaVariantsMap( dnaSeq, dnaToProtein); @@ -895,12 +927,28 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient String desc = StringUtils.listToDelimitedString(peptideVariants, ", "); SequenceFeature sf = new SequenceFeature( - SequenceOntology.SEQUENCE_VARIANT, desc, peptidePos, + SequenceOntologyI.SEQUENCE_VARIANT, desc, peptidePos, peptidePos, 0f, null); peptide.addSequenceFeature(sf); count++; } } + + /* + * ugly sort to get sequence features in start position order + * - would be better to store in Sequence as a TreeSet instead? + */ + Arrays.sort(peptide.getSequenceFeatures(), + new Comparator() + { + @Override + public int compare(SequenceFeature o1, SequenceFeature o2) + { + int c = Integer.compare(o1.getBegin(), o2.getBegin()); + return c == 0 ? Integer.compare(o1.getEnd(), o2.getEnd()) + : c; + } + }); return count; } @@ -920,7 +968,7 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient * LinkedHashMap ensures we add the peptide features in sequence order */ LinkedHashMap variants = new LinkedHashMap(); - SequenceOntology so = SequenceOntology.getInstance(); + SequenceOntologyI so = SequenceOntologyFactory.getInstance(); SequenceFeature[] dnaFeatures = dnaSeq.getSequenceFeatures(); if (dnaFeatures == null) @@ -943,7 +991,7 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient // not handling multi-locus variant features continue; } - if (so.isSequenceVariant(sf.getType())) + if (so.isA(sf.getType(), SequenceOntologyI.SEQUENCE_VARIANT)) { int[] mapsTo = dnaToProtein.locateInTo(dnaCol, dnaCol); if (mapsTo == null) @@ -1081,7 +1129,8 @@ public abstract class EnsemblSeqProxy extends EnsemblRestClient */ public static boolean isTranscript(String featureType) { - return NMD_VARIANT.equals(featureType) - || SequenceOntology.getInstance().isA(featureType, SequenceOntology.TRANSCRIPT); + return SequenceOntologyI.NMD_TRANSCRIPT_VARIANT.equals(featureType) + || SequenceOntologyFactory.getInstance().isA(featureType, + SequenceOntologyI.TRANSCRIPT); } }