In [1]:
print(__doc__)
In [2]:
import numpy as np
from bokeh.plotting import HBox, VBox, figure, show, output_file, GridPlot
from bokeh.models.mappers import LinearColorMapper
from bokeh.models import BasicTicker, Grid
from sklearn.datasets import fetch_olivetti_faces
from sklearn.utils.validation import check_random_state
from sklearn.ensemble import ExtraTreesRegressor
from sklearn.neighbors import KNeighborsRegressor
from sklearn.linear_model import LinearRegression
from sklearn.linear_model import Ridge
from sklearn.cross_validation import train_test_split
from sklearn.covariance import GraphLassoCV, ledoit_wolf
from sklearn.grid_search import GridSearchCV
from scipy.spatial import distance
import sklearn
import OVFM.Model as md
import OVFM.FeatureMap as fm
import OVFM.Risk as rsk
import OVFM.LearningRate as lr
import OVFM.DataGeneration as dg
import OVFM.SGD as sgd
import time
import sys
In [3]:
def tuneSGD( D, L, features, targets ):
best_score = 1
last_score = 1
best_params = np.empty( 3 )
for gamma in np.logspace( -5, -1, 5 ):
phi_l = fm.DecomposableFF( gamma, X_train.shape[ 1 ], D, L )
loss = ls.L2( )
for C in np.logspace( -10, 1, 5 ):
re = reg.L1( C )
for eta0 in np.logspace( -2, 2, 5 ):
l = lr.Constant( eta0 )
estimator = sgd.SGD(phi_l, l, loss, re, 1, 5, False, 0)
scores = sgd.scoreCV( estimator, features, targets )
print scores, gamma, C, eta0
if np.mean( scores ) >= last_score:
last_score = 1
break
last_score = np.mean( scores )
if np.mean( scores ) <= best_score:
best_params[ 0 ] = C
best_params[ 1 ] = eta0
best_params[ 2 ] = gamma
best_score = np.mean( scores )
return best_params
In [4]:
# Load the faces datasets
data = fetch_olivetti_faces()
targets = data.target
data = data.images.reshape((len(data.images), -1)).astype( float )
scaler = sklearn.preprocessing.StandardScaler( )
train, test = train_test_split( np.arange( 40 ).reshape( ( 40, 1 ) ), test_size = 0.25 )
train = data[ np.in1d( targets, train.ravel( ) ), : ]
test = data[ np.in1d( targets, test .ravel( ) ), : ] # Test on independent people
In [5]:
# Test on a subset of people
# n_faces = 10
# rng = check_random_state(4)
# face_ids = rng.randint(test.shape[0], size=(n_faces, ))
# test = test[face_ids, :]
n_pixels = data.shape[1]
X_train = scaler.fit_transform( train[:, :np.ceil(0.5 * n_pixels)] ) # Upper half of the faces
y_train = train[:, np.floor(0.5 * n_pixels):] # Lower half of the faces
X_test = scaler.fit_transform( test[:, :np.ceil(0.5 * n_pixels)] )
y_test = test[:, np.floor(0.5 * n_pixels):]
In [6]:
def double2rgba( M ):
low, high = M.min( ), M.max( )
ratio = ( M - low ) / ( high - low )
r = np.asarray( 255 * ratio, np.uint32 )
g = np.asarray( 255 * 0, np.uint32 ) << 8
b = np.asarray( 255 * 0, np.uint32 ) << 16
a = 255 << 24
return r + g + b + a
def double2gray( M ):
low, high = M.min( ), M.max( )
ratio = ( M - low ) / ( high - low )
r = np.asarray( 255 * ratio, np.uint32 )
g = np.asarray( 255 * ratio, np.uint32 ) << 8
b = np.asarray( 255 * ratio, np.uint32 ) << 16
a = 255 << 24
return r + g + b + a
In [7]:
# model = GraphLassoCV( )
# model.fit( y_train )
# M = model.covariance_
# M = model.precision_
# M = np.exp( - 100. * distance.squareform( distance.pdist( y_train.T, 'sqeuclidean' ) ) )
# alpha = 0.9995
# D = alpha * np.eye( M.shape[ 0 ] ) + ( 1- alpha ) * np.ones( M.shape ) #+ ( 1 - alpha ) * np.eye( M.shape[ 0 ] ) * ( M.sum( axis = 1 ).reshape( ( M.shape[ 0 ], 1 ) ) + M )
# L = np.linalg.pinv( D - M )
L = np.eye( y_train.shape[ 1 ] )
In [8]:
#Begin OVFM stuff
D = 250
gamma = 0.00005
C = 1e-4
eta0 = 1.0
# gff = fm.GaussianFF( gamma, X_train.shape[ 1 ], D )
# Kex = gff.kernel_exact( X_train )
# Kap = gff.kernel_approx( X_train )
# output_file( "kernel_approx.html", title="kernel_approx" )
# p1 = figure( title = 'Ground truth', x_range = [ 0, X_train.shape[ 0 ] ], y_range = [ 0, X_train.shape[ 0 ] ] )
# p1.image_rgba( image = [ double2rgba( Kex ) ] , x = [ 0 ], y = [ 0 ], dw = [ X_train.shape[ 0 ] ], dh = [ X_train.shape[ 0 ] ] )
# p2 = figure( title = 'Approximation', x_range = [ 0, X_train.shape[ 0 ] ], y_range = [ 0, X_train.shape[ 0 ] ] )
# p2.image_rgba( image = [ double2rgba( Kap ) ] , x = [ 0 ], y = [ 0 ], dw = [ X_train.shape[ 0 ] ], dh = [ X_train.shape[ 0 ] ] )
# print 'Kernel approximation MSE:', np.max( np.abs( Kex - Kap ) )
# alpha = 1 - 1.0 / X_train.shape[ 1 ]
# L = alpha * np.eye( ( X_train.shape[ 1 ] ) ) + ( 1 - alpha ) * np.ones( ( X_train.shape[ 1 ], X_train.shape[ 1 ] ) )
# L = np.linalg.pinv( np.corrcoef( y_train.T ) )
# L = L / np.linalg.norm( L, 2 )
# best_params = tuneSGD( 300, L, np.vstack( ( X_train, X_test ) ), np.vstack( ( y_train, y_test ) ) )
# C = best_params[ 0 ]
# eta0 = best_params[ 1 ]
# gamma = best_params[ 2 ]
# fs = fm.DecomposableFF( gamma, X_train.shape[ 1 ], D, np.eye( 1 ) )
# X_map = fs.map( X_train )
# Xt_map = fs.map( X_test )
model = md.Model( fm.DecomposableFF( gamma, X_train.shape[ 1 ], D, L ) )
risk = rsk.Ridge( C, 0 )
lc = lr.Constant( 1. * eta0 )
lb = lr.Constant( 0.01 * eta0 )
print L.shape
#End OVFM stuff
In [9]:
# Fit estimators
ESTIMATORS = {
"Extra trees": ExtraTreesRegressor(n_estimators=100, max_features=64,
random_state=0),
"K-nn": KNeighborsRegressor( ),
"Linear regression": LinearRegression(),
"Ridge": Ridge( fit_intercept = True ),
"OVFM": sgd.SGD( risk, 1.0, lc, lb, 10, 0 )
}
y_test_predict = dict()
for name, estimator in ESTIMATORS.items():
if name == 'Extra trees':
estimator.fit( X_train, y_train )
y_test_predict[ name ] = estimator.predict( X_test )
elif name == 'K-nn':
gcvr = GridSearchCV( estimator, { 'n_neighbors': [ 1, 5, 10, 25, 50, 100, 200 ] }, cv = 5 )
gcvr.fit( X_train, y_train )
estimator.set_params( n_neighbors = gcvr.best_params_[ 'n_neighbors' ] )
t = time.clock( )
estimator.fit( X_train, y_train )
y_test_predict[ name ] = estimator.predict( X_test )
elif name == 'Linear regression':
t = time.clock( )
estimator.fit( X_train, y_train )
y_test_predict[ name ] = estimator.predict( X_test )
elif name == 'Ridge':
gcvr = GridSearchCV( estimator, { 'alpha': [ 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7 ] }, cv = 5 )
gcvr.fit( X_train, y_train )
estimator.set_params( alpha = gcvr.best_params_[ 'alpha' ] )
t = time.clock( )
estimator.fit( X_train, y_train )
y_test_predict[ name ] = estimator.predict( X_test )
elif name == 'OVFM':
# for fold in folds:
model.reset( )
t = time.clock( )
estimator.fit( model, X_train, y_train )
y_test_predict[ name ] = model( X_test )
pass
print 'Learning time: ', name, ' ', time.clock( ) - t, np.mean( ( y_test_predict[ name ] - y_test ) ** 2 )
sys.stdout.flush()
In [ ]:
p = []
for i in range( 100 ):
true_face = np.hstack( ( scaler.inverse_transform( X_test[ i ] ), y_test[ i ] ) ).reshape( ( 64, 64 ) )[::-1,:]
if i == 0:
title = 'true face'
else:
title = ''
p.append( [ figure( title = title,
plot_width = 256, plot_height = 256,
min_border = 0, h_symmetry = False, v_symmetry = False,
x_range = [ 0, true_face.shape[ 0 ] ], y_range = [ 0, true_face.shape[ 1 ] ],
tools = "save" ) ] )
p[ i ][ 0 ].image_rgba( image = [ double2gray( true_face ) ] , x = [ 0 ], y = [ 0 ], dw = [ true_face.shape[ 0 ] ], dh = [ true_face.shape[ 1 ] ] )
for j, est in enumerate( sorted( ESTIMATORS ) ):
completed_face = np.hstack( ( scaler.inverse_transform( X_test[ i ] ), y_test_predict[ est ][ i ] ) ).reshape( ( 64, 64 ) )[::-1,:]
if i == 0:
title = est
else:
title = ''
p[ i ].append( figure( title = title,
plot_width = 256, plot_height = 256,
min_border = 0, h_symmetry = False, v_symmetry = False,
x_range = [ 0, true_face.shape[ 0 ] ], y_range = [ 0, true_face.shape[ 1 ] ],
tools = "save" ) )
p[ i ][ j + 1 ].image_rgba( image = [ double2gray( completed_face ) ] , x = [ 0 ], y = [ 0 ], dw = [ true_face.shape[ 0 ] ], dh = [ true_face.shape[ 1 ] ] )
In [ ]:
g = GridPlot( children = p )
show( VBox( HBox( p1, p2 ), g ) )
In [ ]: