Python | R 潜在混合模型

PreviousPython社会经济 | 怀特的异方差一致估计量 NextPython | MATLAB | R 心理认知数学图形模型推断

Last updated 10 months ago

Was this helpful?

Python | R 潜在混合模型

关键词

Python | R | 概率 | 分布 | 统计 | 混合模型 | 高斯 | 数学

📜用例

📜Python | MATLAB | R 心理认知数学图形模型推断 | 📜信用卡消费高斯混合模型 | 📜必修课学业成绩分布异常背景混合模型潜在类别分析

✒️潜在混合模型

本质上，混合模型（或混合分布）是将多个概率分布组合成一个概率分布。

\begin{gathered} P(x)=\pi_0 p_0(x)+\pi_1 p_1(x)+\ldots+\pi_i p_i(x) \\ \text { s.t. } \sum \pi_i=1 \end{gathered}

为了将这些分布组合在一起，我们为每个成分分布分配一个权重，使得该分布下的总概率总和为 1。一个简单的例子是包含 2 个高斯分布的混合分布。我们可以有 2 个具有不同均值和方差的分布，并使用不同的权重将这 2 个分布组合在一起。

具体来说，我们可以认为该分布源自一个两步生成过程。在此过程中，可以从 n 个不同的概率分布中生成一个数据点。首先，我们确定它来自哪个概率分布。这个概率就是权重 π_i。一旦选择了组件概率分布，就可以通过模拟组件概率分布本身来生成数据点。

高斯混合模型本质上是一种混合模型，其中所有分量分布都是高斯分布。

\begin{gathered} f(x)=\pi_0 N\left(\mu_0, \Sigma_0\right)+\pi_1 N\left(\mu_1, \Sigma_1\right)+\ldots+\pi_i N\left(\mu_i, \Sigma_i\right) \\ \text { s.t. } \sum \pi_i=1 \end{gathered}

现在让我们试着理解为什么使用高斯分布来对混合物的成分进行建模。当查看数据集时，我们希望将相似的点聚类在一起。这些聚类通常本质上是球形或椭圆形的，因为我们希望将靠近的点聚类在一起。因此，正态分布是集群的良好模型。分布的均值将是簇的中心，而簇的形状和分布可以通过分布的协方差很好地建模。

集群的第二个变量是不同集群的相对大小。在有机数据集中，我们通常不期望集群的大小相同，这意味着某些集群的点数会比其他集群多。然后，集群的大小将由集群权重 $π_i$ 决定。

在聚类的背景下，我们假设有 k 个影响因素影响数据的生成。每个影响因素都有不同的权重，对应于簇权重π。

💦Python高斯混合模型

基于概率的软聚类方法
每个簇：一个生成模型（高斯或多项式）
参数（例如平均值/协方差未知）

让我们生成一个示例数据集，其中点是从两个高斯过程之一生成的。第一个分布的平均值为 100，第二个分布的平均值为 90；和分布的标准差分别为 5 和 2。

第一个过程我们将获得60,000积分；第二个过程中50,000个点并将它们混合在一起。

 import numpy as np
 np.random.seed(0)
 
 import matplotlib.pyplot as plt
 import seaborn as sns
 import warnings
 warnings.simplefilter(action='ignore', category=FutureWarning)

 Mean1, Mean2  = 100.0, 90.0
 Standard_dev1, Standard_dev2 = 5.0, 2.0
 sample1, sample2 = 60000, 50000
 
 print('Input Normal_distb {:}: μ = {:1}, σ = {:.2}, n = {} '.format("1", Mean1, Standard_dev1, sample1))
 print('Input Normal_distb {:}: μ = {:1}, σ = {:.2}, n = {} '.format("2", Mean2, Standard_dev2, sample2))

 X1 = np.random.normal(loc = Mean1, scale = Standard_dev1, size = sample1)
 X2 = np.random.normal(loc = Mean2, scale = Standard_dev2, size = sample2)
 
 # plot
 fig = plt.figure(figsize=(8, 6), dpi=100)
 sns.distplot(X1, bins=50, kde=True, norm_hist=True, label='Normal distribution 1')
 sns.distplot(X2, bins=50, kde=True, norm_hist=True, label='Normal distribution 2')
 plt.legend()
 plt.show()
 # save
 fig.savefig('norm_distrib12.jpeg')

 
 # mix two distrib together
 X = np.hstack((X1, X2))
 
 # plot
 fig = plt.figure(figsize=(8, 6), dpi=100)
 sns.distplot(X, bins=50, kde=True, norm_hist=True, label='gaussian mixture')
 plt.legend()
 plt.show()
 # save
 fig.savefig('final_dt.jpeg')

