第四节:线性回归——让机器从数据中重新发现胡克定律
本节摘要
上一节我们使用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(index=False)) 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", s=70, label="训练样本" ) plt.scatter( X_test["应变"], y_test, color="#00a7a0", marker="D", s=70, label="测试样本" ) plt.plot( strain_line, stress_line, color="#f28e2b", linewidth=2.8, label="回归直线" ) plt.xlabel("应变 ε") plt.ylabel("应力 σ(MPa)") plt.title("应力—应变数据与线性回归结果") plt.grid(alpha=0.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( index=False, float_format=lambda 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(f"MAE:{mae:.2f} MPa") print(f"RMSE:{rmse:.2f} MPa") print(f"R²:{r2:.4f}")预期输出
MAE:2.45 MPa RMSE:2.76 MPa R²:0.9996这些结果表明:
测试集平均预测误差约为
;
较大的误差也得到较好控制;
模型解释了约99.96%的应力变化。
但是,不能只因为很高,就立即认为模型绝对可靠。还需要结合物理规律、数据范围和误差分布进行判断。
十三、绘制预测值与真实值对比图
plt.scatter( y_test, y_pred, color="#1769aa", s=85 ) minimum = min(y_test.min(), y_pred.min()) maximum = max(y_test.max(), y_pred.max()) plt.plot( [minimum, maximum], [minimum, maximum], "--", color="#f28e2b", linewidth=2.5, label="理想预测 y=x" ) plt.xlabel("真实应力(MPa)") plt.ylabel("预测应力(MPa)") plt.title("测试集预测值与真实值对比") plt.grid(alpha=0.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_intercept=False ) physical_model.fit(X_train, y_train) print( f"无截距模型系数:" f"{physical_model.coef_[0]:.2f} MPa" )这样会强制回归直线通过原点。
但是,是否强制截距为0需要根据数据来源决定:
理论数据或高质量有限元数据,可以考虑强制通过原点;
实验数据存在零点漂移时,应保留截距;
截距明显偏大时,应检查传感器、预载荷或数据处理流程。
这体现了物理知识与数据建模的结合。
十八、常见误区
误区一:
高就说明模型一定正确
只能描述数据拟合程度,不能保证:
单位正确;
数据没有泄漏;
外推可靠;
模型满足物理规律;
数据覆盖了真实工程范围。
误区二:训练误差越小越好
如果只关心训练误差,模型可能记住已有样本,却无法预测新工况。这就是过拟合。
误区三:相关性等于因果关系
模型发现应变与应力存在关系,不代表所有相关变量都具有直接物理因果关系。
误区四:可以忽略单位
如果应力同时包含Pa和MPa,模型可能仍然能够运行,但结果将失去工程意义。
误区五:任何力学关系都能直接使用线性回归
例如杆件位移:
它对和
是线性的,但对
和
并不是简单线性关系。
这时可以构造具有物理意义的新特征:
再建立:
这称为特征工程。
十九、本节练习
基础练习
把噪声全部改为0:
noise = np.zeros_like(strain)重新训练模型,观察:
回归系数是否等于
;
截距是否接近0;
是否等于1。
进阶练习
将噪声幅度扩大5倍:
noise = noise * 5重新计算MAE、RMSE和,观察模型性能如何变化。
思考练习
如果训练数据只覆盖:
却要求模型预测:
这个预测属于插值还是外推?它可能存在哪些风险?
工程挑战
将第三节生成的三种材料数据分别建立线性回归模型,观察回归系数是否接近各材料的弹性模量:
二十、本节小结
本节完成了第一个真正的力学机器学习模型:
我们不仅得到了预测结果,还从回归系数中识别出了材料弹性模量:
这说明机器学习与经典力学并不是相互替代的关系。力学知识可以帮助我们:
选择正确的输入和输出;
解释模型参数;
判断结果是否合理;
发现数据和模型中的异常;
限制不可靠的外推。
下一节将学习岭回归、Lasso回归和弹性网络,进一步解决多参数力学数据中的特征共线性、过拟合和关键参数筛选问题。