孟詳喻，靳英輝，程小柯，張永剛，曹越，曾憲濤

? 循證理論與實踐 ?

生存資料的二次研究系列之一：應(yīng)用R軟件實現(xiàn)生存曲線中數(shù)據(jù)提取與生存曲線的合并

孟詳喻1，靳英輝2，程小柯3，張永剛4，曹越1，曾憲濤1

涉及生存資料的Meta分析最常用的方法是對兩組比較的風(fēng)險比（hazard ratio，HR）進(jìn)行合并，其數(shù)據(jù)基礎(chǔ)是原始研究的HR及其95%可信區(qū)間（95%CI）。然而，原始研究中HR及其95%CI未報道或不完整報道的情況并不少見，這限制了基于HR的Meta分析的開展。絕大多數(shù)涉及生存結(jié)局的原始研究均提供了生存曲線，通過對生存曲線進(jìn)行合并能對人群的生存結(jié)局進(jìn)行全面評估。本文介紹利用R軟件實現(xiàn)生存曲線數(shù)據(jù)提取并對其進(jìn)行合并以生成匯總生存曲線的方法。

生存資料；生存曲線；R軟件

涉及生存資料的Meta分析最常用的方法是以倒方差法對（hazard ratio，HR）進(jìn)行合并［1，2］。該類型的Meta分析主要著眼于兩組比較的評價，能區(qū)分優(yōu)劣但很難反映全局問題。此外，原始研究中

HR及其95%可信區(qū)間（95%CI）未報道或不完整報道的情況并不少見，這限制了基于HR的Meta分析的開展［3，4］。另一方面，對特定人群的生存曲線進(jìn)行合并能在擴(kuò)大樣本量的基礎(chǔ)上全面反映出研究人群的特定生存結(jié)局的總體狀況。絕大多數(shù)涉及生存結(jié)局的原始研究均提供了生存曲線，從而為此類研究的開展提供了便利［5］。國外已有方法學(xué)文獻(xiàn)報道了通過乘積限估計法計算合并生存率的方法，并開發(fā)了相關(guān)的R軟件包以支持其應(yīng)用［6］。本文擬通過實例介紹利用R軟件digitize程序包對生存曲線進(jìn)行標(biāo)準(zhǔn)化提取的方法，并介紹利用MetaSurv程序包對自曲線提取出的生存數(shù)據(jù)進(jìn)行合并以生成匯總后的生存曲線估計。所使用R軟件的版本為3.1.3。

1 程序包簡介

因程序包研發(fā)者及維護(hù)者未能及時回復(fù)CRAN（The Comprehensive R Archive Network）工作組的聯(lián)系，digitize程序包現(xiàn)已自CRAN存儲空間移除。我們自GitHub獲取該程序包最近一次更新的源代碼并使用RStudio重新打包后經(jīng)檢測可在最新版本的R軟件成功安裝并加載運行。MetaSurv程序包目前也無法通過R官網(wǎng)途徑獲取。讀者可自武漢大學(xué)中南醫(yī)院循證與轉(zhuǎn)化醫(yī)學(xué)中心的官網(wǎng)（http://www.whuzncebtm.com/）中的“科學(xué)研究→資料下載”欄中下載這兩個程序包的zip格式的安裝文件，安裝方法詳見曾憲濤主編的著作［7］。

由于digitize程序包的運行需調(diào)用jpeg程序包，因此在安裝digitize程序包之前，首先安裝及加載jpeg程序包，命令如下：

該程序包安裝完成后即可使用下述命令加載digitize程序包和MetaSurv程序包：

2 利用digitize程序包自K-M曲線提取數(shù)據(jù)

