Skip to content

Commit

Permalink
book chapter synced version
Browse files Browse the repository at this point in the history
  • Loading branch information
Andrew Dowsey committed Apr 21, 2021
1 parent da5afb3 commit 03d0a57
Show file tree
Hide file tree
Showing 24 changed files with 199 additions and 140 deletions.
2 changes: 2 additions & 0 deletions DESCRIPTION
View file
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,8 @@ Imports:
fitdistrplus (>= 1.0-14),
fst (>= 0.9.0),
ggplot2 (>= 3.3.2),
ggrepel (>= 0.9.1),
gridExtra (>= 2.3),
igraph (>= 1.2.5),
MCMCglmm (>= 2.29),
plotly (>= 4.9.2.1),
Expand Down
6 changes: 3 additions & 3 deletions NAMESPACE
View file
Original file line number Diff line number Diff line change
Expand Up @@ -35,9 +35,6 @@ export(open_delta)
export(open_sigma)
export(open_theta)
export(pinaka)
export(plot_fdr)
export(plot_pr)
export(plot_volcano)
export(plst)
export(qinaka)
export(qlst)
Expand Down Expand Up @@ -102,14 +99,17 @@ exportMethods(plot_component_deviations)
exportMethods(plot_component_means)
exportMethods(plot_component_stdevs)
exportMethods(plot_dists)
exportMethods(plot_fdr)
exportMethods(plot_group_means)
exportMethods(plot_group_quants)
exportMethods(plot_group_quants_de)
exportMethods(plot_group_quants_fdr)
exportMethods(plot_group_standards)
exportMethods(plot_measurement_means)
exportMethods(plot_measurement_stdevs)
exportMethods(plot_pr)
exportMethods(plot_robust_pca)
exportMethods(plot_volcano)
exportMethods(priors)
exportMethods(read)
exportMethods(render_markdown)
Expand Down
89 changes: 50 additions & 39 deletions R/delta_plots.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -22,12 +22,17 @@ setMethod("plots", "seaMass_delta", function(object, batch, job.id) {
# markdown folder
report.index1 <- list()
root1 <- file.path(filepath(object), "markdown", paste0("group.", as.integer(item)))
dir.create(root1, showWarnings = F)
dir.create(root1, recursive = T)

if ("group.quants.de" %in% ctrl@plot) {
group <- control(root(object))@group[1]

fig <- plot_group_quants_fdr(object, item, value.limits = lims$group.quants, summary = T)
fig <- plot_group_quants_fdr(
object, item, value.limits = lims$group.quants, summary = T,
variable.summary.cols = c("Batch", "Effect", "Contrast", "Baseline", "Group", "Cont.uS", "Base.uS", "Cont.qS", "Base.qS",
"Cont.qC", "Base.qC", "Cont.qM", "Base.qM", "lfdr", "lfsr", "qvalue", "svalue", "NegativeProb", "PositiveProb"),
variable.label.cols = c("Group", "Batch", "qvalue")
)
text1 <- paste0(group, " differential expression", ifelse(name(object) == "default", "", paste0(" (", name(object), ")")))
text2 <- paste0(group, " differential expression", ifelse(name(object) == "default", "", paste0(" (", name(object), ") ")), " for ", item)
report.index1$group.quant.de <- data.table(
Expand All @@ -41,19 +46,14 @@ setMethod("plots", "seaMass_delta", function(object, batch, job.id) {
)
}

if (length(report.index1) > 0) {
# zip
render_markdown(object, root1)
if (!("markdown" %in% ctrl@keep)) unlink(root1, recursive = T)
}
# zip
if (length(report.index1) > 0) render_markdown(object, root1)

return(report.index1)
}, nthread = ctrl@nthread))

# save index
if (length(report.index) > 0) {
fst::write.fst(rbindlist(report.index), file.path(filepath(object), "report", paste0("groups.", batch, ".report.fst")))
}
if (length(report.index) > 0) fst::write.fst(rbindlist(report.index), file.path(filepath(object), "report", paste0("groups.", batch, ".report.fst")))

return(invisible(NULL))
})
Expand All @@ -66,14 +66,17 @@ setMethod("plots", "seaMass_delta", function(object, batch, job.id) {
#' @import data.table
#' @export
#' @include generics.R
plot_volcano <- function(
data.fdr,
#' @include seaMass_delta.R
setMethod("plot_volcano", "seaMass_delta", function(
object,
contours = NULL,
error.bars = TRUE,
labels = 25,
stdev.col = "PosteriorSD",
x.col = "PosteriorMean",
y.col = "qvalue"
y.col = "qvalue",
data.fdr = group_quants_fdr(object),
output = "plotly"
) {
DT.fdr <- as.data.table(data.fdr)
DT.fdr[, s := get(stdev.col)]
Expand Down Expand Up @@ -174,20 +177,20 @@ plot_volcano <- function(
g <- g + ggplot2::geom_vline(xintercept = 0)
g <- g + ggplot2::geom_hline(yintercept = ylim.plot[1])
if (x.col == "m") {
g <- g + ggplot2::xlab("log2(Fold Change)")
g <- g + ggplot2::xlab("log2 fold change")
} else if (x.col == "PosteriorMean") {
g <- g + ggplot2::xlab("log2(Shrunk Fold Change)")
g <- g + ggplot2::xlab("log2 moderated fold change")
} else {
g <- g + ggplot2::xlab(paste0("log2(", x.col, ")"))
g <- g + ggplot2::xlab(paste0("log2 ", x.col))
}
if (y.col == "s") {
g <- g + ggplot2::ylab(paste0("-log2(Fold Change) Posterior Standard Deviation"))
g <- g + ggplot2::ylab(paste0("-log2 fold change Posterior Standard Deviation"))
} else if (y.col == "PosteriorSD") {
g <- g + ggplot2::ylab(paste0("-log2(Shrunk Fold Change) Posterior Standard Deviation"))
g <- g + ggplot2::ylab(paste0("-log2 moderated fold change Posterior Standard Deviation"))
} else {
g <- g + ggplot2::ylab(paste0(paste0("-log10(", y.col, ")")))
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.frame(yintercept = -log10(0.01)), linetype = "dashed")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.frame(yintercept = -log10(0.05)), linetype = "dashed")
g <- g + ggplot2::ylab(paste0(paste0("-log10 ", y.col)))
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.table(yintercept = -log10(0.01)), linetype = "dashed")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.table(yintercept = -log10(0.05)), linetype = "dashed")
}
if ("truth" %in% colnames(data.fdr)) {
g <- g + ggplot2::geom_vline(ggplot2::aes(color = Truth, xintercept = truth), DT.meta[N >= 5])
Expand All @@ -196,12 +199,12 @@ plot_volcano <- function(
} else {
g <- g + ggplot2::theme(legend.position = "none")
}
g <- g + ggrepel::geom_label_repel(ggplot2::aes(label = label), size = 2.5, na.rm = T)
g <- g + ggrepel::geom_label_repel(ggplot2::aes(label = label), size = 2.5, na.rm = T, max.overlaps = Inf)
g <- g + ggplot2::coord_cartesian(xlim = xlim.plot, ylim = ylim.plot, expand = F)
g <- g + ggplot2::scale_colour_hue(l = 50)
g <- g + ggplot2::scale_fill_discrete(guide = NULL)
g
}
})


