Python OpenCV PCA计算特征值

Question

在使用Python 2.7.5和OpenCV（在OSX系统上）时，我对一系列图像进行了主成分分析（PCA），这里的列是像素，行是帧，具体可以参考这个回答。

我想知道如何获取与特征向量对应的特征值？看起来在C++的PCA对象中有这个属性，但在Python中对应的PCACompute()只是一个简单的函数。

感觉省略了PCA中这么重要的部分有点奇怪。

Answer 1

matmul.cpp 这个文件确认了 PCA::Operator() 是被 PCACompute() 使用的，但特征值被丢弃了。所以我做了这个：

# The following mimics PCA::operator() implementation from OpenCV's
# matmul.cpp() which is wrapped by Python cv2.PCACompute(). We can't
# use PCACompute() though as it discards the eigenvalues.

# Scrambled is faster for nVariables >> nObservations. Bitmask is 0 and
# therefore default / redundant, but included to abide by online docs.
covar, mean = cv2.calcCovarMatrix(PCAInput, cv2.cv.CV_COVAR_SCALE |
                                            cv2.cv.CV_COVAR_ROWS  |
                                            cv2.cv.CV_COVAR_SCRAMBLED)

eVal, eVec = cv2.eigen(covar, computeEigenvectors=True)[1:]

# Conversion + normalisation required due to 'scrambled' mode
eVec = cv2.gemm(eVec, PCAInput - mean, 1, None, 0)
# apply_along_axis() slices 1D rows, but normalize() returns 4x1 vectors
eVec = numpy.apply_along_axis(lambda n: cv2.normalize(n).flat, 1, eVec)

(简单假设：行代表观察值，列代表变量；而且变量的数量远多于观察值。在我的情况下，这两点都是成立的。)

这个方法 基本上 是有效的。在下面的内容中，old_eVec 是来自 cv2.PCACompute() 的结果：

In [101]: eVec
Out[101]: 
array([[  3.69396088e-05,   1.66745325e-05,   4.97117583e-05, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       [ -7.23531536e-06,  -3.07411122e-06,  -9.58259793e-06, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       [  1.01496237e-05,   4.60048715e-06,   1.33919606e-05, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       ..., 
       [ -1.42024751e-04,   5.21386198e-05,   3.59923394e-04, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       [ -5.28685812e-05,   8.50139472e-05,  -3.13278542e-04, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       [  2.96546917e-04,   1.23437674e-04,   4.98598461e-04, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00]])

In [102]: old_eVec
Out[102]: 
array([[  3.69395821e-05,   1.66745194e-05,   4.97117981e-05, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       [ -7.23533140e-06,  -3.07411415e-06,  -9.58260534e-06, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       [  1.01496662e-05,   4.60050160e-06,   1.33920075e-05, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       ..., 
       [ -1.42029530e-04,   5.21366564e-05,   3.60067672e-04, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       [ -5.29163444e-05,   8.50261567e-05,  -3.13150231e-04, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00],
       [ -7.13724992e-04,  -8.52700090e-04,   1.57953508e-03, ...,
          0.00000000e+00,   0.00000000e+00,   0.00000000e+00]], dtype=float32)

在输出的最后部分，有一些精度损失是可见的（虽然实际上快速绘制绝对差异并没有显示出精度不高的规律）。

57%的元素有非零的绝对差异。
在这些元素中，95%的差异小于 2e-16，平均绝对差异是 5.3e-4 - 不过，绝对差异最大可以达到 0.059，这在考虑到所有特征向量值都在 -0.048 到 0.045 之间时，算是比较大的。

在 PCA::Operator() 中有代码会转换为最大的 ctype；另一方面，old_eVec 是 float32，而我自己的代码生成的是 float64。值得一提的是，当我编译 numpy 时，遇到了一些与精度相关的错误。

总体来看，精度损失似乎与低特征值的特征向量有关，这又指向了舍入误差等问题。上述实现的结果与 PCACompute() 类似，成功复制了其行为。