在 MATLAB 中使用回归函数时出现秩不足警告
Getting rank deficient warning when using regress function in MATLAB
我有一个包含 30 个自变量的数据集,我尝试使用 regress
函数在 MATLAB R2010b 中执行线性回归。
我收到一条警告,指出我的矩阵 X is rank deficient to within machine precision
.
现在执行这个函数得到的系数和实验的不匹配
谁能建议我如何对这个由 30 个变量组成的方程进行回归分析?
继续我们的讨论,您收到该警告的原因是因为您拥有所谓的 underdetermined system。基本上,您有一组约束,其中您想要解决的变量多于可用数据。欠定系统的一个例子是这样的:
x + y + z = 1
x + y + 2z = 3
(x,y,z)
有无数种组合可以解出上述系统。例如,(x, y, z) = (1, −2, 2), (2, −3, 2), and (3, −4, 2)
。在您的情况下,等级不足意味着 多于 一组回归系数可以满足描述输入变量和输出观察值之间关系的控制方程。这可能是 regress
的输出与您的真实回归系数不匹配的原因。虽然答案不同,但要知道输出是 一个 可能的答案。通过 运行 到 regress
与您的数据,如果我将您的观察矩阵定义为 X
并且您的输出向量为 Y
:[=42=,这就是我得到的]
>> format long g;
>> B = regress(Y, X);
>> B
B =
0
0
28321.7264417536
0
35241.9719076362
899.386999172398
-95491.6154990829
-2879.96318251964
-31375.7038251919
5993.52959752106
0
18312.6649115112
0
0
8031.4391233753
27923.2569044728
7716.51932560781
-13621.1638587172
36721.8387047613
80622.0849069525
-114048.707780113
-70838.6034825939
-22843.7931997405
5345.06937207617
0
106542.307940305
-14178.0346010715
-20506.8096166108
-2498.51437396558
6783.3107243113
您可以看到有 7 个回归系数等于 0,对应于 30 - 23 = 7。我们有 30 个变量和 23 个约束条件可供使用。请注意,这不是唯一可能的解决方案。 regress
本质上是计算最小平方误差解,使得 Y - X*B
的残差和具有最少的误差。这基本上简化为:
B = X^(*)*Y
X^(*)
就是所谓的pseudo-inverse of the matrix. MATLAB has this available, and it is called pinv
。因此,如果我们这样做:
B = pinv(X)*Y
我们得到:
B =
44741.6923363563
32972.479220139
-31055.2846404536
-22897.9685877566
28888.7558524005
1146.70695371731
-4002.86163441217
9161.6908044046
-22704.9986509788
5526.10730457192
9161.69080479427
2607.08283489226
2591.21062004404
-31631.9969765197
-5357.85253691504
6025.47661106009
5519.89341411127
-7356.00479046122
-15411.5144034056
49827.6984426955
-26352.0537850382
-11144.2988973666
-14835.9087945295
-121.889618144655
-32355.2405829636
53712.1245333841
-1941.40823106236
-10929.3953469692
-3817.40117809984
2732.64066796307
您看到没有零系数,因为 pinv
使用 L2 范数找到解决方案,这促进了 "spreading" 错误(因为缺少更好的术语)。您可以通过执行以下操作来验证这些是正确的回归系数:
>> Y2 = X*B
Y2 =
16.1491563400241
16.1264219600856
16.525331600049
17.3170318001845
16.7481541301999
17.3266932502295
16.5465094100486
16.5184456100487
16.8428701100165
17.0749421099829
16.7393450000517
17.2993993099419
17.3925811702017
17.3347117202356
17.3362798302375
17.3184486799219
17.1169638102517
17.2813552099096
16.8792925100727
17.2557945601102
17.501873690151
17.6490477001922
17.7733493802508
类似地,如果我们使用 regress
的回归系数,那么 B = regress(Y,X);
然后做 Y2 = X*B
,我们得到:
Y2 =
16.1491563399927
16.1264219599996
16.5253315999987
17.3170317999969
16.7481541299967
17.3266932499992
16.5465094099978
16.5184456099983
16.8428701099975
17.0749421099985
16.7393449999981
17.2993993099983
17.3925811699993
17.3347117199991
17.3362798299967
17.3184486799987
17.1169638100025
17.281355209999
16.8792925099983
17.2557945599979
17.5018736899983
17.6490476999977
17.7733493799981
存在一些细微的计算差异,这是可以预料的。同样,我们也可以使用mldivide
:
来求答案
B = X \ Y
B =
0
0
28321.726441712
0
35241.9719075889
899.386999170666
-95491.6154989513
-2879.96318251572
-31375.7038251485
5993.52959751295
0
18312.6649114859
0
0
8031.43912336425
27923.2569044349
7716.51932559712
-13621.1638586983
36721.8387047123
80622.0849068411
-114048.707779954
-70838.6034824987
-22843.7931997086
5345.06937206919
0
106542.307940158
-14178.0346010521
-20506.8096165825
-2498.51437396236
6783.31072430201
你可以看到这与 regress
给你的很奇怪。那是因为 \
是一个更聪明的运算符。根据矩阵的结构,它会通过不同的方法找到系统的解决方案。我想让你看一下 Amro 的 post,它讨论了 mldivide
在检查正在操作的输入矩阵的属性时使用的算法:
How to implement Matlab's mldivide (a.k.a. the backslash operator "\")
你应该从这个答案中得到的是,你当然可以继续使用那些回归系数,它们或多或少会为你提供每组输入的每个 Y
值的预期输出对于 X
。但是,请注意这些系数 不是唯一的 。正如您所说,您的地面真值系数与 regress
的输出不匹配,这一点很明显。它不匹配,因为它生成了 另一个 满足您提供的限制的答案。
如果你有一个欠定的系统,有不止一个答案可以描述这种关系,正如你在我上面的实验中看到的那样。
我有一个包含 30 个自变量的数据集,我尝试使用 regress
函数在 MATLAB R2010b 中执行线性回归。
我收到一条警告,指出我的矩阵 X is rank deficient to within machine precision
.
现在执行这个函数得到的系数和实验的不匹配
谁能建议我如何对这个由 30 个变量组成的方程进行回归分析?
继续我们的讨论,您收到该警告的原因是因为您拥有所谓的 underdetermined system。基本上,您有一组约束,其中您想要解决的变量多于可用数据。欠定系统的一个例子是这样的:
x + y + z = 1
x + y + 2z = 3
(x,y,z)
有无数种组合可以解出上述系统。例如,(x, y, z) = (1, −2, 2), (2, −3, 2), and (3, −4, 2)
。在您的情况下,等级不足意味着 多于 一组回归系数可以满足描述输入变量和输出观察值之间关系的控制方程。这可能是 regress
的输出与您的真实回归系数不匹配的原因。虽然答案不同,但要知道输出是 一个 可能的答案。通过 运行 到 regress
与您的数据,如果我将您的观察矩阵定义为 X
并且您的输出向量为 Y
:[=42=,这就是我得到的]
>> format long g;
>> B = regress(Y, X);
>> B
B =
0
0
28321.7264417536
0
35241.9719076362
899.386999172398
-95491.6154990829
-2879.96318251964
-31375.7038251919
5993.52959752106
0
18312.6649115112
0
0
8031.4391233753
27923.2569044728
7716.51932560781
-13621.1638587172
36721.8387047613
80622.0849069525
-114048.707780113
-70838.6034825939
-22843.7931997405
5345.06937207617
0
106542.307940305
-14178.0346010715
-20506.8096166108
-2498.51437396558
6783.3107243113
您可以看到有 7 个回归系数等于 0,对应于 30 - 23 = 7。我们有 30 个变量和 23 个约束条件可供使用。请注意,这不是唯一可能的解决方案。 regress
本质上是计算最小平方误差解,使得 Y - X*B
的残差和具有最少的误差。这基本上简化为:
B = X^(*)*Y
X^(*)
就是所谓的pseudo-inverse of the matrix. MATLAB has this available, and it is called pinv
。因此,如果我们这样做:
B = pinv(X)*Y
我们得到:
B =
44741.6923363563
32972.479220139
-31055.2846404536
-22897.9685877566
28888.7558524005
1146.70695371731
-4002.86163441217
9161.6908044046
-22704.9986509788
5526.10730457192
9161.69080479427
2607.08283489226
2591.21062004404
-31631.9969765197
-5357.85253691504
6025.47661106009
5519.89341411127
-7356.00479046122
-15411.5144034056
49827.6984426955
-26352.0537850382
-11144.2988973666
-14835.9087945295
-121.889618144655
-32355.2405829636
53712.1245333841
-1941.40823106236
-10929.3953469692
-3817.40117809984
2732.64066796307
您看到没有零系数,因为 pinv
使用 L2 范数找到解决方案,这促进了 "spreading" 错误(因为缺少更好的术语)。您可以通过执行以下操作来验证这些是正确的回归系数:
>> Y2 = X*B
Y2 =
16.1491563400241
16.1264219600856
16.525331600049
17.3170318001845
16.7481541301999
17.3266932502295
16.5465094100486
16.5184456100487
16.8428701100165
17.0749421099829
16.7393450000517
17.2993993099419
17.3925811702017
17.3347117202356
17.3362798302375
17.3184486799219
17.1169638102517
17.2813552099096
16.8792925100727
17.2557945601102
17.501873690151
17.6490477001922
17.7733493802508
类似地,如果我们使用 regress
的回归系数,那么 B = regress(Y,X);
然后做 Y2 = X*B
,我们得到:
Y2 =
16.1491563399927
16.1264219599996
16.5253315999987
17.3170317999969
16.7481541299967
17.3266932499992
16.5465094099978
16.5184456099983
16.8428701099975
17.0749421099985
16.7393449999981
17.2993993099983
17.3925811699993
17.3347117199991
17.3362798299967
17.3184486799987
17.1169638100025
17.281355209999
16.8792925099983
17.2557945599979
17.5018736899983
17.6490476999977
17.7733493799981
存在一些细微的计算差异,这是可以预料的。同样,我们也可以使用mldivide
:
B = X \ Y
B =
0
0
28321.726441712
0
35241.9719075889
899.386999170666
-95491.6154989513
-2879.96318251572
-31375.7038251485
5993.52959751295
0
18312.6649114859
0
0
8031.43912336425
27923.2569044349
7716.51932559712
-13621.1638586983
36721.8387047123
80622.0849068411
-114048.707779954
-70838.6034824987
-22843.7931997086
5345.06937206919
0
106542.307940158
-14178.0346010521
-20506.8096165825
-2498.51437396236
6783.31072430201
你可以看到这与 regress
给你的很奇怪。那是因为 \
是一个更聪明的运算符。根据矩阵的结构,它会通过不同的方法找到系统的解决方案。我想让你看一下 Amro 的 post,它讨论了 mldivide
在检查正在操作的输入矩阵的属性时使用的算法:
How to implement Matlab's mldivide (a.k.a. the backslash operator "\")
你应该从这个答案中得到的是,你当然可以继续使用那些回归系数,它们或多或少会为你提供每组输入的每个 Y
值的预期输出对于 X
。但是,请注意这些系数 不是唯一的 。正如您所说,您的地面真值系数与 regress
的输出不匹配,这一点很明显。它不匹配,因为它生成了 另一个 满足您提供的限制的答案。
如果你有一个欠定的系统,有不止一个答案可以描述这种关系,正如你在我上面的实验中看到的那样。