r-statistics-fanの日記

統計好き人間の覚書のようなもの

メタアナリシスじゃない信頼区間をforest plotで書いてみたい

メタアナリシスじゃない信頼区間をforest plotで書いてみたい

同僚に、信頼区間を出すよう頼まれた。 比率なので、二項分布の正確な信頼区間binom.confintを使うことにする。 折角なので、計算したCIをforest plotで書いてみたいと思った。 しかし、低水準関数で一から書くのは何となくハードルが高い。 rmetaの書式を流用して書いてみる。 いずれ、低水準関数で一から書くのもやってみたい。

なお、メタアナリシスではない(別々の薬)ので、totalはそのまま合計して計算した値である。値の架空のもの。

## forest plotを作成する 比率の信頼区間を描画してみる
library(rmeta)
## Loading required package: grid
library(binom)

total.no <- c(1596, 1346, 852, 576, 178, 46)
event.no <- c(195, 160, 100, 50, 30, 7)
drug.name <- c("drug A", "drug B", "drug C", "drug D", "drug E", "drug F")

tabletext <- cbind(c("", "Drug", drug.name, NA, "ALL"), c("side effect", "(No.)", 
    event.no, NA, sum(event.no)), c("Total", "(No.)", total.no, NA, sum(total.no)), 
    c("", "Rate", format(event.no/total.no, digits = 2), NA, format(sum(event.no)/sum(total.no), 
        digits = 2)))

m <- c(NA, NA, binom.confint(event.no, total.no, conf.int = 0.95, methods = "exact")$mean, 
    NA, binom.confint(sum(event.no), sum(total.no), conf.int = 0.95, methods = "exact")$mean)
l <- c(NA, NA, binom.confint(event.no, total.no, conf.int = 0.95, methods = "exact")$lower, 
    NA, binom.confint(sum(event.no), sum(total.no), conf.int = 0.95, methods = "exact")$lower)
u <- c(NA, NA, binom.confint(event.no, total.no, conf.int = 0.95, methods = "exact")$upper, 
    NA, binom.confint(sum(event.no), sum(total.no), conf.int = 0.95, methods = "exact")$upper)

forestplot(tabletext, m, l, u, zero = sum(event.no)/sum(total.no), is.summary = c(rep(FALSE, 
    9), TRUE), clip = c(-Inf, Inf), xlog = FALSE, xlim = c(0.1, 0.2), col = meta.colors(box = "royalblue", 
    line = "darkblue", summary = "royalblue"))

plot of chunk unnamed-chunk-1

うむ、カッコイイのができた。