王雨萌 孫瑞華 黃 傲 徐 凱 李 歡
?
·計算機應用·
單因素多水平臨床試驗定量指標統(tǒng)計分析報表的SAS宏實現(xiàn)
王雨萌1孫瑞華2△黃傲1徐凱1李歡1
SAS軟件輸出的統(tǒng)計結(jié)果多而復雜,常需要人為將統(tǒng)計結(jié)果復制粘貼到統(tǒng)計報告相應位置上,這一過程耗時耗力且容易出錯。為此,可根據(jù)需要編制相應的SAS宏程序,直接生成針對各種類型資料的統(tǒng)計分析報表。以往對于成組計量或計數(shù)資料的統(tǒng)計分析報表的SAS宏程序已經(jīng)有了一定的探討[1-2],本文主要介紹單因素多水平臨床試驗定量指標統(tǒng)計分析報表自動實現(xiàn)的SAS宏程序。
根據(jù)統(tǒng)計學原理,對于符合正態(tài)分布的定量指標,通常采用均數(shù)、標準差、中位數(shù)、最小值、最大值以及95%置信區(qū)間來描述其集中和離散趨勢,而非正態(tài)分布的定量指標,則采用均數(shù)、中位數(shù)、最小值、最大值以及四分位間距來進行統(tǒng)計描述。
對于多水平定量指標的組間比較,若滿足方差分析的前提條件,則采用方差分析處理資料;若指標滿足正態(tài)性但不滿足方差齊性要求,則選用Welch近似方差分析,此時Welch近似方差分析的結(jié)果較方差分析的結(jié)果更穩(wěn)定;當指標既不滿足正態(tài)性也不滿足方差齊性要求時,則選用Kruskal-Wallis秩和檢驗。方差分析與Welch近似方差分析的統(tǒng)計量用F_value表示,Kruskal-Wallis秩和檢驗的統(tǒng)計量用H表示。為直接得到如表1所示的統(tǒng)計報表,編寫如下三組定量指標統(tǒng)計分析的SAS宏程序。
表1 三組基線一般資料的組間比較(FAS)
在該宏程序中設(shè)置四個宏變量&database、&var、&index和&datasty,分別表示數(shù)據(jù)庫、統(tǒng)計量、指標項和統(tǒng)計分析數(shù)據(jù)集,并假定分組變量為group。
%macro sanzujl(database=,var=,index=,datasty=);
利用output語句將描述性統(tǒng)計量和正態(tài)檢驗的結(jié)果輸出到SAS指定的數(shù)據(jù)集中,對該數(shù)據(jù)集進行拆分后再橫向合并,使三組的描述性統(tǒng)計量和正態(tài)檢驗的結(jié)果顯示在一行,方便后面的調(diào)用和輸出。
/*利用univariate程序輸出部分描述性結(jié)果到result1數(shù)據(jù)集*/
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%置信區(qū)間到result2數(shù)據(jù)集*/
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為數(shù)據(jù)集result*/
data result;
merge result1 result2;
by group;
run;
/*將三組描述性結(jié)果分別輸出到獨立的數(shù)據(jù)集中,并定義變量的長度*/
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;
/*將三個獨立的數(shù)據(jù)集分別重新命變量名,再合并到result數(shù)據(jù)集中*/
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近似方差分析以及方差齊性檢驗的結(jié)果分別輸出到指定的數(shù)據(jù)集中,將所有有用變量橫向合并到一個數(shù)據(jù)集。
/*利用anova程序輸出方差分析結(jié)果到ModelANOVA數(shù)據(jù)集、輸出方差齊性檢驗結(jié)果到HOVFTest1數(shù)據(jù)集、輸出welch近似方差分析結(jié)果到Welch1數(shù)據(jù)集*/
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近似方差分析、方差齊性檢驗的結(jié)果合并到數(shù)據(jù)集test中,保留有用變量,添加方差分析與Welch近似方差分析統(tǒng)計量標簽*/
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功能將秩和檢驗的結(jié)果輸出到指定的數(shù)據(jù)集,將統(tǒng)計量和P值拆分后橫向合并到一個數(shù)據(jù)集。
/*利用npar1way wilcoxon程序輸出秩和結(jié)果到KruskalWallisTest數(shù)據(jù)集*/
ods listing close;
Ods Output KruskalWallisTest=KruskalWallisTest;
proc npar1way wilcoxon data=&database;
where &datasty=1;
var &var;
class group;
run;
ods listing;
/*只保留統(tǒng)計量和P值*/
data KruskalWallisTest1;
set KruskalWallisTest;
if Label1=“DF” then delete;
run;
/*將統(tǒng)計量和P值分別獨立保存*/
data d(where=(Name1=“_KW_”))e(where=(Name1=“P_KW”));
set KruskalWallisTest1;
run;
/*將秩和檢驗的統(tǒng)計量和P值合并到kwttest數(shù)據(jù)集,并定義變量的長度和統(tǒng)計量的標簽*/
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;
將描述性統(tǒng)計量、正態(tài)檢驗結(jié)果、方差分析、Welch近似方差分析以及秩和檢驗結(jié)果的數(shù)據(jù)集橫向合并,生成打印的數(shù)據(jù)集,從該數(shù)據(jù)集中輸出特定條件下的描述性統(tǒng)計量和P值。
/*將數(shù)據(jù)集result、test、kwttest合并*/
data _null_;
merge result test kwttest;
/*如果三組均正態(tài)且方差齊,用方差分析*/
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 “例數(shù)(缺失)” @22 NA‘(‘MA’)’ @38 NB‘(‘MB’)’ @55 NC‘(‘MC’)’ @70 stat’=‘FValue @88 ProbF
#3 @5 “均數(shù)±標準差” @22 MeanA‘±’StdA @38 MeanB‘±’StdB @55 MeanC‘±’StdC
#4 @5 “中位數(shù)” @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;
/*如果三組均正態(tài)但方差不齊,用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 “例數(shù)(缺失)”@22 NA‘(‘MA’)’ @38 NB‘(‘MB’)’ @55 NC‘(‘MC’)’ @70 stat’=‘FValueW @88 ProbFW
#3 @5 “均數(shù)±標準差” @22 MeanA‘±’StdA @37 MeanB‘±’StdB @55 MeanC‘±’StdC
#4 @5 “中位數(shù)” @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;
/*如果三組有一組非正態(tài),用秩和檢驗*/
else if PnorA<0.05 or PnorB<0.05 or PnorC<0.05 then do;
file print notitle;
put #1 @3 “&index”
#2 @5 “例數(shù)(缺失)”@22 NA‘(‘MA’)’ @38 NB‘(‘MB’)’ @55 NC‘(‘MC’)’ @70 statkw’=‘cValue1d @88 cValue1e
#3 @5 “均數(shù)” @22 MeanA @38 MeanB @55 MeanC
#4 @5 “中位數(shù)” @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 ‘統(tǒng)計量’ @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系統(tǒng)運行,即可完成如表1所示的統(tǒng)計分析表格。
/*其中,knfps為數(shù)據(jù)庫,設(shè)定數(shù)據(jù)集為FAS數(shù)據(jù)集,要分析比較的變量是weigh和high,分別表示體重和身高。*/
%ctit(表1.三組基線一般資料組間比較(FAS));
%sanzujl(database=knfps,var=weigh,datasty=FAS,index=體重);
%sanzujl(database=knfps,var=high,datasty=FAS,index=身高);
%cleg;
目前無論美國FDA和國內(nèi)SFDA都要求采用世界公認的統(tǒng)計分析軟件SAS來進行統(tǒng)計分析,但是SAS軟件自身的統(tǒng)計分析結(jié)果比較豐富,通常只需要將其中一小部分關(guān)鍵的結(jié)果按照一定格式制作成表格[3]。盡管SAS提供了如proc report、proc template等制表的功能,但定義過程較為復雜,需要對SAS軟件有一定的掌握程度。本文介紹了一段簡單實用且可操作性非常強的SAS宏程序供廣大讀者調(diào)用,其運行結(jié)果經(jīng)驗證準確而可靠,適合應用于單因素多水平定量指標的組間比較。
[1]鄒建東,熊寧寧,卜擎燕,等.正態(tài)分布定量指標統(tǒng)計分析報表的SAS宏程序.中國臨床藥理學與治療學,2004,9(7):838-840.
[2]鄒建東,熊寧寧,卜擎燕,等.四格表指標統(tǒng)計分析報表的SAS宏程序.中國臨床藥理學與治療學,2005,10(3):357-360.
[3]殷紅.臨床試驗中統(tǒng)計分析報告自動化生成的研究與應用.復旦大學,2009.
(責任編輯:劉壯)
孫瑞華,E-mail:sunruihua@263.net
1.北京中醫(yī)藥大學管理學院(100029)
2.中日友好醫(yī)院科研處