Network dynamics visualizer?

[galen]

Hi,

Has anyone made a visualizer, using this as a backend, that shows the network and states of the nodes over time?

thanks in advance

[adam]

Hi @galen,
This is something I’ve been daydreaming about for a while, but haven’t been able to work on much recently. I’d start by checking out this post, which is perhaps the closest to what you mentioned:

Brian running locally in the browser via Javascript using pyodide

Check out this cool demo of Brian running in the browser using Javascript. It uses the Pyodide package. For the moment it’s a bit clunky but we think this could be a cool low maintenance tool for teaching Brian or computational neuroscience.

Followed perhaps by the interactive jupyter notebooks in the Brian docs:
https://brian2.readthedocs.io/en/stable/resources/tutorials/index.html?highlight=interactive#interactive-notebooks-and-files

Outside of Brian, there have been several other interesting network visualizations, including this nice “scrollytelling” blogpost / mini textbook chapter page by Patrick Mineault (again, not Brain, but simple LIF code, and a good template for what could be done)

Balanced excitation and inhibition

and then I put together an interactive network / directed graph / adjacency matrix visualizer which has simple linear gaussian dynamics

• demo: https://awillats.github.io/circuit-visualizer-p5/ (seems to work in safari but not chrome for me at the moment. Hit the w key to start activity, see the README.md for other commands)
• github page: GitHub - awillats/circuit-visualizer-p5: (Experimental) visualizer for understanding networks and properties of their adjacency matrices

I’d love to use this as a interface / visualization for Brian networks some day – and I think something like pyodide / pyscript would be a great way to do that

• introduction to pyscript - includes comparison with pyodide
• A First Look at PyScript: Python in the Web Browser

[galen]

Hey thanks Adam,

All sounds interesting!

Here’s approximately why I’m interested by the way:

https://www.nature.com/articles/s41598-022-19218-0

https://youtu.be/WyAspVjo6VI

Above links are to the paper and a talk about this research. (Starts at 1m35s)

We show how random threshold networks can compute complex Boolean functions in cascades or avalanches. This has many implications for neuroscience and other domains, and may help in discovering more efficient methods for learning in artificial networks.

Best,
Galen

[galen]

Today did a bit of digging around.

Seems these days the best option for coding this oneself is d3 JavaScript.

I’d like to be able to click and drag nodes and edges, then perturb them in various ways and watch it happen, perhaps in slow motion.

More generally it would be nice if the activation rules could even be kuramoto or other, and if a variety of types of information transfer along the edges could happen.

You can learn a lot from very simple toy models.

Then on the side you can plot other state aspects of the network.

It’s 2022! Seems long overdue

Best,
Galen

[adam]

That does sound like an interesting application!

You can get just about anything you want to do with plotting and interactive visualization done with d3, it’s very flexible and powerful, but (in my opinion) it’s quite difficult to learn and put together something with d3 from the ground up.

• thorough tutorial from 2015: Introduction to D3

Depending on your prior experience with javascript / web development, I would recommend instead building a prototype with a “higher level” set of tools. The circuit visualizer I linked above uses Processing / p5.js which is very beginner-friendly and handles most of the hassle of setting up and managing state for the drawing canvas and mouse interactions, and can be quickly prototyped in a live web editor. I was able to set up dragging nodes around and adding and removing connections in a network with relatively little code!

In fact, here’s a live sandbox example of what code for dragging a circle around in p5js looks like:

p5.js Web Editor

It’s also worth considering solutions which stay inside python. Your options for interactivity are more limited, but you don’t have to keep passing data back and forth between javascript and python. There are several interactive visualization libraries in python, including some which are essentially python interfaces to d3.js behind the scenes:

IPython Cookbook - 6.4. Visualizing a NetworkX graph in the Notebook with D3.js

Here are some more resources on various combinations of python and javascript for interactive data science:

• Dynamic scientific visualizations in the browser for Python users also by Patrick Mineault
• Interactive Visualization in Science: Resources and advice for Python + JavaScript by me

