概率统计随机过程之相关分析与因果推断

概率统计随机过程之相关分析与因果推断

因果关系是人类不断探寻的深刻议题,不过想要探究因果联系并不是那么容易的,因此很多学者都会退一步从更弱的关联性分析入手,尤其在大数据时代,关联性的作用也不容小觑。因此,统计学上的关系分析是非常重要的一个环节。本文主要讲解分类数据和数值数据的列联分析和方差分析等内容。

因果分析的复杂性

统计学上如果想要进行因果分析,通常会有如下图的阶段:

统计学探索变量间关系

统计学探索变量间关系

首先,查看变量间是否具有关联性,没有关联性的就是相互独立的变量,一个变量的变化并不会对另一个变量产生影响;当发现两个变量具备关联性时,我们还得查看关联性的强弱,是强相关还是弱相关;在之后,需要检查这种相关性是不是有什么其他隐含的因素,比如二者都是同一个原因的结果,二者本身不具备因果性;最后才能进一步确认因果性。在确认因果性时,一般通过以下模型实现:

因果关系理论抽象.drawio.svg

因果关系理论抽象.drawio.svg

在复杂的统计模型中,其中上边的每一步也需要仔细、系统的研究。当然上述只是因果分析的简要流程。

在因果分析中,对于自变量、因变量不同的类型,也有着不同的分析方法,对于本科生水平大概如下:

统计分析一般方法

统计分析一般方法

其中,数据类型一般可以分成定性的分类数据(品质数据)、顺序数据和定量的数值数据,数值数据还可分为离散数据和连续数据。这些数据的级别是由低到高的,高阶数据可以转换为低阶数据,例如连续数据可以归并成离散数据,数值数据可以按照大小排成顺序数据,顺序数据也可以分成几类形成分类数据。但是低阶数据无法转换成高阶数据。

对于硕士研究生可能需要掌握到下面的成程度:

分类数据统计分析方法:

分类数据统计分析方法

分类数据统计分析方法

数值数据统计分析方法:

数值数据统计分析方法

数值数据统计分析方法

对于博士生大概是这样的:

组间差异检验.png

组间差异检验.png

接下来我将有选择的挑几个阐明。首先,先从本科阶段的内容说起吧。<(  ̄^ ̄)>

分类数据的\(\chi^2\)拟合优度检验

对于离散分布,比如对于二项分布\(B(n,p)\),我们希望验证其服从二项分布,进行了\(M\)\(n\)重伯努利实验得到\(M\)个值,其中实验成功次数为\(0,1,2,3,\dotsb,n\)的频数分别为\(m_0,m_1,m_2,\dotsb,m_n,M=\sum\limits_{i=0}^n m_i\)。那么如果想要验证该n重伯努利分布得到的随机变量\(X\)是服从二项分布\(B(n,p)\)的,需要怎么做呢?一个简单可用的方法是\(\chi^2\)拟合优度检验。

分布的拟合检验是在随机变量\(X\)分布未知时的检验(因为我们要验证的即是其分布类型),因此不同于参数的假设检验问题,属于非参数检验。一般而言分类数据的结果是频数,\(\chi^2\)检验是对分类数据的频数进行分析的统计方法。

之所以叫\(\chi^2\)拟合优度检验,是因为在1900年,统计学四大天王之一卡尔-皮尔逊证明提出对于实验统计出来的频数\(f_i,i\in \mathrel{\Theta}\),它和理论期望的频数\(e_i=M×p_i\),(\(M\)为总数,\(p_i\)为对应概率)存在以下关系: \[ X^2=\sum_{i\in\mathrel{\Theta}} \frac{(f_i-e_i)^2}{e_i}\sim \chi^2(|\mathrel{\Theta}|-1)\tag{1} \] 即构造的统计量\(X^2=\sum\limits_{i\in\mathrel{\Theta}} \frac{(f_i-e_i)^2}{e_i}\)应该服从自由度为\(|\mathrel{\Theta}|-1\)的卡方分布,\(|\mathrel{\Theta}|\)为实验可能出现结果的样数。并且期望频数越大,该分布与卡方分布越接近。当期望频数大于5时,与卡方分布符合比较好。此外,卡方分布只适用于观测数均不小于5的大样本场合。

