侵權投訴
訂閱
糾錯
加入自媒體

森林圖 (forestplot包)

2020-12-28 10:23
科研菌
關注

在Meta分析中森林圖比較常見,但其實掌握了用R語言中的forestplot包繪制森林圖的各個用法,森林圖可以用于表示其他數(shù)據(jù)類型各組間的指標的中值和四分位距的范圍。它是在平面直角坐標系中,以一條垂直的無效線(橫坐標刻度為1或0)為中心,平行于橫軸的多條線段描述了每個組的指標的中值和可信區(qū)間,最后一行(Summary)則用棱形(或其它圖形)描述了多個組別合并的效應量及可信區(qū)間。

首先安裝forestplot 包:

install.packages("forestplot")

主函數(shù) forestplot :

用法:

forestplot(labeltext, mean, lower, upper, align,           is.summary = FALSE, graph.pos = "right", hrzl_lines, clip = c(-Inf,           Inf), xlab = "", zero = ifelse(xlog, 1, 0), graphwidth = "auto",           colgap, lineheight = "auto", line.margin, col = fpColors(),           txt_gp = fpTxtGp(), xlog = FALSE, xticks, xticks.digits = 2,           grid = FALSE, lwd.xaxis, lwd.zero, lwd.ci, lty.ci = 1, ci.vertices,           ci.vertices.height = 0.1, boxsize, mar = unit(rep(5, times = 4),           "mm"), title, legend, legend_args = fpLegend(),           new_page = getOption("forestplot_new_page", TRUE),           fn.ci_norm = fpDrawNormalCI, fn.ci_sum = fpDrawSummaryCI, fn.legend,           ...)

參數(shù):這里只列出了大部分參數(shù),還有一些比較不常用的可以自行探索

labeltext主要是以矩陣或者list形式將數(shù)據(jù)導入函數(shù),最好以矩陣,因為數(shù)據(jù)一般都是矩陣的。mean誤差條的均值lower誤差條 95%置信區(qū)間下限upper誤差條 95%置信區(qū)間上限align每列文字的對齊方式,偶爾會用到。如:align=c("l","c","c")l:左對齊r:右對齊c:居中對齊is.summary主要的功能是讓表格的每一行字體出現(xiàn)差異,從而區(qū)分表頭。其值主要用TRUE/FALSE進行差異化分配。graph.pos定位森林圖所在的位置。通過數(shù)字來確定為第幾列。hrzl_lines以list形式設置表中線條的類型、影響范圍。Eg:“3”=gpar(lwd=1,columns=1:4,col=’red’)意思就是第3行的線條,寬度為1,線段延伸至第四列。Col指的顏色。clipx軸的最大最小范圍xlabx軸的標題zero森林圖中基準線的位置(無效線的橫坐標)graphwidth森林圖在表中的寬度如:graphwidth = unit(.4,"npc")colgap

列與列之間的間隙寬度,默認是 6 mm,需要用 unit 的形式

lineheight行的高度,可以是數(shù)字,也可以是 unit 的形式line.margin行與行之間的間隙的寬度col森林圖橫線以及點的顏色。box:box(點估計值)的顏色line:穿過方塊的橫線的顏色zero:中間那條基準線的顏色summary:summary中菱形的顏色hrz_lines:表中第一條橫線的顏色eg:col=fpcolors(box=’royblue’,line=’darkblue’, summary=’royblue’, hrz_lines=’red’)txt_gp設置表格中文本的格式:用gpar進行賦值,其中cex為文本字體大小,ticks為坐標軸大小,xlab為坐標軸文字字體大小。label:表格主體文字的格式ticks:森林圖下方的坐標軸的刻度文字格式xlab:定義的x軸標題格式title:標題文字的格式eg:txt_gp=fpTxtGp(label=gpar(cex=1.25),                                  ticks=gpar(cex=1.1),                                  xlab=gpar(cex = 1.2),                                  title=gpar(cex = 1.2))xticks橫坐標刻度根據(jù)需要可隨意設置,如:xticks = c(0.5, 1,1.5, 2)lwd.xaxisX軸線寬lwd.zero無效線的寬度lwd.ci置信區(qū)間線條的寬度(粗細)lty.ci置信區(qū)間的線條類型ci.vertices森林圖可信區(qū)間兩端添加小豎線(TRUE)ci.vertices.height設置森林圖可信區(qū)間兩端的小豎線高度,默認是10%行高boxsizebox(點估計值)的大小mar圖形頁邊距,如:mar=unit(rep(1.25, times = 4), "cm")title添加標題legend當同時顯示多個置信區(qū)間時,需要添加圖例

