贝叶斯模型平均法的基本原理及其在logistic回归中的应用实例.pdf

3 国家自然科学基金重大项目 30590374 十五 国家科技攻关项目 2004BA718B04 通讯作者 姜庆五 qwjiang 021 54237435 贝叶斯模型平均法的基本原理及其在 logistic回归中的应用实例 3 复旦大学公共卫生学院流行病学教研室 200032 张志杰 彭文祥 周艺彪 庄建林 姜庆五 提 要 目的 介绍Bayesian模型平均法的基本原理 并对实际资料进行分析 指出该方法的优越性 方法 以 Hosmer和Lemeshow研究低出生体重婴儿影响因素的队列研究为例 分别以Bayesian模型平均法和逐步logistic回归法选 择最佳模型 并分析比较二者的差异及原因 结果 Bayesian模型平均法确定的10个后验概率最大的模型的累积后验概 率仅为0159 模型本身的不确定性是很大的 而逐步logistic回归法确定的最佳模型的后验概率 P k 0 D 01 032 要 远低于Bayesian模型平均法确定的最佳模型的后验概率 P k 0 D 01 12 从回归系数的估计值 标准误和P值比 较两种方法的结果发现 Bayesian模型平均法估计的精度较高 而逐步logistic回归法由于没有考虑模型本身的不确定性 偏向于高估结果 结论 Bayesian模型平均法考虑了模型本身的不确定性 其分析结果更可靠 在统计建模中具有较好的 应用前景 关键词 logistic回归模型 Bayesian模型平均法 模型不确定性 后验概率 通常的数据分析一般是先描述资料的基本特征 然后指定某一个模型 如 logistic回归模型 通过两 阶段变量筛选法 two stage method 1 或逐步回归 法 2 确定最佳模型 并在此基础上进行参数检验和可 信区间的计算及预测 然而由于数据的部分信息被用 于人为地指定某个特定模型而导致所建立模型的效能 与预测能力均减弱 3 4 这是频率统计 即我们通常所 说的经典统计学 建立模型的一个弱点 其主要的原因 是未考虑模型本身的不确定性 model uncertainty 当用回归模型或广义线性模型对所收集的数据建立模 型时 模型本身的不确定性可能很大 5 忽略模型本 身的不确定性而仅仅以单一模型的结果进行 感兴趣 量 quantities of interest 的推断 一方面会低估 感兴 趣量 的不确定性 高估结果 导致结果的适用范围 coverge 要比研究中指定的范围 stated coverage lev2 el 要窄 6 7 另一方面会比正常情况下更趋向于拒绝 无效假设产生误导性结果 8 本文以二值应变量为 例介绍考虑了模型本身不确定性的具有较高预测能力 的Bayesian模型平均法 Bayesian model averaging 简 称BMA 然后对实际资料进行了分析 阐明了其优越 性 原理与方法 11Bayesian模型平均法的基本原理 对于应变量为0 1变量的资料 通常是使用logis2 tic回归模型进行分析 logistic回归模型通常可以表 示为 log P Y 1 P Y 0 0 n i 1 iXi Y是0 1变量 Xi i 1 2 n 是自变量 i i 1 2 n 是回归系 数 exp i 即OR值 频率统计中logistic回归模型建模方法 如 向前 变量筛选法 向后变量筛选法 逐步变量筛选法等 都 忽略了模型本身的不确定性 而Bayesian模型平均法 则是一种基于Bayesian基本原理的将模型本身的不确 定性考虑在内的统计分析方法 6 9 10 它综合考虑了 模型及参数的先验信息和样本提供的信息 其分析更 具有科学性 假设 是我们研究感兴趣的量 本文指的是OR 值或 D是我们调查所得的数据 M M1 Mn 代表所有可能的模型组成的模型空间 而哪一个模型 是最佳模型事先并不知道 即模型本身的不确定性 通常情况下 可以考虑的模型数量是很大的 例如 lo2 gistic回归模型中如果有p个自变量 那么在不考虑 交互作用的情况下可以考虑的模型数量为2 p 根据 Bayesian模型平均法的理论 的后验分布为 P D n i 1 P Mi D P Mi D i 1 2 3 n 1 其中 P Mi D 是在给定数据D和模型Mi的条 件下 的后验分布 P Mi D 是在给定数据D的情 况下Mi为最优模型的概率 从方程 1 可以看出 的 后验分布实际上是以后验模型概率P Mi D 为权 重 对所有模型的后验分布进行加权的一个平均值 在 1 中 假定Mi为最优模型的情况下 的预 测分布为 P Mi D P Mi D P Mi D d 0 1 n 为模型Mi回归系数的向 量 764 中国卫生统计2007年10月第24卷第5期 在 1 中 根据Bayesian原理知模型的后验概率P Mi D 为 P Mi D Pl D Mi P Mi n k 1 Pl D Mk P Mk Pl D Mi P Mi 2 其中 P Mi 是模型Mi为最优模型的先验概率 Pl D Mi 是模型Mi的边际似然 marginal likelihood of model Mi 它是通过对未知参数的积分来获得 Pl D Mi Pl D Mi P d Mi d 3 0 1 n 为模型Mi回归 系数的向量 在 3 中 是模型Mi的参数向量 Pd Mi 是 给定模型Mi的条件下 的先验密度 prior density 即 参数的先验信息 Pl D Mi 是给定模型Mi和参 数向量 的条件下数据D的似然 通过 1 2 和 3 的计算即可得到 本文指OR 值 的后验分布 然后进行相应的假设检验 点估计以 及可信区间等的计算 21Bayesian模型平均法的实现 Bayesian模型平均法的真正实现涉及了模型空间 中所有模型的和式计算 1 和积分运算 3 模型空间 中所有模型的数量是相当大的 如 自变量个数为15 个 那么即使在不考虑交互作用的情况下其模型的数 量可达215 32 768 本文采用Occam窗 10 的方法来 适当缩减模型的数量 当一个模型的后验概率小于最 佳模型后验概率的5 时则从模型空间中去掉该模 型 该方法经验证与全模型空间的运算具有相似的效 能 具体算法见文献 10 积分运算 3 由于常常不能 得到可分析形式的表达式导致其计算通常无法实现 本文采用Laplace近似法来得到其近似值 4 10 log P l D Mi log P l D Mi pilogn O 1 4 其中 是模型Mk的参数向量 的后验均数 pk是模 型Mk中的参数个数 n是样本量 这就是所谓的 Bayesian信息准则 BIC 4 10 Bayesian模型平均法还需要指定模型先验概率和 参数的先验信息 为了不对模型的优劣性加入主观色 彩 本文指定相等的模型先验概率即均匀分布的模型 空间 对于参数的先验信息 我们指定单位信息先验 unit information prior 如 多变量正态先验 以最大 似然法估计其均数 以观察值的期望信息矩阵作为方 差的估计 这种先验信息的指定方法由Raftery 4 10 提 出并证明了其合理性 它可以看作是一种包含了与原 始观察值等量信息的一种先验分布 通过R软件实 现 5 Bayesian模型平均法 31Bayesian模型平均法的统计推断 像经典统计分析一样 1 的后验概率提供了一种 在考虑了模型不确定性的基础上对校正的OR值进行 假设检验 点估计以及可信区间的计算方法 11 从 Bayesian的观点看Bayesian模型平均法的假设检验 其对应的备择假设H1为 Xk为某疾病的危险因素的 后验概率有多大 即 k 即logORk 不等于0的后验概 率有多大 Bayesian模型平均法以包含Xk的所有模 型的后验概率的和作为 k不等于0的后验概率的估 计 P k 0 D Mi A P Mi D Ii k 0 5 其中 A表示经Occam窗方法缩减后的模型空间 Ii为 0 1指示变量 当 k在模型Mi中时 Ii 1 当 k不在 模型Mi中时 Ii 0 本文以 P k 0 D 表示Bayesian模型平均法的 k的后验概率 其与频率统计的P值推断方法不同 对于Bayesian模型平均法后验概率推断的一般规 则 12 为 P k 0 D 015表示没有证据表明Xk是某疾 病的危险因素 015 P k 0 D 0175表示有弱的证据表明 Xk是某疾病的危险因素 0175 P k 0 D 0195表示有强的证据表明 Xk是某疾病的危险因素 P k 0 D 0195表示有非常强的证据表明 Xk是某疾病的危险因素 如果认为Xk是某疾病的危险因素 那么就需要 对其效应的大小进行相应的估计 Bayesian模型平均 法中Xk的点估计的后验均值和后验方差 13 14 分别 为 E k D n i 1 kP Mi D 6 Var k D n i 1 Var k Mi D 2 k P Mi D E k D 2 7 其中 k E k Mi D 关于 k即Xk效应估计的推断 如 假设检验 点 估计和区间估计 都来源于 5 6 和 7 的计算 对于 某一特定协变量效应大小的估计类似于对经典统计的 方法 Bayesian模型平均法采用的是参数的后验均数 不同的是它考虑了模型本身的不确定性 所有的计算 结果可以在R中自动获取 具体操作见参考文献 5 本文之所以要引进Bayesian模型平均法的模型选 择方法 主要是因为Bayesian模型平均法的两个很重 要的特性 1 考虑了模型本身的不确定性 因此在模 型的预测能力上要优于任何一个频率统计所选择的单 一模型 9 而很多资料的分析往往需要预测能力较高 864 Chinese Journal of Health Statistics Oct 2007 Vol 24 No 5 的模型 2 统计推断得到了很好地校正 如 可信区间 和 值 15 都考虑了模型的不确定性 使得结果精度更 高 统计推断的风险性更低 实例分析 我们对Hosmer和Lemeshow队列研究的数据进 行了再分析 比较了Bayesian模型平均法与逐步logis2 tic回归模型在分析结果中的差异 指出Bayesian模型 平均法的优越性 11 数据集简介 Hosmer和Lemeshow 16 研究了低出生体重婴儿 的危险因素并希望对低出生体重婴儿的概率进行预 测 他们收集了美国医院189名孕妇及其新出生的婴 儿资料 共收集了9个变量的信息 Venables和Ripley 17 对该数据集通过逐步logis2 tic回归模型进行了再分析 为了与其结果比较 以上 数据与Venables和Ripley的分析数据集保持一致 同 时按照Venables和Ripley的分析方法将孕妇种族和 孕妇前三个月体检次数分别转变成2个哑变量孕妇种 族1 孕妇种族2和孕妇前三个月体检次数1 孕妇前 三个月体检次数2 增加了3个一阶交互作用项作为 新的变量 孕妇年龄与体检次数1的交互作用 孕妇年 龄与体检次数2的交互作用及孕妇怀孕期间的吸烟史 和子宫过敏的交互作用 相当于数据集中共13个自变 量 变量的处理方法同文献 17 Bayesian模型平均 法首先考虑自变量所有可能的组合 产生初始的213 8 192个模型的模型空间 本文通过 leaps and bounds 算法 18 19 和Occam窗的方法适当减少模型数量后进 行假设检验 点估计和可信区间的计算 21Bayesian模型平均法与逐步logistic回归模型 分析结果的比较 1 两种方法确定的最佳模型后验概率的比较 在Occam窗中共筛选了42个模型用来估计 Bayesian模型平均法的参数估计值 表1列出了 Bayesian模型平均法获得的模型后验概率最大的10 个模型及所选择的变量和Bayesian信息准则 Bayesian Information Criterion 简称BIC BIC值越小 模型的拟 合效果越好 从表1可见最佳模型的后验概率仅为 0112 前10个模型的累计后验概率为0159 由此可 见 模型的不确定性 model uncertainty 在该数据集中 是相当大的 而逐步logistic回归模型确定的最佳模型 不在前10个后验概率最大的模型中 其模型的后验概 率小于01032 从模型的后验概率角度看 Bayesian模 型平均法确定的最佳模型要优于逐步logistic回归模 型 表1 Bayesian模型平均法得到的前10个后验概率最大的模型 模型 孕妇 年龄 孕妇最后一次 月经时的体重 孕妇 生育史 孕妇高 血压史 孕妇子宫 过敏史 孕妇前三个月 体检次数1 孕妇前三个月 体检次数2 孕妇年龄与体检次 数1的交互作用 Bayesian 信息准则 后验 概率 模型1 759160001120 模型2 758198501088 模型3 758191701085 模型4 758131001063 模型5 757166701046 模型6 757164401045 模型7 757146801041 模型8 757130301038 模型9 757104801034 模型10 756104801032 3 11 表示该变量在相应的模型中 21 孕妇年龄与孕妇前三个月体检次数2的交互作用 孕妇怀孕期间是否吸烟与子宫是否过敏的交互作用及孕 妇怀孕期间是否吸烟三项未出现在前10个后验概率最大的模型中 故未在表中列出 2 两种方法确定的参数估计值 标准误和P 值的比较 Bayesian模型平均法与逐步logistic回归模型确 定的参数估计值 标准误和P值的结果见表2 从参数的点估计看 Bayesian模型平均法估计的 回归系数与逐步logistic回归模型估计的回归系数相 比 Bayesian模型平均法估计的回归系数更趋向于0 这是Bayesian统计分析的一个特点 主要是因为 Bayesian模型平均法考虑了模型本身的不确定性 在 本例的分析中甚至发现孕妇年龄 孕妇前三个月体检 次数1 孕妇怀孕期间是否吸烟和子宫是否过敏的交 互作用项的系数估计是反向的 二者的P值分别为 011616 012013和010489 而后验概率分别为01171 01127和01018 对于孕妇年龄和孕妇前三个月体检次 数1的推断两种方法一致 可以认为它们与新生儿的 低出生体重无关 其反向的系数是由抽样误差引起 而 对于孕妇怀孕期间是否吸烟和子宫是否过敏的交互作 用项则不一致 逐步logistic回归模型认为该交互作用 项有统计学意义但由于P值接近0105而需要进一步 的研究 而Bayesian模型平均法则认为该交互作用项 没有统计学意义 P i 0 D 01018 关于它们的 交互作用一方面在近年来的研究表明二者间不存在交 964 中国卫生统计2007年10月第24卷第5期 互作用 24 25 另一方面从理论上也支持二者间不存在 交互作用 因为子宫过敏是人体的 固有特性 在无干 预的情况下 不会随着吸烟状态的变化而变化 所以可 以认为Bayesian模型平均法的结论比逐步logistic回 归模型要准确 更能发现细微的差别 特别是当因素的 效应处于 边缘状态 的情况 从参数估计的标准误来看 Bayesian模型平均法 估计的绝大多数的标准误均小于逐步logistic回归模 型估计的标准误 其相应的可信区间更窄 精度更高 这主要是因为Bayesian模型平均法在估计标准误时考 虑了模型本身的不确定性 因而结果更精确 从P值的角度分析 逐步logistic回归模型认为 孕妇怀孕前最后一次月经时的体重 孕妇生育史 孕妇 高血压史 孕妇子宫过敏史 孕妇前三个月体检次数 2 孕妇年龄和孕妇前三个月体检次数2的交互作用六 项与新生儿的低出生体重有关 而Bayesian模型平均 法仅认为有很强的证据表明孕妇生育史与新生儿的低 出生体重有关 而孕妇怀孕前最后一次月经时的体重 和孕妇高血压史的后验概率分别为01485和01530 仅有比较弱的证据表明其可能与新生儿的低出生体重 有关 如果这两个变量是研究的重要变量 那么 Bayesian模型平均法的结果则提示我们需要对这两个 变量做进一步地研究 其余变量的后验概率均小于 20 Bayesian模型平均法认为没有证据认为其与新 生儿的低出生体重有关 逐步logistic回归模型中P 值在0105附近的有孕妇怀孕期间的吸烟史以及孕妇 怀孕期间是否吸烟与子宫是否过敏的交互作用两项 认为它们可能与新生儿的低出生体重有关 但需要进 一步的研究 而Bayesian模型平均法则表明它们与新 生儿的低出生体重无关 后验概率均低于011 逐步 logistic回归模型中P值 大 于0105的 结 果 则 与 Bayesian模型平均法的结果一致 而近年来的研究结 果更支持Bayesian模型平均法的结论 26 27 我们认为 逐步logistic回归模型由于没有考虑模型本身的不确 定性 偏向于高估结果 容易导致误导性的结论 表2 Bayesian模型平均法与逐步logistic回归模型确定的参数估计值 标准误和P值 变量 逐步logistic回归模型 回归系数标准误P值 Bayesian模型平均法 后验均数后验标准误后验概率 孕妇年龄010755010540011616 010111010284011710 孕妇最后一次月经时的体重 010204010075010066 010077010093014850 孕妇种族1 010072010970010070 孕妇种族2 010030010482010070 孕妇怀孕期间的吸烟史017800014204010635010451011764010910 孕妇生育史115603014970010017115145014745019930 孕妇高血压史210657017487010058019192110213015300 孕妇子宫过敏史118185016676010064010958012962011310 孕妇前三个月体检次数1219211212858012013 010265015193011270 孕妇前三个月体检次数2912449216615010005017695211554011270 孕妇年龄与体检次数1的交互作用 011618010968010946 010130010278013170 孕妇年龄与体检次数2的交互作用 014110011191010006 010349010963011530 孕妇怀孕期间的吸烟史和子宫过敏的交互作用 119167019731010489010092011095010180 讨 论 在资料的处理过程中 我们经常需要通过统计技 术来控制混杂因素的影响 然后在此基础上建立模型 并进行统计推断 然而频率统计的模型筛选方法 如 向前法 向后法 逐步法等 一方面由于涉及了模型间 的多重检验 模糊地讲是非嵌套模型间的比较 因而导 致了其抽样属性非常复杂 20 另一方面 P值的使用 近年来引起了很大的争议 10 如 的确定没有一 定的标准方法可以参考 带有很大的主观性 当样本 量很大时 P值甚至明显的拒绝拟合效果非常好的模 型 而当样本量较小时 P值基于的统计量的分布近似 效果很差 结果不可靠 有模拟研究表明 11 频率统计 的模型选择方法 如 逐步法 所确定的最佳模型及通 过P值判断某变量是否有统计学意义具有较大的偏 差 假阳性率 高达51 认为其结果具有较大地误 导性 Madigan和Raftery 7 9 以对数评分为评价准则 logarithm scoring rule 研究并得出了Bayesian模型 平均法获得的最佳模型在预测能力要优于任何一个频 率统计所确定的单一模型 已经通过统计模拟技术得 到了初步验证 本文使用Bayesian模型平均法和逐步logistic回 归模型法 17 对同一数据集进行分析比较 发现 Bayesian模型平均法不仅可以像频率统计中逐步lo2 gistic回归模型那样提供变量是否具有统计学意义的 信息 而且可以准确地提供其后验概率的大小 这比 一刀切 如 01 05 的频率统计方法能提供更多更 实用的信息 而且由于其考虑了模型本身的不确定性 校正了相应参数的点估计值和可信区间 使得结果更 加可靠 统计分析的结论多是从概率的角度来进行解 074 Chinese Journal of Health Statistics Oct 2007 Vol 24 No 5 释的 而Bayesian模型平均法的基本原理正是基于后 验概率的思想而发展起来的 其结果更加符合其概率 思想 从文献追溯来看 最早涉及Bayesian模型平均法 雏形思想的开创性工作是在1965年 Robert建议使用 一种类似于Bayesian模型平均法思想的加权平均的后 验分布来进行研究工作 Leamer在1978年详述了 Robert的思想 给出了Bayesian模型平均法的基本范 例 并从理论上指出了将模型的不确定性考虑在内的 模型筛选方法 遗憾地是 Bayesian模型平均法在算 法的实现上却很难 因此没有引起足够的重视 而忽略 模型的不确定性所导致的缺点却不断的被发现 21 随着新的统计计算方法 如 Markov chain Monte Car2 lo 的出现 Bayesian模型平均法在社会学领域首先得 到了一定地发展 并在流行病学 药学 社会学 经济学 以及 生 态 学 等 领 域 都 开 始 了 初 步 的 应 用 研 究 9 19 22 23 Bayesian模型平均法不仅可以应用在lo2 gistic回归模型中 还可以应用在线性回归模型 广义 线性模型和生存分析等方法中 并且其算法已经在 R 5 和S Plus软件中得到了实现 我们相信 随着对 Bayesian模型平均法方法的深入研究 其必将在统计 方法学的研究中发挥越来越重要的作用 Bayesian Model Averaging and its Application in logistic Regres2 sion Model Zhang Zhijie Peng Wenxiang Zhou Yibiao et at Department of Epidemiology School of Public Health Fudan U2 niversity 200032 Shanghai Abstract Objective To introduce the basic theory of Bayesian model averaging and explain its superiority by analysing an actual example Methods To take the cohort study of Hosmer and Lemeshow as an exam2 ple we analyzed it with stepwise logistic regression and Bayesian model av2 eraging respectively compared their differences and explained their causes Results Bayesian model averaging s cumulative posterior probability of best 10 models is only 0 59 which means the model uncertainty is very large The posterior probability of best model acquired by stepwise logistic refression P k 0 D 0 032 is less than that of best model acquired by BMA P k 0 D 0 12 From their comparison of point estima2 tion standard error and P values We find that BMA has a narrower confi2 dence interval a higher precision of parameter estimation while stepwise lo2 gistic regression model tends to underestimate model uncertainty leading to overconfident inferences and decisions that are more risky than one thinks they are Conclusion Bayesian model averaging takes into account model uncertainty and it shows a promising prospect in statistical modelling Key words logistic regression model Bayesian model aver2 aging Model uncertainty Posterior probability 参 考 文 献 11Mickey RM Greenland S The impact of confounder selection criteria on effect estimation American Journal of Epidemiology 1989 129 1 125 137 21Goodman LA The analysis of multidimensional contingency tables Step2 wise procedures and direct estimation methods for building models for multiple classifications Technometrics 1971 13 1 33 61 31Lohman T Advance in Body Composition Assessment Current Issues in Exercise Science Human Kinetics Publishers Champaign IL 1992 91 41Raftery AE Approximate Bayes factors and accounting for model uncer2 tainty in generalized linear models Biometrika 1996 83 2 251 266 51Raftery AE Painter IS Volinsky CT BMA An R package for bayesian model averaging R News 5 2 2 8 November 2005 61Hoeting JA Raftery MD Volinsky AE Bayesian Model Averaging A Tutorial Statistical Science 1999 14 4 382 401 71Raftery AE Madigan D Hoeting JA Bayesian model averaging for linear regression models Journal of the American Statistical Association 1991 92 437 179 191 81Freedman DA A note on screening regression equations American Statis2 tician 1983 37 2 152 155 91Madigan D Raftery AE Model selection and accounting for model uncer2 tainty in graphical models using Occam s window Journal of the Ameri2 can Statistical Association 1994 89 428 1535 1546 101Raftery AE Bayesian model selection in social research with Discus2 sion In Sociological Methodology 1995 Marsden PV ed Blackwell Publishers Cambridge Mass 111 196 111Viallefont V Raftery AE Richardson S Variable selection and Bayesian model averaging in case control studies Statistics in Medicine 2001 20 3215 3230 121Kass RE Raftery AE Bayes factors Journal of the American Statistical Association 1995 90 430 773 795 131Raftery AE Bayesian model selection in structural equation models In Testing Structural Equation Models K Bollen and J Long eds 1993 163 180 141Draper D Assessment an