我们从数学角度单纯地看式(1),\((f_i-e_i)^2\)实际是指实验做出的实际结果与理论分布差值的平方,体现的是实际与理论的差异,这个值越大,说明二者越不相符,分母的\(e_i\)更像是正则系数,降低绝对差值的比例影响。所以,从式(1)也可以看出随机变量\(X^2\)越大,我们越倾向于不认同实验分布服从理论分布,从数学角度讲是拒绝域在右侧的\(\chi^2\)检验。这种方式最早由卡尔-皮尔逊提出,通常也称之为皮尔逊\(\chi^2\)拟合优度检验。

下面我们用泰坦尼克号的存活率与性别是否相关举一个简单的例子:

卡方优度检验

卡方优度检验

拟合优度检验只针对一个分类变量进行检验,如果需要对两个或多个分类变量进行出里就需要列联分析。

列联分析

如果我们希望分析两个或多个分类变量之间的是否独立,可以使用列联表。列联分析是一种独立性检验,通常列联表常用于分类数据的两两分析,多维数据的多维列联表不太直观,一般用的较少。

原理也是使用卡方统计量。

列联表中,若两个分类变量\(A,B\),其中\(A\)\(r\)个可取值,记为\(A_1,\dotsb,A_r\)\(B\)\(c\)个可取值,记为\(B_1,\dotsb,B_c\),从总体中抽取样本容量为\(n\)的样本,设其中有\(n_{ij}\)个个体,其属性为\(A_i,B_j\),称其为频数,我们根据上述信息可制作频数列联表:

\(A\setminus B\) \(B_1\) \(\dotsb\) \(B_j\) \(\dotsb\) \(B_c\) 行和
\(A_1\) \(n_{11}\) \(\dotsb\) \(n_{1j}\) \(\dotsb\) \(n_{1c}\) \(n_{1\cdot}\)
\(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\)
\(A_i\) \(n_{i1}\) \(\dotsb\) \(n_{ij}\) \(\dotsb\) \(n_{ic}\) \(n_{i\cdot}\)
\(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\)
\(A_r\) \(n_{r1}\) \(\dotsb\) \(n_{rj}\) \(\dotsb\) \(n_{rc}\) \(n_{r\cdot}\)
列和 \(n_{\cdot 1}\) \(\dotsb\) \(n_{\cdot j}\) \(\dotsb\) \(n_{\cdot c}\) \(n\)

以上列联表是会根据实验或数据集给出的数据制作而成,都是已知的数据。我们将上述表中数据都除以总数\(n\)得到频率/概率表:

\(A\setminus B\) \(B_1\) \(\dotsb\) \(B_j\) \(\dotsb\) \(B_c\) 行和
\(A_1\) \(p_{11}\) \(\dotsb\) \(p_{1j}\) \(\dotsb\) \(p_{1c}\) \(p_{1\cdot}\)
\(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\)
\(A_i\) \(p_{i1}\) \(\dotsb\) \(p_{ij}\) \(\dotsb\) \(p_{ic}\) \(p_{i\cdot}\)
\(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\) \(\vdots\)
\(A_r\) \(p_{r1}\) \(\dotsb\) \(p_{rj}\) \(\dotsb\) \(p_{rc}\) \(p_{r\cdot}\)
列和 \(p_{\cdot 1}\) \(\dotsb\) \(p_{\cdot j}\) \(\dotsb\) \(p_{\cdot c}\) \(1\)

根据频率/概率关系有\(\sum_i\sum_j p_{ij}=1, \sum_j p_{\cdot j}=1, \sum_i p_{i\cdot}=1\)。这个表格的目的就是计算\(p_{i\cdot}\)\(p_{\cdot j}\)。如果变量\(A,B\)是独立的,那么会有\(p_{ij}=p_{i\cdot}p_{\cdot j}\),但是实际统计频率\(p_{ij}\)必然和理论值有所偏差,我们用计算得到的\(p_{i\cdot}\)\(p_{\cdot j}\)相乘,得到独立假设下的理论概率\(\hat{p}_{ij}=p_{i\cdot}×p_{\cdot j}\),再乘以总数\(n\)得到期望频数\(n\hat{p}_{ij}\),那么这就可以看成有\(r×c\)个可选值的卡方拟合优度检验,其自由度为\((r-1)×(c-1)\)。据此,检验统计量为: \[ X^2=\sum_{i=1}^r\sum_{j=1}^c \frac{(n_{ij}-n\hat{p}_{ij})^2}{n\hat{p}_{ij}}\tag{2} \] 其中,\(\hat{p}_{ij}=p_{i\cdot}×p_{\cdot j}=\frac{n_{i\cdot}}{n}×\frac{n_{\cdot j}}{n}\),同样的\(n_{ij}\)\(n\hat{p}_{ij}\)差别越大,统计量值越大,概率分布服从性也越差。又因为理论概率为独立假设下的概率分布,概率服从性差意味着两个分类变量独立性也越差。

