// Simple neural nets
// (c) Alasdair Turner 2009
 
// Free software: you can redistribute this program and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
 
Network neuralnet;
 
PFont font;
 
void setup()
{
  size(400,400);
 
  font = loadFont("LucidaSans-20.vlw");
  textFont(font);
 
  setupSigmoid();
  loadData();
 
  neuralnet = new Network(196,49,10);
 
  background(220,204,255);
  noStroke();
  smooth();
  pushMatrix();
  neuralnet.draw();
  popMatrix();
  fill(0);
}
 
boolean b_dataloaded = false;
 
// left click to test, right click (or ctrl+click on a Mac) to train
boolean b_train = false, b_test = false;
 
void draw()
{
  int response = -1, actual = -1;
  if (!b_dataloaded) {
    loadData();
    b_dataloaded = true;
    b_test = true;
  }
  if (b_train) {
    // this allows some fast training without displaying:
    for (int i = 0; i < 500; i++) {
      // select a random training input and train
      int row = (int) floor(random(0,training_set.length));
      response = neuralnet.respond(training_set[row].inputs);
      actual = training_set[row].output;
      neuralnet.train(training_set[row].outputs); 
    }
  }
  else if (b_test) {
    int row = (int) floor(random(0,testing_set.length));
    response = neuralnet.respond(testing_set[row].inputs);
    actual = testing_set[row].output;
  }
  if (b_train || b_test) {
    // draw
    background(220,204,255);
    noStroke();
    smooth();
    pushMatrix();
    neuralnet.draw();
    popMatrix();
    b_train = b_test = false;
    fill(0);
    text(str(response),350,27);
    text(str(actual),350,275);
  }
  fill(0);
}
 
void mousePressed() 
{
  if (mouseButton == LEFT) {
    b_test = true;
  }
  else {
    b_train = true;
  }
}
 
 
 
 

// Simple neural nets: load data
// (c) Alasdair Turner 2009
 
// Free software: you can redistribute this program and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
 
// this uses the MNIST database of handwritten digits
// http://yann.lecun.com/exdb/mnist/ (accessed 04.06.09)
// Yann LeCun and Corinna Cortes
 
// note: I have reduced the originals to 14 x 14 from 28 x 28
 
Datum [] training_set;
Datum [] testing_set;
 
class Datum
{
  float [] inputs;
  float [] outputs;
  int output;
  Datum()
  {
    inputs = new float [196];
    outputs = new float[10];
  }
  void imageLoad(byte [] images, int offset)
  {
    for (int i = 0; i < 196; i++) {
      // note, you must use int() to convert rather than (int) to cast:
      inputs[i] = int(images[i+offset]) / 128.0 - 1.0;
    }
  }
  void labelLoad(byte [] labels, int offset)
  {
    output = int(labels[offset]);
    for (int i = 0; i < 10; i++) {
      if (i == output) {
        outputs[i] = 1.0;
      }
      else {
        outputs[i] = -1.0;
      }
    }
  }
}
 
void loadData()
{
  byte [] images = loadBytes("t10k-images-14x14.idx3-ubyte");
  byte [] labels = loadBytes("t10k-labels.idx1-ubyte");
 
  training_set = new Datum [8000];
  int tr_pos = 0;
  testing_set = new Datum [2000];
  int te_pos = 0;
 
  for (int i = 0; i < 10000; i++) {
    if (i % 5 != 0) {
      training_set[tr_pos] = new Datum();
      training_set[tr_pos].imageLoad(images,16 + i * 196);
      training_set[tr_pos].labelLoad(labels,8 + i);
      tr_pos++;
    }
    else {
      testing_set[te_pos] = new Datum();
      testing_set[te_pos].imageLoad(images,16 + i * 196);
      testing_set[te_pos].labelLoad(labels,8 + i);
      te_pos++;
    }
  }
}
 
 
 
 
 

// Simple neural nets: network
// (c) Alasdair Turner 2009
 
// Free software: you can redistribute this program and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
 
// This class is for the neural network,
// which is hard coded with three layers:
// input, hidden and output
 
