基于MATLAB的回归分析实践：建模农村居民家庭人均纯收入与消费支出间的关系

1. 数据说明

本文使用的数据来自于国家统计局1997年-2012年的年度数据。选取的数据指标为农村居民家庭平均每人纯收入（元），农村居民家庭平均每人消费支出（元）。

2. 模型构建

模型：$Y={\beta }_0+{\beta }_1\ast x+\sigma$
自变量$x$：农村居民家庭平均每人纯收入（元）；
因变量$y$：农村居民家庭平均每人消费支出（元）。
利用MATLAB中的regress函数实现回归分析，具体的代码如下：

x = [2090.1, 2162, 2210.3, 2253.4, 2366.4, 2475.6, 2622.2, 2936.4, 3254.9, 3587, 4140.4, 4760.6, 5153.2, 5919, 6977.3, 7916.6];
y = [1617.2, 1590.3, 1577.4, 1670.1, 1741.1, 1834.3, 1943.3, 2184.7, 2555.4, 2829, 3223.9, 3660.7, 3993.5, 4381.8, 5221.1, 5908];
X = [ones(length(y), 1), x'];
Y = y';
[b, bint, r, rint, stats] = regress(Y, X);
b, bint, format short g;stats

3. 运行结果及评估分析

对于regress函数解读如下：
[b, bint, r, rint, stats] = regress(Y, X, alpha) %求回归系数的点估计和区间估计、并检验回归模型

1. bint表示回归系数的区间估计；
2. r表示残差；
3. rint表示置信区间；
4. stats表示用于检验回归模型的统计量，有三个数值：相关系数$R^2$、$F$值、与$F$对应的概率$p$。
   其中相关系数$R^2$越接近 1，说明回归方程越显著；
   当$F>F_{1-\alpha }(1,n-2)$时，拒绝$H_0$；$F$越大，说明回归方程越显著；
   与$F$对应的概率$p<\alpha$时，拒绝$H_0$。
5. alpha表示显著性水平(缺省时为 0.05)

因此，根据运行结果可进行评估与分析。

3.1 评估

3.1.1 回归方程的显著性

（F检验） 回归方程为Y = 6.6638 + 0.75339 * x， $F=4997.3$， $p=2.7855e-19<0.05$，故回归模型成立；

（t检验） 假设检验$H_0:{\beta }_1=0\leftrightarrow H_1:{\beta }_1\ne 0$，计算${\hat{\sigma }}^2=\frac{1}{n-2}{\sum }_{i=1}^n(y_i-\widehat{b_0}-\widehat{b_1}{x}_i{)}^2$，代码如下：

n=length(x);
yh=beta_0+beta_1*x;
s2=sum((y-yh).^2)/(n-2);
s1=sqrt(s2);

上述代码的运行结果如下图所示，得到${\hat{\sigma }}^2=5816.7080$，$\sigma =76.2673$。

继续计算$T=\frac{\widehat{b_1}\sqrt{l_{xx}}}{\hat{\sigma }}$，代码如下：

sb1=s1/sqrt(lxx);
T=beta_1/sb1;
alpha=0.05;
t0=tinv(1-alpha/2,n-2);
Tt=[T,t0]

运行结果如下图所示，由于$T=70.69152>t_0=2.14479$，因此拒绝$H_0$，认为回归方程具有显著性。

3.1.2 拟合优度

由stats的第一个数得$R^2=0.99721$。

3.1.3 相关系数

$R^2$开方得相关系数$R=0.9986$。

3.2 分析

3.2.1 各回归参数的点估计

从b里可以得到回归系数的参数估计（也为无偏估计），因而点估计${\beta }_0$为6.6638，${\beta }_1$为0.75339。

3.2.2 各回归参数的0.95区间估计

从bint里可以得到回归系数区间估计（0.95）：${\beta }_0$：[-89.374 , 102.7]， ${\beta }_1$：[0.73053, 0.77625]。

3.2.3 方差分析表

方差分析表计算代码如下：

sx = sum(x);
sy = sum(y);
sxy = sum(x.*y);
sxx = sum(x.^2);
syy = sum(y.^2);

mx = 1/16* sum(x);
my = 1/16* sum(y);

lxy = sxy - 16 * mx * my;
lxx = sxx - 16 * mx^2;
lyy = syy - 16 * my^2;

beta_1 = lxy / lxx;
beta_0 = my - beta_1 * mx;

%计算方差分析表
ssr = beta_1 * lxy;
sse = lyy - beta_1 * lxy;

运行结果如下图所示

方差来源	平方和	自由度	均方	F值	显著性
回归	29067783.3249	1	29067783.3249	4997.2911	*
剩余	81433.9126	14	5816.70804
总和	29149217.2375	15

由于计算出的$F$值为$4997.2911$，比$F_{0.05}(1,14)=4.60$大很多，因此显著性很大。

3.2.4 残差分析图

采用rcoplot(r, rint)函数利用之前回归分析计算的残差矩阵进行残差图绘制，残差图如下所示：

由于实际中影响支出的因素不只有收入一个，仅以收入作为自变量进行一元线性回归确实误差较大，因此图中残差还可以接受。对于异方差性，它是指随机误差项具有不同的方差，即从图上看出方差具有明显的不同的分布特性，而从残差图上可以看出，残差的分布没有明显不同的分布特性，故没有太大的异方差性。

4.结论及预测

对$x_0=10000$进行预测，得$y=7540.5524$，预测区间为[7320.31665539802, 7760.78818825828]，具体代码如下：

%x0=10000的预测
x0 = 10000;
y0hat = beta_0 + beta_1 * x0;
sig2hat = 1 / (n - 2)*(lyy - lxy^2/lxx);
sighat = sqrt(sig2hat)
deltax0 = sighat*tinv(0.975, n-2)*sqrt(1+1/n+(x0-mean(x))^2/lxx)

y0hat - deltax0
y0hat + deltax0
y0hat

运行结果如下图所示：

利用fitlm函数对预测图像进行绘制，代码及图像如下：

tb = table(x', y', 'VariableNames', {'income', 'expend'});
model = fitlm(tb, 'expend~income');
plot(model);
model