下面举一个例子:

列联表分析.png

列联表分析.png

此外,对于数值型数据,我们也可以通过将其分割归类成几段,降维成分类数据,从而使用列联表分析。不过,我们还有一种对其有更佳的处理方式,即方差分析。

方差分析(Analysis of variance, ANOVA)

当自变量是分类变量,因变量是数值变量时的相关性分析,可以使用方差分析。比如探究不同教学方式是否对成绩有影响、不同专业毕业之后的薪资是否有区别等等。区别于列联分析,方差分析的因变量都是数值。

如果只是想知道是否数值型因变量是否受到分类型自变量影响,那么使用假设检验也是可以的。但是需要研究的目标变多时,例如设4个总体的均值分别为\(\mu_1,\mu_2,\mu_3,\mu_4\),如果用一般假设检验方法,如t检验,一次只能研究两个样本,要检验4个均值是否相等,就需要检验6次:

  • 检验1:\(H_0:\mu_1=\mu_2\)
  • 检验2:\(H_0:\mu_1=\mu_3\)
  • 检验3:\(H_0:\mu_1=\mu_4\)
  • 检验4:\(H_0:\mu_2=\mu_3\)
  • 检验5:\(H_0:\mu_2=\mu_4\)
  • 检验6:\(H_0:\mu_3=\mu_4\)

很显然,这样做十分的繁琐,并且多次检验会导致出错概率增加,如果设拒绝域\(\alpha=0.05\),即每次检验犯第一类错误的概率为0.05,做6次检验会使犯第一类(至少一次)的概率变成\(1-(1-\alpha)^6\approx 0.265\),相应置信水平会降低到0.735。因此使用方差分析一是可以提升检验的效率,二是可以增加分析的可靠性,避免多次检验造成的误差累积。

方差分析:在共同的显著性水平\(\alpha\)下,同时考虑多个平均值的差异。通常以F分布来进行检验,称为方差分析。

方差分析由统计学四天王之一的Fisher于1923年提出。我们在进行方差分析之前还要注意其需要满足以下三个条件:

  1. 正太总体。每个组的总体应服从正态分布,对于因素的每一个水平,其观测值是来自正太分布总体的简单随机样本。
  2. 方差齐性。各个总体的方差\(\sigma^2\)必须相同。
  3. 独立性。每个观测值必须是独立的。

在上述假设成立的前提下,要分析自变量对因变量是否有影响,在形式上也就转化成为检验自变量的各个水平(总体)的均值是否相等。而以上三种假设都有对应的检验方法。如正态性检验、方差齐性检验以及独立性检验。

方差分析的原理是对数据误差的来源进行分类分析。对于同一因素的不同处理水平,产生的结果可能有不同。根据误差来源,结果的不同可能确实是不同处理水平导致的,也有可能仅仅是因为随机误差。但是这种不同到底是不是不同因素水平导致的呢?方差分析就是通过将误差分解成组间误差和组内误差,用二者的比值的偏离程度来进行分析。而数据中,误差的体现可以由方差透露,因此对误差的分析就能够变成对方差的分析。

  • 组内误差:在因素同一水平处理下,数据的差异,这种差异只可能是随机性导致的。
  • 组间误差:在因素不同水平处理下,数据的差异,这种差异虽然也包含随机性,但也可能是不同处理水平造成的系统性差异。

单因素方差分析

