假设$X$是个随机变量服从标准正态分布，另一个随机变量$Y$满足$Y=X^2$，那么它们的协方差

$$\text{cov}(X,Y)=\text{E}(XY)-\text{E}(X)\text{E}(Y)=\text{E}(X^3)-\text{E}(X)\text{E}(X^2)=0-0\times\text{E}(X^2)=0$$

协方差为0，说明$X$和$Y$不相关。但是显然$X$和$Y$不独立。

这样的例子其实不少，相关系数其实描述的是两个变量的线性相关性。如果线性关系强，相关系数的绝对值就大；线性关系弱，相关系数就少，甚至为0。

一个简单的例子是，$X\in [-1, 1]$，$Y=|X|$，显然$X$和$Y$的相关系数为0，但是又显然不是独立的。类似地，可以得到很多这样的例子，比如$X\in[-1, 1]$，$Y=X^m$，$m$是任意一个偶数。

这个问题挂这么久了，一直没有人回答，那我就试试吧。

先说说MCMC是什么。第一个MC是Markov Chain，第二个MC是Monte Carlo。MCMC就是两者的结合，顾名思义，就是带有马尔可夫链性质的蒙特卡洛模拟方法。

-----------什么是马尔可夫链-----------

假设随机变量$X_t$表示$t$时刻发生的事件。一个随机过程$X_0,X_1,X_2,\cdots,X_T$，如果满足

$$P(X_{n+1}|X_n)=P(X_{n+1}|X_n, X_{n-1}, X_{n-2}, \cdots, X_0),$$

就称这个过程是一个马尔可夫分链。换句话说，在一个马尔可夫链当中，下个时刻的事件状态只和当前状态有关。

-----------什么是蒙特卡洛模拟-----------

蒙特卡洛模拟是基于大数定律的随机重复抽样方法。比如说，为了估计抛某个有偏差硬币落在正面的概率，我们可以重复抛$m$次，得到$k$次正面，那么$p=\frac{k}{m}$。比如说，为了估计圆周率，我们可以在正方形中画一个内切圆，然后对这个正方形随机重复投点$m$次，如果有$k$次落在圆内，那么可以估计$\pi$为$\frac{4k}{m}$。

-----------到底什么是MCMC-----------

举个简单的例子：假设如果今天晴天，明天下雨的概率是0.1；如果今天晴天，明天晴天的概率是0.9；如果今天下雨，明天下雨(rainy)的概率是0.5；如果今天下雨，明天晴天的概率也是0.5。问题来了，如果我们不知道今天的天气如何，怎么通过随机抽样来模拟未来10天的天气呢？

步骤1: 随机选定初始点$X_0$，比如可选择为晴天。

步骤2: 根据上述的概率，随机产生一个天气$X_1$（0.1概率为雨天，0.9概率为晴天）

步骤3: 根据上述的概率和$X_2$，随机生成$X_3$。

步骤4: 如此反复100次（越多越好）。

步骤5: 取出$X_{101},X_{102},\cdots, X_{110}$即可。

这样得到的10个天气就是随机抽取出来的。为什么步骤4中次数越多越好呢？因为$i$越小，$X_i$越容易被初始值$X_0$影响。当$i$变得很大时，就趋于稳定（依照转移矩阵稳定态概率的的随机）。

-----------在贝叶斯网络和隐式马尔可夫模型中的应用-----------

MCMC经常被用来估计复杂的贝叶斯网络中的后验概率分布。

类似地，在隐式马尔可夫模型中，需要计算似然估计，而求这些似然估计，需要对很多隐藏状态求和，计算量很大，所以可以通过MCMC来模拟求解。

似然函数

\begin{eqnarray*}&&L(\theta|x_1=0.32,x_2=0.12,x_3=0.81,\cdots,x_{10}=0.44)\\&=&P(x_1=0.32|X=\theta)P(x_2=0.12|X=\theta)\times\cdots\times P(x_{10}=0.44|X=\theta)\\&=&\frac{1}{\theta^{10}}\end{eqnarray*}

显然$\theta$越小，似然函数越大。$\theta$必须大于样本中的最大数。所以这个均匀分布的上界$X$的最大似然估计就是这个样本中的最大值1.72。

X的最大似然估计是最大值。但是最大似然太激进了，只用了一个sample，并不是统计意义上的估计。因为均匀分布是对称的，两个参数0和X应该有相同的统计意义，比如所有samples到这两点的L1distance的和相等，所得到的X是均值的两倍。为了防止xi>X,最后X=Max(2E(x),max(x)).

$X\sim N(0,1)$，$Y\sim N(0,1)$并且$X$和$Y$独立。要确定它们的乘积$XY$的概率分布，我们可以通过正态分布的密度函数来计算出$XY$的密度函数。

\begin{eqnarray*}f_{XY}(z)&=&\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}\frac{e^{-x^2/2}}{\sqrt{2\pi}}\frac{e^{-y^2/2}}{\sqrt{2\pi}}\delta(xy-z)dxdy\\&=&\frac{K_0(|z|)}{\pi}\end{eqnarray*}

这个$K_0(z)$是第二类修正贝塞尔函数。以上是正规解法，具体可参考本链接。

以下是粗暴解法，我们可以直接做similation看看。下图是根据10000个点描绘出来的概率密度曲线。

其实，我们还有另外一种聪明的方式来理解$XY$的分布。

$$XY = \frac{(X+Y)^2 - X^2 - Y^2}{2}.$$

$X$和$Y$都是正态分布，所以$X+Y$是正态分布。所以可见$XY$其实是一个卡方分布减去两个卡方分布的结果。