new_page

是否新頁fn.ci_normbox(點估計值)的形狀,默認是方塊。如:fn.ci_norm="fpDrawDiamondCI":box 類型選擇鉆石

示例代碼①:先從構建的最簡單的數(shù)據(jù)開始

# 構建示例數(shù)據(jù)

library(forestplot)# Cochrane data from the 'rmeta'-packagecochrane_from_rmeta <- structure(list(                    mean =c(NA, NA, 0.578, 0.165, 0.246, 0.700, 0.348, 0.139, 1.017, NA, 0.531),                    lower =c(NA, NA, 0.372, 0.018, 0.072, 0.333, 0.083, 0.016, 0.365, NA, 0.386),                    upper =c(NA, NA, 0.898, 1.517, 0.833, 1.474, 1.455, 1.209, 2.831, NA, 0.731)),                    .Names =c("mean", "lower", "upper"),                    row.names =c(NA, -11L),                    class ="data.frame")

tabletext<-cbind(c("", "Study", "Auckland", "Block","Doran", "Gamsu", "Morrison", "Papageorgiou","Tauesch", NA, "Summary"),c("Deaths", "(steroid)", "36", "1","4", "14", "3", "1","8", NA, NA),c("Deaths", "(placebo)", "60", "5","11", "20", "7", "7","10", NA, NA),c("", "OR", "0.58", "0.16","0.25", "0.70", "0.35", "0.14","1.02", NA, "0.53"))

以下使用 forestplot 函數(shù)畫森林圖,注意查看每個代碼發(fā)生變化的參數(shù)以及對應圖片中明顯變化的地方。    

forestplot(tabletext, graph.pos =4,            hrzl_lines =list("3"=gpar(lty=2),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="#000044"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="#000044")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(TRUE,TRUE,rep(FALSE,8),TRUE),           clip=c(0.1,2.5),           xlog=TRUE,           col=fpColors(box="royalblue",line="darkblue",summary="royalblue",hrz_lines ="#444444"))

forestplot(tabletext, graph.pos =4,            hrzl_lines =list("3"=gpar(lty=2,col ="red"),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="#000044"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="#000044")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(TRUE,TRUE,rep(FALSE,8),TRUE),           clip=c(0.1,2.5),           xlog=TRUE,           col=fpColors(box="royalblue",line="darkblue",summary="royalblue",hrz_lines ="#444444"))

forestplot(tabletext, graph.pos =4,            hrzl_lines =list("3"=gpar(lty=2),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="red"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="#000044")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(TRUE,TRUE,rep(FALSE,8),TRUE),           clip=c(0.1,2.5),           xlog=TRUE,           col=fpColors(box="royalblue",line="darkblue",summary="royalblue",hrz_lines ="#444444"))

forestplot(tabletext, graph.pos =4,            hrzl_lines =list("3"=gpar(lty=2),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="#000044"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="red")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(TRUE,TRUE,rep(FALSE,8),TRUE),           clip=c(0.1,2.5),           xlog=TRUE,           col=fpColors(box="royalblue",line="darkblue",summary="royalblue",hrz_lines ="#444444"))

forestplot(tabletext, graph.pos =4,            hrzl_lines =list("3"=gpar(lty=2),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="#000044"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="#000044")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(TRUE,TRUE,rep(FALSE,8),TRUE),           clip=c(0.1,2.5),           xlog=TRUE,           col=fpColors(box="red",line="darkblue",summary="royalblue",hrz_lines ="#444444"))

