Zhifu、JonyKai、汤圆、夏雅薇编译

不想去健身房的我,最后被贝叶斯分析说服了...

可能经常你会听到一些很主观的评价比如“你太瘦了”或者“你怎么那么高”,但这里瘦或者高都是基于评价者的主观判断和视觉记忆做出的评述,并没有严格的参照。

作者从小被人说体型瘦小,于是他用了贝叶斯分析最终得出了自己的体重确实低于本国平均水平的结论。没有比直白的数字更有说服力了,想要说服自己健身的小伙伴不妨也试试作者的统计学暴击法!

这篇文章将叙述一个在线性回归理论中应用贝叶斯分析的有趣试验(一个小秘密:我在这篇文章中使用了公制计量单位呦)。

如文章标题所述,我将会对自己的体格进行一番科学的研究。

在开始之前,我希望大家多了解我一些:我在越南出生,在新加坡上高中,现在在美国上大学。我经常因为身形瘦小被人们取笑,说我应该多去健身房锻炼增肌,拥有更强壮的体魄。这些评价我一般都一笑置之。对于一个身高169厘米(5尺6寸),体重58公斤(127磅)的人来说,我的BMI指数(20.3)几近完美。

仔细一想,大家可能没说错:我比一般的越南男性要高,但却只有平均体重(维基百科里越南男性的平均体重是58公斤,平均身高是162厘米),“看起来”可能是要稍微瘦一些。

这里“看起来”是关键:背后的逻辑很清楚,不是吗?如果你把自己抻长一些,体重不变,那确实应该看起来苗条一些。我把这个看作是严肃的科学问题,并准备深入研究。

对于一个169厘米高的越南男性来说,我到底轻了多少呢?

我们需要一种有理有据的方式来研究这个问题。有个好方法是尽可能多地找到越南男子身高和体重的数据,来判断我在这个样本中的位置。

越南人口数据

在浏览各种网页后,我找到了一份调查研究数据,包含超过10,000名越南人的人口统计信息。将抽样条件设置为年龄18-29岁的越南男性,从而得到数量为383的样本,这个样本足以用来进行接下来的分析啦。

首先,通过人口体重的直方图,看看我在年轻越南男性中的体重分布位置。

红线表示样本的中位数,而橙色线表示平均值

这张直方图表明我的体重略低于383名年轻越南男子的体重平均值和中位数。看起来是与我们要研究的相关呀!然而问题并不在于我的体重与这个样本本身的比较,是假设这383个人可以代表越南男性,在身高169厘米的情况下,我的体重与整个越南人口相比处于一个什么位置。为此,我们需要进行回归分析

首先绘制一个身高和体重的二维散点图

好吧,看起来我处在平均水平。但是如果我们只看身高169厘米的数据(想象一条垂直x轴于169厘米这个刻度并穿过红点的直线),我的体重在他们之中处于下游。

另一个重要的发现是越南男性身高和体重呈正相关。我们将进行定量分析来进一步了解这种关系。

首先,让我们快速添加“普通最小二乘”线。我稍后会详细介绍这一点,现在先在图上展示出来。

最小二乘线可以表示为y = -86.32 + 0.889x,这表明通常情况下,我这个年龄的越南男性,每增加1厘米的身高,体重会增加0.88千克。

但是,这并没有解决我们的问题;身高169厘米,体重58公斤到底是太沉,太轻还是刚刚好呢?要以定量的方式进一步解释这个问题,如果有所有身高1米68的人的体重分布,那么我的体重排在前25%,50%或75%的几率是多少?要弄清这一点,我们需要深入学习并理解回归背后的理论。

线性回归理论

线性回归模型中,Y变量(在我们的例子中,是人的体重)是x(身高)的线性函数。在这个线性关系中截距和斜率分别为β0和β1;也就是说,假设E(Y | X = x)=β0+β1x。我们不知道β0和β1是多少,所以将它们视为未知参数

在大多数标准线性回归模型中,我们进一步假设给定X = x的情况下,Y的条件分布是正态分布的。

这就是基本的线性回归模型:

可以被改写成:

注意,在许多模型中,我们可以用精度参数τ替换方差参数σ,其中τ= 1 /σ。

总结:因变量Y遵循由平均数μi和精度参数τ决定的正态分布。μi是由β0和β1决定的X的线性函数。

最后,我们还需假设未知方差不依赖于x;这种假设称为同方差性。

涉及的内容很多,都涵括在下面这张图里啦。

图像来源:Joseph Chang(耶鲁大学)

在实际的数据分析问题中,我们掌握的只是黑点 - 数据。使用这些数据,我们的目标是推断不知道的事情,包括β0,β1(在图片中的蓝色虚线)和σ(它决定了在给定一个y值的时候,红色正态分布密度的宽度)。注意,每个黑点周围的正态分布看起来完全相同。这是同方差性的本质。

