多维度的NumPy PolyFit和PolyVal？

Assume an n-dimensional array of observations that are reshaped to be a 2d-array with each row being one observation set. Using this reshape approach, np.polyfit can compute 2nd order fit coefficients for the entire ndarray (vectorized):

假设n维的观察阵列被重新整形为2d阵列,每行是一个观察集。使用这种重塑方法,np.polyfit可以计算整个ndarray(矢量化)的二阶拟合系数:

fit = np.polynomial.polynomialpolyfit(X, Y, 2)

where Y is shape (304000, 21) and X is a vector. This results in a (304000,3) array of coefficients, fit.

其中Y是形状(304000,21),X是矢量。这导致(304000,3)系数数组合适。

Using an iterator it is possible to call np.polyval(fit, X) for each row. This is inefficient when a vectorized approach may exist. Could the fit result be applied to the entire observation array without iterating? If so, how?

使用迭代器可以为每一行调用np.polyval(fit,X)。当存在矢量化方法时,这是低效的。拟合结果是否可以应用于整个观察阵列而无需迭代?如果是这样,怎么样?

This is along the lines of this SO question.

这与这个SO问题一致。

2 个解决方案

#1

np.polynomial.polynomial.polyval takes multidimensional coefficient arrays:

np.polynomial.polynomial.polyval采用多维系数数组:

>>> x = np.random.rand(100)
>>> y = np.random.rand(100, 25)
>>> fit = np.polynomial.polynomial.polyfit(x, y, 2)
>>> fit.shape # 25 columns of 3 polynomial coefficients
(3L, 25L)
>>> xx = np.random.rand(50)
>>> interpol = np.polynomial.polynomial.polyval(xx, fit)
>>> interpol.shape # 25 rows, each with 50 evaluations of the polynomial
(25L, 50L)

And of course:

而且当然:

>>> np.all([np.allclose(np.polynomial.polynomial.polyval(xx, fit[:, j]),
...                     interpol[j]) for j in range(25)])
True

#2

np.polynomial.polynomial.polyval is a perfectly fine (and convenient) approach to efficient evaluation of polynomial fittings.

np.polynomial.polynomial.polyval是一种完美的(方便)方法,可有效评估多项式拟合。

However, if 'speediest' is what you are looking for, simply constructing the polynomial inputs and using the rudimentary numpy matrix multiplication functions results in slightly faster ( roughly 4x faster) computational speeds.

但是,如果您正在寻找'最快',那么简单地构造多项式输入并使用基本的numpy矩阵乘法函数会导致计算速度稍快(大约快4倍)。

Setup

Using the same setup as above, we'll create 25 different line fittings.

使用与上面相同的设置,我们将创建25个不同的线路配件。

>>> num_samples = 100000
>>> num_lines = 100
>>> x = np.random.randint(0,100,num_samples)
>>> y = np.random.randint(0,100,(num_samples, num_lines))
>>> fit = np.polyfit(x,y,deg=2)
>>> xx = np.random.randint(0,100,num_samples*10)

Numpy's polyval Function

res1 = np.polynomial.polynomial.polyval(xx, fit)

Basic Matrix Multiplication

inputs = np.array([np.power(xx,d) for d in range(len(fit))])
res2 = fit.T.dot(inputs)

Timing the functions

Using the same parameters above...

使用上面相同的参数......

%timeit _ = np.polynomial.polynomial.polyval(xx, fit)
1 loop, best of 3: 247 ms per loop

%timeit inputs = np.array([np.power(xx, d) for d in range(len(fit))]);_ = fit.T.dot(inputs)
10 loops, best of 3: 72.8 ms per loop

To beat a dead horse...

击败死马......

mean Efficiency bump of ~3.61x faster. Speed fluctuations probably come from random computer processes in background.

意味着效率提升约3.61倍。速度波动可能来自后台的随机计算机进程。