pcvr (version 1.2.0)

pcv.emd: Earth Mover's Distance between spectral histograms

Description

pcv.emd can be used to calculate Earth Mover's Distance between pairwise histograms in a wide dataframe of multi value traits. The is expected to be used with output from mv_ag. See also pcv.euc for euclidean distance between histograms.

Usage

pcv.emd(
  df,
  cols = NULL,
  reorder = NULL,
  include = reorder,
  mat = FALSE,
  plot = TRUE,
  parallel = getOption("mc.cores", 1),
  trait = "trait",
  id = "image",
  value = "value",
  raiseError = TRUE,
  method = "emd"
)
pcv.euc(
  df,
  cols = NULL,
  reorder = NULL,
  include = reorder,
  mat = FALSE,
  plot = TRUE,
  parallel = getOption("mc.cores", 1),
  trait = "trait",
  id = "image",
  value = "value",
  raiseError = TRUE,
  method = "euc"
)

Value

A dataframe/matrix (if plot=FALSE) or a list with a dataframe/matrix and\ a ggplot (if plot=TRUE). The returned data contains pairwise EMD values.

Arguments

df

Data frame to use with multi value traits in wide format or long format

cols

Columns to use. Defaults to NULL in which case all columns are used. Single strings will be used to regex a pattern in column names (see examples). A vector of names, positions, or booleans will also work. For long data this is taken as a regex pattern (or full name) to use in filtering the trait column.

reorder

Should data be reordered to put similar rows together in the resulting plot? This takes a vector of column names of length 1 or more (see examples).

include

if a long dataframe is returned then these columns will be added to the dataframe, labelled for i and j (the row positions for compared histograms). If a matrix is returned then this information is stored in the row names. This defaults to reorder.

mat

Logical, should data be returned as an nrow x nrow matrix or as a long dataframe? By Default this is FALSE and a long dataframe is returned. Both options are comparable in terms of speed, although for large datasets the matrix version may be slightly faster.

plot

Logical, should a plot be returned? For a matrix this is made with heatmap(), for a dataframe this uses ggplot.

parallel

Number of cores to use. Defaults to 1 unless the "mc.cores" option is set.

trait

Column name for long data to identify traits. This defaults to "trait". If this and value are in the column names of the data then it is assumed to be in long format, otherwise it is assumed to be in wide format.

id

A vector of column names that uniquely identifies observations if the data is in long format. Defaults to "image".

value

A column name for the values to be drawn from in long data. Defaults to "value".

raiseError

Logical, should warnings/errors be raised for potentially large output? It is easy to ask for very many comparisons with this function so the goal of this argument is to catch a few of those and give estimates of how much time something may take. If the function is expected to take very long then a warning or an error is raised. If this is set to FALSE then no time estimates are made.

method

Which method to use (one of "emd" or "euc"). Defaults to "emd".

Examples

set.seed(123)
test <- mvSim(
  dists = list(
    runif = list(min = 0, max = 100),
    rnorm = list(mean = 90, sd = 20)
  ),
  n_samples = 10
)
test$meta1 <- rep(LETTERS[1:3], length.out = nrow(test))
test$meta2 <- rep(LETTERS[4:5], length.out = nrow(test))

x <- pcv.emd(
  df = test, cols = "sim", reorder = "group",
  include = c("meta1", "meta2"), mat = FALSE,
  plot = FALSE, parallel = 1
)
head(x)
x2 <- pcv.emd(
  df = test, cols = "sim", reorder = "group",
  include = c("meta1", "meta2"), mat = FALSE,
  plot = FALSE, parallel = 1, method = "euc"
)
head(x2)

# \donttest{
tryCatch(
  {
    library(data.table)
    file <- paste0(
      "https://media.githubusercontent.com/media/joshqsumner/",
      "pcvrTestData/main/pcv4-multi-value-traits.csv"
    )
    df1 <- read.pcv(file, "wide", reader = "fread")

    df1$genotype <- substr(df1$barcode, 3, 5)
    df1$genotype <- ifelse(df1$genotype == "002", "B73",
      ifelse(df1$genotype == "003", "W605S",
        ifelse(df1$genotype == "004", "MM", "Mo17")
      )
    )
    df1$fertilizer <- substr(df1$barcode, 8, 8)
    df1$fertilizer <- ifelse(df1$fertilizer == "A", "100",
      ifelse(df1$fertilizer == "B", "50", "0")
    )

    w <- pcv.emd(df1,
      cols = "hue_frequencies", reorder = c("fertilizer", "genotype"),
      mat = FALSE, plot = TRUE, parallel = 1
    )
  },
  error = function(err) {
    message(err)
  }
)

# Note on computational complexity
# This scales as O^2, see the plot below for some idea
# of the time for different input data sizes.
emdTime <- function(x, n = 1) {
  x^2 / n * 0.0023
}
plot(
  x = c(18, 36, 54, 72, 108, 135), y = c(0.74, 2.89, 6.86, 10.99, 26.25, 42.44),
  xlab = "N Input Images", ylab = "time (seconds)"
) # benchmarked test data
lines(x = 1:150, y = emdTime(1:150)) # exponential function

plot(
  x = 1:1000, y = emdTime(1:1000), type = "l",
  xlab = "N Input Images", ylab = "time (seconds)"
)
# }