1、数学建模作业8多元分析实验一、 基本实验1.回归分析为估计山上积雪融化后对下游灌溉的影响，在山上建立一个观测站，测量最大积雪深度 X(米)与当年灌溉面积 Y(公顷) ，测得连续 10 年的数据如表 8.1 所示。表 8.1 10 年中最大积雪深度与当年灌溉面积的数据X Y X Y1 5.1 1907 6 7.8 30002 3.5 1287 7 4.5 19473 7.1 2700 8 5.6 22734 6.2 2373 9 8.0 31135 8.8 3260 10 6.4 2493（1）建立一元线性回归模型，求解，并验证系数、方程或相关系数是否通过检验；（2）观测得今年的数据是 X=7
2、米，给出今年灌溉面积的预测值、预测区间和置信区间（=0.05） ；（3）将数据散点、回归预测值、回归的预测区间和置信区间均画在一张图上，分析线性回归的拟合情况。解：（1）利用 R 软件中的 lm()求回归参数并作相应的检验。写出相应的 R 程序，如下：x|t|) (Intercept) 140.95 125.11 1.127 0.293 x 364.18 19.26 18.908 6.33e-08 *-Signif. codes: 0 * 0.001 * 0.01 * 0.05 . 0.1 1Residual standard error: 96.42 on 8 degrees of free
3、domMultiple R-squared: 0.9781, Adjusted R-squared: 0.9754 F-statistic: 357.5 on 1 and 8 DF, p-value: 6.33e-08系数检验：常数项截距 Intercept 的 p 值0.05,接受 0 假设，说明不显著。x 的系数的 p 值 predict(lm.sol, new, interval=“confidence“)fit lwr upr1 2690.227 2613.35 2767.105有分析结果可知：今年灌溉渠的预测值为 2690.227 公顷；预测区间是2613.35, 2767.105
4、(单位：公顷)；置信区间（=0.05）为2454.971, 2925.484 (单位：公顷)。（3）参考 exam0806.R,编写 R 程序：x|t|) (Intercept) 64.97 118.67 0.547 0.601 x 373.75 17.40 21.485 1.19e-07 *-Signif. codes: 0 * 0.001 * 0.01 * 0.05 . 0.1 1Residual standard error: 71.08 on 7 degrees of freedomMultiple R-squared: 0.9851, Adjusted R-squared: 0.98
5、29 F-statistic: 461.6 on 1 and 7 DF, p-value: 1.192e-07与修正之前相比，R 2 与调整的 R2 有所提高。回归直线为： 64.93.5YX检验：修正后是否有效，看一下回归诊断的结果，编写 R 程序：op |t|) (Intercept) 43.65007 18.05442 2.418 0.02984 * X1 1.78534 0.53977 3.308 0.00518 *X2 -0.08329 0.42037 -0.198 0.84579 X3 0.16102 0.11158 1.443 0.17098 -Signif. codes: 0
6、* 0.001 * 0.01 * 0.05 . 0.1 1Residual standard error: 19.97 on 14 degrees of freedomMultiple R-squared: 0.5493, Adjusted R-squared: 0.4527 F-statistic: 5.688 on 3 and 14 DF, p-value: 0.009227由运行结果可知：常数项一般相关，X1 系数高度相关，但是 X2,X3 系数的 P 值0.05,不相关，没有通过检验。方程 F 检验 P 值 5583.7 111.27- X3 1 830.6 6414.4 111.77
7、- X1 1 4363.4 9947.2 119.66Step: AIC=109.32Y X1 + X3Df Sum of Sq RSS AIC5599.4 109.32- X3 1 833.2 6432.6 109.82- X1 1 5169.5 10768.9 119.09从程序的运行结果可见，当用全部变量作回归方程时，AIC 值为 111.27，接下来的数据表显示，如果去掉 X2,则相应的 ACI 值为 109.32，此时 ACI 值达到最小，因此 R 软件自动去掉变量 X2,进行下一轮计算。 R 软件认为此时得到“最优”回归方程。下面我们分析一下回归结果，用函数 summary()提取
8、相关信息，看看是否是最优的回归方程:summary(lm.step)Call:lm(formula = Y X1 + X3, data = cement)Residuals:Min 1Q Median 3Q Max -29.713 -11.324 -2.953 11.286 48.679 Coefficients:Estimate Std. Error t value Pr(|t|) (Intercept) 41.4794 13.8834 2.988 0.00920 *X1 1.7374 0.4669 3.721 0.00205 *X3 0.1548 0.1036 1.494 0.15592
9、-Signif. codes: 0 * 0.001 * 0.01 * 0.05 . 0.1 1Residual standard error: 19.32 on 15 degrees of freedomMultiple R-squared: 0.5481, Adjusted R-squared: 0.4878 F-statistic: 9.095 on 2 and 15 DF, p-value: 0.002589有结果可见，回归系数检验的显著性水平有提高，但是 X3 系数的检验水平仍不理想。下面用 drop1 来做逐步回归：drop1(lm.step)Single term deletion
10、sModel:Y X1 + X3Df Sum of Sq RSS AIC5599.4 109.32X1 1 5169.5 10768.9 119.09X3 1 833.2 6432.6 109.82从运算结果可见如果去掉 X3,则 AIC 值从 109,32 增加到 109.82，增加的最少。去掉 X3,残差的平方和上升的也最少，拟合越好的方程，残差的平方和应越小，从这两项指标，应该再去掉变量 X3.lm.opt|t|) (Intercept) 59.2590 7.4200 7.986 5.67e-07 *X1 1.8434 0.4789 3.849 0.00142 * -Signif. co
11、des: 0 * 0.001 * 0.01 * 0.05 . 0.1 1Residual standard error: 20.05 on 16 degrees of freedomMultiple R-squared: 0.4808, Adjusted R-squared: 0.4484 F-statistic: 14.82 on 1 and 16 DF, p-value: 0.001417最终结果还是满意的，尽管 R2 还是 0.5 左右。但是系数检验的显著性水平已经很好了。最后得到的“ 最优” 回归方程为：59.21.843YX4. 方差分析 I(单因素方差分析)24 只小鼠随机地分成
12、3 组，每组 8 只，每组喂养不同的饲料组，喂养一定时间后，测得小鼠肝中含铁量，结果如表 8.3 所示。表 8.3 不同饲料组小鼠肝脏中铁含量（单位： ug/g）饲料 1 2 3 4 5 6 7 8A 1.00 1.01 1.13 1.14 1.70 2.01 2.23 2.63B 0.96 1.23 1.54 1.96 2.94 3.68 5.59 6.96C 2.07 3.72 4.50 4.90 6.00 6.84 8.23 10.33（1）试用单因素方差分析方法分析不同饲料的小鼠肝中的铁含量是否有显著差异？（2）如果有显著差异，那些水平之间有显著差异？（3）数据是否满足正态性和方差齐性
13、的要求，如果不满足请选择合适的分析方法，并给出你的最终结论。解：（1）利用方差分析表的方法进行分析：设小白鼠喂养的饲料为因素，三种不同的饲料为三个水平，喂养饲料后的小白鼠肝脏中铁含量视为来自三个正态分布总体 （i=1,2,3）的样本观2(,)iN测值。问题归结为：。0123123:,:,H不 全 相 等参考 exam0813.R 编写 R 程序：mouseF) A 2 7 3.118 36.559 9.104 0.001422 *Residuals 21 84.329 4.016 -Signif. codes: 0 * 0.001 * 0.01 * 0.05 . 0.1 1由计算结果显示，P
14、值F) B 2 3.974 1.987 26.69 0.000164 *A 2 4.441 2.221 29.83 0.000107 *B:A 4 21.159 5.290 71.06 8.34e-07 *Residuals 9 0.670 0.074 -Signif. codes: 0 * 0.001 * 0.01 * 0.05 . 0.1 1可见在显著性水平=0.05下，反应压力（因素B ）对提高产量的效应是高度显著的；反应温度（因素A ）对提高产量的效应是高度显著的；反应压力（因素B）和反应温度 (因素A) 交互作用对提高产量的效应是显著的。（2）计算反应压力（因素B） 、反应温度（因素
15、A ）及交互效应的均值：编写R 程序：fanying.lm - lm(Y B+A+B:A, data=fanying)anova(fanying.lm)attach(fanying);tapply(Y, B, mean) # 反应压力的均值tapply(Y, A, mean) # 反应温度的均值运行结果为：A1 A2 A3 5.150000 4.533333 5.750000即第三种反应温度对提高产量最有利。tapply(Y, B, mean)B1 B2 B3 5.783333 4.666667 4.983333即第一种反应压力对提高产量最有利。matrix(tapply(Y, A:B, me
16、an), nr=3, nc=3, byrow=T, dimnames=list(levels(A),levels(B) #计算交互效应的均值运行结果为：B1 B2 B3A1 4.45 6.5 4.50A2 6.30 3.6 3.70A3 6.60 3.9 6.75交互效应的均值最大为A3B3交互作用时的6.75，即交互作用时第三种温度和第三种压力共同作用对提高产量效果最好。（3）通过（1）与（2）的计算可知反应压力、反应温度以及反应温度与反应压力交互作用对提高产量的效应都是显著的。从反应温度角度看，第三种反应温度对提高产量效果最好，均值为5.75；从反应压力角度看，地一种反应压力对提高产量效果最好，均值为5.78；从反应压力和反应温度交互作用看，第三种反应温度和第三种反应压力共同作用对提高产量的效果最好，均值为6.75.比较上述三种情况的最大均值，可知6.75最大，此种情况对提高产量效果最优。即同时选择第三种反应压力和第三种反应温度共同作用时，对生产最有利。
