Hierarchical Parallelization of Gene Differential Association Analysis
 Mark Needham^{1},
 Rui Hu^{2}Email author,
 Sandhya Dwarkadas^{1} and
 Xing Qiu^{2}
DOI: 10.1186/1471210512374
© Needham et al; licensee BioMed Central Ltd. 2011
Received: 28 January 2011
Accepted: 21 September 2011
Published: 21 September 2011
Abstract
Background
Microarray gene differential expression analysis is a widely used technique that deals with high dimensional data and is computationally intensive for permutationbased procedures. Microarray gene differential association analysis is even more computationally demanding and must take advantage of multicore computing technology, which is the driving force behind increasing compute power in recent years. In this paper, we present a twolayer hierarchical parallel implementation of gene differential association analysis. It takes advantage of both fine and coarsegrain (with granularity defined by the frequency of communication) parallelism in order to effectively leverage the nonuniform nature of parallel processing available in the cuttingedge systems of today.
Results
Our results show that this hierarchical strategy matches data sharing behavior to the properties of the underlying hardware, thereby reducing the memory and bandwidth needs of the application. The resulting improved efficiency reduces computation time and allows the gene differential association analysis code to scale its execution with the number of processors. The code and biological data used in this study are downloadable from http://www.urmc.rochester.edu/biostat/people/faculty/hu.cfm
Conclusions
The performance sweet spot occurs when using a number of threads per MPI process that allows the working sets of the corresponding MPI processes running on the multicore to fit within the machine cache. Hence, we suggest that practitioners follow this principle in selecting the appropriate number of MPI processes and threads within each MPI process for their cluster configurations. We believe that the principles of this hierarchical approach to parallelization can be utilized in the parallelization of other computationally demanding kernels.
Background
Microarray gene differential expression analysis has been widely used to uncover the underlying biological mechanism. Researchers utilize this technology to identify potentially "interesting" genes. More specifically, a statistical test is applied to each individual gene to detect whether the mean expression level of this gene is the same or not across different biological conditions or phenotypes studied in an experiment. A chosen multiple testing procedure (MTP) is then employed to control certain perfamily Type I errors. Genes work together to fulfill certain biological functions and they are known to be strongly correlated [1, 2]. The structure of intergene correlation contains rich information that cannot be extracted from mean expression levels. Recent years have seen more and more research focusing on gene dependence structures. For example, some procedures, such as gene set enrichment analysis [3, 4], incorporate existing biological gene sets information into statistical procedures. Gene cluster analysis uses gene dependence and similarity to group genes [5–11]. Gene network analysis, such as method based on Gaussian or Bayesian networks, employs gene dependence to study gene dynamics and reasoning [12–14]. Another approach is to directly select genes based on the phenotypic differences of their dependence structure [15–20]. In this paper, we consider the very last approach and focus on a gene differential association analysis (henceforth denoted as GDAA) procedure proposed in [19]. Unlike traditional differential gene expression analysis, GDAA is designed to select genes that have different dependence structures with other genes in two phenotypes. It complements the analysis of differentially expressed genes. Combining both gene differential association analysis and gene differential expression analysis provides a more comprehensive functional interpretation of the experimental results. As an example, GDAA was applied in [20] to two sets of Childhood Leukemia data (HYPERDIP and TEL) [21] and selected differentially associated (DA) genes that could not be detected by differential gene expression analysis. Furthermore, the TEL group is differentiated from other leukemia subtypes by the presence of t(12;21)(p13;q22) translocation, generating the TELAML1 fusion gene. Through the overrepresentation of DA genes, the chromosomal band 21q22.3 containing the TELAML1 fusion gene was identified. This chromosomal band was not identified by differential gene expression analysis.
A typical microarray data set reports expression levels for tens of thousands of genes. For example, both sets of Childhood Leukemia data HYPERDIP and TEL [21] have expression levels for m = 7, 084 genes updated from the original expression levels by using a custom CDF file to produce values of gene expressions. The CDF files can be found at http://brainarray.mbni.med.umich.edu. Please see [19] for more details. Each slide is then represented by an array reporting the logarithm (base 2) of expression level on the set of 7,084 genes. For convenience, the words "gene" and "gene expression" are used interchangeably to refer to these gene expressions in this paper. Due to such a high dimensionality, the computation of traditional gene differential expression analysis is considered to be more time consuming than many traditional statistical analyses in medical research. A gene selection procedure based on gene dependence structures has to be even more computationally intensive. This is because the dependence structure is typically measured by a pertinent association score, such as the Pearson correlation coefficient for all gene pairs, of which the multiplicity (dimensionality) is $\frac{m\left(m1\right)}{2}$ instead of m. It is therefore more computationally intensive to detect the differences hidden in the correlation matrix. In particular, for the procedure proposed in [19], the length of the computation is O(m × m × n × K), where m = 7, 084 is the number of genes, n = 79 is the number of subjects in each phenotypic group, and K = 10, 000 is the number of permutations for approximating the statistical null distribution. Such large number of permutations is necessary because statistical inference for microarray analysis is based on multiple testing adjusted pvalues, which demands much finer estimation of unadjusted pvalues compared to regular permutation tests. With a large number of genes and a medium sample size, running GDAA can take several days or even a month. For example, a sequential implementation of the procedure in [19] took nearly two months to complete the calculation on a computer with a 2 GHz AMD Opteron processor and 2GB SDRAM. Until about 2003, processor designers were able to leverage technology advances that allowed increasing numbers of smaller and faster transistors on a single chip in order to improve the performance of sequential computation. Hence, it was possible for computational scientists who wanted their codes to run faster to simply wait for the next generation of machines. However, the reality is that around 2003, chipmakers discovered that they were no longer able to sustain faster sequential execution due to the inability to dissipate the heat generated by the computation [22]. Consequently, designers turned to using the increasing transistor counts to add more processors, each of which execute independent sequential computation. The processors typically share access to the memory subsystem and offchip bandwidth. These multicore chips now dominate the desktop market and are used to put together multiprocessor servers consisting of multiple processor chips, as well as networked clusters of such servers for highend computation. Parallel computing (utilizing multiple compute resources simultaneously for the same application) that effectively leverages these increasingly multicore clusters of multiprocessors is thus even more critical than in the past in order to obtain results in a timely manner.
In this paper, we propose a new parallel design for the gene differential association analysis procedure in [19]. The key to our parallelization strategy is that it takes advantage of both fine and coarsegrain parallelism (the granularity representing the frequency of sharing/communication in the concurrent computation). The hardwarebased memory sharing within a multicore is utilized for the finegrain parallelism (with higher need for sharing/communication). Sharing memory in hardware avoids the need for data replication. Since GDAA utilizes a multivariate nonparametric test, it has more memory needs than a comparable gene differential expression analysis. Therefore, the memory sharing feature in our strategy is also critical to reducing the bandwidth demands of the GDAA procedure. The results show that our strategy leverages GDAA's characteristics to reduce the memory and bandwidth needs of the application, thereby improving computational efficiency.
Implementation
Gene Differential Association Analysis Procedure
We outline the related GDAA procedure below. More details can be found in [19].
Statistical Hypothesis Testing
Assume there are two biological conditions or phenotypes A and B. Under each condition n subjects are sampled, each measured with m gene expression levels.
where k = 1, ⋯, i  1, i + 1, ⋯, m. We denote the correlation vectors (w_{i 1}, ⋯, w_{i,i1}, w_{i,i+1}, ⋯, w_{ im } ) by w _{i}. This vector represents the relationship between the i th gene and all other genes.
where ${F}_{{w}_{i}\left(A\right)}\left(x\right)$ and ${F}_{{w}_{i}\left(B\right)}\left(x\right)$ are the joint distribution functions of w _{i}(A) and w _{i}(B), respectively. If H_{ i } is rejected, we declare the i th gene to be a differentially associated gene.
The Nstatistic
In order to test H _{i}, we need to create samples of correlation vectors to mimic the joint distributions ${F}_{{w}_{i}\left(A\right)}\left(x\right)$ and ${F}_{{w}_{i}\left(B\right)}\left(x\right)$, respectively. We divide the dataset under condition A into$G\left(1\le G\le \frac{n}{2}\right)$ subgroups, each subgroup containing $\frac{n}{G}$ subjects. In order to compute correlation coefficients, every subgroup must contain at least two subjects. Sample sizes of subgroups do not have to be equal. When G does not divide n, the last few subgroups can have a slightly larger or smaller sample size. That being said, an approximately even partition of subgroups is still desirable because it leads to better statistical power than unbalanced partitions.
From these subgroups, we compute a sample of size G correlation vectors for the i th gene, denoted by w _{i}(A, k), 1 ≤ k ≤ G. Similarly, we have a sample of size G correlation vectors for the i th gene under condition B, denoted by w _{i}(B, k), 1 ≤ k ≤ G.
The Nstatistic can serve as a measurement of how much the intergene correlation structure of the i th gene has changed from condition A to condition B.
Permutationbased Null Distribution and pvalue
Denote ${N}_{i}^{*}$ as the Nstatistic associated with the i th gene. To determine the statistical significance of ${N}_{i}^{*}$, which is represented by a pvalue, we need to model the null distribution of this statistic. This can be done by the following resampling method. First, we combine the gene expression data under both conditions and randomly permute subjects. Then we divide them into two groups of equal size, mimicking two biological conditions without differentially associated genes. By applying formula (1), we get a permutation based Nstatistic for the i th gene, which can be considered as an observation from the null distribution of N_{ i } , i.e., the distribution of N_{ i } when H_{ i } holds. Repeating this permutation process K times produces K permutation based Nstatistics for the i th gene, denoted by N_{ ik } , 1 ≤ k ≤ K.
The smaller ${\stackrel{\u0303}{p}}_{i}$ is, the more likely w _{i}(A) is different from w _{i}(B), i.e., the i th gene changes its relationship with all other genes across conditions A and B. If ${\stackrel{\u0303}{p}}_{i}$ is less than a predefined threshold, we reject H_{ i } and declare the i th gene to be a differentially associated gene.
Summary of the GDAA Procedure
 1.
