基于指數分布的生存資料盲態下樣本量再估計方法*

朱　典　王　陵　李嬋娟　蔣志偉　張　威　王治東　吳克堅　李　晨　夏結來△

【提　要】　目的　針對生存分析樣本量再估計的適應性設計臨床試驗,本文基于指數分布,提出盲態下生存資料的樣本量再估計方法。方法　采用蒙特卡洛模擬方法,預設4個期中分析時間點,以參數初始估計值M10為均數暫定搜索范圍為M10±0.5M10,在該搜索范圍內產生指數分布隨機數對期中分析時的截尾生存數據進行填補,對填補后的數據進行極大似然估計,得到試驗參數的再估計值若落在搜索范圍外,則需更改搜索范圍以包含并重新進行搜索直至落在搜索范圍內,此時的即為試驗參數的再估計值,在此參數的基礎上重新估計樣本量,并比較4個期中分析點的樣本量再估計結果,確定1個最合適的期中分析時間點。 結果　經過期中分析調整后的樣本量接近真實值,且期中分析時間點越向后移,樣本量估計結果越接近真實值,變異越小。 結論　建議在入組結束并完成1/4最短隨訪時間時進行一次期中分析重新估計樣本量,根據估計結果考慮是否增加樣本量。

傳統的臨床試驗在試驗開展前,研究者根據前期研究結果并結合相關文獻設定試驗關鍵參數,對樣本量進行估計,以期望達到設定的檢驗效能。但是,試驗開始前對參數的設定可能存在較大偏差,以此為基礎估計的樣本量有可能不足以達到預設的檢驗效能。而在適應性設計中允許在試驗開始之后根據已得到的部分試驗數據進行期中分析,重新估計樣本量,以使試驗達到預設的檢驗效能[1]。

對于不同類型的數據,需采用不同的樣本量再估計方法。對于正態資料,Gould[2]等人提出一種基于EM算法的盲態下樣本量再估計方法,在期中分析時對試驗數據的方差進行重新估計。對于二分類變量,Shih[3]等人提出一種基于分層理念的盲態下樣本量再估計方法,在期中分析時對兩組的事件發生率重新估計;Friede[4]等人提出一種基于負二項分布的盲態下樣本量再估計方法,在期中分析時以負二項分布為基礎建模,對總體率進行重新估計。對于生存資料,Togo[5]等人提出一種基于偏似然估計的樣本量再估計方法,在期中分析時對數據揭盲,并根據構建的模型對兩組之間的風險比(HR)進行重新估計;Todd[6]等人提出一種基于外推法的盲態下生存數據樣本量再估計方法,假設兩組之間的風險比(HR)不變,試驗進展過程中在多個預設時間點上進行多次期中分析,并根據多次期中分析結果求效應的平均差別,推斷試驗結束時的生存率。

期中分析時對試驗揭盲會導致Ⅰ型錯誤的膨脹,從而影響試驗的完整性。而在盲態下進行期中分析,對Ⅰ型錯誤的影響可以忽略,無需對檢驗水準進行校正[1]。本研究提出的生存資料樣本量再估計方法是在盲態下只進行一次期中分析,該方法以極大似然估計為基礎,結合Todd[6]等人的研究方法,假定生存數據服從指數分布,兩組之間的風險比(HR)不變,在期中分析時對其他試驗參數進行重新估計,提出生存資料樣本量再估計的新方法。

原理與方法

對于服從指數分布的生存資料,期中分析在盲態下進行,其分組未知,相當于求解混合指數分布的參數問題。EM算法是Dempster[7]等人提出一種對于不完全數據進行參數估計的方法,以極大似然估計為基礎進行迭代,最終得到一個局部最優參數估計的迭代算法。朱利平[8]等人將EM算法應用到電子元件的壽命試驗中。研究表明,樣本的截尾比例越大,EM算法參數估計的效率越差,且估計結果往往依賴于初始值和收斂標準的設定[8-9]。在對生存資料進行期中分析時,往往產生大量的截尾數據。因此,本研究提出填補數據的方法,在期中分析時,用參數的初始估計值產生模擬數據,并用模擬數據中大于期中分析時間點的數據填補截尾數據,經過填補后的數據包含一定比例的真實參數信息和一定比例的填補參數信息。本研究結合EM算法的思想,但考慮到其參數估計結果依賴于初始值的設定,因此,提出在期中分析時以初始估計值為均值確定搜索范圍(hunting zone,HZ)對實際數據進行填補,然后,對填補后的數據進行極大似然估計。若初始估計值非常接近參數真實值,則此時參數的初始估計值與極大似然估計值也非常接近,則在搜索范圍內篩選出參數的極大似然估計值與初始估計值之差最小的參數極大似然估計值作為期中分析時的參數估計值。

