1. 生存分析在轉錄組流程中的定位與核心價值1.1 為什么表達矩陣跑完差異分析之后要接生存分析很多做轉錄組的朋友把差異表達基因DEG列表跑出來之后第一反應是去做GO/KEGG富集畫完氣泡圖就覺得這個項目差不多收尾了。但如果你的樣本來自臨床隊列、隨訪記錄齊全那DEG列表其實只是半成品——真正能決定這批數(shù)據(jù)能不能發(fā)好刊、能不能落到臨床轉化上的是生存分析這一步。生存分析要回答的問題非常直接某個基因的表達高低和樣本的生存時間總生存OS、無進展生存PFS、無復發(fā)生存RFS等之間到底有沒有統(tǒng)計學上的關聯(lián)。如果高表達組的中位生存期比低表達組短了一半且log-rank檢驗p值小于0.05那這個基因就具備成為預后標志物prognostic biomarker的潛力。這跟單純看“腫瘤vs正常”的差異表達完全是兩個維度的事情——差異表達告訴你“這個基因在疾病里變了”生存分析告訴你“這個基因的變化對結局有沒有影響”。我在實際項目里見過太多這樣的情況一個基因在腫瘤組織里上調了5倍差異極顯著但一做生存分析p值0.6毫無預后價值反過來有些基因差異倍數(shù)只有1.3倍看著不起眼但生存曲線分得干干凈凈p值0.001這種基因反而更值得深挖。所以生存分析本質上是給DEG列表做一次“臨床價值篩選”把那些真正和結局掛鉤的基因挑出來。1.2 生存分析適合什么樣的數(shù)據(jù)場景不是所有轉錄組數(shù)據(jù)都能做生存分析這一點必須先說清楚否則后面全是白費功夫。能跑生存分析的數(shù)據(jù)至少要滿足三個條件第一樣本必須有隨訪信息。每個樣本要對應一個生存時間time和一個生存狀態(tài)status通常1代表死亡/事件發(fā)生0代表刪失/censored。這兩個字段缺一不可而且時間單位要統(tǒng)一天、月、年都行但全表必須一致。第二樣本量要夠。生存分析對樣本量的要求比差異分析更苛刻。差異分析30vs30就能跑但生存分析如果總樣本低于50分組之后每組可能就十幾例log-rank檢驗的效能會非常低很容易出現(xiàn)假陰性。我個人的經(jīng)驗是總樣本至少60例以上每組不少于20例結果才相對可靠。如果樣本量在40-60之間可以做但結論要謹慎最好用連續(xù)變量而非二分高低組的Cox回歸來補充驗證。第三隨訪時間要足夠長。如果中位隨訪時間只有3個月而事件發(fā)生率又很低那大部分樣本都是刪失的生存曲線會非常平幾乎分不開。理想情況下事件數(shù)死亡或復發(fā)例數(shù)至少要達到總樣本的30%以上。提示如果你的數(shù)據(jù)來自公共數(shù)據(jù)庫比如某類腫瘤的公開隊列下載時一定要看清楚隨訪字段的定義有些數(shù)據(jù)集里“status1”代表存活而不是死亡搞反了整篇文章的結論就全錯了。1.3 生存分析在整個轉錄組流程中的位置把生存分析放進完整的轉錄組分析流程里看它的位置大概是這樣原始數(shù)據(jù)質控與比對表達定量count矩陣或TPM/FPKM矩陣差異表達分析得到DEG列表功能富集分析GO/KEGG/GSEA生存分析本篇重點預后模型構建多基因Cox模型、風險評分獨立預后因素驗證單因素/多因素Cox可以看到生存分析是連接“分子層面發(fā)現(xiàn)”和“臨床預后價值”的橋梁。它往下可以延伸到多基因預后模型的構建往上則依賴前面幾步產出的干凈表達矩陣。所以這一步做得好不好直接決定了后續(xù)模型構建的天花板。2. 數(shù)據(jù)準備與預處理的關鍵細節(jié)2.1 表達矩陣的格式要求與常見坑生存分析的輸入數(shù)據(jù)通常是一張表達矩陣行是基因列是樣本外加一張臨床信息表含生存時間和狀態(tài)。表達矩陣的格式有幾種常見類型處理方式不太一樣Raw count矩陣直接來自featureCounts、HTSeq等定量工具整數(shù)。做生存分析前建議做log2(count1)轉換或voom/TPM標準化否則極端高表達基因會主導分組。TPM/FPKM矩陣已經(jīng)是標準化后的連續(xù)值通常再做一次log2(x1)轉換即可。微陣列數(shù)據(jù)可能是已經(jīng)log2轉換過的信號值需要確認是否還需要進一步標準化。我踩過最多的坑是樣本ID不匹配。表達矩陣的列名可能是“TCGA-XX-XXXX-01A”而臨床表的樣本ID可能是“TCGA-XX-XXXX”少了后面的“-01A”后綴直接merge就會得到空表。解決辦法是用字符串截取統(tǒng)一ID長度比如取前12位或前15位。另一個坑是重復樣本。同一個患者可能有多個樣本比如癌和癌旁或者不同時間點取樣做生存分析時一個患者只能保留一個樣本否則會違反生存分析的獨立性假設。通常保留腫瘤樣本去掉癌旁如果有多個腫瘤樣本取表達值的中位數(shù)或第一個時間點的樣本。2.2 生存時間與狀態(tài)字段的清洗生存時間字段最常見的三個問題問題一時間為0或負數(shù)。有些數(shù)據(jù)集里生存時間記錄為0這在生存分析里是沒有意義的log-rank檢驗和Cox模型都無法處理。處理方式是直接剔除這些樣本或者如果數(shù)量很少3例可以賦值為0.5最小時間單位的一半。問題二狀態(tài)字段編碼不統(tǒng)一。有的數(shù)據(jù)集用“Dead/Alive”有的用“1/0”有的用“DECEASED/LIVING”。統(tǒng)一轉成數(shù)值型事件發(fā)生死亡/復發(fā)1刪失存活/失訪0。問題三時間單位混亂。有的數(shù)據(jù)集時間是“天”有的是“月”有的是“年”。建議全部統(tǒng)一成“月”因為臨床文獻里中位生存期通常以月為單位報告方便對比。下面是一段我常用的R代碼用來清洗和合并表達矩陣與臨床數(shù)據(jù)# 假設 expr 是表達矩陣clin 是臨床信息表 # 統(tǒng)一樣本ID colnames(expr) - substr(colnames(expr), 1, 12) clin$sample_id - substr(clin$sample_id, 1, 12) # 去重同一患者保留一個樣本 clin - clin[!duplicated(clin$sample_id), ] expr - expr[, colnames(expr) %in% clin$sample_id] # 對齊順序 common_samples - intersect(colnames(expr), clin$sample_id) expr - expr[, common_samples] clin - clin[match(common_samples, clin$sample_id), ] # 清洗生存時間和狀態(tài) clin - clin[clin$OS_time 0, ] clin$OS_status - ifelse(clin$OS_status %in% c(Dead, DECEASED, 1), 1, 0) # log2轉換 expr - log2(expr 1)這段代碼看起來簡單但每一步都對應一個實際踩過的坑。尤其是match()那一步如果不做后面分組的時候樣本順序錯亂結果會完全不可復現(xiàn)。2.3 基因篩選策略用全部基因還是只用DEG這是一個經(jīng)常被問到的問題做生存分析時是把所有基因都跑一遍還是只跑差異表達基因兩種策略各有道理全基因掃描適合探索性分析可以發(fā)現(xiàn)一些差異不顯著但有預后價值的基因。缺點是多重檢驗校正后可能沒有基因顯著尤其是樣本量小的時候。僅用DEG適合驗證性分析把范圍縮小到已經(jīng)證明有差異的基因上多重檢驗負擔小更容易得到顯著結果。缺點是可能漏掉那些“差異不大但預后很強”的基因。我通常的做法是兩步走先對全基因做單因素Cox回歸篩選不校正或只做寬松校正把p0.05的基因挑出來然后在這些基因里再看哪些是DEG。這樣既不會漏又能和差異分析結果互相印證。注意如果基因數(shù)量超過2萬做單因素Cox回歸時一定要用并行計算否則R會跑到你懷疑人生。survival包配合future.apply可以輕松并行。3. 生存分析的核心方法拆解3.1 Kaplan-Meier曲線與log-rank檢驗最直觀的分組比較KM曲線是生存分析里最直觀的工具沒有之一。它的核心思想是把樣本按某個基因的表達值分成高/低兩組然后畫出兩組的生存概率隨時間變化的階梯曲線最后用log-rank檢驗判斷兩條曲線是否有顯著差異。分組方式是這一步的關鍵。常見的有三種中位數(shù)分組表達值大于中位數(shù)的為高表達組小于等于的為低表達組。這是最常用的方法簡單、穩(wěn)健每組樣本量相等。最佳截斷值optimal cutpoint用survminer包的surv_cutpoint()函數(shù)自動尋找使log-rank統(tǒng)計量最大的截斷值。優(yōu)點是能最大化組間差異缺點是容易過擬合尤其是在小樣本里。四分位數(shù)分組取上四分位和下四分位中間部分丟棄。適合樣本量很大的情況但會損失信息。我個人的建議是主分析用中位數(shù)分組敏感性分析用最佳截斷值。如果兩種分組方式結論一致那結果就很穩(wěn)如果不一致就要謹慎解讀說明這個基因的預后價值可能不夠穩(wěn)健。下面是一段完整的KM分析代碼library(survival) library(survminer) # 按中位數(shù)分組 expr_gene - expr[GENE_OF_INTEREST, ] group - ifelse(expr_gene median(expr_gene), High, Low) group - factor(group, levels c(Low, High)) # 構建生存對象 surv_obj - Surv(time clin$OS_time, event clin$OS_status) # KM擬合 fit - survfit(surv_obj ~ group) # 畫圖 ggsurvplot( fit, data clin, pval TRUE, risk.table TRUE, palette c(#2E9FDF, #E7B800), xlab Time (months), ylab Overall Survival Probability, title Kaplan-Meier Curve for GENE_OF_INTEREST, legend.title Expression, legend.labs c(Low, High) )這段代碼跑出來的圖可以直接放進論文。risk.table TRUE會在圖下方顯示每個時間點的風險人數(shù)這是很多期刊的硬性要求。3.2 單因素與多因素Cox回歸從分組到連續(xù)變量KM曲線只能告訴你“高表達組和低表達組有沒有差異”但沒法告訴你“表達值每升高一個單位風險增加多少”。Cox比例風險模型就是干這個的。單因素Cox回歸每次只放一個變量基因表達值得到HR風險比、95%置信區(qū)間和p值。HR1表示高表達是風險因素表達越高死亡風險越大HR1表示保護因素。多因素Cox回歸同時放入多個變量比如基因表達年齡性別分期校正混雜因素后看基因是否還是獨立預后因素。這一步非常重要因為很多基因的單因素顯著性其實是被臨床分期帶出來的——晚期患者表達高、生存差基因本身可能沒有獨立價值。# 單因素Cox cox_univ - coxph(Surv(OS_time, OS_status) ~ expr_gene, data clin) summary(cox_univ) # 多因素Cox校正年齡、性別、分期 cox_multi - coxph( Surv(OS_time, OS_status) ~ expr_gene age gender stage, data clin ) summary(cox_multi)解讀Cox結果時要注意幾點HR的置信區(qū)間不能跨1否則即使p0.05也不可靠。比例風險假設PH assumption必須檢驗。用cox.zph()函數(shù)如果p0.05說明違反了PH假設需要考慮時依系數(shù)或分層Cox。連續(xù)變量vs二分變量連續(xù)變量的HR解釋為“每單位增加的風險變化”二分變量的HR解釋為“高組相對低組的風險變化”。兩者可以都做互相驗證。3.3 森林圖與多基因預后模型從單基因到組合當你篩出多個有預后價值的基因后下一步自然是把它們組合起來構建預后模型。最常見的是LASSO-Cox或多因素Cox逐步回歸得到一個風險評分公式riskScore β1×expr1 β2×expr2 ... βn×exprn然后按風險評分的中位數(shù)把樣本分成高風險組和低風險組再做KM曲線和ROC曲線驗證模型效能。森林圖是展示多基因Cox結果的標配用forestplot包或survminer的ggforest()都能畫。我一般會把單因素和多因素的結果放在一張圖里對比這樣審稿人一眼就能看出哪些基因在校正后仍然顯著。library(forestplot) # 假設已經(jīng)有一個包含HR、lower、upper、pvalue的數(shù)據(jù)框 forestplot( labeltext c(Gene, rownames(result_table)), mean c(NA, result_table$HR), lower c(NA, result_table$lower), upper c(NA, result_table$upper), zero 1, xlog TRUE )4. 實操全流程與關鍵參數(shù)記錄4.1 從表達矩陣到生存分析結果的完整步驟把前面幾節(jié)的內容串起來一個完整的生存分析實操流程大概是這樣第一步數(shù)據(jù)加載與清洗。讀入表達矩陣和臨床表統(tǒng)一樣本ID去重剔除時間為0的樣本log2轉換。第二步基因篩選。如果基因太多先用方差或中位絕對偏差MAD過濾掉低變異基因保留前5000-10000個高變異基因做后續(xù)分析。這一步能大幅減少計算量而且不會損失關鍵信息。第三步批量單因素Cox。對每個基因跑單因素Cox記錄HR、p值、置信區(qū)間。用BH方法做多重檢驗校正得到FDR。第四步篩選顯著基因。通常取FDR0.05或p0.01作為閾值。如果顯著基因太多100個可以進一步收緊到FDR0.01。第五步KM曲線可視化。對top顯著基因畫KM曲線確認分組差異方向與Cox結果一致。第六步多因素Cox驗證。把top基因和臨床變量一起放入多因素模型篩選獨立預后因素。第七步模型構建與驗證。用LASSO或逐步回歸構建風險評分畫ROC曲線、校準曲線做內部驗證如Bootstrap或外部驗證如果有獨立隊列。4.2 參數(shù)選擇與計算過程實錄這里記錄幾個關鍵參數(shù)的選擇過程方便你直接抄作業(yè)多重檢驗校正方法用BHBenjamini-Hochberg而不是Bonferroni。Bonferroni太保守在2萬個基因的背景下幾乎不可能有基因通過。BH控制的是FDR更適合高通量篩選場景。Cox回歸的ties處理方法默認是efron這也是大多數(shù)情況下的推薦值。如果事件時間有大量重復比如按天記錄很多樣本同一天死亡可以考慮breslow但efron通常更準確。LASSO的lambda選擇用cv.glmnet()做10折交叉驗證取lambda.min或lambda.1se。lambda.1se更保守選出的基因更少模型更簡潔lambda.min選出的基因更多預測性能可能略好。我一般兩個都跑對比一下。ROC曲線的 cutoff 選擇用timeROC包做時間依賴ROC通常報告1年、3年、5年的AUC。AUC0.7算不錯0.8算很好0.9要警惕過擬合。library(glmnet) library(timeROC) # LASSO-Cox x - t(expr[top_genes, ]) y - Surv(clin$OS_time, clin$OS_status) cv_fit - cv.glmnet(x, y, family cox, nfolds 10) lasso_coef - coef(cv_fit, s lambda.min) selected_genes - rownames(lasso_coef)[lasso_coef[, 1] ! 0] # 風險評分 risk_score - predict(cv_fit, newx x, s lambda.min)[, 1] # 時間依賴ROC roc - timeROC( T clin$OS_time, delta clin$OS_status, marker risk_score, cause 1, times c(12, 36, 60), iid TRUE ) roc$AUC4.3 結果解讀與圖表輸出規(guī)范生存分析的結果通常以三類圖呈現(xiàn)KM曲線必須包含p值、風險表、分組圖例。線條顏色建議用對比明顯的配色如藍色和紅色不要用相近的綠色和青色。森林圖橫軸是HR對數(shù)刻度每個基因一個點加誤差線垂直參考線在HR1處。顯著基因用不同顏色標注。ROC曲線多條曲線畫在一起圖例標注AUC值。對角線是隨機猜測的參考線。表格方面單因素和多因素Cox的結果通常合并成一張表列包括基因名、HR、95%CI下限、95%CI上限、p值、FDR。這張表是論文附錄的標配。提示所有圖表的分辨率至少300dpiKM曲線建議保存為PDF矢量圖方便后期調整。R里用ggsave()時注意設置width和height默認尺寸往往偏小。5. 常見問題與排查技巧實錄5.1 生存曲線交叉或分不開怎么辦KM曲線交叉是生存分析里最讓人頭疼的問題之一。交叉意味著兩組之間的風險關系隨時間變化——早期高表達組風險低后期反而風險高。這種情況違反了Cox模型的PH假設log-rank檢驗的結果也不可靠。遇到交叉我的處理思路是檢查分組方式。中位數(shù)分組可能把真正的高表達和中等表達混在一起了。試試用上三分位vs下三分位或者最佳截斷值。檢查時間分段。如果交叉發(fā)生在某個時間點之后可以考慮做時間分段Cox比如前12個月和后12個月分別建模?？紤]非線性關系。有些基因的表達和風險不是線性關系而是U型的——太高太低都不好中間最安全。這種情況可以用限制性立方樣條restricted cubic spline來建模。接受它。如果交叉很輕微且log-rank p值仍然顯著可以在論文里如實報告并說明可能的原因。5.2 樣本量不足時的替代方案樣本量不足是臨床轉錄組研究的常態(tài)。如果總樣本只有40例直接做KM分組分析基本沒戲。這時候可以考慮用連續(xù)變量做Cox回歸而不是二分分組。連續(xù)變量利用了每個樣本的完整信息統(tǒng)計效能比二分變量高。用GSEA風格的生存分析把基因按表達值排序看預后相關基因是否富集在某個通路里而不是逐個基因檢驗。合并多個數(shù)據(jù)集。如果有多個小隊列可以用ComBat或removeBatchEffect做批次校正后合并但要注意平臺差異。降低顯著性閾值。在探索性分析里p0.05不校正也可以接受但要在論文里明確說明這是探索性結果需要獨立隊列驗證。5.3 多重檢驗校正后沒有顯著基因怎么辦這是高通量生存分析里最常見的結果之一。2萬個基因跑完FDR最小的可能是0.08就是過不了0.05的線。這時候不要慌幾個應對策略策略一放寬到FDR0.1或p0.01不校正。很多發(fā)表的預后模型文章用的就是p0.05不校正雖然不夠嚴謹?shù)谔剿餍噪A段可以接受。策略二先做通路水平的生存分析。把基因按通路分組計算通路評分如GSVA分數(shù)然后對通路評分做生存分析。通路數(shù)量少幾百個多重檢驗負擔小更容易得到顯著結果。策略三用機器學習方法降維。隨機生存森林random survival forest或深度學習方法可以在不依賴單基因顯著性的情況下構建預后模型。策略四換一個終點。OS可能太難分了試試PFS、RFS或DSS疾病特異性生存。有時候OS因為競爭風險非癌癥死亡的干擾信號被稀釋了。5.4 常見報錯與快速排查表報錯信息可能原因解決方法Error in Surv(time, event) : time must be numeric生存時間列是字符型用as.numeric()轉換Error in coxph: NA/NaN/Inf in foreign function call表達矩陣有缺失值或無窮值檢查并剔除含NA的基因/樣本Error in survfit: object group not found分組變量沒有正確創(chuàng)建檢查ifelse條件確認沒有NAWarning: Logrank test p-value is NA兩組樣本量差異過大或事件數(shù)太少檢查分組比例考慮換分組方式Error in ggsurvplot: could not find function包沒有加載library(survminer)Error in cv.glmnet: x should be a matrix輸入不是矩陣用as.matrix()轉換注意coxph報錯“NA/NaN/Inf”是最常見的90%的情況是表達矩陣里有NA。建議在分析前統(tǒng)一做一次expr - na.omit(expr)或者用expr[is.na(expr)] - 0填充但填充要謹慎只適合少量缺失。5.5 獨家避坑經(jīng)驗分享最后分享幾個我在實際項目里踩過的坑都是文檔里不會寫的坑一基因名大小寫不一致。表達矩陣里是“TP53”臨床注釋里是“Tp53”直接匹配就找不到。統(tǒng)一轉成大寫再匹配。坑二log轉換的底數(shù)。有的數(shù)據(jù)集已經(jīng)log2轉換過了你又轉了一次等于做了log2(log2(x1)1)表達值的分布會完全變形。做之前一定要看數(shù)據(jù)范圍——如果最大值在20以內大概率已經(jīng)log過了如果最大值上萬那就是原始值需要轉換?？尤S訪時間的單位。有一次我拿到一個數(shù)據(jù)集時間字段寫的是“days”但實際記錄的是“months”結果中位生存期算出來是3000多“天”差點鬧笑話。做之前一定要和臨床合作者確認單位?？铀淖罴呀財嘀档倪^擬合。surv_cutpoint()在小樣本里很容易找到“完美”的截斷值p值能到0.0001但換一個隊列就完全不顯著。所以最佳截斷值的結果只能作為敏感性分析不能作為主要結論?？游錍ox模型的PH假設。很多人跑完Cox直接看p值從來不檢驗PH假設。如果PH假設不成立HR的解釋就是錯的。養(yǎng)成習慣每次coxph()之后都跑一次cox.zph()?？恿L險表的人數(shù)對不上。KM曲線下方的風險表有時候數(shù)字會跳變比如從20直接跳到18中間沒有19這是因為那個時間點沒有刪失也沒有事件屬于正?，F(xiàn)象不用慌?？悠叨嘁蛩谻ox的變量選擇。不要把所有臨床變量一股腦全放進去樣本量不夠會導致過擬合。經(jīng)驗法則是每個變量至少對應10個事件。如果有30個事件最多放3個變量?？影薘OC曲線的時間點選擇。不要隨便選1年、3年、5年要根據(jù)你的中位隨訪時間來定。如果中位隨訪只有2年那5年ROC的AUC會非常不穩(wěn)定因為5年時剩下的樣本太少了?？泳磐獠框炞C的數(shù)據(jù)集選擇。如果要做外部驗證盡量選同一個平臺比如都是RNA-seq或者都是Affymetrix U133 Plus 2.0跨平臺驗證的AUC通常會掉0.1-0.2這是正常的但要在論文里說明。坑十結果的可復現(xiàn)性。R的隨機過程如交叉驗證、Bootstrap一定要設set.seed()否則每次跑的結果都不一樣審稿人復現(xiàn)不了就麻煩了。6. 從生存分析到預后模型的擴展思路6.1 多基因風險評分的構建與驗證單基因生存分析只是起點真正有臨床價值的是多基因風險評分。構建流程通常是單因素Cox篩出候選基因p0.05LASSO-Cox進一步降維得到最終基因集和系數(shù)計算每個樣本的風險評分按評分中位數(shù)分高/低風險組KM曲線和ROC曲線驗證多因素Cox確認風險評分是獨立預后因素構建列線圖nomogram整合風險評分和臨床變量列線圖是臨床醫(yī)生最喜歡的工具因為它可以把復雜的模型轉化成一個直觀的評分系統(tǒng)。用rms包的nomogram()函數(shù)就能畫。6.2 免疫浸潤與生存分析的聯(lián)合分析現(xiàn)在純生存分析的文章越來越難發(fā)審稿人通常會要求加上機制層面的解釋。一個常見的思路是把生存分析和免疫浸潤分析結合起來用CIBERSORT或ESTIMATE計算每個樣本的免疫細胞浸潤比例看預后相關基因的表達和免疫細胞浸潤的相關性看高/低風險組的免疫微環(huán)境差異用TIDE或IPS預測免疫治療響應這樣文章的故事線就完整了基因表達異?！绊懨庖呶h(huán)境→導致預后差異→可能影響免疫治療響應。6.3 單細胞數(shù)據(jù)中的生存相關基因驗證如果你有單細胞數(shù)據(jù)可以把生存分析篩出的關鍵基因放到單細胞圖譜里看它在哪些細胞類型中表達。這一步能大幅提升文章的機制深度——從“這個基因和預后相關”升級到“這個基因在腫瘤上皮細胞中高表達可能通過某種通路影響預后”。單細胞驗證的基本流程降維聚類→注釋細胞類型→看目標基因在各細胞類型的表達→看目標基因在腫瘤細胞vs正常細胞的差異→做擬時序分析看表達動態(tài)。6.4 藥物敏感性預測與臨床轉化最后一步是臨床轉化。用GDSC或CTRP數(shù)據(jù)庫看高風險組對哪些藥物更敏感。常用的方法是pRRophetic包它基于基因表達預測IC50。如果高風險組對某種化療藥或靶向藥更敏感那這個風險評分就有了指導用藥的潛在價值。當然這一步只是預測離真正的臨床應用還有很長的路。但作為文章的最后一塊拼圖它能顯著提升研究的轉化意義。我個人在實際操作中的體會是生存分析看起來簡單——不就是畫個KM曲線、跑個Cox嗎——但真正做好、做扎實需要關注的細節(jié)非常多。從樣本ID匹配到PH假設檢驗從分組方式選擇到多重檢驗校正每一步都可能影響最終結論。我見過太多文章因為生存分析的某個細節(jié)沒處理好被審稿人質疑結果可靠性。所以寧可前面多花時間把數(shù)據(jù)清洗干凈、把假設檢驗做全也不要等到投稿了再回頭補。另外生存分析的結論一定要有外部驗證或至少內部Bootstrap驗證單隊列的結果永遠只是“關聯(lián)”不是“因果”這一點在解讀時務必保持清醒。