#' @title Plot the results of a hierarchical Bayesian model
#' @description Plot the results of a hierarchical Bayesian model
#' @param x An object of class hierarchicalBayesianModel (generated by \code{\link[mpMap2]{fitClusterModel}}) or magicHeuristicHBC (generated by \code{\link[mpMap2]{runHeuristics}}).
#' @param allFounderNames Optional list of the names of the founder lines replicates within the genetic data.
#' @param chainIndex The index of the MCMC chain to show.
#' @param ... Extra arguments to these functions are ignored.
#' @details
#' This function plots the result of using the end-point of one of the Markov chains as a fitted model. Homozygote clusters are shown in red and green, and the heterozygote cluster in dark blue.
#' Light blue represents points falling in the outlier cluster. Black points are unclassified. Ovals represent the covariance matrixces of the homozygote and heterozygote clusters.
#' If Input allfounderNames is input, then the founding lines of the population are given a different shape, and shown in black.
#' @method plot hierarchicalBayesianModel
#' @S3method plot hierarchicalBayesianModel
#' @rdname plot.model
#' @examples
#' data("eightWayExampleData", package="magicCalling")
#' data <- eightWayExampleData[[1]]
#' meanY <- mean(data[,2])
#' startingPoints <- list(
#' rbind(c(0.5, meanY), c(0.5, meanY)),
#' rbind(c(0.5, meanY), c(0.5, meanY)),
#' rbind(c(0.25, meanY), c(0.5, meanY)),
#' rbind(c(0.25, meanY), c(0.5, meanY)),
#' rbind(c(0.75, meanY), c(0.5, meanY)),
#' rbind(c(0.75, meanY), c(0.5, meanY)),
#' rbind(c(0.8, meanY), c(0.2, meanY)),
#' rbind(c(0.8, meanY), c(0.2, meanY))
#' )
#' result <- fitClusterModel(data, startingPoints, n.iter = 200, D_hom = diag(2)*4, V_hom = cbind(c(0.005, 0), c(0, 0.1))/3, n_hom = 30, D_err = diag(2), V_err = diag(2)*10/3, n_err = 300, V_het = diag(2)*0.025/3, n_het = 1500)
#' plot(result, chainIndex = 1)
#' heuristicResults <- runHeuristics(result, minHomozygoteSize = 200)
#' plot(heuristicResults, chainIndex = heuristicResults$chainIndex)
#' @export
plot.hierarchicalBayesianModel <- function(x, allFounderNames, chainIndex, ...)
{
heuristicResults <- runHeuristics(x, minHomozygoteSize = 200)
if(missing(chainIndex) && length(heuristicResults$classification) == 1)
{
chainIndex <- 1
}
if(missing(chainIndex))
{
stop("Input chainIndex is required, if there are multiple Markov chains")
}
data <- heuristicResults$data
ellipse1 <- as.data.frame(ellipse(center = heuristicResults$clusterMeans[[chainIndex]][1, ], shape = heuristicResults$covariances[[chainIndex]][1,,], radius=5, col = 2, center.pch = NULL, draw = FALSE))
ellipse2 <- as.data.frame(ellipse(center = heuristicResults$clusterMeans[[chainIndex]][2, ], shape = heuristicResults$covariances[[chainIndex]][2,,], radius=5, col = 3, center.pch = NULL, draw = FALSE))
ellipse3 <- as.data.frame(ellipse(center = heuristicResults$clusterMeans[[chainIndex]][3, ], shape = heuristicResults$covariances[[chainIndex]][3,,], radius=5, col = 4, center.pch = NULL, draw = FALSE))
if(missing(allFounderNames))
{
data <- data.frame(x = data[,1], y = data[,2], colour = factor(heuristicResults$classification[[chainIndex]], levels = 1:5), row.names = rownames(data))
basePlot <- ggplot2::ggplot(mapping = ggplot2::aes_string(x = "x", y = "y"), data = data) + ggplot2::geom_point(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour"), shape = 16)
built <- ggplot2::ggplot_build(basePlot)
y.range <- built$layout$panel_ranges[[1]]$y.range
x.range <- built$layout$panel_ranges[[1]]$x.range
basePlot + ggplot2::geom_path(data = ellipse1, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 2) + ggplot2::geom_path(data = ellipse2, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 3) + ggplot2::geom_path(data = ellipse3, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 4) + ggplot2::theme_bw() + ggplot2::scale_color_manual(values = palette(), guide = FALSE) + ggplot2::scale_shape_identity() + coord_cartesian(xlim = x.range, ylim = y.range, expand = FALSE)
}
else
{
data <- data.frame(x = data[,1], y = data[,2], colour = factor(heuristicResults$classification[[chainIndex]], levels = 1:5), shape = 16, row.names = rownames(data))
founderIndices <- na.omit(match(allFounderNames, rownames(data)))
data[founderIndices, "shape"] <- 2
data[founderIndices, "colour"] <- 1
#Split out the founder indices subset, because these will have different plot symbols, which we want to appear *on top* of everything else.
founderIndicesSubset <- data[data[,"shape"] == 2, ]
basePlot <- ggplot2::ggplot(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour"), data = data) + ggplot2::geom_point(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour", shape = "shape")) + ggplot2::scale_shape_identity()
built <- ggplot2::ggplot_build(basePlot)
y.range <- built$layout$panel_ranges[[1]]$y.range
x.range <- built$layout$panel_ranges[[1]]$x.range
basePlot + ggplot2::geom_path(data = ellipse1, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 2) + ggplot2::geom_path(data = ellipse2, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 3) + ggplot2::geom_path(data = ellipse3, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 4) + ggplot2::theme_bw() + ggplot2::scale_color_manual(values = palette(), guide = FALSE) + ggplot2::geom_point(data = founderIndicesSubset) + coord_cartesian(xlim = x.range, ylim = y.range, expand = FALSE)
}
}
#' @rdname plot.model
#' @method plot magicHeuristicHBC
#' @S3method plot magicHeuristicHBC
#' @export
plot.magicHeuristicHBC <- function(x, allFounderNames, chainIndex, ...)
{
if(missing(chainIndex) && length(x$classifications) == 1)
{
chainIndex <- 1
}
if(missing(chainIndex) && x$chainIndex == -1)
{
stop("Input chainIndex is required, if there are multiple Markov chains, and no chain was assessed as being good quality.")
}
if(missing(chainIndex) && x$chainIndex != -1) chainIndex <- x$chainIndex
heuristicResults <- x
data <- heuristicResults$data
ellipse1 <- as.data.frame(ellipse(center = heuristicResults$clusterMeans[[chainIndex]][1, ], shape = heuristicResults$covariances[[chainIndex]][1,,], radius=5, col = 2, center.pch = NULL, draw = FALSE))
ellipse2 <- as.data.frame(ellipse(center = heuristicResults$clusterMeans[[chainIndex]][2, ], shape = heuristicResults$covariances[[chainIndex]][2,,], radius=5, col = 3, center.pch = NULL, draw = FALSE))
ellipse3 <- as.data.frame(ellipse(center = heuristicResults$clusterMeans[[chainIndex]][3, ], shape = heuristicResults$covariances[[chainIndex]][3,,], radius=5, col = 4, center.pch = NULL, draw = FALSE))
if(missing(allFounderNames))
{
data <- data.frame(x = data[,1], y = data[,2], colour = factor(heuristicResults$classifications[[chainIndex]], levels = 1:5), row.names = rownames(data))
basePlot <- ggplot2::ggplot(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour"), data = data) + ggplot2::geom_point(shape = 16)
built <- ggplot2::ggplot_build(basePlot)
y.range <- built$layout$panel_ranges[[1]]$y.range
x.range <- built$layout$panel_ranges[[1]]$x.range
basePlot + ggplot2::geom_path(data = ellipse1, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 2) + ggplot2::geom_path(data = ellipse2, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 3) + ggplot2::geom_path(data = ellipse3, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 4) + ggplot2::theme_bw() + ggplot2::scale_color_manual(values = palette(), guide = FALSE) + ggplot2::scale_shape_identity() + coord_cartesian(xlim = x.range, ylim = y.range, expand = FALSE)
}
else
{
data <- data.frame(x = data[,1], y = data[,2], colour = factor(heuristicResults$classifications[[chainIndex]], levels = 1:5), shape = 16, row.names = rownames(data))
founderIndices <- na.omit(match(allFounderNames, rownames(data)))
data[founderIndices, "shape"] <- 2
data[founderIndices, "colour"] <- 1
#Split out the founder indices subset, because these will have different plot symbols, which we want to appear *on top* of everything else.
founderIndicesSubset <- data[data[,"shape"] == 2,]
basePlot <- ggplot2::ggplot(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour"), data = data) + ggplot2::geom_point(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour", shape = "shape")) + ggplot2::scale_shape_identity()
built <- ggplot2::ggplot_build(basePlot)
y.range <- built$layout$panel_ranges[[1]]$y.range
x.range <- built$layout$panel_ranges[[1]]$x.range
basePlot + ggplot2::geom_path(data = ellipse1, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 2) + ggplot2::geom_path(data = ellipse2, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 3) + ggplot2::geom_path(data = ellipse3, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 4) + ggplot2::theme_bw() + ggplot2::scale_color_manual(values = palette(), guide = FALSE) + ggplot2::geom_point(data = founderIndicesSubset) + ggplot2::scale_shape_identity() + coord_cartesian(xlim = x.range, ylim = y.range, expand = FALSE)
}
}
#' @title Plot a manually called marker
#' @description Plot the results of manually calling a marker
#' @param x An object of class markerResult.
#' @param allFounderNames Optional list of the names of the founder lines replicates within the genetic data.
#' @param ... Extra arguments to these functions are ignored.
#' @S3method plot markerResult
#' @details
#' This function plots the results of calling a marker manually using interactiveCall. In these cases a marker might have been called using
#' the hierarchical Bayesia approach, or using DBSCAN.
#' @export
plot.markerResult <- function(x, allFounderNames, ...)
{
data <- x$data
if(x$hasVariability)
{
if("dbscan" %in% names(x) && x$dbscan)
{
if(missing(allFounderNames))
{
data <- data.frame(x = data[,1], y = data[,2], colour = factor(x$classification + 1, levels = 1:(max(x$classification)+1)))
ggplot2::ggplot(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour"), data = data) + ggplot2::geom_point(shape = 16) + ggplot2::theme_bw() + ggplot2::scale_color_manual(values = palette(), guide = FALSE)
}
else
{
founderIndices <- na.omit(match(allFounderNames, rownames(data)))
data <- data.frame(x = data[,1], y = data[,2], colour = factor(x$classification+1, levels = 1:(max(x$classification)+1)), shape = 16, row.names = rownames(data))
data[founderIndices, "shape"] <- 2
data[founderIndices, "colour"] <- 1
#Split out the founder indices subset, because these will have different plot symbols, which we want to appear *on top* of everything else.
founderIndicesSubset <- data[data[,"shape"] == 2,]
ggplot2::ggplot(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour", shape = "shape"), data = data) + ggplot2::geom_point() + ggplot2::theme_bw() + ggplot2::scale_color_manual(values = palette(), guide = FALSE) + ggplot2::geom_point(data = founderIndicesSubset) + ggplot2::scale_shape_identity()
}
}
else
{
ellipse1 <- as.data.frame(ellipse(center = x$clusterMeans[1, ], shape = x$covariances[1,,], radius=5, col = 2, center.pch = NULL, draw = FALSE))
ellipse2 <- as.data.frame(ellipse(center = x$clusterMeans[2, ], shape = x$covariances[2,,], radius=5, col = 3, center.pch = NULL, draw = FALSE))
ellipse3 <- as.data.frame(ellipse(center = x$clusterMeans[3, ], shape = x$covariances[3,,], radius=5, col = 4, center.pch = NULL, draw = FALSE))
if(missing(allFounderNames))
{
data <- data.frame(x = data[,1], y = data[,2], colour = factor(x$classification, levels = 1:5), row.names = rownames(data))
basePlot <- ggplot2::ggplot(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour"), data = data) + ggplot2::geom_point(shape = 16)
built <- ggplot2::ggplot_build(basePlot)
y.range <- built$layout$panel_ranges[[1]]$y.range
x.range <- built$layout$panel_ranges[[1]]$x.range
basePlot + ggplot2::geom_path(data = ellipse1, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 2) + ggplot2::geom_path(data = ellipse2, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 3) + ggplot2::geom_path(data = ellipse3, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 4) + ggplot2::theme_bw() + ggplot2::scale_color_manual(values = palette(), guide = FALSE) + coord_cartesian(xlim = x.range, ylim = y.range, expand = FALSE)
}
else
{
data <- data.frame(x = data[,1], y = data[,2], colour = factor(x$classification, levels = 1:5), shape = 16, row.names = rownames(data))
founderIndices <- na.omit(match(allFounderNames, rownames(data)))
data[founderIndices, "shape"] <- 2
data[founderIndices, "colour"] <- 1
#Split out the founder indices subset, because these will have different plot symbols, which we want to appear *on top* of everything else.
founderIndicesSubset <- data[data[,"shape"] == 2,]
basePlot <- ggplot2::ggplot(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour"), data = data) + ggplot2::geom_point(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour", shape = "shape")) + ggplot2::scale_shape_identity()
built <- ggplot2::ggplot_build(basePlot)
y.range <- built$layout$panel_ranges[[1]]$y.range
x.range <- built$layout$panel_ranges[[1]]$x.range
notFoundersIndicesSubset <- data[data[,"shape"] != 2,]
ggplot2::ggplot(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour"), data = notFoundersIndicesSubset) + ggplot2::geom_point(mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour", shape = "shape")) + ggplot2::scale_shape_identity() + ggplot2::geom_path(data = ellipse1, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 2) + ggplot2::geom_path(data = ellipse2, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 3) + ggplot2::geom_path(data = ellipse3, mapping = ggplot2::aes_string(x = "x", y = "y"), colour = 4) + ggplot2::theme_bw() + ggplot2::scale_color_manual(values = palette(), guide = FALSE) + ggplot2::geom_point(data = founderIndicesSubset, mapping = ggplot2::aes_string(x = "x", y = "y", colour = "colour", shape = "shape")) + coord_cartesian(xlim = x.range, ylim = y.range, expand = FALSE)
}
}
}
else
{
data <- data.frame(x = data[,1], y = data[,2], row.names = rownames(data))
ggplot2::ggplot(mapping = ggplot2::aes_string(x = "x", y = "y"), data = data) + ggplot2::geom_point(colour = 1, shape = 16) + ggplot2::scale_shape_identity()
}
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.