如果是非数值的样本(multinomial)，那么可以用卡方检验。

如果是数值样本，可以用柯尔莫哥洛夫-斯摩洛夫检验(K-S test)。

Categorical的数据：用卡方检验

Numerical的数据：用KS

Binary的数据：用T test

提一个我目前正在看而且比较冷门的吧，Wasserstein distance. 

我觉得挺有帮助的参考文章有 Wasserstein GAN and the Kantorovich-Rubinstein Duality 还有 stackExchange上面这个问题Kullback-Leibler distance for comparing two distribution from sample points (原本以为KL divergence也可以的但是看了这个之后就有些犹豫了就不说出来了hhh)

说一个不用假设检验的想法，对于数值的概率分布，把两个分布的累积分布函数的曲线画出来，然后求两个曲线的距离（L1或者L2）。人为设置一个阈值，小于它就说明它们两条曲线够接近，分布够相似。

$p=0.4$, $np=400$, $np(1-p)=240$，成功次数是符合$N(400,240)$正态分布的。

$$z=\frac{500-400}{\sqrt{240}}=6.455$$

z值大于3的概率是0.0013。大于6的概率就几乎等于0了。

假设$P$是代表p值的随机变量，要证明$P$服从[0, 1]均匀分布，只需要证明$\text{Pr}(P\leq z)=z$。

假定我们假设检验中的统计量是$T$，那么根据p值的定义，我们有$P=F(T)$，其中$F$是统计量$T$的概率累积函数。

接下来，证明很简单

$$\text{Pr}(P\leq z) = \text{Pr}(F(T)\leq z) = \text{Pr}(T\leq F^{-1}(z))=F(F^{-1}(z))=z$$

第二个等号是利用了概率累积函数弱单调增的性质。第三个是根据概率累积函数的定义得到的。

“我在一个参考书里看到说，假设检验的p value是服从[0, 1]均匀分布的。”

------------

准确地说，应该是“在零假设下，假设检验的p value是服从[0, 1]均匀分布的”

分情况讨论吧

第一种：

第一杯大于0.5升的概率

$$P_1=0.5$$

第二种：

假设第一杯倒出$x$升水，那么第二杯大于0.5升的概率

$$P_2=\int_{0}^{0.5}\frac{1-x-0.5}{1-x}dx$$

第三种：

假设第一杯倒出$x$升水，那么第一杯水和第二杯加起来小于等于0.5升的概率

$$P_3=\int_{0}^{0.5}\frac{0.5-x}{1-x}dx$$

其中有一杯大于0.5升的概率是

$$P_1+P_2+P_3=0.5+\int_{0}^{0.5}\frac{1-2x}{1-x}dx=2x+\log(|x-1|)+0.5$$

最后答案就是

$$1.5+\log(0.5)\approx 0.807$$

面试题都这么变态的吗，正在找工作的我感受到了恐惧。。

想了一晚上没想出个具体的解，但是想出个思路感觉也许可行（计算也比较麻烦），希望能得到指点。

假设事件A、B、C分别表示每杯水的水量，皆服从均匀分布，其中有一杯大于0.5升的相反事件就是三杯水都小于0.5升，于是有：

P(至少一杯水>0.5L) = 1 - P(A<0.5)P(B<0.5|A<0.5)P(C<0.5|A<0.5,B<0.5)

P(A<0.5)就是0.5了，接下来要解P(B<0.5|A<0.5)就要考虑条件分布P(B|A)，根据贝叶斯公式有：

$P(B|A) = \frac{P(A|B)P(B)}{P(A)}$

边缘密度p(A)就是1；现求先验分布p(B)：已知B服从均匀分布U(0,1-A)，这里的超参数A是未知的，但是A的分布我们知道。这里可以采用双层贝叶斯，相当于p(A)作为第二层先验，那第二层的贝叶斯分母部分就是B的边缘概率密度：

$p(B)=\int_{0}^{1-B} 1/(1-A)*1\ dA$

把A积掉以后得到一个关于B的概率函数p(B)；然后，易知密度函数p(A|B)=1/(1-B)，与p(B)相乘之后得后验密度p(B|A) ，将其对B和A在0到0.5上进行二重积分可以得到P(B<0.5|A<0.5)。

同样的道理P(C<0.5|A<0.5,B<0.5)我想也是可求的。

不过我想你们的面试单位一定有什么奇技淫巧，这种计算这么麻烦显然不是正确答案。抛砖引玉，等一个更好的解。最后祝愿楼主面试顺利通过。

感觉和这个问题（一米长的绳子，随机剪两刀，最长的一段有多长？）类似

试试Bonferroni方法

單向變異數分析的事後檢定这个也许对你用帮助的

假设$X_i\sim U(0,1)$，$i=1,2,\cdots, k$

那么$$F(x)=P(\min(X_i) \leq x)=1-P(\min(X_i)>x)=1-(1-x)^k$$

对$F(x)$求导，得到概率密度函数$$f(x)=k(1-x)^{k-1}.$$

根据期望的定义，进行积分

$$E=\int_{0}^1xf(x)=\int_{0}^1kx(1-x)^{k-1}=k\int^1_0(1-z)z^{k-1}dz=\frac{1}{k+1}$$

1/(n+1)吧

k个独立随机变量服从[0,1]均匀分布，把它们的采样排序后，第i小变量的分布是$Beta(i,k+1-i)$参考这里。而$Beta$分布的均值是$\frac{i}{i+k+1-i}=\frac{i}{k+1}$。所以$i=1$时，最小数的均值为$\frac{1}{k+1}$。

正巧刚才在看A/B test的维基词条（英语），其中提到了一下E-test。