[galen]

Thanks, well I have a PhD in computer science, so no get out of jail card for me as far as learning what I need to.

It seems d3 automatically lets you drag things around, plus the web interface, so probably worth the trouble.

Also looked at Julia today, which may be a good option due to speed.

Best

[mstimberg]

Hi everyone. As @adam knows, I’m very much interested in this general topic as well
@galen: Just to make sure I understand your requirements: are you only interested in visualizing the results of existing simulations, or do you want to interactively trigger new simulations based on user input?

[galen]

Hi @mstimberg, for starters, a wrapper of the Brian simulator would be nice.

However, I think it would also be very nice to allow simple threshold or Kuramoto activity at each node. I do think for students and researchers, toy models can be incredibly useful, and it would be nice to watch their activity, probably in slow motion.

Personally I know that I can learn a great deal from incredibly simple models, then learn even more by very gradually making them more complex.
(Of course even the simplest activity can have very complex emergent properties both statically and dynamically when scaled to 10^9 or 10^10 interacting ‘units’.

I.e. simple stories can already very compelling. That is, I can learn a great deal from drawing a network of 3 neurons interacting, especially if I can impose a variety of ways of interacting and activation at each ‘node’.

(Going the other direction, based on recent research it seems that individual neurons are much more complex than previously thought (i.e. they are not ‘nodes’).

So eventually that needs to be addressed in these simulators.)

Best

[mstimberg]

@galen

for starters, a wrapper of the Brian simulator would be nice.

Can you explain what you mean by this?

@galen

However, I think it would also be very nice to allow simple threshold or Kuramoto activity at each node. I do think for students and researchers, toy models can be incredibly useful, and it would be nice to watch their activity, probably in slow motion.

I am a bit confused, you mean that Brian should support simpler models? It does already support rate models or coupled oscillators of the Kuramoto style, the only limitation is that you cannot have delayed continuous interactions.

@galen

Personally I know that I can learn a great deal from incredibly simple models, then learn even more by very gradually making them more complex.

I agree with that, but I am not quite sure what you want to say here with respect to changes needed for Brian or with respect to interactive visualizations in general?

[galen]

As a side note, it seems that Julia could be a useful platform for this, since it’s significantly faster than Python and could scale up, and has some web interfaces as well, and seems to be growing as a community.

The good thing is that the network only needs to be laid out once, and the dynamics do not have to grow in the size of the entire network (each node does not activate O(N) nodes, so I think dynamics for pretty large networks can be visualized with today’s platforms, esp. using GPUs, etc.

So I think one can build large enough brain simulations to learn a great deal if one slows down and thinks about what things mean. For example, the combinatorics of the Boolean function space are massive. (Again, see paper link.)

See you at braincriticality.org ?

[galen]

Not the best way to have a discussion, it seems.

I’m new to Brian, a quick search for Kuramoto in the docs didn’t yield anything.

Happy to discuss by more human means.

[mstimberg]

@galen no worries, I did not intend you to drag you into a discussion I guess I somewhat misread your earlier posts as being about something specific that you were missing from Brian or that you wanted to build on top of Brian, but they were probably meant to be more general.

In either case, please feel free to drop a note / open a new post if you want to discuss anything specific that you suggest we should do on the Brian side. Or if you have anything else to share (some interactive visualization, for example) that you think could be relevant to the community here, even if it doesn’t use Brian.

@galen

I’m new to Brian, a quick search for Kuramoto in the docs didn’t yield anything.

That’s true, we don’t have any example of this kind – Brian’s focus is on biological spiking neural networks. I think we only have one example for a network of non-spiking neurons (Example: Tetzlaff_2015 — Brian 2 2.5.4.post90 documentation), but this example is already very complex. I’ll add an example for simple coupled oscillators soon, just to show that you can implement them quite easily in Brian.

[rth]

Examples also have non-spiking Jansen and Rit 1995 model, just in case

[galen]

No worries, great job on building it!

The internet can make brutes of us all, so I try to be aware. :^]

Yes I’m coming from a more basic and abstract set of research questions: how do interactions process information (doesn’t have to be neurons, could just be banging on the table).

In the case, by “wrapper” I meant GUI, especially showing the node states and edges.
Imagine running the network, slowing it down, and even reversing in time, clicking and dragging new connections, interrogating or editing node states and rules, etc etc etc.

Basically network vis. in the 21st century.
Exhausted now and moving countries, so unfollowing for the moment.

Best and keep up the good work!
8’}

