第四节线性回归——让机器从数据中重新发现胡克定律本节摘要上一节我们使用NumPy和Pandas建立了应力、应变与位移数据集。本节将第一次真正训练机器学习模型以带有少量测量误差的应力—应变数据为例通过线性回归自动识别材料弹性模量并使用测试集、MAE、RMSE和决定系数评估预测效果。学习者将理解“训练模型”并不是让计算机记住答案而是让它从样本中寻找可以解释和预测力学响应的规律。一、本节学习目标完成本节后你将能够理解特征、标签、训练集和测试集掌握线性回归的基本原理使用scikit-learn训练力学预测模型从回归系数中识别材料弹性模量使用模型预测新的应力结果使用MAE、RMSE和评估模型判断模型结果是否满足基本物理规律。二、从胡克定律到机器学习在线弹性阶段应力与应变满足胡克定律其中为应力为应变为弹性模量。如果我们已经知道就可以直接计算应力。但在真实实验中材料参数可能未知测量结果也往往带有误差。现在反过来思考如果只给计算机一组应变和应力数据它能不能自己推断出弹性模量这就是一个典型的监督学习问题线性回归模型写成其中为模型预测的应力为回归系数为截距。将它与胡克定律比较可以发现理想情况下因此本节训练的不只是一个数学模型还可以从模型参数中识别材料的物理性质。三、机器学习的基本流程线性回归不是简单地在数据上画一条直线而是一个完整的建模过程。图4-1 线性回归学习流程整个过程可以概括为四、准备应力—应变数据假设某钢材的真实弹性模量为如果数据完全来自理论公式所有样本都会严格落在一条直线上。但真实实验和有限元后处理数据可能受到以下因素影响测量误差网格离散误差数值舍入误差边界条件扰动材料本身的离散性。因此我们在理论应力中加入少量扰动使案例更接近实际数据。import numpy as np import pandas as pd strain np.arange(0.0, 0.0021, 0.0002) elastic_modulus_true 210000.0 noise np.array([ -1.5, 2.0, -2.5, 1.0, 3.5, -2.0, 2.5, -1.0, 1.5, -3.0, 2.0 ]) stress elastic_modulus_true * strain noise data pd.DataFrame({ 应变: strain, 应力_MPa: stress }) print(data.head().to_string(indexFalse)) print(\n样本数量, len(data))预期输出应变 应力_MPa 0.0000 -1.5 0.0002 44.0 0.0004 81.5 0.0006 127.0 0.0008 171.5 样本数量 11当应变为0时数据中的应力为。这并不表示材料在零应变时真的承受压应力而是模拟实验中的零点漂移或测量误差。五、认识特征和标签在监督学习中数据通常分成两部分。1. 特征特征是模型的输入。本案例的特征是应变2. 标签标签是希望模型预测的结果。本案例的标签是应力对应代码为X data[[应变]] y data[应力_MPa] print(X的形状, X.shape) print(y的形状, y.shape)预期输出X的形状 (11, 1) y的形状 (11,)为什么X要写成data[[应变]]而不是data[应变]因为scikit-learn要求输入特征通常是二维结构其中为样本数量为特征数量。本案例有11个样本和1个特征所以的形状为六、划分训练集和测试集不能使用全部数据训练模型后再用同一批数据证明模型很好。因为这样只能说明模型拟合了已经见过的数据。我们需要把数据分成训练集用于寻找模型参数测试集用于检验模型对未见数据的预测能力。train_indices [0, 2, 3, 5, 6, 8, 9] test_indices [1, 4, 7, 10] X_train X.iloc[train_indices] y_train y.iloc[train_indices] X_test X.iloc[test_indices] y_test y.iloc[test_indices] print(训练样本数量, len(X_train)) print(测试样本数量, len(X_test)) print(测试集应变, X_test[应变].to_numpy())预期输出训练样本数量 7 测试样本数量 4 测试集应变 [0.0002 0.0008 0.0014 0.002 ]这里手动指定训练集和测试集是为了让初学者清楚看到哪些数据参与训练、哪些数据只用于测试。实际项目中通常使用from sklearn.model_selection import train_test_split进行自动划分。七、训练第一个线性回归模型1. 导入模型from sklearn.linear_model import LinearRegression model LinearRegression() model.fit(X_train, y_train)fit()表示让模型根据训练数据寻找最合适的回归系数和截距。2. 查看模型参数coefficient model.coef_[0] intercept model.intercept_ print(f回归系数{coefficient:.2f} MPa) print(f回归截距{intercept:.4f} MPa)预期输出回归系数210540.54 MPa 回归截距-1.0811 MPa因此机器学习模型得到的关系为八、模型参数有什么物理意义模型识别出的回归系数为材料的真实弹性模量为二者相对误差为可以用代码验证relative_error abs( coefficient - elastic_modulus_true ) / elastic_modulus_true * 100 print(f弹性模量识别误差{relative_error:.3f}%)预期输出弹性模量识别误差0.257%这说明模型从带有扰动的数据中较准确地识别出了材料弹性模量。截距为理想的胡克定律通过坐标原点因此截距理论上应为0。较小的非零截距可以被理解为测量零点偏差或数据噪声造成的系统误差。九、线性回归究竟在优化什么对于第个样本模型预测值为真实值与预测值之间的差称为残差线性回归通常通过最小化均方误差寻找最优参数也就是说模型在寻找一条直线使所有样本点到直线的误差平方和尽可能小。之所以对误差进行平方是因为正误差和负误差不会相互抵消较大的误差会受到更强的惩罚目标函数可以方便地进行数学优化。十、绘制训练数据与回归直线import matplotlib.pyplot as plt plt.rcParams[font.sans-serif] [ Microsoft YaHei, SimHei ] plt.rcParams[axes.unicode_minus] False strain_line np.linspace( strain.min(), strain.max(), 200 ).reshape(-1, 1) stress_line model.predict(strain_line) plt.scatter( X_train[应变], y_train, color#1769aa, s70, label训练样本 ) plt.scatter( X_test[应变], y_test, color#00a7a0, markerD, s70, label测试样本 ) plt.plot( strain_line, stress_line, color#f28e2b, linewidth2.8, label回归直线 ) plt.xlabel(应变 ε) plt.ylabel(应力 σMPa) plt.title(应力—应变数据与线性回归结果) plt.grid(alpha0.3) plt.legend() plt.tight_layout() plt.show()预期绘图结果蓝色圆点表示训练样本绿色菱形表示模型没有参与训练的测试样本橙色直线表示模型学习到的应力—应变关系数据点应分布在回归直线附近。图4-2 应力—应变回归结果十一、使用模型预测测试集y_pred model.predict(X_test) comparison pd.DataFrame({ 应变: X_test[应变].to_numpy(), 真实应力_MPa: y_test.to_numpy(), 预测应力_MPa: y_pred, 残差_MPa: y_test.to_numpy() - y_pred }) print( comparison.to_string( indexFalse, float_formatlambda value: f{value:.4f} ) )预期输出应变 真实应力_MPa 预测应力_MPa 残差_MPa 0.0002 44.0000 41.0270 2.9730 0.0008 171.5000 167.3514 4.1486 0.0014 293.0000 293.6757 -0.6757 0.0020 422.0000 420.0000 2.0000残差定义为因此残差为正表示模型预测偏小残差为负表示模型预测偏大残差越接近0预测越准确。十二、评估模型预测效果1. 平均绝对误差平均绝对误差为MAE表示预测值平均偏离真实值多少。2. 均方根误差RMSE对大误差更加敏感。3. 决定系数通常越接近1模型解释能力越强接近0模型与直接使用平均值相差不大模型可能比直接使用平均值更差。4. 使用代码计算指标from sklearn.metrics import ( mean_absolute_error, mean_squared_error, r2_score ) mae mean_absolute_error(y_test, y_pred) rmse np.sqrt( mean_squared_error(y_test, y_pred) ) r2 r2_score(y_test, y_pred) print(fMAE{mae:.2f} MPa) print(fRMSE{rmse:.2f} MPa) print(fR²{r2:.4f})预期输出MAE2.45 MPa RMSE2.76 MPa R²0.9996这些结果表明测试集平均预测误差约为较大的误差也得到较好控制模型解释了约99.96%的应力变化。但是不能只因为很高就立即认为模型绝对可靠。还需要结合物理规律、数据范围和误差分布进行判断。十三、绘制预测值与真实值对比图plt.scatter( y_test, y_pred, color#1769aa, s85 ) minimum min(y_test.min(), y_pred.min()) maximum max(y_test.max(), y_pred.max()) plt.plot( [minimum, maximum], [minimum, maximum], --, color#f28e2b, linewidth2.5, label理想预测 yx ) plt.xlabel(真实应力MPa) plt.ylabel(预测应力MPa) plt.title(测试集预测值与真实值对比) plt.grid(alpha0.3) plt.legend() plt.tight_layout() plt.show()预期绘图结果如果模型预测完全准确所有数据点都会落在这条直线上。数据点距离虚线越近模型预测越准确。图4-3 预测应力与真实应力对比十四、预测一个新的应力假设现在输入一个训练数据中没有出现过的应变使用模型进行预测new_strain pd.DataFrame({ 应变: [0.0011] }) predicted_stress model.predict(new_strain)[0] print(f预测应力{predicted_stress:.2f} MPa)预期输出预测应力230.51 MPa使用真实弹性模量计算理论应力机器学习预测结果为二者十分接近。十五、将预测结果转换为结构载荷如果杆件截面积为由可以计算对应载荷area 100.0 predicted_force predicted_stress * area print(f预测载荷{predicted_force:.2f} N) print(f预测载荷{predicted_force / 1000:.2f} kN)预期输出预测载荷23051.35 N 预测载荷23.05 kN这样我们就把一个机器学习预测结果重新转化为了工程上可以理解的载荷。十六、模型为什么不能随意外推当前训练数据的应变范围为如果直接预测模型仍然会输出但这并不代表钢材真的可以在线弹性状态下承受如此高的应力。原因是材料可能已经发生屈服塑性变形损伤颈缩断裂。因此机器学习模型给出数值答案不等于这个答案满足物理规律。必须区分插值在训练数据范围内预测外推超出训练数据范围预测。通常外推的风险远高于插值。十七、是否应该强制截距为0胡克定律理论上满足因此可以设置physical_model LinearRegression( fit_interceptFalse ) physical_model.fit(X_train, y_train) print( f无截距模型系数 f{physical_model.coef_[0]:.2f} MPa )这样会强制回归直线通过原点。但是是否强制截距为0需要根据数据来源决定理论数据或高质量有限元数据可以考虑强制通过原点实验数据存在零点漂移时应保留截距截距明显偏大时应检查传感器、预载荷或数据处理流程。这体现了物理知识与数据建模的结合。十八、常见误区误区一高就说明模型一定正确只能描述数据拟合程度不能保证单位正确数据没有泄漏外推可靠模型满足物理规律数据覆盖了真实工程范围。误区二训练误差越小越好如果只关心训练误差模型可能记住已有样本却无法预测新工况。这就是过拟合。误区三相关性等于因果关系模型发现应变与应力存在关系不代表所有相关变量都具有直接物理因果关系。误区四可以忽略单位如果应力同时包含Pa和MPa模型可能仍然能够运行但结果将失去工程意义。误区五任何力学关系都能直接使用线性回归例如杆件位移它对和是线性的但对和并不是简单线性关系。这时可以构造具有物理意义的新特征再建立这称为特征工程。十九、本节练习基础练习把噪声全部改为0noise np.zeros_like(strain)重新训练模型观察回归系数是否等于截距是否接近0是否等于1。进阶练习将噪声幅度扩大5倍noise noise * 5重新计算MAE、RMSE和观察模型性能如何变化。思考练习如果训练数据只覆盖却要求模型预测这个预测属于插值还是外推它可能存在哪些风险工程挑战将第三节生成的三种材料数据分别建立线性回归模型观察回归系数是否接近各材料的弹性模量二十、本节小结本节完成了第一个真正的力学机器学习模型我们不仅得到了预测结果还从回归系数中识别出了材料弹性模量这说明机器学习与经典力学并不是相互替代的关系。力学知识可以帮助我们选择正确的输入和输出解释模型参数判断结果是否合理发现数据和模型中的异常限制不可靠的外推。下一节将学习岭回归、Lasso回归和弹性网络进一步解决多参数力学数据中的特征共线性、过拟合和关键参数筛选问题。