I’m trying to access the synaptic variables of a specific subset of neurons, but I’m having a hard time identifying all the synapse objects that are efferent from these. How would one find the synapse objects efferent from specific neurons in a population?
assuming you have a single group of synpases S, you could look for those which come from neuron neuron_id with
S[neuron_id, :]
accessing synaptic variables
that approach should also work for a list of neuron_ids.
Another option is to use Brian’s string-based syntax for indexing and the S.i property which stores the index of the presynaptic neuron:
S['i==5'] #simple enough for hardcoded values
S[f'i=={neuron_id}'] #not sure whether this string substitution is necessary
S['i<5'] #string-based syntax allows for additional range selection
Got it, thanks!
print('Synapses efferent from 2 to 4\n\t', S[2:4, :])
print('Synapses afferent to 2 to 4\n\t', S[:, 2:4])
print('Value of w for synapses efferent from neurons 2 to 4\n\t', S.w[2:4, :])
print('Value of w for synapses afferent to neurons 2 to 4\n\t', S.w[:, 2:4])
I was also able to index non-contiguous neurons this way:
spiked_idx = np.where(spiked_neurons)[0]
synapses_spiked = synapse_monitor[spiked_idx, :]
Another question I have is how I could slice the monitors to only return data within a certain timeframe. For example, how would I return the values recorded between 1 and 2 seconds? I can do something like synapse_monitor.t[(synapse_monitor.t>2*second)&(synapse_monitor.t<3*second)] but the synaptic variables themselves don’t have any time information attached, so I wouldn’t be sure of how to filter them by time.
I guess I could do np.where((synapse_monitor.t>2*second)&(synapse_monitor.t<3*second)) to get the indices and then access the variable arrays using those, but is there a cleaner way?
1 Like
i would have to double check this works, but yeah I do something along the lines of what you’ve outlined: (I don’t think the np.where is necessary)
#generally the time vector is identical across monitors, but this assumes they're recording at the same rate:
t = synapse_monitor.t
t_range = (t>t0) & (t<t1) # should be an array of booleans which can be used as an index
voltage_monitor.v[ t_range ]
the detail here I haven’t double checked is you just need to make sure you’re indexing the right dimension (i.e. time and not neuron index).
Great, thanks!
I have another query that I hope you can help me with:
what if I wanted to get the pairs of neurons connected within a subgroup?
Let’s say I want to get all the connected pairs within the first 64 neurons, I’m trying something like:
for i, j in zip( S.i[ :64, : ], S.j[ :, :64 ] ):
print( I,'->' ,j )
but I’m not at all sure that this is correct.
And if I wanted to get all the connected pairs with a pre-synaptic neuron inside the first 64 subgroup? Would this be correct?:
for i, j in zip( S.i[ :64, : ], S.j[ :, : ] ):
print( I,'->' ,j )
I’m trying to build a graph from the network in order to be able to visualise it in something like NetworkX or Dash. Any suggestions of the best way of approaching this? I’d like to be able to interactively see what the neuronal topology looks like.
Hi @Tioz90 . No, this is not quite correct. The expression zip( S.i[ :64, : ], S.j[ :, :64 ] ) will give you the source indices of all synapses originating from the group, and the target indices of all synapses targetting the group. But the two lists of indices will not correspond to each other – the first will also include entries for synapses coming from the group but going outside the group, and the second will include entries for synapses coming from outside the group and going inside. In general, you will have to use the same indexing for i and j to get two arrays of the same size, where individual entries correspond to the same synapses. Here’s a quick example (using a full size of 10, and a subgroup size of 5 for simplicity):
>>> G = NeuronGroup(10, '')
>>> S = Synapses(G, G)
>>> S.connect(0.5)
>>> len(S)
43
>>> list(zip(S.i, S.j)) # all synapses
[(0, 0), (0, 3), (0, 5), (0, 9), (1, 3), (1, 6), (2, 0), (2, 1), (2, 2), (2, 3), (2, 6), (3, 0), (3, 1), (3, 4), (3, 6), (3, 8), (4, 0), (4, 1), (4, 8), (5, 4), (5, 7), (5, 9), (6, 1), (6, 2), (6, 3), (6, 4), (6, 7), (6, 8), (6, 9), (7, 1), (7, 2), (7, 4), (7, 5), (7, 7), (8, 0), (8, 1), (8, 3), (8, 5), (8, 8), (9, 3), (9, 4), (9, 6), (9, 7)]
Now, using your proposed loop would give us the following source and target indices:
>>> S.i[:5, :]
array([0, 0, 0, 0, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 4, 4, 4],
dtype=int32)
>>> S.j[:, :5]
array([0, 3, 3, 0, 1, 2, 3, 0, 1, 4, 0, 1, 4, 1, 2, 3, 4, 1, 2, 4, 0, 1,
3, 3, 4], dtype=int32)
You can immediately see that this does not work, since the arrays do not have the same length. Also, if you look at the first four entries of the two arrays, they seem to imply that there are two connections 0→3, and two connections 0→0. Instead, you get the correct result if you use the same indexing for both S.i and S.j:
>>> S.i[:5, :5]
array([0, 0, 1, 2, 2, 2, 2, 3, 3, 3, 4, 4], dtype=int32)
>>> S.j[:5, :5]
array([0, 3, 3, 0, 1, 2, 3, 0, 1, 4, 0, 1], dtype=int32)
This is identical to the list of all synapses, by filtering out everything that starts or ends outside of the group. To verify, we can compare Brian’s filtering and a manual filtering in Python:
>>> list(zip(S.i[:5, :5], S.j[:5, :5]))
[(0, 0), (0, 3), (1, 3), (2, 0), (2, 1), (2, 2), (2, 3), (3, 0), (3, 1), (3, 4), (4, 0), (4, 1)]
>>> [(i, j) for i, j in zip(S.i, S.j) if i < 5 and j < 5]
[(0, 0), (0, 3), (1, 3), (2, 0), (2, 1), (2, 2), (2, 3), (3, 0), (3, 1), (3, 4), (4, 0), (4, 1)]
The same holds when you want to look at neurons originating from the subgroup, you are interested in the source and target indices from the same synapses, so they have to use the same indexing:
>> list(zip(S.i[:5, :], S.j[:5, :]))
[(0, 0), (0, 3), (0, 5), (0, 9), (1, 3), (1, 6), (2, 0), (2, 1), (2, 2), (2, 3), (2, 6), (3, 0), (3, 1), (3, 4), (3, 6), (3, 8), (4, 0), (4, 1), (4, 8)]
Finally, note that instead of indices, you can also use an actual subgroup as an index:
>>> subgroup = G[:5]
>>> list(zip(S.i[subgroup, subgroup], S.j[subgroup, subgroup])) # within the group
[(0, 0), (0, 3), (1, 3), (2, 0), (2, 1), (2, 2), (2, 3), (3, 0), (3, 1), (3, 4), (4, 0), (4, 1)]
>>> list(zip(S.i[subgroup, :], S.j[subgroup, :])) # originating from group
[(0, 0), (0, 3), (0, 5), (0, 9), (1, 3), (1, 6), (2, 0), (2, 1), (2, 2), (2, 3), (2, 6), (3, 0), (3, 1), (3, 4), (3, 6), (3, 8), (4, 0), (4, 1), (4, 8)]
I think any tool for network visualization should be happy with input like the one above, i.e. a lists of source and target indices or pairs of (source, target) indices. Sometimes a dictionary structure can also be handy, i.e. a mapping from source indices to all its targets (here created for all connections within the subgroup):
>>> import collections
>>> conn_dict = collections.defaultdict(list)
>>> for i, j in zip(S.i[subgroup, subgroup], S.j[subgroup, subgroup]):
... conn_dict[i].append(j)
...
>>> print(dict(conn_dict))
{0: [0, 3], 1: [3], 2: [0, 1, 2, 3], 3: [0, 1, 4], 4: [0, 1]}
For simple visualizations, you can also have a look at the plotting code in the Synapses tutorial, and at the plot functions in the brian2tools package.
Got it, even though it can all be a bit confusing. Basically, S.i[:5, :5] means take all synapses where the pre and post neurons have id in range(0,5) and give me the pre synaptic neurons of those synapses. Correct?
Adding some explicit methods like Synapses.get_pre_neurons() and Synapses.get_post_neurons() - instead of relying on slicing - could go a long way. 
I think the classic adjacency matrix view that brian2tools is a great general-purpose high-level view.
If you’re okay focusing on a small (10-node) subset of nodes you could take a look at my (work in progress) circuit visualizer
[code and list of keyboard shortcuts here]
[live in-browser demo here]
marcel’s code for generating lists of pairs, or a dictionary of all outputs for each input node should work in the textbox entry at the bottom… but unfortunately integer indices will need to be replaced with letters of the alphabet (let me know if you’re interested in a snippet of python to help generate that format)
otherwise you should just need to import that output list into a networkx graph, and they you can display the network interactively with plotly / dash:
if you run into any roadblocks in doing that translation let me know and I’ll try to help
Very nice, thanks! I’ll look into it, but I’m currently trying to look at larger networks than that. I would definitely like to receive some tips on the best way of representing and characterising in a visual way the connectivity of large(-ish) networks (~300 nodes for now).
I’ve currently had some luck with PyVis (see pic)but I’m open to trying other tools. Would you particularly recommend PlotLy?
Plotly is great for general purpose interactive plots from python. [example from delayed gaussian connections in brian]
But PyVis looks like it has all the features you’d want from networkx + plotly anyway, so if PyVis is working well for you, it’s not worth switching in my opinion.
It looks like your network is quite densely connected which can make these plots quite hard to follow. My initial guess as to some tricks which might help include:
- sorting / ordering the node position based on any clustering / subnetwork structure
- color coding edges by node
- trying an “adjacency matrix” view like those in brian2tools
you might also find some more inspiration in the networkx gallery here:
https://networkx.org/documentation/stable/auto_examples/index.html
1 Like
One thing I’m still not clear about is why looking at the synapses between two subpopulations of E and I brings in a ton of neurons which aren’t supposed to be there.
The synapse I_E is connecting populations I to E, and i_neurons and e_neurons are indices of neurons defining a subpopulation within those two.
I believe this operation should give the empty set, but instead it’s returning many neurons in E, that are outside the subpopulation e_neurons:
set(I_E.j[ i_neurons, e_neurons ]).difference(set(e_neurons))
And, in fact, my graph shows many nodes which shouldn’t be there, as they’re not in e_neurons (all those which are light blue should not be there):
It seems like they’re all centered around ‘i_24’, but calling
net.I_E.j[ 24,e_neurons ] gives [109 159 252], which are in e_neurons.
Hi. Could you give more information on how your populations, synapses, and subpopulations are defined? Or ideally a minimal example. If I run what I assume to be a (highly simplified) version of your code, all seems to be working as expected:
I = NeuronGroup(10, '')
E = NeuronGroup(40, '')
I_E = Synapses(I, E)
I_E.connect(p=0.5)
i_neurons = np.arange(3, 8)
e_neurons = np.arange(20, 30)
print(list(zip(I_E.i[i_neurons, e_neurons], I_E.j[i_neurons, e_neurons])))
assert len(set(I_E.j[i_neurons, e_neurons]).difference(set(e_neurons))) == 0
[(3, 21), (3, 22), (3, 23), (3, 25), (3, 28), (4, 20), (4, 21), (4, 23), (4, 24), (4, 27), (5, 20), (5, 21), (5, 23), (5, 25), (5, 26), (5, 28), (6, 20), (6, 21), (6, 23), (6, 24), (6, 26), (6, 29), (7, 21), (7, 25), (7, 26), (7, 27)]
Unfortunately I wouldn’t know how to do that, without sharing a ton of code.
It really seems like it’s including post-synaptic neurons which shouldn’t be there as:
e_neurons = [ 44 47 53 0 3 159 103 139 109 119 185 200 251 234 252]
i_neurons = [24 12 58 1 38 39 23 46 24 17 37 25 13 8 9]
and
net.I_E.i[ i_neurons, e_neurons ] = [ 1 1 1 1 8 8 8 8 9 9 9 12 12 12 17 17 17 17 23 23 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 24 25 25 25 25 37 37 38 38 38 39 39 46 46 46 58 58]
net.I_E.j[ i_neurons, e_neurons ] = [ 47 103 109 200 3 119 234 251 44 109 251 119 200 252 3 139 251 252 185 200 2 7 10 17 18 20 21 27 31 34 36 41 50 54 58 61 64 65 78 79 81 83 88 92 96 100 112 115 124 125 128 131 134 141 148 151 152 157 159 159 160 167 171 176 177 180 188 191 192 194 204 212 219 222 225 226 229 236 238 241 245 247 253 0 3 185 251 47 185 53 109 159 103 159 53 103 200 3 200]
but
set(net.I_E.j[i_neurons, e_neurons]).difference(set(e_neurons)) = {2, 7, 10, 17, 18, 20, 21, 27, 31, 34, 36, 41, 50, 54, 58, 61, 64, 65, 78, 79, 81, 83, 88, 92, 96, 100, 112, 115, 124, 125, 128, 131, 134, 141, 148, 151, 152, 157, 160, 167, 171, 176, 177, 180, 188, 191, 192, 194, 204, 212, 219, 222, 225, 226, 229, 236, 238, 241, 245, 247, 253}
instead of being the empty set.
Thanks for the additional info, I will have a closer look on Monday. I don’t think the pickled synapses would help, but could you give more info on how you created the synapses, in particular whether they use full NeuronGroups as their sourse/target or whether they use subgroups? Thanks!
I managed to reproduce some of your issue, there definitely is something wrong here. I will try to get to the bottom of it this afternoon/evening.
As a workaround, you can get the indices of the synapses you are interested in manually via:
syn_indices = np.arange(len(net.I_E))[np.isin(net.I_E.i, i_neurons) & np.isin(net.I_E.j, e_neurons)]
You can then use 1-dimensional indexing on whatever you are interested in, e.g. to get the all the target neurons:
net.I_E.j[syn_indices]
Thanks!
For now, I worked around the issue by doing:
[ (i, j) for i, j in zip( net.I_E.i, net.I_E.j ) if i in i_neurons and j in e_neurons ]
1 Like
Ah, I finally found out what happened. The algorithm we use for determining the synaptic indices when you index with [source_indices, target_indices] is quite simple/inefficient. We basically do:
correct_pre = []
for index in source_indices:
correct = np.flatnonzero(S.i == index)
correct_pre.extend(correct)
correct_post = []
for index in target_indices:
correct = np.flatnonzero(S.j == index)
correct_post.extend(correct)
And finally:
synapses = np.intersect1d(correct_pre, correct_post, assume_unique=True)
The key here is the assume_unique=True which makes the intersect1d function a bit faster. As it turns out, this assumption is violated here, since your i_neurons list contains the index 24 twice. The documentation says " If True but ar1 or ar2 are not unique, incorrect results and out-of-bounds indices could result." which, well, is what you are seeing. We should clearly handle this better, but I guess having the duplicated index there is not something you did intentionally either.
2 Likes
Great!
I actually had noticed I was getting duplicates in my neuron sampling, so I switched to using np.random.choice().
The Synapse indexing issue had also disappeared around that time, but I hadn’t connected the two!
2 Likes