Analyse der gewichteten Hauptkomponenten

Nach einiger Suche finde ich sehr wenig über die Einbeziehung von Beobachtungsgewichten / Messfehlern in die Hauptkomponentenanalyse. Was ich finde, stützt sich tendenziell auf iterative Ansätze zur Einbeziehung von Gewichtungen (z . B. hier ). Meine Frage ist, warum dieser Ansatz notwendig ist? Warum können wir die Eigenvektoren der gewichteten Kovarianzmatrix nicht verwenden?

Es kommt darauf an, worauf genau sich Ihre Gewichte beziehen.

Reihengewichte

Sei die Datenmatrix mit Variablen in Spalten und Beobachtungen in Zeilen. Wenn jeder Beobachtung ein Gewicht , ist es in der Tat einfach, diese Gewichte in die PCA aufzunehmen.$\mathbf{X}$$n$$\mathbf x_i$$w_i$

Zunächst braucht man zur Berechnung des gewichteten Mittels und subtrahieren sie aus den Daten , um zu zentrieren es.$\boldsymbol \mu = \frac{1}{\sum w_i}\sum w_i \mathbf x_i$

Dann berechnen wir die gewichtete Kovarianzmatrix $\frac{1}{\sum w_i}\mathbf X^\top \mathbf W \mathbf X$ , wobei $\mathbf W = \operatorname{diag}(w_i)$ die diagonale Matrix von Gewichten ist, und wenden Sie Standard-PCA an, um es zu analysieren.

Zellgewichte

Die Arbeit von Tamuz et al., 2013 , die Sie gefunden haben, betrachtet einen komplizierteren Fall, wenn auf jedes Element der Datenmatrix unterschiedliche Gewichte angewendet werden. Dann gibt es in der Tat keine analytische Lösung und man muss eine iterative Methode anwenden. Beachten Sie, dass sie, wie von den Autoren anerkannt, das Rad neu erfunden haben, da solche allgemeinen Gewichte sicherlich zuvor in Betracht gezogen wurden, z. B. in Gabriel und Zamir, 1979, Näherung von Matrizen mit niedrigerem Rang durch kleinste Quadrate mit einer beliebigen Auswahl an Gewichten . Dies wurde auch hier diskutiert .$w_{ij}$

Als zusätzliche Bemerkung: Wenn die Gewichte sowohl mit Variablen als auch mit Beobachtungen variieren, aber symmetrisch sind, so dass , dann ist wieder eine analytische Lösung möglich, siehe Koren und Carmel, 2004, Robust Lineare Dimensionsreduktion .$w_{ij}$$w_{ij}=w_{ji}$

Vielen Dank, Amöbe, für den Einblick in die Zeilengewichte. Ich weiß, dass dies kein Stapelüberlauf ist, aber ich hatte einige Schwierigkeiten, eine Implementierung von zeilengewichteten PCA mit Erklärung zu finden, und da dies eines der ersten Ergebnisse ist, wenn ich nach gewichteten PCA google, dachte ich, es wäre gut, meine Lösung beizufügen Vielleicht kann es anderen in der gleichen Situation helfen. In diesem Python2-Code-Snippet wird eine PCA, die mit einem RBF-Kernel wie oben beschrieben gewichtet ist, verwendet, um die Tangenten eines 2D-Datensatzes zu berechnen. Über Feedback würde ich mich sehr freuen!

def weighted_pca_regression(x_vec, y_vec, weights):
    """
    Given three real-valued vectors of same length, corresponding to the coordinates
    and weight of a 2-dimensional dataset, this function outputs the angle in radians
    of the line that aligns with the (weighted) average and main linear component of
    the data. For that, first a weighted mean and covariance matrix are computed.
    Then u,e,v=svd(cov) is performed, and u * f(x)=0 is solved.
    """
    input_mat = np.stack([x_vec, y_vec])
    weights_sum = weights.sum()
    # Subtract (weighted) mean and compute (weighted) covariance matrix:
    mean_x, mean_y =  weights.dot(x_vec)/weights_sum, weights.dot(y_vec)/weights_sum
    centered_x, centered_y = x_vec-mean_x, y_vec-mean_y
    matrix_centered = np.stack([centered_x, centered_y])
    weighted_cov = matrix_centered.dot(np.diag(weights).dot(matrix_centered.T)) / weights_sum
    # We know that v rotates the data's main component onto the y=0 axis, and
    # that u rotates it back. Solving u.dot([x,0])=[x*u[0,0], x*u[1,0]] gives
    # f(x)=(u[1,0]/u[0,0])x as the reconstructed function.
    u,e,v = np.linalg.svd(weighted_cov)
    return np.arctan2(u[1,0], u[0,0]) # arctan more stable than dividing


# USAGE EXAMPLE:
# Define the kernel and make an ellipse to perform regression on:
rbf = lambda vec, stddev: np.exp(-0.5*np.power(vec/stddev, 2))
x_span = np.linspace(0, 2*np.pi, 31)+0.1
data_x = np.cos(x_span)[:-1]*20-1000
data_y = np.sin(x_span)[:-1]*10+5000
data_xy = np.stack([data_x, data_y])
stddev = 1 # a stddev of 1 in this context is highly local
for center in data_xy.T:
    # weight the  points based on their euclidean distance to the current center
    euclidean_distances = np.linalg.norm(data_xy.T-center, axis=1)
    weights = rbf(euclidean_distances, stddev)
    # get the angle for the regression in radians
    p_grad = weighted_pca_regression(data_x, data_y, weights)
    # plot for illustration purposes
    line_x = np.linspace(-5,5,10)
    line_y = np.tan(p_grad)*line_x
    plt.plot(line_x+center[0], line_y+center[1], c="r")
    plt.scatter(*data_xy)
    plt.show()

Und eine Beispielausgabe (für jeden Punkt gilt das Gleiche):

Prost,
Andres