线性模型出现非正定矩阵的问题解释

常常在Mplus中求解线性结构方程的时候出现 如下警告：

 WARNING:  THE LATENT VARIABLE COVARIANCE MATRIX (PSI) IS NOT POSITIVE

1 背景：

         大约不少人找了很多书籍，要么一笔带过，要么只给出结论，让人终觉疑虑满满。是线性代数的知识，但是你去翻书，人也未必告诉你，那么这类问题在那里呢，这个问题属于数值分析领域。

百科是这样解释的：数值分析(numerical analysis)是研究分析用计算机求解数学计算问题的数值计算方法及其理论的学科，是数学的一个分支，它以数字计算机求解数学问题的理论和方法为研究对象。为计算数学的主体部分。

2 方程分解

现在从矩阵分析与计算的角度来解释

考虑这样一个方程组  A x = b，假设 A是两个被试样本的智商和成绩(两个自变量)

             智商  期中成绩  高考成绩      

A同学      a1          a2          a

B同学     b1           b2          b

我们在收集问卷的时候，一定是把这几个变量都收集了，这里的a和b是已知数，你在收集的时候一定找她们班主任要一份成绩，那个韦氏智力量表加上平时成绩和过往的期末成绩。依据频率学派的观点，我们建模的目的是求出是求出未知的参数，就是暂时不知道，但是是存在这么一个确定的值，建立一个这样的模型。我有变量了，这一堆值在我先验下，知道智商和期中成绩是会影响高考成绩的，（而且正太分布的）。那么，需要求解的过程是什么的呢。就是这自变量前面的系数，定量研究，就是用数值来精确的刻画，到底在多大的程度上起作用。

比如，最后得出这样的一个模型

        在期中成绩不变的情况下，1个单位的智商的引起2个单位的成绩的改变量，因为前面系数为正，所以是正向预测。这个2就是智商这个变量的系数，我们建立模型的目的就是用这些数据得出变量的系数。然后再来一个 学生，我们知道了两个自变量，就把这些系数乘以对应的变量，得出的最后结果就是预测了他的高考成绩。也就是说，我给定的这个智商 期中成绩数值中包含了两个部分，系数和变量的值。

看上面这个方程   A x = b  就是 A就是一堆自变量组成的矩阵 b就是因变量了，这个X矩阵正是我们建立模型要求的。

我们都求解过这样的方程，已知  

这里两个变量 x1   x2，我们很容易求解，用消元法。

这样的方法叫做，解析式求解;

但是当我们有很多的变量和样本的时候，假设是  6个变量100个样本这样的矩阵，我是没有办法用解析式求解的，于是我们只能使用近似的迭代方法来解。

推而广之：

在多元线性模型中，y-因变量       b - 系数矩阵  X- 自变量矩阵      E- 残差矩阵

给出这么一个矩阵

抽象的方程表示如下图（2）

图（2）表示的是最小二乘法就是要求解出， y - （系数和自变量乘积组成的理论）之前最小值

图（3）求解广义逆矩阵

我们关心的是什么时候，这个广义逆矩阵无解，

1，就是如果上图的自变量X矩阵中，（代表样本个数的行数n 小于变量的个数P ） ，代表要求解的未知个数大于方程的个数，（此时求解的模型过拟合）则无解

2 如果求解的矩阵的秩小于p ，就是说不是满秩矩阵，则无解，什么时候会出现这种情况呢，就是当自变量中的某一个变量可以由其他的矩阵表示出来，类似于出现 x1=2 x2 +3

在求解初等行变换的时候，就是这么干的。

我们在建模的时候假定变量之间是独立的，（如果你还记得的话，回归分析的前提是残差正态，变量独立，方差齐性，存在线性关系）比如在取得100个样本的智商和学习成绩这两个自变量中，会不会出现这100个人的  智商=2 倍学习成绩 +3呢，不会的，没有那么巧合的时候。

我们关心的是什么时候，这个逆矩阵无解，

1，就是如果上图的自变量X矩阵中，（代表样本个数的行数n 小于变量的个数P ） ，代表要求解的未知个数大于方程的个数，（此时求解的模型过拟合）则无解

2 如果求解的矩阵的秩小于p ，就是说不是满秩矩阵，则无解，什么时候会出现这种情况呢，就是当自变量中的某一个变量可以由其他的矩阵表示出来，类似于出现 x1=2 x2 +3

