双变量贯序排列python代码 双变量有序probit

第2章 Probability Distribution

• 2.2 Multinomial Variables

• P74 两变量的条件期望与条件方差

• 2.3 The Gaussian Distribution

• P86 高斯分布的参数辨识
• P94 序列估计

• Robbins-Monro 算法

• P99 高斯分布参数的贝叶斯估计

• 一般性序列估计
• 一维高斯分布均值的后验推断
• 一维高斯分布方差的后验推断 Gamma分布
• 一维高斯分布均值和方差联合的后验推断 Gaussian-gamma分布
• 高维高斯分布均值的后验推断
• 高维高斯分布方差的后验推断 Wishart分布
• 高维高斯分布均值和方差联合的后验推断 Gaussian-Wishart分布

• P103 学生t分布
• P107 von Mises 分布

• 2.4 The Exponential Family

• P113 一般形式
• P115 参数估计与充分统计量
• P117 共轭先验
• P117 无信息先验

• 尺度参数的无信息先验分布

• 2.5 无参数概率密度估计

• P122 核密度估计和近邻方法
• P125 kNN算法的一种无参解释

2.2 Multinomial Variables

P74 两变量的条件期望与条件方差

由Exercise2.8：考虑两个变量$x$和$y$，联合概率分布为$p(x,y)$. 那么
$\mathbb{E}[x]=\mathbb{E}_y [\mathbb{E}_x[x|y]], \text{ 这条较为广知} \\ var[x]=\mathbb{E}_y[var_x[x|y]]+var_y[\mathbb{E}_x[x|y]].$
这里$\mathbb{E}_x[x|y]$表示在条件分布$p(x|y)$下，$x$的期望。条件方差记号类似。
所以可知
$\mathbb{E}_\theta[\theta]=\mathbb{E}_\mathcal{D} [\mathbb E_\theta[\theta|\mathcal D]] \\ var_\theta [\theta] = \mathbb E_\mathcal D [var_\theta [\theta | \mathcal D]] + var_\mathcal D [\mathbb E _\theta [\theta | \mathcal D]]$
注意二式的右侧，第一项为$\theta$的后验分布方差的期望，第二项为后验分布期望的方差。
其中，$var_\mathcal D [\mathbb E _\theta [\theta | \mathcal D]] > 0$，所以$var_\theta [\theta] > \mathbb E_\mathcal D [var_\theta [\theta | \mathcal D]]$。也就是说观测到数据后，$\theta$的不确定性会减小。不过这只对平均而言成立。可以构造特殊的数据集，并让$\theta$的后验分布的方差变大。
（疑问：这好像并不能证明$card[\mathcal D]$越大，不确定性越小。这里猜测可以用类似方法证明，写出两个数据集$\mathcal D_1$和$\mathcal D_2$，$\mathcal D_1$放到两侧，$\mathcal D_2$放到右边，构造一个类似上述的式子？有空试一波！）

2.3 The Gaussian Distribution

P86 高斯分布的参数辨识

高斯分布在给定形式后，如何看出参数$\mu$和方差$\Sigma$，直接关注指数表达式即可：
$-\frac{1}{2}(x-\mu)^T \Sigma^{-1} (x - \mu) = -\frac{1}{2}x^T \Sigma^{-1}x + x^T\Sigma^{-1}\mu+const$
只要写成这样的形式，就能直接从二次项中读出$\frac{1}{2}\Sigma^{-1}$，从一次项读出$\Sigma^{-1}\mu$.
用这样的方法，2.3.1节写出当全变量为高斯分布时的条件分布
2.3.2节给出了边缘分布
2.3.3节给出了线性高斯模型的边缘分布和条件分布，即
$\begin{aligned} x &\sim \mathcal{N}(x|\mu, \Lambda^{-1}) \\ y|x &\sim \mathcal{N}(y|Ax+b, L^{-1}) \end{aligned}$其中$\Lambda=\Sigma^{-1}$被称为精度矩阵Precision Matrix）
线性高斯模型的结果$y$仍然是高斯分布
$\begin{aligned} \mathbb E[y] &=A\mu+b \\ cov[y] &=L^{-1}+A\Lambda^{-1}A^T \end{aligned}$

