source: orange-bioinformatics/orangecontrib/bio/obiGene.py @ 1927:f3f90e6b14dc

Revision 1927:f3f90e6b14dc, 32.6 KB checked in by Ales Erjavec <ales.erjavec@…>, 4 months ago (diff)

Optimize the case where the organism identifier is a known common taxid.

Line 
1from __future__ import absolute_import
2
3import os, time
4
5from Orange.orng import orngServerFiles
6
7from . import obiTaxonomy
8from . import obiKEGG
9
10default_database_path = orngServerFiles.localpath("NCBI_geneinfo")
11
12class GeneInfo(object):
13    """ An object representing the NCBI information for a gene.
14    """
15
16    NCBI_GENEINFO_TAGS = ("tax_id", "gene_id", "symbol", "locus_tag", "synonyms",
17                          "dbXrefs", "chromosome", "map_location", "description", "type",
18                          "symbol_from_nomenclature_authority", "full_name_from_nomenclature_authority",
19                          "nomenclature_status", "other_designations", "modification_date")
20    NCBI_MULTIPLE_CARDINALITY_TAGS = ("synonyms", "dbXrefs", "other_designations")
21   
22    __slots__ = NCBI_GENEINFO_TAGS
23    def __init__(self, line):
24        """ Construct the GeneInfo object from a line in the NCBI gene_info file
25        """
26        line = line.split("\t")
27        for attr, value in zip(self.__slots__, line):
28            if value == "-":
29                value = None
30            if attr in GeneInfo.NCBI_MULTIPLE_CARDINALITY_TAGS:
31                value = value.split("|") if value != None else []
32            setattr(self, attr, value)
33
34    def __repr__(self):
35        def format(value):
36            if not value:
37                return "-"
38            elif type(value) == list:
39                return "|".join(value)
40            else:
41                return value
42        return "\t".join(format(getattr(self, slot)) for slot in self.__slots__)
43
44    def __str__(self):
45        return repr(self)
46
47class GeneHistory(object):
48    NCBI_GENE_HISTORY_TAGS = ("tax_id", "gene_id", "discontinued_gene_id", "discontinued_symbol", "discontinue_date")
49    __slots__ = NCBI_GENE_HISTORY_TAGS
50    def __init__(self, line):
51        for attr, value in zip(self.__slots__, line.split("\t")):
52            setattr(self, attr, value)
53           
54           
55class NCBIGeneInfo(dict):
56    TAX_MAP = {
57            "2104": "272634",  # Mycoplasma pneumoniae
58            "4530": "39947",  # Oryza sativa
59            "5833": "36329",  # Plasmodium falciparum
60            "4932": "559292",  # Saccharomyces cerevisiae
61            }
62       
63    def __init__(self, organism, genematcher=None):
64        """ An dictionary like object for accessing NCBI gene info
65        Arguments::
66                - *organism*    Organism id
67
68        Example::
69            >>> info = NCBIGeneInfo("Homo sapiens")
70        """
71       
72        self.taxid = self.organism_name_search(organism)
73
74
75        fname = orngServerFiles.localpath_download("NCBI_geneinfo", "gene_info.%s.db" % self.taxid)
76        file = open(fname, "rb")
77        self.update(dict([(line.split("\t", 3)[1], line) for line in file.read().splitlines() if line.strip() and not line.startswith("#")]))
78
79        # NOTE orig init time for gene matcher: 2.5s, new 4s: investigate the slowdown
80        # NOTE matches are not the same because aliases are build a bit
81        # differently (main name versus old aliases conflict!)
82
83        self.matcher = genematcher
84        if self.matcher == None:
85            if self.taxid == '352472':
86                self.matcher = matcher([GMNCBI(self.taxid), GMDicty(), [GMNCBI(self.taxid), GMDicty()]])
87            else:
88                self.matcher = matcher([GMNCBI(self.taxid)])
89
90        #if this is done with a gene matcher, pool target names
91        self.matcher.set_targets(self.keys())
92       
93    def history(self):
94        if getattr(self, "_history", None) is None:
95            fname = orngServerFiles.localpath_download("NCBI_geneinfo", "gene_history.%s.db" % self.taxid)
96            try:
97                self._history = dict([(line.split("\t")[2], GeneHistory(line)) for line in open(fname, "rb").read().splitlines()])
98               
99            except Exception, ex:
100                print >> sys.srderr, "Loading NCBI gene history failed.", ex
101                self._history = {}
102        return self._history
103       
104    @classmethod
105    def organism_version(cls, name):
106        oname = cls.organism_name_search(name)
107        #FIXME, dirty hack to ensure file id downloaded
108        orngServerFiles.localpath_download("NCBI_geneinfo", "gene_info.%s.db" % oname) 
109        return orngServerFiles.info("NCBI_geneinfo", "gene_info.%s.db" % oname)["datetime"]
110
111    @classmethod
112    def organism_name_search(cls, org):
113        if org in cls.common_taxids():
114            return org
115        elif org in NCBIGeneInfo.TAX_MAP:
116            return NCBIGeneInfo.TAX_MAP[org]
117
118        taxids = obiTaxonomy.to_taxid(org, mapTo=cls.common_taxids())
119        if not taxids:
120            taxids = obiTaxonomy.search(org, onlySpecies=False)
121            taxids = set(cls.common_taxids()).intersection(taxids) #onlySpecies=False needed to find correct dicty
122        if len(taxids) == 0:
123            raise obiTaxonomy.UnknownSpeciesIdentifier, org
124        elif len(taxids) > 1:
125            raise obiTaxonomy.MultipleSpeciesException, ", ".join(["%s: %s" % (id, obiTaxonomy.name(id)) for id in taxids])
126        taxid = taxids.pop()
127        return cls.TAX_MAP.get(taxid, taxid)
128
129    @classmethod   
130    def load(cls, file):
131        """ A class method that loads gene info from file
132        """
133        if type(file) in [str, unicode]:
134            file = open(file, "rb")
135        return cls((line.split("\t", 3)[1], line) for line in file.read().splitlines() if line.strip() and not line.startswith("#"))
136       
137    def get_info(self, gene_id, def_=None):
138        """ Search and return the GeneInfo object for gene_id
139        """
140        try:
141            return self(gene_id)
142        except KeyError:
143            return def_
144       
145    def __call__(self, name):
146        """ Search and return the GeneInfo object for gene_id
147        """
148        #id = self.translate.get(name, name)
149        #print self.matcher.umatch(name), self.matcher.match(name)
150        id = self.matcher.umatch(name)
151        return self[id]
152
153    def __getitem__(self, key):
154#        return self.get(gene_id, self.matcher[gene_id])
155        return GeneInfo(dict.__getitem__(self, key))
156
157    def __setitem__(self, key, value):
158        if type(value) == str:
159            dict.__setitem__(self, key, value)
160        else:
161            dict.__setitem__(self, key, repr(value))
162
163    def get(self, key, def_=None):
164        try:
165            return self[key]
166        except KeyError:
167            return def_
168
169    def itervalues(self):
170        for val in dict.itervalues(self):
171            yield GeneInfo(val)
172
173    def iteritems(self):
174        for key, val in zip(self.iterkeys(), self.itervalues()):
175            yield key, val
176
177    def values(self):
178        return list(self.itervalues())
179   
180    def items(self):
181        return list(self.iteritems())
182
183    @staticmethod
184    def get_geneinfo_from_ncbi(file, progressCallback=None):
185        import urllib2, gzip, shutil, tempfile
186        from cStringIO import StringIO
187        if isinstance(file, basestring):
188            file = open(file, "wb")
189       
190        stream = urllib2.urlopen("ftp://ftp.ncbi.nih.gov/gene/DATA/gene_info.gz")
191        tmpfile = tempfile.TemporaryFile()
192        shutil.copyfileobj(stream, tmpfile)
193        tmpfile.seek(0)
194        stream = gzip.GzipFile(None, "rb", fileobj=tmpfile)
195        shutil.copyfileobj(stream, file)
196       
197    @staticmethod
198    def get_gene_history_from_ncbi(file, progressCallback=None):
199        import urllib2, gzip, shutil, tempfile
200        from cStringIO import StringIO
201        if isinstance(file, basestring):
202            file = open(file, "wb")
203       
204        stream = urllib2.urlopen("ftp://ftp.ncbi.nih.gov/gene/DATA/gene_history.gz")
205        tmpfile = tempfile.TemporaryFile()
206        shutil.copyfileobj(stream, tmpfile)
207        tmpfile.seek(0)
208        stream = gzip.GzipFile(None, "rb", fileobj=tmpfile)
209        shutil.copyfileobj(stream, file)
210       
211    @classmethod
212    def common_taxids(cls):
213        taxids = obiTaxonomy.common_taxids()
214        return [cls.TAX_MAP.get(id, id) for id in taxids if cls.TAX_MAP.get(id, id)]
215   
216    @classmethod
217    def essential_taxids(cls):
218        taxids = obiTaxonomy.essential_taxids()
219        return [cls.TAX_MAP.get(id, id) for id in taxids if cls.TAX_MAP.get(id, id)]
220
221
222class EnsembleGeneInfo(object):
223#    BIO_MART_DATASETS = {"9606": "hsapiens"}
224    DEF_ATTRS = ["ensembl_gene_id", "external_gene_id", "entrezgene"]
225    BIOMART_CONF = {"9606": ("hsapiens_gene_ensembl", DEF_ATTRS) # + ["hgnc_symbol"])
226                    }
227    def __init__(self, organism, gene_matcher=None):
228        self.organism = self.organism_name_search(organism) 
229        self.load()
230        if gene_matcher is None:
231            self.gene_matcher = matcher([GMEnsembl(self.organism), GMNCBI(self.organism)])
232        else:
233            self.gene_matcher = gene_matcher
234       
235    @classmethod
236    def organism_name_search(cls, name):
237        taxids = obiTaxonomy.to_taxid(name, mapTo=cls.common_taxids())
238        if len(taxids) == 1:
239            return taxids[0]
240        else:
241            raise ValueError(name)
242       
243    @classmethod
244    def common_taxids(cls):
245        return ["3702", "9913", "6239", "7955", "9606", "7227",
246                "10090", "10116", "4932", "4896", "8355"]
247   
248    @classmethod
249    def version(cls):
250        return "v1.0"
251       
252    def filename(self):
253        return "ensembl_" + self.organism
254   
255    def create_info(self):
256        from . import obiBioMart
257        if self.organism in self.BIOMART_CONF:
258            dset_name, attrs = self.BIOMART_CONF[self.organism]
259        else:
260            dset_name, attrs = self.default_biomart_conf(self.organism)
261       
262        dataset = obiBioMart.BioMartDataset("ensembl", dset_name)
263        table = dataset.get_example_table(attributes=attrs, unique=True)
264        print len(table)
265        table.save(dset_name + ".tab")
266        from collections import defaultdict
267        names = defaultdict(set)
268        for ex in table:
269            row = [str(v) for v in ex if not v.isSpecial() and str(v)]
270            names[row[0]].update(row)
271        return names
272       
273    def default_biomart_conf(self, taxid):
274        name = obiTaxonomy.name(self.organism).lower()
275        name1, name2 = name.split(" ")[: 2]
276        dset_name = name1[0] + name2 + "_gene_ensembl"
277        return dset_name, self.DEF_ATTRS
278   
279    def load(self):
280        import cPickle
281        dir = orngServerFiles.localpath("EnsembleGeneInfo")
282        if not os.path.exists(dir):
283            os.makedirs(dir)
284       
285        try:
286            filename = orngServerFiles.localpath_download("EnsembleGeneInfo", self.filename())
287            info = cPickle.load(open(filename, "rb"))
288        except Exception, ex:   
289            filename = orngServerFiles.localpath("EnsembleGeneInfo", self.filename())
290            info = self.create_info()
291            cPickle.dump(info, open(filename, "wb"))
292           
293        self.info = info
294       
295    def __getitem__(self, key):
296        return self.info[key]
297    def __contains__(self, key):
298        return key in self.info
299    def __len__(self):
300        return len(self.info)
301    def __iter__(self):
302        return iter(self.info)
303    def keys(self):
304        return self.info.keys()
305    def values(self):
306        return self.info.values()
307    def items(self):
308        return self.info.items()
309    def get(self, key, default=None):
310        return self.info.get(key, default)
311   
312    def genes(self):
313        return self.info.keys()
314   
315    def aliases(self):
316        return [(key,) + tuple(value) for key, value in self.items()]
317   
318    def ensembl_id(self, name):
319        return self.gene_matcher.umatch(name)
320       
321"""
322Gene matcher.
323
324"Database" for each organism is a list of sets of gene aliases.
325"""
326
327from collections import defaultdict
328import os
329
330gene_matcher_path = None
331
332def ignore_case(gs):
333    """ Transform names in sets in list to lower case """
334    return [ set([a.lower() for a in g]) for g in gs ]
335
336def create_mapping(groups, lower=False):
337    """
338    Returns mapping of aliases to the group index. If lower
339    is True, lower case forms of gene aliases are mapped to indices.
340
341    Unpickling the results of this function (binary format)
342    is slower than running it.
343
344    TIMING NOTES:
345    - lower() costs are neglible (< 10%)
346    - building sets instead of lists also costs about 10 percent
347    """
348    togroup = defaultdict(set)
349
350    # code duplicated because running a function in relatively expensive here.
351    if lower: 
352        for i,group in enumerate(groups):
353            for alias in group:
354                togroup[alias.lower()].add(i)
355    else:
356        for i,group in enumerate(groups):
357            for alias in group:
358                togroup[alias].add(i)
359
360    return togroup
361
362def join_sets(set1, set2, lower=False):
363    """
364    Joins two sets of gene set mappings. If lower is True, lower case
365    forms of gene aliases are compared.
366
367    A group g1 from set1 is joined to a group of aliases g2 from set2,
368    if the groups share at least one gene.
369    Returns all joined groups and groups that were not matched, which
370    remain unchanged.
371
372    The operation both commutative and associative.
373    """
374
375    set1 = [ set(a) for a in set1 ]
376    set2 = [ set(a) for a in set2 ]
377
378    currentmap = create_mapping(set1, lower=lower)
379
380    new = [] #new groups
381
382    #remember used to find unused
383    set1used = set() 
384    set2used = set()
385
386    fn = lambda x: x
387    if lower:
388        fn = lambda x: x.lower()
389
390    for i, group in enumerate(set2):
391
392        #find groups of aliases (from set1)  intersecting with a current
393        #group from set2
394        cross = reduce(set.union, 
395            [ currentmap[fn(alias)] for alias in group if fn(alias) in currentmap], set())
396
397        for c in cross:
398            #print c, group & set1[c], group, set1[c]
399            set1used.add(c)
400            set2used.add(i)
401            new.append(group | set1[c]) #add a union
402
403    #add groups without matches (from both sets)
404    set1add = set(range(len(set1))) - set1used
405    set2add = set(range(len(set2))) - set2used
406
407    for a in set1add:
408        new.append(set1[a])
409    for a in set2add:
410        new.append(set2[a])
411
412    return new
413 
414def join_sets_l(lsets, lower=False):
415    """
416    Joins multiple gene set mappings using join_sets function.
417    """
418    current = lsets[0]
419    for b in lsets[1:]:
420        current = join_sets(current, b, lower=lower)
421    return current
422
423class Matcher(object):
424    """
425    Gene matcher tries to match an input gene to some target.
426    """
427
428    def copy(self):
429        return notImplemented()
430
431    def __call__(self, targets):
432        return self.set_targets(targets)
433
434    def set_targets(self, targets):
435        """
436        Set input list of gene names as targets.
437        Abstract function.
438        """
439        notImplemented()
440
441    def match(self, gene):
442        """Returns a list of matching target gene names."""
443        notImplemented()
444
445    def umatch(self, gene):
446        """Returns an unique (only one matching target) target or None"""
447        mat = self.match(gene)
448        return mat[0] if len(mat) == 1 else None
449
450    def explain(self, gene):
451        """
452        Returns an gene matches with explanations as lists of tuples.
453        Each tuple consists of a list of target genes in a set
454        of aliases matched to input gene, returned as a second part
455        of the tuple.
456        """
457        notImplemented()
458
459def buffer_path():
460    """ Returns buffer path from Orange's setting folder if not
461    defined differently (in gene_matcher_path). """
462    if  gene_matcher_path == None:
463        from Orange.orng import orngEnviron
464        pth = os.path.join(orngEnviron.directoryNames["bufferDir"], 
465            "gene_matcher")
466        try:
467            os.makedirs(pth)
468        except:
469            pass
470        return pth
471    else:
472        return gene_matcher_path
473
474def auto_pickle(filename, version, func, *args, **kwargs):
475    """
476    Run function func with given arguments and save the results to
477    a file named filename. If results for a given filename AND
478    version were already saved, just read and return them.
479    """
480
481    import cPickle as pickle
482
483    output = None
484    outputOk = False
485
486    try:
487        f = open(filename,'rb')
488
489        try:
490            versionF = pickle.load(f)
491            if version == None or versionF == version:
492                outputOk = True
493                output = pickle.load(f)
494        except:
495            pass
496        finally:
497            f.close()
498
499    except:
500        pass
501
502    if not outputOk:
503        output = func(*args, **kwargs)
504
505        #save output before returning
506        f = open(filename,'wb')
507        pickle.dump(version, f, -1)
508        pickle.dump(output, f, -1)
509        f.close()
510
511    return output
512
513class MatcherAliases(Matcher):
514    """
515    Genes matcher based on a list of sets of given aliases.
516
517    Target genes belonging to same sets of aliases as the input gene are
518    returned as matching genes.
519
520    """
521    def __init__(self, aliases, ignore_case=True):
522        self.aliases = aliases
523        self.ignore_case = ignore_case
524        self.mdict = create_mapping(self.aliases, self.ignore_case)
525
526    def to_ids(self, gene):
527        """ Return ids of sets of aliases the gene belongs to. """
528        if self.ignore_case:
529            gene = gene.lower()
530        return self.mdict[gene]
531
532    def set_targets(self, targets):
533        """
534        A reverse dictionary is made according to each target's membership
535        in the sets of aliases.
536        """
537        d = defaultdict(list)
538        #d = id: [ targets ], where id is index of the set of aliases
539        for target in targets:
540            ids = self.to_ids(target)
541            if ids != None:
542                for id in ids:
543                    d[id].append(target)
544        mo = MatchAliases(d, self)
545        self.matcho = mo #backward compatibility - default match object
546        return mo
547
548    #this two functions are solely for backward compatibility
549    def match(self, gene):
550        return self.matcho.match(gene)
551    def explain(self, gene):
552        return self.matcho.explain(gene)
553
554class Match(object):
555
556    def umatch(self, gene):
557        """Returns an unique (only one matching target) target or None"""
558        mat = self.match(gene)
559        return mat[0] if len(mat) == 1 else None
560 
561class MatchAliases(Match):
562
563    def __init__(self, to_targets, parent):
564        self.to_targets = to_targets
565        self.parent = parent
566
567    def match(self, gene):
568        """
569        Input gene is first mapped to ids of sets of aliases which contain
570        it. Target genes belonding to the same sets of aliases are returned
571        as input's match.
572        """
573        inputgeneids = self.parent.to_ids(gene)
574        #return target genes with same ids
575        return list(set( \
576            reduce(lambda x,y: x+y, 
577                [ self.to_targets[igid] for igid in inputgeneids ], [])))
578
579    def explain(self, gene):
580        inputgeneids = self.parent.to_ids(gene)
581        return [ (self.to_targets[igid], self.parent.aliases[igid]) for igid in inputgeneids ]
582
583class MatcherAliasesPickled(MatcherAliases):
584    """
585    Gene matchers based on sets of aliases supporting pickling should
586    extend this class. Subclasses must define functions "filename",
587    "create_aliases_version" and "create_aliases". Those are crucial for
588    pickling of gene aliases to work.
589
590    Loading of gene aliases is done lazily: they are loaded when they are
591    needed. Loading of aliases for components of joined matchers is often
592    unnecessary and is therefore avoided.
593    """
594   
595    def set_aliases(self, aliases):
596        self.saved_aliases = aliases
597
598    def get_aliases(self):
599        if not self.saved_aliases: #loads aliases if not loaded
600            self.aliases = self.load_aliases()
601        #print "size of aliases ", len(self.saved_aliases)
602        return self.saved_aliases
603
604    aliases = property(get_aliases, set_aliases)
605
606    def get_mdict(self):
607        """ Creates mdict. Aliases are loaded if needed. """
608        if not self.saved_mdict:
609            self.saved_mdict = create_mapping(self.aliases, self.ignore_case)
610        return self.saved_mdict
611
612    def set_mdict(self, mdict):
613        self.saved_mdict = mdict
614
615    mdict = property(get_mdict, set_mdict)
616
617    def set_targets(self, targets):
618        return MatcherAliases.set_targets(self, targets)
619
620    def filename(self):
621        """ Returns file name for saving aliases. """
622        notImplemented()
623       
624    def create_aliases_version(self):
625        """ Returns version of the source database state. """
626        notImplemented()
627
628    def create_aliases(self):
629        """ Returns gene aliases. """
630        notImplemented()
631
632    def load_aliases(self):
633        fn = self.filename()
634        ver = self.create_aliases_version() #if version == None ignore it
635        if fn != None:
636            if isinstance(fn, tuple): #if you pass tuple, look directly
637               filename = fn[0]
638            else:
639               filename = os.path.join(buffer_path(), fn)
640            return auto_pickle(filename, ver, self.create_aliases)
641        else:
642            #if either file version of version is None, do not pickle
643            return self.create_aliases()
644
645    def __init__(self, ignore_case=True):
646        self.aliases = []
647        self.mdict = {}
648        self.ignore_case = ignore_case
649        self.filename() # test if valid filename can be built
650
651
652class MatcherAliasesKEGG(MatcherAliasesPickled):
653
654    def _organism_name(self, organism):
655        return obiKEGG.organism_name_search(organism)
656
657    def create_aliases(self):
658        org = obiKEGG.KEGGOrganism(self.organism, genematcher=GMDirect())
659        osets = org.gene_aliases()
660        return osets
661
662    def create_aliases_version(self):
663        # use KEGG short release string (e.g. '66.0+')
664        release = obiKEGG.KEGGOrganism.organism_version(self.organism) + ".2"
665        release, _ = release.split("/")
666        return release
667
668    def filename(self):
669        return "kegg_2_" + self._organism_name(self.organism)
670
671    def __init__(self, organism, ignore_case=True):
672        self.organism = organism
673        MatcherAliasesPickled.__init__(self, ignore_case=ignore_case)
674
675
676class MatcherAliasesFile(MatcherAliasesPickled):
677
678    def create_aliases(self):
679        canNotCreateButCanOnlyOpen()
680
681    def create_aliases_version(self):
682        return None
683
684    def filename(self):
685        return (self.filename_,)
686
687    def __init__(self, filename, ignore_case=True):
688        self.filename_ = filename
689        MatcherAliasesPickled.__init__(self, ignore_case=ignore_case)
690
691
692class MatcherAliasesGO(MatcherAliasesPickled):
693
694    def _organism_name(self, organism):
695        """ Returns internal GO organism name. Used to define file name. """
696        from . import obiGO
697        return obiGO.organism_name_search(self.organism)
698
699    def create_aliases(self):
700        from . import obiGO
701        annotations = obiGO.Annotations(self.organism, genematcher=GMDirect())
702        names = annotations.geneNamesDict
703        return map(set, list(set([ \
704            tuple(sorted(set([name]) | set(genes))) \
705            for name,genes in names.items() ])))
706
707    def filename(self):
708        return "go_" + self._organism_name(self.organism)
709
710    def create_aliases_version(self):
711        from . import obiGO
712        return "v2." + obiGO.Annotations.organism_version(self.organism)
713
714    def __init__(self, organism, ignore_case=True):
715        self.organism = organism
716        MatcherAliasesPickled.__init__(self, ignore_case=ignore_case)
717
718class MatcherAliasesDictyBase(MatcherAliasesPickled):
719
720    def create_aliases(self):
721        from . import obiDicty
722        db = obiDicty.DictyBase()
723        #db.info, db.mappings
724        infoa = [ set([id,name]) | set(aliases) for id,(name,aliases,_) in db.info.items() ]
725        mappingsa = [ set(filter(None, a)) for a in db.mappings ]
726        joineda = join_sets(infoa, mappingsa, lower=True)
727        return joineda
728
729    def create_aliases_version(self):
730        from . import obiDicty
731        return "v1." + obiDicty.DictyBase.version()
732
733    def filename(self):
734        return "dictybase" 
735
736    def __init__(self, ignore_case=True):
737        MatcherAliasesPickled.__init__(self, ignore_case=ignore_case)
738
739class MatcherAliasesNCBI(MatcherAliasesPickled):
740
741    def _organism_name(self, organism):
742        return NCBIGeneInfo.organism_name_search(organism)
743
744    def create_aliases(self):
745        ncbi = NCBIGeneInfo(self.organism, genematcher=GMDirect())
746        out = []
747        for k in ncbi.keys():
748            out.append(set(filter(None, [k, ncbi[k].symbol, ncbi[k].locus_tag] + [ s for s in ncbi[k].synonyms ] )))
749        return out
750
751    def filename(self):
752        return "ncbi_" + self._organism_name(self.organism)
753
754    def create_aliases_version(self):
755        return "v2." + NCBIGeneInfo.organism_version(self.organism)
756
757    def __init__(self, organism, ignore_case=True):
758        self.organism = organism
759        MatcherAliasesPickled.__init__(self, ignore_case=ignore_case)
760       
761class MatcherAliasesAffy(MatcherAliasesPickled):
762    def create_aliases(self):
763        filename = orngServerFiles.localpath_download("Affy", self.organism + ".pickle")
764        import cPickle
765        return cPickle.load(open(filename, "rb"))
766   
767    def filename(self):
768        return "affy_" + self.organism
769   
770    def create_aliases_version(self):
771        orngServerFiles.localpath_download("Affy", self.organism + ".pickle")
772        return orngServerFiles.info("Affy", self.organism + ".pickle")["datetime"]
773       
774    def __init__(self, organism, **kwargs):
775        self.organism = organism
776        MatcherAliasesPickled.__init__(self, **kwargs)
777   
778       
779class MatcherAliasesEnsembl(MatcherAliasesPickled):
780    """ A matcher for Ensemble ids
781    """
782    DEF_ATTRS = ["ensembl_gene_id", "external_gene_id", "entrezgene"]
783    # taxid: (dataset_name, [name_attr1, name_attr2 ...])
784    BIOMART_CONF = {}
785    def __init__(self, organism, **kwargs):
786        self.organism = organism
787        MatcherAliasesPickled.__init__(self, **kwargs)
788       
789    def filename(self):
790        return "ensembl_" + self.organism
791   
792    def create_aliases_version(self):
793        return "v1"
794   
795    def create_aliases(self):
796        from . import obiBioMart
797        if self.organism in self.BIOMART_CONF:
798            dset_name, attrs = self.BIOMART_CONF[self.organism]
799        else:
800            dset_name, attrs = self.default_biomart_conf(self.organism)
801       
802        dataset = obiBioMart.BioMartDataset("ensembl", dset_name)
803        table = dataset.get_example_table(attributes=attrs)
804        from collections import defaultdict
805        names = defaultdict(set)
806        for ex in table:
807            row = [str(v) for v in ex if not v.isSpecial() and str(v)]
808            names[row[0]].update(row)
809        return names.values()
810       
811    def default_biomart_conf(self, taxid):
812        name = obiTaxonomy.name(self.organism).lower()
813        name1, name2 = name.split(" ")[: 2]
814        dset_name = name1[0] + name2 + "_gene_ensembl"
815        return dset_name, self.DEF_ATTRS
816       
817
818class MatcherAliasesPickledJoined(MatcherAliasesPickled):
819    """
820    Creates a new matcher by joining gene aliases from different data sets.
821    Sets of aliases are joined if they contain common genes.
822
823    The joined gene matcher can only be pickled if the source gene
824    matchers are picklable.
825    """
826
827    def filename(self):
828        # do not pickle if any is unpicklable
829        try:
830            filenames = [ mat.filename() for mat in self.matchers ]
831            if self.ignore_case:
832                filenames += [ "icj" ]
833            return "__".join(filenames)
834        except:
835            return None
836
837    def create_aliases(self):
838        return join_sets_l([ mat.aliases for mat in self.matchers ], lower=self.ignore_case)
839
840    def create_aliases_version(self):
841        try:
842            return "v4_" + "__".join([ mat.create_aliases_version() for mat in self.matchers ])
843        except:
844            return None
845
846    def __init__(self, matchers):
847        """
848        Join matchers together. Groups of aliases are joined if
849        they share a common name.
850
851        If ignore_case is True, ignores case when joining gene aliases.
852        """
853        #FIXME: sorting of matchers to avoid multipying pickled files for
854        #different orderings.
855        self.matchers = matchers
856        allic = set([ m.ignore_case for m in self.matchers ])
857        if len(allic) > 1:
858            notAllMatchersHaveEqualIgnoreCase()
859        ignore_case = list(allic)[0]
860
861        MatcherAliasesPickled.__init__(self, ignore_case=ignore_case)
862       
863class MatcherSequence(Matcher):
864    """
865    Each gene goes through sequence of gene matchers (in the same order
866    as in the matchers arguments) until a match is found.
867    """
868   
869    def __init__(self, matchers):
870        self.matchers = matchers
871
872    def set_targets(self, targets):
873        ms = []
874        targets = list(targets) #copy targets as multiple use would
875                                #be problematic if a generator was passed
876        for matcher in self.matchers:
877            ms.append(matcher.set_targets(targets))
878        om = MatchSequence(ms)
879        self.matcho = om
880        return om
881
882    #this two functions are solely for backward compatibility
883    def match(self, gene):
884        return self.matcho.match(gene)
885
886    def explain(self, gene):
887        return self.matcho.explain(gene)
888
889class MatchSequence(Match):
890
891    def __init__(self, ms):
892        self.ms = ms
893
894    def match(self, gene):
895        for match in self.ms:
896            m = match.match(gene)
897            if m: 
898                return m
899        return []
900
901    def explain(self, gene):
902        for match in self.ms:
903            m = match.match(gene)
904            if m: 
905                return match.explain(gene)
906        return []
907
908class MatcherDirect(Matcher):
909    """
910    Direct matching to targets.
911    """
912
913    def __init__(self, ignore_case=True):
914        self.ignore_case = ignore_case
915
916    def set_targets(self, targets):
917        targets = list(targets) #would be problematic if generator was passed
918                                #as it is used twice
919        aliases = [ set([a]) for a in targets]
920        self.am = MatcherAliases(aliases, ignore_case=self.ignore_case)
921        self.matcho = self.am.set_targets(targets)
922        return self.matcho
923
924    #this two functions are solely for backward compatibility
925    def match(self, gene):
926        return self.matcho.match(gene)
927    def explain(self, gene):
928        return self.matcho.explain(gene)
929
930               
931GMDirect = MatcherDirect
932GMKEGG = MatcherAliasesKEGG
933GMGO = MatcherAliasesGO
934GMNCBI = MatcherAliasesNCBI
935GMDicty = MatcherAliasesDictyBase
936GMAffy = MatcherAliasesAffy
937GMEnsembl = MatcherAliasesEnsembl
938
939def issequencens(x):
940    return hasattr(x, '__getitem__') and not isinstance(x, basestring)
941
942def matcher(matchers, direct=True, ignore_case=True):
943    """
944    Build a matcher from a sequence of matchers. If a sequence element is a
945    sequence, join matchers in the subsequence.
946
947    direct - if True, add a direct matcher to targets
948    ignore_case - if True, ignores case with optionally added direct matcher
949    """
950    seqmat = []
951    if direct:
952        seqmat.append(MatcherDirect(ignore_case=ignore_case))
953    for mat in matchers:
954        if issequencens(mat):
955            mat = MatcherAliasesPickledJoined(list(mat))
956            seqmat.append(mat)
957        else:
958            seqmat.append(mat)
959    return MatcherSequence(seqmat)
960
961if __name__ == '__main__':
962
963    #m1 = matcher([[GMNCBI('44689'), GMDicty()]])
964    #print m1.matchers[1].aliases[:100]
965
966    #m2 = GMNCBI('Dictyostelium discoideum')
967    #print m2.aliases
968
969
970
971    """
972    gi = info(list(info)[0])
973    print gi.tax_id, gi.synonyms, gi.dbXrefs, gi.symbol_from_nomenclature_authority, gi.full_name_from_nomenclature_authority
974    """
975
976    #dobim z joinom prave stvari?
977
978    import time
979    from . import obiGeneSets
980
981    def testsets():
982        return obiGeneSets.collections([":kegg:hsa", ":go:hsa"])
983
984    def names1():
985        import orange
986        data = orange.ExampleTable("DLBCL.tab")
987        return [ a.name for a in  data.domain.attributes ]
988
989    def namesd():
990        import orange
991        data = orange.ExampleTable("dd_ge_biorep1.tab")
992        names = [ ex["gene"].value for ex in data ]
993        return names
994
995    genesets = auto_pickle("testcol", "3", testsets)
996    names = auto_pickle("testnames", "4", names1)
997    names2 = auto_pickle("testnamesdicty", "4", namesd)
998
999    info = NCBIGeneInfo('Dictyostelium discoideum')
1000    for a in names2:
1001        print a
1002        info.get_info(a)
1003
1004    t = time.time()
1005    mat5 = matcher([[GMKEGG('human'),GMGO('human')]], direct=False, ignore_case=True)
1006    mat7 = GMDicty()
1007    mat8 = GMNCBI('Homo sapiens')
1008    print "initialized all", time.time()-t
1009
1010    print "using targets"
1011
1012    mat5.set_targets(names)
1013    mat7.set_targets(names)
1014    mat8.set_targets(names)
1015
1016    print "before genes"
1017    genes = reduce(set.union, genesets.values()[:1000], set())
1018    genes = list(genes)
1019    print genes[:20]
1020    print len(genes)
1021    print "after genes"
1022
1023    for g in sorted(genes):
1024        print "KGO ", g, mat5.match(g), mat5.explain(g)
1025        print "DICT", g, mat7.match(g)
1026        print "NCBI", g, mat8.match(g)
1027
1028
Note: See TracBrowser for help on using the repository browser.