Commit 924f9e93 by Eric Coissac

Cleaned version with new version of simulation

parent 31363933
...@@ -30,7 +30,6 @@ Collate: ...@@ -30,7 +30,6 @@ Collate:
'mprocuste.R' 'mprocuste.R'
'anova.pm.R' 'anova.pm.R'
'covls.R' 'covls.R'
'mcov.R'
'plot.pm.R' 'plot.pm.R'
'procuste.R' 'procuste.R'
'simnull.R' 'simnull.R'
......
...@@ -47,7 +47,6 @@ export(protate) ...@@ -47,7 +47,6 @@ export(protate)
export(simulate_correlation) export(simulate_correlation)
export(simulate_matrix) export(simulate_matrix)
export(varls) export(varls)
export(varls2)
export(vars.procmod) export(vars.procmod)
export(weighted.residuals) export(weighted.residuals)
import(MASS) import(MASS)
......
#' @include procmod.frame.R
#' @include multivariate.R
#'
#' @import expm
#'
#' @author Christelle Gonindard-Melodelima
#' @author Eric Coissac
NULL
#' Compute the trace of a square matrix.
#'
#'
#' @param X a square matrix
#' @return the trace of X
#'
#' @examples
#' m = matrix(1:16,nrow=4)
#' ProcMod:::.Trace(m)
#'
#' @note Internal function do not use.
#'
#' @rdname internal.getPermuteMatrix
#' @author Eric Coissac
#' @author Christelle Gonindard-Melodelima
#'
.Trace = function(X) sum(diag(X))
#' Compute the variance, covariance matrix of K coordinate matrices.
#'
#' Covariance between two matrices is defined as the sum of the
#' sigular values of the X'Y matrix. All the matrices must have
#' the same number of rows.
#'
#' @param ... the set of matrices
#'
#' @examples
#' # Build Three matrices of 3 rows.
#' A <- matrix(1:9,nrow=3)
#' B <- matrix(10:15,nrow=3)
#' C <- matrix(20:31,nrow=3)
#' # compute the variance covariance matrix
#' mvar(A,B,C)
#' mvar(A=A,B=B,C=C)
#'
#' @author Eric Coissac
#' @author Christelle Gonindard-Melodelima
#' @export
varls = function(...,permutations = how(nperm = 999)) {
Xs <- list(...)
if (length(Xs)==1) {
x = Xs[[1]]
if (is.procmod.frame(x))
Xs=x
else if (is.pm(x))
return(x$cov)
else
Xs=procmod.frame(x)
}
else
Xs=as.procmod.frame(Xs)
Xnames=names(Xs)
Xs <- ortho(Xs)
XXs = as.procmod.frame(mapply(tcrossprod, Xs,
SIMPLIFY = FALSE))
nX = length(Xs)
N = nrow(Xs)-1
if (! is.null(permutations)) {
pmatrix = .getPermuteMatrix(perm=permutations,N=nrow(Xs))
rCovXXs = matrix(0, nrow = nX * 2, ncol = nX * 2)
pval = matrix(0, nrow = nX * 2, ncol = nX * 2)
for (i in 1:(2*nX))
for (j in 1:(2*nX)) {
if (i %% 2 && j %% 2) {
rCovXXs[i,j]=.Trace(sqrtm(XXs[[ceiling(i / 2)]] %*% XXs[[ceiling(j/2)]])) / N
pval[i,j]=-1
}
else if (i ==j) {
vv = numeric(nrow(pmatrix))
for (k in seq_len(nrow(pmatrix))) {
d = Xs[[ceiling(i / 2)]][pmatrix[k,],]
dd = tcrossprod(d)
vv[k] = Re(.Trace(sqrtm(dd %*% dd)))
}
rCovXXs[i,j]=mean(vv) / N
pval[i,j]=shapiro.test(vv)$p.value
}
else if (i <=j) {
vv = numeric(nrow(pmatrix))
for (k in seq_len(nrow(pmatrix))) {
d = Xs[[ceiling(i / 2)]][pmatrix[k,],]
dd = tcrossprod(d)
vv[k] = Re(.Trace(sqrtm(dd %*% XXs[[ceiling(j/2)]])))
}
rCovXXs[i,j]=mean(vv) / N
rCovXXs[j,i]=rCovXXs[i,j]
pval[i,j]=shapiro.test(vv)$p.value
pval[j,i]=shapiro.test(vv)$p.value
}
}
CovXXs=rCovXXs
Xnames = rep(Xnames,rep(2,length(Xnames)))
Xsuff = rep("",length(Xnames))
Xsuff[seq(from=2,to=length(Xnames),by=2)]="r_"
Xnames = mapply(paste0, Xsuff,Xnames,collapse="")
colnames(CovXXs)=Xnames
rownames(CovXXs)=Xnames
colnames(pval)=Xnames
rownames(pval)=Xnames
print(pval)
}
else {
Xx <- rep(1:nX,nX)
Xy <- rep(1:nX,rep(nX,nX))
CovXXs <- mapply(function(x,y) .Trace(sqrtm(XXs[[x]] %*% XXs[[y]])),
Xx,Xy,
SIMPLIFY = TRUE) / N
dim(CovXXs)=c(nX,nX)
colnames(CovXXs)=Xnames
rownames(CovXXs)=Xnames
}
return(Re(CovXXs))
}
#' Compute the variance, covariance matrix of K coordinate matrices.
#'
#' Covariance between two matrices is defined as the sum of the
#' sigular values of the X'Y matrix. All the matrices must have
#' the same number of rows.
#'
#' @param ... the set of matrices
#'
#' @examples
#' # Build Three matrices of 3 rows.
#' A <- matrix(1:9,nrow=3)
#' B <- matrix(10:15,nrow=3)
#' C <- matrix(20:31,nrow=3)
#' # compute the variance covariance matrix
#' varls2(A,B,C)
#' varls2(A=A,B=B,C=C)
#'
#' @author Eric Coissac
#' @author Christelle Gonindard-Melodelima
#' @export
varls2 = function(...,permutations = how(nperm = 999)) {
Xs <- list(...)
if (length(Xs)==1) {
x = Xs[[1]]
if (is.procmod.frame(x))
Xs=x
else if (is.pm(x))
return(x$cov)
else
Xs=procmod.frame(x)
}
else
Xs=as.procmod.frame(Xs)
Xnames=names(Xs)
Xs <- ortho(Xs)
XXs = as.procmod.frame(mapply(tcrossprod, Xs,
SIMPLIFY = FALSE))
nX = length(Xs)
N = nrow(Xs)-1
CovXXs = matrix(0, nrow = nX, ncol = nX)
for (i in seq_len(nX))
for (j in i:nX) {
CovXXs[i,j] = .Trace(sqrtm(XXs[[i]] %*% XXs[[j]]))
}
if (! is.null(permutations)) {
pmatrix = .getPermuteMatrix(perm=permutations,N=nrow(Xs))
nP = nrow(pmatrix)
rCovXXs = array(0,dim=c(nX,nX,nP))
for (k in seq_len(nP)) {
Xp = Xs[pmatrix[k,],]
for (i in seq_len(nX)) {
dd = tcrossprod(Xp[[i]])
for (j in i:nX)
rCovXXs[i,j,k] = Re(.Trace(sqrtm(dd %*% XXs[[j]])))
}
}
for (i in seq_len(nX))
for (j in i:nX)
CovXXs[i,j] = CovXXs[i,j] - mean(rCovXXs[i,j,])
}
for (i in seq_len(nX))
for (j in i:nX)
CovXXs[j,i] = CovXXs[i,j]
CovXXs = CovXXs / N
colnames(CovXXs)=Xnames
rownames(CovXXs)=Xnames
return(Re(CovXXs))
}
#' Compute the person correlation matrix of K coordinate matrices
#'
#' @author Eric Coissac
#' @author Christelle Gonindard-Melodelima
#' @export
corls = function(...,permutations = how(nperm = 999)) {
cov = varls(...,permutations = permutations)
s = sqrt(diag(cov))
vv= outer(s,s)
return(cov/vv)
}
#' Compute the person partial correlation matrix of K coordinate matrices
#'
#' @author Eric Coissac
#' @author Christelle Gonindard-Melodelima
#' @export
corls.partial = function(...,permutations = how(nperm = 999)) {
C = solve(corls(...,permutations = permutations))
D = sqrt(diag(C) %o% diag(C))
return(C/D)
}
% Generated by roxygen2: do not edit by hand % Generated by roxygen2: do not edit by hand
% Please edit documentation in R/covls.R, R/mcov.R % Please edit documentation in R/covls.R
\name{corls} \name{corls}
\alias{corls} \alias{corls}
\title{Compute the person correlation matrix of K coordinate matrices} \title{Compute the person correlation matrix of K coordinate matrices}
\usage{ \usage{
corls(..., permutations = how(nperm = 999)) corls(..., nperm = 100, rcovls = FALSE)
corls(..., permutations = how(nperm = 999))
} }
\description{ \description{
Compute the person correlation matrix of K coordinate matrices Compute the person correlation matrix of K coordinate matrices
Compute the person correlation matrix of K coordinate matrices
} }
\author{ \author{
Eric Coissac Eric Coissac
Christelle Gonindard-Melodelima Christelle Gonindard-Melodelima
Eric Coissac
Christelle Gonindard-Melodelima
} }
% Generated by roxygen2: do not edit by hand % Generated by roxygen2: do not edit by hand
% Please edit documentation in R/covls.R, R/mcov.R % Please edit documentation in R/covls.R
\name{corls.partial} \name{corls.partial}
\alias{corls.partial} \alias{corls.partial}
\title{Compute the person partial correlation matrix of K coordinate matrices} \title{Compute the person partial correlation matrix of K coordinate matrices}
\usage{ \usage{
corls.partial(..., permutations = how(nperm = 999)) corls.partial(..., nperm = 100, rcovls = FALSE)
corls.partial(..., permutations = how(nperm = 999))
} }
\description{ \description{
Compute the person partial correlation matrix of K coordinate matrices Compute the person partial correlation matrix of K coordinate matrices
Compute the person partial correlation matrix of K coordinate matrices
} }
\author{ \author{
Eric Coissac Eric Coissac
Christelle Gonindard-Melodelima Christelle Gonindard-Melodelima
Eric Coissac
Christelle Gonindard-Melodelima
} }
% Generated by roxygen2: do not edit by hand % Generated by roxygen2: do not edit by hand
% Please edit documentation in R/covls.R, R/mcov.R % Please edit documentation in R/covls.R
\name{.Trace} \name{.Trace}
\alias{.Trace} \alias{.Trace}
\title{Compute the trace of a square matrix.} \title{Compute the trace of a square matrix.}
\usage{ \usage{
.Trace(X) .Trace(X)
.Trace(X)
} }
\arguments{ \arguments{
\item{X}{a square matrix} \item{X}{a square matrix}
\item{X}{a square matrix}
} }
\value{ \value{
the trace of X the trace of X
the trace of X
} }
\description{ \description{
Compute the trace of a square matrix. Compute the trace of a square matrix.
Compute the trace of a square matrix.
} }
\note{ \note{
Internal function do not use. Internal function do not use.
Internal function do not use.
} }
\examples{ \examples{
m = matrix(1:16,nrow=4) m = matrix(1:16,nrow=4)
ProcMod:::.Trace(m) ProcMod:::.Trace(m)
m = matrix(1:16,nrow=4)
ProcMod:::.Trace(m)
} }
\author{ \author{
Eric Coissac Eric Coissac
Christelle Gonindard-Melodelima Christelle Gonindard-Melodelima
Eric Coissac
Christelle Gonindard-Melodelima
} }
% Generated by roxygen2: do not edit by hand % Generated by roxygen2: do not edit by hand
% Please edit documentation in R/covls.R, R/mcov.R % Please edit documentation in R/covls.R
\name{varls} \name{varls}
\alias{varls} \alias{varls}
\title{Compute the variance, covariance matrix of K coordinate matrices.} \title{Compute the variance, covariance matrix of K coordinate matrices.}
\usage{ \usage{
varls(..., permutations = how(nperm = 999)) varls(..., nperm = 100, rcovls = FALSE)
varls(..., permutations = how(nperm = 999))
} }
\arguments{ \arguments{
\item{...}{the set of matrices} \item{...}{the set of matrices}
\item{rcovls}{} \item{rcovls}{}
\item{...}{the set of matrices}
} }
\description{ \description{
Covariance between two matrices is defined as the sum of the Covariance between two matrices is defined as the sum of the
sigular values of the X'Y matrix. All the matrices must have sigular values of the X'Y matrix. All the matrices must have
the same number of rows. the same number of rows.
Covariance between two matrices is defined as the sum of the
sigular values of the X'Y matrix. All the matrices must have
the same number of rows.
} }
\examples{ \examples{
# Build Three matrices of 3 rows. # Build Three matrices of 3 rows.
...@@ -33,21 +25,9 @@ the same number of rows. ...@@ -33,21 +25,9 @@ the same number of rows.
varls2(A,B,C) varls2(A,B,C)
varls2(A=A,B=B,C=C) varls2(A=A,B=B,C=C)
# Build Three matrices of 3 rows.
A <- matrix(1:9,nrow=3)
B <- matrix(10:15,nrow=3)
C <- matrix(20:31,nrow=3)
# compute the variance covariance matrix
mvar(A,B,C)
mvar(A=A,B=B,C=C)
} }
\author{ \author{
Eric Coissac Eric Coissac
Christelle Gonindard-Melodelima Christelle Gonindard-Melodelima
Eric Coissac
Christelle Gonindard-Melodelima
} }
% Generated by roxygen2: do not edit by hand
% Please edit documentation in R/mcov.R
\name{varls2}
\alias{varls2}
\title{Compute the variance, covariance matrix of K coordinate matrices.}
\usage{
varls2(..., permutations = how(nperm = 999))
}
\arguments{
\item{...}{the set of matrices}
}
\description{
Covariance between two matrices is defined as the sum of the
sigular values of the X'Y matrix. All the matrices must have
the same number of rows.
}
\examples{
# Build Three matrices of 3 rows.
A <- matrix(1:9,nrow=3)
B <- matrix(10:15,nrow=3)
C <- matrix(20:31,nrow=3)
# compute the variance covariance matrix
varls2(A,B,C)
varls2(A=A,B=B,C=C)
}
\author{
Eric Coissac
Christelle Gonindard-Melodelima
}
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment