For learning purposes, I'm trying to implement a simple neural network with only linear layers followed by logistic activation. As far as I understand, the backpropagation algorithm exploits the structure of the chain rule in order to efficiently decompose the derivative of the loss (with respect to the weights and inputs) into factors corresponding to each layer.
So, I've written some training code corresponding to the following steps:
- Initialize weights for each layer to random values centered around 0
Forward propagate a training example through the network, i.e.:
For each layer in forward order, compute the output of the layer:
y = logistic(w * [1 | x])
wis the matrix of weights of which the first column is the bias,
yis the output of the layer,
xis the input of the layer (i.e. the output of the previous layer or the training example for the first layer)
- Store the output of the last layer
- Backward propagate error gradients, i.e.:
- For each layer in backward order:
- compute the gradient of the error in the layer's output with respect to its input
- compute the same error gradient with respect to the weights of the layer, which, by the chain rule, is the product of the gradient of the error w.r.t. the output and the gradient of the output w.r.t the weights
- adjust weights in the direction of this (second) antigradient (this is effectively stochastic gradient descent). The input for this procedure for the last layer is the gradient of the loss w.r.t. the output of the last layer; for the rest of the layers, it's the gradient w.r.t. the input of the next layer.
- For each layer in backward order:
- Repeat for another randomly-selected training example several times until convergence.
The code for the neural network is as follows. First, some helpers:
from scipy.special import erf import numpy as np # nonlinearity def logistic(x): return 0.5 * (1 + erf(x)) # Computes gradient (possibly, the Jacobian) of y with respect to xs def grad(y, xs): dx = 0.0001 return np.stack( (y(xs + dxs) - y(xs - dxs)) / (2 * dx) for dxs in dx * np.eye(len(xs.T)) ).T # Helper for Layer in the forward phase def forward_pure(weights, inputs): return logistic(weights @ inputs) def loss(actual, predicted): return np.sum((predicted - actual) ** 2) # Gradient of loss w.r.t. predictions def loss_grad(actual, predicted): return grad(lambda xs: loss(actual, xs), predicted)
Then, the class for a single layer:
class Layer: def __init__(self, input_size, output_size, learning_rate=1): # weights column 0 is bias self.weights = np.random.normal(size=(output_size, input_size + 1)) self.learning_rate = learning_rate def forward(self, inputs): # Concatenate leading 1 to inputs for applying the bias too self.inputs = np.array([1, *inputs]) self.outputs = forward_pure(self.weights, self.inputs) return self.outputs def backward(self, error_grad): # Derivative of output w.r.t. weights self.dout_dw = grad(lambda w: forward_pure(w, self.inputs), self.weights) # Derivative of output w.r.t. inputs # (input 0 is constant, so differentiating wrt. it is useless) self.dout_din = grad(lambda i: forward_pure(self.weights, i), self.inputs)[:, 1:] # Adjust weights # We note that according to the chain rule: # d[Error]/d[weights] = d[Error]/d[outputs] @ d[outputs]/d[weights] # self.derror_dw = error_grad @ self.dout_dw self.derror_dw = error_grad @ self.dout_dw self.weights -= self.learning_rate * self.derror_dw # Return gradient propagated backwards return self.dout_din
Finally, the class for a full multi-layer network:
class NN: def __init__(self, layers): self.layers = layers def forward(self, inputs): for layer in self.layers: inputs = layer.forward(inputs) self.outputs = inputs return inputs def backward(self, actual): # derivative of error w.r.t. output de_do = loss_grad(actual, self.outputs) for layer in reversed(self.layers): de_do = layer.backward(de_do) return de_do
Although I'm using linear layers only, each followed by the logistic activation, I'm performing numeric differentiation instead of analytically computing the derivatives because I'd like to abstract away the details of the layers, I'm trying to understand backpropagation in general.
The problem is that I've got an interesting bug. When my output is a single value (i.e. each training example comes with a label that is a one-element array), then the network learns perfectly fine. However, when the output is a vector with more than one element, then it converges to the wrong values. For example, with a single-layer network, I can successfully predict the AND and the OR functions, however, trying to do both at the same time fails:
The training code:
def train(nn, num_epochs, Xs, ys): for i in range(num_epochs): samples = list(zip(Xs, ys)) np.random.shuffle(samples) for X, y in samples: pred = nn.forward(X) err = nn.backward(y) for X, y in zip(Xs, ys): pred = nn.forward(X) print(X, y, np.round(pred, 3))
And its invocation:
print("Learning logical AND") train( nn = NN(layers=[ Layer(input_size=2, output_size=1), ]), num_epochs=10000, Xs = np.array([ [0, 0], [0, 1], [1, 0], [1, 1], ]), ys = np.array([ , , , , ]) ) print("Learning logical OR") train( nn = NN(layers=[ Layer(input_size=2, output_size=1), ]), num_epochs=10000, Xs = np.array([ [0, 0], [0, 1], [1, 0], [1, 1], ]), ys = np.array([ , , , , ]) ) print("Learning both at the same time") train( nn = NN(layers=[ Layer(input_size=2, output_size=2), ]), num_epochs=10000, Xs = np.array([ [0, 0], [0, 1], [1, 0], [1, 1], ]), ys = np.array([ [0, 0], [0, 1], [0, 1], [1, 1], ]) )
The output is:
Learning logical AND [0 0]  [0.] [0 1]  [0.003] [1 0]  [0.003] [1 1]  [0.997] Learning logical OR [0 0]  [0.003] [0 1]  [0.998] [1 0]  [0.998] [1 1]  [1.] Learning both at the same time [0 0] [0 0] [0. 0.] [0 1] [0 1] [0.127 0.279] [1 0] [0 1] [1. 0.996] [1 1] [1 1] [1. 1.]
For your reference, the full code can be found here so you don't have to copy-paste together the above snippets.