class Network
{
  // this network is hard coded to only have one hidden layer
  Neuron [] m_input_layer; 
  Neuron [] m_hidden_layer;
  Neuron [] m_output_layer;
  // create a network specifying numbers of inputs, hidden layer neurons
  // and number of outputs, e.g. Network(4,4,3)
  Network(int inputs, int hidden, int outputs)
  {
    m_input_layer = new Neuron [inputs];
    m_hidden_layer = new Neuron [hidden];
    m_output_layer = new Neuron [outputs];
    // set up the network topology
    for (int i = 0; i < m_input_layer.length; i++) {
      m_input_layer[i] = new Neuron();
    }
    // route the input layer to the hidden layer
    for (int j = 0; j < m_hidden_layer.length; j++) {
      m_hidden_layer[j] = new Neuron(m_input_layer);
    }
    // route the hidden layer to the output layer
    for (int k = 0; k < m_output_layer.length; k++) {
      m_output_layer[k] = new Neuron(m_hidden_layer);
    }  
  }
  int respond(float [] inputs)
  {
    float [] responses = new float [m_output_layer.length];
    // feed forward
    // simply set the input layer to display the inputs
    for (int i = 0; i < m_input_layer.length; i++) {
      m_input_layer[i].m_output = inputs[i];
    }
    // now feed forward through the hidden layer
    for (int j = 0; j < m_hidden_layer.length; j++) {
      m_hidden_layer[j].respond();
    }
    // and finally feed forward to the output layer
    for (int k = 0; k < m_output_layer.length; k++) {
      responses[k] = m_output_layer[k].respond();
    }
    // now check the best response:
    int response = -1;
    float best = max(responses);
    for (int a = 0; a < responses.length; a++) {
      if (responses[a] == best) {
        response = a;
      }
    }
    return response;
  }
  void train(float [] outputs)
  {
    // adjust the output layer
    for (int k = 0; k < m_output_layer.length; k++) {
      m_output_layer[k].finderror(outputs[k]);
      m_output_layer[k].train();
    }
    // propagate back to the hidden layer
    for (int j = 0; j < m_hidden_layer.length; j++) {
      m_hidden_layer[j].train();
    }
    // the input layer doesn't learn:
    // it is simply the inputs
  }
  void draw()
  {
    // note, this draw is hard-coded for Network(196,49,10)
    // which reflects my use of the MNIST database of handwritten digits
    for (int i = 0; i < m_input_layer.length; i++) {
      pushMatrix();
      translate((i%14) * width / 25.0 + width * 0.22, (i/14) * height / 25.0 + height * 0.42);
      m_input_layer[i].draw(true);
      popMatrix();
    }
    for (int j = 0; j < m_hidden_layer.length; j++) {
      pushMatrix();
      translate((j%7) * width / 25.0 + width * 0.36, (j/7) * height / 25.0 + height * 0.12);
      m_hidden_layer[j].draw(false);
      popMatrix();
    }
    // this is slightly tricky -- I've switched the order so the output
    // neurons are arrange 1,2,3...8,9,0 rather than 0,1,2...7,8,9
    // (that's what the (k+9) % 10 is doing)
    for (int k = 0; k < m_output_layer.length; k++) {
      pushMatrix();
      translate(((k+9)%10) * width / 20.0 + width * 0.25, height * 0.05);
      m_output_layer[k].draw(true);
      popMatrix();
    }
  }
}
 
 
 

// Simple neural nets: neuron
// (c) Alasdair Turner 2009
 
// Free software: you can redistribute this program and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
 
// This class is for each neuron.  It works
// as a feed-forward multilayer perceptron, 
// learning by backpropagation
 
float LEARNING_RATE = 0.01;
 
class Neuron
{
  Neuron [] m_inputs;
  float [] m_weights;
  float m_threshold;
  float m_output;
  float m_error;
  // the input layer of neurons have no inputs:
  Neuron()
  {
    m_threshold = 0.0;
    m_error = 0.0;
    // initial random output
    m_output = lookupSigmoid(random(-5.0,5.0));
  }
  // all other layers (hidden and output) have 
  // neural inputs
  Neuron(Neuron [] inputs)
  {
    m_inputs = new Neuron [inputs.length];
    m_weights = new float [inputs.length];
    for (int i = 0; i < inputs.length; i++) {
      m_inputs[i] = inputs[i];
      m_weights[i] = random(-1.0,1.0);
    }
    m_threshold = random(-1.0,1.0);
    m_error = 0.0;
    // initial random output
    m_output = lookupSigmoid(random(-5.0,5.0));
  }
  // respond looks at the layer below, and prepares a response:
  float respond()
  {
    float input = 0.0;
    for (int i = 0; i < m_inputs.length; i++) {
      input += m_inputs[i].m_output * m_weights[i];
    }
    m_output = lookupSigmoid(input + m_threshold);
    // reset our error value ready for training
    m_error = 0.0;
    return m_output;
  }
  // find error is used on the output neurons
  void finderror(float desired)
  {
    m_error = desired - m_output;
  }
  // train adjusts the weights by comparing actual output to correct output
  void train()
  {
    float delta = (1.0 - m_output) * (1.0 + m_output) * m_error * LEARNING_RATE;
    for (int i = 0; i < m_inputs.length; i++) {
      // tell the next layer down what it's doing wrong
      m_inputs[i].m_error += m_weights[i] * m_error;
      // correct our weights
      m_weights[i] += m_inputs[i].m_output * delta;
    }
  }
  void draw(boolean not_hidden)
  {
    float level = (0.5-m_output*0.5); 
    if (not_hidden) {
      fill(255 * level);
    }
    else {
      // merge hidden layer with background color
      fill(110 + 128 * level, 102 + 128 * level, 127 + 128 * level);
    }
    ellipse(0,0,16,16);
  }
}
 
 
 
 
 

// Simple neural nets: sigmoid functions
// (c) Alasdair Turner 2009
 
// Free software: you can redistribute this program and/or modify
// it under the terms of the GNU General Public License as published by
// the Free Software Foundation, either version 3 of the License, or
// (at your option) any later version.
 
// a sigmoid function is the neuron's response to inputs
// the sigmoidal response ranges from -1.0 to 1.0
 
// for example, the weighted sum of inputs might be "2.1"
// our response would be lookupSigmoid(2.1) = 0.970
 
// this is a look up table for sigmoid (neural response) values
// which is valid from -5.0 to 5.0
float [] g_sigmoid = new float [200];
 
// this function precalculate a sigmoid (response) function
void setupSigmoid()
{
  for (int i = 0; i < 200; i++) {
    float x = (i / 20.0) - 5.0;
    g_sigmoid[i] = 2.0 / (1.0 + exp(-2.0 * x)) - 1.0;
  }
}
 
// once the sigmoid has been set up, this function accesses it:
float lookupSigmoid(float x)
{
  return g_sigmoid[constrain((int) floor((x + 5.0) * 20.0),0,199)];
}