PCA احتمالي

احتمالي تحليل المكونات الرئيسية (PCA) هي تقنية للحد من الأبعاد التي بتحليل بيانات عبر الفضاء كامنة أقل الأبعاد ( البقشيش والمطران 1999 ). غالبًا ما يتم استخدامه عند وجود قيم مفقودة في البيانات أو للقياس متعدد الأبعاد.

الواردات

import functools
import warnings

import matplotlib.pyplot as plt
import numpy as np
import seaborn as sns

import tensorflow.compat.v2 as tf
import tensorflow_probability as tfp

from tensorflow_probability import bijectors as tfb
from tensorflow_probability import distributions as tfd

tf.enable_v2_behavior()

plt.style.use("ggplot")
warnings.filterwarnings('ignore')

الموديل

النظر في مجموعة البيانات \(\mathbf{X} = \{\mathbf{x}_n\}\) من \(N\) نقاط البيانات، حيث كل نقطة البيانات هي \(D\)الأبعاد، $ \ mathbf {س} _n \ في \ mathbb {R} ^ D\(. We aim to represent each \)\ mathbf {س} _n $ تحت متغير الكامنة \(\mathbf{z}_n \in \mathbb{R}^K\) مع انخفاض البعد، $ K <D\(. The set of principal axes \)\ mathbf {W} $ يتعلق المتغيرات الكامنة للبيانات.

على وجه التحديد ، نفترض أن كل متغير كامن يتم توزيعه بشكل طبيعي ،

\[ \begin{equation*} \mathbf{z}_n \sim N(\mathbf{0}, \mathbf{I}). \end{equation*} \]

يتم إنشاء نقطة البيانات المقابلة عبر الإسقاط ،

\[ \begin{equation*} \mathbf{x}_n \mid \mathbf{z}_n \sim N(\mathbf{W}\mathbf{z}_n, \sigma^2\mathbf{I}), \end{equation*} \]

حيث مصفوفة \(\mathbf{W}\in\mathbb{R}^{D\times K}\) يعرف المحاور الرئيسية. في PCA احتمالي، نحن عادة المهتمة في تقدير المحاور الرئيسية \(\mathbf{W}\) وعلى المدى الضوضاء\(\sigma^2\).

يعمم PCA الاحتمالي PCA الكلاسيكي. بتهميش المتغير الكامن ، يكون توزيع كل نقطة بيانات

\[ \begin{equation*} \mathbf{x}_n \sim N(\mathbf{0}, \mathbf{W}\mathbf{W}^\top + \sigma^2\mathbf{I}). \end{equation*} \]

PCA الكلاسيكي هو الحال محددة من PCA احتمالي عندما يصبح التغاير من الضوضاء المتناهي الصغر، \(\sigma^2 \to 0\).

أنشأنا نموذجنا أدناه. في تحليلنا، ونحن نفترض \(\sigma\) هو معروف، وبدلا من نقطة تقدير \(\mathbf{W}\) كمعلمة النموذج، نضع قبل أكثر من ذلك من أجل استنتاج التوزيع على المحاور الرئيسية. سنقوم التعبير عن نموذج باعتباره TFP JointDistribution، وتحديدا، سنستخدم JointDistributionCoroutineAutoBatched .

def probabilistic_pca(data_dim, latent_dim, num_datapoints, stddv_datapoints):
  w = yield tfd.Normal(loc=tf.zeros([data_dim, latent_dim]),
                 scale=2.0 * tf.ones([data_dim, latent_dim]),
                 name="w")
  z = yield tfd.Normal(loc=tf.zeros([latent_dim, num_datapoints]),
                 scale=tf.ones([latent_dim, num_datapoints]),
                 name="z")
  x = yield tfd.Normal(loc=tf.matmul(w, z),
                       scale=stddv_datapoints,
                       name="x")

num_datapoints = 5000
data_dim = 2
latent_dim = 1
stddv_datapoints = 0.5

concrete_ppca_model = functools.partial(probabilistic_pca,
    data_dim=data_dim,
    latent_dim=latent_dim,
    num_datapoints=num_datapoints,
    stddv_datapoints=stddv_datapoints)

model = tfd.JointDistributionCoroutineAutoBatched(concrete_ppca_model)

البيانات

يمكننا استخدام النموذج لتوليد البيانات عن طريق أخذ عينات من التوزيع المسبق المشترك.

actual_w, actual_z, x_train = model.sample()

print("Principal axes:")
print(actual_w)

Principal axes:
tf.Tensor(
[[ 2.2801023]
 [-1.1619819]], shape=(2, 1), dtype=float32)

نحن نتخيل مجموعة البيانات.

plt.scatter(x_train[0, :], x_train[1, :], color='blue', alpha=0.1)
plt.axis([-20, 20, -20, 20])
plt.title("Data set")
plt.show()

الحد الأقصى للاستدلال اللاحق

نبحث أولاً عن تقدير النقاط للمتغيرات الكامنة التي تزيد من كثافة الاحتمال اللاحق. ويعرف هذا أقصى البعدية (MAP) الاستدلال، ويتم ذلك عن طريق احتساب قيم \(\mathbf{W}\) و \(\mathbf{Z}\) لتعظيم الخلفي كثافة \(p(\mathbf{W}, \mathbf{Z} \mid \mathbf{X}) \propto p(\mathbf{W}, \mathbf{Z}, \mathbf{X})\).