forestplot(tabletext, graph.pos =4,            hrzl_lines =list("3"=gpar(lty=2),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="#000044"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="#000044")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(TRUE,TRUE,rep(FALSE,8),TRUE),           clip=c(0.1,2.5),           xlog=TRUE,           col=fpColors(box="royalblue",line="red",summary="royalblue",hrz_lines ="#444444"))

forestplot(tabletext, graph.pos =4,            hrzl_lines =list("3"=gpar(lty=2),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="#000044"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="#000044")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(TRUE,TRUE,rep(FALSE,8),TRUE),           clip=c(0.1,2.5),           xlog=TRUE,           col=fpColors(box="royalblue",line="darkblue",summary="red",hrz_lines ="#444444"))

forestplot(tabletext, graph.pos =4,            hrzl_lines =list("3"=gpar(lty=2),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="#000044"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="#000044")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(TRUE,TRUE,rep(FALSE,8),TRUE),           clip=c(0.1,2.5),           xlog=TRUE,           col=fpColors(box="royalblue",line="darkblue",summary="royalblue",hrz_lines ="red"))

# 添加標題forestplot(tabletext, graph.pos =4,           title="Hazard Ratio",           hrzl_lines =list("3"=gpar(lty=2),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="#000044"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="#000044")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(TRUE,TRUE,rep(FALSE,8),TRUE),           clip=c(0.1,2.5),           xlog=TRUE,           col=fpColors(box="royalblue",line="darkblue",summary="royalblue",hrz_lines ="#444444"))

# 定義x軸forestplot(tabletext, graph.pos =4,           title="Hazard Ratio",           hrzl_lines =list("3"=gpar(lty=2),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="#000044"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="#000044")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(TRUE,TRUE,rep(FALSE,8),TRUE),           xlab=" <---PCI Better--- ---Medical Therapy Better--->",           clip=c(0.1,2.5),           xlog=TRUE,           col=fpColors(box="royalblue",line="darkblue",summary="royalblue",hrz_lines ="#444444"))

# 取消對頭2行和最后1行字體的特殊設置forestplot(tabletext, graph.pos =4,           title="Hazard Ratio",           hrzl_lines =list("3"=gpar(lty=2),                            "11" =gpar(lwd=1, columns=c(1:3,5), col ="#000044"),                            "12" =gpar(lwd=1, lty=2, columns=c(1:3,5), col ="#000044")),           cochrane_from_rmeta,new_page =TRUE,           is.summary=c(rep(FALSE,11)),           xlab=" <---PCI Better--- ---Medical Therapy Better--->",           clip=c(0.1,2.5),           xlog=TRUE,           zero=1,           col=fpColors(box="royalblue",line="darkblue",summary="royalblue",hrz_lines ="#444444"))

示例代碼②:需要定義亞組的數(shù)據(jù)

### 準備數(shù)據(jù)

library(forestplot)#數(shù)據(jù)來源:https://www.r-bloggers.com/forest-plot-with-h(huán)orizontal-bands/data <- read.csv("forestplotdata.csv", stringsAsFactors=FALSE)
head(data)#   Variable Count Percent Point.Estimate  Low High PCI.Group Medical.Therapy.Group P.Value# 1  Overall  2166     100            1.3 0.90 1.50      17.2                  15.6      NA# 2             NA      NA             NA   NA   NA        NA                    NA      NA# 3      Age    NA      NA             NA   NA   NA        NA                    NA    0.05# 4    <= 65  1534      71            1.5 1.05 1.90      17.0                  13.2      NA# 5     > 65   632      29            0.8 0.60 1.25      17.8                  21.3      NA# 6             NA      NA             NA   NA   NA        NA                    NA      NA

### 繪制森林圖

## 簡單森林圖

# 構建tabletext,更改列名稱,將 count 和 percent 合并np <- ifelse(!is.na(data$Count), paste(data$Count," (",data$Percent,")",sep=""), NA)View(np)

