陶麗新，劉一松，胡良平

(1.首都醫(yī)科大學(xué)公共衛(wèi)生學(xué)院，北京 100069；2.世界中醫(yī)藥聯(lián)合會(huì)臨床科研統(tǒng)計(jì)學(xué)專業(yè)委員會(huì)，北京 100029；3.北京醫(yī)普科諾科技有限公司，北京 100190；4.軍事醫(yī)學(xué)科學(xué)院生物醫(yī)學(xué)統(tǒng)計(jì)學(xué)咨詢中心，北京 100850*通信作者：胡良平：E-mail: lphu812@sina.com)

?

基于SAS軟件實(shí)現(xiàn)樣本含量估計(jì)及應(yīng)用

陶麗新1,2，劉一松3，胡良平2,4*

(1.首都醫(yī)科大學(xué)公共衛(wèi)生學(xué)院，北京 100069；2.世界中醫(yī)藥聯(lián)合會(huì)臨床科研統(tǒng)計(jì)學(xué)專業(yè)委員會(huì)，北京 100029；3.北京醫(yī)普科諾科技有限公司，北京 100190；4.軍事醫(yī)學(xué)科學(xué)院生物醫(yī)學(xué)統(tǒng)計(jì)學(xué)咨詢中心，北京 100850*通信作者：胡良平：E-mail: lphu812@sina.com)

本文的目的是使讀者能方便快捷地運(yùn)用SAS軟件中的POWER過(guò)程和GLMPOWER過(guò)程實(shí)現(xiàn)樣本含量估計(jì)。在不同的場(chǎng)合下估計(jì)樣本含量需要提供不同的前提條件，即使僅限于假設(shè)檢驗(yàn)時(shí)估計(jì)樣本含量，也必須進(jìn)一步弄清對(duì)應(yīng)的設(shè)計(jì)類型、結(jié)果變量的性質(zhì)、某些先驗(yàn)知識(shí)和對(duì)結(jié)果精確度的要求。本文通過(guò)一些實(shí)例，介紹了估計(jì)樣本含量與檢驗(yàn)效能的SAS實(shí)現(xiàn)方法。讀者只需要修改本文中所呈現(xiàn)的SAS程序中的少量參數(shù)，就可方便地使用SAS軟件實(shí)現(xiàn)樣本含量與檢驗(yàn)效能估計(jì)。事實(shí)說(shuō)明，盡管SAS軟件非常難學(xué)難用，但借助現(xiàn)成的SAS程序，讀者可以輕松自如地解決很多與統(tǒng)計(jì)分析有關(guān)的具體問(wèn)題。

SAS軟件；樣本含量；檢驗(yàn)效能；均值；標(biāo)準(zhǔn)差；率

1 估計(jì)樣本含量的必要性與前提條件

在開展試驗(yàn)性科學(xué)研究中，往往需要在相同試驗(yàn)條件下做多次重復(fù)試驗(yàn)，以便使隨機(jī)變量的變化規(guī)律能充分顯露出來(lái)。文獻(xiàn)[1]就“重復(fù)原則與樣本含量”、“重復(fù)的三層含義”、“估計(jì)樣本含量和檢驗(yàn)效能的意義”、“估計(jì)樣本含量需要的前提條件”等內(nèi)容，作了較為詳細(xì)的介紹，本文將直接運(yùn)用SAS軟件中的POWER過(guò)程和GLMPOWER過(guò)程估計(jì)樣本含量和檢驗(yàn)效能。

2 定量資料假設(shè)檢驗(yàn)時(shí)樣本含量估計(jì)[2-5]

【例1】某研究者欲比較兩種物理療法增加肌肉彈性的效果。希望采用t檢驗(yàn)比較兩總體均值之間差別并希望得出具有統(tǒng)計(jì)學(xué)意義的檢驗(yàn)結(jié)果，設(shè)檢驗(yàn)水準(zhǔn)α為0.05，檢驗(yàn)效能1-β為0.9，假定的均值和標(biāo)準(zhǔn)差的數(shù)值(在實(shí)踐中，這些數(shù)據(jù)應(yīng)基于文獻(xiàn)資料或預(yù)實(shí)驗(yàn)的結(jié)果)將列在下面。擬回答的問(wèn)題是：至少需要檢測(cè)多少例受試對(duì)象？

