X-Git-Url: http://source.jalview.org/gitweb/?a=blobdiff_plain;f=src%2Fjalview%2Fio%2FFlatFile.java;fp=src%2Fjalview%2Fio%2FFlatFile.java;h=0000000000000000000000000000000000000000;hb=2ed4d627d3eb73a49299aad68bf719881ec4538d;hp=9e5c65224b4d8d4490348325b6ea4c742962b3a5;hpb=9ae61a06b581ffba65614a2f08acc36e7f7685e7;p=jalview.git
diff --git a/src/jalview/io/FlatFile.java b/src/jalview/io/FlatFile.java
deleted file mode 100644
index 9e5c652..0000000
--- a/src/jalview/io/FlatFile.java
+++ /dev/null
@@ -1,833 +0,0 @@
-package jalview.io;
-
-import java.io.IOException;
-import java.text.ParseException;
-import java.util.ArrayList;
-import java.util.Arrays;
-import java.util.HashMap;
-import java.util.Hashtable;
-import java.util.List;
-import java.util.Locale;
-import java.util.Map;
-import java.util.Map.Entry;
-import java.util.TreeMap;
-
-import jalview.bin.Cache;
-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.DnaUtils;
-import jalview.util.MapList;
-import jalview.util.MappingUtils;
-
-/**
- * A base class to support parsing of GenBank, EMBL or DDBJ flat file format
- * data. Example files (rather than formal specifications) are provided at
- *
- *
- * https://ena-docs.readthedocs.io/en/latest/submit/fileprep/flat-file-example.html
- * https://www.ncbi.nlm.nih.gov/Sitemap/samplerecord.html
- *
- *
- * or to compare the same entry, see
- *
- *
- * https://www.ebi.ac.uk/ena/browser/api/embl/X81322.1
- * https://www.ncbi.nlm.nih.gov/nuccore/X81322.1
- *
- *
- * The feature table part of the file has a common definition, only the start of
- * each line is formatted differently in GenBank and EMBL. See
- * http://www.insdc.org/files/feature_table.html#7.1.
- */
-public abstract class FlatFile extends AlignFile
-{
- protected static final String LOCATION = "location";
-
- protected static final String QUOTE = "\"";
-
- protected static final String DOUBLED_QUOTE = QUOTE + QUOTE;
-
- protected static final String WHITESPACE = "\\s+";
-
- /**
- * Removes leading or trailing double quotes (") unless doubled, and changes
- * any 'escaped' (doubled) double quotes to single characters. As per the
- * Feature Table specification for Qualifiers, Free Text.
- *
- * @param value
- * @return
- */
- protected static String removeQuotes(String value)
- {
- if (value == null)
- {
- return null;
- }
- if (value.startsWith(QUOTE) && !value.startsWith(DOUBLED_QUOTE))
- {
- value = value.substring(1);
- }
- if (value.endsWith(QUOTE) && !value.endsWith(DOUBLED_QUOTE))
- {
- value = value.substring(0, value.length() - 1);
- }
- value = value.replace(DOUBLED_QUOTE, QUOTE);
- return value;
- }
-
- /**
- * Truncates (if necessary) the exon intervals to match 3 times the length of
- * the protein; also accepts 3 bases longer (for stop codon not included in
- * protein)
- *
- * @param proteinLength
- * @param exon
- * an array of [start, end, start, end...] intervals
- * @return the same array (if unchanged) or a truncated copy
- */
- protected static int[] adjustForProteinLength(int proteinLength,
- int[] exon)
- {
- if (proteinLength <= 0 || exon == null)
- {
- return exon;
- }
- int expectedCdsLength = proteinLength * 3;
- int exonLength = MappingUtils.getLength(Arrays.asList(exon));
-
- /*
- * if exon length matches protein, or is shorter, or longer by the
- * length of a stop codon (3 bases), then leave it unchanged
- */
- if (expectedCdsLength >= exonLength
- || expectedCdsLength == exonLength - 3)
- {
- return exon;
- }
-
- int origxon[];
- int sxpos = -1;
- int endxon = 0;
- origxon = new int[exon.length];
- System.arraycopy(exon, 0, origxon, 0, exon.length);
- int cdspos = 0;
- for (int x = 0; x < exon.length; x += 2)
- {
- cdspos += Math.abs(exon[x + 1] - exon[x]) + 1;
- if (expectedCdsLength <= cdspos)
- {
- // advanced beyond last codon.
- sxpos = x;
- if (expectedCdsLength != cdspos)
- {
- // System.err
- // .println("Truncating final exon interval on region by "
- // + (cdspos - cdslength));
- }
-
- /*
- * shrink the final exon - reduce end position if forward
- * strand, increase it if reverse
- */
- if (exon[x + 1] >= exon[x])
- {
- endxon = exon[x + 1] - cdspos + expectedCdsLength;
- }
- else
- {
- endxon = exon[x + 1] + cdspos - expectedCdsLength;
- }
- 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;
- }
-
- /*
- * when true, interpret the mol_type 'source' feature attribute
- * and generate an RNA sequence from the DNA record
- */
- protected boolean produceRna=true;
-
-
- /*
- * values parsed from the data file
- */
- protected String sourceDb;
-
- protected String accession;
-
- protected String version;
-
- protected String description;
-
- protected int length = 128;
-
- protected List dbrefs;
-
- protected boolean sequenceStringIsRNA=false;
-
- protected String sequenceString;
-
- protected Map cds;
-
- /**
- * Constructor
- *
- * @param fp
- * @param sourceId
- * @throws IOException
- */
- public FlatFile(FileParse fp, String sourceId) throws IOException
- {
- super(false, fp); // don't parse immediately
- this.sourceDb = sourceId;
- dbrefs = new ArrayList<>();
-
- /*
- * using TreeMap gives CDS sequences in alphabetical, so readable, order
- */
- cds = new TreeMap<>(String.CASE_INSENSITIVE_ORDER);
-
- parse();
- }
-
- /**
- * process attributes for 'source' until the next FT feature entry
- * only interested in 'mol_type'
- * @param tokens
- * @return
- * @throws IOException
- */
- private String parseSourceQualifiers(String[] tokens) throws IOException
- {
- if (!"source".equals(tokens[0]))
- {
- throw (new RuntimeException("Not given a 'source' qualifier line"));
- }
- // search for mol_type attribute
-
- StringBuilder sb = new StringBuilder().append(tokens[1]); // extent of
- // sequence
-
- String line = parseFeatureQualifier(sb, false);
- while (line != null)
- {
- if (!line.startsWith("FT ")) // four spaces, end of this feature table
- // entry
- {
- return line;
- }
-
- // case sensitive ?
- int p = line.indexOf("\\mol_type");
- int qs = line.indexOf("\"", p);
- int qe = line.indexOf("\"", qs + 1);
- String qualifier=line.substring(qs,qe).toLowerCase(Locale.ROOT);
- if (qualifier.indexOf("rna") > -1)
- {
- sequenceStringIsRNA = true;
- }
- if (qualifier.indexOf("dna") > -1)
- {
- sequenceStringIsRNA = false;
- }
- line=parseFeatureQualifier(sb, false);
- }
- return line;
- }
-
-
- /**
- * Parses one (GenBank or EMBL format) CDS feature, saves the parsed data, and
- * returns the next line
- *
- * @param location
- * @return
- * @throws IOException
- */
- protected String parseCDSFeature(String location) throws IOException
- {
- String line;
-
- /*
- * parse location, which can be over >1 line e.g. EAW51554
- */
- CdsData data = new CdsData();
- StringBuilder sb = new StringBuilder().append(location);
- line = parseFeatureQualifier(sb, false);
- data.cdsLocation = sb.toString();
-
- while (line != null)
- {
- if (!isFeatureContinuationLine(line))
- {
- // e.g. start of next feature "FT source..."
- break;
- }
-
- /*
- * extract qualifier, e.g. FT /protein_id="CAA37824.1"
- * - the value may extend over more than one line
- * - if the value has enclosing quotes, these are removed
- * - escaped double quotes ("") are reduced to a single character
- */
- int slashPos = line.indexOf('/');
- if (slashPos == -1)
- {
- Cache.log.error("Unexpected EMBL line ignored: " + line);
- line = nextLine();
- continue;
- }
- int eqPos = line.indexOf('=', slashPos + 1);
- if (eqPos == -1)
- {
- // can happen, e.g. /ribosomal_slippage
- line = nextLine();
- continue;
- }
- String qualifier = line.substring(slashPos + 1, eqPos);
- String value = line.substring(eqPos + 1);
- value = removeQuotes(value);
- sb = new StringBuilder().append(value);
- boolean asText = !"translation".equals(qualifier);
- line = parseFeatureQualifier(sb, asText);
- String featureValue = sb.toString();
-
- if ("protein_id".equals(qualifier))
- {
- data.proteinId = featureValue;
- }
- else if ("codon_start".equals(qualifier))
- {
- try
- {
- data.codonStart = Integer.parseInt(featureValue.trim());
- } catch (NumberFormatException e)
- {
- Cache.log.error("Invalid codon_start in XML for " + this.accession
- + ": " + e.getMessage());
- }
- }
- else if ("db_xref".equals(qualifier))
- {
- String[] parts = featureValue.split(":");
- if (parts.length == 2)
- {
- String db = parts[0].trim();
- db = DBRefUtils.getCanonicalName(db);
- DBRefEntry dbref = new DBRefEntry(db, "0", parts[1].trim());
- data.xrefs.add(dbref);
- }
- }
- else if ("product".equals(qualifier))
- {
- data.proteinName = featureValue;
- }
- else if ("translation".equals(qualifier))
- {
- data.translation = featureValue;
- }
- else if (!"".equals(featureValue))
- {
- // throw anything else into the additional properties hash
- data.cdsProps.put(qualifier, featureValue);
- }
- }
-
- if (data.proteinId != null)
- {
- this.cds.put(data.proteinId, data);
- }
- else
- {
- Cache.log.error("Ignoring CDS feature with no protein_id for "
- + sourceDb + ":" + accession);
- }
-
- return line;
- }
-
- protected abstract boolean isFeatureContinuationLine(String line);
-
- /**
- * Output (print) is not (yet) implemented for flat file format
- */
- @Override
- public String print(SequenceI[] seqs, boolean jvsuffix)
- {
- return null;
- }
-
- /**
- * Constructs and saves the sequence from parsed components
- */
- protected void buildSequence()
- {
- if (this.accession == null || this.sequenceString == null)
- {
- Cache.log.error("Failed to parse data from EMBL");
- return;
- }
-
- String name = this.accession;
- if (this.sourceDb != null)
- {
- name = this.sourceDb + "|" + name;
- }
-
- if (produceRna && sequenceStringIsRNA)
- {
- sequenceString = sequenceString.replace('T', 'U').replace('t', 'u');
- }
-
- SequenceI seq = new Sequence(name, this.sequenceString);
- seq.setDescription(this.description);
-
- /*
- * add a DBRef to itself
- */
- DBRefEntry selfRef = new DBRefEntry(sourceDb, version, accession);
- int[] startEnd = new int[] { 1, seq.getLength() };
- selfRef.setMap(new Mapping(null, startEnd, startEnd, 1, 1));
- seq.addDBRef(selfRef);
-
- for (DBRefEntry dbref : this.dbrefs)
- {
- seq.addDBRef(dbref);
- }
-
- processCDSFeatures(seq);
-
- seq.deriveSequence();
-
- addSequence(seq);
- }
-
- /**
- * Process the CDS features, including generation of cross-references and
- * mappings to the protein products (translation)
- *
- * @param seq
- */
- protected void processCDSFeatures(SequenceI seq)
- {
- /*
- * record protein products found to avoid duplication i.e. >1 CDS with
- * the same /protein_id [though not sure I can find an example of this]
- */
- Map proteins = new HashMap<>();
- for (CdsData data : cds.values())
- {
- processCDSFeature(seq, data, proteins);
- }
- }
-
- /**
- * Processes data for one parsed CDS feature to
- *
- * - create a protein product sequence for the translation
- * - create a cross-reference to protein with mapping from dna
- * - add a CDS feature to the sequence for each CDS start-end range
- * - add any CDS dbrefs to the sequence and to the protein product
- *
- *
- * @param SequenceI
- * dna
- * @param proteins
- * map of protein products so far derived from CDS data
- */
- void processCDSFeature(SequenceI dna, CdsData data,
- Map proteins)
- {
- /*
- * parse location into a list of [start, end, start, end] positions
- */
- int[] exons = getCdsRanges(this.accession, data.cdsLocation);
-
- MapList maplist = buildMappingToProtein(dna, exons, data);
-
- int exonNumber = 0;
-
- for (int xint = 0; exons != null && xint < exons.length - 1; xint += 2)
- {
- int exonStart = exons[xint];
- int exonEnd = exons[xint + 1];
- int begin = Math.min(exonStart, exonEnd);
- int end = Math.max(exonStart, exonEnd);
- exonNumber++;
- String desc = String.format("Exon %d for protein EMBLCDS:%s",
- exonNumber, data.proteinId);
-
- SequenceFeature sf = new SequenceFeature("CDS", desc, begin, end,
- this.sourceDb);
- for (Entry val : data.cdsProps.entrySet())
- {
- sf.setValue(val.getKey(), val.getValue());
- }
-
- sf.setEnaLocation(data.cdsLocation);
- boolean forwardStrand = exonStart <= exonEnd;
- sf.setStrand(forwardStrand ? "+" : "-");
- sf.setPhase(String.valueOf(data.codonStart - 1));
- sf.setValue(FeatureProperties.EXONPOS, exonNumber);
- sf.setValue(FeatureProperties.EXONPRODUCT, data.proteinName);
-
- dna.addSequenceFeature(sf);
- }
-
- boolean hasUniprotDbref = false;
- for (DBRefEntry xref : data.xrefs)
- {
- dna.addDBRef(xref);
- if (xref.getSource().equals(DBRefSource.UNIPROT))
- {
- /*
- * construct (or find) the sequence for (data.protein_id, data.translation)
- */
- SequenceI protein = buildProteinProduct(dna, xref, data, proteins);
- Mapping map = new Mapping(protein, maplist);
- map.setMappedFromId(data.proteinId);
- xref.setMap(map);
-
- /*
- * add DBRefs with mappings from dna to protein and the inverse
- */
- DBRefEntry db1 = new DBRefEntry(sourceDb, version, accession);
- db1.setMap(new Mapping(dna, maplist.getInverse()));
- protein.addDBRef(db1);
-
- hasUniprotDbref = true;
- }
- }
-
- /*
- * if we have a product (translation) but no explicit Uniprot dbref
- * (example: EMBL M19487 protein_id AAB02592.1)
- * then construct mappings to an assumed EMBLCDSPROTEIN accession
- */
- if (!hasUniprotDbref)
- {
- SequenceI protein = proteins.get(data.proteinId);
- if (protein == null)
- {
- protein = new Sequence(data.proteinId, data.translation);
- protein.setDescription(data.proteinName);
- proteins.put(data.proteinId, protein);
- }
- // assuming CDSPROTEIN sequence version = dna version (?!)
- DBRefEntry db1 = new DBRefEntry(DBRefSource.EMBLCDSProduct,
- this.version, data.proteinId);
- protein.addDBRef(db1);
-
- DBRefEntry dnaToEmblProteinRef = new DBRefEntry(
- DBRefSource.EMBLCDSProduct, this.version, data.proteinId);
- Mapping map = new Mapping(protein, maplist);
- map.setMappedFromId(data.proteinId);
- dnaToEmblProteinRef.setMap(map);
- dna.addDBRef(dnaToEmblProteinRef);
- }
-
- /*
- * comment brought forward from EmblXmlSource, lines 447-451:
- * TODO: if retrieved from EMBLCDS, 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
- */
- }
-
- /**
- * Computes a mapping from CDS positions in DNA sequence to protein product
- * positions, with allowance for stop codon or incomplete start codon
- *
- * @param dna
- * @param exons
- * @param data
- * @return
- */
- MapList buildMappingToProtein(final SequenceI dna, final int[] exons,
- final CdsData data)
- {
- MapList dnaToProteinMapping = null;
- int peptideLength = data.translation.length();
-
- int[] proteinRange = new int[] { 1, peptideLength };
- if (exons != null && exons.length > 0)
- {
- /*
- * We were able to parse 'location'; do a final
- * product length truncation check
- */
- int[] cdsRanges = adjustForProteinLength(peptideLength, exons);
- dnaToProteinMapping = new MapList(cdsRanges, proteinRange, 3, 1);
- }
- else
- {
- /*
- * workaround until we handle all 'location' formats fully
- * e.g. X53828.1:60..1058 or <123..>289
- */
- Cache.log.error(String.format(
- "Implementation Notice: EMBLCDS location '%s'not properly supported yet"
- + " - Making up the CDNA region of (%s:%s)... may be incorrect",
- data.cdsLocation, sourceDb, this.accession));
-
- int completeCodonsLength = 1 - data.codonStart + dna.getLength();
- int mappedDnaEnd = dna.getEnd();
- if (peptideLength * 3 == completeCodonsLength)
- {
- // this might occur for CDS sequences where no features are marked
- Cache.log.warn("Assuming no stop codon at end of cDNA fragment");
- mappedDnaEnd = dna.getEnd();
- }
- else if ((peptideLength + 1) * 3 == completeCodonsLength)
- {
- Cache.log.warn("Assuming stop codon at end of cDNA fragment");
- mappedDnaEnd = dna.getEnd() - 3;
- }
-
- if (mappedDnaEnd != -1)
- {
- int[] cdsRanges = new int[] {
- dna.getStart() + (data.codonStart - 1), mappedDnaEnd };
- dnaToProteinMapping = new MapList(cdsRanges, proteinRange, 3, 1);
- }
- }
-
- return dnaToProteinMapping;
- }
-
- /**
- * Constructs a sequence for the protein product for the CDS data (if there is
- * one), and dbrefs with mappings from CDS to protein and the reverse
- *
- * @param dna
- * @param xref
- * @param data
- * @param proteins
- * @return
- */
- SequenceI buildProteinProduct(SequenceI dna, DBRefEntry xref,
- CdsData data, Map proteins)
- {
- /*
- * check we have some data to work with
- */
- if (data.proteinId == null || data.translation == null)
- {
- return null;
- }
-
- /*
- * Construct the protein sequence (if not already seen)
- */
- String proteinSeqName = xref.getSource() + "|" + xref.getAccessionId();
- SequenceI protein = proteins.get(proteinSeqName);
- if (protein == null)
- {
- protein = new Sequence(proteinSeqName, data.translation, 1,
- data.translation.length());
- protein.setDescription(data.proteinName != null ? data.proteinName
- : "Protein Product from " + sourceDb);
- proteins.put(proteinSeqName, protein);
- }
-
- return protein;
- }
-
- /**
- * Returns the CDS location as a single array of [start, end, start, end...]
- * positions. If on the reverse strand, these will be in descending order.
- *
- * @param accession
- * @param location
- * @return
- */
- protected int[] getCdsRanges(String accession, String location)
- {
- if (location == null)
- {
- return new int[] {};
- }
-
- try
- {
- List ranges = DnaUtils.parseLocation(location);
- return MappingUtils.rangeListToArray(ranges);
- } catch (ParseException e)
- {
- Cache.log.warn(
- String.format("Not parsing inexact CDS location %s in ENA %s",
- location, accession));
- return new int[] {};
- }
- }
-
- /**
- * Reads the value of a feature (FT) qualifier from one or more lines of the
- * file, and returns the next line after that. Values are appended to the
- * string buffer, which should be already primed with the value read from the
- * first line for the qualifier (with any leading double quote removed).
- * Enclosing double quotes are removed, and escaped (repeated) double quotes
- * reduced to one only. For example for
- *
- *
- * FT /note="gene_id=hCG28070.3
- * FT ""foobar"" isoform=CRA_b"
- * the returned value is
- * gene_id=hCG28070.3 "foobar" isoform=CRA_b
- *
- *
- * Note the side-effect of this method, to advance data reading to the next
- * line after the feature qualifier (which could be another qualifier, a
- * different feature, a non-feature line, or null at end of file).
- *
- * @param sb
- * a string buffer primed with the first line of the value
- * @param asText
- * @return
- * @throws IOException
- */
- String parseFeatureQualifier(StringBuilder sb, boolean asText)
- throws IOException
- {
- String line;
- while ((line = nextLine()) != null)
- {
- if (!isFeatureContinuationLine(line))
- {
- break; // reached next feature or other input line
- }
- String[] tokens = line.split(WHITESPACE);
- if (tokens.length < 2)
- {
- Cache.log.error("Ignoring bad EMBL line for " + this.accession
- + ": " + line);
- break;
- }
- if (tokens[1].startsWith("/"))
- {
- break; // next feature qualifier
- }
-
- /*
- * if text (e.g. /product), add a word separator for a new line,
- * else (e.g. /translation) don't
- */
- if (asText)
- {
- sb.append(" ");
- }
-
- /*
- * remove trailing " and unescape doubled ""
- */
- String data = removeQuotes(tokens[1]);
- sb.append(data);
- }
-
- return line;
- }
-
- /**
- * Reads and saves the sequence, read from the lines following the ORIGIN
- * (GenBank) or SQ (EMBL) line. Whitespace and position counters are
- * discarded. Returns the next line following the sequence data (the next line
- * that doesn't start with whitespace).
- *
- * @throws IOException
- */
- protected String parseSequence() throws IOException
- {
- StringBuilder sb = new StringBuilder(this.length);
- String line = nextLine();
- while (line != null && line.startsWith(" "))
- {
- line = line.trim();
- String[] blocks = line.split(WHITESPACE);
-
- /*
- * the first or last block on each line might be a position count - omit
- */
- for (int i = 0; i < blocks.length; i++)
- {
- try
- {
- Long.parseLong(blocks[i]);
- // position counter - ignore it
- } catch (NumberFormatException e)
- {
- // sequence data - append it
- sb.append(blocks[i]);
- }
- }
- line = nextLine();
- }
- this.sequenceString = sb.toString();
-
- return line;
- }
-
- /**
- * Processes a feature line. If it declares a feature type of interest
- * (currently, only CDS is processed), processes all of the associated lines
- * (feature qualifiers), and returns the next line after that, otherwise
- * simply returns the next line.
- *
- * @param line
- * the first line for the feature (with initial FT omitted for EMBL
- * format)
- * @return
- * @throws IOException
- */
- protected String parseFeature(String line) throws IOException
- {
- String[] tokens = line.trim().split(WHITESPACE);
- if (tokens.length < 2 || (!"CDS".equals(tokens[0]) && (!"source".equals(tokens[0]))))
- {
- return nextLine();
- }
- if (tokens[0].equals("source"))
- {
- return parseSourceQualifiers(tokens);
- }
- return parseCDSFeature(tokens[1]);
- }
-}
-
-/**
- * A data bean class to hold values parsed from one CDS Feature
- */
-class CdsData
-{
- String translation; // from /translation qualifier
-
- String cdsLocation; // the raw value e.g. join(1..1234,2012..2837)
-
- int codonStart = 1; // from /codon_start qualifier
-
- String proteinName; // from /product qualifier; used for protein description
-
- String proteinId; // from /protein_id qualifier
-
- List xrefs = new ArrayList<>(); // from /db_xref qualifiers
-
- Map cdsProps = new Hashtable<>(); // other qualifiers
-}