Принципиальные оси ковариационной матрицы не совпадают с ее углом - PullRequest
Новая
       32

Принципиальные оси ковариационной матрицы не совпадают с ее углом

0 голосов
/ 15 февраля 2019

Я пытаюсь получить большую ось ковариации (градиент и перехват).Я использую отсортированные собственные векторы для вычисления угла эллипса, но когда я строю получившийся эллипс относительно большой оси, они не совпадают.

У кого-нибудь есть идеи? 

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');

covariance

1 Ответ

0 голосов
/ 15 февраля 2019

Проблема была в следующей строке 

# gradients = [v[0] / v[1] for v in vectors] # wrong
gradients = [v[1] / v[0] for v in vectors] # correct

, потому что градиент - это изменение y по сравнению с x.Затем фигура выглядит следующим образом.Ellipse Plot

Я также добавил plt.figure() до начала построения графика и plt.axis("equal") после вызова plot_ellipse.

Я также хотел бы привести numpy.linalg.eigh документация:

w: (…, M) ndarray

Собственные значения в порядке возрастания, каждое повторяется в соответствии со своей кратностью.

и, следовательно, функция eigsorted может быть опущена.

...