[mstimberg]

For future reference: Brian’s documentation now has an example of the Kuramoto model (and the plot is animated )
https://brian2.readthedocs.io/en/latest/examples/coupled_oscillators.html

[galen]

That’s great Marcel.

Suppose I want to

1. create a network in networkx using the Synapse connectivity
2. OR create the network in networkx, then make the Synapses connect that way
3. then build this same Kuramoto model using that network topology?
4. And then display the states of the network as they change?

Here is some code to animate networks states, just needs Brian network and node states:
thanks for any help!
best

# code to make an animation of node states over time and save as mp4

from matplotlib import pyplot as plt, animation
from matplotlib import colormaps
import networkx as nx
import numpy as np

#%matplotlib inline

plt.rcParams["figure.autolayout"] = True

num_steps = 50
N = 10

#the number of possible states of each node
num_states_per_node = 2


####################### REPLACE #######################
# REPLACE A WITH THE SYNAPSE CONNECTIVITY FROM BRIAN 
# create a network from a random adjacency matrix (2-D array actually)
A = np.random.randint(num_states_per_node, size = [N,N])
A = np.triu(A, k = 1)
####################### REPLACE #######################

# create the networkx graph
G = nx.from_numpy_array(A)

##### NOTE, pos (DICTIONARY OF NODE POSITIONS) COULD BE REPLACED WITH ACTUALLY NEURON LOCATIONS, IF KNOWN
pos = nx.spring_layout(G)
fig = plt.figure(figsize=[10,10])

node_states = np.random.randint(num_states_per_node, size = N)

nx.draw(G, node_color = node_states, edgecolors = 'grey', edge_color = 'lightgrey', cmap=colormaps['binary'])

def animate(frame):
    fig.clear()

    ####################### REPLACE #######################
    # REPLACE THIS LINE WITH THE STATES FROM THE BRIAN NODE STATES:
    # an array of 0s and 1s giving node states
    node_states = np.random.randint(num_states_per_node, size = N)
    ####################### REPLACE #######################

    nx.draw(G, pos, node_color = node_states, edgecolors = 'grey', edge_color = 'lightgrey', cmap=colormaps['binary'])

# create the animation
anim = animation.FuncAnimation(fig, animate, frames=num_steps, interval=1000, repeat=False)

# saving to m4 using ffmpeg writer
writervideo = animation.FFMpegWriter(fps=5)
anim.save('percolation_zoom.mp4', writer=writervideo)
plt.close()

[mstimberg]

Hi @galen

Here’s a merge of your code and the one I posted in the example. I give the connectivity matrix to Brian (it could also be a matrix of weights). The other option would be to generate the connectivity in Brian and then extract the connectivity/weight matrix – both approaches are described here: Synapses — Brian 2 2.5.1 documentation
I changed your plotting code slightly, since the Kuramoto model does not have not a binary activation but a continuous phase for each neuron:

# code to make an animation of node states over time and save as mp4

from matplotlib import pyplot as plt, animation
from matplotlib import colormaps
import networkx as nx
import numpy as np

from brian2 import *
defaultclock.dt = 1*ms


