Change the input and update weigths between runs

Support
1

Description of problem

I created spiking neural network and i have train dataLoader from which i want to take input for my network. Also after run i wnat to update weights offline and run next example with updated weights

Minimal code to reproduce problem

from mit_bih import download_mitbih

train_loader, test_loader = download_mitbih()

from brian2 import *
import numpy as np
import matplotlib.pyplot as plt
import torch
from datetime import datetime

prefs.codegen.target = "numpy"
np.random.seed(25)
start_scope()

def get_random_20_percent_indices(n_count, exc_p=0.8, p=0.2, exc_w=1, inh_w=-1):
    main_range = np.arange(n_count)
    n = int(p * len(main_range))
    idx = np.random.choice(main_range, n, replace=False)
    exc_n = int(exc_p * len(idx))
    wegths = []
    for i in range(exc_n):
        wegths.append(exc_w * np.random.rand())

    for i in range(exc_n, len(idx)):
        wegths.append(inh_w * np.random.rand())

    return idx, wegths


dt = 0.1 * ms
defaultclock.dt = dt
single_exeample_time = 500 * ms
stimulus = TimedArray([], dt=dt)
batch_ecg, _ = next(iter(train_loader))
ecg = batch_ecg[0].numpy() if isinstance(batch_ecg, torch.Tensor) else batch_ecg[0]

# Input neuron params
N_input = 1
v_rest_input = 0*volt
v_threshold_input = 3*mV
g_l_input = 0.1
tau_input = 5*ms
input_eqs = '''
dv/dt = (v_rest_input - v + I/g_l_input) / tau_input : volt
I = abs(I_previous - I_current) * volt : volt
I_previous = stimulus(t-dt) : 1
I_current = stimulus(t) : 1
'''


input_group = NeuronGroup(N_input, input_eqs, threshold='v>v_threshold_input', reset='v=v_rest_input', method='euler')
input_group.v = 'v_rest_input'

print(f"[{datetime.now()}] Input neuron group created")

# LSM neuron params
N_lsm = 1000
N_lsm_exc = int(0.8 * N_lsm)
N_lsm_inh = N_lsm - N_lsm_exc
v_th_lsm = -55*mV
v_rs_lsm = -65*mV
tau_lsm = 5*ms
g_l_lsm = 1
refreac_lsm = 4*ms

lsm_eqs = '''
dv/dt = (v_rs_lsm - v) / tau_lsm : volt
'''

lsm_group = NeuronGroup(N_lsm, lsm_eqs, threshold='v>v_th_lsm', reset='v=v_rs_lsm', refractory=refreac_lsm, method='euler')
lsm_group.v = 'v_rs_lsm'

print(f"[{datetime.now()}] LSM neuron group created")

import numpy as np

def get_random_20_percent_indices(n_count = N_lsm, exc_p = 0.8, p = 0.2, exc_w = 1, inh_w = -1):
    main_range = np.arange(n_count)
    n = int(p * len(main_range))
    idx = np.random.choice(main_range, n, replace=False)
    exc_n = int(exc_p * len(idx))
    wegths = []
    for i in range(exc_n):
        wegths.append(exc_w*np.random.rand())

    for i in range(exc_n, len(idx)):
        wegths.append(inh_w*np.random.rand())

    return idx, wegths

w_i_lsm = 0.6
S_input_lsm = Synapses(input_group, lsm_group, 'w : 1', on_pre='v_post += 10*w*mV')
indeces, weigths = get_random_20_percent_indices(exc_w=w_i_lsm, inh_w=w_i_lsm)
for i in indeces:
    S_input_lsm.connect(i=0, j=i)

S_input_lsm.w = weigths

print(f"[{datetime.now()}] Input -> LSM synapses created")

