相關系數是管理學、社會學、心理學及護理學領域常用的統計學指標，基于相關系數數據的Meta分析也日漸增多。R軟件meta程序包與metafor程序包能夠實現多種類型的Meta分析，包括相關系數數據的Meta分析。本文通過實例展示了應用R軟件meta程序包和metafor程序包實現相關系數數據的Meta分析過程，并簡要比較了兩款程序包的結果與繪圖等功能。
引用本文： 魏雪梅, 張永剛, 陶華, 湯紅明, 陳虎, 曾憲濤. 應用R軟件meta程序包與metafor程序包實現相關系數數據的Meta分析. 中國循證醫學雜志, 2015, 15(7): 855-860. doi: 10.7507/1672-2531.20150144 復制
在社會、管理、心理及護理學領域的研究，多報告為相關系數（correlation coefficient，r值）；臨床醫學領域的相關性研究，亦有報告r值的情況[1-5]。對于這類研究所開展的Meta分析即為相關系數的Meta分析。R軟件meta程序包與metafor程序包是兩款功能強大、專用于Meta分析的程序包，能夠實現多種類型數據的Meta分析，包括r值的Meta分析。先前已有文章介紹這兩款程序包實現其他類型數據的Meta分析[6-10]，本文通過實例展示這兩款程序包實現相關系數數據Meta分析的方法。
1 程序包安裝及加載
meta程序包目前最新的版本為v4.2-0，metafor程序包最新版本為v1.9-7，這兩款程序包更加詳細的信息建議參閱官方網站資料。使用前需先行安裝及加載該程序包，本文使用的為R 3.2.0版本。
兩款程序包安裝及加載的命令如下：
install.packages(“meta”)
library(meta)
install.packages("metafor")
library(metafor)
2 數據加載
本文以《基于Pearson相關系數的老年人社會支持與心理健康相關性研究的Meta分析》一文[11]中“客觀支持與軀體化因子”結局指標數據為例進行演示。該文是將相關系數轉換為fisher's Z值與標準誤（standard error，SE）后進行合并，再轉換為r值。本文直接從納入研究中提取了r值與樣本量，以及例文的fisher's Z值和SE值，見表 1所示。
表1
示例“客觀支持與軀體化因子”結局指標數據*
上述數據排列完成后，儲存在桌面的“Rwork”文件夾中的“cor.xlsx”文件中。文件夾及文件的名稱、類型（如可以選擇“cor.txt”格式）及儲存位置均可由操作者自行設置。此外，還可采用R數據框的格式建立數據集。因Excel是常用的Meta分析數據整理軟件，故本文采用直接從Excel中加載數據的方法。
因本文選擇了“.xlsx”格式，故加載前還需安裝及加載xlsx程序包，命令如下：
install.packages("xlsx")
library(xlsx)
完成后就可以進行數據的加載了，具體命令如下：
dat＜-read.xlsx("C:\\Users\\ Administrator \\Desktop\\cor.xlsx",1)
dat
加載后的數據同表 1。
3 meta程序包實現
3.1 分析計算及森林圖繪制
meta程序包依靠“metacor( )”函數可以使用隨機效應模型和固定效應模型進行相關系數數據的Meta分析，基于的合并方法為倒方差法。
在完成數據加載之后，即可開始實現數據的分析。“metacor( )”函數中的相關參數為：cor、n分別為相關系數和總樣本量；data為指定數據集；sm為合并效應量的類型，可選的參數為“zcor”、“cor”分別表示Fisher’ Z值相關系數轉換和直接相關系數合并，一般來說，除非樣本量很大，否則都應該進行Fisher’ Z值相關系數轉換。backtransf的作用為是否顯示Fisher’ Z值，當backtransf=FALSE時將會顯示Fisher’ Z值，否則顯示相關系數數值。
具體的命令展示如下：
datacor＜-metacor(studlab=paste(Study,Year,sep=""),R,N,data=dat,sm="zcor",backtransf=FALSE)
datacor
具體結果見表 2，可以發現，納入研究間無明顯異質性（Q統計量=1.02，相應P=0.796 7；I2=0%），故隨機效應模型和固定效應模型的結果相同，為-0.201 8[95%CI（-0.247 9，-0.155 8）]。
在計算結果顯示完畢后，就可以進行相關圖形的繪制了。以隨機效應模型繪制為例，繪制森林圖的命令如下：
forest(datacor,comb.fixed=FALSE)
繪制的圖形見圖 1。
表2
meta程序包執行命令后的結果
圖1
meta程序包繪制的森林圖
3.2 發表偏倚分析
以隨機效應模型為例，繪制漏斗圖的命令如下（圖 2）：
funnel(datacor,comb.fixed=FALSE)
除開傳統的漏斗圖，該程序包還可以使用Egger線性回歸法進行檢驗并繪制Egger圖，基本的命令如下：
metabias(datacor,method.bias="linreg",plotit=T,k.min=5)
命令運行后的結果顯示因納入研究數目過小，故無法執行Egger檢測。
圖2
meta程序包繪制的漏斗圖
3.3 敏感性分析
該程序包可以實現依次逐項剔除納入研究的敏感性分析。仍以隨機效應模型為例，實現敏感性分析計算及其圖形繪制的命令如下：
datacor＜-metacor(studlab=Study,R,N,data=dat,sm="zcor",backtransf=FALSE,comb.fixed=FALSE)
metainf(datacor)
forest(metainf(datacor))
命令執行后，即可完成敏感性分析，圖形如圖 3所示（計算數據略）。
圖3
meta程序包繪制的敏感性分析森林圖
4 metafor程序包實現
4.1 分析計算及森林圖繪制
metafor程序包依靠“rma( )”函數進行相關系數數據的Meta分析，分析時直接提供原始相關系數及樣本量即可；當然，也可以將相關系數轉化為Fisher’s Z值后進行合并。
為了展示兩款程序包的異同之處，本處使用原始相關系數及樣本量基于直接模型擬合的方法進行合并，具體函數命令如下：
res＜-rma(ri=R,ni=N,measure="COR",method="REML",data=dat)
res
運行上述命令后，計算結果如表 3所示。
表3
metafor程序包計算后的結果
由上述結果可知，采用REML法進行的異質性檢驗結果為P=0.783 5、I2=0%，表明各納入研究間同質性較好，因此應采用固定效應模型的結果。改用固定效應模型后，結果與隨機效應模型相同。合并的總效應量為-0.199 6[95%CI（-0.243 8，
-0.155 5）]。
森林圖繪制的具體命令如下：
forest(res,xlim=c(-6,4) ,order=order(dat$Year),slab=paste(dat$Study,dat$Year,sep="-"),ilab=cbind(dat$R,dat$N),ilab.xpos=c(-3,-2) ,cex=.8)
ref＜-par(cex=.8,font=2)
text(c(-3,-2) ,5.5,c("R","N"))
text(-6,5.5,"Study and Year",pos=4)
text(4,5.5,"R [95%CI]",pos=2)
par(ref)
運行上述命令后得到的森林圖如圖 4所示。
圖4
metafor程序包繪制的森林圖
4.2 發表偏倚分析
漏斗圖法的命令如下：
funnel(res,main="Funnel Plot")
命令運行后，即可得到圖 5所示圖形。
圖5
metafor程序包繪制的漏斗圖
Egger線性回歸檢驗使用“regtest( )”函數實現，具體命令為：
reg＜-regtest(res)
reg
命令運行后結果如下：
Regression Test for Funnel Plot Asymmetry
model: mixed-effects meta-regression model
predictor: standard error
test for funnel plot asymmetry: Z=-0.6562,P= 0.511 7
可以看出，P=0.511 7＞0.05，故表明本Meta分析無發表偏倚存在。
4.3 敏感性分析
如同meta程序包一樣，metafor程序包可以實現依次逐項剔除納入研究的敏感性分析，具體使用的是“influence( )”函數。實現敏感性分析計算及其圖形繪制的命令如下：
influence(res)
plot(influence(res))
命令執行后，即可完成敏感性分析，圖形如圖 6所示（計算數據略）。
圖6
metafor程序包繪制的敏感性分析圖
5 結語
相關系數數據的Meta分析是社會學、管理學、心理學及護理學等領域常見的Meta分析類型。本文使用的實例是通過將r值先行通過公式轉換為fisher's Z值和SE值，再使用RevMan 5軟件進行合并，然后再行使用公式轉換為合并的r值及其可信區間[11]。這種方法需要作者具有較深厚的統計學功底，且在計算時容易出錯。R軟件meta程序包與metafor程序包均可以直接合并r值與樣本量，也可通過命名先行轉換為fisher's Z值和SE值后通過命名合并，再通過命令將合并后的值轉換回r值。這樣對于統計學功底不好的使用者來講簡便易行，且不會出錯。
從分析結果來看，兩款程序包的最終結果與例文使用的RevMan 5軟件結果[2, 12]一致，保留兩位小數的結果均為-0.20[95%CI（-0.25，-0.16）]。兩款程序包在異質性定性檢測與定量檢測方面的結果均一致，也與RevMan 5相同。兩款程序包進行敏感性分析時的方法相同，均是通過依次剔除納入研究以判斷納入研究對整體結果的影響，但展示的圖形卻不相同，metafor程序包提供的信息更為豐富。
此外，兩款程序包均可以實現相關系數數據的亞組分析。meta程序包實現亞組分析僅需要在整體分析的命令中增加變量的定義（通過“byvar”設定），森林圖亦可通過命令一次性展示。相比較之下，metafor程序包實現亞組分析可通過“subset”設定來依次展示亞組分析的結果，亦可使用命令一次性展示；森林圖的展示也可以依次或一次性展示。本文因為納入研究較少，未行展示。
發表偏倚是執行Meta分析時必要的組成部分，當前有多種定性或定量的檢測發表偏倚的方法[1-4]。Cochrane系統評價手冊建議對每個結局指標采用傳統漏斗圖評價發表偏倚，但當納入研究數目過少時則不能采用定量檢測的方法（如Egger’s檢驗）[13]。一般而言，建議納入研究的數目在9項以上時采用定量檢測。本文為了展示meta程序包與metafor程序包的完整功能及其代碼，方便讀者學習，故給出了定性與定量檢測的代碼與結果。
總之，兩款程序包均能很好的滿足相關系數數據的Meta分析需求。從繪圖功能角度而言，metafor程序包的功能更為強大；從數據計算而言，兩款程序包相似；此外，兩款程序包均可繪制輪廓增強漏斗圖[8, 10]。使用者可以根據自身習慣來靈活選用。
在社會、管理、心理及護理學領域的研究，多報告為相關系數（correlation coefficient，r值）；臨床醫學領域的相關性研究，亦有報告r值的情況[1-5]。對于這類研究所開展的Meta分析即為相關系數的Meta分析。R軟件meta程序包與metafor程序包是兩款功能強大、專用于Meta分析的程序包，能夠實現多種類型數據的Meta分析，包括r值的Meta分析。先前已有文章介紹這兩款程序包實現其他類型數據的Meta分析[6-10]，本文通過實例展示這兩款程序包實現相關系數數據Meta分析的方法。
1 程序包安裝及加載
meta程序包目前最新的版本為v4.2-0，metafor程序包最新版本為v1.9-7，這兩款程序包更加詳細的信息建議參閱官方網站資料。使用前需先行安裝及加載該程序包，本文使用的為R 3.2.0版本。
兩款程序包安裝及加載的命令如下：
install.packages(“meta”)
library(meta)
install.packages("metafor")
library(metafor)
2 數據加載
本文以《基于Pearson相關系數的老年人社會支持與心理健康相關性研究的Meta分析》一文[11]中“客觀支持與軀體化因子”結局指標數據為例進行演示。該文是將相關系數轉換為fisher's Z值與標準誤（standard error，SE）后進行合并，再轉換為r值。本文直接從納入研究中提取了r值與樣本量，以及例文的fisher's Z值和SE值，見表 1所示。
表1
示例“客觀支持與軀體化因子”結局指標數據*
上述數據排列完成后，儲存在桌面的“Rwork”文件夾中的“cor.xlsx”文件中。文件夾及文件的名稱、類型（如可以選擇“cor.txt”格式）及儲存位置均可由操作者自行設置。此外，還可采用R數據框的格式建立數據集。因Excel是常用的Meta分析數據整理軟件，故本文采用直接從Excel中加載數據的方法。
因本文選擇了“.xlsx”格式，故加載前還需安裝及加載xlsx程序包，命令如下：
install.packages("xlsx")
library(xlsx)
完成后就可以進行數據的加載了，具體命令如下：
dat＜-read.xlsx("C:\\Users\\ Administrator \\Desktop\\cor.xlsx",1)
dat
加載后的數據同表 1。
3 meta程序包實現
3.1 分析計算及森林圖繪制
meta程序包依靠“metacor( )”函數可以使用隨機效應模型和固定效應模型進行相關系數數據的Meta分析，基于的合并方法為倒方差法。
在完成數據加載之后，即可開始實現數據的分析。“metacor( )”函數中的相關參數為：cor、n分別為相關系數和總樣本量；data為指定數據集；sm為合并效應量的類型，可選的參數為“zcor”、“cor”分別表示Fisher’ Z值相關系數轉換和直接相關系數合并，一般來說，除非樣本量很大，否則都應該進行Fisher’ Z值相關系數轉換。backtransf的作用為是否顯示Fisher’ Z值，當backtransf=FALSE時將會顯示Fisher’ Z值，否則顯示相關系數數值。
具體的命令展示如下：
datacor＜-metacor(studlab=paste(Study,Year,sep=""),R,N,data=dat,sm="zcor",backtransf=FALSE)
datacor
具體結果見表 2，可以發現，納入研究間無明顯異質性（Q統計量=1.02，相應P=0.796 7；I2=0%），故隨機效應模型和固定效應模型的結果相同，為-0.201 8[95%CI（-0.247 9，-0.155 8）]。
在計算結果顯示完畢后，就可以進行相關圖形的繪制了。以隨機效應模型繪制為例，繪制森林圖的命令如下：
forest(datacor,comb.fixed=FALSE)
繪制的圖形見圖 1。
表2
meta程序包執行命令后的結果
圖1
meta程序包繪制的森林圖
3.2 發表偏倚分析
以隨機效應模型為例，繪制漏斗圖的命令如下（圖 2）：
funnel(datacor,comb.fixed=FALSE)
除開傳統的漏斗圖，該程序包還可以使用Egger線性回歸法進行檢驗并繪制Egger圖，基本的命令如下：
metabias(datacor,method.bias="linreg",plotit=T,k.min=5)
命令運行后的結果顯示因納入研究數目過小，故無法執行Egger檢測。
圖2
meta程序包繪制的漏斗圖
3.3 敏感性分析
該程序包可以實現依次逐項剔除納入研究的敏感性分析。仍以隨機效應模型為例，實現敏感性分析計算及其圖形繪制的命令如下：
datacor＜-metacor(studlab=Study,R,N,data=dat,sm="zcor",backtransf=FALSE,comb.fixed=FALSE)
metainf(datacor)
forest(metainf(datacor))
命令執行后，即可完成敏感性分析，圖形如圖 3所示（計算數據略）。
圖3
meta程序包繪制的敏感性分析森林圖
4 metafor程序包實現
4.1 分析計算及森林圖繪制
metafor程序包依靠“rma( )”函數進行相關系數數據的Meta分析，分析時直接提供原始相關系數及樣本量即可；當然，也可以將相關系數轉化為Fisher’s Z值后進行合并。
為了展示兩款程序包的異同之處，本處使用原始相關系數及樣本量基于直接模型擬合的方法進行合并，具體函數命令如下：
res＜-rma(ri=R,ni=N,measure="COR",method="REML",data=dat)
res
運行上述命令后，計算結果如表 3所示。
表3
metafor程序包計算后的結果
由上述結果可知，采用REML法進行的異質性檢驗結果為P=0.783 5、I2=0%，表明各納入研究間同質性較好，因此應采用固定效應模型的結果。改用固定效應模型后，結果與隨機效應模型相同。合并的總效應量為-0.199 6[95%CI（-0.243 8，
-0.155 5）]。
森林圖繪制的具體命令如下：
forest(res,xlim=c(-6,4) ,order=order(dat$Year),slab=paste(dat$Study,dat$Year,sep="-"),ilab=cbind(dat$R,dat$N),ilab.xpos=c(-3,-2) ,cex=.8)
ref＜-par(cex=.8,font=2)
text(c(-3,-2) ,5.5,c("R","N"))
text(-6,5.5,"Study and Year",pos=4)
text(4,5.5,"R [95%CI]",pos=2)
par(ref)
運行上述命令后得到的森林圖如圖 4所示。
圖4
metafor程序包繪制的森林圖
4.2 發表偏倚分析
漏斗圖法的命令如下：
funnel(res,main="Funnel Plot")
命令運行后，即可得到圖 5所示圖形。
圖5
metafor程序包繪制的漏斗圖
Egger線性回歸檢驗使用“regtest( )”函數實現，具體命令為：
reg＜-regtest(res)
reg
命令運行后結果如下：
Regression Test for Funnel Plot Asymmetry
model: mixed-effects meta-regression model
predictor: standard error
test for funnel plot asymmetry: Z=-0.6562,P= 0.511 7
可以看出，P=0.511 7＞0.05，故表明本Meta分析無發表偏倚存在。
4.3 敏感性分析
如同meta程序包一樣，metafor程序包可以實現依次逐項剔除納入研究的敏感性分析，具體使用的是“influence( )”函數。實現敏感性分析計算及其圖形繪制的命令如下：
influence(res)
plot(influence(res))
命令執行后，即可完成敏感性分析，圖形如圖 6所示（計算數據略）。
圖6
metafor程序包繪制的敏感性分析圖
5 結語
相關系數數據的Meta分析是社會學、管理學、心理學及護理學等領域常見的Meta分析類型。本文使用的實例是通過將r值先行通過公式轉換為fisher's Z值和SE值，再使用RevMan 5軟件進行合并，然后再行使用公式轉換為合并的r值及其可信區間[11]。這種方法需要作者具有較深厚的統計學功底，且在計算時容易出錯。R軟件meta程序包與metafor程序包均可以直接合并r值與樣本量，也可通過命名先行轉換為fisher's Z值和SE值后通過命名合并，再通過命令將合并后的值轉換回r值。這樣對于統計學功底不好的使用者來講簡便易行，且不會出錯。
從分析結果來看，兩款程序包的最終結果與例文使用的RevMan 5軟件結果[2, 12]一致，保留兩位小數的結果均為-0.20[95%CI（-0.25，-0.16）]。兩款程序包在異質性定性檢測與定量檢測方面的結果均一致，也與RevMan 5相同。兩款程序包進行敏感性分析時的方法相同，均是通過依次剔除納入研究以判斷納入研究對整體結果的影響，但展示的圖形卻不相同，metafor程序包提供的信息更為豐富。
此外，兩款程序包均可以實現相關系數數據的亞組分析。meta程序包實現亞組分析僅需要在整體分析的命令中增加變量的定義（通過“byvar”設定），森林圖亦可通過命令一次性展示。相比較之下，metafor程序包實現亞組分析可通過“subset”設定來依次展示亞組分析的結果，亦可使用命令一次性展示；森林圖的展示也可以依次或一次性展示。本文因為納入研究較少，未行展示。
發表偏倚是執行Meta分析時必要的組成部分，當前有多種定性或定量的檢測發表偏倚的方法[1-4]。Cochrane系統評價手冊建議對每個結局指標采用傳統漏斗圖評價發表偏倚，但當納入研究數目過少時則不能采用定量檢測的方法（如Egger’s檢驗）[13]。一般而言，建議納入研究的數目在9項以上時采用定量檢測。本文為了展示meta程序包與metafor程序包的完整功能及其代碼，方便讀者學習，故給出了定性與定量檢測的代碼與結果。
總之，兩款程序包均能很好的滿足相關系數數據的Meta分析需求。從繪圖功能角度而言，metafor程序包的功能更為強大；從數據計算而言，兩款程序包相似；此外，兩款程序包均可繪制輪廓增強漏斗圖[8, 10]。使用者可以根據自身習慣來靈活選用。