def run_sim(K, connection_matrix, runtime=10*second):
    assert connection_matrix.ndim == 2, "need 2D matrix"
    assert connection_matrix.shape[0] == connection_matrix.shape[1], "need quadratic matrix"
    N = connection_matrix.shape[0]

    eqs = '''
    dTheta/dt = omega + K/N*coupling : radian
    omega : radian/second (constant) # intrinsic frequency
    coupling : 1
    '''

    oscillators = NeuronGroup(N, eqs, method='euler')
    oscillators.Theta = 'rand()*2*pi'  # random initial phase
    oscillators.omega = 'clip(3 + randn()*1, 0, inf)*radian/second'  # 𝒩(0.5, 0.5)

    connections = Synapses(oscillators, oscillators,
                        'coupling_post = sin(Theta_pre - Theta_post) : 1 (summed)')
    sources, targets = connection_matrix.nonzero()
    connections.connect(i=sources, j=targets)

    mon = StateMonitor(oscillators, 'Theta', record=True)
    run(runtime)
    return mon.Theta[:]

#%matplotlib inline

plt.rcParams["figure.autolayout"] = True

N = 10

#the number of possible states of each node
num_states_per_node = 2

A = np.random.randint(num_states_per_node, size = [N,N])
A = np.triu(A, k = 1)

# create the networkx graph
G = nx.from_numpy_array(A)

##### NOTE, pos (DICTIONARY OF NODE POSITIONS) COULD BE REPLACED WITH ACTUALLY NEURON LOCATIONS, IF KNOWN
pos = nx.spring_layout(G)
fig = plt.figure(figsize=[10,10])

node_states = np.random.randint(num_states_per_node, size = N)

nx.draw(G, node_color = node_states, edgecolors = 'grey', edge_color = 'lightgrey', cmap=colormaps['binary'])

# run simulation with given adjacency matrix
runtime = 10*second
results = run_sim(K=15/second, connection_matrix=A, runtime=runtime)

num_steps = 250
time_between_frames = runtime/num_steps

def animate(frame):
    fig.clear()
    time_step = int(round(frame*time_between_frames/defaultclock.dt))
    node_states = results.T[time_step] % (2*pi)

    nx.draw(G, pos, node_color = node_states, edgecolors = 'grey', edge_color = 'lightgrey', cmap=colormaps['twilight'], vmin=0, vmax=2*np.pi)

# create the animation
anim = animation.FuncAnimation(fig, animate, frames=num_steps, interval=time_between_frames/ms, repeat=False)

# saving to m4 using ffmpeg writer
writervideo = animation.FFMpegWriter(fps=25)
anim.save('percolation_zoom.mp4', writer=writervideo)
plt.close()

I increased the neurons’ intrinsic frequencies a bit, and made the connection quite strong so that one sees a clear phase synchronization between the units :

PS: Discouse classifed your post as “potential spam” for some reason, so I had to manually accept to make it appear here for everyone.

[galen]

Thanks @mstimberg , that’s exactly what I was looking for.

I noticed I neglected to change the output filename. Should be ‘Kuramoto_network.mp4’ or so.

Just a thought, since networks are fundamental to these models, perhaps useful to build this in?

Feel free to use this as a canonical example in the docs or tutorials.

best

[galen]

I’ve just cleaned up the code and commented it a lot for easier understanding. Seems to work fine!

# code to make a network animation of Kuramoto node states over time and save as mp4

from matplotlib import pyplot as plt, animation
from matplotlib import colormaps
import networkx as nx
import numpy as np

from brian2 import *
defaultclock.dt = 1*ms


