Hi @nvar . This topic that is much more complex than one might think when one reads papers that “add a Gaussian noise current” to the model It has been on my list for a while to write a longer document on this topic, but well…

Here are a few remarks:

There are different reasons to model “noise” in a model, and therefore different mathematical descriptions. Two commonly considered sources of noise is synaptic bombardment, where one does not want to model the individual spikes: see @rth 's comment above regarding `PoissonInput`

. If one assumes a very simple synapse model, low individual weights and the limit of “infinitely many incoming spikes” then this leads to an Ornstein-Uhlenbeck process with the time constant \tau_m added to the membrane potential, which is what we use in many examples. There are similar approximations when considering e.g. ion channel noise.

Now, you want to add “noise as a current with a mean of 1nA and a sigma of 10nA”. Taken literally, this is easy: you can include `I_noise`

in the `dv/dt`

equation, as you did above, and then define

```
I_noise = mu_noise + sigma_noise*randn() : amp (constant over dt)
```

Something to this effect is actually done in a number of papers. The problem with this formulation is that this describes an implementation in a step-based simulation with a certain time step, and not an “ideal” mathematical model. You could think of this as similar to describing the membrane potential not by a differential equation, but by something like `v_new = v_old + v_update`

. Or to change the perspective a bit: in electrophysiology, it wouldn’t be enough to say “we inject a \cal{N}(\mu, \sigma) white noise current”, but you’d either say something like “a \cal{N}(\mu, \sigma) white noise current filtered in this way”, or “we inject current, with its strength drawn every 1ms from \cal{N}(\mu, \sigma)”.

Here’s a simple simulation showing that this matters: it uses the equation from above to choose a new current every time step:

## Code example 1

```
from brian2 import *
mu_noise = .1*nA
sigma_noise = 1*nA
E_L = -70*mV
g_L = 10*nS
C_m = 100*pF
def run_sim(dt):
defaultclock.dt = dt
group = NeuronGroup(1,
'''dv/dt = (g_L*(E_L - v) + I_noise)/C_m : volt
I_noise = mu_noise + sigma_noise*randn() : amp (constant over dt)''',
method='exact')
group.v = E_L
mon = StateMonitor(group, ['v', 'I_noise'], record=0)
run(200*ms)
return mon.t[:], mon.v[0], mon.I_noise[0]
t1, v1, I_noise1 = run_sim(0.5*ms)
t2, v2, I_noise2 = run_sim(0.01*ms)
fig, axs = plt.subplots(4, 2, sharey='row',
gridspec_kw={'width_ratios': [3, 1]},
constrained_layout=True)
axs[0, 0].plot(t1/ms, I_noise1/nA)
axs[0, 0].set_ylabel('I_noise (nA)')
axs[0, 1].hist(I_noise1/nA, bins=np.arange(-5, 5, 0.5), orientation='horizontal')
axs[1, 0].plot(t1/ms, v1/mV)
axs[1, 0].set_ylim(-80, -20)
axs[1, 0].set_ylabel('v (mV)')
axs[1, 1].axis('off')
axs[2, 0].plot(t2/ms, I_noise2/nA)
axs[2, 0].set_ylabel('I_noise (nA)')
axs[2, 1].hist(I_noise2/nA, bins=np.arange(-5, 5, 0.5), orientation='horizontal')
axs[3, 0].plot(t2/ms, v2/mV)
axs[3, 0].set_ylim(-80, -20)
axs[3, 0].set_ylabel('v (mV)')
axs[3, 1].axis('off')
axs[3, 0].set_xlabel('t (ms)')
plt.show()
```

For demonstration purposes, I run this with `dt=0.5*ms` (top) and `dt=0.01*ms` (bottom):

As you can see in the histograms on the right, in both cases we have currents distributed according to a Gaussian with the same mean and variance. But as you can clearly see, the membrane potential reacts quite differently in the two cases. An alternative would be to use a `run_regularly`

operation to only update the `I_noise`

variable at regular intervals (e.g. once very ms), independent of the time step.

To abstract from the implementation and use a mathematical description, you can use differential equations and the variable \xi to describe the input current as an OU-process that reverts back to the mean:

```
group = NeuronGroup(1,
'''dv/dt = (g_L*(E_L - v) + I_noise)/C_m : volt
dI_noise/dt = (mu_noise - I_noise)/tau_noise + sigma_noise*sqrt(2/tau_noise)*xi : amp''',
method='euler')
```

Note that `tau_noise`

is a parameter of the noise (how fast the noise changes) and not necessarily `Cm/gL`

– using the membrane time constant for the noise is something that comes out of certain approximations, but it depends on what kind of noise you want to model specifically. This gives you noise that has a comparable effect independent of the time step, e.g. with the same time steps as above and `tau_noise=1*ms`

it gives

The implementation in the Bono & Clopath paper is very similar – it might even be equivalent, I did not look into this in detail.

As you noted in the other thread, using `xi`

does not work with certain algorithms like `exponential_euler`

, so in this case unfortunately you have to again leave the beautiful realm of mathematics and come back into the ugly reality of simulation. Replacing the above by:

```
'''....
dI_noise/dt = (mu_noise - I_noise)/tau_noise + noise : amp
noise = sigma_noise*sqrt(2/tau_noise)*randn()/sqrt(dt) : amp/second (constant over dt)
'''
```

will give you the Euler approximation to the noise term (so in the above code where it uses `euler`

already, this will not change anything), which is then compatible with `exponential_euler`

and all other algorithms that do not support stochastic differential equations.

I hope that made things a bit clearer?