Changeset 1607:1daf61f47beb in orange-bioinformatics


Ignore:
Timestamp:
03/26/12 20:55:03 (2 years ago)
Author:
markotoplak
Branch:
default
Message:

Added ASSESS to obiGeneSetSig.

Files:
2 edited

Legend:

Unmodified
Added
Removed
  • obiAssess.py

    r1600 r1607  
    673673    choosen_cv = ["Iris-setosa", "Iris-versicolor"] 
    674674    #ass = AssessLearner()(data, matcher, gsets, rankingf=AT_loessLearner()) 
    675     #ass = AssessLearner()(data, matcher, gsets, minPart=0.0) 
     675    ass = AssessLearner()(data, matcher, gsets, minPart=0.0) 
    676676    #ass = MeanLearner()(data, matcher, gsets, default=False) 
    677     ass = CORGsLearner()(data, matcher, gsets) 
     677    #ass = CORGsLearner()(data, matcher, gsets) 
    678678    #ass = MedianLearner()(data, matcher, gsets) 
    679679    #ass = PLSLearner()(data, matcher, gsets, classValues=choosen_cv, minPart=0.0) 
  • obiGeneSetSig.py

    r1601 r1607  
    11import Orange 
     2import obiGsea 
    23import Orange.utils 
    34import obiGeneSets 
     
    1011import obiExpression 
    1112import scipy.stats 
    12  
    13 #STILL MISSING: Assess 
    14  
     13import math 
     14import statc 
     15 
     16class GeneSetTrans(object): 
     17 
     18    __new__ = Orange.utils._orange__new__(object) 
     19 
     20    def _mat_ni(self, data): 
     21        """ With cached gene matchers. """ 
     22        if data.domain not in self._cache: 
     23            self._cache[data.domain] = mat_ni(data, self.matcher) 
     24        return self._cache[data.domain] 
     25 
     26    def _match_instance(self, instance, geneset, takegenes=None): 
     27        nm, name_ind = self._mat_ni(instance) 
     28        genes = [ nm.umatch(gene) for gene in geneset ] 
     29        if takegenes: 
     30            genes = [ genes[i] for i in takegenes ] 
     31        return nm, name_ind, genes 
     32 
     33    def _match_data(self, data, geneset, odic=False): 
     34        nm, name_ind = self._mat_ni(data) 
     35        genes = [ nm.umatch(gene) for gene in geneset ] 
     36        if odic: 
     37            to_geneset = dict(zip(genes, geneset)) 
     38        takegenes = [ i for i,a in enumerate(genes) if a != None ] 
     39        genes = [ genes[i] for i in takegenes ] 
     40        if odic: 
     41            return nm, name_ind, genes, takegenes, to_geneset 
     42        else: 
     43            return nm, name_ind, genes, takegenes 
     44 
     45    def __init__(self, matcher=None, gene_sets=None, min_size=3, max_size=1000, min_part=0.1, class_values=None): 
     46        self.matcher = matcher 
     47        self.gene_sets = gene_sets 
     48        self.min_size = min_size 
     49        self.max_size = max_size 
     50        self.min_part = min_part 
     51        self.class_values = class_values 
     52        self._cache = {} 
     53 
     54    def __call__(self, data, weight_id=None): 
     55 
     56        #selection of classes and gene sets 
     57        data = takeClasses(data, classValues=self.class_values) 
     58        nm,_ =  self._mat_ni(data) 
     59        gene_sets = select_genesets(nm, self.gene_sets, self.min_size, self.max_size, self.min_part) 
     60 
     61        #build a new domain 
     62        newfeatures = self.build_features(data, gene_sets) 
     63        newdomain = Orange.data.Domain(newfeatures, data.domain.class_var) 
     64        return Orange.data.Table(newdomain, data) 
     65 
     66    def build_features(self, data, gene_sets): 
     67        return [ self.build_feature(data, gs) for gs in gene_sets ] 
     68 
     69 
     70def normcdf(x, mi, st): 
     71    return 0.5*(2. - stats.erfcc((x - mi)/(st*math.sqrt(2)))) 
     72 
     73class AT_edelmanParametric(object): 
     74 
     75    def __init__(self, **kwargs): 
     76        for a,b in kwargs.items(): 
     77            setattr(self, a, b) 
     78 
     79    def __call__(self, nval): 
     80 
     81        if self.mi1 == None or self.mi2 == None or self.st1 == None or self.st2 == None: 
     82            return 0.0  
     83 
     84        val = nval 
     85 
     86        try: 
     87            if val >= self.mi1: 
     88                p1 = 1 - normcdf(val, self.mi1, self.st1) 
     89            else: 
     90                p1 = normcdf(val, self.mi1, self.st1) 
     91 
     92            if val >= self.mi2: 
     93                p2 = 1 - normcdf(val, self.mi2, self.st2) 
     94            else: 
     95                p2 = normcdf(val, self.mi2, self.st2) 
     96 
     97            #print p1, p2 
     98            return math.log(p1/p2) 
     99        except: 
     100            #print p1, p2, "exception" 
     101            return 0 
     102 
     103class AT_edelmanParametricLearner(object): 
     104    """ 
     105    Returns attribute transfromer for Edelman parametric measure for a 
     106    given attribute in the dataset.  Edelman et al, 06. Modified a bit. 
     107    """ 
     108 
     109    def __init__(self, a=None, b=None): 
     110        """ 
     111        a and b are choosen class values. 
     112        """ 
     113        self.a = a 
     114        self.b = b 
     115 
     116    def __call__(self, i, data): 
     117        cv = data.domain.classVar 
     118        #print data.domain 
     119 
     120        if self.a == None: self.a = cv.values[0] 
     121        if self.b == None: self.b = cv.values[1] 
     122 
     123        def avWCVal(value): 
     124            return [ex[i].value for ex in data if ex[-1].value == value and not ex[i].isSpecial() ] 
     125 
     126        list1 = avWCVal(self.a) 
     127        list2 = avWCVal(self.b) 
     128 
     129        mi1 = mi2 = st1 = st2 = None 
     130 
     131        try: 
     132            mi1 = statc.mean(list1) 
     133            st1 = statc.std(list1) 
     134        except: 
     135            pass 
     136     
     137        try: 
     138            mi2 = statc.mean(list2) 
     139            st2 = statc.std(list2) 
     140        except: 
     141            pass 
     142 
     143        return AT_edelmanParametric(mi1=mi1, mi2=mi2, st1=st1, st2=st2) 
     144 
     145class AT_loess(object): 
     146 
     147    def __init__(self, **kwargs): 
     148        for a,b in kwargs.items(): 
     149            setattr(self, a, b) 
     150 
     151    def __call__(self, nval): 
     152 
     153        val = nval 
     154 
     155        def saveplog(a,b): 
     156            try: 
     157                return math.log(a/b) 
     158            except: 
     159                if a < b: 
     160                    return -10 
     161                else: 
     162                    return +10 
     163 
     164        try: 
     165            ocene = self.condprob(val) 
     166            if sum(ocene) < 0.01: 
     167                return 0.0 
     168            return saveplog(ocene[0], ocene[1]) 
     169 
     170        except: 
     171            return 0.0 
     172 
     173class AT_loessLearner(object): 
     174 
     175    def __call__(self, i, data): 
     176        try: 
     177            ca = Orange.statistics.contingency.VarClass(data.domain.attributes[i], data) 
     178            a =  Orange.statistics.estimate.ConditionalLoess(ca, nPoints=5) 
     179            return AT_loess(condprob=a) 
     180        except: 
     181            return AT_loess(condprob=None) 
     182 
     183def nth(l, n): 
     184    return [a[n] for a in l] 
     185 
     186class Assess(GeneSetTrans): 
     187    """ 
     188    Uses the underlying GSEA code to select genes. 
     189    Takes data and creates attribute transformations. 
     190    """ 
     191 
     192    def __init__(self, rankingf=None, **kwargs): 
     193        self.rankingf = rankingf 
     194        if self.rankingf == None: 
     195            self.rankingf = AT_edelmanParametricLearner() 
     196        super(Assess, self).__init__(**kwargs) 
     197 
     198    def build_features(self, data, gene_sets): 
     199 
     200        attributes = [] 
     201 
     202        #attrans: { i_orig: ranking_function } 
     203        attrans = [ self.rankingf(iat, data) for iat, at in enumerate(data.domain.attributes) ] 
     204 
     205        nm_all, _ =  self._mat_ni(data) 
     206 
     207        for gs in gene_sets: 
     208 
     209            at = Orange.feature.Continuous(name=str(gs)) 
     210 
     211            geneset = list(gs.genes) 
     212            nm, name_ind, genes, takegenes, to_geneset = self._match_data(data, geneset, odic=True) 
     213            genes = set(genes) 
     214             
     215            def t(ex, w, geneset=geneset, takegenes=takegenes, nm=nm, attrans=attrans): 
     216 
     217                nm2, name_ind2, genes2 = self._match_instance(ex, geneset, takegenes) 
     218 
     219                ex_atts = [ at.name for at in ex.domain.attributes ] 
     220                new_atts = [ name_ind[nm.umatch(an)] if nm.umatch(an) != None else None 
     221                    for an in ex_atts ] 
     222                #new_atts: indices of genes in original data for that sample  
     223                #POSSIBLE REVERSE IMPLEMENTATION (slightly different 
     224                #for data from different chips): 
     225                #save pairs together and sort (or equiv. dictionary transformation) 
     226 
     227                indexes = filter(lambda x: x[0] != None, zip(new_atts, range(len(ex_atts)))) 
     228 
     229                lcor = [ attrans[index_in_data](ex[index_in_ex].value)  
     230                    for index_in_data, index_in_ex in indexes if 
     231                    ex[index_in_ex].value != '?' ] 
     232                #indexes in original lcor, sorted from higher to lower values 
     233                ordered = obiGsea.orderedPointersCorr(lcor)  
     234                #subset = list of indices, lcor = correlations, ordered = order 
     235                subset = [ name_ind2[g] for g in genes2 ] 
     236                return obiGsea.enrichmentScoreRanked(subset, lcor, ordered)[0]  
     237 
     238            at.get_value_from = t 
     239            attributes.append(at) 
     240 
     241        return attributes 
     242    
    15243def setSig_example_geneset(ex, data): 
    16244    """ Gets learning data and example with the same domain, both 
     
    69297 
    70298    return filter(ok_sizes, gene_sets)  
    71  
    72 class GeneSetTrans(object): 
    73  
    74     __new__ = Orange.utils._orange__new__(object) 
    75  
    76     def _mat_ni(self, data): 
    77         """ With cached gene matchers. """ 
    78         if data.domain not in self._cache: 
    79             self._cache[data.domain] = mat_ni(data, self.matcher) 
    80         return self._cache[data.domain] 
    81  
    82     def _match_instance(self, instance, geneset, takegenes=None): 
    83         nm, name_ind = self._mat_ni(instance) 
    84         genes = [ nm.umatch(gene) for gene in geneset ] 
    85         if takegenes: 
    86             genes = [ genes[i] for i in takegenes ] 
    87         return nm, name_ind, genes 
    88  
    89     def _match_data(self, data, geneset, odic=False): 
    90         nm, name_ind = self._mat_ni(data) 
    91         genes = [ nm.umatch(gene) for gene in geneset ] 
    92         if odic: 
    93             to_geneset = dict(zip(genes, geneset)) 
    94         takegenes = [ i for i,a in enumerate(genes) if a != None ] 
    95         genes = [ genes[i] for i in takegenes ] 
    96         if odic: 
    97             return nm, name_ind, genes, takegenes, to_geneset 
    98         else: 
    99             return nm, name_ind, genes, takegenes 
    100  
    101     def __init__(self, matcher=None, gene_sets=None, min_size=3, max_size=1000, min_part=0.1, class_values=None): 
    102         self.matcher = matcher 
    103         self.gene_sets = gene_sets 
    104         self.min_size = min_size 
    105         self.max_size = max_size 
    106         self.min_part = min_part 
    107         self.class_values = class_values 
    108         self._cache = {} 
    109  
    110     def __call__(self, data, weight_id=None): 
    111  
    112         #selection of classes and gene sets 
    113         data = takeClasses(data, classValues=self.class_values) 
    114         nm,_ =  self._mat_ni(data) 
    115         gene_sets = select_genesets(nm, self.gene_sets, self.min_size, self.max_size, self.min_part) 
    116  
    117         #build a new domain 
    118         newfeatures = self.build_features(data, gene_sets) 
    119         newdomain = Orange.data.Domain(newfeatures, data.domain.class_var) 
    120         return Orange.data.Table(newdomain, data) 
    121  
    122     def build_features(self, data, gene_sets): 
    123         return [ self.build_feature(data, gs) for gs in gene_sets ] 
    124299 
    125300def vou(ex, gn, indices): 
     
    403578        print '\n'.join([ a + ": " +str(b) for a,b in ol]) 
    404579 
    405     ass = CORGs(data, matcher=matcher, gene_sets=gsets, class_values=choosen_cv, min_part=0.0) 
     580    ass = Assess(data, matcher=matcher, gene_sets=gsets, class_values=choosen_cv, min_part=0.0) 
    406581    ar = to_old_dic(ass.domain, data[:5]) 
    407582    pp2(ar) 
Note: See TracChangeset for help on using the changeset viewer.