標(biāo)準(zhǔn)治療方法的肌肉彈性均值為13，如果新療法的肌肉彈性被定位14.0、14.5和15.0三種情況，且假定分別有兩種可能的標(biāo)準(zhǔn)差，即σ=1.2與σ=1.7；假定研究者還希望兩組的樣本含量具有下列三種比例關(guān)系，即：①每個(gè)組的樣本含量相同(1:1)；②新療法的樣本含量是標(biāo)準(zhǔn)療法的兩倍(1:2)；③新療法的樣本含量是標(biāo)準(zhǔn)療法的三倍(1:3)。

【分析與解答】根據(jù)前面給出的條件可知，均值有3種可能的取值、標(biāo)準(zhǔn)差有兩種可能的取值、兩組樣本含量的比例有三種可能的取值，故這樣一共有3×2×3=18種情況。也就是說(shuō)，要在18種條件下估計(jì)所需要的樣本含量。

使用下面的SAS程序，可以獲得與上面18種條件對(duì)應(yīng)的樣本含量估計(jì)結(jié)果：

proc power;

twosamplemeans

groupmeans = (13 14) (13 14.5) (13 15)

stddev = 1.2 1.7

groupweights = 1 | 1 2 3

power = 0.9

ntotal = .;

run;

【程序說(shuō)明】在power過(guò)程中使用twosamplemeans語(yǔ)句來(lái)決定18種情況下為了達(dá)到90%的檢驗(yàn)效能所需要的樣本含量。將選項(xiàng)“ntotal=.”指定樣本含量為結(jié)果變量；用“groupmenas=”選項(xiàng)指定假設(shè)的均值，對(duì)每種情況的均值使用匹配標(biāo)記，即用圓括號(hào)隔開不同的均值組合條件；使用“stddev=”選項(xiàng)為每種情況指定標(biāo)準(zhǔn)差；使用“groupweights=”選項(xiàng)指定權(quán)重比例，可以再使用匹配標(biāo)記，即使用“|”符號(hào)隔開不同的取值。

【輸出結(jié)果】

Computed N Total

【結(jié)果解讀】最好的設(shè)計(jì)類型是(均值差最大為2，小的標(biāo)準(zhǔn)差為1.2和平衡設(shè)計(jì))，達(dá)到至少0.9的檢驗(yàn)效能需要的樣本含量為18(n1=n2=9)；最差的設(shè)計(jì)為(小的均值差為1，大的標(biāo)準(zhǔn)差為1.7，1:3的非平衡設(shè)計(jì))，達(dá)到至少0.9的檢驗(yàn)效能需要樣本含量為164(n1=41、n2=123)。在“固定設(shè)計(jì)中”指定的檢驗(yàn)效能為0.9，在“計(jì)算的樣本含量”表格中輸出根據(jù)樣本比例調(diào)整的實(shí)際的檢驗(yàn)效能(見上面輸出結(jié)果中倒數(shù)第2列)。

值得一提的是，在上面的SAS程序中，若給定總樣本含量，而將power設(shè)置為“缺省值”，就是在各種不同條件下，估計(jì)檢驗(yàn)效能。如將上面的SAS程序修改如下：

procpower; twosamplemeans groupmeans=(1314)(1314.5)(1315) stddev=1.21.7 groupweights=1|123 power=. ntotal=100;run;

讀者可以自己運(yùn)行上面這段SAS程序，會(huì)得到各種條件下的檢驗(yàn)效能估計(jì)值，因篇幅所限，此輸出結(jié)果從略。在下面的其他SAS程序中，也可進(jìn)行同樣的修改。

【分析與解答】這是一個(gè)結(jié)果變量為定量變量的單因素三水平設(shè)計(jì)試驗(yàn)研究的樣本含量估計(jì)問(wèn)題。所需的SAS程序如下。

odshtml; procpower; onewayanovagroupmeans=290|658|763stddev=174155127 groupweights=(111)alpha=0.05Ntotal=.power=0.8; run;quit; odshtmlclose;

【程序說(shuō)明】選項(xiàng)“groupmeans =290 | 658 | 763”定義各水平組指標(biāo)總體均值的估計(jì)值；選項(xiàng)“stddev = 174 155 127”定義幾種總體標(biāo)準(zhǔn)差的估計(jì)值；選項(xiàng)“groupweights = (1 1 1)”定義各水平組樣本含量的分配比例。對(duì)于類似的問(wèn)題，上述程序中的數(shù)字需要根據(jù)用戶的實(shí)際情況進(jìn)行修改。