因此，将这些过程混合在一起后，我们就得到了上图所示的数据集。我们可以注意到 2 个峰值：大约 90 和 100，但对于峰值中间的许多点，很难确定它们来自哪个分布。那么我们如何解决这个任务呢？我们可以使用高斯混合模型，该模型将使用期望最大化算法来估计分布的参数。

💦潜在变量最大似然估计（算法）：

 X = X.reshape((len(X), 1))

 from sklearn.mixture import GaussianMixture
 
 GMM = GaussianMixture(n_components = 2, init_params = 'random')
 GMM.fit(X)
 print('Converged:', GMM.converged_) # check if the model has converged

 Y = np.array([[105.0]])
 prediction = GMM.predict_proba(Y)
 print('Probability each Gaussian (state) in the model given each sample p = {}'.format(prediction))
 print()
 yhat = GMM.predict(X)
 
 print(yhat[:100])
 print(yhat[-100:])

 print(len(yhat[yhat==0]))
 print(len(yhat[yhat==1]))

多元高斯：d > 1

 from sklearn.datasets.samples_generator import make_blobs
 from scipy.stats import multivariate_normal
 
 X, Y = make_blobs(cluster_std=1.0, random_state=123, n_samples=12000, centers=3)
 
 X = np.dot(X, np.random.RandomState(0).randn(2,2))
 
 x, y = np.meshgrid(np.sort(X[:,0]), np.sort(X[:,1]))
 XY = np.array([x.flatten(), y.flatten()]).T
 
 GMM = GaussianMixture(n_components=3).fit(X) # instantiate and fit the model
 print('Converged:', GMM.converged_) # check if the model has converged
 means = GMM.means_ 
 covariances = GMM.covariances_
 
 Y = np.array([[0.5], [0.5]])
 prediction = GMM.predict_proba(Y.T)
 print('Probability each Gaussian (state) in the model given each sample p = {}'.format(prediction))
  
 fig = plt.figure(figsize = (12,12), dpi = 100)
 ax0 = fig.add_subplot(111)
 ax0.scatter(X[:,0], X[:,1])
 ax0.scatter(Y[0,:], Y[1,:], c = 'orange', zorder = 10, s = 100)
 for m,c in zip(means,covariances):
     multi_normal = multivariate_normal(mean = m, cov = c)
     ax0.contour(np.sort(X[:,0]), np.sort(X[:,1]), multi_normal.pdf(XY).reshape(len(X), len(X)), colors='black', alpha=0.3)
     ax0.scatter(m[0], m[1], c = 'grey', zorder = 10, s = 100)
     
 plt.show()
 
 fig.savefig('2d.jpeg')

PreviousPython社会经济 | 怀特的异方差一致估计量 NextPython | MATLAB | R 心理认知数学图形模型推断

Last updated 10 months ago

Was this helpful?

关键词

Python | R | 概率 | 分布 | 统计 | 混合模型 | 高斯 | 数学

🏈指点迷津 | Brief

📜用例

📜Python | MATLAB | R 心理认知数学图形模型推断 | 📜信用卡消费高斯混合模型 | 📜必修课学业成绩分布异常背景混合模型潜在类别分析

✒️潜在混合模型

本质上，混合模型（或混合分布）是将多个概率分布组合成一个概率分布。

\begin{gathered} P(x)=\pi_0 p_0(x)+\pi_1 p_1(x)+\ldots+\pi_i p_i(x) \\ \text { s.t. } \sum \pi_i=1 \end{gathered}

高斯混合模型本质上是一种混合模型，其中所有分量分布都是高斯分布。

\begin{gathered} f(x)=\pi_0 N\left(\mu_0, \Sigma_0\right)+\pi_1 N\left(\mu_1, \Sigma_1\right)+\ldots+\pi_i N\left(\mu_i, \Sigma_i\right) \\ \text { s.t. } \sum \pi_i=1 \end{gathered}

