Levenberg-Marquardt (LM) 算法进行非线性拟合


目录

  • 1. LM算法
  • 2. 调包实现
  • 3. LM算法实现
  • 4. 源码地址


1. LM算法

LM算法是一种非线性最小二乘优化算法,用于求解非线性最小化问题。LM主要用于解决具有误差函数的非线性最小二乘问题,其中误差函数是参数的非线性函数,需要通过调整参数使误差函数最小化。算法的基本思想是通过迭代的方式逐步调整参数,使得误差函数在参数空间中逐渐收敛到最小值。在每一次迭代中,算法通过求解一个线性方程组来更新参数。这个线性方程组由误差函数的雅可比矩阵和参数更新量构成。

LM算法的优点在于它能够快速收敛到局部最小值,并且对于初始参数的选择不太敏感。此外,算法还能够处理参数个数多于观测数据个数的问题,并且对于存在噪声的数据也比较鲁棒。

2. 调包实现

如图1所示,调用scipy.optimize的least_squares函数实现对测试函数 exp ⁡ ( − a x 2 − b y 2 ) \exp(-ax^2-by^2) exp(ax2by2)的拟合结果。目标参数为 [ 0.5 , 0.5 ] [0.5, 0.5] [0.5,0.5],初始参数设置为 [ 1.0 , 1.0 ] [1.0, 1.0] [1.0,1.0],经过22次迭代,由于观测值暂未添加噪声,所以最终拟合参数与目标参数完全一致。

在这里插入图片描述

Fig. 1. 三维目标拟合: $\exp(-ax^2-by^2)$

3. LM算法实现

使用Python对LM做了简单实现,并对测试函数 exp ⁡ ( a x 2 + b x + c ) \exp(ax^2+bx+c) exp(ax2+bx+c)进行拟合,观测值添加高斯噪声。目标参数为 [ 1.0 , 2.0 , 3.0 ] [1.0, 2.0, 3.0] [1.0,2.0,3.0],初始参数设置为 [ 3.0 , 9.0 , 6.0 ] [3.0, 9.0, 6.0] [3.0,9.0,6.0],经过41次迭代,拟合参数为 [ 2.0 , 0.6 , 3.5 ] [2.0, 0.6, 3.5] [2.0,0.6,3.5],MSE损失小于0.000001,符合拟合误差要求。图2绘制了第12(蓝),13(黄),15(绿)次迭代结果以及最终拟合结果(红)。

在这里插入图片描述

Fig. 2. 二维目标拟合: $\exp(ax^2+bx+c)$
# 部分函数代码:def Func(abc,iput):   # 需要拟合的函数,abc是包含三个参数的一个矩阵[[a],[b],[c]]a = abc[0,0]b = abc[1,0]c = abc[2,0]return np.exp(a*iput**2+b*iput+c)def Deriv(abc,iput,n):  # 对函数求偏导x1 = abc.copy()x2 = abc.copy()x1[n,0] -= 0.000001x2[n,0] += 0.000001p1 = Func(x1,iput)p2 = Func(x2,iput)d = (p2-p1)*1.0/(0.000002)return dxk_l = []  # 用来存放每次迭代的结果
while conve:mse,mse_tmp = 0,0step += 1  fx = Func(xk,h) - ymse += sum(fx**2)for j in range(3): J[:,j] = Deriv(xk,h,j) # 数值求导                                                    mse /= n  # 范围约束H = J.T*J + u*np.eye(3)   # 3*3dx = -H.I * J.T*fx        # xk_tmp = xk.copy()xk_tmp += dxfx_tmp =  Func(xk_tmp,h) - y  mse_tmp = sum(fx_tmp[:,0]**2)mse_tmp /= n#判断是否下降q = float((mse - mse_tmp)/((0.5*dx.T*(u*dx - J.T*fx))[0,0]))if q > 0:s = 1.0/3.0v = 2mse = mse_tmpxk = xk_tmptemp = 1 - pow(2*q-1,3)if s > temp:u = u*selse:u = u*tempelse:u = u*vv = 2*vxk = xk_tmpprint ("step = %d,abs(mse-lase_mse) = %.8f" %(step,abs(mse-lase_mse)))  if abs(mse-lase_mse)<0.000001:breaklase_mse = mse  # 记录上一个 mse 的位置conve -= 1xk_l.append(xk)