在求解初等行变换的时候，就是这么干的。

我们在建模的时候假定变量之间是独立的，（如果你还记得的话，回归分析的前提是残差正态，变量独立，方差齐性，存在线性关系）比如在取得100个样本的智商和学习成绩这两个自变量中，会不会出现这100个人的  智商=2 倍学习成绩 +3呢，不会的，没有那么巧合的时候。为什么？因为有残差的存在，总是存在这我们无法知道的一些东西影响了变量的值，我们无法知晓提出的不重要的变量，一股脑的放入方程外，叫做残差。

3  矩阵求解

再次拿出这个方程  A x = b，先给出这样一个简单的矩阵

行数等于列数是一个方阵，这个方阵大家看着很眼熟，通常我们在表示的时候，把左下角的一拿掉，就是上三角矩阵，这个1很好理解，

我们经常挂在嘴边的固定方差，这里就是固定方差，对角线上的叫做方差，不在对角线上的叫啥？协方差！

所以有协方差矩阵，到底拿协方差矩阵放入liserl模型中还是拿相关矩阵呢，协方差矩阵在在很多软件中计算的都有基于矩阵乘以它的逆矩阵的计算语句，比如在python中，所以协方差矩阵方便计算。同时 ，中心化矩阵在数理统计中非常有用。

这么说 相关系数就是标准化的协方差，如果有疑问，请看知乎，https://www.zhihu.com/question/20852004

这里给出的是一个相关系数组成的2*2矩阵。

0.99的意思是，这两个变量的相关度为0.99

给定的数据中，因变量是b，b在这两个同学的取值为1  0  

我们要求的是这个矩阵

得到   

 

如果我们将这个0.99替换成 0.991则得到的会是 

    二者差是

可以知道仅有输入系数0.001的误差，二者产生的结果却是相差 5000多倍，我们本来的目的是使得当前求解的样本矩阵和理论矩阵差值变小，这里非但没有变小还变大了，

如果在模型的迭代求解过程中出现第n 次的迭代值和低n+ 1次的值相差巨大。是非常不利于模型拟合的。

4 运算中的解释

想象一下，如果软件迭代再运行一步试这给系数为 0.992，如果看书认真的话，或者模型跑的次数多，你会看到迭代出来的值，如果出现模型相邻两步的差越来越大，那么软件没有办法求解的，候杰泰老师的书中解释过迭代的说法，简单来说是这样的，我拿观察的实例矩阵 与 模型求解的理论矩阵，让软件求解二者越来越小的距离，设定一个阈值，如果它们二者的差异小过0.05，那么我们就说，哈，不要迭代了，理论模型也不可能和你给定的矩阵100%一致，有一些出入是正常的。这时候数据验证了你的理论，论文可以发了。

这个理论矩阵是怎么工作的呢，试，第一步我先给它一个初始值（start），带进去估算二者的差异，大了，就在当前迭代减去上一步的迭代值。

这个拉姆达 字母可以是给定的一个值，比如0.05，也可以是让软件在迭代过程中，求偏导给出。那么，如果我第一次 调整的值和第二次的值乘0.05，还是有0.05*50000=250倍。如果软件在求解过程中一某几次出现在250再来一个 250倍，有三次这样的指数级变化，计算机运行速度再快也招架不住指数增长，势必会无法计算，导致这样接近无穷的值变为非数字，叫做上溢（overflow）出现由于输入轻微的扰动而迅速改变的函数对于科学计算来说是有问题的，这样的状况叫做病态条件（ill conditioned）敏感的（sensitive）。

原因是矩阵A的行列式 |A|  等于 

 0.01990000000000003 几乎为0  （故意把小数点拉的这么多，就是为了显示，着就是我的python软件能够计算出来的最大精度）。计算机计算这样的矩阵时候，如果副对角线元素之积的与1相差很小的时候，就是有一行出现0.923这样相关系数高的。就是约等于这个系数矩阵不可逆（奇异矩阵），行列式很小的时候，会给迭代运算带来巨大的麻烦。

是时候再冒出来一个概念了 ：矩阵不可逆，一个n阶方阵A称为可逆的，或非奇异的，如果存在一个n阶方阵B，使得 AB+BA=E

