python中的3D样条插值

Question

我正在搜索等效的 Matlab 命令

Vq = interp3(X,Y,Z,V,Xq,Yq,Zq)

在 Python 中。在 Matlab 中，我可以使用方法“样条”插值，我在 python 中找不到 3D 数据。存在 scipy.interpolate.griddata，但没有 3D 数据样条选项。

我要插值的数据是一个 3D 矩阵 (51x51x51)，它有规律地分布在 3D 网格上。

scipy.interpolate.Rbf 可能是选项，但我不明白它的工作原理：

xi = yi = zi = np.linspace(1, 132651, 132651) interp = scipy.interpolate.Rbf(xi, yi, zi, data, function='cubic')

导致内存错误。

编辑： 我想要的一个最小的例子（没有插值）： Matlab代码

v=rand([51,51,51]);
isosurface (v, 0.3);

为简单起见，我在此示例中使用随机数据。我想制作等值面图（特别是费米曲面图）。由于某些结构非常小，因此需要 51x51x51 的高网格分辨率。

进一步说明：矩阵中的数据集是相互独立的，z（或第 3 个分量）不是 x 和 y 的函数。

Answer 1

可以按照您的描述使用scipy.interpolate.Rbf 完成 3 维以上的样条插值。出于绘图目的，您可以使用较小的分辨率（1000 点是一个很好的经验法则），当您想要评估您的样条时，您可以毫无问题地在远大于 132000 的点上进行插值（参见下面的示例）。

您能否为您在 matlab 中尝试做的事情添加 Minimal, Complete, and Verifiable example？这将解释为什么需要创建分辨率为 132000 点的网格空间。另外，请注意，存在维度灾难。 Matlab 使用cubic spline or a piecewise polynomial，由于过度拟合可能很危险。我建议您使用更理智的方法来训练 51 个数据点并应用于 132000 多个数据点。 This 是多项式曲线拟合和模型选择的一个很好的例子。

示例：

生成数据：

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import axes3d

%matplotlib inline
import random
# set seed to reproducible
random.seed(1)
data_size = 51
max_value_range = 132651
x = np.array([random.random()*max_value_range for p in range(0,data_size)])
y = np.array([random.random()*max_value_range for p in range(0,data_size)])
z = 2*x*x*x + np.sqrt(y)*y + random.random()
fig = plt.figure(figsize=(10,6))
ax = axes3d.Axes3D(fig)
ax.scatter3D(x,y,z, c='r')

拟合样条和插值

x_grid = np.linspace(0, 132651, 1000*len(x))
y_grid = np.linspace(0, 132651, 1000*len(y))
B1, B2 = np.meshgrid(x_grid, y_grid, indexing='xy')
Z = np.zeros((x.size, z.size))

import scipy as sp
import scipy.interpolate
spline = sp.interpolate.Rbf(x,y,z,function='thin_plate',smooth=5, episilon=5)

Z = spline(B1,B2)
fig = plt.figure(figsize=(10,6))
ax = axes3d.Axes3D(fig)
ax.plot_wireframe(B1, B2, Z)
ax.plot_surface(B1, B2, Z,alpha=0.2)
ax.scatter3D(x,y,z, c='r')

在大数据上拟合样条

predict_data_size = 132000
x_predict = np.array([random.random()*max_value_range for p in range(0,predict_data_size)])
y_predict = np.array([random.random()*max_value_range for p in range(0,predict_data_size)])
z_predict = spline(x_predict, y_predict)
fig = plt.figure(figsize=(10,6))
ax = axes3d.Axes3D(fig)
ax.plot_wireframe(B1, B2, Z)
ax.plot_surface(B1, B2, Z,alpha=0.2)
ax.scatter3D(x_predict,y_predict,z_predict, c='r')