scipy插值速度:是否有比RectBivariateSpline更快的方法?

5

我正在开发一个项目,该项目在很大程度上依赖于2D插值。经过一些初步搜索,我发现很多文章都说地图坐标或可能的RegularGridInterpolator可以提供最快的速度,只要我能接受结构化网格。我编写了以下代码试用scipy中的所有多变量插值函数,结果表明RectBivariateSpline是最佳选择。

from __future__ import division
import numpy as np
import time
import math
from scipy.interpolate import RectBivariateSpline, RegularGridInterpolator, interp2d, CloughTocher2DInterpolator, Rbf, griddata, interpn
from scipy.ndimage.interpolation import map_coordinates
import random

grid_size = 100
xrange = np.linspace(0, 1, grid_size)
yrange = np.linspace(0, 1, grid_size)
x2d, y2D = np.meshgrid(xrange, yrange, indexing='ij')

z2D = x2d ** 2.0 + y2D

i2d = interp2d(xrange, yrange, z2D, 'linear')
rbs = RectBivariateSpline(xrange, yrange, z2D)
RGI = RegularGridInterpolator((xrange, yrange), z2D)
CT2d = CloughTocher2DInterpolator(np.array(list(zip(x2d.flatten(), y2D.flatten()))), z2D.flatten())
rbfi = Rbf(x2d, y2D, z2D, function='linear')
grd_inputs1 = np.array(list(zip(x2d.flatten(), y2D.flatten())))
grd_inputs2 = z2D.flatten()

print 'Exact: ',0.5**2.0 + 0.5
print 'interp2d',i2d(0.5, 0.5)[0]
print 'RectBivariateSpline',rbs(0.5, 0.5)[0][0]
print 'RegularGridInterpolator',RGI([0.5, 0.5])[0]

coords = np.asarray(([0.5], [0.5]))
coords = [(c - lo) * (n - 1) / (hi - lo) for (lo, hi), c, n in zip([[0,1],[0,1]], coords, np.shape(z2D.flatten()))]
print 'Map Coordinates:', map_coordinates(z2D.flatten(), coords, order=1)[0]

print 'Rbf', rbfi([0.5, 0.5])[0]
print 'CloughTocher2DInterpolator', CT2d([0.5, 0.5])[0]
print 'griddata', griddata(grd_inputs1, grd_inputs2, [0.5, 0.5], 'linear')[0]
print 'interpn', interpn([xrange, yrange], z2D, [0.5, 0.5])[0]

print ''

samples = 25000

start = time.time()
for n in range(samples):
    i2d(0.5, 0.5)
end = time.time()
print 'interp2d: %0.4f [us]' % (((end-start)/samples)*1e6)

start = time.time()
for n in range(samples):
    rand = random.uniform(0,1)
    rbs(0.5, 0.5)
end = time.time()
print 'RectBivariateSpline: %0.4f [us]' % (((end-start)/samples)*1e6)

start = time.time()
for n in range(samples):
    rand = random.uniform(0,1)
    RGI([0.5, 0.5])
end = time.time()
print 'RegularGridInterpolator: %0.4f [us]' % (((end-start)/samples)*1e6)

start = time.time()
for n in range(samples):
    map_coordinates(z2D.flatten(), coords, order=1)
end = time.time()
print 'Map Coordiantes: %0.4f [us]' % (((end-start)/samples)*1e6)

start = time.time()
for n in range(samples):
    rand = random.uniform(0,1)
    rbfi([0.5, 0.5])[0]
end = time.time()
print 'Rbf: %0.4f [us]' % (((end-start)/samples)*1e6)

start = time.time()
for n in range(samples):
    rand = random.uniform(0,1)
    CT2d([0.5, 0.5])[0]
end = time.time()
print 'CloughTocher2DInterpolator: %0.4f [us]' % (((end-start)/samples)*1e6)

start = time.time()
for n in range(int(samples/100)):
    rand = random.uniform(0,1)
    griddata(grd_inputs1, grd_inputs2, [0.5, 0.5], 'linear')
end = time.time()
print 'griddata: %0.4f [us]' % (((end-start)/int(samples/100))*1e6)

start = time.time()
for n in range(int(samples)):
    rand = random.uniform(0,1)
    #grdd(250e3+rand, bar2Pa(125+rand))[0]
    interpn([xrange, yrange], z2D, [0.5, 0.5])[0]
end = time.time()
print 'interpn: %0.4f [us]' % (((end-start)/int(samples))*1e6)

