Package 'OmicFlow'

Description

Visualizes group differences using a half-boxplot, half-jitter plot with pairwise statistical testing and p-value annotation. This function is built into the class omics with method alpha_diversity().

Usage

boxjitter_test(
  data,
  values,
  col_name = lifecycle::deprecated(),
  groups,
  group_by = lifecycle::deprecated(),
  split_by = NULL,
  palette = NULL,
  method = lifecycle::deprecated(),
  test = "wilcox",
  paired = FALSE,
  p.adjust.method = "fdr",
  step.increase = 0.05
)

Arguments

Value

A ggplot2 object to be further modified

Examples

n_row <- 1000
n_col <- 100
density <- 0.2
num_entries <- n_row * n_col
num_nonzero <- round(num_entries * density)

set.seed(123)
positions <- sample(num_entries, num_nonzero, replace=FALSE)
row_idx <- ((positions - 1) %% n_row) + 1
col_idx <- ((positions - 1) %/% n_row) + 1

values <- runif(num_nonzero, min = 0, max = 1)
sparse_mat <- Matrix::sparseMatrix(
   i = row_idx,
   j = col_idx,
   x = values,
   dims = c(n_row, n_col)
 )

div <- OmicFlow::diversity(
  x = sparse_mat,
  metric = "shannon"
)

dt <- data.table::data.table(
  "shannon" = div,
  "treatment" = c(rep("healthy", n_col / 2), rep("tumor", n_col / 2)),
  "sex" = c(rep("male", n_col / 4), rep("female", n_col / 4))
)

colors <- OmicFlow::colormap(dt, "treatment")

# Comparing two groups
plt <- OmicFlow::boxjitter_test(
 data = dt,
 values = "shannon",
 groups = "treatment",
 palette = colors,
 test = "wilcox",
 paired = FALSE,
 p.adjust.method = "fdr"
)

# Performing a test while stratifying the plot in two groups
plt <- OmicFlow::boxjitter_test(
 data = dt,
 values = "shannon",
 groups = "treatment",
 split_by = "sex",
 palette = colors,
 test = "wilcox",
 paired = FALSE,
 p.adjust.method = "fdr"
)

Description

Calculates the Bray-Curtis dissimilarity of a Matrix pairwise for each column.

Usage

bray(x, weighted = TRUE, threads = 1)

Arguments

Details

The Bray-Curtis dissimilarity between two samples $A$ and $B$, each of length $n$, is defined as:

$d(A,B) = \frac{\sum_{i}^n |A_i - B_i|}{\sum_{i}^n (A_i + B_i)}$

where $A_i$ and $B_i$ are the abundances of the $i$-th feature in sample $A$ and $B$, respectively. When weighted is set to FALSE, counts are replaced by presence/absence data.

Value

A column x column dist object.

References

Bray, J.R. & Curtis, J.T. (1957) An Ordination of the Upland Forest Communities of Southern Wisconsin. Ecological Monographs, 27(4), 325–349.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
    metaData = metadata_file,
    countData = counts_file,
    featureData = features_file,
    treeData = tree_file
)

taxa$feature_subset(Kingdom == "Bacteria")
taxa$scale(method = "tss")

bray(taxa$countData)

Description

Calculates the Canberra dissimilarity of a Matrix pairwise for each column.

Usage

canberra(x, weighted = TRUE, threads = 1)

Arguments

Details

The Canberra dissimilarity between two samples $A$ and $B$, each of length $n$, is defined as:

$d(A,B) = \frac{1 / NZ} \sum_{i}^n \frac{|A_i - B_i|}{|A_i| + |B_i|}$

where $A_i$ and $B_i$ are the abundances of the $i$-th feature in sample $A$ and $B$, respectively. NZ are the number of non-zero entries. When weighted is set to FALSE, counts are replaced by presence/absence data.

Value

A column x column dist object.

References

Lance, G.N. & Williams, W.T. (1967) Mixed-data classificatory programs. I. Agglomerative systems. Australian Computer Journal, 1(1), 15-20.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
    metaData = metadata_file,
    countData = counts_file,
    featureData = features_file,
    treeData = tree_file
)

taxa$feature_subset(Kingdom == "Bacteria")
taxa$scale(method = "tss")

canberra(taxa$countData)

Description

Creates a color map of setNames class with names as unique group labels given a column name and values as hexcode colors. This function is built into the composition method from the abstract class omics and inherited by other omics classes, such as; metagenomics and proteomics.

Usage

colormap(
  data,
  col_name = lifecycle::deprecated(),
  groups,
  Brewer.palID = "Set2"
)

Arguments

Value

A setNames.

Examples

library("data.table")
dt <- data.table(
  "SAMPLE_ID" = c("sample_1", "sample_2", "sample_3"),
  "treatment" = c("healthy", "tumor", NA)
)

colors <- colormap(data = dt,
                   groups = "treatment")

Description

Mainly used within omics and other functions to check if given column name(s) exist in the table and is not completely empty (containing NAs).

Usage

column_exists(column, table)

Arguments

Value

A boolean value.

Description

Creates a stacked barchart of features. It is possible to both show barcharts for each sample or group them by a categorical variable. The function is compatible with the class omics method composition().

Usage

composition_plot(
  data,
  palette,
  feature_rank,
  title_name = NULL,
  group_by = NULL
)

Arguments

Value

A ggplot2 object to be further modified

Examples

# Create mock_data as data.frame (data.table is also supported)
mock_data <- data.frame(
  SAMPLE_ID = rep(paste0("Sample", 1:10), each = 5),
  Genus = rep(c("GenusA","GenusB","GenusC","GenusD","GenusE"), times = 10),
  value = c(
    0.1119, 0.1303, 0.0680, 0.5833, 0.1065,      # Sample1
    0.2080, 0.1179, 0.0211, 0.4578, 0.1951,      # Sample2
    0.4219, 0.1189, 0.2320, 0.1037, 0.1235,      # Sample3
    0.4026, 0.0898, 0.1703, 0.1063, 0.2309,      # Sample4
    0.1211, 0.0478, 0.5721, 0.1973, 0.0618,      # Sample5
    0.2355, 0.0293, 0.2304, 0.1520, 0.3528,      # Sample6
    0.2904, 0.0347, 0.3651, 0.0555, 0.2544,      # Sample7
    0.4138, 0.0299, 0.0223, 0.4996, 0.0345,      # Sample8
    0.4088, 0.0573, 0.0155, 0.2888, 0.2296,      # Sample9
    0.4941, 0.0722, 0.2331, 0.1023, 0.0983       # Sample10
  ),
  Group = rep(c("Group1","Group2","Group1",
                "Group1","Group2","Group2",
                 "Group1","Group1","Group1","Group2"),
               each = 5)
)

# Create a colormap
mock_palette <- c(
  GenusA = "#1f77b4",  # blue
  GenusB = "#ff7f0e",  # orange
  GenusC = "#2ca02c",  # green
  GenusD = "#d62728",  # red
  GenusE = "#9467bd"   # purple
)

# Optionally: Use OmicFlow::colormap()
mock_palette <- colormap(
  data = mock_data,
  col_name = "Genus",
  Brewer.palID = "RdYlBu"
)

composition_plot(
  data = mock_data,
  palette = mock_palette,
  feature_rank = "Genus",
  title_name = "Mock Genus Composition"
)

