Classification

Notebook version: 1.1 (Oct 25, 2017)

Authors: Jesús Cid Sueiro (jcid@tsc.uc3m.es)

Changes: v.1.0 - First version. Python version
         v.1.1 - Updated to sklearn.model_selection. Python 3 compatibility.



In [1]:

    
# Import some libraries that will be necessary for working with data and displaying plots

# To visualize plots in the notebook
%matplotlib inline

#import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import scipy.io       # To read matlab files
from sklearn.preprocessing import PolynomialFeatures
from sklearn import svm

from sklearn import model_selection

import pylab
pylab.rcParams['figure.figsize'] = 9, 7

1. Introduction

In this notebook we will analyze the behavior of logistic regression and support vector machines on the dataset in file Dataset2D.mat. We first load the dataset.



In [2]:

    
# Load dataset
matvar = scipy.io.loadmat('Dataset2D.mat')
Xtrain = matvar['xTrain']
Xtest = matvar['xTest']
Xval = matvar['xVal']

# We must use astype(int) to convert the original target values (which are unsigned integers) to int.
Ytrain = matvar['yTrain'].astype(int)
Ytest = matvar['yTest'].astype(int)
Yval = matvar['yVal'].astype(int)

1.1 Data Preparation.

Normalize the dataset. Remind that the same transformation must be applied to training, validation and test data. Store train, validation and test input data in variables Xtrain, Xval and Xtest



In [3]:

    
# <SOL>
# </SOL>

# Check normalization
print(np.mean(Xtrain, axis=0))
print(np.mean(Xval, axis=0))
print(np.mean(Xtest, axis=0))
print(np.std(Xtrain, axis=0))
print(np.std(Xval, axis=0))
print(np.std(Xtest, axis=0))









    



The number of training samples is 400
The number of validation samples is 400
The number of test samples is 400
The data dimension is 2
[-3.24185123e-16  2.84217094e-16]
[-0.00966556  0.03965134]
[0.00125431 0.02832642]
[1. 1.]
[0.98064167 1.02630613]
[0.9809003  1.03485737]

Visualize the input variables from the training set in a 2-dimensional plot.



In [4]:

    
# Data visualization. This works for dimension 2 only.
if Xtrain.shape[1]==2:
    plt.scatter(Xtrain[:, 0], Xtrain[:, 1], c=Ytrain.flatten(), s=50, cmap='copper')
    plt.xlabel("$x_0$", fontsize=14)
    plt.ylabel("$x_1$", fontsize=14)
plt.show()

2. Linear Classification with Logistic Regression.

First we will analyze the behavior of logistic regression for this dataset.

2.1. MAP estimator.

Implement a function to compute the MAP estimate of the parameters of a linear logistic regression model with Gaussian prior and a given value of the inverse regularization parameter $C$. The method should return the estimated parameter and the negative log-likelihood, $\text{NLL}({\bf w})$. The sintaxis must be w, NLL = logregFitR(Z_tr, Y_tr, rho, C, n_it) where

Z_tr is the input training data matrix (one instance per row)
Y_tr contains the labels of corresponding to each row in the data matrix
rho is the learning step
C is the inverse regularizer
n_it is the number of iterations



In [5]:

    
# <SOL>
#</SOL>

2.2 Log-likelihood

Compute the MAP estimate for a polynomial regression with degree 5, for $C$ ranging from -0.01 to 100. Sample $C$ uniformly in a log scale, an plot using plt.semilogx.

Plot the final value of $\text{NLL}$ as a function of $C$. Can you explain the qualitative behavior of $\text{NLL}$ as $C$ grows?

The plot may show some oscillation because of the random noise introduced by random initializations of the learning algoritm. In order to smooth the results, you can initialize the random seed right before calling the logregFitR method, using

np.random.seed(3)



In [6]:

    
# <SOL>
# </SOL>









    



As C grows, the regularization effect dissapears, and the fit method minimizes NLL.



In [7]:

    
# This is a plot for the last value of C used in the code above.

if Xtrain.shape[1]==2:

    # Create a regtangular grid.
    x_min, x_max = Xtrain[:, 0].min(), Xtrain[:, 0].max() 
    y_min, y_max = Xtrain[:, 1].min(), Xtrain[:, 1].max()
    dx = x_max - x_min
    dy = y_max - y_min
    h = dy /400
    xx, yy = np.meshgrid(np.arange(x_min - 0.1 * dx, x_max + 0.1 * dx, h),
                         np.arange(y_min - 0.1 * dx, y_max + 0.1 * dy, h))
    X_grid = np.array([xx.ravel(), yy.ravel()]).T

    # Compute Z_grid
    Z_grid = poly.fit_transform(X_grid)
    n_grid = Z_grid.shape[0]
    Zn, mz, sz = normalize(Z_grid[:,1:], mz, sz)
    Z_grid = np.concatenate((np.ones((n_grid,1)), Zn), axis=1)

    # Compute the classifier output for all samples in the grid.
    pp, dd = logregPredict(Z_grid, w)
    pp = pp.reshape(xx.shape)

    # Paint output maps
    plt.figure()
    pylab.rcParams['figure.figsize'] = 8, 4  # Set figure size
    for i in [1, 2]:
        ax = plt.subplot(1,2,i)
        ax.set_xlabel('$x_0$')
        ax.set_ylabel('$x_1$')
        ax.axis('equal')
        if i==1:
            ax.contourf(xx, yy, pp, cmap=plt.cm.copper)
        else:
            ax.contourf(xx, yy, np.round(pp), cmap=plt.cm.copper)
        ax.scatter(Xtrain[:, 0], Xtrain[:, 1], c=Ytrain.flatten(), s=4, cmap='summer')
        
    plt.show()

2.3. Training and test errors.

Plot the training and validation error rates as a function of $C$. Compute the value of $C$ minimizing the validation error rate.



In [8]:

    
# <SOL>
# </SOL>









    



The optimal value of C is 0.22773062069049685

3. Non-linear classification with Support Vector Machines

In this section we will train an SVM with Gaussian kernels. In this case, we will select parameter $C$ of the SVM by cross-validation.

3.1. Dataset preparation.

Join the training and validation datasets in a single input matrix X_tr2 and a single label vector Y_tr2



In [9]:

    
# <SOL>
# </SOL>

3.2. Cross validated error estimate

Apply a 10-fold cross validation procedure to estimate the average error rate of the SVM for $C=1$ and $\gamma$ (which is the kernel width) equal to 5.



In [10]:

    
# <SOL>
# </SOL>









    



The average error rate is 0.24225

3.3. Influence of $C$.

Repeate exercise 3.2 for $\gamma=5$ and different values of $C$, ranging from $10^{-3}$ to $10^{4}$, obtained by uniform sampling in a logarithmic scale. Plot the average number of errors as function of $C$.

Note that fitting the SVM may take some time, specially for the largest values of $C$.



In [11]:

    
# <SOL>
# </SOL>









    



C = 0.001
C = 0.005994842503189409
C = 0.03593813663804628
C = 0.21544346900318845
C = 1.291549665014884
C = 7.742636826811269
C = 46.41588833612782
C = 278.2559402207126
C = 1668.100537200059
C = 10000.0

3.3. Hyperparameter optimization.

Compute the value of $C$ minimizing the validation error rate.



In [12]:

    
# <SOL>
# </SOL>









    



The optimal value of C in the explored range is 0.21544346900318845

3.4. Test error

Evaluate the classifier performance using the test data, for the selected hyperparameter values.



In [13]:

    
# <SOL>
# </SOL>









    



The test error for the selected model is 0.20499999999999996