S_lsm_lsm = Synapses(lsm_group, lsm_group, 'w : 1', on_pre='v_post += w*mV')
weigths = []
for i in range(N_lsm):
    post_indices, post_weights = get_random_20_percent_indices(exc_w=0.3, inh_w=-0.5)
    S_lsm_lsm.connect(i=i, j=post_indices)

    for indx in range(len(post_indices)):
        weigths.append(post_weights[indx])
    # for indx in range(len(post_indices)):
    #     if i != post_indices[indx]:
    #         S_lsm_lsm.connect(i=i, j=post_indices[indx])
    #         weigths.append(post_weights[indx])

S_lsm_lsm.w = weigths

print(f"[{datetime.now()}] LSM -> LSM synapses created")

# Output neuron params
N_output = 5
v_th_out = -55*mV
v_rs_out = -65*mV
tau_out = 5*ms
g_l_out = 1

out_eqs = '''
dv/dt = (v_rs_out - v) / tau_out : volt
'''

output_group = NeuronGroup(N_output, out_eqs, threshold='v>v_th_out', reset='v=v_rs_out', method='euler')
output_group.v = 'v_rs_out'

print(f"[{datetime.now()}] Output neuron group created")

S_lsm_out = Synapses(lsm_group, output_group, 'w : 1', on_pre='v_post += w*volt', name = 'S_lsm_out')
weigths = []
for f in range(N_output):
    post_indices, post_weights = get_random_20_percent_indices(exc_w=0.3, inh_w=-1)
    for indx in range(len(post_indices)):
        S_lsm_out.connect(i=post_indices[indx], j=f)
        weigths.append(post_weights[indx])

S_lsm_out.w = weigths

print(f"[{datetime.now()}] LSM -> output synapses created")

input_group_ST = StateMonitor(input_group, variables=True, record=True)
input_group_SP = SpikeMonitor(input_group, record=True)

lsm_group_ST = StateMonitor(lsm_group, variables=True, record=True)
lsm_group_SP = SpikeMonitor(lsm_group, record=True)

output_group_ST = StateMonitor(output_group, variables=True, record=True)
output_group_SP = SpikeMonitor(output_group, record=True)

print(f"[{datetime.now()}] Monitors created created")

network = Network([input_group, lsm_group, output_group, S_input_lsm, S_lsm_lsm, S_lsm_out, input_group_ST, input_group_SP, lsm_group_ST, lsm_group_SP, output_group_ST, output_group_SP])
network.run(0*ms)
print(f"[{datetime.now()}] Network built")

def print_nn():
    fig, axs = plt.subplots(nrows=3, ncols=1)

    axs[0].plot(input_group_SP.t/ms, input_group_SP.i, '.')
    axs[0].set_title('Time (ms)')
    axs[0].set_xlabel('Neuron index')
    axs[0].set_ylabel('Input Neuron Spikes')

    axs[1].plot(output_group_SP.t/ms, output_group_SP.i, '.')
    axs[1].set_title('Time (ms)')
    axs[1].set_xlabel('Neuron index')
    axs[1].set_ylabel('Output Neuron Spikes')

    axs[2].plot(lsm_group_SP.t/ms, lsm_group_SP.i, '.')
    axs[2].set_title('Time (ms)')
    axs[2].set_xlabel('Neuron index')
    axs[2].set_ylabel('LSM Neuron Spikes')

    plt.show()

def train_network(target_idx):
    target_ratio = 0.8
    out_counts = np.bincount(output_group_SP.i, minlength=N_output)
    total_spikes = out_counts.sum()
    ratio = out_counts[target_idx] / total_spikes if total_spikes > 0 else 0
            
    if (ratio >= target_ratio):
        print(f"[{datetime.now()}] Target ratio {ratio} reached, no weight update needed.")
        return
    
    for j in range(N_output):
        for i in range(N_lsm):
            if len(network['S_lsm_out'].w[i, j]) == 0:
                continue

            w = network['S_lsm_out'].w[i, j][0]
            if j == target_idx:
                delta = 0.01 * (- (ratio - target_ratio))
                if w > 0:
                    delta = +abs(delta)
            else:
                suppression = out_counts[j] / total_spikes if total_spikes > 0 else 0
                delta = -0.01 * suppression
                if w < 0:
                    delta = -abs(delta)
            new_w = w + delta
            network['S_lsm_out'].w[i, j] = new_w

    print(f"[{datetime.now()}] Output counts: {out_counts}, total spikes: {total_spikes}, ratio for target {target_idx}: {ratio}")