composition_plot(
  data = mock_data,
  palette = mock_palette,
  feature_rank = "Genus",
  title_name = "Mock Genus Composition by Group",
  group_by = "Group"
)

Description

Calculates the cosine disimilarity of a Marix pairwise for each column.

Usage

cosine(x, weighted = TRUE, threads = 1)

Arguments

Details

The cosine dissimilarity between two samples $A$ and $B$, each of length $n$, is defined as:

$d(A,B) = 1 - \frac{\sum_{i}^n A_i B_i}{\sqrt{\sum_{i}^n A_i^2} \sqrt{\sum_{i}^n B_i^2}}$

where $A_i$ and $B_i$ are the abundances of the $i$-th feature in sample $A$ and $B$, respectively. When weighted is set to FALSE, counts are replaced by presence/absence data.

Value

A column x column dist object.

References

Deza, M. M., & Deza, E. (2009). Encyclopedia of Distances. Springer Science & Business Media., 308.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
    metaData = metadata_file,
    countData = counts_file,
    featureData = features_file,
    treeData = tree_file
)

taxa$feature_subset(Kingdom == "Bacteria")
taxa$scale(method = "tss")

cosine(taxa$countData)

Description

Computes the alpha diversity based on Shannon index, simpson or invsimpson. Code is adapted from diversity and uses sparseMatrix in triplet format over the dense matrix. The code is much faster and memory efficient, while still being mathematically correct. This function is built into the class omics with method alpha_diversity() and inherited by other omics classes, such as; metagenomics and proteomics.

Usage

diversity(x, metric = "shannon", normalize = TRUE, base = exp(1))

Arguments

Value

A numeric vector with type double.

See Also

diversity

Examples

n_row <- 1000
n_col <- 100
density <- 0.2
num_entries <- n_row * n_col
num_nonzero <- round(num_entries * density)

set.seed(123)
positions <- sample(num_entries, num_nonzero, replace=FALSE)
row_idx <- ((positions - 1) %% n_row) + 1
col_idx <- ((positions - 1) %/% n_row) + 1

values <- runif(num_nonzero, min = 0, max = 1)
sparse_mat <- Matrix::sparseMatrix(
  i = row_idx,
  j = col_idx,
  x = values,
  dims = c(n_row, n_col)
)

# Alpha diversity is computed on column level
## Transpose the sparseMatrix if required with `t()` from Matrix R package.
result <- OmicFlow::diversity(
  x = sparse_mat,
  metric = "shannon"
)

Description

Calculates the Euclidean dissimilarity of a Matrix pairwise for each column.

Usage

euclidean(x, weighted = TRUE, threads = 1)

Arguments

Details

The Euclidean dissimilarity between two samples $A$ and $B$, each of length $n$, is defined as:

$d(A,B) = \sqrt{\sum_{i=1}^n (A_i - B_i)^2}$

where $A_i$ and $B_i$ are the abundances of the $i$-th feature in sample $A$ and $B$, respectively. When weighted is set to FALSE, counts are replaced by presence/absence data.

Value

A column x column dist object.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
    metaData = metadata_file,
    countData = counts_file,
    featureData = features_file,
    treeData = tree_file
)

taxa$feature_subset(Kingdom == "Bacteria")
taxa$scale(method = "tss")

euclidean(taxa$countData)

Description

Computes the hill numbers for q is 0, 1 or 2. Code is adapted from hill_taxa and uses sparseMatrix in triplet format over the dense matrix. The code is much faster and memory efficient, while still being mathematical correct.

Usage

hill_taxa(x, q = 0, normalize = TRUE, base = exp(1))

Arguments

Value

A numeric vector with type double.

See Also

hill_taxa

Examples

n_row <- 1000
n_col <- 100
density <- 0.2
num_entries <- n_row * n_col
num_nonzero <- round(num_entries * density)

set.seed(123)
positions <- sample(num_entries, num_nonzero, replace=FALSE)
row_idx <- ((positions - 1) %% n_row) + 1
col_idx <- ((positions - 1) %/% n_row) + 1

values <- runif(num_nonzero, min = 0, max = 1)
sparse_mat <- Matrix::sparseMatrix(
  i = row_idx,
  j = col_idx,
  x = values,
  dims = c(n_row, n_col)
)

result <- OmicFlow::hill_taxa(
 x = sparse_mat,
 q = 2
)

Description

Calculates the Jaccard dissimilarity of a Matrix pairwise for each column.

Usage

jaccard(x, weighted = TRUE, threads = 1)

Arguments

Details

The weighted Jaccard disimilarity between two samples $A$ and $B$, each of length $n$, is defined as:

$d(A,B) = 1 - \frac{ \sum_{i}^{n} \min(A_i, B_i) }{ \sum_{i}^{n} \max(A_i, B_i) }$

where $A_i$ and $B_i$ are the abundances of the $i$-th feature in sample $A$ and $B$, respectively. When weighted is set to FALSE, abundances are changed to 1 (classical Jaccard for binary data).

Value

A column x column dist object.

References

Jaccard, P. (1912) The distribution of the flora in the alpine zone. New Phytologist, 11(2), 37–50. library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow") counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow") features_file <- system.file("extdata", "features.tsv", package = "OmicFlow") tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new( metaData = metadata_file, countData = counts_file, featureData = features_file, treeData = tree_file )

taxa$feature_subset(Kingdom == "Bacteria") taxa$scale(method = "tss")

jaccard(taxa$countData)

Description

Calculates the Jensen-Shannon divergence of a Matrix pairwise for each column.

Usage

jsd(x, weighted = TRUE, threads = 1)

Arguments

Details

The Jensen-Shannon divergence between two probability distributions $A$ and $B$, each of length $n$, is defined as:

$d(A, B) = \frac{1}{2} D_{KL}(A \parallel M) + \frac{1}{2} D_{KL}(B \parallel M)$

where $M = \frac{1}{2} (A + B)$ is the mixture distribution, and $D_{KL}$ is the Kullback-Leibler divergence. When weighted is set to FALSE, counts are replaced by presence/absence data.

Value

A column x column dist object.

References

Lin, J. (1991). Divergence measures based on the Shannon entropy. IEEE Transactions on Information Theory, 37(1), 145-151.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
    metaData = metadata_file,
    countData = counts_file,
    featureData = features_file,
    treeData = tree_file
)

taxa$feature_subset(Kingdom == "Bacteria")
taxa$scale(method = "tss")

jsd(taxa$countData)

Description

Calculates the Manhattan dissimilarity of a Matrix pairwise for each column.

Usage

manhattan(x, weighted = TRUE, threads = 1)

Arguments

Details

The Manhattan dissimilarity between two samples $A$ and $B$, each of length $n$, is defined as:

$d(A, B) = \sum_{i}^n |A_i - B_i|$

where $A_i$ and $B_i$ are the abundances of the $i$-th feature in sample $A$ and $B$, respectively. When weighted is set to FALSE, counts are replaced by presence/absence data.

Value

A column x column dist object.

References

Deza, M. M., & Deza, E. (2009). Encyclopedia of Distances. Springer Science & Business Media., 313.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
    metaData = metadata_file,
    countData = counts_file,
    featureData = features_file,
    treeData = tree_file
)

taxa$feature_subset(Kingdom == "Bacteria")
taxa$scale(method = "tss")

manhattan(taxa$countData)

Description

Wrapper function that converts a sparseMatrix or countData field from omics classes to a data.table

Usage

matrix_to_dtable(x)

