Plot stochastic trajectories deviations from 'real' path using a colormesh in matplotlib (Python)

Question

Hi I created a program that will create deviations from a real trajectory, it is complicated and I do not have a simple example unfortunately.

It calculates a path with stochastic initial conditions from the real path and does this for x iterations, the goal is to show that the deviations become larger at greater times.

The real path and the deviations are showed below.

However I want to show that the deviations become greater the longer in time we are. Ofcourse I could just calculate the variance and plot mean+var and mean-var at each time step but I was wondering if I could plot something like this, using hist2d

You see that the blocks are not as smooth as a like and this is not that great to use.

Then I went and looked at python's kde and created the following.

This is also not preferable as I think it bins more points at the minima and maxima. Also it is 'too smeared out'. Especially in the beginning, all the points are the same so I want there just to be a straight line to really show that the deviations start later on.

I guess my question is; is what I want even possible and what package/command should I use. I haven't found what I am looking for on other questions. Or has anyone a suggestion to nicely show what I want in a any other way?

Answer 1

Here is an idea plotting multiple curves with transparency on top of each other:

import numpy as np
import matplotlib.pyplot as plt

x = np.linspace(0, 10, 200)
for _ in range(1000):
    plt.plot(x, np.sin(x * np.random.normal(1, 0.1)) * np.random.normal(1, 0.1), color='r', alpha=0.02)
plt.plot(x, np.sin(x), color='b')
plt.margins(x=0)
plt.show()

Another option creates a 2d histogram:

import numpy as np
import matplotlib.pyplot as plt

x = np.linspace(0, 10, 200)
all_curves = np.array([np.sin(x * np.random.normal(1, 0.1)) * np.random.normal(1, 0.1) for _ in range(100)])
plt.hist2d(x=np.tile(x, all_curves.shape[0]), y=all_curves.ravel(), bins=(100, 100), cmap='inferno')
plt.show()

Still another approach would use fill_between (as suggested by @bramb) between confidence intervals:

import numpy as np
import matplotlib.pyplot as plt

x = np.linspace(0, 10, 200)
all_curves = np.array([np.sin(x * np.random.normal(1, 0.1)) * np.random.normal(1, 0.1) for _ in range(1000)])

confidence_interval1 = 95
confidence_interval2 = 80
confidence_interval3 = 50
for ci in [confidence_interval1, confidence_interval2, confidence_interval3]:
    low = np.percentile(all_curves, 50 - ci / 2, axis=0)
    high = np.percentile(all_curves, 50 + ci / 2, axis=0)
    plt.fill_between(x, low, high, color='r', alpha=0.2)
plt.plot(x, np.sin(x), color='b')
plt.margins(x=0)
plt.show()