协方差矩阵的主轴与其角度不对齐

我试图获得协方差的主轴(梯度和截距)。我正在使用排序的特征向量来计算椭圆的角度,但是当我针对主轴绘制生成的椭圆时,它们没有对齐。


有没有人有任何想法?


import numpy as np

import matplotlib.pyplot as plt

from matplotlib.patches import Ellipse


def eigsorted(cov):

    vals, vecs = np.linalg.eigh(cov)

    order = vals.argsort()[::-1]

    return vals[order], vecs[:,order]


def orientation_from_covariance(cov, sigma):

    vals, vecs = eigsorted(cov)

    theta = np.degrees(np.arctan2(*vecs[:,0][::-1]))

    w, h = 2 * sigma * np.sqrt(vals)

    return w, h, theta


def plot_ellipse(ax, mu, covariance, color, linewidth=2, alpha=0.5):

    x, y, angle = orientation_from_covariance(covariance, 2)

    e = Ellipse(mu, x, y, angle=angle)

    e.set_alpha(alpha)

    e.set_linewidth(linewidth)

    e.set_edgecolor(color)

    e.set_facecolor(color)

    e.set_fill(False)

    ax.add_artist(e)

    return e



from statsmodels.stats.moment_helpers import corr2cov

corr = np.eye(2)

corr[0, 1] = corr[1, 0] = 0.7

cov = corr2cov(corr, [1, 5])

mu = [1, 1]


vectors = eigsorted(cov)[1].T

gradients = [v[0] / v[1] for v in vectors]

intercepts = [mu[1] - (gradient*mu[0]) for gradient in gradients]


plt.scatter(*np.random.multivariate_normal(mu, cov, size=9000).T, s=1);

plot_ellipse(plt.gca(), mu, cov, 'k')

_x = np.linspace(*plt.xlim())

for i,g in zip(intercepts, gradients):

    plt.plot(_x, i + (_x * g), 'k');

http://img2.mukewang.com/618a67850001ddde03750246.jpg

四季花海
浏览 258回答 1
1回答

蛊毒传说

问题是以下行# gradients = [v[0] / v[1] for v in vectors] # wrong gradients = [v[1] / v[0] for v in vectors] # correct因为梯度是 的变化y超过 的变化x。该图看起来像这样。我还在plt.figure()绘图开始之前和通话plt.axis("equal")之后添加了plot_ellipse。我还想引用numpy.linalg.eigh文档:w : (..., M) ndarray按升序排列的特征值,每个都根据其多重性重复。因此eigsorted可以省略该功能。
打开App,查看更多内容
随时随地看视频慕课网APP

相关分类

Python