Measure timescale separation in time series

[ftavella]

From our preliminary results, we see that networks with more "spiky" signals tend to be more period tunable. It would be nice to have a measure of this spikiness to see if it correlates with period tunability. Therefore, it would be great to implement a function that given all the time series for a particular network (and parameter set), calculates the ratio between the rising part of the signal and the falling part of the signal. This measure could be a proxy for a "biphasic" behavior in oscillations and might correlate with period tunability.

[ftavella]

A good way of quantifying the timescale separation would be through the following formula

time_scale_separation = abs((rising - falling) / (rising + falling))

where rising corresponds to the time it takes to go from the trough to the peak and falling the time it takes to go from the peak to the trough. This quantity can be averaged over nodes to get a single number per topology and parameter set.

[biphy]

I can work on this.

The formula gives ~1 for sawtooth waveforms, which I think are not necessarily related to relaxation oscillators. I think we can instead measure the ratio of reaction rates ($(dP_i/dt)/(dP_j/dt)$ or its absolute value) during one period and take the maximum (possibly after some regularization, or some other statistics) as a quantity directly related to the time scale separation. If this is large, it indicates there are two species whose populations are changing at drastically different rates, which is exactly what we mean by the time scale separation in relaxation oscillations.

[ftavella]

That's a great idea! Thanks for working on this.

[biphy]

I slightly changed the approach and just measured how long phase trajectories stay near the nullcline of each variable in one cycle. The following numbers are that time divided by the respective period, so values closer to 1 indicate larger time scale separation.

It seems this factor can distinguish spiky waveforms (03, 04) from others--in both cases dynamics takes place near the nullcline of the first and fastest variable (blue lines). I'm not sure what happened for the condition 02, where the waveforms look pretty spiky but the factor is ~0.15. But it's still larger than ~0.01 factors for the other sinusoidal cases.

# Goodwin + IP, random oscillatory parameters
# var1, var2, var3
01: 0.00847457627118645, 0.004237288135593225, 0.00847457627118645
02: 0.15420560747663514, 0.13084112149532673, 0.13551401869158894
03: 0.922885572139304, 0.7363184079601993, 0.2114427860696517
04: 0.8602150537634407, 0.4139784946236558, 0.14516129032258082
05: 0.028901734104046263, 0.011560693641618505, 0.011560693641618505

# Goodwin, random oscillatory parameters
06: 0.014925373134328382, 0.014925373134328382, 0.0
07: 0.008130081300813016, 0.012195121951219523, 0.012195121951219523
08: 0.012987012987013005, 0.006493506493506503, 0.006493506493506503
09: 0.013157894736842105, 0.0, 0.013157894736842105
10: 0.006756756756756767, 0.006756756756756767, 0.006756756756756767

[ftavella]

Love this idea! The numbers are looking good too. Let me know when you want me to review a PR!