Tomorrow: Add a third zoom in. Improve again? (Yes! 70%@0, 90% @ .002

Expand + examples to include 2 mod 3 fields. Double - examples accordingly

Best Zoom2Extractor with Daisy: ROC = 'max_features': 'auto', 'min_samples_split': 4, 'n_jobs': 2, 'criterion': 'infogain', 'n_estimators': 800


In [1]:
%pylab

import json
import random
import itertools

import numpy as np
from scipy import stats
from sklearn.metrics import auc_score
from sklearn.ensemble import GradientBoostingClassifier

from bubbly.model import Model
import bubbly.extractors as ext
from bubbly.util import rfp_curve
from bubbly.dr1 import WideLocationGenerator, highest_quality_on_params
from bubbly.hyperopt import rf_space

from skimage.feature import daisy


Welcome to pylab, a matplotlib-based Python environment [backend: module://IPython.kernel.zmq.pylab.backend_inline].
For more information, type 'help(pylab)'.
/Users/beaumont/Library/Python/2.7/lib/python/site-packages/scikits/__init__.py:1: UserWarning: Module argparse was already imported from /opt/local/Library/Frameworks/Python.framework/Versions/2.7/lib/python2.7/argparse.pyc, but /Users/beaumont/Library/Python/2.7/lib/python/site-packages is being added to sys.path
  __import__('pkg_resources').declare_namespace(__name__)

In [22]:
class DaisyExtractor(ext.Extractor):    
    
    def _extract_rgb(self, rgb):
        kwargs = dict(step=rgb.shape[0]/5, radius=rgb.shape[0] / 10, rings=2,
                      histograms=6, orientations=8)
        return np.hstack(daisy(rgb[:, :, i], **kwargs).ravel() for i in [0, 1])
        
class DaisySuper(ext.CompositeExtractor):
    composite_classes = [ext.RingExtractor, ext.MultiWaveletExtractor, ext.CompressionExtractor,
                         DaisyExtractor]


class ZoomedExtractor(ext.Extractor):
    def __init__(self, orig):
        self.orig = orig
        
    def extract(self, lon, l, b, r):
        return np.hstack((self.orig.extract(lon, l, b, r), 
                          self.orig.extract(lon, l, b, r * 2)))
    
class Zoom2Extractor(ext.Extractor):
    def __init__(self, orig):
        self.orig = orig
        
    def extract(self, lon, l, b, r):
        return np.hstack((self.orig.extract(lon, l, b, r), 
                          self.orig.extract(lon, l, b, r / 2),
                          self.orig.extract(lon, l, b, r * 2),
                          self.orig.extract(lon, l, b + r / 2, r)))

In [244]:
data = json.load(open('../models/benchmark_training_data.json'))

lg = WideLocationGenerator(mod3 = 1)
data['pos'] = filter(lambda x: lg._longitude_check(x[0]), highest_quality_on_params())

for k, v in data.items():
    data[k] = sorted(v)
    
ex = Zoom2Extractor(DaisySuper())
ex.shp = (60, 60)

In [245]:
from sklearn.utils import resample

def _xy(ex, on, off):
    x = np.vstack(ex.extract(*o).ravel().astype(np.float32) for o in on + off)
    x = np.nan_to_num(x)
    y = np.hstack((np.ones(len(on), dtype=np.int), np.zeros(len(off), dtype=np.int)))
    return x, y

npos = len(data['pos'])
xtrain, ytrain = _xy(ex, data['pos'], random.sample(data['neg'], npos))
xvalidate, yvalidate = _xy(ex, data['cv_pos'], data['cv_neg'])

In [252]:
from PyWiseRF import WiseRF
import os
os.environ['WISERF_ROOT'] = '/Users/beaumont/WiseRF-1.5.9-macosx-x86_64-rc1'

from sklearn.metrics import auc_score, auc
import sklearn.metrics as met


def decision_function(self, x):
    p = self.predict_proba(x)
    return p[:, 1] - p[:, 0]

WiseRF.decision_function = decision_function

def auc_below_fpos(y, yp, fpos):
    fp, tp, th = met.roc_curve(y, yp)
    good = (fp <= fpos)
    return auc(fp[good], tp[good])


class Choice(object):
    def __init__(self, *choices):
        self._choices = choices
        
    def rvs(self):
        return random.choice(self._choices)
    
class Space(object):
    def __init__(self, **hyperparams):
        self._hyperparams = hyperparams
        
    def __iter__(self):
        while True:
            yield {k: v.rvs() for k, v in self._hyperparams.items()}

            
gb_space = Space(learning_rate = stats.uniform(1e-3, 1 - 1.01e-3),
                 n_estimators = Choice(50, 100, 200),
                 max_depth = Choice(1, 2, 3),
                 subsample = stats.uniform(1e-3, 1 - 1.01e-3))

rf_space = Space(n_estimators = Choice(200, 400, 800, 1600),
                 min_samples_split = Choice(1, 2, 4),                 
                 criterion = Choice('gini', 'gainratio', 'infogain'),
                 max_features = Choice('auto'),
                 n_jobs = Choice(2))
                 
                 
def gb_objective(**params):
    clf = GradientBoostingClassifier(**params)
    clf.fit(xtrain, ytrain)
    
    df = clf.decision_function(xvalidate).ravel()
    return (df[yvalidate == 1] <= df[yvalidate == 0].max()).mean(), clf

def rf_objective(**params):
    clf = WiseRF(**params)
    clf.fit(xtrain, ytrain)

    df = clf.decision_function(xvalidate).ravel()
    return -auc_below_fpos(yvalidate, df, .0005), clf
    
    
def fmin(objective, space, threshold=np.inf):
    best = threshold
    
    try:        
        for p in space:
            f, clf = objective(**p)
            if f < best:
                best = f
                yield best, p, clf
    except KeyboardInterrupt:
        pass

In [253]:
import sys
from bubbly.util import roc_curve

for best, best_params, clf in fmin(rf_objective, rf_space, threshold):
    print best, best_params
    sys.stdout.flush()

    roc_curve(yvalidate, clf.decision_function(xvalidate), label='val')
    roc_curve(ytrain, clf.decision_function(xtrain), label='train')
    xlim(0, .002)
    legend(loc='lower right')
    show()


-0.000282484076433 {'max_features': 'auto', 'min_samples_split': 1, 'n_jobs': 2, 'criterion': 'gainratio', 'n_estimators': 200}
-0.000380573248408 {'max_features': 'auto', 'min_samples_split': 2, 'n_jobs': 2, 'criterion': 'gini', 'n_estimators': 1600}
-0.000385987261146 {'max_features': 'auto', 'min_samples_split': 4, 'n_jobs': 2, 'criterion': 'gini', 'n_estimators': 1600}

In [145]:
from skimage.util.montage import montage2d

def montage(arrs):
    print "montaging %i images" % len(arrs)
    print 'image dim:', arrs[0].shape
    r, g, b = tuple(montage2d(np.array([a[:, :, i] for a in arrs])) 
                    for i in range(3))
    return np.dstack((r, g, b)).astype(np.uint8)

def ex(*params):
    p = list(params)
    r1 = rgb.extract(*p)
    p[-1] *= 2.5
    r2 = rgb.extract(*p)
    return np.hstack((r1, r2))

rgb = ext.RGBExtractor()
rgb.shp = (80, 80)

df = clf.decision_function(xvalidate).ravel()

on_ind = np.argsort(df[yvalidate == 1])
off_ind = np.argsort(df[yvalidate == 0])[::-1]

figure(figsize=(10, 10))
im = montage([ex(*data['cv_pos'][i]) for i in on_ind[:49]])
imshow(im)
title("Hard Positives")
show()

figure(figsize=(10, 10))
im = montage([ex(*data['cv_neg'][i]) for i in off_ind[:49]])

imshow(im)
title("Hard Negatives")
show()


montaging 49 images
image dim: (80, 160, 3)
montaging 49 images
image dim: (80, 160, 3)

In [190]:
xtrain2, ytrain2 = _xy(ex, data['pos'], random.sample(data['neg'], npos))

In [246]:
xtrain_all, ytrain_all = _xy(ex, data['pos'], data['neg'])

In [200]:
class BaggedEstimator(object):
    def __init__(self, estimators):
        self.estimators = estimators
        
    def decision_function(self, X):
        return reduce(np.add, (e.decision_function(X) for e in self.estimators))
    
    def predict(self, X):
        return self.decision_function(X) > 0

In [247]:
from sklearn.utils import resample

def balanced_resample(x, y):
    npos = (y == 1).sum()
    xn = x[y == 0]
    xresample = resample(x, n_samples = npos)
    return np.vstack((x[y == 1], xresample)), np.hstack((np.ones(npos, dtype=np.int), np.zeros(npos, dtype=np.int)))

In [248]:
from sklearn.base import clone
clfs = []
for i in range(5):
    x, y = balanced_resample(xtrain_all, ytrain_all)
    clfs.append(clone(clf).fit(x, y))
    
bag = BaggedEstimator(clfs)

In [249]:
roc_curve(yvalidate, bag.decision_function(xvalidate), label='val')
roc_curve(ytrain, bag.decision_function(xtrain), label='train')
xlim(0, .002)
legend(loc='lower right')
show()



In [216]:
from skimage.util.montage import montage2d

def montage(arrs):
    print "montaging %i images" % len(arrs)
    print 'image dim:', arrs[0].shape
    r, g, b = tuple(montage2d(np.array([a[:, :, i] for a in arrs])) 
                    for i in range(3))
    return np.dstack((r, g, b)).astype(np.uint8)

def ex(*params):
    p = list(params)
    r1 = rgb.extract(*p)
    p[-1] *= 2.5
    r2 = rgb.extract(*p)
    return np.hstack((r1, r2))

rgb = ext.RGBExtractor()
rgb.shp = (80, 80)

df = bag.decision_function(xvalidate).ravel()

on_ind = np.argsort(df[yvalidate == 1])
off_ind = np.argsort(df[yvalidate == 0])[::-1]

figure(figsize=(15, 10))
im = montage([ex(*data['cv_pos'][i]) for i in on_ind[:49]])
imshow(im, origin='upper')
title("Hard Positives")
show()

figure(figsize=(15, 10))
im = montage([ex(*data['cv_neg'][i]) for i in off_ind[:49]])

imshow(im, origin='upper')
title("Hard Negatives")
show()


montaging 49 images
image dim: (80, 160, 3)