NumPy PolyFit and PolyVal in Multiple Dimensions?

arrays, numpy, python, scipy

Solution

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

>>> 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

Problem

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): ``` 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. 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? This is along the lines of this SO question.

Original source

Related problems