如何在MATLAB中没有任何内置函数或循环的情况下计算协方差矩阵？

Question

是否有可能在不使用MATLAB中的任何内置函数或循环的情况下找到矩阵的协方差？我完全不知道解决这个问题的想法.

我想的是:

cov(x,y) = 1/(n-1) .* (x*y)

但是,我不认为这会奏效.有任何想法吗？

Answer 1

这是如何数值计算协方差矩阵的一个很好的例子. http://www.itl.nist.gov/div898/handbook/pmc/section5/pmc541.htm.但是,为了完整起见,让我们把它放在这篇文章中. 我对"内置"函数的含义感到困惑,因为协方差要求你对矩阵的列求和.如果您不能使用任何内置函数来总结这些元素,那么我不知道如何在不使用for循环的情况下执行此操作. 编辑:我想出了如何在不使用内置函数或循环的情况下执行此操作,但您需要使用它size来确定矩阵中有多少行...除非您在函数中将其指定为常量.

在数值上,您可以计算协方差矩阵,如下所示:

从本质上讲,我^个行和第j ^个的协方差矩阵的列是这样的,你拿的列的产品的总和i减去列的平均值i与柱j减柱的平均值j.现在,添加它们,然后除以n - 1.这被称为无偏估计.你还会注意到这个矩阵是对称的,因为即使你翻转顺序(即查看j列i后面的列),答案应该仍然是相同的.我假设你不能使用meanMATLAB,所以让我们从第一原则开始.

首先,计算一个计算每列平均值的行向量.你可以做什么来计算所有列的总和而不使用sum,因为它也是一个内置函数,将这个1s的行向量与你的矩阵相乘A,输出将是一个包含所有列的总和的行向量的列.因此,这样做:

one_vector(1:size(A,1)) = 1;
mu = (one_vector * A) / size(A,1);

第一行代码的技巧是我们动态创建一个与矩阵中行数相同的数组A.我们完全填满了1s.请注意,您可以使用ones,但您说您不能使用任何内置函数. mu将包含我们在所有列上的向量.

现在,让我们通过用均值减去每一列来预处理数据,因为这就是我们所做的定义.要做到这一点没有任何内置函数,你可以做的是用他们各自的方法减去所有列,重复mu多次,因为我们有1s one_vector.因此:

A_mean_subtract = A - mu(one_vector, :);

这里有点棘手(而且很酷).如果我们转置矩阵A,您将看到行成为列,列成为行.如果我们采用这个转置并乘以原始矩阵,我们实际上会得到矩阵的列i和列之间的乘积之和.这是我们协方差计算的第一部分.然后我们除以.因此,我们的协方差很简单:jAn - 1

covA = (A_mean_subtract.' * A_mean_subtract) / (size(A,1) - 1);

这是一个简单的例子,以及我在上面展示的网站上看到的内容.假设A是这样的:

A = [4 2 0.5; 4.2 2.1 0.59; 3.9 2.0 0.58; 4.3 2.1 0.62; 4.1 2.2 0.63]

A =

    4.0000    2.0000    0.5000
    4.2000    2.1000    0.5900
    3.9000    2.0000    0.5800
    4.3000    2.1000    0.6200
    4.1000    2.2000    0.6300

运行上面的代码,这是我们得到的:

covA =

    0.0250    0.0075    0.0042
    0.0075    0.0070    0.0034
    0.0042    0.0034    0.0026

您会看到这也与covMATLAB中的函数匹配:

>> cov(A)

ans =

0.0250    0.0075    0.0042
0.0075    0.0070    0.0034
0.0042    0.0034    0.0026

有点暗示

如果你输入edit cov你的MATLAB命令提示符下,你可以看到他们是如何计算没有任何的协方差矩阵for循环....这基本上是我给你同样的答案:)

如果你想更有效地做到这一点

假设你可以使用sum和bsxfun,我们可以用更少(更有效率)的代码行来做到这一点.首先,像我们上面使用的那样计算你的平均向量sum:

mu = sum(A) / size(A,1);

现在,要A使用每列的相应均值减去矩阵,您可以使用它bsxfun来帮助您减少此减法:

A_mean_subtract = bsxfun(@minus, A, mu);

现在,像以前一样计算你的协方差矩阵:

covA = (A_mean_subtract.' * A_mean_subtract) / (size(A,1) - 1);

您应该得到与我们之前看到的完全相同的结果.

关于稳定性的小注

我们使用定义来计算两列之间的协方差的直接定义.但是,已经证明,如果提供某些类型的数据,使用直线定义可能会导致数值不稳定.请参阅此维基百科页面,该页面通过各种算法计算两个n更稳定的长度向量之间的协方差.