这是一个用来比较两个服从泊松分布随机变量均值是否相同的检验。

楼主若是感兴趣可以阅读一下这篇论文：A more powerful test for comparing two Poisson means

所谓E-test，是wafer test的一种，而且它并不像其他functional wafer testing一样。E-test并不关心芯片的功能情况，而只关心一个工艺下最基本的性能。E-test的结构位于Scribe line里面。

我感觉这是芯片设计相关的内容，和本社区的内容无关吧。

我也没听过。或者题主记错了？

也许是别人笔误了，把F test写成了E test。

0-1分布(伯努利分布)、多项式分布、高斯分布(正态分布)、泊松分布等都属于指数分布族，指数分布族可以构建广义线性模型。

只要是概率密度函数能够写成特定的形式（参考第一页的公式），那么这个分布就可以称作指数族分布。

指数族分布有：

高斯分布

二项分布

伯努利分布

多项分布

泊松分布

指数分布

beta分布

拉普拉斯分布

gamma分布

对于回归来说，如果应变量y服从某个指数族分布，那么我们就可以用广义线性回归来建模。比如说如果y是服从伯努利分布，我们可以使用逻辑回归（也是一种广义线性模型）。

KL divergence又可以叫相对熵

从它的计算过程：

$$KL(p||q)=-\int p(x)\ln q(x)\mathrm{d}x-\left ( -\int p(x)\ln p(x)\mathrm{d}x \right )$$

可以看到不满足对称。

信息论里可以理解右边一项是x分布为p(x)时，传输x信息的最小编码长度。左边一项是x的分布为p(x)，却被误认为分布是q(x)情况下的实际编码长度。所以非负性也可以理解。

从信息论的角度来说，$KL(P||Q)$是指当我们用$Q$去逼近$P$的信息损失量。

从统计角度来说可能更好理解一点。两个概率分布$P$和$Q$，概率密度函数分别为$p(x)$和$q(x)$。

我们根据概率分布$P$去生成一个随机点$x$，那么这个$x$有多大的可能性是属于概率分布$Q$的呢？它们的似然函数分别为$p(x)$和$q(x)$。因为$P$是所谓的原分布，那么$q(x)$越接近$p(x)$，所以$\frac{q(x)}{p(x)}$就该接近1。

顺理成章，就很好理解$KL$为什么不对称了。$KL(P||Q)$指的是$P$为原分布，$KL(Q||P)$指的是$Q$为原分布。

$$KL(P||Q)=\lim_{n\rightarrow\infty} \frac{1}{n}\sum_{i=1}^n\log\frac{{q(x_i)}}{p(x_i)}$$

此外英文好得可以看看这篇文章，写得更具体。

如果题主需要一个具有对称性的度量的话，推荐Jensen–Shannon divergence，它可以看作是对称版的Kullback–Leibler divergence。

这里有两个层次的"独立"。第一是两个随机变量X和Y的独立，是一般讲的independent random  variable。第二是独立取得样本x和y，为了区别于第一个independent，可以叫Unpaired。图片中样本独立是指Independent Samples t Test，又叫Unpaired t Test，X和Y的分布是独立的，并且样本的取得也是独立的。比如分别在A和B两学校找100人称体重。和它相对应的是Paired t Test，等效于对X-Y做One Sample t Test, 看X-Y的均值是否为常数。比如同样100个人依次在两台秤上称体重。

图片一中的注释和图片二中的“简单随机样本”都可以理解为unpaired。可能举个paired的反例更容易理解。

扔两次，如果是【正面】【反面】，算是0；如果是【反面】【正面】，算是1。这两种情况是等概率的，概率都是$p(1-p)$。

如果出现其他情况，则无效，重新抛两次硬币。

扔两次，00、11时无输出，而01输出0，10输出1 。

所谓T检验的自由度实际上就是T分布的自由度。

那么为什么是$n-1$呢？

比如你有$5$个数，$x_1,x_2,\cdots,x_5$，它们的均值是$2$，也就是说

$$x_1+x_2+x_3+x_4+x_5=5\times 2 = 10$$

显然只要确定其中的四个量的值，第五个就自然确定了。所以自由的变量数是样本数减一。

