我想大家肯定都用过森林图,应用比较多的场景可能是展示meta分析
,回归分析
结果的时候。🥳
画森林图的包还是挺多的,今天介绍一下forplo包
的用法。😘
rm(list = ls())
library(tidyverse)
library(forplo)
library(meta)
library(autoReg)
library(survival)
这里我准备了一个随机数据,假设是我们收集的各篇文章中的某个指标。
dat <- read.csv("./forplo_meta.csv",row.names = 1)
DT::datatable(dat)
forplo(dat[,1:3],
groups=dat$groups,
font = "Arial",
grouplabs=c('Low risk of bias',
'Some concerns',
'High risk of bias'))
logORs <- round(log(dat$OR),2)
SE <- round((log(dat$UCI)-logORs)/1.96,2)
meta1 <- meta::metagen(logORs[1:3],SE[1:3])
meta2 <- meta::metagen(logORs[4:9],SE[4:9])
meta3 <- meta::metagen(logORs[10:12],SE[10:12])
metatot <- meta::metagen(logORs,SE)
dat2 <- dat
dat2$logORs <- logORs
dat2$SE <- SE
dat2$weights <- round(metatot$w.random/sum(metatot$w.random)*100,2)
dat2 <- rbind(subset(dat2,groups==1),
c(round(exp(meta1$TE.random),2),
round(exp(meta1$lower.random),2),
round(exp(meta1$upper.random),2),
1,round(meta1$TE.random,2),round(meta1$seTE.random,2),sum(dat2$weights[1:3])),
subset(dat2,groups==2),
c(round(exp(meta2$TE.random),2),
round(exp(meta2$lower.random),2),
round(exp(meta2$upper.random),2),
2,round(meta2$TE.random,2),round(meta2$seTE.random,2),sum(dat2$weights[4:9])),
subset(dat2,groups==3),
c(round(exp(meta3$TE.random),2),
round(exp(meta3$lower.random),2),
round(exp(meta3$upper.random),2),
3,round(meta3$TE.random,2),round(meta3$seTE.random,2),sum(dat2$weights[10:12])),
c(round(exp(metatot$TE.random),2),
round(exp(metatot$lower.random),2),
round(exp(metatot$upper.random),2),
4,round(metatot$TE.random,2),round(metatot$seTE.random,2),100))
dat2 <- data.frame(dat2)
rownames(dat2)[c(4,11,15,16)] <- c('Effect 1','Effect 2','Effect 3','Effect 4')
DT::datatable(dat2)
forplo(dat2[,1:3],
font = "Arial",
groups=dat2$groups,
grouplabs=c('Low risk of bias',
'Some concerns',
'High risk of bias',
'Overall'),
left.align=TRUE,
add.columns=dat2[,5:7],
add.colnames=c('log(OR)','SE','Weights'),
scaledot.by = dat2$weights,
col= '#FEBE8C',
ci.edge= T,
char = 15, # 10为空心圆, 15为实心方格, 20实心圆
diamond=c(4,11,15,16),
diamond.col= '#EB6440',
favorlabs=c('Favours other models','Favours midwife-led'), # x.lab
## 阴影参数
shade.every=1, # 每行加背景底色
)
本期使用survival
里的cancer
数据,类型为lung
。😘
data("cancer")
DT::datatable(lung)
这里我们先做个factor
转换,方便后续操作。🧐
lung$sex <- factor(lung$sex, levels = unique(lung$sex))
lung$ph.ecog <- factor(lung$ph.ecog, levels = c(0,1,2,3))
把所有变量都纳入进去进行回归建模吧。👇
fit1 <- coxph(Surv(time, status) ~ sex + age + ph.ecog + ph.karno + pat.karno + meal.cal + wt.loss,
data = lung)
summary(fit1)
这里我们用一下之前介绍的autoreg
, 看一下回归结果。🤩
Note! 由于这个dataset
存在缺失值,所以使用了mice
包进行填补,想要了解mice
具体教程的小伙伴可以在公众号内搜索ggmice
或者mice
。😏
autoReg(fit1, uni = T, threshold = 0.1,
final = T,imputed = T) %>%
myft()
不知道大家有没有发现, 这个结果和autoreg
出的结果会有一点细微差异, 仁者见仁,智者见智吧。🤒
forplo(fit1,
font='Arial',
sort = T,
#flipbelow1=T,
left.align=T,
ci.edge= T,
scaledot.by=abs(coef(fit1)),
col = '#BAD1C2',
char = 20,
shade.every = 1,
shade.col = '#FDF0E0',
shade.alpha = 0.8
)
lung$status <- factor(lung$status, levels = rev(unique(lung$status)))
fit2 <- glm(status ~ sex + age + ph.karno + pat.karno + meal.cal + wt.loss,
data = lung,
family="binomial")
summary(fit2)
autoReg(fit2, uni = T, threshold = 0.1,
final = T,imputed = T) %>%
myft()
forplo(fit2,
font='Arial',
sort = T,
#flipbelow1=T,
left.align=T,
ci.edge= T,
scaledot.by=abs(coef(fit2)),
col = '#DBA39A',
char = 15,
shade.every = 1,
shade.col = '#FDF0E0',
shade.alpha = 0.8,
right.bar=T,
rightbar.ticks=T,
leftbar.ticks=T
)
最后祝大家早日不卷!~