Python SciPy optimize模块详解:非线性方程、最小二乘与全局优化
本文系统介绍Python SciPy库中optimize模块的核心用法。涵盖非线性方程组求解(fsolve)、最小二乘拟合(leastsq与curve_fit)、局部极小值算法(Nelder-Mead等)及全局优化策略(basinhopping)。通过具体代码示例解析参数传递、雅可比矩阵加速及初始值敏感性等关键问题,为数据拟合与数值优化提供实用参考。
SciPy optimize模块概述
前言
SciPy库中的optimize模块封装了多种高效的数值优化算法,本文将对其中常用功能进行梳理与记录。
非线性方程组求解
非线性方程组求解
在SciPy中,求解非线性方程组主要依赖fsolve()函数。其标准调用格式为fsolve(fun, x0),其中fun是定义方程组误差的计算函数,接收参数x并返回各方程的残差值;x0则是变量x的初始猜测值。
"""
计算非线性方程组:
5x1+3 = 0
4x0^2-2sin(x1x2)=0
x1x2-1.5=0
"""
## 误差函数
def fun(x):
x0,x1,x2 = x.tolist()
return[5*x1+3,4x0^2-2sin(x1x2),x1x2-1.5]
result = optimize.fsolve(fun,[1,1,1])
## result
[-0.70622057 -0.6 -2.5]
在求解非线性方程时,若涉及类似坐标上升算法的场景,往往需要利用未知数的导数信息。同样,scipy的fslove()也提供了fprime参数,用于传入未知数的雅各比矩阵以加速计算过程。该雅各比矩阵的每一行代表对应方程对各未知变量的偏导数。针对前述示例,可构造如下雅各比矩阵并传入:
def j(x):
x0,x1,x2 = x.tolist()
return[[0,5,0],[8*x0,-2*x2*cos(x1*x2],[0,x2,x1]]
result = optimize.fsolve(fun,[1,1,1],fprime=j)
#result
[-0.70622057 -0.6 -2.5]
关于scipy内部实现fslove时具体采用的算法(如是否结合坐标上升或梯度相关优化),目前尚未经过详细考证,感兴趣读者可直接查阅源代码。
最小二乘拟合方法
最小二乘拟合
此处不展开最小二乘算法的理论推导,相关解释资料已十分丰富。在optimize模块中,可利用leastsq()函数对数据进行最小二乘拟合。该函数使用简便,只需传入计算误差的函数及待确定参数的初始值即可。
x = np.array([8.19,2.72,6.39,8.71,4.7,2.66,3.78])
y = np.array([7.01,2.78,6.47,6.71,4.1,4.23,4.05])
def residual(p):
k,b = p
return y-(k*x+b)
r = optimize.leastsq(residual,[1,0])
k,b = r[0]
# print k
.613495349193
# print b
.79409254326
def func(x,p):
"""
计算的正弦波 :A*sin(2*pi*k*x+theta)
"""
A,k,theta = p
return A*sin(2*np.pi*k*x+theta)
def redis(p,y,x):
return y-func(x,p)
x = np.linspace(0,2*np.pi,100)
A,k,theta = 10,0.34,np.pi/6
y0 = func(x,[A,k,theta])
# 加入噪声
np.random.seed(0)
y1 = y0+2*np.random.randn(len(x))
p0 = [7,0.40,0]
# p0是A,k,theta的初始值,y1,x要拟合的数据
plsq = optimize.leastsq(redis, p0,args=(y1,x))
print [A,k,theta] #真是的参数值
print plsq[0] #拟合后的参数值
针对正弦波或余弦波等曲线拟合任务,optimize模块提供了curve_fit()函数。其用法与leastsq()略有差异,该函数直接计算曲线拟合值。例如,上述正弦波拟合可通过curve_fit()实现:
def func2(x,p):
"""
计算的正弦波 :A*sin(2*pi*k*x+theta)
"""
A,k,theta = p
return A*sin(2*np.pi*k*x+theta)
ret,_=optimize.curve_fit(func2,x,y1,p0=p0)
需注意,该函数对初始值较为敏感。若初始频率与真实频率偏差过大,可能导致算法无法收敛至真实频率值。
局部与全局最小值算法
局部最小值算法
optimize模块还提供了多种常用的局部最小值算法,包括Nelder-Mead、Powell、CG、BFGS、Newton-CG等。在计算最小值时,通常可传入一阶导数矩阵(雅可比矩阵)或二阶导数矩阵(黑塞矩阵)以加速收敛。需注意的是,这些算法通常无法保证收敛至全局最小值,多数情况下仅能收敛至局部极小值。其调用方式如下:
optimize.minimize(target_fun,init_val,method,jac,hess)
target_fun:函数的表达式计算;
init_val:初始值;
method:最小化的算法;
jac:雅各比矩阵
hess:黑塞矩阵。
全局最小值算法
全局最小值可通过optimize.basinhopping()函数实现。该函数首先需定义误差计算方式(如平方误差函数),指定迭代次数niter,并通过minimizer_kwargs传入局部极小值优化方法。例如,针对前述正弦函数拟合可编写如下代码:
def func1(x,p):
"""
计算的正弦波 :A*sin(2*pi*k*x+theta)
"""
A,k,theta = p
return A*sin(2*np.pi*k*x+theta)
def func_error(p,y,x):
return np.sum((y-func1(x,p)**2)
result = optimize.basinhopping(func_error,[1,1,1],niter=10,
minimizer_kwargs={"method":"L-BFGS-B",
"args":(y1,x1)})
## [1,1,1]是传入的初始值,args是需要拟合的数据
以上关于Python科学计算中SciPy optimize模块用法的介绍到此结束,希望能为读者提供有价值的参考。
本站声明
本文内容仅供参考,不代表本站观点。转载或引用请注明出处。如发现内容问题,请通过联系我们页面反馈。