Arguments

Value

A data.table class.

Description

OmicFlow is constructed around the main superclass omics that in turn is inherited by metagenomics and proteomics, that offer additional fields or methods. The superclass omics contains both public and private methods, only the public componenents are documented, whereas the private methods can be accessed via metagenomics$private_methods, which is not supported.

The omics class components use Matrix and data.table in the background for fast loading and data wrangling, these data structures return the data by reference, so does the R6 class, therefore classes in OmicFlow do not create copies of the object unless called via $copy().

All classes needs to be first initialised via the $new() method and require the metaData and countData (or biomData in case of metagenomics) field components. Alternatively, you can also initialise with only the metaData and later add the other fields via the active binding, but this will not create a back-up and is not advised! For more hands-on examples see Getting Started with OmicFlow.

Additional arguments in metagenomics$new()

The metagenomics class also offers the loading of a phylogenetic tree as a phylo class (see as.phylo) that can be supplied to treeData. When a treeData is supplied than every class field is arranged according to the tree tip labels.

The biomData argument expects a BIOM v2.1.0. in either HDF5 or JSON format. The biomData must be supplied via a filepath and required contents are checked upfront, such as the existence of "/observation/matrix/data", "/observation/matrix/indptr", "/observation/matrix/indices", "/observation/metadata/taxonomy", "/observation/ids", "/sample/ids" .

Input requirements

Every class in OmicFlow requires at least a metaData and countData field to be specified. The metaData must contain a SAMPLE_ID and optionally a SAMPLEPAIR_ID, see section "Metadata validation" below.

The countData can be both a matrix or data frame. If the original contains rownames, then these will become the FEATURE_ID in featureData (automatically created), otherwise a random feature identifier is assigned. If the countData is a filetype or a data frame with both matching sample identifiers in SAMPLE_ID and other additional columns, then these will be automatically split, columns matching the SAMPLE_ID will be placed in the countData, any other columns will be placed in front of the FEATURE_ID column in featureData.

Metadata validation

Every class in OmicFlow begins by validating the sample metadata. This validation checks that sample identifiers match the sample IDs present in the abundance tables.

The validation rules are defined in JSON format in the omics class. Metadata can be supplied as either:

• A CSV or TSV file.

• A data.table.

In both cases, metadata must include a header row and follow a one-row-per-sample structure. For file inputs, the header must be on the first line. Columns not described below are permitted and are ignored during metadata validation.

Minimum requirement

The metadata must contain a SAMPLE_ID column.

• Each row must have a unique, non-empty sample identifier.

• Sample identifiers cannot contain spaces; use underscores (⁠_⁠) or dashes (-) instead.

Example metadata:

SAMPLE_ID  SAMPLEPAIR_ID  CONTRAST_Treatment  VARIABLE_Age
S1         P1             Drug                42
S2         P1             Placebo             36
S3         P2             Drug                51

Column types and naming rules

Required column

Optional standard columns

Some functions require a SAMEPLEPAIR_ID in the metaData when paired = TRUE is used. The user can add or remove the SAMPLEPAIR_ID at anytime to the metaData via active binding.

Pattern-based columns

Additional variables can be defined with special column-name prefixes:

• CONTRAST_...: Grouping or category labels for differential comparisons. For example, CONTRAST_Treatment may contain Drug and Placebo.

• VARIABLE_...: Numeric or string variables for statistical analysis. For example, VARIABLE_Age may contain values such as 42 and 51.

Pattern-based columns are only used in ⁠$autoFlow()⁠. Currently, only columns with the CONTRAST_ prefix are supported.

Super class

omics -> metagenomics

Active bindings

Methods

Public methods

• metagenomics$new()

• metagenomics$write_biom()

• metagenomics$clone()

metagenomics$new()

Usage

metagenomics$new(
  countData = NULL,
  metaData = NULL,
  featureData = NULL,
  treeData = NULL,
  biomData = NULL,
  feature_names = c("Kingdom", "Phylum", "Class", "Order", "Family", "Genus", "Species")
)

Arguments

Returns

A new metagenomics object.

Examples

library("OmicFlow")

## Method 1: Load from filepath
metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,  # optional
 treeData = tree_file          # optional
)

## Method 2: Load from a biom file
biom_file <- system.file("extdata", "biomv2.biom", package = "OmicFlow")

taxa <- metagenomics$new(
 metaData = metadata_file,
 biomData = biom_file
)

metagenomics$write_biom()

Creates a BIOM file in HDF5 format, which is compatible to the python biom-format v2.1.0. It only saves the featureData and countData to the biom format, no treeData is supported thus far.

Usage

metagenomics$write_biom(filename)

Arguments

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
 treeData = tree_file
)

taxa$write_biom(filename = "output.biom")
file.remove("output.biom")

metagenomics$clone()

The objects of this class are cloneable with this method.

Usage

metagenomics$clone(deep = FALSE)

Arguments

Examples

## ------------------------------------------------
## Method `metagenomics$new()`
## ------------------------------------------------

library("OmicFlow")

## Method 1: Load from filepath
metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,  # optional
 treeData = tree_file          # optional
)

## Method 2: Load from a biom file
biom_file <- system.file("extdata", "biomv2.biom", package = "OmicFlow")

taxa <- metagenomics$new(
 metaData = metadata_file,
 biomData = biom_file
)

## ------------------------------------------------
## Method `metagenomics$write_biom()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
 treeData = tree_file
)

taxa$write_biom(filename = "output.biom")
file.remove("output.biom")

Description

OmicFlow is constructed around the main superclass omics that in turn is inherited by metagenomics and proteomics, that offer additional fields or methods. The superclass omics contains both public and private methods, only the public componenents are documented, whereas the private methods can be accessed via omics$private_methods, which is not supported.

The omics class components use Matrix and data.table in the background for fast loading and data wrangling, these data structures return the data by reference, so does the R6 class, therefore classes in OmicFlow do not create copies of the object unless called via $copy().

All classes needs to be first initialised via the $new() method and require the metaData and countData field components. Alternatively, you can also initialise with only the metaData and later add the other fields via the active binding, but this will not create a back-up and is not advised! For more hands-on examples see Getting Started with OmicFlow.

Input requirements

Every class in OmicFlow requires at least a metaData and countData field to be specified. The metaData must contain a SAMPLE_ID and optionally a SAMPLEPAIR_ID, see section "Metadata validation" below.

The countData can be both a matrix or data frame. If the original contains rownames, then these will become the FEATURE_ID in featureData (automatically created), otherwise a random feature identifier is assigned. If the countData is a filetype or a data frame with both matching sample identifiers in SAMPLE_ID and other additional columns, then these will be automatically split, columns matching the SAMPLE_ID will be placed in the countData, any other columns will be placed in front of the FEATURE_ID column in featureData.

Metadata validation

Every class in OmicFlow begins by validating the sample metadata. This validation checks that sample identifiers match the sample IDs present in the abundance tables.

The validation rules are defined in JSON format in the omics class. Metadata can be supplied as either:

• A CSV or TSV file.

• A data.table.

In both cases, metadata must include a header row and follow a one-row-per-sample structure. For file inputs, the header must be on the first line. Columns not described below are permitted and are ignored during metadata validation.

Minimum requirement

The metadata must contain a SAMPLE_ID column.

• Each row must have a unique, non-empty sample identifier.

• Sample identifiers cannot contain spaces; use underscores (⁠_⁠) or dashes (-) instead.