# 寫出將要在圖中展現(xiàn)出來的文本tabletext <- cbind(c("Subgroup","",data$Variable),                   c("No. of Patients (%)","",np),                   c("4-Yr Cum. Event Rate PCI","",data$PCI.Group),                   c("4-Yr Cum. Event Rate Medical Therapy","",data$Medical.Therapy.Group),                   c("P Value","",data$P.Value)View(tabletext)                  

##繪制森林圖forestplot(labeltext=tabletext, graph.pos=3,           mean=c(NA,NA,data$Point.Estimate),           lower=c(NA,NA,data$Low), upper=c(NA,NA,data$High),           boxsize=0.5)

接下來要對森林圖進行優(yōu)化:

## 定義亞組subgps <- c(4,5,8,9,12,13,16,17,20,21,24,25,28,29,32,33)data$Variable[subgps] <- paste("  ",data$Variable[subgps])View(data)

png(filename = "Forestplot.png",width=960, height=640)forestplot(labeltext=tabletext,           graph.pos=3, #為Pvalue箱線圖所在的位置           mean=c(NA,NA,data$Point.Estimate),           lower=c(NA,NA,data$Low),            upper=c(NA,NA,data$High),           title="Hazard Ratio Plot", #定義標題           xlab="    <---PCI Better---   ---Medical Therapy Better--->", #定義x軸           ##根據(jù)亞組的位置,設置線型,寬度造成“區(qū)塊感”           hrzl_lines=list("3" = gpar(lwd=1, col="black"),                           "7" = gpar(lwd=60, lineend="butt", columns=c(2:6), col="#99999922"),                           "15" = gpar(lwd=60, lineend="butt", columns=c(2:6), col="#99999922"),                           "23" = gpar(lwd=60, lineend="butt", columns=c(2:6), col="#99999922"),                           "31" = gpar(lwd=60, lineend="butt", columns=c(2:6), col="#99999922")),           #fpTxtGp函數(shù)中的cex參數(shù)設置各部分字體大小           txt_gp=fpTxtGp(label=gpar(cex=1.25),                          ticks=gpar(cex=1.1),                          xlab=gpar(cex = 1.2),                          title=gpar(cex = 1.2)),           col=fpColors(box="#1c61b6", lines="#1c61b6", zero = "gray50"), ##fpColors函數(shù)設置顏色           zero=1, #箱線圖中基準線的位置           cex=0.9, lineheight = "auto",           colgap=unit(8,"mm"),           lwd.ci=2, boxsize=0.5, #箱子大小,線的寬度           ci.vertices=TRUE, ci.vertices.height = 0.4) #森林圖可信區(qū)間兩端添加小豎線,設置高度dev.off()

森林圖怎么看:

(1)森林圖中橫短線與中線相交表示無統(tǒng)計學意義;

(2)95% CI上下限均>1,即在森林圖中,其95% CI橫線不與無效豎線相交,且該橫線落在無效線右側(cè)時,說明該指標大于豎線代表的結(jié)局;

(3)95% CI上下限均<1,即在森林圖中,其95% CI橫線不與無效豎線相交,且該橫線落在無效線左側(cè)時,說明該指標小于于豎線代表的結(jié)局。

(4)最后以菱形所在位置代表總體的評價結(jié)果。具體數(shù)據(jù)具體分析哈!

聲明: 本文由入駐維科號的作者撰寫,觀點僅代表作者本人,不代表OFweek立場。如有侵權或其他問題,請聯(lián)系舉報。

發(fā)表評論

0條評論,0人參與

請輸入評論內(nèi)容...

請輸入評論/評論長度6~500個字

您提交的評論過于頻繁,請輸入驗證碼繼續(xù)

暫無評論

暫無評論

醫(yī)療科技 獵頭職位 更多
文章糾錯
x
*文字標題:
*糾錯內(nèi)容:
聯(lián)系郵箱:
*驗 證 碼:

粵公網(wǎng)安備 44030502002758號