【輸出結(jié)果】

ComputedNTotalIndexStdDevActualPowerNTotal11740.8661221550.9331231270.9119

【結(jié)果解讀】分析結(jié)果表明：本研究各類人群只需隨機(jī)抽取4例，3類人群共12例即可滿足要求。

【例3】用三種方法治療腦卒中抑郁患者，觀察神經(jīng)功能康復(fù)狀況。估計(jì)治療后三種方法斯堪地那維亞卒中量表(Scemelinavian Stroke Scale,SSS)評(píng)分均數(shù)分別為11.0、23.0、9.0，標(biāo)準(zhǔn)差分別為3、3、2。希望各組樣本含量的比例或權(quán)重(本例取1:1:1)、各組的平均標(biāo)準(zhǔn)差(保守的做法，可采用所有組的最大標(biāo)準(zhǔn)差)、允許犯第Ⅰ類(或假陽(yáng)性)錯(cuò)誤的概率α(通常取其值為0.05)、給定檢驗(yàn)效能1-β的值[相當(dāng)于給定允許犯第Ⅱ類(或假陰性)錯(cuò)誤的概率β(通常取其值為0.2或0.1)，本例取1-β=0.9]。如果要得到三組之間的差別以及兩兩比較的結(jié)果都有統(tǒng)計(jì)學(xué)意義的結(jié)論，問(wèn)每組各需要多少例患者？

【分析與解答】 這是做單因素三水平設(shè)計(jì)一元定量資料均值比較及兩兩比較之前進(jìn)行樣本含量估計(jì)的問(wèn)題。所需的SAS程序如下。

data d11_20;

input method $ score CellWgt;

datalines;

A 11.0 1

B 23.0 1

C 9.0 1

;

run;

ods html;

PROC GLMPOWER data=d11_20;

class method;

model score=method;

weight CellWgt;

contrast "A vs. B" method 1 -1 0;

contrast "A vs. C" method 1 0 -1;

contrast "B vs. C" method 0 1 -1;

POWER stddev = 3.0 alpha = 0.05 power = 0.90 ntotal = .；

run; quit;

ods html close;

【程序說(shuō)明】 data步定義影響因素各水平組指標(biāo)的總體均值的估計(jì)值及各水平組的樣本含量分配比例；GLMPOWER過(guò)程中的選項(xiàng)“stddev=3.0”中的3.0為各水平組指標(biāo)的總體標(biāo)準(zhǔn)差的估計(jì)值，本例取三個(gè)水平組總體標(biāo)準(zhǔn)差估計(jì)值的最大值帶入計(jì)算，也可以同時(shí)考慮多個(gè)總體標(biāo)準(zhǔn)差的估計(jì)值，將選項(xiàng)寫成“stddev=2.0 2.5 3.0”或“stddev=2.0 to 3.0 by 0.5”的形式；程序中的數(shù)字需要根據(jù)實(shí)際情況進(jìn)行修改。

具體地說(shuō)，在數(shù)據(jù)步中：method實(shí)際上指“試驗(yàn)因素”，若用factor更直觀，它是一個(gè)字符型變量，其具體取值或水平可用一個(gè)字母或一個(gè)單詞、甚至可用一句話(當(dāng)字符個(gè)數(shù)超過(guò)8時(shí)，不能僅用一個(gè)美元符號(hào)，如最長(zhǎng)為20個(gè)字符，可寫成：method $20.)來(lái)表示；“score”指各組的算術(shù)平均值，若用“mean”更直觀；“CellWgt”指各組的樣本含量比例或權(quán)重，本例取1:1:1，根據(jù)實(shí)際需要，也可取其它比例關(guān)系，如1:2:3；

在過(guò)程步中，有三個(gè)Contrast語(yǔ)句，相當(dāng)于單因素三水平之間的兩兩比較，分別表示A與B比較、A與C比較、B與C比較；POWER語(yǔ)句中寫了4個(gè)內(nèi)容，可以寫在同一行上，寫成4行是為了看起來(lái)更清楚一些，前3行代表提供的3個(gè)已知條件，它們分別為平均的或最大的樣本標(biāo)準(zhǔn)差、α值和1-β的值；而“ntotal = .”表示需要估計(jì)總樣本含量，即三組樣本含量之總和。

【輸出結(jié)果】

