Updated XOR spiking neuron example to work with current implementation.

Renamed spiking neuron output to fired.
Report mean genetic distance in DefaultSpeciesSet as an aid to choosing compatibility threshold.

README.md

@@ -7,8 +7,6 @@
 Significant, wide-spread changes have been made recently to the design of the library. 
 Once these are more stable they will be released as version 0.9.  This should happen in the first week or two of 2017.
 
-Currently, Izhikevich (spiking) neurons are not working.  
-
 ## About ##
 
 NEAT (NeuroEvolution of Augmenting Topologies) is a method developed by Kenneth O. Stanley for evolving arbitrary neural 

examples/xor/clean.sh

@@ -1,2 +1,3 @@
 #!/usr/bin/env bash
-rm *.csv *.gv *.svg
+rm *.csv *.gv *.svg
+rm neat-checkpoint-*

examples/xor/config-spiking

@@ -0,0 +1,94 @@
+#--- parameters for the XOR-2 experiment ---#
+
+[NEAT]
+max_fitness_threshold = 0.9
+pop_size              = 500
+reset_on_extinction   = False
+
+[IZGenome]
+# node bias options
+bias_init_mean          = 0.0
+bias_init_stdev         = 20.0
+bias_max_value          = 100.0
+bias_min_value          = -100.0
+bias_mutate_power       = 5.0
+bias_mutate_rate        = 0.7
+bias_replace_rate       = 0.1
+
+# genome compatibility options
+compatibility_disjoint_coefficient = 1.0
+compatibility_weight_coefficient   = 0.5
+
+# connection add/remove rates
+conn_add_prob           = 0.2
+conn_delete_prob        = 0.2
+
+# connection enable options
+enabled_default         = True
+enabled_mutate_rate     = 0.01
+
+feed_forward            = False
+initial_connection      = full
+
+# node add/remove rates
+node_add_prob           = 0.1
+node_delete_prob        = 0.1
+
+# network parameters
+num_hidden              = 0
+num_inputs              = 2
+num_outputs             = 2
+
+# node parameters for regular spiking
+a_init_mean      = 0.02
+a_init_stdev     = 0.0
+a_max_value      = 30.0
+a_min_value      = -30.0
+a_mutate_power   = 0.0
+a_mutate_rate    = 0.0
+a_replace_rate   = 0.0
+
+b_init_mean      = 0.2
+b_init_stdev     = 0.0
+b_max_value      = 30.0
+b_min_value      = -30.0
+b_mutate_power   = 0.0
+b_mutate_rate    = 0.0
+b_replace_rate   = 0.0
+
+c_init_mean      = -65.0
+c_init_stdev     = 0.0
+c_max_value      = 30.0
+c_min_value      = -30.0
+c_mutate_power   = 0.0
+c_mutate_rate    = 0.0
+c_replace_rate   = 0.0
+
+d_init_mean      = 8.0
+d_init_stdev     = 0.0
+d_max_value      = 30.0
+d_min_value      = -30.0
+d_mutate_power   = 0.0
+d_mutate_rate    = 0.0
+d_replace_rate   = 0.0
+
+# connection weight options
+weight_init_mean        = 0.0
+weight_init_stdev       = 10.0
+weight_max_value        = 30
+weight_min_value        = -30
+weight_mutate_power     = 5.0
+weight_mutate_rate      = 0.8
+weight_replace_rate     = 0.1
+
+[DefaultSpeciesSet]
+compatibility_threshold = 3.0
+
+[DefaultStagnation]
+species_fitness_func = mean
+max_stagnation       = 20
+
+[DefaultReproduction]
+elitism            = 2
+survival_threshold = 0.2
+

examples/xor/evolve-spiking.py

@@ -2,24 +2,24 @@
 from __future__ import print_function
 
 import os
+from matplotlib import pylab as plt
+from matplotlib import patches
 import neat
 
+import visualize
+
 # Network inputs and expected outputs.
 xor_inputs = ((0, 0), (0, 1), (1, 0), (1, 1))
 xor_outputs = (0, 1, 1, 0)
 
-# Use fast spiking neurons.
-neuron_params = neat.iznn.FAST_SPIKING_PARAMS
 # Maximum amount of simulated time (in milliseconds) to wait for the network to produce an output.
-max_time = 50.0
-# Use a simulation time step of 0.25 millisecond.
-dt = 0.25
+max_time_msec = 20.0
 
 
 def compute_output(t0, t1):
     '''Compute the network's output based on the "time to first spike" of the two output neurons.'''
     if t0 is None or t1 is None:
-        # If one of the output neurons failed to fire within the allotted time,
+        # If neither of the output neurons fired within the allotted time,
         # give a response which produces a large error.
         return -1.0
     else:
