Python调用实现最小二乘法的方法详解

2026-04-02 21:05:29

最小二乘法是非常经典的算法，而且这个名字我们在高中的时候就已经接触了，属于极其常用的算法。此前曾经写过线性最小二乘法的原理，并用Python实现：最小二乘法及其Python实现；以及scipy中非线性最小二乘法的调用方式：非线性最小二乘法（文末补充内容）；还有稀疏矩阵的最小二乘法：稀疏矩阵最小二乘法。

下面讲对numpy和scipy中实现的线性最小二乘法进行说明，并比较二者的速度。

numpy实现

numpy中便实现了最小二乘法，即lstsq(a,b)用于求解类似于a@x=b中的x，其中，a为M×N的矩阵；则当b为M行的向量时，刚好相当于求解线性方程组。对于Ax=b这样的方程组，如果A是满秩仿真，那么可以表示为x=A⁻¹b，否则可以表示为x=(A^TA)⁻¹A^Tb。

当b为M×K的矩阵时，则对每一列，都会计算一组x。

其返回值共有4个，分别是拟合得到的x、拟合误差、矩阵a的秩、以及矩阵a的单值形式。

import numpy as np
np.random.seed(42)
M = np.random.rand(4,4)
x = np.arange(4)
y = M@x
xhat = np.linalg.lstsq(M,y)
print(xhat[0])
#[0. 1. 2. 3.]

scipy封装

scipy.linalg同样提供了最小二乘法函数，函数名同样是lstsq，其参数列表为

lstsq(a, b, cond=None, overwrite_a=False, overwrite_b=False, check_finite=True, lapack_driver=None)

其中a, b即Ax=b，二者均提供可覆写开关，设为True可以节省运行时间，此外，函数也支持有限性检查，这是linalg中许多函数都具备的选项。其返回值与numpy中的最小二乘函数相同。

cond为浮点型参数，表示奇异值阈值，当奇异值小于cond时将舍弃。

lapack_driver为字符串选项，表示选用何种LAPACK中的算法引擎，可选'gelsd', 'gelsy', 'gelss'。

import scipy.linalg as sl
xhat1 = sl.lstsq(M, y)
print(xhat1[0])
# [0. 1. 2. 3.]

速度对比

最后，对着两组最小二乘函数做一个速度上的对比

from timeit import timeit
N = 100
A = np.random.rand(N,N)
b = np.arange(N)

timeit(lambda:np.linalg.lstsq(A, b), number=10)
# 0.015487500000745058
timeit(lambda:sl.lstsq(A, b), number=10)
# 0.011151800004881807

这一次，二者并没有拉开太大的差距，即使将矩阵维度放大到500，二者也是半斤八两。

N = 500
A = np.random.rand(N,N)
b = np.arange(N)

timeit(lambda:np.linalg.lstsq(A, b), number=10)
0.389679799991427
timeit(lambda:sl.lstsq(A, b), number=10)
0.35642060000100173

补充

Python调用非线性最小二乘法

简介与构造函数

在scipy中，非线性最小二乘法的目的是找到一组函数，使得误差函数的平方和最小，可以表示为如下公式

其中ρ表示损失函数，可以理解为对f_i(x)的一次预处理。

scipy.optimize中封装了非线性最小二乘法函数least_squares，其定义为

least_squares(fun, x0, jac, bounds, method, ftol, xtol, gtol, x_scale, f_scale, loss, jac_sparsity, max_nfev, verbose, args, kwargs)

其中，func和x0为必选参数，func为待求解函数，x0为函数输入的初值，这两者无默认值，为必须输入的参数。

bound为求解区间，默认(−∞,∞)，verbose为1时，会有终止输出，为2时会print更多的运算过程中的信息。此外下面几个参数用于控制误差，比较简单。

	默认值	备注
ftol	10^-8	函数容忍度
xtol	10^-8	自变量容忍度
gtol	10^-8	梯度容忍度
x_scale	1.0	变量的特征尺度
f_scale	1.0	残差边际值

loss为损失函数，就是上面公式中的ρ \rhoρ，默认为linear，可选值包括

迭代策略

上面的公式仅给出了算法的目的，但并未暴露其细节。关于如何找到最小值，则需要确定搜索最小值的方法，method为最小值搜索的方案，共有三种选项，默认为trf

trf：即Trust Region Reflective，信赖域反射算法
dogbox：信赖域狗腿算法
lm：Levenberg-Marquardt算法

这三种方法都是信赖域方法的延申，信赖域的优化思想其实就是从单点的迭代变成了区间的迭代，由于本文的目的是介绍scipy中所封装好的非线性最小二乘函数，故而仅对其原理做简略的介绍。

其中r为置信半径，假设在这个邻域内，目标函数可以近似为线性或二次函数，则可通过二次模型得到区间中的极小值点s_k。然后以这个极小值点为中心，继续优化信赖域所对应的区间。

