Changeset 1583:2b90e7c95085 in orange-bioinformatics


Ignore:
Timestamp:
03/05/12 13:51:06 (2 years ago)
Author:
markotoplak
Branch:
default
Message:

Fixed a bug in obiGene (umatch did not work in Match objects). Gene set signatures method now accept examples from different domains on the input.

Files:
2 edited

Legend:

Unmodified
Added
Removed
  • obiGene.py

    r1578 r1583  
    545545 
    546546class Match(object): 
    547     pass 
     547 
     548    def umatch(self, gene): 
     549        """Returns an unique (only one matching target) target or None""" 
     550        mat = self.match(gene) 
     551        return mat[0] if len(mat) == 1 else None 
    548552  
    549553class MatchAliases(Match): 
  • obiGeneSetSig.py

    r1578 r1583  
    66import numpy 
    77from collections import defaultdict 
    8 import statc 
    98import stats 
     9import obiGsea 
    1010 
    11 def selectGenesetsData(data, matcher, geneSets, minSize=3, maxSize=1000, minPart=0.1, classValues=None): 
    12     """ 
    13     Returns gene sets and data which falling under upper criteria. 
    14     """ 
    15     gso = obiGsea.GSEA(data, matcher=matcher, classValues=classValues, atLeast=0) 
    16     gso.addGenesets(geneSets) 
    17     okgenesets = gso.selectGenesets(minSize=minSize, maxSize=maxSize, minPart=minPart).keys() 
    18     gsetsnum = gso.to_gsetsnum(okgenesets) 
    19     return gso.data, okgenesets, gsetsnum 
     11def setSig_example_geneset(ex, data): 
     12    """ Gets learning data and example with the same domain, both 
     13    containing only genes from the gene set. """ 
     14 
     15    distances = [ [], [] ]     
     16 
     17    def pearsonr(v1, v2): 
     18        return numpy.corrcoef([v1, v2])[0,1] 
     19 
     20    def pearson(ex1, ex2): 
     21        #leaves undefined elements out 
     22 
     23        attrs = range(len(ex1.domain.attributes)) 
     24        vals1 = [ ex1[i].value for i in attrs ] 
     25        vals2 = [ ex2[i].value for i in attrs ] 
     26 
     27        common = [ True if v1 != "?" and v2 != "?" else False \ 
     28            for v1,v2 in zip(vals1,vals2) ] 
     29        vals1 = [ v for v,c in zip(vals1, common) if c ] 
     30        vals2 = [ v for v,c in zip(vals2, common) if c ] 
     31 
     32        return pearsonr(vals1, vals2) 
     33 
     34    def ttest(ex1, ex2): 
     35        try: 
     36            return stats.lttest_ind(ex1, ex2)[0] 
     37        except: 
     38            return 0.0 
     39     
     40    #maps class value to its index 
     41    classValueMap = dict( [ (val,i) for i,val in enumerate(data.domain.class_var.values) ]) 
     42  
     43    #create distances to all learning data - save or other class 
     44    for c in data: 
     45        distances[classValueMap[c[-1].value]].append(pearson(c, ex)) 
     46 
     47    return ttest(distances[0], distances[1]) 
     48 
     49def mat_ni(data, matcher): 
     50    nm = matcher([at.name for at in data.domain.attributes]) 
     51    name_ind = dict((n.name,i) for i,n in enumerate(data.domain.attributes)) 
     52    return nm, name_ind 
     53 
     54def select_genesets(nm, gene_sets, min_size=3, max_size=1000, min_part=0.1): 
     55    """ Returns a list of gene sets that have sizes in limits """ 
     56 
     57    def ok_sizes(gs): 
     58        """compares sizes of genesets to limitations""" 
     59        transl = filter(lambda x: x != None, [ nm.umatch(gene) for gene in gs.genes ]) 
     60        if len(transl) >= min_size \ 
     61            and len(transl) <= max_size \ 
     62            and float(len(transl))/len(gs.genes) >= min_part: 
     63            return True 
     64        return False 
     65 
     66    return filter(ok_sizes, gene_sets)  
    2067 
    2168class SetSig(object): 
     
    3279 
    3380    def __call__(self, data, weight_id=None): 
    34         data, oknames, gsetsnum = obiAssess.selectGenesetsData(data,  
    35             self.matcher, self.gene_sets, 
    36             minSize=self.min_size, maxSize=self.max_size,  
    37             minPart=self.min_part, classValues=self.class_values) 
    3881 
    39         def setSig_example_geneset(ex, data): 
    40             """ ex contains only selected genes """ 
    41  
    42             distances = [ [], [] ]     
    43  
    44             def pearsonr(v1, v2): 
    45                 try: 
    46                     return statc.pearsonr(v1, v2)[0] 
    47                 except: 
    48                     return numpy.corrcoef([v1, v2])[0,1] 
    49  
    50             def pearson(ex1, ex2): 
    51                 attrs = range(len(ex1.domain.attributes)) 
    52                 vals1 = [ ex1[i].value for i in attrs ] 
    53                 vals2 = [ ex2[i].value for i in attrs ] 
    54                 return pearsonr(vals1, vals2) 
    55  
    56             def ttest(ex1, ex2): 
    57                 try: 
    58                     return stats.lttest_ind(ex1, ex2)[0] 
    59                 except: 
    60                     return 0.0 
    61              
    62             #maps class value to its index 
    63             classValueMap = dict( [ (val,i) for i,val in enumerate(data.domain.classVar.values) ]) 
    64           
    65             #create distances to all learning data - save or other class 
    66             for c in data: 
    67                 distances[classValueMap[c[-1].value]].append(pearson(c, ex)) 
    68  
    69             return ttest(distances[0], distances[1]) 
     82        data = obiGsea.takeClasses(data, classValues=self.class_values) 
     83        nm,_ =  mat_ni(data, matcher) 
     84        gene_sets = select_genesets(nm, self.gene_sets, self.min_size, self.max_size, self.min_part) 
    7085 
    7186        attributes = [] 
    7287 
    73         for name, gs in gsetsnum.items(): #for each geneset 
    74             #for each gene set: take the attribute subset and work on the attribute subset only 
    75             #only select the subset of genes from the learning data 
    76             at = Orange.feature.Continuous(name=name.id) 
     88        for gs in gene_sets: 
     89            at = Orange.feature.Continuous(name=gs.id) 
    7790 
    78             def t(ex, w, gs=gs, ldata=data): 
    79                 domain = Orange.data.Domain([ldata.domain.attributes[ai] for ai in gs], ldata.domain.classVar) 
    80                 datao = Orange.data.Table(domain, ldata) 
    81                 example = Orange.data.Instance(domain, ex) #domains need to be the same 
    82                 return setSig_example_geneset(example, datao) 
     91            def t(ex, w, gs=gs, data=data): #copy od the data 
     92                geneset = list(gs.genes) 
     93 
     94                nm, name_ind = mat_ni(data, matcher) 
     95                nm2, name_ind2 = mat_ni(ex, matcher) 
     96 
     97                genes = [ nm.umatch(gene) for gene in geneset ] 
     98                genes2 = [ nm2.umatch(gene) for gene in geneset ] 
     99 
     100                genes, genes2 = zip(*[ (g,g2) for g,g2 in zip(genes, genes2) if g != None]) 
     101 
     102                domain = Orange.data.Domain([data.domain.attributes[name_ind[gene]] for gene in genes], data.domain.class_var) 
     103                datao = Orange.data.Table(domain, data) 
     104 
     105                def vou(ex, gn, indices): 
     106                    """ returns the value or unknown for the given gene name""" 
     107                    if gn not in indices: 
     108                        return "?" 
     109                    else: 
     110                        return ex[indices[gn]].value 
     111                 
     112                #convert the example to the same domain 
     113                exvalues = [ vou(ex, gn, name_ind2) for gn in genes2 ] + [ "?" ] 
     114                example = Orange.data.Instance(domain, exvalues) 
     115 
     116                return setSig_example_geneset(example, datao) #only this one is setsig specific 
    83117          
    84118            at.get_value_from = t 
     
    104138 
    105139    choosen_cv = ["Iris-setosa", "Iris-versicolor"] 
     140 
    106141    def to_old_dic(d, data): 
    107142        ar = defaultdict(list) 
     
    116151        print '\n'.join([ a + ": " +str(b) for a,b in ol]) 
    117152 
    118     ass = SetSig(ldata, matcher=matcher, gene_sets=gsets, class_values=choosen_cv, min_part=0.0) 
    119     print ass.domain 
     153    ass = SetSig(matcher=matcher, gene_sets=gsets, class_values=choosen_cv, min_part=0.0) 
     154    ass = ass(data) 
    120155    ar = to_old_dic(ass.domain, data[:5]) 
    121156    pp2(ar) 
Note: See TracChangeset for help on using the changeset viewer.