本章思維導(dǎo)圖：

近兩年RCS和Nomogram一樣炙手可熱(泛濫成災(zāi))，有必要做一個更專業(yè)更SCI化的專業(yè)梳理。出于此目的，匯總心得如下。

前言

臨床連續(xù)數(shù)據(jù)與結(jié)局的線性關(guān)聯(lián)方式與模型優(yōu)化是一個重要問題。醫(yī)學(xué)常見非線性曲線形態(tài)如下圖，有S L log -log J threshold U ∩ 形態(tài)，數(shù)據(jù)擬合越接近真實關(guān)聯(lián)曲線形態(tài)才越利于提高預(yù)測性能。統(tǒng)計學(xué)家George Box 的名言說到 “all models are wrong, but some are useful”。意思是，不能說某個構(gòu)建的模型100%正確，只可以說它看起來似乎有點用。這幾年熱門的RCS(限制性立方樣條)就是這樣一種方法，它可以為模型擬合提供線索, 讓模型看起來更好用，讓數(shù)據(jù)分析更亮眼。

我們經(jīng)常通過線性回歸、Logistic回歸、Cox回歸探討連續(xù)自變量和因變量的關(guān)系，但上述模型要求合適的鏈接函數(shù)g(y)，保證自變量和因變量服從一般/廣義線性前提條件。其中線性前提條件轉(zhuǎn)換可參考后續(xù)轉(zhuǎn)換章節(jié)，例如正態(tài)轉(zhuǎn)換、logitP轉(zhuǎn)換。這里需要注意，上述假設(shè)不成立時，我們需要采用什么樣的方法探討自變量和因變量是線性還是非線性關(guān)聯(lián)呢？平滑函數(shù) 非線性關(guān)系可以平滑工具來實現(xiàn)。平滑工具主要包括核平滑、多項式平滑和樣條平滑。核平滑是使用一組局部權(quán)重來生成平滑的估計，在實際應(yīng)用中很難實現(xiàn)；多項式是最簡單的函數(shù)平滑結(jié)構(gòu)，但只是一種粗略的估計方法，存在過擬合、共線性、全局性等諸多問題。樣條(spline)本質(zhì)上是一個分段多項式, 但它一般要求每個分段點上連續(xù)并且二階可導(dǎo)，這樣可以保證曲線的平滑性，如下圖雞蛋曲面knot可導(dǎo)。Spline適用條件為①數(shù)據(jù)x與y關(guān)系不符合直線或廣義線性前提 ②數(shù)據(jù)多項式回歸R2低③某個節(jié)點前后趨勢發(fā)生明顯改變。

RCS

RCS(Restricted Cubic Spline或Natural cubic splines )。立方樣條本質(zhì)上是連續(xù)平滑的分段三次多項式。其中“分段pieces”的數(shù)量由使用的“節(jié)點knots”數(shù)量決定。每個節(jié)點內(nèi)實際上是一個三次多項式。經(jīng)典的分段回歸使每個段的內(nèi)部效應(yīng)被強制統(tǒng)一, 在節(jié)點的位置跳躍,“瞬時變化”不合理， 這不但不符合很多實際情況, 而且不能發(fā)現(xiàn)最大值和最小值的點。RCS的數(shù)學(xué)原理實際上比每個節(jié)點擬合三次多項式要復(fù)雜一些，需要施加進一步限制以便spline是連續(xù)平滑沒有間隙。簡而言之，RCS實質(zhì)上是通過選擇節(jié)點的位置和數(shù)量，擬合樣條函數(shù)RCS（X），使得連續(xù)變量X在整個取值范圍內(nèi)呈現(xiàn)光滑的曲線，如圖。Restricted RCS在回歸樣條的基礎(chǔ)上附加Restricted：樣條函數(shù)在自變量數(shù)據(jù)范圍左右兩端的兩個區(qū)間內(nèi)為線性函數(shù)。因為左右兩端預(yù)測的區(qū)間非常寬加一個線性函數(shù)作為約束條件，這樣使得兩邊的預(yù)測更為準(zhǔn)確。由于加了約束條件，正常情況下會導(dǎo)致各段區(qū)間內(nèi)RCS擬合效果稍差于分段回歸，但是總體趨勢上會更加合理。

RCS通過樣條函數(shù)RCS（X）轉(zhuǎn)換自變量X后，然后根據(jù)因變量的分布類型選擇合適的鏈接函數(shù)，進而擬合模型g（Y）=常數(shù)項+RCS（X）+其他自變量，其中g(shù)為鏈接函數(shù)。樣條函數(shù)RCS（X）包括一個線性項X以及K-2個立方項（S），即RCS（X）=β0X+β1S1+…+β（k-2）S（k-2），完整公式見下圖。Ci(x)是落在第i個節(jié)點中的三次分量，g是所謂的鏈接函數(shù)。

其中，節(jié)點的位置對樣條函數(shù)的擬合影響不大，一般根據(jù)連續(xù)變量的百分位數(shù)選擇。節(jié)點數(shù)量對樣條函數(shù)的影響較大，節(jié)點的數(shù)量決定曲線的形狀，當(dāng)節(jié)點的數(shù)量為2時，得到的擬合曲線就是一條直線。研究顯示，節(jié)點數(shù)量為3~5時樣條函數(shù)擬合較好，一般推薦knot=4。