for batch_idx, (inputs, targets) in enumerate(train_loader):
    for idx, input in enumerate(inputs):

        print(f"[{datetime.now()}] Processing batch {batch_idx}, sample {idx}")

        ecg = input.numpy() if isinstance(inputs, torch.Tensor) else input
        stimulus = TimedArray(ecg.copy(order='C'), dt=dt)
        
        network.run(single_exeample_time)

        train_network(targets[idx].item())  
        # print_nn()

    break

What you have aready tried

Expected output (if relevant)

Actual output (if relevant)

Full traceback of error (if relevant)

2

Hi @dmytr0 . Could you be a bit more specific about what you want to do and what you tried?

3

Hi @mstimberg

I want to train synapse weights between LSM and Output layer with my custom offline method. I don’t want to use STDP or something like this. Purpose of training loop: get train example, run it, update weights, and the next example should be run with updated weights from previous step.

Firstly, I tried network.store(‘init’) before the training loop. In this case stimulus was updated. But weights weren’t updated for the next example run. The final code was like this

network = Network([input_group, lsm_group, output_group, S_input_lsm, S_lsm_lsm, S_lsm_out, input_group_ST, input_group_SP, lsm_group_ST, lsm_group_SP, output_group_ST, output_group_SP])
network.run(0*ms)

network.store("init")
for batch_idx, (inputs, targets) in enumerate(train_loader):    
    for idx, input in enumerate(inputs):
        print(f"[{datetime.now()}] Processing batch {batch_idx}, sample {idx}")
        
        ecg = input.numpy() if isinstance(inputs, torch.Tensor) else input
        stimulus = TimedArray(ecg.copy(order='C'), dt=dt)
        network.run(single_exeample_time)
        train_network(targets[idx].item())
        network.restore("init")
             
break

Then I tried with devices. I tried set_device('cuda_standalone', build_on_run=True) , but on building device.build(run=False) all my RAM was used, and the further code wasn’t executed on Mac ( I have MacBook on M1 chip). I also tried to execute it on Colab Pro environment with High RAM (52 Gb) and workspace crashed due to using all RAM.

Here is code of creating train and test data loadres.

import wfdb
import numpy as np
from sklearn.model_selection import train_test_split
import torch
from torch.utils.data import DataLoader, TensorDataset
from scipy.interpolate import interp1d
import os
from datetime import datetime

