571A另一个系列的文章介绍了回归分析的理论，这个系列的文章介绍R语言做回归分析的实践，但不会涉及R语言编程，只是介绍回归分析需要的命令、怎么输入以及怎么解释输出。

NBA球员的draft number与他们的工资之间是有一定的关系的。NBA Draft简单理解就是球队pick心仪的新球员，players selected number 1 overall相当于就是C位，比如2002年的C位就是姚明，他的draft number就是1。所以我们的一个直觉就是新球员的draft number越小（名次越高），他的工资就应该越高。在这个系列的博文中，我们用2017-2018 NBA draft的数据为例，来验证一下我们的直觉，这个数据我上传了的，需要的话可以找来下载。

基础回归分析

首先读取数据，先简单看一下数据在excel里面的样子

setwd("D:/Stat PhD/taking course/summer1/ref/regression")
Data <- read.csv("Salary1.csv", header = TRUE, sep = ",", quote = "\"",
         dec = ".", fill = TRUE, comment.char = "")
X <- as.numeric(Data[,2])
Y <- as.numeric(Data[,3])

自己尝试的时候记得用setwd修改工作目录，读取数据的时候如果数据存在工作目录下read.csv第一个输入的位置可以直接写filename不用写文件的路径，如果数据没有存在工作目录下的，就要把文件的路径写完整。这个文件读进来以后可能是变量类型是char，要做回归的话用as.numeric转成num就可以了。

接下来用lm命令做线性回归就可以了，lm是R语言估计线性模型的函数。第一条命令做工资关于名次的回归，把结果存在ureg01.lm中，第二行输出结果的一些总结信息。

> ureg01.lm <- lm(Y~X)
> summary(ureg01.lm)

Call:
lm(formula = Y ~ X)

Residuals:
     Min       1Q   Median       3Q      Max 
-9901919 -4761884 -1586787  2215989 27551106 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 10285438     440112   23.37   <2e-16 ***
X            -139730      11404  -12.25   <2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 6365000 on 649 degrees of freedom
Multiple R-squared:  0.1879,	Adjusted R-squared:  0.1866 
F-statistic: 150.1 on 1 and 649 DF,  p-value: < 2.2e-16

第二行以下的内容就是回归结果的总结。call后面这个是回归命令，公式是Y~X，意思就是我们执行了Y关于X的回归。residual后面的是残差的描述性统计，从左到右分别是残差的最小值、25%分位点，中位数，75%分位点以及最大值。Coefficients后面是系数的估计，一元回归只有截距项和解释变量X的系数。第一列是系数的估计值，10285438在这个问题下面的解释是如果有第0名存在，那么工资估计就是这么多，-139730表示名次每降一名，工资平均会降139730，也就是说我们之前的直觉是对的，名次越低工资越低。第二列是系数的估计量的方差，第三列是系数估计量的t统计量，第四列是系数估计量的t检验的p值，这里的t检验原假设是系数为0，因为这两个p值都非常小，所以我们可以很自信地拒绝原假设，认同名次越低工资也会越低的假设。倒数第三行是残差的标准误（标准差）以及相应的自由度，这里一共有651个样本，所以总自由度是650，回归模型占一个自由度，所以误差自由度是649。倒数第二行是R方的结果，第一个R方是多元回归的R方，就是根据残差平方和、回归平方和以及自由度调整计算出来的，一般看这个R方就可以了，这个0.1879表示这个一元线性回归模型可以解释18.79%的工资的变化（也就是说这个模型解释力其实很低，名次对于球员工资没有主要的解释力）；第二个R方是考虑到只要我们不断增加解释变量，第一个R方根据定义的话它就会不断变大，但这种变大没有意义，因为模型可能是过拟合的状态，所以第二个R方相对于第一个会把模型的复杂度考虑进去，在第一个R方的基础上，模型越复杂，第二个R方就会越小。最后一行是对模型整体的检验，它的原假设是截距与X的系数都为0，这里F统计量是150.1，自由度是1和649，p值非常小，说明我们可以拒绝原假设，认为并非截距和X的系数都是0，也就是说这个模型还是有意义的。

现在我们有了第一个解释球员工资的模型，它可以用回归方程表示出来：
$\hat{Y} = 10285438-139730X$
用这个模型我们可以做一些简单的拟合与预测。如果某位球员名次是43名，根据这个模型我们可以估计他的工资应该是

> predict(ureg01.lm,newdata = data.frame(X=43),interval = "conf",level = 0.95)
      fit     lwr     upr
1 4277046 3726681 4827410

用来做拟合和预测的都是predict函数，输入第一项是模型对象，我们这个问题的模型对象就是之前估计得到的ureg01.lm，输入第二项newdata = 后面要接的类型数据框，需要用data.frame做一个转换，第三项是选择区间估计的类型，在做predict的时候，区间估计类型有两种，拟合和预测，拟合用conf表示，预测用predict表示，同样的置信水平下拟合的置信区间更短，因为预测会把新样本也看成是随机变量，会让Y的估计值方差更大。第四项是置信水平。输出第一列是fit，也就是拟合值，如果某球员名次是43，那么根据模型估计的工资就是4277046，置信区间是[3726681,4827410]，也就是说我们有95%的把握他的工资会在这个区间内。

再多提一下这个数据框，因为不用数据框或者用的数据框和模型对象的不一致就会报错，如果是多个待拟合对象，我们也要用data.frame把它变成数据框，

扫描二维码关注公众号，回复： 11342923 查看本文章

> predict(ureg01.lm,newdata = data.frame(X=c(43,44)),interval = "conf",level = 0.95)
      fit     lwr     upr
1 4277046 3726681 4827410
2 4137316 3576386 4698245

现在讨论预测，假设有一个球员名次是43名，想要预测他的工资大概是多少，我们也用predict，只是区间估计命令换成pred

> predict(ureg01.lm,newdata = data.frame(X=43),interval = "pred",level = 0.95)
      fit      lwr      upr
1 4277046 -8232809 16786901

可以发现第一列的结果和拟合是一样的，只是区间估计的结果不一样了，相同的置信水平下，预测的置信区间会比拟合的更大，但如果这名球员是没有被估计模型的样本包括在内的，一般我们还是需要用预测的。