極大似然法是混合分布參數估計的常用方法,期中分析時對服從混合指數分布的生存數據應用極大似然法,重新估計混合指數分布的參數進行樣本量再估計。期中分析時的數據是含有缺失數據的觀測數據,若用K表示分組(未知),X表示觀測數據,則包含分組的完整數據用Y表示:Y=(X,K)。λ是指數分布的參數,f(K|X,λ)是給定觀測數據X=Xi,λ=λi的條件下缺失數據K的條件密度。通過對觀測數據X的對數似然函數lnL(λ|X)求極大值得到λ的極大似然估計值。期中分析時的生存時間為觀測數據,用X表示,xi表示每個受試者的生存時間。受試者的分組情況用K表示,Ki=1表示受試者接受對照組治療,Ki=0表示受試者接受試驗組治療,因對照組與試驗組的樣本量比例為1:1,所以Ki服從二項分布,π(Ki=1)=π(Ki=0)=0.5,但是,期中分析是在盲態下進行,因此Ki是不能被觀測到的隨機變量。用λc表示對照組的風險率,λt表示試驗組的風險率,則對照組受試者的生存時間服從指數分布f1i,其概率密度函數為:

f1i=λce-λcxi(λc>0),

試驗組受試者的生存時間服從指數分布f2i,其概率密度為:

f2i=λte-λtxi(λt>0),

則xi服從混合指數分布fi,其概率密度函數為:

fi=0.5f1i+0.5f2i.

在混合指數分布fi中,記θ=(λc,λt),則xi和Ki的聯合分布為:

g(xi,Ki,θ)=(0.5f1i)Ki(0.5f2i)1-Ki,

因此,Ki在xi給定時的條件分布為:

xi給定時θ的似然函數為:

對數似然函數為:

給定參數初始估計值θ0(λc0,λt0)時,該對數似然函數的期望為:

通過極大化Q(θ,θ0)函數可得到參數θ(λc,λt)的極大似然估計值[8]。

在以事件為終點的隨機對照臨床試驗中,研究者根據以往數據估計受試者的生存時間,進而轉換為風險率λ進行樣本量計算。指數分布下風險率λ與中位生存時間M的轉換公式為:

本研究中定義隨訪時間為最后一個受試者入組開始至研究結束所經歷的時間。由于入組時間不同,因此最后一個入組受試者的隨訪時間決定了總研究時長,且所有受試者均隨訪至研究結束時間點,此時樣本量計算公式為(1)～(3)所示[10]:

(1)

(2)

(3)

公式(1)～(3)中,λc和λt分別是對照組和試驗組的風險率,T0是入組時間,T是總研究時長,k是試驗組和對照組的樣本含量比例,nc是對照組樣本量。

采用Monte Carlo模擬構建基于指數分布的生存資料模型,對模擬產生的生存資料進行log-rank檢驗,探討本研究方法對生存資料樣本量再估計的準確性和穩定性。本研究所有數據均由模擬產生,采用SAS 9.1.3統計分析軟件編寫程序,并進行數據分析。

1.參數設置

對于一個基于log-rank檢驗的雙臂臨床試驗,對照組與試驗組樣本比例為1:1,生存時間服從指數分布。試驗的真實參數如下,對照組的中位生存時間為M1=1年,對照組與試驗組的風險比HR為0.65和0.75,試驗入組時間為T0=0.5年,受試者在入組期內均勻入組,最短隨訪時間為t=2年,總研究時長為T=T0+t=2.5年,檢驗水準α為雙側0.05,檢驗效能1-β為0.8,則對照組的風險率λc為ln2/M1,試驗組的中位生存時間為M2=M1/HR,風險率λt為ln2/M2,試驗所需樣本量為N。

