close all
%%%%%%%%%%%%%%%%%%%%%%%%%生成实验数据集
rand('state',0)
sigma_matrix1=eye(2)
sigma_matrix2=50*eye(2)
u1=[0,0]
u2=[30,30]
m1=100
m2=300%样本数
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%sm1数据集
Y1=multivrandn(u1,m1,sigma_matrix1)
Y2=multivrandn(u2,m2,sigma_matrix2)
scatter(Y1(:,1),Y1(:,2),'bo')
hold on
scatter(Y2(:,1),Y2(:,2),'r*')
title('SM1数据集')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%sm2数据集
u11=[0,0]
u22=[5,5]
u33=[10,10]
u44=[15,15]
m=600
sigma_matrix3=2*eye(2)
Y11=multivrandn(u11,m,sigma_matrix3)
Y22=multivrandn(u22,m,sigma_matrix3)
Y33=multivrandn(u33,m,sigma_matrix3)
Y44=multivrandn(u44,m,sigma_matrix3)
figure(2)
scatter(Y11(:,1),Y11(:,2),'bo')
hold on
scatter(Y22(:,1),Y22(:,2),'r*')
scatter(Y33(:,1),Y33(:,2),'go')
scatter(Y44(:,1),Y44(:,2),'c*')
title('SM2数据集')
end
function Y = multivrandn(u,m,sigma_matrix)
%%生成指定均值和协方差矩阵的高斯数据
n=length(u)
c = chol(sigma_matrix)
X=randn(m,n)
Y=X*c+ones(m,1)*u
end
如果随机变量X的所有取值都可以逐个列举出来，则称X为离散型随机变量。相应的概率分布有二项分布，泊松分布。
如果随机变量X的所有取值无法逐个列举出来，而是取数轴上某一区间内的任一点，则称X为连续型随机变量。相应的概率分布有正态分布，均匀分布，指数分布，伽马分布，偏态分布，卡方分布，beta分布等。(真多分布，好恐怖~~)
在离散型随机变量X的一切可能值中，各可能值与其对应概率的乘积之和称为该随机变量X的期望值，记作E(X) 。比如有随机变量，取值依次为：2，2，2，4，5。求其平均值：(2+2+2+4+5)/5 = 3。
期望值也就是该随机变量总体的均值。 推导过程如下：
= (2+2+2+4+5)/5
= 1/5 2 3 + 4/5 + 5/5
= 3/5 2 + 1/5 4 + 1/5 5
= 0.6 2 + 0.2 4 + 0.2 5
= 60% 2 + 20% 4 + 20%*5
= 1.2 + 0.8 + 1
= 3
倒数第三步可以解释为值为2的数字出现的概率为60%，4的概率为20%，5的概率为20%。 所以E(X) = 60% 2 + 20% 4 + 20%*5 = μ = 3。
0-1分布（两点分布），它的随机变量的取值为1或0。即离散型随机变量X的概率分布为：P{X=0} = 1-p, P{X=1} = p，即：
则称随机变量X服从参数为p的0-1分布，记作X~B（1，p)。
在生活中有很多例子服从两点分布，比如投资是否中标，新生婴儿是男孩还是女孩，检查产品是否合格等等。
大家非常熟悉的抛硬币试验对应的分布就是二项分布。抛硬币试验要么出现正面，要么就是反面，只包含这两个结果。出现正面的次数是一个随机变量，这种随机变量所服从的概率分布通常称为 二项分布 。
像抛硬币这类试验所具有的共同性质总结如下：（以抛硬币为例）
通常称具有上述特征的n次重复独立试验为n重伯努利试验。简称伯努利试验或伯努利试验概型。特别地，当试验次数为1时，二项分布服从0-1分布(两点分布)。
举个栗子：抛3次均匀的硬币，求结果出现有2个正面的概率 。
已知p = 0.5 (出现正面的概率) ，n = 3 ，k = 2
所以抛3次均匀的硬币，求结果出现有2个正面的概率为3/8。
二项分布的期望值和方差 分别为：
泊松分布是用来描述在一 指定时间范围内或在指定的面积或体积之内某一事件出现的次数的分布 。生活中服从泊松分布的例子比如有每天房产中介接待的客户数，某微博每月出现服务器瘫痪的次数等等。 泊松分布的公式为 ：
其中 λ 为给定的时间间隔内事件的平均数，λ = np。e为一个数学常数，一个无限不循环小数，其值约为2.71828。
泊松分布的期望值和方差 分别为：
使用Python绘制泊松分布的概率分布图：
因为连续型随机变量可以取某一区间或整个实数轴上的任意一个值，所以通常用一个函数f(x)来表示连续型随机变量，而f(x)就称为 概率密度函数 。
概率密度函数f(x)具有如下性质 ：
需要注意的是，f(x)不是一个概率，即f(x) ≠ P(X = x) 。在连续分布的情况下，随机变量X在a与b之间的概率可以写成：
正态分布（或高斯分布）是连续型随机变量的最重要也是最常见的分布，比如学生的考试成绩就呈现出正态分布的特征，大部分成绩集中在某个范围（比如60-80分），很小一部分往两端倾斜（比如50分以下和90多分以上）。还有人的身高等等。
正态分布的定义 ：
如果随机变量X的概率密度为( -∞<x<+∞)：
则称X服从正态分布，记作X~N(μ,σ²)。其中-∞<μ<+∞，σ>0， μ为随机变量X的均值，σ为随机变量X的标准差。 正态分布的分布函数
正态分布的图形特点 ：
使用Python绘制正态分布的概率分布图：
正态分布有一个3σ准则，即数值分布在(μ-σ,μ+σ)中的概率为0.6827，分布在（μ-2σ,μ+2σ)中的概率为0.9545，分布在(μ-3σ,μ+3σ)中的概率为0.9973，也就是说大部分数值是分布在(μ-3σ,μ+3σ)区间内，超出这个范围的可能性很小很小，仅占不到0.3%，属于极个别的小概率事件，所以3σ准则可以用来检测异常值。
当μ=0，σ=1时，有
此时的正态分布N(0,1) 称为标准正态分布。因为μ，σ都是确定的取值，所以其对应的概率密度曲线是一条 形态固定 的曲线。
对标准正态分布，通常用φ(x)表示概率密度函数，用Φ(x)表示分布函数：
假设有一次物理考试特别难，满分100分，全班只有大概20个人及格。与此同时语文考试很简单，全班绝大部分都考了90分以上。小明的物理和语文分别考了60分和80分，他回家后告诉家长，这时家长能仅仅从两科科目的分值直接判断出这次小明的语文成绩要比物理好很多吗？如果不能，应该如何判断呢？此时Z-score就派上用场了。 Z-Score的计算定义 ：
即 将随机变量X先减去总体样本均值，再除以总体样本标准差就得到标准分数啦。如果X低于平均值，则Z为负数，反之为正数 。通过计算标准分数，可以将任何一个一般的正态分布转化为标准正态分布。
小明家长从老师那得知物理的全班平均成绩为40分，标准差为10，而语文的平均成绩为92分，标准差为4。分别计算两科成绩的标准分数：
物理：标准分数 = (60-40)/10 = 2
语文：标准分数 = (85-95)/4 = -2.5
从计算结果来看，说明这次考试小明的物理成绩在全部同学中算是考得很不错的，而语文考得很差。
指数分布可能容易和前面的泊松分布混淆，泊松分布强调的是某段时间内随机事件发生的次数的概率分布，而指数分布说的是 随机事件发生的时间间隔 的概率分布。比如一班地铁进站的间隔时间。如果随机变量X的概率密度为：
则称X服从指数分布，其中的参数λ>0。 对应的分布函数 为：
均匀分布的期望值和方差 分别为：
使用Python绘制指数分布的概率分布图：
均匀分布有两种，分为 离散型均匀分布和连续型均匀分布 。其中离散型均匀分布最常见的例子就是抛掷骰子啦。抛掷骰子出现的点数就是一个离散型随机变量，点数可能有1，2，3，4，5，6。每个数出现的概率都是1/6。
设连续型随机变量X具有概率密度函数：
则称X服从区间(a,b)上的均匀分布。X在等长度的子区间内取值的概率相同。对应的分布函数为：
f(x)和F(x)的图形分别如下图所示：
均匀分布的期望值和方差 分别为：
本文翻译自 https://jakevdp.github.io/PythonDataScienceHandbook/05.12-gaussian-mixtures.html
上一节中探讨的k-means聚类模型简单易懂，但其简单性导致其应用中存在实际挑战。具体而言，k-means的非概率特性及简单地计算点与类蔟中心的欧式距离来判定归属，会导致其在许多真实的场景中性能较差。本节，我们将探讨高斯混合模型(GMMs)，其可以看成k-means的延伸，更可以看成一个强有力的估计工具，而不仅仅是聚类。
我们将以一个标准的import开始
我们看下k-means的缺陷，思考下如何提高聚类模型。正如上一节所示，给定简单，易于分类的数据，k-means能找到合适的聚类结果。
举例而言，假设我们有些简单的数据点，k-means算法能以某种方式很快地将它们聚类，跟我们肉眼分辨的结果很接近：
从直观的角度来看，我可能期望聚类分配时，某些点比其他的更确定：举例而言，中间两个聚类之间似乎存在非常轻微的重叠，这样我们可能对这些数据点的分配没有完全的信心。不幸的是，k-means模型没有聚类分配的概率或不确定性的内在度量（尽管可能使用bootstrap 的方式来估计这种不确定性）。为此，我们必须考虑泛化这种模型。
k-means模型的一种理解思路是，它在每个类蔟的中心放置了一个圈（或者，更高维度超球面)，其半径由聚类中最远的点确定。该半径充当训练集中聚类分配的一个硬截断：任何圈外的数据点不被视为该类的成员。我们可以使用以下函数可视化这个聚类模型：
观察k-means的一个重要发现，这些聚类模式必须是圆形的。k-means没有内置的方法来计算椭圆形或椭圆形的簇。因此，举例而言，假设我们将相同的数据点作变换，这种聚类分配方式最终变得混乱：
高斯混合模型（GMM）试图找到一个多维高斯概率分布的混合，以模拟任何输入数据集。在最简单的情况下，GMM可用于以与k-means相同的方式聚类。
但因为GMM包含概率模型，因此可以找到聚类分配的概率方式 - 在Scikit-Learn中，通过调用predict_proba方法实现。它将返回一个大小为[n_samples, n_clusters]的矩阵，用于衡量每个点属于给定类别的概率：
我们可以可视化这种不确定性，比如每个点的大小与预测的确定性成比例；如下图，我们可以看到正是群集之间边界处的点反映了群集分配的不确定性：
本质上说，高斯混合模型与k-means非常相似：它使用期望-最大化的方式，定性地执行以下操作：
有了这个，我们可以看看四成分的GMM为我们的初始数据提供了什么：
同样，我们可以使用GMM方法来拟合我们的拉伸数据集；允许full的协方差，该模型甚至可以适应非常椭圆形，伸展的聚类模式：
这清楚地表明GMM解决了以前遇到的k-means的两个主要实际问题。
如果看了之前拟合的细节，你将看到covariance_type选项在每个中都设置不同。该超参数控制每个类簇的形状的自由度；对于任意给定的问题，必须仔细设置。默认值为covariance_type =“diag”，这意味着可以独立设置沿每个维度的类蔟大小，并将得到的椭圆约束为与轴对齐。一个稍微简单和快速的模型是covariance_type =“spherical”，它约束了类簇的形状，使得所有维度都相等。尽管它并不完全等效，其产生的聚类将具有与k均值相似的特征。更复杂且计算量更大的模型（特别是随着维数的增长）是使用covariance_type =“full”，这允许将每个簇建模为具有任意方向的椭圆。
对于一个类蔟，下图我们可以看到这三个选项的可视化表示：
尽管GMM通常被归类为聚类算法，但从根本上说它是一种密度估算算法。也就是说，GMM适合某些数据的结果在技术上不是聚类模型，而是描述数据分布的生成概率模型。
例如，考虑一下Scikit-Learn的make_moons函数生成的一些数据：
如果我们尝试用视为聚类模型的双成分的GMM模拟数据，则结果不是特别有用：
但是如果我们使用更多成分的GMM模型，并忽视聚类的类别，我们会发现更接近输入数据的拟合：
这里，16个高斯分布的混合不是为了找到分离的数据簇，而是为了对输入数据的整体分布进行建模。这是分布的一个生成模型，这意味着GMM为我们提供了生成与我们的输入类似分布的新随机数据的方法。例如，以下是从这个16分量GMM拟合到我们原始数据的400个新点：
GMM非常方便，可以灵活地建模任意多维数据分布。
GMM是一种生成模型这一事实为我们提供了一种确定给定数据集的最佳组件数的自然方法。生成模型本质上是数据集的概率分布，因此我们可以简单地评估模型下数据的可能性，使用交叉验证来避免过度拟合。校正过度拟合的另一种方法是使用一些分析标准来调整模型可能性，例如 Akaike information criterion (AIC) 或 Bayesian information criterion (BIC) 。Scikit-Learn的GMM估计器实际上包含计算这两者的内置方法，因此在这种方法上操作非常容易。
让我们看看在moon数据集中，使用AIC和BIC函数确定GMM组件数量：
最佳的聚类数目是使得AIC或BIC最小化的值，具体取决于我们希望使用的近似值。 AIC告诉我们，我们上面选择的16个组件可能太多了：大约8-12个组件可能是更好的选择。与此类问题一样，BIC建议使用更简单的模型。
注意重点：这个组件数量的选择衡量GMM作为密度估算器的效果，而不是它作为聚类算法的效果。我鼓励您将GMM主要视为密度估算器，并且只有在简单数据集中保证时才将其用于聚类。
我们刚刚看到了一个使用GMM作为数据生成模型的简单示例，以便根据输入数据定义的分布创建新样本。在这里，我们将运行这个想法，并从我们以前使用过的标准数字语料库中生成新的手写数字。
首先，让我们使用Scikit-Learn的数据工具加载数字数据：
接下来让我们绘制前100个，以准确回忆我们正在看的内容：
我们有64个维度的近1,800位数字，我们可以在这些位置上构建GMM以产生更多。 GMM可能难以在如此高维空间中收敛，因此我们将从数据上的可逆维数减少算法开始。在这里，我们将使用一个简单的PCA，要求它保留99％的预测数据方差：
结果是41个维度，减少了近1/3，几乎没有信息丢失。根据这些预测数据，让我们使用AIC来计算我们应该使用的GMM组件的数量：
似乎大约110个components最小化了AIC；我们将使用这个模型。我们迅速将其与数据拟合并确保它已收敛合：
现在我们可以使用GMM作为生成模型在这个41维投影空间内绘制100个新点的样本：
最后，我们可以使用PCA对象的逆变换来构造新的数字：
大部分结果看起来像数据集中合理的数字！
考虑一下我们在这里做了什么：给定一个手写数字的样本，我们已经模拟了数据的分布，这样我们就可以从数据中生成全新的数字样本：这些是“手写数字”，不是单独的出现在原始数据集中，而是捕获混合模型建模的输入数据的一般特征。这种数字生成模型可以证明作为贝叶斯生成分类器的一个组成部分非常有用，我们将在下一节中看到。