#' Add together two numbers.
Expand All @@ -210,10 +213,14 @@ plot_volcano <- function(
#' @return The sum of \code{x} and \code{y}.
#' @import data.table
#' @export
plot_fdr <- function(
data.fdr,
#' @include generics.R
#' @include seaMass_delta.R
setMethod("plot_fdr", "seaMass_delta", function(
object,
y.max = NULL,
y.col = "qvalue"
y.col = "qvalue",
data.fdr = group_quants_fdr(object),
output = "plotly"
) {
DT <- as.data.table(data.fdr)
DT <- DT[!is.na(get(y.col))]
Expand All @@ -227,17 +234,17 @@ plot_fdr <- function(
ylabels <- function() function(x) format(x, digits = 2)

g <- ggplot2::ggplot(DT, ggplot2::aes_string(x = "Discoveries", y = y.col))
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.frame(yintercept = 0.01), linetype = "dotted")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.frame(yintercept = 0.05), linetype = "dotted")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.frame(yintercept = 0.10), linetype = "dotted")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.table(yintercept = 0.01), linetype = "dotted")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.table(yintercept = 0.05), linetype = "dotted")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.table(yintercept = 0.10), linetype = "dotted")
g <- g + ggplot2::geom_step(direction = "vh")
g <- g + ggplot2::scale_x_continuous(expand = c(0, 0))
g <- g + ggplot2::scale_y_reverse(breaks = sort(c(pi, 0.0, 0.01, 0.05, 0.1, 0.2, 0.5, 1.0)), labels = ylabels(), expand = c(0.001, 0.001))
g <- g + ggplot2::coord_cartesian(xlim = c(0, xmax), ylim = c(y.max, 0))
g <- g + ggplot2::xlab("Number of Discoveries")
g <- g + ggplot2::ylab("False Discovery Rate")
g <- g + ggplot2::xlab("number of discoveries")
g <- g + ggplot2::ylab("qvalue")
g
}
})