在聚类的背景下，我们假设有 k 个影响因素影响数据的生成。每个影响因素都有不同的权重，对应于簇权重π。

💦Python高斯混合模型

基于概率的软聚类方法
每个簇：一个生成模型（高斯或多项式）
参数（例如平均值/协方差未知）

第一个过程我们将获得60,000积分；第二个过程中50,000个点并将它们混合在一起。

 import numpy as np
 np.random.seed(0)
 
 import matplotlib.pyplot as plt
 import seaborn as sns
 import warnings
 warnings.simplefilter(action='ignore', category=FutureWarning)

 Mean1, Mean2  = 100.0, 90.0
 Standard_dev1, Standard_dev2 = 5.0, 2.0
 sample1, sample2 = 60000, 50000
 
 print('Input Normal_distb {:}: μ = {:1}, σ = {:.2}, n = {} '.format("1", Mean1, Standard_dev1, sample1))
 print('Input Normal_distb {:}: μ = {:1}, σ = {:.2}, n = {} '.format("2", Mean2, Standard_dev2, sample2))

 X1 = np.random.normal(loc = Mean1, scale = Standard_dev1, size = sample1)
 X2 = np.random.normal(loc = Mean2, scale = Standard_dev2, size = sample2)
 
 # plot
 fig = plt.figure(figsize=(8, 6), dpi=100)
 sns.distplot(X1, bins=50, kde=True, norm_hist=True, label='Normal distribution 1')
 sns.distplot(X2, bins=50, kde=True, norm_hist=True, label='Normal distribution 2')
 plt.legend()
 plt.show()
 # save
 fig.savefig('norm_distrib12.jpeg')

 
 # mix two distrib together
 X = np.hstack((X1, X2))
 
 # plot
 fig = plt.figure(figsize=(8, 6), dpi=100)
 sns.distplot(X, bins=50, kde=True, norm_hist=True, label='gaussian mixture')
 plt.legend()
 plt.show()
 # save
 fig.savefig('final_dt.jpeg')

💦潜在变量最大似然估计（算法）：

 X = X.reshape((len(X), 1))

 from sklearn.mixture import GaussianMixture
 
 GMM = GaussianMixture(n_components = 2, init_params = 'random')
 GMM.fit(X)
 print('Converged:', GMM.converged_) # check if the model has converged

 Y = np.array([[105.0]])
 prediction = GMM.predict_proba(Y)
 print('Probability each Gaussian (state) in the model given each sample p = {}'.format(prediction))
 print()
 yhat = GMM.predict(X)
 
 print(yhat[:100])
 print(yhat[-100:])

 print(len(yhat[yhat==0]))
 print(len(yhat[yhat==1]))

多元高斯：d > 1

 from sklearn.datasets.samples_generator import make_blobs
 from scipy.stats import multivariate_normal
 
 X, Y = make_blobs(cluster_std=1.0, random_state=123, n_samples=12000, centers=3)
 
 X = np.dot(X, np.random.RandomState(0).randn(2,2))
 
 x, y = np.meshgrid(np.sort(X[:,0]), np.sort(X[:,1]))
 XY = np.array([x.flatten(), y.flatten()]).T
 
 GMM = GaussianMixture(n_components=3).fit(X) # instantiate and fit the model
 print('Converged:', GMM.converged_) # check if the model has converged
 means = GMM.means_ 
 covariances = GMM.covariances_
 
 Y = np.array([[0.5], [0.5]])
 prediction = GMM.predict_proba(Y.T)
 print('Probability each Gaussian (state) in the model given each sample p = {}'.format(prediction))
  
 fig = plt.figure(figsize = (12,12), dpi = 100)
 ax0 = fig.add_subplot(111)
 ax0.scatter(X[:,0], X[:,1])
 ax0.scatter(Y[0,:], Y[1,:], c = 'orange', zorder = 10, s = 100)
 for m,c in zip(means,covariances):
     multi_normal = multivariate_normal(mean = m, cov = c)
     ax0.contour(np.sort(X[:,0]), np.sort(X[:,1]), multi_normal.pdf(XY).reshape(len(X), len(X)), colors='black', alpha=0.3)
     ax0.scatter(m[0], m[1], c = 'grey', zorder = 10, s = 100)
     
 plt.show()
 
 fig.savefig('2d.jpeg')