#########################
过度分散
#########################

泊松回归模型是处理计数数据最基本的模型，应用最为广泛。
然而泊松模型也不是万能的，它有一个局限性。
泊松模型没有分散参数，并且它的方差和期望相同，
这导致泊松模型无法处理方差和期望不相同的数据。
本章我们讨论离散数据中一个常见的问题，
即离散数据的方差大于它的期望，
此时泊松模型无法很好的处理，
需要借助其他的手段或者模型才行。

在GLM中，响应变量 :math:`Y` 的方差 :math:`V(Y)` 等于分散函数 :math:`a(\phi)` 和方差函数 :math:`\nu(\mu)`
的乘积。

.. math::

    V(Y) = a(\phi)\nu(\mu)

其中方差函数 :math:`\nu(\mu)` 是一个关于期望 :math:`\mu`
的函数，其代表着方差 :math:`V(Y)` 和期望 :math:`\mu` 的关系。
比如高斯模型，其方差函数为常量 :math:`\nu(\mu)=1` ，
说明高斯变量的方差与其期望是独立不相关的。
反之，对于泊松模型，其方差函数为 :math:`\nu(\mu)=\mu`
，这意味着泊松模型的方差和期望存在关系。

分散函数(参数) :math:`a(\phi)` 相当于一个缩放(scale)因子，
对 :math:`\nu(\mu)` 起到一个缩放的作用，
通过这个缩放因子使得分布可以灵活适应(拟合)任何方差的数据。
下表给出了 ``GLM`` 中常见分布的方差函数和分散函数，
从中可以发现一些特点。



.. csv-table:: GLM中的方差
    :header: "分布","类型","分散函数 :math:`a(\\phi)`", "方差函数 :math:`\\nu(\\mu)`"


    "Normal(Gaussian) :math:`N(\mu,\sigma^2)`", "连续分布",":math:`\sigma^2`", 1
    "Inverse Gaussian","连续分布", ":math:`-\sigma^2`", ":math:`-\mu^3`"
    "Gamma", "连续分布", ":math:`-\phi`", ":math:`\mu^2`"
    "Power","连续分布","xxx", ":math:`\mu^k`"
    "Bernoulli :math:`B(\mu)`", "离散分布",1, ":math:`\mu(1-\mu)`"
    "Binomial :math:`B(n,\mu)`", "离散分布",1, ":math:`\mu(1-\mu/n)`"
    "Poisson :math:`Poisson(\mu)`", "离散分布",1 , ":math:`\mu`"
    "Categorical :math:`Cat(K,\mu)`", "离散分布", "1", ":math:`\mu_k(1-\mu_k)`"
    "Negative binomial :math:`\mathop{NB}(\mu,\alpha)`", "离散分布", "1", ":math:`\mu+\alpha\mu^2`"




从上表中可以发现离散分布的分散函数都是常数，而其方差函数都是期望的函数，
这表明离散分布的方差都是由它的期望决定的。由于没有分散参数，
这就使得离散分布存在一个特有的现象，过度分散(overdispersion)。



什么是过度分散
#################################################

当训练好一个模型后，利用模型的预测值 :math:`\hat{\mu}` 和分散参数 :math:`\phi`
就能计算出模型预测值的方差

.. math::

    V_e = a(\phi)\nu(\hat{\mu}_i)


模型预测值的方差称之为名义方差(nominal variance)，用符号 :math:`V_{e}` 表示。
名义方差 :math:`V_e` 就是模型(概率分布)的理论方差，由方差函数计算得到的。
与之相对应的，是观测数据的真实方差，
用符号 :math:`V_{o}` 表示，代表观测(observed)方差。

.. math::

    V_{o} = V(Y_i)


名义方差 :math:`V_e` 和观测方差 :math:`V_{o}` 通常并不相同，
这和模型对数据的拟合程度相关，模型对数据拟合的越好，名义方差 :math:`V_e` 和观测方差 :math:`V_{o}`
就越接近。


模型的预测值是模型的期望值， :math:`\hat{y}_i=\hat{\mu}_i` ，
方差函数 :math:`\nu(\hat{\mu}_i)` 是期望值
:math:`\hat{\mu}_i` 的函数，反映的是期望值，
分散函数 :math:`a(\phi)` 起到缩放和调节的作用，
显然，对于名义方差 :math:`V_e` 来说，分散函数 :math:`a(\phi)` 的作用至关重要，
分散函数 :math:`a(\phi)` 使得名义方差 :math:`V_e` 尽量接近和观测方差 :math:`V_{o}`
，分散参数是用来控制模型方差的。
然而，对于离散模型，没有分散参数 :math:`\phi`
，这就是使得名义方差 :math:`V_e` 不能自由的适配任何数据的观测方差 :math:`V_{o}` 。

当模型不存在分散参数时，就无法拟合数据的真实方差，就会发生
:math:`V_e` 和 :math:`V_o` 相差较大的现象。
根据两者的关系，有两种异常情况。

:过度分散(overdispersion): 当 :math:`V_{o}` 大于 :math:`V_e` 时，称之为过度分散，这种情况在离散数据中经常发生。