P94 序列估计

假定样本是一个一个序列观测的，记第$N$次观测后，均值估计为$\mu_{ML}^{(N)}$，则易知
$\mu_{ML}^{(N)} = \frac{1}{N}\sum_{n=1}^N x_n = \mu_{ML}^{(N-1)} + \frac{1}{N}(x_N - \mu_{ML}^{(N-1)})$
上式可看作是对$\mu$的不断修正。这里考虑一个一般化的序列学习算法：

Robbins-Monro 算法

对于一对随机变量$\theta$和$z$，并假定$f(\theta)= \mathbb E[z|\theta]$. 希望通过序列数据找到根$\theta^*$满足$f(\theta ^*)=0$.

假定$z$的条件方差有限，即$\mathbb E[(z-f)^2|\theta] < \infty$. 不失一般性，我们认为$\theta > \theta^*$时，$f(\theta) > 0$；$\theta < \theta^*$时，$f(\theta) < 0$. 则

其中$z(\theta^{N})$是给定$\theta^N$下$z$的观测。

$\{\alpha _N\}$表示正数序列满足

$\lim_{N\rightarrow \infty} \alpha_N = 0 \\ \sum_{N=1}^\infty \alpha_N = \infty \\ \sum_{N=1}^\infty \alpha_N^2 < \infty$

上式会以概率为1收敛到根。 第一项确保了修正项会收敛到一个有限值，第二项确保了不会对根欠收敛，第三项确保了累积噪声的方差有限，所以不会破坏收敛。（这个算法在强化学习的摇臂赌博机中也用到了）考虑一般的最大似然问题，参数$\theta_{ML}$是一个驻点，满足

当$N \rightarrow \infty$，上式即

注意这个形式，和Robbins-Monro的要求是一样的，可以得到

$z$可以看作是其中的$-\frac{\partial}{\partial \theta^{(N-1)}}\ln p(x_N| \theta^{(N-1)})$。

对于高斯分布的均值估计$\mu_{ML}^{(N)}$，即$z=-\frac{1}{\sigma^2}(x - \mu_{ML})$，取$\alpha_N=\frac{\sigma^2}{N}$，则得到一致的更新公式。

P99 高斯分布参数的贝叶斯估计

一般性序列估计

$p(\mu|D) \propto \left[ p(\mu) \prod_{n=1}^{N-1}p(x_n|\mu) \right] p(x_N|\mu)$

上式括号中的项可以看作是读入到第$N-1$个数据之后，得到的参数分布，可以看作是第$N$次的先验分布。

一维高斯分布均值的后验推断

如果已知方差，不知道均值，假定$\mu \sim \mathcal{N} (\mu|\mu_0, \sigma_0^2)$，$x|\mu \sim \mathcal{N} (x|\mu,\sigma)$，那么由
$p(\mu|X) \propto p(X|\mu) p(\mu)$

可得$p(\mu|X) = \mathcal {N} (\mu|\mu_N, \sigma^2_N)$，其中

$\begin{aligned} \mu_N &= \frac{\sigma^2}{N\sigma_0^2+\sigma^2}\mu_0 + \frac{N \sigma_0^2}{N\sigma^2_0+\sigma^2}\mu_{ML} \\ \frac{1}{\sigma^2_N} &=\frac{1}{\sigma_0^2} + \frac{N}{\sigma^2} \end{aligned}$
$N$是$X$中样本数，$\mu_{ML}=\frac{1}{N}\sum_{n=1}^N x_n$.

这个式子很有趣

• 当$N=0$时，等同于先验分布
• 当$N=\infty$时，等同于极大似然
• 随着$N$增大时，方差越来越小，$\mu$越来越确定
• 当$\sigma^2_0=\infty$时，等同于最大似然，方差很大意味着先验没有提供什么稳定的信息

一维高斯分布方差的后验推断 Gamma分布

如果已知均值，不知道方差，采用精确度$\lambda=\frac{1}{\sigma^2}$进行表示。高斯分布的方差后验为：

