【问题标题】:scipy.interp2d warning and different result than expectedscipy.interp2d 警告和与预期不同的结果
【发布时间】:2016-01-16 09:42:44
【问题描述】:

我正在尝试将 MATLAB 代码转换为等效的 python。 我有 3 个数组,我想计算 interp2d:

nuA = np.asarray([2.439,2.5,2.6,2.7,2.8,3.0,3.2,3.5,4.0,5.0,6.0,8.0,10,15,25])
nuB = np.asarray([0,0.1,0.2,0.3,0.5,0.7,1])
a, b = np.meshgrid(nuA, nuB)
betaTab  = np.transpose(np.asarray([[0.0,2.16,1.0,1.0,1.0,1.0,1.0],[0.0,1.592,3.39,1.0,1.0,1.0,1.0],[0.0,0.759,1.8,1.0,1.0,1.0,1.0],[0.0,0.482,1.048,1.694,1.0,1.0,1.0],[0.0,0.36,0.76,1.232,2.229,1.0,1.0],[0.0,0.253,0.518,0.823,1.575,1.0,1.0],[0.0,0.203,0.41,0.632,1.244,1.906,1.0],[0.0,0.165,0.332,0.499,0.943,1.56,1.0],[0.0,0.136,0.271,0.404,0.689,1.23,2.195],[0.0,0.109,0.216,0.323,0.539,0.827,1.917],[0.0,0.096,0.19,0.284,0.472,0.693,1.759],[0.0,0.082,0.163,0.243,0.412,0.601,1.596],[0.0,0.074,0.147,0.22,0.377,0.546,1.482],[0.0,0.064,0.128,0.191,0.33,0.478,1.362],[0.0,0.056,0.112,0.167,0.285,0.428,1.274]]))
ip = scipy.interpolate.interp2d(a,b,betaTab)

当我尝试运行它时,会显示以下警告:

/usr/local/lib/python2.7/dist-packages/scipy/interpolate/fitpack.py:981: RuntimeWarning: No more knots can be added because the additional knot would
coincide with an old one. Probable cause: s too small or too large
a weight to an inaccurate data point. (fp>s)
    kx,ky=1,1 nx,ny=4,14 m=105 fp=21.576347 s=0.000000
  warnings.warn(RuntimeWarning(_iermess2[ierm][0] + _mess))

我知道interp2dmatlab interp2 不同,在python 中RectBivariateSpline 函数是首选。但由于我的数据长度,我不能使用后一个函数。此外,ip(xi,yi) 的最终结果与 MATLAB 答案不同。

如何在没有警告的情况下计算 interp2d 并正确计算?

【问题讨论】:

  • MATLAB 代码做了什么样的插值?样条线,分段线性(三角形平面)?
  • 该错误实际上是说对于您的数据,线性样条插值是不好的。例如,它可能太平了。表面是什么样的?

标签: python scipy linear-interpolation


【解决方案1】:

您的输入数据似乎定义不明确。这是您的输入点的表面图:

这不是一个容易插值的问题。顺便说一句,我有recently ran into problems interp2d 甚至无法插入平滑的数据集。所以我建议改为查看scipy.interpolate.griddata

import numpy as np
import scipy.interpolate as interp
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

#define your data as you did in your question: a, b and betaTab

ip = interp.interp2d(a,b,betaTab)           # original interpolator

aplotv = np.linspace(a.min(),a.max(),100)   # to interpolate at
bplotv = np.linspace(b.min(),b.max(),100)   # to interpolate at
aplot,bplot = np.meshgrid(aplotv,bplotv)    # mesh to interpolate at

# actual values from interp2d:
betainterp2d = ip(aplotv,bplotv)

# actual values from griddata:
betagriddata = interp.griddata(np.array([a.ravel(),b.ravel()]).T,betaTab.ravel(),np.array([aplot.ravel(),bplot.ravel()]).T)
# ^ this probably could be written in a less messy way,
# I'll keep thinking about it

#plot results
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(aplot,bplot,betainterp2d,cmap='viridis',cstride=1,rstride=1)

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot_surface(aplot,bplot,betagriddata,cmap='viridis',cstride=1,rstride=1)

结果:(左:interp2d,右:griddata

结论:使用scipy.interpolate.griddata

【讨论】:

  • 如何从特定点获得价值? griddata 的结果是一维 (shape=(10000,)) 数组,而在 ip = interp2d(...) 中,我们可以使用 ip(xi,yi) 获得最终插值。 (xi , yi) 在一种情况下的对应值是(2.439, 0.097),结果必须大约是2.1
  • 已解决:使用print(griddata(np.array([a.ravel(),b.ravel()]).T,betaTab.ravel(),np.array([2.439, 0.097]).T)) 获取2.09
  • 没错。虽然来自interp2d 的插值器函数接受一维网格上的点(并输出从这些一维网格生成的二维网格上的插值点),但griddata 在输入上逐点进行(输入和输出点都应该是@ 987654346@s,其中每个 对应于每个点)。请注意,在您只有 1 个输入点的情况下,您不需要转置:您的单点应该是单行向量,np.array([2.43‌​9, 0.097]) 等价于(实际上是 transposing that 1d array won't do anything)。
猜你喜欢
  • 2016-04-07
  • 1970-01-01
  • 1970-01-01
  • 2017-06-12
  • 2019-12-27
  • 2015-08-13
  • 2015-02-17
  • 2015-09-14
  • 1970-01-01
相关资源
最近更新 更多