單一率Meta分析的Matlab軟件實現*

肖 靜 吳桂云 高月霞 沈 毅 李奕辰 陳德喜 蔡 輝△

Meta分析是對以往的研究結果進行系統定量綜合的統計學方法[1]。近年來，隨著循證醫學的迅速興起，Meta分析研究的結論在有關疾病或醫學健康問題的防治、治療等方面越來越受到重視[2]。其中單一率資料的Meta分析是一種只提供人群總人數和事件發生人數，不同于其它有對照人群的Meta分析，大多為患病率、檢出率、病死率和感染率等[3]，屬于橫斷面研究，是揭示暴露與疾病關系不可或缺的，因此基于計算機的實現是很有必要的[4]。本文運用實例數據基于Matlab軟件編寫了Meta分析以及發表偏倚和異質性檢驗實現的相關程序，供讀者參考。

實例分析

數據來自英國萊切斯特大學Eaden等[5]對大腸潰瘍患者發生結直腸癌風險的Meta分析，摘錄年齡組為11～12歲，數據見表1。

表1 患大腸潰瘍人群發生結直腸癌風險Meta分析數據(部分)

1.Matlab程序實現過程：

(1)建立數據集并計算各研究效應量：

clear all;

clc;

s=[8 54

12 2295

9 429

……];%錄入數據，第一列為各個研究事件發生例數，第二列為各個研究總個體數。

r=s(:,1)./s(:,2);%取出s矩陣中的第一列除以第二列，計算各個研究事件的發生率ri。

a=s(:,1);%定義s矩陣中第一列為ai。

n=s(:,2);%定義s矩陣中第二列數據為ni。

Dr=exp(log(r)-1.96*se);%計算各個研究發生率95%的下限，注意Matlab軟件中log表示ln。

Ur=exp(log(r)+1.96*se);%計算各個研究發生率95%的上限。

w=1./se.^2;%各個研究的權重。

Dhr=exp(log(hr)-1.96*sehr);

Uhr=exp(log(hr)+1.96*sehr);

Q=sum(w.*(log(r)-log(hr)).^2);%Q檢驗。

k=size(s,1);

Qp=1-chi2cdf(Q,k-1);%計算Q檢驗的p值。

I2=(Q-(k-1))/Q*100;%計算I2統計量。

If I2<0

I2=0;

楊琳說，以前朋友的老公。以前朋友是指歐陽橘紅。聽楊琳這一說，一種從來不曾有的溫暖，流入他的心中。老天爺不睜開眼睛看一看，這樣的好人，為什么偏偏得這種惡病？為什么好人就命不長呢？

end

If I2>0

I2=I2;

end

I2

Z=abs(log(hr)/sehr);%I2檢驗統計量值的計算。

Zp=2*(1-normcdf(Z,0,1));%I2檢驗統計量對應的P值。

H=(Q/(k-1))^0.5;%H檢驗統計量的計算。