$p(X|\lambda) = \prod_{n=1}^N \mathcal {N} (x_n| \mu, \lambda^{-1}) \propto \lambda^{N/2} exp \left\{ - \frac{\lambda}{2} \sum_{n=1}^N (x_n - \mu)^2 \right\}$

注意，这种写法下，对应的先验共轭分布其实是Gamma分布！
$Gam(\lambda |a,b) = \frac{1}{\Gamma(a)}b^a \lambda^{a-1} exp(-b\lambda)$

如果记先验为$Gam(\lambda|a_0, b_0)$，则对应的后验为
$p(\lambda|X) \propto \lambda^{a_0 - 1} \lambda^{N/2} exp \left \{ -b_0 \lambda - \frac{\lambda}{2} \sum_{n=1}^{N}(x_n - \mu)^2 \right \}$

从中可以辨识出分布为$Gam(\lambda| a_N, b_N)$
$\begin{aligned} a_N &= a_0 + \frac{N}{2} \\ b_N &= b_0 + \frac{1}{2}\sum_{n=1}^N (x_n - \mu)^2=b_0 + \frac{N}{2} \sigma^2_{ML} \end{aligned}$

• 当$N$增大时，$a_N$增大，实际上，可以把$a_0$解释成是已经有了的$2a_0$个先验伪观测，$b_0$解释成是这$2a_0$个先验观测具有方差$\frac{b_0}{a_0}$
• 如果直接估计$\sigma^2$，而不是$\lambda$，那么得到对应先验分布是Inverse Gamma 分布。

一维高斯分布均值和方差联合的后验推断 Gaussian-gamma分布

如果方差和均值都不知道，那么$p(X|\mu,\lambda)$的连乘可以写成如下形式：
$p(\mu|\lambda)p(\lambda) \propto \mathcal{N} (\mu|\mu_0, (\beta \lambda)^{-1})Gam(\lambda|a, b)$
这也即共轭先验的形式，该分布叫做normal-gamma或Gaussian-gamma分布

高维高斯分布均值的后验推断

如果已知方差，不知道均值，这种情况下，均值仍然是高斯分布。

高维高斯分布方差的后验推断 Wishart分布

如果已知均值，不知道方差，如果记精确度矩阵$\Lambda = \Sigma^{-1}$，那么共轭先验分布为Wishart分布，这种分布可以看作是Gamma分布的高维推广，就类似于Beta分布和Dirichlet分布的关系。表达式为
$\mathcal{W} (\Lambda|W, \nu) = B|\Lambda|^{(\nu - D - 1) / 2}exp\left \{ -\frac{1}{2}Tr(W^{-1}\Lambda) \right\}$
其中$\nu$是自由度，$B$为归一化因子
$B(W, \nu) = |W|^{-\nu / 2} \left(2^{\nu D / 2} \pi ^{D(D-1)/4} \prod_{i=1}^{D} \Gamma (\frac{\nu + 1 - i} {2}) \right) ^{-1}$

如果直接对$\Sigma$估计，而不是$\Lambda$，则得到对应的共轭先验为Inverse Wishart分布

高维高斯分布均值和方差联合的后验推断 Gaussian-Wishart分布

如果方差和均值都不知道，那么共轭先验的形式为：
$p(\mu, \Lambda| \mu_0, \beta, W, \nu) = \mathcal {N} (\mu| \mu_0, (\beta\Lambda)^{-1}) \mathcal{W} (\Lambda|W, \nu)$

称之为Normal-Wishart或Gaussian-Wishart分布。

P103 学生t分布

如果一维高斯分布精确度先验为Gamma分布，均值已知，则$x$的边缘分布为
$\begin{aligned} p(x|\mu, a,b ) &= \int_0^{\infty}\mathcal {N} (x | \mu, \tau^{-1})Gam(\tau|a,b)d\tau \\ &=\frac{b^a}{\Gamma(a)}\left(\frac{1}{2\pi}\right)^{1/2} \left [b+\frac{(x-\mu)^2}{2}\right]^{-a-1/2}\Gamma(a+1/2) \end{aligned}$