在進行期中分析的臨床試驗中,一般進行1～2次期中分析,本模擬研究考慮最大可能性,設定4個期中分析時間點(interim,I),并比較4個期中分析時間點的樣本量再估計結果,最終選定一個合理的期中分析時間點。4個期中分析時間點分別為I1-I4,I1是入組結束時(T0),I2是完成1/4隨訪時間時(T0+0.25t),I3是完成1/2隨訪時間時(T0+0.5t),I4是完成3/4隨訪時間時(T0+0.75t)。

2.模擬步驟

(1)根據預先設定的參數,運用SAS隨機函數ranexp(seed)先產生參數為M1、M2的指數分布數據x,受試者均勻入組,入組間隔為tn=T0/N0,入組順序即為隨機數的序號ob。令int=I-T0,令y=x-(N0-ob)×tn,將其中y>int者設定為y=int。

(2)重復步驟(1)產生參數為HM10、HM20的指數分布數據x,用其中y>int者填補參數為M1,M2數據中y=int的截尾數據,即為期中分析的數據。每個搜索范圍內有多個HM10,以HZ1(0.25～0.75)為例,搜索間隔為0.1年,則有0.25、0.35、0.45、0.55、0.65、0.75共6種填補參數,根據HR求得HM20,并以HM10、HM20為參數產生6組混合指數分布隨機數,用該6組隨機數中y>int的數據分別填補一次截尾數據,產生6組期中分析數據。

模擬結果

表1　不同HR時4個期中分析點不同HZ的樣本量再估計結果

表2　不同HR時4個期中分析點HZ11的樣本量再估計結果

圖1　不同HR時再估計樣本量與期中分析時間點、HZ之間的關系

圖2　不同HR時再估計樣本量與期中分析時間點、HZ(包含參數真實值)之間的關系

圖3　各期中分析時間點上HZ11-HZ5再估計樣本量的α(a)和power(b)

模擬示例

討　　論

藥物臨床試驗的投入巨大,而樣本量是否充足直接關系到試驗的成敗。若在試驗設計階段對試驗參數的估計存在較大的不確定性,建議采用適應性設計,預設期中分析點對樣本量進行再估計,根據估計結果考慮是否增加樣本量。一般期中分析時只考慮是否增加樣本量以維持把握度,而不考慮降低樣本量[1]。Todd[6]等人提出的盲態下樣本量再估計方法,需要設定多個期中分析時間點進行多次期中分析,具體實施有一定的難度。本研究提出的盲態下生存資料樣本量再估計方法,只需要進行一次期中分析。本研究假設對照組與試驗組的樣本量之比為1:1,風險比為HR,α為0.05,power為0.8,并假設5個不同的M10,并以M10為均值暫定搜索范圍為M10±0.5M10,比較4個期中分析點上不同搜索范圍下樣本量再估計結果。模擬研究表明,在期中分析點I2～I4,當搜索范圍包含參數真實值時,經過期中分析調整后的樣本量接近真實值,且樣本量估計結果變異較小。在I1時間點時經填補的期中分析數據包含真實參數信息較少,搜索結果與填補參數相關性較大,因此樣本量再估計結果呈現出與搜索范圍值大小相關的現象,但在I2～ I4時間點時未出現該現象。

綜上所述,并考慮到實際應用中,若樣本量再估計結果大于初始樣本量,需要繼續追加樣本量,并完成相同的隨訪時間t,此時期中分析時間點越晚,則研究時間就越長,研究成本也隨之增加。因此,相似背景下建議在I2(完成1/4隨訪時間)時進行期中分析,根據參數估計結果重新確定搜索范圍,以獲得更加準確的樣本量。若是其他差異較大的背景時,1/4隨訪時間作為期中分析時間點劃分的參考依據,可進行適當的調整。若樣本量再估計結果比試驗設計之初的樣本量少,則維持原定樣本量不變;反之,則根據樣本量再估計結果考慮增加樣本量以維持把握度。本研究的不足之處是在HR恒定不變的假設下進行模擬研究,后續研究將會考慮HR改變時的情況。