public class GeneLoci
{
/*
- * implemented as an adapter over DBRefEntry with
- * source -> species id
- * version -> reference
- * accession -> chromosome
+ * an identifier for the species
*/
- private DBRefEntry loci;
+ public final String species;
- boolean forwardStrand;
-
- /**
- * Constructor
- *
- * @param taxon
- * @param ref
- * @param chrId
- * @param map
- * @param forward
- */
- public GeneLoci(String taxon, String ref, String chrId, MapList map,
- boolean forward)
- {
- loci = new DBRefEntry(taxon, ref, chrId, new Mapping(map));
- forwardStrand = forward;
- }
-
- /**
- * Answers the identifier for the species
- *
- * @return
+ /*
+ * an identifier for a genome assembly, e.g. GRCh38
*/
- public String getSpecies()
- {
- return loci.getSource();
- }
+ public final String assembly;
- /**
- * Answers the identifier for the genomic reference assembly
+ /*
+ * a chromosome identifier, e.g. "5" or "X"
*/
- public String getReference()
- {
- return loci.getVersion();
- }
+ public final String chromosome;
- /**
- * Answers the chromosome identifier
- *
- * @return
+ /*
+ * mapping from sequence positions to chromosome locations;
+ * any regions with start > end are on the reverse strand
*/
- public String getChromosome()
- {
- return loci.getAccessionId();
- }
+ public final MapList mapping;
/**
- * Answers the mapping from sequence positions (in sequence start..end
- * coordinates) to the corresponding loci in the chromosome (in reference
- * assembly coordinates, base 1)
+ * Constructor
*
- * @return
+ * @param taxon
+ * @param ref
+ * @param chrId
+ * @param map
*/
- public MapList getMapping()
- {
- return loci.getMap().getMap();
- }
-
- public boolean isForwardStrand()
+ public GeneLoci(String taxon, String ref, String chrId, MapList map)
{
- return forwardStrand;
+ species = taxon;
+ assembly = ref;
+ chromosome = chrId;
+ mapping = map;
}
}
boolean forwardStrand = "1".equals(tokens[5]);
String species = ""; // dunno yet!
int[] from = new int[] { start, end };
- int[] to = new int[] { chStart, chEnd };
+ int[] to = new int[] { forwardStrand ? chStart : chEnd,
+ forwardStrand ? chEnd : chStart };
MapList map = new MapList(from, to, 1, 1);
- GeneLoci gl = new GeneLoci(species, ref, chrom, map, forwardStrand);
+ GeneLoci gl = new GeneLoci(species, ref, chrom, map);
setGeneLoci(gl);
} catch (NumberFormatException e)
{
* patch to ensure gene to chromosome mapping is complete
* (in case created before gene length was known)
*/
- MapList geneMapping = loci.getMapping();
+ MapList geneMapping = loci.mapping;
if (geneMapping.getFromRanges().get(0)[1] == 0)
{
geneMapping.getFromRanges().get(0)[0] = gene.getStart();
List<int[]> transcriptRange = Arrays.asList(new int[] {
transcript.getStart(), transcript.getEnd() });
MapList mapList = new MapList(transcriptRange, transcriptLoci, 1, 1);
- GeneLoci gl = new GeneLoci(loci.getSpecies(), loci.getReference(),
- loci.getChromosome(), mapList, loci.isForwardStrand());
+ GeneLoci gl = new GeneLoci(loci.species, loci.assembly,
+ loci.chromosome, mapList);
transcript.setGeneLoci(gl);
}
return null; // not used
}
+ /**
+ * Constructs a URL of the format <code>
+ * http://rest.ensembl.org/map/human/GRCh38/17:45051610..45109016:1/GRCh37?content-type=application/json
+ * </code>
+ *
+ * @param species
+ * @param chromosome
+ * @param fromRef
+ * @param toRef
+ * @param startPos
+ * @param endPos
+ * @return
+ * @throws MalformedURLException
+ */
protected URL getUrl(String species, String chromosome, String fromRef,
String toRef, int startPos, int endPos)
throws MalformedURLException
{
- String url = getDomain() + "/map/" + species + "/" + fromRef + "/"
- + chromosome + ":" + startPos + ".." + endPos + ":1/" + toRef
- + "?content-type=application/json";
+ /*
+ * start-end might be reverse strand - present forwards to the service
+ */
+ boolean forward = startPos <= endPos;
+ int start = forward ? startPos : endPos;
+ int end = forward ? endPos : startPos;
+ String strand = forward ? "1" : "-1";
+ String url = String.format(
+ "%s/map/%s/%s/%s:%d..%d:%s/%s?content-type=application/json",
+ getDomain(), species, fromRef, chromosome, start, end, strand,
+ toRef);
try
{
return new URL(url);
{
url = getUrl(species, chromosome, fromRef, toRef, queryRange[0],
queryRange[1]);
+ // System.out.println("Calling " + url);
br = getHttpResponse(url, null);
return (parseResponse(br));
} catch (Throwable t)
// todo check for "mapped"
JSONObject val = (JSONObject) rvals.next();
JSONObject mapped = (JSONObject) val.get("mapped");
- String start = mapped.get("start").toString();
- String end = mapped.get("end").toString();
- result = new int[] { Integer.parseInt(start), Integer.parseInt(end) };
+ int start = Integer.parseInt(mapped.get("start").toString());
+ int end = Integer.parseInt(mapped.get("end").toString());
+ String strand = mapped.get("strand").toString();
+ if ("1".equals(strand))
+ {
+ result = new int[] { start, end };
+ }
+ else
+ {
+ result = new int[] { end, start };
+ }
}
} catch (IOException | ParseException | NumberFormatException e)
{
{
String choice = chooser.getSelectedFile().getPath();
Cache.setProperty("LAST_DIRECTORY", choice);
- new VCFLoader(viewport.getAlignment()).loadVCF(choice);
+ new VCFLoader(viewport.getAlignment()).loadVCF(choice, this);
}
}
import htsjdk.variant.vcf.VCFHeaderLine;
import jalview.analysis.AlignmentUtils;
+import jalview.api.AlignViewControllerGuiI;
import jalview.datamodel.AlignmentI;
import jalview.datamodel.DBRefEntry;
import jalview.datamodel.GeneLoci;
import jalview.io.gff.SequenceOntologyI;
import jalview.util.MapList;
+import java.io.IOException;
import java.util.HashMap;
import java.util.List;
import java.util.Map;
* instances of this class.
*
* @param filePath
+ * @param status
*/
- public void loadVCF(String filePath)
+ public void loadVCF(String filePath, AlignViewControllerGuiI status)
{
VCFReader reader = null;
try
{
- long start = System.currentTimeMillis();
+ // long start = System.currentTimeMillis();
reader = new VCFReader(filePath);
VCFHeader header = reader.getFileHeader();
computePeptideVariants(seq);
}
}
- long elapsed = System.currentTimeMillis() - start;
- System.out.println(String.format(
- "Added %d VCF variants to %d sequence(s) (%dms)", varCount,
- seqCount, elapsed));
-
- reader.close();
+ // long elapsed = System.currentTimeMillis() - start;
+ String msg = String.format("Added %d VCF variants to %d sequence(s)",
+ varCount, seqCount);
+ if (status != null)
+ {
+ status.setStatus(msg);
+ }
} catch (Throwable e)
{
System.err.println("Error processing VCF: " + e.getMessage());
e.printStackTrace();
+ } finally
+ {
+ if (reader != null)
+ {
+ try
+ {
+ reader.close();
+ } catch (IOException e)
+ {
+ // ignore
+ }
+ }
}
}
return 0;
}
- MapList mapping = seqCoords.getMapping();
- List<int[]> seqChromosomalContigs = mapping.getToRanges();
+ List<int[]> seqChromosomalContigs = seqCoords.mapping.getToRanges();
for (int[] range : seqChromosomalContigs)
{
count += addVcfVariants(seq, reader, range, isVcfRefGrch37);
{
GeneLoci seqCoords = ((Sequence) seq).getGeneLoci();
- String chromosome = seqCoords.getChromosome();
- String seqRef = seqCoords.getReference();
- String species = seqCoords.getSpecies();
+ String chromosome = seqCoords.chromosome;
+ String seqRef = seqCoords.assembly;
+ String species = seqCoords.species;
// TODO handle species properly
if ("".equals(species))
/*
* query the VCF for overlaps
+ * (convert a reverse strand range to forwards)
*/
int count = 0;
- MapList mapping = seqCoords.getMapping();
+ MapList mapping = seqCoords.mapping;
+ int fromLocus = Math.min(range[0], range[1]);
+ int toLocus = Math.max(range[0], range[1]);
CloseableIterator<VariantContext> variants = reader.query(chromosome,
- range[0], range[1]);
+ fromLocus, toLocus);
while (variants.hasNext())
{
/*
/*
* save mapping for possible future re-use
*/
- String key = species + EXCL + chromosome + EXCL + fromRef + EXCL
- + toRef;
+ String key = makeRangesKey(chromosome, species, fromRef, toRef);
if (!assemblyMappings.containsKey(key))
{
assemblyMappings.put(key, new HashMap<int[], int[]>());
protected int[] findSubsumedRangeMapping(int[] queryRange, String chromosome,
String species, String fromRef, String toRef)
{
- String key = species + EXCL + chromosome + EXCL + fromRef + EXCL
- + toRef;
+ String key = makeRangesKey(chromosome, species, fromRef, toRef);
if (assemblyMappings.containsKey(key))
{
Map<int[], int[]> mappedRanges = assemblyMappings.get(key);
/*
* mapping is 1:1 in length, so we trust it to have no discontinuities
*/
- if (fromRange[0] <= queryRange[0] && fromRange[1] >= queryRange[1])
+ if (rangeContains(fromRange, queryRange))
{
/*
* fromRange subsumes our query range
}
return null;
}
+
+ /**
+ * Answers true if range's start-end positions include those of queryRange,
+ * where either range might be in reverse direction, else false
+ *
+ * @param range
+ * @param queryRange
+ * @return
+ */
+ protected static boolean rangeContains(int[] range, int[] queryRange)
+ {
+ int min = Math.min(range[0], range[1]);
+ int max = Math.max(range[0], range[1]);
+
+ return (min <= queryRange[0] && max >= queryRange[0]
+ && min <= queryRange[1] && max >= queryRange[1]);
+ }
+
+ /**
+ * Formats a ranges map lookup key
+ *
+ * @param chromosome
+ * @param species
+ * @param fromRef
+ * @param toRef
+ * @return
+ */
+ protected static String makeRangesKey(String chromosome, String species,
+ String fromRef, String toRef)
+ {
+ return species + EXCL + chromosome + EXCL + fromRef + EXCL
+ + toRef;
+ }
}