Boolean indexing demystified: an __array__ case!

Support
1

Description of problem

This issue involves indexing synapses. There seems to be a deviation from numpy conversion. Please look at the script below:

Minimal code to reproduce problem

import brian2 as b2
import numpy as np
import matplotlib.pyplot as plt

G = b2.NeuronGroup(1000, model='v:1')
S = b2.Synapses(G, model='w:1')
S.connect(p=0.1)
S.w = 0 # optional, it is anyway zero-initialized


idxs = S.j== 5
print(idxs.shape, sum(idxs)) # expected values: (N**2 *p, N*p)  --> correct

# Here comes the weird part
S.w[idxs] = 1 # expecting N*p non-zero element

# however
print('\n\nOnly {} elements actually changed their value'.format(sum(S.w))) # always 2 !
plt.figure()
plt.plot(S.w, label='S.w') 
plt.legend()

print('Also, strangely both S.w and S.w[idx] have equal shapes:')
print('Shape of S.w: {}'.format(S.w.shape))
print('Shape of S.w[idxs]: {}'.format(S.w[idxs].shape)) 

#plt.plot(S.w[idxs], label='S.w[idxs]')


# But somehow this can be changed. Let's reset everything
S.w = 0
S.w.__array__()[idxs] = 1  # added __array__()

print('\n\n{} are correctly changed their value'.format(sum(S.w))) 

plt.figure()
plt.plot(S.w, label='S.w') 
plt.legend()
 
print('Yet, the shapes are still equal:')
print('Shape of S.w: {}'.format(S.w.shape))
print('Shape of S.w[idxs]: {}'.format(S.w[idxs].shape))

with the following output:

(99897,) 82 # might be different for you...


Only 2.0 elements actually changed their value
Also, strangely both S.w and S.w[idx] have equal shapes:
Shape of S.w: (99897,)
Shape of S.w[idxs]: (99897,)


82.0 are correctly changed their value
Yet, the shapes are still equal:
Shape of S.w: (99897,)
Shape of S.w[idxs]: (99897,)

which shows that boolean indexing have the following mysterious behaviors:

  1. If idxs is a boolean array of size N with n True elements, and S a Synapse object N synapse containing variable w, then S.w[idx] does not have size n but N. Why?
  2. The operation S.w[idxs]=<whatever_value>, changes only the value of the first two components independent of the idxs. Why?
  3. Why using __array__() on w fixes the issue?
2

Hi @arashgmn . Thanks for this report, we need to handle this better.
I will need to look into the technical reasons in more detail, but we are actually not supporting boolean indexing for synaptic variables in the first place! The documentation is not really clear about this, we are just saying:

Note that it is also possible to index synaptic variables with a single index (integer, slice, or array), but in this case synaptic indices have to be provided.

(Synapses — Brian 2 2.5.4 documentation)

So, why doesn’t your code raise an error and instead does something weird? The reason is that the boolean variables are interpreted as 0 or 1. So your

S.w[idxs] = 1

is interpreted as something like

S.w[[0, 1, 1, 0, 1, 0, 0, 0, 1, 0, 1, 1, 1, …]] = 1

I.e. it will only set the weights with indices 0 and 1 (but repeatedly…).

Using S.w[...] goes through quite a bit of machinery, since we support indexing in many ways that are not supported by numpy, e.g. with string expressions. At the moment, string expressions are the only way to do “boolean” indexing, i.e. your example would work as:

S.w["j == 5"] = 1

Under the hood, the weights are stored as a simple numpy array, though, and for efficiency reasons, we allow direct access to this array without copying. This is why all the following do work:

S.w.__array__()[idxs] = 1
S.w[:][idxs] = 1
np.asarary(S.w)[idxs] = 1

Again, I agree that we should fix this, but hopefully the explanations make things less mysterious :ghost:

3

Thank you, @mstimberg! As always, very helpful response. :pray:

1 Like