Matlab多元线性回归与显著性检验:从原理到完整实践 📅 发布时间:2026/9/18 21:06:03 👁 浏览次数: 简介围绕多元线性回归及显著性检验这份Word文档提供了一套在Matlab中可直接运行的完整程序并基于研究生教材《数理统计》例4.4.1展开完整覆盖从数据读取、模型构建到显著性检验的流程。资源包共1个docx文件约50KB内容涵盖程序说明、示例数据、完整可复制代码与运行结果结构清晰便于对照学习。程序不仅完成回归方程的整体F检验还进一步对每个自变量的回归系数做t检验尤其补充了对x2、x3的检验计算精度更高并可输出拒绝或接受原假设的明确结论用户可根据需要输入显著性水平α只需替换Excel数据即可移植到其他多元回归场景扩展性与可读性都不错。目前已有246人学习/下载适合统计建模、数据分析及Matlab编程学习者尤其适用于互联网数据驱动决策场景能帮助快速实现回归模型拟合与显著性评估。1. 多元线性回归与显著性检验一份 Matlab 程序文档背后的完整交付拿到一份名为《多元线性回归及显著性检验Matlab程序.docx》的文档最常见的使用场景是课程设计、论文附录或者是业务侧需要留存一份“能跑、能复现、能看懂”的算法说明。标题本身已经把任务切成了四块多元线性回归的模型构造、显著性检验的判断流程、Matlab 的实现方式以及最终如何把结果整理成一份供他人阅读的程序文档。多数人会以为难点是拟合本身但在 Matlab 里X\b一行就能得到最小二乘解真正花时间的是之后的工作回归方程整体是否显著、每个自变量是否有统计学意义、残差是否满足假设、共线性有没有让 t 检验失真。这篇文章按我交付这类程序时的常规顺序讲从矩阵写法到regress参数详解再到两类显著性检验和残差诊断最后落到如何用表格输出一份能直接贴进 Word 的结果。适合正在做课程设计或论文回归分析的人也适合写惯了机器学习库、想回到底层矩阵运算的工程师。2. 多元线性回归的矩阵写法与 Matlab 最小二乘实现2.1 把观测写成设计矩阵含截距与不含截距的两种写法多元线性回归的标准形式是y β₀ β₁x₁ β₂x₂ ... β_p x_p ε在 Matlab 里落地第一步不是写公式而是把数据摆成矩阵。约定所有观测按行排列每个自变量占一列所有样本拼成设计矩阵 X因变量 y 是列向量长度必须等于 X 的行数。行是观测、列是变量的这个约定要贯穿始终后面所有检验统计量的计算都建立在这套维度约定上。含截距项时需要在 X 的最左侧补一列全 1对应 β₀。这一步非常容易漏漏掉它的后果是回归会被强制过原点得到的 R² 和 F 值都失去常规意义。有些人用zscore对数据做标准化之后就不再补常数项这种做法只在你明确知道要做无截距回归时才成立。使用corrcoef计算相关系数时本质上是中心化后的协方差结构和带截距回归不是一回事不要混淆。% 生成一份模拟数据x2 与 x1 部分相关便于后面演示共线性 rng(42); n 30; x1 randn(n, 1); x2 0.6 * x1 randn(n, 1) * 0.8; y 2 1.5 * x1 - 0.8 * x2 randn(n, 1) * 0.3; % 带截距的设计矩阵第一列全 1 X_design [ones(n, 1), x1, x2];这里的rng(42)是固定随机数种子保证任何人在任何机器上重跑这段代码都得到相同结果课程作业和论文复现最好都保留这一行。x2用0.6*x1 噪声构造目的是让两个自变量之间存在一定的相关性后面计算 VIF 多重共线性诊断时能直接看到效果。2.2 正规方程、左除与 regress 的取舍最小二乘的目标是让残差平方和最小对 β 求导后得到正规方程(XᵀX)β Xᵀy其解写作 β (XᵀX)⁻¹Xᵀy。这个形式在教科书上非常干净但在 Matlab 里我不建议直接写成inv(X*X) * X*y。原因是当 X 的列之间存在近似共线性时(XᵀX) 接近奇异显式求逆会把数值误差放大矩阵条件数很大时结果甚至符号都会变化。Matlab 的\运算符mldivide会走 QR 分解或 Cholesky 分解路线数值上稳定得多。% 不推荐显式求逆 beta_normal inv(X_design * X_design) * X_design * y; % 推荐左除内部走 QR 分解 beta_left X_design \ y;当 X 是方阵时A\b等价于inv(A)*b但回归场景里 X 几乎总是长方形矩阵此时\实际上是求最小二乘解和pinv(X)*y等价但性能更好。左除返回的是列向量顺序与 X 的列顺序严格对应第一个元素是截距第二个开始是 x1、x2 的系数。在交付代码时我一般会保留beta_normal作为注释对照但实际运行只走左除这样审阅代码的人既能看到数学定义又能看到工程实现。2.3 regress 输出字段拆解b、bint、r、rint、statsregress是 Statistics and Machine Learning Toolbox 里的经典函数它把最小二乘、置信区间、残差和检验统计量一次性打包。典型调用是[b, bint, r, rint, stats] regress(y, X, alpha)其中alpha默认 0.05表示 95% 置信水平。它不要求你手工拼ones(n,1)但如果你不给它这一列它就真的默认模型不含截距。[b, bint, r, rint, stats] regress(y, X_design, 0.05);每个输出都是一个可以在后续代码里直接引用的变量含义和维度如下输出维度含义bp1 行 1 列回归系数顺序与 X 列一致bintp1 行 2 列每个系数的 95% 置信区间下界与上界rn 行 1 列每个样本的残差rintn 行 2 列每个残差的置信区间stats1 行 4 列R²、F 统计量、F 检验 p 值、误差方差估计stats的第 4 个数经常被误读为某个系数的 t 统计量实际上它是均方误差 MSE即残差方差的无偏估计数值等于sum(r.^2) / (n-p-1)。后面第 3 章的单变量 t 检验需要从 MSE 出发手动推导标准误所以这一步的认知很重要。如果你只需要做预测b和设计矩阵相乘就够了如果你需要写“显著性检验”章节bint和stats才是重点。3. 显著性检验回归方程 F 检验与回归系数 t 检验3.1 F 检验先回答整个方程有没有解释力显著性检验分两个层次第一层是回归方程整体的显著性第二层是单个系数的显著性。整体检验的原假设是“所有非截距系数都等于 0”如果这个假设不能被推翻那整个回归就没有统计学意义后面看单个系数也没意义。计算路径是方差分解总平方和 SST Σ(yᵢ − ȳ)²回归平方和 SSR Σ(ŷᵢ − ȳ)²残差平方和 SSE Σ(yᵢ − ŷᵢ)²。F 统计量定义为F (SSR / p) / (SSE / (n − p − 1))其中 p 是自变量个数n − p − 1 是残差自由度。F 服从 F(p, n−p−1) 分布stats的第 2、3 个数就是这个 F 值和对应的 p 值不需要你手动再算一遍。但为了理解它我建议在自己的代码里复算一次。p size(X_design, 2) - 1; % 自变量个数 n length(y); % 样本量 alpha 0.05; SST sum((y - mean(y)).^2); SSE sum(r.^2); SSR SST - SSE; F_manual (SSR / p) / (SSE / (n - p - 1)); pF_manual 1 - fcdf(F_manual, p, n - p - 1); % 与 regress 输出对比 fprintf(手动计算 F %.4f, p %.4g\n, F_manual, pF_manual); fprintf(regress F %.4f, p %.4g\n, stats(2), stats(3));fcdf是 F 分布的累积分布函数1 - fcdf(...)就是右侧尾部概率即 p 值。stats(2)和手动计算结果应该完全一致如果不一致优先检查 X 是否漏了ones列或 y 是否取成了行向量。p 值小于 0.05 时拒绝原假设结论是“回归方程总体显著”但这只说明自变量组合起来对 y 有解释力不代表每个变量都有效。3.2 t 检验与 bint 置信区间判断单个变量是否显著回归方程整体显著之后还要逐个检查每个自变量。第 j 个系数的 t 检验原假设是 βⱼ 0统计量为tⱼ bⱼ / se(bⱼ)标准误 se(bⱼ) 的计算需要用到系数协方差矩阵Cov(β̂) MSE · (XᵀX)⁻¹对角线开方就是各系数的标准误。MSE stats(4); % 均方误差 CovBeta MSE * inv(X_design * X_design); se sqrt(diag(CovBeta)); t_stat b ./ se; df n - p - 1; p_t 2 * (1 - tcdf(abs(t_stat), df)); % 与 bint 的结论对照 fprintf(x1: t %.3f, p %.4g, bint [%.3f, %.3f]\n, ... t_stat(2), p_t(2), bint(2,1), bint(2,2));tcdf是 t 分布的累积分布函数乘以 2 是因为双侧检验。判定规则有两条等价路径看 p 值是否小于 0.05或者看bint的行是否包含 0。bint(2,1)和bint(2,2)是 x1 系数的置信下界和上界如果区间不跨越 0说明 x1 显著。我每次写报告都会把这两种判据同时列上因为审阅者里有习惯看置信区间的人也有只看 p 值的人。3.3 查临界值还是看 p 值Matlab 里的 4 个判断代码有些课程要求使用临界值法这需要查 F 分布表和 t 分布表。Matlab 里可以用finv和tinv直接算临界值比查纸质表更准F_crit finv(1 - alpha, p, n - p - 1); % 单侧上尾 t_crit tinv(1 - alpha / 2, df); % 双侧 fprintf(F 临界值 %.4f, t 临界值 %.4f\n, F_crit, t_crit);整个显著性检验流程可以收敛成 4 个判断步骤先看stats(3)是否小于 0.05决定方程整体是否显著再看stats(1)即 R²判断解释力度是否达到业务期待然后逐个看p_t或bint剔除不显著的自变量最后保留的变量书写回归方程。顺序不能反很多人一上来就盯 R²R² 很高但 F 检验 p 值大这种情况通常发生在样本量极小的时候结论是宁可别用这个模型。4. 用 regress 与 fitlm 跑通一个销售回归案例4.1 构造数据与设计矩阵造一个 x4 干扰变量便于观察不显著把前面的片段拼成一个完整案例。模拟一个销售场景y 是月销售额x1 是电视广告投入x2 是线上广告投入x3 是促销活动次数另外加入 x4 作为与 y 完全无关的干扰变量。这样跑完显著性检验后x1、x2 显著x4 不显著结果正好符合教学预期。% 完整回归案例构造 50 条样本 rng(42); n 50; x1 randn(n, 1); x2 0.5 * x1 randn(n, 1); x3 randn(n, 1); x4 randn(n, 1); % 干扰变量理论上不显著 y 5 2.0 * x1 - 1.2 * x2 0.6 * x3 0.5 * randn(n, 1); % 构造设计矩阵 X [ones(n, 1), x1, x2, x3, x4]; varnames {截距, x1_电视, x2_线上, x3_促销, x4_干扰}; % regress 建模 [b, bint, ~, ~, stats] regress(y, X, 0.05);在交付的 .docx 里这段代码通常就是“数据说明”章节的起点。注意x4是纯随机数它与 y 之间只可能存在抽样误差层面的相关性所以它的系数 p 值大概率大于 0.05。如果你重跑这段代码发现 x4 的 p 值小于 0.05不要紧张0.05 的显著性水平本身就允许 5% 的假阳性率把样本量改到 100 会发现结果稳定。4.2 两种建模路线regress 适合教学、fitlm 适合预测regress的所有输出都是裸矩阵教学文档里容易逐项解释但做预测和模型诊断时fitlm更方便它返回一个LinearModel对象内部已经算好系数表、p 值、调整后 R²、残差属性并且自动处理截距项。% fitlm 路线第一参数是自变量矩阵第二参数是 y mdl fitlm([x1, x2, x3, x4], y); % 查看系数表Estimate / SE / tStat / pValue disp(mdl.Coefficients); % 查看整体指标R²、调整 R²、F 检验 p 值 fprintf(R2 %.4f, adjR2 %.4f, F p-value %.4g\n, ... mdl.Rsquared.Ordinary, mdl.Rsquared.Adjusted, mdl.Coefficients.pValue(1));mdl.Coefficients是一张 table行名分别是截距和各个变量列包含 Estimate、SE、tStat、pValue这个结构可以直接用writetable写出去正好满足 docx 交付需求。两者的边界是regress让你清楚每一步在算什么适合写原理章节fitlm让你少写很多代码适合最终交付程序和做预测。regress与fitlm的核心差异如下维度regressfitlm截距处理手动加ones列自动加截距输出形式矩阵与向量LinearModel 对象系数显著性需手动推导标准误系数表直接给出 p 值预测函数X_new * b手动乘predict(mdl, X_new)适用场景教学、原理说明建模、预测、残差诊断4.3 预测与区间预测新样本必须拼接截距列回归模型落地到业务里的最重要动作就是预测。如果使用regress结果新样本矩阵必须与训练时的 X 结构完全一致第一列是全 1后面列与训练变量顺序一致数量一致否则矩阵乘法直接报错。% 新样本电视广告 0.8线上广告 -0.3促销 1.2干扰变量 0 X_new [1, 0.8, -0.3, 1.2, 0]; y_hat X_new * b; % 使用 fitlm 的 predict 可以同时拿到预测值与 95% 置信区间 [y_hat2, y_ci] predict(mdl, X_new); fprintf(点预测 %.3f, 预测区间 [%.3f, %.3f]\n, y_hat2, y_ci(1), y_ci(2));注意predict的第二个输出是置信区间不是预测区间如果文档里要写“预测区间”应使用predict(mdl, X_new, Prediction, observation)。这两个区间在交付文档时容易混淆置信区间描述的是均值的不确定性预测区间包含个体噪声区间更宽。我在 docx 里通常同时给出点预测和预测区间并标注用的是哪种区间避免业务人员把均值置信区间当成单次预测范围来用。5. 显著性检验通过还不够残差诊断与多重共线性5.1 rcoplot 交互式残差图与异常点处理F 检验和 t 检验都依赖一个前提残差 ε 是独立同分布的正态随机变量均值为 0方差恒定。如果这个前提严重不满足前面所有 p 值都不可信。regress输出的rint残差置信区间就派上了用场。% 画残差区间图超出零线的点标记为红色 figure; rcoplot(r, rint);rcoplot画出每个样本的残差及置信区间。置信区间不跨越 0 的样本在图上表现为红色它们就是可疑的异常点。在图上手动点击红色点可以把它标成异常并返回交互结束后这些点会被标记但rcoplot并不会自动从数据里剔除它需要你记录下索引后自己过滤。我一般只把异常点作为“数据审查线索”而不会直接静默删除因为删点必须记录理由否则论文审稿人或业务审计会质疑结果稳定性。5.2 三行代码检查正态性与同方差性残差正态性最直接的检查是normplot它把残差分位数和理论正态分位数画在一起点大致落在直线上就说明正态性成立。同方差性则看残差与拟合值的散点图理想情况是点在 0 附近水平带状分布不呈喇叭形。% 检查残差正态性 figure; normplot(r); % 检查等方差性残差 vs 拟合值 figure; plot(mdl.Fitted, r, .); xlabel(拟合值); ylabel(残差); refline(0, 0);normplot末端如果出现明显的 S 形弯曲说明残差尾部偏厚这时 p 值的可信度下降可以考虑对 y 做对数变换或 Box-Cox 变换。散点图如果左边紧右边散开呈“漏斗形”说明存在异方差t 检验标准误会偏低导致 p 值虚小变量可能被错误判断为显著。5.3 VIF 大于 10 时 t 检验会失效多重共线性是一个容易被忽略但杀伤力很大的问题。它的典型症状是回归方程 F 检验高度显著R² 也不低但每个系数的 t 检验 p 值都很大出现“整体显著、个体全不显著”的矛盾。原因是自变量之间相关性高导致系数标准误膨胀。方差膨胀因子 VIF 的计算只需要一行corrcoef加一行求逆% 计算三个真实自变量之间的 VIF X_vars [x1, x2, x3]; R corrcoef(X_vars); VIF diag(inv(R)); % 输出结果 for i 1:length(VIF) fprintf(变量 %d 的 VIF %.2f\n, i, VIF(i)); endcorrcoef计算相关系数矩阵inv(R)的对角线元素就是每个变量的 VIF。经验阈值是VIF 大于 10 认为存在严重共线性大于 5 需要警惕。案例中的 x1 和 x2 我故意设计了 0.5 的相关性VIF 大概率在 2~4 之间属于可接受范围如果把 x2 的构造改为0.95*x1 噪声VIF 会迅速飙过 10。处理策略优先选择删变量或做岭回归而不是只报告 p 值。6. 把回归结果整理成 docx 交付表格格式化输出与选型边界6.1 用 table 与 writetable 生成可粘贴的结果表程序交付的最终产物不只是.m文件还要有一份能让不跑代码的人看懂的结果表。Matlab 的table可以组装回归系数、置信区间和 p 值然后writetable输出成 Excel再从 Excel 粘贴到 Word。T table(b, bint(:,1), bint(:,2), p_t, ... VariableNames, {系数, 置信下限, 置信上限, t检验p值}); T.Properties.RowNames varnames; disp(T); writetable(T, regression_results.xlsx, WriteRowNames, true); % p 值格式化避免出现 0.0000 这种让人误以为“绝对不显著”的写法 for i 1:length(p_t) fprintf(%s: b %.4f, p %.4g\n, varnames{i}, b(i), p_t(i)); end%.4g格式会把 1e-6 以下的小数显示成科学计数法比保留四位小数更严谨。如果 docx 里要求纯文本格式直接disp(T)后全选复制到 Word再套用表格样式即可。6.2 线性表达式边界什么时候换逐步回归、什么时候交给非线性拟合显著性和 R² 都达标的线性模型只说明“在线性框架内可用”。当业务数据存在明显指数衰减期或饱和增长期时线性假设本身就不成立此时继续做线性的显著性检验没有意义。常见做法是先用stepwiselm处理高维共线性逐步加入和剔除变量减少人工筛选成本如果 y 与 x 的关系呈 S 型曲线可以用fitnlm拟合非线性最小二乘模型或者交给 BP 神经网络这类拟合曲线能力更强的工具但代价是回归系数的可解释性和显著性检验体系会丢失。线性回归的显著性检验讲到底回答的是“这个线性关系是否可靠”不是“这个模型是否最优”。交付程序时在注释里写明这一点可以省下后续大量沟通成本。本文还有配套的精品资源点击获取