使用numpy进行多变量多项式回归

29

我有许多样本(y_i, (a_i, b_i, c_i)),其中y被假定为关于a,b,c的多项式,达到一定的次数。例如,对于给定的数据集和二次项,我可能会产生以下模型

y = a^2 + 2ab - 3cb + c^2 +.5ac

可以使用最小二乘法来实现这一点,并且是numpy中polyfit例程的轻微扩展。在Python生态系统中是否有标准实现?


2
我在这里发布了代码来解决这个问题https://github.com/mrocklin/multipolyfit。 - MRocklin
3个回答

21
sklearn 提供了一种简单的方法来实现这一点。
在这里发布的示例的基础上构建:
import numpy as np
from sklearn.preprocessing import PolynomialFeatures
from sklearn import linear_model

#X is the independent variable (bivariate in this case)
X = np.array([[0.44, 0.68], [0.99, 0.23]])

#vector is the dependent data
vector = np.array([109.85, 155.72])

#predict is an independent variable for which we'd like to predict the value
predict= np.array([[0.49, 0.18]])

#generate a model of polynomial features
poly = PolynomialFeatures(degree=2)

#transform the x data for proper fitting (for single variable type it returns,[1,x,x**2])
X_ = poly.fit_transform(X)

#transform the prediction to fit the model type
predict_ = poly.fit_transform(predict)

#here we can remove polynomial orders we don't want
#for instance I'm removing the `x` component
X_ = np.delete(X_,(1),axis=1)
predict_ = np.delete(predict_,(1),axis=1)

#generate the regression object
clf = linear_model.LinearRegression()
#preform the actual regression
clf.fit(X_, vector)

print("X_ = ",X_)
print("predict_ = ",predict_)
print("Prediction = ",clf.predict(predict_))

这是输出结果:
>>> X_ =  [[ 0.44    0.68    0.1936  0.2992  0.4624]
>>>  [ 0.99    0.23    0.9801  0.2277  0.0529]]
>>> predict_ =  [[ 0.49    0.18    0.2401  0.0882  0.0324]]
>>> Prediction =  [ 126.84247142]

你能否包含delete函数的实现呢?谢谢! - Shivam Gaur
1
抱歉,这是numpy的文档,https://docs.scipy.org/doc/numpy/reference/generated/numpy.delete.html - David Hoffman
PolynomialFeatures是什么?它的作用是什么?我可以看到代码吗? - Charlie Parker
1
这对我来说没有意义,为什么 fit_transform 既返回多项式特征矩阵(范德蒙矩阵),又返回预测值? :/ - Charlie Parker
这与手动执行如 c_pinv = np.dot(np.linalg.pinv( Kern_train ),Y_train) 相比如何? - Charlie Parker
predict variable should be a 2d array [[0.49, 0.18]] - Luca


2
polyfit可以使用,但是有更好的最小二乘解决方案。我建议使用kmpfit,可以在以下网址找到:http://www.astro.rug.nl/software/kapteyn-beta/kmpfittutorial.html。它比polyfit更稳健,并且他们的页面上有一个示例,展示了如何进行简单的线性拟合,这应该提供了进行二次多项式拟合的基础知识。

def model(p, v, x, w):       
   a,b,c,d,e,f,g,h,i,j,k = p      #coefficients to the polynomials      
   return  a*v**2 + b*x**2 + c*w**2 + d*v*x + e*v*w + f*x*w + g*v + h*x + i*y + k
def residuals(p, data): # Function needed by fit routine v, x, w, z = data # The values for v, x, w and the measured hypersurface z a,b,c,d,e,f,g,h,i,j,k = p #coefficients to the polynomials return (z-model(p,v,x,w)) # Returns an array of residuals. #This should (z-model(p,v,x,w))/err if # there are error bars on the measured z values
#initial guess at parameters. Avoid using 0.0 as initial guess par0 = [1.0, 1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0,1.0]
#create a fitting object. data should be in the form #that the functions above are looking for, i.e. a Nx4 #list of lists/tuples like (v,x,w,z) fitobj = kmpfit.Fitter(residuals=residuals, data=data)
# call the fitter fitobj.fit(params0=par0)

这些事情的成功与拟合的起始值密切相关,因此如果可能的话,请仔细选择。由于有太多的自由参数,因此可能很难得出解决方案。


1
你能否提供一个使用polyfit进行多元回归的示例?我不确定是否支持这个功能。在查看kmpfit文档后,我担心这个库也可能不支持。 - MRocklin
你想拟合什么,y(x) = ax**2 + bx + c吗?无论如何,你可以使用mpfit/kmpfit进行多变量拟合。 - reptilicus
不,y(v, x, w) = av**2 + bx2 + c*w2 + dvx + evw + fxw + gv + hx + i*y + k - MRocklin
2
这个库可以工作,但是它通过迭代方法解决问题。最小二乘多项式拟合可以通过解线性系统一步完成。我在另一个答案中发布了使用NumPy实现此操作的代码。 - MRocklin

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