@@ -30,57 +30,67 @@ def compute_output(t0, t1):
         return max(0.0, min(1.0, response))
 
 
-def simulate(genome):
+def simulate(genome, config):
     # Create a network of "fast spiking" Izhikevich neurons.
-    net = neat.iznn.create_phenotype(genome, **neuron_params)
+    net = neat.iznn.IZNN.create(genome, config)
+    dt = net.get_time_step_msec()
     sum_square_error = 0.0
     simulated = []
-    for inputData, outputData in zip(xor_inputs, xor_outputs):
+    for idata, odata in zip(xor_inputs, xor_outputs):
         neuron_data = {}
         for i, n in net.neurons.items():
             neuron_data[i] = []
 
         # Reset the network, apply the XOR inputs, and run for the maximum allowed time.
         net.reset()
-        net.set_inputs(inputData)
+        net.set_inputs(idata)
         t0 = None
         t1 = None
         v0 = None
         v1 = None
-        num_steps = int(max_time / dt)
+        num_steps = int(max_time_msec / dt)
+        net.set_inputs(idata)
         for j in range(num_steps):
             t = dt * j
             output = net.advance(dt)
 
             # Capture the time and neuron membrane potential for later use if desired.
             for i, n in net.neurons.items():
-                neuron_data[i].append((t, n.v))
+                neuron_data[i].append((t, n.current, n.v, n.u, n.fired))
 
             # Remember time and value of the first output spikes from each neuron.
             if t0 is None and output[0] > 0:
-                t0, v0 = neuron_data[net.outputs[0]][-2]
+                t0, I0, v0, u0, f0 = neuron_data[net.outputs[0]][-2]
 
             if t1 is None and output[1] > 0:
-                t1, v1 = neuron_data[net.outputs[1]][-2]
+                t1, I1, v1, u1, f0 = neuron_data[net.outputs[1]][-2]
 
         response = compute_output(t0, t1)
-        sum_square_error += (response - outputData) ** 2
+        sum_square_error += (response - odata) ** 2
+
+        #print(genome)
+        #visualize.plot_spikes(neuron_data[net.outputs[0]], False)
+        #visualize.plot_spikes(neuron_data[net.outputs[1]], True)
 
-        simulated.append((inputData, outputData, t0, t1, v0, v1, neuron_data))
+        simulated.append((idata, odata, t0, t1, v0, v1, neuron_data))
 
     return sum_square_error, simulated
 
 
-def eval_fitness(genomes):
-    for genome in genomes:
-        sum_square_error, simulated = simulate(genome)
-        genome.fitness = 1 - sum_square_error
+def eval_genome(genome, config):
+    sum_square_error, simulated = simulate(genome, config)
+    return 1.0 - sum_square_error
+
+
+def eval_genomes(genomes, config):
+    for genome_id, genome in genomes:
+        genome.fitness = eval_genome(genome, config)
 
 
 def run(config_path):
     # Load the config file, which is assumed to live in
     # the same directory as this script.
