Sunday, October 4, 2015

Biểu đồ rừng (forest plot)

Biểu đồ rừng là một dạng biểu đồ hay được dùng trong phân tích tổng hợp để thể hiện [thường là] effect size và khoảng tin cậy 95%. Effect size có thể là odds ratio, relative risk, hay hiệu số giữa hai nhóm cho biến liên tục. Bài này sẽ hướng dẫn cách sử dụng ggplot2 để vẽ biểu đồ rừng qua một ví dụ thực tế dưới đây.


Chúng ta có số liệu về odds ratio và khoảng tin cậy 95% của 10 SNP (một dạng biến thể gen). Mỗi SNP, có 3 giá trị: odds ratio (OR), và hai giá trị phản ảnh khoảng tin cậy 95% (LOR và HOR).


SNP = c("rs1", "rs2", "rs3", "rs4", "rs5", "rs6", "rs7", "rs8", "rs9", "rs10")

OR = c(1.17, 0.82, 1.12, 1.02, 1.01, 1.01, 0.95, 0.99, 1.10, 1.11)

LOR = c(1.07, 0.73, 1.02, 0.92, 0.91, 0.91, 0.85, 0.89, 1.01, 0.97)

HOR = c(1.27, 0.93, 1.22, 1.22, 1.10, 1.11, 1.05, 1.09, 1.21, 1.17)

d = data.frame(SNP, OR, LOR, HOR)


Dùng ggplot2 để vẽ:

p = ggplot(data=d, aes(x=OR, y=SNP))
p = p + geom_point(aes(size=OR, fill=SNP), colour="black", shape=16)
p = p + geom_errorbarh(aes(xmin=LOR,xmax=HOR), height=0.2)
p = p + geom_vline(xintercept=1, linetype="dashed")
p = p + scale_x_continuous(name="Odds ratio")
p = p + theme(legend.position=" ")
p + theme_bw() + ylab(" ")




Sắp xếp theo thứ tự và thêm màu:

ddtransform(d, SNP=reorder(SNP, OR))

p = ggplot(data=dd, aes(x=OR, y=SNP))
p = p + geom_point(aes(size=OR, fill=SNP), colour="black", shape=21)
p = p + geom_errorbarh(aes(xmin=LOR,xmax=HOR), height=0.2)
p = p + geom_vline(xintercept=1, linetype="dashed")
p = p + scale_x_continuous(name="Odds ratio")
p = p + theme(legend.position=" ")
p + theme_bw() + ylab(" ")



# Position of snp
# p = p + geom_text(aes(x=0.5,label=rs),size=4) 

Thêm một biểu đồ khác:

Group = factor(c("Infliximab","Etanercept","Adalimumab","Golimumab","Certolizumab","Abatacept","Rituximab","Tocilizumab","Total"),

# actual order
levels=c("Total","Certolizumab","Golimumab","Adalimumab","Etanercept","Infliximab","Tocilizumab","Rituximab","Abatacept"))

ES = c(1.88,2.67,2.35,2.14,5.08,1.68,2.11,2.01,2.16)
LOW = c(1.01,1.44,1.52,1.59,3.46,1.47,1.64,1.57,1.83)
HIGH = c(3.51,4.94,3.65,2.89,7.48,1.90,2.72,2.57,2.55)

dat = data.frame(group, ES, LOW, HIGH)

p = ggplot(dat,aes(ES, Group))
p = p + geom_point(size=5, shape=18)
p = p + geom_errorbarh(aes(xmax=HIGH, xmin=LOW), height=0.15)
p = p + geom_vline(xintercept=1, linetype="longdash")
p = + scale_x_continuous(breaks = seq(0,14,1), labels = seq(0,14,1)) + labs(x="Favours Experimental", y="")
p = p + theme(text = element_text(size=20)) + ggtitle("Risk Ratio for ACR20 response at 6 months")

p





No comments:

Post a Comment