source: orange-bioinformatics/obiGene.py @ 1468:13f7c83ffda0

Revision 1468:13f7c83ffda0, 30.8 KB checked in by ales_erjavec <ales.erjavec@…>, 3 years ago (diff)

Added EnsembleGeneInfo class.

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