View source: R/vc_score_h_perm.R
vc_score_h_perm | R Documentation |
This function computes an approximation of the variance component test for homogeneous trajectory
based on the asymptotic distribution of a mixture of χ^{2}s using Davies method
from davies
vc_score_h_perm( y, x, indiv, phi, w, Sigma_xi = diag(ncol(phi)), na_rm = FALSE, n_perm = 1000, progressbar = TRUE, parallel_comp = TRUE, nb_cores = parallel::detectCores() - 1 )
y |
a numeric matrix of dim |
x |
a numeric design matrix of dim |
indiv |
a vector of length |
phi |
a numeric design matrix of size |
w |
a vector of length |
Sigma_xi |
a matrix of size |
na_rm |
logical: should missing values (including |
n_perm |
the number of permutation to perform. Default is |
progressbar |
logical indicating whether a progress bar should be displayed when computing permutations (only in interactive mode). |
parallel_comp |
a logical flag indicating whether parallel computation
should be enabled. Only Linux and MacOS are supported, this is ignored on Windows.
Default is |
nb_cores |
an integer indicating the number of cores to be used when
|
A list with the following elements:
score
: an approximation of the observed set score
scores_perm
: a vector containing the permuted set scores
gene_scores_unscaled
: approximation of the individual gene scores
gene_scores_unscaled_perm
: a list of approximation of the permuted individual gene scores
davies
rm(list=ls()) set.seed(123) ##generate some fake data ######################## ng <- 100 nindiv <- 30 nt <- 5 nsample <- nindiv*nt tim <- matrix(rep(1:nt), nindiv, ncol=1, nrow=nsample) tim2 <- tim^2 sigma <- 5 b0 <- 10 #under the null: beta1 <- rnorm(n=ng, 0, sd=0) #under the (heterogen) alternative: beta1 <- rnorm(n=ng, 0, sd=0.1) #under the (homogen) alternative: beta1 <- rnorm(n=ng, 0.06, sd=0) y.tilde <- b0 + rnorm(ng, sd = sigma) y <- t(matrix(rep(y.tilde, nsample), ncol=ng, nrow=nsample, byrow=TRUE) + matrix(rep(beta1, each=nsample), ncol=ng, nrow=nsample, byrow=FALSE) * matrix(rep(tim, ng), ncol=ng, nrow=nsample, byrow=FALSE) + #matrix(rep(beta1, each=nsample), ncol=ng, nrow=nsample, byrow=FALSE) * # matrix(rep(tim2, ng), ncol=ng, nrow=nsample, byrow=FALSE) + matrix(rnorm(ng*nsample, sd = sigma), ncol=ng, nrow=nsample, byrow=FALSE) ) myindiv <- rep(1:nindiv, each=nt) x <- cbind(1, myindiv/2==floor(myindiv/2)) myw <- matrix(rnorm(nsample*ng, sd=0.1), ncol=nsample, nrow=ng) #run test #We only use few permutations (10) to keep example running time low #Otherwise one can use n_perm = 1000 score_homogen <- vc_score_h_perm(y, x, phi=tim, indiv=myindiv, w=myw, Sigma_xi=cov(tim), n_perm = 10, parallel_comp = FALSE) score_homogen$score score_heterogen <- vc_score_perm(y, x, phi=tim, indiv=myindiv, w=myw, Sigma_xi=cov(tim), n_perm = 10, parallel_comp = FALSE) score_heterogen$score scoreTest_homogen <- vc_test_asym(y, x, phi=tim, indiv=rep(1:nindiv, each=nt), w=matrix(1, ncol=ncol(y), nrow=nrow(y)), Sigma_xi=cov(tim), homogen_traj = TRUE) scoreTest_homogen$set_pval scoreTest_heterogen <- vc_test_asym(y, x, phi=tim, indiv=rep(1:nindiv, each=nt), w=matrix(1, ncol=ncol(y), nrow=nrow(y)), Sigma_xi=cov(tim), homogen_traj = FALSE) scoreTest_heterogen$set_pval
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.