王雨萌 孙瑞华 黄 傲 徐 凯 李 欢
·计算机应用·
单因素多水平临床试验定量指标统计分析报表的SAS宏实现
王雨萌1孙瑞华2△黄傲1徐凯1李欢1
SAS软件输出的统计结果多而复杂,常需要人为将统计结果复制粘贴到统计报告相应位置上,这一过程耗时耗力且容易出错。为此,可根据需要编制相应的SAS宏程序,直接生成针对各种类型资料的统计分析报表。以往对于成组计量或计数资料的统计分析报表的SAS宏程序已经有了一定的探讨[1-2],本文主要介绍单因素多水平临床试验定量指标统计分析报表自动实现的SAS宏程序。
根据统计学原理,对于符合正态分布的定量指标,通常采用均数、标准差、中位数、最小值、最大值以及95%置信区间来描述其集中和离散趋势,而非正态分布的定量指标,则采用均数、中位数、最小值、最大值以及四分位间距来进行统计描述。
对于多水平定量指标的组间比较,若满足方差分析的前提条件,则采用方差分析处理资料;若指标满足正态性但不满足方差齐性要求,则选用Welch近似方差分析,此时Welch近似方差分析的结果较方差分析的结果更稳定;当指标既不满足正态性也不满足方差齐性要求时,则选用Kruskal-Wallis秩和检验。方差分析与Welch近似方差分析的统计量用F_value表示,Kruskal-Wallis秩和检验的统计量用H表示。为直接得到如表1所示的统计报表,编写如下三组定量指标统计分析的SAS宏程序。
表1 三组基线一般资料的组间比较(FAS)
在该宏程序中设置四个宏变量&database、&var、&index和&datasty,分别表示数据库、统计量、指标项和统计分析数据集,并假定分组变量为group。
%macro sanzujl(database=,var=,index=,datasty=);
利用output语句将描述性统计量和正态检验的结果输出到SAS指定的数据集中,对该数据集进行拆分后再横向合并,使三组的描述性统计量和正态检验的结果显示在一行,方便后面的调用和输出。
/*利用univariate程序输出部分描述性结果到result1数据集*/
proc univariate noprint;
var &var;
by group;
where &datasty=1;
output out=result1 n=n mean=mean median=med std=std min=min max=max nmiss=nmiss qrange=qrange PROBn=Pnor;
/*利用means程序输出95%置信区间到result2数据集*/
proc means clm noprint data=&database;
var &var;
by group;
where &datasty=1;
output out=result2(drop=_TYPE_ _FREQ_)uclm=uclm lclm=lclm;
run;
/*合并result1和result2为数据集result*/
data result;
merge result1 result2;
by group;
run;
/*将三组描述性结果分别输出到独立的数据集中,并定义变量的长度*/
data a(where=(group=“A”))b(where=(group=“B”))c(where=(group=“C”));
set result;
mean=round(mean,0.01);
med=round(med,0.01);
std=round(std,0.01);
Pnor=round(Pnor,0.01);
lclm=round(lclm,0.01);
uclm=round(uclm,0.01);
qrange=round(qrange,0.01);
min=round(min,0.01);
max=round(max,0.01);
run;
/*将三个独立的数据集分别重新命变量名,再合并到result数据集中*/
data result(drop=group);
merge a(rename=(n=NA nmiss=MA mean=MeanA med=medA std=StdA qrange=qrangeA lclm=lclmA uclm=uclmA Min=MinA Max=MaxA Pnor=PnorA))
b(rename=(n=NB nmiss=MB mean=MeanB med=medB std=StdB qrange=qrangeB lclm=lclmB uclm=uclmB Min=MinB Max=MaxB Pnor=PnorB))
c(rename=(n=NC nmiss=MC mean=MeanC med=medC std=StdC qrange=qrangeC lclm=lclmC uclm=uclmC Min=MinC Max=MaxC Pnor=PnorC));
run;
利用SAS的ODS功能将方差分析、Welch近似方差分析以及方差齐性检验的结果分别输出到指定的数据集中,将所有有用变量横向合并到一个数据集。
/*利用anova程序输出方差分析结果到ModelANOVA数据集、输出方差齐性检验结果到HOVFTest1数据集、输出welch近似方差分析结果到Welch1数据集*/
ods listing close;
ods output ModelANOVA=ModelANOVA HOVFTest=HOVFTest1 Welch=Welch1;
proc anova data=&database;
where &datasty=1;
class group;
model &var=group;
means group/hovtest welch;
run;
ods listing;
/*方差齐性检验保留有用变量*/
data HOVFTest(keep=FValue ProbF);
set HOVFTest1;
if _n_ ^=1 then delete;
run;
/*welch近似方差分析保留有用变量*/
data Welch(keep=FValue ProbF);
set Welch1;
if _n_ ^=1 then delete;
run;
/*将方差分析、Welch近似方差分析、方差齐性检验的结果合并到数据集test中,保留有用变量,添加方差分析与Welch近似方差分析统计量标签*/
data test(keep=FValue ProbF stat FValueH ProbFH FValueW ProbFW);
merge ModelANOVA(keep=FValue ProbF)HOVFTest(rename=(FValue=FValueH ProbF=ProbFH))Welch(rename=(FValue=FValueW ProbF=ProbFW));
stat=“F_value”;
run;
利用SAS的ODS功能将秩和检验的结果输出到指定的数据集,将统计量和P值拆分后横向合并到一个数据集。
/*利用npar1way wilcoxon程序输出秩和结果到KruskalWallisTest数据集*/
ods listing close;
Ods Output KruskalWallisTest=KruskalWallisTest;
proc npar1way wilcoxon data=&database;
where &datasty=1;
var &var;
class group;
run;
ods listing;
/*只保留统计量和P值*/
data KruskalWallisTest1;
set KruskalWallisTest;
if Label1=“DF” then delete;
run;
/*将统计量和P值分别独立保存*/
data d(where=(Name1=“_KW_”))e(where=(Name1=“P_KW”));
set KruskalWallisTest1;
run;
/*将秩和检验的统计量和P值合并到kwttest数据集,并定义变量的长度和统计量的标签*/
data kwttest(keep=cValue1d cValue1e statkw);
merge d(rename=(Name1=Name1d cValue1=cValue1d nValue1=nValue1d))e(rename=(Name1=Name1e cValue1=cValue1e nValue1=nValue1e));
cValue1d=round(cValue1d,0.01);
label cValue1d=“Chi-Square”;
label cValue1e=“Pr > Chi-Square”;
statkw=“H”;
run;
将描述性统计量、正态检验结果、方差分析、Welch近似方差分析以及秩和检验结果的数据集横向合并,生成打印的数据集,从该数据集中输出特定条件下的描述性统计量和P值。
/*将数据集result、test、kwttest合并*/
data _null_;
merge result test kwttest;
/*如果三组均正态且方差齐,用方差分析*/
if PnorA>0.05 and PnorB>0.05 and PnorC>0.05 and ProbFH>0.05 then do;
file print notitle;
put #1 @3 “&index”
#2 @5 “例数(缺失)” @22 NA‘(‘MA’)’ @38 NB‘(‘MB’)’ @55 NC‘(‘MC’)’ @70 stat’=‘FValue @88 ProbF
#3 @5 “均数±标准差” @22 MeanA‘±’StdA @38 MeanB‘±’StdB @55 MeanC‘±’StdC
#4 @5 “中位数” @22 MedA @38 MedB @55 MedC
#5 @5 “最小值-最大值” @22 MinA‘-’MaxA @38 MinB‘-’MaxB @55 MinC‘-’MaxC
#6 @5 ‘95%CI’ @22 lclmA‘-’uclmA @38 lclmB‘-’uclmB @55 lclmC‘-’uclmC;
end;
/*如果三组均正态但方差不齐,用welch方差分析*/
else if PnorA>0.05 and PnorB>0.05 and PnorC>0.05 and ProbFH<0.05 then do;
file print notitle;
put #1 @3 “&index”
#2 @5 “例数(缺失)”@22 NA‘(‘MA’)’ @38 NB‘(‘MB’)’ @55 NC‘(‘MC’)’ @70 stat’=‘FValueW @88 ProbFW
#3 @5 “均数±标准差” @22 MeanA‘±’StdA @37 MeanB‘±’StdB @55 MeanC‘±’StdC
#4 @5 “中位数” @22 MedA @38 MedB @55 MedC
#5 @5 “最小值-最大值” @22 MinA‘-’MaxA @38 MinB‘-’MaxB @55 MinC‘-’MaxC
#6 @5 ‘95%CI’ @22 lclmA‘-’uclmA @38 lclmB‘-’uclmB @55 lclmC‘-’uclmC;
end;
/*如果三组有一组非正态,用秩和检验*/
else if PnorA<0.05 or PnorB<0.05 or PnorC<0.05 then do;
file print notitle;
put #1 @3 “&index”
#2 @5 “例数(缺失)”@22 NA‘(‘MA’)’ @38 NB‘(‘MB’)’ @55 NC‘(‘MC’)’ @70 statkw’=‘cValue1d @88 cValue1e
#3 @5 “均数” @22 MeanA @38 MeanB @55 MeanC
#4 @5 “中位数” @22 MedA @38 MedB @55 MedC
#5 @5 “最小值-最大值” @22 MinA‘-’MaxA @38 MinB‘-’MaxB @55 MinC‘-’MaxC
#6 @5 “四分位间距” @22 qrangeA @38 qrangeB @55 qrangeC;
end;
run;
%mend sanzujl;
最后,打印出表头及表格线。
%macro ctit(tit);
data _null_;
file print n=ps notitles;put #1 @5 “&tit”
#2 @2 94*‘_’
#3 @5 ‘指标项’ @22 ‘A组’ @38 ‘B组’ @55 ‘C组’ @70 ‘统计量’ @88 ‘P值’
#4 @2 94*‘_’;
run;
%mend ctit;
%macro cleg;
data _null_;
file print n=ps notitles;
put #1 @2 94*‘_’;
run;
%mend cleg;
将以上所有SAS宏程序提交SAS系统运行,即可完成如表1所示的统计分析表格。
/*其中,knfps为数据库,设定数据集为FAS数据集,要分析比较的变量是weigh和high,分别表示体重和身高。*/
%ctit(表1.三组基线一般资料组间比较(FAS));
%sanzujl(database=knfps,var=weigh,datasty=FAS,index=体重);
%sanzujl(database=knfps,var=high,datasty=FAS,index=身高);
%cleg;
目前无论美国FDA和国内SFDA都要求采用世界公认的统计分析软件SAS来进行统计分析,但是SAS软件自身的统计分析结果比较丰富,通常只需要将其中一小部分关键的结果按照一定格式制作成表格[3]。尽管SAS提供了如proc report、proc template等制表的功能,但定义过程较为复杂,需要对SAS软件有一定的掌握程度。本文介绍了一段简单实用且可操作性非常强的SAS宏程序供广大读者调用,其运行结果经验证准确而可靠,适合应用于单因素多水平定量指标的组间比较。
[1]邹建东,熊宁宁,卜擎燕,等.正态分布定量指标统计分析报表的SAS宏程序.中国临床药理学与治疗学,2004,9(7):838-840.
[2]邹建东,熊宁宁,卜擎燕,等.四格表指标统计分析报表的SAS宏程序.中国临床药理学与治疗学,2005,10(3):357-360.
[3]殷红.临床试验中统计分析报告自动化生成的研究与应用.复旦大学,2009.
(责任编辑:刘壮)
孙瑞华,E-mail:sunruihua@263.net
1.北京中医药大学管理学院(100029)
2.中日友好医院科研处