4. 源码地址

如果对您有用的话可以点点star哦~

https://github.com/Jurio0304/cs-math/blob/main/hw4_LM.ipynb


创作不易,麻烦点点赞和关注咯!

本文来自互联网用户投稿,该文观点仅代表作者本人,不代表本站立场。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如若转载,请注明出处:http://www.mzph.cn/news/829825.shtml

如若内容造成侵权/违法违规/事实不符,请联系多彩编程网进行投诉反馈email:809451989@qq.com,一经查实,立即删除!

相关文章

Vue Canvas图片水印的绘制 图片加水印

效果 定义画布 <canvas width"800" height"800" ref"cn" ></canvas>绘制水印 draw(){const img new Image()img.srchttps://img1.baidu.com/it/u3035183739,1826404114&fm253&fmtauto&app138&fJPEGimg.onload(()…

pyqt 动态更换表头和数据

目录 pyqt 动态更换表头和数据代码 效果图&#xff1a; pyqt 动态更换表头和数据代码 from PyQt5.QtGui import QColor, QBrush from PyQt5.QtWidgets import QApplication, QTableWidget, QVBoxLayout, QWidget, QPushButton, QTableWidgetItemclass Example(QWidget):def _…

Redis底层数据结构之ZSkipList

目录 一、概述二、ZSkipList结构三、和平衡树和哈希表的对比 redis底层数据结构已完结&#x1f44f;&#x1f44f;&#x1f44f;&#xff1a; ☑️redis底层数据结构之SDS☑️redis底层数据结构之ziplist☑️redis底层数据结构之quicklist☑️redis底层数据结构之Dict☑️redis…

机器学习和深度学习-- 李宏毅(笔记与个人理解)Day22

Day 22 Transformer seqence to seqence 有什么用呢&#xff1f; Encoder how Block work 仔细讲讲Residual 的过程&#xff1f; 重构 Decoder - AutoRegressive Mask 由于是文字接龙&#xff0c;所以无法考虑右边的 info 另一种decoder Encoder to Decoder – Cross Attend…

llama3本地部署

目录 II.下载 II.验证ollama安装 II.安装llama3 和启动 II.命令行调用 II.api调用 II.参考文献 II.下载 https://ollama.com/download/windows OllamaSetup.exe https://github.com/meta-llama/llama3 II.验证ollama安装 cmd ollama II.安装llama3 和启动 ollama run …

【问题分析】TaskDisplayArea被隐藏导致的黑屏以及无焦点窗口问题【Android 14】

1 问题描述 用户操作出的偶现的黑屏以及无焦点窗口问题。 直接原因是&#xff0c;TaskDisplayArea被添加了eLayerHidden标志位&#xff0c;导致所有App的窗口不可见&#xff0c;从而出现黑屏和无焦点窗口问题&#xff0c;相关log为&#xff1a; 这个log是MTK添加的&#xff0…

Django模型继承之多表继承

在Django模型继承中&#xff0c;支持的第二种模型继承方式是层次结构中的每个模型都是一个单独的模型。每个模型都指向分离的数据表&#xff0c;并且可以被独立查询和创建。在继承关系中&#xff0c;子类和父类之间通过一个自动创建的OneToOneField进行连接。示例代码如下&…

C语言入门课程学习笔记-6

C语言入门课程学习笔记-6 第27课 - 字符数组与字符串&#xff08;上&#xff09;第28课 - 字符数组与字符串&#xff08;下&#xff09;第29课 - 数组专题练习&#xff08;上&#xff09;第30课 - 数组专题练习&#xff08;下&#xff09; 本文学习自狄泰软件学院 唐佐林老师的…

不只有 Spring,这四款Java 基础开发框架同样值得关注!

Java 开发不只有 Spring &#xff0c;今天给大家推荐几个同样优秀的 Java 基础开发框架&#xff0c;为日常项目开发提供更多的选择。答应我&#xff0c;请不要再叫我 Spring 小子了&#xff0c;​好吗&#xff1f; 项目概览&#xff1a; Guice&#xff1a;轻量级依赖注入框架 …

2024Mac系统热门游戏排行榜 Mac支持的网络游戏有哪些?mac能玩哪些大型网游 苹果电脑Mac游戏资源推荐 Mac玩Windows游戏

“游戏是这个世界上唯一能和女性争夺男朋友的东西&#xff08;/滑稽&#xff0c;有不少女生也喜欢玩游戏&#xff09;。” 虽然只是一句玩笑话&#xff0c;不过也可以看出游戏对大多数男生来说是必不可少的一项娱乐活动了。而网络游戏是游戏中的一大分支&#xff0c;能让玩家们…

科技“冷”战:NIST刷新制冷效率,中国实力逆境崛起!

4月23日&#xff0c;美国国家标准与技术研究院&#xff08;NIST&#xff09;的研究人员报道称&#xff0c;他们通过对常用于科研和工业领域的制冷机进行改装&#xff0c;显著降低了将材料冷却至略高于绝对零度所需的时间和能量。 科学家们指出&#xff0c;他们的原型设备每年能…

Linux 学习之路 -- 进程篇 -- 进程控制

目录 一、进程终止 <1>使用语言和系统自带的方法&#xff0c;进行转换 <2>自定义错误码 <3>小结&#xff1a; <2>两个接口exit / _exit 二、进程等待 <1>简单了解 <2>wait调用 <3>waitpid调用 <4>status <1>W…

复杂的字符串算法——KMP算法

字符串算法 模式匹配&#xff08;Pattern Matching&#xff09;&#xff1a;在一篇长度为 &#x1d45b; 的文本 &#x1d446; 中&#xff0c;找某个长度为 &#x1d45a; 的关键词 &#x1d443;。&#x1d443; 可能多次出现&#xff0c;都需要找到。 最优的模式匹配算法复…

AHB传输---突发操作

突发操作 在本协议中定义了4拍、8拍和16拍的突发&#xff0c;以及未定义长度的突发和单次传输。它支持增量和包装突发&#xff1a; 增量突发访问连续位置&#xff0c;每个传输的地址是前一个地址的增量。包装突发在跨越地址边界时会包装。地址边界的计算方法是突发中拍数与传…

Android—统一依赖版本管理

依赖版本管理有多种方式 config.gradle 用于Groovy DSL&#xff0c;新建一个 config.gradle 文件&#xff0c;然后将项目中所有依赖写在里面&#xff0c;更新只需修改 config.gradle文件内容&#xff0c;作用于所有module buildSrc 可用于Kotlin DSL或Groovy DSL&#xff0c;…

MATLAB冒号表示法

MATLAB 冒号表示法 colon(:)是在MATLAB中最有用的运算符之一。它用于创建向量&#xff0c;下标数组和指定迭代。 如果要创建包含1到10的整数的行向量&#xff0c;请编写- 示例 1:10 MATLAB执行该语句并返回包含1到10的整数的行向量- ans 1 2 3 4 5 6 7 8 9 10 如果要指定一…

github Copilot的使用总结

1. 代码建议和补全 GitHub Copilot 的基本使用涉及编写代码时的实时代码建议和补全。一旦你已经安装并配置好 GitHub Copilot 插件&#xff0c;你可以在支持的编辑器&#xff08;如 Visual Studio Code&#xff09;中开始使用 Copilot。以下是一些基本的使用步骤&#xff1a; …

VBA技术资料MF146:发出多次Beep提示声

我给VBA的定义&#xff1a;VBA是个人小型自动化处理的有效工具。利用好了&#xff0c;可以大大提高自己的工作效率&#xff0c;而且可以提高数据的准确度。“VBA语言専攻”提供的教程一共九套&#xff0c;分为初级、中级、高级三大部分&#xff0c;教程是对VBA的系统讲解&#…

K8S 部署和访问 Kubernetes 仪表板(Dashboard)

文章目录 部署 Dashboard UI浏览器访问登陆系统 Dashboard 是基于网页的 Kubernetes 用户界面。 你可以使用 Dashboard 将容器应用部署到 Kubernetes 集群中&#xff0c;也可以对容器应用排错&#xff0c;还能管理集群资源。 你可以使用 Dashboard 获取运行在集群中的应用的概览…

unit4.web服务的部署及高级优化方案

搭建web服务器要求如下&#xff1a; 1.web服务器的主机ip&#xff1a;172.25.254.100 [rootserver101 桌面]# vmset.sh 100 连接已成功激活&#xff08;D-Bus 活动路径&#xff1a;/org/freedesktop/NetworkManager/ActiveConnection/3&#xff09; [rootserver101 桌面]# ifc…