优化meta分析中的森林图 (forest plot)

2020-06-01 11:59:16 浏览数 (1)

一般meta中的森林图是这样的:

见到过高端文章里比较复杂的是这样的,每个物种有很多个效应量:

那么就来尝试一下实现这种图。

读入数据:

代码语言:javascript复制
>library(ggplot2)
>eff_size <- read.table("E:/forestdata.txt", header = T,sep = "t")
>eff_size
                               species          A          B          C      Lower
1              Fusobacterium nucleatum  0.7011937  0.8211937  0.5011937  0.5308094
2                     Parvimonas micra  0.6503721  0.7703721  0.4503721  0.4806099
3                      Parvimonas spp.  0.6292963  0.7492963  0.4292963  0.4534209
4                  Gemella morbillorum  0.6005114  0.7205114  0.4005114  0.4391758
5          Peptostreptococcus stomatis  0.5803963  0.7003963  0.3803963  0.4192075
6                 Solobacterium moorei  0.5323909  0.6523909  0.3323909  0.3720123
7                Clostridium symbiosum  0.4791824  0.5991824  0.2791824  0.3112856
8               Anaerococcus vaginalis  0.4176721  0.5376721  0.2176721  0.2423912
9        Porphyromonas asaccharolytica  0.4153169  0.5353169  0.2153169  0.2558456
10               Prevotella intermedia  0.3888456  0.5088456  0.1888456  0.2138550
11             Streptococcus australis -0.2203435 -0.1003435 -0.4203435 -0.3828304
12 Lachnospiraceae bacterium 1 1 57FAA -0.2314120 -0.1114120 -0.4314120 -0.3896990
13           Ruminococcussp 5 1 39BFAA -0.2384873 -0.1184873 -0.4384873 -0.3968436
14                    Pseudomonas spp. -0.2550211 -0.1350211 -0.4550211 -0.4389134
15              Bifidobacterium longum -0.2570670 -0.1370670 -0.4570670 -0.4156960
16                  Eubacterium hallii -0.2631902 -0.1431902 -0.4631902 -0.4500349
17         Adlercreutzia equolifaciens -0.2660466 -0.1460466 -0.4660466 -0.4244658
18          Streptococcus thermophilus -0.2794899 -0.1594899 -0.4794899 -0.4381088
19        Faecalibacterium prausnitzii -0.2985913 -0.1785913 -0.4985913 -0.4573489
20              Roseburia intestinalis -0.3151002 -0.1951002 -0.5151002 -0.4784761
         Upper  y
1   0.87157807  1
2   0.82013431  2
3   0.80517175  3
4   0.76184690  4
5   0.74158514  5
6   0.69276956  6
7   0.64707925  7
8   0.59295310  8
9   0.57478820  9
10  0.56383624 10
11 -0.05785654 11
12 -0.07312512 12
13 -0.08013106 13
14 -0.07112876 14
15 -0.09843809 15
16 -0.07634546 16
17 -0.10762733 17
18 -0.12087110 18
19 -0.13983379 19
20 -0.15172427 20

#长宽转化
>library(reshape2)
>eff.ggplot = melt(eff_size,
             id.vars = c("species","Lower","Upper","y"),    
             measure.vars = c("A","B","C"),                     
             variable.name='group',
             value.name='values')  
>head(eff.ggplot)
                      species     Lower     Upper y group    values
1     Fusobacterium nucleatum 0.5308094 0.8715781 1     A 0.7011937
2            Parvimonas micra 0.4806099 0.8201343 2     A 0.6503721
3             Parvimonas spp. 0.4534209 0.8051718 3     A 0.6292963
4         Gemella morbillorum 0.4391758 0.7618469 4     A 0.6005114
5 Peptostreptococcus stomatis 0.4192075 0.7415851 5     A 0.5803963
6        Solobacterium moorei 0.3720123 0.6927696 6     A 0.5323909

#设定主题
>plot_theme = theme(panel.background=element_blank(),
                   panel.grid=element_blank(),
                   axis.line.x=element_line(size=.5, colour="black"),
                   axis.line.y=element_line(size=.5, colour="black"),
                   axis.ticks=element_line(color="black"),
                   axis.text=element_text(color="black", size=20),
                   legend.position="right",
                   legend.background=element_blank(),
                   legend.key=element_blank(),
                   legend.text= element_text(size=20),
                   text=element_text(family="sans", size=20))

#散点和置信区间连线
>p = ggplot(data=eff.ggplot,aes(values,y,col=group)) 
    geom_point(shape=group,size=3)
>for (i in 1:nrow(eff.ggplot)){ 
  p <- p   geom_segment(aes_string(
                x = eff.ggplot$Lower[i], y = eff.ggplot$y[i], 
                xend = eff.ggplot$Upper[i], yend = eff.ggplot$y[i]),
                color="tomato3",size=1)
};
>p = p   geom_vline(xintercept=0,color="tomato3")

加物种注释:

代码语言:javascript复制
>left = eff_size[eff_size$Lower>0,]
>right = eff_size[eff_size$Upper<0,]

>p = p   annotate("text",x=0.15,y=right$y,label=right$species,hjust=0) 
        annotate("text",x=-0.1,y=left$y,label=left$species,hjust=1)
>p = p   labs(title="Forest plot", 
              subtitle="Multi- effect size for each species",
              caption="By 水岸风堤",
              x="Effect size",
              y="Species") 
              plot_theme;p

0 人点赞