-    config = neat.Config(neat.DefaultGenome, neat.DefaultReproduction,
+    config = neat.Config(neat.iznn.IZGenome, neat.DefaultReproduction,
                          neat.DefaultSpeciesSet, neat.DefaultStagnation,
                          config_path)
 
@@ -91,52 +101,59 @@ def run(config_path):
     config.output_nodes = 2
 
     pop = neat.population.Population(config)
-    pop.run(eval_fitness, 200)
 
-    print('Number of evaluations: {0}'.format(pop.total_evaluations))
+    # Add a stdout reporter to show progress in the terminal.
+    pop.add_reporter(neat.StdOutReporter())
+    stats = neat.StatisticsReporter()
+    pop.add_reporter(stats)
 
-    # Display the winning genome.
-    winner = pop.statistics.best_genome()
-    print('\nBest genome:\n{!s}'.format(winner))
+    if 0:
+        winner = pop.run(eval_genomes, 300)
+    else:
+        pe = neat.ParallelEvaluator(6, eval_genome)
+        winner = pop.run(pe.evaluate, 300)
 
-    # Show output of the most fit genome against training data.
-    winner = pop.statistics.best_genome()
+    # Display the winning genome.
     print('\nBest genome:\n{!s}'.format(winner))
-    print('\nBest network output:')
 
-    # Create a network of "fast spiking" Izhikevich neurons, simulation time step 0.25 millisecond.
-    net = neat.iznn.create_phenotype(winner, **neuron_params)
-    for inputData, outputData in zip(xor_inputs, xor_outputs):
-        neuron_data = {}
-        for i, n in net.neurons.items():
-            neuron_data[i] = []
+    node_names = {-1:'A', -2: 'B', 0:'A XOR B'}
+    visualize.draw_net(config, winner, True, node_names=node_names)
+    visualize.plot_stats(stats, ylog=False, view=True)
+    visualize.plot_species(stats, view=True)
 
-        # Reset the network, apply the XOR inputs, and run for the maximum allowed time.
-        net.reset()
-        net.set_inputs(inputData)
-        t0 = None
-        t1 = None
-        num_steps = int(max_time / dt)
-        for j in range(num_steps):
-            t = dt * j
-            output = net.advance(dt)
+    # Show output of the most fit genome against training data, and create
+    # a plot of the traces out to the max time for each set of inputs.
+    print('\nBest network output:')
+    plt.figure(figsize=(12, 12))
+    sum_square_error, simulated = simulate(winner, config)
+    for r, (inputData, outputData, t0, t1, v0, v1, neuron_data) in enumerate(simulated):
+        response = compute_output(t0, t1)
+        print("{0!r} expected {1:.3f} got {2:.3f}".format(inputData, outputData, response))
 
-            # Capture the time and neuron membrane potential for later use if desired.
-            for i, n in net.neurons.items():
-                neuron_data[i].append((t, n.v))
+        axes = plt.subplot(4, 1, r + 1)
+        plt.title("Traces for XOR input {{{0:.1f}, {1:.1f}}}".format(*inputData), fontsize=12)
+        for i, s in neuron_data.items():
+            if i in [0, 1]:
+                t, I, v, u, fired = zip(*s)
+                plt.plot(t, v, "-", label="neuron {0:d}".format(i))
 
-            # Remember time and value of the first output spikes from each neuron.
-            if t0 is None and output[0] > 0:
-                t0, v0 = neuron_data[net.outputs[0]][-2]
+        # Circle the first peak of each output.
+        circle0 = patches.Ellipse((t0, v0), 1.0, 10.0, color='r', fill=False)
+        circle1 = patches.Ellipse((t1, v1), 1.0, 10.0, color='r', fill=False)
+        axes.add_artist(circle0)
+        axes.add_artist(circle1)
 
-            if t1 is None and output[1] > 0:
-                t1, v1 = neuron_data[net.outputs[1]][-2]
+        plt.ylabel("Potential (mv)", fontsize=10)
+        plt.ylim(-100, 50)
+        plt.tick_params(labelsize=8)
+        plt.grid()
 
-        response = compute_output(t0, t1)
-        print(inputData, response)
+    plt.xlabel("Time (in ms)", fontsize=10)
+    plt.tight_layout(pad=0.4, w_pad=0.5, h_pad=1.0)
+    plt.savefig("traces.png", dpi=90)
+    plt.show()
 
 
 if __name__ == '__main__':
     local_dir = os.path.dirname(__file__)
-    config_path = os.path.join(local_dir, 'config-feedforward')
-    run(config_path)
+    run(os.path.join(local_dir, 'config-spiking'))

examples/xor/visualize.py

@@ -41,17 +41,14 @@ def plot_stats(statistics, ylog=False, view=False, filename='avg_fitness.svg'):
 
 def plot_spikes(spikes, view=False, filename=None, title=None):
     """ Plots the trains for a single spiking neuron. """
-    if plt is None:
-        warnings.warn("This display is not available due to a missing optional dependency (matplotlib)")
-        return
-
-    t_values = [t for t, I, v, u in spikes]
-    v_values = [v for t, I, v, u in spikes]
-    u_values = [u for t, I, v, u in spikes]
-    I_values = [I for t, I, v, u in spikes]
+    t_values = [t for t, I, v, u, f in spikes]
+    v_values = [v for t, I, v, u, f in spikes]
+    u_values = [u for t, I, v, u, f in spikes]
+    I_values = [I for t, I, v, u, f in spikes]
+    f_values = [f for t, I, v, u, f in spikes]
 
     fig = plt.figure()
-    plt.subplot(3, 1, 1)
+    plt.subplot(4, 1, 1)
     plt.ylabel("Potential (mv)")
     plt.xlabel("Time (in ms)")
     plt.grid()
@@ -62,13 +59,19 @@ def plot_spikes(spikes, view=False, filename=None, title=None):
     else:
         plt.title("Izhikevich's spiking neuron model ({0!s})".format(title))
 
-    plt.subplot(3, 1, 2)
+    plt.subplot(4, 1, 2)
+    plt.ylabel("Fired")
+    plt.xlabel("Time (in ms)")
+    plt.grid()
+    plt.plot(t_values, f_values, "r-")
+
+    plt.subplot(4, 1, 3)
     plt.ylabel("Recovery (u)")
     plt.xlabel("Time (in ms)")
     plt.grid()
     plt.plot(t_values, u_values, "r-")
 
-    plt.subplot(3, 1, 3)
+    plt.subplot(4, 1, 4)
     plt.ylabel("Current (I)")
     plt.xlabel("Time (in ms)")
     plt.grid()