2.1提取前準(zhǔn)備 提取前的主要準(zhǔn)備工作包括自原始文獻(xiàn)電子文檔截取生存曲線，并利用Photoshop軟件添加等分輔助線，建議每條輔助線在時間軸的位置對應(yīng)至單月，具體方法在此不做贅述。本處所選用的原始文獻(xiàn)為2013年發(fā)表于《Lancet Oncology》雜志上的一項比較泊馬度胺聯(lián)合低劑量地塞米松與高劑量地塞米松單藥治療復(fù)發(fā)及難治性多發(fā)性骨髓瘤患者的Ⅲ期隨機(jī)對照試驗，所選取的生存曲線為比較兩種藥物無進(jìn)展生存期（PFS）的K-M曲線［8］。

3 生存曲線的合并

圖3　對生存曲線進(jìn)行取點

表1　提取出的生存數(shù)據(jù)列表

我們通過一實例操作介紹使用MetaSurv程序包合并自不同單個研究的提取出的生存數(shù)據(jù)并生成匯總估計的生存曲線。具體方法學(xué)模型原理有興趣的讀者請參閱相關(guān)方法學(xué)文獻(xiàn)［6］。

本實例所納入的研究為5項以復(fù)發(fā)與難治性多發(fā)性骨髓瘤患者為對象的對照組為地塞米松或地塞米松聯(lián)合安慰劑的隨機(jī)對照試驗［8-12］。按照前述方法提取各對照組各月的總生存（OS）K-M曲線生存數(shù)據(jù)。整理好的原始數(shù)據(jù)表見表2。數(shù)據(jù)表表頭“Study”為研究編號，研究1至5分別對應(yīng)2005 Richardson、2007 Dimopoulos、2007 Weber、2011 Kropff與2013 San-Miguel（發(fā)表年份+第一作者英文名）；“Time”為時間點（月份）；“Survival”為各時間點所對應(yīng)的生存率；“NbRisk”為各時間段事件風(fēng)險人數(shù)，該數(shù)據(jù)自文獻(xiàn)中提取或根據(jù)Thiery等的方法學(xué)假定截尾率在隨訪期內(nèi)恒定而計算得出［13］。

表2　示例5項研究的數(shù)據(jù)

接下來使用下述代碼進(jìn)行合并分析：

命令運行后，結(jié)果如表3所示。

結(jié)果中包含了固定效應(yīng)模型與隨機(jī)效應(yīng)模型下匯總后的平均生存時間與中位生存時間及其95%可信區(qū)間，還顯示反應(yīng)異質(zhì)性的統(tǒng)計量Cochran Q、H2與I2值。

從表3還可以看出，Cochran Q=186.248778，H2= 1.565116，I2= 36.106963%。根據(jù)Higgins等的方法學(xué)還可對異質(zhì)性的統(tǒng)計顯著性進(jìn)行檢測，命令代碼如下：

命令運行后，結(jié)果如下：

由結(jié)果可知，H值可信區(qū)間下限大于1 （LL=1.222318），檢測結(jié)果有統(tǒng)計學(xué)顯著性，故后續(xù)合并后生存曲線的繪制采用隨機(jī)效應(yīng)模型的數(shù)據(jù)。

輸入以下命令代碼進(jìn)行后續(xù)操作：

表3　示例的計算結(jié)果

圖4　根據(jù)各研究所提取數(shù)據(jù)生成的相應(yīng)曲線及合并后的OS曲線與95%可信區(qū)間

命令運行后，最終生成的圖形如圖4所示。圖正中藍(lán)色粗實線為合并后的OS曲線，該曲線兩側(cè)的紅色點線為其95%可信區(qū)間上下限。其余線條為根據(jù)單個研究所提取的生存數(shù)據(jù)得出的近似曲線，各曲線末端標(biāo)注了其隨訪終點。

4 小結(jié)

HR值等關(guān)鍵數(shù)據(jù)的報道缺失是生存資料Meta分析過程中常常遇見的問題，自K-M曲線提取生存數(shù)據(jù)后根據(jù)相關(guān)方法學(xué)文獻(xiàn)進(jìn)行估算是解決這類問題較為可靠的方法之一［14］。本文詳細(xì)介紹了利用R軟件digitize程序包進(jìn)行生存曲線數(shù)據(jù)提取的操作方法。該方法通過添加等分輔助線以單位時間間隔提取數(shù)據(jù)保證了數(shù)據(jù)質(zhì)量及后續(xù)計算結(jié)果的穩(wěn)定可靠。但該方法需逐一手動取點，相對費時。