def download_mitbih() -> tuple[DataLoader, DataLoader]:
    """Download the MIT-BIH Arrhythmia Database if not already present."""
    # Define the class mapping (AAMI classes: 0-N, 1-S, 2-V, 3-F, 4-Q)
    class_map = {
        'N': 0, 'L': 0, 'R': 0, 'e': 0, 'j': 0,  # Normal
        'A': 1, 'a': 1, 'J': 1, 'S': 1,  # Supraventricular
        'V': 2, 'E': 2,  # Ventricular
        'F': 3,  # Fusion
        '/': 4, 'f': 4, 'Q': 4  # Unknown/Paced/Q
    }

    # List of MIT-BIH records (excluding paced ones: 102,104,107,217)
    records = [
        100, 101, 103, 105, 106, 108, 109, 111, 112, 113, 114, 115, 116, 117, 118, 119,
        121, 122, 123, 124, 200, 201, 202, 203, 205, 207, 208, 209, 210, 212, 213, 214,
        215, 219, 220, 221, 222, 223, 228, 230, 231, 232, 233, 234
    ]

    # Download the MIT-BIH database if not already present
    db_path = 'mitdb'
    if not os.path.exists(db_path):
        wfdb.dl_database('mitdb', db_path)

    # Parameters for beat extraction (window around R-peak: 90 samples before, 90 after for 180-sample beats at 360Hz)
    left_pad = 90
    right_pad = 90

    # Collect all beats and labels
    all_beats = []
    all_labels = []

    print(f"[{datetime.now()}] Processing records...")

    for rec in records:
        # Load signal and annotations
        sig, fields = wfdb.rdsamp(f'{db_path}/{rec}')
        ann = wfdb.rdann(f'{db_path}/{rec}', 'atr')

        # Use the first lead (usually MLII)
        signal = sig[:, 0]

        for i in range(len(ann.symbol)):
            sym = ann.symbol[i]
            if sym in class_map:
                label = class_map[sym]
                pos = ann.sample[i]
                # Extract beat segment if within bounds
                if pos - left_pad >= 0 and pos + right_pad < len(signal):
                    beat = signal[pos - left_pad: pos + right_pad]
                    resampled_beat = resample_to_dt(beat)
                    all_beats.append(resampled_beat)
                    all_labels.append(label)

    # Convert to numpy arrays
    all_beats = np.array(all_beats)
    all_labels = np.array(all_labels)

    print(f"[{datetime.now()}] Finished processing records.")
    
    # Stratified split into train/test (80/20)
    X_train, X_test, y_train, y_test = train_test_split(
        all_beats, all_labels, test_size=0.2, stratify=all_labels, random_state=42
    )

    # Create PyTorch datasets
    train_dataset = TensorDataset(
        torch.tensor(X_train, dtype=torch.float32),
        torch.tensor(y_train, dtype=torch.long)
    )
    test_dataset = TensorDataset(
        torch.tensor(X_test, dtype=torch.float32),
        torch.tensor(y_test, dtype=torch.long)
    )

    # Create DataLoaders with shuffle (shuffle=True for train, optional for test but added as per request)
    batch_size = 32  # Adjustable
    train_loader = DataLoader(train_dataset, batch_size=batch_size, shuffle=True)
    test_loader = DataLoader(test_dataset, batch_size=batch_size, shuffle=True)

    # Example usage: print shapes
    print(f"[{datetime.now()}] Train beats: {X_train.shape}, Test beats: {X_test.shape}")
    print(f"[{datetime.now()}] Train labels distribution: {np.bincount(y_train)}")
    print(f"[{datetime.now()}] Test labels distribution: {np.bincount(y_test)}")

    return train_loader, test_loader

def resample_to_dt(data, fs=360.0, dt=0.0001):
        t_orig = np.arange(len(data)) / fs
        interp = interp1d(t_orig, data, kind='linear', fill_value='extrapolate')
        duration = t_orig[-1]
        t_new = np.arange(0, duration, dt)
        return np.array(interp(t_new))

And I tried now set_device('cpp_standalone', directory='brian2_output', build_on_run=False) and

for batch_idx, (inputs, targets) in enumerate(train_loader):
    for idx, input in enumerate(inputs):

        print(f"[{datetime.now()}] Processing batch {batch_idx}, sample {idx}")
        ecg = input.numpy() if isinstance(inputs, torch.Tensor) else input
        stim = TimedArray(ecg.copy(order='C'), dt=dt)
        device.run(run_args={stimulus : stim})

        train_network(targets[idx].item())
        # print_nn()

    break

But got an error on device.run

Exception has occurred: TypeError

stat: path should be string, bytes, os.PathLike or integer, not NoneType

  File "/Users/dmyloserdov/SNN/snn.py", line 217, in <module>
    device.run(run_args={stimulus : stim})
TypeError: stat: path should be string, bytes, os.PathLike or integer, not NoneType
4

Hi @dmytr0 . Thanks for the details, I think I understand the issue now. In your first code block you posted, you run:

        train_network(targets[idx].item())
        network.restore("init")

