setwd("d:/data")
# 读取数据
d2.4 <- read.csv("C:\\Users\\86167\\Desktop\\ex2.4.csv",header = T)
model <- lm(y ~ x1 + x2 + x3 + x4, data = d2.4)
summary(model)##
## Call:
## lm(formula = y ~ x1 + x2 + x3 + x4, data = d2.4)
##
## Residuals:
## Min 1Q Median 3Q Max
## -709.52 -271.18 -35.46 212.92 1026.99
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -1.593 155.914 -0.010 0.9919
## x1 190.588 33.514 5.687 9.08e-07 ***
## x2 1.588 2.674 0.594 0.5555
## x3 1.013 1.315 0.770 0.4453
## x4 304.716 139.869 2.179 0.0346 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 422.4 on 45 degrees of freedom
## (因为不存在，1个观察量被删除了)
## Multiple R-squared: 0.7485, Adjusted R-squared: 0.7261
## F-statistic: 33.48 on 4 and 45 DF, p-value: 5.8e-13# 进行逐步回归
model.step <- step(model, direction = "both") ## Start: AIC=609.33
## y ~ x1 + x2 + x3 + x4
##
## Df Sum of Sq RSS AIC
## - x2 1 62949 8092032 607.72
## - x3 1 105799 8134882 607.98
## <none> 8029083 609.33
## - x4 1 846846 8875929 612.34
## - x1 1 5770382 13799465 634.41
##
## Step: AIC=607.72
## y ~ x1 + x3 + x4
##
## Df Sum of Sq RSS AIC
## - x3 1 141084 8233116 606.58
## <none> 8092032 607.72
## + x2 1 62949 8029083 609.33
## - x4 1 867198 8959230 610.81
## - x1 1 14810227 22902258 657.74
##
## Step: AIC=606.58
## y ~ x1 + x4
##
## Df Sum of Sq RSS AIC
## <none> 8233116 606.58
## + x3 1 141084 8092032 607.72
## + x2 1 98234 8134882 607.98
## - x4 1 862282 9095398 609.56
## - x1 1 15298607 23531723 657.09# 输出逐步回归结果
summary(model.step) ##
## Call:
## lm(formula = y ~ x1 + x4, data = d2.4)
##
## Residuals:
## Min 1Q Median 3Q Max
## -815.95 -279.23 -59.84 245.14 993.05
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -45.58 130.66 -0.349 0.7288
## x1 207.47 22.20 9.345 2.72e-12 ***
## x4 307.22 138.47 2.219 0.0314 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 418.5 on 47 degrees of freedom
## (因为不存在，1个观察量被删除了)
## Multiple R-squared: 0.7421, Adjusted R-squared: 0.7311
## F-statistic: 67.62 on 2 and 47 DF, p-value: 1.474e-14##数据中有50个观测值和5个变量，分别为No（编号）、y（CEO年薪）、x1（在目前职位年数）、x2（前一年股票价格的变化）、x3（前一年公司销售额的变化）和x4（是否有 MBA 学位0表示无,1表示有）。 ##F - 统计量为 33.48，对应的 p 值为 5.8e -13，远小于常见的显著性水平（如 0.05）。这表明整个回归模型是高度显著的，即至少有一个自变量对CEO年薪（y）有显著的线性影响，整体模型具有统计学意义。这意味着我们所建立的模型在解释 CEO 年薪的变化方面是有效的，变量之间存在线性关系的假设是合理的。 ##常数项（Intercept）常数项的估计值为-1.593，在实际意义中，当所有自变量（x1 - x4）都为 0 时,CEO年薪的预期值（在实际情况中，这种情况可能并不具有实际意义，但在模型中是一个基准值）。其t值为-0.010，p值为0.9919，不显著，说明常数项对模型的贡献较小，在解释CEO年薪时，其单独的影响可以忽略不计。 ##在目前职位年数（x1）系数估计值为190.588，意味着在其他自变量保持不变的情况下，在目前职位年数每增加 1 年，CEO年薪预计平均增加 190.588。t 值为 5.687，p 值为 9.08e-07，在0.05的显著性水平下高度显著，说明在目前职位年数是影响CEO年薪的重要因素之一，与年薪呈正相关关系。这表明CEO在公司任职时间越长，其年薪越高，符合一般的经验认知，可能是因为随着任职时间增加，CEO 积累了更多的经验和管理能力，对公司的贡献也可能更大，从而获得更高的薪酬。 ##前一年股票价格的变化（x2）系数为1.588，表明在其他条件不变时，前一年股票价格每增加 1 单位，CEO年薪预计平均增加 1.588。但t值为 0.594，p 值为 0.5555，不显著，说明前一年股票价格的变化在该模型中对CEO年薪的单独影响不明显。可能存在其他因素掩盖了其对年薪的作用，或者股票价格变化与其他自变量之间存在复杂关系影响了其显著性。例如，股票价格变化可能需要与其他公司业绩指标结合起来才能更好地影响CEO薪酬，或者公司可能有其他薪酬政策因素使得股票价格变化对年薪的直接影响不显著。 ##前一年公司销售额的变化（x3）系数为 1.013，即前一年公司销售额每增加 1 单位，CEO年薪预计平均增加1.013。t值为0.770，p值为0.4453，不显著，说明前一年公司销售额的变化单独对CEO年薪的影响不显著。这可能是因为销售额变化可能不是直接决定CEO年薪的关键因素，或者公司在制定薪酬时考虑了更多其他综合因素，如长期战略目标、市场竞争地位等，而不仅仅是短期销售额变化。也有可能是销售额变化与其他自变量存在多重共线性关系，导致其在模型中的显著性被削弱。 ##是否有 MBA 学位（x4）系数为 304.716，意味着有 MBA 学位（x4 = 1）的 CEO 相比没有 MBA 学位（x4 = 0）的CEO，年薪预计平均高出304.716。t值为2.179，p 值为 0.0346，在 0.05的显著性水平下显著，说明是否有MBA学位是影响CEO年薪的一个重要因素，拥有 MBA学位对年薪有正向影响。这可能是因为MBA学位通常代表着更高的管理知识和技能水平，公司认为具有MBA学位的CEO能够为公司带来更多价值，从而给予更高的薪酬。
##从最终模型来看，在目前职位年数（x1）和是否有MBA学位（x4）是影响CEO年薪的两个关键因素。在目前职位年数的系数较大，说明其对CEO年薪的影响更为显著，可能是因为随着任职时间的增加，CEO积累了更多的经验、人脉和管理能力，能够为公司创造更多价值，从而获得更高的薪酬。是否有MBA学位的系数也显著，这表明MBA学位在CEO薪酬决策中具有重要作用。拥有MBA学位可能代表着更高的管理知识水平和专业素养，公司认为具有MBA学位的CEO能够更好地应对复杂的商业环境和管理挑战，因此给予更高的薪酬。 ##给定新的数据（假设为新的公司相关变量值），使用predict函数进行点预测和区间预测（如 95% 置信区间）。 ##计算残差并绘制残差图，以检查模型的基本假定是否成立。
# 计算普通残差
y.res<-residuals(model)
# 计算标准化残差
y.rst<-rstandard(model.step)
print(y.rst)## 1 2 3 4 5 6
## -0.446568862 -0.945218013 0.062004702 -0.816593881 -0.281164377 0.997634957
## 7 8 9 10 11 12
## -0.171591006 -0.016449598 0.743372084 0.189530905 2.410737483 0.446008336
## 13 14 15 16 17 18
## -0.123934084 1.909831850 -0.852775518 -1.980823416 -0.705232906 -0.295270768
## 19 20 21 22 23 24
## 0.725379919 -0.975727136 -1.672437645 -0.268752285 -0.398016503 -0.394519084
## 25 26 27 28 29 30
## -0.611646884 -0.626212592 0.152394172 -1.406352767 0.745391559 2.312642857
## 31 32 33 34 35 36
## -0.043006286 1.686007110 -0.229683682 -0.977540717 -0.004239526 -0.295270768
## 37 38 39 40 41 42
## -1.413701844 -0.799868777 -0.040006782 2.394244916 1.308598927 0.753432011
## 43 44 45 46 47 48
## -0.093586794 0.491488133 1.085409594 -0.908800697 0.629820894 -0.059490423
## 49 50
## -0.345530197 -0.815410137# 计算预测值
y.fit<-predict(model.step)
# 绘制普通残差与预测值的散点图
plot(y.res~y.fit)
# 绘制标准化残差与预测值的散点图
plot(y.rst~y.fit)
##数据点 - 3.116 以下和 3.304 以上的数据很少，所以序号为 11 的数据
2.410737483 相对比较大可能是异常点
lm.step_new<-update(model.step,log(.)~.) #对模型进行对数变换
y.rst<-rstandard(lm.step_new) #计算lm.step_new的标准化残差
y.fit<-predict(lm.step_new) #计算lm.step_new的预测值
plot(y.rst~ y.fit) #绘制以标准化残差为纵坐标,预测值为横坐标的散点图
##去掉异常点在回归
lm.exam<-lm(log(y)~x1+x2+x3+x4,data= d2.4[-c(11),])
lm.step<-step(lm.exam,direction="both") ## Start: AIC=-103.92
## log(y) ~ x1 + x2 + x3 + x4
##
## Df Sum of Sq RSS AIC
## - x3 1 0.0469 4.8388 -105.443
## - x2 1 0.1112 4.9031 -104.796
## <none> 4.7919 -103.920
## - x4 1 0.8659 5.6578 -97.780
## - x1 1 5.9462 10.7381 -66.383
##
## Step: AIC=-105.44
## log(y) ~ x1 + x2 + x4
##
## Df Sum of Sq RSS AIC
## - x2 1 0.1425 4.9812 -106.021
## <none> 4.8388 -105.443
## + x3 1 0.0469 4.7919 -103.920
## - x4 1 0.8598 5.6985 -99.429
## - x1 1 5.9044 10.7431 -68.360
##
## Step: AIC=-106.02
## log(y) ~ x1 + x4
##
## Df Sum of Sq RSS AIC
## <none> 4.9812 -106.021
## + x2 1 0.1425 4.8388 -105.443
## + x3 1 0.0781 4.9031 -104.796
## - x4 1 0.8859 5.8672 -100.000
## - x1 1 16.3562 21.3375 -36.736y.rst<-rstandard(lm.step) #计算回归模型lm.step的标准化残差
y.fit<-predict(lm.step) #计算回归模型lm.step的预测值
plot(y.rst~ y.fit) #绘制以标准化残差为纵坐标，预测值为横坐标的散点图
##残值几乎都落在[-2.2]区域
par(mfrow=c(2,2)) #在一个2×2网格中创建4个绘图区
plot(lm.step_new) #绘制模型诊断图
influence.measures(lm.step_new) ## Influence measures of
## lm(formula = log(y) ~ x1 + x4, data = d2.4) :
##
## dfb.1_ dfb.x1 dfb.x4 dffit cov.r cook.d hat inf
## 1 0.00091 -0.001445 -0.003494 -0.007957 1.101 2.16e-05 0.0313
## 2 -0.00456 0.007236 -0.034835 -0.060711 1.092 1.25e-03 0.0307
## 3 0.16182 -0.028581 -0.123398 0.187323 1.095 1.18e-02 0.0602
## 4 -0.00270 0.004284 -0.020625 -0.035947 1.098 4.40e-04 0.0307
## 5 0.06608 0.071924 -0.136086 0.160383 1.127 8.70e-03 0.0736
## 6 -0.01600 0.025416 0.005887 0.040870 1.120 5.69e-04 0.0494
## 7 -0.18611 0.068627 0.105269 -0.196328 1.104 1.30e-02 0.0670
## 8 0.11286 -0.179246 0.147692 0.223723 1.125 1.68e-02 0.0846
## 9 0.01031 -0.006491 -0.003074 0.010309 1.182 3.62e-05 0.0975
## 10 0.07868 -0.124961 -0.005982 -0.168961 1.114 9.64e-03 0.0669
## 11 -0.03529 0.056054 0.135546 0.308657 0.914 3.05e-02 0.0313
## 12 0.09371 -0.048247 -0.038973 0.094733 1.151 3.05e-03 0.0794
## 13 0.12628 0.011031 -0.130474 0.171849 1.100 9.96e-03 0.0591
## 14 0.08205 0.415308 -0.503170 0.606661 0.995 1.18e-01 0.1107
## 15 -0.02231 0.035436 -0.040831 -0.060128 1.113 1.23e-03 0.0464
## 16 0.03474 -0.055179 -0.133431 -0.303841 0.919 2.96e-02 0.0313
## 17 -0.09151 -0.007994 0.094548 -0.124530 1.116 5.25e-03 0.0591
## 18 -0.33232 0.122541 0.187970 -0.350566 1.025 4.04e-02 0.0670
## 19 0.13598 -0.215954 0.248838 0.366439 0.939 4.31e-02 0.0464
## 20 -0.01658 0.026330 -0.043964 -0.067349 1.097 1.54e-03 0.0358
## 21 -0.17900 0.284278 -0.327565 -0.482373 0.832 7.18e-02 0.0464
## 22 0.06773 0.073721 -0.139486 0.164389 1.126 9.14e-03 0.0736
## 23 0.00011 -0.000175 -0.000424 -0.000965 1.101 3.17e-07 0.0313
## 24 0.06150 0.029076 -0.087842 0.106767 1.127 3.87e-03 0.0635
## 25 0.00371 -0.005892 -0.014248 -0.032444 1.099 3.58e-04 0.0313
## 26 0.00396 -0.006295 -0.015221 -0.034661 1.098 4.09e-04 0.0313
## 27 -0.00484 0.007681 0.004948 0.017473 1.108 1.04e-04 0.0376
## 28 0.02820 -0.431926 0.416171 -0.533866 1.176 9.42e-02 0.1703
## 29 -0.01634 0.025943 0.062735 0.142857 1.057 6.86e-03 0.0313
## 30 -0.06656 0.105709 0.068094 0.240482 1.008 1.91e-02 0.0376
## 31 -0.25141 0.399277 -0.279907 -0.448799 1.156 6.69e-02 0.1453
## 32 0.06690 -0.106242 0.005959 -0.130450 1.159 5.77e-03 0.0900
## 33 0.04092 0.091210 -0.132120 0.156037 1.153 8.25e-03 0.0894
## 34 0.08809 -0.139904 -0.032403 -0.224971 1.054 1.69e-02 0.0494
## 35 -0.20625 0.327566 -0.229635 -0.368193 1.185 4.54e-02 0.1453
## 36 -0.33232 0.122541 0.187970 -0.350566 1.025 4.04e-02 0.0670
## 37 -0.04893 0.077716 -0.129764 -0.198787 1.033 1.32e-02 0.0358
## 38 0.03132 -0.049735 -0.032037 -0.113143 1.085 4.33e-03 0.0376
## 39 0.00184 -0.002916 -0.001878 -0.006634 1.108 1.50e-05 0.0376
## 40 0.12458 -0.197858 0.330369 0.506097 0.727 7.58e-02 0.0358 *
## 41 -0.02323 0.036894 0.089214 0.203153 1.014 1.37e-02 0.0313
## 42 0.02649 -0.016677 -0.007897 0.026486 1.181 2.39e-04 0.0975
## 43 0.00374 -0.005945 -0.003830 -0.013525 1.108 6.23e-05 0.0376
## 44 0.05741 -0.091178 0.152242 0.233222 1.007 1.80e-02 0.0358
## 45 0.46261 -0.291288 -0.137927 0.462615 1.041 6.99e-02 0.0975
## 46 0.26804 -0.425700 0.023878 -0.522700 0.984 8.78e-02 0.0900
## 47 0.01393 -0.022123 0.106508 0.185625 1.026 1.15e-02 0.0307
## 48 0.00253 -0.004014 -0.002585 -0.009131 1.108 2.84e-05 0.0376
## 49 0.01953 -0.031011 0.051780 0.079323 1.094 2.14e-03 0.0358
## 50 -0.45372 0.080137 0.345995 -0.525233 0.868 8.59e-02 0.0602##50为强影响点
# 假设新的数据
preds <- data.frame(x1 = 5, x2 = 0.1, x3 = 0.2, x4 = 1)
# 点预测和区间预测
predict(model.step,newdata=preds,interval="prediction", level=0.95)## fit lwr upr
## 1 1299.007 442.0937 2155.92##计算结果y的点预测为6.90,预测区间为[6.23, 7.57]