Nothing
R/plot.dispRity_fun.R
In dispRity: Measuring Disparity
Defines functions
plot.projection
plot.axes
plot.test.metric
plot.model.sim
plot.model.test
plot.dtt
plot.randtest
plot.preview
plot.rarefaction
plot.continuous
plot.discrete
plot.elements
plot.observed
get.quantile.col
get.shift
get.dots
get.plot.params
get.data.params
## Get the disparity characteristics
get.data.params <- function(data) {
return(list(
"distribution" = length(data$disparity[[1]]$elements) != 1,
"bootstrap" = !is.null(data$call$bootstrap) && data$call$bootstrap[[2]] != "covar",
"rarefaction" = data$call$bootstrap[[3]],
"between.groups" = data$call$disparity$metrics$between.groups,
"elements" = names(data$disparity)
))
}
## Setting up all the default plot parameters:
# $disparity: summarised_data
# $helpers: n_quantiles, npoints
# $options: ylim, ylab, xlab, col, ...
# $observed: data, ylim, ylab, xlab, col, ...
# $elements: elements, ylim, ylab, xlab, col, ...
get.plot.params <- function(data, data_params, cent.tend, quantiles, rarefaction_level, elements_args, type, observed_args, ...) {
## Set up the plotting data
## Summarise the data
if(is.null(data_params$model.sim)) {
# if((!is.na(data$call$subsets["trees"])) && (as.numeric(data$call$subsets["trees"]) > 1)) {
# summarised_data <- summary.dispRity(data, quantiles = quantiles, cent.tend = cent.tend, digits = 5, na.rm = TRUE)
# } else {
summarised_data <- summary.dispRity(data, quantiles = quantiles, cent.tend = cent.tend, digits = 5, na.rm = TRUE)
# }
} else {
summarised_data <- summary.dispRity(data, quantiles = quantiles, cent.tend = cent.tend, digits = 5)
## Making a fake obs column
colnames(summarised_data)[3] <- "obs"
}
## Find the observed data rows (that are not rarefactions)
observed_data <- !is.na(summarised_data[,(2+ifelse(data_params$between.groups, 2, 1))])
if(any(!observed_data)) {
## Find NAs in the last column
last_col_na <- is.na(summarised_data[, ncol(summarised_data)])
## If there are any NA in the last column, the observed NA is a true NA (nor rarefaction)
observed_data <- last_col_na | observed_data
}
## Get a specific rarefaction level
if(!is.null(rarefaction_level)) {
## Rarefaction level
if(data_params$between.groups) {
## Find the matching rarefaction level
rarefaction_match <- which(summarised_data$n_1 == rarefaction_level & summarised_data$n_2 == rarefaction_level)
} else {
rarefaction_match <- which(summarised_data$n == rarefaction_level)
}
## Get both the observed and rarefied data
selected_data <- summarised_data[rarefaction_match ,]
} else {
## Just get the observed data
selected_data <- summarised_data[observed_data,]
}
## Separate the data elements
disparity <- list()
name_part <- c(1:ifelse(data_params$between.groups, 3, 2))
disparity$names <- selected_data[, name_part, drop = FALSE]
disparity$data <- selected_data[, -name_part, drop = FALSE]
## Update the data to be in boxplot format
if(type == "box") {
## Extracting boxplot data
if(is.null(rarefaction_level)) {
## Getting the observed or bootstrapped data
if(data_params$bootstrap && !observed_args$observed) {
## Getting the bootstrapped data
box_data <- unlist(get.disparity(data, observed = FALSE), recursive = FALSE)
} else {
if(data_params$distribution || data_params$between.groups) {
box_data <- lapply(get.disparity(data, observed = TRUE, concatenate = FALSE), c)
} else {
if(ncol(disparity$data) == 1) {
box_data <- get.disparity(data, observed = TRUE)
} else {
box_data <- unlist(get.disparity(data, observed = FALSE), recursive = FALSE)
}
}
}
} else {
## Find the correct rarefaction level
box_data <- unlist(get.disparity(data, observed = FALSE, rarefaction = rarefaction_level), recursive = FALSE)
}
## Updating the disparity data part
if(data_params$between.groups) {
names(box_data) <- data_params$elements
}
disparity$data <- box_data
}
## Set up the helpers options
helpers <- list()
## Detect the number of quantiles
helpers$n_quantiles <- length(quantiles)
## Detect the number of points
if(is.null(data_params$model.sim)) {
helpers$n_points <- length(data$disparity)
} else {
## Selecting helpers for model.sim
helpers$n_points <- nrow(disparity$data)
}
## Set up the plotting options
dots <- list(...)
options <- list()
## Set the xlabel
if(is.null(dots$xlab)) {
## Default is subsets
options$xlab <- "Subsets"
## Default chrono.subset label
if(type == "continuous" && "continuous" %in% data$call$subsets) {
options$xlab <- "Time (Mya)"
}
## Default rarefaction label
if(type == "rarefaction") {
options$xlab <- "Elements"
}
} else {
## User input
check.class(dots$xlab, "character", " must be a character string.")
check.length(dots$xlab, 1, " must be a character string.")
options$xlab <- dots$xlab
dots$xlab <- NULL
}
## Set the ylabel
if(is.null(dots$ylab)) {
## Default is the metric name
options$ylab <- as.character(data$call$disparity$metrics$name)
} else {
## User input
check.class(dots$ylab, "character", " must be a character string.")
if(length(dots$ylab) > 2) stop.call("", "ylab can have maximum of two elements.")
options$ylab <- dots$ylab
dots$ylab <- NULL
}
## Set the y limits
if(is.null(dots$ylim)) {
if(type != "rarefaction") {
## Get the range of the data
options$ylim <- range(disparity$data, na.rm = TRUE)
## Add 2% on each side
percent_change <- 0.02
if(options$ylim[1] > 0) {
options$ylim[1] <- options$ylim[1] - options$ylim[1]*percent_change
} else {
options$ylim[1] <- options$ylim[1] + options$ylim[1]*percent_change
}
options$ylim[2] <- options$ylim[2] + options$ylim[2]*percent_change
} else {
options$ylim <- "rarefaction"
}
} else {
## User input
check.class(dots$ylim, c("numeric", "integer"))
check.length(dots$ylim, 2, " must be a vector of two elements.")
options$ylim <- dots$ylim
dots$ylim <- NULL
}
## Set the colours
if(is.null(dots$col)) {
## For boxplots, the default colour is white
if(type == "box" && type != "rarefaction") {
options$col <- "white"
} else {
options$col <- "black"
}
} else {
## User input
check.class(dots$col, "character", " must be a character string.")
options$col <- dots$col
dots$col <- NULL
}
## Add potential missing colours
if(helpers$n_quantiles > 0 && type != "box") {
## Count if they are any missing colours
cols_missing <- (helpers$n_quantiles + 1) - length(options$col)
## Adding grey scales if quantiles
if(cols_missing > 0) {
colfun <- grDevices::colorRampPalette(c("grey", "lightgrey"))
options$col <- c(options$col, colfun(cols_missing))
}
}
## Add additional options
if(length(dots) > 0) {
options <- c(options, dots)
}
## Observed data
if(observed_args$observed) {
## Adding the observed data
selected_data <- summarised_data[observed_data,]
observed_args$names <- selected_data[, name_part]
observed_args$data <- selected_data[, -name_part]
## Default observed arguments
if(is.null(observed_args$col)) {
if(type == "box" && options$col[[1]] == "white") {
observed_args$col <- "black"
} else {
observed_args$col <- options$col[[1]]
}
}
if(is.null(observed_args$pch)) {
observed_args$pch <- 4
}
if(is.null(observed_args$cex)) {
observed_args$cex <- 1
}
}
## Observed data
if(elements_args$elements) {
## Default observed arguments
if(is.null(elements_args$col)) {
if(type == "box" && options$col[[1]] == "white") {
elements_args$col <- c("black", "darkgrey")
} else {
elements_args$col <- options$col
}
}
if(length(elements_args$col) < 2) {
elements_args$col <- rep(elements_args$col, 2)
}
if(is.null(elements_args$pch)) {
elements_args$pch <- c(15, 15)
} else {
if(length(elements_args$pch) < 2) {
elements_args$pch <- rep(elements_args$pch, 2)
}
}
if(is.null(elements_args$lty)) {
elements_args$lty <- c(2, 2)
} else {
if(length(elements_args$lty) < 2) {
elements_args$lty <- rep(elements_args$lty, 2)
}
}
}
## Output the finalised list
return(list("disparity" = disparity, "helpers" = helpers, "options" = options, "observed_args" = observed_args, "elements_args" = elements_args))
}
## Get dots options
# @param dots ...
# @param args the list of arguments to be affected by dots
# @param name the argument name
# @param default an optional default argument (default is NULL)
# @param fun an optional suffix argument for the function to be used (default is NULL)
get.dots <- function(dots, args, name, default = NULL, fun = NULL) {
## Override the name without fun in dots
if(!is.null(fun) && !is.null(names(dots))) {
name_plus <- paste0(fun, ".", name)
if(name_plus %in% names(dots)) {
dots[[name]] <- dots[[name_plus]]
}
}
## Fill in the dots in the arguments list
if(!is.null(default)) {
if(is.null(dots[[name]])) {
## Set the default
args[[name]] <- default
} else {
## Use the dots
args[[name]] <- dots[[name]]
}
} else {
## Use the dots
if(!is.null(dots[[name]])) {
args[[name]] <- dots[[name]]
}
}
return(args)
}
## Handle the shift argument
get.shift <- function(add, plot_params) {
## Set the shift parameter (for add)
shift = 0
if(add) {
## Is the previous plot the same size?
prev_axis <- par("xaxp")
if(prev_axis[2] == plot_params$helpers$n_points) {
shift = 0
} else {
shift = 0.5
}
}
return(shift)
}
## Getting the columns for the right quantiles
get.quantile.col <- function(cent_tend_col, cis, n_quantiles) {
return(c(
## Lower quantile
cent_tend_col + cis,
## Higher quantile
(cent_tend_col + n_quantiles*2) - (cis-1)
))
}
## Observed points
plot.observed <- function(plot_params) {
if(plot_params$observed_args$observed) {
## Set the observed arguments
points_args <- plot_params$observed_args
points_args$x <- 1:plot_params$helpers$n_points
points_args$y <- points_args$data[, "obs"]
points_args$observed <- NULL
points_args$names <- NULL
points_args$data <- NULL
## Calling the points
do.call(points, points_args)
} else {
return(invisible())
}
}
## Plot elements
plot.elements <- function(plot_params, data_params, type) {
## Elements plots
add_args <- plot_params$elements_args
add_args$elements <- NULL
## Getting the scaling value
percent_change <- 0.02
rescaling <- plot_params$options$ylim
if(rescaling[1] >= 0) {
rescaling[1] <- rescaling[1] + percent_change * rescaling[1]
} else {
rescaling[1] <- rescaling[1] - percent_change * rescaling[1]
}
rescaling[2] <- rescaling[2] - percent_change * rescaling[2]
## Rescale the y values
rescaled_elements <- plot_params$disparity$names[,-1, drop = FALSE]
rescaled_elements <- matrix(scales::rescale(unlist(rescaled_elements), to = rescaling), ncol = ncol(rescaled_elements), byrow = FALSE)
## Add the values
add_args$x <- 1:plot_params$helpers$n_points
add_args$y <- rescaled_elements[,1]
## Select the first colour
add_args$col <- plot_params$elements_args$col[1]
## Check if ylab exists
if(is.null(plot_params$elements_args$ylab)) {
axis_label <- ifelse(length(plot_params$options$ylab) > 1, plot_params$options$ylab[2], "Elements")
} else {
axis_label <- plot_params$elements_args$ylab
}
## Add the different types
if(type == "continuous") {
add_args$lty <- plot_params$elements_args$lty[1]
do.call(lines, add_args)
} else {
add_args$pch <- plot_params$elements_args$pch[1]
do.call(points, add_args)
}
## Add the ones for the second group (if two groups)
if(data_params$between.groups) {
## Get the n_2 column
add_args$y <- rescaled_elements[,2]
## Select the second colour
add_args$col <- plot_params$elements_args$col[2]
## Add the different types
if(type == "continuous") {
add_args$lty <- plot_params$elements_args$lty[2]
do.call(lines, add_args)
} else {
## Slightly shift
add_args$x <- add_args$x + (mode.val(diff(add_args$x))*0.01)
add_args$pch <- plot_params$elements_args$pch[2]
do.call(points, add_args)
}
}
## Add the axis label
seq.along.range <- function(range, out = 5) {
seq(from = range(range)[1],
to = range(range)[2],
length.out = out)
}
## Add the second y axis
axis(4, at = seq.along.range(rescaling), labels = seq.along.range(plot_params$disparity$names[,-1]),
lty = plot_params$elements_args$lty)
## Add the second axis label
mtext(axis_label, side = 4, line = 2)
}
## discrete plotting
plot.discrete <- function(plot_params, data_params, add, density, type) {
## Get the shifting argument
shift <- get.shift(add, plot_params)
## Select the central tendency column
cent_tend_col <- ifelse(data_params$bootstrap, 2, 1)
## Creating the matrix frame (for neat plotting)
frame_matrix <- matrix(1:plot_params$helpers$n_points,
ncol = plot_params$helpers$n_point,
dimnames = list(c(), plot_params$disparity$names[,1]))
if(!add) {
## Set the frame boxplot arguments
box_args <- plot_params$options
box_args$x <- frame_matrix
box_args$col <- box_args$border <- "white"
box_args$ylab <- box_args$ylab[[1]]
box_args$boxwex <- 0.00001
box_args$type <- "n"
## Plot the boxplot
do.call(boxplot, box_args)
}
## Check if bootstrapped
if(data_params$bootstrap || data_params$distribution) {
## Set the width (default)
width <- 0.5
## Initialise the plot arguments
plot_args <- plot_params$options
## Add the quantiles
if(type == "polygon") {
## Add the polygon options
plot_args$border <- plot_params$options$col[[1]]
plot_args$density <- density
## Loop through each point
for(point in 1:plot_params$helpers$n_points) {
## Loop through each quantile
for(cis in 1:plot_params$helpers$n_quantiles) {
## Add the colour option
plot_args$col <- plot_params$options$col[[cis+1]]
## Set the bars width
point_width <- width/(plot_params$helpers$n_quantiles - cis + 1.5)
## Set the x coordinates
plot_args$x <- c(point - point_width,
point + point_width,
point + point_width,
point - point_width) + shift
## Select the quantiles columns
quantiles_col <- get.quantile.col(cent_tend_col, cis, plot_params$helpers$n_quantiles)
## Set the y coordinates
plot_args$y <- c(rep(plot_params$disparity$data[point, quantiles_col[1]], 2),
rep(plot_params$disparity$data[point, quantiles_col[2]], 2))
## Add the polygon
do.call(polygon, plot_args)
}
}
}
if(type == "line") {
## Loop through each point
for(point in 1:plot_params$helpers$n_points) {
## Loop through the quantiles
for(cis in 1:plot_params$helpers$n_quantiles) {
## The the line type
plot_args$lty <- plot_params$helpers$n_quantiles - cis + 1
plot_args$lwd <- 1 + cis * 1.5
plot_args$col <- plot_params$options$col[[cis+1]]
## Set the x values
plot_args$x <- rep(point, 2)
## Select the quantiles columns
quantiles_col <- get.quantile.col(cent_tend_col, cis, plot_params$helpers$n_quantiles)
## Set the y value
plot_args$y <- c(plot_params$disparity$data[point, quantiles_col[1]],
plot_params$disparity$data[point, quantiles_col[2]])
## Plotting the line
do.call(lines, plot_args)
}
}
}
}
## Add the points estimates
point_args <- list()
## Get the points coordinates
point_args$x <- 1:plot_params$helpers$n_points + shift
point_args$y <- plot_params$disparity$data[, cent_tend_col]
## Get the options
point_args$col <- plot_params$options$col[[1]]
if(is.null(plot_params$options$pch)) {
point_args$pch <- 19
} else {
point_args$pch <- plot_params$options$pch
}
## Add the points
do.call(points, point_args)
## Save parameters
return(par())
}
## continuous plotting
plot.continuous <- function(plot_params, data_params, add, density) {
## Get the shifting argument
shift <- get.shift(add, plot_params)
## Select the central tendency column
cent_tend_col <- ifelse(data_params$bootstrap, 2, 1)
## Set up the cur plot options
plot_args <- plot_params$options
plot_args$x <- seq(from = 1, to = plot_params$helpers$n_points)-shift
plot_args$y <- plot_params$disparity$data[, cent_tend_col]
plot_args$type <- "l"
plot_args$col <- plot_params$options$col[[1]]
plot_args$xaxt <- "n"
## Plot the central tendency
if(!add) {
## Plot the thingy
do.call(plot, plot_args)
## Add the axis labels
options(warn = -1)
try_round <- as.numeric(plot_params$disparity$names[,"subsets"])
options(warn = 0)
if(all(is.na(try_round))) {
## Keep the ticks as they are
x_ticks <- plot_params$disparity$names[,"subsets"]
} else {
rounding <- 3-max(nchar(round(try_round)))
x_ticks <- round(try_round, digits = ifelse(rounding < 0, 0, rounding))
}
axis(1, 1:plot_params$helpers$n_points, x_ticks)
} else {
## Plot the line
do.call(lines, plot_args)
}
## Check if bootstrapped
if(data_params$bootstrap || data_params$distribution) {
## Setting the plot parameters
poly_args <- plot_params$options
poly_args$border <- "NA"
poly_args$density <- density
poly_args$x <- numeric()
poly_args$y <- numeric()
## Add the polygons
for (cis in 1:plot_params$helpers$n_quantiles) {
## Set the x values
poly_args$x <- c(1:plot_params$helpers$n_points)
poly_args$x <- c(poly_args$x, rev(poly_args$x))
## Select the quantiles columns
quantiles_col <- get.quantile.col(cent_tend_col, cis, plot_params$helpers$n_quantiles)
## Set the y values
poly_args$y <- plot_params$disparity$data[, quantiles_col[1]]
poly_args$y <- c(poly_args$y, rev(plot_params$disparity$data[, quantiles_col[2]]))
## Set up the colour
poly_args$col <- plot_params$options$col[cis+1]
## Dividing the polygon if NAs
if(any(is_nas <- is.na(poly_args$y[1:plot_params$helpers$n_points]))) {
## Selecting the groups of applicable data
groups <- numeric()
group_label <- 1
## (attributing a group as soon as it jumps an NA)
for(point in seq_along(is_nas)) {
if(is_nas[point]) {
## Increment the group
group_label <- group_label + 1
groups[point] <- NA
} else {
groups[point] <- group_label
}
}
## splitting the data into the groups
split.combine.data <- function(var, groups, n_points) {
vals_1 <- split(var[1:n_points], groups)
vals_2 <- split(var[-c(1:n_points)], rev(groups))
return(mapply(c, vals_1, vals_2, SIMPLIFY = FALSE))
}
y_vals <- split.combine.data(poly_args$y, groups, plot_params$helpers$n_points)
x_vals <- split.combine.data(poly_args$x, groups, plot_params$helpers$n_points)
## Combine args
combined.poly.args <- function(x, y, poly_args) {
poly_args$x <- x
poly_args$y <- y
return(poly_args)
}
list_poly_args <- mapply(combined.poly.args, x_vals, y_vals, SIMPLIFY = FALSE, MoreArgs = list("poly_args" = poly_args))
## Plot all the polygons
lapply(list_poly_args, function(args) do.call(polygon, args))
} else {
## Plot the polygon
do.call(polygon, poly_args)
}
}
## Add the central tendency
do.call(lines, plot_args)
}
## Save parameters
return(par())
}
## Plot the rarefaction
plot.rarefaction <- function(plot_params, data_params, data) {
## How many rarefaction plots?
n_plots <- length(data$subsets)
## Open the multiple plots
plot_size <- ifelse(n_plots == 3, 4, n_plots)
op_tmp <- par(mfrow = c(ceiling(sqrt(plot_size)),round(sqrt(plot_size))))
## Get the list of subsets
subsets_levels <- names(data$subsets)
## Get all the rarefaction values
rarefied_data <- summary(data)
## Get were the central tendency value is
cent_tend_col <- ifelse(data_params$between.groups, 5, 4)
## Setting the plotting arguments
all_plot_args <- plot_params$options
all_plot_args$ylim <- NULL
all_plot_args$lty <- 1
all_plot_args$type <- "l"
## Plot the different curves
for(one_subset in subsets_levels) {
## get the subset data
subset_data <- rarefied_data[rarefied_data[,1] == one_subset, ]
## Setting the plotting args for the specific subset
one_plot_args <- all_plot_args
one_plot_args$x <- subset_data[,2]
one_plot_args$y <- subset_data[,cent_tend_col]
one_plot_args$ylim <- range(subset_data[, -c(1:(cent_tend_col-1))])
one_plot_args$main <- one_subset
## Plot the central tendency
do.call(plot, one_plot_args)
## Add the quantiles
for(cis in 1:plot_params$helpers$n_quantiles) {
## Get the quantile columns
ci_cols <- get.quantile.col(cent_tend_col, cis, plot_params$helpers$n_quantiles)
## Set the plotting arguments
lines_args <- one_plot_args
lines_args$col <- all_plot_args$col[length(all_plot_args$col) - (cis-1)]
lines_args$y <- subset_data[,ci_cols[1]]
## plot the quantiles
do.call(lines, lines_args)
lines_args$y <- subset_data[,ci_cols[2]]
do.call(lines, lines_args)
}
}
## Done!
par(op_tmp)
}
## Plotting a space preview
plot.preview <- function(data, specific.args, ...) {
dots <- list(...)
## The "ggplot" colours
gg.color.hue <- function(n) {
grDevices::hcl(h = seq(15, 375, length = n + 1), l = 65, c = 100)[1:n]
}
make.transparent <- function(colour, levels) {
grDevices::adjustcolor(colour, alpha.f = (1/levels) + 1/3)
}
## Set up the specific args
specific.args <- get.dots(specific.args, specific.args, "matrix", c(1:length(data$matrix)))
specific.args <- get.dots(specific.args, specific.args, "dimensions", c(1:2))
## Plot the legend?
if(is.null(dots$legend)) {
plot_legend <- TRUE
dots$legend <- NULL
} else {
if(is.logical(dots$legend)) {
plot_legend <-dots$legend
dots$legend <- NULL
}
}
if(is.null(specific.args$tree) || is.null(data$tree[[1]])) {
## Don't plot trees
plot_trees <- FALSE
} else {
if(is(specific.args$tree, "logical")) {
## Toggle plot_trees
plot_trees <- specific.args$tree
if(plot_trees) {
## Select the trees to plot
specific.args$tree <- 1:length(data$tree)
}
} else {
## Specific args must be the specific trees to plot
plot_trees <- TRUE
}
}
## Capturing the dots options
plot_args <- c(list(x = NULL, y = NULL), dots)
## Removing specific args from dots
remove <- c(grep(c("legend"), names(plot_args)), grep(c("lines"), names(plot_args)), grep(c("points"), names(plot_args)))
plot_args[remove] <- NULL
## Getting the loadings
loading <- apply(do.call(rbind, lapply(data$matrix[specific.args$matrix], function(matrix) apply(matrix, 2, var, na.rm = TRUE))), 2, mean)
loading <- round(loading/sum(loading)*100, 2)
## Setting the labels
plot_args <- get.dots(dots, plot_args, "xlab", paste0("Dimension ", specific.args$dimensions[1], " (", loading[specific.args$dimensions[1]], "%)"))
plot_args <- get.dots(dots, plot_args, "ylab", paste0("Dimension ", specific.args$dimensions[2], " (", loading[specific.args$dimensions[2]], "%)"))
## Setting plot limits
plot_lim <- range(unlist(lapply(data$matrix[specific.args$matrix], function(matrix, dim) c(matrix[, dim]), dim = specific.args$dimensions)))
plot_args <- get.dots(dots, plot_args, "xlim", plot_lim)
plot_args <- get.dots(dots, plot_args, "ylim", plot_lim)
## Get the number of colour groups
n_groups <- length(data$subsets)
n_groups <- ifelse(n_groups == 0, 1, n_groups)
## Setting the colours
plot_args <- get.dots(dots, plot_args, "col",
default =
## Toggle between different defaults
if(n_groups == 1) {
"black"
} else {
if(data$call$subsets[[1]] %in% c("customised", "covar")) {
gg.color.hue(n_groups)
} else {
grDevices::heat.colors(n_groups+2)[1:n_groups]
}
}
, fun = "points")
## Setting the pch
plot_args <- get.dots(dots, plot_args, "pch", 19, fun = "points")
if(length(plot_args$pch) != n_groups) {
plot_args$pch <- rep(plot_args$pch, n_groups)
}
## Make a colour and pch classifier
col_order <- plot_args$col
pch_order <- plot_args$pch
if(n_groups > 1) {
## Get the list of subsets
subsets <- lapply(data$subsets,`[[`, "elements")
## Make an empty classifier
classifier <- rep(NA, nrow(data$matrix[[1]]))
for(class in match(sort(unlist(lapply(subsets, length)), decreasing = TRUE), unlist(lapply(subsets, length)))) {
## Store the selected subsets in the classifier (potentially overriding)
if(dim(subsets[[1]])[2] == 1) {
classifier[c(subsets[[class]])] <- class
} else {
classifier[c(elements.sampler(subsets[[class]]))] <- class
}
}
plot_args$col <- plot_args$col[classifier]
plot_args$pch <- plot_args$pch[classifier]
}
## Plot the empty plot
do.call(plot, plot_args)
## Plot the trees
if(plot_trees) {
## Plotting one edge
plot.edge <- function(one_edge, points_data, params) {
params$x <- points_data[one_edge, 1]
params$y <- points_data[one_edge, 2]
do.call(lines, params)
}
## Get the lines arguments
lines_args <- plot_args
lines_args <- get.dots(dots, lines_args, "col", "grey", fun = "lines")
if(length(specific.args$tree) > 1) {
lines_args$col <- make.transparent(lines_args$col, levels = length(specific.args$tree))
}
lines_args <- get.dots(dots, lines_args, "lwd", 1, fun = "lines")
lines_args <- get.dots(dots, lines_args, "lty", 1, fun = "lines")
## Plotting each tree
for(one_tree in specific.args$tree) {
## Select the right data
if(length(data$matrix) == length(specific.args$tree)) {
data_matrix <- data$matrix[[one_tree]]
} else {
data_matrix <- data$matrix[[1]]
}
## Selecting the origin points for the tree
points_data <- data_matrix[, specific.args$dimensions][c(data$tree[[one_tree]]$tip.label, data$tree[[one_tree]]$node.label), ]
## Plotting all the edges
silent <- apply(data$tree[[one_tree]]$edge, 1, plot.edge,
points_data = points_data,
params = lines_args)
}
}
## Plot the points
for(matrix in specific.args$matrix) {
point_args <- plot_args
## Transparentize the colour
if(length(specific.args$matrix) > 1) {
point_args$col <- make.transparent(point_args$col, levels = length(specific.args$matrix)*2)
}
## Add the points per matrix
point_args$x <- data$matrix[[specific.args$matrix[matrix]]][, specific.args$dimensions[1]]
point_args$y <- data$matrix[[specific.args$matrix[matrix]]][, specific.args$dimensions[2]]
## Call the points
do.call(points, point_args)
}
## Plot the legend
if(n_groups > 1 && plot_legend) {
## Set the legend arguments
legend_args <- list()
legend_args <- get.dots(dots, legend_args, "x", "topright", fun = "legend")
legend_args <- get.dots(dots, legend_args, "y", fun = "legend")
legend_args <- get.dots(dots, legend_args, "legend", names(data$subsets), fun = "legend")
legend_args <- get.dots(dots, legend_args, "col", col_order, fun = "legend")
legend_args <- get.dots(dots, legend_args, "pch", pch_order, fun = "legend")
legend_args <- get.dots(dots, legend_args, "cex", 2/3, fun = "legend")
legend_args <- get.dots(dots, legend_args, "bg", NULL, fun = "legend")
## Add the legend
do.call(legend, legend_args)
}
## Return invisible
return(invisible())
}
## The following is a modified version of plot.randtest from ade4 v1.4-3
plot.randtest <- function(plot_data) {
## Extracting the elements
dots <- plot_data$dots
data_sub <- plot_data$data_sub
## Plot the legend or not
if(is.null(dots$legend)) {
plot_legend <- TRUE
dots$legend <- NULL
} else {
if(is.logical(dots$legend)) {
plot_legend <-dots$legend
dots$legend <- NULL
}
}
## Add the histogram data
plot_args <- c(list(x = data_sub$plot$hist), dots)
## Removing specific args from dots
remove <- c(grep(c("legend"), names(plot_args)), grep(c("lines"), names(plot_args)), grep(c("points"), names(plot_args)))
plot_args[remove] <- NULL
## Plot arguments
plot_args <- get.dots(dots, plot_args, "xlim", data_sub$plot$xlim)
plot_args <- get.dots(dots, plot_args, "ylim", c(0, max(data_sub$plot$hist$count)))
plot_args <- get.dots(dots, plot_args, "col", "grey")
plot_args <- get.dots(dots, plot_args, "main")
## Plotting the simulated data
do.call(plot, plot_args)
## Observed data
observed <- data_sub$obs
## Adding the observed data
lines_args <- list(x = c(observed, observed), y = c(plot_args$ylim[2]/2, 0))
lines_args <- get.dots(dots, lines_args, "col", "black", fun = "lines")
lines_args <- get.dots(dots, lines_args, "lty", 1, fun = "lines")
lines_args <- get.dots(dots, lines_args, "lwd", 1, fun = "lines")
do.call(lines, lines_args)
point_args <- list(x = observed, y = plot_args$ylim[2]/2)
point_args <- get.dots(dots, point_args, "col", "black", fun = "points")
point_args <- get.dots(dots, point_args, "pch", 18, fun = "points")
point_args <- get.dots(dots, point_args, "cex", 2, fun = "points")
do.call(points, point_args)
## Adding the legend (test results)
if(plot_legend) {
## Set the legend arguments
legend_args <- list()
legend_args <- get.dots(dots, legend_args, "x", "topleft", fun = "legend")
legend_args <- get.dots(dots, legend_args, "y", fun = "legend")
legend_args <- get.dots(dots, legend_args, "bty", "n", fun = "legend")
legend_args <- get.dots(dots, legend_args, "legend", c("p-value", round(data_sub$pvalue, 5)), fun = "legend")
legend_args <- get.dots(dots, legend_args, "cex", 2/3, fun = "legend")
legend_args <- get.dots(dots, legend_args, "adj", 0.2, fun = "legend")
## Add the legend
do.call(legend, legend_args)
}
}
## The following is a modified version of dtt plots (from https://github.com/mwpennell/geiger-v2/blob/master/R/disparity.R)
plot.dtt <- function(data, quantiles, cent.tend, density, ...) {
dots <- list(...)
plot_args <- dots
remove <- c(grep(c("polygon"), names(plot_args)), grep(c("lines"), names(plot_args)))
plot_args[remove] <- NULL
## Set the default options
plot_args <- get.dots(plot_args, plot_args, "ylim",
if(!is.null(data$sim)) {
c(range(pretty(c(data$dtt, data$sim))))
} else {
c(range(pretty(data$dtt)))
}
)
plot_args <- get.dots(plot_args, plot_args, "xlab", "scaled time")
plot_args <- get.dots(plot_args, plot_args, "ylab", paste0("scaled ", paste(as.character(data$call[[3]]), collapse = " ")))
plot_args <- get.dots(plot_args, plot_args, "col", c("black", grDevices::colorRampPalette(c("lightgrey", "grey"))(length(quantiles))))
## Add the data for the plot args
plot_args$x <- data$times
plot_args$y <- data$dtt
plot_args$type <- "n"
## Empty plot
do.call(plot, plot_args)
## Add the simulated data
if(!is.null(data$sim)) {
## Check the quantiles
check.class(quantiles, "numeric", " must be any value between 1 and 100.")
## Are quantiles probabilities or proportions ?
if(any(quantiles < 1)) {
## Transform into proportion
quantiles <- quantiles*100
}
## Are quantiles proper proportions
if(any(quantiles < 0) | any(quantiles > 100)) {
stop.call("", "quantiles(s) must be any value between 0 and 100.")
}
n_quantiles <- length(quantiles)
## Check the central tendency
check.class(cent.tend, "function")
## The function must work
if(make.metric(cent.tend, silent = TRUE)$type != "level1") {
stop.call("", "cent.tend argument must be a function that outputs a single numeric value.")
}
## Summarised data
quantiles_values <- apply(data$sim, 1, quantile, probs = CI.converter(quantiles), na.rm = TRUE)
cent_tend_values <- apply(data$sim, 1, cent.tend)
## Plotting the polygons for each quantile
for(cis in 1:n_quantiles) {
poly_args <- list(x = c(data$times, rev(data$times)),
y = c(quantiles_values[(n_quantiles*2) - (cis-1), ], rev(quantiles_values[cis ,])))
poly_args <- get.dots(dots, poly_args, "border", FALSE, "polygon")
poly_args <- get.dots(dots, poly_args, "density", density, "polygon")
poly_args$col <- plot_args$col[cis+1]
do.call(polygon, poly_args)
}
## Add the central tendency
lines_args <- list(x = data$times,
y = cent_tend_values)
lines_args <- get.dots(dots, lines_args, "lty", 2, "lines")
lines_args <- get.dots(dots, lines_args, "lwd", 1, "lines")
lines_args <- get.dots(dots, lines_args, "col", plot_args$col[1], "lines")
do.call(lines, lines_args)
}
## Add the observed disparity
lines_args <- plot_args
lines_args$type <- NULL
lines_args <- get.dots(dots, lines_args, "lty", 1, "lines")
lines_args <- get.dots(dots, lines_args, "lwd", 1.5, "lines")
lines_args <- get.dots(dots, lines_args, "col", plot_args$col[1], "lines")
do.call(lines, lines_args)
}
## Plotting model tests results
plot.model.test <- function(data, ...) {
plot_args <- list(...)
## Set the default plotting arguments
plot_args <- get.dots(plot_args, plot_args, "ylab", "weighted AIC", "barplot")
plot_args <- get.dots(plot_args, plot_args, "col", "grey", "barplot")
## Extracting the weighted aicc
aic_values <- data$aic.models[, 3]
## Ordering the weighted aicc
plot_args$height <- aic_values[order(aic_values, decreasing = TRUE)]
## Plot
do.call(barplot, plot_args)
}
## Plotting the model simulation results
plot.model.sim <- function(data, add, density, quantiles, cent.tend, ...) {
dots <- list(...)
## Get inherited subsets (always from simulations to avoid NAs)
subset_names <- rev(data$simulation.data$fix$subsets)
## Set up the data parameters
data_params <- list("distribution" = TRUE,
"bootstrap" = TRUE,
"rarefaction" = NULL,
"between.groups" = FALSE,
"elements" = subset_names,
"model.sim" = TRUE)
## Model names
if(is.null(dots$main)) {
## Get the title
if(is(data$model, "character")) {
plot_main <- data$model
} else {
plot_main <- paste0(rownames(data$model)[1], " model\nAICc: ", data$model[1, "aicc"], "; log.lik: ", data$model[1, "log.lik"])
}
} else {
plot_main <- dots$main
}
## Get the plotting parameters
options(warn = -1)
plot_params <- get.plot.params(data, data_params, cent.tend, quantiles,
rarefaction_level = NULL,
type = "continuous",
elements_args = list(elements = FALSE),
observed_args = list(observed = FALSE),
xlab = ifelse(is.null(dots$xlab), "Time", dots$xlab),
ylab = ifelse(is.null(dots$ylab), "Disparity", dots$ylab),
main = plot_main)
options(warn = 0)
## Plot the simulated data
plot.continuous(plot_params, data_params, add = add, density = density)
}
## Plotting test metrics
plot.test.metric <- function(data, specific.args, ...) {
## Adding slopes
add.slope <- function(model, col) {
## Get the slope parameters
slope_param <- try.get.from.model(model, "Estimate")
## Plot the model
if(!any(is.na(slope_param)) || !is.null(slope_param)) {
## Add the slope
abline(a = slope_param[1],
b = slope_param[2],
col = col)
}
}
## Adding fits
add.fit <- function(model) {
fit_param <- try.get.from.model(model, "r.squared")
if(!is.null(fit_param) || length(fit_param) != 0) {
if(any(names(fit_param) == "adj.r.squared")) {
fit_param <- fit_param$adj.r.squared
is_adjusted <- TRUE
} else {
is_adjusted <- FALSE
}
return(paste0(ifelse(is_adjusted, "Adj. R^2: ", "R^2: "), unlist(round(fit_param, 3))))
} else {
return(NA)
}
}
## Detect whether to plot the shift steps or not
shift_plots <- !is.null(data$saved_steps)
if(shift_plots) {
## Find which steps to plot
if(!is.null(specific.args$visualise.steps)) {
check.class(specific.args$visualise.steps, c("numeric", "integer"))
steps_to_visualise <- specific.args$visualise.steps
if(any(check <- steps_to_visualise > n.subsets(data$saved_steps[[1]]))) {
stop(paste0("Impossible to display the step", ifelse(sum(check) > 1, "s ", " "), paste0(steps_to_visualise[check], collapse = ", "), " because the test only contains ", n.subsets(data$saved_steps[[1]]), " steps."), call. = FALSE)
}
} else {
steps_to_visualise <- floor(seq(from = 1, to = n.subsets(data$saved_steps[[1]]) - 1, length.out = 4))
}
n_steps <- length(steps_to_visualise)
}
## Setting plot window size
group_plot <- sapply(names(data$results), function(x) strsplit(x, "\\.")[[1]][[1]])
## Separating the plots in different groups (per plot windows)
n_plots <- length(unique(group_plot))
if(!shift_plots) {
if(n_plots > 1){
## Correct the number of plots if only 3
plot_size <- ifelse(n_plots == 3, 4, n_plots)
## Setting up plotting window
op_tmp <- par(mfrow = c(ceiling(sqrt(plot_size)),floor(sqrt(plot_size))))
}
} else {
## Create the template of step plots
## How many steps in each block?
steps_to_consider <- n_steps + n_steps %% 2
## Create the base matrix for plotting
base_matrix <- matrix(0, ceiling(sqrt(steps_to_consider)), floor(sqrt(steps_to_consider)), byrow = TRUE)
## Get the column for the plot results
results_plots <- matrix(rep(1:n_plots, each = steps_to_consider), ncol = ncol(base_matrix), byrow = TRUE)
## Create the vector of empty steps plots
steps_plots <- rep(c(base_matrix), n_plots)
## Fill the vector of step plots
for(i in 1:n_plots) {
select <- (1+(steps_to_consider*(i-1))):(steps_to_consider*i)
start_val <- ifelse(max(steps_plots) == 0, max(results_plots), max(steps_plots))
steps_plots[select[1:length(steps_to_visualise)]] <- (start_val+1):(start_val+length(steps_to_visualise))
}
## Create the layout
set_layout <- layout(cbind(results_plots, matrix(steps_plots, ncol = ncol(base_matrix), byrow = TRUE)))
# layout.show(set_layout)
## Parameter backup
op_tmp <- NULL
}
## Get the plotting arguments
plot_args <- list(...)
## Get the ylim
if(is.null(plot_args$ylim)) {
plot_args$ylim <- range(unlist(lapply(data$results, function(x) x[,2])), na.rm = TRUE)
}
## Get the colours
if(is.null(plot_args$col)) {
plot_args$col <- c("black", "orange", "blue")
} else {
if((missing_cols <- (3 - length(plot_args$col))) > 0) {
plot_args$col <- c(plot_args$col, rep(plot_args$col, missing_cols))
}
}
## get the title
if(is.null(plot_args$xlab)) {
plot_args$xlab <- "Amount of data considered (%)"
}
if(is.null(plot_args$ylab)) {
plot_args$ylab <- as.character(as.expression(data$call$metric))
}
## get the point size
if(is.null(plot_args$pch)) {
plot_args$pch <- 19
}
## Figuring out the legend
if(is.null(specific.args$legend)) {
plot_legend <- TRUE
legend_args <- list()
} else {
if(is(specific.args$legend, "logical")) {
plot_legend <- specific.args$legend
legend_args <- list()
} else {
plot_legend <- TRUE
legend_args <- specific.args$legend
}
}
## Checking out if it's two lists
if(length(legend_args == 2) && !is.null(names(legend_args[[1]])) && !is.null(names(legend_args[[2]]))) {
if(names(legend_args[[1]]) == names(legend_args[[2]])) {
legend_args_1 <- legend_args[[1]]
legend_args_2 <- legend_args_2_base <- legend_args[[2]]
}
} else {
legend_args_1 <- legend_args_2 <- legend_args_2_base <- legend_args
}
## Separating the data
plot_groups <- list()
for(group in 1:length(unique(group_plot))) {
plot_groups[[group]] <- data$results[grep(unique(group_plot)[[group]], names(data$results))]
}
## Separating the models
if(!is.null(data$models)) {
model_groups <- list()
for(group in 1:length(unique(group_plot))) {
model_groups[[group]] <- data$models[grep(unique(group_plot)[[group]], names(data$models))]
}
}
## Plot all the results
for(one_plot in 1:n_plots) {
## Get the data to plot
plot_data <- plot_groups[[one_plot]]
## Setting the specific plot arguments
one_plot_args <- plot_args
## get the xlim
if(is.null(plot_args$xlim)) {
one_plot_args$xlim <- range(as.numeric(plot_data[[1]][,1]))
}
## Get the title
if(is.null(plot_args$main)) {
one_plot_args$main <- unique(group_plot)[[one_plot]]
}
## Data to plot
one_plot_args$x <- plot_data[[1]][,1]
one_plot_args$y <- plot_data[[1]][,2]
## Which kind of plot
if(length(plot_data) == 1) {
## Set the colours
col_vector <- plot_args$col[1]
one_plot_args$col <- col_vector
do.call(plot, one_plot_args)
} else {
## Set the coulours
col_vector <- plot_args$col[-1]
## First part
one_plot_args$col <- col_vector[1]
do.call(plot, one_plot_args)
## Second part
one_plot_args$x <- plot_data[[2]][,1]
## Add a small shift
step <- abs(mode.val(diff(one_plot_args$x)))
one_plot_args$x <- one_plot_args$x + (step*0.05)
one_plot_args$y <- plot_data[[2]][,2]
one_plot_args$col <- col_vector[2]
do.call(points, one_plot_args)
## Set up the legend arguments
if(plot_legend) {
if(is.null(legend_args_1$x)) {
legend_args_1$x <- "bottomright"
}
if(is.null(legend_args_1$legend)) {
legend_args_1$legend <- names(plot_data)
}
if(is.null(legend_args_1$pch)) {
legend_args_1$pch <- plot_args$pch
}
if(is.null(legend_args_1$col)) {
legend_args_1$col <- col_vector
}
## Plot the legend
do.call(legend, legend_args_1)
}
}
## Adding the fit and models
if(!is.null(data$models)) {
## Get the fit of the first model
fit <- add.fit(model_groups[[one_plot]][[1]])
## Slope for the first model
add.slope(model_groups[[one_plot]][[1]], col = col_vector[1])
## Get the eventual second fit
if(!is.null(model_groups[[one_plot]][[1]])) {
## Fit for the second model
fit <- c(fit, add.fit(model_groups[[one_plot]][[2]]))
## Slope for the second model
add.slope(model_groups[[one_plot]][[2]], col = col_vector[2])
}
## Add the fits as a legend
if(!all(na_fit <- is.na(fit))) {
## Set up the legend arguments
if(plot_legend) {
if(is.null(legend_args_2$x)) {
legend_args_2$x <- "topright"
}
if(is.null(legend_args_2$legend)) {
legend_args_2$legend <- fit[!na_fit]
}
if(is.null(legend_args_2$lty)) {
legend_args_2$lty <- c(1,1)[!na_fit]
}
if(is.null(legend_args_2$col)) {
legend_args_2$col <- col_vector[!na_fit]
}
## Plot the legend
do.call(legend, legend_args_2)
## Reinitialise the legend
legend_args_2 <- legend_args_2_base
}
}
}
}
## Add the steps visualisation
if(shift_plots) {
## List of plot margins
mar_base <- c(2,2,2,1)
xaxts <- yaxts <- rep("n", n_steps)
yaxts[c(nrow(base_matrix), length(base_matrix))-1] <- "s"
xaxts[c(ncol(base_matrix)+1, length(base_matrix))] <- "s"
## Select the shifts
if(length(data$saved_steps)/2 != round(length(data$saved_steps)/2)) {
## There is the random one
shift_list <- c(1, seq(from = 2, to = length(data$saved_steps), by = 2))
} else {
## There is no random one
shift_list <- c(seq(from = 1, to = length(data$saved_steps), by = 2))
}
## Loop through the shifts
for(shift in shift_list) {
## Loop through the different steps
for(step in 1:n_steps){
## Selecting the subsets
step_preview <- get.subsets(data$saved_steps[[shift]],steps_to_visualise[step])
## Making a new subset without the negatives
step_preview$subsets$negatives$elements <- matrix((1:nrow(step_preview$matrix[[1]]))[-c(step_preview$subsets[[1]]$elements)], ncol = 1)
## Rename the subsets
names(step_preview$subsets)[2] <- as.character(100 - as.numeric(names(step_preview$subsets)[1]))
## Default title
step_name <- names(step_preview$subsets)[1]
## Selecting the colours and the title
if(names(data$saved_steps)[[shift]] == "random") {
select_col <- c(plot_args$col[1], "grey")
step_name <- paste0(step_name, "%")
} else {
select_col <- plot_args$col[-1]
step_name <- paste0(step_name, "% - ", as.character(100-as.numeric(step_name)), "%")
}
## Plotting the shift preview
par(mar = mar_base)
plot.dispRity(step_preview, type = "preview", legend = TRUE, col = select_col, main = step_name,xaxt = xaxts[step], yaxt = yaxts[step], legend.bg = "white")
}
}
## Reset the default plot margins
par(mar = c(5, 4, 4, 2) + 0.1)
}
## Restoring the parameters
if(n_plots > 1 && !is.null(op_tmp)) {
par(op_tmp)
}
}
## Plot axes
plot.axes <- function(data, ...) {
## Magic value for below
transparency <- 0.5
## Get the number of groups
n_groups <- length(data$dim.list)
## Get the columns to plot
if(is.null(data$call$colnames)) {
names.arg <- paste0("dim.", 1:length(data$scaled.var[[1]]))
} else {
names.arg <- data$call$colnames
}
## Plotting windows
op_tmp <- par(mfrow = c(ceiling(sqrt(n_groups)), round(sqrt(n_groups))))
## Getting the plotting parameters
plot_params <- list(...)
## Colour
if(is.null(plot_params$col)) {
plot_params$col <- rep("grey", n_groups)
} else {
if(length(plot_params$col) < n_groups) {
plot_params$col <- rep(plot_params$col, n_groups)
}
}
## Main
if(is.null(plot_params$main)) {
plot_params$main <- names(data$dim.list)
} else {
if(length(plot_params$main) < n_groups) {
plot_params$main <- rep(plot_params$main, n_groups)
}
}
## Ylab
if(is.null(plot_params$ylab)) {
plot_params$ylab <- "Scaled variance"
}
## Get the max number of dimensions
max_dim <- length(data$scaled.var[[1]])
## Xlim (selecting the maximum number of bars)
if(is.null(plot_params$xlim)) {
if(max(data$dimensions) != max_dim) {
## Only display the selected dimensions + 15%
display <- ceiling(max(data$dimensions)*1.1)
if(display <= max_dim) {
max_dim <- display
}
}
} else {
display <- max(plot_params$xlim)
if(display <= max_dim) {
max_dim <- display
}
plot_params$xlim <- NULL
}
## Plotting the bars
for(one_group in 1:n_groups) {
## Cum bars
bar_params <- plot_params
bar_params$col <- grDevices::adjustcolor(bar_params$col[one_group], alpha.f = transparency)
bar_params$main <- bar_params$main[one_group]
bar_params$height <- data$cumsum.var[[one_group]][1:max_dim]
bar_params$names.arg <- names.arg[1:max_dim]
bars_coords <- do.call(barplot, bar_params)
## Scaled bars
bar_params <- plot_params
bar_params$col <- bar_params$col[one_group]
bar_params$height <- data$scaled.var[[one_group]][1:max_dim]
bar_params$add <- TRUE
do.call(barplot, bar_params)
## Add the threshold line
abline(h = data$call$threshold, lty = 2)
abline(v = bars_coords[max(data$dim.list[[one_group]])], lty = 2)
}
par(op_tmp)
return(invisible())
}
## Plot projections
plot.projection <- function(data, specific.args, cent.tend, ...) {
dots <- list(...)
## Check whether to do a correlation plot or not
do_correlation_plot <- !is.null(specific.args$correlation.plot)
## Don't do a correlation plot
if(!do_correlation_plot) {
## Find the number of elements to plot
n_plots <- length(data)
## Open the plot window
par(mfrow = c(1,n_plots))
## Shmart plot
shmart.plot <- function(one_data, ylabs, ...) {
dots <- list(...)
if(is.null(dots$ylab)) {
plot.dispRity(one_data, ylab = ylabs, ...)
} else {
plot.dispRity(one_data, ...)
}
}
## Plot all that
if(is.null(dots$ylab)) {
ylab <- names(data)
} else {
ylab <- dots$ylab
dots$ylab <- NULL
}
silent <- mapply(shmart.plot, data, ylabs = ylab, MoreArgs = list(...))
} else {
## Get the central tendencies
if(!all(specific.args$correlation.plot %in% names(data)) && length(specific.args$correlation.plot) != 2) {
stop(paste0("correlation.plot argument must contain 2 elements from data (data contains: ", paste(names(data), collapse = ", "), ")."))
}
## Remove the phylogeny part (if exists)
if(n.subsets(data[[1]]) > 1) {
## Check if any subset is the size of the entire matrix
global <- which(size.subsets(data[[1]]) == nrow(data[[1]]$matrix[[1]]))
if(length(global) != 0) {
nsets <- 1:n.subsets(data[[1]])
data <- lapply(data, function(x, subsets) get.subsets(x, subsets = subsets), subsets = nsets[-global])
}
}
## Get the central tendencies
var1 <- unlist(lapply(lapply(data[[specific.args$correlation.plot[1]]]$disparity, function(X, type) return(X$elements)), function(X, fun) apply(X, 1, fun), fun = cent.tend))
var2 <- unlist(lapply(lapply(data[[specific.args$correlation.plot[2]]]$disparity, function(X, type) return(X$elements)), function(X, fun) apply(X, 1, fun), fun = cent.tend))
plot_data <- cbind(var1, var2)
colnames(plot_data) <- specific.args$correlation.plot
## Plot the results
plot(plot_data, ...)
}
}
Try the dispRity package in your browser
Any scripts or data that you put into this service are public.
dispRity documentation built on Aug. 9, 2022, 5:11 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Defines functions plot.projection plot.axes plot.test.metric plot.model.sim plot.model.test plot.dtt plot.randtest plot.preview plot.rarefaction plot.continuous plot.discrete plot.elements plot.observed get.quantile.col get.shift get.dots get.plot.params get.data.params
## Get the disparity characteristics
get.data.params <- function(data) {
return(list(
"distribution" = length(data$disparity[[1]]$elements) != 1,
"bootstrap" = !is.null(data$call$bootstrap) && data$call$bootstrap[[2]] != "covar",
"rarefaction" = data$call$bootstrap[[3]],
"between.groups" = data$call$disparity$metrics$between.groups,
"elements" = names(data$disparity)
))
}
## Setting up all the default plot parameters:
# $disparity: summarised_data
# $helpers: n_quantiles, npoints
# $options: ylim, ylab, xlab, col, ...
# $observed: data, ylim, ylab, xlab, col, ...
# $elements: elements, ylim, ylab, xlab, col, ...
get.plot.params <- function(data, data_params, cent.tend, quantiles, rarefaction_level, elements_args, type, observed_args, ...) {
## Set up the plotting data
## Summarise the data
if(is.null(data_params$model.sim)) {
# if((!is.na(data$call$subsets["trees"])) && (as.numeric(data$call$subsets["trees"]) > 1)) {
# summarised_data <- summary.dispRity(data, quantiles = quantiles, cent.tend = cent.tend, digits = 5, na.rm = TRUE)
# } else {
summarised_data <- summary.dispRity(data, quantiles = quantiles, cent.tend = cent.tend, digits = 5, na.rm = TRUE)
# }
} else {
summarised_data <- summary.dispRity(data, quantiles = quantiles, cent.tend = cent.tend, digits = 5)
## Making a fake obs column
colnames(summarised_data)[3] <- "obs"
}
## Find the observed data rows (that are not rarefactions)
observed_data <- !is.na(summarised_data[,(2+ifelse(data_params$between.groups, 2, 1))])
if(any(!observed_data)) {
## Find NAs in the last column
last_col_na <- is.na(summarised_data[, ncol(summarised_data)])
## If there are any NA in the last column, the observed NA is a true NA (nor rarefaction)
observed_data <- last_col_na | observed_data
}
## Get a specific rarefaction level
if(!is.null(rarefaction_level)) {
## Rarefaction level
if(data_params$between.groups) {
## Find the matching rarefaction level
rarefaction_match <- which(summarised_data$n_1 == rarefaction_level & summarised_data$n_2 == rarefaction_level)
} else {
rarefaction_match <- which(summarised_data$n == rarefaction_level)
}
## Get both the observed and rarefied data
selected_data <- summarised_data[rarefaction_match ,]
} else {
## Just get the observed data
selected_data <- summarised_data[observed_data,]
}
## Separate the data elements
disparity <- list()
name_part <- c(1:ifelse(data_params$between.groups, 3, 2))
disparity$names <- selected_data[, name_part, drop = FALSE]
disparity$data <- selected_data[, -name_part, drop = FALSE]
## Update the data to be in boxplot format
if(type == "box") {
## Extracting boxplot data
if(is.null(rarefaction_level)) {
## Getting the observed or bootstrapped data
if(data_params$bootstrap && !observed_args$observed) {
## Getting the bootstrapped data
box_data <- unlist(get.disparity(data, observed = FALSE), recursive = FALSE)
} else {
if(data_params$distribution || data_params$between.groups) {
box_data <- lapply(get.disparity(data, observed = TRUE, concatenate = FALSE), c)
} else {
if(ncol(disparity$data) == 1) {
box_data <- get.disparity(data, observed = TRUE)
} else {
box_data <- unlist(get.disparity(data, observed = FALSE), recursive = FALSE)
}
}
}
} else {
## Find the correct rarefaction level
box_data <- unlist(get.disparity(data, observed = FALSE, rarefaction = rarefaction_level), recursive = FALSE)
}
## Updating the disparity data part
if(data_params$between.groups) {
names(box_data) <- data_params$elements
}
disparity$data <- box_data
}
## Set up the helpers options
helpers <- list()
## Detect the number of quantiles
helpers$n_quantiles <- length(quantiles)
## Detect the number of points
if(is.null(data_params$model.sim)) {
helpers$n_points <- length(data$disparity)
} else {
## Selecting helpers for model.sim
helpers$n_points <- nrow(disparity$data)
}
## Set up the plotting options
dots <- list(...)
options <- list()
## Set the xlabel
if(is.null(dots$xlab)) {
## Default is subsets
options$xlab <- "Subsets"
## Default chrono.subset label
if(type == "continuous" && "continuous" %in% data$call$subsets) {
options$xlab <- "Time (Mya)"
}
## Default rarefaction label
if(type == "rarefaction") {
options$xlab <- "Elements"
}
} else {
## User input
check.class(dots$xlab, "character", " must be a character string.")
check.length(dots$xlab, 1, " must be a character string.")
options$xlab <- dots$xlab
dots$xlab <- NULL
}
## Set the ylabel
if(is.null(dots$ylab)) {
## Default is the metric name
options$ylab <- as.character(data$call$disparity$metrics$name)
} else {
## User input
check.class(dots$ylab, "character", " must be a character string.")
if(length(dots$ylab) > 2) stop.call("", "ylab can have maximum of two elements.")
options$ylab <- dots$ylab
dots$ylab <- NULL
}
## Set the y limits
if(is.null(dots$ylim)) {
if(type != "rarefaction") {
## Get the range of the data
options$ylim <- range(disparity$data, na.rm = TRUE)
## Add 2% on each side
percent_change <- 0.02
if(options$ylim[1] > 0) {
options$ylim[1] <- options$ylim[1] - options$ylim[1]*percent_change
} else {
options$ylim[1] <- options$ylim[1] + options$ylim[1]*percent_change
}
options$ylim[2] <- options$ylim[2] + options$ylim[2]*percent_change
} else {
options$ylim <- "rarefaction"
}
} else {
## User input
check.class(dots$ylim, c("numeric", "integer"))
check.length(dots$ylim, 2, " must be a vector of two elements.")
options$ylim <- dots$ylim
dots$ylim <- NULL
}
## Set the colours
if(is.null(dots$col)) {
## For boxplots, the default colour is white
if(type == "box" && type != "rarefaction") {
options$col <- "white"
} else {
options$col <- "black"
}
} else {
## User input
check.class(dots$col, "character", " must be a character string.")
options$col <- dots$col
dots$col <- NULL
}
## Add potential missing colours
if(helpers$n_quantiles > 0 && type != "box") {
## Count if they are any missing colours
cols_missing <- (helpers$n_quantiles + 1) - length(options$col)
## Adding grey scales if quantiles
if(cols_missing > 0) {
colfun <- grDevices::colorRampPalette(c("grey", "lightgrey"))
options$col <- c(options$col, colfun(cols_missing))
}
}
## Add additional options
if(length(dots) > 0) {
options <- c(options, dots)
}
## Observed data
if(observed_args$observed) {
## Adding the observed data
selected_data <- summarised_data[observed_data,]
observed_args$names <- selected_data[, name_part]
observed_args$data <- selected_data[, -name_part]
## Default observed arguments
if(is.null(observed_args$col)) {
if(type == "box" && options$col[[1]] == "white") {
observed_args$col <- "black"
} else {
observed_args$col <- options$col[[1]]
}
}
if(is.null(observed_args$pch)) {
observed_args$pch <- 4
}
if(is.null(observed_args$cex)) {
observed_args$cex <- 1
}
}
## Observed data
if(elements_args$elements) {
## Default observed arguments
if(is.null(elements_args$col)) {
if(type == "box" && options$col[[1]] == "white") {
elements_args$col <- c("black", "darkgrey")
} else {
elements_args$col <- options$col
}
}
if(length(elements_args$col) < 2) {
elements_args$col <- rep(elements_args$col, 2)
}
if(is.null(elements_args$pch)) {
elements_args$pch <- c(15, 15)
} else {
if(length(elements_args$pch) < 2) {
elements_args$pch <- rep(elements_args$pch, 2)
}
}
if(is.null(elements_args$lty)) {
elements_args$lty <- c(2, 2)
} else {
if(length(elements_args$lty) < 2) {
elements_args$lty <- rep(elements_args$lty, 2)
}
}
}
## Output the finalised list
return(list("disparity" = disparity, "helpers" = helpers, "options" = options, "observed_args" = observed_args, "elements_args" = elements_args))
}
## Get dots options
# @param dots ...
# @param args the list of arguments to be affected by dots
# @param name the argument name
# @param default an optional default argument (default is NULL)
# @param fun an optional suffix argument for the function to be used (default is NULL)
get.dots <- function(dots, args, name, default = NULL, fun = NULL) {
## Override the name without fun in dots
if(!is.null(fun) && !is.null(names(dots))) {
name_plus <- paste0(fun, ".", name)
if(name_plus %in% names(dots)) {
dots[[name]] <- dots[[name_plus]]
}
}
## Fill in the dots in the arguments list
if(!is.null(default)) {
if(is.null(dots[[name]])) {
## Set the default
args[[name]] <- default
} else {
## Use the dots
args[[name]] <- dots[[name]]
}
} else {
## Use the dots
if(!is.null(dots[[name]])) {
args[[name]] <- dots[[name]]
}
}
return(args)
}
## Handle the shift argument
get.shift <- function(add, plot_params) {
## Set the shift parameter (for add)
shift = 0
if(add) {
## Is the previous plot the same size?
prev_axis <- par("xaxp")
if(prev_axis[2] == plot_params$helpers$n_points) {
shift = 0
} else {
shift = 0.5
}
}
return(shift)
}
## Getting the columns for the right quantiles
get.quantile.col <- function(cent_tend_col, cis, n_quantiles) {
return(c(
## Lower quantile
cent_tend_col + cis,
## Higher quantile
(cent_tend_col + n_quantiles*2) - (cis-1)
))
}
## Observed points
plot.observed <- function(plot_params) {
if(plot_params$observed_args$observed) {
## Set the observed arguments
points_args <- plot_params$observed_args
points_args$x <- 1:plot_params$helpers$n_points
points_args$y <- points_args$data[, "obs"]
points_args$observed <- NULL
points_args$names <- NULL
points_args$data <- NULL
## Calling the points
do.call(points, points_args)
} else {
return(invisible())
}
}
## Plot elements
plot.elements <- function(plot_params, data_params, type) {
## Elements plots
add_args <- plot_params$elements_args
add_args$elements <- NULL
## Getting the scaling value
percent_change <- 0.02
rescaling <- plot_params$options$ylim
if(rescaling[1] >= 0) {
rescaling[1] <- rescaling[1] + percent_change * rescaling[1]
} else {
rescaling[1] <- rescaling[1] - percent_change * rescaling[1]
}
rescaling[2] <- rescaling[2] - percent_change * rescaling[2]
## Rescale the y values
rescaled_elements <- plot_params$disparity$names[,-1, drop = FALSE]
rescaled_elements <- matrix(scales::rescale(unlist(rescaled_elements), to = rescaling), ncol = ncol(rescaled_elements), byrow = FALSE)
## Add the values
add_args$x <- 1:plot_params$helpers$n_points
add_args$y <- rescaled_elements[,1]
## Select the first colour
add_args$col <- plot_params$elements_args$col[1]
## Check if ylab exists
if(is.null(plot_params$elements_args$ylab)) {
axis_label <- ifelse(length(plot_params$options$ylab) > 1, plot_params$options$ylab[2], "Elements")
} else {
axis_label <- plot_params$elements_args$ylab
}
## Add the different types
if(type == "continuous") {
add_args$lty <- plot_params$elements_args$lty[1]
do.call(lines, add_args)
} else {
add_args$pch <- plot_params$elements_args$pch[1]
do.call(points, add_args)
}
## Add the ones for the second group (if two groups)
if(data_params$between.groups) {
## Get the n_2 column
add_args$y <- rescaled_elements[,2]
## Select the second colour
add_args$col <- plot_params$elements_args$col[2]
## Add the different types
if(type == "continuous") {
add_args$lty <- plot_params$elements_args$lty[2]
do.call(lines, add_args)
} else {
## Slightly shift
add_args$x <- add_args$x + (mode.val(diff(add_args$x))*0.01)
add_args$pch <- plot_params$elements_args$pch[2]
do.call(points, add_args)
}
}
## Add the axis label
seq.along.range <- function(range, out = 5) {
seq(from = range(range)[1],
to = range(range)[2],
length.out = out)
}
## Add the second y axis
axis(4, at = seq.along.range(rescaling), labels = seq.along.range(plot_params$disparity$names[,-1]),
lty = plot_params$elements_args$lty)
## Add the second axis label
mtext(axis_label, side = 4, line = 2)
}
## discrete plotting
plot.discrete <- function(plot_params, data_params, add, density, type) {
## Get the shifting argument
shift <- get.shift(add, plot_params)
## Select the central tendency column
cent_tend_col <- ifelse(data_params$bootstrap, 2, 1)
## Creating the matrix frame (for neat plotting)
frame_matrix <- matrix(1:plot_params$helpers$n_points,
ncol = plot_params$helpers$n_point,
dimnames = list(c(), plot_params$disparity$names[,1]))
if(!add) {
## Set the frame boxplot arguments
box_args <- plot_params$options
box_args$x <- frame_matrix
box_args$col <- box_args$border <- "white"
box_args$ylab <- box_args$ylab[[1]]
box_args$boxwex <- 0.00001
box_args$type <- "n"
## Plot the boxplot
do.call(boxplot, box_args)
}
## Check if bootstrapped
if(data_params$bootstrap || data_params$distribution) {
## Set the width (default)
width <- 0.5
## Initialise the plot arguments
plot_args <- plot_params$options
## Add the quantiles
if(type == "polygon") {
## Add the polygon options
plot_args$border <- plot_params$options$col[[1]]
plot_args$density <- density
## Loop through each point
for(point in 1:plot_params$helpers$n_points) {
## Loop through each quantile
for(cis in 1:plot_params$helpers$n_quantiles) {
## Add the colour option
plot_args$col <- plot_params$options$col[[cis+1]]
## Set the bars width
point_width <- width/(plot_params$helpers$n_quantiles - cis + 1.5)
## Set the x coordinates
plot_args$x <- c(point - point_width,
point + point_width,
point + point_width,
point - point_width) + shift
## Select the quantiles columns
quantiles_col <- get.quantile.col(cent_tend_col, cis, plot_params$helpers$n_quantiles)
## Set the y coordinates
plot_args$y <- c(rep(plot_params$disparity$data[point, quantiles_col[1]], 2),
rep(plot_params$disparity$data[point, quantiles_col[2]], 2))
## Add the polygon
do.call(polygon, plot_args)
}
}
}
if(type == "line") {
## Loop through each point
for(point in 1:plot_params$helpers$n_points) {
## Loop through the quantiles
for(cis in 1:plot_params$helpers$n_quantiles) {
## The the line type
plot_args$lty <- plot_params$helpers$n_quantiles - cis + 1
plot_args$lwd <- 1 + cis * 1.5
plot_args$col <- plot_params$options$col[[cis+1]]
## Set the x values
plot_args$x <- rep(point, 2)
## Select the quantiles columns
quantiles_col <- get.quantile.col(cent_tend_col, cis, plot_params$helpers$n_quantiles)
## Set the y value
plot_args$y <- c(plot_params$disparity$data[point, quantiles_col[1]],
plot_params$disparity$data[point, quantiles_col[2]])
## Plotting the line
do.call(lines, plot_args)
}
}
}
}
## Add the points estimates
point_args <- list()
## Get the points coordinates
point_args$x <- 1:plot_params$helpers$n_points + shift
point_args$y <- plot_params$disparity$data[, cent_tend_col]
## Get the options
point_args$col <- plot_params$options$col[[1]]
if(is.null(plot_params$options$pch)) {
point_args$pch <- 19
} else {
point_args$pch <- plot_params$options$pch
}
## Add the points
do.call(points, point_args)
## Save parameters
return(par())
}
## continuous plotting
plot.continuous <- function(plot_params, data_params, add, density) {
## Get the shifting argument
shift <- get.shift(add, plot_params)
## Select the central tendency column
cent_tend_col <- ifelse(data_params$bootstrap, 2, 1)
## Set up the cur plot options
plot_args <- plot_params$options
plot_args$x <- seq(from = 1, to = plot_params$helpers$n_points)-shift
plot_args$y <- plot_params$disparity$data[, cent_tend_col]
plot_args$type <- "l"
plot_args$col <- plot_params$options$col[[1]]
plot_args$xaxt <- "n"
## Plot the central tendency
if(!add) {
## Plot the thingy
do.call(plot, plot_args)
## Add the axis labels
options(warn = -1)
try_round <- as.numeric(plot_params$disparity$names[,"subsets"])
options(warn = 0)
if(all(is.na(try_round))) {
## Keep the ticks as they are
x_ticks <- plot_params$disparity$names[,"subsets"]
} else {
rounding <- 3-max(nchar(round(try_round)))
x_ticks <- round(try_round, digits = ifelse(rounding < 0, 0, rounding))
}
axis(1, 1:plot_params$helpers$n_points, x_ticks)
} else {
## Plot the line
do.call(lines, plot_args)
}
## Check if bootstrapped
if(data_params$bootstrap || data_params$distribution) {
## Setting the plot parameters
poly_args <- plot_params$options
poly_args$border <- "NA"
poly_args$density <- density
poly_args$x <- numeric()
poly_args$y <- numeric()
## Add the polygons
for (cis in 1:plot_params$helpers$n_quantiles) {
## Set the x values
poly_args$x <- c(1:plot_params$helpers$n_points)
poly_args$x <- c(poly_args$x, rev(poly_args$x))
## Select the quantiles columns
quantiles_col <- get.quantile.col(cent_tend_col, cis, plot_params$helpers$n_quantiles)
## Set the y values
poly_args$y <- plot_params$disparity$data[, quantiles_col[1]]
poly_args$y <- c(poly_args$y, rev(plot_params$disparity$data[, quantiles_col[2]]))
## Set up the colour
poly_args$col <- plot_params$options$col[cis+1]
## Dividing the polygon if NAs
if(any(is_nas <- is.na(poly_args$y[1:plot_params$helpers$n_points]))) {
## Selecting the groups of applicable data
groups <- numeric()
group_label <- 1
## (attributing a group as soon as it jumps an NA)
for(point in seq_along(is_nas)) {
if(is_nas[point]) {
## Increment the group
group_label <- group_label + 1
groups[point] <- NA
} else {
groups[point] <- group_label
}
}
## splitting the data into the groups
split.combine.data <- function(var, groups, n_points) {
vals_1 <- split(var[1:n_points], groups)
vals_2 <- split(var[-c(1:n_points)], rev(groups))
return(mapply(c, vals_1, vals_2, SIMPLIFY = FALSE))
}
y_vals <- split.combine.data(poly_args$y, groups, plot_params$helpers$n_points)
x_vals <- split.combine.data(poly_args$x, groups, plot_params$helpers$n_points)
## Combine args
combined.poly.args <- function(x, y, poly_args) {
poly_args$x <- x
poly_args$y <- y
return(poly_args)
}
list_poly_args <- mapply(combined.poly.args, x_vals, y_vals, SIMPLIFY = FALSE, MoreArgs = list("poly_args" = poly_args))
## Plot all the polygons
lapply(list_poly_args, function(args) do.call(polygon, args))
} else {
## Plot the polygon
do.call(polygon, poly_args)
}
}
## Add the central tendency
do.call(lines, plot_args)
}
## Save parameters
return(par())
}
## Plot the rarefaction
plot.rarefaction <- function(plot_params, data_params, data) {
## How many rarefaction plots?
n_plots <- length(data$subsets)
## Open the multiple plots
plot_size <- ifelse(n_plots == 3, 4, n_plots)
op_tmp <- par(mfrow = c(ceiling(sqrt(plot_size)),round(sqrt(plot_size))))
## Get the list of subsets
subsets_levels <- names(data$subsets)
## Get all the rarefaction values
rarefied_data <- summary(data)
## Get were the central tendency value is
cent_tend_col <- ifelse(data_params$between.groups, 5, 4)
## Setting the plotting arguments
all_plot_args <- plot_params$options
all_plot_args$ylim <- NULL
all_plot_args$lty <- 1
all_plot_args$type <- "l"
## Plot the different curves
for(one_subset in subsets_levels) {
## get the subset data
subset_data <- rarefied_data[rarefied_data[,1] == one_subset, ]
## Setting the plotting args for the specific subset
one_plot_args <- all_plot_args
one_plot_args$x <- subset_data[,2]
one_plot_args$y <- subset_data[,cent_tend_col]
one_plot_args$ylim <- range(subset_data[, -c(1:(cent_tend_col-1))])
one_plot_args$main <- one_subset
## Plot the central tendency
do.call(plot, one_plot_args)
## Add the quantiles
for(cis in 1:plot_params$helpers$n_quantiles) {
## Get the quantile columns
ci_cols <- get.quantile.col(cent_tend_col, cis, plot_params$helpers$n_quantiles)
## Set the plotting arguments
lines_args <- one_plot_args
lines_args$col <- all_plot_args$col[length(all_plot_args$col) - (cis-1)]
lines_args$y <- subset_data[,ci_cols[1]]
## plot the quantiles
do.call(lines, lines_args)
lines_args$y <- subset_data[,ci_cols[2]]
do.call(lines, lines_args)
}
}
## Done!
par(op_tmp)
}
## Plotting a space preview
plot.preview <- function(data, specific.args, ...) {
dots <- list(...)
## The "ggplot" colours
gg.color.hue <- function(n) {
grDevices::hcl(h = seq(15, 375, length = n + 1), l = 65, c = 100)[1:n]
}
make.transparent <- function(colour, levels) {
grDevices::adjustcolor(colour, alpha.f = (1/levels) + 1/3)
}
## Set up the specific args
specific.args <- get.dots(specific.args, specific.args, "matrix", c(1:length(data$matrix)))
specific.args <- get.dots(specific.args, specific.args, "dimensions", c(1:2))
## Plot the legend?
if(is.null(dots$legend)) {
plot_legend <- TRUE
dots$legend <- NULL
} else {
if(is.logical(dots$legend)) {
plot_legend <-dots$legend
dots$legend <- NULL
}
}
if(is.null(specific.args$tree) || is.null(data$tree[[1]])) {
## Don't plot trees
plot_trees <- FALSE
} else {
if(is(specific.args$tree, "logical")) {
## Toggle plot_trees
plot_trees <- specific.args$tree
if(plot_trees) {
## Select the trees to plot
specific.args$tree <- 1:length(data$tree)
}
} else {
## Specific args must be the specific trees to plot
plot_trees <- TRUE
}
}
## Capturing the dots options
plot_args <- c(list(x = NULL, y = NULL), dots)
## Removing specific args from dots
remove <- c(grep(c("legend"), names(plot_args)), grep(c("lines"), names(plot_args)), grep(c("points"), names(plot_args)))
plot_args[remove] <- NULL
## Getting the loadings
loading <- apply(do.call(rbind, lapply(data$matrix[specific.args$matrix], function(matrix) apply(matrix, 2, var, na.rm = TRUE))), 2, mean)
loading <- round(loading/sum(loading)*100, 2)
## Setting the labels
plot_args <- get.dots(dots, plot_args, "xlab", paste0("Dimension ", specific.args$dimensions[1], " (", loading[specific.args$dimensions[1]], "%)"))
plot_args <- get.dots(dots, plot_args, "ylab", paste0("Dimension ", specific.args$dimensions[2], " (", loading[specific.args$dimensions[2]], "%)"))
## Setting plot limits
plot_lim <- range(unlist(lapply(data$matrix[specific.args$matrix], function(matrix, dim) c(matrix[, dim]), dim = specific.args$dimensions)))
plot_args <- get.dots(dots, plot_args, "xlim", plot_lim)
plot_args <- get.dots(dots, plot_args, "ylim", plot_lim)
## Get the number of colour groups
n_groups <- length(data$subsets)
n_groups <- ifelse(n_groups == 0, 1, n_groups)
## Setting the colours
plot_args <- get.dots(dots, plot_args, "col",
default =
## Toggle between different defaults
if(n_groups == 1) {
"black"
} else {
if(data$call$subsets[[1]] %in% c("customised", "covar")) {
gg.color.hue(n_groups)
} else {
grDevices::heat.colors(n_groups+2)[1:n_groups]
}
}
, fun = "points")
## Setting the pch
plot_args <- get.dots(dots, plot_args, "pch", 19, fun = "points")
if(length(plot_args$pch) != n_groups) {
plot_args$pch <- rep(plot_args$pch, n_groups)
}
## Make a colour and pch classifier
col_order <- plot_args$col
pch_order <- plot_args$pch
if(n_groups > 1) {
## Get the list of subsets
subsets <- lapply(data$subsets,`[[`, "elements")
## Make an empty classifier
classifier <- rep(NA, nrow(data$matrix[[1]]))
for(class in match(sort(unlist(lapply(subsets, length)), decreasing = TRUE), unlist(lapply(subsets, length)))) {
## Store the selected subsets in the classifier (potentially overriding)
if(dim(subsets[[1]])[2] == 1) {
classifier[c(subsets[[class]])] <- class
} else {
classifier[c(elements.sampler(subsets[[class]]))] <- class
}
}
plot_args$col <- plot_args$col[classifier]
plot_args$pch <- plot_args$pch[classifier]
}
## Plot the empty plot
do.call(plot, plot_args)
## Plot the trees
if(plot_trees) {
## Plotting one edge
plot.edge <- function(one_edge, points_data, params) {
params$x <- points_data[one_edge, 1]
params$y <- points_data[one_edge, 2]
do.call(lines, params)
}
## Get the lines arguments
lines_args <- plot_args
lines_args <- get.dots(dots, lines_args, "col", "grey", fun = "lines")
if(length(specific.args$tree) > 1) {
lines_args$col <- make.transparent(lines_args$col, levels = length(specific.args$tree))
}
lines_args <- get.dots(dots, lines_args, "lwd", 1, fun = "lines")
lines_args <- get.dots(dots, lines_args, "lty", 1, fun = "lines")
## Plotting each tree
for(one_tree in specific.args$tree) {
## Select the right data
if(length(data$matrix) == length(specific.args$tree)) {
data_matrix <- data$matrix[[one_tree]]
} else {
data_matrix <- data$matrix[[1]]
}
## Selecting the origin points for the tree
points_data <- data_matrix[, specific.args$dimensions][c(data$tree[[one_tree]]$tip.label, data$tree[[one_tree]]$node.label), ]
## Plotting all the edges
silent <- apply(data$tree[[one_tree]]$edge, 1, plot.edge,
points_data = points_data,
params = lines_args)
}
}
## Plot the points
for(matrix in specific.args$matrix) {
point_args <- plot_args
## Transparentize the colour
if(length(specific.args$matrix) > 1) {
point_args$col <- make.transparent(point_args$col, levels = length(specific.args$matrix)*2)
}
## Add the points per matrix
point_args$x <- data$matrix[[specific.args$matrix[matrix]]][, specific.args$dimensions[1]]
point_args$y <- data$matrix[[specific.args$matrix[matrix]]][, specific.args$dimensions[2]]
## Call the points
do.call(points, point_args)
}
## Plot the legend
if(n_groups > 1 && plot_legend) {
## Set the legend arguments
legend_args <- list()
legend_args <- get.dots(dots, legend_args, "x", "topright", fun = "legend")
legend_args <- get.dots(dots, legend_args, "y", fun = "legend")
legend_args <- get.dots(dots, legend_args, "legend", names(data$subsets), fun = "legend")
legend_args <- get.dots(dots, legend_args, "col", col_order, fun = "legend")
legend_args <- get.dots(dots, legend_args, "pch", pch_order, fun = "legend")
legend_args <- get.dots(dots, legend_args, "cex", 2/3, fun = "legend")
legend_args <- get.dots(dots, legend_args, "bg", NULL, fun = "legend")
## Add the legend
do.call(legend, legend_args)
}
## Return invisible
return(invisible())
}
## The following is a modified version of plot.randtest from ade4 v1.4-3
plot.randtest <- function(plot_data) {
## Extracting the elements
dots <- plot_data$dots
data_sub <- plot_data$data_sub
## Plot the legend or not
if(is.null(dots$legend)) {
plot_legend <- TRUE
dots$legend <- NULL
} else {
if(is.logical(dots$legend)) {
plot_legend <-dots$legend
dots$legend <- NULL
}
}
## Add the histogram data
plot_args <- c(list(x = data_sub$plot$hist), dots)
## Removing specific args from dots
remove <- c(grep(c("legend"), names(plot_args)), grep(c("lines"), names(plot_args)), grep(c("points"), names(plot_args)))
plot_args[remove] <- NULL
## Plot arguments
plot_args <- get.dots(dots, plot_args, "xlim", data_sub$plot$xlim)
plot_args <- get.dots(dots, plot_args, "ylim", c(0, max(data_sub$plot$hist$count)))
plot_args <- get.dots(dots, plot_args, "col", "grey")
plot_args <- get.dots(dots, plot_args, "main")
## Plotting the simulated data
do.call(plot, plot_args)
## Observed data
observed <- data_sub$obs
## Adding the observed data
lines_args <- list(x = c(observed, observed), y = c(plot_args$ylim[2]/2, 0))
lines_args <- get.dots(dots, lines_args, "col", "black", fun = "lines")
lines_args <- get.dots(dots, lines_args, "lty", 1, fun = "lines")
lines_args <- get.dots(dots, lines_args, "lwd", 1, fun = "lines")
do.call(lines, lines_args)
point_args <- list(x = observed, y = plot_args$ylim[2]/2)
point_args <- get.dots(dots, point_args, "col", "black", fun = "points")
point_args <- get.dots(dots, point_args, "pch", 18, fun = "points")
point_args <- get.dots(dots, point_args, "cex", 2, fun = "points")
do.call(points, point_args)
## Adding the legend (test results)
if(plot_legend) {
## Set the legend arguments
legend_args <- list()
legend_args <- get.dots(dots, legend_args, "x", "topleft", fun = "legend")
legend_args <- get.dots(dots, legend_args, "y", fun = "legend")
legend_args <- get.dots(dots, legend_args, "bty", "n", fun = "legend")
legend_args <- get.dots(dots, legend_args, "legend", c("p-value", round(data_sub$pvalue, 5)), fun = "legend")
legend_args <- get.dots(dots, legend_args, "cex", 2/3, fun = "legend")
legend_args <- get.dots(dots, legend_args, "adj", 0.2, fun = "legend")
## Add the legend
do.call(legend, legend_args)
}
}
## The following is a modified version of dtt plots (from https://github.com/mwpennell/geiger-v2/blob/master/R/disparity.R)
plot.dtt <- function(data, quantiles, cent.tend, density, ...) {
dots <- list(...)
plot_args <- dots
remove <- c(grep(c("polygon"), names(plot_args)), grep(c("lines"), names(plot_args)))
plot_args[remove] <- NULL
## Set the default options
plot_args <- get.dots(plot_args, plot_args, "ylim",
if(!is.null(data$sim)) {
c(range(pretty(c(data$dtt, data$sim))))
} else {
c(range(pretty(data$dtt)))
}
)
plot_args <- get.dots(plot_args, plot_args, "xlab", "scaled time")
plot_args <- get.dots(plot_args, plot_args, "ylab", paste0("scaled ", paste(as.character(data$call[[3]]), collapse = " ")))
plot_args <- get.dots(plot_args, plot_args, "col", c("black", grDevices::colorRampPalette(c("lightgrey", "grey"))(length(quantiles))))
## Add the data for the plot args
plot_args$x <- data$times
plot_args$y <- data$dtt
plot_args$type <- "n"
## Empty plot
do.call(plot, plot_args)
## Add the simulated data
if(!is.null(data$sim)) {
## Check the quantiles
check.class(quantiles, "numeric", " must be any value between 1 and 100.")
## Are quantiles probabilities or proportions ?
if(any(quantiles < 1)) {
## Transform into proportion
quantiles <- quantiles*100
}
## Are quantiles proper proportions
if(any(quantiles < 0) | any(quantiles > 100)) {
stop.call("", "quantiles(s) must be any value between 0 and 100.")
}
n_quantiles <- length(quantiles)
## Check the central tendency
check.class(cent.tend, "function")
## The function must work
if(make.metric(cent.tend, silent = TRUE)$type != "level1") {
stop.call("", "cent.tend argument must be a function that outputs a single numeric value.")
}
## Summarised data
quantiles_values <- apply(data$sim, 1, quantile, probs = CI.converter(quantiles), na.rm = TRUE)
cent_tend_values <- apply(data$sim, 1, cent.tend)
## Plotting the polygons for each quantile
for(cis in 1:n_quantiles) {
poly_args <- list(x = c(data$times, rev(data$times)),
y = c(quantiles_values[(n_quantiles*2) - (cis-1), ], rev(quantiles_values[cis ,])))
poly_args <- get.dots(dots, poly_args, "border", FALSE, "polygon")
poly_args <- get.dots(dots, poly_args, "density", density, "polygon")
poly_args$col <- plot_args$col[cis+1]
do.call(polygon, poly_args)
}
## Add the central tendency
lines_args <- list(x = data$times,
y = cent_tend_values)
lines_args <- get.dots(dots, lines_args, "lty", 2, "lines")
lines_args <- get.dots(dots, lines_args, "lwd", 1, "lines")
lines_args <- get.dots(dots, lines_args, "col", plot_args$col[1], "lines")
do.call(lines, lines_args)
}
## Add the observed disparity
lines_args <- plot_args
lines_args$type <- NULL
lines_args <- get.dots(dots, lines_args, "lty", 1, "lines")
lines_args <- get.dots(dots, lines_args, "lwd", 1.5, "lines")
lines_args <- get.dots(dots, lines_args, "col", plot_args$col[1], "lines")
do.call(lines, lines_args)
}
## Plotting model tests results
plot.model.test <- function(data, ...) {
plot_args <- list(...)
## Set the default plotting arguments
plot_args <- get.dots(plot_args, plot_args, "ylab", "weighted AIC", "barplot")
plot_args <- get.dots(plot_args, plot_args, "col", "grey", "barplot")
## Extracting the weighted aicc
aic_values <- data$aic.models[, 3]
## Ordering the weighted aicc
plot_args$height <- aic_values[order(aic_values, decreasing = TRUE)]
## Plot
do.call(barplot, plot_args)
}
## Plotting the model simulation results
plot.model.sim <- function(data, add, density, quantiles, cent.tend, ...) {
dots <- list(...)
## Get inherited subsets (always from simulations to avoid NAs)
subset_names <- rev(data$simulation.data$fix$subsets)
## Set up the data parameters
data_params <- list("distribution" = TRUE,
"bootstrap" = TRUE,
"rarefaction" = NULL,
"between.groups" = FALSE,
"elements" = subset_names,
"model.sim" = TRUE)
## Model names
if(is.null(dots$main)) {
## Get the title
if(is(data$model, "character")) {
plot_main <- data$model
} else {
plot_main <- paste0(rownames(data$model)[1], " model\nAICc: ", data$model[1, "aicc"], "; log.lik: ", data$model[1, "log.lik"])
}
} else {
plot_main <- dots$main
}
## Get the plotting parameters
options(warn = -1)
plot_params <- get.plot.params(data, data_params, cent.tend, quantiles,
rarefaction_level = NULL,
type = "continuous",
elements_args = list(elements = FALSE),
observed_args = list(observed = FALSE),
xlab = ifelse(is.null(dots$xlab), "Time", dots$xlab),
ylab = ifelse(is.null(dots$ylab), "Disparity", dots$ylab),
main = plot_main)
options(warn = 0)
## Plot the simulated data
plot.continuous(plot_params, data_params, add = add, density = density)
}
## Plotting test metrics
plot.test.metric <- function(data, specific.args, ...) {
## Adding slopes
add.slope <- function(model, col) {
## Get the slope parameters
slope_param <- try.get.from.model(model, "Estimate")
## Plot the model
if(!any(is.na(slope_param)) || !is.null(slope_param)) {
## Add the slope
abline(a = slope_param[1],
b = slope_param[2],
col = col)
}
}
## Adding fits
add.fit <- function(model) {
fit_param <- try.get.from.model(model, "r.squared")
if(!is.null(fit_param) || length(fit_param) != 0) {
if(any(names(fit_param) == "adj.r.squared")) {
fit_param <- fit_param$adj.r.squared
is_adjusted <- TRUE
} else {
is_adjusted <- FALSE
}
return(paste0(ifelse(is_adjusted, "Adj. R^2: ", "R^2: "), unlist(round(fit_param, 3))))
} else {
return(NA)
}
}
## Detect whether to plot the shift steps or not
shift_plots <- !is.null(data$saved_steps)
if(shift_plots) {
## Find which steps to plot
if(!is.null(specific.args$visualise.steps)) {
check.class(specific.args$visualise.steps, c("numeric", "integer"))
steps_to_visualise <- specific.args$visualise.steps
if(any(check <- steps_to_visualise > n.subsets(data$saved_steps[[1]]))) {
stop(paste0("Impossible to display the step", ifelse(sum(check) > 1, "s ", " "), paste0(steps_to_visualise[check], collapse = ", "), " because the test only contains ", n.subsets(data$saved_steps[[1]]), " steps."), call. = FALSE)
}
} else {
steps_to_visualise <- floor(seq(from = 1, to = n.subsets(data$saved_steps[[1]]) - 1, length.out = 4))
}
n_steps <- length(steps_to_visualise)
}
## Setting plot window size
group_plot <- sapply(names(data$results), function(x) strsplit(x, "\\.")[[1]][[1]])
## Separating the plots in different groups (per plot windows)
n_plots <- length(unique(group_plot))
if(!shift_plots) {
if(n_plots > 1){
## Correct the number of plots if only 3
plot_size <- ifelse(n_plots == 3, 4, n_plots)
## Setting up plotting window
op_tmp <- par(mfrow = c(ceiling(sqrt(plot_size)),floor(sqrt(plot_size))))
}
} else {
## Create the template of step plots
## How many steps in each block?
steps_to_consider <- n_steps + n_steps %% 2
## Create the base matrix for plotting
base_matrix <- matrix(0, ceiling(sqrt(steps_to_consider)), floor(sqrt(steps_to_consider)), byrow = TRUE)
## Get the column for the plot results
results_plots <- matrix(rep(1:n_plots, each = steps_to_consider), ncol = ncol(base_matrix), byrow = TRUE)
## Create the vector of empty steps plots
steps_plots <- rep(c(base_matrix), n_plots)
## Fill the vector of step plots
for(i in 1:n_plots) {
select <- (1+(steps_to_consider*(i-1))):(steps_to_consider*i)
start_val <- ifelse(max(steps_plots) == 0, max(results_plots), max(steps_plots))
steps_plots[select[1:length(steps_to_visualise)]] <- (start_val+1):(start_val+length(steps_to_visualise))
}
## Create the layout
set_layout <- layout(cbind(results_plots, matrix(steps_plots, ncol = ncol(base_matrix), byrow = TRUE)))
# layout.show(set_layout)
## Parameter backup
op_tmp <- NULL
}
## Get the plotting arguments
plot_args <- list(...)
## Get the ylim
if(is.null(plot_args$ylim)) {
plot_args$ylim <- range(unlist(lapply(data$results, function(x) x[,2])), na.rm = TRUE)
}
## Get the colours
if(is.null(plot_args$col)) {
plot_args$col <- c("black", "orange", "blue")
} else {
if((missing_cols <- (3 - length(plot_args$col))) > 0) {
plot_args$col <- c(plot_args$col, rep(plot_args$col, missing_cols))
}
}
## get the title
if(is.null(plot_args$xlab)) {
plot_args$xlab <- "Amount of data considered (%)"
}
if(is.null(plot_args$ylab)) {
plot_args$ylab <- as.character(as.expression(data$call$metric))
}
## get the point size
if(is.null(plot_args$pch)) {
plot_args$pch <- 19
}
## Figuring out the legend
if(is.null(specific.args$legend)) {
plot_legend <- TRUE
legend_args <- list()
} else {
if(is(specific.args$legend, "logical")) {
plot_legend <- specific.args$legend
legend_args <- list()
} else {
plot_legend <- TRUE
legend_args <- specific.args$legend
}
}
## Checking out if it's two lists
if(length(legend_args == 2) && !is.null(names(legend_args[[1]])) && !is.null(names(legend_args[[2]]))) {
if(names(legend_args[[1]]) == names(legend_args[[2]])) {
legend_args_1 <- legend_args[[1]]
legend_args_2 <- legend_args_2_base <- legend_args[[2]]
}
} else {
legend_args_1 <- legend_args_2 <- legend_args_2_base <- legend_args
}
## Separating the data
plot_groups <- list()
for(group in 1:length(unique(group_plot))) {
plot_groups[[group]] <- data$results[grep(unique(group_plot)[[group]], names(data$results))]
}
## Separating the models
if(!is.null(data$models)) {
model_groups <- list()
for(group in 1:length(unique(group_plot))) {
model_groups[[group]] <- data$models[grep(unique(group_plot)[[group]], names(data$models))]
}
}
## Plot all the results
for(one_plot in 1:n_plots) {
## Get the data to plot
plot_data <- plot_groups[[one_plot]]
## Setting the specific plot arguments
one_plot_args <- plot_args
## get the xlim
if(is.null(plot_args$xlim)) {
one_plot_args$xlim <- range(as.numeric(plot_data[[1]][,1]))
}
## Get the title
if(is.null(plot_args$main)) {
one_plot_args$main <- unique(group_plot)[[one_plot]]
}
## Data to plot
one_plot_args$x <- plot_data[[1]][,1]
one_plot_args$y <- plot_data[[1]][,2]
## Which kind of plot
if(length(plot_data) == 1) {
## Set the colours
col_vector <- plot_args$col[1]
one_plot_args$col <- col_vector
do.call(plot, one_plot_args)
} else {
## Set the coulours
col_vector <- plot_args$col[-1]
## First part
one_plot_args$col <- col_vector[1]
do.call(plot, one_plot_args)
## Second part
one_plot_args$x <- plot_data[[2]][,1]
## Add a small shift
step <- abs(mode.val(diff(one_plot_args$x)))
one_plot_args$x <- one_plot_args$x + (step*0.05)
one_plot_args$y <- plot_data[[2]][,2]
one_plot_args$col <- col_vector[2]
do.call(points, one_plot_args)
## Set up the legend arguments
if(plot_legend) {
if(is.null(legend_args_1$x)) {
legend_args_1$x <- "bottomright"
}
if(is.null(legend_args_1$legend)) {
legend_args_1$legend <- names(plot_data)
}
if(is.null(legend_args_1$pch)) {
legend_args_1$pch <- plot_args$pch
}
if(is.null(legend_args_1$col)) {
legend_args_1$col <- col_vector
}
## Plot the legend
do.call(legend, legend_args_1)
}
}
## Adding the fit and models
if(!is.null(data$models)) {
## Get the fit of the first model
fit <- add.fit(model_groups[[one_plot]][[1]])
## Slope for the first model
add.slope(model_groups[[one_plot]][[1]], col = col_vector[1])
## Get the eventual second fit
if(!is.null(model_groups[[one_plot]][[1]])) {
## Fit for the second model
fit <- c(fit, add.fit(model_groups[[one_plot]][[2]]))
## Slope for the second model
add.slope(model_groups[[one_plot]][[2]], col = col_vector[2])
}
## Add the fits as a legend
if(!all(na_fit <- is.na(fit))) {
## Set up the legend arguments
if(plot_legend) {
if(is.null(legend_args_2$x)) {
legend_args_2$x <- "topright"
}
if(is.null(legend_args_2$legend)) {
legend_args_2$legend <- fit[!na_fit]
}
if(is.null(legend_args_2$lty)) {
legend_args_2$lty <- c(1,1)[!na_fit]
}
if(is.null(legend_args_2$col)) {
legend_args_2$col <- col_vector[!na_fit]
}
## Plot the legend
do.call(legend, legend_args_2)
## Reinitialise the legend
legend_args_2 <- legend_args_2_base
}
}
}
}
## Add the steps visualisation
if(shift_plots) {
## List of plot margins
mar_base <- c(2,2,2,1)
xaxts <- yaxts <- rep("n", n_steps)
yaxts[c(nrow(base_matrix), length(base_matrix))-1] <- "s"
xaxts[c(ncol(base_matrix)+1, length(base_matrix))] <- "s"
## Select the shifts
if(length(data$saved_steps)/2 != round(length(data$saved_steps)/2)) {
## There is the random one
shift_list <- c(1, seq(from = 2, to = length(data$saved_steps), by = 2))
} else {
## There is no random one
shift_list <- c(seq(from = 1, to = length(data$saved_steps), by = 2))
}
## Loop through the shifts
for(shift in shift_list) {
## Loop through the different steps
for(step in 1:n_steps){
## Selecting the subsets
step_preview <- get.subsets(data$saved_steps[[shift]],steps_to_visualise[step])
## Making a new subset without the negatives
step_preview$subsets$negatives$elements <- matrix((1:nrow(step_preview$matrix[[1]]))[-c(step_preview$subsets[[1]]$elements)], ncol = 1)
## Rename the subsets
names(step_preview$subsets)[2] <- as.character(100 - as.numeric(names(step_preview$subsets)[1]))
## Default title
step_name <- names(step_preview$subsets)[1]
## Selecting the colours and the title
if(names(data$saved_steps)[[shift]] == "random") {
select_col <- c(plot_args$col[1], "grey")
step_name <- paste0(step_name, "%")
} else {
select_col <- plot_args$col[-1]
step_name <- paste0(step_name, "% - ", as.character(100-as.numeric(step_name)), "%")
}
## Plotting the shift preview
par(mar = mar_base)
plot.dispRity(step_preview, type = "preview", legend = TRUE, col = select_col, main = step_name,xaxt = xaxts[step], yaxt = yaxts[step], legend.bg = "white")
}
}
## Reset the default plot margins
par(mar = c(5, 4, 4, 2) + 0.1)
}
## Restoring the parameters
if(n_plots > 1 && !is.null(op_tmp)) {
par(op_tmp)
}
}
## Plot axes
plot.axes <- function(data, ...) {
## Magic value for below
transparency <- 0.5
## Get the number of groups
n_groups <- length(data$dim.list)
## Get the columns to plot
if(is.null(data$call$colnames)) {
names.arg <- paste0("dim.", 1:length(data$scaled.var[[1]]))
} else {
names.arg <- data$call$colnames
}
## Plotting windows
op_tmp <- par(mfrow = c(ceiling(sqrt(n_groups)), round(sqrt(n_groups))))
## Getting the plotting parameters
plot_params <- list(...)
## Colour
if(is.null(plot_params$col)) {
plot_params$col <- rep("grey", n_groups)
} else {
if(length(plot_params$col) < n_groups) {
plot_params$col <- rep(plot_params$col, n_groups)
}
}
## Main
if(is.null(plot_params$main)) {
plot_params$main <- names(data$dim.list)
} else {
if(length(plot_params$main) < n_groups) {
plot_params$main <- rep(plot_params$main, n_groups)
}
}
## Ylab
if(is.null(plot_params$ylab)) {
plot_params$ylab <- "Scaled variance"
}
## Get the max number of dimensions
max_dim <- length(data$scaled.var[[1]])
## Xlim (selecting the maximum number of bars)
if(is.null(plot_params$xlim)) {
if(max(data$dimensions) != max_dim) {
## Only display the selected dimensions + 15%
display <- ceiling(max(data$dimensions)*1.1)
if(display <= max_dim) {
max_dim <- display
}
}
} else {
display <- max(plot_params$xlim)
if(display <= max_dim) {
max_dim <- display
}
plot_params$xlim <- NULL
}
## Plotting the bars
for(one_group in 1:n_groups) {
## Cum bars
bar_params <- plot_params
bar_params$col <- grDevices::adjustcolor(bar_params$col[one_group], alpha.f = transparency)
bar_params$main <- bar_params$main[one_group]
bar_params$height <- data$cumsum.var[[one_group]][1:max_dim]
bar_params$names.arg <- names.arg[1:max_dim]
bars_coords <- do.call(barplot, bar_params)
## Scaled bars
bar_params <- plot_params
bar_params$col <- bar_params$col[one_group]
bar_params$height <- data$scaled.var[[one_group]][1:max_dim]
bar_params$add <- TRUE
do.call(barplot, bar_params)
## Add the threshold line
abline(h = data$call$threshold, lty = 2)
abline(v = bars_coords[max(data$dim.list[[one_group]])], lty = 2)
}
par(op_tmp)
return(invisible())
}
## Plot projections
plot.projection <- function(data, specific.args, cent.tend, ...) {
dots <- list(...)
## Check whether to do a correlation plot or not
do_correlation_plot <- !is.null(specific.args$correlation.plot)
## Don't do a correlation plot
if(!do_correlation_plot) {
## Find the number of elements to plot
n_plots <- length(data)
## Open the plot window
par(mfrow = c(1,n_plots))
## Shmart plot
shmart.plot <- function(one_data, ylabs, ...) {
dots <- list(...)
if(is.null(dots$ylab)) {
plot.dispRity(one_data, ylab = ylabs, ...)
} else {
plot.dispRity(one_data, ...)
}
}
## Plot all that
if(is.null(dots$ylab)) {
ylab <- names(data)
} else {
ylab <- dots$ylab
dots$ylab <- NULL
}
silent <- mapply(shmart.plot, data, ylabs = ylab, MoreArgs = list(...))
} else {
## Get the central tendencies
if(!all(specific.args$correlation.plot %in% names(data)) && length(specific.args$correlation.plot) != 2) {
stop(paste0("correlation.plot argument must contain 2 elements from data (data contains: ", paste(names(data), collapse = ", "), ")."))
}
## Remove the phylogeny part (if exists)
if(n.subsets(data[[1]]) > 1) {
## Check if any subset is the size of the entire matrix
global <- which(size.subsets(data[[1]]) == nrow(data[[1]]$matrix[[1]]))
if(length(global) != 0) {
nsets <- 1:n.subsets(data[[1]])
data <- lapply(data, function(x, subsets) get.subsets(x, subsets = subsets), subsets = nsets[-global])
}
}
## Get the central tendencies
var1 <- unlist(lapply(lapply(data[[specific.args$correlation.plot[1]]]$disparity, function(X, type) return(X$elements)), function(X, fun) apply(X, 1, fun), fun = cent.tend))
var2 <- unlist(lapply(lapply(data[[specific.args$correlation.plot[2]]]$disparity, function(X, type) return(X$elements)), function(X, fun) apply(X, 1, fun), fun = cent.tend))
plot_data <- cbind(var1, var2)
colnames(plot_data) <- specific.args$correlation.plot
## Plot the results
plot(plot_data, ...)
}
}
Try the dispRity package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.