Interpolation Methods in Python: Lagrange, Newton, and Piecewise Linear
Introduction
Interpolation methods are essential in data analysis, where we need to estimate function values at specific points based on known data points. In this article, we will discuss three commonly used interpolation methods: Lagrange, Newton, and piecewise linear interpolation. We will implement these methods in Python and provide examples of their usage.
Lagrange Interpolation
Lagrange interpolation is a method of constructing a polynomial function that passes through a set of given points. The basic idea is to find the n-th polynomial interpolation function pn(x) and rewrite it into another representation, then determine the undetermined function of interpolation conditions to obtain the interpolation polynomial.
Given n + 1 mutually different points (xi, fi), we can consider a configuration over these n + 1 points, where n does not exceed the number of polynomial. To estimate any point ξ, where ξ ≠ xi, i = 0, 1, 2, …, n, we can use the values Pn(ξ) as an exact value f(ξ) approximation.
The Lagrange interpolation polynomial is defined as:
Pn(ξ) = ∑[f(i) * L(i, ξ)]
where L(i, ξ) is the Lagrange basis polynomial.
Code Implementation
from matplotlib import pyplot as plt
def Lg(data, testdata):
"""
Lagrange interpolation function
"""
data_x = [data[i][0] for i in range(len(data))]
data_y = [data[i][1] for i in range(len(data))]
if testdata in data_x:
return data_y[data_x.index(testdata)]
predict = 0
for i in range(len(data_x)):
af = 1
for j in range(len(data_x)):
if j != i:
af *= (1.0 * (testdata - data_x[j]) / (data_x[i] - data_x[j]))
predict += data_y[i] * af
return predict
def plot(data, nums):
"""
Plot the interpolated function
"""
data_x = [data[i][0] for i in range(len(data))]
data_y = [data[i][1] for i in range(len(data))]
Area = [min(data_x), max(data_x)]
X = [Area[0] + 1.0 * i * (Area[1] - Area[0]) / nums for i in range(nums)]
X[len(X) - 1] = Area[1]
Y = [Lg(data, x) for x in X]
plt.plot(X, Y, label='result')
for i in range(len(data_x)):
plt.plot(data_x[i], data_y[i], 'ro', label="point")
plt.savefig('Lg.jpg')
plt.show()
data = [[0, 0], [1, 2], [2, 3], [3, 8], [4, 2]]
print(Lg(data, 1.5))
plot(data, 100)
Newton Interpolation
Newton interpolation is a method of constructing a polynomial function that passes through a set of given points. The basic idea is to be solved n-th interpolation polynomial Pn(x) and rewrite it into another representation, then use the interpolation conditions to determine Pn(x) is determined coefficients in order interpolation function a desired.
Newton’s difference is a concept that makes it easy to calculate the difference between the node increases. We define:
Δf(x0) = f(x0) - f(x0)
Δf(x1) = f(x1) - f(x0)
…
Δf(xk) = f(xk) - f(xk-1)
…
The k-order difference quotient is:
Δk f(xk) = Δk-1 f(xk) / Δk-1 xk
…
We can express the k-order difference quotient as a function of the value:
Δk f(xk) = ∑[(-1)^(k-i) * (k choose i) * f(xk-i)]
…
Code Implementation
from matplotlib import pyplot as plt
def calF(data):
"""
Calculate the difference quotient
"""
data_x = [data[i][0] for i in range(len(data))]
data_y = [data[i][1] for i in range(len(data))]
F = [1 for i in range(len(data))]
FM = []
for i in range(len(data)):
FME = []
if i == 0:
FME = data_y
else:
for j in range(len(FM[len(FM) - 1]) - 1):
delta = data_x[i + j] - data_x[j]
value = 1.0 * (FM[len(FM) - 1][j + 1] - FM[len(FM) - 1][j]) / delta
FME.append(value)
FM.append(FME)
F = [fme[0] for fme in FM]
return F
def NT(data, testdata, F):
"""
Newton interpolation function
"""
data_x = [data[i][0] for i in range(len(data))]
data_y = [data[i][1] for i in range(len(data))]
if testdata in data_x:
return data_y[data_x.index(testdata)]
else:
predict = 0
for i in range(len(data_x)):
Eq = 1
if i == 0:
for j in range(i):
Eq *= (testdata - data_x[j])
predict += (F[i] * Eq)
return predict
def plot(data, nums):
"""
Plot the interpolated function
"""
data_x = [data[i][0] for i in range(len(data))]
data_y = [data[i][1] for i in range(len(data))]
Area = [min(data_x), max(data_x)]
X = [Area[0] + 1.0 * i * (Area[1] - Area[0]) / nums for i in range(nums)]
X[len(X) - 1] = Area[1]
F = calF(data)
Y = [NT(data, x, F) for x in X]
plt.plot(X, Y, label='result')
for i in range(len(data_x)):
plt.plot(data_x[i], data_y[i], 'ro', label="point")
plt.savefig('Newton.jpg')
plt.show()
data = [[0, 0], [1, 2], [2, 3], [3, 8], [4, 2]]
plot(data, 100)
Piecewise Linear Interpolation
Piecewise linear interpolation is a method of constructing a linear function that passes through a set of given points. The basic idea is to find the nearest point to the interpolation point and use the linear function to estimate the value.
Code Implementation
from matplotlib import pyplot as plt
def DivideLine(data, testdata):
"""
Piecewise linear interpolation function
"""
data_x = [data[i][0] for i in range(len(data))]
data_y = [data[i][1] for i in range(len(data))]
if testdata in data_x:
return data_y[data_x.index(testdata)]
else:
index = 0
for j in range(len(data_x)):
if data_x[j] < testdata and data_x[j + 1] > testdata:
index = j
break
predict = 1.0 * (testdata - data_x[index]) * (data_y[index + 1] - data_y[index]) / (data_x[index + 1] - data_x[index]) + data_y[index]
return predict
def plot(data, nums):
"""
Plot the interpolated function
"""
data_x = [data[i][0] for i in range(len(data))]
data_y = [data[i][1] for i in range(len(data))]
Area = [min(data_x), max(data_x)]
X = [Area[0] + 1.0 * i * (Area[1] - Area[0]) / nums for i in range(nums)]
X[len(X) - 1] = Area[1]
Y = [DivideLine(data, x) for x in X]
plt.plot(X, Y, label='result')
for i in range(len(data_x)):
plt.plot(data_x[i], data_y[i], 'ro', label="point")
plt.savefig('DivLine.jpg')
plt.show()
data = [[0, 0], [1, 2], [2, 3], [3, 8], [4, 2]]
print(DivideLine(data, 1.5))
plot(data, 100)