我可以从奇异值分解中得到数据传播(噪声)吗?

Can I get data spread (noise) from singular value decomposition?

我希望使用奇异值分解来估计椭圆数据的标准差。我不确定这是否是最好的方法,我可能对整个过程想得太多了,所以我需要一些帮助。

我使用以下脚本模拟了一些数据...

from matplotlib import pyplot as plt
import numpy


def svd_example():
    # simulate some data...
    # x values have standard deviation 3000
    xdata = numpy.random.normal(0, 3000, 5000).reshape(-1, 1)
    # y values standard deviation 300
    ydata = numpy.random.normal(0, 300, 5000).reshape(-1, 1)
    # apply some rotation
    ydata_rotated = ydata + (xdata * 0.5)
    data = numpy.hstack((xdata, ydata_rotated))

    # get singular values
    left_singular_matrix, singular_values, right_singular_matrix = numpy.linalg.svd(data)
    print 'singular values', singular_values

    # plot data....
    plt.scatter(data[:, 0], data[:, 1], s=5)
    plt.ylim(-15000, 15000)
    plt.show()

svd_example()

我得到...的奇异值...

>>> singular values [ 234001.71228678   18850.45155942]

我的数据是这样的...

我假设奇异值会给我一些数据传播的指示,而不管它是旋转的,对吗?但是这些值 [234001.71228678 18850.45155942] 对我来说毫无意义。我的标准偏差是 3000 和 300。这些奇异值是否代表方差?我如何转换它们?

奇异值确实给出了传播的一些迹象。事实上,它们与这些方向的标准偏差有关。但是,它们没有标准化。如果除以样本数的平方根,您将得到与用于创建数据的标准差非常相似的值:

singular_values / np.sqrt(5000)
# array([ 3398.61320614,   264.00975837])

为什么得到的是3400和264,而不是3000和300?那是因为ydata + (xdata * 0.5)不是旋转而是剪切操作。真正的旋转将保留原始标准偏差。

例如,以下代码会将数据旋转 40 度:

# apply some rotation
s = numpy.sin(40 * numpy.pi / 180)
c = numpy.cos(40 * numpy.pi / 180)
data = numpy.hstack((xdata, ydata)).dot([[c, s], [-s, c]])

通过这样的旋转,您将获得非常接近原始标准偏差的归一化奇异值。


编辑: 规范化

我不得不承认,规范化可能不是适用于此的正确术语。它并不一定意味着将值缩放到某个范围。归一化,正如我的意思,是将值置于定义的范围内,与样本数量无关。

要了解除以 sqrt(5000) 的来源,让我们谈谈标准偏差。设 xn 个零均值样本的数据向量。然后标准偏差计算为 sqrt(sum(x**2)/n)sqrt(sum(x**2)) / sqrt(n)。现在,你可以认为奇异值分解只计算 sqrt(sum(x**2)) 部分,所以我们必须自己除以 sqrt(n)

恐怕这不是一个非常数学化的解释,但希望它传达了这个想法。