Details

Type: Bug

Status: Closed

Priority: Major

Resolution: Fixed

Affects Version/s: None

Fix Version/s: 1.1

Component/s: TVBSimulator

Epic Link:
Description
Due to the initial method in models. The amplitude of the diffusion process,
nsig, is rescaled by the length and the integration time step size. So, while
the scaling makes sense in order to respect the time scale at which the
diffusion process occurs, the speed affects the initial conditions, which is a
problem.
Discussions with Stuart:
In the current implementation initial conditions are set so the total amount
of diffusion is bounded to the statevariableranges.
Even without the scaling by length, the change in length will mean a different
time over which the basic diffusion process operates, so a change in
conduction speed will still produce different initial conditions – possibly
more noticeable due to the different amplitude achieved by the diffusion as a
result of the time over which it operates.
Solutions: Read more about setting the initial history for delayed systems. h
ttp://reference.wolfram.com/mathematica/tutorial/NDSolveDelayDifferentialEquations.html
0) Form the literature, very often it is seen that the initial history is a
constant vector with the initial values for the ODE system (dummy solution).
1a) Separate the noise source's nsig and the statevariablerange bounding.
1b) Set nsig using an actual integration so it's properly scaled by dt, and
including a simple forcing term set based on statevariablerange,
ie: dx/dt =a x where "a" is automatically set inversely proportional statevariable
range.
2) TVB's integrators may as well be used, with their noise, which probably
means we could drop the separate Model noise entirely.
+ Make sure to reverse the time after the integration, so that the changes in
speed, and thus length, only changes history furthest in the past.
3) To stick to the diffusion approach: the desired result should be achievable by
using a constant scaling factor for nsig, instead of the current history
length dependent one. For a single constant value to achieve the current
purpose of the scaling (ie, keeping it within sensible bounds), the constant
would need to be set by replacing the current history length by the longest
possible history length, that is, that caused by the slowest permitted
speed... This would produce less broadly diffused initial conditions than
currently for most speeds.
There was also the issue that for some models whose state variables
represent firing rates, the range had to be set very narrow, otherwise in the
initial history some negative values appeared, which is not strictly
correct for that type of models.
4) Using a constant scaling, as mentioned above, that considers the longest permitted
history for the rate of diffusion, combined with a truncated normal distribution for the
random variates, scipy has one:
scipy.stats.truncnorm
should make it possible to turn the "probably bounded" above into a "strictly bounded".
That is, something like setting set the truncation to:
((max_svr  min_svr)/2)/max_steps
should work, shouldn't it??? It should mean that, starting from the mid point, taking the
extremely unlikely worst case of every random variate returned being a maximum, ie at
the cutoff, would over the history length (max_steps) leave you at the boundary of the
statevariablerange.
Implementing this neatly would probably require modifying TVB's noise so
that there is an option to use truncnorm instead of normal as the form of the
random stream, and then making use of this for the Model's noise definition...
And the longest possible history is based on the steps required to represent
the maximum time delay, which depends on dt, the speed and the longest fiber
distance – on average 100 mm according to:
Schüz A, Braitenberg V. 2002. The human cortical white matter: quantitative
aspects of corticocortical longrange connectivity. In: Schüz A, Miller R,
editors. Cortical areas: unity and diversity. London: CRC Press. p 377–385.
Schuz and Braitenberg
However, in some datasets (eg http://heliosphan.org/pittsburghbraincomp.html)
the longest white matter fiber is said to be over 200mm always speaking of
corticocortical connections.
How the delays operate in TVB:

The delayed_state gives the state of the process at a time prior to the
current time. It takes the coupling variable, "cvar", given by the model to
know which of the history variables to report. Then, a delay index gives
the time delay (state further in the history array) at which to obtain the
value/state of the selected coupling variable.