Example metadata:

SAMPLE_ID  SAMPLEPAIR_ID  CONTRAST_Treatment  VARIABLE_Age
S1         P1             Drug                42
S2         P1             Placebo             36
S3         P2             Drug                51

Column types and naming rules

Required column

Optional standard columns

Some functions require a SAMEPLEPAIR_ID in the metaData when paired = TRUE is used. The user can add or remove the SAMPLEPAIR_ID at anytime to the metaData via active binding.

Pattern-based columns

Additional variables can be defined with special column-name prefixes:

• CONTRAST_...: Grouping or category labels for differential comparisons. For example, CONTRAST_Treatment may contain Drug and Placebo.

• VARIABLE_...: Numeric or string variables for statistical analysis. For example, VARIABLE_Age may contain values such as 42 and 51.

Pattern-based columns are only used in ⁠$autoFlow()⁠. Currently, only columns with the CONTRAST_ prefix are supported.

Active bindings

Methods

Public methods

• omics$new()

• omics$copy()

• omics$print()

• omics$reset()

• omics$removeNAs()

• omics$feature_subset()

• omics$sample_subset()

• omics$samplepair_subset()

• omics$feature_merge()

• omics$scale()

• omics$rankstat()

• omics$alpha_diversity()

• omics$composition()

• omics$distance()

• omics$ordination()

• omics$foldchange()

• omics$autoFlow()

• omics$clone()

omics$new()

Usage

omics$new(metaData = NULL, countData = NULL, featureData = NULL)

Arguments

Returns

A new omics object.

Examples

library("OmicFlow")

## Method 1: load from filepath
metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
obj <- omics$new(
 metaData = metadata_file,
 countData = counts_file
)

## Method 2: Load from data.frame or matrix
n_cols <- 5
n_rows <- 100
n_vals <- n_cols * n_rows
metadata <- data.frame("SAMPLE_ID" = paste0("Sample_", 1:n_cols))
features <- data.frame("FEATURE_ID" = paste0("protein_", 1:n_rows))
counts <- Matrix::Matrix(
 1:n_vals, nrow = n_rows, ncol = n_cols, 
 dimnames = list(features$FEATURE_ID, metadata$SAMPLE_ID)
)

obj <- omics$new(
 metaData = metadata,
 featureData = features, # optional
 countData = counts
)

## Method 3: you have features and counts in a single data.frame
metadata <- data.frame(SAMPLE_ID = c("S1", "S2", "S3"))
counts <- data.frame(
 S1 = c(2,3,0), 
 S2 = c(2,"",1), 
 S3 = c(2,1,NA), 
 proteins = c("ZEB1", "ZEB2", "MAPK")
)

obj <- omics$new(
 metaData = metadata,
 countData = counts
)

omics$copy()

This method is very similar to the existing ⁠$clone()⁠ function, except it also creates a new save of the back-up of each class that is created when ⁠$new()⁠ is called and can be reverted to via ⁠$reset()⁠

Usage

omics$copy(deep = FALSE)

Arguments

Returns

A copy of omics object.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")

obj <- omics$new(
 metaData = metadata_file,
 countData = counts_file
)

# Perform a modification and copy
obj$scale()

cloned <- obj$copy(deep=TRUE)
cloned$scale(method = "clr")
cloned$reset() # resets to data after clone creation.

omics$print()

Displays field parameters of the class via stdout.

Usage

omics$print()

Returns

A cli formatted message

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")

obj <- omics$new(
 metaData = metadata_file,
 countData = counts_file
)

# method 1 to call print function
obj

# method 2 to call print function
obj$print()

omics$reset()

Upon creation of a new class object (i.e. omics, metagenomics or proteomics) a backup of the original data is created. Since modification of the object is done by reference and duplicates are not made, it is possible to revert back to the unmodified data via $reset(). The methods from the super class omics also contains a private method to prevent any changes to the original object when using methods, such as ordination, alpha_diversity or foldchange.

Usage

omics$reset()

Returns

The object is modified in place and returned (invisibly). If a copy is required, please use ⁠$copy()⁠ beforehand with obj2 <- obj$copy().

Examples

library(ggplot2)
library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

taxa <- omics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file
)

# Performs modifications
taxa$scale(transform = log2)

# resets
taxa$reset()

# An inbuilt reset function prevents unwanted modification to the taxa object.
taxa$rankstat(feature_ranks = c("Kingdom", "Phylum", "Family", "Genus", "Species"))

omics$removeNAs()

Subset metaData via the complete.cases given a single or vector of column names, then synchronize the object.

Usage

omics$removeNAs(column)

Arguments

Returns

The object is modified in place and returned (invisibly). If a copy is required, please use ⁠$copy()⁠ beforehand with obj2 <- obj$copy().

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

obj$removeNAs(column = "treatment")

omics$feature_subset()

Subset features using a logical expression evaluated in featureData, then synchronize the object.

Usage

omics$feature_subset(expr)

Arguments

Returns

The object is modified in place and returned (invisibly). If a copy is required, please use ⁠$copy()⁠ beforehand with obj2 <- obj$copy().

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

obj$feature_subset(Genus == "Pseudomonas")

omics$sample_subset()

Subset samples using a logical expression evaluated in metaData, then synchronize the object.

Usage

omics$sample_subset(expr)

Arguments

Returns

The object is modified in place and returned (invisibly). If a copy is required, please use ⁠$copy()⁠ beforehand with obj2 <- obj$copy().

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

obj$sample_subset(treatment == "tumor")

omics$samplepair_subset()

Subset samplepairs by choosing the number of unique pairs in metaData, then synchronize the object.

Usage

omics$samplepair_subset(num_unique_pairs = NULL)

Arguments

Returns

The object is modified in place and returned (invisibly). If a copy is required, please use ⁠$copy()⁠ beforehand with obj2 <- obj$copy().

omics$feature_merge()

Agglomerates features by column, then synchronize the object.

Usage

omics$feature_merge(feature_rank, feature_filter = NULL)

Arguments

Returns

The object is modified in place and returned (invisibly). If a copy is required, please use ⁠$copy()⁠ beforehand with obj2 <- obj$copy().

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

obj$feature_merge(feature_rank = c("Kingdom", "Phylum"))
obj$feature_merge(feature_rank = "Genus", feature_filter = c("uncultured", "metagenome"))

omics$scale()

Feature scaling on the countData. The $scale function is able to apply transformations element-wise on the positive values, (optional: add pseudocounts) and perform normalisation or standardisation methods.

Usage

omics$scale(
  method = "tss",
  transform = NULL,
  base = exp(1),
  pseudocount = NULL
)

Arguments

Returns

The object is modified in place and returned (invisibly). If a copy is required, please use ⁠$copy()⁠ beforehand with obj2 <- obj$copy().

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)
# standard relative abundance computation
obj$scale()

# CLR
obj$reset()
obj$scale(method = "clr")

# transform
obj$reset()
obj$scale(method = "none", transform = log2)

omics$rankstat()

Rank statistics based on featureData. Counts the number of (unique) features identified for each column of interest from the featureData (excluding the FEATURE_ID).

Usage

omics$rankstat(feature_ranks, unique = FALSE)

Arguments

Returns

A ggplot object.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

plt <- obj$rankstat(feature_ranks = c("Kingdom", "Phylum", "Family", "Genus", "Species"))
plt

omics$alpha_diversity()

