TcGSA.LR.parallel: Parallel computing the Likelihood Ratios for the Gene Sets under Scrutiny

Description

A parallel version of the function TcGSA.LR to be used on a cluster of computing processors. This function computes the Likelihood Ratios for the gene sets under scrutiny, as well as estimations of genes dynamics inside those gene sets through mixed models.

Usage

TcGSA.LR.parallel(Ncpus, type_connec, expr, gmt, design, subject_name = "Patient_ID", time_name = "TimePoint", crossedRandom = FALSE, covariates_fixed = "", time_covariates = "", time_func = "linear", group_name = "", separateSubjects = FALSE, minGSsize = 10, maxGSsize = 500, monitorfile = "")

Arguments

Ncpus

The number of processors available on the cluster.

type_connec

The type of connection between the processors. Supported cluster types are "SOCK", "PVM", "MPI", and "NWS". See also makeCluster.

expr

a matrix or dataframe of gene expression. Its dimension are $n$x$p$, with the $p$ samples in column and the $n$ genes in row.

gmt

a gmt object containing the gene sets definition. See GSA.read.gmt and definition on www.broadinstitute.org.

design

a matrix or dataframe containing the experimental variables that used in the model, namely subject_name, time_name, and covariates_fixed and time_covariates if applicable. Its dimension are $p$x$m$ and its row are is in the same order as the columns of expr.

subject_name

the name of the factor variable from design that contains the information on the repetition units used in the mixed model, such as the patient identifiers for instance. Default is 'Patient_ID'. See Details.

time_name

the name of the numeric or factor variable from design contains the information on the time replicates (the time points at which gene expression was measured). Default is 'TimePoint'. See Details.

crossedRandom

logical flag indicating wether the random effects of the subjects and of the time points should be modeled as one crossed random effect or as two separated random effects. Default is FALSE. See details.

covariates_fixed

a character vector with the names of numeric or factor variables from the design matrix that should appear as fixed effects in the model. See details. Default is "", which corresponds to no covariates in the model.

time_covariates

the name of a numeric variable from design that contains the information on the time replicates (the time points at which gene expression was measured). Default is 'TimePoint'. See Details.

time_func

the form of the time trend. Can be either one of "linear", "cubic", "splines" or specified by the user, or the column name of a factor variable from design. If specified by the user, it must be as an expression using only names of variables from the design matrix with only the three following operators: +, *, / . The "splines" form corresponds to the natural cubic B-splines (see also ns). If there are only a few timepoints, a "linear" form should be sufficient. Otherwise, the "cubic" form is more parsimonious than the "splines" form, and should be sufficiently flexible. If the column name of a factor variable from design is supplied, then time is considered as discrete in the analysis. If the user specify a formula using column names from design, both factor and numeric variables can be used.

group_name

in the case of several treatment groups, the name of a factor variable from the design matrix. It indicates to which treatment group each sample belongs to. Default is "", which means that there is only one treatment group. See Details.

separateSubjects

logical flag indicating that the analysis identifies gene sets that discriminates patients rather than gene sets than have a significant trend over time. Default is FALSE. See Details.

minGSsize

the minimum number of genes in a gene set. If there are less genes than this number in one of the gene sets under scrutinity, the Likelihood Ratio of this gene set is not computed (the mixed model are not fitted). Default is 10 genes as the minimum.

maxGSsize

the maximum number of genes in a gene set. If there are more genes than this number in one of the gene sets under scrutinity, the Likelihood Ratio of this gene set is not computed (the mixed model are not fitted). This is to avoid very long computation times. Default is 500 genes as the maximum.

monitorfile

a writable connections or a character string naming a file to write into, to monitor the progress of the analysis. Default is "" which is no monitoring. See Details.

Value

TcGSA.LR returns a tcgsa object, which is a list with the 5 following elements:

fit a data frame that contains the 3 following variables:
- LR: the likelihood ratio between the model under the null hypothesis and the model under the alternative hypothesis.
- CVG_H0: convergence status of the model under the null hypothesis.
- CVG_H1: convergence status of the model under the alternative hypothesis.
time_func: a character string passing along the value of the time_func argument used in the call.
GeneSets_gmt: a gmt object passing along the value of the gmt argument used in the call.
group.var: a factor passing along the group_name variable from the design matrix.
separateSubjects: a logical flag passing along the value of the separateSubjects argument used in the call.
Estimations: a list of 3 dimensions arrays. Each element of the list (i.e. each array) corresponds to the estimations of gene expression dynamics for each of the gene sets under scrutiny (obtained from mixed models). The first dimension of those arrays is the genes included in the concerned gene set, the second dimension is the Patient_ID, and the third dimension is the TimePoint. The values inside those arrays are estimated gene expressions.
time_DF: the degree of freedom of the natural splines functions

Details

This Time-course Gene Set Analysis aims at identifying gene sets that are not stable over time, either homogeneously or heterogeneously (see Hejblum et al, 2012) in terms of their probes. And when the argument separatePatients is TRUE, instead of identifying gene sets that have a significant trend over time (possibly with probes heterogeneity of this trend), TcGSA identifies gene sets that have significantly different trends over time depending on the patient.

If the monitorfile argument is a character string naming a file to write into, in the case of a new file that does not exist yet, such a new file will be created. A line is written each time one of the gene sets under scrutiny has been analysed (i.e. the two mixed models have been fitted, see TcGSA.LR) by one of the parallelized processors.

References

Hejblum BP, Skinner J, Thiebaut R, (2015) Time-Course Gene Set Analysis for Longitudinal Gene Expression Data. PLoS Computat Biol 11(6): e1004310. doi: 10.1371/journal.pcbi.1004310

Examples

Run this code


data(data_simu_TcGSA)

tcgsa_sim_1grp <- TcGSA.LR(expr=expr_1grp, gmt=gmt_sim, design=design, 
                          subject_name="Patient_ID", time_name="TimePoint",
                          time_func="linear", crossedRandom=FALSE)
                          
## Not run:  
# require(doParallel)
# tcgsa_sim_1grp <- TcGSA.LR.parallel(Ncpus = 2, type_connec = 'SOCK',
#                             expr=expr_1grp, gmt=gmt_sim, design=design, 
#                             subject_name="Patient_ID", time_name="TimePoint",
#                             time_func="linear", crossedRandom=FALSE, 
#                             separateSubjects=TRUE)
# ## End(Not run)
tcgsa_sim_1grp
summary(tcgsa_sim_1grp)

Run the code above in your browser using DataLab