类似问题 怎么理解库克距离(Cook's distance)?

先验分布是均匀分布的话就叫做flat了。

对于连续变量，就是pdf就是常数；对于离散变量，pmf就是1/n。

当我们没有什么信息的时候，就会使用flat prior，因为它比较“公平”，不失一般性。

Flat prior就是prior pdf是常数，画prior pdf时是一条水平线段。

主要的作用是对没有先验信息的问题，假设一个均匀的先验分布，然后套入贝叶斯的框架。具体说就是把一个最大似然函数（ML，maximum likelihood）套入一个最大后验概率（MAP，Maximum a posteriori）。

flat prior就是没有什么信息量的先验，或者常数先验。

假设转移概率矩阵为$P$，平稳概率分布为$\pi$，很显然

$$\pi P=\pi$$

也就是说

$$P^T\pi^T=\pi^T$$

所以只要求$P^T$的特征值为1的特征向量就可以了。

参考这里

一般是用查表，凑出pdf。如果不查表，直接用公式的话：

Moment generating function: $M_x(t)=E(e^{tx})=\int_{-\infty}^{\infty}e^{tx}f_X(x)dx$

把$t$换成$it$，得到Characteristic function：$\varphi_x(t)=E(e^{itx})=\int_{-\infty}^{\infty}e^{itx}f_X(x)dx$

再由Inversion formulae $f_X(x)=\frac{1}{2\pi}\int_{-\infty}^{\infty}e^{itx}\varphi_x(t)dt$

注意最后一个公式是傅里叶反变换。一般到这步还是查表。

如果是离散的情况，应该好办，就是求一组线性方程吧。

如果是连续的话，就看你能不能求出每个moment的通项了，然后和一些常用分布的moment进行比较，然后确定分布，然后就有pdf了。

如果是直接公式从mgf到pdf的话，我就不知道了。

期待有人能够回答。

这是个经典的问题，叫做经典矩问题(Classical Moments Problem)。

对于连续分布，如果只知道有限个矩函数，那么是无法确定概率密度分布的，这时可以考虑最大熵方法来估计这个概率分布。如果想知道细节的话，可以看这篇论文Maximum entropy and the problem of moments: A stable algorithm

正面朝上的概率是$p$，是一个0到1之间的随机变量

在开始扔硬币前，$p$的flat prior是$U(0,1)$，也就是说$p$服从0到1的均匀分布，其实也就等价于$\text{Beta}(1,1)$分布。

然后扔了10次硬币，6次朝上，4次朝下，此时$p$的后验分布为$\text{Beta}(1+6, 1+4)=\text{Beta}(7, 5)$。

如果$p$是正面朝上的概率的话，那么这个问题的flat prior是$p\sim U(0,1)$，而不是$p=0.5$。

指数分布是用来描述两次事件发生的时间间隔，指数分布的参数$\lambda$是单位时间内事件发生的平均次数。

所以对于一个指数分布$\text{Exp}(\lambda)$，它的期望是$\frac{1}{\lambda}$

比如说，你用一个指数分布来表示一个公交站台，每班车之间的时间间隔（分钟数），那么这个分布的参数$\lambda$就是平均每分钟有多少车进这个公交站台。

对于supervised learning，分布是指关于特征$X$和结果$Y$的联合分布$F(X,Y)$或者条件$F(Y|X)$。

我们说训练集和测试集服从同分布的意思是训练集和测试集都是由服从同一个分布的随机样本组成的，也就是

$$(X_{train},Y_{train}), (X_{test},Y_{test}) \text{ i.i.d. }~\sim F(X,Y)$$

对于unsupervised learning，分布是指特征$X$的分布$F(X)$，也就是

$$X_{train}, X_{test} \text{ i.i.d. }~\sim F(X)$$

补充几句：

但是现实中比较难做到这点，特别是当训练集是过去的数据，测试集是当下的数据，由于时间的因素，它们很可能不是完全同分布的，这就增加了预测难度。

这也是为什么一般交叉验证的误差往往小于实际的测试误差。因为交叉验证中每折数据都是来自训练集，它们肯定是同分布的。

如果训练集和测试集的分布风马牛不相及，那么根据训练集学习得到的模型在测试集上就几乎没有什么用了。所以我们训练模型和应用模型时一个重要的前提假设就是训练集和测试集是同分布的。

另外一个方面是牵涉到过拟合问题，即使训练集和测试集是同分布的，由于数据量的问题，训练集的分布可能无法完整体现真实分布，当我们过分去学习训练集分布的时候，我们反而会远离真实分布（以及测试集的分布），造成预测不准确，这就造成过拟合。

几何分布和泊松分布都是属于离散分布，指数分布是属于连续分布。

一个人不停地掷骰子，两个连续的6之间一共掷了多少次骰子。这个次数就是服从几何分布。

指数分布是几何分布的连续版本。如果平均每小时有三辆公交车经过，两辆公交车的间隔的时间就是服从指数分布。

如果平均每小时有三辆公交车经过，每小时经过的公交车的个数就付出泊松分布。

指数分布和泊松分布有点互逆的感觉。指数分布是事件每发生一次需要的时间，泊松分布是单位时间内事件发生的次数。

如果是想做假设检验的话，可以把数据集$S$从中位数$m$一分为二，得到两份，$S_1$和$S_2$

先对$S_1$里的元素$x$，做$x=m-x$的操作，得到$S'_1$

再对$S_2$里的元素$x$，做$x=x-m$的操作，得到$S'_2$

如果$S$是对称的，那么$S'_1$和$S'_2$一定是同分布的。所以检验$S$是否对称，就转化成了检验$S'_1$和$S'_2$是否同分布。我们可以通过KS-test来完成同分布检验，并且得到p值

应该用偏度来衡量。

$skew[x]=\gamma_1=E[(\frac{x-\mu}{\sigma})^3]=\frac{\mu_3}{\sigma^3}$

就等价于判断 平均数=中位数？

该用Mann-Whitney_U_test（和Wilcoxon rank sum test还不一样）。是比较两个sample的median/mean是否有区别。

这里有个例子。

我也不懂，学习了下，体会它的思路是：即使原数据（KPI）非正太分布，但是它们的rank $U$大概是正太分布，然后对$U$做normalization后再查表检验是否能否定null假设。

--------------------------------------------------------------

做了个实验，对于非正态分布，t test的independent pair test会有很大偏差，比如期望的p接近0，但是得到的p分布较广；Mann-Whitney_U_test才能得到正确的p值。而对于接近正态的分布，t test 更准确。

结论：如果接近正态分布，用t-test更准确，否则用Mann-Whitney_U_test。

实验如下：

产生左图里的两个gaussian mixture的sample x1,x2。

共生成K=1000次，sample的个数比较多，$n1=10000,n2=20000$。可

以看到$mean(x1)\approx mean(x2)$, $median(x1)>median(x2)$。所以如果检测中比较的是x1,x2的mean，p值应该接近0.

分别进行t test和Mann-Whitney_U_test，p值的分布如下。可以看到Mann-Whitney_U_test的p值才接近于期望。因为sample size$n1,n2$足够大，说明此时中心极限定理并不能让t test满足正态分布的前提。

如果x1,x2来源于同一个gaussian，t test 的p值更接近0。

--------------------------------------------------------------

https://en.wikipedia.org/wiki/Location_test

注意第一个表是针对ordinal measurement，也就是median，不是针对mean的。

#Compare t_test and Mann-Whitney_U_test(or rank sum test) for Gaussian mixture case
import numpy as np
import matplotlib.pyplot as plt
from numpy import random
import scipy as sp
from scipy.stats import ttest_ind,ranksums,mannwhitneyu
n1=10000 # sample size of group1
n2=2*n1 # sample size of group2
K=1000 # number of Monte Carlo tests
random.seed(0)
x1=np.zeros([K,n1])
x2=np.zeros([K,n2])

# 2 component Gaussian mixture
sigma1=.5
sigma2=.5
mu1=np.array([-2.0,1.]*10)
mu2=np.array([-1.,2.]*10)

pi1b=2./3 # prob(sample x1 from Gaussian1b)
pi2b=1./3 # prob(sample x2 from Gaussian2b)
n1_b=np.int(n1*pi1b)
n1_a=n1-n1_b
n2_b=np.int(n2*pi2b)
n2_a=n2-n2_b


mean1=np.zeros([K,1])
mean2=np.zeros([K,1])
median1=np.zeros([K,1])
median2=np.zeros([K,1])
var1=np.zeros([K,1])
var2=np.zeros([K,1])

t_value=np.zeros([K,1])
p_value=np.zeros([K,1])

ranksum_value=np.zeros([K,1])
p_ranksum_value=np.zeros([K,1])

for k in range(K):
    x1_temp0=random.normal(loc=mu1[0],scale=sigma1,size=[1,n1_a])
    x1_temp1=random.normal(loc=mu1[1],scale=sigma1,size=[1,n1_b])
    x1_temp=np.concatenate([x1_temp0,x1_temp1],axis=-1)
    x1[k,:]=x1_temp
    mean1[k]=x1_temp.mean()
    median1[k]=np.median(x1_temp)
    var1[k]=np.var(x1_temp)
    
    x2_temp0=random.normal(loc=mu2[0],scale=sigma2,size=[1,n2_a])
    x2_temp1=random.normal(loc=mu2[1],scale=sigma2,size=[1,n2_b])
    x2_temp=np.concatenate([x2_temp0,x2_temp1],axis=-1)
    x2[k,:]=x2_temp
    mean2[k]=x2_temp.mean()
    median2[k]=np.median(x2_temp)
    var2[k]=np.var(x2_temp)
    
    result=ttest_ind(x1_temp.flatten(),x2_temp.flatten(),equal_var=False)
    t_value[k]=result[0]
    p_value[k]=1.-result[1]
#     result=ranksums(x1_temp.flatten(),x2_temp.flatten())
    result=mannwhitneyu(x1_temp.flatten(),x2_temp.flatten())
    ranksum_value[k]=result[0]
    p_ranksum_value[k]=result[1]
    
    
bins = np.linspace(-3, 3, 100)

plt.figure(figsize=[15,5])
plt.subplot(131)
plt.hist(x1.flatten(),bins,label='x1',alpha=0.5,normed=True)
plt.hist(x2.flatten(),bins,label='x2',alpha=0.5,normed=True)
plt.title('True pdf of x1 and x2')
plt.legend()

plt.subplot(132)
plt.hist(mean1.flatten(),bins,label='x1',alpha=0.5,normed=True)
plt.hist(mean2.flatten(),bins,label='x2',alpha=0.5,normed=True)
plt.title('hist of mean')
plt.legend()

plt.subplot(133)
plt.hist(median1.flatten(),bins,label='x1',alpha=0.5,normed=True)
plt.hist(median2.flatten(),bins,label='x2',alpha=0.5,normed=True)
plt.title('hist of median')
plt.legend()
plt.show()

bins1 = np.linspace(0, 1, 100)

plt.figure()
plt.hist(p_value.flatten(),bins1,label='t_test',alpha=.5,normed=True,color='red')
plt.hist(p_ranksum_value.flatten(),bins1,label='Mann-Whitney_U_test',alpha=.5,normed=True,color='green')
plt.title('p_value of t_test and Mann-Whitney_U_test')
plt.legend()
plt.show()

1.中位数

2.数值变换

可以考虑用一些nonparametric bootstrap的方法

当数据量很大的时候，t-test或者z-test应该不太需要担心正态性了吧（中心极限定理）

感觉并没有什么联系。

几何分布是描述一个伯努力实验第一次成功时尝试的总次数。

超几何分布是描述从$N$个球中，无放回地抽出$k$个球，其中红球的个数。

$r$个独立几何分布rv的和是negative binomial,和超几何分布没有关系

用一句话概括：Shapiro–Wilk test是利用顺序统计量(order statistics)来检验样本是否服从正态分布。

如果想知道更详细的，可以看huluobo发的链接。

https://en.wikipedia.org/wiki/Shapiro%E2%80%93Wilk_test

不要伸手

可以用Markov Chain hitting time求解：

甲、乙状态转移图为

初始状态$x_0=[1,0,0,0]$，$x_{t+1}=x_{t}P$,$P$是状态转移矩阵

$K_i^A $是状态i到状态A的步数的期望

对于甲：

$$P_A=\begin{bmatrix}0.5 & 0.5 & 0\\0.5 & 0 & 0.5\\0 & 0 & 1\end{bmatrix}$$

$K_1^3 = 1 + \sum\limits_{j = 1,2} {{P_{1j}}K_j^3} =1+0.5K_1^3+0.5K_2^3$

$K_2^3 = 1 + \sum\limits_{j = 1,2} {{P_{2j}}K_j^3} =1+0.5K_1^3$

有$K_1^3=1+0.5K_1^3+0.5(1+0.5K_1^3）$

最后$K_1^3=6$

对于乙：

$$P_B=\begin{bmatrix}0.5 & 0.5 & 0\\0 & 0.5 & 0.5\\0 & 0 & 1\end{bmatrix}$$

$K_1^3 = 1 + \sum\limits_{j = 1,2} {{P_{1j}}K_j^3} =1+0.5K_1^3+0.5K_2^3$

$K_2^3 = 1 + \sum\limits_{j = 1,2} {{P_{2j}}K_j^3} =1+0.5K_2^3$

有$K_2^3=2$,$K_1^3=4$

由上可知，甲平均6步能达到HH，乙平均4步能到达HT，所以乙胜利的概率高。如果还要考虑平局的情况，要计算出甲乙胜出的概率会更复杂。改了下代码，加入了最终步数的计算。

import numpy as np
n=10000
a_wins, b_wins, draw = 0, 0, 0
a_step_all,b_step_all=[],[]
for j in range(n):
    coin_flips1, coin_flips2, flag = [],[], True
    a_step=0
    while flag:
        coin_flips1.append(np.random.randint(0, 2))
        a_step+=1
        if coin_flips1[-2:] == [1, 1]:
            a_step_all.append(a_step)
            flag = False
    b_step=0
    flag=True
    while flag:
        coin_flips2.append(np.random.randint(0, 2))
        b_step+=1
        if coin_flips2[-2:] == [1, 0]:
            b_step_all.append(b_step)
            flag = False 
    if a_step<b_step: a_wins+=1 elif a_step>b_step:
        b_wins+=1
    else:
        draw+=1
print ("P(A wins)=%.3f,P(B wins)=%.3f, P(draw)=%.3f" % (a_wins/n, b_wins/n, draw/n))
print ("P(A win+draw)=%.3f, P(B win+draw)=%.3f"%((a_wins+draw)/(n+draw), (b_wins+draw)/(n+draw)))
print("Mean of steps to reach target, A=%.3f, b=%.3f"%(np.mean(a_step_all),np.mean(b_step_all)))

结果：

P(A wins)=0.329,P(B wins)=0.527, P(draw)=0.144
P(A win+draw)=0.413, P(B win+draw)=0.587
Mean of steps to reach target, A=5.941, b=3.999

根据的xwemin解答，找到一个解法。

重贴下状态转移图。

#https://en.wikipedia.org/wiki/Absorbing_Markov_chain

import numpy as np

#transitioning matrixes

#states of(A,B) [11,12,21,22,13,23,31,32,33]

P=np.matrix([

            [1,1,1,1,0,0,0,0,0],

            [0,1,0,1,1,1,0,0,0],

            [1,1,0,0,0,0,1,1,0],

            [0,1,0,0,1,0,0,1,1],

            [0,0,0,0,4,0,0,0,0],

            [0,0,0,0,0,4,0,0,0],

            [0,0,0,0,0,0,4,0,0],

            [0,0,0,0,0,0,0,4,0],

            [0,0,0,0,0,0,0,0,4]])/4.

print('Transitioning matrix:')

print(P)



Q=P[:4,:4]

R=P[:4,4:]

Ir=P[4:,4:]



print('Fundamental matrix:')

N=np.linalg.pinv(np.eye(4)-Q)

print(N)



print('Initial state')

s0=np.array([1,0,0,0])

print('Absorbing probabilities')

B=N*R

#Absorbing probabilities from state 11

B0=s0*B

print(B0)

print('P(A wins)=%.6f'%np.sum(B0[0,2:4]))

print('P(B wins)=%.6f'%np.sum(B0[0,:2]))

print('P(draws)=%.6f'%np.sum(B0[0,4]))

答案是

P(A wins)=0.322314

P(B wins)=0.537190

P(draws)=0.140496

甲胜：P甲； 乙胜：P乙； 平：P平

甲H乙H

甲胜：1/4

乙胜：1/4

平： 1/4

1/4*(甲T乙H)

甲H乙T

甲胜：1/2

乙胜：0 

平 ：0 

1/4*（甲T乙H） 

1/4*（甲T乙T）

甲T乙H

甲胜：0

乙胜：1/2

1/4（甲T乙H）

1/4（甲H乙H）

甲T乙T

甲胜：P甲

乙胜：P乙

平 ：P平

列方程求解

P甲 = 39/121

P乙 = 65/121

P平 = 17/121

To: 得得得

可能是我说的不太清楚，甲乙之间的比赛是有可能出现平局的，也就是甲扔出了HH乙扔出了HT， 你的code里是甲乙同时看一个硬币，而我的问题时甲乙两个人同时扔两个硬币。我把你的code修改了一下结果好像确实是乙胜出的概率大。具体解释我还没有想好。

import json
import numpy as np


a_wins, b_wins, draw = 0, 0, 0
for j in range(100000):
    coin_flips1, coin_flips2, flag = [],[], True
    while flag:
        coin_flips1.append(np.random.randint(0, 2))
        coin_flips2.append(np.random.randint(0, 2))
        if coin_flips1[-2:] == [1, 1] and coin_flips2[-2:] != [1, 0]:
            a_wins += 1
            flag = False
        elif coin_flips1[-2:] != [1, 1] and coin_flips2[-2:] == [1, 0]:
            b_wins += 1
            flag = False
        elif coin_flips1[-2:] == [1, 1] and coin_flips2[-2:] == [1, 0]:
            draw += 1
            flag = False
print (a_wins, b_wins, draw)

>> 32116    53763     14121

一样大，python写个simulation

a_wins, b_wins = 0, 0
for j in range(100000):
    coin_flips, flag = [], True
    while flag:
        coin_flips.append(np.random.randint(0, 2))
        if coin_flips[-2:] == [1, 1]:
            a_wins += 1
            flag = False
        elif coin_flips[-2:] == [1, 0]:
            b_wins += 1
            flag = False

最后得到的a_wins和b_wins几乎是相等的

具体分析的话，就是扔到T对于甲和乙都是等价的没用；扔到H对甲和乙一样，下一轮获胜的概率各50%。所以甲乙一样。

肯定是乙啊

并没有严格的定义。

尾部是指远离整个分布的中心（比如中位数，均值）的那片区域。如果概率累积函数(cdf)尾部在整个cdf所占的比例过高，这就是重尾分布(heavy-tailed)。

重尾分布有两种情形，一种是尾巴非常长，所以叫长尾分布(long-tailed)；另一种是尾巴相对没有那么长，但是比较大，所以叫肥尾分布(fat-tailed)。

都差不多吧

顺序统计量：

把样本$X_1,X_2,\cdots,X_n$按照观察值从小到大的顺序把它们重新排列成

$$X_{(1)}\leq X_{(2)}\leq \cdots \leq X_{(n)},$$

称统计量$X_{(1)},X_{(2)},\cdots, X_{(n)}$为顺序统计量。

$R=X_{(n)}-X_{(1)}$称为样本极差。

Tweedie分布是一种泊松分布和伽马分布的复合分布。

有三个参数，一个是$p$，当$p=1$，Tweedie就是泊松分布，当$p=2$，Tweedie就是伽马分布。

第二个参数是$\mu$，是Tweedie分布的期望。

第三个参数是$\phi$，控制Tweedie分布的方差。方差$Var=\phi \mu^p$。

直白一点理解就是，Tweedie随机变量是$X$个伽马随机变量的和。其中$X$是服从期望为$\mu^{(2-p)/(\phi(2-p))}$的泊松分布，其中每个伽马随机变量是独立同分布的，服从$\Gamma((2-p)/(p-1), \phi(p-1)\mu^{p-1})$。

Tweedie分布最明显的一个特点是以一定的概率生成数值为0的样本。Tweedie分布在商业场景中有很多例子，比如人们进入某个商城后的消费额。一部分人只是进去随便逛逛（比如蹭空调），所以消费额是0，另一部分客户是有消费额的，是随机的连续数值。

直方图是如下效果

用python可以简单实现一下生成tweedie分布随机变量

import numpy as np

import matplotlib.pyplot as plt



def tweedie(n,p,mu,phi):

    #checking the value of variance power between 1-2

    result = np.full(n, np.nan)

    # calculating mean of poisson distribution

    lambdaa = mu ** (2 - p) / (phi * (2 - p))

    # shape parameter of gamma distribution

    alpha = (2 - p) / (1 - p)

    # scale parameter of gamma distribution

    gam = phi * (p - 1) * ( mu ** (p - 1))

    # Generating Poisson random sample

    rs = np.random.poisson(lambdaa, n)

    for i in range(n):

    # Generate single data point of gamma distribution using poisson random variable

        result[i] = np.random.gamma(rs[i] * np.abs(alpha), gam, 1)

    return result



x = tweedie(100000, 1.2, 3, 2)

plt.hist(x, bins=40)

plt.show()

上面随机产生了100000个随机点，分布如下。可以看出有很多0。

弱大数定律：样本均值依概率收敛于期望值。

强大数定律：样本均值以概率1收敛于期望值。

一个事件的odds就是这个事件发生的概率比上不发生的概率，也就是$\frac{p}{1-p}$。

odds ratio就是两个事件odds的比率。就等于

$$\frac{p_1(1-p_2)}{(1-p_1)p_2}$$

准确来说不是$\chi^2$分布，而是$\chi^2$ test。

$\chi^2$ test是用来检测两个categorical特征的独立性。每个$\chi^2$ test会对应一个$\chi^2$ statistic，这个数值越大说明两个特征的相关性越大。

对于分类问题，我们可以用$\chi^2$ test得到每个categorical特征和y的$\chi^2$ statistic，我们最终只保留$\chi^2$ statistic比较大的特征。这个就是特征选择的过程。

也不一定只能是分类问题，对于回归和数值特征也可以，但是我们要先做分箱处理，把连续特征和y进行离散化，然后再用卡方检测进行特征选择。

威布尔分布是指数分布的一个推广。

指数分布的例子就包括灯泡的寿命，而且是无记忆性；也就是说一个使用了一年的灯泡和一个新灯泡在未来某一个时间点坏掉的概率是一样的。但是这个似乎又不是太合理。

威布尔分布对此进行了修改，引进了一个变量$\beta$，当$\beta=1$它就是指数分布；当$\beta<1$，此时威布尔分布就是“越用越耐用”，在商业上的应用的例子，注册越久的会员越不容易流失；当$\beta>1$，此时威布尔分布就是“越老越不行”，比如说越老的车越容易报废，越老的零件越不可靠。

令$Z=max(X,Y)$, 作图如下。

可以看到红色区域的mean=1.5，紫色mean=2.5，蓝色mean=2.5，

绿色mean= $\int_{1}^{2}\int_{x}^{2}ydydx=5/3$

橙色mean= $\int_{1}^{2}\int_{y}^{2}xdxdy=5/3$

整体mean=$(1.5+2.5+2.5+5/3*0.5+5/3*0.5)/4=2.04$

import numpy as np
n=100000
x = np.random.uniform(0, 2,n)
y = np.random.uniform(1, 3,n)
z=np.maximum(x,y)
print(np.mean(z))

2.0405704490035905

偷懒做个模拟

import numpy as np
x = np.random.uniform(0, 2, 1000000)
y = np.random.uniform(1, 3, 1000000)
print(np.mean(np.max([x, y], axis=0)))

得到的结果是2.041519543743509

log-normal分布是指对一个正值的随机变量取自然对数后，这个随机变量服从$\mathcal{N}(\mu,\sigma^2)$的正态分布。

log会对大的数值压缩更大，log-normal分布是右斜的，例如

也很明显是长尾的分布，生活中很多例子都是长尾的分布，比如每个家庭的收入，产品的销量等等。

log-normal适用于各种因素相乘得到结果的应用。其原理是乘法在log scale上是加法，并且基于中心极限定理，很多随机变量相加会趋近于正态分布，从而有log-normal。比如说身高，有基因，营养，运动，睡眠等等很多因素影响生长率，这些生长率的乘积是总生长率，本身生长率也是时间的函数。这些乘法性的随机变量在log scale变加法性变量，其和又满足正太分布。

Black–Scholes model

https://en.wikipedia.org/wiki/Black%E2%80%93Scholes_model

假设有两个骰子$X$和$Y$，分别扔了$n_X$次和$n_Y$次。

对于$s=1,2,\cdots, 6$，$N_X(s)$是骰子$X$扔出$s$的次数，$N_Y(s)$是骰子$Y$扔出$s$的次数。

我们的零假设是$X$和$Y$同分布，在零假设下扔出$s$的概率的估计为

$$\hat p_s=\frac{N_X(s)+N_Y(s)}{n_X+n_Y}$$

根据这个概率估计，骰子$X,Y$得到$s$的次数的期望分别为

$$\hat{\mathbb{E}}(N_X(S))=n_X\hat p_s$$

$$\hat{\mathbb{E}}(N_Y(S))=n_Y\hat p_s$$

下面计算卡方检验的统计量$D$，

$$D=\sum_{i=1}^6\left(\frac{(N_X(s)-n_X\hat p_s)^2}{n_X\hat p_s}+\frac{(N_Y(s)-n_Y\hat p_s)^2}{n_Y\hat p_s}\right)$$

自由度就是$6-1=5$

当$D$大的时候，我们拒绝零假设。具体的p值需要查卡方分布的表。

可以先对pdf进行Kernel smoothing，然后计算平滑后的曲线的极大值点，如果算出来是两个，那就是双峰分布。这样可以避免画出pdf再肉眼观测判断。

这个文章讲得很清楚，推荐

对高斯分布的贝叶斯推断

假设两个骰子都是公平的，每抛一次，a是第一骰子的数值，b是第二个骰子的数值，a是从1到6等概率取值，b同理。重复多次进行计数，如果两个骰子是独立的，那么每一对(a, b)出现的概率都是1/36。

下面就按照卡方检验，类似于检测一个骰子是否是公平的，检验这36对数值是否是等概率出现的。

可以用峰度（Kurtosis）来表示长尾得程度，正态分布得峰度是3，如果一个分布的峰度大于3，你可以认为它是比正态分布更长尾的。

峰度的计算公式是根据四阶中心距和标准差得到的

$$\text{Kurtosis}=\frac{\frac{1}{n}\sum_{i=1}^n (x-\bar x)^4}{\sigma^4}$$

$\bar x$是样本中的均值

用二阶矩、三阶矩、四阶矩是一个思路

我在网上还找到一篇论文，基本思想是用去除离群点的方法来定义长尾的程度，用了分位数和IQR，论文Measuring heavy-tailedness of distributions

在贝叶斯统计中，某一不确定量p的先验概率（Prior probability）分布是在考虑"观测数据"前，能表达p不确定性的概率分布。它旨在描述这个不确定量的不确定程度，而不是这个不确定量的随机性。

在贝叶斯统计中，一个随机事件或者一个不确定事件的后验概率（Posterior probability）是在考虑和给出相关证据或数据后所得到的条件概率。同样，后验概率分布是一个未知量（视为随机变量）基于试验和调查后得到的概率分布。“后验”在本文中代表考虑了被测试事件的相关证据。

比如，当我们得到一枚硬币，由于没有任何额外信息，我们只能猜测这个硬币是公平的，随机一抛正面朝上的概率$p$服从$\text{Beta}(1,1)$，这个分布就是先验分布。现在我们抛了硬币两次，并且两次都正面朝上，我们就获得了额外的信息，可以更新分布，$p\sim \text{Beta}(3,1)$。这个分布就是后验分布。

PMF is for concrete variable

PDF is for continuous variable

PMF = probability mass function（概率质点函数）

比如骰子：每个面的概率就是PMF

pdf = probability dense function（概率密度函数）

比如均匀分布

一个是分类变量，一个是连续变量

$k$阶样本原点矩

$$A_k=\frac{1}{n}\sum_{i=1}^n X_{i}^k$$

原点矩指的是和原点$0$相关的，中心矩就是和中心（均值）相关的。

$k$阶样本中心矩

$$B_k=\frac{1}{n}\sum_{i=1}^n (X_{i}-\bar X)^k$$

其中$\bar X$是样本的均值，也就是$A_1$，1阶原点矩