Alpha diversity computation, statistical testing and visualisation.

It computes the alpha diversity via diversity, performs a pairwise wilcoxon statistical test given groups via the boxjitter_test.

Usage

omics$alpha_diversity(
  col_name = lifecycle::deprecated(),
  groups,
  metric = "shannon",
  Brewer.palID = "Set2",
  group_by = lifecycle::deprecated(),
  split_by = NULL,
  evenness = FALSE,
  paired = FALSE,
  p.adjust.method = "fdr"
)

Arguments

Returns

A list of components:

data

A data.frame from diversity.

stats

A pairwise statistics from pairwise_wilcox_test.

plot

A ggplot object.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

plt <- obj$alpha_diversity(groups = "treatment",
                           metric = "shannon")

omics$composition()

Creates a long data.table of SAMPLE_IDs, the top 10-15 most abundant bacteria based on the feature_rank with any additional columns specified via add_cols. The $composition() method works with aggregated values by running first ⁠$feature_merge()⁠, then sorting for the most abundant bacteria based on the row sums. The output of $composition() can be used as input to composition_plot.

Usage

omics$composition(
  feature_rank,
  feature_filter = NULL,
  col_name = lifecycle::deprecated(),
  add_cols = NULL,
  feature_top = 10,
  Brewer.palID = "RdYlBu"
)

Arguments

Returns

A list of components:

data

A data.table of feature compositions.

palette

A setNames palette from colormap matching the top features.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

result <- obj$composition(feature_rank = "Genus",
                          feature_filter = c("uncultured"),
                          feature_top = 10)

plt <- composition_plot(data = result$data,
                        palette = result$palette,
                        feature_rank = "Genus")

omics$distance()

Compute a distance metric from countData

Usage

omics$distance(
  metric,
  weighted = TRUE,
  threads = 1,
  normalize = TRUE,
  base = exp(1)
)

Arguments

Returns

A column x column dist object.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
    metaData = metadata_file,
    countData = counts_file,
    featureData = features_file
)

obj$feature_subset(Kingdom == "Bacteria")
dist <- obj$distance(metric = "bray")

omics$ordination()

The $ordination() method either computes a distance matrix by specification of the metric or takes a custom distance matrix via distmat. It then computes either the pairwise PERMANOVA test on the defined contrast via group_by via pairwise_adonis or the pairwise ANOSIM test via pairwise_anosim.

Which pairwise test to compute is decided via the method argument, where "pcoa" results in the pairwise_adonis and "nmds" to pairwise_anosim.

It also visualises these results as a scree-plot (only when "pcoa" is used), scores-plot (see ordination_plot) and an anova-plot (see plot_pairwise_stats). These plots, alongisde the distance matrix, principal components are returned as a list.

Usage

omics$ordination(
  metric = "bray",
  method = "pcoa",
  group_by,
  distmat = NULL,
  weighted = TRUE,
  normalize = TRUE,
  threads = 1,
  perm_design = NULL,
  perm = 999
)

Arguments

Returns

A list of components:

distmat

A distance dissimilarity in matrix format.

stats

A statistical test as a data.frame.

pcs

principal components as a data.frame.

scree_plot

A ggplot object.

anova_plot

A ggplot object.

scores_plot

A ggplot object.

Examples

library("ggplot2")
library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

pcoa_plots <- obj$ordination(metric = "bray",
                             method = "pcoa",
                             group_by = "treatment",
                             weighted = TRUE)
pcoa_plots

omics$foldchange()

Differential feature expression for both paired and non-paired data.

The function performs feature agglomeration via ⁠$feature_merge()⁠ (which can be disabled via feature_merge=FALSE), subsetting to remove NAs in condition.group via ⁠$removeNAs()⁠ and finding samplepairs when paired = TRUE is supplied via ⁠$subset_samplepairs()⁠.

The fold-changes can be computed differently based on the method and aggregate_method options. Transformations of the data should be done beforehand via ⁠$scale()⁠. Finally, homogeneity of variance is computed based on the selected aggregate_method option. If aggregate_method = "mean" then the row_levene is applied, and if aggregate_method = "median" is used then the row_brownforsythe is applied. Any filtering of the results is left to the end-user, volcano plots are automatically added to the list name volcano_plot via volcano_plot.

Usage

omics$foldchange(
  condition.group,
  condition_A,
  condition_B,
  method = "identity",
  aggregate_method = "median",
  group_by = lifecycle::deprecated(),
  split_by = NULL,
  feature_merge = FALSE,
  feature_rank = "FEATURE_ID",
  feature_filter = NULL,
  paired = FALSE,
  pvalue.threshold = 0.05,
  logfold.threshold = 0.06,
  abundance.threshold = 0
)

Arguments

Returns

A list of components:

split_by subsets

A matrix subset from each split_by separated by condition_A and condition_B.

data

A data.table main output of fold-changes between conditions, contains abundance, fold-change, and homogeneity tests.

volcano_plot

A list of ggplot plots for each contrast of A vs B, number of plots depend on split_by and condition_A input.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- omics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file
)
## Method 1: using non-log transformation
obj$scale(method = "tss")

dfe <- obj$foldchange(
 feature_rank = "Genus",
 feature_merge = TRUE,
 method = "identity",
 aggregate_method = "median",
 condition.group = "treatment",
 condition_A = "tumor",
 condition_B = "healthy"
)

## Method 2: using log transformation
obj$reset()
obj$scale(method = "clr")

dfe <- obj$foldchange(
 feature_rank = "Genus",
 feature_merge = TRUE,
 method = "log",
 aggregate_method = "mean",
 condition.group = "treatment",
 condition_A = "tumor",
 condition_B = "healthy"
)

omics$autoFlow()

Automated Omics Analysis only works if there are column names in the metaData with the prefix "CONTRAST_", which is detected via sub-string pattern matching. The transformation of the data is left to the user via the ⁠$scale()⁠ method. The $autoFlow() incorporates the following functions:

$rankstat()

This function is only applied when feature_ranks is specified, useful in metagenomics, see ⁠$rankstat()⁠ for more information.

$alpha_diversity()

This function only works with the metagenomics class, uses "shannon" metric, and works for both paired and non-paired data, applied on each "CONTRAST_" column, see ⁠$alpha_diversity()⁠ for more information.

$ordination()

This function is applied on all classes, multiple distance_metrics can be supplied, and only the Principal Component method ($ordination(method = "pcoa")) is used. Additional arguments can be controlled with distance_metrics, distmat, weighted, perm, threads, see ⁠$ordination()⁠ for more information. NOTE: perm_design is not supported in ⁠$autoFlow()⁠.

$composition()

This function is applied on the significant contrasts that are obtained from either $alpha_diversity() or $ordination(). Multiple feature_contrasts can be supplied, by default it only uses the feature_contrasts = "FEATURE_ID". The composition plots are created by default via the composition_plot(group_by) with feature_top = 15. Additional arguments can be controlled with feature_filter, see ⁠$composition()⁠ for more information.

$foldchange()

This function is applied on the significant contrasts that are obtained from either $alpha_diversity() or $ordination(). Multiple feature_contrasts can be supplied, by default it only uses the feature_contrasts = "FEATURE_ID". Additional arguments can be controlled with fc_method, aggregate_method, pvalue.threshold, logfold.threshold, abundance.threshold, see ⁠$foldchange()⁠ for more information.