在生存資料匯總分析過程中，單純通過HR合并雖然能區(qū)分優(yōu)劣但很難反映全局問題。對某些治療組的生存曲線進(jìn)行合并能提供很多其他有意義的信息。本文介紹了使用R軟件MetaSurv程序包實現(xiàn)生存曲線數(shù)據(jù)合并的方法，并通過一實例進(jìn)行了具體介紹。我們希望這一方法能為今后更多的研究所采用，既是對HR分析很好的補(bǔ)充，更有助于對生存資料進(jìn)行更為全面的反映。

［1］ Spruance SL，Reid JE，Grace M，et al. Hazard ratio in clinical trials［J］. Antimicrob Agents Chemother，2004，48（8）: 2787-92.

［2］ 曾憲濤. 系統(tǒng)評價/Meta分析理論與實踐［M］. 北京: 軍事醫(yī)學(xué)科學(xué)出版社，2012:117-9.

［3］ Michels S，Piedois P，Burdett S，et al. Meta-analysis when only the median survival times are known: a comparison with individual patient data results［J］. Int J Technol Assess Health Care，2005，21（1）: 119-25.

［4］ Parmar MKB，Torri V，Stewart L. Extracting summary statistics to perform meta-analyses of the published literature for survival endpoints［J］. Stat Med，1998，17（24）: 2815-34.

［5］ Arends L，Hunink M，Stijnen T. Meta-analysis of summary survival curve［J］. Stat Med，2008，27（22）:4381-96.

［6］ Combescure C，F(xiàn)oucher Y，and Jackson D. Meta-analysis of singlearm survival studies: a distribution-free approach for estimating summary survival curves with random effects［J］. Stat Med，2014，33（15）:2521-37.

［7］ 曾憲濤，張超. R與Meta分析［M］. 北京: 軍事醫(yī)學(xué)科學(xué)出版社，2015:26.

［8］ San Miguel J，Weisel K，Moreau P，et al. Pomalidomide plus lowdose dexamethasone versus high-dose dexamethasone alone for patients with relapsed and refractory multiple myeloma （MM-003）: a randomised， open-label， phase 3 trial［J］. Lancet Oncol， 2013， 14（11）: 1055-66.

［9］ Kropff M，Baylon HG，Hillengass J，et al. Thalidomide versus dexamethasone for the treatment of relapsed and/or refractory multiple myeloma: results from OPTIMUM，a randomized trial［J］. Haematologica，2012，97（5）: 784-91.

［10］ Weber DM，Chen C，Niesvizky R，et al. Lenalidomide plus dexamethasone for relapsed multiple myeloma in North America［J］. N Engl J Med，2007，357（21）: 2133-42.

［11］ Dimopoulos M，Spencer A，Attal M，et al. Lenalidomide plus dexamethasone for relapsed or refractory multiple myeloma［J］. N Engl J Med，2007，357（21）:2123-32.

［12］ Richardson PG，Sonneveld P，Schuster MW，et al. Bortezomib or highdose dexamethasone for relapsed multiple myeloma［J］. N Engl J Med，2005，352（24）:2487-98.

［13］ Tierney JF，Stewart LA，Ghersi D，et al. Practical methods for incorporating summary time-to-event data into meta-analysis［J］. Trials， 2007， 8: 16.

［14］ Earle C，Wells G. An Assessment of Methods to Combine Published Survival Curves［J］. Medical Decision Making，2000，20:104-11.

本文編輯：田國祥

Extracting data from survival curve and pooled survival curves by using R software

MENG Xiang-yu*，JIN Ying-hui， CHENG Xiao-ke， ZHANG Yong-gang， CAO Yue， ZENG Xian-tao. *Center for Evidence-Based and Translational Medicine， Zhongnan Hospital， Wuhan University， Wuhan 430071， China.

