/* * Jalview - A Sequence Alignment Editor and Viewer ($$Version-Rel$$) * Copyright (C) $$Year-Rel$$ The Jalview Authors * * This file is part of Jalview. * * Jalview is free software: you can redistribute it and/or * modify it under the terms of the GNU General Public License * as published by the Free Software Foundation, either version 3 * of the License, or (at your option) any later version. * * Jalview is distributed in the hope that it will be useful, but * WITHOUT ANY WARRANTY; without even the implied warranty * of MERCHANTABILITY or FITNESS FOR A PARTICULAR * PURPOSE. See the GNU General Public License for more details. * * You should have received a copy of the GNU General Public License * along with Jalview. If not, see . * The Jalview Authors are detailed in the 'AUTHORS' file. */ package jalview.datamodel.xdb.embl; import jalview.analysis.SequenceIdMatcher; import jalview.datamodel.DBRefEntry; import jalview.datamodel.DBRefSource; import jalview.datamodel.FeatureProperties; import jalview.datamodel.Mapping; import jalview.datamodel.Sequence; import jalview.datamodel.SequenceFeature; import jalview.datamodel.SequenceI; import jalview.util.DBRefUtils; import jalview.util.MapList; import jalview.util.MappingUtils; import jalview.util.StringUtils; import java.util.Hashtable; import java.util.List; import java.util.Map; import java.util.Map.Entry; import java.util.Vector; import java.util.regex.Pattern; /** * Data model for one entry returned from an EMBL query, as marshalled by a * Castor binding file * * For example: * http://www.ebi.ac.uk/Tools/dbfetch/dbfetch?db=ena_sequence&id=J03321 * &format=emblxml * * @see embl_mapping.xml */ public class EmblEntry { private static final Pattern SPACE_PATTERN = Pattern.compile(" "); String accession; String version; String taxDivision; String desc; String rCreated; String rLastUpdated; String lastUpdated; Vector keywords; Vector dbRefs; Vector features; EmblSequence sequence; /** * @return the accession */ public String getAccession() { return accession; } /** * @param accession * the accession to set */ public void setAccession(String accession) { this.accession = accession; } /** * @return the dbRefs */ public Vector getDbRefs() { return dbRefs; } /** * @param dbRefs * the dbRefs to set */ public void setDbRefs(Vector dbRefs) { this.dbRefs = dbRefs; } /** * @return the desc */ public String getDesc() { return desc; } /** * @param desc * the desc to set */ public void setDesc(String desc) { this.desc = desc; } /** * @return the features */ public Vector getFeatures() { return features; } /** * @param features * the features to set */ public void setFeatures(Vector features) { this.features = features; } /** * @return the keywords */ public Vector getKeywords() { return keywords; } /** * @param keywords * the keywords to set */ public void setKeywords(Vector keywords) { this.keywords = keywords; } /** * @return the lastUpdated */ public String getLastUpdated() { return lastUpdated; } /** * @param lastUpdated * the lastUpdated to set */ public void setLastUpdated(String lastUpdated) { this.lastUpdated = lastUpdated; } /** * @return the releaseCreated */ public String getRCreated() { return rCreated; } /** * @param releaseCreated * the releaseCreated to set */ public void setRCreated(String releaseCreated) { this.rCreated = releaseCreated; } /** * @return the releaseLastUpdated */ public String getRLastUpdated() { return rLastUpdated; } /** * @param releaseLastUpdated * the releaseLastUpdated to set */ public void setRLastUpdated(String releaseLastUpdated) { this.rLastUpdated = releaseLastUpdated; } /** * @return the sequence */ public EmblSequence getSequence() { return sequence; } /** * @param sequence * the sequence to set */ public void setSequence(EmblSequence sequence) { this.sequence = sequence; } /** * @return the taxDivision */ public String getTaxDivision() { return taxDivision; } /** * @param taxDivision * the taxDivision to set */ public void setTaxDivision(String taxDivision) { this.taxDivision = taxDivision; } /** * @return the version */ public String getVersion() { return version; } /** * @param version * the version to set */ public void setVersion(String version) { this.version = version; } /** * Recover annotated sequences from EMBL file * * @param sourceDb * @param peptides * a list of protein products found so far (to add to) * @return dna dataset sequence with DBRefs and features */ public SequenceI getSequence(String sourceDb, List peptides) { SequenceI dna = new Sequence(sourceDb + "|" + accession, sequence.getSequence()); dna.setDescription(desc); DBRefEntry retrievedref = new DBRefEntry(sourceDb, version, accession); dna.addDBRef(retrievedref); // add map to indicate the sequence is a valid coordinate frame for the // dbref retrievedref.setMap(new Mapping(null, new int[] { 1, dna.getLength() }, new int[] { 1, dna.getLength() }, 1, 1)); // TODO: transform EMBL Database refs to canonical form if (dbRefs != null) { for (DBRefEntry dbref : dbRefs) { dna.addDBRef(dbref); } } try { for (EmblFeature feature : features) { if (feature.dbRefs != null) { for (DBRefEntry dbref : feature.dbRefs) { dna.addDBRef(dbref); } } if (FeatureProperties.isCodingFeature(sourceDb, feature.getName())) { parseCodingFeature(feature, sourceDb, dna, peptides); } } } catch (Exception e) { System.err.println("EMBL Record Features parsing error!"); System.err .println("Please report the following to help@jalview.org :"); System.err.println("EMBL Record " + accession); System.err.println("Resulted in exception: " + e.getMessage()); e.printStackTrace(System.err); } return dna; } /** * Extracts coding region and product from a CDS feature and properly decorate * it with annotations. * * @param feature * coding feature * @param sourceDb * source database for the EMBLXML * @param dna * parent dna sequence for this record * @param peptides * list of protein product sequences for Embl entry */ void parseCodingFeature(EmblFeature feature, String sourceDb, SequenceI dna, List peptides) { boolean isEmblCdna = sourceDb.equals(DBRefSource.EMBLCDS); int[] exon = getCdsRanges(feature); String prseq = null; String prname = ""; String prid = null; Map vals = new Hashtable(); SequenceIdMatcher matcher = new SequenceIdMatcher(peptides); /* * codon_start 1/2/3 in EMBL corresponds to phase 0/1/2 in CDS * (phase is required for CDS features in GFF3 format) */ int codonStart = 1; /* * parse qualifiers, saving protein translation, protein id, * codon start position, product (name), and 'other values' */ if (feature.getQualifiers() != null) { for (Qualifier q : feature.getQualifiers()) { String qname = q.getName(); if (qname.equals("translation")) { // remove all spaces (precompiled String.replaceAll(" ", "")) prseq = SPACE_PATTERN.matcher(q.getValues()[0]).replaceAll(""); } else if (qname.equals("protein_id")) { prid = q.getValues()[0]; } else if (qname.equals("codon_start")) { try { codonStart = Integer.parseInt(q.getValues()[0]); } catch (NumberFormatException e) { System.err.println("Invalid codon_start in XML for " + accession + ": " + e.getMessage()); } } else if (qname.equals("product")) { // sometimes name is returned e.g. for V00488 prname = q.getValues()[0]; } else { // throw anything else into the additional properties hash String[] qvals = q.getValues(); if (qvals != null) { String commaSeparated = StringUtils.arrayToSeparatorList(qvals, ","); vals.put(qname, commaSeparated); } } } } // SequenceI product = null; DBRefEntry protEMBLCDS = null; exon = MappingUtils.removeStartPositions(codonStart - 1, exon); boolean noProteinDbref = true; SequenceI product = null; Mapping map = null; if (prseq != null && prname != null && prid != null) { /* * look for product in peptides list, if not found, add it */ product = matcher.findIdMatch(prid); if (product == null) { product = new Sequence(prid, prseq, 1, prseq.length()); product.setDescription(((prname.length() == 0) ? "Protein Product from " + sourceDb : prname)); peptides.add(product); matcher.add(product); } // we have everything - create the mapping and perhaps the protein // sequence if (exon == null || exon.length == 0) { System.err .println("Implementation Notice: EMBLCDS records not properly supported yet - Making up the CDNA region of this sequence... may be incorrect (" + sourceDb + ":" + getAccession() + ")"); if (prseq.length() * 3 == (1 - codonStart + dna.getSequence().length)) { System.err .println("Not allowing for additional stop codon at end of cDNA fragment... !"); // this might occur for CDS sequences where no features are // marked. exon = new int[] { dna.getStart() + (codonStart - 1), dna.getEnd() }; map = new Mapping(product, exon, new int[] { 1, prseq.length() }, 3, 1); } if ((prseq.length() + 1) * 3 == (1 - codonStart + dna.getSequence().length)) { System.err .println("Allowing for additional stop codon at end of cDNA fragment... will probably cause an error in VAMSAs!"); exon = new int[] { dna.getStart() + (codonStart - 1), dna.getEnd() - 3 }; map = new Mapping(product, exon, new int[] { 1, prseq.length() }, 3, 1); } } else { // Trim the exon mapping if necessary - the given product may only be a // fragment of a larger protein. (EMBL:AY043181 is an example) if (isEmblCdna) { // TODO: Add a DbRef back to the parent EMBL sequence with the exon // map // if given a dataset reference, search dataset for parent EMBL // sequence if it exists and set its map // make a new feature annotating the coding contig } else { // final product length truncation check // TODO should from range include stop codon even if not in protein // in order to include stop codon in CDS sequence (as done for // Ensembl)? int[] cdsRanges = adjustForProteinLength(prseq.length(), exon); map = new Mapping(product, cdsRanges, new int[] { 1, prseq.length() }, 3, 1); // reconstruct the EMBLCDS entry // TODO: this is only necessary when there codon annotation is // complete (I think JBPNote) DBRefEntry pcdnaref = new DBRefEntry(); pcdnaref.setAccessionId(prid); pcdnaref.setSource(DBRefSource.EMBLCDS); pcdnaref.setVersion(getVersion()); // same as parent EMBL version. MapList mp = new MapList(new int[] { 1, prseq.length() }, new int[] { 1 + (codonStart - 1), (codonStart - 1) + 3 * prseq.length() }, 1, 3); pcdnaref.setMap(new Mapping(mp)); if (product != null) { product.addDBRef(pcdnaref); protEMBLCDS = new DBRefEntry(pcdnaref); protEMBLCDS.setSource(DBRefSource.EMBLCDSProduct); product.addDBRef(protEMBLCDS); } } } // add cds feature to dna seq - this may include the stop codon for (int xint = 0; exon != null && xint < exon.length; xint += 2) { SequenceFeature sf = makeCdsFeature(exon, xint, prname, prid, vals, codonStart); sf.setType(feature.getName()); // "CDS" sf.setFeatureGroup(sourceDb); dna.addSequenceFeature(sf); } } /* * add dbRefs to sequence, and mappings for Uniprot xrefs */ if (feature.dbRefs != null) { boolean mappingUsed = false; for (DBRefEntry ref : feature.dbRefs) { ref.setSource(DBRefUtils.getCanonicalName(ref.getSource())); if (ref.getSource().equals(DBRefSource.UNIPROT)) { String proteinSeqName = DBRefSource.UNIPROT + "|" + ref.getAccessionId(); if (map != null && map.getTo() != null) { if (mappingUsed) { /* * two or more Uniprot xrefs for the same CDS - * each needs a distinct Mapping (as to a different sequence) */ map = new Mapping(map); } mappingUsed = true; /* * try to locate the protein mapped to (possibly by a * previous CDS feature) */ SequenceI proteinSeq = matcher.findIdMatch(proteinSeqName); if (proteinSeq == null) { proteinSeq = new Sequence(proteinSeqName, product .getSequenceAsString()); matcher.add(proteinSeq); peptides.add(proteinSeq); } map.setTo(proteinSeq); map.getTo().addDBRef( new DBRefEntry(ref.getSource(), ref.getVersion(), ref .getAccessionId())); ref.setMap(map); } noProteinDbref = false; } if (product != null) { DBRefEntry pref = new DBRefEntry(ref.getSource(), ref.getVersion(), ref.getAccessionId()); pref.setMap(null); // reference is direct product.addDBRef(pref); // Add converse mapping reference if (map != null) { Mapping pmap = new Mapping(dna, map.getMap().getInverse()); pref = new DBRefEntry(sourceDb, getVersion(), this.getAccession()); pref.setMap(pmap); if (map.getTo() != null) { map.getTo().addDBRef(pref); } } } dna.addDBRef(ref); } if (noProteinDbref && product != null) { // add protein coding reference to dna sequence so xref matches if (protEMBLCDS == null) { protEMBLCDS = new DBRefEntry(); protEMBLCDS.setAccessionId(prid); protEMBLCDS.setSource(DBRefSource.EMBLCDSProduct); protEMBLCDS.setVersion(getVersion()); protEMBLCDS .setMap(new Mapping(product, map.getMap().getInverse())); } product.addDBRef(protEMBLCDS); // Add converse mapping reference if (map != null) { Mapping pmap = new Mapping(product, protEMBLCDS.getMap().getMap() .getInverse()); DBRefEntry ncMap = new DBRefEntry(protEMBLCDS); ncMap.setMap(pmap); if (map.getTo() != null) { dna.addDBRef(ncMap); } } } } } /** * Helper method to construct a SequenceFeature for one cds range * * @param exons * array of cds [start, end, ...] positions * @param exonStartIndex * offset into the exons array * @param proteinName * @param proteinAccessionId * @param vals * map of 'miscellaneous values' for feature * @param codonStart * codon start position for CDS (1/2/3, normally 1) * @return */ protected SequenceFeature makeCdsFeature(int[] exons, int exonStartIndex, String proteinName, String proteinAccessionId, Map vals, int codonStart) { int exonNumber = exonStartIndex / 2 + 1; SequenceFeature sf = new SequenceFeature(); sf.setBegin(Math.min(exons[exonStartIndex], exons[exonStartIndex + 1])); sf.setEnd(Math.max(exons[exonStartIndex], exons[exonStartIndex + 1])); sf.setDescription(String.format( "Exon %d for protein '%s' EMBLCDS:%s", exonNumber, proteinName, proteinAccessionId)); sf.setPhase(String.valueOf(codonStart - 1)); sf.setStrand(exons[exonStartIndex] <= exons[exonStartIndex + 1] ? "+" : "-"); sf.setValue(FeatureProperties.EXONPOS, exonNumber); sf.setValue(FeatureProperties.EXONPRODUCT, proteinName); if (!vals.isEmpty()) { StringBuilder sb = new StringBuilder(); boolean first = true; for (Entry val : vals.entrySet()) { if (!first) { sb.append(";"); } sb.append(val.getKey()).append("=").append(val.getValue()); first = false; sf.setValue(val.getKey(), val.getValue()); } sf.setAttributes(sb.toString()); } return sf; } /** * Returns the CDS positions as a list of [start, end, start, end...] * positions. If on the reverse strand, these will be in descending order. * * @param feature * @return */ protected int[] getCdsRanges(EmblFeature feature) { if (feature.locations == null) { return new int[] {}; } int cdsBoundaryCount = 0; // count of all start/stop locations int[][] cdsLocations = new int[feature.locations.size()][]; int locationNumber = 0; for (EmblFeatureLocations loc : feature.locations) { int[] locationRanges = loc.getElementRanges(accession); cdsLocations[locationNumber++] = locationRanges; cdsBoundaryCount += locationRanges.length; } int[] cdsRanges = new int[cdsBoundaryCount]; int copyTo = 0; for (int[] ranges : cdsLocations) { System.arraycopy(ranges, 0, cdsRanges, copyTo, ranges.length); copyTo += ranges.length; } return cdsRanges; } /** * truncate the last exon interval to the prlength'th codon * * @param prlength * @param exon * @return new exon */ private int[] adjustForProteinLength(int prlength, int[] exon) { int origxon[], sxpos = -1, endxon = 0, cdslength = prlength * 3; // first adjust range for codon start attribute if (prlength >= 1 && exon != null) { origxon = new int[exon.length]; System.arraycopy(exon, 0, origxon, 0, exon.length); int cdspos = 0; for (int x = 0; x < exon.length && sxpos == -1; x += 2) { cdspos += Math.abs(exon[x + 1] - exon[x]) + 1; if (cdslength <= cdspos) { // advanced beyond last codon. sxpos = x; if (cdslength != cdspos) { System.err .println("Truncating final exon interval on region by " + (cdspos - cdslength)); } // locate the new end boundary of final exon as endxon endxon = exon[x + 1] - cdspos + cdslength; break; } } if (sxpos != -1) { // and trim the exon interval set if necessary int[] nxon = new int[sxpos + 2]; System.arraycopy(exon, 0, nxon, 0, sxpos + 2); nxon[sxpos + 1] = endxon; // update the end boundary for the new exon // set exon = nxon; } } return exon; } }