《基于多水平模型的工具变量方法研究及应用》由会员上传分享,免费在线阅读,更多相关内容在学术论文-天天文库。
分类号论文编号._________________密级________________^^fi^?博士学位论文基于多水平模型的工具变量方法研究及应用TheStudyandApplicationofMultilevelInstrumentalVariableModelinDatawithHierarchicalStructure研究生姓名:秦婴逸学号:20122148指导教师:贺佳教授卫生统计学教研室吴骋副教授卫生统计学教研室学科、专业:流行病与卫生统计学学位类型:科学学位答辩日期:2015年5月二〇一五年五月 第工苹*太摩博士学位论文基于多水平模型的工具变量方法研究及应用TheStudyandApplicationofMultilevelInstrumentalVariableModelinDatawithHierarchicalStructure研究生姓名:秦婴逸指导教师:贺佳教授吴骋副教授学科、专业:流行病与卫生统计学培养单位:第二军医大学卫生统计教研室二〇一五年五月 独创性声明本人声明所呈交的学位论文是我个人在导师指导下进行的研究工作。除了文中特别加以标注和致谢的地方外,论文中不包含其他人已经发表或撰写过的研究成果。与我一同工作的同志对本研究所做的任何贡献均已在论文中作了明确的说明并表示谢意。本人承担本声明的法律责任。学位论文作者签名日期:年r月曰学位论文版权使用授权声明本人完全了解第二军医大学有关保留、使用学位论文的规定,第二军医大学有权保留并向国家有关部门或机构送交论文的复印件和电子版,允许论文被查阅和借阅。本人授权第二军医大学可以将学位论文全文或部分内容编入《中国学位论文全文数据库》、《中国优秀博硕士学位论文全文数据库》等数据库进行检索。可以采用影印、缩印或扫描等复制手段保存、汇编学位论文。(保密的学位论文在解密后适用本授权书)导师签名:日期:如&年r月日日期:年r月日 目录摘要.........................................................................................................-1-Abstract.......................................................................................................-5-缩略词表...................................................................................................-10-第一部分概述.........................................................................................-12-一、研究背景.....................................................................................-12-二、研究目的与意义.........................................................................-18-三、研究内容.....................................................................................-19-四、研究步骤.....................................................................................-25-五、资料来源、分析工具及研究平台.............................................-26-第二部分多水平工具变量模型在模拟数据中的拟合和验证..............-28-一、处理/暴露因素和结局变量为连续型变量................................-29-二、处理/暴露因素和结局变量为连续型变量(均存在层级效应)-38-三、处理/暴露因素和结局变量为分类变量....................................-46-四、讨论.............................................................................................-54-第三部分实例应用.................................................................................-58-一、概况.............................................................................................-58-二、应用实例一:60岁以上老人每周运动时间和身体健康状况的关系 .............................................................................................................-59-三、应用实例二:60岁以上男性吸烟对其患高血压的影响.........-63-四、讨论.............................................................................................-68-第四部分研究结论与展望.....................................................................-70-一、研究结论.................................................................................-70-二、研究特色和创新点.................................................................-71-三、尚待开展的研究.....................................................................-71-附录:核心程序.......................................................................................-73-文献综述...................................................................................................-87-参考文献.................................................................................................-100- 基于多水平模型的工具变量方法研究及应用摘要研究背景:随着医疗卫生信息化的不断发展,对分析方法的需求不断增加,并且“真实世界的研究”在目前越来越受到关注,随着数据集收集范围的不断扩大,数据来源常常包括不同的地区、不同的医院,如全市医院信息数据、全国卫生服务调查数据等。这些数据具有层次结构特征,对于这样的数据进行分析,首先需要考虑数据中不同水平单位对结果可能产生的影响,对于此类问题,多水平分析模型可以很好地进行处理。多水平模型将方差成分模型和多元回归模型相结合,把广义线性模型中的差异拆分为固定效应和随机效应两部分,从而更加准确地估计处理/暴露因素的效应值。在利用多水平模型对具有层次结构特征的数据进行分析时,不仅能够很好地控制不同的水平因素对结果所产生的影响,而且通过纳入多个已知观测的混杂因素,能较好地控制这些已知观测混杂因素对结果产生的影响。也有研究者将倾向性评分法(PropensityScoreAnalysis,PSA)引入多水平模型,采取倾向性评分匹配法、分层法和加权法对数据集中的已知观测混杂因素进行控制,从而更好地获得准确的结果。但是,目前大部分卫生服务方面的调查是关于人群健康方面的普查,当研究者利用这样的数据进行某专项疾病或健康方面的研究时,通常所需的变量并不能完全满足研究要求,研究结果通常会受到未知观测混杂因素的影响,如分析每周运动时间对自身健康状况的影响,数据中已包括了一些已知观测混杂因素(年龄、BMI、患病情况、吸烟、喝酒等),但对于本人的心情、家庭关系、病情轻重程度等因素调查数据中可能未包含或难以测量,这些因素同样可能会影响到分析结果的准确性,而目前这些常用于具有层次结构特征数据的分析方法并不能控制这方面的影响。在普通数据分析时,可以利用工具变量方法对未知观测混杂因素进行控制,在本研究中,我们将工具变量的思想引入多水平模型数据分析中,用以处理未知观测混杂因素所产生的影响。研究目的:目前,对于层次结构特征数据中未观测混杂偏倚的控制鲜有研究进行探索,本研究针对此问题,将构建出多水平工具变量模型(MultilevelInstrumentalVariable,MIV),从而较为全面地控制层次结构特征数据中水平因素、已知观测混杂因素和未知观测混杂因素对结果所产生的偏倚,并且基于资料中数据类型的不同(连续性变量和分类变量),将构建出不同的多水平工具变量模型,以分别适用于连续型变量和分类变量资料分析中。本研究还将对所构建的模型进行准确性和精确性方面的评价,探索各种数据条件下所应当采用的最佳参数估计模型,在模型构建的基础上,笔者引入自助法-1- 基于多水平模型的工具变量方法研究及应用(Bootstrap),使模型估计得到的结果更加可靠。研究方法:研究首先进行数据模拟,数据的模拟过程主要根据数据类型的不同分为3部分,在数据模拟过程中将考虑到不同强度的未知观测混杂因素和不同强度的工具变量,从而较为全面对模型进行评价。(1)模型构建在处理/暴露因素和结局变量为连续型变量情况下,将构建出两阶段最小二乘多水平工具变量模型与两阶段残差纳入多水平工具变量模型;在处理/暴露因素和结局变量为连续型变量且均存在层次效应情况下,将构建出两阶段多水平回归工具变量模型和两阶段多水平回归残差纳入工具变量模型;在处理/暴露因素和结局变量为分类变量情况下,将构建出两阶段logistic回归多水平工具变量模型和线性回归+logistic回归多水平工具变量模型。在模型构建的过程中我们还引入了自助法(Bootstrap),在本研究中自助法采用的是分层个例重复抽样法,根据原始样本量的大小进行等样本重复抽样,每次抽500次,然后用所构建的模型对500个复样本进行分析。(2)模型评价模型评价部分将所构建的多水平工具变量模型和普通多水平回归模型所得的结果用四个指标进行客观科学的评价,分别为绝对偏倚、置信区间宽度、标准误、置信区间覆盖率。根据这四个指标可以反映模型在不同数据情况下的准确性和精确性,为后续模型的调整和应用提供了科学根据。(3)实例分析最后将构建的多水平工具变量模型应用于实例分析中。实例分析数据来源于第五次全国卫生服务调查数据(上海)。针对结局变量和处理/暴露因素为连续型变量,本研究所选的实例为分析上海60岁以上老人每周体育锻炼时间对其健康状况的影响,男性和女性分别进行分析,其中可能存在的未知观测混杂包括本人的心情、家庭关系、病情轻重程度等,结局变量采用欧洲五维健康量表(Europeanqualityoflife5-dimensions,EQ-5D)评分,工具变量选择为其爱人每周运动的次数。针对结局变量和处理/暴露因素为分类变量,本研究所选的实例为分析上海市60以上岁男性是否吸烟对其是否患有高血压的影响,其中可能存在的未知混杂包括基因特征、周围环境因素等,工具变量选择为其家人是否吸烟。实例分析中应用普通多水平模型和模拟中所获得的最优多水平工具变量模型进行分析,并对不同方法所获得的结果进行比较。研究结果:数据模拟的结果显示,研究发现在资料中存在未知观测混杂因素时,所构建的多水平工具变量模型有较好的效果,具体如下:-2- 基于多水平模型的工具变量方法研究及应用(1)处理/暴露因素和结局变量为连续型变量当不存在未知观测混杂因素时,所有模型均能获得较为理想的结果,但当研究中存在未知观测混杂因素时,普通多水平线性回归模型和自助法多水平线性回归模型会获得偏倚较大的结果,偏差最大的出现在β୳=6、ߙ=1时的普通多水平线性回归模型中,其绝对误差为-2.8219,但多水平工具变量的结果较为稳定,当β୳=6、ߙ=5时两阶段最小二乘多水平工具变量模型、两阶段残差纳入多水平工具变量模型、自助法两阶段最小二乘多水平工具变量模型和自助法两阶段残差纳入多水平工具变量模型结果的绝对偏倚分别为-0.0004、-0.0009、0.0012和0.0006。在四种多水平工具变量模型中,自助法引入的模型结果的置信区间较宽,提示结果更为保守,当工具变量的强度增加时,其区间会相应的变窄。(2)处理/暴露因素和结局变量为连续型变量且均存在层次效应普通多水平线性回归模型仅适用于无未知观测混杂因素的数据中,该模型在数据中存在未知观测混杂因素时所得的结果偏离金标准较大。虽然两阶段自助法两阶段最小二乘多水平工具变量模型可以基本准确估计得出处理/暴露因素的效应值,但其置信区间过宽。两阶段多水平回归工具变量模型、两阶段多水平回归残差纳入工具变量模型、自助法两阶段多水平回归工具变量模型和自助法两阶段多水平回归残差纳入工具变量模型在不同的数据情况下均能得到理想的结果,其中自助法两阶段多水平回归工具变量模型的准确度和精确度总体最佳,当β୳=6、ߙ=5时该模型的绝对偏倚仅为0.0009。(3)处理/暴露因素和结局变量为分类变量结果展示当数据中无未知观测混杂因素存在的情况下普通多水平logistic回归模型所获得结果最佳,但当混杂因素存在时,普通多水平logistic回归模型所得的结果将偏离金标准较远,并且置信区间覆盖率较低,自助法两阶段logistic回归多水平工具变量模型和自助法线性回归+logistic回归多水平工具变量模型两种模型在有未知观测混杂因素数据情况下表现较好,点估计最接近所设的金标准,但此两个模型的置信区间受到工具变量强度影响较大,在弱工具变量时,模型估计的结果过于保守,区间过宽,两阶段logistic回归多水平工具变量模型和线性回归+logistic回归多水平工具变量模型两种模型在各种数据情况下均未表现出很好的效果。实例分析的结果显示,在分析每周运动时间同自身健康状况间的关系时,普通多水平回归模型和多水平工具变量模型均提示在60岁以上人群中,男性和女性每周运动时间同健康评分间存在正相关关系,但在男性中普通多水平回归模型所得回归系数为0.42(0.41-0.43),多水平工具变量模型所得回归系数为0.70(0.53-0.86),两者相差约0.3,女性中,普通多水平回归模型所得回归系数为0.49(0.48-0.50),多-3- 基于多水平模型的工具变量方法研究及应用水平模型所得结果为0.37(0.21-0.53),两者相差约0.1。说明在该实例中,两者都能较好地识别出感兴趣的影响因素与应变量之间的关系,但对关系大小的衡量存在差别。在分析60岁以上男性吸烟和患高血压间关系时,普通多水平logistic回归提示吸烟为保护因素,OR值为0.74(0.65-0.83),此与目前所公认的结论相违背,但多水平工具变量模型提示吸烟是患高血压的危险因素,OR值为5.05(1.40-18.26)。此项研究中仅纳入五项协变量,很多高血压的危险因素在卫生服务调查中未收集,如家族史、血液生化指标等,普通分析方法无法控制这些未知观测混杂因素对结果产生的影响,从而得到了错误的结论,当利用多水平工具变量模型对这些因素进行控制后,所得结果将更为可靠。研究结论:本研究通过模拟研究和实例分析发现多水平工具变量模型均能很好地获得较为准确的结果。当研究资料收集较全,均不存在十分重要的未知观测混杂因素时,普通的分析模型即可获得较好的结果。当资料并非为专项研究调查,层次结构特征数据中遗漏了较多或一些较为重要的影响因素时,普通分析模型将不再适用,可以采用本研究所构建的多水平工具变量模型。当针对连续性变量数据时,首先需要看数据中处理/暴露因素在水平2单位上是否存在异质性,当处理/暴露因素存在层级效应时,建议采用自助法两阶段多水平回归工具变量模型;当针对分类数据时,建议采用自助法两阶段logistic回归多水平工具变量模型和自助法线性回归+logistic回归多水平工具变量模型两种模型。分析过程中可以采用多个工具变量模型进行分析,当结果一致时,可以更加肯定研究的结论。在模型使用过程中,应当尽可能地寻找强度较高的工具变量,从而可以获得更为准确的结果。关键词:层次结构特征数据,多水平模型,工具变量,未知观测混杂因素。-4- 基于多水平模型的工具变量方法研究及应用AbstractBackground:Withthedevelopmentofmedicalandhealthinformation,weneedmoreeffectiveanalysismethods.Furthermore,the“realworldstudy”hadbeenpaidmoreattention.Recently,manyresearcheshavecollecteddatafrommanysources,suchasdifferenthospitals,differentdistricts,sothedatahaveahierarchicalstructure.Toanalysisthesedata,weshouldtakeconsideroftheinfluenceofmultilevelfactorsonresults.Themultilevelregressionmodelcouldcontroltheseinfluence,andthismodelcombinestheprinciplesofvariancecomponentmodelandmultipleregressionmodel.Themultilevelregressionmodeldividesthedifferenceintofixedeffectsandrandomeffectstoassessthemoreaccurateeffectoftreatment/exposurefactors.Thismodelnotonlycancontroltheinfluenceofmultilevelfactors,butalsocanincludemanymeasuredconfoundingfactorstocontroltheconfoundingbiasfromthesefactors.Recently,someresearchershadintroducedthePropensityScoreAnalysis(PSA)intothemultilevelregressionmodel.TheycouldcontrolthemeasuredconfoundingbiasmuchbetterusingPSmatching,PSstratification,andPSweighting.Themainpurposesofthemajorityofthehealthservicedataaretoevaluatethehealthstatusamongthewholepopulation.Whenweusethesedatatoconducttheresearchesonspecificdiseasesorhealth,thevariablesinthedatacouldnotfullymeettherequirementsofourresearchesusually.Therefore,theresultswouldbeaffectedbysomeunmeasuredconfoundingbias.Forexample,whenwewanttofindtheassossationbetweenthemovementtimeofweekandphysicalcondition,therearesomeknownconfoundingfactors(suchasage,BMI,disease,smoke,drinkandetal.)whichcouldbecontroled.However,someunmeasuredconfoundingfactors(suchasmood,familyrelationships,andetal.)alsoexistinthedata,andtheseconfoundingfactorsmayaffecttheaccuracyofanalysisresults.Thereislittleofresearcheswhichcontroltheunmeasuredconfoundingbiasinthedatawithhierarchicalstructure.Weusuallyusetheinstrumentalvariable(IV)methodtocontrolthiskindofbiasinthedatawithouthierarchicalstructure.Inourresearch,wewillintroducetheprincipleofIVintothemultilevelregressionmodeltobuildthemultilevelinstrumentalvariablemodel(MIV).Aim:-5- 基于多水平模型的工具变量方法研究及应用Recently,therewaslittleresearchwhichfocusonthecontrolofunmeasuredconfoundingfactorsinthedatawithhierarchicalstructure.Tosolvethisproblem,wewillestablishtheMultilevelInstrumentalVariableModel(MIV)tocontroltheconfoundingfactorsinthedatawithhierarchicalstructure,suchaslevelfactors,knownconfoundingfactorsandunknownunmeasuredfactors.Furthermore,wewillbuilddifferentkindsofMIVsaccordingtothedifferenttypesofdata(continuousandbinaryvariables).WewillmakeobjectiveandscientificevaluationforthemultilevelregressionmodelandMIVmodelsthroughapplyingthemodelsintothesimulateddata.Andwealsointroducethebootstrapintothemodeltomaketheresultsmorereliable.Methods:First,weconductasimulationstudy.Asimulationstudyhadbeendividedintothreepartsaccordingtothedifferenttypesofdata.Andwewillconsiderthedifferentstrengthofunmeasuredconfounderandinstumentalvariableintheprocessofsimulationtoprovideacomprehensiveevaluationforthemodels.(1)modelbuildingWhentreatment/exposurefactorsandoutcomearecontinuousvariables,webuildMultilevelTwoStageLeastSquare(M2SLS)andMultilevelTwoStageResidualInclusion(M2SRI).Whentreatment/exposurefactorsandoutcomearecontinuousvariablesandbothhaveeffectofhierarchy,webuildTwoStageLinearMultilevelModel(2SLM)andTwoStageLinearMultilevelResidualInclusionModel(2SLMRI).Whentreatment/exposurefactorsandoutcomearebinaryvariables,webuildMultilevelTwoStagePredictionSubstitute(M2SPS)andMultilevelLeastSquarePredictionSubstitute(MLSPS).Wewillintroducebootstrapintothemodelsandusehierarchicalcase-resampling.Thesamplesizeofresampleisequaltothatoftheoriginalsample.WeapplytheMIVmodelstothe500resamples.(2)modelevaluationWeprovidetheobjectiveandscientificevaluationforthemultilevelregressionmodelandMIVmodelsbasedonfourindexes:bias,confidenceintervals,standarderror,andconfidenceintervalscoverage.Thesefourindexescouldreflecttheaccuracyandprecisionofthemodelsindifferentdatasituation.(3)instanceanalysisFinally,weusetheseMIVmodelsintherealworlddata,whichfromthefifthnational-6- 基于多水平模型的工具变量方法研究及应用healthservicessurvey(Shanghai).Fortreatment/exposurefactorsandoutcomeascontinuousvariables,wewillexploretherelationshipbetweentheweeklyphysicalexercisetimeandhealthstatusamongthepeopleover60yearsold(menandwomen).Thepotentialunmeasuredconfoundersmightbemood,familyrelationship,severityofdisease,andetc.TheoutcomeisEuropeanqualityoflife5-dimensions(EQ-5D).Theinstrumentalvariableishusbandorwife’sfrequenceofweeklyphysicalexercise.Fortreatment/exposurefactorsandoutcomeasbinaryvariables,wewillexploretheassociationbetweenthesmokeandhypertensionamongmenover60yearsold.Thepotentialunmeasuredconfoundersmightbethecharacteristicsofgene,environmentandetc.Theinstrumentalvariableisthesmokestatusoftheirfamilymembers.WewillapplyboththemultilevelregressionmodelandMIVmodelsintheseanalyses,andcomparetheresultsofthesemodels.Results:Inthesimulationstudy,wefoundthattheMIVmodelshadobtainedbetterresultswhenthereareunmeasuredconfoundersinthedata.(1)treatment/exposurefactorsandoutcomearecontinuousvariables.Whentherewasnounmeasuredconfounderinthedata,allthesesixmodelshadobtainedgoodresults.Howerver,whenunmeasuredconfounderswereaddedinthedata,thepointestimationsofmultilevelregressionmodelandbootstrapmultilevelregressionmodelwerefarawayfromthegoldstandard,andthemaxbiaswas-2.8219(multilevelregressionmodel,β୳=6,ߙ=1).TheresultsoffourMIVmodelsweremuchmorestable,andthebiasforM2SLS,M2SRI,boot-M2SLSandboot-M2SRIwere-0.0004、-0.0009、0.0012and0.0006(β୳=6,ߙ=1).TheconfidenceintervalsofBoot-M2SLSistoowide.Withtheincreaseofthestrengthoftheinstrumentalvariable,theconfidenceintervalswouldbenarrow.(2)treatment/exposurefactorsandoutcomearecontinuousvariablesandbothhaveeffectofhierarchy.Themultilevelregressionmodelisonlyapplicabletodatawithoutunmeasuredconfounder.However,thismodelwouldobtainbiasresultswhenweaddunmeasuredconfoundersinthedata.AlthoughBoot-M2SLScouldbasicallyestimatetheeffectoftreatment/exposurefactor,thewidthofCIsinthismodelwastoowide.Wefoundthat2SLM,2SLMRI,Boot-2SLM,andBoot-2SLMRIcouldobtainaccurateresultsunderany-7- 基于多水平模型的工具变量方法研究及应用datasituation.TheaccuracyandprecisionofBoot-2SLMwerebestamongthesefourmodels,andthebiasforthismodelwas0.0009(β୳=6,ߙ=5).(3)treatment/exposurefactorsandoutcomearebinaryvariables.Whentherewasnounmeasuredconfounderinthedata,theresultsofthemultilevellogisticregressionmodelweremuchbetter.Thepointestimationsofthismodelwerefarawayfromthegoldstandardandtheconfidenceintervalscoveragewastoolowwhentherewereunmeasuredconfoundersinthedata.Boot-M2SPSandBoot-MLSPScouldobtainbetterresultsunderthisdatasituation.However,theCIsofthesetwomodelswereaffectedbythestrengthoftheinstrumentalvariable.TheCIsweretoowideusingtheweakinstrumentalvariable.AndwealsofoundtheM2SPSandMLSPSdidnotshowgoodresultsunderanydatasituation.Toexploretherelationshipbetweentheweeklyphysicalexercisetimeandhealthstatus,allofthesesixmodelshadindicatedtherewerepositiverelationshipsamongmenandwomenpeopleover60yearsold.Amongthemen,theregressioncoefficientofthemultilevelregressionmodelwas0.42(0.41-0.43),andtheregressioncoefficientsofallofMIVmodelswereabout0.70(0.53-0.86).Thedifferenceofthesetwokindsofmodelswas0.3.Amongthewomen,theregressioncoefficientofthemultilevelregressionmodelwas0.49(0.48-0.50),andtheregressioncoefficientsofMIVmodelswereabout0.37(0.21-0.53).Thedifferenceofthesetwokindsofmodelswas0.1.AlthoughtbothtraditionalregressionmodelandMIVcouldfindthesignificantrelationshipbetweenthesetwothings,thesewasdifferenceonthethesizeofeffect.Toexploretherelationshipbetweenthesmokeandhypertensionamongmenover60yearsold,themultilevellogisticregressionmodelhadindicatedsmokewasaprotectivefactor(OR=0.74,0.65-0.83).However,theMIVmodelshadshownthatsmokewasariskfactor(OR=5.05,1.40-18.26).WeconsideredthattheresultsofMIVmodelswereconsistentwiththegenerallyacceptedconclusion.Weincludedonly5covariantsinthisstudy,andmanyfactorsforhypertensionwerenotcollectedinthenationalhealthservicessurvey,suchasfamilyhistory,bloodbiochemicalindex,andetal.Thetraditionalregressionmodelcouldnotcontroltheeffectoftheseunmeasuredconfoundingfactors.However,theMIVcouldcontrolthesefactorsandprovidemorereliableresults.Conclusions:TheresultsofsimulationstudyhadindicatedthattheMIVcouldobtainmoreaccurate-8- 基于多水平模型的工具变量方法研究及应用results.Thetraditionalregressionmodelcouldbeappliedinthestudywhichcollectmostofimpotantfactors.However,whenthereweresomeimpotantunmeasuredconfoundingfactorsinthedatawithhierarchicalstructure,weshouldusetheMIVtofindthetruerelationship.Whentreatment/exposurefactorsandoutcomearecontinuousvariablesandbothhaveeffectofhierarchy,wesuggestedtouseBootstrapTwoStageLinearMultilevelModel(boot-2SLM).Whentreatment/exposurefactorsandoutcomearebinaryvariables,wesuggestedtouseBootstrapMultilevelTwoStagePredictionSubstitute(boot-M2SPS)andBootstrapMultilevelLeastSquarePredictionSubstitute(boot-MLSPS).Furthermore,wesuggestedthatweshouldapplyseveralMIVsinthestudy.Andwhentheresultswereconsistent,theresultwillbemorereliable.Moreover,thestronginstrumentalvariableswouldprovidemuchmoreaccurateresults.KEYWORDS:datawithhierarchicalstructure;multilevelmodel;instrumentalvariables;unmeasuredconfounder.-9- 基于多水平模型的工具变量方法研究及应用缩略词表缩略词英文全称中文名称RCTRandomizedcontrolledtrial随机对照试验PSAPropensityScoreAnalysis倾向性评分法ANNArtificialNeuralNetwork人工神经网络SVMsSupportVectorMachines支持向量机CARTCatigoricalandRegressionTree分类与回归树MSMMarginalStructuralModels边缘结构模型IPWInverseProbabilityWeighting逆概率加权法IVInstrumentalVariable工具变量2SLS2-StageLeastSquares两阶段最小二乘法RDRiskDifference风险差异OROddsRatio优势比RRRelativeRisk相对危险度2SPSTwo-StagePredictorSubstitution两阶段logistic回归模型MIVMultilevelInstrumentalVariable多水平工具变量模型M2SLSMultilevelTwoStageLeastSquare两阶段最小二乘多水平工具变量模型M2SRIMultilevelTwoStageResidualInclusion两阶段残差纳入多水平工具变量模型2SLMTwoStageLinearMultilevelModel两阶段多水平回归工具变量模型2SLMRITwoStageLinearMultilevelResidual两阶段多水平回归残差纳InclusionModel入工具变量模型M2SPSMultilevelTwoStagePredictionSubstitute两阶段logistic回归多水平工具变量模型MLSPSMultilevelLeastSquarePredictionSubstitute线性回归+logistic回归多水平工具变量模型CIConfidenceIntervals置信区间SEStandardError标准误SDStandardDifference标准差-10- 基于多水平模型的工具变量方法研究及应用EQ-5DEuropeanqualityoflife5-dimensions欧洲五维健康量表2SPStwo-stagepredictorsubstitution线性回归+logistic回归模型MLMMultilevelLinearModel普通多水平线性回归模型MLRMMultilevelLogisticRegressionModel普通多水平logistic回归模型BMIBodyMassIndex身体质量指数-11- 基于多水平模型的工具变量方法研究及应用第一部分概述一、研究背景随着医疗卫生事业的不断进步和信息化产业的发展,以及广大群众对健康的需求1-3度不断提高,目前医疗卫生服务行业已逐渐从传统的一致性治疗方案模式向个体4,5化治疗方案模式发展,从最初简单的就医问药向三级预防的最高层病因预防发展。在这样的医疗卫生大背景下,循证医学思想已经深入人心,2001牛津证据分级将循6-8证医学证据等级从高到底分为五类等级,分别为:1a同质RCT(randomizedcontrolledtrial)的系统评价,1b单个RCT(可信区间窄),1c全或无病案系列;2a同质队列研究的系统评价,2b单个队列研究(包括低质量RCT,如随访率<80%),2c结果研究,生态学研究;3a同质病例对照研究的系统评价,3b单个病例对照;4病例系列研究(包括低质量队列和病例对照研究);5基于经验未经严格论证的专家意见。处于证据等级顶端的为随机对照试验,虽然此种研究方法可以通过对患者或受试对象进行严格选取和干预措施,从而可以很好地控制研究中混杂因素对结果的影响,但该类试验方法本身存在一定的局限性,由于受试对象和干预措施进过严格的限定,导致所得试验结果的外延性受到一定影响,当干预措施的创伤较大或是针对某些罕见病种的研究时,医学伦理和病例数较少都会限制随机对照试验的应用。由于随机对照试验的试验周期一般不会过长,所以对一些药物不良反应的研究,随机对照试验也是9-11不适用的。针对这样的问题,目前对“真实世界研究”(realworldstudy)的关注越来越深,并且,医疗卫生信息化的高速发展将医疗卫生行业研究推向“大数据”(bigdata)时代,对于“真实世界研究”和医疗卫生“大数据”的研究,随机对照试验研12,13究模式并不十分适合,因而观察性研究的方法越来越被关注。观察性研究是在真实世界情况下所实施的研究,并未对观察对象实施人为干预,其与随机对照试验是相互补充的,观察性研究的结果既能进一步验证随机对照试验所得到的结果,又能反映随机对照试验应用于真实世界中的效果,并且观察性研究也能10,14-17在一定程度上情况下弥补随机对照试验的局限性:(1)观察性研究的人群来源于真实世界情况,并未对其施加刻意的人为干预措施,观察人群所处的环境同现实中患者所处的环境类似,其所获得结果的外延性较好;(2)当研究目的为发现某些疾病的危险因素,以完成病因预防的目的时,随机对照试验设计并不能很好实现,观察性研究可以对这样的问题进行很好的分析;(3)由于随机对照试验的周期有一定的局限性,有些药物或治疗措施的不良反应很难在短时间内被发现,此时,我们可以-12- 基于多水平模型的工具变量方法研究及应用通过观察性研究发现一些药物或治疗措施的不良反应,从而对患者的健康进行很好的保护。虽然观察性研究目前越来越被广泛的重视,但是由于观察性研究仅为观察真实世界人群中的病情或疾病发展情况,并未对其刻意的施加人为干预措施,也没有刻意规范被观测人群在就医或生活过程中的一些行为,并且一些病情或基本情况并非是通过实验室检查或严格医疗机构进行的检测所获得,而是通过调查问卷或面对面访谈所获得,故其研究结果通常会受到一些混杂因素的干扰,尤其是对医疗卫生“大数据”而言,由于其海量和多样性特征决定了其中存在大量的背景噪声,因而影响了所获得结果的可靠性,这些混杂因素中一部分可以通过规范观察性研究的设计,在患者的选取和信息的收集过程中减少混杂因素的产生,另外一部分就需要通过统计学的算法对一些混杂因素进行校正和去除,从而获得更加准确的结果。(一)混杂偏倚控制方法的研究情况在研究分析中常见的混杂因素包括:已知观测混杂因素、未知观测混杂因素、时依性混杂因素、水平间差异造成的混杂等。随着观察性研究在实际应用中越来越多,18混杂偏倚控制方法的研究也得到了很大的发展,基于最初始的多元线性回归模型和logistic回归模型,目前常用的混杂偏倚控制方法有了较大的进步,各种方法都有各自适用的范围和针对不同类型的混杂偏倚。1.倾向性评分法(PropensityScoreAnalysis,PSA)倾向性评分法主要用于均衡处理/暴露组间已观测混杂因素所产生的偏倚,该方19法是由Rosenbaum和Rubin在1983年首次提出,其基本原理为通过将多个协变量所产生的影响综合成为一个倾向性评分值,从而对协变量降维,然后利用所获得的倾向性评分值对不同对照组进行分层、匹配或者加权等处理,以达到协变量所产生的影响在组间均衡,最后利用均衡好的数据进行处理/暴露因素影响的分析。最初估计倾向性评分值的算法为logistic回归模型,该方法较为直观和常用,其公式如下:ഁబశഁభభశഁమమశ⋯శశഁೖೖP(Y=1|X)=(1-1)ଵାഁబశഁభభశഁమమశ⋯శశഁೖೖ式(1-1)中P(Y=1|X)表示为在K个协变量(X1,X2,……,XK)条件存在下事件Y发生的概率值,此为常规的logistic回归,其中的Y事件为我们所规定的处理20/暴露因素,P(Y=1|X)即为我们所需的倾向性评分值,然后我们可以通过分层,212223匹配(如贪婪匹配,马氏距离匹配等)或加权方法根据此倾向性评分值进行进一步的处理。此方法也是目前应用最广泛的倾向性评分方法,如丹麦的一项关于在二型糖尿病患者中分析不同胰岛素促泌剂和二甲双胍死亡率和心血管风险比较的研究24分析中,作者通过倾向性评分法对研究中的已观测混杂因素进行均衡,发现了格列-13- 基于多水平模型的工具变量方法研究及应用美脲、格列苯脲、格列吡嗪和甲苯磺丁脲与二甲双胍相比心血管保护作用稍差。目前,倾向性评分值的计算方法还有很多种,如人工神经网络(artificialneuralnetwork,ANN)、支持向量机(supportvectormachines,SVMs)、分类与回归树(catigoricaland25-28regressiontree,CART)、Boosting算法等,还有研究将倾向性评分法同贝叶斯思29,30想进行结合,这些算法各有自己的优缺点和适用范围,在很大程度上扩充了倾向性评分的适用范围。2.边缘结构模型(MarginalStructuralModels,MSM)31,32边缘结构模型是Robins于1997年首次提出的,其主要是用来控制观察性研究中所存在的时依性混杂因素,所谓的时依性混杂因素就是其随时间变化,是结局的影响因素,并且会影响随后的处理/暴露,同时又会受到前次处理/暴露的影响。边缘结构模型的核心思想主要是通过逆概率加权法(InverseProbabilityWeighting,IPW)为每个观察对象都赋予相应的权重,类似于构建出一个虚拟的人群,虚拟人群中的每个个体都相当于接受了处理/暴露因素的所有水平,然后对这个虚拟人群进行回归模型分析,以获得最终处理/暴露因素的效应值。3.多水平模型(MultilevelModel)由于观察性研究中所抽取的人群样本往往存在一定的聚集性,这样的数据具有层次结构特征(hierarchicalstructure),所谓的层次结构特征即为同一群体内的个体由于社会背景、生活习俗等相同而具有相似性,即反应变量的分布在个体间不具备独立性,存在地理距离内、某行政划区内或特定空间范围内的聚集性,对于这样的数据采用传统的回归模型分析时,由于模型拟合为考虑残差的分层,最终获得的结果可能同真实结果存在一定的偏差,所以多水平模型更加适合这样数据的分析。该类模型是20世纪70年代发展起来的,是将方差成分模型和多元回归模型结合为一体的分析方法,主要思想是将回归模型的差异拆分为固定效应和随机效应两部分,从而更加准确地估33,34计处理/暴露因素的效应值。该方法的优势在于考虑了数据分析时误差的层次性,还可以通过多元分析,有效地校正已观测混杂因素对处理/暴露因素效应值所产生的影响。以上所介绍的校正混杂因素的方法多为对已观测混杂因素进行的控制,但在现实的观察性研究中,由于试验设计和现实情况的局限,我们很可能会遗漏或无法观测一些重要的混杂因素,即使我们采用了一些混杂因素校正的方法,但对这些遗漏或无法观测的混杂因素所产生的影响没有进行控制,以至于所获得的结果将会同真实结果存在差异,对于未观测混杂因素所产生的影响,工具变量(InstrumentalVariable,IV)方法可以较好的进行解决。-14- 基于多水平模型的工具变量方法研究及应用(二)工具变量方法的研究现状与进展工具变量方法是Brookhart等于2006年首先将该方法从计量经济学中引入到观察35性研究中来,该方法通过所选取的工具变量采用回归方法将处理/暴露分解为与混杂因素相关和不相关的两部分,利用分解出的与同混杂因素不相关的部分同结局变量进行分析,消除处理/暴露因素与混杂因素间的关系,最终获得处理/暴露因素真实的效应值。工具变量的特点是a.工具变量同处理/暴露因素存在相关性,相关性越强,说明工具变量强度越强;b.工具变量同结局变量无直接相关关系;c.工具变量同其他混杂因素无相关关系。工具变量的选取需要以上三个条件。混杂因素工具变量处理/结局暴露变量注:实线表示有相关关系,虚线表示无相关关系图1-1工具变量所满足的三个特点根据数据类型的不同,目前工具变量的算法有多种:1.暴露因素和结局变量均为连续型变量。对于此类数据目前常用的工具变量方法为两阶段最小二乘法(2-StageLeastSquares,2SLS),该方法第一步是以处理/暴露因素为因变量,以工具变量和已知混杂因素为自变量进行最小二乘回归,求得干预措施的估计值;第二步是以结局变量为因变量,用第一步估计获得的处理/暴露因素估计值和已知混杂因素为自变量进行最小二乘回归,获得处理/暴露因素的效应值。目前针对于连续型变量,除了采用两阶36段最小二乘法外,还有学者采用两阶段残差纳入模型进行分析,此方法就是在第一步是计算出处理/暴露因素估计值同原始值之间的差值,即残差,再将残差替代第二步中的估计值进行运算。2.暴露因素为连续变量,而结局变量为分类变量目前对于这种类型的数据,第一步同两阶段最小二乘法相同,第二步利用logistic回归取代OLS回归,获得处理/暴露因素的效应值。-15- 基于多水平模型的工具变量方法研究及应用3.暴露因素和结局变量均为分类变量对于这样的数据类型,有学者认为可以采用两阶段最小二乘法进行估计,但是无法得出OR(OddsRatio)或RR(RelativeRisk),但是可以计算得到RD(RiskDifference)37值,如下式所示:ா[|ୀଵ]ିா[|ୀ]ߚ୍=(1-2)ா[|ୀଵ]ିா[|ୀ]式(1-2)其中β୍就是RD值,很多研究采用RD*100表示两组间的危险差异。目前还有研究采用两阶段logistic回归模型(Two-StagePredictorSubstitution,36,382SPS)进行分析,此方法可以得到通常研究中所需的OR或RR值。两阶段logistic回归模型第一步利用工具变量和协变量同处理/暴露因素做logistic回归,得到处理/暴露因素的概率预测值,第二步用概率预测值同结局变量进行第二步logistic回归。针对二分类数据目前还有一些工具变量算法可以计算,如GMM模型、两阶段Probit37,39,40模型、三阶段回归模型等方法。工具变量方法目前被应用的越来越广泛,其相关的文献发表数量也在逐年的递增,Pubmed数据库中工具变量相关文献发表情况从2006的50篇增长到2014年的202篇(如图1-2所示)。250202200文161献150141发表数(1009190篇)69565050500200620072008200920102011201220132014年份图1-2Pubmed数据库中工具变量相关文献发表情况(2006-2014)(三)常用混杂偏倚控制方法的优缺点以上所介绍的四种常用的混杂偏倚控制方法针对观察性研究数据中不同的问题,这些方法具有各自的适用范围和缺点,如倾向性评分法主要控制已观测混杂因素对最终分析结果产生偏倚,但对未知观测混杂控制不够理想;多水平模型主要是控制层次结构特征数据中不同水平对结果产生的偏倚,但其还是基于普通回归模型,对其他混-16- 基于多水平模型的工具变量方法研究及应用杂所产生的偏倚控制不足;边缘结构模型则是控制纵向观测数据中所存在的时依性混杂因素所产生的偏倚,但模型要求较高,当存在交互作用时模型应用需要十分注意;工具变量方法目的是控制未知观测混杂因素,此方法工具变量的选取较为严格(如表1-1所示)。表1-1四种控制混杂偏倚数据挖掘方法间的比较控制的混杂是否适用于层次结构方法优点局限性类型特征数据对未观测混杂控制不理想;可能存在残余实现较为方便、直观,混杂;匹配法的匹配已知观测混其中PS匹配法,可以倾向性评分部分适用规则有待进一步研杂直观地均衡组别间的究,样本量减少可能基线。会影响到最终结果。可以很好的控制未知强工具变量的选择是观测混杂因素对结果本方法的重点和难未知观测混工具变量暂不适用的影响,尤其在自变量点,弱工具变量会造杂和应变量均为连续型成最终结果的偏倚和变量的情况下。置信区间过宽。模型研究较为成熟,可以很好的应用于层次对未知混杂、共线性水平间差异结构特征数据中;可以多水平模型适用问题等方面的控制力造成的混杂在一定程度上控制已度不足。观测混杂因素所产生的偏倚对于模型的设定较为敏感;不能适用于在很好地控制时依性混适宜性混杂某一混杂条件下暴露边缘结构模型暂不适用杂因素;对多分类或连因素水平条件概率为1的续变量也同样适用;不情况下;协变量同暴存在残余混杂。露因素间应当避免关联或交互作用。针对这些控制混杂偏倚方法的缺点,研究者常通过对模型进行合理的修改来弥补,也有研究者通过多种方法的结合使用,对于其相互的优缺点进行弥补,如BoLu和41SueMarcus将倾向性评分法和工具变量进行结合,以便于同时控制数据中存在的已观测混杂偏倚和未观测混杂因素,从而得到更加可靠的分析结果。以上这四种控制混杂因素方法中,仅工具变量方法可以较好的控制未知观测混杂因素对结果产生的影响,但由于模型的局限性,目前该方法仅应用于普通结构数据中,-17- 基于多水平模型的工具变量方法研究及应用鲜有研究将该方法进行改进并应用于层次结构特征数据中。但是具有层次结构特征的数据越来越常见,其中较为典型的为全国卫生服务调查数据。全国卫生服务调查平均每五年进行一次,从1993年开始至今,共进行了五次全国卫生服务调查,此项调查设计了家庭基本情况、家庭成员的基本情况(如基本身体情况、医保情况、生活习惯、慢性病史等)、就医住院情况、儿童妇女相关健康情况,根据调查的结果,政府机构可以较为全面地掌握城乡居民健康状况、卫生服务利用、医疗保健费用及负担等信息,42-44这已经成为中国卫生调查制度的重要组成部分。卫生服务调查的抽样方法主要为分层随机抽样,其数据的基本结构为:省市—区县—街道—居委会—家庭—个人。抽样方法决定了数据具有明显的层次性,由于街道或居委会的地理位置、医疗机构设施或相关医疗健康政策存在一定特异性,这决定了被调查的个体人员的健康数据可能在某个街道或居委会这个群体中存在一定的相关聚集性。另外,此为调查数据,由于所覆盖的信息范围较为广泛,当进行某一项研究时,可能其所包含的信息不能完全满足研究所需的信息量,即可能存在一些遗漏或者无法观测的混杂因素。由此可见,卫生服务调查数据中除了存在一些已观测混杂因素的影响,还存在水平因素和未知观测混杂因素对结果产生的影响,基于此种情况,采用传统的分析方法可能会得出偏差较大的结果,故需要对工具变量模型进行研究和改进,使其更好地适用于层次结构特征数据中。二、研究目的与意义(一)研究目的目前对于层次结构特征数据中未观测混杂偏倚的控制鲜有研究,在这种类型的数据中常存在三种方面的混杂偏倚:已知观测混杂、未知观测混杂,以及水平因素产生的偏倚。如果采用传统的统计模型分析,可能无法全面控制此类数据中的混杂偏倚。针对这样的问题,本研究基于多水平模型和工具变量模型中控制不同混杂偏倚的思想,拟构建多水平工具变量模型(MultilevelInstrumentalVariable,MIV),以实现同时控制层次结构特征数据中存在的三种主要混杂偏倚。针对数据类型的不同,将构建不同的多水平工具变量模型,以分别适用于连续型变量和分类变量,并且探索比较不同工具变量方法同多水平模型结合时的稳定性和准确性,从而获得最佳的参数估计模型,在此模型构建的基础上,将引入自助法(Bootstrap),使模型估计得到的结果更加可靠。(二)研究意义本研究拟构建多水平工具变量模型,解决层次结构特征数据中未知观测混杂因素-18- 基于多水平模型的工具变量方法研究及应用难以控制的问题,为观察性研究中层次结构特征数据混杂因素的控制提供更为全面有效的手段,同时控制了已知观测混杂、未知观测混杂和水平因素三方面因素对结果产生的影响,提高了分析结果的准确性,增强分析结果的可靠性。最后将模型应用于全国卫生服务调查数据(上海)中,解决实际问题。三、研究内容本研究首先对目前多水平模型和常用的工具变量方法进行学习和进行相关的系统文献回顾,掌握这些方法的原理和应用条件;根据数据类型的不同,对连续型变量和分类变量分别进行模拟,并且设置不同强弱的混杂因素和工具变量;将构建的多个多水平工具变量模型和普通分析方法进行比较,从估计结果的准确度和精确度对多种方法进行评价,发现各方法存在的优缺点;最后将模型应用于实际数据中,解决实际问题。(一)含有未观测混杂因素数据的模拟在模拟分层数据时,首先生成已知观测混杂因素(C)、未知观测混杂因素(U)和工具变量(Z),然后根据以上三类变量生成处理/暴露因素(X),最后通过已知观测混杂、未知观测混杂和处理/暴露因素生成最终的结局变量(Y)。在数据模拟过程中考虑:(1)不同数据类型的处理/暴露因素和结局变量(连续型数据和二分类数据);(2)不同强弱的未知混杂因素;(3)不同强弱的工具变量。生成已知混杂C生成N个层,并确定各设定C、U、Z与干预措施X的生成未知混杂U层的截距ߚ关系生成X生成工具变量Z在各层中,设定C、U、X和结局将生成的数变量Y的关系生成Y据集中的未结局变量Y为连续型变量X为连续型变量知混杂U去除,得到最终分析的数结局变量Y为分类变量X为分类变量据集图1-3数据模拟流程(二)模型构建当处理/暴露因素和结局变量为连续型变量时,构建出四种多水平工具变量模型:-19- 基于多水平模型的工具变量方法研究及应用两阶段最小二乘多水平工具变量模型(MultilevelTwoStageLeastSquare,M2SLS)、两阶段残差纳入多水平工具变量模型(MultilevelTwoStageResidualInclusion,M2SRI)、两阶段多水平回归工具变量模型(TwoStageLinearMultilevelModel,2SLM)和两阶段多水平回归残差纳入工具变量模型(TwoStageLinearMultilevelResidualInclusionModel,2SLMRI)。当处理/暴露因素和结局变量为分类变量时,构建出两种多水平工具变量模型:两阶段logistic回归多水平工具变量模型(MultilevelTwoStagePredictionSubstitute,M2SPS)和线性回归+logistic回归多水平工具变量模型(MultilevelLeastSquarePredictionSubstitute,MLSPS)。对于这些模型算法,均将引入Bootstrap算法,分别构建出Boot-M2SLS(BootstrapMultilevelTwostageleastsquare)、Boot-M2SRI(BootstrapMultilevelTwoStageResidualInclusion)、Boot-2SLM(BootstrapTwoStageLinearMultilevelModel)、Boot-2SLMRI(BootstrapTwoStageLinearMultilevelResidualInclusionModel)、Boot-M2SPS(BootstrapMultilevelTwoStagePredictionSubstitute)和Boot-MLSPS(BootstrapMultilevelLeastSquarePredictionSubstitute)。具体模型构建思路如下:1.处理/暴露因素和结局变量为连续型变量(1)两阶段最小二乘多水平工具变量模型(M2SLS)该模型采用传统的两阶段最小二乘估计工具变量模型(2SLS)的思想,进行模型的构建。第一步:利用已知观测混杂C、工具变量Z同处理/暴露因素X做多元线性回归模型,通过多元线性回归模型对处理/暴露因素X进行重新估计,获得新的处理/暴露因素ܺ,其公式如下:X=ߙ+ߙܼ+ߙܥ+ߝଵ(1-3)X=ܧ[ܺ|ܼ,ܥ]第二步:由于各因素对结局变量的影响存在数据层次因素的干扰,所以在第二步运算时将进行多水平回归模型。第二步运算利用已知观测混杂因素C、重新估计得到的新的处理/暴露因素ܺ和结局变量Y做多水平线性回归模型,其公式如下:Y=ߚ+ߚଵܺ+ߚଶܥ+݁(1-4)式(1-4)中可以对ߚ进行分解,ߚ=ߚ+ݑݑ~ܰ(0,ߪଶ)݁~ܰ(0,ߪଶ)(1-5)௨బబ将(1-4)式和(1-5)式进行结合得到下式:Y=ߚ+ߚଵܺ+ߚଶܥ+(ݑ+݁)(1-6)式中(1-6)i=1,2,……,nj,表示为水平1单位(即所观测的个体数据),j=1,-20- 基于多水平模型的工具变量方法研究及应用2,……,J,表示为水平2单位(如社区、街道)。式(1-6)中ߚ,ߚଵ和ߚଶ为固定系数,ݑ和݁为两个随机参数,其中ߪଶ表示水平2单位的方差成分,݁表示水௨బ平1单位的方差成分。此为多水平模型中较为简单的方差成分模型表示公式。两阶段最小二乘多水平工具变量模型构建的基本思想就是用一步的多元回归模型重新估计处理/暴露因素,再用新得到的处理/暴露因素和结局变量进行水平线性回归模型拟合,得到处理/暴露因素的效应值。(2)两阶段残差纳入多水平工具变量模型(M2SRI)45两阶段残差纳入模型首先由Hausman在1978年应用在线性模型中,Joseph36V.Terza等人也认为在线性模型中2SRI和2SLS的结果是一致的。两阶段残差纳入模型在第一步估计时,计算出处理/暴露因素的残差,应用残差进行第二步运算。第一步:利用已知观测混杂C、工具变量Z同处理/暴露因素X做多元线性回归模型,重新估计得到新的处理/暴露因素ܺ,如式(1-3)。用新得到的处理/暴露因素ܺ减去原始的处理/暴露因素X,得到相应的残差ε。ε=E[X|ܼ,ܥ]−X(1-7)第二步:利用已知观测混杂因素C、处理/暴露因素X、处理/暴露因素的残差ε和结局变量Y做多水平线性回归模型,其公式如下:Y=ߚ+ߚଵX+ߚఌε+ߚଶܥ+(ݑ+݁)(1-8)式(1-8)中ߚଵ即为处理/暴露因素的效应值。2.处理/暴露因素和结局变量为连续变量并且均存在层级效应(1)两阶段多水平回归工具变量模型(2SLM)第一步:利用已知观测混杂C、工具变量Z同处理/暴露因素X做多水平线性回归模型,通过多元线性回归模型对处理/暴露因素X进行重新估计,获得新的处理/暴露因素ܺ,其公式如下:X=ߙ+ߙܼ+ߙܥ+ߝ+ߤ(1-9)X=ܧ[ܺ|ܼ,ܥ]第二步:由于各因素对结局变量的影响存在数据层次因素的干扰,所以在第二步运算时将进行多水平回归模型。第二步运算利用已知观测混杂因素C、重新估计得到的新的处理/暴露因素ܺ和结局变量Y做多水平线性回归模型,其公式如下:Y=ߚ+ߚଵܺ+ߚଶܥ+݁(1-10)式(1-10)中可以对ߚ进行分解,ߚ=ߚ+ݑݑ~ܰ(0,ߪଶ)݁~ܰ(0,ߪଶ)(1-11)௨బబ将(1-10)式和(1-11)式进行结合得到下式:-21- 基于多水平模型的工具变量方法研究及应用Y=ߚ+ߚଵܺ+ߚଶܥ+(ݑ+݁)(1-12)ߚଵ为最终的效应值。(2)两阶段多水平回归残差纳入工具变量模型(2SLMRI)第一步:利用已知观测混杂C、工具变量Z同处理/暴露因素X做多水平线性回归模型,重新估计得到新的处理/暴露因素ܺ,如式(1-9)。用新得到的处理/暴露因素ܺ减去原始的处理/暴露因素X,得到相应的残差ε。ε=E[X|ܼ,ܥ]−X(1-13)第二步:同两阶段残差纳入多水平工具变量模型的第二步运算。3.处理/暴露因素和结局变量为分类变量(1)两阶段logistic回归多水平工具变量模型(M2SPS)针对处理/暴露因素和结局变量均为分类变量时,我们常常需要获得OR或RR值,而通过先前应用于连续型变量的两种多水平工具变量模型是无法获取的,为了得到所需的OR或RR值,本研究将用普通logistic回归和多水平logistic回归模型对先前的多水平工具变量模型进行替换,构建出两阶段logistic回归多水平工具变量模型。普通logistic回归模型中,应变量Y为一个二分类变量,在m个自变量(X1,X2,……,Xm)作用下阳性结果发生的概率为:ଵP=P[Y=1|ܺଵ,ܺଶ,…,ܺ]=(1-14)ଵା୶୮[ି(ఉబାఉభభାఉమమା⋯ାఉ]对式(1-14)作对数转化,模型可化为如下线性形式:Lnቀቁ=ߚ+ߚଵܺଵ+ߚଶܺଶ+⋯++ߚܺ(1-15)ଵି估计logistic回归模型的参数时,目前通常采用的是最大似然估计(maximumlikelihoodestimate,MLE),即为构建出一个样本似然函数:L=∏୬ܲ(1−ܲ)ଵିୀଵ(1-16)其中ܲ表示第i例观测对象在暴露条件下阳性结果发生的概率,似然函数L达到最大值时即为模型的参数估计值。第一步:利用已知观测混杂C、工具变量Z同处理/暴露因素X做logistic回归模型,以得到处理/暴露因素X的概率预测值P(X)。ଵP(X)=PX=1หC,Z൧=(1-17)ଵା୶୮[ି(ఈబାఈୋఈೋ]第二步:利用所得到处理/暴露因素X的概率预测值P(X)、已知观测混杂C和结局变量Y做多水平logistic回归模型,解决数据层次因素对结果的干扰。其公式如下:ܻ=݈݃݅ݐ൫ܲ൯=ߚ+ߚଵܲ(ܺ)+ߚܥ(1-18)-22- 基于多水平模型的工具变量方法研究及应用ߚ=ߚ+ݑݑ~ܰ(0,ߪଶ)(1-19)௨బ由式(1-18)和式(1-19)进行合并得到下式:݈݃݅ݑ+ߚ=൯ܲ൫ݐ+ߚଵܲ(ܺ)+ߚܥ(1-20)式(1-20)中i表示水平1单位(即所观测的个体数据),j表示水平2单位(如社区、街道),ߚଵ即为处理/暴露因素的效应值,又称固定效应参数,ݑ为随机效应,又称高水平(水平2)的残差。(2)线性回归+logistic回归多水平工具变量模型(MLSPS)有研究者认为在第一步采用logistic回归时,由于并不清楚两变量logistic的误差分布,并且第一步中logistic回归估计所得的处理/暴露因素预测概率值P(x)并不能保证其与其余已观测混杂无相关关系,故JosephV.Terza和BingCai提出了线性回归+logistic回归模型(two-stagepredictorsubstitution,2SPS)应用于非层次结构特征数36,38据中。线性回归+logistic回归多水平工具变量模型的思想是将两阶段logistic回归多水平工具变量模型中第一步重新估计处理/暴露因素所采用的logistic回归用多元线性回归进行代替。第一步:利用已知观测混杂C、工具变量Z同处理/暴露因素X做多元线性回归模型,原始的处理/暴露因素X是二分类数据,通过多元线性回归模型重新估计的处理/暴露因素ܺ变为连续性变量,应用新估计所得的处理/暴露因素ܺ代入第二步运算中。X=ߙ+ߙܼ+ߙܥ+ߝଵ(1-21)X=ܧ[ܺ|ܼ,ܥ]第二步:利用所得到处理/暴露因素X的概率预测值ܺ、已知观测混杂C和结局变量Y做多水平logistic回归模型,公式如下:ܻ=݈݃݅ݐ൫ܲ൯=ߚ+ߚଵܺ+ߚܥ(1-22)ߚ=ߚ+ݑݑ~ܰ(0,ߪଶ)(1-23)௨బ式(1-22)中的ߚଵ即为处理/暴露因素的效应值。4.自助法同多水平工具变量模型结合46-51自助法是一种为统计推断所建立的重复抽样方法,利用计算机对样本进行复位抽样或放回抽样(samplingwithreplacement),该方法能够有效地校正参数估计偏倚,提高统计推断的准确性。其基本原理就是从原始样本中进行随机重复抽样,一般所抽的复样本的样本量与原始样本量一致,对反复所抽取的N个复样本用同一模型分别进行拟合分析,则模型参数被估计N次,从而可以从参数估计的自助法样本分49布来计算出参数估计的偏倚、偏倚校正估计值、标准误和置信区间。目前应用于层级结构特征数据的自助法有三种:(1)个例重复抽样法-23- 基于多水平模型的工具变量方法研究及应用(case-resampling);(2)参数残差自助法(parametricresidualbootstrap);(3)52非参数残差自助法(nonparametricresidualbootstrap)。为了节省模拟和分析时间,本研究在模拟阶段和实例分析阶段所用的自助法均为分层个例重复抽样法,重抽的复样本例数与原样本相同。在模型构建时将自助法分别同普通多水平回归模型和两种多水平工具变量模型进行结合,以增加模型的稳定性。(三)模型比较与评价针对普通多水平模型和所构建的多水平工具变量模型,根据模拟数据所获得的结果对各种模型从准确性、精确性、置信区间宽度和置信区间覆盖率四个方面进行比较,从而发现各种方法所存在的优缺点。评价指标详见表1-2。表1-2模型评价指标评估标准计算方法准确性̅绝对偏倚(Bias)ߜ=ߚመ−ߚ置信区间宽度所有ߚመ置信区间上限和下限的平均值精确性标准误(StandardError,SE)ܵܧ൫ߚመ൯置信区间覆盖率(Confidenceintervals所有ߚመ的置信区间中包含ߚ所占的比例coverage)̅注:表中ߚ为模拟数据中所设参数估计的金标准,ߚመ为每次分析是所获得的参数估计,ߚመ为多次模拟所获得参数估计的平均值,SE(ߚመ)为多次模拟所获得参数估计的标准误。(四)实例应用53本研究实例分析部分研究采用的数据为第五次全国卫生服务调查数据(上海),研究内容分为两部分。针对结局变量和处理/暴露因素为连续型变量,所选的实例为分析上海60岁以上老人每周体育锻炼时间对其健康状况的影响,男性和女性分别进行分析,其中可能存在的未知观测混杂包括本人的心情、家庭关系、病情轻重程度等,结局变量采用欧洲五维健康量表(Europeanqualityoflife5-dimensions,EQ-5D)评分54-56,工具变量选择为其爱人每周运动的次数。针对结局变量和处理/暴露因素为分类变量,本研究所选的实例为分析上海市60岁以上男性是否吸烟对其是否患有高血压的影响,其中可能存在的未知混杂包括基因特征、周围环境因素等,工具变量选择为其家人是否吸烟。实例分析中应用普通多水平模型和模拟中所获得的最优多水平工具变量模型进行分析,并对不同方法所获得的结果进行比较。-24- 基于多水平模型的工具变量方法研究及应用四、研究步骤本研究首先需要进行多水平模型和普通工具变量模型的文献回顾和方法学习,掌握这些方法的基本原理和适用范围。本研究主要分为两大部分:方法学研究部分和实例应用部分。方法学研究,即模型构建和数据模拟,首先需要进行含有未观测混杂因素的层次结构特征数据的模拟,在模拟数据的过程中需要考虑到强工具变量或弱工具变量和强未观测混杂因素或弱未观测混杂因素等情况,将前面所介绍的普通多水平模型和多水平工具变量模型应用于模拟数据中进行分析,根据所规定的模型比较评价指标对各种模型进行评价,发现模型存在的优缺点,然后对模型进行进一步的调整,从而确定最优模型。实例应用部分主要利用第五次全国卫生服务调查数据(上海),首先进行数据清理,然后进行已知观测混杂因素和工具变量的筛选,最后应用所建的最优多水平工具变量模型进行分析和结果的解释。本研究的技术路线如图1-4所示。-25- 基于多水平模型的工具变量方法研究及应用方法学研究实例应用数据模拟卫卫生服务调查查数据已知混杂强工具具变量强未观测混杂杂C数据清理理弱工具具变量弱未观测混杂杂已知混杂因素素筛选生生成暴露因素X和结局变量Y工具变量筛筛选暴暴露因素为连续续型变量的暴露露因素为分类变变量的多多水平工具变变量模型水水平工具变量模模型多水平工具具变量模型选择择方法的的评价比较结论解释释模型的的优化及确定图1-4研究技术路线图五、资料来源、分析工具及研究平台本研究模拟数据中包括了已知观测混杂因素、未知观测混杂因素,处理/暴露因素和结局变量,在进行模型运算时将不纳入未知观测混杂因素;实例数据分析采用的数据为第五次全国卫生服务调查数据(上海)。模拟数据的产生、实例分析数据的清理,模型构建运算和模型结果评价比较均采用的是StatisticsAnalysisSystem9.2(SAS9.2,64位)软件实现。最终模型结果展示采用MicrosoftExcel2010和Stata10.0软件进行实现。模拟数据的生成和模型运算的SAS宏程序附于论文中,并可重现。本研究所使用的计算机配置为InterCorei7CPU3.40GHz,内存16G,硬盘空间1TB,操作系统为Windows7(64位)系统。本研究受到第二军医大学2012年研究生苗子培育基金“医疗卫生数据挖掘方法研究及应用”和第二军医大学2014年博士研究生创新基金“基于多水平模型的工具变量方法研究及应用”两项基金资助,还受到上海市软科学研究重点项目“大数据与新兴业态关键技术研究——医疗卫生大数据-26- 基于多水平模型的工具变量方法研究及应用分析技术研究”(编号:14692101700)。教研室在数据挖掘方法和工具变量模型研究的经验,并有专家教授进行指导,也为本研究的完成奠定了坚实的基础。-27- 基于多水平模型的工具变量方法研究及应用第二部分多水平工具变量模型在模拟数据中的拟合和验证数据模拟,即数据仿真,是模拟流行病学研究中所遇的不同情况,也是流行病学研究中常用于模型构建评价的一种研究手段。在数据模拟中可以应用计算机对各种协变量、处理/暴露因素和结局变量进行模拟,并根据实际情况对参数间的关系进行调整,在模拟的过程中也可以很好地掌握模拟数据的真实情况,并且可以进行多次模拟,故可以很好地掌握各种模型的优劣。本研究主要是根据不同的处理/暴露因素和结局变量类型,将数据模拟过程分为两个主要部分:数据为连续型变量和分类变量。在模拟这两类数据时均含有未知观测混杂因素,并考虑到强未知观测混杂因素、弱未知观测混杂因素和无未知观测混杂因素(未知观测混杂因素和结局变量间的回归系数为0),根据不同强度的未知观测混杂因素存在的情况下,以观测模型的稳定性和准确性。在模拟数据集中还需要考虑到强工具变量、中度工具变量和弱工具变量三种情况,分别应用这三种强弱的工具变量来验证模型,以观测不同强弱工具变量对模型的影响。针对连续型数据,将分别应用普通多水平线性回归模型(MultilevelLinearModel,MLM)、两阶段最小二乘多水平工具变量模型(M2SLS)、两阶段残差纳入多水平工具变量模型(M2SRI)、自助法多水平线性回归模型(Boot-MLM)、自助法两阶段最小二乘多水平工具变量模型(Boot-M2SLS)和自助法两阶段残差纳入多水平工具变量模型(Boot-M2SRI)六个模型进行拟合分析;当处理/暴露因素在水平2单位上存在异质性,将分别应用普通多水平线性回归模型、两阶段多水平回归工具变量模型(2SLM)、两阶段多水平回归残差纳入工具变量模型(2SLMRI)、自助法两阶段多水平回归工具变量模型(Boot-2SLM)、自助法两阶段多水平回归残差纳入工具变量模型(Boot-2SLMRI)和自助法两阶段最小二乘多水平工具变量模型六个模型进行拟合分析;针对分类数据,将分别应用普通多水平logistic回归模型(MultilevelLogisticRegressionModel,MLRM)、两阶段logistic回归多水平工具变量模型(M2SPS)、线性回归+logistic回归多水平工具变量模型(MLSPS)、自助法多水平logistic回归模型(Boot-MLRM)、自助法两阶段logistic回归多水平工具变量模型(Boot-M2SPS)和自助法线性回归+logistic回归多水平工具变量模型(Boot-MLSPS)六个模型进行拟合分析。最后将这六种模型的结果进行评价比较,从而获得这六种模型的优劣,以便后续在实例中进行应用。-28- 基于多水平模型的工具变量方法研究及应用一、处理/暴露因素和结局变量为连续型变量(一)模拟数据集构建1.已知观测混杂因素构建模拟数据集中设两个已知观测混杂因素C1和C2,这两个已知观测混杂因素均为正态分布数据,如下式:C1~N(4,1)(2-1)C2~N(2,0.5)(2-2)2.未知观测混杂因素构建模拟数据中设一个未知观测混杂因素U,此未知观测混杂因素为正态分布数据,如下式:U~N(5,1)(2-3)3.工具变量构建模拟数据中设一个工具变量Z,此工具变量为正态分布数据,如下式:Z~N(3,0.5)(2-4)4.处理/暴露因素构建模拟数据中处理/暴露因素由已知观测混杂因素C1、C2、工具变量Z和未知观测混杂因素U构建回归方程计算得到,回归方程如下:X=ߙ+ߙଵܥଵ+ߙଶܥଶ+ߙ௨ܷ+ߙ௭ܼ(2-5)式(2-5)中ߙଵ、ߙଶ分别设为8、-6,ߙ௨设为-2。ߙ௭表示工具变量的强弱,在模拟过程中分别设为1、3、5,分别表示弱、中、强三种类型的工具变量。5.结局变量构建模拟数据中结局变量由已知观测混杂因素C1、C2、未知观测混杂因素U和处理/暴露因素X构建回归方程计算得到,回归方程如下:Y=ߚ+ߚଵܥଵ+ߚଶܥଶ+ߚ௨ܷ+ߚ௫ܺ(2-6)式(2-6)中已知观测混杂因素的回归系数ߚଵ、ߚଶ分别设为3、-5,未知观测混杂因素的回归系数ߚ௨在模拟的过程中分别设为0、2、4、6,分别表示无、弱、中、强四种类型的未知观测混杂因素。ߚ௫为处理/暴露因的实际效应值,在本研究中设为4,此也为本研究的金标准。详见表2-1,本研究共模拟12种情况下数据集。-29- 基于多水平模型的工具变量方法研究及应用表2-1处理/暴露因素和结局变量为连续型变量参数取值数据集未知观测混杂因素强度(ߚ௨)工具变量强度(ߙ௭)模拟101模拟203模拟305模拟421模拟523模拟625模拟741模拟843模拟945模拟1061模拟1163模拟1265由于本研究针对的是层次结构特征数据,在模拟数据中设第二水平共有200层,每层中含有5个个体数据,每次模拟的数据集大小为1000例,在模拟结局变量Y时,每层的常数项ߚ(j=1,2,3,……,200)服从均匀分布,如下式:ߚ~U(−10,10)(2-7)数据模拟过程中,对每一种情况下的模拟数据进行模拟500次。模拟完成的数据集如表2-2所示。表2-2处理/暴露因素和结局变量为连续型变量数据模拟示例sitec1c2uzxy14.02221.32354.88711.921934.4187153.190711.82522.23846.18463.115514.726765.852314.74571.90384.45112.823642.1056182.335013.53242.08045.88103.759933.1609144.892813.19762.08854.24002.859029.2116124.768923.93252.31695.27763.750322.910397.332521.55012.04503.65422.71243.539510.814922.02662.19064.83132.79534.537917.863924.24452.33364.98403.260223.441899.723523.94512.13255.07822.659019.059282.489334.97711.74936.46172.710714.934770.999435.10582.37697.68112.08836.648437.540634.17271.69385.78192.22397.758038.796834.20032.63145.23092.59485.309023.293134.95562.79475.41172.928311.677650.5789-30- 基于多水平模型的工具变量方法研究及应用sitec1c2uzxy44.52592.30794.64572.496538.6354167.072144.22701.97023.52002.437340.2269171.978642.82702.31074.62423.131828.2483120.370443.80292.10586.32732.288629.6625133.385644.74610.98774.98902.953849.9186220.153255.42201.41265.13722.419319.323592.474550.91901.58174.63692.1679-17.9708-72.057355.32931.54444.34633.462024.5862111.006554.60812.57385.92972.91266.726835.4251……………………………………注:表中site表示水平2单位。(二)模型在模拟数据中的应用模拟数据在分析的时候,首先将数据集中未知观测混杂因素U去除,并不纳入模型分析中,故分析数据集中存在未知观测混杂因素。1.普通多水平线性回归模型在普通多水平线性回归模型中纳入C1、C2、X和结局变量Y,其公式如下:Y=ߚ+ߚ௫ܺ+ߚଵܥ1୧୨+ߚଶܥ2୧୨+(ݑ+݁)(2-8)式(2-8)中i=1,2,3,4,5,表示为水平1单位;式中j=1,2,……,200,表示为水平2单位(下同);ߚ௫为处理/暴露因素的效应值。2.两阶段最小二乘多水平工具变量模型第一步:利用已观测混杂因素C1、C2、工具变量Z同处理/暴露因素X做多元线性回归,重新估计得到新的处理/暴露因素X:X=ߙ+ߙܼ+ߙଵܥ1+ߙଶܥ2(2-9)X=ܧ[ܺ|ܼ,ܥ1,C2]第二步:利用C1、C2、X和结局变量Y做多水平线性回归模型:Y=ߚ+ߚ௫X+ߚଵܥ1୧୨+ߚଶܥ2୧୨+(ݑ+݁)(2-10)式(2-10)中ߚ௫为处理/暴露因素的效应值。3.两阶段残差纳入多水平工具变量模型第一步:利用C1、C2、Z和处理/暴露因素X做多元线性回归并估计得到X,再计算相应的残差ε:ε=E[X|ܼ,ܥ]−X(2-11)第二步:利用C1、C2、ε、X和结局变量Y做多水平线性回归模型:Y=ߚ+ߚ௫X+ߚఌε+ߚଵܥ1୧୨+ߚଶܥ2୧୨+(ݑ+݁)(2-12)-31- 基于多水平模型的工具变量方法研究及应用式(2-12)中ߚ௫为处理/暴露因素的效应值。4.自助法与模型结合在模拟数据中所采用的自助法为分层个例重复抽样法,分别在各层中(即在水平2单位内)进行有放回的重复抽样,每层中有放回地重抽出5例样本,共在模拟数据集中的200层内进行,故所得的复样本和原始样本大小相同。自助法共反复抽500次,共得500个复样本,再分别利用普通多水平线性回归模型、两阶段最小二乘多水平工具变量模型和两阶段残差纳入多水平工具变量模型三种模型对反复重抽的复样本进行运算,每个模型将估计共得到500个样本参数,再根据这500个样本参数重新计算出最终的效应值、标准误和置信区间。以自助法多水平线性回归模型为例:首先通过自助法重新获取了k个复样本,然后利用普通多水平线性回归模型对这k个复样本分别进行分析:Y=ߚ+ߚ௫ܺ+ߚଵܥ1୧୨୩+ߚଶܥ2୧୨୩+(ݑ+݁)(2-13)式中k=1,2,3,……,500,故得到k个ߚ௫,其均值തߚത௫ത为最终处理/暴露因素的效应值,效应值的置信区间的计算为തߚത௫ത±Z1−0.05/2SD൫ߚx൯。两阶段最小二乘多水平工具变量模型和两阶段残差纳入多水平工具变量模型的自助法计算方法同普通多水平线性回归模型相同。(三)结果比较本研究针对处理/暴露因素和结局变量均为连续型变量的模型拟合共构建和纳入了6个模型。1.模型准确度的比较表2-3展示了这6种模型的准确度。当无未知观测混杂因素(β୳=0)时,我们发现这6个模型均能获得较为准确的点估计结果,绝对偏倚均在0.01以下,其中以普通多水平线性回归模型和自助法多水平线性回归模型所得的结果最为准确,可以准确地获得点估计。多水平工具变量模型随着工具变量强度的增加,其点估计的准确度也有相应的提高。模型置信区间宽度随着工具变量强度增加而变窄。当模型中的未知混杂因素强度逐渐增强时,普通多水平线性回归模型和自助法多水平线性回归模型的估计偏差越来越大,其偏差最大的出现在β୳=6、ߙ=1时的普通多水平线性回归模型中,其绝对误差为-2.8219。对于其余四种多水平工具变量而言,我们发现当工具变量的强度较弱时,其模型估计的绝对偏倚较大,最大的偏差为β୳=4、ߙ=1时自助法两阶段最小二乘多水平工具变量模型所得的结果(绝对偏倚为0.0727)。当工具变量强度增加时,四种多水平工具变量模型所得结果也越准确。-32- 基于多水平模型的工具变量方法研究及应用当β୳=6、ߙ=5时两阶段最小二乘多水平工具变量模型、两阶段残差纳入多水平工具变量模型、两阶段最小二乘多水平工具变量模型和两阶段残差纳入多水平工具变量模型结果的绝对偏倚分别为-0.0004、-0.0009、0.0012和0.0006。模型中自助法所得结果的置信区间较未引入自助法模型所获得的置信区间宽,其中置信区间最宽的为自助法两阶段最小二乘多水平工具变量模型。同无未知观测混杂因素情况下所显示的趋势类似,模型的置信区间的宽度受到工具变量强度影响较大,区间宽度在强工具变量时较窄。2.模型精确度的比较表2-4展示了这6种模型的精确度。当无未知观测混杂因素时,我们发现多水平线性回归模型、两阶段残差纳入多水平工具变量模型、自助法多水平线性回归模型和自助法两阶段残差纳入多水平工具变量模型四种方法的标准误较小。针对模型覆盖率,水平线性回归模型和两阶段残差纳入多水平工具变量模型的参数覆盖率均为100%,自助法两阶段最小二乘多水平工具变量和自助法两阶段残差纳入多水平工具变量模型的覆盖率也较好,但两阶段最小二乘多水平工具变量模型的参数估计覆盖率较低,均在50%以下。当模型中存在未知观测混杂因素时,普通多水平线性回归模型和自助法多水平线性回归模型的参数估计覆盖率均为0%,说明这两种方法在存在未知观测混杂因素时,模型估计的精确度很差。两阶段最小二乘多水平工具变量模型和两阶段残差纳入多水平工具变量模型的参数估计覆盖率随着未知观测混杂因素强度的增加而下降,也说明此两种方法的参数估计的精确度也较差。模型参数估计覆盖率较好的为自助法两阶段最小二乘工具变量模型和自助法两阶段残差纳入多水平工具变量模型,均在90%以上,结合前面模型准确度评价的结果可以发现,其模型参数估计覆盖率较好的原因可能是因为其所估计的参数置信区间较宽,从而导致其包含金标准(ߚ୶=4)的可能性较大。这两个模型的参数估计标准误SE(β୶)较大,最大的为β୳=6、ߙ=1时自助法两阶段最小二乘工具变量模型(SE(β୶)=0.4109)。-33- 区间上限4.00614.00214.00124.28644.07604.04404.57124.15754.09384.77784.23144.1346区间下限3.99333.99803.99883.78263.92693.95533.56763.85833.91503.32993.78313.8666-34-00绝对偏倚-0.00030.03450.0015-0.00030.06940.00790.00440.05380.00720.0006区间上限4.85164.27054.16144.79644.23634.14144.81384.23574.13924.86214.25784.1512区间下限)3.16453.73103.84103.30473.77083.86293.33163.78173.86853.24943.75633.8511Boot-M2SLSBoot-M2SRI=4୶ߚ绝对偏倚0.00800.00080.00120.05060.00350.00220.07270.00870.00380.05570.00710.0012基于多水平模型的工具变量方法研究及应用区间上限4.00144.00124.00093.07413.39153.64062.14732.78003.28181.22312.17052.9232区间下限3.99853.99883.99913.04393.33023.57862.08702.65713.15681.13311.98622.7364000绝对偏倚-0.9410-0.6392-0.3904-1.8829-1.2814-0.7807-2.8219-1.9216-1.1702区间上限4.14214.04654.02774.16124.04644.02694.18114.05124.03134.15084.04924.0269区间下限的置信区间上限和下限的平均值。୶3.85743.95363.97233.87733.95373.97143.89773.95833.97563.87233.95653.9713መߚ00绝对偏倚-0.00030.01920.0001-0.00090.03940.00480.00350.01150.0028-0.0009区间上限4.13674.04534.02844.16784.04774.02914.17804.05134.03064.14934.04864.0274;区间上下限是由所有结局变量和混杂因素均为连续型变量时各方法准确度的比较(区间下限௫3.85193.95243.97303.88393.95493.97363.89463.95843.97493.87083.95603.9718−ߚ2-3௫̅መ表绝对偏倚-0.0057-0.00110.00070.02580.00130.00130.03630.00490.00280.01010.0023-0.0004ߜ=ߚ区间上限4.03384.02784.02163.09273.38873.63132.15082.74633.24101.21172.10592.8516区间下限MLMBoot-MLMM2SLSM2SRI3.96623.97223.97843.02523.33303.58802.08352.69083.19751.14452.05042.8082000绝对偏倚注:表中绝对偏倚表示为-0.9410-0.6392-0.3904-1.8828-1.2814-0.7808-2.8219-1.9218-1.1701)ߙ135135135135工具变量强度()β୳0246混杂因素强度( (%)覆盖率)Boot-M2SRI୶-35-β(SE(%)覆盖率)Boot-M2SLS୶β(SE基于多水平模型的工具变量方法研究及应用(%)覆盖率)୶β(SE(%)覆盖率)୶β(SE(%)覆盖率结局变量和混杂因素均为连续型变量时各方法精确度的比较)2-4୶β表(SE(%)覆盖率MLMBoot-MLMM2SLSM2SRI)୶β(0.00071000.261141.40.00291000.000893.20.264999.80.003793.80.008600.239648.60.1416690.008700.24771000.148692.00.017100.316435.20.285339.20.017100.334198.20.300693.20.024400.390231.20.386530.20.024400.410994.60.405793.0SE)(ߙ130.00061000.082343.60.00101000.000789.60.08281000.001291.250.00041000.048545.40.00051000.000593.40.048999.80.000794.4130.016400.073948.20.041572.80.016400.07431000.041893.250.017400.048242.60.025176.80.017400.048599.20.025191.4130.034700.089637.80.083141.60.034800.089999.40.083393.450.036700.056738.20.052642.00.036800.056698.60.052489.8130.053300.126228.60.125226.80.053300.126995.60.125993.650.051900.079823.60.078626.00.052000.080192.80.078990.8强度工具变量)୳(β0246混杂因素强度 基于多水平模型的工具变量方法研究及应用图2-1~图2-4直观地展示了这六种模型应用不同强度工具变量模型在不同强度未知观测混杂因素存在的情况下参数估计值和置信区间。总体来看,四种多水平工具变量所得结果较为准确,基本很接近所设的金标准,但对于普通的多水平线性回归模型和自助法多水平线性回归模型当资料中存在未知观测混杂因素时,其获得的结果明显偏离所设的金标准。自助法多水平工具变量模型所得结果较为保守,其区间较宽。从图中我们可以发现多水平工具变量模型的参数估计随着工具变量强度的不断增加也越来越接近我们所设的金标准,并且其置信区间也在逐渐缩短。图2-1连续型数据无未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=0)-36- 基于多水平模型的工具变量方法研究及应用图2-2连续型数据弱未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=2)图2-3连续型数据中度未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=4)-37- 基于多水平模型的工具变量方法研究及应用图2-4连续型数据强未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=6)二、处理/暴露因素和结局变量为连续型变量(均存在层级效应)(一)模拟数据集构建该部分的数据集建立过程同前面连续型变量的模拟数据集建立类似,仅在处理/暴露因素构建过程中,将原先的ߙ改为ߙ,回归方程如下:X=ߙ+ߙଵܥଵ+ߙଶܥଶ+ߙ௨ܷ+ߙ௭ܼ(2-14)常数项ߙ(j=1,2,3,……,200)服从均匀分布,如下式:ߙ~U(−20,20)(2-15)整体数据共模拟500次。数据集中所设的未知观测混杂因素强度和工具变量强度的取值共有12种情况,同前面连续型变量模拟中相同,详见表2-5。研究中金标准ߚ௫同样设为4。-38- 基于多水平模型的工具变量方法研究及应用表2-5处理/暴露因素和结局变量为连续型变量并且均存在层级效应的数据模拟示例sitec1c2u1za0b0xy14.21512.10875.43962.717210.34690.292334.1225149.763213.29641.73194.03212.238510.34690.292329.4546127.404812.92681.86464.38992.665910.34690.292327.1233117.022713.53081.69124.82452.525210.34690.292331.4232137.770414.87842.04905.00133.646110.34690.292345.3078195.916423.35872.66525.23133.2850-2.8451-5.077013.995358.116522.81451.84253.48482.7213-2.8451-5.077015.253062.135925.93512.46514.18043.3995-2.8451-5.077038.4818162.690824.47673.17417.06423.9407-2.8451-5.077019.498884.606024.83091.33476.65862.7614-2.8451-5.077028.2837129.194133.78731.50395.10113.8188-15.0167-7.848215.149566.794235.19412.12684.52911.9715-15.0167-7.848214.575364.459834.14222.46186.06133.0212-15.0167-7.84826.334029.728235.50312.11433.99702.2160-15.0167-7.848219.408583.717834.83651.55054.25782.2136-15.0167-7.848216.924675.122443.64841.76514.27293.478413.08521.201140.5278173.977742.55331.79116.11443.442713.08521.201127.7499123.134043.14912.08674.50063.591013.08521.201134.7119148.063944.14791.60836.03662.207913.08521.201135.5849160.016144.62661.96296.09942.730913.08521.201139.7755176.567054.32092.04213.91921.8901-17.3986-4.29666.528332.407452.89862.36944.48333.0333-17.3986-4.2966-2.2268-7.388955.36182.29575.36233.3014-17.3986-4.296617.503781.049352.64392.08223.63342.7876-17.3986-4.2966-2.0694-7.7865注:表中a0和b0分别表示在模拟生成处理/暴露因素和结局变量时所取的截距。(二)模型在模拟数据中的应用该部分的构型构建过程已在前面叙述,主要是将两阶段最小二乘多水平工具变量模型和两阶段残差纳入多水平工具变量模型中第一步重新估计处理/暴露因素和估计残差中所用的最小二乘回归用多水平线性回归模型进行替代,最后形成两阶段多水平回归工具变量模型、两阶段多水平回归残差纳入工具变量模型、自助法两阶段多水平回归工具变量模型和自助法两阶段多水平回归残差纳入工具变量模型,再将这四个模型所获得的结果同普通多水平线性回归模型和自助法两阶段最小二乘多水平工具变量模型进行比较。在引入自助法时,同前面分析时类似,也采用分层个例重复抽样法进行,分200层进行有放回的重复抽样,共抽500次,然后利用构建的两阶段最小二乘多水平工具变量模型、两阶段多水平回归工具变量模型和两阶段多水平回归残差纳入工具变量模型对这500个复样本进行拟合,共得500个样本参数,再利用这500个样本参数来估-39- 基于多水平模型的工具变量方法研究及应用计处理/暴露因素的效应值、标准误和置信区间。(三)结果比较本研究针对处理/暴露因素和结局变量均为连续型变量并且均存在层级效应的模型拟合共构建和纳入了6个模型。1.模型准确度的比较表2-6展示了6种模型的准确度。当无未知观测混杂因素(β୳=0)时,我们发现除自助法两阶段最小二乘多水平工具变量模型外,其余5种模型均能较为准确地估计处理/暴露因素的效应值,绝对偏倚均在0.01以下,其中以两阶段多水平回归工具变量模型和自助法两阶段多水平回归工具变量模型所得结果最为准确,尤其在使用强工具变量时(ߙ=5),这两个模型的绝对偏倚为0。自助法两阶段最小二乘多水平工具变量模型的置信区间较宽,其余模型的置信区间均很窄,说明模型估计较为准确。当模型中存在未知观测混杂因素时,普通多水平线性回归的结果出现了偏倚,并且随着未知观测混杂因素强度的不断增强模型结果的偏倚增大,当β୳=6、α=1时,其绝对误差为-2.8038。虽然自助法两阶段最小二乘多水平工具变量模型也可以大致估计出最终的效应值,但结果依然存在一定的偏倚,最大的偏倚为0.3681(β୳=6、α=3),并且该模型所得的置信区间宽度较宽,结果较为保守。其余四种多水平工具变量模型均能获得准确的效应值点估计,最大的偏倚仅为-0.0446(β୳=6、α=1,自助法两阶段多水平回归残差纳入工具变量模型),并且这四种模型的置信区间比较窄,说明模型估计的较为准确,并且模型结果受到工具变量强度强度较小。2.模型精确度的比较表2-4展示了这6种模型的精确度。当数据中无未知观测混杂因素时,我们发现多水平线性回归模型、两阶段多水平回归工具变量模型和两阶段多水平回归残差纳入工具变量模型三个模型的置信区间覆盖率均为100%,但是自助法两阶段多水平回归残差纳入工具变量模型的置信区间覆盖率较低,最高仅为58.4%。其余两个模型的覆盖率均在90%左右。参数估计标准误SE(β୶)方面除了自助法两阶段最小二乘多水平工具变量模型外其余模型的标准误均在0.003以下。当模型中存在未知观测混杂因素时,普通多水平线性回归模型的参数估计覆盖率均为0%,该结果类似于前面连续型变量但处理/暴露因素无层次效应时的结果,进一步说明该方法并不适用于数据中存在未知观测混杂因素的情况,尽管此模型在参数估计标准误方面表现较好。两阶段多水平回归工具变量模型和两阶段多水平回归残差纳入工具变量模型的置信区间覆盖率随着未知观测混杂因素强度的增加而下降,说明这两种模型在数据中存在强未知观测混杂因素情况下的结果精确度并不佳。但这两种模-40- 基于多水平模型的工具变量方法研究及应用型的参数估计标准误SE(β୶)也较小。三种自助法多水平工具变量模型在置信区间覆盖率方面表现较好,基本均在90%以上,自助法两阶段多水平回归工具变量模型和自助法两阶段多水平回归残差纳入工具变量模型的参数估计标准误SE(β୶)较小,并且受到工具变量强度的影响较小。自助法两阶段最小二乘多水平工具变量模型的参数估计标准误较大,均在0.5以上,最大为0.9803。自助法两阶段最小二乘多水平工具变量模型的置信区间覆盖率较高的原因可能为其模型参数估计的置信区间宽度较宽,包含金标准的概率高。图2-5~图2-8直观地展示了这六种模型应用不同强度工具变量模型在不同强度未知观测混杂因素存在的情况下参数估计值和置信区间。图中可见普通多水平线性回归模型在数据中存在未知观测混杂因素时所得的结果很不理想,但是两阶段多水平回归工具变量模型、两阶段多水平回归残差纳入工具变量模型、自助法两阶段多水平回归工具变量模型和自助法两阶段多水平回归残差纳入工具变量模型四个模型在不同的数据情况下均能获得较为理想的结果。自助法两阶段最小二乘多水平工具变量模型虽然能够基本准确地估计出处理/暴露因素的效应值,但模型的置信区间宽度过宽,提示模型估计准确性较差。-41- 区间上限区间下限-42-绝对偏倚区间上限)=4୶ߚ区间下限绝对偏倚基于多水平模型的工具变量方法研究及应用区间上限M2SLSMBoot-2SLMRIBoot-区间下限绝对偏倚区间上限区间下限绝对偏倚区间上限区间下限结局变量和混杂因素为连续型变量且均存在层级效应时各方法准确度的比较(绝对偏倚2-6表区间上限区间下限Boot-2SLMLM2SLM2SLMRI绝对偏倚-0.00143.96794.0294-0.00123.9325-0.00103.97294.06524.0251-0.00043.9601-0.00063.9786-0.00634.03924.02033.9274-0.89313.0741-0.00224.060103.1398-0.00933.92373.9582-0.0012-0.60443.36863.97394.03734.05773.99363.4226-0.00103.95934.0260-0.0004-0.36863.6103-0.01324.00404.03873.9959-0.00093.6525-0.00143.97253.9197-0.0065-1.85242.1142-0.00244.00323.97304.05383.98894.02462.1811-0.01893.91223.9579-0.00264.0251-0.01573.9981-1.25312.7194-0.00204.03733.99514.05003.9356-0.25151.33686.16022.77450.00153.96133.9719-0.00263.9997-0.75903.2196-0.021404.03304.02403.93884.04173.2624-0.00513.96873.90970.32052.32996.3110-0.0198-0.00193.9974-2.80381.16264.05594.04750.00063.93184.02113.95444.00251.2297-0.02273.9052-0.0042-0.03174.0287-0.0055-1.89872.07363.96054.0419-0.00113.93734.04933.8596-0.20631.29446.29303.96832.1290-0.00653.95294.0408-0.00263.99764.0543-1.15322.82534.0770-0.02394.0207-0.00093.95364.00024.03410.36662.35836.37482.86840.00133.9749-0.03433.9041-0.00563.87594.04120.18523.04785.3225-0.00693.85724.04813.90524.12230.15893.03575.28214.02773.95254.0741-0.04344.0836-0.00194.03370.0011-0.31921.26026.10153.7622-0.00613.8749-0.00773.97474.15103.90474.12143.79544.0275-0.04464.08320.36102.34306.37904.18913.76120.00090.18783.03945.3361-0.00814.14973.86293.7950-0.34461.22496.08594.13894.18870.00080.36812.35146.38483.86274.13880.10873.01075.2067)ߙ135135135135工具变量强度()β୳0246混杂因素强度( (%)覆盖率)Boot-M2SLS୶β(-43-SE(%)覆盖率)୶β(SE(%)基于多水平模型的工具变量方法研究及应用MBoot-2SLMRI覆盖率)୶β(SE(%)覆盖率)୶β(SE(%)覆盖率)୶β(结局变量和混杂因素均为连续型变量且均存在层级效应时各方法精确度的比较SE2-7(%)表覆盖率Boot-2SLMLM2SLM2SLMRI)୶β(SE)(ߙ10.00051000.00281000.00251000.002987.60.002620.80.959693.230.00051000.00181000.00111000.002190.40.001338.00.951093.650.00041000.00111000.00061000.001395.80.000758.40.567996.610.009200.033394.00.033193.40.027790.60.027589.40.980392.630.018000.034174.60.034174.80.032093.60.031993.40.940696.050.016800.024372.20.024372.60.023892.60.023792.80.623094.610.017400.071669.20.071470.20.062091.00.061890.60.918492.630.033600.070244.00.070244.00.067793.20.067793.20.944294.250.033900.049638.40.049638.00.048792.00.048792.20.637696.410.025100.116848.60.1167490.108395.40.108295.60.954592.630.051000.114127.00.114127.00.112793.20.112793.20.949394.850.051700.076428.20.076428.40.076092.40.076092.40.599592.8强度工具变量)୳(β0246混杂因素强度 基于多水平模型的工具变量方法研究及应用图2-5连续型数据且均存在层级效应无未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=0)图2-6连续型数据且均存在层级效应弱未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=2)-44- 基于多水平模型的工具变量方法研究及应用图2-7连续型数据且均存在层级效应中度未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=4)图2-8连续型数据且均存在层级效应弱未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=6)-45- 基于多水平模型的工具变量方法研究及应用三、处理/暴露因素和结局变量为分类变量(一)模拟数据集构建1.已知观测混杂因素构建模拟数据集中设两个已知观测混杂因素C1和C2,已知观测混杂因素C1为正态分布,已知观测混杂因素C2为二分类分布:C1~N(2,0.5)(2-16)C2~Bernoulli(0.4)(2-17)2.未知观测混杂因素构建模拟数据集中设有一个未知观测混杂因素U,未知观测混杂因素U设定为正态分布数据,如下式:U~N(1,0.5)(2-18)3.工具变量构建模拟数据集中设有一个工具变量Z,此工具变量为二分类变量,如下式:Z~Bernoulli(0.5)(2-19)4.处理/暴露因素构建模拟数据中处理/暴露因素X是由已知观测混杂因素C1、C2、工具变量Z和未知观测混杂因素U构建回归方程计算得到,处理/暴露因素X设定为二分类分布数据,计算公式如下:ଵP(X)=PX=1หC1,C2,U,Z൧=(2-20)ଵା୶୮[ି(ఈబାఈభେଵାఈమେଶାఈೠାఈ]X=Bernoulli(P(X))式(2-20)中已知观测混杂因素回归系数ߙଵ、ߙଶ分别设为-1、1,未知观测混杂因素回归系数ߙ௨设为-1。ߙ௭表示工具变量的强弱,在模拟过程中分别设为1、2、3,分别表示弱、中、强三种类型的工具变量。5.结局变量构建模拟数据集中结局变量Y由已知观测混杂因素C1、C2、未知观测混杂因素U和处理/暴露因素X构建回归方程计算得到,结局变量Y设定为二分类分布数据,计算公式如下:ଵP(Y)=PY=1หC1,C2,U,X൧=(2-21)ଵା୶୮[ି(ఉబౠାఉభେଵାఉమେଶାఉೠାఉଡ଼]Y=Bernoulli(P(Y))式(2-21)中已知观测混杂因素的回归系数ߚଵ、ߚଶ分别设为-1、1,未知观测混杂因素的回归系数ߚ௨在模拟过程中分别设定为0、1、2,分别表示为无、弱、强三种-46- 基于多水平模型的工具变量方法研究及应用类型的未知观测混杂因素。ߚ௫为处理/暴露因素的回归系数,其OR值表示为exp(ߚ௫),在本研究中回归系数设为2,OR值为exp(2),为7.3891,此为本研究的金标准,后面进行模型比较时以回归系数为准。本研究共模拟9种情况的数据集,详见表2-8。表2-8处理/暴露因素和结局变量为分类变量参数取值数据集未知观测混杂因素强度(ߚ௨)工具变量强度(ߙ௭)模拟101模拟202模拟303模拟411模拟512模拟613模拟721模拟822模拟923在本研究中,所模拟的数据为层次结构特征数据,在模拟数据中所设的第二水平共有200层,每层中含有10个个体数据,每次模拟的数据集大小为2000例,在模拟结局变量Y时,每层的常数项ߚ(j=1,2,3,……,200)服从正态分布,如下式:ߚ~N(0,0.5)(2-22)在数据模拟时,对每一种情况下的模拟数据进行模拟500次,模拟完成的数据集如表2-9所示。表2-9处理/暴露因素和结局变量为连续型变量数据模拟示例sitec1c2uzpixxpiy12.894710.074210.187100.0560010.934611.009310.390810.9528112.064411.644210.099000.7586111.521601.163800.039700.4322011.119112.224000.055000.9627111.859811.196200.072000.6115111.786200.595510.132200.1579012.395800.332000.038100.0568011.679910.884100.112700.5024011.318010.937800.147300.6175122.282000.941310.061600.2182022.159102.059000.008900.7469120.702400.562300.146210.8243121.854801.302000.025200.4682021.823100.986210.090400.32581-47- 基于多水平模型的工具变量方法研究及应用sitec1c2uzpixxpiy21.441410.762000.154000.5513122.044011.352600.052300.6867122.209601.562900.013800.5099122.058501.215300.022400.3765021.430611.015100.125000.6733132.196300.224210.127810.3833131.196110.748610.390600.6396031.538110.994410.262600.6733132.129300.492710.107000.1334031.850401.458210.056900.5838032.081600.235410.139800.0880131.299810.510810.423010.8802132.314111.436510.095300.69660………………………………………………注:表中site表示水平2单位;pix为X的概率分布,即P[X=1│C1,C2,U,Z];pi为Y的概率分布,P[Y=1│C1,C2,U,X]。(二)模型在模拟数据中的应用同处理/暴露因素和结局变量为连续型变量类似,在对模拟数据分析时,首先将数据中的未知观测混杂因素U去除,并不纳入后续的模型分析中,故所分析的模拟数据中就存在了未知观测混杂因素。1.普通多水平logistic回归模型在普通多水平logistic回归模型中纳入已知观测混杂因素C1、C2、处理/暴露因素X和结局变量Y,公式如下:ܻ=݈݃݅ݑ+ߚ=൯ܲ൫ݐ+ߚ௫ܺ+ߚଵܥ1+ߚଶܥ2(2-23)ݑ~ܰ(0,ߪଶ)௨బ式中i表示水平1单位,i=1,2,3,……,10;式中j为水平2单位,j=1,2,……,200(下同);ߚ௫为处理/暴露因素的回归系数,可根据ORX=exp(ߚ௫)计算出处理/暴露因素的OR值。2.两阶段logistic回归多水平工具变量模型第一步:利用已观测混杂因素C1、C2、工具变量Z同处理/暴露因素X做logistic回归,以得到处理/暴露因素X的概率预测值P(X)。ଵP(X)=PX=1หC1,C2,Z൧=(2-24)ଵା୶୮[ି(ఈబାఈభେଵାఈమେଶାఈೋ]第二步:利用已观测混杂因素C1、C2、处理/暴露因素X的概率预测值P(X)和结局变量Y做多水平logistic回归模型,其公式如下:-48- 基于多水平模型的工具变量方法研究及应用ܻ=݈݃݅ݑ+ߚ=൯ܲ൫ݐ+ߚ௫ܲ(ܺ)+ߚଵܥ1+ߚଶܥ2(2-25)ݑ~ܰ(0,ߪଶ)௨బ式中ߚ௫为处理/暴露因素的回归系数。3.线性回归+logistic回归多水平工具变量模型第一步:利用工具变量Z、已观测混杂因素C1、C2和处理/暴露因素X做多元线性回归,获得重新估计的处理/暴露因素X,如2-8式所示。第二步:利用已观测混杂因素C1、C2、重新估计得到的处理/暴露因素X预测值X和结局变量Y做多水平logistic回归模型,其公式如下:ܻ=݈݃݅ݑ+ߚ=൯ܲ൫ݐ+ߚ௫X+ߚଵܥ1+ߚଶܥ2(2-26)ݑ~ܰ(0,ߪଶ)௨బ式中ߚ௫为处理/暴露因素的回归系数。4.自助法与模型结合在此部分的模拟研究中,与处理/暴露因素和结局变量同为连续型变量部分类似,采用分层个例重复抽样法,分别在200层中(即在水平2单位内)进行有放回的重复抽样,每层中有放回地重抽出10例样本,每次所得复样本为2000例,同原始样本量一致,自助法共反复抽取500次,再分别利用普通多水平logistic回归模型、两阶段logistic回归多水平工具变量模型和线性回归+logistic回归多水平工具变量模型三种模型对500个复样本进行分析,用所得的500个样本参数重新计算出最终的效应值、标准误和置信区间,自助法的具体实施过程同自助法普通多水平线性回归模型类似。(三)结果比较本研究针对处理/暴露因素和结局变量均为分类变量的模型拟合共构建和纳入了6个模型。1.模型准确度的比较表2-10展示了这6种模型的准确度。当无未知观测混杂因素(β୳=0)存在时,从表中可以发现多水平logistic回归模型参数估计的绝对偏倚最小,其他五种模型均存在不同程度的偏倚,其中偏倚最大的是自助法线性回归+logistic回归多水平工具变量模型,模型参数估计的置信区间也以多水平logistic回归模型和自助法多水平logistic回归模型较窄,对于四种多水平工具变量模型,当使用弱工具变量(ߙ௭=1)时,模型参数估计的置信区间很宽,说明模型参数估计的准确度较差,随着工具变量强度的不断增加,模型参数估计准确度也在不断增加。-49- 基于多水平模型的工具变量方法研究及应用在模拟数据中加入不同强度的未知观测混杂因素后,多水平logistic回归模型和自助法多水平logistic回归模型随着未知观测混杂因素强度的增加,虽然参数估计的置信区间仍然是最窄的,但其所得的参数估计值离金标准越来越远,说明这两种模型在存在未知观测混杂因素的情况下,模型参数估计的准确度较差。当数据中存在弱未知观测混杂时,我们发现自助法线性回归+logistic回归多水平工具变量模型的参数估计绝对偏倚较小,均在0.1以下;自助法两阶段logistic回归多水平工具变量模型的参数估计绝对偏倚也较小,尤其是在利用中度和强度工具变量进行估计时,其绝对偏倚也未超过0.1。在强未知观测混杂因素存在的情况下,同样是自助法两阶段logistic回归多水平工具变量模型和自助法线性回归+logistic回归多水平工具变量模型的参数估计绝对偏倚较小。其模型参数估计的置信区间随着工具变量的强度增加区间宽度也会逐渐变窄,说明工具变量强度越强,多水平工具变量模型参数估计的准确度越好。2.模型精确度的比较表2-11展示了6中模型的精确度。当数据集中不存在未知观测混杂因素时,我们发现普通多水平logistic回归模型、两阶段logistic回归多水平工具变量模型和线性回归+logistic回归多水平工具变量模型在模型参数估计的覆盖率方面较为稳定,当使用自助法后,模型参数估计的覆盖率有着一定的变化,尤其是在使用强工具变量时,覆盖率下降很多,这可能由于参数估计的置信区间变窄所造成。多次模拟所得参数的标准误以普通多水平logistic回归模型最小,从而说明当数据中无未知观测混杂因素时,普通logistic回归模型较其他模型在精确度方面更加稳定。在模拟数据集中加入不同强度的未知观测混杂因素后,普通多水平logistic回归模型的参数估计覆盖率均在65%以下,尤其是当存在强未知观测混杂因素情况下,该模型的参数估计覆盖率未超过5%。在六个模型中,参数估计覆盖率表现最好的是自助法线性回归+logistic回归多水平工具变量模型,均在85%以上,其次为自助法两阶段logistic回归多水平工具变量模型。多水平工具变量模型多次模拟所得估计参数的标准误均随着工具变量强度的不断增加而减少,说明强工具变量可以得到更为精确的结果。图2-9至图2-11直观地展示了这六种模型在不同强度未知观测混杂因素存在的情况下,利用不同强度工具变量后模型参数估计情况。当数据集中无未知观测混杂因素时,普通多水平logistic回归模型的参数点估计最接近金标准,并且区间最窄,但当数据集中存在了未知观测混杂因素时,自助法两阶段logistic回归多水平工具变量模型和自助法线性回归+logistic回归多水平工具变量模型两种模型的参数估计较为准确,但是当所采用的工具变量较弱时,该两种模型参数估计的置信区间较宽,说明模型估计的准确性较差,采用较强的工具变量进行模型估计时,该两种模型所得参数-50- 基于多水平模型的工具变量方法研究及应用的置信区间将逐渐变窄,说明在自助法多水平工具变量模型中采用强工具变量进行运算时,模型所得的结果会更加准确。综上所述,对于处理/暴露因素和结局变量均为分类变量时,模型方法的选择需事先仔细考虑数据中是否存在未知观测混杂因素,当数据无未知观测混杂因素时,需要考虑应用普通多水平logistic回归模型;当存在未知观测混杂时,需考虑采用自助法两阶段logistic回归多水平工具变量模型和自助法线性回归+logistic回归多水平工具变量模型两种模型。对于两阶段logistic回归多水平工具变量模型和线性回归+logistic回归多水平工具变量模型在各种数据条件的情况下并未表现出很好的效果。-51- 61区间下限(%)90.869.246.297.094.494.896.090.685.8覆盖率区间上限)Boot-MLSPS୶β-52-(绝对偏倚1.75890.56070.29101.44970.50120.25901.45210.50730.2680(%)SE区间下限95.486.271.695.293.293.694.290.284.4覆盖率区间上限)୶)β(=21.25190.47250.26901.20000.45280.25461.21410.45440.2579୶绝对偏倚ߚ(%)SE基于多水平模型的工具变量方法研究及应用区间下限MBoot-M2SPS61.853.237.689.091.091.619.49.410.4覆盖率)区间上限୶β(0.22250.17740.15070.23830.16490.14910.19830.16310.1372绝对偏倚(%)SE94.688.684.496.893.889.294.887.865.6限区间下覆盖率)୶β区间上限(1.49880.49250.25621.24300.44450.22691.25600.44730.2338绝对偏倚(%)SE95.896.695.694.290.482.894.082.660.2覆盖率区间下限结局变量和混杂因素均为分类变量时各方法精确度的比较)୶β结局变量和混杂因素均为分类变量时各方法准确度的比较(2-11(区间上限0.40970.23531.05360.40340.22371.07430.40260.2254表2-10表(%)SE绝对偏倚95.895.661.250.249.43.60.80.0覆盖率区间下限Boot-MLRMLRMM2SPSMLSPS)୶β(SE区间上限Boot-MLSPSBoot-M2SPSBoot-MLRMMLRMM2SPSMLSPS)௭(ߙ10.185694.41.066520.146830.126010.209920.144130.128210.174720.144030.1198绝对偏倚工具变量强度)ߙ௭10.01431.652022.3767-0.1313-0.31410.00461.71364.051632.29560.03071.17912.88220.00901.75800.5372-0.42801-0.30461.30012.26000.04631.57312.51945.50240.42211.40743.43692-0.29601.39682.0906-0.4624-0.58930.33851.94050.26111.73832.78400.32041.99833-0.26431.48063.66452.73652.0113-0.29280.88432.53002.64260.31502.0345-0.2628-0.87091-0.70800.92840.2078-0.1430-0.19840.90602.69731.9908-0.23811.30332.22040.32651.45843.19472.59554.34534.55852-0.67641.0368-0.08051.5975-0.19581.32952.27891.6556-0.4634-0.56930.72561.69623.75490.32941.85022.8085-0.09331.49430.8598-0.32716.04682.24153-0.60731.15253.6424-0.04221.68681.6103-0.43010.75302.38682.31910.54902.01893.0791-0.08551.07582.7531-0.3504-0.90512.2288-0.2579-0.5365-0.37340.73752.51561.6329-0.38681.16162.06474.20430.01591.10152.9302-0.02471.51022.44044.0208-0.51241.1880-0.36751.16502.1001-0.55081.0705-0.0535-0.88664.7790.02191.54042.50341.7871)-0.43721.30971.8279-0.23640.92262.60461.8158-0.2550-0.5223-0.17540.90812.7411-0.19671.34042.26624.0123-0.17651.34422.3027-0.1321-0.92544.661工具变量强度(୳(β012)β୳012混杂因素强度混杂因素强度( 基于多水平模型的工具变量方法研究及应用图2-9分类数据无未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=0)图2-10分类数据弱未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=1)-53- 基于多水平模型的工具变量方法研究及应用图2-11分类数据强未知观测混杂因素情况下各模型参数点估计和置信区间比较(β୳=2)四、讨论在本研究的模型构建和模拟研究验证阶段分别采用了普通多水平回归模型和构建出两种多水平工具变量模型,数据分析时将数据主要分为连续型变量数据(其中分为处理/暴露因素存在层级效应和无层级效应)和分类数据两种类型,在数据模拟过程中设立了不同强度的未知观测混杂因素和不同强度的工具变量,在这样数据情况下分别用各种分析模型进行拟合,并根据初设的数据金标准对各种模型进行评价比较。(一)普通多水平回归模型模拟研究结果在数据模拟研究分析过程中,本研究发现当数据中无未知观测混杂因素存在的情况下,普通多水平线性回归模型和普通多水平logistic回归模型分别在处理/暴露因素和结局变量为连续型变量和分类变量的情况下得到很好的分析效果,此结果在参数的准确度和精确度方面都较好,这说明经典的多水平回归模型在数据中影响因素较全的情况下可以分析得到准确的结果。当研究中存在了未知观测混杂因素时,传统的多水平回归模型所得的结果将会产生不同程度的偏倚,在未知观测混杂因素强度增加时模型估计所得的参数偏离真值越远,在处理/暴露因素和结局变量为连续型变量时,普通多水平线性回归模型离金标准最大的为2.8219,在处理/暴露因素和结局变量为分-54- 基于多水平模型的工具变量方法研究及应用类变量时,普通多水平logistic回归模型离金标准最大的为-0.7080。结果说明在实际研究中,当我们无法收集到较全的影响因素(即遗漏了较为重要的影响因素)时,传统的多水平回归模型将不再适用,否则将会获得错误的分析结果。(二)普通多水平工具变量模型模拟研究结果根据参数模拟的结果,处理/暴露因素和结局变量为连续型变量时发现多水平工具变量在参数估计的准确度和精确度方面均有较普通多水平较好,普通多水平工具变量模型在置信区间覆盖率方面逊色于自助法多水平工具变量模型,其原因在于自助法多水平所得结果较为保守,其置信区间较宽。在处理/暴露因素和结局变量为连续型变量且均存在层次效应的数据中构建的四种多水平工具变量在结果的准确度和精确度方面均有较好的表现,较普通多水平线性回归模型和自助法两阶段最小二乘多水平工具变量模型更为准确和稳定。对于处理/暴露因素和结局变量为分类变量的数据,当工具变量强度较弱的时候,这两种模型的参数点估计离真值较远,但当工具变量强度为中度或强度时,参数点估计明显准确很多,普通多水平回归模型在此部分并未获得较为理想的结果。(三)自助法引入模型后模拟研究效果在模拟研究中,根据模拟数据的数据特征,我们均采取了有放回的分层个例重复抽样法,共重复抽出500个复样本,再用各个模型对这些复样本进行拟合,从而每个模型将会获得500个参数,这500个参数的均数即为原始样本的参数估计的点估计值,其标准差即为原始样本参数估计的标准误。本研究希望通过将自助法引入到多水平工具变量模型中,增加模型的准确度和精确度。模拟研究的结果显示自助法多水平工具变量在连续型和分类数据中都有较好的表现。在处理/暴露因素和结局变量均为连续型变量模拟数据中,我们发现自助法两阶段最小二乘多水平工具变量模型和自助法两阶段残差纳入多水平工具变量模型在模型参数估计的点估计和区间的覆盖率方面都有着较好的表现。两种模型中自助法两阶段最小二乘多水平工具变量模型的参数估计较为稳定,受到工具变量强度的影响并不十分明显,但是该模型也较为保守,所估计参数的置信区间在这六种模型中也是最宽的。当我们考虑资料中处理/暴露因素也存在层级效应时,自助法两阶段多水平回归工具变量模型在模型的点估计、置信区间宽度、标准误和置信区间覆盖率方面均有较好的表现,其他三个多水平工具变量模型所估计的参数点估计也很接近金标准。从此可见,对于处理/暴露因素和结局变量均为连续型变量的实例分析中的模型选择,首先我们需要确定是否数据中有遗漏的重要混杂因素,如果没有,我们可以用普通多水-55- 基于多水平模型的工具变量方法研究及应用平回归模型进行分析,当存在重要的未知观测混杂因素时,我们首先需要找到一个合适的工具变量,然后检测处理/暴露因素是否在水平单位上存在异质性,如果不存在可以利用两阶段最小二乘多水平工具变量模型和两阶段残差纳入多水平工具变量模型(普通+自助法)进行模型拟合,增加模型结果的稳定性;如果在水平单位上存在异质性可以采用两阶段多水平回归工具变量模型和两阶段多水平回归残差纳入工具变量模型(普通+自助法)进行分析。在处理/暴露因素和结局变量均为分类变量模拟数据中,我们发现自助法两阶段logistic回归多水平工具变量模型和自助法线性回归+logistic回归多水平工具变量模型两种模型所获得的结果较为相似,在无未知观测混杂因素存在时,这两种自助法多水平工具变量模型所获得的结果均有一定程度的偏倚,但当未知观测混杂因素存在时,这两种自助法多水平工具变量模型均能获得较为准确的参数估计,但是这两种模型参数估计的置信区间均受到工具变量强度影响较大,弱工具变量时,其所获得的置信区间很宽,当随着工具变量强度增加时,其区间也有明显缩短,但参数点估计受工具变量强度的影响较小。由此可见,当我们对处理/暴露因素和结局变量均为分类变量的实际数据分析时模型选择需要考虑未知观测混杂因素是否存在,无未知观测混杂因素存在时,普通多水平logistic回归模型较为适合,但当未知观测混杂因素存在时,可以同时引用上述两种自助法多水平工具变量模型进行分析,从而增加结果的可靠性,若需要获得更为准确的结果,即需要选择强度较强的工具变量进行运算。(四)局限性在模拟研究中所考虑的数据集情况均为最理想的数据环境,数据的分布和数据间的关系均在最理想的环境下设定的,在此条件下所构建的模型将更加合适,但应用于实际数据集分析中时,模型的可靠性我们将继续进行探索。在数据分析时我们需要考虑到各种情况,故我们在设置模拟数据参数的时候通过设置不同强度的未知观测混杂和不同强度的工具变量,从而全面的评估多水平工具变量模型的适用范围和模型的稳定性。在考虑工具变量的强弱时,本研究首先基于的前提是工具变量同结局变量间是相互独立的(即无相关关系),其强弱是由工具变量同处理/暴露因素间回归系数的大小来确定,而在实际分析中,工具变量的强弱不仅通过其与处理/暴露因素间关系的大小来确定,还要根据其与结局变量间是否存在相关关系而确定,与结局变量的相关性对最后参数估计值得影响方面将在后续的研究中进行重点探索。工具变量模型应用于实际数据中,其重点和难点就是寻找一个合理的强工具变量,对于工具变量强度有文献是根据第一步工具变量和处理/暴露因素间分析时所得到的-56- 基于多水平模型的工具变量方法研究及应用57F值统计量确定,但此统计量会受到样本量的影响,也有研究用相关系数和OR值58,59来确定工具变量的强度。在模拟研究中我们通过设定不同强度的工具变量来评价模型,根据模拟研究的结果,我们针对不同数据类型和不同强度的工具变量给出了最佳的模型来对结果进行分析,在弱工具变量存在时,模型研究所展示的问题主要存在于其参数估计的置信区间较宽,但点估计的偏差并不大,在实际应用时其模型获得的结果可能较为保守,如何增进弱工具变量估计的准确度将进一步进行研究。在多水平工具变量模型的构建过程中,我们主要考虑了在第二步的时候融入多水平回归模型思想,在第一步的时候应用的是普通回归模型进行估计,在处理/暴露因素和结局变量为连续型变量时我们在第一步运算时考虑了在处理/暴露因素上存在水平因素的影响,故在第一步运算时采用多水平线性回归对处理/暴露因素进行重估计和残差的计算,然后在将进行第二步的多水平回归运算,通过两步均为多水平模型的运算,针对处理/暴露因素和结局变量均存在其余水平因素影响情况下的数据可以获得更为准确的参数估计,针对处理/暴露因素和结局变量为分类变量的数据,模型尚在构建中,将考虑到二分类变量数据中变量重新估计如何更加合理,如何克服弱工具变量造成参数置信区间较宽的弊端,将进行后续研究。目前有研究为了更好的控制已知观测混杂因素和未知观测混杂因素,研究者先用倾向性评分法对模型进行倾向性评分匹配,从而控制了已知观测混杂因素,再利用工具变量方法来控制未知观测混杂因素,从而获得更加准确的结果,在多水平工具变量中也可以引入这一思想,在后续的研究中考虑加入倾向性匹配思想,从而更加全面的对各种混杂因素进行控制。-57- 基于多水平模型的工具变量方法研究及应用第三部分实例应用本研究实例分析部分应用的是第五次全国卫生服务调查数据(上海),其中所针对的人群为60岁以上的老人,分别关心每周运动时间对其身体健康状况的影响和是否吸烟与其是否患有高血压间的关系。一、概况卫生服务调查问卷中共设计了6个部分的调查表:1.家庭一般情况调查表,此部分主要是关心家庭的人口、收入、消费等一般情况;2.家庭成员个人情况调查表,此部分主要调查家庭中个人的基本情况、生活习惯和健康情况;3.调查前两周内病伤情况调查表,此部分调查家庭成员两周内病伤和门诊就诊情况;4.调查前一年内住院情况调查表,此部分主要关心个人一年内住院情况;5.5岁以下下儿童调查表,反映了了5岁以下小孩的基本情况和疫苗接种情况;6.15-64岁妇女调查表,反映了妇女生育情况。通过这6部分的调查表进行调查,可以全面的反映居民的的一般身体健康状况。在调查过程中采用的是分层随机抽样,决定了卫生服务调查数据为层次结构特征,其数据结构为:省市—区县—街道—居委会—家庭—个人。其具体格式如图3-1。图3-1卫生服务调查数据(上海)的数据结构构图本次的上海市卫生服务调查覆盖了全市18个区,共200个居委,每个居委各调查60户家庭,共调查12002户,共计31532人,其中60岁以上的男性共5415人,已婚并配偶健在的有4781人,60岁以上已婚并配偶健在的女性4285人。卫生服务调查的主要目的是进行一次人口健康状况的普查,从而全面了解和掌握居民的健康状况,但此调查并非为某一项专病的专项调查,如果利用这样的数据进行某一项专项的分析,如分析每周运动时间对健康状况的影响,数据中所提供的信息信(即所存在的协变量)可能不足已全面的反映这样的问题,如果我们采用传统的回归模型进行分析时,由于存在遗漏的混杂因素,会对最终的结果造成较为严重的偏倚,从而-58- 基于多水平模型的工具变量方法研究及应用影响到结果的准确性和真实性。二、应用实例一:60岁以上老人每周运动时间和身体健康状况的关系(一)人群选择和变量选择在此项研究中所选的人群为60岁以上已婚并配偶健在的老人,分析过程中男性和女性将分别进行分析,其中男性共4781人,女性4285人。处理/暴露因素:定为每周运动时间,调查问卷中有平均每周锻炼次数变量“近6个月,您平均每周体育锻炼几次”,对变量进行转化,如果从不锻炼,此项设为0,如果锻炼不到1次,设为1,锻炼1-2次,设为2,锻炼3-5次设为4,锻炼6次及以上则设为6,调查问卷中还平均每次锻炼时间(分钟),将次变量转化为小时,最终每周运动时间为上述两项(每周锻炼次数和平均每次锻炼时间)的乘积。结局变量:定为自评身体健康评分,在调查问卷中自评健康情况采用的是欧洲五维健康量表(Europeanqualityoflife5-dimensions,EQ-5D)评分,其中包括了心动、自我照顾、平常活动(工作、读书或做家务)、疼痛不舒服和焦虑和抑郁五个方面的调查,每项调查中各包含三个选项,分别为没有问题、有些问题和很有问题三个等级,除了这五方面的调查外,该量表中还包含了一项自评健康的评分,其分值为0-100分,分值越大反映的健康状况越好,分值越小反映的健康状况越差。在此项研究中,我们采用EQ-5D自评健康评分为最终的结局变量。已知观测混杂因素:共从调查问卷中选出10个已知观测混杂因素,年龄、户口性质(农村和城市)、家庭一年总收入、bodymassindex(BMI,BMI小于18.5定为过轻,BMI在18.5-24定为正常,BMI在24-28定为超重,BMI等于28及以上定为肥胖)、吸烟(每天吸、非每天吸和不吸)、饮酒(饮酒和不饮)、文化程度(没上过学和小学设为1、初中设为2,高中和技工学校设为3,中专和大专设为4,本科及以上设为5)、慢性病情况(所患有慢性病的种数)、两周有误病伤情况和一年内是否有过住院。以上十个已知观测混杂因素都可能会影响到所调查居民的自评健康状况。工具变量选择:此研究所选择的工具变量为其爱人每周运动的次数,选择的理由是爱人每周运动的次数会影响自己每周运动时间,但是与自身的健康评分并无关系。此研究由于所包含的变量较少,很难较为全面地对自身健康评分的影响因素进行探索,如本人心情状况、疾病轻重程度、本人体质等因素,这些因素为未知观测混杂因素,故对于这样的问题,普通回归模型可能无法准确的分析影响因素的影响系数。-59- 基于多水平模型的工具变量方法研究及应用(二)人群基本情况本研究对研究人员的基本情况进行分析和展示,详情见表3-1。通过分析我们发现在男性中自评健康评分为77.67、女性为76.88,两者相近,每周的运动时间也是相近的,分别为2.55小时和2.34小时。在10个已观测混杂中,男性和女性在吸烟和饮酒这两项生活习惯中分布有很大的不同,男性每天吸烟共1684人,占35.22%,在女性中每天吸烟共27人,占0.63%,男性中饮酒共1651人,占34.66%,女性中饮酒仅122人,占2.87%,在其他八个已知观测混杂因素中,男性和女性的基本情况较为一致。表3-1资料中各变量的基本情况男性女性变量(N=4781)(N=4285)结局变量EQ-5D(分)(Mean±SD)77.67±13.1176.88±13.28处理/暴露因素每周运动时间(小时)(Mean±SD)2.55±3.352.34±3.13已知观测混杂因素年龄(Mean±SD)64.41±7.5768.28±7.05家庭年收入(万元)(Mean±SD)7.20±5.427.14±5.61慢性病种数(Mean±SD)0.91±0.930.89±0.92城市户口3920(81.99)3398(79.30)文化程度没上过学和小学1485(31.06)1973(46.04)初中1535(32.11)1255(29.29)高中和技工学校974(20.37)703(16.41)中专和大专390(8.16)202(4.71)本科及以上397(8.30)152(3.55)BMI过轻(<18.5)267(5.58)319(7.44)正常(18.5-24)2682(56.10)2361(55.10)超重(24-28)1571(32.86)1287(30.04)-60- 基于多水平模型的工具变量方法研究及应用男性女性变量(N=4781)(N=4285)肥胖(≥28)261(5.46)318(7.42)吸烟每天吸1684(35.22)27(0.63)非每天吸188(3.93)49(1.14)不吸2909(60.85)4209(98.23)饮酒1651(34.66)122(2.87)两周有病伤2705(56.79)2363(55.35)一年内有住院461(9.84)336(8.02)(三)工具变量的强度分析在研究分析中,我们将男性和女性分开进行分析,在分析过程中分别以其配偶的每周运动次数作为工具变量,当进行男性老年人分析时,利用工具变量和其他已知混杂因素同处理/暴露因素男性每周运动时间进行回归分析时,得出工具变量配偶每周运动次数的回归系数是0.4829,F值为765.69。当进行女性老年人分析时,也利用工具变量和其他已知混杂因素同处理/暴露因素女性每周运动时间进行回归分析时,得出工具变量配偶每周运动次数的回归系数是0.4925,F值为892.79。根据这两个数据集工具变量分析的结果,我们发现工具变量的F值比较大,均在750以上,虽然回归系数较小,仅为0.49左右,但总体看来,本研究所选的工具变量强度较好。(四)实例分析结果1.男性老年人每周运动时间对身体健康状况的影响在进行实例分析时,数据分为两个水平,水平一单位定为所调查的居民个体,水平二单位定为居委会。我们首先建立空模型,检测数据在第二水平(居委会)层面的残差方差为25.2443,标准误为3.2082,Z值为7.87,p<0.0001,说明此研究的结局变量男性自评健康评分在水平间存在异质性,故应该采用多水平模型进行分析。我们应用普通多水平回归模型和八个多水平工具变量模型分别进行分析,分析结果如表3-2所示。我们发现普通多水平线性模型的效应值为0.4185(0.4090-0.42810),但其他八个多水平工具变量的点估计分别为0.6433、0.6878、0.7216、0.7246、0.6447、0.6903、0.7258、0.7282,所有多水平工具变量模型所得结果均在0.7左右,但普通多水平回归模型所得的结果仅在0.4,即使所有结果均有统计学意义,提示每周运动时间越长,自身健康状况越好,但是多水平工具变量所得的回归系数与普通多水平工具变量所得-61- 基于多水平模型的工具变量方法研究及应用的点估计相差约0.3。我们又对处理/暴露因素进行空模型的建立,检测出处理/暴露因素在第二水平(居委会)层面的残差方差为1.5721,标准误为0.2024,Z值为7.77,p<0.0001,说明处理/暴露因素在第二水平存在异质性,所以我们在工具变量模型的第一步中应该应用多水平回归模型,故两阶段多水平回归工具变量模型、两阶段多水平回归残差纳入工具变量模型、自助法两阶段多水平回归工具变量模型和自助法两阶段多水平回归残差纳入工具变量模型更加适合本数据。这四个多水平工具变量模型所得的结果较为一致,均约为0.72。此工具变量的强度较好,由于普通多水平工具变量模型和加入自助法后的多水平工具变量模型所获得的结果相近,但普通多水平工具变量所得结果的置信区间较窄,自助法所获得结果的置信区间比较,故我们可以采用两阶段多水平回归工具变量模型和两阶段多水平回归残差纳入工具变量模型的结果来解释本研究。表3-2男性老年人每周运动时间对身体健康状况影响的结果分析方法效应值点估计置信区间普通多水平线性模型0.41850.4090-0.4281两阶段最小二乘多水平工具变量模型0.64330.6179-0.6686两阶段残差纳入多水平工具变量模型0.68780.6625-0.7132两阶段多水平回归工具变量模型0.72160.6931-0.7501两阶段多水平回归残差纳入工具变量模型0.72460.6962-0.7531自助法两阶段最小二乘多水平工具变量模型0.64470.3659-0.9234自助法两阶段残差纳入多水平工具变量模型0.69030.4117-0.9690自助法两阶段多水平回归工具变量模型0.72580.4098-1.0419自助法两阶段多水平回归残差纳入工具变量模型0.72820.4121-1.04422.女性老年人每周运动时间对身体健康状况的影响在进行实例分析时,我们参照了分析男性老年人的分析过程,首先对结局变量女性身体健康状况评分建立空模型,检测数据在第二水平(居委会)层面的残差方差为25.5211,标准误为3.2924,Z值为7.75,p<0.0001,说明该结局变量在第二水平间存在一定的异质性,建议采用多水平模型进行分析。我们随后对处理/暴露因素女性每周运动时间进行空模型的建立,检测出数据在第二水平(居委会)层面的残差方差为1.5300,标准误为0.1939,Z值为7.89,p<0.0001,说明处理/暴露因素在第二水平间也存在一定的异质性,同样建议在第一步回归时应用多水平模型进行分析。我们同时应用普通多水平线性回归模型和八种多水平工具变量模型对数据进行分析,结果同男性老年人分析时所得结果类似,八种工具变量所得结果较为一致,与-62- 基于多水平模型的工具变量方法研究及应用普通多水平线性模型所得结果有一定差异,普通多水平模型所得效应值为0.4897,多水平工具变量模型所得结果分别为0.3251、0.3929、0.3705、0.3761、0.3373、0.4063、0.3857和0.3901。两者结果相差近0.1左右,此结果并无分析男性老年人分析时普通多水平线性模型同多水平工具变量模型间的差异大,说明在分析女性时,已有的变量可以在较大程度上解释自变量与结局变量间的关系。在工具变量第一步回归时,我们检测出处理/暴露因素在第二水平上存在异质性,说明有水平效应的存在,故我们认为第一步也应该采用多水平回归模型,故女性每周运动时间与自身健康状况的关系为正相关,两阶段多水平回归工具变量模型、两阶段多水平回归残差纳入工具变量模型、自助法两阶段多水平回归工具变量模型和自助法两阶段多水平回归残差纳入工具变量模型所得系数均为0.38左右,此反映了女性每周运动时间和自身健康状况间较为真实准确的关系。表3-3女性老年人每周运动时间对身体健康状况影响的结果分析方法效应值点估计置信区间普通多水平线性模型0.48970.4787-0.5006两阶段最小二乘多水平工具变量模型0.32510.2995-0.3507两阶段残差纳入多水平工具变量模型0.39290.3672-0.4185两阶段多水平回归工具变量模型0.37050.3413-0.3997两阶段多水平回归残差纳入工具变量模型0.37610.3469-0.4053自助法两阶段最小二乘多水平工具变量模型0.33730.0604-0.6141自助法两阶段残差纳入多水平工具变量模型0.40630.1305-0.6821自助法两阶段多水平回归工具变量模型0.38570.0694-0.7019自助法两阶段多水平回归残差纳入工具变量模型0.39010.0740-0.7063三、应用实例二:60岁以上男性吸烟对其患高血压的影响(一)人群选择和变量选择在此项研究中我们所选择的人群为60岁以上的男性老年人,从总体人群中共选择5415人。处理/暴露因素:根据此项研究的研究目的,处理/暴露因素我们定为男性老年人是否吸烟,在调查问卷表中有“您现在的吸烟状况”这项调查条目,其中包括了三个选项:每天吸、非每天吸和不吸。在分析过程中我们将变量进行处理,将每天吸和非每天吸的居民归为吸烟组,并将不吸的居民归为非吸烟组。故在分析时,原来的三分-63- 基于多水平模型的工具变量方法研究及应用类变量就转化为现在的二分类变量。结局变量:在此项研究中,结局变量我们定为居民是否患有高血压,在调查问卷中有“您是否被医生确诊患有高血压病”这项调查项目,此项目有两个选项:是或否,故结局变量为二分类变量。已知观测混杂因素:由于问卷调查中包含的疾病危险因素相关的变量较少,我们仅从调查表中找出5个已观测混杂因素,分别为:1.户口性质(城市和农村);2.家庭一年总收入;3.BMI系数(其分类同分析连续型变量时类似,共分为四类,过轻、正常、超重和肥胖);4.饮酒(饮酒和不饮);5.文化程度(其分类同分析连续型变量时类似,共分为5类)。这五个已知观测混杂因素中属于连续型变量的为家庭一年总收入,属于二分类变量的为户口性质和饮酒,属于有序多分类变量的为BMI和文化程度。(二)人群基本情况在全部60岁以上男性老年人群(共有5415人)中,吸烟的人群(每天吸烟和非每天吸烟)共有2093人,占38.65%,总体人群中患有高血压的人群共有2637人,占48.70%。本研究共纳入五个已知观测混杂因素混杂因素,这些已知观测混杂因素在两组间的比较中均存在统计学的差异,如饮酒中吸烟的人群更加倾向于饮酒占了51.58%。两组体型的比较发现,吸烟组中提醒过轻的和正常的较非吸烟组所占的比例较高,非吸烟组中超重和肥胖较吸烟组所占的比例高(超重:32.84%vs30.24%;肥胖:5.78%vs4.97%)。所调查的人群中吸烟组较非吸烟组的受教育程度低,年收入也较低(6.63万vs7.22万),吸烟人群中城市户口所占的比例也较低一些。详见表3-4。表3-4资料中两组人群基线的基本情况已知观测混杂因素吸烟组非吸烟组p值总人数20933322家庭一年总收入6.63±5.067.22±5.76<0.0001户口性质城市1617(77.26)2782(83.74)农村476(22.74)540(16.26)<0.0001饮酒饮酒1077(51.58)755(22.83)不饮1011(48.42)2552(77.17)<0.0001BMI过轻(<18.5)129(6.16)202(6.08)-64- 基于多水平模型的工具变量方法研究及应用已知观测混杂因素吸烟组非吸烟组p值正常(18.5-24)1227(58.62)1837(55.30)超重(24-28)633(30.24)1091(32.84)肥胖(≥28)104(4.97)192(5.78)0.0184文化程度没上过学和小学772(36.88)1037(31.22)初中787(37.60)932(28.06)高中和技工学校350(16.72)695(20.92)中专和大专109(5.21)308(9.27)本科及以上75(3.58)350(10.54)<0.0001注:在分析过程中连续型变量采用t检验,二分类变量采用߯ଶ检验,有序多分类变量采用非参数检验。(三)工具变量强度在研究分析中,我们利用家人是否吸烟作为工具变量进行分析,首先我们利用纳入的五个已知观测混杂因素和工具变量家人是否有人吸烟同处理/暴露因素自身是否吸烟做logistic回归和多元线性回归,在logistic回归中工具变量的OR值为1.81,在多元线性回归分析中,工具变量的回归系数仅为0.0087,但F值为49.14,从以上的结果我们可以发现,研究分析中我们所选取的工具变量强度并不强,仅为中等强度工具变量。由于研究所调查内容的局限性,我们并没有找到强工具变量的存在,故在后续分析中我们将利用此工具变量。(四)实例分析结果此项研究中我们共设有了两个水平,水平1单位代表的是个体(居民个人),水平2单位代表的是居委会,即生活在同一居委会的居民为一个群体。在进行实例模型拟合分析前,我们首先对结局变量是否还有高血压进行了空模型的构建,检测结局变量在水平2单位上是否存在异质性,空模型拟合的结果为ICC值为1.7741,标准误为0.3269,t值为5.43,p<0.0001,说明数据存在相当程度的组间异质性,即结局变量在水平2单位上存在异质性,故需要采用多水平回归模型。在进行分析时,我们分别采用了普通多水平logistic回归模型、两阶段logistic回归多水平工具变量模型、线性回归+logistic回归多水平工具变量模型、自助法两阶段logistic回归多水平工具变量模型和自助法线性回归+logistic回归多水平工具变量模型五个模型分别对数据集进行分析。普通的回归模型所得结果的回归系数为-0.3057,-65- 基于多水平模型的工具变量方法研究及应用OR值为0.7366(95%CI,0.6502-0.8344),结果表明吸烟相对不吸烟患有高血压的风险较低,为不吸烟组的0.7366倍,此结果与我们平时所认为的情况和文献中所报道60-63的结果有所差异,造成这样的结果的原因很可能是由于数据分析时纳入已知观测混杂因素过少,仅有五个。当我们进行多水平工具变量进行数据拟合时,发现四种分类数据多水平工具变量模型所得的结果较为一致,均提示吸烟为患有高血压的危险因素,回归系数在1.5左右,所得的OR值为5左右,但是多水平工具变量所得结果的置信区间较宽,前面我们分析了所找的工具变量的强度,发现工具变量的强度不是很强,故所得置信区间较宽可能是由于工具变量强度不足造成的。总体看来,多水平工具变量所分析的结果更加符合目前公认两者间的真实关系,普通多水平logistic回归模型所得结果在公共卫生方面可能会引起一个错误的导向,即吸烟更不容易患有高血压,从而鼓励大家吸烟,这样一个有悖于真实情况的结果会误导大众,从而造成较坏的影响。-66- -67-0.834414.293518.255717.935224.8208~~~~~值的置信区间0.65021.16101.40071.11681.26610.73664.07365.05684.47565.6058基于多水平模型的工具变量方法研究及应用OROR-0.18102.65982.90452.88683.2117~~~~~值的置信区间ߚ-0.43040.14920.33700.11050.2359值ߚ-0.30571.40451.62071.49861.7238男性老年人是否吸烟对患有高血压危险分析的结果3-5表回归多水平工具变量模型分析方法回归多水平工具变量模型回归模型回归多水平工具变量模型+logistic回归多水平工具变量模型logisticlogistic+logisticlogistic普通多水平两阶段线性回归自助法两阶段自助法线性回归 基于多水平模型的工具变量方法研究及应用四、讨论此部分主要是将之前所构建的多水平工具变量模型和普通多水平回归模型应用于实际问题分析中,针对先前所构建模型时考虑的数据类型,笔者分别进行了实例应用,根据所获得的结果,可以发现在数据类型是分类变量和连续型变量中,不同多水平工具变量模型之间的结果较为一致,但其同普通多水平回归模型所得的结果存在一定的差异。在分析60岁老年人每周运动时间同自身健康状况间的关系时,我们将男性和女性分别进行分析,男性分析时,多水平工具变量模型所得的结果同普通多水平回归模型的差异较女性分析时大,男性和女性分析时所采用的工具变量强度都很强,所得的结果的置信区间也比较理想,但是此部分结果仅仅是在回归系数上工具变量模型进行了一定成都的校正,同普通多水平回归模型其最终并未获得不同的结论,均为每周运动时间越长,人群的身体健康状况也越好。但当我们进行男性老年人吸烟和高血压间关系分析时,此部分针对的时结局变量和处理/暴露因素均为二分类变量的情况,多水平工具变量模型所得的结果同普通多水平logistic回归模型所得结果存在较大的差异,甚至于逆反了普通多水平logistic回归模型所得的结论,在此部分由于所纳入的已知观测混杂因素数量有限,不能较好的反应研究的真实情况,从而获得了吸烟会降低患高血压风险的结论,但当我们引入了工具变量构建多水平工具变量模型,控制未知观测混杂因素对结果的影响,最终得出吸烟会增高患高血压的风险,此结果更加接近真实情况。从分析结果我们可以发现,之前所构建的多水平工具变量模型可以很好地控制未知观测混杂因素,获得的结果和真实情况下的结果较为接近,这同我们先前进行的模拟研究的结论较为一致。工具变量强度对多水平工具变量模型所得结果具有一定的影响。在进行每周运动时间和自身健康状况关系的研究中我们所寻找的两个工具变量为强工具变量,所得的结果点估计较为一致,置信区间也较窄;但当我们研究吸烟和高血压间的关系时,我们所寻找的工具变量家人是否吸烟仅为中等强度,故在分析的时候多水平工具变量所得结果的置信区间较宽,虽然在此研究中并未影响到最终结论,但也反映出工具变量越弱,其所得的结果越保守(即越容易获得阴性结果)。在今后我们利用多水平工具变量进行实例研究时,需要尽量寻找强度较强的工具变量对数据进行分析,从而可以获得更为准确可靠的结果。我们将自助法引入工具变量模型中时,同模拟研究一样,采用的是分层个例重复抽样法,在引入自助法前,首先获取数据资料的水平2单位各层的样本量从而定义出一个样本框,然后按照原来样本大小进行有放回的重抽样,所得的最终样本量大小不-68- 基于多水平模型的工具变量方法研究及应用变,自助法共进行500次。分析的结果展示,自助法引入后所得的置信区间较宽,在模拟研究中发现自助法所得的结果点估计较为准确,区间也会相应的变宽。层次结构特征数据的自助法还有两种:参数残差自助法和非参数自助法。有研究认为这两种自助法所得结果更加准确,在后续的研究中我们将进一步研究引入不同自助法,研究多水平工具变量模型结果的变化。在连续型变量分析时我们检测了处理/暴露因素在水平2单位上是否存在异质性,从而建议在第一步的回归估计过程中采用多水平回归模型,但对于分类变量,我们并未考虑到处理/暴露因素在水平2单位上的异质性,此部分研究我们也将列入下一步研究计划中,如何更加准确地对处理/暴露因素进行重估计是此部分的难点。综上所述,在实例分析过程中,我们发现所构建的多水平工具变量可以较好地对数据进行分析,尤其是当我们寻找到较为合适的工具变量时,其所得的结果较为可靠。-69- 基于多水平模型的工具变量方法研究及应用第四部分研究结论与展望一、研究结论本研究的主要目的是研究控制层次结构特征数据中未知观测混杂因素对结果的影响,研究针对数据类型的不同构建出不同的多水平工具变量模型,通过数据模拟研究对所构建的模型进行评价,最终将多水平工具变量模型应用于实际数据集中,解决实际问题。在模拟研究过程中我们设定了不同强度的工具变量和未知观测混杂因素,并且从绝对偏倚、标准误、置信区间宽度和置信区间覆盖率四个方面对模型的准确性和精确性进行了评价。当处理/暴露因素和结局变量均为连续型变量并且处理/暴露因素在水平2单位上无异质性时,在无未知观测混杂因素存在的情况下,六种模型均能获得较为准确的点估计结果,但当未知观测混杂因素存在时,普通多水平线性回归模型和自助法多水平线性回归模型所获得结果偏离金标准较远,多水平工具变量所得的结果较为一致,均接近于金标准,但自助法多水平工具变量结果的置信区间较宽,提示结果更加保守。当处理/暴露因素和结局变量均为连续型变量并且处理/暴露因素在水平2单位上存在异质性时,普通多水平线性回归模型仅适用于无未知观测混杂因素的数据中,但两阶段多水平回归工具变量模型、两阶段多水平回归残差纳入工具变量模型、自助法两阶段多水平回归工具变量模型和自助法两阶段多水平回归残差纳入工具变量模型在不同的数据情况下均能得到较为理想的结果,其中自助法两阶段多水平回归工具变量模型的准确度和精确度总体最佳。自助法两阶段最小二乘多水平工具变量模型虽然能基本准确估计出效应值,但置信区间过宽。当处理/暴露因素和结局变量均为分类变量时,普通多水平logistic回归模型仅在无未知观测混杂因素情况比较适用,当存在未知观测混杂因素时,普通多水平logistic回归模型将会产生偏倚较大的结果。对于分类变量的多水平工具变量模型受工具变量强度的影响较大,采用强度较弱的工具变量进行多水平工具变量拟合所获得结果的置信区间较宽,结果较为保守,容易犯二类错误,从而影响到最终结果的准确性,但采用强工具变量,模型结果的置信区间将会明显变窄。四个多水平工具变量模型中以自助法两阶段logistic回归多水平工具变量模型和自助法线性回归+logistic回归多水平工具变量模型两种模型所得结果较好。实例分析过程中我们采用了多水平工具变量模型解决实际问题,发现60岁以上老年人每周运动时间越长,其身体健康状况越好,还发现60岁以上男性老年人吸烟-70- 基于多水平模型的工具变量方法研究及应用会增加患有高血压的风险。在分析过程中我们发现所有多水平工具变量所得结果较为一致,都与普通多水平模型所得结果存在一定的差异,尤其是在分析吸烟和高血压关系的研究中,多水平工具变量模型和普通多水平模型所得结论相反,并且多水平工具变量所得结果更加接近真实情况,但在此项研究中由于我们所寻找的工具变量强度不足,故所得结果的置信区间较宽,但未影响最终结论的获得。二、研究特色和创新点本研究主要是从实际研究中存在的问题出发,根据实际分析的需求进行分析模型的研究和构建,并且通过大量的模拟,参照实际情况,模拟出不同条件下的数据集,如数据类型、混杂强度、工具变量强度等,在模型构建过程中,不仅构建了不同的多水平工具变量模型,并且还将自助法引入模型中,将不同模型对模拟数据进行拟合,对所得到的结果进行客观的评价,获得在不同数据条件下最优模型的选择,最终又将构建的多水平工具变量模型应用于实际问题中,解决实际问题。多水平工具变量模型目前鲜有研究进行报道,在此项研究中构建出了多水平工具变量模型,用来解决层次结构特征数据中未知观测混杂因素对结果产生的影响,并提示在各种数据情况下应当选择的最优模型。研究还针对不同多水平工具变量模型编写了相应的SAS程序,可以直接用于实际数据的分析。三、尚待开展的研究(一)工具变量强弱对模型的影响在模型构建和评价过程中我们发现,工具变量的强度对多水平工具变量模型所分析出来的结果有较大的影响,在弱工具变量时,模型分析出来的结果有一定的偏倚,并且其置信区间很宽,此现象在进行分类变量多水平工具变量模型拟合的过程中尤为明显,模拟过程中在此部分模型中弱工具变量所得结果虽然点估计较为接近于我们所设定的金标准,但由于置信区间较宽,所获得的结论为阴性结果,说明弱工具变量所得的结果过于保守,犯二类错误的可能性较大。在实例分析过程中,连续型变量研究过程中所选取的工具变量强度较强,故最终分析所得的结果较为理想,但当我们进行分类变量实例分析研究时,由于所选取的工具变量的强度仅为中等强度,即使在分析的过程中多水平工具变量模型所得的结果较普通多水平logistic回归模型更为合理,但是其OR值的可信区间过宽,这也提示了分析的准确度并不是很高。在后续的研究中需要针对此问题进行深入研究,通过对模型第一步中重新估计处理/暴露因素的方法进行修改和完善,从而可以在弱工具变量存在的情况下更加准确地进行模型拟合。-71- 基于多水平模型的工具变量方法研究及应用(二)工具变量同结局变量间的关系在本研究中工具变量的设定前提是其与结局变量并无直接关系,在模拟过程中,我们设定工具变量同结局变量间的相关关系为0,但此种情况在现实数据分析过程中是很少见的,在实例分析中,我们仅从理论上认定工具变量同结局变量并无直接相关关系,即使存在相关关系也是间接通过处理/暴露因素而造成的,故在下一步研究过程中我们将要在模拟数据中加上工具变量同结局变量的关系(分为不同强弱),如此便可以更加深入的探讨不同强度的工具变量对多水平工具变量模型的稳定性的影响,若工具变量模型受到其影响过大,我们将会对模型进行进一步的调整。(三)自助法的引入和改进如前所述,目前应用于层次结构特征数据中的自助法有个例重复抽样法、参数残差自助法、非参数残差自助法三种。在本研究中为了节约分析时间,主要是对个例重复抽样法进行了改进,采用了分层个例重复抽样法,在水平2单位的各个层内按照原始样本量大小进行有放回的重复抽样,有学者认为个例重复抽样法在层次结构特征数据中应用时会引起设计矩阵的变异,从而影响到结果的稳定性,但也提示数据样本量46,52较大时这种变异可以被忽略。我们应用了分层个例重复抽样的自助法时发现,其在弱工具变量存在的情况下,所得结果的置信区间比较宽,故在后续的研究中我们将引入参数残差自助法和非参数残差自助法对模型进行拟合,再对模型进行评价。(四)考虑处理/暴露在水平2单位上的异质性我们在连续型变量的多水平工具变量模型构建过程中我们考虑了处理/暴露因素在水平2单位上的异质性,并在模型构建过程中构建出两阶段多水平回归工具变量模型和两阶段多水平回归残差纳入工具变量模型,应用这两种模型对模拟数据进行分析,从而获得了较好的结果。由于分类变量的数据模拟较为复杂,并且第一步用多水平模型重估计处理/暴露因素的方法还有待研究,基于目前我们所构建的分类变量多水平工具变量模型的分析结果,发现其置信区间较宽,对第一步采用多水平logistic回归模型是否会影响结果准确性的问题将继续进行深入研究。。-72- 基于多水平模型的工具变量方法研究及应用附录:核心程序现简要介绍本研究数据模拟和匹配方法实现的SAS程序:一、处理/暴露因素和结局变量为连续型变量libnameresult"E:2014工具变量结果连续性变量iv";/*数据集生成的宏,其中data为所生成数据集的名称;rnum为生成均匀分布时的随机数;jbz为研究中所设的金标准;hz为混杂因素的强度;inst为工具变量同结局变量间的关系,在此研究中默认设为0;gj为工具变量同处理因素间的关系,即为工具变量强度,size为水平2单位上每层中所含有的样本量*/%macroshujuji_new(data,rnum,jbz,hz,inst,gj,size);/*模拟生成结局变量中截距by0,此研究中水平2单位上共有200层*/dataa;doi=1to200;by0=-10+(10+10)*uniform(&rnum);output;end;run;data_null_;setaend=end;callsymput("by0"||left(put(_n_,5.)),left(by0));ifendthencallsymput("num",left(put(_n_,5.)));run;data&data;%doi=1%to#doid=1to&size;site=&i;c1=rand('normal',4,1);c2=rand('normal',2,0.5);u1=rand('normal',5,1);z=rand('normal',3,0.5);x=10+8*c1-6*c2-2*u1+&gj*z;y=&&by0&i+&jbz*x+3*c1-5*c2+&hz*u1+&inst*z;output;end;%end;run;%mend;%macrotiaoxuan(str=,prefix=,num=,split=%str());%localijword;%leti=1;%do%while(%qscan(&str,&i,&split)ne);%letj=&i;%leti=%eval(&i+1);%if%length(&prefix)%then%do;%global&prefix&j;%letword=%scan(&str,&j,&split);%let&prefix&j=&word;-73- 基于多水平模型的工具变量方法研究及应用%put&&&prefix&j;%end;%end;%global#%let&num=%eval(&i-1);%mend;/*此为普通多水平线性回归的宏程序,其中ds为数据集名称;cov为模型中所包含的协变量;newset为多水平线性回归模型回归系数导出的数据集名称*/%macromulreg(ds,cov,newset);odslistingclose;odsoutputSolutionF=&newset;procmixedcovtestdata=&ds;classsite;modely=&cov/solutionCL;randomintercept/subject=site;parms(0.001)(1)/hold=2;run;odsoutputclose;odslisting;%mend;/*此为自助法多水平线性回归的宏程序,其中ds为数据集名称;cov为模型中所包含的协变量;newset为多自助法水平线性回归模型回归系数导出的数据集名称*/%macroboot_mulreg(ds,cov,newset);odslistingclose;odsoutputSolutionF=&newset;procmixedcovtestdata=&ds;byReplicate;classsite;modely=&cov/solution;randomintercept/subject=site;parms(0.001)(1)/hold=2;run;odsoutputclose;odslisting;%mend;/*此为模型拟合的核心程序,其中n为模拟程序的循环次数,在本研究中设为500;hzxs为混杂因素的强度;gjbl为工具变量的强度;iv_xs为工具变量同结局变量间的关系,在此项研究中设为0*/%macroxunhuan(n,hzxs,gjbl,iv_xs);%tiaoxuan(str=&hzxs,prefix=hunza,num=hz_n,split=%str());%tiaoxuan(str=&iv_xs,prefix=ivnew,num=iv_n,split=%str());%tiaoxuan(str=&gjbl,prefix=gjnew,num=gj_n,split=%str());%dom=1%to&gj_n;dataresult.iv_&m;run;%dok=1%to&iv_n;%doj=1%to&hz_n;datajieguo;run;%doxh=1%to&n;%shujuji_new(yuanshi,2014,4,&&hunza&j,&&ivnew&k,&&gjnew&m,5);datayuanshi;setyuanshi;dropid;run;/*此步为自助法的重抽样过程,根据水平2单位site进行,每层5例,共进行500次,最后生成500个复样本数据集boot*/procsurveyselectdata=yuanshimethod=ursn=5outhitsrep=500seed=20140305out=bootnoprint;-74- 基于多水平模型的工具变量方法研究及应用STRATAsite;run;procsortdata=boot;byReplicate;run;/*普通多水平线性回归模型*/%mulreg(yuanshi,xc1c2,xishu);dataxishu1;setxishu(where=(Effect="x"));x_mean=Estimate;x_se=StdErr;x_low=Lower;x_up=Upper;keepx_meanx_sex_lowx_up;run;/*两阶段最小二乘多水平工具变量模型*//*第一步重新得出新的处理因素xx*/procregdata=yuanshinoprint;modelx=c1c2z;outputout=yuanshi_lxPREDICTED=xx;run;quit;/*第二步估计处理因素的效应值*/%mulreg(yuanshi_lx,xxc1c2,xishu);dataxishu2;setxishu(where=(Effect="xx"));iv_x_mean=Estimate;iv_x_se=StdErr;iv_x_low=Lower;iv_x_up=Upper;keepiv_x_meaniv_x_seiv_x_lowiv_x_up;run;/*两阶段残差纳入多水平工具变量模型*//*第一步重新得出残差*/datayuanshi_lx;setyuanshi_lx;ca=xx-x;run;/*第二步估计处理因素的效应值*/%mulreg(yuanshi_lx,xc1c2ca,xishu);dataxishu3;setxishu(where=(Effect="x"));ca_x_mean=Estimate;ca_x_se=StdErr;ca_x_low=Lower;ca_x_up=Upper;keepca_x_meanca_x_seca_x_lowca_x_up;run;/*自助法普通多水平线性回归模型*/%boot_mulreg(boot,xc1c2,xishu);dataxishu;setxishu;ifEffect="x"thenoutput;keepEstimate;run;procunivariatedata=xishunoprint;varEstimate;outputout=xishu4pctlpts=2.55097.5pctlpre=boot_xmean=boot_x_meanSTD=boot_x_sd;run;/*自助法两阶段最小二乘多水平工具变量模型*//*第一步重新得出新的处理因素xx*/procregdata=bootnoprint;byReplicate;modelx=c1c2z/selection=stepwisesle=0.10sls=0.15stb;outputout=boot_lxPREDICTED=xx;run;quit;/*第二步估计获得500个处理因素的效应值*/%boot_mulreg(boot_lx,xxc1c2,xishu);dataxishu;setxishu;ifEffect="xx"thenoutput;keepEstimate;run;/*计算处理因素效应值的均数和标准差,即为最终效应值和相应的标准误*/procunivariatedata=xishunoprint;varEstimate;outputout=xishu5pctlpts=2.55097.5pctlpre=boot_iv_xmean=boot_iv_x_meanSTD=boot_iv_x_sd;run;/*自助法两阶段残差纳入多水平工具变量模型*//*第一步重新得出残差*/-75- 基于多水平模型的工具变量方法研究及应用databoot_lx;setboot_lx;ca=xx-x;run;/*第二步估计获得500个处理因素的效应值*/%boot_mulreg(boot_lx,xc1c2ca,xishu);dataxishu;setxishu;ifEffect="x"thenoutput;keepEstimate;run;procunivariatedata=xishunoprint;varEstimate;outputout=xishu6pctlpts=2.55097.5pctlpre=boot_ca_xmean=boot_ca_x_meanSTD=boot_ca_x_sd;run;dataxishu;mergexishu1xishu2xishu3xishu4xishu5xishu6;run;datajieguo;setjieguoxishu;hz=&&hunza&j;iv=&&ivnew&k;iv_qd=&&gjnew&m;run;dmlog'clear';dmoutput'clear';%end;dataresult.iv_&m;setresult.iv_&mjieguo;run;%end;%end;%end;%mend;/*强工具变量(gjbl=5)情况下的模拟*/%xunhuan(500,6,5,0);dataresult.iv_6_5;setresult.iv_1;run;%xunhuan(500,4,5,0);dataresult.iv_4_5;setresult.iv_1;run;%xunhuan(500,2,5,0);dataresult.iv_2_5;setresult.iv_1;run;%xunhuan(500,0,5,0);dataresult.iv_0_5;setresult.iv_1;run;/*中等强度工具变量(gjbl=3)情况下的模拟*/%xunhuan(500,6,3,0);dataresult.iv_6_3;setresult.iv_1;run;%xunhuan(500,4,3,0);dataresult.iv_4_3;setresult.iv_1;run;%xunhuan(500,2,3,0);dataresult.iv_2_3;setresult.iv_1;run;%xunhuan(500,0,3,0);dataresult.iv_0_3;setresult.iv_1;run;/*弱工具变量(gjbl=1)情况下的模拟*/%xunhuan(500,6,1,0);dataresult.iv_6_1;setresult.iv_1;run;%xunhuan(500,4,1,0);dataresult.iv_4_1;setresult.iv_1;run;%xunhuan(500,2,1,0);dataresult.iv_2_1;setresult.iv_1;run;%xunhuan(500,0,1,0);dataresult.iv_0_1;setresult.iv_1;run;-76- 基于多水平模型的工具变量方法研究及应用二、处理/暴露因素和结局变量为连续型变量并且均存在层级效应libnameresult"E:2014工具变量结果新连续性变量iv";/*数据集生成的宏,其中data为所生成数据集的名称;rnum为生成均匀分布时的随机数;jbz为研究中所设的金标准;hz为混杂因素的强度;inst为工具变量同结局变量间的关系,在此研究中默认设为0;gj为工具变量同处理因素间的关系,即为工具变量强度,size为水平2单位上每层中所含有的样本量*/%macroshujuji_new(data,rnum,jbz,hz,inst,gj,size);dataa;doi=1to200;by0=-10+(10+10)*uniform(&rnum);bx0=-20+(20+20)*uniform(&rnum);/*增加了生成暴露因素中截距bx0,为-20到20均匀分布*/output;end;run;data_null_;setaend=end;callsymput("by0"||left(put(_n_,5.)),left(by0));callsymput("bx0"||left(put(_n_,5.)),left(bx0));ifendthencallsymput("num",left(put(_n_,5.)));run;data&data;%doi=1%to#doid=1to&size;site=&i;c1=rand('normal',4,1);c2=rand('normal',2,0.5);u1=rand('normal',5,1);z=rand('normal',3,0.5);x=&&bx0&i+8*c1-6*c2-2*u1+&gj*z;y=&&by0&i+&jbz*x+3*c1-5*c2+&hz*u1+&inst*z;output;end;%end;run;%mend;%macrotiaoxuan(str=,prefix=,num=,split=%str());%localijword;%leti=1;%do%while(%qscan(&str,&i,&split)ne);%letj=&i;%leti=%eval(&i+1);%if%length(&prefix)%then%do;%global&prefix&j;%letword=%scan(&str,&j,&split);%let&prefix&j=&word;%put&&&prefix&j;%end;%end;%global#%let&num=%eval(&i-1);%mend;-77- 基于多水平模型的工具变量方法研究及应用%macromulreg(ds,cov,newset);odslistingclose;odsoutputSolutionF=&newset;procmixedcovtestdata=&ds;classsite;modely=&cov/solutionCL;randomintercept/subject=site;parms(0.001)(1)/hold=2;run;odsoutputclose;odslisting;%mend;%macroboot_mulreg(ds,cov,newset);odslistingclose;odsoutputSolutionF=&newset;procmixedcovtestdata=&ds;byReplicate;classsite;modely=&cov/solution;randomintercept/subject=site;parms(0.001)(1)/hold=2;run;odsoutputclose;odslisting;%mend;/*此为普通多水平线性回归的宏程序,其中ds为数据集名称;cov为模型中所包含的协变量;newset为多水平线性回归模型回归系数导出的数据集名称*/%macroxunhuan(n,hzxs,gjbl,iv_xs);%tiaoxuan(str=&hzxs,prefix=hunza,num=hz_n,split=%str());%tiaoxuan(str=&iv_xs,prefix=ivnew,num=iv_n,split=%str());%tiaoxuan(str=&gjbl,prefix=gjnew,num=gj_n,split=%str());%dom=1%to&gj_n;dataresult.iv_&m;run;%dok=1%to&iv_n;%doj=1%to&hz_n;datajieguo;run;%doxh=1%to&n;%shujuji_new(yuanshi,2014,4,&&hunza&j,&&ivnew&k,&&gjnew&m,5);datayuanshi;setyuanshi;dropid;run;/*此步为自助法的重抽样过程,根据水平2单位site进行,每层5例,共进行500次,最后生成500个复样本数据集boot*/procsurveyselectdata=yuanshimethod=ursn=5outhitsrep=500seed=20140305out=bootnoprint;STRATAsite;run;procsortdata=boot;byReplicate;run;/*普通多水平线性回归模型*/%mulreg(yuanshi,xc1c2,xishu);dataxishu1;setxishu(where=(Effect="x"));x_mean=Estimate;x_se=StdErr;x_low=Lower;x_up=Upper;keepx_meanx_sex_lowx_up;run;-78- 基于多水平模型的工具变量方法研究及应用/*两阶段多水平回归工具变量模型*//*第一步利用多水平线性回归模型重新得出新的处理因素Pred和残差Resid*/procmixedcovtestdata=yuanshi;classsite;modelx=c1c2z/solutionCLoutp=yuanshi_lx;randomintercept/subject=site;parms(0.001)(1)/hold=2;run;/*第二步估计处理因素的效应值*/%mulreg(yuanshi_lx,Predc1c2,xishu);dataxishu2;setxishu(where=(Effect="Pred"));m_iv_x_mean=Estimate;m_iv_x_se=StdErr;m_iv_x_low=Lower;m_iv_x_up=Upper;keepm_iv_x_meanm_iv_x_sem_iv_x_lowm_iv_x_up;run;/*两阶段多水平回归残差纳入工具变量模型*//*第二步带入前面所获得的残差,估计处理因素的效应值*/%mulreg(yuanshi_lx,xc1c2Resid,xishu);dataxishu3;setxishu(where=(Effect="x"));m_ca_x_mean=Estimate;m_ca_x_se=StdErr;m_ca_x_low=Lower;m_ca_x_up=Upper;keepm_ca_x_meanm_ca_x_sem_ca_x_lowm_ca_x_up;run;/*利用多水平线性回归模型重新得出500个复样本中新的处理因素Pred和残差Resid*/procmixedcovtestdata=boot;byReplicate;classsite;modelx=c1c2z/solutionCLoutp=boot_lx;randomintercept/subject=site;parms(0.001)(1)/hold=2;run;/*自助法两阶段多水平回归工具变量模型*/%boot_mulreg(boot_lx,Predc1c2,xishu);dataxishu;setxishu;ifEffect="Pred"thenoutput;keepEstimate;run;procunivariatedata=xishunoprint;varEstimate;outputout=xishu4pctlpts=2.55097.5pctlpre=boot_m_iv_xmean=boot_m_iv_x_meanSTD=boot_m_iv_x_sd;run;/*自助法两阶段多水平回归残差纳入工具变量模型*/%boot_mulreg(boot_lx,xc1c2Resid,xishu);dataxishu;setxishu;ifEffect="x"thenoutput;keepEstimate;run;procunivariatedata=xishunoprint;varEstimate;outputout=xishu5pctlpts=2.55097.5pctlpre=boot_m_ca_xmean=boot_m_ca_x_meanSTD=boot_m_ca_x_sd;run;/*自助法两阶段最小二乘多水平工具变量模型*//*第一步重新得出新的处理因素xx*/procregdata=bootnoprint;byReplicate;modelx=c1c2z/selection=stepwisesle=0.10sls=0.15stb;outputout=boot_lxPREDICTED=xx;run;quit;/*第二步估计获得500个处理因素的效应值*/%boot_mulreg(boot_lx,xxc1c2,xishu);-79- 基于多水平模型的工具变量方法研究及应用dataxishu;setxishu;ifEffect="xx"thenoutput;keepEstimate;run;procunivariatedata=xishunoprint;varEstimate;outputout=xishu6pctlpts=2.55097.5pctlpre=boot_iv_xmean=boot_iv_x_meanSTD=boot_iv_x_sd;run;dataxishu;mergexishu1xishu2xishu3xishu4xishu5xishu6;run;datajieguo;setjieguoxishu;hz=&&hunza&j;iv=&&ivnew&k;iv_qd=&&gjnew&m;run;dmlog'clear';dmoutput'clear';%end;dataresult.iv_&m;setresult.iv_&mjieguo;run;%end;%end;%end;%mend;/*强工具变量(gjbl=5)情况下的模拟*/%xunhuan(500,6,5,0);dataresult.iv_6_5;setresult.iv_1;run;%xunhuan(500,4,5,0);dataresult.iv_4_5;setresult.iv_1;run;%xunhuan(500,2,5,0);dataresult.iv_2_5;setresult.iv_1;run;%xunhuan(500,0,5,0);dataresult.iv_0_5;setresult.iv_1;run;/*中等强度工具变量(gjbl=3)情况下的模拟*/%xunhuan(500,6,3,0);dataresult.iv_6_3;setresult.iv_1;run;%xunhuan(500,4,3,0);dataresult.iv_4_3;setresult.iv_1;run;%xunhuan(500,2,3,0);dataresult.iv_2_3;setresult.iv_1;run;%xunhuan(500,0,3,0);dataresult.iv_0_3;setresult.iv_1;run;/*弱工具变量(gjbl=1)情况下的模拟*/%xunhuan(500,6,1,0);dataresult.iv_6_1;setresult.iv_1;run;%xunhuan(500,4,1,0);dataresult.iv_4_1;setresult.iv_1;run;%xunhuan(500,2,1,0);dataresult.iv_2_1;setresult.iv_1;run;%xunhuan(500,0,1,0);dataresult.iv_0_1;setresult.iv_1;run;三、处理/暴露因素和结局变量为分类变量libnameresult"E:2014工具变量结果分类变量iv";/*数据集生成的宏,其中data为所生成数据集的名称;rnum为生成均匀分布时的随机数;jbz为研究中所设的金标准,本研究设为2;hz为混杂因素的强度;inst为调整生成结局变量中0和1比例的参-80- 基于多水平模型的工具变量方法研究及应用数;gj为工具变量同处理因素间的关系,即为工具变量强度,size为水平2单位上每层中所含有的样本量*/%macroshujuji_new(data,jbz,hz,inst,gj,size);/*模拟生成结局变量中截距by0,此研究中水平2单位上共有200层*/dataa;doi=1to200;by0=rand('normal',0,0.5);output;end;run;data_null_;setaend=end;callsymput("by0"||left(put(_n_,5.)),left(by0));ifendthencallsymput("num",left(put(_n_,5.)));run;data&data;%doi=1%to#doid=1to&size;site=&i;c1=rand('normal',2,0.5);c2=rand('bernoulli',0.4);u=rand('normal',1,0.5);z=rand('bernoulli',0.5);ptx=-0.5-1*c1+1*c2-1*u+&gj*z;pix=1/(1+exp(-ptx));x=rand('bernoulli',pix);pt=&&by0&i+&jbz*x-1*c1+1*c2+&hz*u-&inst;pi=1/(1+exp(-pt));y=rand('bernoulli',pi);output;end;%end;run;%mend;%macrotiaoxuan(str=,prefix=,num=,split=%str());%localijword;%leti=1;%do%while(%qscan(&str,&i,&split)ne);%letj=&i;%leti=%eval(&i+1);%if%length(&prefix)%then%do;%global&prefix&j;%letword=%scan(&str,&j,&split);%let&prefix&j=&word;%put&&&prefix&j;%end;%end;%global#%let&num=%eval(&i-1);%mend;/*此为普通多水平线性回归的宏程序,其中data为数据集名称;&ybl为模型中的结局变量;group-81- 基于多水平模型的工具变量方法研究及应用为模型中处理因素;var为其他已观测混杂因素。在此过程中首先应用glimmix过程初步估计出各变量的系数,然后应用nlmixed过程进一步进行估计,回归系数导入数据集xishu中*/%macromutlogistic(data,ybl,group,var);odslistingclose;odsoutputParameterEstimates=xishu;procglimmixdata=&datamethod=rspl;classsite;model&ybl(event='1')=&group&var/sdist=binarylink=logitddfm=bw;randomint/subject=site;nloptionstech=nrridg;run;odsoutputclose;odslisting;dataxishu;setxishu;parm=compress("b"||Effect||"="||Estimate);run;data_null_;setxishuend=end;callsymput("effect_new"||left(put(_n_,2.)),left(Effect));callsymput("parm_new"||left(put(_n_,2.)),left(parm));ifendthencallsymput("num_parm",left(put(_n_,2.)));run;%letparm=&parm_new1;%doi=2%to&num_parm;%letparm=&parm&&parm_new&i;%end;%doi=2%to&num_parm;%letbeffect_new&i=b&&effect_new&i;%end;%letgongshi=b&effect_new1+u0j;%doi=2%to&num_parm;%letgongshi=&gongshi+&&beffect_new&i*&&effect_new&i;%end;odslistingclose;odsoutputParameterEstimates=xishu;procnlmixeddata=&data;parms&parmv_u0=0.5;z=&gongshi;if(&ybl=1)thenp=1/(1+exp(-z));elsep=1-(1/(1+exp(-z)));ll=log(p);model&ybl~general(ll);randomu0j~normal(0,v_u0)subject=site;run;odsoutputclose;odslisting;%mend;/*此为自助法普通多水平线性回归的宏程序,其中data_old为原始样本的数据集;data为自助法所受样本的数据集名称;&ybl为模型中的结局变量;group为模型中处理因素;var为其他已观测混杂因素。回归系数导入数据集xishu中*/%macroboot_mutlogistic(data_old,data,ybl,group,var);-82- 基于多水平模型的工具变量方法研究及应用odslistingclose;odsoutputParameterEstimates=xishu;procglimmixdata=&data_oldmethod=rspl;classsite;model&ybl(event='1')=&group&var/sdist=binarylink=logitddfm=bw;randomint/subject=site;nloptionstech=nrridg;run;odsoutputclose;odslisting;dataxishu;setxishu;parm=compress("b"||Effect||"="||Estimate);run;data_null_;setxishuend=end;callsymput("effect_new"||left(put(_n_,2.)),left(Effect));callsymput("parm_new"||left(put(_n_,2.)),left(parm));ifendthencallsymput("num_parm",left(put(_n_,2.)));run;%letparm=&parm_new1;%doi=2%to&num_parm;%letparm=&parm&&parm_new&i;%end;%doi=2%to&num_parm;%letbeffect_new&i=b&&effect_new&i;%end;%letgongshi=b&effect_new1+u0j;%doi=2%to&num_parm;%letgongshi=&gongshi+&&beffect_new&i*&&effect_new&i;%end;odslistingclose;odsoutputParameterEstimates=xishu;procnlmixeddata=&data;byReplicate;parms&parmv_u0=0.5;z=&gongshi;if(&ybl=1)thenp=1/(1+exp(-z));elsep=1-(1/(1+exp(-z)));ll=log(p);model&ybl~general(ll);randomu0j~normal(0,v_u0)subject=site;run;odsoutputclose;odslisting;%mend;/*此为模型拟合的核心程序,其中n为模拟程序的循环次数,在本研究中设为500;hzxs为混杂因素的强度;gjbl为工具变量的强度;iv_xs为工具变量同结局变量间的关系,在此项研究中设为0*/%macroxunhuan(n,hzxs,gjbl,iv_xs);%tiaoxuan(str=&hzxs,prefix=hunza,num=hz_n,split=%str());%tiaoxuan(str=&iv_xs,prefix=ivnew,num=iv_n,split=%str());%tiaoxuan(str=&gjbl,prefix=gjnew,num=gj_n,split=%str());%dom=1%to&gj_n;-83- 基于多水平模型的工具变量方法研究及应用dataresult.iv_&m;run;%dok=1%to&iv_n;%doj=1%to&hz_n;datajieguo;run;%doxh=1%to&n;%shujuji_new(yuanshi,2,&&hunza&j,&&ivnew&k,&&gjnew&m,10);datayuanshi;setyuanshi;dropiduptxpixptpi;run;/*此步为自助法的重抽样过程,根据水平2单位site进行,每层10例,共进行500次,最后生成500个复样本数据集boot*/procsurveyselectdata=yuanshimethod=ursn=10outhitsrep=500seed=20140305out=bootnoprint;STRATAsite;run;procsortdata=boot;byReplicate;run;/*普通多水平logistic回归模型*/%mutlogistic(yuanshi,y,x,c1c2);dataxishu1;setxishu;mutlog_x=Estimate;mutlog_x_low=Lower;mutlog_x_up=Upper;mutlog_x_se=StandardError;ifParameter="bx"thenoutput;keepmutlog_xmutlog_x_semutlog_x_lowmutlog_x_up;run;/*两阶段logistic回归多水平工具变量模型*//*第一步应用logistic回归重新得出新的处理因素xx*/proclogisticdata=yuanshiDESCENDINGnoprint;modelx=c1c2z;outputout=yuanshi_newpred=xx;run;/*第二步估计处理因素的效应值*/%mutlogistic(yuanshi_new,y,xx,c1c2);dataxishu2;setxishu;iv_x=Estimate;iv_x_low=Lower;iv_x_up=Upper;iv_x_se=StandardError;ifParameter="bxx"thenoutput;keepiv_xiv_x_seiv_x_lowiv_x_up;run;/*线性回归+logistic回归多水平工具变量模型*//*第一步应用线性回归重新得出新的处理因素xx*/procregdata=yuanshinoprint;modelx=c1c2z;outputout=yuanshi_lxPREDICTED=xx;run;quit;/*第二步估计处理因素的效应值*/%mutlogistic(yuanshi_lx,y,xx,c1c2);dataxishu3;setxishu;lx_x=Estimate;lx_x_low=Lower;lx_x_up=Upper;lx_x_se=StandardError;ifParameter="bxx"thenoutput;keeplx_xlx_x_selx_x_lowlx_x_up;run;/*自助法多水平logistic回归模型*/%boot_mutlogistic(yuanshi,boot,y,x,c1c2);dataxishu;setxishu;ifParameter="bx"thenoutput;keepEstimate;run;procunivariatedata=xishunoprint;varEstimate;outputout=xishu4pctlpts=2.55097.5pctlpre=boot_xmean=boot_x_meanSTD=boot_x_sd;run;-84- 基于多水平模型的工具变量方法研究及应用/*自助法两阶段logistic回归多水平工具变量模型*//*第一步应用logistic回归重新得出新的处理因素xx*/proclogisticdata=bootDESCENDINGnoprint;byReplicate;modelx=c1c2z;outputout=boot_ivpred=xx;run;/*第二步估计处理因素的效应值*/%boot_mutlogistic(yuanshi_new,boot_iv,y,xx,c1c2);dataxishu;setxishu;ifParameter="bxx"thenoutput;keepEstimate;run;procunivariatedata=xishunoprint;varEstimate;outputout=xishu5pctlpts=2.55097.5pctlpre=boot_iv_xmean=boot_iv_x_meanSTD=boot_iv_x_sd;run;/*自助法线性回归+logistic回归多水平工具变量模型*//*第一步应用线性回归重新得出新的处理因素xx*/procregdata=bootnoprint;byReplicate;modelx=c1c2z/selection=stepwisesle=0.10sls=0.15stb;outputout=boot_lxPREDICTED=xx;run;quit;/*第二步估计处理因素的效应值*/%boot_mutlogistic(yuanshi_lx,boot_lx,y,xx,c1c2);dataxishu;setxishu;ifParameter="bxx"thenoutput;keepEstimate;run;procunivariatedata=xishunoprint;varEstimate;outputout=xishu6pctlpts=2.55097.5pctlpre=boot_lx_xmean=boot_lx_x_meanSTD=boot_lx_x_sd;run;dataxishu;mergexishu1xishu2xishu3xishu4xishu5xishu6;run;datajieguo;setjieguoxishu;hz=&&hunza&j;iv=&&ivnew&k;iv_qd=&&gjnew&m;run;dmlog'clear';dmoutput'clear';%end;dataresult.iv_&m;setresult.iv_&mjieguo;run;%end;%end;%end;%mend;/*强工具变量(gjbl=3)情况下的模拟*/%xunhuan(500,2,3,1);dataresult.iv_2_3;setresult.iv_1;run;%xunhuan(500,1,3,0);dataresult.iv_1_3;setresult.iv_1;run;%xunhuan(500,0,3,0);dataresult.iv_0_3;setresult.iv_1;run;/*中等强度工具变量(gjbl=2)情况下的模拟*/%xunhuan(500,2,2,1);-85- 基于多水平模型的工具变量方法研究及应用dataresult.iv_2_2;setresult.iv_1;run;%xunhuan(500,1,2,0);dataresult.iv_1_2;setresult.iv_1;run;%xunhuan(500,0,2,0);dataresult.iv_0_2;setresult.iv_1;run;/*弱工具变量(gjbl=1)情况下的模拟*/%xunhuan(500,2,1,1);dataresult.iv_2_1;setresult.iv_1;run;%xunhuan(500,1,1,0);dataresult.iv_1_1;setresult.iv_1;run;%xunhuan(500,0,1,0);dataresult.iv_0_1;setresult.iv_1;run;-86- 基于多水平模型的工具变量方法研究及应用文献综述-87- 基于多水平模型的工具变量方法研究及应用-88- 基于多水平模型的工具变量方法研究及应用-89- 基于多水平模型的工具变量方法研究及应用-90- 基于多水平模型的工具变量方法研究及应用-91- 基于多水平模型的工具变量方法研究及应用-92- 基于多水平模型的工具变量方法研究及应用-93- 基于多水平模型的工具变量方法研究及应用-94- 基于多水平模型的工具变量方法研究及应用-95- 基于多水平模型的工具变量方法研究及应用-96- 基于多水平模型的工具变量方法研究及应用-97- 基于多水平模型的工具变量方法研究及应用-98- 基于多水平模型的工具变量方法研究及应用-99- 基于多水平模型的工具变量方法研究及应用参考文献1.姜辉.新医改背景下的卫生信息化建设探讨.中国医院管理2011;31(8):71-2.2.徐庐生.近十年来医院信息化的发展.中国医疗器械信息2010;16(3):1-8.3.许健,查佳凌,尤超,唐静雯,李先锋,吴韬.医疗信息化集成平台在医院的建设与思考.中国医院2012;(2):5-8.4.李刚,杨宇军.社区医院信息化建设过程中的问题与对策.牡丹江医学院学报2012;33(1):79-80.5.曲晨.中国健康教育信息化平台建设研究进展.江苏预防医学2014;25(5):93-4.6.http://www.cebm.net/levels_of_evidence.asp.7.陈耀龙,李幼平,杜亮,王莉,文进,杨晓妍.医学研究中证据分级和推荐强度的演进.中国循证医学杂志2008;8(2):127-33.8.史楠楠,王思成,韩学杰,王丽颖,赵静,吕爱平.证据分级体系的演进及其对中医临床实践指南的启示.北京中医药大学学报2011;34(2):87-91.9.BrittonA,McKeeM,BlackN,McPhersonK,SandersonC,BainC.Choosingbetweenrandomisedandnon-randomisedstudies:asystematicreview.Healthtechnologyassessment1998;2(13):i-iv,1-124.10.ConcatoJ,ShahN,HorwitzRI.Randomized,controlledtrials,observationalstudies,andthehierarchyofresearchdesigns.TheNewEnglandjournalofmedicine2000;342(25):1887-92.11.高妍,杨文英.观察性研究在临床实践中的重要意义.中华内科杂志2009;48(6):449-51.12.FoxKA.Anintroductiontotheglobalregistryofacutecoronaryevents:GRACE.2000.13.SharpeN.Clinicaltrialsandtherealworld:selectionbiasandgeneralisabilityoftrialresults.Cardiovasculardrugsandtherapy/sponsoredbytheInternationalSocietyofCardiovascularPharmacotherapy2002;16(1):75-7.14.BlackN.Whyweneedobservationalstudiestoevaluatetheeffectivenessofhealthcare.Bmj1996;312(7040):1215-8.-100- 基于多水平模型的工具变量方法研究及应用15.GrzeskowiakLE,GilbertAL,MorrisonJL.Investigatingoutcomesassociatedwithmedicationuseduringpregnancy:areviewofmethodologicalchallengesandobservationalstudydesigns.Reproductivetoxicology2012;33(3):280-9.16.SheikhA,SmeethL,AshcroftR.Randomisedcontrolledtrialsinprimarycare:scopeandapplication.TheBritishjournalofgeneralpractice:thejournaloftheRoyalCollegeofGeneralPractitioners2002;52(482):746-51.17.VandenbrouckeJP.Whenareobservationalstudiesascredibleasrandomisedtrials?Lancet2004;363(9422):1728-31.18.GrimesDA,SchulzKF.Biasandcausalassociationsinobservationalresearch.TheLancet2002;359(9302):248-52.19.RosenbaumPR,RubinDB.Thecentralroleofthepropensityscoreinobservationalstudiesforcausaleffects.Biometrika1983;70(1):41-55.20.RosenbaumPR,RubinDB.Reducingbiasinobservationalstudiesusingsubclassificationonthepropensityscore.JournaloftheAmericanStatisticalAssociation1984;79(387):516-24.21.KatajainenJ,RaitaT.Ananalysisofthelongestmatchandthegreedyheuristicsintextencoding.JournaloftheACM(JACM)1992;39(2):281-94.22.RubinDB.BiasreductionusingMahalanobis-metricmatching.Biometrics1980:293-8.23.HiranoK,ImbensGW.Estimationofcausaleffectsusingpropensityscoreweighting:Anapplicationtodataonrightheartcatheterization.HealthServicesandOutcomesresearchmethodology2001;2(3-4):259-78.24.SchrammTK,GislasonGH,VaagA,etal.Mortalityandcardiovascularriskassociatedwithdifferentinsulinsecretagoguescomparedwithmetforminintype2diabetes,withorwithoutapreviousmyocardialinfarction:anationwidestudy.Europeanheartjournal2011;32(15):1900-8.25.SweredoskiMJ,BaldiP.COBEpro:anovelsystemforpredictingcontinuousB-cellepitopes.Proteinengineering,design&selection:PEDS2009;22(3):113-20.26.SetoguchiS,SchneeweissS,BrookhartMA,GlynnRJ,CookEF.Evaluatingusesofdataminingtechniquesinpropensityscoreestimation:asimulationstudy.Pharmacoepidemiologyanddrugsafety2008;17(6):546-55.-101- 基于多水平模型的工具变量方法研究及应用27.HarderVS,MorralAR,ArkesJ.Marijuanauseanddepressionamongadults:Testingforcausalassociations.Addiction2006;101(10):1463-72.28.StoneRA,ObroskyDS,SingerDE,KapoorWN,FineMJ.Propensityscoreadjustmentforpretreatmentdifferencesbetweenhospitalizedandambulatorypatientswithcommunity-acquiredpneumonia.PneumoniaPatientOutcomesResearchTeam(PORT)Investigators.Medicalcare1995;33(4Suppl):AS56-66.29.McCandlessLC,GustafsonP,AustinPC.Bayesianpropensityscoreanalysisforobservationaldata.Statisticsinmedicine2009;28(1):94-112.30.KaplanD,ChenJ.ATwo-StepBayesianApproachforPropensityScoreAnalysis:SimulationsandCaseStudy.Psychometrika2012;77(3):581-609.31.JMR.Marginalstructuralmodels.In:1997ProceedingsoftheSectiononBayesianStatisticalScience,Alexandria,VA:AmericanStatisticalAssociation.1998;(1–10).32.RobinsJM,HernanMA,BrumbackB.Marginalstructuralmodelsandcausalinferenceinepidemiology.Epidemiology2000;11(5):550-60.33.王济川,谢海义,姜宝法.多层统计分析模型:方法与应用:高等教育出版社;2008.34.杨珉,李晓松.医学和公共卫生研究常用多水平统计模型:北京大学医学出版社;2007.35.BrookhartMA,WangPS,SolomonDH,SchneeweissS.Evaluatingshort-termdrugeffectsusingaphysician-specificprescribingpreferenceasaninstrumentalvariable.Epidemiology2006;17(3):268-75.36.TerzaJV,BasuA,RathouzPJ.Two-stageresidualinclusionestimation:addressingendogeneityinhealtheconometricmodeling.Journalofhealtheconomics2008;27(3):531-43.37.RassenJA,SchneeweissS,GlynnRJ,MittlemanMA,BrookhartMA.Instrumentalvariableanalysisforestimationoftreatmenteffectswithdichotomousoutcomes.Americanjournalofepidemiology2009;169(3):273-84.38.CaiB,HennessyS,FloryJH,ShaD,TenHaveTR,SmallDS.Simulationstudyofinstrumentalvariableapproacheswithanapplicationtoastudyoftheantidiabeticeffectofbezafibrate.Pharmacoepidemiologyanddrugsafety2012;21Suppl2:114-20.39.AmemiyaT.Qualitativeresponsemodels:Asurvey.Journalofeconomic-102- 基于多水平模型的工具变量方法研究及应用literature1981:1483-536.40.JohnstonK,GustafsonP,LevyA,GrootendorstP.Useofinstrumentalvariablesintheanalysisofgeneralizedlinearmodelsinthepresenceofunmeasuredconfoundingwithapplicationstoepidemiologicalresearch.Statisticsinmedicine2008;27(9):1539-56.41.LuB,MarcusS.Evaluatinglong-termeffectsofapsychiatrictreatmentusinginstrumentalvariableandmatchingapproaches.HealthServicesandOutcomesResearchMethodology2012;12(4):288-301.42.蔡敏,饶克勤,钱军程,高军,徐玲.3次国家卫生服务调查居民主要疾病患病情况变化与分析.中国医院统计2005;12(2):107-11.43.石光,贡森.改革开放以来中国卫生投入及其绩效分析.中国发展评论2005;1:32-42.44.陈心广,魏晟,饶克勤.中国城市基本医疗服务需求弹性经济学模型研究.中国卫生经济1996;2:55-7.45.HausmanJA.Specificationtestsineconometrics.Econometrica:JournaloftheEconometricSociety1978:1251-71.46.DavisonAC.Bootstrapmethodsandtheirapplication:Cambridgeuniversitypress;1997.47.EfronB.Bootstrapmethods:anotherlookatthejackknife.TheannalsofStatistics1979:1-26.48.MooneyCZ,DuvalRD.Bootstrapping:Anonparametricapproachtostatisticalinference:Sage;1993.49.EfronB,EfronB.Thejackknife,thebootstrapandotherresamplingplans:SIAM;1982.50.MillerRG.Thejackknife-areview.Biometrika1974;61(1):1-15.51.EfronB.Nonparametricstandarderrorsandconfidenceintervals.CanadianJournalofStatistics1981;9(2):139-58.52.CarpenterJR,GoldsteinH,RasbashJ.Anovelbootstrapprocedureforassessingtherelationshipbetweenclasssizeandachievement.JournaloftheRoyalStatisticalSociety:SeriesC(AppliedStatistics)2003;52(4):431-43.53.第五次国家卫生服务调查方案,http://www.nhfpc.gov.cn/mohwsbwstjxxzx/s8211/201308/cecaaee775f849cea0186cd23a4fbcba.shtml.-103- 基于多水平模型的工具变量方法研究及应用54.SavoiaE,FantiniMP,PandolfiPP,DallolioL,CollinaN.AssessingtheconstructvalidityoftheItalianversionoftheEQ-5D:preliminaryresultsfromacross-sectionalstudyinNorthItaly.Healthandqualityoflifeoutcomes2006;4:47.55.MathewsWC,MayS.EuroQol(EQ-5D)measureofqualityoflifepredictsmortality,emergencydepartmentutilization,andhospitaldischargeratesinHIV-infectedadultsundercare.Healthandqualityoflifeoutcomes2007;5:5.56.BrooksR.EuroQol:thecurrentstateofplay.Healthpolicy1996;37(1):53-72.57.StockJH,WrightJH,YogoM.Asurveyofweakinstrumentsandweakidentificationingeneralizedmethodofmoments.JournalofBusiness&EconomicStatistics2002;20(4).58.PrattN,RougheadEE,RyanP,SalterA.Antipsychoticsandtheriskofdeathintheelderly:aninstrumentalvariableanalysisusingtwopreferencebasedinstruments.Pharmacoepidemiologyanddrugsafety2010;19(7):699-707.59.Ionescu‐IttuR,DelaneyJA,AbrahamowiczM.Bias–variancetrade‐offinpharmacoepidemiologicalstudiesusingphysician‐preference‐basedinstrumentalvariables:asimulationstudy.Pharmacoepidemiologyanddrugsafety2009;18(7):562-71.60.TalukderMH,JohnsonWM,VaradharajS,etal.Chroniccigarettesmokingcauseshypertension,increasedoxidativestress,impairedNObioavailability,endothelialdysfunction,andcardiacremodelinginmice.AmericanJournalofPhysiology-HeartandCirculatoryPhysiology2011;300(1):H388-H96.61.BullenC.Impactoftobaccosmokingandsmokingcessationoncardiovascularriskanddisease.2008.62.李素清,房婉波.农村局部地区高血压患病率与吸烟之间的关系.临床医药实践2013;22(7):529-31.63.向全永,潘晓群,吕淑荣,武鸣.吸烟与高血压相关关系研究.中国预防医学杂志2010;(11):1129-31.-104-
此文档下载收益归作者所有