The "MVN" Class
MVN class is a tool used to generate gene expressions that
follow multivariate normal distribution.
1 2 3 4 5 6 7 8 9 10 11 12 13
MVN(mu, Sigma, tol = 1e-06) covar(object) correl(object) ## S4 method for signature 'MVN' alterMean(object, TRANSFORM, ...) ## S4 method for signature 'MVN' alterSD(object, TRANSFORM, ...) ## S4 method for signature 'MVN' nrow(x) ## S4 method for signature 'MVN' rand(object, n, ...) ## S4 method for signature 'MVN' summary(object, ...)
numeric vector representing k-dimensional means
numeric k-by-k covariance matrix containing the measurement of the linear coupling between every pair of random vectors
numeric scalar roundoff error that will be tolerated when assessing the singularity of the covariance matrix
object of class
function that takes a vector of mean expression or standard deviation and returns a transformed vector that can be used to alter the appropriate slot of the object.
numeric scalar representing number of samples to be simulated
extra arguments for generic or plotting routines
The implementation of
MVN class is designed for efficiency when
generating new samples, since we expect to do this several times.
Basically, this class separates the
mvrnorm function from the
MASS package into several steps. The computationally expensive step
(when the dimension is large) is the eigenvector decomposition of the
covariance matrix. This step is performed at construction and the
pieces are stored in the object. The
rand method for
objects contains the second half of the
Note that we typically work on expression value with its logarithm to some appropriate base. That is, the multivariate normal should be used on the logarithmic scale in order to contruct engine.
alterMean for an
MVN simply replaces the appropriate slot by
the transformed vector.
alterSD for an
MVN is trickier,
because of the way the data is stored. In order to have some hope of getting
this correct, we work in the space of the covariance matrix, Sigma.
If we let R denote the correlation matrix and let Delta be the
diagonal matrix whose entries are the individual standard deviations, then
Sigma = Delta \%*\% R \%*\% Delta.
So, we can change the standard deviations by replacing Delta in this
product. We then construct a new
MVN object with the old mean vector
and the new covariance matrix.
correl functions calculate the covariance matrix
and correlation matrix underlies the covariance matrix for the objects of
MVN, respectively. We have four assertions as shown below,
and will be tested in the examples section:
covarshould return the same matrix that was used in the function call to construct the
After applying an
alterMeanmethod, the covariance matrix is unchanged.
The diagonal of correlation matrix consists of all ones.
After applying an
alterSDmethod, the correlation matrix is unchanged.
Objects from the Class
Although objects of the class can be created by a direct call to
new, the preferred method is to use the
MVN generator function.
numeric vector containing the k-dimentional means
numeric vector containing the square roots of eigenvalues of the covariance matrix
numeric matrix with k*k dimensions whose columns contain the eigenvectors of the covariance matrix
- alterMean(object, TRANSFORM, ...)
Takes an object of class
MVN, loops over the
muslot, alters the mean as defined by TRANSFORM function, and returns an object of class
- alterSD(object, TRANSFORM, ...)
Takes an object of class
MVN, works on the diagonal matrix of the covariance matrix, alters the standard deviation as defined by TRANSFORM function, and reconstructs an object of class
MVNwith the old
muand reconstructed covariance matrix.
Returns the number of genes (i.e, the length of the
- rand(object, n, ...)
Generates nrow(MVN)*n matrix representing gene expressions of
nsamples following the multivariate normal distribution captured in the object of
- summary(object, ...)
Prints out the number of multivariate normal random variables in the object of
Returns the covariance matrix of the object of class
Returns the correlation matrix of the object of class
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28
showClass("MVN") ## Not run: tolerance <- 1e-10 ## Create a random orthogonal 2x2 matrix a <- runif(1) b <- sqrt(1-a^2) X <- matrix(c(a, b, -b, a), 2, 2) ## Now choose random positive squared-eigenvalues Lambda2 <- diag(rev(sort(rexp(2))), 2) ## Construct a covariance matrix Y <- t(X) ## Use the MVN constructor marvin <- MVN(c(0,0), Y) ## Check the four assertions print(paste('Tolerance for assertion checking:', tolerance)) print(paste('Covar assertion 1:', all(abs(covar(marvin) - Y) < tolerance))) mar2 <- alterMean(marvin, normalOffset, delta=3) print(paste('Covar assertion 2:', all(abs(covar(marvin) - covar(mar2)) < tolerance))) print(paste('Correl assertion 1:', all(abs(diag(correl(marvin)) - 1) < tolerance))) mar3 <- alterSD(marvin, function(x) 2*x) print(paste('Correl assertion 1:', all(abs(correl(marvin) - correl(mar2)) < tolerance))) rm(a, b, X, Lambda2, Y, marvin, mar2, mar3) ## End(Not run)