In [1]:
from matplotlib import pylab
import nengo
import random
import numpy as np
import gzip as gz
import cPickle
from cPickle import load
try:
    import Image
except ImportError:
    from PIL import Image
from scipy.sparse.linalg import svds
import scipy
from scipy import ndimage
import matplotlib.pyplot as plt
import matplotlib.animation as animation
from nengo.utils.compat import is_iterable
import collections

#%matplotlib inline #Makes visualizations appar inline (Commented out because animation popup as new window)

In [2]:
#The code in this cell is for reading the images from the MNIST database and not part of the brain model.
def load_img(path, dims):
    """Load the image at path and return an array representing the raster.
    Flattens image. Shifts pixel activations such that 0 represents gray,
    normalizes the output array.
    Keyword arguments:
    path -- str, path of the image to be loaded.
    dims -- (w, h), where w,h are ints indicating dimensions of the image (in
        px)."""

    img = Image.open(path).resize(dims).getdata()
    img.convert('L')
    img = subtract(array(img).flatten(), 127.5)
    return img/norm(img)


def load_data(filename):
    """Uncompress, unpickle and return a .pkl.gz file.
    Keyword arguments:
    filename -- str, a valid file path"""

    return load(gz.open(filename))

def load_mini_mnist(option=None):
    """Load and return the first \%10 of the images in the mnist dataset.
    Does not return labels. Pass in 'train', 'valid' or 'test' if you want to
    load a specific subset of the dataset.
    Keyword arguments:
    option -- str (default=None)."""

    mini_mnist = load(gz.open('./mini_mnist.pkl.gz', 'rb'))
    if option == 'train':
        return mini_mnist[0]
    elif option == 'valid':
        return mini_mnist[1]
    elif option == 'test':
        return mini_mnist[2]
    else:
        return mini_mnist

In [3]:
def rotate_img(img, degrees):
    '''
    img is the dim**2 by 1 vector representing the pixel values.
    Rotates image the degrees passed in counterclockwise
    Returns the Reshaped image (to original shape which is the one dimensional vector)
    dim is a global variable
    '''
    original = img.shape
    newImg = scipy.ndimage.interpolation.rotate(np.reshape(img, (dim,dim), 'F'),degrees,reshape=False)
    newImg = np.reshape(newImg, original, 'F')
    return newImg

In [26]:
conn_synapse = 0.1 #post synaptic time constant to use for filtering (pstc) - what does changing this do?
probe_synapse = 0.01 #pstc
multiplier = 2 #not used
n_neurons = 5000
direct = True #Direct - function computed explicitly instead of in neurons 
stop_time = 3.0
run_time = 3.0 #in seconds

In [5]:
dim = 28 #size of the image
mnist = load_mini_mnist()
#train = mnist[0] #collection of training images
img = mnist[1][0] #image to be used for testing
#compress_size = 400 #?
#basis, S, V = svds(train.T, k=compress_size) #Used for encoding and decoding information 
    #a set of vectors representing what a hand drawn number should look like?

In [6]:
#Need same number of vectors in basis as number of neurons (randomly sample from basis)
#expanded_basis = np.array([random.choice(basis.T) for _ in range(n_neurons)])

In [7]:
def stim_func(t):
    '''returns the image for first 0.1s'''
    if t < 0.01:
        return img
    else:
        return [0 for _ in range(len(img))]

In [8]:
def connection_func(x):
    '''takes the output from the first ensemble and rotates it 1 degree'''
    return rotate_img(x,1)

In [57]:
'''def iter_func(i, x):
    matrix = np.zeros((dim,dim))
    matrix[i//dim][i%dim] = x

    matrix = scipy.ndimage.interpolation.rotate(matrix,1,reshape=False)

    out = np.reshape(matrix, (dim**2,), 'F')

def iter_rotate(x):
    for i in range(dim**2):
            
        yield iter_func(i,x)


#generator_func = iter_rotate(x)

E = []
for i in range(dim**2):
    E += iter_rotate(x) 
print(is_iterable(iter_rotate))
print(isinstance(iter_rotate, collections.Iterable))

def coordinates():
    for i in range(10):
        x=i
        y=0
        yield(x,y)
        
print(list(coordinates()))


list = range(dim**2)

list = iter(list)


def iter_rotate(x):
    matrix = np.zeros((dim,dim))
    i = list.next()
    print(i)
    if(i>=dim**2):
        return matrix
    matrix[i//dim][i%dim] = x

    matrix = scipy.ndimage.interpolation.rotate(matrix,1,reshape=False)

    
    out = np.reshape(matrix, (784L,), 'F')
    return out
    
'''