如果只考虑一个因素(自变量)对结果(因变量)的影响,那么只需要单因素方差分析。我们将数据集用\(X\)表示,每一个样品\(x_{ij}\)表示在自变量第\(i,i\in\{1,2,\dotsb,r\}\)个处理水平下,获取到的第\(j,j\in\{1,2,\dotsb,m\}\)个结果。那么可以用\(\overline{\overline{X}}\)表示数据集整体的均值,\(\overline{X}_i\)表示第\(i\)个因素水平的组内均值。\(x_{ij}\)还可进行如下分解: \[ \begin{aligned} x_{ij}-\overline{\overline{X}}&=(x_{ij}-\overline{X}_i)+(\overline{X}_i-\overline{\overline{X}})\\ x_{ij}&=\overline{\overline{X}}+\underbrace{(x_{ij}-\overline{X}_i)}_{\text{组内误差}}+\underbrace{(\overline{X}_i-\overline{\overline{X}})}_{组间误差} \end{aligned}\tag{3} \] 为了进一步分析误差来源,我们假设\(x_{ij}\)所在的处理组其服从分布为\(N(\mu_i,\sigma^2)\)的正态分布(正态总体假设,且由于方差齐性,各个处理组方差都是\(\sigma^2\)),那么\(x_{ij}\)也可以表示成\(x_{ij}=\mu_i+\varepsilon_{ij}\),其中\(\varepsilon_{ij}\)是数据\(x_{ij}\)与组内真实均值\(\mu_i\)的离差,该离差的来源就是随机性。

而根据中心极限定理,组内数据平均值\(\overline{X}_i\)应服从\(N(\mu_i,\frac{\sigma^2}{m})\)的正态分布,\(m\)为组内数据数量,那么\(\overline{X}_i\)可以表示成\(\overline{X}_i=\mu_i+\varepsilon_i\),其中\(\varepsilon_i\)是数据的组内平均值\(\overline{X}_i\)与组内真实均值\(\mu_i\)的离差,不难证明\(\varepsilon_i=\frac{1}{m}\sum_{j=1}^m \varepsilon_{ij}\),即组内均值的离差等于组内数据离差的均值。式(3)中的组内误差等价于: \[ x_{ij}-\overline{X}_i=(\mu_i+\varepsilon_{ij})-(\mu_i+\varepsilon_i)=\varepsilon_{ij}-\varepsilon_i\tag{4} \] 所以组内误差来源是随机性

同理,若我们将总体平均表示成\(\overline{\overline{X}}=\mu+\varepsilon\)的形式,其中\(\mu\)是所有分布真实的均值(各组都服从正态分布,和也服从正态分布),\(\varepsilon\)是数据均值与真正均值的离差,式(3)中的组间误差等价于 \[ \overline{X}_i-\overline{\overline{X}}=(\mu_i+\varepsilon_i)-(\mu+\varepsilon)\tag{5} \] 如果组间没有系统性差异,那么组内真实均值\(\mu_i\)应该和总体均值\(\mu\)相同\(\mu_i=\mu\),此时造成差异的只有离差\(\varepsilon_i-\varepsilon\),就只有随机性造成的偏差。但当不同处理水平确实有影响时,那么某些组的均值就不会等同于整体均值\(\mu_i\neq \mu\),此时式(5)中就会存在系统误差\(\mu_i-\mu\)随机误差\(\varepsilon_i-\varepsilon\)。系统误差越大,组间误差就会越大。

为了去除误差正负号的影响以及方便计算,我们将总误差进行平方求和\(S_T=\sum_{i=1}^r\sum_{j=1}^m(x_{ij}-\overline{\overline{X}})^2\),结合式(3)有: \[ \begin{aligned} S_T&=\sum_{i=1}^r\sum_{j=1}^m(x_{ij}-\overline{\overline{X}})^2=\sum_{i=1}^r\sum_{j=1}^m[(x_{ij}-\overline{X}_i)+(\overline{X}_i-\overline{\overline{X}})]^2\\ S_T&=\sum_{i=1}^r\sum_{j=1}^m[(x_{ij}-\overline{X}_i)^2+(\overline{X}_i-\overline{\overline{X}})^2]+\underbrace{\sum_{i=1}^r\sum_{j=1}^m 2(x_{ij}-\overline{X}_i)(\overline{X}_i-\overline{\overline{X}})}_{=0}\\ S_T&=\underbrace{\sum_{i=1}^r\sum_{j=1}^m (x_{ij}-\overline{X}_i)^2}_{\text{组内误差平方和}}+\underbrace{\sum_{i=1}^r\sum_{j=1}^m (\overline{X}_i-\overline{\overline{X}})^2}_{组间误差平方和} \end{aligned} \] 我们令组内误差平方和\(S_e=\sum_{i=1}^r\sum_{j=1}^m (x_{ij}-\overline{X}_i)^2\),组间误差平方和\(S_A=\sum_{i=1}^r\sum_{j=1}^m (\overline{X}_i-\overline{\overline{X}})^2=\sum_{i=1}^r m (\overline{X}_i-\overline{\overline{X}})^2\),所以有 \[ S_T=S_E+S_A\tag{6} \] 这就是误差的平方和分解。下面我们不加证明地给出: \[ \frac{S_e}{\sigma^2}\sim \chi^2(n-r)\tag{7} \] 当组间误差没有系统误差,即\(\mu_i=\mu\)时: \[ \frac{S_A}{\sigma^2}\sim \chi^2(r-1)\tag{8} \]\(S_e, S_A\)二者独立。