The above functions are called in order, only significant contrasts are further used to create composition and volcano plots. The autoFlow is intended as a first screening of certain contrasts and extra analysis, such as the use of a permutation design are not supported and should be applied by the user.

The report argument controls if an rmarkdown document is created, and with filename the output filename can be specified.

Usage

omics$autoFlow(
  feature_contrast = "FEATURE_ID",
  feature_filter = NULL,
  feature_ranks = NULL,
  distance_metrics = c("bray"),
  distmat = NULL,
  weighted = TRUE,
  fc_method = "log",
  aggregate_method = "median",
  pvalue.threshold = 0.05,
  logfold.threshold = 0.06,
  abundance.threshold = 0,
  perm = 999,
  threads = 1,
  report = TRUE,
  filename = paste0(getwd(), "/report.html")
)

Arguments

Returns

if report = TRUE a HTML report is generated, otherwise a list containing:

plots

• rankstat_plot as a ggplot.

• alpha_div_plots as a list of ggplot.

• pcoa_plots as a nested matrix of lists with ggplot (contrasts x distance_metrics).

• composition_plots as a nested matrix of lists with ggplot (contrasts x feature_contrasts).

• Log2FC_plots as a nested matrix of lists with ggplot (contrasts x feature_contrasts).

data

• alpha_div_data as a list of lists of data and stats.

• pcoa_data as a nested matrix of lists of stats, dist_mat and pcs (contrasts x distance_metrics).

• composition_data as a nested matrix of lists of data (contrasts x feature_contrasts).

• Log2FC_data as a nested matrix of lists of data (contrasts x feature_contrasts).

omics$clone()

The objects of this class are cloneable with this method.

Usage

omics$clone(deep = FALSE)

Arguments

Examples

## ------------------------------------------------
## Method `omics$new()`
## ------------------------------------------------

library("OmicFlow")

## Method 1: load from filepath
metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
obj <- omics$new(
 metaData = metadata_file,
 countData = counts_file
)

## Method 2: Load from data.frame or matrix
n_cols <- 5
n_rows <- 100
n_vals <- n_cols * n_rows
metadata <- data.frame("SAMPLE_ID" = paste0("Sample_", 1:n_cols))
features <- data.frame("FEATURE_ID" = paste0("protein_", 1:n_rows))
counts <- Matrix::Matrix(
 1:n_vals, nrow = n_rows, ncol = n_cols, 
 dimnames = list(features$FEATURE_ID, metadata$SAMPLE_ID)
)

obj <- omics$new(
 metaData = metadata,
 featureData = features, # optional
 countData = counts
)

## Method 3: you have features and counts in a single data.frame
metadata <- data.frame(SAMPLE_ID = c("S1", "S2", "S3"))
counts <- data.frame(
 S1 = c(2,3,0), 
 S2 = c(2,"",1), 
 S3 = c(2,1,NA), 
 proteins = c("ZEB1", "ZEB2", "MAPK")
)

obj <- omics$new(
 metaData = metadata,
 countData = counts
)

## ------------------------------------------------
## Method `omics$copy()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")

obj <- omics$new(
 metaData = metadata_file,
 countData = counts_file
)

# Perform a modification and copy
obj$scale()

cloned <- obj$copy(deep=TRUE)
cloned$scale(method = "clr")
cloned$reset() # resets to data after clone creation.


## ------------------------------------------------
## Method `omics$print()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")

obj <- omics$new(
 metaData = metadata_file,
 countData = counts_file
)

# method 1 to call print function
obj

# method 2 to call print function
obj$print()


## ------------------------------------------------
## Method `omics$reset()`
## ------------------------------------------------

library(ggplot2)
library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

taxa <- omics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file
)

# Performs modifications
taxa$scale(transform = log2)

# resets
taxa$reset()

# An inbuilt reset function prevents unwanted modification to the taxa object.
taxa$rankstat(feature_ranks = c("Kingdom", "Phylum", "Family", "Genus", "Species"))


## ------------------------------------------------
## Method `omics$removeNAs()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

obj$removeNAs(column = "treatment")


## ------------------------------------------------
## Method `omics$feature_subset()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

obj$feature_subset(Genus == "Pseudomonas")


## ------------------------------------------------
## Method `omics$sample_subset()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

obj$sample_subset(treatment == "tumor")


## ------------------------------------------------
## Method `omics$feature_merge()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

obj$feature_merge(feature_rank = c("Kingdom", "Phylum"))
obj$feature_merge(feature_rank = "Genus", feature_filter = c("uncultured", "metagenome"))


## ------------------------------------------------
## Method `omics$scale()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)
# standard relative abundance computation
obj$scale()

# CLR
obj$reset()
obj$scale(method = "clr")

# transform
obj$reset()
obj$scale(method = "none", transform = log2)


## ------------------------------------------------
## Method `omics$rankstat()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

plt <- obj$rankstat(feature_ranks = c("Kingdom", "Phylum", "Family", "Genus", "Species"))
plt

## ------------------------------------------------
## Method `omics$alpha_diversity()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

plt <- obj$alpha_diversity(groups = "treatment",
                           metric = "shannon")


## ------------------------------------------------
## Method `omics$composition()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

result <- obj$composition(feature_rank = "Genus",
                          feature_filter = c("uncultured"),
                          feature_top = 10)

plt <- composition_plot(data = result$data,
                        palette = result$palette,
                        feature_rank = "Genus")


## ------------------------------------------------
## Method `omics$distance()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
    metaData = metadata_file,
    countData = counts_file,
    featureData = features_file
)

obj$feature_subset(Kingdom == "Bacteria")
dist <- obj$distance(metric = "bray")

## ------------------------------------------------
## Method `omics$ordination()`
## ------------------------------------------------

library("ggplot2")
library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- metagenomics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file,
)

pcoa_plots <- obj$ordination(metric = "bray",
                             method = "pcoa",
                             group_by = "treatment",
                             weighted = TRUE)
pcoa_plots


## ------------------------------------------------
## Method `omics$foldchange()`
## ------------------------------------------------

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")

obj <- omics$new(
 metaData = metadata_file,
 countData = counts_file,
 featureData = features_file
)
## Method 1: using non-log transformation
obj$scale(method = "tss")

dfe <- obj$foldchange(
 feature_rank = "Genus",
 feature_merge = TRUE,
 method = "identity",
 aggregate_method = "median",
 condition.group = "treatment",
 condition_A = "tumor",
 condition_B = "healthy"
)

## Method 2: using log transformation
obj$reset()
obj$scale(method = "clr")

dfe <- obj$foldchange(
 feature_rank = "Genus",
 feature_merge = TRUE,
 method = "log",
 aggregate_method = "mean",
 condition.group = "treatment",
 condition_A = "tumor",
 condition_B = "healthy"
)

Description

Creates an ordination plot from pre-computed principal components. This function is built into the class omics with method ordination() and inherited by other omics classes, such as; metagenomics and proteomics.

Usage

ordination_plot(
  data,
  col_name = lifecycle::deprecated(),
  groups = lifecycle::deprecated(),
  pair,
  dist_explained = NULL,
  dist_metric = NULL
)

Arguments

Value

A ggplot2 object to be further modified

Examples

# Mock principal component scores
set.seed(123)
mock_data <- data.frame(
  SampleID = paste0("Sample", 1:10),
  PC1 = rnorm(10, mean = 0, sd = 1),
  PC2 = rnorm(10, mean = 0, sd = 1),
  groups = rep(c("Group1", "Group2"), each = 5)
)

