Say I am given data as follows:
x = [1, 2.5, 3.4, 5.8, 6] y = [2, 4, 5.8, 4.3, 4] I want to design a function that will interpolate linearly between 1 and 2.5, 2.5 to 3.4, and so on using Python.
I have tried looking through this Python tutorial, but I am still unable to get my head around it.
import scipy.interpolate y_interp = scipy.interpolate.interp1d(x, y) print y_interp(5.0) scipy.interpolate.interp1d does linear interpolation by and can be customized to handle error conditions.
As I understand your question, you want to write some function y = interpolate(x_values, y_values, x), which will give you the y value at some x? The basic idea then follows these steps:
x_values which define an interval containing x. For instance, for x=3 with your example lists, the containing interval would be [x1,x2]=[2.5,3.4], and the indices would be i1=1, i2=2 (y_values[i2]-y_values[i1])/(x_values[i2]-x_values[i1]) (ie dy/dx).x is now the value at x1 plus the slope multiplied by the distance from x1.You will additionally need to decide what happens if x is outside the interval of x_values, either it's an error, or you could interpolate "backwards", assuming the slope is the same as the first/last interval.
Did this help, or did you need more specific advice?
def interpolate(x1: float, x2: float, y1: float, y2: float, x: float): """Perform linear interpolation for x between (x1,y1) and (x2,y2) """ return ((y2 - y1) * x + x2 * y1 - x1 * y2) / (x2 - x1) 
lerp function everyone who searches for wants, in pure Python, even has a docstring—a perfect answer. And I had to scroll to see it.I thought up a rather elegant solution (IMHO), so I can't resist posting it:
from bisect import bisect_left class Interpolate(object): def __init__(self, x_list, y_list): if any(y - x <= 0 for x, y in zip(x_list, x_list[1:])): raise ValueError("x_list must be in strictly ascending order!") x_list = self.x_list = map(float, x_list) y_list = self.y_list = map(float, y_list) intervals = zip(x_list, x_list[1:], y_list, y_list[1:]) self.slopes = [(y2 - y1)/(x2 - x1) for x1, x2, y1, y2 in intervals] def __getitem__(self, x): i = bisect_left(self.x_list, x) - 1 return self.y_list[i] + self.slopes[i] * (x - self.x_list[i]) I map to float so that integer division (python <= 2.7) won't kick in and ruin things if x1, x2, y1 and y2 are all integers for some iterval.
In __getitem__ I'm taking advantage of the fact that self.x_list is sorted in ascending order by using bisect_left to (very) quickly find the index of the largest element smaller than x in self.x_list.
Use the class like this:
i = Interpolate([1, 2.5, 3.4, 5.8, 6], [2, 4, 5.8, 4.3, 4]) # Get the interpolated value at x = 4: y = i[4] I've not dealt with the border conditions at all here, for simplicity. As it is, i[x] for x < 1 will work as if the line from (2.5, 4) to (1, 2) had been extended to minus infinity, while i[x] for x == 1 or x > 6 will raise an IndexError. Better would be to raise an IndexError in all cases, but this is left as an exercise for the reader. :)
__call__ instead of __getitem__ to be preferrable in general, its usually an interpolation function.
Building on Lauritz` answer, here's a version with the following changes
__call__ instead of __getitem__from bisect import bisect_right class Interpolate: def __init__(self, x_list, y_list): if any(y - x <= 0 for x, y in zip(x_list, x_list[1:])): raise ValueError("x_list must be in strictly ascending order!") self.x_list = x_list self.y_list = y_list intervals = zip(x_list, x_list[1:], y_list, y_list[1:]) self.slopes = [(y2 - y1) / (x2 - x1) for x1, x2, y1, y2 in intervals] def __call__(self, x): if not (self.x_list[0] <= x <= self.x_list[-1]): raise ValueError("x out of bounds!") if x == self.x_list[-1]: return self.y_list[-1] i = bisect_right(self.x_list, x) - 1 return self.y_list[i] + self.slopes[i] * (x - self.x_list[i]) Example usage:
>>> interp = Interpolate([1, 2.5, 3.4, 5.8, 6], [2, 4, 5.8, 4.3, 4]) >>> interp(4) 5.425 Instead of extrapolating off the ends, you could return the extents of the y_list. Most of the time your application is well behaved, and the Interpolate[x] will be in the x_list. The (presumably) linear affects of extrapolating off the ends may mislead you to believe that your data is well behaved.
Returning a non-linear result (bounded by the contents of x_list and y_list) your program's behavior may alert you to an issue for values greatly outside x_list. (Linear behavior goes bananas when given non-linear inputs!)
Returning the extents of the y_list for Interpolate[x] outside of x_list also means you know the range of your output value. If you extrapolate based on x much, much less than x_list[0] or x much, much greater than x_list[-1], your return result could be outside of the range of values you expected.
def __getitem__(self, x): if x <= self.x_list[0]: return self.y_list[0] elif x >= self.x_list[-1]: return self.y_list[-1] else: i = bisect_left(self.x_list, x) - 1 return self.y_list[i] + self.slopes[i] * (x - self.x_list[i]) 
__call__ instead of __getitem__ to be preferrable in general, its usually an interpolation function.Your solution did not work in Python 2.7. There was an error while checking for the order of the x elements. I had to change to code to this to get it to work:
from bisect import bisect_left class Interpolate(object): def __init__(self, x_list, y_list): if any([y - x <= 0 for x, y in zip(x_list, x_list[1:])]): raise ValueError("x_list must be in strictly ascending order!") x_list = self.x_list = map(float, x_list) y_list = self.y_list = map(float, y_list) intervals = zip(x_list, x_list[1:], y_list, y_list[1:]) self.slopes = [(y2 - y1)/(x2 - x1) for x1, x2, y1, y2 in intervals] def __getitem__(self, x): i = bisect_left(self.x_list, x) - 1 return self.y_list[i] + self.slopes[i] * (x - self.x_list[i]) 