Divide the subjects (slides) under each condition (A or B) into G subgroups such that there are approximately $\frac{n}{G}$ subjects for each subgroup.
 2.
For each gene, compute its correlation vectors from all subgroups. This step produces G correlation vectors for one gene in each condition.
 3.
Compute the Nstatistic for the i th gene from these 2 × G samples using Equation(1) and record it as ${N}_{i}^{*}$.
 4.
Pool the subjects in both conditions together. Randomly shuffle the subjects, and then split them into two groups of equal size.
 5.
Divide the subjects in each group into G subgroups and compute G correlation vectors in each subgroup for each gene.
 6.
Compute the Nstatistics for each gene based on these 2 × G correlation vectors.
 7.
Repeat steps 4 to 6 K times and record the permutationbased Nstatistics as N_{ ik } , i = 1, ..., m, k = 1, ..., K.
 8.
Obtain the permutationbased pvalue, p_{ i } , using Equation(2).
 9.
Adjust pvalue by using Equation(3). Select differentially associated genes based on the adjusted pvalues and a prespecified PFER level.
Computation
 1.
Read in and initialize data (performed in Python on the master process).
 2.
Calculate ${N}_{i}^{*}$, 1 ≤ i ≤ m, for the unpermuted dataset using a single core (C++ code) on the master process.
 3.