如果记$\nu=2a, \lambda=a/b$，则上式化学生t分布

$St(x|\mu, \lambda, \nu) = \frac{\Gamma(\nu/2 + 1/2)}{\Gamma(\nu/2)} \left( \frac{\lambda}{\pi \nu}\right)^{1/2} \left[1+\frac{\lambda(x-\mu)^2}{\nu}\right] ^{-\nu/2 - 1/2}$

$\lambda$有时称为t分布的precision，$\nu$称为自由度。$\nu=1$时，退化为Cauchy distribution; $\nu \rightarrow \infty$时，成为高斯分布$\mathcal{N} (x| \mu, \lambda^{-1})$.

• 相比于高斯分布，学生t分布的一个优点抗离群点robust，学生t分布的尾巴比较厚，没有高斯分布那么敏感。另外，如果一组数据，高斯分布拟合得好，学生t分也能拟合好，因为高斯分布是学生t分布的一个特例。如图所示

如果再另$\eta=\tau b/a$，则学生t分布又可写为
$St(x|\mu, \lambda, \nu)=\int _0^\infty \mathcal{N} (x|\mu, (\eta \lambda)^{-1})Gam(\eta| \nu/2, \nu/2)d\eta$
通过该形式，可以扩展出高维学生t分布
$\begin{aligned} St( x|\mu, \Lambda, \nu) &=\int _0^\infty \mathcal{N} (x|\mu, (\eta \Lambda)^{-1})Gam(\eta| \nu/2, \nu/2)d\eta \\ &= \frac{\Gamma(\nu/2 + D/2)}{\Gamma(\nu/2)} \frac{|\Lambda|^{1/2}}{(\pi \nu) ^{D/2}} \left[1+\frac{\Delta^2}{\nu}\right] ^{-\nu/2 - D/2} \end{aligned}$

其中$D$是维度， $\Delta^2=(x-\mu)^T\Lambda(x-\mu)$

P107 von Mises 分布

一个二维高斯分布，关注其在以原点为圆心的单位圆下的条件概率分布，角度的分布为von Mises分布（循环正态分布）
$p(\theta|\theta_0, m)=\frac{1}{2\pi I_0(m)} exp\{ m cos(\theta - \theta_0)\}$
其中$m = r_0/\sigma^2, r_0 = \left \| \mu \right \|_2,\theta_0=tan^{-1}(\mu_y/\mu_x)$，而
$I_0(m)=\frac{1}{2\pi} \int _0^{2\pi} exp \{m cos\theta\} d\theta$
是归一化因子。

• 当$m$变大时，von Mises分布近似高斯分布

2.4 The Exponential Family

P113 一般形式

$p(x|\eta) = h(x)g(\eta)exp\{ \eta^T u(x)\}$

其中$x$可以是一维或多维，也可以是离散或连续。$g(\eta)$叫做natural parameters，可看作归一化因子
实际上，本章中上述讨论过的概率分布都是指数族分布的特例。

P115 参数估计与充分统计量

考虑一般参数$\eta$估计问题，最大似然得到
$p(X|\eta) \propto g(\eta)^Nexp \left\{ \eta^T \sum_{n=1}^N u(x_n) \right\}$

对数求导后得到
$-\nabla \ln g(\eta_{ML}) = \frac{1}{N}\sum_{n=1}^N u(x_n)$

• 注意这里$\sum_n u(x_n)$足够计算$\eta$，所以被称为充分统计量。例如对于Bernoulli分布，仅需要保存$\{x_n\}$的和，对于高斯分布，需要保存$\{x_n\}, \{x_n^2\}$各自的和。
• 当$N \rightarrow \infty$时，右侧变为$\mathbb E _x[u(x)]$.

P117 共轭先验