:分散不足(underdispersion): 当 :math:`V_{o}` 小于 :math:`V_e` 时，称之为分散不足，这种情况虽然理论上存在，但实际中并不常见。
    因此本书我们重点讨论过度分散的问题，分散不足的问题暂不讨论。

在实际应用中，分散不足的情况很少见，过度分散却是非常常见的，尤其是对于计数数据。
计数数据最基础的模型就是泊松模型，而泊松模型的方差等于期望，实际应用很少有数据满足这个假设。
在实际的计数数据中通常都是方差大于期望的，因此用传统的泊松模型处理计数数据必然会面临过度分散的问题，
而一旦发生过度分散，就意味着模型无法拟合数据的实际方差。



通常在实际应用中遇到的真实数据，并不是某个单一的、标准的概率分布产生的，生成数据的真实分布是无从得知的。
我们做的工作就是寻找一个尽可能接近数据真实分布的模型（分布）去拟合数据，
比如，遇到连续值数据就用高斯分布；遇到二值离散数据就用二项分布；
遇到计数数据就用泊松分布。
**过度分散(overdispersion)或者分散不足(underdispersion)现象的本质就是：**
**我们选取的模型的理论方差和数据的真实方差是不相符的，并且相差比较大。**


.. hint:: 为什么只有离散数据会有过度分散的问题，而连续值分布就没有呢？

    这是因为连续值分布通常都有一个额外的分散参数 :math:`\phi`
    ，这个参数的作用就是缩放(scale)方差函数 :math:`\nu(\mu)` ，
    使得模型的名义方差 :math:`V_e` 可以更好的拟合观测样本的方差 :math:`V_{o}` ，
    模型的分散参数 :math:`\phi` 就是专门用来拟合数据方差的参数，然而离散分布模型缺少这样一个参数，
    导致模型不能拟合数据的实际方差。

    然而，虽然连续值分布有分散参数 :math:`\phi` ，但这并不意味着连续值模型就不会发生过度分散现象，
    如果使用不当，同样也会导致过度分散。比如，在应用传统的高斯模型（线性回归）时，
    通常都是假设分散参数 :math:`\phi` 是常数 :math:`1`，这就相当于分散参数 :math:`\phi` 失去了作用，
    如果响应数据的方差和这个假设并不相符，同样会导致过度分散现象。

    此外，前面章节中我们介绍过，可以通过皮尔逊卡方统计量 :math:`\chi^2` 除以其自由度的方式估计出分散参数 :math:`\phi`。
    这种方式建立在分散参数与样本无关的假设之上，
    即假设所有样本拥有相同的 :math:`\phi` 。
    如果样本方差和这个假设严重不符，同样也会存在过度分散的可能，
    只是这种严重不相符的情况很少见，多数情况下，这个方式足够用了。


过度分散的检测
##################################

离散分布之所以存在过度分散问题，是因为离散分布的模型缺少一个分散参数 :math:`\phi`
去拟合数据的方差。
那么我们可以先假设模型存在分散参数 :math:`\phi` ，
此时观测方差 :math:`V_o` 就等于 :math:`\phi` 乘上名义方差 :math:`V_e`
。

.. math::
    V_o = \phi V_e

然后估计出 :math:`\phi` 的值，
如果 :math:`\phi \approx 1` ，意味着名义方差 :math:`V_e` 和观测方差 :math:`V_o`
是近似相等的，不存在过度分散，或者说问题不严重；反之，
如果 :math:`\phi \gg 1` ，
则意味着存在过度分散现象。
在 :numref:`ch_glm_estimate_phi` 和 :numref:`ch_glm_gof_chi` 讨论过，
可以通过皮尔逊卡方统计量估计分散参数 :math:`\phi`
。


.. math::


   \hat{\phi}=\frac{\chi^2}{N-p}

皮尔逊卡方统计量 :math:`\chi^2` 与其自由度(:math:`N-p`)
的商称为皮尔逊分散统计量。
皮尔逊分散统计量也是一个估计量，并不是十分准确的，
因此一般情况下，当皮尔逊分散统计量的值大于 :math:`2`
时才认为发生了明显的过度分散。
在早期的一些资料中，检测过度分散也可以用偏差统计量，
然而后续的一些研究中发现，皮尔逊卡方统计量更准确一些。




皮尔逊卡方统计量 :math:`\chi^2` 是模型拟合优度统计量，
本身可以用来评估模型拟合程度的。
理论上模型对数据拟合的足够好时， :math:`\chi^2` 渐近服从自由度为 :math:`N-p`
的卡方分布，它的期望值是 :math:`N-p` 。
分散统计量的估计是建立在 :math:`\chi^2` 服从卡方分布并且期望为 :math:`N-p` 的假设基础上。
如果模型对数据拟合的不够好，那么计算出的 :math:`\hat{\phi}` 的也会远大于 :math:`1`
，因此利用分散统计量判断是否存在过度分散并不是十分可靠的。
通常我们把由于模型欠拟合导致的过度分散称为表象过度分散(apparent overdispersion)，
真正的数据过度分散称为真实过度分散(real overdispersion)。