#' Precision-Recall plot
Expand All @@ -247,12 +254,16 @@ plot_fdr <- function(
#' @return A ggplot2 object.
#' @import data.table
#' @export
plot_pr <- function(
data.fdr,
#' @include generics.R
#' @include seaMass_delta.R
setMethod("plot_pr", "seaMass_delta", function(
object,
plot.fdr = TRUE,
y.max = NULL,
legend.nrow = 1,
y.col = "qvalue"
y.col = "qvalue",
data.fdr = group_quants_fdr(object),
output = "plotly"
) {
if (is.data.frame(data.fdr)) {
DTs.pr <- list(unknown = data.fdr)
Expand Down Expand Up @@ -285,9 +296,9 @@ plot_pr <- function(
if (is.null(y.max)) y.max <- pi

g <- ggplot2::ggplot(DTs.pr, ggplot2::aes(x = TrueDiscoveries, y = FDP, colour = Method, fill = Method, linetype = Method))
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.frame(yintercept = 0.01), linetype = "dotted")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.frame(yintercept = 0.05), linetype = "dotted")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.frame(yintercept = 0.10), linetype = "dotted")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.table(yintercept = 0.01), linetype = "dotted")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.table(yintercept = 0.05), linetype = "dotted")
g <- g + ggplot2::geom_hline(ggplot2::aes(yintercept=yintercept), data.table(yintercept = 0.10), linetype = "dotted")
g <- g + ggplot2::geom_ribbon(ggplot2::aes(ymin = lower, ymax = upper), colour = NA, alpha = 0.3)
if (plot.fdr) g <- g + ggplot2::geom_line(ggplot2::aes(y = y), lty = "dashed")
g <- g + ggplot2::geom_step(direction = "vh")
Expand All @@ -303,4 +314,4 @@ plot_pr <- function(
} else {
g + ggplot2::theme(legend.position = "top") + ggplot2::guides(lty = ggplot2::guide_legend(nrow = legend.nrow))
}
}
})
13 changes: 6 additions & 7 deletions R/delta_process.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -34,8 +34,8 @@ setMethod("process", "seaMass_delta", function(object, chain, job.id) {

# markdown folder
report.index <- list()
root <- file.path(filepath(fit.theta), "markdown", "study")
dir.create(root, recursive = T, showWarnings = F)
root <- file.path(filepath(object), "markdown", "study")
dir.create(root, recursive = T)
group <- control(root(object))@group[1]

# calculate plot limits
Expand All @@ -62,11 +62,11 @@ setMethod("process", "seaMass_delta", function(object, chain, job.id) {

# plot
if ("group.quants.de.batch" %in% ctrl@plot) {
cat(paste0("[", Sys.time(), "] generating group quants differential expression plot...\n"))
group_quants_de(object, summary = T, as.data.table = T)
#cat(paste0("[", Sys.time(), "] generating group quants differential expression summary...\n"))
#group_quants_de(object, summary = T, as.data.table = T)

cat(paste0("[", Sys.time(), "] generating group quants differential expression plot...\n"))
text <- paste0(group, "differential expression for '", gsub("\\.", "' effect, comparison '", batch), "'", name)
text <- paste0(group, " differential expression for '", gsub("\\.", "' effect, comparison '", batch), "'", name)
report.index$assay.stdevs <- data.table(
section = "Study-level", section.order = 0, item = text, item.order = 75000,
item.href = generate_markdown(
Expand All @@ -82,10 +82,9 @@ setMethod("process", "seaMass_delta", function(object, chain, job.id) {

# zip
render_markdown(object, root)
if (!("markdown" %in% ctrl@keep)) unlink(root, recursive = T)

# save index
fst::write.fst(rbindlist(report.index), file.path(filepath(fit.theta), "report", "study.report.fst"))
fst::write.fst(rbindlist(report.index), file.path(filepath(object), "report", "study.report.fst"))

increment_completed(filepath(object))
}
Expand Down
3 changes: 3 additions & 0 deletions R/generics.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -52,10 +52,13 @@ if (!exists("plot_group_quants")) setGeneric("plot_group_quants", function(objec
if (!exists("plot_group_quants_de")) setGeneric("plot_group_quants_de", function(object, ...) standardGeneric("plot_group_quants_de"))
if (!exists("plot_group_quants_fdr")) setGeneric("plot_group_quants_fdr", function(object, ...) standardGeneric("plot_group_quants_fdr"))
if (!exists("plot_group_standards")) setGeneric("plot_group_standards", function(object, ...) standardGeneric("plot_group_standards"))
if (!exists("plot_fdr")) setGeneric("plot_fdr", function(object, ...) standardGeneric("plot_fdr"))
if (!exists("plot_measurement_means")) setGeneric("plot_measurement_means", function(object, ...) standardGeneric("plot_measurement_means"))
if (!exists("plot_measurement_stdevs")) setGeneric("plot_measurement_stdevs", function(object, ...) standardGeneric("plot_measurement_stdevs"))
if (!exists("plot_pr")) setGeneric("plot_pr", function(object, ...) standardGeneric("plot_pr"))
if (!exists("plot_priors")) setGeneric("plot_priors", function(object, ...) standardGeneric("plot_priors"))
if (!exists("plot_robust_pca")) setGeneric("plot_robust_pca", function(object, ...) standardGeneric("plot_robust_pca"))
if (!exists("plot_volcano")) setGeneric("plot_volcano", function(object, ...) standardGeneric("plot_volcano"))
if (!exists("plots")) setGeneric("plots", function(object, ...) standardGeneric("plots"))
if (!exists("prepare_delta")) setGeneric("prepare_delta", function(object, ...) standardGeneric("prepare_delta"))
if (!exists("prepare_sigma")) setGeneric("prepare_sigma", function(object, ...) standardGeneric("prepare_sigma"))
Expand Down
40 changes: 29 additions & 11 deletions R/report.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -34,7 +34,7 @@ setMethod("render_markdown", "seaMass", function(object, root) {
# generate _site.yml
cat(paste0(
'navbar:\n',
' title: ', name(seaMass::root(object)), '\n',
' title: ', name(root(object)), '\n',
' left:\n',
' - text: "Index"\n',
' href: index.html\n',
Expand Down Expand Up @@ -65,7 +65,11 @@ setMethod("render_markdown", "seaMass", function(object, root) {
zip::zip(file.path(path, paste0(basename(root), ".report.zip")), files, compression_level = 6, include_directories = F, mode = "cherry-pick")

# tidy up
unlink(file.path(root, "_site"), recursive = T)
if (!("markdown" %in% ctrl@keep)) {
unlink(root, recursive = T)
} else {
unlink(file.path(root, "_site"), recursive = T)
}

return(invisible(NULL))
})
Expand All @@ -76,17 +80,17 @@ setMethod("render_markdown", "seaMass", function(object, root) {
#' @import data.table
#' @export
#' @include generics.R
setMethod("assemble_report", "seaMass_sigma", function(object, filepath = "report", title = seaMass::name(object)) {
setMethod("assemble_report", "seaMass_sigma", function(object, filename = paste0(name(object), "_seaMass_report.zip"), title = seaMass::name(object)) {
ctrl <- control(object)

# create markdown directory for index
root <- file.path(dirname(filepath(object)), "markdown")
if (file.exists(root)) unlink(root, recursive = T)
dir.create(root)
dir.create(file.path(root, "_site", name(object)), recursive = T)

# generate _site.yml
cat(paste0(
#'name: ', name, '\n',
'output_dir: ', file.path("_site", name(object)), '\n',
'navbar:\n',
' title: ', title, '\n',
' left:\n',
Expand Down Expand Up @@ -139,24 +143,38 @@ setMethod("assemble_report", "seaMass_sigma", function(object, filepath = "repor

# render site
rmarkdown::render_site(root, quiet = T)
unlink(file.path(root, "_site", name(object), "fig.html"))

# root html index
cat(paste0(
'<!DOCTYPE HTML>\n',
'<meta charset="UTF-8">\n',
'<meta http-equiv="refresh" content="1; url=', name(object), '/index.html">\n',
'<script>\n',
' window.location.href = "', name(object), '/index.html"\n',
'</script>\n',
'<title>Page Redirection</title>\n',
'If you are not redirected automatically, follow the <a href="', name(object), '/index.html">link to ', filename, '</a>\n'
), file = file.path(root, "_site", "index.html"))

# zip
zipfile <- file.path(dirname(filepath(object)), paste0(filepath, ".zip"))
files <- file.path(root, "_site", setdiff(list.files(file.path(root, "_site")), c("fig.html")))
zip::zip(zipfile, files, compression_level = 6, include_directories = F, mode = "cherry-pick")
zipfile <- file.path(dirname(filepath(object)), filename)
files <- c(file.path(root, "_site", "index.html"), file.path(root, "_site", name(object)))
zip::zip(zipfile, files, compression_level = 6, mode = "cherry-pick")

# tidy up
unlink(file.path(root, "_site"), recursive = T)

# append zips
files <- list.files(dirname(filepath(object)), pattern = "^.*\\.report\\.zip$", recursive = T, full.names = T)
for (zipfile1 in files) {
root1 <- file.path(dirname(filepath(object)), "markdown", basename(tools::file_path_sans_ext(zipfile1)))
root1 <- file.path(dirname(filepath(object)), "markdown", "_site", name(object))
zip::unzip(zipfile1, exdir = root1)
files1 <- list.files(root1, full.names = T)
zip::zip_append(zipfile, files1, compression_level = 6, include_directories = F, mode = "cherry-pick")
zip::zip_append(zipfile, root1, compression_level = 6, include_directories = F, mode = "cherry-pick")
unlink(root1, recursive = T)
}

if (!("markdown" %in% ctrl@keep)) unlink(root, recursive = T)

return(invisible(NULL))
})
Loading

0 comments on commit 03d0a57

Please sign in to comment.