def run_sim(K, connection_matrix, runtime=10*second):
    '''
    run the Kuramoto simulation using a pre-determined network connectivity
    
    K is the coupling constant
    connection_matrix is the standard adjacency matrix
    runtime is the total runtime in seconds
    
    see: https://en.wikipedia.org/wiki/Kuramoto_model and the brian documentation
    '''
    
    # check for existence of 2-d square graph adjacency matrix
    assert connection_matrix.ndim == 2, "need 2D matrix"
    assert connection_matrix.shape[0] == connection_matrix.shape[1], "need quadratic matrix"
    N = connection_matrix.shape[0]

    # the Kuramoto dynamics
    eqs = '''
    dTheta/dt = omega + K/N*coupling : radian
    omega : radian/second (constant) # intrinsic frequency
    coupling : 1
    '''

    # create the oscillators
    oscillators = NeuronGroup(N, eqs, method='euler')
    oscillators.Theta = 'rand()*2*pi'  # random initial phase
    oscillators.omega = 'clip(3 + randn()*1, 0, inf)*radian/second'  # 𝒩(0.5, 0.5)

    # build the connections from the adjacency matrix
    connections = Synapses(oscillators, oscillators,
                        'coupling_post = sin(Theta_pre - Theta_post) : 1 (summed)')
    sources, targets = connection_matrix.nonzero()
    connections.connect(i=sources, j=targets)

    # monitor the node states
    mon = StateMonitor(oscillators, 'Theta', record=True)
    
    # actually run the simulation
    run(runtime)
    
    # return node states over time
    return mon.Theta[:]

def animate(frame):
    '''
    animate each frame of the simulation
    
    frame is one frame
    see matplotlib animation documentation
    '''
    
    # reset the figure
    fig.clear()
    
    # update the node states at each timestep
    time_step = int(round(frame*time_between_frames/defaultclock.dt))
    node_states = results.T[time_step] % (2*pi)

    # draw the network with node states
    nx.draw(G, pos, node_color = node_states, edgecolors = 'grey', 
            edge_color = 'lightgrey', cmap=colormaps['twilight'], vmin=0, vmax=2*np.pi)


    
    
#%matplotlib inline  

plt.rcParams["figure.autolayout"] = True

# number of nodes ('neurons')
N = 10

# random unweighted adjacency matrix
num_weights = 2
A = np.random.randint(num_weights, size = [N,N])
A = np.triu(A, k = 1)  # keep only upper triangle

# create the networkx graph
G = nx.from_numpy_array(A)

##### NOTE, pos (DICTIONARY OF NODE POSITIONS) COULD BE REPLACED WITH ACTUALLY NEURON LOCATIONS, IF KNOWN
pos = nx.spring_layout(G)
fig = plt.figure(figsize=[10,10])

# run simulation with given adjacency matrix (network connectivity)
runtime = 10*second
results = run_sim(K=15/second, connection_matrix=A, runtime=runtime)

# length and speed
num_steps = 250
time_between_frames = runtime/num_steps

# create the animation
anim = animation.FuncAnimation(fig, animate, frames=num_steps, interval=time_between_frames/ms, repeat=False)

# saving to m4 using ffmpeg writer
writervideo = animation.FFMpegWriter(fps=25)
anim.save('Kuramoto_network.mp4', writer=writervideo)
plt.close()

Here is a jupyter notebook: (change the filename to .ipynb)
Kuramoto network brian.py (106.6 KB)

[galen]

This has been great. I’ve also made code that plots the angle (phase) of each node as a vector. Will post it soon if I remember.

Suppose I want to implement something even simpler: The linear threshold model (LTM).

Here is the LTM:

• Some graph topology (usually Erdos-Renyi)
• Each node has a threshold (phi) uniformly distributed in [0,1]
• Nodes can be activated or unactivated (1 or 0)
• We set the state of all nodes to 0 unactivated

Now, to run the simulation:

• set one or a few seed nodes to 1 (to perturb the network)
• examine each unactivated node u randomly
• check the fraction of u’s graph neighbors that are labeled 1
• if that fraction >= phi, u becomes labeled
• stop when nothing more changes

Of course it is possible to implement this in matrix vector form.

The LTM has the attraction of being one of the simplest models to exhibit cascade criticality, and can also be connected theoretically to graph percolation. (See Duncan Watts “A simple model of global cascades on random networks”)

https://www.pnas.org/doi/10.1073/pnas.082090499

Similarly, I wonder if it is possible to implement the Kuramoto model in such a way that the nodes are not oscillating in the beginning, then we ‘pluck’ one or a few seed nodes so that they start oscillating, then see how this propagates through the network.

thanks for any starting points.

Best