草庐IT

线性代数Python计算:线性方程组的最小二乘解

戌崂石 2023-04-08 原文

给定ℝ上无解线性方程组 A x = b \boldsymbol{Ax}=\boldsymbol{b} Ax=b,构造 A T A \boldsymbol{A}^\text{T}\boldsymbol{A} ATA A T b \boldsymbol{A}^\text{T}\boldsymbol{b} ATb,然后调用博文《线性方程组的通解》定义的mySolve函数,解方程组 A T A x = A T b \boldsymbol{A}^\text{T}\boldsymbol{Ax}=\boldsymbol{A}^\text{T}\boldsymbol{b} ATAx=ATb。取任一特解 x 0 \boldsymbol{x}_0 x0即为解线性方程组 A x = b \boldsymbol{Ax}=\boldsymbol{b} Ax=b的一个最小二乘解(即 A \boldsymbol{A} A的列向量生成空间中距 b \boldsymbol{b} b最近的向量)。
例1 用Python计算无解方程组 { 4 x 1 + 2 x 2 − x 3 = 2 3 x 1 − x 2 + 2 x 3 = 10 11 x 1 + 3 x 2 = 8 \begin{cases}4x_1+2x_2-x_3=2\\3x_1-x_2+2x_3=10\\11x_1+3x_2\quad\quad=8\end{cases} 4x1+2x2x3=23x1x2+2x3=1011x1+3x2=8的最小二乘解。

import numpy as np                              #导入numpy
from fractions import Fraction as F             #导入Fraction
np.set_printoptions(formatter=                  #设置输出数据格式
                    {'all':lambda x:str(F(x).limit_denominator())})
A=np.array([[4,2,-1],                           #设置系数矩阵A
            [3,-1,2],
            [11,3,0]],dtype='float')
b=np.array([2,10,8])                            #常数项量b
B=np.matmul(A.T,A)                              #A的转置与A的积
c=np.matmul(A.T,b.reshape(3,1))                 #A的转置与b的积
X=mySolve(B,c)                                  #解最小二乘方程组
print(X[:,0])

程序的第5~7行设置原方程组的系数矩阵A,第8行设置原方程组的常数向量b。第9行调用numpy的matmul函数计算 A T A \boldsymbol{A}^\text{T}\boldsymbol{A} ATA,存于B。第10行计算 A T A b \boldsymbol{A}^\text{T}\boldsymbol{Ab} ATAb,存于c。第11行调用函数mySolve(见博文《线性方程组的通解》)解方程组 A T A x = A T b \boldsymbol{A}^\text{T}\boldsymbol{Ax}=\boldsymbol{A}^\text{T}\boldsymbol{b} ATAx=ATb,解集记为X。注意,X中第1列(X[:,0])存储的是方程组的特解。运行程序,输出

[9/5 -18/5 0]

即原方程组 A x = b \boldsymbol{Ax}=\boldsymbol{b} Ax=b的最小二乘解为 x 0 = ( 9 5 − 18 5 0 ) \boldsymbol{x}_0=\begin{pmatrix}\frac{9}{5}\\-\frac{18}{5}\\0\end{pmatrix} x0= 595180 。它是 α 1 = ( 4 3 1 ) , α 2 = ( 2 − 1 3 ) , α 3 = ( 1 2 0 ) \boldsymbol{\alpha}_1=\begin{pmatrix}4\\3\\1\end{pmatrix},\boldsymbol{\alpha}_2=\begin{pmatrix}2\\-1\\3\end{pmatrix},\boldsymbol{\alpha}_3=\begin{pmatrix}1\\2\\0\end{pmatrix} α1= 431 ,α2= 213 ,α3= 120 的生成子空间 W W W中距 b = ( 2 10 8 ) \boldsymbol{b}=\begin{pmatrix}2\\10\\8\end{pmatrix} b= 2108 最近的向量。
写博不易,敬请支持:
如果阅读本文于您有所获,敬请点赞、评论、收藏,谢谢大家的支持!

有关线性代数Python计算:线性方程组的最小二乘解的更多相关文章

  1. python - 如何使用 Ruby 或 Python 创建一系列高音调和低音调的蜂鸣声? - 2

    关闭。这个问题是opinion-based.它目前不接受答案。想要改进这个问题?更新问题,以便editingthispost可以用事实和引用来回答它.关闭4年前。Improvethisquestion我想在固定时间创建一系列低音和高音调的哔哔声。例如:在150毫秒时发出高音调的蜂鸣声在151毫秒时发出低音调的蜂鸣声200毫秒时发出低音调的蜂鸣声250毫秒的高音调蜂鸣声有没有办法在Ruby或Python中做到这一点?我真的不在乎输出编码是什么(.wav、.mp3、.ogg等等),但我确实想创建一个输出文件。

  2. ruby-on-rails - 使用一系列等级计算字母等级 - 2

    这里是Ruby新手。完成一些练习后碰壁了。练习:计算一系列成绩的字母等级创建一个方法get_grade来接受测试分数数组。数组中的每个分数应介于0和100之间,其中100是最大分数。计算平均分并将字母等级作为字符串返回,即“A”、“B”、“C”、“D”、“E”或“F”。我一直返回错误:avg.rb:1:syntaxerror,unexpectedtLBRACK,expecting')'defget_grade([100,90,80])^avg.rb:1:syntaxerror,unexpected')',expecting$end这是我目前所拥有的。我想坚持使用下面的方法或.join,

  3. Python 相当于 Perl/Ruby ||= - 2

    这个问题在这里已经有了答案:关闭10年前。PossibleDuplicate:Pythonconditionalassignmentoperator对于这样一个简单的问题表示歉意,但是谷歌搜索||=并不是很有帮助;)Python中是否有与Ruby和Perl中的||=语句等效的语句?例如:foo="hey"foo||="what"#assignfooifit'sundefined#fooisstill"hey"bar||="yeah"#baris"yeah"另外,类似这样的东西的通用术语是什么?条件分配是我的第一个猜测,但Wikipediapage跟我想的不太一样。

  4. java - 什么相当于 ruby​​ 的 rack 或 python 的 Java wsgi? - 2

    什么是ruby​​的rack或python的Java的wsgi?还有一个路由库。 最佳答案 来自Python标准PEP333:Bycontrast,althoughJavahasjustasmanywebapplicationframeworksavailable,Java's"servlet"APImakesitpossibleforapplicationswrittenwithanyJavawebapplicationframeworktoruninanywebserverthatsupportstheservletAPI.ht

  5. 华为OD机试用Python实现 -【明明的随机数】 2023Q1A - 2

    华为OD机试题本篇题目:明明的随机数题目输入描述输出描述:示例1输入输出说明代码编写思路最近更新的博客华为od2023|什么是华为od,od薪资待遇,od机试题清单华为OD机试真题大全,用Python解华为机试题|机试宝典【华为OD机试】全流程解析+经验分享,题型分享,防作弊指南华为o

  6. python - 如何读取 MIDI 文件、更改其乐器并将其写回? - 2

    我想解析一个已经存在的.mid文件,改变它的乐器,例如从“acousticgrandpiano”到“violin”,然后将它保存回去或作为另一个.mid文件。根据我在文档中看到的内容,该乐器通过program_change或patch_change指令进行了更改,但我找不到任何在已经存在的MIDI文件中执行此操作的库.他们似乎都只支持从头开始创建的MIDI文件。 最佳答案 MIDIpackage会为您完成此操作,但具体方法取决于midi文件的原始内容。一个MIDI文件由一个或多个音轨组成,每个音轨是十六个channel中任何一个上的

  7. 「Python|Selenium|场景案例」如何定位iframe中的元素? - 2

    本文主要介绍在使用Selenium进行自动化测试或者任务时,对于使用了iframe的页面,如何定位iframe中的元素文章目录场景描述解决方案具体代码场景描述当我们在使用Selenium进行自动化测试的时候,可能会遇到一些界面或者窗体是使用HTML的iframe标签进行承载的。对于iframe中的标签,如果直接查找是无法找到的,会抛出没有找到元素的异常。比如近在咫尺的例子就是,CSDN的登录窗体就是使用的iframe,大家可以尝试通过F12开发者模式查看到的tag_name,class_name,id或者xpath来定位中的页面元素,会抛出NoSuchElementException异常。解决

  8. 计算机毕业设计ssm+vue基本微信小程序的小学生兴趣延时班预约小程序 - 2

    项目介绍随着我国经济迅速发展,人们对手机的需求越来越大,各种手机软件也都在被广泛应用,但是对于手机进行数据信息管理,对于手机的各种软件也是备受用户的喜爱小学生兴趣延时班预约小程序的设计与开发被用户普遍使用,为方便用户能够可以随时进行小学生兴趣延时班预约小程序的设计与开发的数据信息管理,特开发了小程序的设计与开发的管理系统。小学生兴趣延时班预约小程序的设计与开发的开发利用现有的成熟技术参考,以源代码为模板,分析功能调整与小学生兴趣延时班预约小程序的设计与开发的实际需求相结合,讨论了小学生兴趣延时班预约小程序的设计与开发的使用。开发环境开发说明:前端使用微信微信小程序开发工具:后端使用ssm:VU

  9. python ffmpeg 使用 pyav 转换 一组图像 到 视频 - 2

    2022/8/4更新支持加入水印水印必须包含透明图像,并且水印图像大小要等于原图像的大小pythonconvert_image_to_video.py-f30-mwatermark.pngim_dirout.mkv2022/6/21更新让命令行参数更加易用新的命令行使用方法pythonconvert_image_to_video.py-f30im_dirout.mkvFFMPEG命令行转换一组JPG图像到视频时,是将这组图像视为MJPG流。我需要转换一组PNG图像到视频,FFMPEG就不认了。pyav内置了ffmpeg库,不需要系统带有ffmpeg工具因此我使用ffmpeg的python包装p

  10. Python 刷Leetcode题库,顺带学英语单词(31) - 2

    ValidPalindromeGivenastring,determineifitisapalindrome,consideringonlyalphanumericcharactersandignoringcases. [#125]Example:"Aman,aplan,acanal:Panama"isapalindrome."raceacar"isnotapalindrome.Haveyouconsiderthatthestringmightbeempty?Thisisagoodquestiontoaskduringaninterview.Forthepurposeofthisproblem

随机推荐