以上就是信赖域方法的基本原理。

雅可比矩阵

在了解了信赖域方法之后，就会明白雅可比矩阵在数值求解时的重要作用，而如何计算雅可比矩阵，则是接下来需要考虑的问题。jac参数为计算雅可比矩阵的方法，主要提供了三种方案，分别是基于两点的2-point；基于三点的3-point；以及基于复数步长的cs。一般来说，三点的精度高于两点，但速度也慢一倍。

此外，可以输入自定义函数来计算雅可比矩阵。

测试

最后，测试一下非线性最小二乘法

import numpy as np
from scipy.optimize import least_squares

def test(xs):
    _sum = 0.0
    for i in range(len(xs)):
        _sum = _sum + (1-np.cos((xs[i]*i)/5)*(i+1))
    return _sum

x0 = np.random.rand(5)
ret = least_squares(test, x0)
msg = f"最小值" + ", ".join([f"{x:.4f}" for x in ret.x])
msg += f"\nf(x)={ret.fun[0]:.4f}"
print(msg)
'''
最小值0.9557, 0.5371, 1.5714, 1.6931, 5.2294
f(x)=0.0000
'''

到此这篇关于Python调用实现最小二乘法的方法详解的文章就介绍到这了,更多相关Python最小二乘法内容请搜索我们以前的文章或继续浏览下面的相关文章希望大家以后多多支持我们！

Python实现曲线拟合的最小二乘法

