天天看点

numpy和MATLAB中协方差矩阵的计算

前言:最近重新看了下PCA算法的实现,在里面用到了协方差矩阵的计算。为了今后方便回顾,特此记录。

这里只针对协方差矩阵的计算做说明,并与numpy和MATLAB中的协方差矩阵计算做对比。关于协方差矩阵的定义这里不再赘述。

混淆点

关于矩阵中数据是按照列排列和按照行排列的出来的协方差矩阵是不同的。

如果数据按照列排列,例如一个3*2的矩阵,那么该矩阵有2个数据(或者叫观测值),3个变量。

如果数据按照行排列,同样以3*2的矩阵为例,那么该矩阵有3个数据(或者叫观察值),2个变量。

不过不管矩阵中数据是按照行排列还是列排列,最后的协方差矩阵的维度只与变量有关。

计算

下面将分别解释数据按行排列和列排列的协方差矩阵的计算。

一、数据按照行排列

矩阵A包含2个观察值,每个观察值有3个变量。

numpy和MATLAB中协方差矩阵的计算

①计算变量的均值(也就是每列的均值)

numpy和MATLAB中协方差矩阵的计算

 ②变量减去均值

numpy和MATLAB中协方差矩阵的计算

③计算协方差矩阵(m代表数据的个数)

numpy和MATLAB中协方差矩阵的计算

 ④代码验证

numpy

In [39]: A = np.array([[1,5,3],[6,7,2]])

In [40]: Amean = np.mean(A,axis=0)

In [41]: U = A - Amean

In [42]: U.T @ U /(2-1)
Out[42]:
array([[12.5,  5. , -2.5],
       [ 5. ,  2. , -1. ],
       [-2.5, -1. ,  0.5]])

In [43]: np.cov(A,rowvar=0)#numpy中默认是行是变量,列是数据。所有将rowvar置0,代表行为数据
Out[43]:
array([[12.5,  5. , -2.5],
       [ 5. ,  2. , -1. ],
       [-2.5, -1. ,  0.5]])
           

MATLAB

>> A = [1 5 3;6 7 2]

A =

     1     5     3
     6     7     2

>> cov(A)

ans =

   12.5000    5.0000   -2.5000
    5.0000    2.0000   -1.0000
   -2.5000   -1.0000    0.5000
           

二、数据按照列排列

矩阵A包含3个观察值,每个观察值有2个变量。

numpy和MATLAB中协方差矩阵的计算

①计算变量的均值(也就是每行的均值)

numpy和MATLAB中协方差矩阵的计算

 ②变量减去均值

numpy和MATLAB中协方差矩阵的计算

 ③计算协方差矩阵(m代表数据的个数)

numpy和MATLAB中协方差矩阵的计算

 ④代码验证

In [53]: A = np.array([[1,5,3],[6,7,2]])

In [54]: Amean = np.mean(A,axis=1,keepdims=True)

In [55]: U = A - Amean

In [56]: U @ U.T / (3 - 1)
Out[56]:
array([[4., 1.],
       [1., 7.]])

In [57]: np.cov(A)#numpy默认行是变量
Out[57]:
array([[4., 1.],
       [1., 7.]])
           

总结

在计算矩阵的协方差矩阵时关键点是要弄清楚矩阵是一行对应一组数据还是一列对应一组数据,计算出来的协方差矩阵都是n维的方阵(n代表变量的个数)。在计算时,第一步都是用原始矩阵减去变量的均值得到变量和为0的矩阵U;第二步根据数据的排列方式进行计算,若一行代表一个数据则选用公式U.T*U/(m-1),若一列代表一个数据则U*U/(m-1)。