Try Theano with Kaggle's MNIST Data ~ Logistic Regression ~

Overview

Theano tutorial (http://deeplearning.net/tutorial/contents.html) Based on, I implemented logistic regression with Kaggle's MNIST data.

Implementation

from __future__ import print_function

__docformat__ = 'restructedtext en'

import six.moves.cPickle as pickle
import timeit
import numpy 
import theano
import theano.tensor as T
import numpy as np

rng = numpy.random.RandomState(8000)

#Definition of shared variables
#The data created by this expression is

#-------------------------------------------------------------------------
#Logistic Regression class
#-------------------------------------------------------------------------
class LogisticRegression(object):
    """Multi-class Logistic Regression Class
    """
    def __init__(self, input, n_in, n_out):


        #MODEL
        #Shared variable representation
        #About W
        self.W = theano.shared(
            value=numpy.zeros(
                (n_in, n_out),
                dtype=theano.config.floatX
            ),
            name='W',
            borrow=True
        )
        #borrow determines whether the entity of data defined in Python space is also shared by shared variables.
        #About b
        self.b = theano.shared(
            value=numpy.zeros(
                (n_out,),
                dtype=theano.config.floatX
            ),
            name='b',
            borrow=True
        )


        #Calculation of softmax function (calculate the probability of being assigned to a class)
        #p(y=i|x,W,b)
        self.p_y_given_x = T.nnet.softmax(T.dot(input, self.W) + self.b)
        #Take argmax and output the predicted value
        #ypred
        self.y_pred = T.argmax(self.p_y_given_x, axis=1)
        #.Stored in paramsn
        self.params = [self.W, self.b]
        #.Store in input
        self.input = input

    #Loss function
    def negative_log_likelihood(self, y):
        #Definition of loss function
        #self.p_y_given_x is the softmax function (output yn)
        #T.arange(y.shape[0])Is the row and column of Numpy[]
        #-Ln(p(t|w)).PRML4.Type 90
        return -T.mean(T.log(self.p_y_given_x)[T.arange(y.shape[0]), y])

    #Error rate
    def errors(self, y):

        #Checking the dimensions of the output
        if y.ndim != self.y_pred.ndim:
            raise TypeError(
                'y should have the same shape as self.y_pred',
                ('y', y.type, 'y_pred', self.y_pred.type)
            )
        #Check if the data type is correct
        if y.dtype.startswith('int'):
            #neq outputs not equal output answer is correct and correct probability
            #First argument: argument with maximum output probability image
            #Second argument: Correct label
            return T.mean(T.neq(self.y_pred, y))
        else:
            raise NotImplementedError()

#Mini-batch stochastic gradient descent
def sgd_optimization_mnist(learning_rate=0.10, n_epochs=100,
                           dataset='mnist.pkl.gz',
                           batch_size=500):
                               
    dtype = np.float32
    data = np.loadtxt("../input/train.csv", dtype=dtype,
                           delimiter=',', skiprows=1)
    test_data = np.loadtxt("../input/test.csv", dtype=dtype,
                          delimiter=',', skiprows=1)
    test_datasets=[test_data]
    
    print (data.shape)

    labels = data[:,0]
    data = data[:, 1:]                               

    index = T.lscalar()  #minibatch index
    NUM_TRAIN = len(data)
    NUM_TEST = len(test_data)
    if NUM_TRAIN % batch_size != 0: 
        whole = (NUM_TRAIN // batch_size) * batch_size
        data = data[:whole]
        NUM_TRAIN = len(data) 
        
    #random sort
    indices = rng.permutation(NUM_TRAIN)
    data, labels = data[indices, :], labels[indices]
    # batch_size is 500, (480, 20)And train on 96 percent of the dataset and validate the rest
    is_train = numpy.array( ([0]* (batch_size - 20) + [1] * 20) * (NUM_TRAIN // batch_size))
    
    test_indices = rng.permutation(NUM_TEST)
    test_data, test_labels = test_data[test_indices, :], labels[test_indices]
    is_test = numpy.array( ([0]* (batch_size - 20) + [1] * 20) * (NUM_TEST // batch_size))
    # trai_set,valid_set,test_Set preparation
    train_set_x, train_set_y = numpy.array(data[is_train==0]), labels[is_train==0]
    valid_set_x, valid_set_y = numpy.array(data[is_train==1]), labels[is_train==1]
    test_set_x,test_set_y    = numpy.array(data[is_test==0]), labels[is_test==0]
    #minibatch calculation
    n_train_batches = len(train_set_y) // batch_size
    n_valid_batches = len(valid_set_y) // batch_size
    n_test_batches = len(test_set_y)   // batch_size
    
    ##############
    #Model building#
    ##############

    print('... building the model')

    #long type scalar
    index = T.lscalar()  

    #matrix type
    x = T.matrix('x')
    #int type vector
    y = T.ivector('y') 

    #Building logistic regression class
    #MNIST data is 28X28
    classifier = LogisticRegression(input=x, n_in=28 * 28, n_out=10)

    #cost function is negative log-likelihood
    cost = classifier.negative_log_likelihood(y)

    #Mini batch type theano.Summarize with function
    train_set_x = theano.shared(numpy.asarray(train_set_x, dtype=theano.config.floatX))
    train_set_y = T.cast(theano.shared(numpy.asarray(train_set_y, dtype=theano.config.floatX)), 'int32')
    test_set_x  = theano.shared(numpy.asarray(test_set_x, dtype=theano.config.floatX))
    test_set_y   = T.cast(theano.shared(numpy.asarray(test_set_y, dtype=theano.config.floatX)), 'int32')
    valid_set_x = theano.shared(numpy.asarray(valid_set_x, dtype=theano.config.floatX)) 
    valid_set_y = T.cast(theano.shared(numpy.asarray(valid_set_y, dtype=theano.config.floatX)), 'int32')

    #validate model
    validate_model = theano.function(
        inputs=[index],
        outputs=classifier.errors(y),
        givens={
            x: valid_set_x[index * batch_size: (index + 1) * batch_size],
            y: valid_set_y[index * batch_size: (index + 1) * batch_size]
        }
    )

    #test model
    test_model = theano.function(
        inputs=[index],
        outputs=classifier.errors(y),
        givens={
            x: valid_set_x[index * batch_size:(index + 1) * batch_size],
            y: valid_set_y[index * batch_size:(index + 1) * batch_size]
        }
    )

    #Differentiate cost by W
    g_W = T.grad(cost=cost, wrt=classifier.W)
    #Differentiate cost by b
    g_b = T.grad(cost=cost, wrt=classifier.b)
    #cost function is negative log-likelihood
    
    #Parameter (W,b) Update
    updates = [(classifier.W, classifier.W - learning_rate * g_W),
               (classifier.b, classifier.b - learning_rate * g_b)]

    #train model
    #Parameters defined in updates (W,b) Update
    train_model = theano.function(
        inputs=[index],
        outputs=cost,
        updates=updates,
        givens={
            x: train_set_x[index * batch_size: (index + 1) * batch_size],
            y: train_set_y[index * batch_size: (index + 1) * batch_size]
        }
    )

    #################
    #Model training#
    #################
    print('... training the model')
    #Maximum number of LOOPs
    patience = 50000  
    #If it doesn't increase more than 2
    patience_increase = 5  
    #Upper limit of grades
    improvement_threshold = 0.995         
    #Frequency of meeting patience. In this case, check with every epoch
    validation_frequency = min(n_train_batches, patience // 2)

    best_validation_loss = numpy.inf
    test_score = 0.
    start_time = timeit.default_timer()

    #The beginning of LOOP
    done_looping = False
    epoch = 0
    minibatch_index = 1
    while (epoch < n_epochs) and (not done_looping):
        epoch = epoch + 1
        for minibatch_index in range(n_train_batches):

            minibatch_avg_cost = train_model(minibatch_index)
            #Iteration calculation
            iter = (epoch - 1) * n_train_batches + minibatch_index

            if (iter + 1) % validation_frequency == 0:
                # zero-One loss calculated by Validation
                validation_losses = [validate_model(i)
                                     for i in range(n_valid_batches)]
                this_validation_loss = numpy.mean(validation_losses)

                print(
                    'epoch %i, minibatch %i/%i, validation error %f %%' %
                    (
                        epoch,
                        minibatch_index + 1,
                        n_train_batches,
                        this_validation_loss * 100.
                    )
                )

                #If the validation loss is the lowest
                if this_validation_loss < best_validation_loss:

                    if this_validation_loss < best_validation_loss *  \
                       improvement_threshold:
                        patience = max(patience, iter * patience_increase)

                    best_validation_loss = this_validation_loss
                    #Test with test set

                    test_losses = [test_model(i)
                                   for it in range(len(test_datasets))]
                    test_score = numpy.mean(test_losses)

                    print(
                        (
                            '     epoch %i, minibatch %i/%i, test error of'
                            ' best model %f %%'
                        ) %
                        (
                            epoch,
                            minibatch_index + 1,
                            n_train_batches,
                            test_score * 100.
                        )
                    )

                    #Save best model
                    with open('best_model.pkl', 'wb') as f:
                        pickle.dump(classifier, f)

            if patience <= iter:
                done_looping = True
                break

    end_time = timeit.default_timer()
    print(
        (
            'Optimization complete with best validation score of %f %%,'
            'with test performance %f %%'
        )
        % (best_validation_loss * 100., test_score * 100.)
    )
    print('The code run for %d epochs, with %f epochs/sec' % (
        epoch, 1. * epoch / (end_time - start_time)))
          
def predict():
    """
    An example of how to load a trained model and use it
    to predict labels.
    """
    #Call the trained model
    classifier = pickle.load(open('best_model.pkl'))

    #predict model
    predict_model = theano.function(
    inputs=[classifier.input],
        outputs=classifier.y_pred) 

    #Call the test set
    test_data = np.loadtxt("../input/test.csv", dtype=dtype,
                          delimiter=',', skiprows=1)
    #predict
    predicted_values = predict_model(test_data)
    print("Predicted values for the first 10 examples in test set:")
    print(predicted_values)

    np.savetxt('Submit_TheanoLR3.csv', np.c_[range(1, len(test_data) + 1), predicted_values], delimiter=',', comments = '', header = 'ImageId,Label', fmt='%d')

if __name__ == '__main__':    
    sgd_optimization_mnist()

result

The result is 0.90757. It was 0.90900 in scikit-learn, and it seems like this.

Recommended Posts

Try Theano with Kaggle's MNIST Data ~ Logistic Regression ~
Try regression with TensorFlow
Try TensorFlow MNIST with RNN
Implementing logistic regression with NumPy
Try data parallelism with Distributed TensorFlow
Logistic regression analysis Self-made with python
[Logistic regression] Implement k-validation with stats models
Logistic regression
Check raw data with Kaggle's Titanic (kaggle ⑥)
Try converting to tidy data with pandas
Try to aggregate doujin music data with pandas
Implement a discrete-time logistic regression model with stan
Logistic regression implementation with particle swarm optimization method
[Logistic regression] Implement holdout verification with stats models
Points to note when performing logistic regression with Statsmodels
Try MNIST with VAT (Virtual Adversarial Training) in Keras
Challenge image classification with TensorFlow2 + Keras 3 ~ Visualize MNIST data ~
Solving the iris problem with scikit-learn ver1.0 (logistic regression)
Predicting Kaggle's Hello World, Titanic Survivors with Logistic Regression-Modeling-