一、基本思想
假设,在长期大量的消费支出记录中,每个月最终的消费金额,与这个月中,任意一天的消费金额存在着一定的关系。同时我们假设这个关系可以用多元回归和均值回归任意种方式来描述,记作:
- 多元回归:
f m l r ( x 1 , x 2 , . . . , x n ) = ξ + ρ 1 x 1 + ρ 2 x 2 + . . . + ρ n x n ( n < 31 ) f_{mlr}(x_1,x_2,...,x_n) = \xi + \rho_1x_1 + \rho_2x_2 + ... + \rho_nx_n \space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space (n<31) fmlr(x1,x2,...,xn)=ξ+ρ1x1+ρ2x2+...+ρnxn (n<31)
- 均值回归:
f m r ( x 1 , x 2 , . . . , x n ) = T ∗ x 1 + x 2 + . . . + x n n ( T 表示当月天数, n < 31 ) f_{mr}(x_1,x_2,...,x_n) = T * \frac{x_1 + x_2 + ... + x_n}{n} \space\space\space\space\space (T表示当月天数,n<31) fmr(x1,x2,...,xn)=T∗nx1+x2+...+xn (T表示当月天数,n<31)
接下来将对两种模型进行探究。
二、步骤
1、获取数据
我从一份网上的账目清单中获取了一份支出报告。其中包含近八年共5572条消费记录:
由于我们只需要研究一个月里面,每天的支出和月总支出的关系,所以我们把数据整理成下面的样式:
使用Excel的透视表功能可以完成上述转换。
对于上图,行索引为[年-月]
组合,比如202101
表示2021年1月份。每行的数据表示每天的实际支出,最特殊的是最后一个值,为当月的总支出。为了抵消通胀影响,我们只取最近四年的数据,月份取30天,空值替换为0。由此我们获得了由30个列向量组成的30*55
日消费矩阵,以及一个单独的向量组成的1*55
总金额矩阵。以下代称X
矩阵和Y
矩阵。
将这两个矩阵写入Matlab。
2、观察数据
刚拿到的数据可能存在脏数据,所以第一步一定是将离群值处理成一个合适的值,我们将清洗后的列向量矩阵命名为x_zscore
:
x_zscore = filloutliers(X, 'clip', 'mean'); % 离群值清洗,以均值的三倍为离群值,将离群值替换为均值
% -- 在先期实验中,已通过对30个自变量与因变量共3288条数据的模拟,确定了不同参数对应的绝对值累计误差:
% 'clip' 'mean' 98,674
% 'center' 'mean' 286,966
% 'clip' 'quartiles' 637,770
% 'center' 'quartiles' 695,930
% 'clip' 'median' 945,326
% 'center' 'median' 1,358,290
% 所以可得最好的离群值处理办法为'clip','mean'组合
清洗完成之后,我们需要观察x_zscore 中的每个列向量,是否与Y存在着线性关系,所以我们做散点图:
[m, ~] = size(X);
for i = 1:LENGTH
subplot(5, 6, i), plot(x_zscore(:, i), Y, '+'), xlabel(i)
end
从图中可以看出,大部分的天数与最终的消费金额是存在一定关系的,这就意味着我们可以进行下一步试探。
3、建立多元回归
3.1求参数向量
回顾一下多元线性回归公式:
f m l r ( x 1 , x 2 , . . . , x n ) = ξ + ρ 1 x 1 + ρ 2 x 2 + . . . + ρ n x n ( n < 31 ) f_{mlr}(x_1,x_2,...,x_n) = \xi + \rho_1x_1 + \rho_2x_2 + ... + \rho_nx_n \space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space\space (n<31) fmlr(x1,x2,...,xn)=ξ+ρ1x1+ρ2x2+...+ρnxn (n<31)
我们需要求出每个自变量对应的参数。我们定义参数向量:
b ( n ) → = [ ξ ρ 1 ρ 2 . . . ρ n ] \overrightarrow{b_{\left(n\right)}} = \begin{bmatrix} \xi &\rho_1& \rho_2 &... &\rho_n \end{bmatrix} b(n)=[ξρ1ρ2...ρn]
受益于计算机技术的发展,我们可以使用计算机计算出公式里面的参数:
% 获取多元回归模型的所有参数
b = ones(30, 1);
for i = 1:30
b(1:i + 1, i + 1) = funcx(i, x_zscore , Y);
end
b(:, 1) = [];
function result = funcx(LENGTH, X, Y) % LENGTH表示给与的自变量个数,也就是给多少个天数作为预测的自变量
x = [ones(m, 1), X(:, 1:LENGTH)];
y = Y;
[b, ~, ~, ~, ~] = regress(y, x);
result = b;
end
最终输出的b是一个31*30
列向量矩阵,我们称之为参数向量矩阵b
。其中第n
个参数向量(列向量)表示自变量个数为n时,所对应的回归系数:
b ( n ) → = [ ξ ρ 1 ρ 2 . . . ρ n ] \overrightarrow{b_{\left(n\right)}} = \begin{bmatrix} \xi &\rho_1& \rho_2 &... &\rho_n \end{bmatrix} b(n)=[ξρ1ρ2...ρn]
如果不方便理解请看下面的实际例子:
比如我们用2021年1月份的前三天来预测整个月份的开支:
2021年1月份的前三天,对应的自变量为:
X ( 202101 , 3 ) = [ 0.000 100.380 336.000 ] X_{(202101,3)} = \begin{bmatrix} 0.000& 100.380 &336.000 \end{bmatrix} X(202101,3)=[0.000100.380336.000]
对应的参数向量为:
b ( 3 ) → = [ 4615.743 5.455 − 0.893 1.049 ] \overrightarrow{b_{\left(3\right)}} = \begin{bmatrix} 4615.743 &5.455& -0.893 &1.049 \end{bmatrix} b(3)=[4615.7435.455−0.8931.049]
那么基于这三天所做出的预测就是
Y ( 202101 , 3 ) = 4615.743 + 5.455 ∗ 0.000 − 0.893 ∗ 100.380 + 1.049 ∗ 336.000 = 4878.568 Y_{(202101,3)} = 4615.743 +5.455 * 0.000 - 0.893 * 100.380 + 1.049 * 336.000 = 4878.568 Y(202101,3)=4615.743+5.455∗0.000−0.893∗100.380+1.049∗336.000=4878.568
3.2 定义多元线性回归方程
有了上面的参数向量矩阵b
,我们就可以定义一个用于预测的函数:
% 构建多元线性回归函数
% 输入自变量组成的行向量a 以及参数向量矩阵b,输出预测的总开支result:
% Example: MLR([100 200 300] ,b)
function result = MLR(a, b)
[~, LENGTH] = size(a);
result = a * b(2:1 + LENGTH, LENGTH) + b(1, LENGTH);
end
4、建立均值回归
回顾一下均值回归:
f m r ( x 1 , x 2 , . . . , x n ) = T ∗ x 1 + x 2 + . . . + x n n ( T 表示当月天数, n < 31 ) f_{mr}(x_1,x_2,...,x_n) = T * \frac{x_1 + x_2 + ... + x_n}{n} \space\space\space\space\space (T表示当月天数,n<31) fmr(x1,x2,...,xn)=T∗nx1+x2+...+xn (T表示当月天数,n<31)
均值回归中,T 表示当月的天数。均值回归比较简单,我们可以直接根据公式得出均值回归方程:
% 构建均值回归函数
% 输入自变量组成的行向量a 以及天数 T,输出预测的总开支result:
% Example: MR([100 200 300] ,31)
function result = MR(a, T)
[~, LENGTH] = size(a);
if LENGTH == 1
result = (a / LENGTH) * T;
else
result = transpose((sum(transpose(a)) / LENGTH) * T);
end
end
三、模型搭建
为了让预测模型更加精确,我们将采用两边结合的方法,将多元模型与均值回归有机的结合在一起,最终形成一个相对可靠的预测模型。
1、思路说明
对于上文的两个模型,其最大的问题在于其精准度与提供的自变量个数高度相关。对于这一影响,我们将采取相关系数作为衡量标准,求出不同自变量数量下,预测值与实际值的相关系数分布。
2、相关系数评估
2.1 获取预测值向量矩阵
从提供1个变量开始,一直到提供30个变量,我们分别求出不同自变量个数下面,所获得的预测值
Predicted_BY_MLR = ones(55, 1);
Predicted_BY_MR = ones(55, 1);
for i = 1:30
Predicted_BY_MLR(:, i + 1) = MLR(X(:, 1:i), b);
Predicted_BY_MR(:, i + 1) = MR(X(:, 1:i), 30);
end
Predicted_BY_MLR(:, 1) = []; % 存储多元回归预测值
Predicted_BY_MR(:, 1) = []; % 存储均值回归预测值
2.2 获取相关系数矩阵
将每一组预测值与实际值做协方差,并求出对应的相关系数
Correlation_By_MLR = ones(1,30); % 存储多元回归相关系数
Correlation_By_MR = ones(1,30); % 存储均值回归相关系数
for i = 1:30
Cov_By_MLR = cov(Predicted_BY_MLR(:, i), Y);
Cov_By_MR = cov(Predicted_BY_MR(:, i), Y);
Correlation_By_MLR(i) = Cov_By_MLR(1, 2) / sqrt(Cov_By_MLR(1, 1) * Cov_By_MLR(2, 2));
Correlation_By_MR(i) = Cov_By_MR(1, 2) / sqrt(Cov_By_MR(1, 1) * Cov_By_MR(2, 2));
end
2.3 绘图
% 绘图
x1_lab = 1:30;
x2_lab = 1:30;
plot(x1_lab,Correlation_By_MLR,'+',Color='red'),hold on,plot(x2_lab,Correlation_By_MR,'*',Color='blue')
% 红线表示多元回归,蓝线表示均值回归
从图上可以很明显的得出如下结论:
- 1、提供的天数低于5天时,两个模型预测精度都不高
- 2、在提供的天数没超过23天之前,多元回归会优于均值回归
- 3、在提供的天数超过23天之后,均值回归会优于多元回归
2.4 模型构建
基于2.3的结论,实际上我们已经得出了一个比较合理的模型方法,如果自变量的个数小于23天,则采用多元回归,如果超过23天,则采用均值回归。具体的预测值可以使用上文中对于的回归函数进行计算。