基于右刪失數據下加速失效跡回歸模型的估計

中圖分類號：O212 文獻標志碼：A 文章編號：1671-5489（2025）04-1059-09

Estimation of Accelerated Failure Trace Regression Model Based on Right Censored Data

FAN Yifan，XU Ping，XIAO Nannan，WANG Chunjie （School of Mathematics and Statistics，Changchun University of Technology， Changchun 130ol2，China）

Abstract： Aiming at the challenges of high-dimensional medical image data in survival analysis，we proposed an accelerated failure trace regression model. The regression parameters were estimated by using Kaplan-Meier weighting and the Peaceman-Rachford algorithm. The numerical simulation results show that the estimation performance of the accelerated failure trace regression model is better than that of the traditional Lasso regression model. We apply the model to Alzheimer's disease image data to further verify its effectiveness and practical value.

Keywords：accelerated failure trace regression model； high-dimensional right censored data； Kaplan-Meier weighting; Peaceman-Rachford spliting algorithm； Alzheimer's disease image data

隨著信息采集技術的進步和大數據時代的到來，醫療圖像數據在生存分析中應用日益廣泛.然而，與常規數據結構不同，醫療圖像數據具有高維特性，并反映了像素點之間的復雜關系，為數據分析提出了新挑戰．近年來，針對高維數據背景下的生存分析模型研究備受關注，其中較常見的是加速失效時間（accelerated failure time，AFT）模型．Huang 等[1]提出了一種穩健的加權最小絕對偏差法估計高維AFT模型；Liu等2提出了一種快速預測高維變量AFT模型的變量選擇和收縮估計的方法，并使用乳腺癌研究的基因表達條形碼數據進行模擬和驗證；Schmid等[3]在高維非參AFT模型中，利用 Hilbert 范數懲罰進行估計，提出了一種新的增強算法；Khan等[4在具有高維預測變量的 AFT模型中提出了變量選擇方法．盡管這些研究為高維數據的生存分析提供了有效的方法，但傳統線性模型很難有效處理圖像數據中的像素關系，可能導致偏差，因此，有必要對基于高維圖像數據的生存分析模型開展進一步研究.

圖像數據是很多模型研究的核心對象，其中跡回歸（trace regresson）模型因其能有效處理高維數據而備受關注，該模型主要用于建立圖像數據預測變量與響應變量之間的關系．Koltchinskii等[5]研究了跡回歸模型，并將其應用到矩陣補全中；Negahban 等[6]對高維度有噪聲的低秩矩陣進行參數估計，提出了一種基于矩陣核范數或跡范數正則化的M-估計量，為跡回歸模型的參數估計奠定了基礎．在跡回歸模型研究中，許多研究者考慮在高維數據背景下研究稀疏約束或低秩約束的矩陣回歸問題，例如：Zhou等[7提出了低秩矩陣回歸模型，建立了一種高效且可擴展的估計算法；Elsener等[8]考慮到系數矩陣的低秩性，提出了帶有核范數懲罰的正則化Huber 矩陣回歸方法；此外，為研究COVID-19 疾病產生的二維圖像，Zhang 等[9提出了一種協變量為高維的情況下新的潛在矩陣因子回歸模型.

在完整數據條件下，跡回歸模型已得到廣泛研究，但在刪失數據條件下的跡回歸模型研究相對較少．因此，本文基于高維右刪失數據，考慮帶有圖像矩陣結構變量的加速失效時間模型的建模問題，不僅擴展了跡回歸模型在刪失數據條件下的應用范圍，還為高維數據分析提供了一種新方法.

1 數據與模型

假設試驗中有 n 個個體，對第 i 個個體，記感興趣的事件失效時間為 T_i ，刪失時間為 C_i ，存在可觀測的協變量 X_i ，當數據為右刪失時，得到的觀測數據結構為

O={O_i=（Y_i^*，X_i，δ_i），i=1，2，…，n}，

其中： Y_i^*=min{T_i，C_i} ; δ_i=I（T_i?C_i） ， i=1，2，…，n ， I（?） 為示性變量.

失效時間 T_i 滿足如下加速失效跡回歸（accelerated failure time trace regression）模型：

log（T_i）=?X_i，Θ^*?+ε_i，

其中 為 d₁×d₂ 維協變量矩陣， 是系數矩陣， 表示 d₁×d₂ 的實矩陣空間， ε_i∈R 是誤差項．對任意兩個矩陣 ，定義內積 ?A，B?：=tr（A^TB） ，其中 tr（?） 表示跡.

本文對系數矩陣 θ^* 施加秩約束，矩陣的秩被其奇異值的 l_q 范數控制，所以 θ^* 應滿足：

其中： d₁∧d₂=min{d₁，d₂} ; σ_i（θ^?） 是 θ^* 的第 j 個最大奇異值； ρ 是一個非負常數， ρ 隨維數和樣本量增長而變化．注意 0?q?1 ，當 q=0 時，約束式（2）是一個精確秩約束．對 B_q（θ^?） 的限制保證了奇異值衰減足夠快，且比精確的低階假設更普遍[10].

2參數估計

定義失效時間 T 的分布函數為 F ，其Kaplan-Meier 估計量[1]為 ： 可寫成如下形式： ，其中 ω_i 是 Kaplan-Meier 估計量中的跳躍值，可表示為

ω_i 也稱為Kaplan-Meier權重， Y_?（1）?…?Y_?（n） 是 Y 的次序統計量， X_（1）?…?X_（n） 是 X 的次序統計量，（20 δ_（1）?…?δ_（n） 是示性變量.

根據文獻[12]的思想，本文提出加權最小二乘估計方法，通過最小化以下目標函數，估計回歸系數 θ^* ：

通過 ω_i 權重均值分別中心化 X_（i） 和 Y_?（i） ，定義

令 ， ．所以，加權最小二乘目標函數（5）可改寫為

對式（6）求期望，有

Π_ω（i）-?Θ，X_ω（i）?）²=EY_ω（i）²-2?Θ，EY_ω（i）X_ω（i）?+vec（Θ）^?E（vec（X_ω（i））vec（X_ω（i））^?）vec（X_ω（i）），

由于 Y_ω（i）² 與 Θ 無關，所以可忽略 EY_ω（i）² ，并用樣本協方差代替 Σ_Yω（i）X_ω（i） 和 Σ_Xω（i）X_ω（i） ，得到以下公式：

其中 和 分別是 EY_ω（i）X_ω（i） 和 E（vec（X_ω（i））vec（X_ω（i））^T） 的估計量.

為估計 θ^* ，將跡回歸模型與正則化相結合，通過核范數懲罰加強稀疏和低秩結構，得到以下廣義l₂ 損失的M-估計量方法：

其中 是一凸集， λ_N 表示懲罰系數，是 Θ 的核范數．類似 l₁ 范數正則化產生稀疏估計，對核范數正則化使得解具有稀疏奇異值．對于響應變量是單變量的加速失效跡回歸模型， 和 采用如下形式：

3 PRSM算法

利用 Peaceman-Rachford 分裂算法（PRSM）[13]在線性約束下最小化兩個凸函數的和：

min{f₁（θ_x）+f₂（θ_y）}，

將 C₂=-I ， c=0 ， ， 代入式（11），得到加速失效跡回歸模型的PRSM算法．vec （X）= ， X_j 是 X 的第 j 列．對于向量 ，用 mat（x） 表示由 x 構造的 d₁×d₂ 矩陣，其中 （x_{（j-1）d1+1}，…，x_jd1）^T 是 mat（x） 的第 j 列， 表示 d₁d₂ 維的實向量空間.

利用Peaceman-Rachford分裂算法，通過算法中定義的交叉驗證方法，選取最優懲罰系數，從而得到參數 θ^* 的估計值．設 X 是一個 n×d₁d₂ 矩陣，其行為相互獨立的向量 {vec（X_i）}_i=1ⁿ ， Y 是 n 維響應向量，加速失效跡回歸模型的PRSM算法如下：

其中 是Lagrange乘子， β 是懲罰參數， α 是松弛因子，并根據Eckstein等[14]和He等[15]的研究結果，取 α∈{0.1，0.2，…，0.9，1}，β∈{0.2，0.3，…，0.9，1，2，4，8，16}， 再通過BIC準則對其進行選取，當右刪失數據刪失率為 20% 時， α=0，1 ， β=0.2 ；當右刪失數據刪失率為 40% 時， α=0，3 ，

（2 β=0.2 ： S （z）是矩陣 的奇異值軟閾值函數，設 ， Z=UAV^T= Udiag（λ₁，…，λ_r）V^T 為其奇異值分解， S（z）=vec（Udiag（（λ₁）₊，（λ₂）₊，…，（λ_r）₊）V^T），（λ_r）₊ （x）₊= max{x，0} .當 |θ_x-θ_y| 小于某個預定的閾值時，返回的 mat（θ_s） 作為 θ^* 的最終估計量.

4數值模擬

下面進行數值模擬以驗證本文方法的有效性．考慮維度 d₁=d₂=d=20，40，60 ，回歸系數矩陣的Frobenius 范數 ，回歸系數矩陣的秩 rank（Θ^*）=5 ，構造 ，其中 ν_i 是100個獨立同分布且中心化的高斯隨機向量的樣本協方差矩陣的第 i 個特征向量，這些向量的協方差為 I_d ．對于式（1），協變量矩陣 X_i 服從標準多元正態分布，刪失時間 C_i 服從均勻分布 U（0，u） ，其中常數 u 控制刪失率為 20% 和 40% ，失效時間 T_i 由式（1）生成，參考文獻[10]，對隨機誤差 ε_i 設置如下3種分布類型

1）隨機誤差服從對數正態分布： s_i=（Z-EZ）/50 ，其中 ， σ²=6.25 ：

2）隨機誤差服從Cauchy分布： ε_i=min{∣Z∣，10³}/10， Z 服從Cauchy分布；

3）隨機誤差服從高斯分布： ，其中 σ²=0.25

模擬計算分別產生 300×k 個樣本 （k=1，2，…，7） ，為檢驗參數 θ^* 的估計效果，參考文獻[10]，選取" 統計量進行評估， 統計量的值越小，估計效果越好．此外，還引入Lasso 回歸，與本文的加速失效跡回歸進行比較．表1～表3分別列出了不同模型、不同維度和不同刪失率下回歸系數 θ^* 的 統計誤差.

由表 1～ 表3可見：加速失效跡回歸模型對參數估計的效果較好；在同一維度下，隨著樣本的增大，3種誤差分布下的 統計誤差逐漸減小；在不同刪失率下，3種誤差分布下的 F統計誤差結果相似，說明模型不依賴于刪失率的變化；而在不同誤差分布下， F統計誤差結果有差異，當誤差分布為高斯分布和對數正態分布時，模型的估計效果較好.

在不同刪失率及不同誤差分布下 ：與對數樣本量 的統計誤差分別如圖 1～ 圖6所示，其中實線表示加速失效跡回歸模型下的估計效果，虛線表示Lasso回歸下的估計效果．由圖 1～

圖6可見：在不同的刪失率下，3種誤差分布的 F統計誤差隨樣本量的增大而減小；不同維度下，3種誤差分布的 ：統計誤差隨樣本增大均出現下降趨勢；加速失效跡回歸模型下的估計效果較好.

5 實例分析

阿爾茲海默病神經影像研究（the Alzheimer's disease neuroimaging initiative，ADNI）是一個多中心、縱向隨訪的研究項目，通過深人分析臨床數據、影像學資料、遺傳信息以及生化標志物等多維度信息[16-17]，可實現阿爾茲海默病（AD）的早期檢測與持續追蹤．在 ADNI中，受試者會按一定時間間隔接受定期的隨訪觀察和健康狀態評估．研究人員會基于他們的認知能力，將受試者分為3個不同群體：認知功能正常（CN）、輕度認知障礙（MCI）和阿爾茲海默病（AD）．MCI受試者中的一部分個體可能會逐漸發展到AD階段，而另一部分則可能保持穩定或恢復到CN狀態，其中需特別關注的一個變量是受試者確切患上AD的時間．將AD轉化的時間視為失效時間，而最后一次監測時間視為刪失時間，該數據集可視為右刪失數據.

本文將加速失效跡回歸模型應用于阿爾茲海默病神經影像數據中進行實證分析．采用的實例數據集調查時間為2014年5月至2022年5月，包含了120名受試者及其對應的磁共振圖像（magneticresonanceimaging，MRI，該數據集中有59名受試者為右刪失，占比約 49.2% ．磁共振圖像采集自通用電氣醫療系統（GEMedicalSystems） 1.5T 的 MRI掃描儀的標準T1加權像，圖像數據參數包括：TR=8.9ms ： TI=1000.0ms ： FA=8. 0 ；三維采集矩陣為 256×256×166 ；體素大小為 （0.9×0.9×

1.2）mm³ ．MRI圖像數據遵循Rician分布[18]．然而，在MRI進行預處理過程中，降采樣過程會改變誤差特性，使誤差分布可能更趨近于高斯分布[19-20]．MRI通過標準步驟進行預處理，涵蓋了前聯合與后聯合（AC-PC）線的校正、顱骨及小腦組織的去除以及圖像的空間歸一化與配準等關鍵環節.

對MRI圖像進行預處理，流程概述如下：

1）將原始的DICOM格式數據轉換成NIFTI格式，以便于后續處理；2）對每個受試者的T1加權結構圖像進行手動AC-PC校正[21]；3）進行頭骨剝離和小腦切除操作，并進一步分割出白質、灰質和腦脊液區域[22]；4）將分割出的灰質圖像與Montreal神經學研究所（Montreal neurological institute，MNI）的標準

空間進行配準[23]；5）應用放射變換和非線性扭曲調制，保持空間配準前的組織體積[24」；6）使用 8mm 半高全寬（full width at half maxima，FWHM）的高斯核對圖像進行空間平滑處理，

以減少誤差；7）對圖像降采樣處理，降至 100×100×100

預處理步驟如圖7所示．經過圖像預處理，選擇每個受試者切片中像素點之和最大的腦結構冠狀面 100×100 大小的MRI圖像，如圖8所示．基于加速失效跡回歸模型，對MRI圖像數據進行系數估計，結合MRI圖像數據和估計的回歸系數，生成熱力圖，結果如圖9和圖10所示.

圖9和圖10分別為高維圖像加速失效跡回歸模型和高維圖像Lasso回歸模型的實驗結果，其中紅色深淺表示不同像素對AD發病時間的影響程度，顏色越深表示影響越大．腦島在阿爾茲海默病發展過程中具有重要作用，這主要歸因于其與感知能力、自我意識以及認知功能的緊密聯系，與Hu等[25]的研究結果，腦島是阿爾茲海默病生物標志物的核心價值，結論相似．此外，阿爾茲海默病同樣受顳葉區域的顯著影響，該區域涵蓋了杏仁核、海馬體等重要結構．另一個重要影響位置是海馬體，它負責將短期記憶轉化為皮層長期記憶進行存儲，該區域在學習與鞏固過程中起重要作用，特別易受阿爾茲海默病病理的影響，且在臨床出現阿爾茲海默病癥狀時，該區域會受相當大的損傷[26]．本文提出的加速失效跡回歸模型，找到的海馬體和腦島等區域對阿爾茲海默病有顯著影響；而在Lasso 回歸模型中，只找到了腦島區域，估計效果較差.同時兩種模型也找到了一些未驗證的腦部區域作為影響因素.

綜上所述，跡回歸模型能有效處理高維圖像數據，并在完整數據條件下已得到廣泛研究，但在刪失數據條件下的跡回歸模型研究相對有限．本文考慮帶有圖像矩陣結構的高維右刪失數據下加速失效跡回歸模型的估計建模問題，采用Kaplan-Meier 加權的方法對高維右刪失數據進行處理，使用Peaceman-Rachford分裂算法對回歸參數進行估計，并與Lasso 回歸進行了對比研究．模擬研究結果表明：在同一維度條件下，隨樣本數量的增加，3種誤差分布的 統計誤差呈現出逐漸降低的趨勢；在不同刪失率下，3種誤差分布的 統計誤差結果具有相似性，模型不依賴于刪失率的變化；即使誤差分布各異， F統計誤差的結果也相近，當誤差分布為高斯分布和對數正態分布時，模型的估計效果較好；對比不同模型的處理結果可見，加速失效跡回歸模型的估計效果比Lasso 回歸估計效果更好．將該模型應用于阿爾茲海默病圖像數據，進一步驗證了其有效性和實用價值.

參考文獻

[1]HUANG J，MA SG，XIE HL. Least Absolute Deviations Estimation for the Acelerated Failure Time Model [J].Statistica Sinica，2007，17（4）：1533-1548.

[2]LIU F，DUNSON D， ZOU F. High-Dimensional Variable Selection in Meta-Analysis for Censored Data[J]. Biometrics，2011，67（2）：504-512.

[3]SCHMID M，HOTHORN T. Flexible Boosting of Accelerated Failure Time Models [J]. BMC Bioinformatics， 2008，9： 269-1-269-13.

[4]KHAN M H R，SHAWJE H. Variable Selection for Accelerated Lifetime Models with Synthesized Estimation Techniques [J]. Statistical Methods in Medical Research，2019，28（3）： 937-952.

[5]KOLTCHINSKII V，LOUNICI K，TSYBAKOV A B. Nuclear-Norm Penalization and Optimal Rates for Noisy Low-Rank Matrix Completion [J]. The Annals of Statistics， 2011，39（5）： 2302-2329.

[6]NEGAHBAN S，WAINWRIGHT M J. Estimation of （Near） Low-Rank Matrices with Noise and HighDimensional Scaling [J]. The Annals of Statistics，2011，39（2）：1069-1097.

[7]ZHOU H，LI L X. Regularized Matrix Regresson[J]. Journal of the Royal Statistical Society，Series B： Statistical Methodology，2014，76（2）： 463-483.

[8] ELSENER A，VAN DE GEER S. Robust Low-Rank Matrix Estimation [J]. The Annals of Statistics，2018， 46（6B）：3481-3509.

[9]ZHANG Y Z， ZHANG X， ZHANG H，et al. Low-Rank Latent Matrix-Factor Prediction Modeling for Generalized High-Dimensional Matrix-Variate Regression[J]. Statistics in Medicine，2023，42（2O）：3616-3635.

[10]FAN JQ，WANG W， ZHU Z W. A Shrinkage Principle for Heavy-Tailed Data：High-Dimensional Robust LowRank Matrix Recovery [J]. Annals of Statistics，2021，49（3）：1239-1266.

[11]HU J W，CHAI H. Adjusted Regularized Estimation in the Accelerated Failure Time Model with High Dimensional Covariates [J]. Journal of Multivariate Analysis，2013，122： 96-114.

[12] STUTE W. Distributional Convergence under Random Censorship When Covariables Are Present [J]. Scandinavian Journal of Statistics，1996，23（4）：461-471.

[13] PEACEMAN D W，RACHFORD H H，Jr. The Numerical Solution of Parabolic and Elliptic Differential Equations [J]. Journal of the Society for Industrial and Applied Mathematics，1955，3（1）： 28-41.

[14]ECKSTEINJ，BERTSEKAS D P. On the Douglas-Rachford Spliting Method and the Proximal Point Algorithm for Maximal Monotone Operators [J]. Mathematical Programming，1992，55（3）： 293-318.

[15] HE B S，LIU H，WANG Z R，et al. A Strictly Contractive Peaceman-Rachford Splitting Method for Convex Programming [J]. SIAM Journal on Optimization，2014，24（3）：1011-1040.

[16]LI K，CHAN W，DOODY R S，et al. Prediction of Conversion to Alzheimer’s Disease with Longitudinal Measuresand Time-to-Event Data [J]. Journal of Alzheimer's Disease， 2017，58（2）： 361-371.

[17] WU Q W， ZHAO H， ZHU L，et al. Variable Selection for High-Dimensional Partly Linear Aditive Cox Model with Application to Alzheimer's Disease [J]. Statistics in Medicine，2020，39（23）： 3120-3134.

[18] GUDBJARTSSON H，PATZ S. The Rician Distribution of Noisy MRI Data [J]. Magnetic Resonance in Medicine，1995，34（6）：910-914.

[19]MANJON JV， COUPE P，MARTi-BONMATi L，et al. Adaptive Non-local Means Denoising of MR Images with Spatially Varying Noise Levels [J]. Journal of Magnetic Resonance Imaging，20lo，31（1）：192-203.

[20]KOAY C G，BASSER PJ. Analytically Exact Correction Scheme for Signal Extraction from Noisy Magnitude MR Signals [J]. Journal of Magnetic Resonance，2006，179（2）：317-322.

[21]BOUTTELGIER R，VANDAMME S，VERVERKEN F，et al. Deep Brain Stimulation for Essential Tremor in Patients with Ventriculomegaly [J/OL]. Surgical Neurology International，（2024-07-19）[2024-08-2O]. doi： 10.25259/SNI_979_2023.

[22] SMITH S M. Fast Robust Automated Brain Extraction [J]. Human Brain Mapping，2002，17（3）： 143-155.

[23] ZENNADI M M， PTITO M， REDOUTE J，et al. MRI Atlas of the Pituitary Gland in Young Female Adults [J]. Brain Structure and Function，2024，229（4）：1001-1010.

[24]STROTHER S，LA CONTE S， HANSEN L K，et al. Optimizing the fMRI Data-Processing Pipeline Using Prediction and Reproducibility Performance Metrics：I.A Preliminary Group Analysis [J].Neuroimage，2004， 23（Suppl 1） ： S196-S207.

[25] HU X，MEIBERTH D，NEWPORT B，et al. Anatomical Correlates of the Neuropsychiatric Symptoms in Alzheimer's Disease [J].Current Alzheimer Research，2015，12（3）：266-277.

[26]SCHUFF N，WOERNER N，BORETA L，et al. MRI of Hippocampal Volume Loss in Early Alzheimer’s Disease in Relation to ApoE Genotype and Biomarkers [J]. Brain，2009，132（Pt4）： 1067-1077.

（責任編輯：李琦）