单一率Meta分析的Matlab软件实现*
2014-04-03吴桂云高月霞李奕辰陈德喜
肖 静 吴桂云 高月霞 沈 毅 李奕辰 陈德喜 蔡 辉△
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.