Mercurial > repos > yufei-luo > s_mart
view SMART/Java/Python/structure/Interval.py @ 67:f4de72c80eac draft
Uploaded
author | m-zytnicki |
---|---|
date | Mon, 16 Nov 2015 11:59:35 -0500 |
parents | 169d364ddd91 |
children |
line wrap: on
line source
# # Copyright INRA-URGI 2009-2010 # # This software is governed by the CeCILL license under French law and # abiding by the rules of distribution of free software. You can use, # modify and/ or redistribute the software under the terms of the CeCILL # license as circulated by CEA, CNRS and INRIA at the following URL # "http://www.cecill.info". # # As a counterpart to the access to the source code and rights to copy, # modify and redistribute granted by the license, users are provided only # with a limited warranty and the software's author, the holder of the # economic rights, and the successive licensors have only limited # liability. # # In this respect, the user's attention is drawn to the risks associated # with loading, using, modifying and/or developing or reproducing the # software by the user in light of its specific status of free software, # that may mean that it is complicated to manipulate, and that also # therefore means that it is reserved for developers and experienced # professionals having in-depth computer knowledge. Users are therefore # encouraged to load and test the software's suitability as regards their # requirements in conditions enabling the security of their systems and/or # data to be ensured and, more generally, to use and operate it in the # same conditions as regards security. # # The fact that you are presently reading this means that you have had # knowledge of the CeCILL license and that you accept its terms. # from SMART.Java.Python.structure.Bins import * from commons.core.coord.Range import Range class Interval(Range): """ Store a genomic interval @ivar name: name of the interval [optional] @type name: string @ivar id: id of the interval [optional] @type id: int @ivar bin: bin in which the interval should be if stored in a database [computed] @type bin: int @ival tags: information about the transcript [optional] @type tags: dict @ivar verbosity: verbosity @type verbosity: int [default: 0] """ def __init__(self, interval = None, verbosity = 0): """ Constructor @param interval: interval to be copied @type interval: class L{Interval<Interval>} @param verbosity: verbosity @type verbosity: int """ Range.__init__(self) self.name = None self.id = None self.bin = None self.verbosity = verbosity self.tags = {} if interval != None: self.copy(interval) #!!!! Warning: two methods getStart() and getEnd() give the information maximum and minimum in interval.!!!!# #In case strand = "+", start < end; strand = "-", start > end def getStart(self): if self.start == -1: return -1 if self.end == -1: return self.start return self.getMin() def getEnd(self): if self.end == -1: return -1 if self.start == -1: return self.end return self.getMax() def getChromosome(self): return self.getSeqname() def getDirection(self): return 1 if self.getStrand() == "+" else -1 def getName(self): return self.name def isSet(self): """ Check if the interval is set """ return self.getStart() == None and self.getEnd() == None def copy(self, interval): """ Copy method @param interval: interval to be copied @type interval: class L{Interval<Interval>} """ self.setStart(interval.getStart()) self.setEnd(interval.getEnd()) self.setChromosome(interval.getChromosome()) self.setDirection(interval.getDirection()) self.name = interval.name self.id = interval.id self.bin = interval.bin self.tags = {} for tag in interval.tags: self.tags[tag] = interval.tags[tag] self.verbosity = interval.verbosity def setName(self, name): """ Set the name @param name: name of the interval @type name: string """ if len(name) > 100: name = name[:100] self.name = name def setChromosome(self, chromosome=""): """ Set the chromosome @param chromosome: chromosome on which the interval is @type chromosome: string """ if not chromosome: self.seqname = None else: self.seqname = chromosome.replace("|", "_") def setStart(self, start): """ Set the start point Possibly reset bin @param start: start point of the interval @type start: int """ self.bin = None direction = self.getDirection() if self.start == -1: self.start = start elif self.end == -1: self.end = start else: if direction == 1: self.start = start else: self.end = start if direction == 1: self.start, self.end = min(self.start, self.end), max(self.start, self.end) else: self.start, self.end = max(self.start, self.end), min(self.start, self.end) def setEnd(self, end): """ Set the end point Possibly reset bin @param end: end point of the interval of the interval @type end: int """ self.bin = None direction = self.getDirection() if self.end == -1: self.end = end elif self.start == -1: self.start = end else: if direction == 1: self.end = end else: self.start = end if direction == 1: self.start, self.end = min(self.start, self.end), max(self.start, self.end) else: self.start, self.end = max(self.start, self.end), min(self.start, self.end) def setSize(self, size): """ Possibly modify the end point @param size: size of the transcript @type size: int """ if self.end == None and self.start != None: self.setEnd(self.start + self.getSize() - 1) elif self.start == None and self.end != None: self.setStart(self.end - self.getSize() + 1) def getSize(self): """ Get the size """ return self.getEnd() - self.getStart() + 1 def _setDirection(self, direction): """ Set the direction of the interval (connection to Range) @param direction: direction of the transcript (+ / -) @type direction: int (1 or -1) """ if direction * self.getDirection() < 0: self.reverse() def setDirection(self, direction): """ Set the direction of the interval Possibly parse different formats @param direction: direction of the transcript (+ / -) @type direction: int or string """ if type(direction).__name__ == 'int': self._setDirection(direction / abs(direction)) elif type(direction).__name__ == 'str': if direction == "+": self._setDirection(1) elif direction == "-": self._setDirection(-1) elif direction == "1" or direction == "-1": self._setDirection(int(direction)) elif direction.lower() == "plus": self._setDirection(1) elif direction.lower() == "minus": self._setDirection(-1) else: raise Exception("Cannot understand direction %s" % (direction)) else: raise Exception("Cannot understand direction %s" % (direction)) def extendStart(self, size): """ Extend the interval by the 5' end @param size: the size to be exended @type size: int """ if self.getDirection() == 1: self.setStart(max(0, self.getStart() - size)) else: self.setEnd(self.getEnd() + size) self.bin = None def extendEnd(self, size): """ Extend the interval by the 3' end @param size: the size to be exended @type size: int """ if self.getDirection() == 1: self.setEnd(self.getEnd() + size) else: self.setStart(max(0, self.getStart() - size)) self.bin = None def restrictStart(self, size = 1): """ Restrict the interval by some nucleotides, start from its start position Remove the exons @param size: the size to be restricted to @type size: int """ if self.getDirection() == 1: self.setEnd(min(self.getEnd(), self.getStart() + size - 1)) else: self.setStart(max(self.getStart(), self.getEnd() - size + 1)) self.bin = None def restrictEnd(self, size = 1): """ Restrict the interval by some nucleotides, end from its end position Remove the exons @param size: the size to be restricted to @type size: int """ if self.getDirection() == 1: self.setStart(max(self.getStart(), self.getEnd() - size + 1)) else: self.setEnd(min(self.getEnd(), self.getStart() + size - 1)) self.bin = None def setTagValue(self, name, value): """ Set a tag @param name: name of the tag @type name: string @param value: value of the tag @type value: int or string """ self.tags[name] = value def getTagNames(self): """ Get all the names of the tags """ return self.tags.keys() def getTagValue(self, tag): """ Get the value of a tag @param tag: name of a tag @type tag: string """ if tag not in self.tags: return None return self.tags[tag] def getTagValues(self, tagSep = "; ", fieldSep = " ", surrounder = ""): """ Get the formatted tag values @param tagSep: separator between tags @type tagSep: string @param fieldSep: separator between tag name and tag value @type fieldSep: string @param surrounder: string which optionally surround values @type surrounder: string """ tags = [] for name, value in self.tags.iteritems(): if value == None: continue if isinstance(value, basestring): tags.append("%s%s%s%s%s" % (name, fieldSep, surrounder, value.replace("'", "\\'"), surrounder)) elif type(value) is int: tags.append("%s%s%s%i%s" % (name, fieldSep, surrounder, value, surrounder)) elif type(value) is float: tags.append("%s%s%s%f%s" % (name, fieldSep, surrounder, value, surrounder)) else: raise Exception("Do not know how to print '" + value + "'.") if self.getName() != None: tags.append("%s%s%s%s%s" % ("Name", fieldSep, surrounder, self.getName(), surrounder)) return tagSep.join(tags) def setTagValues(self, tags, tagSep = "; ", fieldSep = " "): """ Set the tag values using given string @param tags: the tags, concatenated @type tags: string @param tagSep: separator between tags @type tagSep: string @param fieldSep: separator between tag name and tag value @type fieldSep: string """ if tags == "": self.tags = {} return for splittedTag in tags.split(tagSep): if fieldSep not in splittedTag: raise Exception("Weird field '%s' in tags '%s'" % (splittedTag, tags)) tag, value = splittedTag.split(fieldSep, 1) if tag == "Name": self.setName(value) continue try: intValue = int(value) self.tags[tag] = intValue except ValueError: try: floatValue = float(value) self.tags[tag] = floatValue except ValueError: self.tags[tag] = value def deleteTag(self, tag): """ Remove a tag @param tag: the tag to be removed @type tag: string """ if tag in self.tags: del self.tags[tag] def setNbOccurrences(self, nbOccurrences): """ Set the number of occurrences of the interval @param nbOccurrences: number of occurrences of the interval @type nbOccurrences: int """ self.setTagValue("nbOccurrences", nbOccurrences) def setOccurrence(self, occurrence): """ Set the occurrence of this interval @param occurrence: an occurrence for this transcript @type occurrence: int """ self.setTagValue("occurrence", occurrence) def __eq__(self, interval): """ Whether two intervals are equal (start and end at same position) @param interval: object to be compared to @type interval: class L{Interval<Interval>} """ if not interval: return False return self.getChromosome() == interval.getChromosome() and self.getStart() == interval.getStart() and self.getEnd() == interval.getEnd() and self.getDirection() == interval.getDirection() def overlapWith(self, interval, nbNucleotides = 1): """ Whether two intervals overlap @param interval: object to be compared to @type interval: class L{Interval<Interval>} @param nbNucleotides: minimum number of nucleotides to declare and overlap @type nbNucleotides: int """ if self.getChromosome() != interval.getChromosome(): return False return (min(self.getEnd(), interval.getEnd()) - max(self.getStart(), interval.getStart()) + 1 >= nbNucleotides) def isIncludeIn(self, interval): return interval.include(self) def include(self, interval): """ Whether this interval includes the other one @param interval: object to be compared to @type interval: class L{Interval<Interval>} """ if self.getChromosome() != interval.getChromosome(): return False return ((self.getStart() <= interval.getStart()) and (self.getEnd() >= interval.getEnd())) def getDifference(self, interval, sameStrand = False): """ Get the difference between this cluster and another one @param interval: object to be compared to @type interval: class L{Interval<Interval>} @param sameStrand: do the comparison iff the intervals are on the same strand @type sameStrand: boolean @return: a (possibly empty) list of intervals """ newInterval = Interval() newInterval.copy(self) if self.getChromosome() != interval.getChromosome(): return [newInterval] if not self.overlapWith(interval): return [newInterval] if sameStrand and self.getDirection() != interval.getDirection(): return [newInterval] intervals = [] if self.getStart() < interval.getStart(): newInterval = Interval() newInterval.copy(self) newInterval.setEnd(min(self.getEnd(), interval.getStart() - 1)) intervals.append(newInterval) if self.getEnd() > interval.getEnd(): newInterval = Interval() newInterval.copy(self) newInterval.setStart(max(self.getStart(), interval.getEnd() + 1)) intervals.append(newInterval) return intervals def getIntersection(self, interval): """ Get the intersection between this interval and another one @param interval: object to be compared to @type interval: class L{Interval<Interval>} @return: an other interval """ if not self.overlapWith(interval): return None newInterval = Interval() newInterval.setChromosome(self.getChromosome()) newInterval.setDirection(self.getDirection()) newInterval.setName("%s_intersect_%s" % (self.getName(), interval.getName())) newInterval.setStart(max(self.getStart(), interval.getStart())) newInterval.setEnd(min(self.getEnd(), interval.getEnd())) return newInterval def getDistance(self, interval): """ Get the distance between two intervals (a non-negative value) @param interval: another interval @type interval: class L{Interval<Interval>} """ if self.overlapWith(interval): return 0 if self.getChromosome() != interval.getChromosome(): raise Exception("Cannot get the distance between %s and %s" % (str(self), str(interval))) return min(abs(self.getStart() - interval.getEnd()), abs(self.getEnd() - interval.getStart())) def getRelativeDistance(self, interval): """ Get the distance between two intervals (negative if first interval is before) @param interval: another interval @type interval: class L{Interval<Interval>} """ if self.overlapWith(interval): return 0 if self.getChromosome() != interval.getChromosome(): raise Exception("Cannot get the distance between %s and %s" % (str(self), str(interval))) if self.getEnd() < interval.getStart(): distance = interval.getStart() - self.getEnd() else: distance = interval.getEnd() - self.getStart() distance *= self.getDirection() return distance def merge(self, interval, normalization = False): """ Merge two intervals @param interval: another interval @type interval: class L{Interval<Interval>} @param normalization: whether the sum of the merge should be normalized wrt the number of mappings of each elements @type normalization: boolean """ if self.getChromosome() != interval.getChromosome(): raise Exception("Cannot merge '%s' and '%s' for they are on different chromosomes." % (str(self), str(interval))) direction = None if self.getStart() == self.getEnd(): direction = interval.getDirection() elif interval.getStart() == interval.getEnd(): direction = self.getDirection() elif self.getDirection() != interval.getDirection(): raise Exception("Cannot merge '%s' and '%s' for they are on different strands." % (str(self), str(interval))) self.setStart(min(self.getStart(), interval.getStart())) self.setEnd(max(self.getEnd(), interval.getEnd())) if direction != None: self.setDirection(direction) nbElements = 0.0 for element in (self, interval): for tagName in ("nbElements", "nbOccurrences"): if tagName not in element.getTagNames(): element.setTagValue(tagName, 1) nbElements += float(element.getTagValue("nbElements")) / float(element.getTagValue("nbOccurrences")) if normalization else float(element.getTagValue("nbElements")) self.setTagValue("nbElements", nbElements) self.bin = None for tagName in ("identity", "nbOccurrences", "occurrence", "nbMismatches", "nbGaps", "rank", "evalue", "bestRegion"): if tagName in self.getTagNames(): del self.tags[tagName] def getBin(self): """ Get the bin of the interval Computed on the fly """ if self.bin == None: self.bin = getBin(self.getStart(), self.getEnd()) return self.bin def getBins(self): """ Get all the bin this interval could fall into """ return getOverlappingBins(self.getStart(), self.getEnd()) def getSqlVariables(cls): """ Get the properties of the object that should be saved in a database """ variables = ["name", "chromosome", "start", "end", "direction", "tags", "bin"] return variables getSqlVariables = classmethod(getSqlVariables) def setSqlValues(self, array): """ Set the values of the properties of this object as given by a results line of a SQL query """ self.id = array[0] self.name = array[1].strip("'") self.setChromosome(array[2].strip("'")) self.setStart(array[3]) self.setEnd(array[4]) self.setDirection(array[5]) self.setTagValues(array[6].strip("'"), ";", "=") self.bin = array[7] def getSqlValues(self): """ Get the values of the properties that should be saved in a database """ values = dict() values["name"] = self.name values["chromosome"] = self.getChromosome() values["start"] = self.getStart() values["end"] = self.getEnd() values["direction"] = self.getDirection() values["tags"] = self.getTagValues(";", "=") values["bin"] = self.getBin() return values def getSqlTypes(cls): """ Get the values of the properties that should be saved in a database """ types = dict() types["name"] = "varchar" types["chromosome"] = "varchar" types["start"] = "int" types["end"] = "int" types["direction"] = "tinyint" types["tags"] = "varchar" types["bin"] = "int" return types getSqlTypes = classmethod(getSqlTypes) def getSqlSizes(cls): """ Get the sizes of the properties that should be saved in a database """ sizes = dict() sizes["name"] = 255 sizes["chromosome"] = 255 sizes["start"] = 11 sizes["end"] = 11 sizes["direction"] = 4 sizes["tags"] = 1023 sizes["bin"] = 11 return sizes getSqlSizes = classmethod(getSqlSizes) def printCoordinates(self): """ Print start and end positions (depending on the direction of the interval) """ if self.getDirection() == 1: return "%d-%d" % (self.getStart(), self.getEnd()) else: return "%d-%d" % (self.getEnd(), self.getStart()) def extractSequence(self, parser): """ Get the sequence corresponding to this interval @param parser: a parser to a FASTA file @type parser: class L{SequenceListParser<SequenceListParser>} @return : a instance of L{Sequence<Sequence>} """ return parser.getSubSequence(self.getChromosome(), self.getStart(), self.getEnd(), self.getDirection(), self.name) def extractWigData(self, parser): """ Get the data retrieved from a wig file @param parser: a parser class to a WIG file @type parser: class L{WigParser<WigParser>} """ data = parser.getRange(self.getChromosome(), self.getStart(), self.getEnd()) if self.getDirection() == -1: if parser.strands: newData = {} for strand in data: data[strand].reverse() newData[-strand] = data[strand] data = newData else: data.reverse() return data def __str__(self): """ Output a simple representation of this interval """ direction = "+" if self.getDirection() == -1: direction = "-" string = "%s:%d-%d (%s)" % (self.getChromosome(), self.getStart(), self.getEnd(), direction) if self.name != "": string = "(%s) %s" % (self.name, string) return string