简介这份资源是广东工业大学数值分析课程第四次上机实验的完整报告面向正在学习数值分析、需要完成Python上机实验的高校学生。内容围绕超定方程组的最小二乘解与最小二乘拟合展开涵盖正规方程ATAxATb的推导与求解、线性组合模型拟合以及石油产量预测等实际应用场景并配有numpy与sklearn的完整源代码。压缩包共1个docx文件约112KB内含题目描述、算法介绍、结果分析、源代码与函数图像结构完整可直接对照复现。目前已有213人学习下载。读者可借此掌握最小二乘法的数学原理与Python实现路径理解从建立正规方程到调用线性代数求解器的完整流程并参考拟合与预测案例的建模思路适合作为实验报告模板或复习数值计算方法的参考资料。1. 从一份实验报告说起这套 Python 数值分析代码能解决什么如果你手头正好有一份数值分析上机实验的源码包打开一看是四个.py文件加一份报告别急着关掉。这套代码解决的是数值计算里最常被问到的三类问题超定方程组的最小二乘解、非线性模型转线性组合的曲线拟合、以及用多项式拟合做趋势预测。它不依赖任何私有库numpy、scipy、sklearn、matplotlib四个常见包就能跑起来对刚接触数值分析、想照着复现实验结果的开发者来说是一份能直接对照公式看代码的参考实现。整份代码围绕正规方程A^T A x A^T b展开把矩阵运算、拟合系数求解、残差评估串成了一条完整链路适合拿来练手也适合当作二次开发的起点。2. 环境与目录把四个实验脚本跑起来的最小配置2.1 依赖清单与版本建议这套代码用到的库不多但版本差异会直接影响结果。numpy负责矩阵运算和polyfit多项式拟合scipy.optimize.leastsq做非线性最小二乘sklearn.linear_model.LinearRegression在拟合题里做线性回归matplotlib负责出图。常见做法是建一个干净的虚拟环境避免和系统里已有的科学计算包打架。# 建虚拟环境并激活Windows 用 venv\Scripts\activate python -m venv na_env source na_env/bin/activate # 安装依赖numpy 建议 1.24 以上scipy 与 numpy 版本要匹配 pip install numpy scipy scikit-learn matplotlib这里有个血泪经验scipy.optimize.leastsq在新版本里仍然可用但它对初始参数param0很敏感初始值给得离谱会直接报TypeError或返回未收敛的结果。我一般先把初始参数全设成 0 或 1跑通再调。matplotlib的中文字体设置在 Windows 和 Linux 上表现不同代码里用SimHei指定字体Linux 上大概率找不到这个字体需要换成系统里实际存在的中文字体否则标题和坐标轴会显示成方块。2.2 目录结构与文件职责源码包按题目拆成独立脚本结构大致如下numerical_analysis_lab4/ ├── problem3_overdetermined.py # 超定方程组 a、b 的最小二乘解 ├── problem6_fitting.py # 1/y c0 c1/x 型拟合 ├── experiment4_oil.py # 石油产量预测直线/抛物线/立方曲线对比 └── report.md # 实验报告含公式推导与结果表三个脚本各自独立没有共享模块复制到任何目录都能单独运行。problem3和problem6是纯计算加打印experiment4_oil会弹图跑之前确认有图形界面或者把plt.show()换成plt.savefig()。报告文件里记录了每个题目的算法说明和结果分析对照代码看能省不少推导时间。2.3 运行方式与预期输出每个脚本直接python 文件名.py即可。problem3会打印两组解problem6会打印拟合系数并弹出一张散点加拟合线的图experiment4_oil会打印四种拟合的多项式系数、2010 年预测值以及残差统计最后弹出一张多曲线对比图。第一次跑建议先注释掉plt.show()确认数值输出正确再开图避免图形窗口卡住终端。提示如果experiment4_oil.py里leastsq报参数维度错误检查param0到param3的长度是否和对应拟合函数的待定系数个数一致线性是 2 个、抛物线 3 个、立方 4 个、四次 5 个少一个都会翻车。3. 核心模块拆解正规方程、拟合变换与残差评估3.1 超定方程组用正规方程把无解变有解超定方程组的方程个数多于未知数一般没有精确解最小二乘的思路是找让残差平方和最小的那组解。定理很直接x*是最小二乘解的充要条件是它满足正规方程A^T A x A^T b。代码就是把这条定理翻译成两行矩阵乘法。import numpy as np # 方程组 a4 个方程 2 个未知数 A np.array([[2, 4], [3, -5], [1, 2], [2, 1]]) b np.array([[11], [3], [6], [7]]) # 构造正规方程 A^T A x A^T b ATA np.dot(A.T, A) # 系数矩阵2x2 ATb np.dot(A.T, b) # 右端向量2x1 x np.linalg.solve(ATA, ATb) print(x) # 预期 x1≈3.0403, x2≈1.2418A.T是转置np.dot做矩阵乘法np.linalg.solve解线性方程组。这里没有用lstsq而是显式构造正规方程目的是让读者看清定理的每一步。参数上唯一要注意的是b的维度写成(4,1)的列向量和A.T相乘后才是(2,1)如果写成(4,)一维数组solve也能跑但返回形状会变后续拼接容易出错。方程组 b 是 5 个方程 3 个未知数套路完全一样把A换成[[1,-2,4],[1,-1,1],[1,0,0],[1,1,1],[1,2,4]]、b换成[[0],[1],[2],[1],[0]]即可解出来是x1≈1.6571, x20, x3≈-0.4286。3.2 非线性拟合先取倒数再套线性回归题目六要拟合y c0 c1/x这不是线性组合形式但取倒数后变成1/y c0 c1·(1/x)基函数是1和1/x就回到了线性组合模型。代码用sklearn的LinearRegression做这件事关键在特征构造。import numpy as np from sklearn.linear_model import LinearRegression # 原始观测数据 X np.array([[1, 0.33], [2, 0.40], [4, 0.44], [5, 0.45]]) X_n 1 / X # 对 x 和 y 同时取倒数 x_train X_n[:, :1] # 特征1/x y_train X_n[:, 1:2] # 目标1/y model LinearRegression() model.fit(x_train, y_train) c1 model.coef_[0][0] # 1/x 的系数 c0 model.intercept_[0] # 截距 print(1/y %.4f * 1/x %.4f % (c1, c0))X_n[:, :1]保留二维形状LinearRegression要求特征矩阵是二维的直接取一维会报错。coef_是二维数组取[0][0]才是标量系数。拟合出来c0≈2.0143, c1≈1.0093代回原式就是y 2.0143 1.0093/x。这里有个容易忽略的点取倒数会改变误差结构原来在y上的最小二乘和在1/y上的最小二乘不是一回事报告里按题目要求做变换拟合实际工程中如果对y的绝对误差更敏感应该直接用非线性最小二乘拟合原式这是二次开发时值得改的地方。3.3 石油产量预测三种多项式拟合与残差对比实验四用 1994 到 2003 年共 10 个数据点分别做直线、抛物线、立方曲线拟合再用残差指标判断哪种更合理。代码用np.polyfit求系数用leastsq做非线性拟合对比两条路径都保留了下来。import numpy as np x np.arange(0, 10, 1, dtypefloat) # 年份偏移1994 对应 0 y np.array([67.052, 68.008, 69.803, 72.024, 73.400, 72.063, 74.669, 74.487, 74.065, 76.777]) for deg, name in [(1, 直线), (2, 抛物线), (3, 立方曲线)]: coef np.polyfit(x, y, deg) # 最小二乘多项式拟合 p np.poly1d(coef) pred_2010 p(2010 - 1994) # 年份偏移后代入 residual y - p(x) print(name, 系数:, coef) print(2010 预测: %.4f % pred_2010) print(二范数平方: %.4f, 残差标准差: %.4f, 绝对残差均值: %.4f % ((residual**2).sum(), residual.std(), abs(residual).mean()))np.polyfit(x, y, deg)返回从高次到低次的系数np.poly1d把它包装成可调用对象代入2010-199416得到预测值。残差三个指标各有侧重二范数平方是残差平方和直接对应最小二乘目标标准差反映离散程度绝对残差均值对异常值没那么敏感。报告里立方曲线的二范数平方 5.7434 最小直线 9.5320 最大所以结论是立方拟合最好。但这里有个坑三次多项式在数据范围外推时增长很快2010 年预测 99.0420 明显高于抛物线的 74.4106报告里也讨论了这一点认为立方曲线更符合增长趋势。实际做预测时多项式次数不是越高越好外推区间离样本越远高次项越容易失控。3.4 残差评估三个指标怎么读残差评估这块代码把leastsq和polyfit的结果放在一起对比leastsq需要自己写目标函数返回残差向量polyfit内部已经做了最小二乘。两种方式得到的系数在数值上应该接近如果差很多多半是leastsq的初始参数没给好或者迭代次数不够。读残差指标时二范数平方适合比较同一组数据上不同模型的拟合优度标准差适合看残差分布的离散程度绝对残差均值适合快速判断整体偏差量级。三个指标同向变化时结论比较可靠如果出现一个指标好另一个指标差就要回头看数据里有没有离群点。4. 避坑与排查跑这份代码时最容易翻车的五处4.1 中文字体显示成方块现象是图标题、坐标轴标签全是方框。原因是matplotlib默认字体不含中文代码里mpl.rcParams[font.sans-serif][SimHei]只在 Windows 上大概率有效。解决方法是先查系统可用字体把SimHei换成实际存在的比如 Linux 上常用WenQuanYi Micro Hei或者干脆把标签改成英文避免字体依赖。4.2 leastsq 报参数形状错误现象是leastsq抛TypeError或返回的var[0]长度不对。原因是初始参数列表长度和拟合函数里解包的系数个数不一致比如立方拟合函数解包 4 个系数param2却只给了 3 个。解决方法是逐个核对param0到param3的长度线性 2、抛物线 3、立方 4、四次 5一个都不能少。4.3 2010 年预测值代入错误现象是预测值大得离谱或小得离谱。原因是polyfit用的是偏移后的x0 到 9预测时却直接代入了2010。解决方法是统一用偏移量2010-199416所有预测都代 16不要混用原始年份和偏移年份。4.4 取倒数后出现除零或极端值现象是拟合系数异常或报inf。原因是原始数据里x或y含 0取倒数直接爆炸。这份数据里x从 1 开始、y都大于 0暂时安全但二次开发换数据时很容易踩到。解决方法是在取倒数前加一层检查过滤或平移含 0 的样本。4.5 图形窗口阻塞终端现象是脚本跑完不退出终端卡住。原因是plt.show()会阻塞直到关闭窗口。解决方法是在批量跑或服务器环境里把plt.show()换成plt.savefig(output.png)需要看图时再单独开交互环境。5. 二次开发把实验代码改成可复用的拟合工具这份代码目前是“一题一脚本”的形态适合交作业不适合反复调用。我一般会做三件事把它变成能复用的工具。第一把超定方程组求解抽成一个函数输入A和b返回解和残差内部用np.linalg.lstsq替代显式正规方程数值稳定性更好尤其是矩阵条件数大的时候显式构造A^T A会放大误差。第二把拟合部分做成配置驱动用字典描述“变换方式 多项式次数”比如{transform: reciprocal, degree: 1}对应题目六{transform: none, degree: 3}对应立方拟合这样换数据不用改代码。第三把残差评估统一成一个函数输入真实值和预测值输出二范数平方、标准差、绝对均值三个指标方便横向比较不同模型。import numpy as np def fit_overdetermined(A, b): 超定方程组最小二乘解返回解向量和残差 x, residuals, rank, sv np.linalg.lstsq(A, b, rcondNone) return x, residuals def evaluate(y_true, y_pred): 统一残差评估三个指标一次返回 r np.asarray(y_true) - np.asarray(y_pred) return { sse: float((r ** 2).sum()), std: float(r.std()), mae: float(abs(r).mean()) }验证方法很直接拿题目三的A、b跑fit_overdetermined结果应该和原来np.linalg.solve得到的x1≈3.0403, x2≈1.2418一致拿题目六的数据跑变换加拟合系数应该落在c0≈2.0143, c1≈1.0093附近。如果偏差超过小数点后两位先查变换是否对称、特征维度是否对得上。改造时还有一个技巧把年份偏移量做成参数而不是硬编码polyfit前先减基准年预测时再加回来这样换数据集不用改预测逻辑。改造点原实现建议做法收益超定方程求解显式A^T Anp.linalg.lstsq条件数大时更稳拟合配置每题硬编码字典驱动变换与次数换数据不改代码残差评估分散在各脚本统一evaluate函数模型对比口径一致年份处理硬编码偏移基准年参数化可迁移到其他年份数据从那以后我每次拿到这种实验型代码都强制先跑一遍数值输出、再开图、最后才动结构顺序反了很容易在图形和数值之间来回折腾。希望这份拆解能帮你把这份实验代码顺利跑起来也能顺着上面的改造点改成自己项目里能用的拟合小工具。本文还有配套的精品资源点击获取