numpy: how interpolate between two arrays for various timesteps?

Question

I'm looking for a way to do a simple linear interpolation between two numpy arrays that represent a start and endpoint in time.

The two arrays have the same length:

fst = np.random.random_integers(5, size=(10.))
>>> array([4, 4, 1, 3, 1, 4, 3, 2, 5, 2])
snd = np.random.random_integers(5, size=(10.))
>>> array([1, 1, 3, 4, 1, 5, 5, 5, 4, 3])

Between my start and endpoint there are 3 timesteps. How can I interpolate between fst and snd? I want to be able, taking the first entry of fst and snd as an example, to retrieve the value of each timestep like

np.interp(1, [1,5], [4,1])    
np.interp(2, [1,5], [4,1])
...
# that is
np.interp([1,2,3,4,5], [1,5], [4,1])
>>> array([ 4.  ,  3.25,  2.5 ,  1.75,  1.  ])

But than not just for the first entry but over the whole array.

Obviously, this won't do it:

np.interp(1, [1,5], [fst,snd])

Well I know I get there in a loop, e.g.

[np.interp(2, [1,5], [item,snd[idx]]) for idx,item in enumerate(fst)]
>>> [3.25, 3.25, 1.5, 3.25, 1.0, 4.25, 3.5, 2.75, 4.75, 2.25]

but I believe when you are lopping over numpy arrays you are doing something fundamentally wrong.

Answer 1

The facilities in scipy.interpolate.interp1d allow this to be done quite easily if you form your samples into a 2D matrix. In your case, you can construct a 2xN array, and construct an interpolation function that operates down the columns:

from scipy.interpolate import interp1d
fst = np.array([4, 4, 1, 3, 1, 4, 3, 2, 5, 2])
snd = np.array([1, 1, 3, 4, 1, 5, 5, 5, 4, 3])
linfit = interp1d([1,5], np.vstack([fst, snd]), axis=0)

You can then generate an interpolated vector at any time of interest. For example linfit(2) produces:

array([ 3.25,  3.25,  1.5 ,  3.25,  1.  ,  4.25,  3.5 ,  2.75,  4.75,  2.25])

or you can invoke linfit() with a vector of time values, e.g. linfit([1,2,3]) gives:

array([[ 4.  ,  4.  ,  1.  ,  3.  ,  1.  ,  4.  ,  3.  ,  2.  ,  5.  ,  2.  ],
       [ 3.25,  3.25,  1.5 ,  3.25,  1.  ,  4.25,  3.5 ,  2.75,  4.75,           2.25],
       [ 2.5 ,  2.5 ,  2.  ,  3.5 ,  1.  ,  4.5 ,  4.  ,  3.5 ,  4.5 , 2.5 ]])

If you're only doing linear interpolation, you could also just do something like:

((5-t)/(5-1)) * fst + ((t-1)/(5-1)) * snd

to directly compute the interpolated vector at any time t.