因此,当不同处理组没有区别,即不存在系统误差时,式(7)(8)都服从卡方分布,那么他们的商(还需除以自由度)就应该服从F分布: \[ \frac{\frac{S_A}{\sigma^2×(r-1)}}{\frac{S_e}{\sigma^2×(n-r)}}=\frac{S_A/(r-1)}{S_e/(n-r)}\sim F(r-1,n-r)\tag{9} \] 考虑到系统误差\(\mu_i-\mu\)越大,\(S_A\)也就越大,那么式(9)的F统计量就越大。因此,该检验的拒绝域应为: \[ W=\{F≥F_{1-\alpha}(r-1,n-r)\}\tag{10} \] 通常我们会将上述计算过程的结果汇总成单因子方差分析表

来源 平方和 自由度 均方 F比 p值
因子 \(S_A\) \(f_A=r-1\) \(MS_A=S_A/f_A\) \(F=MS_A/MS_e\) p
误差 \(S_e\) \(f_e=n-r\) \(MS_e=S_E/f_e\) - -
总和 \(S_T\) \(f_t=n-1\) - - -

对于给定的\(\alpha\),若\(F≥F_{1-\alpha}(f_A,f_e)\),则认为因子显著,自变量的不同水平会对因变量有影响。

如果认为因子的影响是显著的,想要进一步检验到底是哪个水平影响比较显著,还可继续进行多重检验(也叫事后检验),常用的一种多重检验是利用了t分布的LSD-t检验,其他还有SNK-q检验、Turkey、Duncan、Scheffe检验等。

此外,还可以用组间方差占总误差的比例来衡量关系的强度,记为\(R^2=SSA/SST\Rightarrow R=\sqrt{SSA/SST}\),当\(R>0.5\),可认为是中等相关,若\(R>0.8\)可认为是强相关。

双因素方差分析

根据名字我们就知道,双因素方差分析是分析两个分类变量(常称为行因素和列因素)对试验结果的影响。根据两个因素对试验结果的影响是否独立,还可以分为无交互作用的(无重复)双因素方差分析和有交互作用的(可重复)双因素方差分析

双因素方差分析也需要满足方差分析的三个假设:正态性、方差齐性、独立性。

双因素方差分析的基本方法和单因素类似,区别是:

无交互作用时:分为行因素误差平方和、列因素误差平方和、随机误差平方和三类。可以将行因素、列因素当成两个单因素误差分析来看。

有交互作用时:还要添加第四类交互效应误差平方和。有交互作用时数据量会比较大。

正态性检验

正态分布时统计学中最重要的分布之一,判断一组数据是否来自正态总体是很多分析步骤的前置要求,甚至还有国标GB/T4882-2001专门设计了正态概率图来辅助我们判断数据是否服从正态分布。

正态概率图是一种特殊绘制的坐标图,我们将所给数据绘制在坐标图上,如果这些数据大概处于一条直线上,那么可以认为服从正态分布。

当然,数理统计学还有更技术化的方法,比例最常用的W检验和EP检验。

W检验

W检验全称sharpiro-wilk检验,是由二人于1965年提出的正态检验方法。其适用范围为样本容量\(8\leq n \leq 2000\),非常适合小样本的正态性检验,但是当样本容量小于7时,对偏离正态分布的检验不太有效,同时当数据量大于5000时也不适用。

正态性检验方法总结图:

正态性检验.png

正态性检验.png

贝叶斯检验主要利用了KL散度,这是衡量两种分布偏差程度的一种度量,也叫相对熵。

正态性检验的计算往往十分复杂且常常需要查表,最好使用计算机程序辅助计算。

方差齐性检验

总结下这几种方法的利弊及适用条件:方差比、Hartley检验、Bartlett检验都需要原始数据是正态分布,Levene检验和BF法对正态分布不是很依赖。比较常用的是Levene检验,适用于多组方差的比较,且对正态性没要求。 https://zhuanlan.zhihu.com/p/313397172