|
选自即将面世新书《疯狂统计学》第6章 1 引言 2016年底，针刺临床研究论文《针刺治疗慢性严重功能性便秘的随机对照试验》见刊于内科学年鉴，为针刺治疗慢性便秘的疗效提供了国际认可的高质量循证医学证据。该研究共纳入1075例患者，历时3年完成，共有国内15家三甲医院参与。这项随机对照研究通过中心分层区组随机化的方法，设立电针组与假电针组两组进行平行对照，分别接受8周电针针刺穴位天枢、腹结、上巨虚与假电针浅刺双侧天枢旁、腹结旁、上巨虚旁非穴位点。结果发现8周电针治疗可以安全有效地增加慢性严重便秘患者的完全自主排便次数，治疗结束后疗效可持续12周。 随访期的设置对于评价针刺的远期疗效具有重要的意义，该研究的随访期为12周，但与之相矛盾的是随着随访期长度的增加，参与试验受试者的依从性则会下降，即脱落率与随访时间的长度呈正相关关系。所以，如何利用恰当的统计分析方法处理缺失数据对于一项随机对照试验十分重要。本文主要介绍内容为在这项随机对照试验中，针对缺失值处理的具体方法及基于SAS软件的操作步骤进行说明。本文首先介绍如何使用SAS软件对试验基线数据进行统计描述。 2 研究项目的基线描述 2.1 问题与数据 本研究共纳入1 075例患者，人口统计学资料及基线特征按组别采用描述性统计方法展示，有关数据资料如下图6-1所示。 2.2 数据结构的分析 对研究资料图6-1而言，试验因素是针灸治疗(电针组和假电针组)和医院(15家医院)，具体变量信息如表6-1所示。
图6-1 患者数据整理 2.3 分析的目的与统计分析方法的选择 分析目的及方法：数据从excel导入SAS数据库，按组别计算年龄、体重指数等连续性变量的平均值、中位数、四分位数、标准差、最小值和最大值；性别、民族、等分类变量提供各类别的例数和百分比。 2.4 SAS程序及说明 proc import datafile='D:\test.xls' 表6-1 变量信息定义表
out=test dbms=excel replace; getnames=yes; run; 程序说明：语句“proc import”表示调用import过程，进行外部文件导入SAS过程。该语句中的“datafile=”规定要读入外部文件的地址及名称；“out=”规定要输出的SAS数据集，并命名为test；“dbms”选项规定外部数据文件格式的标识名，excel文件可以统一设置为dbms=excel；replace规定替换已存在的文件，便于我们重复导入。“getnames=Yes”规定了导入数据第一行为变量名。 proc means data=test N NMISS MEAN STD MIN MAX Median Q1 Q3 CLM; var BMI Course; class group; run; 程序说明：以计算连续性变量BMI为例，语句“ proc means”表示调用means过程，计算变量的基本描述统计量。该语句中的“class”定义观测组并分别计算观测组的描述统计量；“VAR”规定要求计算简单描述统计量的数值变量和次序；“N NMISS MEAN STD MIN MAX Median Q1 Q3”规定了描述统计量关键词分别是样本量、含缺失值样本数量、均数、标准差、最小值、最大值、中位数、1/4分位数和3/4分位数。 proc freq data=test; tables sex*group/nopercent nocol chisq; run; 程序说明：该程序以分类变量性别举例，关键选项table语句表示输出4个表，同时不显示行百分比(nocol)、以及百分数(nopercent)；chisq目的是输出集中常用的卡方检验统计量，列联表为4个表时，还可以输出Fisher精确检验结果。 2.5 主要分析结果及解释 图6-2为在基线描述过程中proc import、proc means以及proc freq语句在SAS软件中的基本格式。 图6-3为以连续性变量BMI以及age为例，计算变量的基本描述统计量后的 输出结果。 图6-4为以分类变量sex为例，计算变量的基本描述统计量后的输出结果。
图6-2 基线描述SAS程序代码
图6-3 两组基线连续性指标分析结果
图6-4 两组基线分类指标分析结果 3 主要指标的敏感性分析 3.1 问题与数据 在实际临床试验中，常常会出现缺失数据。数据缺失有各种形式：譬如在临床试验中患者中途退出试验研究，或者在某些时间点没有进行检查，或者患者不情愿回答某些项目。这些缺失的原因也许与患者的病情有关，也许无关。 本研究对主要指标的缺失值采用多重填补方法(Multiple imputation，MI)，该方法通常假设数据缺失机制为随机缺失(Missing at random，MAR)。随机缺失是目前较常见处理缺失值的方法，指反应变量的缺失只依赖于已观测到的反应变量値，而与未观测到的反应变量値无关。然而，缺失数据无法被观测，导致随机缺失假设无法被证实。对于采用假定随机缺失的临床试验，国家研究委员会建议对违背MAR的统计假设进行统计推断，认为主要指标缺失机制的敏感性分析是统计报告中的必要组成部分。 本研究采用两组等比例多中心随机区组设计，治疗组为接受电针试验组，对照组为接受假电针安慰剂组。变量Group为分组变量，Y0为主要指标CSBM基线得分，Y1~Y4为治疗期和随访期CSBM不同时间点测量值，有关资料如表6-2所示。 3.2 数据结构分析 对研究资料表6-2而言，试验因素是针灸治疗(电针组和假电针组)，试验效应是治疗效果(Y0-Y4)，具体变量信息如表6-3所示。 3.3 分析目的与统计分析方法选择 据本研究统计分析计划，假定缺失数据机制为随机缺失，对研究主要指标的缺失值处理采用了多重填补方法。本章节分析目的是，抽取378例患者数据，针对多重填补方式的随机缺失假定进行敏感性分析，采用模式混合模型(Pattern mixture model，PMM)(Little 1993；Molenberghs and Kenward 2007，pp.30，34–37)，该方法不仅适用于多种缺失机制的数据，且可在机制下得到更为准确的参数估计值；并不要求一定知道缺失机制的具体分布形式，模型参数估计稳健，该方法被国际医学杂志如JAMA，内科学年鉴等杂志的统计委员一致推荐(http:///aim/pages/AuthorInformationStatisticsOnly)。然而，该方法需要估计较多的特定模式参数，故一般形式下的PMM通常可识别性较低。为了增强模型的识别能力，需设定一些约束条件。根据本研究中观察到的患者CSBM最大为10次，模型约束条件限定值选取10，该模型的基本原理及参数估计方法详见。 表6-2 主要指标数据结构
表6-3 变量信息定义表
3.4 SAS程序及说明 proc mi data=Mono3 seed=1423741; class group; var y0 y1 y2 y3 y4; run; 程序说明：首先检查主要指标缺失类型。语句“proc mi”表示调用MI过程，该语句中的“data=test”是指分析的数据集名称为“test”；“seed=”是指种子数，可以随机输入；“Class”选项要求输入分组变量；Var选项要求Y0~Y4进行缺失值分析。 %midata( data=Mono2, smin=-2, smax=0, sinc=0.1, out=out1); 程序说明：%midata是sas批处理程序，详见附录1。主要功能是对PMM模型进行条件约束，并输入结果。Smin选项要求约束条件最小值；smax选项要求约束条件最大值；sinc表示最小值到最大值的步长；out选项输入最终数据集为out1。 proc reg data=out1; model y1= Trt y0; by Shift _Imputation_; ods output parameterestimates=regparms; run; 程序说明：语句“proc reg”表示调用reg过程，该语句中的“data=test”是指分析的数据集名称为“test”；Model选项中等号前放入因变量，等号后放入可能影响疗效指标的因素；by选项要求数据分析按照模型不同数据集及不同约束条件进行分析，其中“shift”为PMM模型约束条件；“_Imputation_”为多重填补数据集标识变量；“ods output parameterestimates=regparms”选项要求回归模型结果输出到同一数据集中，并命名该数据集为“regparms”。 proc mianalyze parms=regparms; modeleffects Trt; by Shift; ods output parameterestimates=miparm1; run; 程序说明：语句“proc mianalyze”表示调用mianalyze过程，该语句中的parms选项制定数据来源集为regparms；modeleffects选项要求整合研究中组间结果，by选项要求数据分析按照不同约束条件进行分析；“ods output parameterestimates=regparms”选项要求回归模型结果输出到同一数据集中，并命名该数据集为“miparm1”。 proc print label data=miparm1; var Shift Probt; format Probt 8.4; run; 程序说明：语句“ proc print”表示调用print过程，该语句中的“data=”要求调用数据集miparm1的数据，进行后续分析；“var”选项要求对shift和Prrobt两个变量进行分析；“format”规定变量输出格式，此处要求变量Probt小数点后4位。 3.5 主要分析结果及解释 图6-5显示数据中数据缺失模式，各个变量意义：“Group”显示原始数据的缺失类型包括三种形式；“Freq”显示每个类型的的频次；“Prercent”规定了每个类型占样本总数的百分比；COL1、COL2、COL3分别指示不同缺失类型下，即将填补的三个变量平均数。 图6-6显示为在结果模式数据中不同shift下，每次填补中每组治疗效应(“估计”变量)以及可信区间(“下限”、“上限”变量)。 上表结果显示在规定的shift(−10，10)的范围内，两组间的差值及可信区间，结果显示表明在规定的shift内，两组间依然存在统计学差异，P值均<0.05，结果稳定。
图6-5 数据缺失类型结果
图6-6 不同模式下两组多重填补结果 表6-4 不同shift情况下两组差值及可信区间
附录1 PMM SAS macro code /*----------------------------------------------------------------*/ /*--- Generate imputed data set for specified shift parameters ---*/ /*--- data= input data set ---*/ /*--- smin= min shift parameter ---*/ /*--- smax= max shift parameter ---*/ /*--- sinc= increment of the shift parameter ---*/ /*--- out= output imputed data set ---*/ /*----------------------------------------------------------------*/ %macro midata( data=, smin=, smax=, sinc=, out=); data &out; set _null_; run; /*------------ # of shift values ------------*/ %let ncase= %sysevalf( (&smax-&smin)/&sinc, ceil ); /*------- Imputed data for each shift -------*/ %do jc=0 %to &ncase; %let sj= %sysevalf( &smin + &jc * &sinc); proc mi data=&data seed=14823 nimpute=10 out=outmi; class Trt; monotone reg; mnar adjust( y1 / shift=&sj adjustobs=(Trt='1') ); var Trt y0 y1; run; data outmi; set outmi; Shift= &sj; run; data &out; set &out outmi; run; %end; %mend midata; 活动 8月26日（周六）上午将迎来第一期网络学术直播，届时《傻瓜统计学》、《聪明统计学》人气作者胡志德、周支瑞分别空降广州、上海两地授课并直播。其中，AME广州办、上海办、长沙办、台北办更有现场收看直播活动，名额有限，报名从速。无法到场的读者亦可通过关注【AME直播间】在线轻松观看直播，敬请期待！欢迎点击《8月26日，我们相约云端》了解直播详情。 |
|
|
