Excel軟件結(jié)合Visual Basic編程工具探索MicroRNA-501對(duì)肺鱗癌患者生存情況的影響

黃致遠(yuǎn)，江文豪，楊福貴，朱 磊，董 林，汪令偉，馮明昊，魏 亮，王光學(xué)＊，李欽傳,＊

1.同濟(jì)大學(xué)附屬東方醫(yī)院心胸外科（上海 200120）

2.同濟(jì)大學(xué)附屬東方醫(yī)院轉(zhuǎn)化醫(yī)學(xué)研究中心（上海 200120）

非小細(xì)胞肺癌（Non-small cell lung cancer, NSCLC）是肺癌的主要類型，約占所有病例的80%[1]。NSCLC包括幾個(gè)組織學(xué)亞型，如腺癌、鱗狀細(xì)胞癌和大細(xì)胞癌。microRNA（miRNAs）是19-25個(gè)核苷酸長(zhǎng)的單鏈非編碼RNA，通常通過(guò)與3’非翻譯區(qū)域(3’-UTR)結(jié)合在轉(zhuǎn)錄后水平負(fù)調(diào)控基因表達(dá)[2-3]。人們發(fā)現(xiàn)在肝細(xì)胞癌和前列腺癌中miR-501-3p可能起抑癌作用[4-5]。此外，miRNA微陣列研究表明，miR-501-3p與肺癌骨轉(zhuǎn)移相關(guān)[6]。生存分析是探究microRNA表達(dá)情況對(duì)腫瘤患者生存時(shí)間影響大小的常用方法。目前科研工作者在需要使用生存分析相關(guān)操作的時(shí)候多使用SPSS等專業(yè)的統(tǒng)計(jì)軟件[7-8]，其缺點(diǎn)在于專業(yè)性高、操作性差。Excel軟件由美國(guó)的Microsoft公司開(kāi)發(fā)，在全世界普遍使用，可以實(shí)現(xiàn)一般的t檢驗(yàn)，F(xiàn)檢驗(yàn)，線性回歸和方差分析等常見(jiàn)統(tǒng)計(jì)方法[9-12]。本文旨在利用Excel軟件中自帶的Visual Basic（VBA）編程功能，分析肺鱗癌患者組織中的miR-501表達(dá)水平與患者生存時(shí)間之間的關(guān)系。

1 資料與方法

1.1 獲取數(shù)據(jù)

登錄https://kmplot.com/analysis/網(wǎng)站，點(diǎn)擊start miRPower for pan-cancer按鈕，在新頁(yè)面的Gene symbol框中輸入hsa-miR-501。在Split patients by選項(xiàng)中選擇median，使用腫瘤組織中has-miR-501表達(dá)量的中位數(shù)區(qū)分高表達(dá)和低表達(dá)，在plot options中選擇Export plot data as text并在癌癥類型中勾選Lung squamous cell carcinoma。其余選項(xiàng)保持默認(rèn)即可。點(diǎn)擊Draw Kaplan-meier plot按鈕。在新頁(yè)面中點(diǎn)擊Export plot data as text即可獲得原始數(shù)據(jù)，見(jiàn)圖1。

圖1 數(shù)據(jù)格式Figure 1.Data Format

1.2 Kaplan-Meier法生存分析的介紹

Kaplan-Meier法又稱為乘積極限法，用于未分組資料的生存分析，是一個(gè)充分考慮了在不同的分組中可能有患者中途退出實(shí)驗(yàn)時(shí)候的統(tǒng)計(jì)學(xué)處理方法。它可以充分考慮每個(gè)時(shí)間點(diǎn)上患者的生存率，因此比較適合用于腫瘤等慢性病的統(tǒng)計(jì)學(xué)分析[13-14]。生存曲線，即Kaplan-Meier曲線。它是以生存時(shí)間tk為橫軸，生存率為縱軸，繪制而成的連續(xù)型的階梯形曲線。中位生存時(shí)間，是指生存率為0.5時(shí)對(duì)應(yīng)的生存時(shí)間，本文提供的案例采用內(nèi)插法估計(jì)。

1.3 對(duì)數(shù)秩檢驗(yàn)法的介紹

對(duì)數(shù)秩檢驗(yàn)法用于對(duì)比兩個(gè)或者多個(gè)生存曲線的位置是否有差異。它的基本思想是假定無(wú)效假設(shè)成立，即兩總體的生存曲線位置相同，那么據(jù)此應(yīng)該得出兩總體中不同生存時(shí)間的初期觀察人數(shù)和通過(guò)理論死亡概率計(jì)算獲得的理論死亡數(shù)與實(shí)際死亡數(shù)應(yīng)當(dāng)相差不大，否則無(wú)效假設(shè)不成立，認(rèn)為兩生存曲線的位置不相同。該方法不要求生存時(shí)間服從某特定的分布，且是對(duì)整個(gè)生存曲線的比較，而不是某個(gè)特定時(shí)間生存率的比較，適合于本案例中按照has-mir-501表達(dá)量高低分組的腫瘤患者的生存時(shí)間分析。本文提供的數(shù)據(jù)表明兩個(gè)分組的生存曲線沒(méi)有相交，粗略預(yù)測(cè)可以通過(guò)ph檢驗(yàn)，可使用對(duì)數(shù)秩檢驗(yàn)法。

1.4 Excel與VBA在生存分析中的使用

1.4.1 用戶界面設(shè)計(jì)

在Excel自帶的VBE開(kāi)發(fā)環(huán)境中點(diǎn)擊菜單欄＞插入，插入3個(gè)用戶窗體并使用VBE中的工具箱插入控件，如圖2所示。

1.4.2 模塊設(shè)計(jì)與程序邏輯

向該工程插入三個(gè)模塊?！癒M法統(tǒng)計(jì)和作圖”模塊中主要包含5個(gè)VB過(guò)程，單擊窗體上的“統(tǒng)計(jì)描述并繪制生存曲線”按鈕將會(huì)調(diào)用他們，其過(guò)程名字和參數(shù)及返回值如下所示：

圖2 Excel插件設(shè)計(jì)界面Figure 2.Excel plug-in design int

“錄入數(shù)據(jù)”模塊包含2個(gè)VB過(guò)程，用于實(shí)現(xiàn)“刪除數(shù)據(jù)并初始化表格”和“新增加案例”兩個(gè)按鈕的功能，其過(guò)程名字和參數(shù)及返回值如下所示:

“對(duì)數(shù)秩檢驗(yàn)”模塊包含 5個(gè)VB過(guò)程，用于實(shí)現(xiàn)“對(duì)數(shù)秩檢驗(yàn)”窗體中提供的功能。其中LogRankTest過(guò)程用于調(diào)用其余4個(gè)過(guò)程。代碼如下：

對(duì)數(shù)秩檢驗(yàn)的算法步驟較多，操作步驟即圖3所示。

以上諸多步驟和最后的畫圖功能由代碼調(diào)用下面的其余四個(gè)vb過(guò)程實(shí)現(xiàn)：

以上各個(gè)模塊各個(gè)函數(shù)的代碼見(jiàn)附件1，讀者可自《醫(yī)學(xué)新知》官網(wǎng)（http://www.jnewmed.com/）相對(duì)應(yīng)文章中數(shù)據(jù)與資料欄下載使用，同時(shí)包含具體操作視頻。

1.4.3 插件的使用

圖3 對(duì)數(shù)秩檢驗(yàn)程序設(shè)計(jì)思路Figure 3.Design ideas of log-rank test program

在生存分析窗體上點(diǎn)擊“刪除數(shù)據(jù)并初始化表格”，輸入案例數(shù)量2并點(diǎn)擊確定，把附件2（下載查看方式同附件1）中兩個(gè)案例的數(shù)據(jù)復(fù)制到Excel中，點(diǎn)擊統(tǒng)計(jì)描述并繪制生存曲線，Excel會(huì)自動(dòng)進(jìn)行所有案例的生存分析并繪圖，如果用戶需要添加誤差線可以選中生存曲線，點(diǎn)擊圖表工具＞設(shè)計(jì)＞添加圖表元素＞誤差線＞其他誤差線選項(xiàng)，在設(shè)置誤差線格式中設(shè)置垂直誤差線的誤差量為自定義，點(diǎn)擊指定值，在正錯(cuò)誤值和負(fù)錯(cuò)誤值中選擇誤差線繪制數(shù)據(jù)一表當(dāng)中的數(shù)據(jù)，可以根據(jù)需要選擇方法一或者方法二產(chǎn)生的誤差線長(zhǎng)度?；氐皆紨?shù)據(jù)表格，點(diǎn)擊“生存分析窗體中的對(duì)數(shù)秩檢驗(yàn)按鈕”。在彈出的對(duì)數(shù)秩檢驗(yàn)窗體中選擇兩個(gè)案例然后點(diǎn)擊“開(kāi)始按鈕”即可獲得對(duì)數(shù)秩檢驗(yàn)結(jié)果和生存曲線對(duì)比圖。

2 結(jié)果

2.1 生存率及其標(biāo)準(zhǔn)誤的計(jì)算

