如何实现线性插值？

Posted 2023-02-16

【中文标题】如何实现线性插值？

【英文标题】：How to implement linear interpolation?

【发布时间】：2011-11-12 17:01:03

【问题描述】：

假设我得到如下数据：

x = [1, 2.5, 3.4, 5.8, 6]
y = [2, 4, 5.8, 4.3, 4]

我想设计一个函数，使用 Python 在 1 和 2.5、2.5 到 3.4 等之间进行线性插值。

我尝试过查看this Python tutorial，但我仍然无法理解它。

【问题讨论】：

这……不容易。你试过什么？

-1 太笼统了。你不懂怎么编程，或者怎么用python做算法？？

作为一个新学习者，可以这么说，我已经把自己投入了深渊。我正在考虑在算法中使用“for”或“if”语句。所以在 x 的众多范围之间。

en.wikipedia.org/wiki/Linear_least_squares_(mathematics)

IMO，作为一个新手并不是不尝试的好借口——事实上恰恰相反（即一个很好的理由这样做）。

【参考方案1】：

import scipy.interpolate
y_interp = scipy.interpolate.interp1d(x, y)
print y_interp(5.0)

scipy.interpolate.interp1d 进行线性插值，可以自定义处理错误情况。

【讨论】：

问题问的是如何实现这个功能，而不是什么库提供的。

【参考方案2】：

据我了解你的问题，你想写一些函数y = interpolate(x_values, y_values, x)，它会给你y的值x？然后基本思想遵循以下步骤：

在x_values 中查找值的索引，这些索引定义了一个包含x 的区间。例如，对于带有示例列表的x=3，包含间隔为[x1,x2]=[2.5,3.4]，索引为i1=1、i2=2 通过(y_values[i2]-y_values[i1])/(x_values[i2]-x_values[i1])（即dy/dx）计算此区间的斜率。 x 的值现在是x1 的值加上斜率乘以与x1 的距离。

您还需要确定如果x 在x_values 的区间之外会发生什么，或者这是一个错误，或者您可以“向后”插值，假设斜率与第一个/最后一个区间相同。

这有帮助吗，还是您需要更具体的建议？

【参考方案3】：

我想出了一个相当优雅的解决方案（恕我直言），所以我忍不住发布它：

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])

我映射到float，因此如果x1、x2、y1 和y2 都是某个迭代的整数，则整数除法（python

在__getitem__ 中，我利用 self.x_list 以升序排序这一事实，使用bisect_left 来（非常）快速找到小于x 在self.x_list 中的最大元素的索引.

像这样使用类：

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[x] for x < 1 将像从 (2.5, 4) 到 (1, 2) 的线已扩展到负无穷大一样工作，而 i[x] for x == 1 或 x > 6 将提高一个IndexError。最好在所有情况下都提出 IndexError ，但这留给读者作为练习。 :)

【讨论】：

我发现使用__call__ 而不是__getitem__ 通常更可取，它通常是一个插值函数。

不适用于 Python3，因为不再支持索引地图

【参考方案4】：

以 Lauritz 的回答为基础，这是一个具有以下更改的版本

更新到 python3（地图给我带来了问题，没有必要） 修正了边缘值的行为 当 x 超出范围时引发异常 使用__call__ 而不是__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])

示例用法：

>>> interp = Interpolate([1, 2.5, 3.4, 5.8, 6], [2, 4, 5.8, 4.3, 4])
>>> interp(4)
5.425

【参考方案5】：

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 函数，每个搜索的人都想要，在纯 Python 中，甚至有一个文档字符串——一个完美的答案。我必须滚动才能看到它。

【参考方案6】：

您可以返回y_list 的范围，而不是从末端推断。大多数情况下，您的应用程序表现良好，Interpolate[x] 将位于x_list 中。从末端外推的（大概）线性影响可能会误导您相信您的数据表现良好。

返回非线性结果（以x_list 和y_list 的内容为界）您的程序行为可能会提醒您注意与x_list 之外的值存在很大问题。 （当给定非线性输入时，线性行为会变得很糟糕！）

在x_list 之外返回y_list 的Interpolate[x] 的范围也意味着您知道输出值的范围。如果您根据x 多远小于x_list[0] 或x 远多于x_list[-1] 进行推断，则您的返回结果可能超出您预期的值范围。

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__ 而不是__getitem__ 通常更可取，它通常是一个插值函数。

【参考方案7】：

您的解决方案在 Python 2.7 中不起作用。检查 x 元素的顺序时出错。我必须更改代码才能让它工作：

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])

【讨论】：

我得到一个错误.... TypeError 'Interpolate' object is not callable ??解决办法是什么？