估算参数

现在,有几种方法可以估算β0和β1。如果你使用普通最小二乘估计这样的模型,你不必担心概率公式,因为你正在寻找β0和β1的最佳值,而这是通过最小化拟合值与预测值的平方误差得到的。

另一方面,你可以使用最大似然估计(MLE)来估计此类模型:通过最大化似然函数来寻找参数的最佳值。


注意:一个有趣的结果(我没有放上具体的数学证明)是,如果我们假设误差项也属于正态分布的话,那么最小二乘估计的结果也是最大似然估计的结果。

贝叶斯方法的线性回归

不同于最大化似然函数,贝叶斯方法会假设参数服从一个先验分布。根据贝叶斯公式计算出参数后验概率

贝叶斯方法的似然函数同上面的类似,不同之处在于加入了对估计参数β0,β1,τ的先验分布:

等等,什么是先验分布,为什么这会让公式看上去更加复杂?

请相信我,先验分布虽然看上去很奇怪,但实际上很直观。事实上,我们对未知的参数(例子中的β0,β1,τ)分配一个看上去很随意的先验分布,这里存在着很强的哲学缘由。

先验分布能够反映出在没看见数据之前我们对数据的假设理解。在观察过一些数据之后,我们应用贝叶斯公式,就得到了同时考虑到了先验和数据的未知参数后验分布。根据后验分布,我们就能预测出未来的数据的分布。

最终的参数估计虽然取决于数据和先验分布,但是如果数据中包含的信息越多,那先验的影响就越小。

那么我该如何选择先验分布

这是个好问题,因为这里存在着无数种可能。理论上只有存在一个正确的先验能够能够反应出你的先验假设。但是在实际中,先验分布的选取是相当的主观,甚至有时可以是任意的。

我们可以选择一个有着较大标准方差(意味着精度很低)的正态分布先验。比如,我们假设参数β0和β1是服从均值为0标准方差为10000的正态分布。这种分布是毫无信息的分布,因为分布十分平坦(这意味着,参数在任意区间的取值概率几乎相同)。

如果选取了这种类型的先验分布,那么我们就不用考虑在这类分布中哪种分布更好,因为分布几乎都很平坦,在每个地方的概率都可以忽略不计。此外,后验分布不会受这种分布的影响。

同样,对于精度τ,因为其必须是非负的,所以需要选取一个取值限定在非负范围的分布。比如,在这里可以选取一个带有较小形状和尺度参数的伽马分布。

另外一种很有用的不附带信息的分布是均匀分布。如果对 σ 或 τ选择了均匀分布,那么你最终将会得到这样的模型。如下图所示,由John K. Kruschke绘制。

John K. Kruschke绘制的模型

用R语言和JAGS模拟数据

到目前为止,我们仍只停留在理论阶段。大多数情况下。后验分布并不能直接得到(想想正态分布和伽马分布有多复杂,然后还要再将他们乘起来)。Markov Chain Monte Carlo方法常常用来估计模型的参数。利用JAGS就能帮我们完成。

“等一下!!!这个工具真能够帮我们解决这些复杂的公式么?”

JAGS模型是一个基于Markov Chain Monte Carlo(MCMC)的仿真过程,它能够生成出参数空间θ=(β0;β1;τ)的许多次迭代。由参数空间中每个参数生成的样本分布会估计出参数最有可能的分布

为什们是这样呢?这个解释起来十分复杂,已经超出了本篇介绍的范围。直接来说就是:MCMC通过构建一个马尔可夫链产生了服从后验分布的样本,这个马尔可夫链有着同样的目标后验分布!?

这个过程很没意思。最快的方法就是:不去解等式(2),因为通常不可能得到解析解。我们能做的就是聪明的采样,而在数学上证明了这些样本确实服从以β0,β1,τ为参数的分布。

那么挥别了数学之后,我该怎么使用JAGS?

我们按下面的步骤在R语言中运行JAGS

首先以文本的形式写下模型

然后,我们让JAGS执行仿真模拟。这里我使用JAGs对参数空间θ进行10000次模拟

抽样之后,我们就得到了θ=(β0;β1;τ)的抽样数据,如下表所示:

“看上去好酷,那又怎样呢?”

现在我们对参数空间θ进行10000次迭代,根据等式

这就意味着,如果用x=169cm替代每个迭代值,我们将会得到10000个体重值,也就是身高169cm情况下体重的分布。

我们都对以我身高为条件下体重的百分比分布很感兴趣。为了达到这个目的,需要找到基于我身高的体重分布。