w = tf.Variable(tf.random.normal([data_dim, latent_dim]))
z = tf.Variable(tf.random.normal([latent_dim, num_datapoints]))

target_log_prob_fn = lambda w, z: model.log_prob((w, z, x_train))
losses = tfp.math.minimize(
    lambda: -target_log_prob_fn(w, z),
    optimizer=tf.optimizers.Adam(learning_rate=0.05),
    num_steps=200)

plt.plot(losses)

[<matplotlib.lines.Line2D at 0x7f19897a42e8>]

يمكننا استخدام نموذج لبيانات عينة للقيم الاستدلال ل \(\mathbf{W}\) و \(\mathbf{Z}\)، ومقارنتها مع بيانات الفعلي نحن مشروطة.

print("MAP-estimated axes:")
print(w)

_, _, x_generated = model.sample(value=(w, z, None))

plt.scatter(x_train[0, :], x_train[1, :], color='blue', alpha=0.1, label='Actual data')
plt.scatter(x_generated[0, :], x_generated[1, :], color='red', alpha=0.1, label='Simulated data (MAP)')
plt.legend()
plt.axis([-20, 20, -20, 20])
plt.show()

MAP-estimated axes:
<tf.Variable 'Variable:0' shape=(2, 1) dtype=float32, numpy=
array([[ 2.9135954],
       [-1.4826864]], dtype=float32)>

الاستدلال المتغير

يمكن استخدام MAP للعثور على الوضع (أو أحد الأوضاع) للتوزيع اللاحق ، لكنه لا يوفر أي رؤى أخرى حوله. نحن بجانب استخدام الاستدلال التغاير، حيث الخلفي distribtion \(p(\mathbf{W}, \mathbf{Z} \mid \mathbf{X})\) ويقترب باستخدام توزيع التغاير \(q(\mathbf{W}, \mathbf{Z})\) parametrised التي كتبها \(\boldsymbol{\lambda}\). والهدف من ذلك هو العثور على المعلمات التغييرية \(\boldsymbol{\lambda}\) التي تقلل من الاختلاف بين KL ف والخلفي، \(\mathrm{KL}(q(\mathbf{W}, \mathbf{Z}) \mid\mid p(\mathbf{W}, \mathbf{Z} \mid \mathbf{X}))\)، أو مكافئ، التي تزيد من الأدلة ملزمة أقل، \(\mathbb{E}_{q(\mathbf{W},\mathbf{Z};\boldsymbol{\lambda})}\left[ \log p(\mathbf{W},\mathbf{Z},\mathbf{X}) - \log q(\mathbf{W},\mathbf{Z}; \boldsymbol{\lambda}) \right]\).

qw_mean = tf.Variable(tf.random.normal([data_dim, latent_dim]))
qz_mean = tf.Variable(tf.random.normal([latent_dim, num_datapoints]))
qw_stddv = tfp.util.TransformedVariable(1e-4 * tf.ones([data_dim, latent_dim]),
                                        bijector=tfb.Softplus())
qz_stddv = tfp.util.TransformedVariable(
    1e-4 * tf.ones([latent_dim, num_datapoints]),
    bijector=tfb.Softplus())
def factored_normal_variational_model():
  qw = yield tfd.Normal(loc=qw_mean, scale=qw_stddv, name="qw")
  qz = yield tfd.Normal(loc=qz_mean, scale=qz_stddv, name="qz")

surrogate_posterior = tfd.JointDistributionCoroutineAutoBatched(
    factored_normal_variational_model)

losses = tfp.vi.fit_surrogate_posterior(
    target_log_prob_fn,
    surrogate_posterior=surrogate_posterior,
    optimizer=tf.optimizers.Adam(learning_rate=0.05),
    num_steps=200)

print("Inferred axes:")
print(qw_mean)
print("Standard Deviation:")
print(qw_stddv)

plt.plot(losses)
plt.show()

Inferred axes:
<tf.Variable 'Variable:0' shape=(2, 1) dtype=float32, numpy=
array([[ 2.4168603],
       [-1.2236133]], dtype=float32)>
Standard Deviation:
<TransformedVariable: dtype=float32, shape=[2, 1], fn="softplus", numpy=
array([[0.0042499 ],
       [0.00598824]], dtype=float32)>

posterior_samples = surrogate_posterior.sample(50)
_, _, x_generated = model.sample(value=(posterior_samples))

# It's a pain to plot all 5000 points for each of our 50 posterior samples, so
# let's subsample to get the gist of the distribution.
x_generated = tf.reshape(tf.transpose(x_generated, [1, 0, 2]), (2, -1))[:, ::47]

plt.scatter(x_train[0, :], x_train[1, :], color='blue', alpha=0.1, label='Actual data')
plt.scatter(x_generated[0, :], x_generated[1, :], color='red', alpha=0.1, label='Simulated data (VI)')
plt.legend()
plt.axis([-20, 20, -20, 20])
plt.show()

شكر وتقدير

وقد كتب هذا البرنامج التعليمي أصلا في إدوارد 1.0 ( مصدر ). نشكر جميع المساهمين في كتابة ومراجعة هذا الإصدار.

مراجع

[1]: مايكل إي تيبينج وكريستوفر إم بيشوب. تحليل المكون الرئيسي الاحتمالي. مجلة الجمعية الإحصائية الملكية: سلسلة B (الإحصائية المنهجية)، 61 (3): 611-622، 1999.