/*
* 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.analysis;
import jalview.datamodel.AlignedCodonFrame;
import jalview.datamodel.Alignment;
import jalview.datamodel.AlignmentI;
import jalview.datamodel.DBRefEntry;
import jalview.datamodel.Mapping;
import jalview.datamodel.Sequence;
import jalview.datamodel.SequenceFeature;
import jalview.datamodel.SequenceI;
import jalview.util.Comparison;
import jalview.util.DBRefUtils;
import jalview.util.MapList;
import jalview.ws.SequenceFetcherFactory;
import jalview.ws.seqfetcher.ASequenceFetcher;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.List;
/**
* Functions for cross-referencing sequence databases. user must first specify
* if cross-referencing from protein or dna (set dna==true)
*
* @author JimP
*
*/
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);
}
}
/**
* Returns a list of distinct database sources for which sequences have either
*
* - a (dna-to-protein or protein-to-dna) cross-reference
* - an indirect cross-reference - a (dna-to-protein or protein-to-dna)
* reference from another sequence in the dataset which has a cross-reference
* to a direct DBRefEntry on the given sequence
*
*
* @param dna
* true if seqs are nucleotide
* @param seqs
* sequences whose xrefs we are seeking
* @param dataset
* an alignment to search for indirect references
* @return
*/
public static List findXrefSourcesForSequences(boolean dna,
SequenceI[] seqs, AlignmentI dataset)
{
List sources = new ArrayList();
for (SequenceI seq : seqs)
{
if (seq != null)
{
findXrefSourcesForSequence(seq, dna, dataset, sources);
}
}
return sources;
}
/**
* Returns a list of distinct database sources for which a sequence has either
*
* - a (dna-to-protein or protein-to-dna) cross-reference
* - an indirect cross-reference - a (dna-to-protein or protein-to-dna)
* reference from another sequence in the dataset which has a cross-reference
* to a direct DBRefEntry on the given sequence
*
*
* @param seq
* the sequence whose dbrefs we are searching against
* @param dna
* true if the sequence is nucleotide
* @param dataset
* an alignment to search for indirect references
* @param sources
* a list of sources to add matches to
*/
static void findXrefSourcesForSequence(SequenceI seq, boolean dna,
AlignmentI dataset, List sources)
{
/*
* first find seq's xrefs (dna-to-peptide or peptide-to-dna)
*/
DBRefEntry[] rfs = DBRefUtils.selectDbRefs(!dna, seq.getDBRefs());
addXrefsToSources(rfs, sources);
if (dataset != null)
{
/*
* find sequence's direct (dna-to-dna, peptide-to-peptide) xrefs
*/
DBRefEntry[] lrfs = DBRefUtils.selectDbRefs(dna, seq.getDBRefs());
List rseqs = new ArrayList();
/*
* find sequences in the alignment which xref one of these DBRefs
* i.e. is xref-ed to a common sequence identifier
*/
CrossRef.searchDatasetXrefs(seq, !dna, lrfs, dataset, rseqs, null);
/*
* add those sequences' (dna-to-peptide or peptide-to-dna) dbref sources
*/
for (SequenceI rs : rseqs)
{
DBRefEntry[] xrs = DBRefUtils.selectDbRefs(!dna, rs.getDBRefs());
addXrefsToSources(xrs, sources);
}
}
}
/**
* Helper method that adds the source identifiers of some cross-references to
* a (non-redundant) list of database sources
*
* @param xrefs
* @param sources
*/
static void addXrefsToSources(DBRefEntry[] xrefs, List sources)
{
if (xrefs != null)
{
for (DBRefEntry ref : xrefs)
{
String source = ref.getSource();
if (!sources.contains(source))
{
sources.add(source);
}
}
}
}
/**
*
* @param seqs
* sequences whose xrefs are being retrieved
* @param dna
* true if sequences are nucleotide
* @param source
* @param al
* alignment to search for cross-referenced sequences (and possibly
* add to)
* @return products (as dataset sequences)
*/
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)
{
SequenceI dss = seq;
while (dss.getDatasetSequence() != null)
{
dss = dss.getDatasetSequence();
}
boolean found = false;
DBRefEntry[] xrfs = DBRefUtils.selectDbRefs(!dna, dss.getDBRefs());
if ((xrfs == null || xrfs.length == 0) && dataset != null)
{
/*
* found no suitable dbrefs on sequence - look for sequences in the
* alignment which share a dbref with this one
*/
DBRefEntry[] lrfs = DBRefUtils.selectDbRefs(dna, seq.getDBRefs());
/*
* find sequences (except this one!), of complementary type,
* which have a dbref to an accession id for this sequence,
* and add them to the results
*/
found = CrossRef.searchDatasetXrefs(dss, !dna, lrfs, dataset,
rseqs, cf);
}
for (int r = 0; xrfs != null && r < xrfs.length; r++)
{
DBRefEntry xref = xrfs[r];
if (source != null && !source.equals(xref.getSource()))
{
continue;
}
if (xref.hasMap())
{
if (xref.getMap().getTo() != null)
{
found = true;
SequenceI rsq = new Sequence(xref.getMap().getTo());
rseqs.add(rsq);
if (xref.getMap().getMap().getFromRatio() != xref
.getMap().getMap().getToRatio())
{
// get sense of map correct for adding to product alignment.
if (dna)
{
// map is from dna seq to a protein product
cf.addMap(dss, rsq, xref.getMap().getMap());
}
else
{
// map should be from protein seq to its coding dna
cf.addMap(rsq, dss, xref.getMap().getMap().getInverse());
}
}
}
}
if (!found)
{
// do a bit more work - search for sequences with references matching
// xrefs on this sequence.
if (dataset != null)
{
found = searchDataset(dss, xref, dataset, rseqs, cf, false,/*true?*/
!dna);
if (found)
{
xrfs[r] = null; // we've recovered seqs for this one.
}
}
}
}
if (!found)
{
if (xrfs != null && xrfs.length > 0)
{
ASequenceFetcher sftch = SequenceFetcherFactory
.getSequenceFetcher();
SequenceI[] retrieved = null;
int l = xrfs.length;
for (int r = 0; r < xrfs.length; r++)
{
// filter out any irrelevant or irretrievable references
if (xrfs[r] == null
|| ((source != null && !source.equals(xrfs[r]
.getSource())) || !sftch.isFetchable(xrfs[r]
.getSource())))
{
l--;
xrfs[r] = null;
}
}
if (l > 0)
{
// System.out
// .println("Attempting to retrieve cross referenced sequences.");
DBRefEntry[] t = new DBRefEntry[l];
l = 0;
for (int r = 0; r < xrfs.length; r++)
{
if (xrfs[r] != null)
{
t[l++] = xrfs[r];
}
}
xrfs = t;
try
{
retrieved = sftch.getSequences(Arrays.asList(xrfs), !dna);
// problem here is we don't know which of xrfs resulted in which
// retrieved element
} catch (Exception e)
{
System.err
.println("Problem whilst retrieving cross references for Sequence : "
+ 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].getDBRefs();
if (dbr != null && dbr.length > 0)
{
for (int di = 0; di < dbr.length; di++)
{
// find any entry where we should put in the sequence being
// cross-referenced into the map
Mapping map = dbr[di].getMap();
if (map != null)
{
if (map.getTo() != null && map.getMap() != null)
{
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
// if map is congruent
SequenceI ms = map.getTo();
int sf = map.getMap().getToLowest();
int st = map.getMap().getToHighest();
SequenceI mappedrg = ms.getSubSequence(sf, st);
// SequenceI loc = dss.getSubSequence(sf, st);
if (mappedrg.getLength() > 0
&& ms.getSequenceAsString().equals(
dss.getSequenceAsString()))
// && mappedrg.getSequenceAsString().equals(
// loc.getSequenceAsString()))
{
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());
}
// TODO remove this 'else' and the cf.addMap above?
else
{
cf.addMap(retrieved[rs].getDatasetSequence(),
map.getTo(), map.getMap());
}
} catch (Exception e)
{
System.err
.println("Exception when consolidating Mapped sequence set...");
e.printStackTrace(System.err);
}
}
}
}
}
retrieved[rs].updatePDBIds();
rseqs.add(retrieved[rs]);
}
}
}
}
}
}
Alignment ral = null;
if (rseqs.size() > 0)
{
ral = new Alignment(rseqs.toArray(new SequenceI[rseqs.size()]));
if (cf != null && !cf.isEmpty())
{
ral.addCodonFrame(cf);
}
}
return ral;
}
/**
* 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.
*
* @param sequenceI
* @param lrfs
* @param dataset
* @param rseqs
* @return true if matches were found.
*/
private static boolean searchDatasetXrefs(SequenceI sequenceI,
boolean dna, DBRefEntry[] lrfs, AlignmentI dataset,
List rseqs, AlignedCodonFrame cf)
{
boolean found = false;
if (lrfs == null)
{
return false;
}
for (int i = 0; i < lrfs.length; i++)
{
DBRefEntry xref = new DBRefEntry(lrfs[i]);
// add in wildcards
xref.setVersion(null);
xref.setMap(null);
found |= searchDataset(sequenceI, xref, dataset, rseqs, cf, false,
dna);
}
return found;
}
/**
* Searches dataset for DBRefEntrys matching the given one (xrf) and adds the
* associated sequence to rseqs
*
* @param sequenceI
* a sequence to ignore (start point of search)
* @param xrf
* a cross-reference to try to match
* @param dataset
* sequences to search in
* @param rseqs
* result list to add to
* @param cf
* a set of sequence mappings to add to
* @param direct
* - search all references or only subset
* @param dna
* search dna or protein xrefs (if direct=false)
* @return true if relationship found and sequence added.
*/
public static boolean searchDataset(SequenceI sequenceI, DBRefEntry xrf,
AlignmentI dataset, List rseqs, AlignedCodonFrame cf,
boolean direct, boolean dna)
{
boolean found = false;
if (dataset == null)
{
return false;
}
if (dataset.getSequences() == null)
{
System.err.println("Empty dataset sequence set - NO VECTOR");
return false;
}
List ds;
synchronized (ds = dataset.getSequences())
{
for (SequenceI nxt : ds)
{
if (nxt != null)
{
if (nxt.getDatasetSequence() != null)
{
System.err
.println("Implementation warning: getProducts passed a dataset alignment without dataset sequences in it!");
}
if (nxt == sequenceI || nxt == sequenceI.getDatasetSequence())
{
continue;
}
// check if this is the correct sequence type
{
// TODO 'direct' is always set to false - remove?
// or should it be 'true' from findXrefSequences?
// also its Javadoc conflicts with its use:
// test below implies 'direct' means find complementary sequences,
// !direct means select same molecule type
boolean isDna = Comparison
.isNucleotide(new SequenceI[] { nxt });
if ((direct && isDna == dna) || (!direct && isDna != dna))
{
// skip this sequence because it is wrong molecule type
continue;
}
}
// look for direct or indirect references in common
DBRefEntry[] poss = nxt.getDBRefs();
List cands = null;
/*
* TODO does this make any sense?
* if 'direct', search the dbrefs for xrf
* else, filter the dbrefs by type and then search for xrf
* - the result is the same isn't it?
*/
if (direct)
{
cands = DBRefUtils.searchRefs(poss, xrf);
}
else
{
poss = DBRefUtils.selectDbRefs(!dna, poss);
cands = DBRefUtils.searchRefs(poss, xrf);
}
if (!cands.isEmpty())
{
if (!rseqs.contains(nxt))
{
found = true;
rseqs.add(nxt);
if (cf != null)
{
// don't search if we aren't given a codon map object
for (DBRefEntry candidate : cands)
{
Mapping mapping = candidate.getMap();
if (mapping != null)
{
MapList map = mapping.getMap();
if (mapping.getTo() != null
&& map.getFromRatio() != map.getToRatio())
{
// get sense of map correct for adding to product
// alignment.
if (dna)
{
// map is from dna seq to a protein product
cf.addMap(sequenceI, nxt, map);
}
else
{
// map should be from protein seq to its coding dna
cf.addMap(nxt, sequenceI, map.getInverse());
}
}
}
}
}
// TODO: add mapping between sequences if necessary
}
}
}
}
}
return found;
}
}