（概率论与数理统计专业论文）gplsi模型的bayes估计.pdf

中文摘要 摘要 1 作为广义线性模型和非参数回归模型的深入推广，广义部分线性单指数 模型在现代统计中有重要的作用本文利用自由节点的b a y e s 样条技术对广义 部分线性单指数模型进行统计分析我们用截幂样条函数来逼近未知函数，但 样条节点的个数和位置均视为随机变量，而由数据去自动寻找建议的算法通 过可逆跳m c m c 方法来实现本文还设计了一个有效的随机游动m e t r o p o l i s 算法以更新模型中的单指数参数本文方法的有效性通过仿真的例子被验证， 并被应用于分析一个实际数据例子 关键词：截幂样条；m e t r o p o l i s - h a s t i n g s ；可逆跳m c m c a b s t r a c t a sag e n e r a l i z a t i o no fg e n e r a l i z e dl i n e a rm o d e la n dn o n p a r a m e t r i cr e g r e s - s i o nm o d e l ，g e n e r a l i z e dp a r t i a l l yl i n e a rs i n g l e - i n d e x ( g p l s i ) m o d e lp l a y sa n i m p o r t a n tr o l ei na d v a n c e ds t a t i s t i c s b a y e s i a ns p l i n et e c h n i q u ew i t hf r e e - k n o t si su s e dt oa n a l y z eg p l s im o d e li nt h ep r e s e n tp a p e r 。w ea p p r o x i m a t e t h eu n k n o w nf u n c t i o nb yt h et r u n c a t e dp o w e rb a s i ss p l i n e ，a n dt r e a tt h en u i i l - b e ra n dl o c a t i o n so fk n o t sa sr a n d o mv a r i a b l e st h a tw i l lb es e a r c h e da u t o - m a t i c a l l yb yt h ed a t a t h ep r o p o s e da l g o r i t h mi sc a r r i e do u tb yt h er e v e r s i b l e i u m pm a r k o vc h a i nm o n t ec a r l os a m p l e r w 色a l s op r o v i d e a ne f f i c i e n tr a n d o m w a l km e t r o p o l i sa l g o r i t h mt ou p d a t et h es i n g l e - i n d e xv e c t o r t h ep r o p o s e d m e t h o di sv e r i f i e db ys i m u l a t i o n ，a n di sa p p l i e dt oar e a ld a t a k e y w o r d s ：g p l s i ；b a y e s i a ns p l i n e ；r j m c m c 厦门大学学位论文原创性声明 本人呈交的学位论文是本人在导师指导下，独立完成的研究成果。 本人在论文写作中参考其他个人或集体己经发表的研究成果，均在文 中以适当方式明确标明，并符合法律规范和厦门大学研究生学术活 动规范( 试行) 。 另外，该学位论文为() 课题( 组) 的研究成果，获得() 课题( 组) 经费或实验室的 资助，在() 实验室完成。( 请在以上括号内填写课 题或课题组负责人或实验室名称，未有此项声明内容的，可以不作特 别声明。) 声明人( 签名) ： 年月日 厦门大学学位论文著作权使用声明 本人同意厦门大学根据中华人民共和国学位条例暂行实施办法 等规定保留和使用此学位论文，并向主管部门或其指定机构送交学位 论文( 包括纸质版和电子版) ，允许学位论文进入厦门大学图书馆及 其数据库被查阅、借阅。本人同意厦门大学将学位论文加入全国博士、 硕士学位论文共建单位数据库进行检索，将学位论文的标题和摘要汇 编出版，采用影印、缩印或者其它方式合理复制学位论文。 本学位论文属于： () 1 经厦门大学保密委员会审查核定的保密学位论文， 于年月日解密，解密后适用上述授权。 () 2 不保密，适用上述授权。 ( 请在以上相应括号内打“ 或填上相应内容。保密学位论文 应是已经厦门大学保密委员会审定过的学位论文，未经厦门大学保密 委员会审定的学位论文均为公开学位论文。此声明栏不填写的，默认 为公开学位论文，均适用上述授权。) 声明人( 签名) ： 年月日 第一章引言 假设y l ，y n 是独立且均服从自然指数族分布 f气 p ( 统l 晚，妒) = e x p a 7 1 ) ( 坎或一b ( e d ) + e ( 坎，) ， t = 1 ，纯。 ( 1 1 1 ) lj 其中0 i 为自然参数，妒为一个不依赖于i 的公共尺度参数，锄( ) ，6 ( ) 和c ( ) 为给定的函数，且倪( 钟一触，峨代表权重本文取峨= 1 ，i = 1 。， 假设满足 、 统= g 萝8 f 霪1 ) 罗? z i ， i l ，站 ( 1 1 。2 ) 其中。l ，z i 分剐为d 维和蟊维预测向量，8 = ( a l ，& 罐) 芏和卢= 够l 。，勉) ? 分别是d 维单指数向量和也维回归系数向量，g 是未知的一元实值函数， 而g 是已知的单调可微函数为了满足函数的可识别性，通常假设单指数变 量盘是单位南量且8 l 0 + 由式子( 1 1 1 ) 和( 1 1 。2 ) 所确立的模型嘲广义部 分线性单指数模型【l 】，简称为g p l s i 模型特别地，当g 为恒等函数，即 0 i = g ( c l t x t ) + z t z i 时，称为自然或典则形式 g p l s i 模缝是一类非常广泛的模垄，它既包含了经典的广义线性模型， 如线性避归模型，l o g i t 模型，p r o b i t 模型，p o i s s o n 回归模型，g a m m ai i i i 归模型和一些比例风险模氆【2 】等；叉包含了许多非参数模型，如部分线性模 型，单指数模型等 众所周知，单指数模型在统计学和经济学上有广泛的应用，其本质是用 未知的一元函数来代替已知的联系函数，因而比广义线性模型【2 】更具有灵活 性当卢= 0 ，即没有协变量z ，模型( 1 1 2 ) 只是一般的具有未知函数的广义线 性模型。丽w e i s b e r g 穰w e l s h ( 1 9 9 4 ) 跫经解决了广义线性模型中“未知联系蘧 数”的问题1 3 1 ；这些模型引进指数q t x ，解决了在拟合多元非参数回归函数时 的“维数泛滥”问题( 尽管降低了灵活性) 。 对革指数8 需要进行简单的解释，这墨引用l i ( 1 9 9 1 ) 的三点论点f 4 】： g p l s l 模型的b a y e s 估计 4 ( a ) 显然的，在实际问题中，在拟舍数据之前降低维数是很重要的( b 的推论 1 2 1 4 】给出了更完整的的论述且说明在很多情况下这是很关键的一步) ，而 单指数模型的优点为降低维数提供了一个容易解释的方法 ( b ) 如果箩是单调的，那么单指数毡可以跟一般线性模型串的参数一样解 释为“影响”因子。 ( e ) 给定了a 这个单指数向量后，模型的优劣评判就更容易了。 对广义部分线性单指数模型，我们需要估计单指数a ，回归系数移和拟合 未知函数g 在经典的统计中，对于单指数和回归系数，可以用平均导数方法 5 6 7 ( s t o k e r ，1 9 8 6 ；h i i r d l e 和s t o c k e r ，1 9 8 9 ；h o r o w i t z 和h s r d l e ，1 9 9 6 ) 或半参数矩 估计方法来估计；褫对于未知瑟数，可以爰局部多项式圆归 1 8 9 ( i c h i m u r a ，1 9 9 3 ； h s r d l ee ta 1 ，1 9 9 3 ；c a r r o l le ta 1 ，1 9 9 7 ) 或是样条回归方法f l u 】( y u 和r u p p e r t ，2 0 0 2 ) 来拟合然而，局部多项式回归和样条回归的关键是选择一个合适的带宽或平 滑参数。虽然已经有一些方法可以鳃决此类问题，徨对手小样本，所得到的估 计并不稳定因此，我们可以用b a y e s 方法进行估计 对于单指数模型，a n t o n i a d i s ，g r 占g o i r e 和m c k e a g u e ( 2 0 0 4 ) 提出了一种半 炙叶斯分析方法f l l 】，鄢用随机游动m e t r o p o l i s 算法估计单指数，丽焉惩罚祥条 方法估计未知联系函数，其中的平滑参数根据广义交互有效性准则来确定。 w a n g 1 2 】利用自由节点样条技术提供了一个完全的b a y e s 估计方法 本文用囊由节点的b a y e s 样条技术分析广义部分线性单指数模型，包括估 计单指数a 和回归系数p ，拟合未知函数g 。 具体地，第二章首先构造广义部分线性单指数模型的b a y e s 框架在节点 给定的条件下，用截幂祥条逼近未知函数，导出了模型的似然函数但由于实 际问题中，样条节点个数和位置通常是未知的，也需要被估计因此我们通过 设置所有未知量的先验分布，获得未知量的联合后验分布第三章具体介绍了 用自由节点的b a y e s 祥条技术分析广义部分线性单指数模型的步骤在第四章 中，我们用i a i l - 仿真例子验证了所建议方法的有效性第五章应用建议的方法 分析一个实际的例子第六章对本文的研究内容进行简单的总结 g p l s i 模型的b a y e s 估计 第二章g p l s i 模型的b a y e s 分析 5 2 1模型描述 一般地，未知的一元函数夕( ) 可以用截幂样条 1 3 1 1 4 来逼近假设 g ( u ) = 5 0 + 6 1 t 正+ + 妨矿+ 叁1 + t ( t 正一n ) 华， t 正 鼬，6 h 】( 2 1 1 ) 其中p 是多项式的最高阶，让+ = m a x ( u ，o ) ，七代表样条节点个数，r = ( r l ，仉) t 代表节点位置且满足r 1 他，和k 分别代表集合 口t ，i = 1 ，礼) 的最小值和最大值 定义样条系数向量 6 = ( 6 0 ，5 l ，k ) t ， 和样条基向量 b ( 札) = ( 1 ，仳，矿，( 仳一r 1 ) 军，( 札一。，p ) t ， 则函数可以表示为 g ( u ) = b t ( u ) 6u a o t ，6 h 】， 进而自然参数仇表示为 o i = gb t ( c r t x i ) 6 + 矿旎) ，江1 n ( 2 1 2 ) 从而，模型的似然函数( 1 1 1 ) 可以表示为 p ( 可i q ，p ，疋k ，r ) = 唧 0 是一个标量超参数，y o 一( y 0 1 ，) t 是礼维超参数向量 2 节点个数k 的先验密度 对节点个数k ，我们用超参数为a 的p o i s s o n 分布作为先验分布，即 啡) = 筹e ，k k m 黻 其中七m 觚为最大节点个数 g p l s i 模型的b a y e s 估计 7 3 节点位置r 的先验密度 由于k 个节点 r l 他 6 口将区间( a a ，k ) 分解成k + 1 个 子区间( n 一1 ，r 1 ) ( f = l ，2 ，k + 1 ) ，每一个子区间占整个区间长度的比值为 而1 1 - - r l - - 1 ，通常假定该比值服从d i r i c h l e t ( 1 ，l ；1 ) 分布，即( r l ，r 2 ，他) 为 区间( a a ，k ) 上均匀分布的次序统计量 f ( r = 面 其中a k = i l a a r 0 ，因此我们用d 维半单位球面上 的均匀分布作为先验分布，也就是 巾) = 丌由( 罢) ， 其中r ( ) 是r 函数 2 3 联合后验密度 由b a y e s 定理，所有未知量的联合后验密度为 p ( q ，p ，6 ，k ，ry ，a o ，y o ) o ( p ( y l q ，卢，5 ，k ，r ) 丌( 6 ，p i q ，k ，r ) 7 r ( r i a ，k ) r ( a l k ) r ( k ) e x p 丝竺竖生盟竿塑幽型生型 l痧i 唧 詈k g ( 巩酬6 + 刎- - j t b ( g ( 巩酬6 + 刎) ) d 每a kx ，r - 瓤争蔷e 以 唧 芋 ( 等署) t g ( 巩x ? ) 6 + 刎- j t b ( g ( 巩酬巧+ 刎) 灯讯罢，( 击) 七e 以札， 仁3 1 ) g p l s i 模型的b a y e s 估计8 由上面的式子可以发现，联合后验密度的形式极其复杂，直接从联合后验 分布推断参数的全条件后验分布是非常困难的，因此我们采用r j m c m c 算法 产生参数后验分布的样本，以达到推断参数的磊的。 2 4 超参数n o ，珈的选取 本节介绍共轭先验分布( 2 2 1 ) 中的超参数a 0 ，y o 的取值方法 如d i a c o n i s 和y l v i s a k e r 1 6 ( 1 9 7 9 ) 所描述的，y o e ( 印) ) ，其中占p ) 是备矽) 的切向量，且期望由先验分布( 2 2 。1 ) 获得。由予在g p l s i 模型中，e ( y 1 8 ) = 6 汐) ， 我们有 e ( y ) 一e e e ( y l o ) j e ( b ( o ) ) 一y o ( 2 4 1 ) 因此裳望( 2 4 。1 ) 表瞬y o 是y 的边缘均值，丽曼可以看做是e ( y ) 酶先验预测 参数a o 可以看成先验样本大小在本文中，可以表示成罂，其中n o 是先验样 本量，即先验信息量其实，只要a o 足够小，系数向量( j ，p ) 的共轭先验分布 ( 2 2 1 ) 相当于无信息先验 g p l s i 模型的b a y e s 估计 9 第三章模型估计 由于无法完全导出参数的全条件后验密度，我们用可逆跳m a x k o vc h a i n m o n t ec a r l o ( 筒记为鲇m c m c ) 算法来分析模型，傣计参数( q ，反鑫毒，尹) 这个 过程需要经过以下三个步骤： ( a ) 对所有参数o ，熊正k ，r 赋予初值； ( b ) 固定p ，疋七，r 的值，通过m e t r o p o l i s - h a s t i n g s 算法更新指数向量血； ( c ) 固定q ，通过r j m c m c 算法更新p ，正k ，r 反复进行( b ) ，( c ) 步，直至收敛；最后，我们用后验样本均值来估计所有的未 知量 3 1 可逆跳m c m c 算法的描述 首先简单介绍一下r j m c m c 算法r j m c m c 抽样器是m e t r o p o l i s - h a s t i n g s 算法【1 7 l ( t i e r n e y , 1 9 9 4 ) 的推广，它允许改变参数向量的维数我们用p ( 惫，护( 忍) ) 来装示参数k 和拶鳃睽合后验分南，其中k 代表模型指标，拶是铭惫维参 数向量假设当前的m a x k o v 链状态是( k ，p ( 知) ) ，在不同的维数空间之间的转变 有不同的移动类型r j m c m c 抽样器按如下的步骤进行抽样： ( a ) 以概率分布五( r ) 从模型k 到模型中选取移动类型r ( b ) 从建议分布p k k ，( 乱1 分砷) 抽取牡 ( c ) 令p ，“) = h k ，j ，p ( 斛，乱) ，其中h 是已知的可逆函数，且n k + l u f = n k ，+ l u ，i 代表囊量程的维数 ( d ) 为了维持m a x k o v 链蛉平衡，我们以概率r a i n ( 1 ，a ) 接收( k l , 8 ( ) ) ，其中 a = 糟筹撩| 帮 涵n ， = p o s t e r i o rr a t i o p r o p o s a lr a t i o j a c o b i a n g p l s i 模型的b a y e s 估计 1 0 对于模型之间的转换，最简单的方法是从现有的模型中加入或减少参 数而得到新的模型不失一般性，假设n k ， n k 和i u l = n k ，一n 南，那么从高 维模型到低维模型的转换就完全确定了，丽接受率a ( 3 1 1 ) 可以简化为 a = 黼顽淼l 帮1 3 。2 更新单指数巍量o t 假定o t 憋维数不低于2 。我们通过对现有的单指数向量进行随枧歪交变换 得到候选单指数向量i 蚴首先，随机的选择a 的两个不同元素，记为o t t 和 ，然后从区间卜7 r ，州上的截断正态分布0 9 n ( 0 ，口 ) + o 1 n ( 0 ，1 ) 中随机抽取值 霹。通过以下鳇变换来分别得到候选晦量的第i 稚第j 个元素 8 ；= a t c o s t + s t a r 和 = 一a ts i n + c o s ，7 ， 并保持其它元素不变这是典型的随机游动m e t r o p o l i s 算法，即与q ：是 瓴，a j 通过澍绕以零为中心、以任意角度旋转丽得到其中，超参数( 7 1 可以 调节，以产生一个合理的接受概率最后，以概率r a i n 1 ，a l l 接受候选向量 其中 p ( a 4 f 声，最惫，r ，a 0 ，y o ，y ) a 1 2 p ( c v l 3 , 5 , k , r , a o , y o , y ) 唧 学 ( 帮) 丁即t ( x a v + 刎一，6 ( g ( 巩榭) 6 + 刎) ) e x p 学l ( 黹) t c ( b t ( x a ) 6 十z p ) 一j t b ( g ( b t ( x a ) 6 + 邵) ) l = 1 p _ r 1 r 一， 虽然在我们更新单指数向量的时候固定了节点个数，但是候选节点位置可能 随着候选单指数向量的变化而变化我们通过以下方式来重新确定候选节点 位置广= ( r 4 1 ，幢) t t 疗一鼬 r i n a b a * - a a * 5 石= i ，魄一鼬 g p l s i 模型的b a y e s 估计 l l 另一方面，如果接受了候选向量矿，但n 1 0 的假设不再满足，那么我们分 别对q + ，矿，a o 和b 乘以1 ，并对节点位鬣进行重新排序 3 3 更新节点个数惫和位置7 以及参数鼠p 在现有的节点个数和位置基础上，以概率k c m i n l ，掣铲) 增加个 新节点；以概率呶= c - 幽 1 ，帮) 删除一个原有节点；以概率叼k = l “奄一毗 移动现有节点中的一个显然，有 b k 7 r ( 奄) = d k + l 7 r ( 奄+ 1 ) 3 3 1 建议分布的选取 首先介绍参数艿和多的建议分布的基本选取方法+ 由于鑫和声的全条件 后验分布并非为一个标准的分布，因此，我们利用l a p l a c e 近似【1 9 】来得到近 似的分布具体地，由参数民p ，惫，r 的联合后验 l l o p 珂，i “，耵，o o ，y o ) 鼍p ( 剪l a ，p ，正k ，r ) 丌( 瓦p j o ，k ，r ，a o ，y o ) 竹( r i 貔，k ) 万( 惫) 一c 唧 吾p r g ( 巩圳艿+ 钠一户粥( b ? 胁) 6 十刎) ) 唧 詈陋g ( 矾酬艿+ 钠- - j t b ( g ( 巩酬6 + 刎) 一c 唧 芋 ( 背) 瑚b r 叫十刎- j t b ( g ( b 丁降叫+ 刎) ， i ( 最厣，毫，r l 醴，箩，a o ，y o ) 。一l n l ( 5 ，反奄，r 睁，y 。a o ，y o ) 一a o + l 加妒( x 学+ 刎- - j t b ( g 妒( x 叫+ z 3 ) ) 十c ， 其中雾= ( 帮) 。 等式两边关于0 l = 酽，多t ) r 求导碍到 g p l s i 模型的b a y e s 估计 1 2 丽o l = 丁a o + l ( e l ，e 疹d i a g ( 酥孤) ) 雪- ( e 小，址d i a g ( g ( e 融) ) 蛾鸱( ( 积) ) ) 刁 = 竽 e t w i 雪i e t w l v l j ， 式中e t = 磐? 鬈办万) r ，矿= 慨。，) ，溉= 焱a g ( 文p ，) ) ，辑= 、 d i a g 鑫( g ( 印1 ) ) ) ，雾= 慨t 。，g j ) r ， 0 2 z 8 8 1 8 蜓 一璺号卫p t d i a g ( g ( e 秽t ) ) d t 粥( e f 反) 一刀t d t 鲳( 6 ( g ( e p t ) ) ) d t a s ( 酋( e 驷，) ) e 一内t a g ( g ( e 孤) ) 岫a g ( b ( g ( e 孤) ) ) 岫a g ( c ( 印t ) ) 司 一等笋 矿u e t v l w 2 e - - e t m 吼司， 其中溉= d ；a g ( 。( e 参t ) ) ，矽一d ；8 9 ( e 磊) ，= 盛蝤( 若( g ( t ) ) ) 。 = - 弼0 2 1r 3 3 2r j m c m c 算法 1 。增加节点： 假设现有的模型包含k 个节点，且节点位置( r l ，。，强) 把区阌( 妃，魄) 分割成k + 1 个子区间为了增加一个节点，我们首先从k + 1 个子区间里 随机的取一个区间，设为i ，则这个子区间记为m 一1 ，心) 然后，从均匀分 布v 哦一l ，哟中随机抽取值努佟为新增节点的位置。受| j 候选节点位置为 哥= ( r l ，镌一1 ，露，r i ，魄) ， g p l s i 模型的b a y e s 估计 1 3 同时，从正态分布n ( ( 曝l ，p t ) t ，( 嚷。，p r ) r ) 中抽取候选参数向量( 舌r ，声t ) t ， 这里 氍l = ( 南。，一1 ，0 ，锄。，k 老) 2 减少节点： 从k 个现有节点 r l ，强 中随机的抽取一个作为删除节点，则 新的候选节点位置为 庐= ( r l ，n 一1 ，r i + 1 ，r k ) ， 同时，从正态分布n ( ( 占0 1 ，z t ) t , e p 乏，芦t p ) 中抽取候选参数向量( 护，声t ) t ， 这里 疆l = ( 南，瞻1 i + 1 ，k 毒) 。 3 移动节点： 首先，从现有节点 r 1 ，r k 中随机抽取一个节点r i ；其次，从截断 正态分布n ( r i ，( r i + 1 一r i 一1 ) 2 盯) i ( n 一1 ，r i + 1 ) 中抽取候选节点噶，其中t 7 2 为 调节参数，则新的候选节点位置为 = ( r l ，。，毪一1 ，r i + l 。，r k ) ， 同时，从正态分布n ? ，矿) t ，r ，z t ) 掌) 中抽取候选参数向量( 善? ，萨) t 因增力嚣一个节点和删除一个节点为可逆的过程，因此薪增一个节点的 接受率为r a i n 1 ，a 2 ，其中 a 2 =翌l i ! 壁：墨! ：! ! ! ：翌：竺q ! 麴! 兰墨! 壶 p ( 6 ，声，k ，r q ，y ，a o ，蜘) b k 南f 石1 d n ( ( z ) t ，( 反t 伊) t ，甲) e x p 学p g ( 妒岱n ) 方+ 彩) 一户6 ( g ( 雪丁( x q ) 占十邪) ) ) 2 司可赢鬲再而瓦面丙瓦讯 印 幽乒l 雪7 g ( 嚣r ( x 位) 6 + z p ) 一户6 ( g ( b t ( x 口) 6 + z p ) ) i g p l s i 模型的b a y e s 估计 ! 生；! f 垒二竺= ! ! ，! ! 三! ! ：! 兰! i 至：! ：乏l 二兰。 如一她 i ( 。，胪) r f、 e x p 卜妒1 巧m ) ) ? - 簖1 阳r ( 0 1 - - ) ， 司i i 丙磊焉硐 1 4 式中瓦穿= ( 晶，5 。p + i - 1 + ，+ 知+ 1 ) ，七) = ( 讶，厣t ) ，口轰七+ 1 ) = ( 鼹l ，z r ) ，d g 表示歪态分布的密度两数 相应地，删除一个节点的接受率为m i n 1 ，a 3 ) ，其中 a s 一 丛墨垒：查= ! ：! ! 竺：望：垫迎! 塾二! ：i ：壶 p ( 磊疹，k ，叫q ，y ，a o ，y o ) d k 妻 觑( ( 只矿一( 磊t ，伊只e 声) ? ) = 书譬箸篙一 坠二塑! 兰鲤：竺兰 !。 南( n + l n 一1 ) ( 2 7 r ) 女i ( d 叠，卢t ) t 广 唧 一撕也) ? ( - 删1t ) 羽，以) ) 式中耀= ( 晶，k t l ，0 ，k 南一1 ) ，口聂船一1 ) = ( 6 乏l ，芦t ) 由于移动一个节点并不改变节点的个数，因此移动一个节点的接受率 为赫n 圭，a 4 ，其中 a | = p 币( 5 两, 3 , 百k , 币乒l a 瓦, y , 瓦a o , 谪y o ) 。銎：i ：塑坠! 堕：堕 讯 d n ( r ：，n ，盯g ) g p l s i 模型的b a y e s 估计 唧f 鼍尹雪t g ( 雪t ( x q ) 否+ z 矽) 一，6 ( g ( 亩r ( x a ) 舌+ z 后) ) i 唧 学i 雪t g ( b 平q ) 刎一j r b ( a ( 矿( ) 郅) ) j ， i z ( 占t ，卢t ) t i i e x p 一( p 1 一百1 ) t - ( 5 1 r ，口t ) 丁( 口1 一痧1 ) i z ( 6 t ，p t ) r | - e x p 一 ( 痧1 一e 1 ) t 历 ，p r ) t ( 痧1 一e 1 ) 式中印= ( 护，声t ) g p l s i 模型的b a y e s 估计 第四章模拟 1 6 我们用= 次截幂样条来拟合未知函数取超参数a 0 = 0 0 1 ，y 0 = ( 0 1 ，o 1 ) r ， a = 2 0 c = 0 4 对于m a r k o v 链的初值选取如下：从标准均匀分布产生( i = 1 ，d ) 并 标准化得到a ，卢= 0 0 ，或= o 0 1 ，( i = 0 ，1 ，1 + k + 2 ) ，k 一2 ，而节点位置r 的 初僮壹接从穗应的先验分布( d i r i c h l e t ( 1 ，王；圭) 分布) 中掬取 对单指数向量q 的估计值矗，用r o t 丘来评价估计的优劣；用8 7 x 取 值范围内的1 0 0 个点的口( ) 函数值与估计值的均方误差 ，1 0 0 ，、2 a s e = 1 - - 0 王- 6 沪夕( u d ) 来判断拟合函数的好坏 本文中的所有程序都用c + + 语言编写，并在c p u 为a m da t h l o n ( t m ) 6 4 p r o c e s s o r3 0 0 0 + 的计算机上运行。 例1 考虑具有典则形式的l o 舭模型其似然函数为 ， 一 、 p ( 纨l 纯，砂) = e x p l y i i n 惫仙( 1 一鼽) ， 则自然参数为吼一i n 啬 本镄子中，已知童然参数p s i n 万( a t x a ) ( b 一脚 + 黟互其中 a = ( 王，1 ，王) 丁娟，卢= 0 5 ，a = 乎一訾，b = 譬+ 器每个协变量 茂= ( 弼1 ，篪2 ，五3 ) t 包禽有三个元索，均独立同均匀分布u ( - 1 ，1 ) ，而当i 为 奇数时，取z i = 0 ；i 为偶数时，取z i = 1 利用随机方法产生容量为8 0 0 的样本，取超参数c r l = 0 2 ，c r 2 = o 。0 1 ，k m a x = 1 0 按照本文的算法迭代7 0 0 0 次，取前面2 0 0 0 次作为预热，取后面的5 0 0 0 个 观察值的平均值作为后验均值的估计这种试验重复进行了2 0 0 次，每次试验 大约费时5 分5 6 秒 表( 4 1 ) 分别列出了l o g i t 模型参数的真值、后验均值、后验中值、标准差 和9 5 最高后验密度区间( h p d ) ，得到的模拟结果如图( 4 1 ) 可见，各个参数 的估计值都饕常接近于真值，丽相应鳍置信区间也包含了冀值。 表4 1 例1 的参数估计结果 a s er c o l u m nn u m b e r 图4 1 例l 的模拟结果( a ) 参数q ，芦估计的箱形图( b ) a s e ，r 的箱形图 为检验本文算法的运行情况，我们从2 0 0 次试验中随机抽出一次试验图 ( 4 2 ) 刻画了参数q ，p 的送行状况图( 4 3 ) 刻画了函数9 ( ) 的拟合状况单指 数估计值& 与真实值的偏差r = 0 9 9 7 6 8 ，函数夕( ) 与真实值之闻的均方误差为 a s e = 0 。0 9 4 5 8 7 7 ，可见拨合程度较高。 g p l s i 模型的b a y e s 估计 图4 2 饲1 中参数搜，声的运行状况 匿4 3 铡王中溺数窖( ) 的拟合图，其中实线为实际函数，星线代表拟合糕数， 点线代表9 5 置信限曲线 例2 考虑具有典则形式的广义p o i s s o n 模型其似然函数为 吲壤，) = 唧 汕a i - a i - l n 坎! 剡自然参数为纯= i n 气。 本例子中，已知自然参数0 = s i n 7 r ( a 丁x a ) ( b a ) + p z ，其中 n = ( 1 ，一1 ，1 ) r 锈，p = 0 3 ，a = 警一糌，b = 乎十糟每个协变量 拖= ( 篾l ，置2 ，整3 ) r 包含有三个元素，且均独立同标准均匀分布u ( 0 ，1 ) ，而 当i 为奇数时，取z i = o ；i 为偶数时，取z i = ，l 。 利用随机方法产生容量为8 0 0 的样本，取超参数o 1 = 0 0 8 ，盯2 0 1 ，k m a x 一 5 按照本文的算法迭代6 0 0 0 次，取前面1 0 0 0 次作为预热，取后面的5 0 0 0 伞 观察值盼平均值作为后验均值的估计这种试验重复进行了2 0 0 次，每次试验 大约费时2 分3 7 秒 g p l s i 模型的b a y e s 估计 表( 4 2 ) 分别列出了广义p o i s s o n 模型参数的真值、后验均僖、后验中值、 标准差和9 5 最高后验密度区间( h p d ) ，得到的模拟结果如图( 4 4 ) 可见，各 个参数的估计值都接近于真德，丽相应的置借区间也包含了真值 为检验本文算法的运行情况，我们从2 0 0 次试验中随枧抽出一次试验。图 ( 4 5 ) 刻画了参数q ，p 的运行状况图( 4 6 ) 刻画了函数9 ( ) 的拟合状况单指 数估计值a 与真实值的偏差r = 0 9 9 8 7 3 9 ，函数9 ( ) 与真实值之间的均方误差 为a s e = 0 。0 1 3 1 6 4 3 ，可觅拟合程度较高 表4 2 例2 的参数估计结果 图4 4 倒2 的模拟结果( a ) 参数a ，声估计的箱形圈( b ) a s e ，r 的箱形图 g p l s i 模型的b a y e s 估计 o - 9 躺 o 3 o 咱3 _ o 6 峨 图4 5 例2 中参数理，的运行状况 g p l s i 模型的b a y e s 估计 图4 6 例2 中函数9 ( ) 的拟合图，其中实线为实际函数，星线代表拟合函数， 点线代表9 5 置信限蓝线。 g p l s i 模型的b a y e s 估计 第五章实例分析 我们用本文的算法来分析瑞士钞票数据，此组数据共2 0 0 个观察值，有6 个骞变量，溉。，拖分别表示钞票的长度，左边高度，右边高度，内框到下 缘的距离，内框到上缘的距离，对角线长度；因变量y 取值0 或1 ，分别表示 真钞和假钞这组数据中，前1 0 0 个为真钞的数据，后1 0 0 个为假钞的数据 对予这缀数据，e f a n g 和x i a 已经用单指数模型进行了变量选择【黝，可 知选择变量x 4 ，甄，x 6 作为单指数模型的自变量可以较好的拟合模型。因此， 我们选择五，托，蕊作为自变量来建立l o g i t 模型 f、 p ( 轨i 吼) = e x p y i o i l i l ( 1 + ) i = 1 ，2 0 0 ， k， 并取吼为其典则形式，即吼= g ( a l 咒4 + a 2 蕊5 + 俚3 五6 ) 我们用二次截幂样条来拟合未知函数g ( ) ，取超参数a o 一0 0 1 ，y o = ( 0 1 ，o 1 ) r ， 天一2 0 ，c = 0 。4 ，o 1 = 。0 8 ，0 2 = 0 1 ，k m a x = 4 。 对于m a r k o v 链的初值选取如下；从标准均匀分布产生( i 一1 ，3 ) 并 标准化得到o l ，屯= 0 1 ( i 一0 ，1 ，1 + k + 2 ) ，k 一1 ，而节点r 的初值直接从相 应的先验分布( d i r i c h l e t ( 1 ，。，王；薹) 分布) 中抽取。 按照本文的算法迭代2 0 0 0 0 次，取前面1 0 0 0 0 次作为预热，取后面的1 0 0 0 0 个观察值的平均值作为后验均值的估计 图( 5 。7 ) 刻画了单指数8 的抽样状况表( 5 3 ) 舞高了参数的的后验均值、 后验中值、标准差和9 5 最高后验密度区阍图( 5 8 ) 为丞数窖( ) 的估计曲线。 表5 3 实例中参数o 的估计结采 g p l s i 模型的b a y e s 估计 x 1 矿 x 1 0 4 雷5 7 实铡中参数貔的信计睦线图 x 矿 图5 8 实例中函数9 ( ) 的估计曲线 g p l s i 模型的b a y e s 估计 第六章总结 作为广义线性模型和非参数回归模型的深入推广，广义部分线性单指数 模銎在现代统计中有重要的作用本文利用自由节点的b a y e s 样条技术对广义 部分线性单指数模型进行统计分析。本文震截幂样条函数来逼近未知函数， 也可以用有更好数值特征的b 样条函数本文把囱变量x 的所有元素都认 为是有用的，但舀变量中可能包含不相关的元素，因此，还可以用自由节点的 b a y e s 榉条技术对广义部分线性单指数模型盼变量选择做连分析。 参考文献 f l lc a r r o l l ，r j 。，f a n ，j ，g i j b e l s ，i ，w a n d ，m p g e n e r a l i z e dp a r t i a l l yl i n e a rs i n g l e - i n d e xm o d e l s 【j 】j o u r n a lo ft h ea m e r i c a ns t a t i s t i c a la s s o c i a t i o n ，1 9 9 7 ，9 2 ( 4 3 8 ) ：4 7 0 4 8 9 【2 】m c c u l l a g h ，p ，n e l d e r ，j a g e n e r a l i z e dl i n e a rm o d e l s ( 2 n de d 。) f m 】l o n d o n ：c h a p - m a na n dh a l l ，1 9 8 9 。 f 3 】w e i s b e r g ，s ，w e l s h ，a h e s t i m a t i n gt h em i s s i n gl i n kf u n c t i o n j t h ea n n a l so f s t a t i s t i c s ，1 9 9 4 ，2 2 ：1 6 7 4 - 1 7 0 0 c 4 】l i ，k 。c s l i c e di n v e r s er e g r e s s i o nf o rd i m e n s i o nr e d u c t i o n j j o u r n a lo ft h ea m e r - i c a ns t a t i s t i c a la s s o c i a t i o n ，1 9 9 1 ，8 6 ：3 1 6 3 4 2 。 f 5 1s t o k e r ，t m c o n s k 8 t e n te s t i m a t i o no fs c a l e dc o e f f i c i e n t s j e c o n o m e t r i c a ，1 9 8 6 ， 5 4 ( 6 ) ：1 4 6 1 1 4 8 1 【6 】h s r d l e ，w ，s t o k e r ，t m i n v e s t i n gs m o o t hm u l t i p l er e g r e s s i o nb yt h em e t h o do f a v e r a g ed e r i v a t i v e s j j o u r n a lo ft h ea m e r i c a ns t a t i s t i c a la s s o c i a t i o n ，1 9 8 9 ，8 4 ： 9 8 6 9 9 5 【7 】h o r o w i t z ，j l ，h i i r d l e ，w d i r e c ts e m i p a r a m e t r i ce s t i m a t i o no fs i n g l e - i n d e xm o d e l s w i t hd i s c r e t ec o v a r i a t e s j 。j o u r n a lo ft h ea m e r i c a ns t a t i s t i c a la s s o c i a t i o n ，1 9 9 6 ，9 1 ： 1 6 3 2 - 1 6 4 0 i s ii c h i m u r a ，h s e m i p a r a m e t r i cl e a s ts q u a r e s ( s l s ) a n dw e i g h t e ds l se s t i m a t i o no f s i n g l e - i n d e xm o d e l s j ，j o u r n a lo fe c o n o m e t r i c s ，1 9 9 3 ，5 8 ：7 1 1 2 0 。 同h h r d l e ，w ，h a l l ，p ，i c h i m u r a ，h o p t i m a ls m o o t h i n gi ns i n g l e - i n d e xm o d e l s j t h ea n n a l so fs t a t i s t i c s ，1 9 9 3 ，2 1 ( 1 ) ：1 5 7 - 1 7 8 1 0 】y u ，y ，r u p p e r t ，d p e n a l i z e ds p l i n ee s t i m a t i o nf o rp a r t i a l l yl i n e a rs i n g l e - i n d e x m o d e l s j j o