mattja / nsim

Large scale simulation of ODEs or SDEs, analyze time series.
GNU General Public License v3.0
24 stars 5 forks source link

plot time scale are wrong #3

Closed fccoelho closed 8 years ago

fccoelho commented 8 years ago

When you plot a simulation with timeseries.plot(), the time axis is not scaled correctly.

if you run a simulation like this:

sim = nsim.RepeatedSim(Model, T=100, repeat=2)
sim.timeseries.plot()

The time scale of the plot should go from 0 to 100, right?.

mattja commented 8 years ago

It should plot from time 0 to 100 only if there is data up to time 100. Please check whether your time series has only NaN (Not a Number) after a short time. Because that particular system reaches extinction after a finite time, when the noise process pushes I above 1 or below 0. After that time those equations will not give a real number to plot.

fccoelho commented 8 years ago

Ok, I Noticed that. I supposed your Runge-kutta impementation does not use adaptive step-sizes, as we have in scipy's ODEINT, right? I don't know if that is easy of even feasible. But since in theory the variance of the wiener process is proportinal to h (or dt), shinking h with the magnitude of Dx/dt might help the system to remain within the boundaries in which is is defined. It might even be interesting to be able to declare the boundaries of expected dynamics, and perhepos reject steps which fall beyond those boundaries...

I am no expert on these SDE solving algorithms. So I apologize in advance if I am proposing something absurd.

mattja commented 8 years ago

For that particular equation we discussed, at the boundary I=1 the noise intensity for the second equation does not go to zero. So hitting the boundary is a true behavior of the sample paths, in this case it's not due to numerical approximation.

I think you're right, it could be useful feature to automatically halt when it exits the boundary. I can think of some other situations (with repeated simulations) where it is useful to still return arrays of the fixed length, so you can easily plot them together. Probably can do both.

By the way, yes there do exist proper stochastic Runge-Kutta algorithms with adaptive step length, here's one: http://epubs.siam.org/doi/abs/10.1137/S1064827500376922 That's one that I haven't got around to implementing in the library yet. I'm no expert either. Thanks for reporting the two recent bugs.