Title: | Mean-Variance Regularization |
---|---|
Description: | Implements a non-parametric method for joint adaptive mean-variance regularization and variance stabilization of high-dimensional data. It is suited for handling difficult problems posed by high-dimensional multivariate datasets (p >> n paradigm). Among those are that the variance is often a function of the mean, variable-specific estimators of variances are not reliable, and tests statistics have low powers due to a lack of degrees of freedom. Key features include: (i) Normalization and/or variance stabilization of the data, (ii) Computation of mean-variance-regularized t-statistics (F-statistics to follow), (iii) Generation of diverse diagnostic plots, (iv) Computationally efficient implementation using C/C++ interfacing and an option for parallel computing to enjoy a faster and easier experience in the R environment. |
Authors: | Jean-Eudes Dazard [aut, cre], Hua Xu [ctb], Alberto Santana [ctb] |
Maintainer: | Jean-Eudes Dazard <[email protected]> |
License: | GPL (>= 3) | file LICENSE |
Version: | 1.34.0 |
Built: | 2025-02-05 03:24:23 UTC |
Source: | https://github.com/jedazard/mvr |
Implements a non-parametric method for joint adaptive mean-variance regularization and variance stabilization of high-dimensional data (Dazard and Rao, 2012).
It is suited for handling difficult problems posed by high-dimensional multivariate datasets
( paradigm), such as in omics-type data, among which are that the variance is
often a function of the mean, variable-specific estimators of variances are not reliable,
and tests statistics have low powers due to a lack of degrees of freedom.
Key features include:
Normalization and/or variance stabilization of the data
Computation of mean-variance-regularized t-statistics (F-statistics to come)
Generation of diverse diagnostic plots
Computationally efficient implementation using C/C++ interfacing
and an option for parallel computing to enjoy a fast and easy experience
in the R environment
The following describes all the end-user functions, and internal R subroutines needed for running a complete MVR procedure. Other internal subroutines are not to be called by the end-user at any time. For computational efficiency, end-user regularization functions offer the option to configure a cluster. This is indicated by an asterisk (* = optionally involving cluster usage). The R functions are categorized as follows:
END-USER REGULARIZATION & VARIANCE STABILIZATION FUNCTION mvr
(*)
Function for Mean-Variance Regularization and Variance Stabilization.
End-user function for Mean-Variance Regularization (MVR) and Variance Stabilization
by similarity statistic under sample group homoscedasticity or heteroscedasticity assumption.
The function takes advantage of the R package parallel, which allows users
to create a cluster of workstations on a local and/or remote machine(s),
enabling parallel execution of this function and scaling up with the number of CPU
cores available.
END-USER REGULARIZED TESTS-STATISTICS FUNCTIONS mvrt.test
(*)
Function for Computing Mean-Variance Regularized T-test Statistic and Its Significance.
End-user function for computing MVR t-test statistic and its significance (p-value)
under sample group homoscedasticity or heteroscedasticity assumption.
The function takes advantage of the R package parallel, which allows users
to create a cluster of workstations on a local and/or remote machine(s),
enabling parallel execution of this function and scaling up with the number of CPU
cores available.
END-USER DIAGNOSTIC PLOTS FOR QUALITY CONTROL cluster.diagnostic
Function for Plotting Summary Cluster Diagnostic Plots.
Plot similarity statistic profiles and the optimal joint clustering configuration for the means and the variances by group.
Plot quantile profiles of means and standard deviations by group and for each clustering configuration, to check that the distributions
of first and second moments of the MVR-transformed data approach their respective null distributions under the optimal configuration found,
assuming independence and normality of all the variables.
target.diagnostic
Function for Plotting Summary Target Moments Diagnostic Plots.
Plot comparative distribution densities of means and standard deviations of the data before and after Mean-Variance Regularization
to check for location shifts between observed first and second moments and their expected target values under a target centered homoscedastic model.
Plot comparative QQ scatterplots to look at departures between observed distributions of first and second moments of the MVR-transformed data
and their theoretical distributions assuming independence and normality of all the variables.
stabilization.diagnostic
Function for Plotting Summary Variance Stabilization Diagnostic Plots.
Plot comparative variance-mean plots to check the variance stabilization across variables before and after Mean-Variance Regularization.
normalization.diagnostic
Function for Plotting Summary Normalization Diagnostic Plots.
Plot comparative Box-Whisker and Heatmap plots of variables across samples check the effectiveness of normalization before and after
Mean-Variance Regularization.
OTHER END-USER FUNCTIONS MVR.news
Display the MVR Package News
Function to display the log file NEWS
of updates of the MVR package.
END-USER DATASETS
A Real
dataset coming from a quantitative proteomics experiment, consisting of
samples split into a control ("M") and a treated group ("S") with
unique peptides or predictor variables. This is a balanced design with two sample groups (
),
under unequal sample group variance.
A Synthetic
dataset with observations (samples) and
variables,
where
of them are significantly different between the two sample groups.
This is a balanced design with two sample groups (
), under unequal sample group variance.
Known Bugs/Problems : None at this time.
This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University. This project was partially funded by the National Institutes of Health (P30-CA043703).
"Jean-Eudes Dazard, Ph.D." [email protected]
"Hua Xu, Ph.D." [email protected]
"Alberto Santana, MBA." [email protected]
Maintainer: "Jean-Eudes Dazard, Ph.D." [email protected]
Dazard J-E. and J. S. Rao (2010). "Regularized Variance Estimation and Variance Stabilization of High-Dimensional Data." In JSM Proceedings, Section for High-Dimensional Data Analysis and Variable Selection. Vancouver, BC, Canada: American Statistical Association IMS - JSM, 5295-5309.
Dazard J-E., Hua Xu and J. S. Rao (2011). "R package MVR for Joint Adaptive Mean-Variance Regularization and Variance Stabilization." In JSM Proceedings, Section for Statistical Programmers and Analysts. Miami Beach, FL, USA: American Statistical Association IMS - JSM, 3849-3863.
Dazard J-E. and J. S. Rao (2012). "Joint Adaptive Mean-Variance Regularization and Variance Stabilization of High Dimensional Data." Comput. Statist. Data Anal. 56(7):2317-2333.
makeCluster
(R package parallel)
justvsn
(R package vsn) Variance stabilization and calibration for microarray data Huber, 2002
eBayes
(R package limma) Bayesian Regularized t-test statistic Smyth, 2004
samr
(R package samr) SAM Regularized t-test statistic Tusher et al., 2001, Storey, 2003
matest
(R package maanova) James-Stein shrinkage estimator-based Regularized t-test statistic Cui et al., 2005
ebam
(R package siggenes) Empirical Bayes Regularized z-test statistic Efron, 2001
bayesT
Hierarchical Bayesian Regularized t-test statistic Baldi et al., 2001
Plot similarity statistic profiles and the optimal joint clustering configuration for the means and the variances by group.
Plot quantile profiles of means and standard deviations by group and for each clustering configuration, to check that the distributions of first and second moments of the MVR-transformed data approach their respective null distributions under the optimal configuration found, assuming independence and normality of all the variables.
cluster.diagnostic(obj, span = 0.75, degree = 2, family = "gaussian", title = "Cluster Diagnostic Plots", device = NULL, file = "Cluster Diagnostic Plots", path = getwd(), horizontal = FALSE, width = 8.5, height = 11, ...)
cluster.diagnostic(obj, span = 0.75, degree = 2, family = "gaussian", title = "Cluster Diagnostic Plots", device = NULL, file = "Cluster Diagnostic Plots", path = getwd(), horizontal = FALSE, width = 8.5, height = 11, ...)
obj |
Object of class " |
title |
Title of the plot. Defaults to "Cluster Diagnostic Plots". |
span |
Span parameter of the |
degree |
Degree parameter of the |
family |
Family distribution in "gaussian", "symmetric" of the |
device |
Graphic display device in {NULL, "PS", "PDF"}. Defaults to NULL (standard output screen). Currently implemented graphic display devices are "PS" (Postscript) or "PDF" (Portable Document Format). |
file |
File name for output graphic. Defaults to "Cluster Diagnostic Plots". |
path |
Absolute path (without final (back)slash separator). Defaults to working directory path. |
horizontal |
|
width |
|
height |
|
... |
Generic arguments passed to other plotting functions. |
In a plot of a similarity statistic profile, one checks the goodness of fit of the transformed data relative to the hypothesized underlying reference
distribution with mean-0 and standard deviation-1 (e.g. ). The red dashed line depicts the LOESS scatterplot smoother estimator.
The subroutine internally generates reference null distributions for computing the similarity statistic under each cluster configuration.
The optimal cluster configuration (indicated by the vertical red arrow) is found where the similarity statistic reaches its minimum plus/minus
one standard deviation (applying the conventional one-standard deviation rule). A smaller cluster number configuration indicates under-regularization,
while over-regularization starts to occur at larger numbers. This over/under-regularization must be viewed as a form of over/under-fitting
(see Dazard, J-E. and J. S. Rao (2012) for more details).
The quantile diagnostic plots uses empirical quantiles of the transformed means and standard deviations to check how
closely they are approximated by theoretical quantiles derived from a standard normal equal-mean/homoscedastic
model (solid green lines) under a given cluster configuration. To assess this goodness of fit of the transformed data, theoretical null distributions
of the mean and variance are derived from a standard normal equal-mean/homoscedastic model with independence of the first two moments,
i.e. assuming i.i.d. normality of the raw data. However, we do not require i.i.d. normality of the data in general: these theoretical null distributions are
just used here as convenient ones to draw from. Note that under the assumptions that the raw data is i.i.d. standard normal ($N(0, 1)$)
with independence of first two moments, the theoretical null distributions of means and standard deviations for each variable
are respectively:
and
, where
denotes the number of sample groups.
The optimal cluster configuration found is indicated by the most horizontal red curve. The single cluster configuration, corresponding to no transformation,
is the most vertical curve, while the largest cluster number configuration reaches horizontality. Notice how empirical quantiles of transformed
pooled means and standard deviations converge (from red to black) to the theoretical null distributions (solid green lines) for the optimal
configuration. One should see a convergence towards the target null, after which overfitting starts to occur (see Dazard, J-E. and J. S. Rao (2012)
for more details).
Both cluster diagnostic plots help determine (i) whether the minimum of the Similarity Statistic is observed within the range of clusters
(i.e. a large enough number of clusters has been accommodated), and (ii) whether the corresponding cluster configuration is a good fit.
If necessary, run the procedure again with larger value of the
nc.max
parameter in the mvr
as well as
in mvrt.test
functions until the minimum of the similarity statistic profile is reached.
Option file
is used only if device is specified (i.e. non NULL
).
None. Displays the plots on the chosen device
.
This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University. This project was partially funded by the National Institutes of Health (P30-CA043703).
End-user function.
"Jean-Eudes Dazard, Ph.D." [email protected]
"Hua Xu, Ph.D." [email protected]
"Alberto Santana, MBA." [email protected]
Maintainer: "Jean-Eudes Dazard, Ph.D." [email protected]
Dazard J-E. and J. S. Rao (2010). "Regularized Variance Estimation and Variance Stabilization of High-Dimensional Data." In JSM Proceedings, Section for High-Dimensional Data Analysis and Variable Selection. Vancouver, BC, Canada: American Statistical Association IMS - JSM, 5295-5309.
Dazard J-E., Hua Xu and J. S. Rao (2011). "R package MVR for Joint Adaptive Mean-Variance Regularization and Variance Stabilization." In JSM Proceedings, Section for Statistical Programmers and Analysts. Miami Beach, FL, USA: American Statistical Association IMS - JSM, 3849-3863.
Dazard J-E. and J. S. Rao (2012). "Joint Adaptive Mean-Variance Regularization and Variance Stabilization of High Dimensional Data." Comput. Statist. Data Anal. 56(7):2317-2333.
loess
(R package stats) Fit a polynomial surface determined by one or more numerical predictors, using local fitting.
## Not run: #=================================================== # Loading the library and its dependencies #=================================================== library("MVR") #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # (see description of datasets) #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups # Without cluster usage =================================================== nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, log = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) #=================================================== # Summary Cluster Diagnostic Plots (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== cluster.diagnostic(obj = mvr.obj, title = "Cluster Diagnostic Plots (Real - Multi-Group Assumption)", span = 0.75, degree = 2, family = "gaussian", device = NULL, horizontal = FALSE, width = 8.5, height = 11) ## End(Not run)
## Not run: #=================================================== # Loading the library and its dependencies #=================================================== library("MVR") #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # (see description of datasets) #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups # Without cluster usage =================================================== nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, log = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) #=================================================== # Summary Cluster Diagnostic Plots (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== cluster.diagnostic(obj = mvr.obj, title = "Cluster Diagnostic Plots (Real - Multi-Group Assumption)", span = 0.75, degree = 2, family = "gaussian", device = NULL, horizontal = FALSE, width = 8.5, height = 11) ## End(Not run)
End-user function for Mean-Variance Regularization (MVR) and Variance Stabilization by similarity statistic under sample group homoscedasticity or heteroscedasticity assumptions.
Return an object of class "MVR
". Offers the option of parallel computation for improved efficiency.
mvr(data, block = rep(1,nrow(data)), tolog = FALSE, nc.min = 1, nc.max = 30, probs = seq(0, 1, 0.01), B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = NULL)
mvr(data, block = rep(1,nrow(data)), tolog = FALSE, nc.min = 1, nc.max = 30, probs = seq(0, 1, 0.01), B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = NULL)
data |
|
block |
|
tolog |
|
nc.min |
Positive |
nc.max |
Positive |
probs |
|
B |
Positive |
parallel |
|
conf |
|
verbose |
|
seed |
Positive |
Argument block
will be converted to a factor, whose levels will match the data groups. It defaults to a single group situation,
that is, under the assumption of equal variance between sample groups. All group sample sizes must be greater than 1,
otherwise the program will stop.
Argument nc.max
currently defaults to 30. Empirically, we found that this is enough for most datasets tested.
This depends on (i) the dimensionality/sample size ratio , (ii) the signal/noise ratio, and
(iii) whether a pre-transformation has been applied (see Dazard, J-E. and J. S. Rao (2012) for more details).
See the cluster diagnostic function
cluster.diagnostic
for more details, whether larger values of nc.max
may be required.
The function mvr
relies on the R package parallel to create a parallel backend within an R session. This enables access to a cluster
of compute cores and/or nodes on a local and/or remote machine(s) and scaling-up with the number of CPU cores available and efficient parallel
execution. To run a procedure in parallel (with parallel RNG), argument parallel
is to be set to TRUE
and argument conf
is to be specified (i.e. non NULL
). Argument conf
uses the options described in function makeCluster
of the R packages
parallel and snow. PRIMsrc supports two types of communication mechanisms between master and worker processes:
'Socket' or 'Message-Passing Interface' ('MPI'). In PRIMsrc, parallel 'Socket' clusters use sockets communication mechanisms only
(no forking) and are therefore available on all platforms, including Windows, while parallel 'MPI' clusters use high-speed interconnects
mechanism in networks of computers (with distributed memory) and are therefore available only in these architectures. A parallel 'MPI'
cluster also requires R package Rmpi to be installed first. Value type
is used to setup a cluster of type 'Socket' ("SOCKET")
or 'MPI' ("MPI"), respectively. Depending on this type, values of spec
are to be used alternatively:
For 'Socket' clusters (conf$type="SOCKET"
), spec
should be a character
vector
naming the hosts on which
to run the job; it can default to a unique local machine, in which case, one may use the unique host name "localhost".
Each host name can potentially be repeated to the number of CPU cores available on the local machine.
It can also be an integer
scalar specifying the number of processes to spawn on the local machine;
or a list of machine specifications (a character value named host specifying the name or address of the host to use).
For 'MPI' clusters (conf$type="MPI"
), spec
should be an integer
scalar
specifying the total number of processes to be spawned across the network of available nodes, counting the workernodes and masternode.
The actual creation of the cluster, its initialization, and closing are all done internally. For more details, see the reference manual of R package snow and examples below.
When random number generation is needed, the creation of separate streams of parallel RNG per node is done internally by distributing the stream states to the nodes. For more details, see the vignette of R package parallel. The use of a seed allows to reproduce the results within the same type of session: the same seed will reproduce the same results within a non-parallel session or within a parallel session, but it will not necessarily give the exact same results (up to sampling variability) between a non-parallelized and parallelized session due to the difference of management of the seed between the two (see parallel RNG and value of returned seed below).
Xraw |
|
Xmvr |
|
centering |
|
scaling |
|
MVR |
|
block |
Value of argument |
tolog |
Value of argument |
nc.min |
Value of argument |
nc.max |
Value of argument |
probs |
Value of argument |
seed |
User seed(s) used: |
This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University. This project was partially funded by the National Institutes of Health (P30-CA043703).
End-user function.
"Jean-Eudes Dazard, Ph.D." [email protected]
"Hua Xu, Ph.D." [email protected]
"Alberto Santana, MBA." [email protected]
Maintainer: "Jean-Eudes Dazard, Ph.D." [email protected]
Dazard J-E. and J. S. Rao (2010). "Regularized Variance Estimation and Variance Stabilization of High-Dimensional Data." In JSM Proceedings, Section for High-Dimensional Data Analysis and Variable Selection. Vancouver, BC, Canada: American Statistical Association IMS - JSM, 5295-5309.
Dazard J-E., Hua Xu and J. S. Rao (2011). "R package MVR for Joint Adaptive Mean-Variance Regularization and Variance Stabilization." In JSM Proceedings, Section for Statistical Programmers and Analysts. Miami Beach, FL, USA: American Statistical Association IMS - JSM, 3849-3863.
Dazard J-E. and J. S. Rao (2012). "Joint Adaptive Mean-Variance Regularization and Variance Stabilization of High Dimensional Data." Comput. Statist. Data Anal. 56(7):2317-2333.
makeCluster
(R package parallel).
justvsn
(R package vsn) Variance stabilization and calibration for microarray data Huber, 2002
#=================================================== # Loading the library and its dependencies #=================================================== library("MVR") ## Not run: #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # Use help for descriptions #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real ## End(Not run) #=================================================== # Mean-Variance Regularization (Synthetic dataset) # Single-Group Assumption # Assuming equal variance between groups # Without cluster usage #=================================================== data("Synthetic", package="MVR") nc.min <- 1 nc.max <- 10 probs <- seq(0, 1, 0.01) n <- 10 mvr.obj <- mvr(data = Synthetic, block = rep(1,n), tolog = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) ## Not run: #=================================================== # Examples of parallel backend parametrization #=================================================== if (require("parallel")) { cat("'parallel' is attached correctly \n") } else { stop("'parallel' must be attached first \n") } #=================================================== # Ex. #1 - Multicore PC # Running WINDOWS # SOCKET communication cluster # Shared memory parallelization #=================================================== cpus <- parallel::detectCores(logical = TRUE) conf <- list("spec" = rep("localhost", cpus), "type" = "SOCKET", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Ex. #2 - Master node + 3 Worker nodes cluster # All nodes equipped with identical setups of multicores # (8 core CPUs per machine for a total of 32) # SOCKET communication cluster # Distributed memory parallelization #=================================================== masterhost <- Sys.getenv("HOSTNAME") slavehosts <- c("compute-0-0", "compute-0-1", "compute-0-2") nodes <- length(slavehosts) + 1 cpus <- 8 conf <- list("spec" = c(rep(masterhost, cpus), rep(slavehosts, cpus)), "type" = "SOCKET", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Ex. #3 - Enterprise Multinode Cluster w/ multicore/node # Running LINUX with SLURM scheduler # MPI communication cluster # Distributed memory parallelization # Below, variable 'cpus' is the total number of requested # taks (threads/CPUs), which is specified from within a # SLURM script. #================================================== if (require("Rmpi")) { print("Rmpi is loaded correctly \n") } else { stop("Rmpi must be installed first to use MPI\n") } cpus <- as.numeric(Sys.getenv("SLURM_NTASKS")) conf <- list("spec" = cpus, "type" = "MPI", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== data("Real", package="MVR") nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, tolog = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = TRUE, conf = conf, verbose = TRUE, seed = 1234) ## End(Not run)
#=================================================== # Loading the library and its dependencies #=================================================== library("MVR") ## Not run: #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # Use help for descriptions #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real ## End(Not run) #=================================================== # Mean-Variance Regularization (Synthetic dataset) # Single-Group Assumption # Assuming equal variance between groups # Without cluster usage #=================================================== data("Synthetic", package="MVR") nc.min <- 1 nc.max <- 10 probs <- seq(0, 1, 0.01) n <- 10 mvr.obj <- mvr(data = Synthetic, block = rep(1,n), tolog = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) ## Not run: #=================================================== # Examples of parallel backend parametrization #=================================================== if (require("parallel")) { cat("'parallel' is attached correctly \n") } else { stop("'parallel' must be attached first \n") } #=================================================== # Ex. #1 - Multicore PC # Running WINDOWS # SOCKET communication cluster # Shared memory parallelization #=================================================== cpus <- parallel::detectCores(logical = TRUE) conf <- list("spec" = rep("localhost", cpus), "type" = "SOCKET", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Ex. #2 - Master node + 3 Worker nodes cluster # All nodes equipped with identical setups of multicores # (8 core CPUs per machine for a total of 32) # SOCKET communication cluster # Distributed memory parallelization #=================================================== masterhost <- Sys.getenv("HOSTNAME") slavehosts <- c("compute-0-0", "compute-0-1", "compute-0-2") nodes <- length(slavehosts) + 1 cpus <- 8 conf <- list("spec" = c(rep(masterhost, cpus), rep(slavehosts, cpus)), "type" = "SOCKET", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Ex. #3 - Enterprise Multinode Cluster w/ multicore/node # Running LINUX with SLURM scheduler # MPI communication cluster # Distributed memory parallelization # Below, variable 'cpus' is the total number of requested # taks (threads/CPUs), which is specified from within a # SLURM script. #================================================== if (require("Rmpi")) { print("Rmpi is loaded correctly \n") } else { stop("Rmpi must be installed first to use MPI\n") } cpus <- as.numeric(Sys.getenv("SLURM_NTASKS")) conf <- list("spec" = cpus, "type" = "MPI", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== data("Real", package="MVR") nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, tolog = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = TRUE, conf = conf, verbose = TRUE, seed = 1234) ## End(Not run)
Function to display the NEWS file of the MVR package.
MVR.news(...)
MVR.news(...)
... |
Further arguments passed to or from other methods. |
None.
This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University. This project was partially funded by the National Institutes of Health (P30-CA043703).
End-user function.
"Jean-Eudes Dazard, Ph.D." [email protected]
"Hua Xu, Ph.D." [email protected]
"Alberto Santana, MBA." [email protected]
Maintainer: "Jean-Eudes Dazard, Ph.D." [email protected]
Dazard J-E. and J. S. Rao (2010). "Regularized Variance Estimation and Variance Stabilization of High-Dimensional Data." In JSM Proceedings, Section for High-Dimensional Data Analysis and Variable Selection. Vancouver, BC, Canada: American Statistical Association IMS - JSM, 5295-5309.
Dazard J-E., Hua Xu and J. S. Rao (2011). "R package MVR for Joint Adaptive Mean-Variance Regularization and Variance Stabilization." In JSM Proceedings, Section for Statistical Programmers and Analysts. Miami Beach, FL, USA: American Statistical Association IMS - JSM, 3849-3863.
Dazard J-E. and J. S. Rao (2012). "Joint Adaptive Mean-Variance Regularization and Variance Stabilization of High Dimensional Data." Comput. Statist. Data Anal. 56(7):2317-2333.
End-user function for computing MVR t-test statistic and its significance (p-value) under sample group homoscedasticity or heteroscedasticity assumption.
Return an object of class "MVRT
". Offers the option of parallel computation for
improved efficiency.
mvrt.test(data, obj=NULL, block, tolog = FALSE, nc.min = 1, nc.max = 30, pval = FALSE, replace = FALSE, n.resamp = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = NULL)
mvrt.test(data, obj=NULL, block, tolog = FALSE, nc.min = 1, nc.max = 30, pval = FALSE, replace = FALSE, n.resamp = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = NULL)
data |
|
obj |
Object of class " |
block |
|
tolog |
|
nc.min |
Positive |
nc.max |
Positive |
pval |
|
replace |
|
n.resamp |
Positive |
parallel |
|
conf |
|
verbose |
|
seed |
Positive |
Argument block
will be converted to a factor, whose levels will match the data groups. It defaults to a single group situation,
that is, under the assumption of equal variance between sample groups. All group sample sizes must be greater than 1,
otherwise the program will stop.
Argument nc.max
currently defaults to 30. Empirically, we found that this is enough for most datasets tested.
This depends on (i) the dimensionality/sample size ratio , (ii) the signal/noise ratio, and
(iii) whether a pre-transformation has been applied (see Dazard, J-E. and J. S. Rao (2012) for more details).
See the cluster diagnostic function
cluster.diagnostic
for more details, whether larger values of nc.max
may be required.
To save un-necessary computations, previously computed MVR clustering can be provided through option obj
(i.e. obj
is fully
specified as a mvr
object). In this case, arguments data
, block
, tolog
, nc.min
, nc.max
are ignored. If obj
is fully specified (i.e. an object of class "MVR
" returned by mvr
), the
the MVR clustering provided by obj
will be used for the computation of the regularized t-test statistics.
If obj
=NULL
, a MVR clustering computation for the regularized t-test statistics and/or p-values will be performed.
The function mvrt.test
relies on the R package parallel to create a parallel backend within an R session, enabling access to a cluster
of compute cores and/or nodes on a local and/or remote machine(s) and scaling-up with the number of CPU cores available and efficient parallel
execution. To run a procedure in parallel (with parallel RNG), argument parallel
is to be set to TRUE
and argument conf
is to be specified (i.e. non NULL
). Argument conf
uses the options described in function makeCluster
of the R packages
parallel and snow. PRIMsrc supports two types of communication mechanisms between master and worker processes:
'Socket' or 'Message-Passing Interface' ('MPI'). In PRIMsrc, parallel 'Socket' clusters use sockets communication mechanisms only
(no forking) and are therefore available on all platforms, including Windows, while parallel 'MPI' clusters use high-speed interconnects
mechanism in networks of computers (with distributed memory) and are therefore available only in these architectures. A parallel 'MPI'
cluster also requires R package Rmpi to be installed first. Value type
is used to setup a cluster of type 'Socket' ("SOCKET")
or 'MPI' ("MPI"), respectively. Depending on this type, values of spec
are to be used alternatively:
For 'Socket' clusters (conf$type="SOCKET"
), spec
should be a character
vector
naming the hosts on which
to run the job; it can default to a unique local machine, in which case, one may use the unique host name "localhost".
Each host name can potentially be repeated to the number of CPU cores available on the local machine.
It can also be an integer
scalar specifying the number of processes to spawn on the local machine;
or a list of machine specifications (a character value named host specifying the name or address of the host to use).
For 'MPI' clusters (conf$type="MPI"
), spec
should be an integer
scalar
specifying the total number of processes to be spawned across the network of available nodes, counting the workernodes and masternode.
The actual creation of the cluster, its initialization, and closing are all done internally. For more details, see the reference manual of R package snow and examples below.
When random number generation is needed, the creation of separate streams of parallel RNG per node is done internally by distributing the stream states to the nodes. For more details, see the vignette of R package parallel. The use of a seed allows to reproduce the results within the same type of session: the same seed will reproduce the same results within a non-parallel session or within a parallel session, but it will not necessarily give the exact same results (up to sampling variability) between a non-parallelized and parallelized session due to the difference of management of the seed between the two (see parallel RNG and value of returned seed below).
In case p-values are desired (pval
=TRUE
), the use of a cluster is highly recommended.
It is ideal for computing embarassingly parallel tasks such as permutation or bootstrap resamplings.
Note that in case both regularized t-test statistics and p-values are desired, in order to maximize computational efficiency and
avoid multiple configurations (since a cluster can only be configured and used one session at a time, which otherwise would result in a run stop),
the cluster configuration will only be used for the parallel computation of p-values, but not for the MVR clustering computation
of the regularized t-test statistics.
statistic |
|
p.value |
|
seed |
User seed(s) used: |
This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University. This project was partially funded by the National Institutes of Health (P30-CA043703).
End-user function.
"Jean-Eudes Dazard, Ph.D." [email protected]
"Hua Xu, Ph.D." [email protected]
"Alberto Santana, MBA." [email protected]
Maintainer: "Jean-Eudes Dazard, Ph.D." [email protected]
Dazard J-E. and J. S. Rao (2010). "Regularized Variance Estimation and Variance Stabilization of High-Dimensional Data." In JSM Proceedings, Section for High-Dimensional Data Analysis and Variable Selection. Vancouver, BC, Canada: American Statistical Association IMS - JSM, 5295-5309.
Dazard J-E., Hua Xu and J. S. Rao (2011). "R package MVR for Joint Adaptive Mean-Variance Regularization and Variance Stabilization." In JSM Proceedings, Section for Statistical Programmers and Analysts. Miami Beach, FL, USA: American Statistical Association IMS - JSM, 3849-3863.
Dazard J-E. and J. S. Rao (2012). "Joint Adaptive Mean-Variance Regularization and Variance Stabilization of High Dimensional Data." Comput. Statist. Data Anal. 56(7):2317-2333.
makeCluster
(R package parallel)
eBayes
(R package limma) Bayesian Regularized t-test statistic Smyth, 2004
samr
(R package samr) SAM Regularized t-test statistic Tusher et al., 2001, Storey, 2003
matest
(R package maanova) James-Stein shrinkage estimator-based Regularized t-test statistic Cui et al., 2005
ebam
(R package siggenes) Empirical Bayes Regularized z-test statistic Efron, 2001
bayesT
Hierarchical Bayesian Regularized t-test statistic Baldi et al., 2001
#================================================ # Loading the library and its dependencies #================================================ library("MVR") ## Not run: #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # Use help for descriptions #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real ## End(Not run) #================================================ # Regularized t-test statistics (Synthetic dataset) # Multi-Group Assumption # Assuming unequal variance between groups # With option to use prior MVR clustering results # Without computation of p-values # Without cluster usage #================================================ data("Synthetic", package="MVR") nc.min <- 1 nc.max <- 10 probs <- seq(0, 1, 0.01) n <- 10 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("G1", "G2")) mvr.obj <- mvr(data = Synthetic, block = GF, tolog = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) mvrt.obj <- mvrt.test(data = NULL, obj = mvr.obj, block = NULL, pval = FALSE, replace = FALSE, n.resamp = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) ## Not run: #=================================================== # Examples of parallel backend parametrization #=================================================== if (require("parallel")) { cat("'parallel' is attached correctly \n") } else { stop("'parallel' must be attached first \n") } #=================================================== # Ex. #1 - Multicore PC # Running WINDOWS # SOCKET communication cluster # Shared memory parallelization #=================================================== cpus <- parallel::detectCores(logical = TRUE) conf <- list("spec" = rep("localhost", cpus), "type" = "SOCKET", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Ex. #2 - Master node + 3 Worker nodes cluster # All nodes equipped with identical setups of multicores # (8 core CPUs per machine for a total of 32) # SOCKET communication cluster # Distributed memory parallelization #=================================================== masterhost <- Sys.getenv("HOSTNAME") slavehosts <- c("compute-0-0", "compute-0-1", "compute-0-2") nodes <- length(slavehosts) + 1 cpus <- 8 conf <- list("spec" = c(rep(masterhost, cpus), rep(slavehosts, cpus)), "type" = "SOCKET", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Ex. #3 - Enterprise Multinode Cluster w/ multicore/node # Running LINUX with SLURM scheduler # MPI communication cluster # Distributed memory parallelization # Below, variable 'cpus' is the total number of requested # taks (threads/CPUs), which is specified from within a # SLURM script. #================================================== if (require("Rmpi")) { print("Rmpi is loaded correctly \n") } else { stop("Rmpi must be installed first to use MPI\n") } cpus <- as.numeric(Sys.getenv("SLURM_NTASKS")) conf <- list("spec" = cpus, "type" = "MPI", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== data("Real", package="MVR") nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, tolog = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = TRUE, conf = conf, verbose = TRUE, seed = 1234) #=================================================== # Regularized t-test statistics (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups # With option to use prior MVR clustering results # With computation of p-values #=================================================== mvrt.obj <- mvrt.test(data = NULL, obj = mvr.obj, block = NULL, pval = TRUE, replace = FALSE, n.resamp = 100, parallel = TRUE, conf = conf, verbose = TRUE, seed = 1234) ## End(Not run)
#================================================ # Loading the library and its dependencies #================================================ library("MVR") ## Not run: #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # Use help for descriptions #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real ## End(Not run) #================================================ # Regularized t-test statistics (Synthetic dataset) # Multi-Group Assumption # Assuming unequal variance between groups # With option to use prior MVR clustering results # Without computation of p-values # Without cluster usage #================================================ data("Synthetic", package="MVR") nc.min <- 1 nc.max <- 10 probs <- seq(0, 1, 0.01) n <- 10 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("G1", "G2")) mvr.obj <- mvr(data = Synthetic, block = GF, tolog = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) mvrt.obj <- mvrt.test(data = NULL, obj = mvr.obj, block = NULL, pval = FALSE, replace = FALSE, n.resamp = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) ## Not run: #=================================================== # Examples of parallel backend parametrization #=================================================== if (require("parallel")) { cat("'parallel' is attached correctly \n") } else { stop("'parallel' must be attached first \n") } #=================================================== # Ex. #1 - Multicore PC # Running WINDOWS # SOCKET communication cluster # Shared memory parallelization #=================================================== cpus <- parallel::detectCores(logical = TRUE) conf <- list("spec" = rep("localhost", cpus), "type" = "SOCKET", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Ex. #2 - Master node + 3 Worker nodes cluster # All nodes equipped with identical setups of multicores # (8 core CPUs per machine for a total of 32) # SOCKET communication cluster # Distributed memory parallelization #=================================================== masterhost <- Sys.getenv("HOSTNAME") slavehosts <- c("compute-0-0", "compute-0-1", "compute-0-2") nodes <- length(slavehosts) + 1 cpus <- 8 conf <- list("spec" = c(rep(masterhost, cpus), rep(slavehosts, cpus)), "type" = "SOCKET", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Ex. #3 - Enterprise Multinode Cluster w/ multicore/node # Running LINUX with SLURM scheduler # MPI communication cluster # Distributed memory parallelization # Below, variable 'cpus' is the total number of requested # taks (threads/CPUs), which is specified from within a # SLURM script. #================================================== if (require("Rmpi")) { print("Rmpi is loaded correctly \n") } else { stop("Rmpi must be installed first to use MPI\n") } cpus <- as.numeric(Sys.getenv("SLURM_NTASKS")) conf <- list("spec" = cpus, "type" = "MPI", "homo" = TRUE, "verbose" = TRUE, "outfile" = "") #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== data("Real", package="MVR") nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, tolog = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = TRUE, conf = conf, verbose = TRUE, seed = 1234) #=================================================== # Regularized t-test statistics (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups # With option to use prior MVR clustering results # With computation of p-values #=================================================== mvrt.obj <- mvrt.test(data = NULL, obj = mvr.obj, block = NULL, pval = TRUE, replace = FALSE, n.resamp = 100, parallel = TRUE, conf = conf, verbose = TRUE, seed = 1234) ## End(Not run)
Plot comparative Box-Whisker and Heatmap plots of variables across samples check the effectiveness of normalization before and after Mean-Variance Regularization.
normalization.diagnostic(obj, pal, title = "Normalization Diagnostic Plots", device = NULL, file = "Normalization Diagnostic Plots", path = getwd(), horizontal = FALSE, width = 7, height = 8, ...)
normalization.diagnostic(obj, pal, title = "Normalization Diagnostic Plots", device = NULL, file = "Normalization Diagnostic Plots", path = getwd(), horizontal = FALSE, width = 7, height = 8, ...)
obj |
Object of class" |
title |
Title of the plot. Defaults to "Normalization Diagnostic Plots". |
pal |
Color palette. |
device |
Graphic display device in {NULL, "PS", "PDF"}. Defaults to NULL (standard output screen). Currently implemented graphic display devices are "PS" (Postscript) or "PDF" (Portable Document Format). |
file |
File name for output graphic. Defaults to "Normalization Diagnostic Plots". |
path |
Absolute path (without final (back)slash separator). Defaults to working directory path. |
horizontal |
|
width |
|
height |
|
... |
Generic arguments passed to other plotting functions. |
Option file
is used only if device is specified (i.e. non NULL
).
The argument pal
can be any color palette, e.g. as provided by R package RColorBrewer.
None. Displays the plots on the chosen device
.
This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University. This project was partially funded by the National Institutes of Health (P30-CA043703).
End-user function.
"Jean-Eudes Dazard, Ph.D." [email protected]
"Hua Xu, Ph.D." [email protected]
"Alberto Santana, MBA." [email protected]
Maintainer: "Jean-Eudes Dazard, Ph.D." [email protected]
Dazard J-E. and J. S. Rao (2010). "Regularized Variance Estimation and Variance Stabilization of High-Dimensional Data." In JSM Proceedings, Section for High-Dimensional Data Analysis and Variable Selection. Vancouver, BC, Canada: American Statistical Association IMS - JSM, 5295-5309.
Dazard J-E., Hua Xu and J. S. Rao (2011). "R package MVR for Joint Adaptive Mean-Variance Regularization and Variance Stabilization." In JSM Proceedings, Section for Statistical Programmers and Analysts. Miami Beach, FL, USA: American Statistical Association IMS - JSM, 3849-3863.
Dazard J-E. and J. S. Rao (2012). "Joint Adaptive Mean-Variance Regularization and Variance Stabilization of High Dimensional Data." Comput. Statist. Data Anal. 56(7):2317-2333.
justvsn
(R package vsn) Variance stabilization and calibration for microarray data.
loess
(R package stats) Fit a polynomial surface determined by one or more numerical predictors, using local fitting.
## Not run: #=================================================== # Loading the library and its dependencies #=================================================== library("MVR") library("RColorBrewer") #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # (see description of datasets) #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups # Without cluster usage #=================================================== nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, log = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) #=================================================== # Summary Normalization Diagnostic Plots (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== normalization.diagnostic(obj = mvr.obj, title = "Normalization Diagnostic Plots (Real - Multi-Group Assumption)", pal = brewer.pal(n=11, name="RdYlGn"), device = NULL, horizontal = FALSE, width = 7, height = 8) ## End(Not run)
## Not run: #=================================================== # Loading the library and its dependencies #=================================================== library("MVR") library("RColorBrewer") #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # (see description of datasets) #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups # Without cluster usage #=================================================== nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, log = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) #=================================================== # Summary Normalization Diagnostic Plots (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== normalization.diagnostic(obj = mvr.obj, title = "Normalization Diagnostic Plots (Real - Multi-Group Assumption)", pal = brewer.pal(n=11, name="RdYlGn"), device = NULL, horizontal = FALSE, width = 7, height = 8) ## End(Not run)
The dataset comes from a quantitative Liquid Chromatography/Mass-Spectrometry (LC/MS)
shotgun (bottom-up) proteomics experiment. It consists of independent cell cultures
of human of Myeloid Dendritic Cells (MDCs) from normal subjects. Samples were split
into a control ("M") and a treated group ("S"), stimulated with either media alone
or a Toll-Like receptor-3 Ligand respectively. The goal was to identify differentially expressed
peptides (or proteins) between the two groups involved in the immune response
of human MDCs upon TLR-3 Ligand binding.
The dataset is assumed to have been pre-processed for non-ignorable missing values,
leaving a complete dataset with unique peptides or predictor variables.
This is a balanced design with two sample groups (), under unequal sample group variance.
data("Real", package="MVR")
data("Real", package="MVR")
A numeric matrix containing observations (samples) by rows and
variables by columns, named after peptide names (
).
Samples are balanced (
,
) between the two groups ("M", "S").
Compressed Rda data file.
This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University. This project was partially funded by the National Institutes of Health (P30-CA043703).
"Jean-Eudes Dazard, Ph.D." [email protected]
"Hua Xu, Ph.D." [email protected]
"Alberto Santana, MBA." [email protected]
Maintainer: "Jean-Eudes Dazard, Ph.D." [email protected]
See real proteomics data application in Dazard et al., 2011, 2012.
Dazard J-E. and J. S. Rao (2010). "Regularized Variance Estimation and Variance Stabilization of High-Dimensional Data." In JSM Proceedings, Section for High-Dimensional Data Analysis and Variable Selection. Vancouver, BC, Canada: American Statistical Association IMS - JSM, 5295-5309.
Dazard J-E., Hua Xu and J. S. Rao (2011). "R package MVR for Joint Adaptive Mean-Variance Regularization and Variance Stabilization." In JSM Proceedings, Section for Statistical Programmers and Analysts. Miami Beach, FL, USA: American Statistical Association IMS - JSM, 3849-3863.
Dazard J-E. and J. S. Rao (2012). "Joint Adaptive Mean-Variance Regularization and Variance Stabilization of High Dimensional Data." Comput. Statist. Data Anal. 56(7):2317-2333.
Plot comparative variance-mean plots to check the variance stabilization across variables before and after Mean-Variance Regularization.
stabilization.diagnostic(obj, span = 0.5, degree = 2, family = "gaussian", title = "Stabilization Diagnostic Plots", device = NULL, file = "Stabilization Diagnostic Plots", path = getwd(), horizontal = FALSE, width = 7, height = 5, ...)
stabilization.diagnostic(obj, span = 0.5, degree = 2, family = "gaussian", title = "Stabilization Diagnostic Plots", device = NULL, file = "Stabilization Diagnostic Plots", path = getwd(), horizontal = FALSE, width = 7, height = 5, ...)
obj |
Object of class " |
title |
Title of the plot. Defaults to "Stabilization Diagnostic Plots". |
span |
Span parameter of the |
degree |
Degree parameter of the |
family |
Family distribution in "gaussian", "symmetric" of the |
device |
Graphic display device in {NULL, "PS", "PDF"}. Defaults to NULL (standard output screen). Currently implemented graphic display devices are "PS" (Postscript) or "PDF" (Portable Document Format). |
file |
File name for output graphic. Defaults to "Stabilization Diagnostic Plots". |
path |
Absolute path (without final (back)slash separator). Defaults to working directory path. |
horizontal |
|
width |
|
height |
|
... |
Generic arguments passed to other plotting functions. |
In the plots of standard deviations vs. means, standard deviations and means are calculated in a feature-wise manner from the expression matrix. The scatterplot allows to visually verify whether there is a dependence of the standard deviation (or variance) on the mean. The black dotted line depicts the LOESS scatterplot smoother estimator. If there is no variance-mean dependence, then this line should be approximately horizontal.
Option file
is used only if device is specified (i.e. non NULL
).
None. Displays the plots on the chosen device
.
This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University. This project was partially funded by the National Institutes of Health (P30-CA043703).
End-user function.
"Jean-Eudes Dazard, Ph.D." [email protected]
"Hua Xu, Ph.D." [email protected]
"Alberto Santana, MBA." [email protected]
Maintainer: "Jean-Eudes Dazard, Ph.D." [email protected]
Dazard J-E. and J. S. Rao (2010). "Regularized Variance Estimation and Variance Stabilization of High-Dimensional Data." In JSM Proceedings, Section for High-Dimensional Data Analysis and Variable Selection. Vancouver, BC, Canada: American Statistical Association IMS - JSM, 5295-5309.
Dazard J-E., Hua Xu and J. S. Rao (2011). "R package MVR for Joint Adaptive Mean-Variance Regularization and Variance Stabilization." In JSM Proceedings, Section for Statistical Programmers and Analysts. Miami Beach, FL, USA: American Statistical Association IMS - JSM, 3849-3863.
Dazard J-E. and J. S. Rao (2012). "Joint Adaptive Mean-Variance Regularization and Variance Stabilization of High Dimensional Data." Comput. Statist. Data Anal. 56(7):2317-2333.
justvsn
(R package vsn) Variance stabilization and calibration for microarray data.
loess
(R package stats) Fit a polynomial surface determined by one or more numerical predictors, using local fitting.
## Not run: #=================================================== # Loading the library and its dependencies #=================================================== library("MVR") #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # (see description of datasets) #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups # Without cluster usage #=================================================== nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, log = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) #=================================================== # Summary Stabilization Diagnostic Plots (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== stabilization.diagnostic(obj = mvr.obj, title = "Stabilization Diagnostic Plots (Real - Multi-Group Assumption)", span = 0.75, degree = 2, family = "gaussian", device = NULL, horizontal = FALSE, width = 7, height = 5) ## End(Not run)
## Not run: #=================================================== # Loading the library and its dependencies #=================================================== library("MVR") #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # (see description of datasets) #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups # Without cluster usage #=================================================== nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, log = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) #=================================================== # Summary Stabilization Diagnostic Plots (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== stabilization.diagnostic(obj = mvr.obj, title = "Stabilization Diagnostic Plots (Real - Multi-Group Assumption)", span = 0.75, degree = 2, family = "gaussian", device = NULL, horizontal = FALSE, width = 7, height = 5) ## End(Not run)
Generation of a synthetic dataset with n=10 observations (samples) and variables,
where
of them are significantly different between the two sample groups.
This is a balanced design with two sample groups (), under unequal sample group variance.
data("Synthetic", package="MVR")
data("Synthetic", package="MVR")
A numeric matrix containing observations (samples) by rows
and
variables by columns, named
.
Samples are balanced (
,
) between the two groups (
).
Compressed Rda data file.
This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University. This project was partially funded by the National Institutes of Health (P30-CA043703).
"Jean-Eudes Dazard, Ph.D." [email protected]
"Hua Xu, Ph.D." [email protected]
"Alberto Santana, MBA." [email protected]
Maintainer: "Jean-Eudes Dazard, Ph.D." [email protected]
See model #2 in Dazard et al., 2011, 2012.
Dazard J-E. and J. S. Rao (2010). "Regularized Variance Estimation and Variance Stabilization of High-Dimensional Data." In JSM Proceedings, Section for High-Dimensional Data Analysis and Variable Selection. Vancouver, BC, Canada: American Statistical Association IMS - JSM, 5295-5309.
Dazard J-E., Hua Xu and J. S. Rao (2011). "R package MVR for Joint Adaptive Mean-Variance Regularization and Variance Stabilization." In JSM Proceedings, Section for Statistical Programmers and Analysts. Miami Beach, FL, USA: American Statistical Association IMS - JSM, 3849-3863.
Dazard J-E. and J. S. Rao (2012). "Joint Adaptive Mean-Variance Regularization and Variance Stabilization of High Dimensional Data." Comput. Statist. Data Anal. 56(7):2317-2333.
Plot comparative distribution densities of means and standard deviations of the data before and after Mean-Variance Regularization to check for location shifts between observed first and second moments and their expected target values under a target centered homoscedastic model.
Plot comparative QQ scatterplots to look at departures between observed distributions of first and second moments of the MVR-transformed data and their theoretical distributions assuming independence and normality of all the variables.
target.diagnostic(obj, title = "Target Moments Diagnostic Plots", device = NULL, file = "Target Moments Diagnostic Plots", path = getwd(), horizontal = FALSE, width = 8.5, height = 6.5, ...)
target.diagnostic(obj, title = "Target Moments Diagnostic Plots", device = NULL, file = "Target Moments Diagnostic Plots", path = getwd(), horizontal = FALSE, width = 8.5, height = 6.5, ...)
obj |
Object of class " |
title |
Title of the plot. Defaults to "Target Moments Diagnostic Plots". |
device |
Graphic display device in {NULL, "PS", "PDF"}. Defaults to NULL (standard output screen). Currently implemented graphic display devices are "PS" (Postscript) or "PDF" (Portable Document Format). |
file |
File name for output graphic. Defaults to "Target Moments Diagnostic Plots". |
path |
Absolute path (without final (back)slash separator). Defaults to working directory path. |
horizontal |
|
width |
|
height |
|
... |
Generic arguments passed to other plotting functions. |
The plots of the density distribution of means and standard deviations checks that the distributions of means and standard deviations
of the MVR-transformed data have correct target first moments, i.e. with mean ~ 0 and mean ~ 1. The expected target mean and standard
deviation are shown in red (before and) after MVR-transformation. Caption shows the p-values from the parametric two-sample two-sided t-tests
for the equality of parameters to their expectations (assuming normality since usually sample sizes are large : ,
or a relative robustness to moderate violations of the normality assumption).
In the general case, the variables are not normally distributed and not even independent and identically distributed before
and after MVR-transformation. Therefore, the distributions of untransformed first and second moments usually differ from their respective
theoretical null distributions, i.e., from for the means and from
for the standard
deviations, where
denotes the number of sample groups (see Dazard, J-E. and J. S. Rao (2012) for more details).
Also, the observed distributions of transformed first and second moments are unknown. This is reflected in the QQ plots,
where theoretical and empirical quantiles do not necessarily align with each other. Caption shows the p-values from the nonparametric two-sample
two-sided Kolmogorov-Smirnov tests of the null hypothesis that a parameter distribution differs from its theoretical distribution.
Each black dot represents a variable. The red solid line depicts the interquartile line, which passes through the first and third quartiles.
Option file
is used only if device is specified (i.e. non NULL
).
None. Displays the plots on the chosen device
.
This work made use of the High Performance Computing Resource in the Core Facility for Advanced Research Computing at Case Western Reserve University. This project was partially funded by the National Institutes of Health (P30-CA043703).
End-user function.
"Jean-Eudes Dazard, Ph.D." [email protected]
"Hua Xu, Ph.D." [email protected]
"Alberto Santana, MBA." [email protected]
Maintainer: "Jean-Eudes Dazard, Ph.D." [email protected]
Dazard J-E. and J. S. Rao (2010). "Regularized Variance Estimation and Variance Stabilization of High-Dimensional Data." In JSM Proceedings, Section for High-Dimensional Data Analysis and Variable Selection. Vancouver, BC, Canada: American Statistical Association IMS - JSM, 5295-5309.
Dazard J-E., Hua Xu and J. S. Rao (2011). "R package MVR for Joint Adaptive Mean-Variance Regularization and Variance Stabilization." In JSM Proceedings, Section for Statistical Programmers and Analysts. Miami Beach, FL, USA: American Statistical Association IMS - JSM, 3849-3863.
Dazard J-E. and J. S. Rao (2012). "Joint Adaptive Mean-Variance Regularization and Variance Stabilization of High Dimensional Data." Comput. Statist. Data Anal. 56(7):2317-2333.
justvsn
(R package vsn) Variance stabilization and calibration for microarray data.
loess
(R package stats) Fit a polynomial surface determined by one or more numerical predictors, using local fitting.
## Not run: #=================================================== # Loading the library and its dependencies #=================================================== library("MVR") library("RColorBrewer") #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # (see description of datasets) #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups # Without cluster usage #=================================================== nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, log = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) #=================================================== # Summary Target Moments Diagnostic Plots (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== target.diagnostic(obj = mvr.obj, title = "Target Moments Diagnostic Plots (Real - Multi-Group Assumption)", device = NULL, horizontal = FALSE, width = 8.5, height = 6.5) ## End(Not run)
## Not run: #=================================================== # Loading the library and its dependencies #=================================================== library("MVR") library("RColorBrewer") #=================================================== # MVR package news #=================================================== MVR.news() #================================================ # MVR package citation #================================================ citation("MVR") #=================================================== # Loading of the Synthetic and Real datasets # (see description of datasets) #=================================================== data("Synthetic", "Real", package="MVR") ?Synthetic ?Real #=================================================== # Mean-Variance Regularization (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups # Without cluster usage #=================================================== nc.min <- 1 nc.max <- 30 probs <- seq(0, 1, 0.01) n <- 6 GF <- factor(gl(n = 2, k = n/2, length = n), ordered = FALSE, labels = c("M", "S")) mvr.obj <- mvr(data = Real, block = GF, log = FALSE, nc.min = nc.min, nc.max = nc.max, probs = probs, B = 100, parallel = FALSE, conf = NULL, verbose = TRUE, seed = 1234) #=================================================== # Summary Target Moments Diagnostic Plots (Real dataset) # Multi-Group Assumption # Assuming unequal variance between groups #=================================================== target.diagnostic(obj = mvr.obj, title = "Target Moments Diagnostic Plots (Real - Multi-Group Assumption)", device = NULL, horizontal = FALSE, width = 8.5, height = 6.5) ## End(Not run)