Wie platziere ich mit ggplot2 eine Risikotabelle unter einem Kaplan-Meier-Diagramm?
Ich möchte einen Kaplan-Meier-Plot mit ggplot2 mit einer Risikonummer-Tabelle erstellen, die die Risikonummer für jede Gruppe zu jedem Zeitpunkt angibt (d. H. X-Achsen-Häkchen). Die Anzahl der gefährdeten Personen sollte mit dem entsprechenden Häkchen übereinstimmen. Links von der Tabelle mit den Risikonummern sollten Zeilennamen stehen, die die Gruppe angeben, zu der die Risikonummern gehören.
Ich habe das folgende Beispiel geschrieben. Ich lerne aus diesem @ die gefährdeten Zahlen zu bestimmFrag. Ich weiß jedoch nicht, wie ich eine schöne, gut ausgerichtete Zahl am Risikotisch unter dem Kaplan-Meier-Plot erstellen soll. Ein Freund hat mir geholfen, im folgenden Beispiel die Risikotabelle mit der Nummer zu erstellen. Die resultierende Zahl meines Beispiels ist jedoch unzureichend.
library(survival)
library(reshape2)
data(colon)
library(Hmisc)
d <- colon[, Cs(time, status, rx)]
rm(colon)
names(d) <- c("days", "event", "group")
d$group <- ifelse(d$group == "Obs", 1, 2)
fit <- survfit(Surv(days,event)~group, data=d)
diff <- survdiff(Surv(days,event)~group, data=d)
risksets <- with(na.omit(d[, Cs(days, event, group)]), table(group, cut(days, seq(0, max(days), by=365) ) ))
number.at.risk <- sapply(1:nrow(risksets), function(i) Reduce("-", risksets[i,], init=rowSums(risksets)[i], accumulate=TRUE))
number.at.risk <- data.frame(number.at.risk)
names(number.at.risk) <- c("Group.A", "Group.B")
number.at.risk
###
p.value <- round(1 - pchisq(diff$chisq, 1), digits=4)
p.value <- ifelse(p.value < 0.001, "<0.001", paste("= ", p.value))
d.mortality <- data.frame(time=fit$time, surv=fit$surv, strata=summary(fit, censored=T)$strata)
zeros <- data.frame(time=0, surv=1, strata=unique(d.mortality$strata))
d.mortality <- rbind(d.mortality, zeros)
levels(d.mortality$strata) <- c("Group A", "Group B")
d.mortality$surv <- (1-d.mortality$surv)*100 # event free to events and in %
###
g <- ggplot(d.mortality, aes(time, surv, group=strata)) +
geom_step(aes(colour=strata), size=1) +
theme_bw() + # white background
theme(
plot.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
legend.position="none",
axis.line = element_line(color = 'black'),
axis.text.x = element_text(size=15),
axis.text.y = element_text(size=15),
axis.title.x = element_text(size=17, hjust=.5, vjust=.25, face="bold"),
axis.title.y = element_text(size=17, hjust=.5, vjust=1.5, face="bold"),
plot.title = element_text(size=20, hjust=-.1, vjust=1, face="bold")
) +
scale_y_continuous("Cumulative event rate [%]", limits=c(0, 60)) +
scale_x_continuous("Time [years]", limits=c(0, 1825), breaks=seq(0, 1825, 365), labels=c(0, 1, 2, 3, 4, 5)) +
annotate("text", x = 1000, y = 45, label = "Group A") +
annotate("text", x = 1000, y = 30, label = "Group B") +
annotate("text", x = 1000, y = 55, label = paste("P ", p.value, "by log-rank test", collapse=""))
number.at.risk = number.at.risk[1:6,]
df_nums = melt(number.at.risk)
df_nums$year = 1:6
tbl = ggplot(df_nums, aes(x = year, y = factor(variable), colour = variable,label=value)) +
geom_text(size = 3.5) + theme(panel.grid.major = element_blank(), legend.position = "none") + theme_bw() +
theme(
plot.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
legend.position="none",
axis.line = element_blank(),
axis.text.x = element_blank(),
axis.ticks=element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
plot.title = element_blank()
) + scale_y_discrete(breaks=c("Group.B","Group.A"), labels=c("Number at Risk\nGroup B", "Group A"))
Layout <- grid.layout(nrow = 2, ncol = 1, heights = unit(c(2, 0.55), c("null", "null")))
grid.show.layout(Layout)
vplayout <- function(...) {
grid.newpage()
pushViewport(viewport(layout = Layout))
}
subplot <- function(x, y) viewport(layout.pos.row = x, layout.pos.col = y)
mmplot <- function(a, b) {
vplayout()
print(a, vp = subplot(1, 1))
print(b, vp = subplot(2, 1))
}
dev.new()
mmplot(g, tbl)
Wie vorgeschlagen habe ich gtable mit der resultierenden Figur benutzt. Ich war mit dem Layout der Variante a (Beispielcode von baptiste) nicht zufrieden und habe etwas anderes ausprobiert. Version B hat jedoch einen weiteren Nachteil: Die Beschriftungen befinden sich innerhalb der x-Dimensionen des Plot-Layers des Hauptplots.
a) Wie kann ich eine vernünftig gestaltete Figur mit gut ausgerichteten Risikozahlen erstellen?
b) Wie kann ich außerdem einen Titel "Numbers at Risk" zwischen dem Hauptplot und der Tabelle platzieren? Der Titel "Numbers at Risk" sollte am linken Ende der Bezeichnungen "Group A" und "Group B" von @ ausgerichtet seitbl.
c) Die Schriftgröße der Risikozahlen in tbl und der entsprechenden Bezeichnungen "Gruppe A" und "Gruppe B" sollte mit den Häkchenbezeichnungen in der Hauptdarstellung übereinstimmen. Wie kann ich das machen
library(survival)
library(reshape2)
data(colon)
library(Hmisc)
d <- colon[, Cs(time, status, rx)]
rm(colon)
names(d) <- c("days", "event", "group")
d$group <- ifelse(d$group == "Obs", 1, 2)
fit <- survfit(Surv(days,event)~group, data=d)
diff <- survdiff(Surv(days,event)~group, data=d)
risksets <- with(na.omit(d[, Cs(days, event, group)]), table(group, cut(days, seq(0, max(days), by=365) ) ))
number.at.risk <- sapply(1:nrow(risksets), function(i) Reduce("-", risksets[i,], init=rowSums(risksets)[i], accumulate=TRUE))
number.at.risk <- data.frame(number.at.risk)
names(number.at.risk) <- c("Group.A", "Group.B")
number.at.risk
###
p.value <- round(1 - pchisq(diff$chisq, 1), digits=4)
p.value <- ifelse(p.value < 0.001, "<0.001", paste("= ", p.value))
d.mortality <- data.frame(time=fit$time, surv=fit$surv, strata=summary(fit, censored=T)$strata)
zeros <- data.frame(time=0, surv=1, strata=unique(d.mortality$strata))
d.mortality <- rbind(d.mortality, zeros)
levels(d.mortality$strata) <- c("Group A", "Group B")
d.mortality$surv <- (1-d.mortality$surv)*100 # event free to events and in %
###
g <- ggplot(d.mortality, aes(time, surv, group=strata)) +
geom_step(aes(colour=strata), size=1) +
# theme_bw() + # white background
theme(
plot.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
legend.position="none",
axis.line = element_line(color = 'black'),
axis.text.x = element_text(size=15),
axis.text.y = element_text(size=15),
axis.title.x = element_text(size=17, hjust=.5, vjust=.25, face="bold"),
axis.title.y = element_text(size=17, hjust=.5, vjust=4, face="bold"),
plot.title = element_text(size=20, hjust=-.1, vjust=1, face="bold")
) +
scale_y_continuous("Cumulative event rate [%]", limits=c(0, 60)) +
scale_x_continuous("Time [years]", limits=c(0, 1825), breaks=seq(0, 1825, 365), labels=c(0, 1, 2, 3, 4, 5)) +
annotate("text", x = 1000, y = 45, label = "Group A") +
annotate("text", x = 1000, y = 30, label = "Group B") +
annotate("text", x = 1000, y = 55, label = paste("P ", p.value, "by log-rank test", collapse=""))
number.at.risk = number.at.risk[1:6,]
df_nums = melt(number.at.risk)
df_nums$year = 1:6
str(df_nums)
tbl <- ggplot(df_nums, aes(x = year, y = factor(variable), colour = variable, label=value)) +
geom_text() +
# theme_bw() +
theme(
panel.grid.major = element_blank(),
legend.position = "none",
plot.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
legend.position="none",
axis.line = element_blank(),
axis.text.x = element_blank(),
axis.ticks=element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
plot.title = element_blank()
) +
scale_y_discrete(breaks=c("Group.B","Group.A"), labels=c("Group B", "Group A"))
library(gtable)
# Version A
both = rbind(ggplotGrob(g), ggplotGrob(tbl), size="last")
grid.newpage()
grid.draw(both)
# Version B
a <- gtable(unit(15, c("cm")), unit(c(10,3), "cm"))
a <- gtable_add_grob(a, ggplotGrob(g), 1, 1)
a <- gtable_add_grob(a, ggplotGrob(tbl), 2, 1)
grid.newpage()
grid.draw(a)
etzt ist es fast perfekt. Zwei kleine Dinge:
a) Wie kann ich dem Plot den Titel "Nummer in Gefahr" (bekannt mit GIMP) hinzufügen, wie in der folgenden Abbildung gezeigt?
b) Warum steht Gruppe B in der Tabelle über Gruppe A? Die Bezeichnung in df_nums für Gruppe A ist 1 und für Gruppe B 2. Wie kann ich Gruppe A in der Risikonummer-Tabelle über Gruppe B setzen?
> str(df_nums$variable)
Factor w/ 2 levels "Group.A","Group.B": 1 1 1 1 1 1 2 2 2 2 ...
Hier der aktualisierte Code:
library(survival)
library(reshape2)
data(colon)
library(Hmisc)
d <- colon[, Cs(time, status, rx)]
rm(colon)
names(d) <- c("days", "event", "group")
d$group <- ifelse(d$group == "Obs", 1, 2)
fit <- survfit(Surv(days,event)~group, data=d)
diff <- survdiff(Surv(days,event)~group, data=d)
risksets <- with(na.omit(d[, Cs(days, event, group)]), table(group, cut(days, seq(0, max(days), by=365) ) ))
number.at.risk <- sapply(1:nrow(risksets), function(i) Reduce("-", risksets[i,], init=rowSums(risksets)[i], accumulate=TRUE))
number.at.risk <- data.frame(number.at.risk)
names(number.at.risk) <- c("Group.A", "Group.B")
number.at.risk
###
p.value <- round(1 - pchisq(diff$chisq, 1), digits=4)
p.value <- ifelse(p.value < 0.001, "<0.001", paste("= ", p.value))
d.mortality <- data.frame(time=fit$time, surv=fit$surv, strata=summary(fit, censored=T)$strata)
zeros <- data.frame(time=0, surv=1, strata=unique(d.mortality$strata))
d.mortality <- rbind(d.mortality, zeros)
levels(d.mortality$strata) <- c("Group A", "Group B")
d.mortality$surv <- (1-d.mortality$surv)*100 # event free to events and in %
###
g <- ggplot(d.mortality, aes(time, surv, group=strata)) +
geom_step(aes(colour=strata), size=1) +
# theme_bw() + # white background
theme(
plot.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
legend.position="none",
axis.line = element_line(color = 'black'),
axis.text.x = element_text(size=15),
axis.text.y = element_text(size=15),
axis.title.x = element_text(size=17, hjust=.5, vjust=.25, face="bold"),
axis.title.y = element_text(size=17, hjust=.5, vjust=4, face="bold"),
plot.title = element_text(size=20, hjust=-.1, vjust=1, face="bold")
) +
scale_y_continuous("Cumulative event rate [%]", limits=c(0, 60)) +
scale_x_continuous("Time [years]", limits=c(0, 1825), breaks=seq(0, 1825, 365), labels=c(0, 1, 2, 3, 4, 5)) +
annotate("text", x = 1000, y = 45, label = "Group A") +
annotate("text", x = 1000, y = 30, label = "Group B") +
annotate("text", x = 1000, y = 55, label = paste("P ", p.value, "by log-rank test", collapse=""))
number.at.risk = number.at.risk[1:6,]
df_nums = melt(number.at.risk)
str(df_nums$variable)
df_nums
df_nums$year = 1:6
str(df_nums)
tbl <- ggplot(df_nums, aes(x = year, y = factor(variable), colour = variable, label=value)) +
geom_text() +
# theme_bw() +
theme(
panel.grid.major = element_blank(),
legend.position = "none",
plot.background = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
legend.position="none",
axis.line = element_blank(),
axis.text.x = element_blank(),
axis.text.y = element_text(size=15, face="bold", color = 'black'),
axis.ticks=element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
plot.title = element_blank()
) +
scale_y_discrete(breaks=c("Group.A", "Group.B"), labels=c("Group A", "Group B"))
library(gtable)
# Version C
both = rbind(ggplotGrob(g), ggplotGrob(tbl), size="last")
panels <- both$layout$t[grep("panel", both$layout$name)]
both$heights[panels] <- list(unit(1,"null"), unit(2, "lines"))
grid.newpage()
grid.draw(both)