本例生存時(shí)間以月為單位，并將t月當(dāng)作一個(gè)時(shí)點(diǎn)看待。對(duì)分析結(jié)果表格中各欄含義解釋如下：第(1)列為序號(hào)k。第(2)列是將生存時(shí)間t由小到大依次排列,如某時(shí)間點(diǎn)既有完全數(shù)據(jù)又有截尾數(shù)據(jù),將截尾數(shù)據(jù)排在后面。第(3)列為t月的死亡人數(shù)dk。截尾患者即便是已死亡,也非死于研究變量,所以相應(yīng)的d=0。第(4)列為初期病例數(shù)nk,即恰好在t時(shí)點(diǎn)以前尚存活的病例數(shù)。第(5)列計(jì)算各時(shí)點(diǎn)死亡概率qk,即在t時(shí)點(diǎn)以前尚存活的患者恰好在t時(shí)點(diǎn)上(第t個(gè)月)死亡的概率。第(6)列計(jì)算各時(shí)點(diǎn)生存概率pk,即在t時(shí)點(diǎn)以前尚存活的患者在t時(shí)點(diǎn)上(第t個(gè)月)繼續(xù)存活的概率。第(7)列計(jì)算各時(shí)點(diǎn)生存率,即在t時(shí)點(diǎn)以前尚存活的患者活過(guò)t時(shí)點(diǎn)的概率。第(8)列為各時(shí)間點(diǎn)生存率的標(biāo)準(zhǔn)誤。見(jiàn)表1。

2.2 估計(jì)各時(shí)點(diǎn)總體生存率的95%置信區(qū)間

可以采用兩種方法估計(jì)某時(shí)點(diǎn)總體生存率的置信區(qū)間。方法一采用正態(tài)近似原理，使用生存率±1.96×標(biāo)準(zhǔn)誤的方法計(jì)算每個(gè)生存率的置信區(qū)間，不適合于曲線尾部或接近尾部總體生存率的置信區(qū)間估計(jì)。因?yàn)榇颂幍恼龖B(tài)性較差，所估計(jì)的置信區(qū)間的上、下限值可能小于0或大于1。此時(shí)可以計(jì)算經(jīng)過(guò)對(duì)數(shù)變換后的生存率值以及相應(yīng)的標(biāo)準(zhǔn)誤，據(jù)此來(lái)估計(jì)其置信區(qū)間，即方法二。見(jiàn)表2。

表1 極限乘積法估計(jì)miR-501不同表達(dá)量水平的生存率及標(biāo)準(zhǔn)誤Table 1.The limit product method to estimate the survival rate and standard error of miR-501 at different expression levels

2.3 中位生存時(shí)間的計(jì)算

Excel把使用內(nèi)插法計(jì)算獲得的中位生存時(shí)間寫入表格尾部。低表達(dá)miR-501的患者中位生存時(shí)間為37.484830。高表達(dá)miR-501的患者中位生存時(shí)間為74.371976。

2.4 繪制生存曲線

生存曲線中水平橫線的長(zhǎng)短代表一個(gè)t時(shí)點(diǎn)到下一個(gè)t時(shí)點(diǎn)的距離,當(dāng)最后一個(gè)時(shí)點(diǎn)的觀察對(duì)象全部死亡時(shí),曲線與橫軸相交。生存曲線圖可直觀地比較各組觀察對(duì)象的生存過(guò)程,也可對(duì)任意時(shí)刻的生存率作出粗略估計(jì)。繪制該圖首先需要生成繪圖專用的數(shù)據(jù)區(qū)域，該數(shù)據(jù)表格給出了每個(gè)時(shí)間開(kāi)始和結(jié)束時(shí)的生存率。然后根據(jù)這些數(shù)據(jù)繪制xy散點(diǎn)圖，將每個(gè)散點(diǎn)連線。對(duì)于截尾數(shù)據(jù)而言，其生存率沒(méi)有變化，代表其時(shí)間頭和時(shí)間尾的散點(diǎn)在圖上重合為一點(diǎn)，根據(jù)這一特點(diǎn)添加刪失記號(hào)即可完成繪圖。見(jiàn)圖4、圖5。

2.5 對(duì)數(shù)秩檢驗(yàn)

最終獲得結(jié)果見(jiàn)表3。

這兩組生存曲線的P值小于0.05，表明兩組數(shù)據(jù)的差異具有統(tǒng)計(jì)學(xué)意義。miR-501在肺鱗癌組織中表達(dá)較高的病人生存時(shí)間較長(zhǎng)。生存曲線圖如圖6。

3 討論