上面这张图表明我的体重(给定169的身高)最有可能在模拟越南人口中的后30%左右。

比如,我们能发现我的体重在前40%甚至更少的位置

因此大量证据表明,在身高169的条件下,体重58kg会让我处于越南人口的较低百分比处。我确实需要去健身房锻炼并增加些体重了。毕竟如果你不信详尽的贝叶斯统计分析,还能相信什么呢?

我有一份包含了美国8169名年轻男性和女性身高体重的数据(国家寿命调查1997),你能做同样的分析么,看看你会得到什么样的结论?

相关报道:

https://towardsdatascience.com/how-bayesian-statistics-convinced-me-to-hit-the-gym-fa737b0a7ac

大数据文摘
大数据文摘

秉承“普及数据思维,传播数据文化,助⼒产业发展”的企业⽂化,我们专注于数据领域的资讯、案例、技术,形成了“媒体+教育+⼈才服务”的良性⽣态,致⼒于打造精准数据科学社区。

理论贝叶斯分析
1
相关数据
参数技术

在数学和统计学裡,参数(英语:parameter)是使用通用变量来建立函数和变量之间关系(当这种关系很难用方程来阐述时)的一个数量。

最大似然估计技术

极大似然估计是统计学中用来估计概率模型参数的一种方法

后验概率技术

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

线性回归技术

在现实世界中,存在着大量这样的情况:两个变量例如X和Y有一些依赖关系。由X可以部分地决定Y的值,但这种决定往往不很确切。常常用来说明这种依赖关系的最简单、直观的例子是体重与身高,用Y表示他的体重。众所周知,一般说来,当X大时,Y也倾向于大,但由X不能严格地决定Y。又如,城市生活用电量Y与气温X有很大的关系。在夏天气温很高或冬天气温很低时,由于室内空调、冰箱等家用电器的使用,可能用电就高,相反,在春秋季节气温不高也不低,用电量就可能少。但我们不能由气温X准确地决定用电量Y。类似的例子还很多,变量之间的这种关系称为“相关关系”,回归模型就是研究相关关系的一个有力工具。

逻辑技术

人工智能领域用逻辑来理解智能推理问题;它可以提供用于分析编程语言的技术,也可用作分析、表征知识或编程的工具。目前人们常用的逻辑分支有命题逻辑(Propositional Logic )以及一阶逻辑(FOL)等谓词逻辑。

似然函数技术

在数理统计学中,似然函数是一种关于统计模型中的参数的函数,表示模型参数中的似然性。 似然函数在统计推断中有重大作用,如在最大似然估计和费雪信息之中的应用等等。“ 似然性”与“或然性”或“概率”意思相近,都是指某种事件发生的可能性,但是在统计学中,“似然性”和“或然性”或“概率”又有明确的区分。

马尔可夫链蒙特卡尔理论技术

马尔可夫链蒙特卡洛方法(MCMC)本质上是将马尔可夫链用于对蒙特卡洛方法的积分过程中,贝叶斯学派(有时亦包括频率派学者)常常需要对高维概率分布进行积分,从而推断出其模型的参数,或利用模型作出预测。传统蒙特卡洛积分的理念是从所要求的概率分布中进行抽样、得到样品平均值、并估算出期望值;而MCMC则是在抽样分布中构造出一条马尔可夫链,即建立一个随机抽样过程,弥补了传统蒙特卡洛积分过去只能静态模拟的缺陷,并为贝叶斯理论在高维分布的应用开辟了新的道路。目前有许多种算法来构造这样一条马尔可夫链,但所有的算法,甚至包括*吉布斯采样*,都是*梅特罗波利斯-黑斯廷斯算法(Metropolis–Hastings algorithm)*总体框架的特例。

回归分析技术

回归分析是一种用于估计变量之间的关系(当一个自变量变化而其它变量固定时,因变量会如何变化)的统计过程,在预测任务中有广泛的应用。回归分析模型有不同的种类,其中最流行的是线性回归和 逻辑回归(Logistic Regression)。另外还有多变量回归、泊松回归、逐步回归、脊回归(Ridge Regression)、套索回归(Lasso Regression)和多项式回归等等。随机梯度下降(SGD)就是一种起源于回归分析的常用方法,可用于控制复杂度。

马尔可夫链技术

马尔可夫链,又称离散时间马尔可夫链,因俄国数学家安德烈·马尔可夫得名,为状态空间中经过从一个状态到另一个状态的转换的随机过程。该过程要求具备“无记忆”的性质:下一状态的概率分布只能由当前状态决定,在时间序列中它前面的事件均与之无关。这种特定类型的“无记忆性”称作马尔可夫性质。

推荐文章
暂无评论
暂无评论~