Which means that you adapt the weights of your network, but directly afterwards you reset its state including the weights. Instead, you’d have to rewrite your train_network function so tha it does not change the weights in the network, but returns them. Your code would then set these new weights after resetting to the initial state:

        new_weights  = train_network(targets[idx].item())
        network.restore("init")
        network['S_lsm_out'].w[:] =  new_weights

Regarding your second comment about running it with cuda_standalone: this is indeed annoying, but hopefully easy to fix. The reason why it uses so much RAM is most likely the compilation, where we run with as many processes as possible by default (calling make -j). You can change this by setting the preference, e.g. to use 8 instead of “as many as possible” parallel processes:

prefs.devices.cpp_standalone.extra_make_args_unix = ['-j8']

This might just be a typo, but if you use build_on_run=True, then the run will trigger a build – if you want to call device.build yourself, you need build_on_run=False.

5

I’ve just seen that you updated your post

This is correct, but means that you first have to call device.build (to build and compile the simulation), and then device.run (to run the compiled binary). If you only call device.run, it will fail because there is no code to run.

6

I’ve tried fix for cuda_standalone and got errors on device.build(directory='brian2_output', compile=True, run=False)

code_objects/S_lsm_out_synapses_create_array_codeobject_104.cu(302): error: identifier "S_lsm_out_multiple_pre_post" is undefined
          S_lsm_out_multiple_pre_post = true;
          ^

code_objects/S_lsm_out_synapses_create_array_codeobject_101.cu(302): error: identifier "S_lsm_out_multiple_pre_post" is undefined
          S_lsm_out_multiple_pre_post = true;
          ^

code_objects/S_lsm_out_synapses_create_array_codeobject_1.cu(302): error: identifier "S_lsm_out_multiple_pre_post" is undefined
          S_lsm_out_multiple_pre_post = true;
          ^

code_objects/S_lsm_out_synapses_create_array_codeobject_100.cu(302): error: identifier "S_lsm_out_multiple_pre_post" is undefined
          S_lsm_out_multiple_pre_post = true;
          ^

code_objects/S_lsm_out_synapses_create_array_codeobject_10.cu(302): error: identifier "S_lsm_out_multiple_pre_post" is undefined
          S_lsm_out_multiple_pre_post = true;
          ^

code_objects/S_lsm_out_synapses_create_array_codeobject.cu(302): error: identifier "S_lsm_out_multiple_pre_post" is undefined
          S_lsm_out_multiple_pre_post = true;
          ^

code_objects/S_lsm_out_synapses_create_array_codeobject_103.cu(302): error: identifier "S_lsm_out_multiple_pre_post" is undefined
          S_lsm_out_multiple_pre_post = true;
          ^

1 error detected in the compilation of "code_objects/S_lsm_out_synapses_create_array_codeobject_104.cu".
code_objects/S_lsm_out_synapses_create_array_codeobject_102.cu(302): error: identifier "S_lsm_out_multiple_pre_post" is undefined
          S_lsm_out_multiple_pre_post = true;

For cpp_standalone got such error on updating stimulus

Exception has occurred: TypeError
Incorrect size for variable '_timedarray.values'. Shape (0,) ≠ (4973,).
  File "/Users/dmyloserdov/SNN/snn.py", line 217, in <module>
    device.run(run_args={stimulus : stim})
TypeError: Incorrect size for variable '_timedarray.values'. Shape (0,) ≠ (4973,).


for batch_idx, (inputs, targets) in enumerate(train_loader):
    for idx, input in enumerate(inputs):

        print(f"[{datetime.now()}] Processing batch {batch_idx}, sample {idx}")
        ecg = input.numpy() if isinstance(inputs, torch.Tensor) else input
        stim = ecg.copy(order='C')
        device.run(run_args={stimulus : stim})

        train_network(targets[idx].item())
        # print_nn()

    break