ComputedNTotalIndexTypeSourceTestDFErrorDFActualPowerNTotal1Effectmethod260.98992ContrastAvs.B160.98193ContrastAvs.C11410.9001444ContrastBvs.C160.9979

【結(jié)果解讀】由計(jì)算結(jié)果可知：若要求通過(guò)單因素三水平設(shè)計(jì)一元定量資料的方差分析及兩兩比較有90%的把握發(fā)現(xiàn)三個(gè)組SSS評(píng)分的總體均值不完全相等，則每組只需3例患者，三組共需9例患者；但若要求通過(guò)單因素三水平設(shè)計(jì)一元定量資料的方差分析及兩兩比較有90%的把握發(fā)現(xiàn)三個(gè)組中任何兩組的SSS評(píng)分的總體均值不相等，則每組需要48例患者，三組共需144例患者。

為什么會(huì)出現(xiàn)這樣大的差異呢？問(wèn)題的癥結(jié)在于A組與C組的均值過(guò)于接近，它們之間的數(shù)值之差量為11-9=2，要想把總體上存在如此小的“差距”都能發(fā)現(xiàn)出來(lái)，必須具有非常大的樣本含量。

試想：若A組的均值不是11，而是16，可以計(jì)算得出：每組只需要5例，三組共需要15例即可。程序修改之后的輸出結(jié)果如下：

ComputedNTotalIndexTypeSourceTestDFErrorDFActualPowerNTotal1Effectmethod260.97692ContrastAvs.B1120.922153ContrastAvs.C1120.922154ContrastBvs.C160.9979

【分析與解答】 這是做具有一個(gè)重復(fù)測(cè)量因素的兩因素設(shè)計(jì)一元定量資料均值比較之前進(jìn)行樣本含量估計(jì)的問(wèn)題。所需的SAS程序如下。

%MACROcanshu; alpha=0.05;beta=0.10;k=5;pc=0.75;delta=10;deltae=190.641;deltau=457.333;%MENDcanshu;%MACROrepeatmesure(a,b,c,d);FILEPRINT;u005=put((-probit(alpha/2)),8.3);u02=put((-probit(beta)),8.3);n=int((2*(deltau+(1+(k-1)*pc)*deltae/k)*(u005+u02)*(u005+u02)/(delta*delta)))+1;PUT#2@15‘兩組各需’n‘例。’;RUN;%MENDrepeatmesure;DATA;%canshu;%repeatmesure(a,b,c,d);RUN;

【程序說(shuō)明】第2行中“alpha=0.05”中的0.05為可能犯第Ⅰ類錯(cuò)誤的概率，“beta=0.10”中的0.10為可能犯第Ⅱ類錯(cuò)誤的概率，“k=5”中的5為重復(fù)測(cè)量次數(shù)，“pc=0.75”中的0.75為條件相關(guān)系數(shù)，“delta=10”中的10為兩組間可以識(shí)別的最小差異，“deltae=190.641”中的190.641為重復(fù)測(cè)量誤差，“deltau=457.333”中的457.333為個(gè)體間差異的方差。

【輸出結(jié)果】每組需要129例。

3 生存分析時(shí)樣本含量估計(jì)

【例5】某臨床試驗(yàn)擬比較TIPS手術(shù)與脊髓改道手術(shù)延長(zhǎng)出血性食管曲張患者的生存時(shí)間有無(wú)差別，試驗(yàn)組為TIPS手術(shù)組，對(duì)照組為脊髓改道手術(shù)組，隨訪期至少1年。根據(jù)既往研究可知對(duì)照組總體1年生存率是45%，預(yù)期試驗(yàn)組總體1年生存率將達(dá)到65%，采用Log-Rank檢驗(yàn)比較兩組生存率，試估計(jì)樣本含量(采用雙側(cè)檢驗(yàn)，檢驗(yàn)水準(zhǔn)α=0.05，檢驗(yàn)效能1-β=0.85)。

【分析與解答】 對(duì)照組生存率為45%、試驗(yàn)組生存率為65%，對(duì)照組和試驗(yàn)組的終點(diǎn)事件發(fā)生率分別為：π1=1-0.65=0.35，πc=1-0.45=0.55。這是對(duì)生存資料分析之前的樣本含量估計(jì)問(wèn)題。所需的SAS程序如下。

