1
$\begingroup$

I've implemented a neural network and am using numerical gradient checking to validate the back-propagation algorithm is working correctly.

I'm using the standard method to calculate the numerical gradient:

$f'(x) \approx \frac{J(\theta + \epsilon) - J(\theta - \epsilon)}{2 \epsilon}$ where $\epsilon = 10^{-4}$, and

norm(gradients - numericalGradients)/norm(gradients + numericalGradients) < 10e-8

to check back-propagation is operating correctly.

I'm currently checking throughout training for testing purposes and I think that it is functioning correctly, however, as the network approaches its solution the gradients become very small and I think this is why the check fails.

So my question is what are the best practices for when to perform gradient checking?

Once at the start only, only when the minimum gradients are N times bigger than epsilon?

enter image description here

#!/usr/bin/python

import numpy as np
import matplotlib.pyplot as plt

class NeuralNet:

    GRAD_CHECK_THRESH = 10e-8
    GRAD_CHECK_ITER = 9999999
    EPSILSON = 10e-4

    def __init__(self, _maxIter=250, _nHidden=3):
        # Set hyperparameters
        self.maxIter = _maxIter
        self.nHidden = _nHidden
        self.learningRate = np.linspace(0.5, 0.05, self.maxIter)
        self.momentum = 0.5
        self.enNesterov = False  # Nesterov accelerated gradient (https://arxiv.org/pdf/1212.0901v2.pdf)
        self.cost = []
        self.gradDiff = []
        self.minGrad = []
        self.dW1 = 0
        self.dW2 = 0

    def initNet(self, _nInput, _nOutput):
        self.nInput = _nInput
        self.nOutput = _nOutput
        self.W1 = np.random.rand(self.nInput, self.nHidden)/100
        self.W2 = np.random.rand(self.nHidden, self.nOutput)/100

    def train(self, X, y):
        # Initialise network structure
        self.initNet(X.shape[1], y.shape[1])

        # Train - full batch
        for m in range(self.maxIter):

            # Check gradient descent is operating as expected
            if m < self.GRAD_CHECK_ITER:
                self.checkGradients(X, y)

            # Train net
            yHat = self.feedforward(X)
            J = self.costFunction(y, yHat)
            self.cost.append(J)
            dJdW1, dJdW2 = self.backprop(X, y, yHat)
            self.updateWeights(dJdW1, dJdW2, self.learningRate[m])

    def feedforward(self, X):
        self.z2 = np.dot(X, self.W1)
        self.a2 = self.sigmoid(self.z2)
        self.z3 = np.dot(self.a2, self.W2)
        yHat = self.sigmoid(self.z3)
        return yHat

    def costFunction(self, y, yHat):
        J = 0.5 * np.sum((y - yHat)**2)
        return J

    def backprop(self, X, y, yHat):
        # NB: delta 3 depends on the cost function applied
        delta3 = np.multiply(-(y - yHat), self.sigmoidPrime(self.z3))
        dJdW2 = np.dot(self.a2.T, delta3)
        delta2 = np.dot(delta3, self.W2.T) * self.sigmoidPrime(self.z2)
        dJdW1 = np.dot(X.T, delta2)
        return dJdW1, dJdW2

    def updateWeights(self, dJdW1, dJdW2, learnRate):
        if self.enNesterov:
            dW1_prev = self.dW1
            dW2_prev = self.dW2
            self.dW1 = learnRate*dJdW1 + self.momentum*self.dW1
            self.dW2 = learnRate*dJdW2 + self.momentum*self.dW2
            self.W1 = self.W1 - (1+self.momentum)*self.dW1 - self.momentum*dW1_prev
            self.W2 = self.W2 - (1+self.momentum)*self.dW2 - self.momentum*dW2_prev
        else:
            self.dW1 = learnRate*dJdW1 + self.momentum*self.dW1
            self.dW2 = learnRate*dJdW2 + self.momentum*self.dW2
            self.W1 = self.W1 - self.dW1
            self.W2 = self.W2 - self.dW2

    def getWeights(self):
        return np.concatenate((self.W1.ravel(), self.W2.ravel()))

    def setWeights(self, weights):
        W1_start = 0
        W1_end = self.nInput * self.nHidden
        self.W1 = np.reshape(weights[W1_start:W1_end], (self.nInput, self.nHidden))
        W2_end = W1_end + self.nHidden * self.nOutput
        self.W2 = np.reshape(weights[W1_end:W2_end], (self.nHidden, self.nOutput))

    def checkGradients(self, X, y):
        # Numerical gradient
        numGrad = self.computeNumericalGradient(X, y)

        # Backprop gradient
        yHat = self.feedforward(X)
        dJdW1, dJdW2 = self.backprop(X, y, yHat)
        grad = np.concatenate((dJdW1.ravel(), dJdW2.ravel()))
        self.minGrad.append(np.min(grad))

        # Compare
        diff = np.linalg.norm(grad - numGrad) / np.linalg.norm(grad + numGrad)
        if diff < self.GRAD_CHECK_THRESH:
            str = 'PASS'
        else:
            str = 'FAIL'
        print('[{0}] Gradient checking. diff = {1}'.format(str, diff))
        self.gradDiff.append(diff)

    def computeNumericalGradient(self, X, y):
        weights = self.getWeights()
        perturb = np.zeros(weights.shape)
        numGrad = np.zeros(weights.shape)

        for p in range(len(weights)):
            # Set pertubation for this weight only
            perturb[p] = self.EPSILSON
            # Positive perturbation
            self.setWeights(weights + perturb)
            yHat = self.feedforward(X)
            Jpos = self.costFunction(y, yHat)
            # Negative perturbation
            self.setWeights(weights - perturb)
            yHat = self.feedforward(X)
            Jneg = self.costFunction(y, yHat)
            # Compute Numerical Gradient
            numGrad[p] = (Jpos - Jneg) / (2 * self.EPSILSON)
            # Reset perturbation for next iteration
            perturb[p] = 0

        # Reset weights
        self.setWeights(weights)
        return numGrad

    def sigmoid(self, z):
        return 1 / (1 + np.exp(-z))

    def sigmoidPrime(self, z):
        x = np.exp(-z)
        return (x / ((1 + x)**2))


