source: orange/orange/orngLR_Jakulin.py @ 6538:a5f65d7f0b2c

Revision 6538:a5f65d7f0b2c, 25.1 KB checked in by Mitar <Mitar@…>, 4 years ago (diff)

Made XPM version of the icon 32x32.

Line 
1# ORANGE Logistic Regression
2#    by Alex Jakulin (jakulin@acm.org)
3#
4#       based on:
5#           Miller, A.J. (1992):
6#           Algorithm AS 274: Least squares routines to supplement
7#                those of Gentleman.  Appl. Statist., vol.41(2), 458-478.
8#
9#       and Alan Miller's F90 logistic regression code
10# # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
11#
12# CVS Status: $Id$
13#
14# Version 1.8 (14/05/2004)
15#   - Normalized calibration
16#
17# Version 1.7 (11/08/2002)
18#   - Support for new Orange (RandomIndices)
19#   - Assertion error was resulting because Orange incorrectly compared the
20#     attribute values as returned from .getclass() and the values array in
21#     the attribute definition. Cast both values to integer before comparison
22#   - Extended error messages.
23#   - Support for setting the preferred ordinal attribute transformation.
24#   - Prevent divide-by-0 errors in discriminant code.
25#
26# Version 1.6 (31/10/2001)
27#
28# The procedure used by this implementation of LR is minimization of log-likelihood
29# via deviance testing.
30#
31#
32# To Do:
33#   - Trivial imputation is used for missing attribute values. You should use
34#     something better. Eventually, this code might do something automatically.
35#
36#   - MarginMetaLearner is not doing anything when the classifier outputs a
37#     probability distribution. But even in such cases, the PD could be corrected.
38#
39#   - CalibrationMetaLearner  -> calibrates the *PD* estimates via CV
40
41import orange
42import orng2Array
43import orngCRS
44import math
45
46MAX_EXP = 50
47
48# BEWARE: these routines do not work with orange tables and are not orange-compatible
49class BLogisticLearner(orange.Learner):
50    def __init__(self, regularization = 0.0):
51        self.regul = regularization
52       
53    def getmodel(self, examples):
54      errors = ["LogReg: ngroups < 2, ndf < 0 -- not enough examples with so many attributes",
55                "LogReg: n[i]<0",
56                "LogReg: r[i]<0",
57                "LogReg: r[i]>n[i]",
58                "LogReg: constant variable",
59                "LogReg: singularity",
60                "LogReg: infinity in beta",
61                "LogReg: no convergence"]
62      model = orngCRS.LogReg(examples,self.regul)
63      errorno = model[8]
64      #print errors[errorno-1]
65      if errorno == 5 or errorno == 6:
66        # dependencies between variables, remove them
67        raise RedundanceException(model[9])
68      elif errorno == 1:
69        raise TooManyAttributes()
70      else:
71        if errorno != 0 and errorno != 7:
72            # unhandled exception (0=all ok, 7=perfect separation)
73            raise errors[errorno-1]
74      return (model,errorno)
75       
76    def __call__(self, examples):
77      (model,errorno) = self.getmodel(examples)
78      return BLogisticClassifier(model,examples)
79
80
81class BLogisticClassifier(orange.Classifier):
82    def __init__(self, model,examples):
83        (self.chisq,self.devnce,self.ndf,self.beta,
84        self.se_beta,self.fit,self.covbeta,
85        self.stdres,errorno,masking) = model
86
87        # set up the parameters for discrimination
88        sum = 1.0
89        for i in self.beta[1:]:
90            if abs(i) > 1e-6:
91                sum *= abs(i)
92        if sum > 1e100:
93            sum = max(self.beta[1:])
94        if sum < 1e-6:
95            sum = 1e-6
96        scale = 1.0/math.sqrt(sum)
97        self.nbeta = [x*scale for x in self.beta]
98        self.prior = 1.0/(len(examples)+2)
99        self.iprior = 1-(2*self.prior)
100               
101    def getmargin(self,example):
102        sum = self.nbeta[0]
103        for i in xrange(len(self.nbeta)-1):
104            sum = sum + example[i]*self.nbeta[i+1]
105        return sum
106
107    def description(self,attnames,classname,classv):
108        print 'Logistic Regression Report'
109        print '=========================='
110        print '\chi^2',self.chisq
111        print 'deviance',self.devnce
112        print 'NDF',self.ndf
113        print
114        print 'Base outcome:',classname[0],'=',classv
115        assert(len(attnames)==len(self.beta)-1)
116        print 'beta_0:',self.beta[0],'+-',self.se_beta[0]
117        for i in xrange(len(attnames)):
118            print attnames[i],self.beta[i+1],'+-',self.se_beta[i+1]
119
120    def geteffects(self,example):
121        # logistic regression
122        sum = self.beta[0]
123        for i in xrange(len(self.beta)-1):
124            sum = sum + example[i]*self.beta[i+1]
125        return sum
126
127    def __call__(self,example):
128        sum = self.geteffects(example)
129        # print sum, example
130        if sum > MAX_EXP:
131            return (1,1-self.prior)
132        elif sum < -MAX_EXP:
133            return (0,1-self.prior)
134        else:
135            sum = math.exp(sum)
136            p = sum/(1.0+sum) # probability that the class is 1
137            if p < 0.5:
138                return (0,1-self.prior-self.iprior*p)
139            else:
140                return (1,self.prior+self.iprior*p)
141
142class MajorityLogClassifier(orange.Classifier):
143    def __init__(self,ratio,lent):
144        self.beta = [math.log(ratio/(1-ratio))]
145        self.nbeta = [1]
146        self.ratio = ratio
147        self.prior = 1.0/(lent+2)
148        self.iprior = 1-(2*self.prior)
149        if ratio>0.5:
150            self.o = 1
151        else:
152            self.o = 0
153            self.ratio = 1-self.ratio
154        self.ratio = self.prior+self.iprior*self.ratio
155
156    def geteffects(self,example):
157        return self.beta[0]
158
159    def getmargin(self,example):
160        return 0
161       
162    def __call__(self,example):
163        return (self.o,self.ratio)
164       
165
166class TooManyAttributes:
167  def __init__(self):
168    pass
169
170  def __str__(self):
171    return "Too many variables."
172
173class RedundanceException:
174  def __init__(self,redundant_vars):
175    self.redundant_vars = redundant_vars
176
177  def __str__(self):
178    return "Logistic regression cannot work with constant or linearly dependent variables."
179
180
181
182#
183# Logistic regression throws an exception upon constant or linearly
184# dependent attributes. RobustBLogisticLearner remembers to ignore
185# such attributes.
186#
187# returns None, if all attributes singular
188#
189class RobustBLogisticLearner(BLogisticLearner):
190    def __init__(self, regularization=0.0):
191        self.regularization = regularization
192       
193    def __call__(self, examples, translate, importances, classfreq):
194        assert(len(classfreq)==2)
195        skipping = 0
196        na = len(examples[0])
197        mask = [0]*na
198        last_importance = 0
199        assert(na > 0)
200        # while there are any unmasked variables
201        blearner = BLogisticLearner(self.regularization)
202        while skipping < na-1: 
203            try:
204                if skipping != 0:
205                    # remove some variables
206                    data = []
207                    for ex in examples:
208                        maskv = []
209                        for i in xrange(len(mask)):
210                            if mask[i] == 0:
211                                maskv.append(ex[i])
212                        data.append(maskv)
213                else:
214                    data = examples
215                classifier = blearner(data)
216                return RobustBLogisticClassifierWrap(classifier,mask)
217            except RedundanceException, exp:
218                ext_offs = 0 # offset in the existing mask
219                for i in exp.redundant_vars:
220                    # skip what's already masked
221                    while mask[ext_offs] == 1:
222                        ext_offs += 1
223                    if i != 0:
224                        # new masking
225                        mask[ext_offs] = 1
226                        skipping += 1
227                    ext_offs += 1
228                #print 'r.skipping',skipping
229            except TooManyAttributes:
230                # remove something
231                removed = 0
232                while removed == 0 or skipping < na-1:
233                    tokill = importances[last_importance][1]
234                    i = translate.trans[tokill].idx
235                    while i < translate.trans[tokill].nidx:
236                        if mask[i] == 0: # if not deleted already
237                            removed += 1
238                            mask[i] = 1
239                            skipping += 1
240                        i += 1
241                    last_importance += 1
242                #print 'i.skipping',skipping
243        return RobustBLogisticClassifierWrap(MajorityLogClassifier(classfreq[1],len(examples)),mask)
244
245# this wrapper transforms the example
246#
247# it is a wrapper, because it has to work with both
248# the discriminant and the LR
249class RobustBLogisticClassifierWrap(orange.Classifier):
250    def __init__(self, classifier, mask):
251        self.classifier = classifier
252        self.mask = mask
253
254    def translate(self,example):
255        #assert(len(example) == len(self.mask) or len(example) == len(self.mask)-1) # note that for classification, the class isn't defined
256        maskv = []
257        for i in xrange(len(example)):
258            if self.mask[i] == 0:
259                maskv.append(example[i])
260        return maskv
261
262    def description(self,variablenames,n):
263        maskv = []
264        for i in xrange(len(variablenames[0])):
265            if self.mask[i] == 0:
266                maskv.append(variablenames[0][i])
267        self.classifier.description(maskv,variablenames[1],n)
268
269    def geteffects(self, example):
270        return self.classifier.geteffects(self.translate(example))
271
272    def getmargin(self, example):
273        return self.classifier.getmargin(self.translate(example))
274
275    def __call__(self, example):
276        return self.classifier(self.translate(example))
277
278
279#
280# Logistic regression works with arrays and not Orange domains
281# This wrapper performs the domain translation
282#
283class BasicLogisticLearner(RobustBLogisticLearner):
284    def __init__(self,regularization = 0.01):
285        self.translation_mode_d = 0 # dummy
286        self.translation_mode_c = 1 # standardize
287        self.regularization = regularization
288
289    def __call__(self, examples, weight = 0,fulldata=0):
290        if examples.domain.classVar.varType != 1:
291            raise "Logistic learner only works with discrete class."
292        translate = orng2Array.DomainTranslation(self.translation_mode_d,self.translation_mode_c)
293        if fulldata != 0:
294            translate.analyse(fulldata, weight, warning=0)
295        else:
296            translate.analyse(examples, weight, warning=0)
297        translate.prepareLR()
298        mdata = translate.transform(examples)
299
300        # get the attribute importances
301        t = examples
302        importance = []
303        for i in xrange(len(t.domain.attributes)):
304            qi = orange.MeasureAttribute_relief(t.domain.attributes[i],t)
305            importance.append((qi,i))
306        importance.sort()
307        freqs = list(orange.Distribution(examples.domain.classVar,examples))
308        s = 1.0/sum(freqs)
309        freqs = [x*s for x in freqs] # normalize
310
311        rl = RobustBLogisticLearner(regularization=self.regularization)
312        if len(examples.domain.classVar.values) > 2:
313            ## form several experiments:
314            # identify the most frequent class value
315            tfreqs = [(freqs[i],i) for i in xrange(len(freqs))]
316            tfreqs.sort()
317            base = tfreqs[-1][1] # the most frequent class
318            classifiers = []
319            for i in xrange(len(tfreqs)-1):
320                # edit the translation
321                alter = tfreqs[i][1]
322                cfreqs = [tfreqs[-1][0],tfreqs[i][0]] # 0=base,1=alternative
323                # edit all the examples
324                for j in xrange(len(mdata)):
325                    c = int(examples[j].getclass())
326                    if c==alter:
327                        mdata[j][-1] = 1
328                    else:
329                        mdata[j][-1] = 0
330                r = rl(mdata,translate,importance,cfreqs)
331                classifiers.append(r)
332            return ArrayLogisticClassifier(classifiers,translate,tfreqs,examples.domain.classVar,len(mdata))
333        else:
334            r = rl(mdata,translate,importance,freqs)
335            return BasicLogisticClassifier(r,translate)
336
337
338class ArrayLogisticClassifier(orange.Classifier):
339    def __init__(self, classifiers, translator, tfreqs, classVar, allex):
340        self.classifiers = classifiers
341        self.translator = translator
342        self.tfreqs = tfreqs
343        self.classVar = classVar
344        self.nc = len(classifiers)+1
345        self.prior = 1.0/(self.nc+allex)
346        self.iprior = 1.0-self.nc*self.prior
347        assert(self.nc == len(classVar.values) and self.nc == len(tfreqs))
348        self._name = 'Multiclass Logistic Classifier'
349
350    def description(self):
351        for x in self.classifiers:
352            x.description(self.translator.description(),self.translator.cv.attr.values[1])
353
354    def __call__(self, example, format = orange.GetValue):
355        tex = self.translator.extransform(example)
356        effects = []
357        for i in xrange(self.nc-1):
358            idx = self.tfreqs[i][1]
359            effect = self.classifiers[i].geteffects(tex)
360            if effect > MAX_EXP:
361                effect = MAX_EXP
362            elif effect < -MAX_EXP:
363                effect = -MAX_EXP
364            effects.append((idx,math.exp(effect)))
365        tfreqs = self.tfreqs
366
367        # aggregate the predictions
368        p = [0.0 for i in xrange(self.nc)]
369        sum = 0.0
370        for (idx,effect) in effects:
371            sum += effect
372        sum += 1.0
373        q = (self.iprior/sum)
374       
375        # prepare PDF
376        p = [0.0 for i in xrange(self.nc)]
377        psum = 0.0
378        maxp = -1
379        maxi = -1
380        for (idx,effect) in effects:
381            rr = self.prior + effect*q
382            p[idx] = rr
383            psum += rr
384            if rr >= maxp: # most likely class
385                maxi = idx
386                maxp = rr
387        p[tfreqs[-1][1]] = 1.0-psum # base class
388        if 1-psum >= maxp:
389            maxi = tfreqs[-1][1] # most likely class
390
391        v = self.classVar(maxi)
392        if format == orange.GetValue:
393            return v
394        if format == orange.GetBoth:
395            return (v,p)
396        if format == orange.GetProbabilities:
397            return p
398
399class BasicLogisticClassifier(orange.Classifier):
400    def __init__(self, classifier, translator):
401        self.classifier = classifier
402        self.translator = translator
403        self._name = 'Basic Logistic Classifier'       
404
405    def getmargin(self,example):
406        tex = self.translator.extransform(example)
407        r = self.classifier.getmargin(tex)
408        return r
409
410    def description(self):
411        self.classifier.description(self.translator.description(),self.translator.cv.attr.values[1])
412
413    def __call__(self, example, format = orange.GetValue):
414        tex = self.translator.extransform(example)
415        r = self.classifier(tex)
416        #print example, tex, r
417        v = self.translator.getClass(r[0])
418        p = [0.0,0.0]
419        for i in xrange(2):
420            if int(v) == i:
421                p[i] = r[1]
422                p[1-i] = 1-r[1]
423                break
424        if format == orange.GetValue:
425            return v
426        if format == orange.GetBoth:
427            return (v,p)
428        if format == orange.GetProbabilities:
429            return p
430
431
432#
433# A margin-based Bayesian learner
434#
435class BasicBayesLearner(orange.Learner):
436    def _safeRatio(self,a,b):
437        if a*10000.0 < b:
438            return -10
439        elif b*10000.0 < a:
440            return 10
441        else:
442            return math.log(a)-math.log(b)
443
444    def _process(self,classifier,examples):
445        # todo todo - support for loess
446        beta = self._safeRatio(classifier.distribution[1],classifier.distribution[0])
447        coeffs = []
448        for i in xrange(len(examples.domain.attributes)):
449            for j in xrange(len(examples.domain.attributes[i].values)):
450                p1 = classifier.conditionalDistributions[i][j][1]
451                p0 = classifier.conditionalDistributions[i][j][0]
452                coeffs.append(self._safeRatio(p1,p0)-beta)
453        return (beta, coeffs)
454
455   
456    def __init__(self):
457        self.translation_mode_d = 1 # binarization
458        self.translation_mode_c = 1 # standardization
459
460    def __call__(self, examples, weight = 0,fulldata=0):
461        if not(examples.domain.classVar.varType == 1 and len(examples.domain.classVar.values)==2):
462            raise "BasicBayes learner only works with binary discrete class."
463        for attr in examples.domain.attributes:
464            if not(attr.varType == 1):
465                raise "BasicBayes learner does not work with continuous attributes."
466        translate = orng2Array.DomainTranslation(self.translation_mode_d,self.translation_mode_c)
467        if fulldata != 0:
468            translate.analyse(fulldata, weight)
469        else:
470            translate.analyse(examples, weight)
471        translate.prepareLR()
472        (beta, coeffs) = self._process(orange.BayesLearner(examples), examples)
473        return BasicBayesClassifier(beta,coeffs,translate)
474
475
476class BasicBayesClassifier(orange.Classifier):
477    def __init__(self, beta, coeffs, translator):
478        self.beta = beta
479        self.coeffs = coeffs
480        self.translator = translator
481        self._name = 'Basic Bayes Classifier'       
482
483    def getmargin(self,example):
484        tex = self.translator.extransform(example)
485        sum = self.beta
486        for i in xrange(len(self.coeffs)):
487            sum += tex[i]*self.coeffs[i]
488        return -sum
489
490    def __call__(self, example, format = orange.GetValue):
491        sum = -self.getmargin(example)
492
493        # print sum, example
494        if sum > MAX_EXP:
495            r = (1,1.0)
496        elif sum < -MAX_EXP:
497            r = (0,1.0)
498        else:
499            sum = math.exp(sum)
500            p = sum/(1.0+sum) # probability that the class is 1
501            if p < 0.5:
502                r = (0,1-p)
503            else:
504                r = (1,p)
505
506        v = self.translator.getClass(r[0])
507        p = [0.0,0.0]
508        for i in xrange(2):
509            if int(v) == i:
510                p[i] = r[1]
511                p[1-i] = 1-r[1]
512                break
513        if format == orange.GetValue:
514            return v
515        if format == orange.GetBoth:
516            return (v,p)
517        if format == orange.GetProbabilities:
518            return p
519
520
521class BasicCalibrationLearner(orange.Learner):
522    def __init__(self, discr = orange.EntropyDiscretization(), learnr = orange.BayesLearner()):
523        self.disc = discr
524        self.learner = learnr
525
526    def __call__(self, examples):
527        if not(examples.domain.classVar.varType == 1 and len(examples.domain.classVar.values)==2):
528            raise "BasicCalibration learner only works with binary discrete class."
529        if len(examples.domain.attributes) > 1 or examples.domain.attributes[0].varType != 2:
530            raise "BasicCalibration learner only works with a single numerical attribute."
531
532        new_a = self.disc(examples.domain.attributes[0],examples)
533        data = examples.select([new_a, examples.domain.classVar])
534        c = self.learner(data)
535
536        return BasicCalibrationClassifier(c)
537
538class BasicCalibrationClassifier(orange.Classifier):
539    def __init__(self, classifier):
540        self.classifier = classifier
541
542    def __call__(self, example, format = orange.GetValue):
543        return self.classifier(example,format)
544       
545#
546# Margin Probability Wrap
547#
548# Margin metalearner attempts to use the margin-based classifiers, such as linear
549# discriminants and SVM to return the class probability distribution. Thie metalearner
550# only works with binary classes.
551#
552# Margin classifiers output the distance from the separating hyperplane, not
553# the probability distribution. However, the distance from the hyperplane can
554# be associated with the probability. This is a regression problem, for which
555# we can apply logistic regression.
556#
557# However, one must note that perfect separating hyperplanes generate trivial
558# class distributions.
559#
560class MarginMetaLearner(orange.Learner):
561    def __init__(self, learner, folds = 10, replications = 1, normalization=0, fulldata=0, metalearner = BasicLogisticLearner()):
562        self.learner = learner
563        self.folds = folds
564        self.metalearner = metalearner
565        self.replications = replications
566        self.normalization = normalization
567        self.fulldata = fulldata
568       
569    def __call__(self, examples, weight = 0):
570        if not(examples.domain.classVar.varType == 1 and len(examples.domain.classVar.values)==2):
571            # failing the assumptions of margin-metalearner...
572            return MarginMetaClassifierWrap(self.learner(examples))
573
574        mv = orange.FloatVariable(name="margin")
575        estdomain = orange.Domain([mv,examples.domain.classVar])
576        mistakes = orange.ExampleTable(estdomain)
577        if weight != 0:
578            mistakes.addMetaAttribute(1)
579
580        for replication in xrange(self.replications):
581            # perform 10 fold CV, and create a new dataset
582            try:
583                selection = orange.MakeRandomIndicesCV(examples, self.folds, stratified=0, randomGenerator = orange.globalRandom) # orange 2.2
584            except:
585                selection = orange.RandomIndicesCVGen(examples, self.folds) # orange 2.1
586            for fold in xrange(self.folds):
587              if self.folds != 1: # no folds
588                  learn_data = examples.selectref(selection, fold, negate=1)
589                  test_data  = examples.selectref(selection, fold)
590              else:
591                  learn_data = examples
592                  test_data  = examples
593
594              # fulldata removes the influence of scaling on the distance dispersion.                 
595              if weight!=0:
596                  if self.fulldata:
597                      classifier = self.learner(learn_data, weight=weight, fulldata=examples)
598                  else:
599                      classifier = self.learner(learn_data, weight=weight)
600              else:
601                  if self.fulldata:
602                      classifier = self.learner(learn_data, fulldata=examples)
603                  else:
604                      classifier = self.learner(learn_data)
605              # normalize the range
606              if self.normalization:
607                  mi = 1e100
608                  ma = -1e100
609                  for ex in learn_data:
610                      margin = classifier.getmargin(ex)
611                      mi = min(mi,margin)
612                      ma = max(ma,margin)
613                  coeff = 1.0/max(ma-mi,1e-16)
614              else:
615                  coeff = 1.0 
616              for ex in test_data:
617                  margin = coeff*classifier.getmargin(ex)
618                  if type(margin)==type(1.0) or type(margin)==type(1):
619                      # ignore those examples which are handled with
620                      # the actual probability distribution
621                      mistake = orange.Example(estdomain,[float(margin), ex.getclass()])
622                      if weight!=0:
623                          mistake.setmeta(ex.getMetaAttribute(weight),1)
624                      mistakes.append(mistake)
625
626        if len(mistakes) < 1:
627            # nothing to learn from
628            if weight == 0:
629                return self.learner(examples)
630            else:
631                return self.learner(examples,weight)
632        if weight != 0:
633            # learn a classifier to estimate the probabilities from margins
634            # learn a classifier for the whole training set
635            estimate = self.metalearner(mistakes, weight = 1)
636            classifier = self.learner(examples, weight)
637        else:
638            estimate = self.metalearner(mistakes)
639            classifier = self.learner(examples)
640
641        # normalize the range
642        if self.normalization:
643            mi = 1e100
644            ma = -1e100
645            for ex in examples:
646                margin = classifier.getmargin(ex)
647                mi = min(mi,margin)
648                ma = max(ma,margin)
649            coeff = 1.0/max(ma-mi,1e-16)
650        else:
651            coeff = 1.0
652        #print estimate.classifier.classifier
653        #for x in mistakes:
654        #    print x,estimate(x,orange.GetBoth)
655
656        return MarginMetaClassifier(classifier, estimate, examples.domain, estdomain, coeff)
657
658
659class MarginMetaClassifierWrap(orange.Classifier):
660    def __init__(self, classifier):
661        self.classifier = classifier
662        self._name = 'MarginMetaClassifier'
663
664    def __call__(self, example, format = orange.GetValue):
665    try:
666        (v,p) = self.classifier(example,orange.GetBoth)
667    except:
668        v = self.classifier(example)
669        if format == orange.GetValue:
670            return v
671        if format == orange.GetBoth:
672            return (v,p)
673        if format == orange.GetProbabilities:
674            return p
675
676class MarginMetaClassifier(orange.Classifier):
677    def __init__(self, classifier, estimator, domain, estdomain, coeff):
678        self.classifier = classifier
679        self.coeff = coeff
680        self.estimator = estimator
681        self.domain = domain
682        self.estdomain = estdomain
683        self.cv = self.estdomain.classVar(0)
684        self._name = 'MarginMetaClassifier'
685
686    def __call__(self, example, format = orange.GetValue):
687        r = self.coeff*self.classifier.getmargin(example)
688        if type(r) == type(1.0) or type(r) == type(1):
689            # got a margin
690            ex = orange.Example(self.estdomain,[r,self.cv]) # need a dummy class value
691            (v,p) = self.estimator(ex,orange.GetBoth)
692        else:
693            # got a probability distribution, which can happen with LR... easy
694            (v,p) = r
695
696        if format == orange.GetValue:
697            return v
698        if format == orange.GetBoth:
699            return (v,p)
700        if format == orange.GetProbabilities:
701            return p
702       
Note: See TracBrowser for help on using the repository browser.