Improving the performance of interpolating between two arrays

Question

The Problem

I am currently writing a script that converts images into numerical array representation and then calculates "in-between" images based on linear interpolation between the start and end array.

My code does exactly what I want but goes over many nested loops which strikes me as something that will lead to very high computation times for many interpolation steps or big images.

The Code

The code is in python

import numpy as np

# Helper function that calculates the interpolation between two points
def interpolate_points(p1, p2, n_steps=3):
    # interpolate ratios between the points
    ratios = np.linspace(0, 1, num=n_steps)
    # linear interpolate vectors
    vectors = list()
    for ratio in ratios:
        v = (1.0 - ratio) * p1 + ratio * p2
        vectors.append(v)
    return np.asarray(vectors)

# final function that interpolates arrays
def interpolate_arrays(start_array,end_array,n_steps=10):
    n = 0
    array_interpolation = []
    while n < n_steps:
        i = 0
        x = []
        while i < len(start_array):
            e = interpolate_points(start_array[i],end_array[i],n_steps)[n]
            x.append(e)
            i += 1
        array_interpolation += [x]
        n += 1
    return array_interpolation

This results in:

#Test
X1 = [1,1]
X2 = [3,3]

interpolate_arrays(X1,X2,n_steps=3)
#[[1.0, 1.0], [2.0, 2.0], [3.0, 3.0]]

Answer 1

There are some easy wins here. Your interpolate_points doesn't need a loop:

def interpolate_points(p1, p2, n_steps=3):
    """Helper function that calculates the interpolation between two points"""
    # interpolate ratios between the points
    ratios = np.linspace(0, 1, num=n_steps)
    # linear interpolate vectors
    vectors = (1.0 - ratios) * p1 + ratios * p2
    return vectors

Also, even without further vectorization, you should be making use of range in your main function:

def interpolate_arrays(start_array, end_array, n_steps=10):
    """final function that interpolates arrays"""
    array_interpolation = []
    for n in range(n_steps):
        x = []
        for i in range(len(start_array)):
            e = interpolate_points(start_array[i], end_array[i], n_steps)[n]
            x.append(e)
        array_interpolation += [x]
    return array_interpolation

However, all of that can be replaced with a call to interp1d:

import numpy as np
from scipy.interpolate import interp1d


def interpolate_arrays(bounds, n_steps=10):
    """final function that interpolates arrays"""
    bounds = np.array(bounds)

    fun = interp1d(
        x=[0, 1],
        y=bounds.T,
    )
    y = fun(np.linspace(0, 1, n_steps))

    return y


def test():
    X1 = [1.5, 1]
    X2 = [5.5, 3]

    y = interpolate_arrays([X1, X2], n_steps=3)
    assert y.T.tolist() == [[1.5, 1.0], [3.5, 2.0], [5.5, 3.0]]

Even easier:

def interpolate_arrays(X1, X2, n_steps=10):
    """final function that interpolates arrays"""
    return np.linspace(X1, X2, n_steps)


def test():
    X1 = [1.5, 1]
    X2 = [5.5, 3]

    y = interpolate_arrays(X1, X2, n_steps=3)
    assert y.tolist() == [[1.5, 1.0], [3.5, 2.0], [5.5, 3.0]]

Notes:

• If you use interp1d, it would be better if your inputs and outputs are both two-dimensional np.ndarray; in their current form they need a transposition
• Write some unit tests such as the one shown, although it would be a better idea to call isclose since this is floating-point math
• If you want, it is trivial to make this extrapolate as well as interpolate

Basically: if there is a math thing in your head, before even thinking about what it would take to implement it yourself, do a search through scipy/numpy to see if it has already been done for you.