# Basic usage
ordination_plot(
  data = mock_data,
  groups = "groups",
  pair = c("PC1", "PC2")
)

# Adding variance/dissimilarity explained.
ordination_plot(
  data = mock_data,
  groups = "groups",
  pair = c("PC1", "PC2"),
  dist_explained = c(45, 22),
  dist_metric = "bray-curtis"
)

Description

Computes pairwise adonis2, given a distance matrix and a vector of labels. This function is built into the class omics with method ordination() and inherited by other omics classes, such as; metagenomics and proteomics.

Usage

pairwise_adonis(
  x,
  groups,
  metadata = NULL,
  perm_design = NULL,
  p.adjust.method = "bonferroni",
  perm = 999,
  ...
)

Arguments

Value

A data.frame containing:

pairs

combinations of group comparisons

Df

Degrees of Freedom

SumsOfSqs

The sums of squares (centroid) of the null hypothesis $H_0$

F.Model

The F-test of the null hypothesis $H_0$

R2

The R squared of the null hypothesis $H_0$

p.value

The number of $F^p$ higher than the null hypothesis divided by the total number of permutations.

p.adj

The adjusted P-value based on the used p.adjust.method

See Also

adonis2

Examples

# Create random data
set.seed(42)
mock_data <- matrix(rnorm(15 * 10), nrow = 15, ncol = 10)

# Create euclidean dissimilarity matrix
mock_dist <- dist(mock_data, method = "euclidean")

# Define group labels, should be equal to number of columns and rows to dist
mock_groups <- rep(c("A", "B", "C"), each = 5)

# Compute pairwise adonis (PERMANOVA)
result <- pairwise_adonis(x = mock_dist, 
                          groups = mock_groups, 
                          p.adjust.method = "bonferroni", 
                          perm = 99)

Description

Computes pairwise anosim, given a distance matrix and a vector of labels. This function is built into the class omics with method ordination() and inherited by other omics classes, such as; metagenomics and proteomics.

Usage

pairwise_anosim(
  x,
  groups,
  metadata = NULL,
  perm_design = NULL,
  p.adjust.method = "bonferroni",
  perm = 999,
  ...
)

Arguments

Value

A data.frame containing:

pairs

combinations of group comparisons

Df

Degrees of Freedom

R2

The R squared of the null hypothesis $H_0$

p.value

The number of $F^p$ higher than the null hypothesis divided by the total number of permutations.

p.adj

The adjusted P-value based on the used p.adjust.method

See Also

anosim

Examples

# Create random data
set.seed(42)
mock_data <- matrix(rnorm(15 * 10), nrow = 15, ncol = 10)

# Create euclidean dissimilarity matrix
mock_dist <- dist(mock_data, method = "euclidean")

# Define group labels, should be equal to number of columns and rows to dist
mock_groups <- rep(c("A", "B", "C"), each = 5)

# Compute pairwise anosim
result <- pairwise_anosim(x = mock_dist, 
                          groups = mock_groups, 
                          p.adjust.method = "bonferroni", 
                          perm = 99)

Description

Creates a pairwise barplot of stats and F-test statistics plot from pairwise_adonis or pairwise_anosim results. This function is built into the class omics with method ordination() and inherited by other omics classes, such as; metagenomics and proteomics.

Usage

plot_pairwise_stats(
  data,
  stats_col,
  group_col,
  label_col,
  y_axis_title = NULL,
  plot_title = NULL
)

Arguments

Value

A ggplot2 object to be further modified.

Examples

# Create random data
set.seed(42)
mock_data <- matrix(rnorm(15 * 10), nrow = 15, ncol = 10)

# Create euclidean dissimilarity matrix
mock_dist <- dist(mock_data, method = "euclidean")

# Define group labels, should be equal to number of columns and rows to dist
mock_groups <- rep(c("A", "B", "C"), each = 5)

# Compute pairwise adonis
adonis_res <- pairwise_adonis(x = mock_dist, 
                              groups = mock_groups, 
                              p.adjust.method = "bonferroni", 
                              perm = 99)
# Compute pairwise anosim
anosim_res <- pairwise_anosim(x = mock_dist, 
                              groups = mock_groups, 
                              p.adjust.method = "bonferroni", 
                              perm = 99)

# Visualize PERMANOVA pairwise stats
plot_pairwise_stats(data = adonis_res,
                    group_col = "pairs",
                    stats_col = "F.Model",
                    label_col = "p.adj",
                    y_axis_title = "Pseudo F test statistic",
                    plot_title = "PERMANOVA")

# Visualize ANOSIM pairwise stats
plot_pairwise_stats(data = anosim_res,
                    group_col = "pairs",
                    stats_col = "anosimR",
                    label_col = "p.adj",
                    y_axis_title = "ANOSIM R statistic",
                    plot_title = "ANOSIM")

Description

OmicFlow is constructed around the main superclass omics that in turn is inherited by metagenomics and proteomics, that offer additional fields or methods. The superclass omics contains both public and private methods, only the public componenents are documented, whereas the private methods can be accessed via proteomics$private_methods, which is not supported.

The omics class components use Matrix and data.table in the background for fast loading and data wrangling, these data structures return the data by reference, so does the R6 class, therefore classes in OmicFlow do not create copies of the object unless called via $copy().

All classes needs to be first initialised via the $new() method and require the metaData and countData field components. Alternatively, you can also initialise with only the metaData and later add the other fields via the active binding, but this will not create a back-up and is not advised! For more hands-on examples see Getting Started with OmicFlow.

Additional arguments in proteomics$new()

The proteomics class also offers the loading of a phylogenetic tree as a phylo class (see as.phylo) that can be supplied to treeData. When a treeData is supplied than every class field is arranged according to the tree tip labels.

Input requirements

Every class in OmicFlow requires at least a metaData and countData field to be specified. The metaData must contain a SAMPLE_ID and optionally a SAMPLEPAIR_ID, see section "Metadata validation" below.

The countData can be both a matrix or data frame. If the original contains rownames, then these will become the FEATURE_ID in featureData (automatically created), otherwise a random feature identifier is assigned. If the countData is a filetype or a data frame with both matching sample identifiers in SAMPLE_ID and other additional columns, then these will be automatically split, columns matching the SAMPLE_ID will be placed in the countData, any other columns will be placed in front of the FEATURE_ID column in featureData.

Metadata validation

Every class in OmicFlow begins by validating the sample metadata. This validation checks that sample identifiers match the sample IDs present in the abundance tables.

The validation rules are defined in JSON format in the omics class. Metadata can be supplied as either:

• A CSV or TSV file.

• A data.table.

In both cases, metadata must include a header row and follow a one-row-per-sample structure. For file inputs, the header must be on the first line. Columns not described below are permitted and are ignored during metadata validation.

Minimum requirement

The metadata must contain a SAMPLE_ID column.

• Each row must have a unique, non-empty sample identifier.

• Sample identifiers cannot contain spaces; use underscores (⁠_⁠) or dashes (-) instead.

Example metadata:

SAMPLE_ID  SAMPLEPAIR_ID  CONTRAST_Treatment  VARIABLE_Age
S1         P1             Drug                42
S2         P1             Placebo             36
S3         P2             Drug                51

Column types and naming rules

Required column

Optional standard columns

Some functions require a SAMEPLEPAIR_ID in the metaData when paired = TRUE is used. The user can add or remove the SAMPLEPAIR_ID at anytime to the metaData via active binding.