本文实例为大家分享了Python曲线拟合的最小二乘法,供大家参考,具体内容如下模块导入 import numpy as np import gaosi as gs 代码 """ 本函数通过创建增广矩阵,并调用高斯列主元消去法模块进行求解. """ import numpy as np import gaosi as gs shape = int(input('请输入拟合函数的次数:')) x = np.array([0.6,1.3,1.64,1
Python最小二乘法矩阵

最小二乘法矩阵 #! /usr/bin/env python # -*- coding: utf-8 -*- import numpy as np def calc_left_k_mat(k): """ 获得左侧k矩阵 :param k: :return: """ k_mat = [] for i in range(k + 1): now_line = [] for j in range(k + 1): now_line.append(j + i
最小二乘法及其python实现详解

最小二乘法Least Square Method,做为分类回归算法的基础,有着悠久的历史(由马里·勒让德于1806年提出).它通过最小化误差的平方和寻找数据的最佳函数匹配.利用最小二乘法可以简便地求得未知的数据,并使得这些求得的数据与实际数据之间误差的平方和为最小.最小二乘法还可用于曲线拟合.其他一些优化问题也可通过最小化能量或最大化熵用最小二乘法来表达. 那什么是最小二乘法呢?别着急,我们先从几个简单的概念说起. 假设我们现在有一系列的数据点 ,那么由我们给出的拟合函数h(x)得到的估计量就是
如何使用Python最小二乘法拟合曲线代码详解

目录一.背景描述二.前期准备最小二乘法 Python类库 NumPy Matplotlib 三.代码实现总结一.背景描述在一个普通的摸鱼早晨,群里居然出现了一个不合时宜颇为突兀的正经问题,原来是一个博士同学需要対实验数据进行曲线拟合并且批量计算出多项式方程一般来说,这种问题对于经常做实验的同学来说并不陌生,通常使用MATLAB或者Origin Pro这类专业的数据计算软件,甚至Excel也可以实现. 但是作为程序员肯定第一想到的还是使用强大的Python来实现,但是因为本人主要做j
Python实现两种稀疏矩阵的最小二乘法

目录最小二乘法返回值测试最小二乘法 scipy.sparse.linalg实现了两种稀疏矩阵最小二乘法lsqr和lsmr,前者是经典算法,后者来自斯坦福优化实验室,据称可以比lsqr更快收敛. 这两个函数可以求解Ax=b,或arg minx ∥Ax−b∥2,或arg minx ∥Ax−b∥2 +d2∥x−x0∥2,其中A必须是方阵或三角阵,可以有任意秩. 通过设置容忍度at ,bt,可以控制算法精度,记r=b-Ax 为残差向量,如果Ax=b是相容的,lsqr在∥r∥⩽at∗∥A∥⋅∥x∥
Python调用C++程序的方法详解

前言大家都知道Python的优点是开发效率高,使用方便,C++则是运行效率高,这两者可以相辅相成,不管是在Python项目中嵌入C++代码,或是在C++项目中用Python实现外围功能,都可能遇到Python调用C++模块的需求,下面列举出集中c++代码导出成Python接口的几种基本方法,一起来学习学习吧. 原生态导出 Python解释器就是用C实现,因此只要我们的C++的数据结构能让Python认识,理论上就是可以被直接调用的.我们实现test1.cpp如下 #include <Pytho
对python 调用类属性的方法详解

测试时候类的调用是经常会用到的.简单看下类的调用使用的方法吧. 来看例子: 目录结构: 我们现在要在do_class.py这个文件里调用class_learn.py里的类代码(do_class.py): #!/usr/bin/env python3 #coding=utf-8 '''@Author:Jock''' from all_python_learn.class_and_funcation.class_learn import * b = Learn(1,2) b.get() print
Python对象类型及其运算方法(详解)

基本要点: 程序中储存的所有数据都是对象(可变对象:值可以修改不可变对象:值不可修改) 每个对象都有一个身份.一个类型.一个值例: >>> a1 = 'abc' >>> type(a1) str 创建一个字符串对象,其身份是指向它在内存中所处的指针(在内存中的位置) a1就是引用这个具体位置的名称使用type()函数查看其类型其值就是'abc' 自定义类型使用class 对象的类型用于描述对象的内部表示及其支持的方法和操作创建特定类型的对象,也将该对象称为该类
Python 常用模块 re 使用方法详解

一.re模块的查找方法: 1.findall 匹配所有每一项都是列表中的一个元素 import re ret = re.findall('\d+','asd鲁班七号21313') # 正则表达式,待匹配的字符串,flag # ret = re.findall('\d','asd鲁班七号21313') # 正则表达式,待匹配的字符串,flag # print(ret) 2.search 只匹配从左到右的第一个,等到的不是直接的结果,而是一个变量,通过这个变量的group方法来获取结果 impo
对python调用RPC接口的实例详解

要调用RPC接口,python提供了一个框架grpc,这是google开源的 rpc相关文档: https://grpc.io/docs/tutorials/basic/python.html 需要安装的python包如下: 1.grpc安装 pip install grpcio 2.grpc的python protobuf相关的编译工具 pip install grpcio-tools 3.protobuf相关python依赖库 pip install protobuf 4.一些常见原型的生成
对Python捕获控制台输出流的方法详解

有时候我们的代码里可能要调用控制台命令,比如我想用Python写一个批量编译 .java 文件的脚本,用到如下代码常规用法 os.system import os,traceback try: p = os.system("javac Test.java") print p except: print "\nexcept:\n" print traceback.format_exc() 如然编译成功会返回一个0,如果错误会返回一个非0的值给p,这种方法可以知道执行
Python安装依赖(包)模块方法详解

Python模块,简单说就是一个.py文件,其中可以包含我们需要的任意Python代码.迄今为止,我们所编写的所有程序都包含在单独的.py文件中,因此,它们既是程序,同时也是模块.关键的区别在于,程序的设计目标是运行,而模块的设计目标是由其他程序导入并使用. 不是所有程序都有相关联的.py文件-比如说,sys模块就内置于Python中,还有些模块是使用其他语言(最常见的是C语言)实现的.不过,Python的大多数库文件都是使用Python实现的,因此,比如说,我们使用了语句import coll
一个Python优雅的数据分块方法详解

目录 1.背景 2.islice 2.1示例 2.2只指定步长 3.iter 3.1常规使用 3.2进阶使用 4.islice 和 iter 组合使用 5.总结 1.背景看到这个标题你可能想一个分块能有什么难度?还值得细说吗,最近确实遇到一个有意思的分块函数,写法比较巧妙优雅,所以写一个分享. 日前在做需求过程中有一个对大量数据分块处理的场景,具体来说就是几十万量级的数据,分批处理,每次处理100个.这时就需要一个分块功能的代码,刚好项目的工具库中就有一个分块的函数.拿过函数来用,发现还挺好用
Python远程控制Windows服务器的方法详解

目录 1. 被控端 windows 启动 winrm 服务检查 winrm 服务监听状态查看 winrm 配置信息(可选) 配置 winrm client 配置 winrm service 2. 控制端 3. 实战一下 4. 总结在很多企业会使用闲置的 Windows 机器作为临时服务器,有时候我们想远程调用里面的程序或查看日志文件 Windows 内置的服务「 winrm 」可以满足我们的需求它是一种基于标准简单对象访问协议( SOAP )的防火墙友好协议,允许来自不同供应商的硬件和操
Python实现创建模块的方法详解

目录楔子 __import__ importlib.machinery 通过 module 类创建模块将一个类的实例变成一个模块小结楔子导入一个模块,我们一般都会使用 import 关键字,但有些场景下 import 难以满足我们的需要.所以除了 import 之外还有很多其它导入模块的方式,下面就来介绍一下. __import__ 这是一个内置函数,解释器在 import 的时候,实际上就执行了这个函数. # import os 等价于如下方式 os = __import__("os

Python调用实现最小二乘法的方法详解

目录

numpy实现

scipy封装

速度对比

补充

相关推荐

随机推荐