表1 harrell 推薦knot 3-7

RCS形狀受 knots影響很大，knot越多曲線約復(fù)雜。由此可見，RCS 只能用來判斷非線性趨勢，為切點選擇提供線索，而不能精確提供切點！

SCI實戰(zhàn)中，BMJ、EI文章對于RCS非線性關(guān)系描述的非常好，值得參考。例如，BMJ文章截取了其中的脂肪量和死亡風(fēng)險的曲線圖，作者采用中位數(shù)設(shè)置為拐點，從圖中可以看出脂肪量在<21Kg，死亡風(fēng)險基本保持不變，而>21Kg以后，死亡風(fēng)險隨著脂肪量增大而逐漸增高。通過切點前后的HR95%CI,側(cè)面說明研究結(jié)論。通過繪制RCS圖，既直觀又能輸出p for non-linearity值用來全面描述x與y之間的非線性關(guān)系，還可發(fā)現(xiàn)潛在的有趣切點，可以為文章增色不少。RCS對切點的線索，多是采用knot 3-8個且視覺上95%CI最狹窄處，但是很多文章都沒有交代，似乎很多都默認(rèn)了。R包會自動給出對應(yīng)的切點位置。通過SCI分析，我們也發(fā)現(xiàn)knot對應(yīng)的分位數(shù)也并不是一成不變的，如BMJ文獻knot=4 (reference value= median，percentage= P5 P35 P65 P95)，EI文獻knot=4 ?(reference value= P10，percentage= P20 P40 P60 P80)。主流RCS表達套路如下：

① RCS+ HR左右開弓，切點median

② RCS+ 直方圖，切點自定義

③ RCS+ 密度圖，切點自定義

發(fā)散思維1：選擇分位數(shù)還是RCS

連續(xù)數(shù)據(jù)，我們往往第一選擇會依據(jù)臨床共識進行分類，常見的如二分類法（如高血壓 SBP≥140mmhg。有時關(guān)注變量并沒有公認(rèn)的分類切點，此時我們會選擇中位數(shù)或多分類/等級分類(三分位、四分位)。這些分類雖然嚴(yán)格遵循臨床意義或分位數(shù)原則，但是可能削弱模型預(yù)測能力。忽略連續(xù)數(shù)據(jù)非線性形態(tài)，而直接引入分類進行模型預(yù)測，顯然容易丟失許多信息，如下圖，具體測試可在https://drjgauthier.shinyapps.io/spliny/ 中自行感受。簡而言之，能RCS就先RCS, RCS后再根據(jù)圖形定性形態(tài)，選擇對應(yīng)的knots點、臨床界值、ROC界值、中位數(shù)、分位數(shù)、segmented包分段對應(yīng)切點，進行分析。

發(fā)散思維2：其他非線性擬合

RCS很強，但是非線性擬合過程也有其他方法，https://www.r-bloggers.com/2014/09/an-exercise-in-non-linear-modeling/ 對主流的幾類非線性方法進行了比較。為了避免過度擬合，根據(jù)AIC最小值標(biāo)準(zhǔn)選擇模型。結(jié)論是RCS的AIC最小，擬合效果最佳。所以一句話，看不明白RCS準(zhǔn)沒錯。B樣條，別名基本樣條，是RCS的常見替代品，后者也使用節(jié)點來控制靈活性。與RCS區(qū)別在于左右末端不受限制，因此它們比受限制的三次樣條具有更靈活的尾部。廣義加性模型 (GAM) 是廣義線性模型的擴展，專門研究非線性關(guān)系，GAM用途很廣。所用懲罰參數(shù)，薄板回歸樣條、三次回歸樣條、P 樣條(懲罰B樣條)。GAM不光用來為非線性關(guān)系提供線索，更重要的是直接參與GAM模型構(gòu)建，尤其是環(huán)境研究領(lǐng)域常用，詳細(xì)查看相關(guān)章節(jié)。

表 spline 平滑相關(guān)R包

R軟件代碼

SAS在unconditional logistic, pooled logistic, conditional logistic, 甚至是 proportional hazards regression models里均可實現(xiàn)RCS，但是繪圖結(jié)果較丑陋，不推薦！這里如果感興趣的同學(xué)深入研究，也會發(fā)現(xiàn)SAS 提供的 p for non-linearity 值定義是 P >0.05 為非線性關(guān)聯(lián)，與R提供的P<0.05 為非線性關(guān)聯(lián)相反，一定要注意這個坑。STATA也可以跑RCS，也很便捷，但是依舊存在圖形需要轉(zhuǎn)換再編輯問題，也就不推薦了。由于日程使用，臨床數(shù)據(jù)大多結(jié)局為Y=01, 考慮OR或HR為效應(yīng)指標(biāo)。本程序僅提供R語言logistic 與 cox主流代碼。

原創(chuàng)不易，轉(zhuǎn)載請說明來自本公眾號。

本文使用 文章同步助手 同步