导致表象过度分散的直接原因是模型对数据拟合的不好，而导致模型拟合不好的原因常见的有如下几种情况：


- 数据包含异常点。
- 缺少有效的特征。
- 缺少高阶的组合特征。
- 特征数据没有进行归一化处理。
- 连接函数不正确。


表象过度分散不同于真实过度分散，表象过度分散仅仅是由于模型对数据拟合的不好导致，
只要针对性的解决上述问题，并提高模型的拟合程度后就会消除过度分散现象。
而真实过度分散是无法通过提高模型拟合程度解决的，
真实过度分散是数据的真实方差和模型的方差假设相违背，
显然要想解决真实过度分散就要改变模型的方差假设，
也就是要更换另一种新的分布模型才行。
下一章要讨论的负二项式模型就是一种全新的计数数据分布模型
，其方差函数是 :math:`\mu+\alpha \mu^2` ，
不再像泊松分布一样限定方差和期望相等，而是多了一个辅助参数 :math:`\alpha`
来适配数据的方差变化。




过度分散的影响
#################################################

在讲 ``GLM`` 的最大似然估计时提到过，分散参数 :math:`\phi`
并不影响协变量参数 :math:`\beta` 的迭代，
这里从 ``IRLS`` 算法的过程再验证一下。
``IRLS`` 算法参数的迭代公式为

.. math::
    :label: eq_od_010

    \beta &= (X^T W X)^{-1} X^T W Z

    &= \frac{X^T W Z}{X^T W X}

其中权重矩阵 :math:`W` 是一个对角矩阵，其值为

.. math::
    :label: eq_od_011

    W = [a(\phi)\nu(\mu) g'^2]^{-1} = [V_e(\mu) g'^2]^{-1}

其中 :math:`a(\phi)=\phi` 是与样本无关的，因此可以单独提取出来，则有

.. math::

     W = [a(\phi)\nu(\mu) g'^2]^{-1} = \phi^{-1}  W_o

.. note::

    在 :math:`a(\phi)=\phi/w_i` 的情形下，同样可以单独提出分散参数 :math:`\phi` ，
    而把样本权重 :math:`w_i` 留在 :math:`W_o` 中，结论是一样的。


代入到 :eq:`eq_od_010` 可得

.. math::

   \beta &=  \frac{X^T W Z}{X^T W X}

    &= \frac{\phi^{-1} X^T W_o Z}{\phi^{-1} X^T W_o X}

    &= \frac{X^T W_o Z}{X^T W_o X}

可以看到在 :math:`\beta` 的迭代过程中，分散参数 :math:`\phi` 并没有起到作用，
有没有分散参数 :math:`\phi` 不影响参数 :math:`\beta` 的更新。
因此过度分散问题并不影响协变量参数 :math:`\beta` 的估计。




虽然过度分散不影响 :math:`\beta` 的迭代，但是会影响它的标准误。
参数估计量 :math:`\beta` 的标准误的计算

.. math::


    SE(\hat{\beta}) &= \sqrt{\mathcal{J}^{-1}}

    &= \sqrt{ {(X^T W X)}^{-1} }

    &=  \sqrt{ {(X^T [a(\phi)\nu(\hat{\mu}) g'^2]^{-1} X)}^{-1} }

    &=  \sqrt{ {(X^T [V_e(\hat{\mu}) g'^2]^{-1} X)}^{-1} }



可以看出协变量参数估计量 :math:`\hat{\beta}`
的标准误差的计算是依赖模型方差的，如果模型的名义方差比观测方差小，
就会导致计算出的标准误 :math:`SE(\hat{\beta})` 比真实值偏小，
这会使得模型的假设检验失效。




标准误差的修正
######################################

既然清楚了过度分散的影响，就可以针对性的去解决它，
显然最直接到的方法就是人为给模型添加一个分散参数 :math:`\phi`
。因为分散参数不影响 ``IRLS`` 算法的迭代结果，
所以可以不用修改 ``IRLS`` 算法过程，
只需要在 ``IRLS`` 算法结束后修正一下估计量的标准误差。
分散参数 :math:`\phi` 的值就用皮尔逊卡分散统计量来表示，
修正后的标准误差为

.. math::


    SE(\hat{\beta}) =  \sqrt{ \hat{\phi}_{\chi^2} {(X^T W X)}^{-1} }



在对标准误差进行缩放后，该模型不再是似然模型，
因为标准误差既不是基于预期的信息矩阵，也不是基于观测的信息矩阵。
修正标准误差的方法还有另外几种的方法，比如，``William`` 过程、稳健方差估计等，
这里不在详细介绍，有兴趣的读者可以参考其它资料。

事实上，仅仅修正标准误并没有从本质上解决过度分散问题，过度分散的本质原因是模型假设和观测数据不匹配，
因此最根本的解决方法是更换更适合数据的模型，
下两章介绍的负二项式模型、零截断模型、零膨胀模型就是从模型根本上解决过度分散的问题。