#  Data
X = np.array(([3, 5], [5, 1], [10, 2]), dtype=float)
y = np.array(([75], [82], [93]), dtype=float)

# Normalize
X = X/np.amax(X, axis=0)
y = y/100

# Train a network
NN = NeuralNet()
NN.train(X, y)

h = NN.feedforward(X)


# PLotting
f, (ax1, ax2, ax3) = plt.subplots(1, 3, sharey=False)
ax1.plot(NN.cost)
ax1.set(xlabel='Iteration', ylabel='Cost', title='')
ax1.grid()
ax2.plot(NN.gradDiff, 'r')
ax2.set(xlabel='', ylabel='norm(grad-numGrad)/norm(grad+numGrad)', title='')
ax2.grid()
ax3.plot(NN.minGrad, 'k')
ax3.set(xlabel='', ylabel='minimum gradient', title='')
ax3.grid()
plt.show(block=True)
$\endgroup$
  • 1
    $\begingroup$ Have you tried smaller values for $\epsilon$ (e.g. $10^{-8}$)? Near a minimum the gradient usually changes more quickly so a large step size might yield too coarse of a numerical approximation. $\endgroup$ – Ruben van Bergen Mar 6 '18 at 20:54
  • $\begingroup$ @RubenvanBergen Thanks, yes I have played with the epsilon value, and indeed $10^{-8}$ produces a pass to a tolerance of $10^{-8}$ throughout training as the gradients tend towards their minima. Is the test just as valid with a much smaller epsilon (within numerical tolerances)? I read $10^{-4}$ throughout the literature and am wondering why this value? $\endgroup$ – CatsLoveJazz Mar 6 '18 at 21:04
  • 1
    $\begingroup$ Personally, as a general rule I always try to use as small a step-size as possible (i.e. without causing numerical underflow issues) to get the best approximation. I don't know why $10^{-4}$ would be special - perhaps if you use some kind of standardization for your variables/cost function this is a step-size that is usually good enough? But yes, as far as I know the test is just as valid, if not more so, with smaller $\epsilon$. $\endgroup$ – Ruben van Bergen Mar 7 '18 at 7:57

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.