$p(\bm {\eta}| \bm \chi, \nu) = f(\bm\chi ,\nu)g(\bm\eta)^\nu exp\{ \nu \bm\eta^T\bm \chi\}$
其中$f$是一个归一化因子，$g$和$p(X|\eta)$中形式一样。易得后验
$p(\eta |\bf X, \bf \chi, \nu) \propto g(\eta)^{\nu + N} exp \left\{ \eta^T \left ( \sum_{n=1}^N \bf u(\bf x_n) + \nu \bf \chi \right )\right \}$
其中$\nu$被看作是先验伪观测数，每一次观测的统计量$u(x)$为$\chi$

P117 无信息先验

无信息先验这个东西稍微抽象，偏贝叶斯思维。解决的问题是在无先验时如何选择先验，选择的思想是先验要对后验的影响最小。
如果没有什么信息，我们假定先验是均匀分布，这么做存在两个困难：

• 在无限连续数域上发散。称之为反常先验分布。但如果后验分布是正常的，那么可以使用这样的分布（称之为广义先验分布）。例如高斯分布，如果假定均值先验是均匀分布，只要观测到一个数据点，那么后验就正常。
• 如果另一个参数是该参数的非线性变换，那么将不再是均匀分布

尺度参数的无信息先验分布

如果一个分布形式为
$p(x|\sigma) = \frac{1}{\sigma} f (\frac{x}{\sigma})$
其中$\sigma > 0$，$f(x)$已经归一化。
考虑$y=cx, \eta = c\sigma$其中$c> 0$. 那么
$p(y|\eta) = \frac{1}{\eta} f (\frac{y}{\eta})$
$x$和$y$的函数形式相同，所以$\eta$和$\sigma$应该有相同的先验分布，如果$\sigma$的先验分布为$\pi_\sigma(\sigma)$，那么
$\begin{aligned} \pi_\eta(\eta) &= \pi_\sigma(\sigma) \left |\frac{d\sigma}{d\eta} \right|=\frac{1}{c}\pi_{\sigma} (\frac{\eta}{c}) \\ \pi_\eta &=\pi_\sigma \end{aligned}$
取$\eta=c$，解得$\pi_\eta(\eta)=\frac{\pi_\eta(1)}{\eta}$，取$\pi_\eta(1)=1$，则先验分布为$1/\eta$.
这样的一个例子是高斯分布中的标准差
$p(x|\sigma) = \sigma^{-1}exp \left \{ -\left (\frac{x}{\sigma} \right )^{2} \right \}$
还有一种位置参数的无信息先验分布，可以看原书，推导出的结果是均匀分布。

2.5 无参数概率密度估计

P122 核密度估计和近邻方法

这种估计方法不明确给出概率分布的表达式，而是通过数据进行感知。柱状图其实就是一种无参数的概率密度估计方法。

另外一种常用的$p(x)$估计方法是观察$x$的小邻域。记$N$为总样本数，$K$为小邻域内样本数，如果小邻域足够小，认为小邻域内概率不变，则有
$p(x) = \frac{K}{NV}$

这里如果固定$V$，则该方法为核密度估计；如果固定$N$，则为近邻估计，即找以$x$为中心包含$K$个点的最小超球，当作$V$.

P125 kNN算法的一种无参解释

在近邻方法当中，如果有多个类，则对于第$\mathcal C_k$类，记样本数为$N_k$，小邻域内有样本数$K_k$，则
$\begin{aligned} \\ p(x|\mathcal C_k) &= \frac{K_k}{N_kV} \\ p(x) &=\frac{K}{NV} \\ p(\mathcal C_k) &= \frac{N_k}{N} \end{aligned}$
则后验为
$p(\mathcal C_k|x) = \frac{p(x|\mathcal C_k)p(C_k)}{p(x)} = \frac{K_k}{K}$
这样，kNN分类就可以解释为是近邻方法中，后验概率最大的类别。

• 1-NN分类器有一个很有趣的性质：当$N \rightarrow \infty$时，分类错误率不会超过贝叶斯最优分类器错误率的两倍

• 最优分类器可以理解为是看到了真实后验分布
• （我记得这个性质是要求概率连续的）
• 可以参考西瓜书P226

参考文献：
[1] Christopher M. Bishop. Pattern Recognition and Machine Learning. 2006