ZENG Xian-tao， E-mail: zengxiantao1128@163.com

The most common method in Meta-analysis involving survival data is pooled hazard ratios （HR）， and the data base is HR in original studies and 95% confidence interval （95%CI）. However， non-reports or incomplete reports of HR in original studies and 95%CI are not uncommon， which limits the conducting of Meta-analysis based on HR. The most of original studies related to survival outcomes can provide survival curves， and the survival outcomes can be completely reviewed through pooling survival curve in the study population. This paper introduces the approach of extracting and pooling survival curve data for generating a pooled survival curve by using R software.

survival data； survival curve； R software

至此生存數(shù)據(jù)提取過程結(jié)束，可將所提取的數(shù)據(jù)用于后續(xù)分析處理。

R4

A

1674-4055（2016）01-0002-05

1430071 武漢，武漢大學(xué)中南醫(yī)院循證與轉(zhuǎn)化醫(yī)學(xué)中心；2300193 天津，天津中醫(yī)藥大學(xué)護(hù)理學(xué)院循證護(hù)理中心；3475000開封，河南大學(xué)淮河醫(yī)院循證醫(yī)學(xué)中心；4610041 成都，四川大學(xué)華西醫(yī)院期刊社

10.3969/j.issn.1674-4055.2016.01.01

曾憲濤，E-mail:zengxiantao1128@163.com.

準(zhǔn)備好的圖像文件如圖1所示。將其命名為“PFS_2013Miguel.jpg”，并存儲于桌面“Rwork”文件夾中（注：可通過代碼“＞setwd （ ）”設(shè)定工作文件夾，括號中為其路徑）。

圖1提取前預(yù)處理好的圖像文件

2.2數(shù)據(jù)提取 運行R軟件，輸入以下命令代碼進(jìn)行數(shù)據(jù)處理：

library（digitize）

cal＜-ReadAndCal（'PFS_2013Miguel.jpg'）

第二行命令為讀取目標(biāo)圖像文件“PFS_2013Miguel.jpg”并在圖上選定橫軸X與縱軸Y坐標(biāo)標(biāo)定點，下圖中橫軸兩標(biāo)記點與縱軸兩標(biāo)記點分別對應(yīng)Xmin、Xmax、Ymin及Ymax，標(biāo)定點的信息存入對象“cal”中。程序運行界面如圖2。

圖2讀取目標(biāo)圖像文件并在橫軸與縱軸選取標(biāo)定點

該命令是對目標(biāo)曲線（本實例中對紅色的地塞米松組曲線提取數(shù)據(jù)）取點提取各月所對應(yīng)的生存率數(shù)據(jù)，“col='blue'”意指以藍(lán)色空心圓點標(biāo)記各數(shù)據(jù)點，結(jié)果數(shù)據(jù)記入對象curv中，如圖3所示。取點完畢后單擊右鍵，出現(xiàn)彈窗選擇“stop”結(jié)束取點。

該行命令是對數(shù)據(jù)進(jìn)行標(biāo)定，curv代表上一步取點所得的數(shù)據(jù)集，cal即標(biāo)定點集合，最后四個數(shù)字分別代表標(biāo)定點Xmin、Xmax、Ymin、Ymax的坐標(biāo)值，將數(shù)據(jù)標(biāo)定后的結(jié)果記為“data”。

該行命令為將數(shù)據(jù)集data中的x欄數(shù)值取整，記為“Time”即時間（月）。

該行命令為將d a t a中的y欄數(shù)值記為“Survival”即生存率。

該行命令為將“Time”與“Survival”組合為新的數(shù)據(jù)集，命名為“PFS_2013Miguel”。

該行命令為將“PFS_2013Miguel”導(dǎo)出為csv格式的文件，命名為“PFS_2013Miguel.csv”，默認(rèn)存儲于Rwork文件夾。該文件可用Excel打開，如表1所示。