并称B是A的一个逆矩阵。不可逆的矩阵称为奇异矩阵。A的逆矩阵记作A-1。

若|A|≠0，则矩阵A可逆。且

如果 |A| 等于0，则这个方程就不成立，叫 A矩阵不可逆矩阵。

         一个n阶的实对称矩阵M是正定的的条件是当且仅当对于所有的非零实系数向量z，都有 

表示Z的转置。

如果实在不明白，看这个性质，正定矩阵有以下性质：正定矩阵的行列式恒为正；

 如上我们的行列式在计算的过程中几乎为0了，又因为在运算过程中精度等问题，导致在求解过程中出现的病态状态，最后计算机精确有限的情况下给出异常，为啥在Mplus中运算出来了呢，只是给你一个警告，而不是直接报错，就是在迭代求解的过程中设置的语句有，如果有几次的迭代求解，几乎指数级变化，那就在output中给一句 warning，这时候很可能迭代求出的解不可信。本来是要求样本矩阵和理论矩阵二者的差很小，但是由于在计算中二者的距离越来越大，在计算中途可能迭代停止，所以警告我们的是，这时候我们求解的系数可能是有问题的。

回到问题本身，剔除高度相关的变量，防止自变量间共线性，需要剔除冗余变量。（就是可以有其他变量生成的变量，拿一个变量乘以2倍再放入模型中，新生成的变量并没有提供多少有效信息），是一定会出现问题的，同时也是求解计算的要求。

 WARNING:  THE LATENT VARIABLE COVARIANCE MATRIX (PSI) IS NOT POSITIVE

这句警告的意思就是，模型非正定。

5 实例验证

我们拿mplus 7.4 跑一下这样一组数据

自变量  x1    x2    因变量 y

          10     20        1

          5       11          0

这样一组程序

title : this is test for matrix 1;
data: file is te.csv;
variable: names are x1 x2 y;
model:x1 x2 on y;
output:modindices;

产生这样的警告

WARNING:  THE SAMPLE COVARIANCE IS SINGULAR. 
WARNING:  THE SAMPLE CORRELATION OF X2 AND X1 IS  1.000.
WARNING:  THE SAMPLE CORRELATION OF Y AND X1 IS  1.000.

意思就是说，发现矩阵是奇异矩阵（也叫不可逆矩阵）

x1 和 x2 的相关系数是1 事实上也很容易看出来，几乎为1

总的来说不管奇异矩阵，非正定矩阵，就是一句话，如果所求的系数（协方差）矩阵的行列式几乎为0 ，那么，计算机运行求解会出现问题。

行列式出现几乎为0的情况是什么呢，就是存在某；两个或者多个变量的系数接近1，依据行列式求解的过程，见上文。

所以在求解过程中如果出现非正定矩阵，需要做的就是

（1）检查变量是不是存在高相关，导致计算过程病态，逆矩阵无法得出，如果在相关矩阵出检查出来有变量存在高相关，（关于检查高相关的几个指标后文有写）要么删除。（逐步回归）要么主成分去掉冗余变量。，可以看作是预处理，降维。

（2）增加样本量，这个可以从两方面来理解，第一如果是两个样本 两个变量取值空间在01的时候，很容易出现相关矩阵，如果有更多有代表性的样本进来，那么，这100个样本组成的两行出现相关矩阵的可能性几乎没有，第二，为了防止模型的过拟合。 

（3） 采用收缩估计（ShrinkageMethods ）常见的有 ridge  regression，lasso（套索） 以及最小角（Least Angle Regression ）对系数矩阵进行约束，使得计算出的解，不再出现不稳定。或者依据最小角的结果删除变量。这是从算法的选择角度，因为单纯采用梯度下降的ML求解遇到困难。

存在几个小问题:

1 是为了手工计算的方便，矩阵选择的得是2*2，但是也为了照顾直观一目了

2  给出的概念百科略微又些多，但也不得已

3  依然有几个点，是省略了的，因为在有限的篇幅不可能说明白

4  基本上采用夹叙夹议的方式来阐述，未免有些不严谨，和小瑕疵

5 突然觉得，讲结构方程的书中没有矩阵分析和高级代数的文献引用和讲解真真是一件很奇怪的事情

参考书籍 ：

1  矩阵分析与计算 p 31-p 33页

2 矩阵分析与应用