Create K permutations of the original dataset; distribute the permutations k (1 ≤ k ≤ K) to independent slave processes (performed in Python) using MPI [28]. Work is distributed at the granularity of a single permutation  when a process completes the computation for one permutation, it requests the next permutation.
 4.
Each worker/slave receives a permutation (using Python), permutes its local copy of the data, and then computes the vector of Nstatistics using C++, parallelized using the Pthreads [29] package. A total number of P threads are created and the pergene computation is distributed among threads so that each thread performs the Nstatistic computation for $\frac{m}{P}$ genes. When m is not divisible by P, each thread receives a slightly different number of genes.
 5.
Once an MPI process has determined that its threads have computed all Nstatistics (N_{ ik } , 1 ≤ i ≤ m) for the k th permutation it was assigned, it then returns them to the master MPI process.
 6.
The master MPI process collects all the N_{ ik } to calculate pvalues p_{ i } (performed in Python).
Steps 1 and 2 of the algorithm are performed sequentially. To parallelize the remaining steps, we use a twotiered approach. At the first level, we distribute the work by spawning processes from Python using MPI. One MPI process, the "master" process, is responsible for distributing different permutations to the other "slave" processes. Each slave independently permutes the gene data according to the permutation indices received from the master process and computes the Nstatistics for the permuted data (this code is optimized C++ code). The computed values (the vector of Nstatistics) are then returned to the master using an MPI message, where the Python code calculates the pvalue. The key to this implementation is that the core computation is performed in optimized C++ code.
Data Sharing
Gene expression level data in each biological condition is represented using a dynamically allocated (m × n)dimensional array, where n is the number of subjects and m is the number of genes. This twodimensional array is readshared within each MPI process and its size grows as a product of m and n. There are two other dynamically allocated twodimensional arrays created for each MPI process. These two arrays with sizes proportional to (m × G) are used for temporary storage of the correlation computation. One stores the sums of the expression levels within subgroups so that its entry in row i and column j is ${\sum}_{l\in \mathsf{\text{Subgroup}}\phantom{\rule{1em}{0ex}}j}{x}_{il}^{\left(k\right)}$. Here ${x}_{il}^{\left(k\right)}$ is the expression level for gene i and subject l in the k th permutation. The other stores the sums of squares of the expression levels within subgroups so that its entry in row i and column j is ${\sum}_{l\in Subgroupj}({x}_{il}^{(k)}}{)}^{2$. They are also read shared within the MPI process. Another twodimensional dynamically allocated array with size proportional to O(m × G) is created for each thread, storing the correlation vectors for each gene. This array is both read and written by the thread. The Nstatistic, which is a vector with the size of the number of genes m, is also shared across all threads within each MPI process. Each thread writes to independent regions of this vector based on the genes allocated to it.
Results
Our evaluation is conducted on a cluster of five machines, each with 16 GBytes of memory and two 3.0 GHz quadcore Intel Xeon 5450 processors, for a total of 40 processors. The machines are interconnected using Gigabit ethernet. Each quadcore processor chip has 6 MBytes of lastlevel cache per pair of cores. Each machine runs Linux version 2.6.16.600.21smp. The application was compiled using Python version 2.4.2, Pypar version 2.1.0, and gcc version 4.1.2 at the O2 optimization level.
Simulation Data
To gain better insight into the effects of different configurations on performance, we simulate several sets of data. Each set has two groups of n = 100 slides representing 100 subjects in each biological conditions.
Each array has m genes, where m takes on values in {1000, 2000, 3000, 4000, 5000, 6000, 7000, 8000, 9000, 10000}. The slides in each group are divided into G = 10 subgroups to calculate correlation vector samples. K = 100 permutations of the subjects in the two groups are created in order to generate the null distributions of Nstatistics.
Performance Analysis
In Figures 4 and 5, we notice that the speedup curves are not very smooth. This step function can be attributed to several causes. The first is load imbalance due to the granularity of workload distribution  permutations at the MPI parallelization level and genes at the Pthread parallelization level. When using 1 thread per MPI process to conduct 100 permutations, as one example, with 5 processors (20 cores), each core runs 5 permutations (⌈100/20⌉). If we increase the number of processors to 6 (24 cores), some cores will still execute 5 permutations while others execute 4, so that execution time remains proportional to ⌈100/24⌉ = 5, resulting in practically no increase in speedup. As the number of permutations executed per MPI process decreases (with an increasing number of cores), the fraction of idle/wasted time on the cores with one less permutation to execute increases, resulting in lower efficiency. In the case of Figure 4, the increased scheduling variability and poor choice of scheduling when adding a quadcore processor within a machine also contributes to the step function in the 2 and 4 threaded curves. Once the scheduling is made both deterministic and ensures appropriate cache sharing, the step function is less pronounced in the multithreaded runs in Figure 5 due to their reduced memory bandwidth demands and smoother load function at the MPI level.
Conclusions
Microarray technology has made it possible for medical researchers to measure and study the behavior of thousands of genes at once. Technology advances have been on a fast track in recent years, making it possible to conduct microarray experiments much faster and less expensive than in the past. This trend has been leveraged with the availability of larger and larger datasets. Turning so much raw information into knowledge presents a major challenge for both statistical analysis and computation. As of now, microarray data are used for crude screening of differentially expressed genes. Exploiting the rich information contained in the intergene dependence structure has not become a routine, despite the availability of several gene association analysis procedures. This is largely due to the computing bottleneck.
In this paper, we present a parallelized implementation of gene differential association analysis that is designed to leverage the features of today's multicore platforms in which resources are shared among processors at a much finer granularity than in the past. We apply the conventional wisdom of parallelizing at the coarsest granularity to distribute permutations among the nodes in a cluster, using MPI for communication. In addition, we parallelize at the finer granularity of pergene computation within a single dual quadcore machine using shared memory (Pthreads). Sharing memory across threads helps reduce demand for the shared lastlevel cache capacity on the chip by allowing independent threads to share a single copy of the gene expression data. Our results show that this strategy utilizes the multicore cluster platform much more effectively. In general, the performance sweet spot occurs when using a number of threads that allows the working sets of the corresponding MPI processes to fit within the machine's shared cache. We suggest that practitioners follow the principle of determining what resources are shared when making decisions on how to allocate compute resources among MPI processes and threads for their cluster machines. We believe that the principles of this hierarchical approach to parallelization can be utilized in the parallelization of other computationally demanding kernels.
Availability and Requirements

Project name: Hierarchical Parallelization of Gene Differential Association Analysis;

Project home page: http://www.urmc.rochester.edu/biostat/people/faculty/hu.cfm;

Operating system: Linux;

Programming language: Python and C++;

Other requirements: MPI (MPICH2 or Open MPI), Python, C++ Compilation tools, SWIG, Numpy, Pypar;

Licence: GNU GENERAL PUBLIC LICENSE, Version 2, June 1991;

No restrictions to use by nonacademics.
Abbreviations
 GDAA:

Gene Differential Association Analysis
 MTP:

Multiple Testing Procedure
 MPI:

Message Passing Interface
 Pthreads:

POSIX Threads.
Declarations
Acknowledgements
This work was supported in part by NSF grants CCF0702505, CNS0411127, CNS0615139, CNS0834451, CNS0509270, and CCF1016902; and NIH grants 5 R21 GM07925902, 1 R21 HG00464801, and NCRR UL1 RR024160. Any opinions, findings, and conclusions or recommendations expressed in this material are those of the author(s) and do not necessarily reflect the views of the abovenamed organizations. In addition, we would like to thank Ms. Christine Brower for her technical assistance with computing and Ms. Malora Zavaglia for her proofreading effort. Finally we are grateful to the associated editor and two anonymous reviewers for their constructive comments which helped us improve the manuscript.
Authors’ Affiliations
References
 Klebanov L, Jordan C, Yakovlev A: A new type of stochastic dependence revealed in gene expression data. Stat Appl Genet Mol Biol 2006, 5: Article7. [http://dx.doi.org/10.2202/1544–6115.1189]PubMed
 Bhardwaj N, Lu H: Correlation between gene expression profiles and proteinprotein interactions within and across genomes. Bioinformatics 2005, 21(11):2730–2738. [http://dx.doi.org/10.1093/bioinformatics/bti398] 10.1093/bioinformatics/bti398View ArticlePubMed
 Mootha V, Lindgren C, Eriksson K, Subramanian A, Sihag S, Lehar J, Puigserver P, Carlsson E, Ridderstråle M, Laurila E, et al.: PGC1 αresponsive genes involved in oxidative phosphorylation are coordinately downregulated in human diabetes. Nature genetics 2003, 34(3):267–273. 10.1038/ng1180View ArticlePubMed
 Subramanian A, Tamayo P, Mootha V, Mukherjee S, Ebert B, Gillette M, Paulovich A, Pomeroy S, Golub T, Lander E, et al.: Gene set enrichment analysis: a knowledgebased approach for interpreting genomewide expression profiles. Proceedings of the National Academy of Sciences 2005, 102(43):15545–15550. 10.1073/pnas.0506580102View Article
 Raychaudhuri S, Stuart J, Altman R: Principal components analysis to summarize microarray experiments: application to sporulation time series. Pac Symp Biocomput 2000, 5: 455–466.
 Liu A, Zhang Y, Gehan E, Clarke R: Block principal component analysis with application to gene microarray data classification. Statistics in medicine 2002., 21(22):
 Wang A, Gehan E: Gene selection for microarray data analysis using principal component analysis. Statistics in medicine 2005., 24(13):
 Eisen M, Spellman P, Brown P, Botstein D: Cluster analysis and display of genomewide expression patterns. Proceedings of the National Academy of Sciences 1998, 95(25):14863–14868. 10.1073/pnas.95.25.14863View Article
 Törönen P, Kolehmainen M, Wong G, Castrén E: Analysis of gene expression data using selforganizing maps. FEBS letters 1999, 451(2):142–146. 10.1016/S00145793(99)005244View ArticlePubMed
 Furey T, Cristianini N, Duffy N, Bednarski D, Schummer M, Haussler D: Support vector machine classification and validation of cancer tissue samples using microarray expression data. 2000.
 Brown M, Grundy W, Lin D, Cristianini N, Sugnet C, Furey T, Ares M, Haussler D: Knowledgebased analysis of microarray gene expression data by using support vector machines. Proceedings of the National Academy of Sciences 2000, 97: 262–267. 10.1073/pnas.97.1.262View Article
 Bahar I, Atilgan AR, Erman B: Direct evaluation of thermal fluctuations in proteins using a singleparameter harmonic potential. Fold Des 1997, 2(3):173–181. 10.1016/S13590278(97)000242View ArticlePubMed
 Friedman N: Inferring cellular networks using probabilistic graphical models. Science 2004, 303(5659):799–805. [http://dx.doi.org/10.1126/science.1094068] 10.1126/science.1094068View ArticlePubMed
 OpgenRhein R, Strimmer K: From correlation to causation networks: a simple approximate learning algorithm and its application to highdimensional plant gene expression data. BMC Syst Biol 2007, 1: 37. [http://dx.doi.org/10.1186/1752–0509–137] 10.1186/17520509137PubMed CentralView ArticlePubMed
 Li K: Genomewide coexpression dynamics: theory and application. Proceedings of the National Academy of Sciences 2002, 99(26):16875–16880. 10.1073/pnas.252466999View Article
 Lai Y, Wu B, Chen L, Zhao H: A statistical method for identifying differential genegene coexpression patterns. Bioinformatics 2004, 20(17):3146–3155. [http://dx.doi.org/10.1093/bioinformatics/bth379] 10.1093/bioinformatics/bth379View ArticlePubMed
 Shedden K, Taylor J: Differential correlation detects complex associations between gene expression and clinical outcomes in lung adenocarcinomas. Methods of Microarray Data Analysis IV 2005, 121–131.View Article
 Choi J, Yu U, Yoo O, Kim S: Differential coexpression analysis using microarray data and its application to human cancer. Bioinformatics 2005, 21(24):4348–4355. 10.1093/bioinformatics/bti722View ArticlePubMed
 Hu R, Qiu X, Glazko G, Klebanov L, Yakovlev A: Detecting intergene correlation changes in microarray analysis: a new approach to gene selection. BMC Bioinformatics 2009, 10: 20. [http://dx.doi.org/10.1186/1471–2105–10–20] 10.1186/147121051020PubMed CentralView ArticlePubMed
 Hu R, Qiu X, Glazko G: A new gene selection procedure based on the covariance distance. Bioinformatics 2010, 26(3):348–354. [http://dx.doi.org/10.1093/bioinformatics/btp672] 10.1093/bioinformatics/btp672PubMed CentralView ArticlePubMed
 Yeoh EJ, Ross ME, Shurtleff SA, Williams WK, Patel D, Mahfouz R, Behm FG, Raimondi SC, Relling MV, Patel A, Cheng C, Campana D, Wilkins D, Zhou X, Li J, Liu H, Pui CH, Evans WE, Naeve C, Wong L, Downing JR: Classification, subtype discovery, and prediction of outcome in pediatric acute lymphoblastic leukemia by gene expression profiling. Cancer Cell 2002, 1(2):133–143. 10.1016/S15356108(02)000326View ArticlePubMed
 Patterson D: The Trouble with Multicore Microprocessors. IEEE Spectrum 2010, 28–32.
 Szabo A, Boucher K, Carroll W, Klebanov L, Tsodikov A, Yakovlev A: Variable selection and pattern recognition with gene expression data generated by the microarray technology. Mathematical Biosciences 2002, 176: 71–98. 10.1016/S00255564(01)001031View ArticlePubMed
 Szabo A, Boucher K, Jones D, Tsodikov AD, Klebanov LB, Yakovlev AY: Multivariate exploratory tools for microarray data analysis. Biostatistics 2003, 4(4):555–567. [http://dx.doi.org/10.1093/biostatistics/4.4.555] 10.1093/biostatistics/4.4.555View ArticlePubMed
 Xiao Y, Frisina R, Gordon A, Klebanov L, Yakovlev A: Multivariate search for differentially expressed gene combinations. BMC Bioinformatics 2004, 5: 164. [http://dx.doi.org/10.1186/1471–2105–5164] 10.1186/147121055164PubMed CentralView ArticlePubMed
 Klebanov L, Gordon A, Xiao Y, Land H, Yakovlev A: A permutation test motivated by microarray data analysis. Computational Statistics and Data Analysis 2005.
 Gordon A, Glazko G, Qiu X, Yakovlev A: Control of the Mean Number of False Discoveries, Bonferroni, and Stability of Multiple Testing. The Annals of Applied Statistics 2007, 1: 179–190. [http://projecteuclid.org/euclid.aoas/1183143734] 10.1214/07AOAS102View Article
 Message Passing Interface Forum: MPI: A MessagePassing Interface Standard, Version 2.2.2009. [http://www.mpiforum.org/docs/]
 Barney B: POSIX Threads Programming.2011. [https://computing.llnl.gov/tutorials/pthreads/]
Copyright
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.