%letpt=0.35;%letpc=0.55;%letalpha=0.05;%letbeta=0.15; dataa6_1; theta=log10(1-&pt)/log10(1-&pc); theta1=round(theta,0.001); e=((theta+1)/(theta-1))**2*(probit(1-&alpha/2)+probit(1-&beta))**2; e1=int(e)+1; n=2/(&pt+&pc)*e; n1=int(n)+1; fileprint; PUT#3@15‘風(fēng)險(xiǎn)比為‘theta1’;兩組終點(diǎn)事件發(fā)生總例數(shù)不小于‘e1’;兩組所需總的樣本含量為‘n1’?！? run;

【程序說(shuō)明】程序當(dāng)中pt、pc、alpha、beta分別代表試驗(yàn)組終點(diǎn)事件發(fā)生率πt、對(duì)照組終點(diǎn)事件發(fā)生率πc、檢驗(yàn)水準(zhǔn)α、犯假陰性錯(cuò)誤的概率β。

【輸出結(jié)果】風(fēng)險(xiǎn)比為0.539；兩組終點(diǎn)事件發(fā)生總例數(shù)不小于101；兩組所需總的樣本含量為223。

說(shuō)明：樣本含量與檢驗(yàn)效能估計(jì)的內(nèi)容非常豐富，本文的目的僅僅希望起到拋磚引玉的作用。這方面的統(tǒng)計(jì)軟件和書籍[6-7]很多，讀者需要時(shí)，還可參閱有關(guān)文獻(xiàn)。

[1] 張效嘉，胡良平. 精神衛(wèi)生科研如何嚴(yán)格遵守試驗(yàn)設(shè)計(jì)四原則之重復(fù)原則[J]. 四川精神衛(wèi)生，2016, 29(4): 303-306.

[2] SAS Institute Inc. SAS/STAT 9.3 User’s Guide[M]. Cary, NC: SAS Institute Inc, 2011: 3361-3400,5729-5962.

[3] 胡良平. 統(tǒng)計(jì)學(xué)三型理論在實(shí)驗(yàn)設(shè)計(jì)中的應(yīng)用[M]. 北京: 人民軍醫(yī)出版社, 2006: 215-265.

[4] 胡良平. 科研設(shè)計(jì)與統(tǒng)計(jì)分析[M]. 北京: 軍事醫(yī)學(xué)科學(xué)出版社, 2012: 136-227.

[5] 胡良平. SAS常用統(tǒng)計(jì)分析教程[M]. 北京: 電子工業(yè)出版社, 2015: 235-252.

[6] Shein-Chung Chow, Jun Shao, Hansheng Wang. Sample Size Calculations in Clinical Research[M]. 2nd Edition. Chapman & Hall/CRC, 2008: 49-408.

[7] Thomas P. Ryan. Sample Size Determination and Power[M]. Wiley, 2013: 17-352.

(本文編輯：吳俊林)

Sample size estimation and its application by using SAS software

TAOLi-Xin1,2,LIUYi-song3,HULiang-ping2,4*

(1.PublicHealthAcademy,CapitalUniversityofMedicalSciences,Beijing100069,China;2.SpecialtyCommitteeofClinicalScientificResearchStatisticsofWorldFederationofChineseMedicineSocieties,Beijing100029,China;3.BeijingImprove-QualityTechCo.,Ltd.,Beijing100190,China;4.ConsultingCenterofBiomedicalStatistics,AcademyofMilitaryMedicalSciences,Beijing100850,China*Correspondingauthor:HULiang-ping,E-mail:lphu812@sina.com)

The paper aims to help the readers to estimate the sample size by using the POWER and GLMPOWER procedures in SAS software conveniently and easily. It is necessary for a user to supply the prerequisites of sample size estimation under the different situations, even if limited under the case of the hypothesis testing. Usually the premise conditions of sample size estimation, which should be given, are as follows: the design type, the characters of outcomes, the prior knowledge, and the accuracy of the results. The methods of estimating the sample size and power by SAS were introduced through several real examples in this article. It is convenient for the readers to realize the estimation of sample size and power by using SAS. The only thing for them to do is to modify a few parameters in the given SAS programs. Although SAS software is very difficult for people to learn and use, as a matter of fact, the users can solve the concrete problems concerned with statistical analysis by means of revisions of ready-made SAS programs easily.

SAS software; Sample size; Power; Mean value; Standard deviation; Rate

R195.1

A doi:10.11886/j.issn.1007-3256.2016.05.004

2016-10-11)