def iterable_fn(n):
    for i in range(n):
        def rotate(x):
            matrix = np.zeros((dim,dim))
            matrix[i//dim][i%dim] = x

            matrix = scipy.ndimage.interpolation.rotate(matrix,1,reshape=False)
            out = np.reshape(matrix, (784L,), 'F')
            return out
        yield rotate

rotator = iterable_fn(784)
#rotator.next()

rotator2 = list(rotator)
#print len(rotator2)


def comb_func(t,x):
    #print("Here")
    out = np.zeros((dim,dim))
    out = np.reshape(out, (784L,), 'F')
    for i in range(dim**2):
        out = np.add(x[(dim**2) * i : (dim**2) *(i+1)], out)
    
    return out

#pylab.imshow(np.reshape(iter_rotate(1), (dim,dim), 'F'), cmap='Greys_r')
#plt.show()
#print(np.size(iter_rotate(1)))

In [58]:
#A network is primarily used for grouping together related objects and connections for visualization purposes
with nengo.Network() as net:
    
    if direct:
        neuron_type = nengo.Direct() #function computed explicitly, instead of in neurons
    else:
        neuron_type = nengo.LIF() #spiking version of the leaky integrate-and-fire neuron model

    #Input stimulus - provide data to the ensemble
    ipt = nengo.Node(stim_func)
    
    '''An array of ensembles. This acts, in some ways, like a single high-dimensional ensemble,
    but actually consists of many sub-ensembles, each one representing a separate dimension. 
    This tends to be much faster to create and can be more accurate than having one huge 
    high-dimensional ensemble. However, since the neurons represent different dimensions separately,
    we cannot compute nonlinear interactions between those dimensions.'''
    ensArr = nengo.networks.EnsembleArray(100, dim**2, ens_dimensions=1,neuron_type=neuron_type)
    #incresing num neurons has smaller effect on run time
    
    #Connect each pixel of the input to its own ensemble
    nengo.Connection(ipt,ensArr.input)
    
    '''When connecting nodes, threw error: 
    Validation error when setting 'Connection.function_info': Cannot apply functions to passthrough nodes
    This is a workaround (https://github.com/nengo/nengo/issues/805)'''
    #ensArr.output.output=lambda t, x: x

    #ensArr.add_output('rotate_output', lambda x: rotate_img(x,1))
    #a.add_output('square', lambda x: x**2)
    
    ensArr.add_output('rotate_output', function=rotator2)
    
    print(len(ensArr.rotate_output))
    
    n = nengo.Node(comb_func,size_in = (28**2)**2, size_out=28**2)
    
    nengo.Connection(ensArr.rotate_output,n)
    
    nengo.Connection(n,ensArr.input)
    
    #output node of ensArr brings all pixels together, connection performs the rotation and feeds into input node of ensArr
    #nengo.Connection(ensArr.output, ensArr.input, function=connection_func)
    #nengo.Connection(ensArr.rotate_output, ensArr.input)

    #Gathering output of ensArr
    probe = nengo.Probe(ensArr.output,# attr='decoded_output',#sample_every=0.001,
                       synapse=probe_synapse)


614656

In [59]:
sim = nengo.Simulator(net)

In [60]:
sim.run(run_time) #LIF 2:25 #Direct 0:25


Simulation finished in 0:10:14.                                                 

In [ ]:
#Original image
pylab.imshow(np.reshape(img, (dim,dim), 'F'), cmap='Greys_r')

In [ ]:
'''Image at stop time'''
pylab.imshow(np.reshape([0. if x < 0.00001 else x for x in sim.data[probe][int(stop_time*1000)-1]], 
                             (dim, dim), 'F'), cmap=plt.get_cmap('Greys_r'),animated=True)

In [ ]:
'''Image at start time'''
pylab.imshow(np.reshape([0. if x < 0.00001 else x for x in sim.data[probe][1]], 
                             (dim, dim), 'F'), cmap=plt.get_cmap('Greys_r'),animated=True)

In [61]:
'''Animation for Probe output'''
fig = plt.figure()

def updatefig(i):
    im = pylab.imshow(np.reshape([0. if x < 0.00001 else x for x in sim.data[probe][i]],
                                 (dim, dim), 'F'), cmap=plt.get_cmap('Greys_r'),animated=True)
    
    return im,

ani = animation.FuncAnimation(fig, updatefig, interval=1, blit=True)
plt.show()

In [ ]:
# save the output
#cPickle.dump(sim.data[probe], open( "Buffer_rotations_in_connection_ensemble_array_direct.p", "wb" ) )
#cPickle.dump(sim.data[probe], open( "Buffer_rotations_in_connection_ensemble_array_LIF_100_stim0.01.p", "wb" ) )