Pattern-based columns

Additional variables can be defined with special column-name prefixes:

• CONTRAST_...: Grouping or category labels for differential comparisons. For example, CONTRAST_Treatment may contain Drug and Placebo.

• VARIABLE_...: Numeric or string variables for statistical analysis. For example, VARIABLE_Age may contain values such as 42 and 51.

Pattern-based columns are only used in ⁠$autoFlow()⁠. Currently, only columns with the CONTRAST_ prefix are supported.

Super class

omics -> proteomics

Active bindings

Methods

Public methods

• proteomics$new()

proteomics$new()

Usage

proteomics$new(
  countData = NULL,
  metaData = NULL,
  featureData = NULL,
  treeData = NULL
)

Arguments

Returns

A new proteomics object.

Examples

library("OmicFlow")

## Method 1: load from filepath
metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
obj <- proteomics$new(
 metaData = metadata_file,
 countData = counts_file
)

## Method 2: Load from data.frame or matrix
n_cols <- 5
n_rows <- 100
n_vals <- n_cols * n_rows
metadata <- data.frame("SAMPLE_ID" = paste0("Sample_", 1:n_cols))
features <- data.frame("FEATURE_ID" = paste0("protein_", 1:n_rows))
counts <- Matrix::Matrix(
 1:n_vals, nrow = n_rows, ncol = n_cols, 
 dimnames = list(features$FEATURE_ID, metadata$SAMPLE_ID)
)

obj <- proteomics$new(
 metaData = metadata,
 featureData = features, # optional
 countData = counts
)

## Method 3: you have features and counts in a single data.frame
metadata <- data.frame(SAMPLE_ID = c("S1", "S2", "S3"))
counts <- data.frame(
 S1 = c(2,3,0), 
 S2 = c(2,"",1), 
 S3 = c(2,1,NA), 
 proteins = c("ZEB1", "ZEB2", "MAPK")
)

obj <- proteomics$new(
 metaData = metadata,
 countData = counts
)

Examples

## ------------------------------------------------
## Method `proteomics$new()`
## ------------------------------------------------

library("OmicFlow")

## Method 1: load from filepath
metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
obj <- proteomics$new(
 metaData = metadata_file,
 countData = counts_file
)

## Method 2: Load from data.frame or matrix
n_cols <- 5
n_rows <- 100
n_vals <- n_cols * n_rows
metadata <- data.frame("SAMPLE_ID" = paste0("Sample_", 1:n_cols))
features <- data.frame("FEATURE_ID" = paste0("protein_", 1:n_rows))
counts <- Matrix::Matrix(
 1:n_vals, nrow = n_rows, ncol = n_cols, 
 dimnames = list(features$FEATURE_ID, metadata$SAMPLE_ID)
)

obj <- proteomics$new(
 metaData = metadata,
 featureData = features, # optional
 countData = counts
)

## Method 3: you have features and counts in a single data.frame
metadata <- data.frame(SAMPLE_ID = c("S1", "S2", "S3"))
counts <- data.frame(
 S1 = c(2,3,0), 
 S2 = c(2,"",1), 
 S3 = c(2,1,NA), 
 proteins = c("ZEB1", "ZEB2", "MAPK")
)

obj <- proteomics$new(
 metaData = metadata,
 countData = counts
)

Description

Calculates the UniFrac dissimilarity between samples based on phylogenetic branch lengths and abundance or presence/absence data.

Usage

unifrac(x, tree, weighted = TRUE, normalize = TRUE, threads = 1)

Arguments

Details

The UniFrac distance between two samples $A$ and $B$, with phylogenetic tree edges $i = 1 \ldots n$ of lengths $L_i$, is computed differently depending on the weighted and normalized flags. When weighted = FALSE, input counts are first converted to presence/absence data.

Weighted UniFrac (normalized = FALSE and weighted = TRUE):

$d(A,B) = \frac{\sum_{i}^n L_i |A_i - B_i|}{\sum_{i}^n L_i (A_i + B_i)}$

Normalized Weighted UniFrac (normalized = TRUE and weighted = TRUE):

$d(A,B) = \sum_{i}^n L_i |A_i - B_i|$

Unweighted UniFrac (weighted = FALSE, unweighted is always normalized):

$d(A,B) = \frac{\sum_{i}^n L_i |A_i - B_i|}{\sum_{i}^n L_i \max(A_i, B_i)}$

Value

A column x column dist object.

References

Lozupone, C., & Knight, R. (2005). UniFrac: a new phylogenetic method for comparing microbial communities. Applied and Environmental Microbiology, 71(12), 8228–8235.

Examples

library("OmicFlow")

metadata_file <- system.file("extdata", "metadata.tsv", package = "OmicFlow")
counts_file <- system.file("extdata", "counts.tsv", package = "OmicFlow")
features_file <- system.file("extdata", "features.tsv", package = "OmicFlow")
tree_file <- system.file("extdata", "tree.newick", package = "OmicFlow")

taxa <- metagenomics$new(
    metaData = metadata_file,
    countData = counts_file,
    featureData = features_file,
    treeData = tree_file
)

taxa$feature_subset(Kingdom == "Bacteria")
taxa$scale(method = "tss")

# Weighted UniFrac
unifrac(x = taxa$countData, tree = taxa$treeData, weighted=TRUE, normalize=FALSE)

# Weighted Normalized UniFrac
unifrac(x = taxa$countData, tree = taxa$treeData, weighted=TRUE, normalize=TRUE)

# Unweighted UniFrac
unifrac(x = taxa$countData, tree = taxa$treeData, weighted=FALSE)

Description

Creates a Volcano plot from the output of foldchange method from class omics, it plots the foldchanges on the x-axis, log10 trasnformed p-values on the y-axis and adjusts the scatter size based on the percentage abundance of the features. This function is built into the class omics with method DFE() and inherited by other omics classes, such as; metagenomics and proteomics.

Usage

volcano_plot(
  data,
  logfold_col,
  pvalue_col,
  abundance_col,
  feature_rank,
  pvalue.threshold = 0.05,
  logfold.threshold = 0.6,
  abundance.threshold = 0.01,
  label_A = "A",
  label_B = "B"
)

Arguments

Value

A ggplot2 object to be further modified.

Examples

# Create mock data frame
mock_volcano_data <- data.table::data.table(

  # Feature names (feature_rank)
  Feature = paste0("Gene", 1:20),   

  # Log2 fold changes (X)        
  log2FC = c(1.2, -1.5, 0.3, -0.7, 2.3,    
             -2.0, 0.1, 0.5, -1.0, 1.8,
             -0.4, 0.7, -1.4, 1.5, 0.9,
             -2.1, 0.2, 1.0, -0.3, -1.8),
  
  # P-values (Y)
  pvalue = c(0.001, 0.02, 0.3, 0.04, 0.0005, 
             0.01, 0.7, 0.5, 0.02, 0.0008,
             0.15, 0.06, 0.01, 0.005, 0.3,
             0.02, 0.8, 0.04, 0.12, 0.03),
  
  # Mean (relative) abundance for point sizing
  rel_abun = stats::runif(20, 0.01, 0.1)            
)

volcano_plot(
  data = mock_volcano_data,
  logfold_col = "log2FC",
  pvalue_col = "pvalue",
  abundance_col = "rel_abun",
  feature_rank = "Feature",
)