本文采用了乘積極限法和對(duì)數(shù)秩檢驗(yàn)等經(jīng)典生存分析方法對(duì)mir-501的表達(dá)情況進(jìn)行分析，發(fā)現(xiàn)mir-501表達(dá)量較高的分組中腫瘤患者的生存時(shí)間較長(zhǎng)。目前已有較多的學(xué)者研究了microRNA對(duì)腫瘤患者生存時(shí)間的影響，Shi ming 等人通過(guò)薈萃分析和生存分析研究了miR-30d-5p表達(dá)量對(duì)卵巢癌患者生存時(shí)間的影響，發(fā)現(xiàn)mir-30d-5p 表達(dá)較高的患者生存時(shí)間較長(zhǎng)[15]。Hui Wang等人通過(guò)TCGA數(shù)據(jù)庫(kù)獲取生存數(shù)據(jù)對(duì)mir-21調(diào)控的 LZTFL1基因進(jìn)行表達(dá)量與生存時(shí)間相關(guān)性的研究[16]。他們發(fā)現(xiàn)LZTFL1表達(dá)較高的患者生存時(shí)間較長(zhǎng)。Rui-Sheng Zhou等人通過(guò)TCGA數(shù)據(jù)庫(kù)的挖掘，探究了多個(gè)microRNA, lncRNA, mRNA對(duì)舌鱗狀細(xì)胞癌患者生存時(shí)間的影響[17]。可見(jiàn)乘積極限法已經(jīng)成為一種廣泛使用的用于生存時(shí)間分析的統(tǒng)計(jì)學(xué)工具，配合對(duì)數(shù)秩檢驗(yàn)這樣的經(jīng)典檢驗(yàn)方法可以良好的區(qū)別某個(gè)基因表達(dá)情況對(duì)腫瘤患者預(yù)后的影響。

表2 miR-501不同表達(dá)量水平各個(gè)時(shí)間生存率的95%置信區(qū)間計(jì)算Table 2.Calculation of 95% confidence interval for survival rate of miR-501 at different expression levels at each time

圖4 癌組織低表達(dá)miR-501的肺鱗癌患者生存曲線Figure 4.Survival curve of lung squamous cell carcinoma patients with low expression of miR-501 in cancer tissue

圖5 癌組織高表達(dá)miR-501的肺鱗癌患者生存曲線Figure 5.Survival curve of lung squamous cell carcinoma patients with high expression of miR-501 in cancer tissue

表3 低表達(dá)miR-501與高表達(dá)miR-501的對(duì)數(shù)秩檢驗(yàn)表Table 3.Log-rank test table of low expression miR-501 and high expression miR-501

圖6 兩組生存曲線對(duì)比繪圖Figure 6.Comparison of survival curves between the two groups

本文利用了Excel軟件中自帶的Visual Basic（VBA）編程功能，以miR-501為研究目標(biāo)分析肺鱗癌患者組織中的miR-501表達(dá)水平與患者生存時(shí)間之間的關(guān)系，拓展了Excel的使用范圍，把生存分析這類醫(yī)學(xué)特有的統(tǒng)計(jì)方式帶入Excel等大眾辦公軟件中，降低了生存分析的軟件門檻，有利于生存分析方法的普及,極大地方便了一般科研人員對(duì)統(tǒng)計(jì)工具的使用。本程序在極端情況下可能會(huì)計(jì)算出錯(cuò)，例如在office安裝不完整的情況下Excel 可能無(wú)法提供VBA二次開(kāi)發(fā)功能。如果用戶提供的數(shù)據(jù)中出現(xiàn)大量的刪失數(shù)據(jù)，此時(shí)即使K-M法統(tǒng)計(jì)已經(jīng)沒(méi)有了統(tǒng)計(jì)學(xué)意義，但是本程序依然會(huì)按照算法得出不可靠的結(jié)果。經(jīng)過(guò)測(cè)試，此程序需要在Excel2013以上版本才能正常執(zhí)行。在完成Excel的對(duì)數(shù)秩檢驗(yàn)后在兩組生存曲線的繪制中可以手動(dòng)在圖片上加上P值，方便讀者閱讀。目前統(tǒng)計(jì)學(xué)已經(jīng)在許多慢性病的科研項(xiàng)目中發(fā)揮了巨大的威力[18]，許多專業(yè)的統(tǒng)計(jì)學(xué)軟件已經(jīng)大量運(yùn)用到了醫(yī)學(xué)科研中。

本文為統(tǒng)計(jì)軟件的使用者提供了更多的選擇，旨在起到拋磚引玉的作用。生存分析中還有一些其他的檢驗(yàn)方法如Breslow檢驗(yàn)等，本程序沒(méi)有探索，讀者可以自行查閱相關(guān)檢驗(yàn)的原理并編程實(shí)現(xiàn)。

綜上所述，使用Excel 結(jié)合Visual Basic語(yǔ)言編程可以較好地分析miR-501 表達(dá)高低對(duì)肺鱗癌患者生存時(shí)間的影響，其結(jié)果可信可靠。