knitr::opts_chunk$set( collapse = TRUE, comment = "#>" )
library(rAMOPLS)
The rAMOPLS
package is a tool for helping you to extract relevant information from multivariate experiments with multiple factors. It provides simple functions to perform and assess an Anova Multiblock Orthogonal Partial Least Square analysis (AMOPLS) as described by @Boccard2016 with the specific calculation of Variable Important in the Projection (VIP²) described by @Gonzalez-Ruiz2017. Since most biological experiments are unbalanced due to a wide range of reasons, the method of stratified subsampling described by @Boccard2019 has also been implemented here.
Most biologist need to extract relevant information (treatment, tissular, organs, genotype, ...) from noisy and highly uncontrollable experiments. Most of us will try to use a design-of-experiment to assess any factors which can influence the response (wanted and unwanted variations). But then the tools to assess the effect of each factor and their interactions in multivariate experiments are not widely available. The AMOPLS methods combine ANOVA and k-OPLS (@Bylesjo2008) statistical methods to address these problems and provide easy-to-use diagnostics parameters to evaluate the effects.
rAMOPLS
can be installed using one of the following methods:
remotes
or devtools
packages utils::install.packages()
) The AMOPLS method depends on the kopls package which is available here. A copy of the package sources is implemented inside the rAMOPLS
and can be installed with rAMOPLS::install_kopls()
.
The single run_AMOPLS()
function is a wrapper to perform AMOPLS model with permutations and return loadings, scores, variances, saliences and statistics. It take as input two table (as matrix, data.frame or data.table) with the first column being unique identifier for rows and the column names from samplemetadata to study (factors):
datamatrix
: with observations as rows and measurments as columns samplemetadata
: with observations as rows and descriptor as columns (sample groups) factor_names
: Column names of samplemetadata with factor to consider Exemple on data_Ruiz2017
dataset on "Exposure time" and "Dose" with level 1 interaction and 100 permutations:
require(rAMOPLS) result <- run_AMOPLS(datamatrix = data_Ruiz2017$datamatrix, samplemetadata = data_Ruiz2017$samplemetadata, factor_names = c("Exposure time", "Dose"))
The following arguments are optional (default value):
interaction_level
(1): Order of interactions to consider (0 for main effect and 1 for first order interaction: Fac_A, Fac_B and Fac_A x Fac_B) scaling
(T): Logical to perform a variance scaling of the datamatrix nb_perm
(100): Number of permutation to test for each effect (default to 100) nb_compo_orthos
(1:3): Number of orthogonal component to consider (1:5 to test 1 to 5) parallel
(FALSE): If TRUE or a number (y), run permutations in parallel one y core using future and furrr for progress barThe optimal numnber parsimony principle is applied in the case of equally performing models with a different number of orthogonal components (@Boccard2016): The significative model with the smallest number of orthogonal component is to be choosen.
In this example the first model is already significant and should be choosen.
fun_plot_ortho(result) result_optimal <- result$orthoNb_1
All graphical output use ggplot2
and can be customized according to ggplot2 nomenclature.
fun_get_summary(result_optimal) %>% knitr::kable(digits = 3, align = 'c')
fun_plot_optimal_scores(result_optimal)
fun_plot_optimal_loadings(result_optimal)
fun_plot_VIPs(result_optimal, "Dose")
Unbalanced experimental designs involve non-orthogonal ANOVA submatrices which can alter effects interpretation (@Drotleff2019). To cope with unbalanced design, the most strategy is to resample the dataset to get a balanced design by either:
Since oversampling may alter variance decomposition by incorporating identical samples to the model, @Boccard2019 used a strategy based on undersampling strategy. To cope with minimized dataset (and therfore deleted observation), the subsampling is randomly performed a high number of times (1000s) and the models are combined using ensemble-based estimates from balanced models (involving median calculation of each parameters obtained from AMOPLS).
This strategy as been implemented in rAMOPLS and is automatically triggered if the input dataset are unbalanced. The number of subsampling steps can be customized with the subsampling
parameter.
Example on the same dataset as @Boccard2019:
require(rAMOPLS) result_unbalanced <- run_AMOPLS( datamatrix = data_Boccard2019$datamatrix, samplemetadata = data_Boccard2019$samplemetadata, factor_names = c("Chemical", "Batch"), nb_perm = 100, subsampling = 10, parallel = 3) result_optimal <- result_unbalanced[[1]]
require(rAMOPLS) result_unbalanced <- invisible( run_AMOPLS( datamatrix = data_Boccard2019$datamatrix, samplemetadata = data_Boccard2019$samplemetadata, factor_names = c("Chemical", "Batch"), nb_perm = 100, subsampling = 10, parallel = 3) ) result_optimal <- result_unbalanced[[1]]
fun_get_summary(result_optimal) %>% knitr::kable(digits = 3, align = 'c')
fun_plot_optimal_scores(result_optimal)
fun_plot_optimal_loadings(result_optimal)
fun_plot_VIPs(result_optimal, "Chemical")
The first step is to ensure the consitency of the ANOVA partitionning, the number of orthogonal components to include in the model and the reliability of the observed effects using the following indices:
Random experimental designs are expected to produce low R2Y values (no model) and RSR indices close to one (no effect). When both the AMOPLS model R2Y and individual RSR indices are deemed significant, biochemical information can be extracted from the loadings of the predictive components (@Boccard2016).
You can find more detailled explanation on AMOPLS in the references provided at the end of this page, starting by the original article from @Boccard2016.
The package is preloaded with the following datasets to test and compare the functionalities with published articles:
liver.toxicity
: Microarray of gene expression from Fisher rats after acetaminophen exposition @Bushel2007 and used in @Boccard2016 data_Boccard2016
: Metabolomic profiles of A. thaliana leaves from @Boccard2010 used in @Boccard2016 data_Ruiz2017
: Metabolomic profiles of human neural cells from @Gonzalez-Ruiz2017 data_Boccard2019
: Metabolomic profiles used to demonstrate the subsampling strategy in @Boccard2019 Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.