Skip to content

Commit f924678

Browse files
committed
air format R/scatterplot.R
1 parent a0fd038 commit f924678

File tree

1 file changed

+95
-43
lines changed

1 file changed

+95
-43
lines changed

R/scatterplot.R

Lines changed: 95 additions & 43 deletions
Original file line numberDiff line numberDiff line change
@@ -6,52 +6,104 @@
66
#' @param dat Output from [harmonise_data()].
77
#' @export
88
#' @return List of plots
9-
mr_scatter_plot <- function(mr_results, dat)
10-
{
11-
# dat <- subset(dat, paste(id.outcome, id.exposure) %in% paste(mr_results$id.outcome, mr_results$id.exposure))
12-
mrres <- plyr::dlply(dat, c("id.exposure", "id.outcome"), function(d)
13-
{
14-
d <- plyr::mutate(d)
15-
if(nrow(d) < 2 | sum(d$mr_keep) == 0)
16-
{
17-
return(blank_plot("Insufficient number of SNPs"))
18-
}
19-
d <- subset(d, mr_keep)
20-
index <- d$beta.exposure < 0
21-
d$beta.exposure[index] <- d$beta.exposure[index] * -1
22-
d$beta.outcome[index] <- d$beta.outcome[index] * -1
23-
mrres <- subset(mr_results, id.exposure == d$id.exposure[1] & id.outcome == d$id.outcome[1])
24-
mrres$a <- 0
25-
if("MR Egger" %in% mrres$method)
26-
{
27-
temp <- mr_egger_regression(d$beta.exposure, d$beta.outcome, d$se.exposure, d$se.outcome, default_parameters())
28-
mrres$a[mrres$method == "MR Egger"] <- temp$b_i
29-
}
9+
mr_scatter_plot <- function(mr_results, dat) {
10+
# dat <- subset(dat, paste(id.outcome, id.exposure) %in% paste(mr_results$id.outcome, mr_results$id.exposure))
11+
mrres <- plyr::dlply(dat, c("id.exposure", "id.outcome"), function(d) {
12+
d <- plyr::mutate(d)
13+
if (nrow(d) < 2 | sum(d$mr_keep) == 0) {
14+
return(blank_plot("Insufficient number of SNPs"))
15+
}
16+
d <- subset(d, mr_keep)
17+
index <- d$beta.exposure < 0
18+
d$beta.exposure[index] <- d$beta.exposure[index] * -1
19+
d$beta.outcome[index] <- d$beta.outcome[index] * -1
20+
mrres <- subset(
21+
mr_results,
22+
id.exposure == d$id.exposure[1] & id.outcome == d$id.outcome[1]
23+
)
24+
mrres$a <- 0
25+
if ("MR Egger" %in% mrres$method) {
26+
temp <- mr_egger_regression(
27+
d$beta.exposure,
28+
d$beta.outcome,
29+
d$se.exposure,
30+
d$se.outcome,
31+
default_parameters()
32+
)
33+
mrres$a[mrres$method == "MR Egger"] <- temp$b_i
34+
}
3035

31-
if("MR Egger (bootstrap)" %in% mrres$method)
32-
{
33-
temp <- mr_egger_regression_bootstrap(d$beta.exposure, d$beta.outcome, d$se.exposure, d$se.outcome, default_parameters())
34-
mrres$a[mrres$method == "MR Egger (bootstrap)"] <- temp$b_i
35-
}
36+
if ("MR Egger (bootstrap)" %in% mrres$method) {
37+
temp <- mr_egger_regression_bootstrap(
38+
d$beta.exposure,
39+
d$beta.outcome,
40+
d$se.exposure,
41+
d$se.outcome,
42+
default_parameters()
43+
)
44+
mrres$a[mrres$method == "MR Egger (bootstrap)"] <- temp$b_i
45+
}
3646

37-
ggplot2::ggplot(data=d, ggplot2::aes(x=beta.exposure, y=beta.outcome)) +
38-
ggplot2::geom_errorbar(ggplot2::aes(ymin=beta.outcome-se.outcome, ymax=beta.outcome+se.outcome), colour="grey", width=0) +
39-
ggplot2::geom_errorbarh(ggplot2::aes(xmin=beta.exposure-se.exposure, xmax=beta.exposure+se.exposure), colour="grey", height=0) +
40-
ggplot2::geom_point() +
41-
ggplot2::geom_abline(data=mrres, ggplot2::aes(intercept=a, slope=b, colour=method), show.legend=TRUE) +
42-
ggplot2::scale_colour_manual(values=c("#a6cee3", "#1f78b4", "#b2df8a", "#33a02c", "#fb9a99", "#e31a1c", "#fdbf6f", "#ff7f00", "#cab2d6", "#6a3d9a", "#ffff99", "#b15928")) +
43-
ggplot2::labs(colour="MR Test", x=paste("SNP effect on", d$exposure[1]), y=paste("SNP effect on", d$outcome[1])) +
44-
ggplot2::theme(legend.position="top", legend.direction="vertical") +
45-
ggplot2::guides(colour=ggplot2::guide_legend(ncol=2))
46-
})
47-
mrres
47+
ggplot2::ggplot(
48+
data = d,
49+
ggplot2::aes(x = beta.exposure, y = beta.outcome)
50+
) +
51+
ggplot2::geom_errorbar(
52+
ggplot2::aes(
53+
ymin = beta.outcome - se.outcome,
54+
ymax = beta.outcome + se.outcome
55+
),
56+
colour = "grey",
57+
width = 0
58+
) +
59+
ggplot2::geom_errorbarh(
60+
ggplot2::aes(
61+
xmin = beta.exposure - se.exposure,
62+
xmax = beta.exposure + se.exposure
63+
),
64+
colour = "grey",
65+
height = 0
66+
) +
67+
ggplot2::geom_point() +
68+
ggplot2::geom_abline(
69+
data = mrres,
70+
ggplot2::aes(intercept = a, slope = b, colour = method),
71+
show.legend = TRUE
72+
) +
73+
ggplot2::scale_colour_manual(
74+
values = c(
75+
"#a6cee3",
76+
"#1f78b4",
77+
"#b2df8a",
78+
"#33a02c",
79+
"#fb9a99",
80+
"#e31a1c",
81+
"#fdbf6f",
82+
"#ff7f00",
83+
"#cab2d6",
84+
"#6a3d9a",
85+
"#ffff99",
86+
"#b15928"
87+
)
88+
) +
89+
ggplot2::labs(
90+
colour = "MR Test",
91+
x = paste("SNP effect on", d$exposure[1]),
92+
y = paste("SNP effect on", d$outcome[1])
93+
) +
94+
ggplot2::theme(legend.position = "top", legend.direction = "vertical") +
95+
ggplot2::guides(colour = ggplot2::guide_legend(ncol = 2))
96+
})
97+
mrres
4898
}
4999

50100

51-
blank_plot <- function(message)
52-
{
53-
ggplot2::ggplot(data.frame(a=0,b=0,n=message)) +
54-
ggplot2::geom_text(ggplot2::aes(x=a,y=b,label=n)) +
55-
ggplot2::labs(x=NULL,y=NULL) +
56-
ggplot2::theme(axis.text=ggplot2::element_blank(), axis.ticks=ggplot2::element_blank())
101+
blank_plot <- function(message) {
102+
ggplot2::ggplot(data.frame(a = 0, b = 0, n = message)) +
103+
ggplot2::geom_text(ggplot2::aes(x = a, y = b, label = n)) +
104+
ggplot2::labs(x = NULL, y = NULL) +
105+
ggplot2::theme(
106+
axis.text = ggplot2::element_blank(),
107+
axis.ticks = ggplot2::element_blank()
108+
)
57109
}

0 commit comments

Comments
 (0)