产生以下结果:
Exact:  0.75
interp2d 0.7500255076012651
RectBivariateSpline 0.7499999999999997
RegularGridInterpolator 0.750025507601265
Map Coordinates: 0.7500255076012652
Rbf 0.7500000434270742
CloughTocher2DInterpolator 0.750003857158422
griddata 0.7500255076012653
interpn 0.750025507601265

interp2d: 14.4000 [us]
RectBivariateSpline: 2.8400 [us]
RegularGridInterpolator: 90.9200 [us]
Map Coordiantes: 8.9600 [us]
Rbf: 105.8000 [us]
CloughTocher2DInterpolator: 18.0000 [us]
griddata: 138824.0004 [us]
interpn: 158.1600 [us]

我的实现有什么问题吗?这似乎难以置信,因为RectBivariateSpline是高阶拟合,并且相对于地图坐标允许非均匀网格间距。还有其他更快的二维插值方法吗?请原谅这段代码中可能存在的任何显著问题,我和Python刚刚开始接触。

谢谢,

1个回答

1

这种行为的解释似乎是因为每次调用不同的函数时,它们在Python中的开销不同。我们可以通过向量化调用来减少这种开销,如下所示的代码:

import numpy as np
import time
from scipy.interpolate import RectBivariateSpline, RegularGridInterpolator, CloughTocher2DInterpolator, Rbf, interpn
from scipy.ndimage import map_coordinates

grid_size = 100
xrange = np.linspace(0, 1, grid_size)
yrange = np.linspace(0, 1, grid_size)
x2D, y2D = np.meshgrid(xrange, yrange, indexing='ij')

z2D = x2D ** 2.0 + y2D

rbs = RectBivariateSpline(xrange, yrange, z2D,kx=1,ky=1)
RGI = RegularGridInterpolator((xrange, yrange), z2D)
CT2d = CloughTocher2DInterpolator(np.array(list(zip(x2D.flatten(), y2D.flatten()))), z2D.flatten())
rbfi = Rbf(x2D, y2D, z2D, function='linear')

samples = 25000
xpt = np.random.rand(samples)
ypt = np.random.rand(samples)

start = time.time()
rbs_vals = rbs.ev(xpt, ypt)
end = time.time()
print('RectBivariateSpline: %0.4f [us]' % (((end-start)/samples)*1e6))

start = time.time()
rgi_vals = RGI(np.asarray([xpt, ypt]).T)
end = time.time()
print('RegularGridInterpolator: %0.4f [us]' % (((end-start)/samples)*1e6))

start = time.time()
mc_vals = map_coordinates(z2D, np.asarray([xpt, ypt])*(grid_size-1), order=1)
end = time.time()
print('map_coordinates: %0.4f [us]' % (((end-start)/samples)*1e6))

start = time.time()
rbf_vals = rbfi(xpt,ypt)
end = time.time()
print('Rbf: %0.4f [us]' % (((end-start)/samples)*1e6))

start = time.time()
ct_vals = CT2d(np.asarray([xpt, ypt]).T)
end = time.time()
print('CloughTocher2DInterpolator: %0.4f [us]' % (((end-start)/samples)*1e6))

start = time.time()
in_vals = interpn([xrange, yrange], z2D, np.asarray([xpt, ypt]).T)
end = time.time()
print('interpn: %0.4f [us]' % (((end-start)/int(samples))*1e6))

vals = xpt**2 + ypt
print('\nRMS error:')
print('  RectBivariateSpline: %g' % np.sqrt(np.mean((rbs_vals-vals)**2)))
print('  RegularGridInterpolator: %g' % np.sqrt(np.mean((rgi_vals-vals)**2)))
print('  map_coordinates: %g' % np.sqrt(np.mean((mc_vals-vals)**2)))
print('  Rbf: %g' % np.sqrt(np.mean((rbf_vals-vals)**2)))
print('  CloughTocher2DInterpolator: %g' % np.sqrt(np.mean((ct_vals-vals)**2)))
print('  interpn: %g' % np.sqrt(np.mean((in_vals-vals)**2)))

这是一个例子的输出:
RectBivariateSpline: 0.1998 [us]
RegularGridInterpolator: 0.3595 [us]
map_coordinates: 0.0799 [us]
Rbf: 114.0414 [us]
CloughTocher2DInterpolator: 7.0727 [us]
interpn: 0.2797 [us]

RMS error:
  RectBivariateSpline: 1.85823e-05
  RegularGridInterpolator: 1.85823e-05
  map_coordinates: 1.85823e-05
  Rbf: 2.00083e-05
  CloughTocher2DInterpolator: 1.63274e-06
  interpn: 1.85823e-05

map_coordinates是最快的,正如你所期望的那样,尽管RectBivariateSpline也不落后(请注意,我将其阶数降低到1)。


网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接