HSE=0.5.*(log(Q)-log(k-1))./((2*Q)^0.5-(2*k-3)^0.5;%H檢驗統計量的標準誤。

Hs=exp(log(H)+1.96*HSE);

Hx=exp(log(H)-1.96*HSE);

If Qp<=0.1;%當Q檢驗p<=0.1時進行隨機效應模型。

If Q<=k-1

tau=0;

end

If Q>k-1

tau=(Q-(k-1))./(sum(w)-(sum(w.^2)./(sum(w))));%隨機效應模型中τ2的計算。

end

tau

wj=1./((se).^2+tau);%w′的計算。

Zj=abs(log(hrj)/sehrj);%z′的計算。

Zpj=2*(1-normcdf(Zj,0,1));%z′對應的P值的計算。

End

(3)森林圖：此處由于森林圖的程序繁雜略去，感興趣的讀者可以聯系本文作者索取，結果見圖1所示。

(4)發表偏倚的檢驗：繪制漏斗圖以及對漏斗圖進行檢驗(Egger′s檢驗)，相關結果見圖2、圖3。

Figure;%開始畫漏斗圖。

y=se;%以精度的倒數為縱坐標。

plot(x,y,′o′);%繪制圖。

hold on

plot([log(hrj)log(hrj)],[0 1.2],′k-′)

plot([log(hrj)log(hrj)-1.96*1.2],[0 1.2],′k′)

plot([log(hrj)log(hrj)+1.96*1.2],[0,1.2],′k:′)

set(gca,′tickdir′,′out′,′ydir′,′reverse′,′fontsize′,15,′linewidth′,1.5,′xlim′,[-8-2],′ylim′,[0,1.2],′xtick′,[-8,-6,-4,-2],′xticklabel′,[{num2str(exp(-8))},{num2str(exp(-6))},{num2str(exp(-4))},{num2str(exp(-2))}]);%設置坐標軸。

xlabel(′log(r)′,′fontsize′,16);%給X軸加標簽。

ylabel(′se′,′fontsize′,16);%給Y軸加標簽。

box off

%漏斗圖開始檢驗。

Y=log(r)./se;%定義Y數值。

X=1./se;%定義X數值。

b=sum((X-mean(X)).*(Y-mean(Y)))./sum((X-mean(X)).^2);%計算斜率b。

a=mean(Y)-b.*mean(X);%計算截距a。

figure

plot(X,Y,′O′,′markersize′,5);

hold on

plot([0,max(X)],[a,a+b.*max(X)],′k-′);

Yi=a+b.*X;%建立回歸方程。

Syx=(sum((Y-Yi).^2)./(size(X,1)-2)).^0.5;

Sa=Syx.*(sum((X).^2)./(size(X,1).*sum((X-mean(X)).^2))).^0.5;

t=abs((a-0)/Sa;%對a進行t檢驗，此處計算t檢驗統計量。

P=2*(1-tcdf(t,size(X,1)-2));%計算t檢驗統計量對應的自由度k-2下的P值。

tinv(0.975,size(X,1)-2);

ax=a-tinv(0.975,size(X,1)-2)*Sa;

as=a+tinv(0.975,size(X,1)-2)*Sa;

set(gca,′tickdir′,′out′,′fontsize′,15,′linewidth′,1.5,′xlim′,[-0.18],′ylim′,[-328],′xtick′,[1:2:8],′ytick′,[-32:4:8]);%定義坐標軸。

xlabel(′1/se′,′fontsize′,16);

ylabel(′log(r)/se′,′fontsize′,16);

box off

plot([0 0],[axas],′b-′,′linewidth′,1);%畫靠近y的線。

plot(0,a,′rs′,′markersize′,6,′markerfacecolor′,′w′);%畫區間線中間的正方形框，為了直觀表現。

plot([0-0.05 0+0.05],[axax],′k-′,′linewidth′,1);

plot([0-0.05 0+0.05],[asas],′k-′,′linewidth′,1);

(5)異質性檢驗：繪制Galbraith圖，結果見圖4。

figure

plot(X,Y,′.′,′markersize′,8);

hold on

plot([0,max(X)],[0,log(hr).*max(X)],′k-′);

hold on

plot([0,max(X)],[2,log(hr).*max(X)+2],′k-′);%在上兩個單位處畫等斜率直線。

plot([0,max(X)],[-2,log(hr).*max(X)-2],′k-′);

set(gca,′tickdir′,′out′,′fontsize′,15,′linewidth′,1.5,′xlim′,[0 6.5],′ylim′,[-35 2]);

xlabel(′1/se′,′fontsize′,16);

ylabel(′log(r)/se′,′fontsize′,16);

box off

for i=1:size(X,1)

text(X(i)+0.02,Y(i),num2str(i));%此處為給每個點加名稱。

end

注：以上程序為通用程序，可以根據異質性檢驗結果自行選擇計算模型；程序中加冒號表示結果中隱藏不顯示，百分號后面的為注釋部分，在程序中不運行。

圖1 森林圖

圖2 漏斗圖

圖3 Egger′s檢驗圖

圖4 Galbraith圖

討 論

隨著循證醫學的發展，Meta分析已被公認為客觀評價和合成針對某一特定問題研究證據的最佳手段[3]。對于目前國內通用的Meta分析軟件RevMan不能進行Meta回歸分析、累積Meta分析和Egger′s檢驗等，Stata軟件雖為目前Meta分析最受推崇的軟件，但無單一率Meta分析的固定模塊。本文基于Matlab軟件編寫了單一率資料Meta分析的實現程序，非專業人員只需替換數據即可方便實現。Matlab人際交互性強，具有強大的數據處理功能[6]，且分析結果的計算精度可以通過程序自由控制，同時Matlab具有強大的作圖功能，所作圖形比同類軟件更為美觀。本文通過實例分析驗證了利用Matlab程序實現單一率Meta分析結果的可行性和有效性，并完成了異質性和發表偏倚的檢驗。

參 考 文 獻

1.Egger M,Smith G D,Altman D G.Systematic reviews in health care.Meta-analysis in context.2nd edition.London:BMJ Publishing Group,2001.

2.俞慧強,鄭輝烈,李悅,等.Meta分析發表偏倚診斷方法研究.中國衛生統計,2011,28(4):402-405.

3.曾憲濤,冷衛東,李勝,等.如何正確理解及使用GRADE系統.中國循證醫學雜志,2011,11(9):985-990.

4.王佩鑫,李宏田,劉建蒙.無對照二分類資料的Meta分析方法及Stata實現.循證醫學,2012,12(1):52-64.

5.Eaden J A,Abrams K R,Mayberry J F.The risk of colorectal cancer in ulcerative colitis:a meta-analysis.Gut,2001,48(4):526-535.

6.胡小剛,陳劍鴻,孫鳳軍,等.基于Matlab的Kruskalk-Wallis和Nemenyi檢驗的界面實現.中國衛生統計,2011,28(4):466-473.