I also want to use this approach with the device. If I update the weights, should I recompile the network to run with the updated weights?

7

For cuda_standalone mode, this certainly looks like a bug, but I’d have to look into it more in detail.

Regarding the error for C++ standalone mode with TimedArray: you can provide a new TimedArray for each run with run_args, but the size of the TimedArray needs to be the same each time. To set the size correctly, the TimedArray that you define in your code before running device.build needs to have that same size, even if its content gets overwritten and is never used. So you’d need to replace the empty array in your first TimedArray definition by a list/array of the correct size.

You can provide the weights via the run_args as well, there is no need for recompilation. Please try something along the lines of:

new_weights = initial_weights
for batch...:
    for input...:
        device.run(run_args={stimulus : stim, network['S_lsm_out'].w: new_weights})

        new_weights = train_network(targets[idx].item())
8

hi @mstimberg

Thank you for suggestions, cpp_standalone version is working now. The issue was that I initialised an empty stimulus array and tried to fill it with a non-empty array.

But I have another question, can I somehow pass a different run time for device.run and stimulus with different sizes? One of the features that I want to use is to process input timed array of different sizes. Now I have error

Incorrect size for variable '_timedarray.values'. Shape (704,) ≠ (893,).
  File "/Users/dmyloserdov/SNN/snn.py", line 242, in <module>
    device.run(run_args={stimulus : stim, network['S_lsm_out'].w: new_weights})
TypeError: Incorrect size for variable '_timedarray.values'. Shape (704,) ≠ (893,).

code

# Previous processing and NN registration

network = Network([input_group, lsm_group, output_group, S_input_lsm, S_lsm_lsm, S_lsm_out, input_group_ST, input_group_SP, lsm_group_ST, lsm_group_SP, output_group_ST, output_group_SP])
network.run(700*ms)

device.build(directory='brian2_output', compile=True, run=False)

new_weights = initial_weights
for batch_idx, (inputs, targets) in enumerate(train_loader):
    for idx, input in enumerate(inputs):

        print(f"[{datetime.now()}] Processing batch {batch_idx}, sample {idx}")
        ecg = input.numpy() if isinstance(inputs, torch.Tensor) else input
        stim = ecg.copy(order='C')
        
        device.run(run_args={stimulus : stim, network['S_lsm_out'].w: new_weights})            
        new_weights = train_network(targets[idx].item())

        print_nn()
        break

    break
9

I figured out how to pass different sizes of stimulus. I tried this, and it works

new_weights = initial_weights
for batch_idx, (inputs, targets) in enumerate(train_loader):
    for idx, input in enumerate(inputs):

        print(f"[{datetime.now()}] Processing batch {batch_idx}, sample {idx}")
        ecg = input.numpy() if isinstance(inputs, torch.Tensor) else input
        stimulus = TimedArray(ecg.copy(order='C'), dt=dt)
        network.run(len(ecg)*2*ms)
        device.run(run_args={network['S_lsm_out'].w: new_weights})            
        new_weights = train_network(targets[idx].item())

        print_nn()

    break
10

Hi @mstimberg

Thank you for the support. For now, I will run cpp_standalone or runtime versions. Please let me know when you have time to look into cuda_standalone bug.

11

Unfortunately, this is currently a limitation of C++ standalone mode: you cannot change the run time of a simulation, since it is hardcoded into the generated code (changing the runtime therefore always requires a re-compilation). The only workaround is to choose the length of your longest stimulus as the runtime, and pad the stimulus array with zeros for the shorter stimuli (and then of course ignore the results after the actual stimulus). If the stimuli are not very different in size, it shouldn’t slow down things that much. There is a way to stop simulations early, but it is a bit clumsy to write and needs custom C++ code.

12

I don’t think this does what you want to do: everything you do after device.build will not affect the generated code. This means that in your case, the definition of the stimulus and the length of the network.run call will be ignored, since it will use the previously compiled code.