257 KiB
257 KiB
In [1]:
%matplotlib inline
"""
Simple code that tests XOR, OR and AND gates with linear regression
"""
# import necessary packages
import numpy as np
import matplotlib.pyplot as plt
from sklearn import datasets
def sigmoid(x):
return 1/(1 + np.exp(-x))
def feed_forward(X):
# weighted sum of inputs to the hidden layer
z_h = np.matmul(X, hidden_weights) + hidden_bias
# activation in the hidden layer
a_h = sigmoid(z_h)
# weighted sum of inputs to the output layer
z_o = np.matmul(a_h, output_weights) + output_bias
# softmax output
# axis 0 holds each input and axis 1 the probabilities of each category
probabilities = sigmoid(z_o)
return probabilities
# ensure the same random numbers appear every time
np.random.seed(0)
# Design matrix
X = np.array([ [0, 0], [0, 1], [1, 0],[1, 1]],dtype=np.float64)
# The XOR gate
yXOR = np.array( [ 0, 1 ,1, 0])
# The OR gate
yOR = np.array( [ 0, 1 ,1, 1])
# The AND gate
yAND = np.array( [ 0, 0 ,0, 1])
# Defining the neural network
n_inputs, n_features = X.shape
n_hidden_neurons = 2
n_categories = 1
n_features = 2
# we make the weights normally distributed using numpy.random.randn
# weights and bias in the hidden layer
hidden_weights = np.random.randn(n_features, n_hidden_neurons)
hidden_bias = np.zeros(n_hidden_neurons) + 0.01
# weights and bias in the output layer
output_weights = np.random.randn(n_hidden_neurons, n_categories)
output_bias = np.zeros(n_categories) + 0.01
probabilities = feed_forward(X)
print(probabilities)In [2]:
# import necessary packages
import numpy as np
import matplotlib.pyplot as plt
from sklearn.neural_network import MLPClassifier
from sklearn.metrics import accuracy_score
import seaborn as sns
# ensure the same random numbers appear every time
np.random.seed(0)
# Design matrix
X = np.array([ [0, 0], [0, 1], [1, 0],[1, 1]],dtype=np.float64)
# The XOR gate
yXOR = np.array( [ 0, 1 ,1, 0])
# The OR gate
yOR = np.array( [ 0, 1 ,1, 1])
# The AND gate
yAND = np.array( [ 0, 0 ,0, 1])
# Defining the neural network
n_hidden_neurons = 2
eta_vals = np.logspace(-5, 1, 7)
lmbd_vals = np.logspace(-5, 1, 7)
# store models for later use
DNN_scikit = np.zeros((len(eta_vals), len(lmbd_vals)), dtype=object)
epochs = 100
for i, eta in enumerate(eta_vals):
for j, lmbd in enumerate(lmbd_vals):
dnn = MLPClassifier(hidden_layer_sizes=(n_hidden_neurons), activation='logistic',
alpha=lmbd, learning_rate_init=eta, max_iter=epochs)
dnn.fit(X, yXOR)
DNN_scikit[i][j] = dnn
print("Learning rate = ", eta)
print("Lambda = ", lmbd)
print("Accuracy score on data set: ", dnn.score(X, yXOR))
print()
sns.set()
test_accuracy = np.zeros((len(eta_vals), len(lmbd_vals)))
for i in range(len(eta_vals)):
for j in range(len(lmbd_vals)):
dnn = DNN_scikit[i][j]
test_pred = dnn.predict(X)
test_accuracy[i][j] = accuracy_score(yXOR, test_pred)
fig, ax = plt.subplots(figsize = (10, 10))
sns.heatmap(test_accuracy, annot=True, ax=ax, cmap="viridis")
ax.set_title("Test Accuracy")
ax.set_ylabel("$\eta$")
ax.set_xlabel("$\lambda$")
plt.show()In [3]:
import autograd.numpy as np
class Scheduler:
"""
Abstract class for Schedulers
"""
def __init__(self, eta):
self.eta = eta
# should be overwritten
def update_change(self, gradient):
raise NotImplementedError
# overwritten if needed
def reset(self):
pass
class Constant(Scheduler):
def __init__(self, eta):
super().__init__(eta)
def update_change(self, gradient):
return self.eta * gradient
def reset(self):
pass
class Momentum(Scheduler):
def __init__(self, eta: float, momentum: float):
super().__init__(eta)
self.momentum = momentum
self.change = 0
def update_change(self, gradient):
self.change = self.momentum * self.change + self.eta * gradient
return self.change
def reset(self):
pass
class Adagrad(Scheduler):
def __init__(self, eta):
super().__init__(eta)
self.G_t = None
def update_change(self, gradient):
delta = 1e-8 # avoid division ny zero
if self.G_t is None:
self.G_t = np.zeros((gradient.shape[0], gradient.shape[0]))
self.G_t += gradient @ gradient.T
G_t_inverse = 1 / (
delta + np.sqrt(np.reshape(np.diagonal(self.G_t), (self.G_t.shape[0], 1)))
)
return self.eta * gradient * G_t_inverse
def reset(self):
self.G_t = None
class AdagradMomentum(Scheduler):
def __init__(self, eta, momentum):
super().__init__(eta)
self.G_t = None
self.momentum = momentum
self.change = 0
def update_change(self, gradient):
delta = 1e-8 # avoid division ny zero
if self.G_t is None:
self.G_t = np.zeros((gradient.shape[0], gradient.shape[0]))
self.G_t += gradient @ gradient.T
G_t_inverse = 1 / (
delta + np.sqrt(np.reshape(np.diagonal(self.G_t), (self.G_t.shape[0], 1)))
)
self.change = self.change * self.momentum + self.eta * gradient * G_t_inverse
return self.change
def reset(self):
self.G_t = None
class RMS_prop(Scheduler):
def __init__(self, eta, rho):
super().__init__(eta)
self.rho = rho
self.second = 0.0
def update_change(self, gradient):
delta = 1e-8 # avoid division ny zero
self.second = self.rho * self.second + (1 - self.rho) * gradient * gradient
return self.eta * gradient / (np.sqrt(self.second + delta))
def reset(self):
self.second = 0.0
class Adam(Scheduler):
def __init__(self, eta, rho, rho2):
super().__init__(eta)
self.rho = rho
self.rho2 = rho2
self.moment = 0
self.second = 0
self.n_epochs = 1
def update_change(self, gradient):
delta = 1e-8 # avoid division ny zero
self.moment = self.rho * self.moment + (1 - self.rho) * gradient
self.second = self.rho2 * self.second + (1 - self.rho2) * gradient * gradient
moment_corrected = self.moment / (1 - self.rho**self.n_epochs)
second_corrected = self.second / (1 - self.rho2**self.n_epochs)
return self.eta * moment_corrected / (np.sqrt(second_corrected + delta))
def reset(self):
self.n_epochs += 1
self.moment = 0
self.second = 0In [4]:
momentum_scheduler = Momentum(eta=1e-3, momentum=0.9)
adam_scheduler = Adam(eta=1e-3, rho=0.9, rho2=0.999)In [5]:
weights = np.ones((3,3))
print(f"Before scheduler:\n{weights=}")
epochs = 10
for e in range(epochs):
gradient = np.random.rand(3, 3)
change = adam_scheduler.update_change(gradient)
weights = weights - change
adam_scheduler.reset()
print(f"\nAfter scheduler:\n{weights=}")In [6]:
import autograd.numpy as np
def CostOLS(target):
def func(X):
return (1.0 / target.shape[0]) * np.sum((target - X) ** 2)
return func
def CostLogReg(target):
def func(X):
return -(1.0 / target.shape[0]) * np.sum(
(target * np.log(X + 10e-10)) + ((1 - target) * np.log(1 - X + 10e-10))
)
return func
def CostCrossEntropy(target):
def func(X):
return -(1.0 / target.size) * np.sum(target * np.log(X + 10e-10))
return funcIn [7]:
from autograd import grad
target = np.array([[1, 2, 3]]).T
a = np.array([[4, 5, 6]]).T
cost_func = CostCrossEntropy
cost_func_derivative = grad(cost_func(target))
valued_at_a = cost_func_derivative(a)
print(f"Derivative of cost function {cost_func.__name__} valued at a:\n{valued_at_a}")In [8]:
import autograd.numpy as np
from autograd import elementwise_grad
def identity(X):
return X
def sigmoid(X):
try:
return 1.0 / (1 + np.exp(-X))
except FloatingPointError:
return np.where(X > np.zeros(X.shape), np.ones(X.shape), np.zeros(X.shape))
def softmax(X):
X = X - np.max(X, axis=-1, keepdims=True)
delta = 10e-10
return np.exp(X) / (np.sum(np.exp(X), axis=-1, keepdims=True) + delta)
def RELU(X):
return np.where(X > np.zeros(X.shape), X, np.zeros(X.shape))
def LRELU(X):
delta = 10e-4
return np.where(X > np.zeros(X.shape), X, delta * X)
def derivate(func):
if func.__name__ == "RELU":
def func(X):
return np.where(X > 0, 1, 0)
return func
elif func.__name__ == "LRELU":
def func(X):
delta = 10e-4
return np.where(X > 0, 1, delta)
return func
else:
return elementwise_grad(func)In [9]:
z = np.array([[4, 5, 6]]).T
print(f"Input to activation function:\n{z}")
act_func = sigmoid
a = act_func(z)
print(f"\nOutput from {act_func.__name__} activation function:\n{a}")
act_func_derivative = derivate(act_func)
valued_at_z = act_func_derivative(a)
print(f"\nDerivative of {act_func.__name__} activation function valued at z:\n{valued_at_z}")In [10]:
import math
import autograd.numpy as np
import sys
import warnings
from autograd import grad, elementwise_grad
from random import random, seed
from copy import deepcopy, copy
from typing import Tuple, Callable
from sklearn.utils import resample
warnings.simplefilter("error")
class FFNN:
"""
Description:
------------
Feed Forward Neural Network with interface enabling flexible design of a
nerual networks architecture and the specification of activation function
in the hidden layers and output layer respectively. This model can be used
for both regression and classification problems, depending on the output function.
Attributes:
------------
I dimensions (tuple[int]): A list of positive integers, which specifies the
number of nodes in each of the networks layers. The first integer in the array
defines the number of nodes in the input layer, the second integer defines number
of nodes in the first hidden layer and so on until the last number, which
specifies the number of nodes in the output layer.
II hidden_func (Callable): The activation function for the hidden layers
III output_func (Callable): The activation function for the output layer
IV cost_func (Callable): Our cost function
V seed (int): Sets random seed, makes results reproducible
"""
def __init__(
self,
dimensions: tuple[int],
hidden_func: Callable = sigmoid,
output_func: Callable = lambda x: x,
cost_func: Callable = CostOLS,
seed: int = None,
):
self.dimensions = dimensions
self.hidden_func = hidden_func
self.output_func = output_func
self.cost_func = cost_func
self.seed = seed
self.weights = list()
self.schedulers_weight = list()
self.schedulers_bias = list()
self.a_matrices = list()
self.z_matrices = list()
self.classification = None
self.reset_weights()
self._set_classification()
def fit(
self,
X: np.ndarray,
t: np.ndarray,
scheduler: Scheduler,
batches: int = 1,
epochs: int = 100,
lam: float = 0,
X_val: np.ndarray = None,
t_val: np.ndarray = None,
):
"""
Description:
------------
This function performs the training the neural network by performing the feedforward and backpropagation
algorithm to update the networks weights.
Parameters:
------------
I X (np.ndarray) : training data
II t (np.ndarray) : target data
III scheduler (Scheduler) : specified scheduler (algorithm for optimization of gradient descent)
IV scheduler_args (list[int]) : list of all arguments necessary for scheduler
Optional Parameters:
------------
V batches (int) : number of batches the datasets are split into, default equal to 1
VI epochs (int) : number of iterations used to train the network, default equal to 100
VII lam (float) : regularization hyperparameter lambda
VIII X_val (np.ndarray) : validation set
IX t_val (np.ndarray) : validation target set
Returns:
------------
I scores (dict) : A dictionary containing the performance metrics of the model.
The number of the metrics depends on the parameters passed to the fit-function.
"""
# setup
if self.seed is not None:
np.random.seed(self.seed)
val_set = False
if X_val is not None and t_val is not None:
val_set = True
# creating arrays for score metrics
train_errors = np.empty(epochs)
train_errors.fill(np.nan)
val_errors = np.empty(epochs)
val_errors.fill(np.nan)
train_accs = np.empty(epochs)
train_accs.fill(np.nan)
val_accs = np.empty(epochs)
val_accs.fill(np.nan)
self.schedulers_weight = list()
self.schedulers_bias = list()
batch_size = X.shape[0] // batches
X, t = resample(X, t)
# this function returns a function valued only at X
cost_function_train = self.cost_func(t)
if val_set:
cost_function_val = self.cost_func(t_val)
# create schedulers for each weight matrix
for i in range(len(self.weights)):
self.schedulers_weight.append(copy(scheduler))
self.schedulers_bias.append(copy(scheduler))
print(f"{scheduler.__class__.__name__}: Eta={scheduler.eta}, Lambda={lam}")
try:
for e in range(epochs):
for i in range(batches):
# allows for minibatch gradient descent
if i == batches - 1:
# If the for loop has reached the last batch, take all thats left
X_batch = X[i * batch_size :, :]
t_batch = t[i * batch_size :, :]
else:
X_batch = X[i * batch_size : (i + 1) * batch_size, :]
t_batch = t[i * batch_size : (i + 1) * batch_size, :]
self._feedforward(X_batch)
self._backpropagate(X_batch, t_batch, lam)
# reset schedulers for each epoch (some schedulers pass in this call)
for scheduler in self.schedulers_weight:
scheduler.reset()
for scheduler in self.schedulers_bias:
scheduler.reset()
# computing performance metrics
pred_train = self.predict(X)
train_error = cost_function_train(pred_train)
train_errors[e] = train_error
if val_set:
pred_val = self.predict(X_val)
val_error = cost_function_val(pred_val)
val_errors[e] = val_error
if self.classification:
train_acc = self._accuracy(self.predict(X), t)
train_accs[e] = train_acc
if val_set:
val_acc = self._accuracy(pred_val, t_val)
val_accs[e] = val_acc
# printing progress bar
progression = e / epochs
print_length = self._progress_bar(
progression,
train_error=train_errors[e],
train_acc=train_accs[e],
val_error=val_errors[e],
val_acc=val_accs[e],
)
except KeyboardInterrupt:
# allows for stopping training at any point and seeing the result
pass
# visualization of training progression (similiar to tensorflow progression bar)
sys.stdout.write("\r" + " " * print_length)
sys.stdout.flush()
self._progress_bar(
1,
train_error=train_errors[e],
train_acc=train_accs[e],
val_error=val_errors[e],
val_acc=val_accs[e],
)
sys.stdout.write("")
# return performance metrics for the entire run
scores = dict()
scores["train_errors"] = train_errors
if val_set:
scores["val_errors"] = val_errors
if self.classification:
scores["train_accs"] = train_accs
if val_set:
scores["val_accs"] = val_accs
return scores
def predict(self, X: np.ndarray, *, threshold=0.5):
"""
Description:
------------
Performs prediction after training of the network has been finished.
Parameters:
------------
I X (np.ndarray): The design matrix, with n rows of p features each
Optional Parameters:
------------
II threshold (float) : sets minimal value for a prediction to be predicted as the positive class
in classification problems
Returns:
------------
I z (np.ndarray): A prediction vector (row) for each row in our design matrix
This vector is thresholded if regression=False, meaning that classification results
in a vector of 1s and 0s, while regressions in an array of decimal numbers
"""
predict = self._feedforward(X)
if self.classification:
return np.where(predict > threshold, 1, 0)
else:
return predict
def reset_weights(self):
"""
Description:
------------
Resets/Reinitializes the weights in order to train the network for a new problem.
"""
if self.seed is not None:
np.random.seed(self.seed)
self.weights = list()
for i in range(len(self.dimensions) - 1):
weight_array = np.random.randn(
self.dimensions[i] + 1, self.dimensions[i + 1]
)
weight_array[0, :] = np.random.randn(self.dimensions[i + 1]) * 0.01
self.weights.append(weight_array)
def _feedforward(self, X: np.ndarray):
"""
Description:
------------
Calculates the activation of each layer starting at the input and ending at the output.
Each following activation is calculated from a weighted sum of each of the preceeding
activations (except in the case of the input layer).
Parameters:
------------
I X (np.ndarray): The design matrix, with n rows of p features each
Returns:
------------
I z (np.ndarray): A prediction vector (row) for each row in our design matrix
"""
# reset matrices
self.a_matrices = list()
self.z_matrices = list()
# if X is just a vector, make it into a matrix
if len(X.shape) == 1:
X = X.reshape((1, X.shape[0]))
# Add a coloumn of zeros as the first coloumn of the design matrix, in order
# to add bias to our data
bias = np.ones((X.shape[0], 1)) * 0.01
X = np.hstack([bias, X])
# a^0, the nodes in the input layer (one a^0 for each row in X - where the
# exponent indicates layer number).
a = X
self.a_matrices.append(a)
self.z_matrices.append(a)
# The feed forward algorithm
for i in range(len(self.weights)):
if i < len(self.weights) - 1:
z = a @ self.weights[i]
self.z_matrices.append(z)
a = self.hidden_func(z)
# bias column again added to the data here
bias = np.ones((a.shape[0], 1)) * 0.01
a = np.hstack([bias, a])
self.a_matrices.append(a)
else:
try:
# a^L, the nodes in our output layers
z = a @ self.weights[i]
a = self.output_func(z)
self.a_matrices.append(a)
self.z_matrices.append(z)
except Exception as OverflowError:
print(
"OverflowError in fit() in FFNN\nHOW TO DEBUG ERROR: Consider lowering your learning rate or scheduler specific parameters such as momentum, or check if your input values need scaling"
)
# this will be a^L
return a
def _backpropagate(self, X, t, lam):
"""
Description:
------------
Performs the backpropagation algorithm. In other words, this method
calculates the gradient of all the layers starting at the
output layer, and moving from right to left accumulates the gradient until
the input layer is reached. Each layers respective weights are updated while
the algorithm propagates backwards from the output layer (auto-differentation in reverse mode).
Parameters:
------------
I X (np.ndarray): The design matrix, with n rows of p features each.
II t (np.ndarray): The target vector, with n rows of p targets.
III lam (float32): regularization parameter used to punish the weights in case of overfitting
Returns:
------------
No return value.
"""
out_derivative = derivate(self.output_func)
hidden_derivative = derivate(self.hidden_func)
for i in range(len(self.weights) - 1, -1, -1):
# delta terms for output
if i == len(self.weights) - 1:
# for multi-class classification
if (
self.output_func.__name__ == "softmax"
):
delta_matrix = self.a_matrices[i + 1] - t
# for single class classification
else:
cost_func_derivative = grad(self.cost_func(t))
delta_matrix = out_derivative(
self.z_matrices[i + 1]
) * cost_func_derivative(self.a_matrices[i + 1])
# delta terms for hidden layer
else:
delta_matrix = (
self.weights[i + 1][1:, :] @ delta_matrix.T
).T * hidden_derivative(self.z_matrices[i + 1])
# calculate gradient
gradient_weights = self.a_matrices[i][:, 1:].T @ delta_matrix
gradient_bias = np.sum(delta_matrix, axis=0).reshape(
1, delta_matrix.shape[1]
)
# regularization term
gradient_weights += self.weights[i][1:, :] * lam
# use scheduler
update_matrix = np.vstack(
[
self.schedulers_bias[i].update_change(gradient_bias),
self.schedulers_weight[i].update_change(gradient_weights),
]
)
# update weights and bias
self.weights[i] -= update_matrix
def _accuracy(self, prediction: np.ndarray, target: np.ndarray):
"""
Description:
------------
Calculates accuracy of given prediction to target
Parameters:
------------
I prediction (np.ndarray): vector of predicitons output network
(1s and 0s in case of classification, and real numbers in case of regression)
II target (np.ndarray): vector of true values (What the network ideally should predict)
Returns:
------------
A floating point number representing the percentage of correctly classified instances.
"""
assert prediction.size == target.size
return np.average((target == prediction))
def _set_classification(self):
"""
Description:
------------
Decides if FFNN acts as classifier (True) og regressor (False),
sets self.classification during init()
"""
self.classification = False
if (
self.cost_func.__name__ == "CostLogReg"
or self.cost_func.__name__ == "CostCrossEntropy"
):
self.classification = True
def _progress_bar(self, progression, **kwargs):
"""
Description:
------------
Displays progress of training
"""
print_length = 40
num_equals = int(progression * print_length)
num_not = print_length - num_equals
arrow = ">" if num_equals > 0 else ""
bar = "[" + "=" * (num_equals - 1) + arrow + "-" * num_not + "]"
perc_print = self._format(progression * 100, decimals=5)
line = f" {bar} {perc_print}% "
for key in kwargs:
if not np.isnan(kwargs[key]):
value = self._format(kwargs[key], decimals=4)
line += f"| {key}: {value} "
sys.stdout.write("\r" + line)
sys.stdout.flush()
return len(line)
def _format(self, value, decimals=4):
"""
Description:
------------
Formats decimal numbers for progress bar
"""
if value > 0:
v = value
elif value < 0:
v = -10 * value
else:
v = 1
n = 1 + math.floor(math.log10(v))
if n >= decimals - 1:
return str(round(value))
return f"{value:.{decimals-n-1}f}"In [11]:
import autograd.numpy as np
from sklearn.model_selection import train_test_split
def SkrankeFunction(x, y):
return np.ravel(0 + 1*x + 2*y + 3*x**2 + 4*x*y + 5*y**2)
def create_X(x, y, n):
if len(x.shape) > 1:
x = np.ravel(x)
y = np.ravel(y)
N = len(x)
l = int((n + 1) * (n + 2) / 2) # Number of elements in beta
X = np.ones((N, l))
for i in range(1, n + 1):
q = int((i) * (i + 1) / 2)
for k in range(i + 1):
X[:, q + k] = (x ** (i - k)) * (y**k)
return X
step=0.5
x = np.arange(0, 1, step)
y = np.arange(0, 1, step)
x, y = np.meshgrid(x, y)
target = SkrankeFunction(x, y)
target = target.reshape(target.shape[0], 1)
poly_degree=3
X = create_X(x, y, poly_degree)
X_train, X_test, t_train, t_test = train_test_split(X, target)In [12]:
input_nodes = X_train.shape[1]
output_nodes = 1
linear_regression = FFNN((input_nodes, output_nodes), output_func=identity, cost_func=CostOLS, seed=2023)In [13]:
linear_regression.reset_weights() # reset weights such that previous runs or reruns don't affect the weights
scheduler = Constant(eta=1e-3)
scores = linear_regression.fit(X_train, t_train, scheduler)In [14]:
linear_regression.reset_weights() # reset weights such that previous runs or reruns don't affect the weights
scores = linear_regression.fit(X_train, t_train, scheduler, lam=1e-4, epochs=1000)In [15]:
from sklearn.datasets import load_breast_cancer
from sklearn.preprocessing import MinMaxScaler
wisconsin = load_breast_cancer()
X = wisconsin.data
target = wisconsin.target
target = target.reshape(target.shape[0], 1)
X_train, X_val, t_train, t_val = train_test_split(X, target)
scaler = MinMaxScaler()
scaler.fit(X_train)
X_train = scaler.transform(X_train)
X_val = scaler.transform(X_val)In [16]:
input_nodes = X_train.shape[1]
output_nodes = 1
logistic_regression = FFNN((input_nodes, output_nodes), output_func=sigmoid, cost_func=CostLogReg, seed=2023)In [17]:
logistic_regression.reset_weights() # reset weights such that previous runs or reruns don't affect the weights
scheduler = Adam(eta=1e-3, rho=0.9, rho2=0.999)
scores = logistic_regression.fit(X_train, t_train, scheduler, epochs=1000, X_val=X_val, t_val=t_val)In [18]:
input_nodes = X_train.shape[1]
hidden_nodes1 = 100
hidden_nodes2 = 30
output_nodes = 1
dims = (input_nodes, hidden_nodes1, hidden_nodes2, output_nodes)
neural_network = FFNN(dims, hidden_func=RELU, output_func=sigmoid, cost_func=CostLogReg, seed=2023)In [19]:
neural_network.reset_weights() # reset weights such that previous runs or reruns don't affect the weights
scheduler = Adam(eta=1e-4, rho=0.9, rho2=0.999)
scores = neural_network.fit(X_train, t_train, scheduler, epochs=1000, X_val=X_val, t_val=t_val)In [20]:
from sklearn.datasets import load_digits
def onehot(target: np.ndarray):
onehot = np.zeros((target.size, target.max() + 1))
onehot[np.arange(target.size), target] = 1
return onehot
digits = load_digits()
X = digits.data
target = digits.target
target = onehot(target)
input_nodes = 64
hidden_nodes1 = 100
hidden_nodes2 = 30
output_nodes = 10
dims = (input_nodes, hidden_nodes1, hidden_nodes2, output_nodes)
multiclass = FFNN(dims, hidden_func=LRELU, output_func=softmax, cost_func=CostCrossEntropy)
multiclass.reset_weights() # reset weights such that previous runs or reruns don't affect the weights
scheduler = Adam(eta=1e-4, rho=0.9, rho2=0.999)
scores = multiclass.fit(X, target, scheduler, epochs=1000)In [21]:
X = np.array([ [0, 0], [0, 1], [1, 0],[1, 1]],dtype=np.float64)
# The XOR gate
yXOR = np.array( [[ 0], [1] ,[1], [0]])
input_nodes = X.shape[1]
output_nodes = 1
logistic_regression = FFNN((input_nodes, output_nodes), output_func=sigmoid, cost_func=CostLogReg, seed=2023)
logistic_regression.reset_weights() # reset weights such that previous runs or reruns don't affect the weights
scheduler = Adam(eta=1e-1, rho=0.9, rho2=0.999)
scores = logistic_regression.fit(X, yXOR, scheduler, epochs=1000)In [22]:
# import necessary packages
import numpy as np
import matplotlib.pyplot as plt
from sklearn import datasets
# ensure the same random numbers appear every time
np.random.seed(0)
# display images in notebook
%matplotlib inline
plt.rcParams['figure.figsize'] = (12,12)
# download MNIST dataset
digits = datasets.load_digits()
# define inputs and labels
inputs = digits.images
labels = digits.target
print("inputs = (n_inputs, pixel_width, pixel_height) = " + str(inputs.shape))
print("labels = (n_inputs) = " + str(labels.shape))
# flatten the image
# the value -1 means dimension is inferred from the remaining dimensions: 8x8 = 64
n_inputs = len(inputs)
inputs = inputs.reshape(n_inputs, -1)
print("X = (n_inputs, n_features) = " + str(inputs.shape))
# choose some random images to display
indices = np.arange(n_inputs)
random_indices = np.random.choice(indices, size=5)
for i, image in enumerate(digits.images[random_indices]):
plt.subplot(1, 5, i+1)
plt.axis('off')
plt.imshow(image, cmap=plt.cm.gray_r, interpolation='nearest')
plt.title("Label: %d" % digits.target[random_indices[i]])
plt.show()In [23]:
from sklearn.model_selection import train_test_split
# one-liner from scikit-learn library
train_size = 0.8
test_size = 1 - train_size
X_train, X_test, Y_train, Y_test = train_test_split(inputs, labels, train_size=train_size,
test_size=test_size)
# equivalently in numpy
def train_test_split_numpy(inputs, labels, train_size, test_size):
n_inputs = len(inputs)
inputs_shuffled = inputs.copy()
labels_shuffled = labels.copy()
np.random.shuffle(inputs_shuffled)
np.random.shuffle(labels_shuffled)
train_end = int(n_inputs*train_size)
X_train, X_test = inputs_shuffled[:train_end], inputs_shuffled[train_end:]
Y_train, Y_test = labels_shuffled[:train_end], labels_shuffled[train_end:]
return X_train, X_test, Y_train, Y_test
#X_train, X_test, Y_train, Y_test = train_test_split_numpy(inputs, labels, train_size, test_size)
print("Number of training images: " + str(len(X_train)))
print("Number of test images: " + str(len(X_test)))In [24]:
# building our neural network
n_inputs, n_features = X_train.shape
n_hidden_neurons = 50
n_categories = 10
# we make the weights normally distributed using numpy.random.randn
# weights and bias in the hidden layer
hidden_weights = np.random.randn(n_features, n_hidden_neurons)
hidden_bias = np.zeros(n_hidden_neurons) + 0.01
# weights and bias in the output layer
output_weights = np.random.randn(n_hidden_neurons, n_categories)
output_bias = np.zeros(n_categories) + 0.01In [25]:
# setup the feed-forward pass, subscript h = hidden layer
def sigmoid(x):
return 1/(1 + np.exp(-x))
def feed_forward(X):
# weighted sum of inputs to the hidden layer
z_h = np.matmul(X, hidden_weights) + hidden_bias
# activation in the hidden layer
a_h = sigmoid(z_h)
# weighted sum of inputs to the output layer
z_o = np.matmul(a_h, output_weights) + output_bias
# softmax output
# axis 0 holds each input and axis 1 the probabilities of each category
exp_term = np.exp(z_o)
probabilities = exp_term / np.sum(exp_term, axis=1, keepdims=True)
return probabilities
probabilities = feed_forward(X_train)
print("probabilities = (n_inputs, n_categories) = " + str(probabilities.shape))
print("probability that image 0 is in category 0,1,2,...,9 = \n" + str(probabilities[0]))
print("probabilities sum up to: " + str(probabilities[0].sum()))
print()
# we obtain a prediction by taking the class with the highest likelihood
def predict(X):
probabilities = feed_forward(X)
return np.argmax(probabilities, axis=1)
predictions = predict(X_train)
print("predictions = (n_inputs) = " + str(predictions.shape))
print("prediction for image 0: " + str(predictions[0]))
print("correct label for image 0: " + str(Y_train[0]))In [26]:
# to categorical turns our integer vector into a onehot representation
from sklearn.metrics import accuracy_score
# one-hot in numpy
def to_categorical_numpy(integer_vector):
n_inputs = len(integer_vector)
n_categories = np.max(integer_vector) + 1
onehot_vector = np.zeros((n_inputs, n_categories))
onehot_vector[range(n_inputs), integer_vector] = 1
return onehot_vector
#Y_train_onehot, Y_test_onehot = to_categorical(Y_train), to_categorical(Y_test)
Y_train_onehot, Y_test_onehot = to_categorical_numpy(Y_train), to_categorical_numpy(Y_test)
def feed_forward_train(X):
# weighted sum of inputs to the hidden layer
z_h = np.matmul(X, hidden_weights) + hidden_bias
# activation in the hidden layer
a_h = sigmoid(z_h)
# weighted sum of inputs to the output layer
z_o = np.matmul(a_h, output_weights) + output_bias
# softmax output
# axis 0 holds each input and axis 1 the probabilities of each category
exp_term = np.exp(z_o)
probabilities = exp_term / np.sum(exp_term, axis=1, keepdims=True)
# for backpropagation need activations in hidden and output layers
return a_h, probabilities
def backpropagation(X, Y):
a_h, probabilities = feed_forward_train(X)
# error in the output layer
error_output = probabilities - Y
# error in the hidden layer
error_hidden = np.matmul(error_output, output_weights.T) * a_h * (1 - a_h)
# gradients for the output layer
output_weights_gradient = np.matmul(a_h.T, error_output)
output_bias_gradient = np.sum(error_output, axis=0)
# gradient for the hidden layer
hidden_weights_gradient = np.matmul(X.T, error_hidden)
hidden_bias_gradient = np.sum(error_hidden, axis=0)
return output_weights_gradient, output_bias_gradient, hidden_weights_gradient, hidden_bias_gradient
print("Old accuracy on training data: " + str(accuracy_score(predict(X_train), Y_train)))
eta = 0.01
lmbd = 0.01
for i in range(1000):
# calculate gradients
dWo, dBo, dWh, dBh = backpropagation(X_train, Y_train_onehot)
# regularization term gradients
dWo += lmbd * output_weights
dWh += lmbd * hidden_weights
# update weights and biases
output_weights -= eta * dWo
output_bias -= eta * dBo
hidden_weights -= eta * dWh
hidden_bias -= eta * dBh
print("New accuracy on training data: " + str(accuracy_score(predict(X_train), Y_train)))In [27]:
class NeuralNetwork:
def __init__(
self,
X_data,
Y_data,
n_hidden_neurons=50,
n_categories=10,
epochs=10,
batch_size=100,
eta=0.1,
lmbd=0.0):
self.X_data_full = X_data
self.Y_data_full = Y_data
self.n_inputs = X_data.shape[0]
self.n_features = X_data.shape[1]
self.n_hidden_neurons = n_hidden_neurons
self.n_categories = n_categories
self.epochs = epochs
self.batch_size = batch_size
self.iterations = self.n_inputs // self.batch_size
self.eta = eta
self.lmbd = lmbd
self.create_biases_and_weights()
def create_biases_and_weights(self):
self.hidden_weights = np.random.randn(self.n_features, self.n_hidden_neurons)
self.hidden_bias = np.zeros(self.n_hidden_neurons) + 0.01
self.output_weights = np.random.randn(self.n_hidden_neurons, self.n_categories)
self.output_bias = np.zeros(self.n_categories) + 0.01
def feed_forward(self):
# feed-forward for training
self.z_h = np.matmul(self.X_data, self.hidden_weights) + self.hidden_bias
self.a_h = sigmoid(self.z_h)
self.z_o = np.matmul(self.a_h, self.output_weights) + self.output_bias
exp_term = np.exp(self.z_o)
self.probabilities = exp_term / np.sum(exp_term, axis=1, keepdims=True)
def feed_forward_out(self, X):
# feed-forward for output
z_h = np.matmul(X, self.hidden_weights) + self.hidden_bias
a_h = sigmoid(z_h)
z_o = np.matmul(a_h, self.output_weights) + self.output_bias
exp_term = np.exp(z_o)
probabilities = exp_term / np.sum(exp_term, axis=1, keepdims=True)
return probabilities
def backpropagation(self):
error_output = self.probabilities - self.Y_data
error_hidden = np.matmul(error_output, self.output_weights.T) * self.a_h * (1 - self.a_h)
self.output_weights_gradient = np.matmul(self.a_h.T, error_output)
self.output_bias_gradient = np.sum(error_output, axis=0)
self.hidden_weights_gradient = np.matmul(self.X_data.T, error_hidden)
self.hidden_bias_gradient = np.sum(error_hidden, axis=0)
if self.lmbd > 0.0:
self.output_weights_gradient += self.lmbd * self.output_weights
self.hidden_weights_gradient += self.lmbd * self.hidden_weights
self.output_weights -= self.eta * self.output_weights_gradient
self.output_bias -= self.eta * self.output_bias_gradient
self.hidden_weights -= self.eta * self.hidden_weights_gradient
self.hidden_bias -= self.eta * self.hidden_bias_gradient
def predict(self, X):
probabilities = self.feed_forward_out(X)
return np.argmax(probabilities, axis=1)
def predict_probabilities(self, X):
probabilities = self.feed_forward_out(X)
return probabilities
def train(self):
data_indices = np.arange(self.n_inputs)
for i in range(self.epochs):
for j in range(self.iterations):
# pick datapoints with replacement
chosen_datapoints = np.random.choice(
data_indices, size=self.batch_size, replace=False
)
# minibatch training data
self.X_data = self.X_data_full[chosen_datapoints]
self.Y_data = self.Y_data_full[chosen_datapoints]
self.feed_forward()
self.backpropagation()In [28]:
epochs = 100
batch_size = 100
dnn = NeuralNetwork(X_train, Y_train_onehot, eta=eta, lmbd=lmbd, epochs=epochs, batch_size=batch_size,
n_hidden_neurons=n_hidden_neurons, n_categories=n_categories)
dnn.train()
test_predict = dnn.predict(X_test)
# accuracy score from scikit library
print("Accuracy score on test set: ", accuracy_score(Y_test, test_predict))
# equivalent in numpy
def accuracy_score_numpy(Y_test, Y_pred):
return np.sum(Y_test == Y_pred) / len(Y_test)
#print("Accuracy score on test set: ", accuracy_score_numpy(Y_test, test_predict))In [29]:
eta_vals = np.logspace(-5, 1, 7)
lmbd_vals = np.logspace(-5, 1, 7)
# store the models for later use
DNN_numpy = np.zeros((len(eta_vals), len(lmbd_vals)), dtype=object)
# grid search
for i, eta in enumerate(eta_vals):
for j, lmbd in enumerate(lmbd_vals):
dnn = NeuralNetwork(X_train, Y_train_onehot, eta=eta, lmbd=lmbd, epochs=epochs, batch_size=batch_size,
n_hidden_neurons=n_hidden_neurons, n_categories=n_categories)
dnn.train()
DNN_numpy[i][j] = dnn
test_predict = dnn.predict(X_test)
print("Learning rate = ", eta)
print("Lambda = ", lmbd)
print("Accuracy score on test set: ", accuracy_score(Y_test, test_predict))
print()In [30]:
# visual representation of grid search
# uses seaborn heatmap, you can also do this with matplotlib imshow
import seaborn as sns
sns.set()
train_accuracy = np.zeros((len(eta_vals), len(lmbd_vals)))
test_accuracy = np.zeros((len(eta_vals), len(lmbd_vals)))
for i in range(len(eta_vals)):
for j in range(len(lmbd_vals)):
dnn = DNN_numpy[i][j]
train_pred = dnn.predict(X_train)
test_pred = dnn.predict(X_test)
train_accuracy[i][j] = accuracy_score(Y_train, train_pred)
test_accuracy[i][j] = accuracy_score(Y_test, test_pred)
fig, ax = plt.subplots(figsize = (10, 10))
sns.heatmap(train_accuracy, annot=True, ax=ax, cmap="viridis")
ax.set_title("Training Accuracy")
ax.set_ylabel("$\eta$")
ax.set_xlabel("$\lambda$")
plt.show()
fig, ax = plt.subplots(figsize = (10, 10))
sns.heatmap(test_accuracy, annot=True, ax=ax, cmap="viridis")
ax.set_title("Test Accuracy")
ax.set_ylabel("$\eta$")
ax.set_xlabel("$\lambda$")
plt.show()In [31]:
from sklearn.neural_network import MLPClassifier
# store models for later use
DNN_scikit = np.zeros((len(eta_vals), len(lmbd_vals)), dtype=object)
for i, eta in enumerate(eta_vals):
for j, lmbd in enumerate(lmbd_vals):
dnn = MLPClassifier(hidden_layer_sizes=(n_hidden_neurons), activation='logistic',
alpha=lmbd, learning_rate_init=eta, max_iter=epochs)
dnn.fit(X_train, Y_train)
DNN_scikit[i][j] = dnn
print("Learning rate = ", eta)
print("Lambda = ", lmbd)
print("Accuracy score on test set: ", dnn.score(X_test, Y_test))
print()In [32]:
# optional
# visual representation of grid search
# uses seaborn heatmap, could probably do this in matplotlib
import seaborn as sns
sns.set()
train_accuracy = np.zeros((len(eta_vals), len(lmbd_vals)))
test_accuracy = np.zeros((len(eta_vals), len(lmbd_vals)))
for i in range(len(eta_vals)):
for j in range(len(lmbd_vals)):
dnn = DNN_scikit[i][j]
train_pred = dnn.predict(X_train)
test_pred = dnn.predict(X_test)
train_accuracy[i][j] = accuracy_score(Y_train, train_pred)
test_accuracy[i][j] = accuracy_score(Y_test, test_pred)
fig, ax = plt.subplots(figsize = (10, 10))
sns.heatmap(train_accuracy, annot=True, ax=ax, cmap="viridis")
ax.set_title("Training Accuracy")
ax.set_ylabel("$\eta$")
ax.set_xlabel("$\lambda$")
plt.show()
fig, ax = plt.subplots(figsize = (10, 10))
sns.heatmap(test_accuracy, annot=True, ax=ax, cmap="viridis")
ax.set_title("Test Accuracy")
ax.set_ylabel("$\eta$")
ax.set_xlabel("$\lambda$")
plt.show()Warning:
Output truncated. This notebook contains too many cells to display efficiently.