臨床預測模型也是大家比較感興趣的，今天就帶著大家看一篇臨床預測模型的文章，並且用一個例子給大家過一遍做法。
這篇文章來自護理領域頂級期刊的文章，文章名在下面
Ballesta-Castillejos A, Gómez-Salgado J, Rodríguez-Almagro J, Hernández-Martínez A. Development and validation of a predictive model of exclusive breastfeeding at hospital discharge: Retrospective cohort study. Int J Nurs Stud. 2021 May;117:103898. doi: 10.1016/j.ijnurstu.2021.103898. Epub 2021 Feb 7. PMID: 33636452.
文章作者做了個出院時單純母乳餵養的預測模型，資料來自兩個佇列，一個佇列做模型，另外一個用來驗證。樣本量的計算用的是“10個樣本一個變數”的標準，預測結局變數是個二分類，變數篩選的方法依然是單因素有意義的都納入預測模型中，具體的建模方法是logistic迴歸模型。然後用AUC進行模型的評估。
在結果報告上，作者報告了所有的有意義的預測因素，每個因素會展示OR和OR的置信區間，還有整體模型的評價指標，包括展示了模型的ROC曲線，和AUC，作者還展示了不同概率階截斷值下的模型的Sensitivity， Specificity， PPV， NPV， LR+， LR-。如圖：
作者整個文章是用SPSS做出來的，今天給大家寫寫論文中的各種指標都是什麼意義以及如何用R語言做出來論文中需要報告的各種指標。
理論鋪墊
首先給大家寫靈敏度（Sensitivity）與特異度（Specificity），這兩個東西都是針對二分類結局來講的，大家先瞅瞅下面的圖：
我們真實的結果有兩種可能，模型的預測也有兩種可能，上圖的AD表示模型預測對的個案數量，那麼靈敏度就是：在真的有病了，你的模型有多大可能檢驗出來，表示為
Sensitivity: A/(A+C) × 100
在論文中就是這個母親真的是單純母乳餵養的，模型有多大可能識別為真的單純母乳餵養。
特異度就是：我是真的真的沒病，你的模型有多大可能說我真的沒病，表示為：
Specificity: D/(D+B) × 100
在論文中就是這個母親真的不會去單純母乳餵養的，模型有多大可能識別為真的不單純母乳餵養。
有些同學說，我知道個模型預測準確率不就好了嗎，用(A+D)/(A+B+C+D)來評估模型不就好了嗎？搞這麼麻煩。。
不能這麼想的，比如你現在有一個傻瓜模型，這個模型傻到它只會將所有的人都預測為沒病，剛好這個模型被用在了一個正常人群中，然後我們發現這個傻瓜模型的正確率也是100%，這個就很離譜，所以並模型預測準確性是不能全面評估模型表現的，需要藉助Sensitivity， Specificity。
我們再看 PPV， NPV這兩個指標：
Positive Predictive Value: A/(A+B) × 100
Negative Predictive Value: D/(D+C) × 100
看上面的公式，相信大家都看得出這兩個其實就是模型的陽性預測準確性和陰性預測準確性，也可以從特定角度說明模型的表現。
再看LR+和 LR-，這兩個就是陽性似然比 (positive likelihood ratio, LR+)和 陰性似然比
(positive likelihood ratio, LR+)，似然比的概念請參考下一段英文描述：
Likelihood ratio (LR) is the ratio of two probabilities: (i) probability that a given test result may be expected in a diseased individual and (ii) probability that the same result will occur in a healthy subject.
那麼：(LR+) = sensitivity / (1 - specificity)，意思就是真陽性率與假陽性率之比。說明模型正確判斷陽性的可能性是錯誤判斷陽性可能性的倍數。比值越大，試驗結果陽性時為真陽性的概率越大。
(LR-) = (1 - sensitivity) / specificity，意思就是假陰性率與真陰性率之比。表示錯誤判斷陰性的可能性是正確判斷陰性可能性的倍數。其比值越小，試驗結果陰性時為真陰性的可能性越大。
所以大家記住：陽性似然比越大越好，陰性似然比越小越好。
最後再把上面的所有的內容總結一個表獻給可愛的粉絲們，嘿嘿。下面就是一個分類結局預測變數需要報告的一些模型評估指標：
再回過頭想想我們所謂的陽性或者陰性，如果用logistics迴歸做的話本身這個陽性陰性的判別都是可以設定的，因為我們的模型擬合出來的是響應概率，就是Logit公式裡面的那個p值，你可以以p=0.5為我們判別陰陽性的cutoff，當然你還可以以0.9或者0.1為cutoff，cutoff不同自然我們模型靈敏度和特異度就不同了，就是說靈敏度和特異度是隨著cutoff不同而變化著的，所以要穩定地評估模型表現還需要另外的指標，這個時候我們就引出了一個很重要的概念：ROC曲線和曲線下面積AUC。
在實戰中理解ROC曲線
我現在手上有資料如下：
我要做一個default的預測模型，default是個二分類變數，取值為“No” 和“Yes”，為了簡單我預測因素只考慮一個balance。於是我建立一個logistics模型:
model_glm = glm(default ~ balance, data = default_trn, family = "binomial")我們將預測為Yes的概率為0.1的時候作為cutoff值，劃分預測結局（p<0.1的時候為No，p>0.1的時候為Yes），同樣地我們還可以將cutoff設定為0.5，0.9。然後我們分別看一看模型的靈敏度和特異度：
test_pred_10 = get_logistic_pred(model_glm, data = default_tst, res = "default",
pos = "Yes", neg = "No", cut = 0.1)
test_pred_50 = get_logistic_pred(model_glm, data = default_tst, res = "default",
pos = "Yes", neg = "No", cut = 0.5)
test_pred_90 = get_logistic_pred(model_glm, data = default_tst, res = "default",
pos = "Yes", neg = "No", cut = 0.9)
test_tab_10 = table(predicted = test_pred_10, actual = default_tst$default)
test_tab_50 = table(predicted = test_pred_50, actual = default_tst$default)
test_tab_90 = table(predicted = test_pred_90, actual = default_tst$default)
test_con_mat_10 = confusionMatrix(test_tab_10, positive = "Yes")
test_con_mat_50 = confusionMatrix(test_tab_50, positive = "Yes")
test_con_mat_90 = confusionMatrix(test_tab_90, positive = "Yes")
metrics = rbind(
c(test_con_mat_10$overall["Accuracy"],
test_con_mat_10$byClass["Sensitivity"],
test_con_mat_10$byClass["Specificity"]),
c(test_con_mat_50$overall["Accuracy"],
test_con_mat_50$byClass["Sensitivity"],
test_con_mat_50$byClass["Specificity"]),
c(test_con_mat_90$overall["Accuracy"],
test_con_mat_90$byClass["Sensitivity"],
test_con_mat_90$byClass["Specificity"])
)
rownames(metrics) = c("c = 0.10", "c = 0.50", "c = 0.90")
metrics執行程式碼後得到響應概率不同cutoff值的情況下模型的靈敏度和特異度，如下圖：
可以看到我們設定的響應概率的cutoff值不同，就是判斷陰陽的標準不同，我們得到的模型的靈敏度和特異度就是不同的。
以上只是為了再次給大家直觀地說明我們的模型的靈敏度和特異度是取決於我們的響應概率界值的，你判斷陰陽的標準會直接影響模型的靈敏度和特異度，於是，我們換個想法，對於一個模型，我們將靈敏度作為橫座標，特異度作為縱座標，然後cutoff隨便取，我們形成一條曲線，這就考慮了所有的cutoff情況了，就可以穩定地評估模型的表現了，這條曲線就是ROC曲線。
那麼我們期望的是一個模型它的靈敏度高的時候，特異度也能高，體現到曲線上就是曲線下的面積能夠越大越好。
預測模型R語言實操
此部分給大家寫如何做出論文中的各種指標，以及如何繪出ROC曲線。
依然是用上一部分的資料，依然是做balance預測default的模型：
model_glm = glm(default ~ balance, data = default_trn, family = "binomial")模型輸出結果如下：
結果中有輸出模型的截距和balance的β值，我們可以用如下程式碼得到balance的OR值以及置信區間：
exp(cbind(OR = coef(model_glm), confint(model_glm)))執行上面的程式碼就可以得到balance的OR值和OR的置信區間：
同時我們有原始的真實值，我們模型擬合好了之後可以用該模型進行預測，得到預測值，形成混淆矩陣：
model_glm_pred = ifelse(predict(model_glm, type = "response") > 0.5, "Yes", "No")在矩陣中就可以得到哪些是原始資料真實的No和Yes，哪些是模型預測出來的No和Yes：
上面就是我們自己資料做出來的混淆矩陣，然後大家就可以直接帶公式計算出需要報告的模型的Sensitivity， Specificity， PPV， NPV， LR+， LR-了。
同時大家可以用pROC包中的roc函式一行程式碼繪製出ROC曲線並得到曲線下面積，比如我做的模型，寫出程式碼如下：
test_roc = roc(default_tst$default ~ test_prob, plot = TRUE, print.auc = TRUE)就可以得到ROC曲線和曲線下面積了：
上面的所有操作都是可以在SPSS中完成的，論文作者也是用SPSS做的。大家感興趣去閱讀原論文哈。
小結
今天結合發表文章給大家寫了分類結局的預測模型需要報告哪些指標，指標的意義以及如何用R做一個分類結局的預測模型，感謝大家耐心看完，自己的文章都寫的很細，重要程式碼都在原文中，希望大家都可以自己做一做，請轉發本文到朋友圈後私信回覆“資料連結”獲取所有資料和本人收集的學習資料。如果對您有用請先記得收藏，再點贊分享。
也歡迎大家的意見和建議，大家想了解什麼統計方法都可以在文章下留言，說不定我看見了就會給你寫教程哦，有疑問歡迎私信。R資料分析：跟隨top期刊手把手教你做一個臨床預測模型