Volume 6, Issue 1

Differentiating strains of a pathogen is often central to investigating its epidemiological aspects. The genetic similarity of a group of strains can be assessed by calculating a matrix of dissimilarities from their DNA fingerprinting profiles. The mean dissimilarity for each strain across other strains within the group is then used as an observation in a statistical analysis. These observations are not independent of each other, and so standard analysis techniques such as the t-test are inappropriate, because they underestimate the variance of the group means, and hence overstate the statistical significance of any differences. By examining the correlation between elements of the dissimilarity matrix, it is shown that the variance is underestimated by a factor of between about 2 and 4. Permutation tests are proposed as a way of addressing the problem of dependence, and are applied to a study of fluconazole resistance in Candida albicans.

This paper discusses characteristics of dye biases in microarray data that the conventional normalization methods do not handle, and proposes a new normalization method involving a mixture of splines model. We also develop a test for between-group comparisons of each gene that is designed to be used with our proposed method.

For testing multiple null hypotheses, the classical approach to dealing with the multiplicity problem is to restrict attention to procedures that control the familywise error rate (FWER), the probability of even one false rejection. In many applications, one might be willing to tolerate more than one false rejection provided the number of such cases is controlled, thereby increasing the ability of the procedure to detect false null hypotheses. This suggests replacing control of the FWER by controlling the probability of k or more false rejections, which is called the k-FWER. In Hommel and Hoffmann (1987) and Lehmann and Romano (2005a), single step and stepdown procedures are derived that control the k-FWER, without making any assumptions concerning the dependence structure of the p-values of the individual tests. However, if the p-values are mutually independent, one can improve the procedures. In fact, Sarkar (2005) provided such an improvement. However, we show other improvements are possible which appear to be generally much better, and are sometimes unimprovable. When k=1, the procedure reduces to the classical method of Sidak, and the stepdown procedure is unimprovable and strictly dominates that of Sarkar. Under a monotonicity condition, an unimprovable procedure is obtained. In the case k=2, the monotonicity condition is satisfied, and the condition can be checked numerically in general. We then develop a stepdown method that controls the false discovery proportion. Except for the case of k-FWER control with k=1, the gains are surprisingly dramatic, and theoretical and numerical evidence is given.

Although researchers use duplicate genotyped data to calculate an inconsistency rate, there is no power analysis to assess the value of the duplicate data. In this paper, we present a model in which the genotyping error rate is related to the inconsistency rate. We extend the g genotype by h phenotype chi-squared test to incorporate the duplicate genotyped data. When a subject is inconsistently genotyped (that is, has two observed genotypes), our procedure is to allocate 0.5 units to each of the two genotypes. We specify the multivariate analysis of variance (MANOVA) test comparing these extended counts. We provide freely available software for this test and also for a permutation test used on small samples. A simulation study shows that the asymptotic null distribution of the MANOVA test holds when the total number of subjects, N, is at least 300. We also document with a simulation study that the asymptotic distribution of this test under various alternative hypotheses is a satisfactory approximation to the simulated power. In all cases, the power of the MANOVA test using the duplicate genotyped data is greater than the power of the chi-squared test ignoring the duplicate data. Power increases ranged from 0.776% to 4.652% for 80% powered tests and 0.292% to 2.028% for 95% powered tests. Researchers now can compute the value of the duplicate genotyped data as part of the design of the study.

Multi-color optical mapping is a new technique being developed to obtain detailed physical maps (indicating relative positions of various recognition sites) of DNA molecules. We consider a study design in which the data consist of noisy observations of multiple copies of a DNA molecule marked with colors at recognition sites. The primary goal is to estimate a physical map. A secondary goal is to estimate error rates associated with the experiment, which are potentially useful for analysis and refinement of the biochemical steps in the mapping procedure. We propose statistical models for various sources of error and use maximum likelihood estimation (MLE) to construct a physical map and estimate error rates. To overcome difficulties arising in the maximization process, a latent-variable Markov chain version of the model is proposed, and the EM algorithm is used for maximization. In addition, a simulated annealing procedure is applied to maximize the profile likelihood over the discrete space of sequences of colors. We apply the methods to simulated data on the bacteriophage lambda genome.

In this paper, we propose a novel method for sparse logistic regression with non-convex regularization Lp (p <1). Based on smooth approximation, we develop several fast algorithms for learning the classifier that is applicable to high dimensional dataset such as gene expression. To the best of our knowledge, these are the first algorithms to perform sparse logistic regression with an Lp and elastic net (Le) penalty. The regularization parameters are decided through maximizing the area under the ROC curve (AUC) of the test data. Experimental results on methylation and microarray data attest the accuracy, sparsity, and efficiency of the proposed algorithms. Biomarkers identified with our methods are compared with that in the literature. Our computational results show that Lp Logistic regression (p <1) outperforms the L1 logistic regression and SCAD SVM. Software is available upon request from the first author.

Many alternative data-adaptive algorithms can be used to learn a predictor based on observed data. Examples of such learners include decision trees, neural networks, support vector regression, least angle regression, logic regression, and the Deletion/Substitution/Addition algorithm. The optimal learner for prediction will vary depending on the underlying data-generating distribution. In this article we introduce the "super learner", a prediction algorithm that applies any set of candidate learners and uses cross-validation to select between them. Theory shows that asymptotically the super learner performs essentially as well as or better than any of the candidate learners. In this article we present the theory behind the super learner, and illustrate its performance using simulations. We further apply the super learner to a data example, in which we predict the phenotypic antiretroviral susceptibility of HIV based on viral genotype. Specifically, we apply the super learner to predict susceptibility to a specific protease inhibitor, nelfinavir, using a set of database-derived non-polymorphic treatment-selected mutations.

A number of computational methods have been proposed for identifying transcription factor binding sites from a set of unaligned sequences that are thought to share the motif in question. We here introduce an algorithm, called cosmo, that allows this search to be supervised by specifying a set of constraints that the position weight matrix of the unknown motif must satisfy. Such constraints may be formulated, for example, on the basis of prior knowledge about the structure of the transcription factor in question. The algorithm is based on the same two-component multinomial mixture model used by MEME, with stronger reliance, however, on the likelihood principle instead of more ad-hoc criteria like the E-value. The intensity parameter in the ZOOPS and TCM models, for instance, is estimated based on a profile-likelihood approach, and the width of the unknown motif is selected based on BIC. These changes allow cosmo to outperform MEME even in the absence of any constraints, as evidenced by 2- to 3-fold greater sensitivity in some simulation studies. Additional improvements in performance can be achieved by selecting the model type (OOPS, ZOOPS, or TCM) data-adaptively or by supplying correctly specified constraints, especially if the motif appears only as a weak signal in the data. The algorithm can data-adaptively choose between working in a given constrained model or in the completely unconstrained model, guarding against the risk of supplying mis-specified constraints. Simulation studies suggest that this approach can offer 3 to 3.5 times greater sensitivity than MEME. The algorithm has been implemented in the form of a stand-alone C program as well as a web application that can be accessed at http://cosmoweb.berkeley.edu. An R package is available through Bioconductor (http://bioconductor.org).

High-dimensional case-control analysis is encountered in many different settings in genomics. In order to rank genes accordingly, many different scores have been proposed, ranging from ad hoc modifications of the ordinary t statistic to complicated hierarchical Bayesian models.Here, we introduce the “shrinkage t” statistic that is based on a novel and model-free shrinkage estimate of the variance vector across genes. This is derived in a quasi-empirical Bayes setting. The new rank score is fully automatic and requires no specification of parameters or distributions. It is computationally inexpensive and can be written analytically in closed form.Using a series of synthetic and three real expression data we studied the quality of gene rankings produced by the “shrinkage t” statistic. The new score consistently leads to highly accurate rankings for the complete range of investigated data sets and all considered scenarios for across-gene variance structures.

Quantitative Polymerase Chain Reaction (Q-PCR) aims at determining the initial quantity of specific nucleic acids from the observation of the number of amplified DNA molecules. The most widely used technology to monitor the number of DNA molecules as they replicate is based on fluorescence chemistry. Considering this measurement technique, the observation of DNA amplification by PCR contains intrinsically two kinds of variability. On the one hand, the number of replicated DNA molecules is random, and on the other hand, the measurement of the fluorescence emitted by the DNA molecules is collected with some random error. Relying on a stochastic model of these two types of variability, we aim at providing estimators of the parameters arising in the proposed model, and, more specifically, of the initial amount of molecules. The theory of branching processes is classically used to model the evolution of the number of DNA molecules at each replication cycle. The model is a binary splitting Galton-Watson branching process. Its unknown parameters are the initial number of DNA molecules and the reaction efficiency of PCR, which is defined as the probability of replication of a DNA molecule. The number of DNA molecules is indirectly observed through noisy fluorescence measurements resulting in a so-called Hidden Markov Model. We aim at inference of the parameters of the underlying branching process, and the parameters of the noise from the fluorescence measurements in a Bayesian framework. Using simulations and experimental data, we investigate the performance of the Bayesian estimators obtained by Markov Chain Monte Carlo methods.

Amplified Fragment Length Polymorphism (AFLP) markers are formed by selective amplification of DNA fragments from digested total genomic DNA. The technique is popular because it is a relatively inexpensive way to produce large numbers of reproducible genetic markers. In this paper, we describe a Bayesian approach to modeling AFLP marker evolution by nucleotide substitution and an MCMC approach to estimate phylogeny from AFLP marker data. We demonstrate the method on species in Carex section Ovales, a group of sedges common in North America. We compare the results of our analysis with a clustering method based on Nei and Li's restriction-site distance and a two-state Bayesian analysis using MrBayes.

The etiology of complex diseases is heterogeneous. The presence of risk alleles in one or more genetic loci affects the function of a variety of intermediate biological pathways, resulting in the overt expression of disease. Hence, there is an increasing focus on identifying the genetic basis of disease by systematically studying phenotypic traits pertaining to the underlying biological functions. In this paper we focus on identifying genetic loci linked to quantitative phenotypic traits in experimental crosses. Such genetic mapping methods often use a one stage design by genotyping all the markers of interest on the available subjects. A genome scan based on single locus or multi-locus models is used to identify the putative loci. Since the number of quantitative trait loci (QTLs) is very likely to be small relative to the number of markers genotyped, a one-stage selective genotyping approach is commonly used to reduce the genotyping burden, whereby markers are genotyped solely on individuals with extreme trait values. This approach is powerful in the presence of a single quantitative trait locus (QTL) but may result in substantial loss of information in the presence of multiple QTLs. Here we investigate the efficiency of sequential two stage designs to identify QTLs in experimental populations. Our investigations for backcross and F2 crosses suggest that genotyping all the markers on 60% of the subjects in Stage 1 and genotyping the chromosomes significant at 20% level using additional subjects in Stage 2 and testing using all the subjects provides an efficient approach to identify the QTLs and utilizes only 70% of the genotyping burden relative to a one stage design, regardless of the heritability and genotyping density. Complex traits are a consequence of multiple QTLs conferring main effects as well as epistatic interactions. We propose a two-stage analytic approach where a single-locus genome scan is conducted in Stage 1 to identify promising chromosomes, and interactions are examined using the loci on these chromosomes in Stage 2. We examine settings under which the two-stage analytic approach provides sufficient power to detect the putative QTLs.

In genetic association studies, there is increasing interest in understanding the joint effects of genetic and nongenetic factors. For rare diseases, the case-control study is a practical design, and logistic regression is the standard method of inference. However, the power to detect statistical interaction is a concern, even with relatively large samples. Under independence of genetic and nongenetic covariates, improved precision of interaction estimators is possible, but logistic regression does not make use of this assumption and consequently is not statistically efficient. In recent work to improve efficiency, profile likelihood methods have been used to develop semi-parametric inference that incorporates the independence assumption. We describe an alternate derivation of these estimators for rare diseases that is based on classic arguments from case-control inference. These arguments lead to a simplification in the variance estimator. We also describe a strategy for relaxing the independence assumption. Under either independence or the proposed dependence model, inference for association parameters is conveniently obtained by fitting a conditional logistic regression. The statistical properties of the proposed methodology are investigated by simulation.

An objective of microarray data analysis is to identify gene expressions that are associated with a disease related outcome. For each gene, a test statistic is computed to determine if an association exists, and this statistic generates a marginal p-value. In an effort to pool this information across genes, a p-value density function is derived. The p-value density is modeled as a mixture of a uniform (0,1) density and a scaled ratio of normal densities derived from the asymptotic normality of the test statistic. The p-values are assumed to be weakly dependent and a quasi-likelihood is used to estimate the parameters in the mixture density. The quasi-likelihood and the weak dependence assumption enables estimation and asymptotic inference on the false discovery rate for a given rejection region, and its inverse, the p-value threshold parameter for a fixed false discovery rate. A false discovery rate analysis on a localized prostate cancer data set is used to illustrate the methodology. Simulations are performed to assess the performance of this methodology.

There have been various attempts to reconstruct gene regulatory networks from microarray expression data in the past. However, owing to the limited amount of independent experimental conditions and noise inherent in the measurements, the results have been rather modest so far. For this reason it seems advisable to include biological prior knowledge, related, for instance, to transcription factor binding locations in promoter regions or partially known signalling pathways from the literature. In the present paper, we consider a Bayesian approach to systematically integrate expression data with multiple sources of prior knowledge. Each source is encoded via a separate energy function, from which a prior distribution over network structures in the form of a Gibbs distribution is constructed. The hyperparameters associated with the different sources of prior knowledge, which measure the influence of the respective prior relative to the data, are sampled from the posterior distribution with MCMC. We have evaluated the proposed scheme on the yeast cell cycle and the Raf signalling pathway. Our findings quantify to what extent the inclusion of independent prior knowledge improves the network reconstruction accuracy, and the values of the hyperparameters inferred with the proposed scheme were found to be close to optimal with respect to minimizing the reconstruction error.

Statistical analysis of microarray gene expression data has recently attracted a great deal of attention. One problem of interest is to relate genes to survival outcomes of patients with the purpose of building regression models for the prediction of future patients' survival based on their gene expression data. For this, several authors have discussed the use of the proportional hazards or Cox model after reducing the dimension of the gene expression data. This paper presents a new approach to conduct the Cox survival analysis of microarray gene expression data with the focus on models' predictive ability. The method modifies the correlation principal component regression (Sun, 1995) to handle the censoring problem of survival data. The results based on simulated data and a set of publicly available data on diffuse large B-cell lymphoma show that the proposed method works well in terms of models' robustness and predictive ability in comparison with some existing partial least squares approaches. Also, the new approach is simpler and easy to implement.

We propose here a simple and robust approach for meta-analysis of molecular association studies. Making use of the binary structure of the data, and by treating the genotypes as independent variables in a logistic regression, we apply a simple and commonly used methodology that performs satisfactorily, being at the same time very flexible. We present simple tests for detecting heterogeneity and we describe a random effects extension of the method in order to allow for between studies heterogeneity. We derive also simple tests for assessing the most plausible genetic model of inheritance, and its between-studies heterogeneity as well as adjusting for covariates. The methodology introduced here is easily extended in cases with polytomous or continuous outcomes as well as in cases with more than two alleles. We apply the methodology in several published meta-analyses of genetic association studies with very encouraging results. The main advantages of the proposed methodology is its flexibility and the ease of use, while at the same time covers almost every aspect of a meta-analysis providing overall estimates without the need of multiple comparisons. We anticipate that this simple method would be used in the future in meta-analyses of genetic association studies. A STATA command performing all the available computations is available at http://bioinformatics.biol.uoa.gr/~pbagos/metagen/.

In a hidden Markov model, one "estimates" the state of the hidden Markov chain at t by computing via the forwards-backwards algorithm the conditional distribution of the state vector given the observed data. The covariance matrix of this conditional distribution measures the information lost by failure to observe directly the state of the hidden process. In the case where changes of state occur slowly relative to the speed at which information about the underlying state accumulates in the observed data, we compute approximately these covariances in terms of functionals of Brownian motion that arise in change-point analysis. Applications in gene mapping, where these covariances play a role in standardizing the score statistic and in evaluating the loss of noncentrality due to incomplete information, are discussed. Numerical examples illustrate the range of validity and limitations of our results.

Microarrays are a tool for measuring the expression levels of a large number of genes simultaneously. In the microarray experiment, however, many undesirable systematic variations are observed. Correct identification and removal of these variations is essential to allow the comparison of expression levels across experiments. We describe the use of linear mixed models for the normalization of two-color spotted microarrays for various sources of variation including printtip variation. Normalization with linear mixed models provides a parametric model of which results compare favorably to intensity dependent normalization LOWESS methods. We illustrate the use of this technique on two datasets. The first dataset contains 24 arrays, each with approximately 600 genes, replicated 3 times per array. A second dataset, coming from the apo AI experiment, was used to further illustrate the methods. Finally, a simulation study was done to compare between methods.

Microarray applications for the study of gene expression are becoming accessible for researchers in more and more systems. Applications from field or laboratory experiments are often complicated by the need to superimpose sample pairing for two-color arrays on experimental designs that may already be complex. For example, split-plot designs are commonly used in biological systems where experiments involve two types of treatments that are not readily applied at the same scale. We demonstrate how effects that are confounded with arrays can still be estimated when there is sufficient replication. To illustrate, we evaluate three methods of sample pairing superimposed on a split-plot design with two treatments, deriving the variance associated with parameter estimates for each. Design A has levels of the whole plot treatment paired on the same microarray within a level of the subplot treatment. Design B has crossed levels paired on the same microarray. Design C has levels of the treatment applied to subplots paired on the same microarray within a whole plot. Designs A and B have lower variance than design C for comparing the levels of the whole plot treatment. Designs B and C have lower variance for comparing the levels of the subplot treatment and design C has lower variance for comparing the levels of the subplot treatment within each level of the whole plot treatment. We provide SAS code for the analyses of variance discussed.

The problem of finding periodically expressed genes from time course microarray experiments is at the center of numerous efforts to identify the molecular components of biological clocks. We present a new approach to this problem based on the cyclohedron test, which is a rank test inspired by recent advances in algebraic combinatorics. The test has the advantage of being robust to measurement errors, and can be used to ascertain the significance of top-ranked genes. We apply the test to recently published measurements of gene expression during mouse somitogenesis and find 32 genes that collectively are significant. Among these are previously identified periodic genes involved in the Notch/FGF and Wnt signaling pathways, as well as novel candidate genes that may play a role in regulating the segmentation clock. These results confirm that there are an abundance of exceptionally periodic genes expressed during somitogenesis. The emphasis of this paper is on the statistics and combinatorics that underlie the cyclohedron test and its implementation within a multiple testing framework.

We propose an approximate distribution for the gapped local score of a two sequence comparison. Our method stands on combining an adapted scoring scheme that includes the gaps and an approximate distribution of the ungapped local score of two independent sequences of i.i.d. random variables. The new scoring scheme is defined on h-tuples of the sequences, using the gapped global score. The influence of h and the accuracy of the p-value are numerically studied and compared with obtained p-value of BLAST. The numerical experiments emphasize that our approximate p-values outperform the BLAST ones, particularly for both simulated and real short sequences.

Liquid Chromatography - Mass Spectrometry (LC-MS) is a powerful method for sensitive detection and quantification of proteins and peptides in complex biological fluids like serum. LC-MS produces complex data sets, consisting of some hundreds of millions of data points per sample at a resolution of 0.1 amu in the m/z domain and 7000 data points in the time domain. However, the detection of the lower abundance proteins from this data is hampered by the presence of artefacts, such as high frequency noise and spikes. Moreover, not all of the tens of thousands of the chromatograms produced per sample are relevant for the pursuit of the biomarkers. Thus in analysing the LC-MS data, two critical pre-processing issues arise. Which of the thousands of the: 1. chromatograms per sample are relevant for the detection of the biomarkers?, and 2. signals per chromatogram are truly compound-related? Each of these issues involves assessing the significance (deviation from noise) of multiple observations and the issue of multiple comparisons arises. Current methods disregard the multiplicity and provide no concrete threshold for significance. However, with such procedures, the probability of one or more false-positives is high as the number of tests to be performed is large, and must be controlled. Realizing that the cut-offs for declaring a chromatogram (or a signal) to be compound-related can hugely influence which proteins are detected, it seems natural to define thresholds that are neither arbitrary nor subjective. We suggest the choice of thresholds guided by the critical aim of controlling the False Discovery Rate (FDR) in multiple hypotheses testing for significance over a large set of features produced per sample. This involves the use of the regression diagnostics to characterize the signals of a chromatogram (e.g. as outliers or influential) and to suggest suitable tests statistics for the multiple testing procedures (MTP) for discriminating noise and spikes from true signals. The role of the Generalized Linear Models (GLM) in this MTP is investigated. The method is applied to LC-MS datasets from trypsin-digested serum spiked with varying levels of horse heart cytochrome C (cytoc).

The objective of the present paper is to develop a truly functional Bayesian method specifically designed for time series microarray data. The method allows one to identify differentially expressed genes in a time-course microarray experiment, to rank them and to estimate their expression profiles. Each gene expression profile is modeled as an expansion over some orthonormal basis, where the coefficients and the number of basis functions are estimated from the data. The proposed procedure deals successfully with various technical difficulties that arise in typical microarray experiments such as a small number of observations, non-uniform sampling intervals and missing or replicated data. The procedure allows one to account for various types of errors and offers a good compromise between nonparametric techniques and techniques based on normality assumptions. In addition, all evaluations are performed using analytic expressions, so the entire procedure requires very small computational effort. The procedure is studied using both simulated and real data, and is compared with competitive recent approaches. Finally, the procedure is applied to a case study of a human breast cancer cell line stimulated with estrogen. We succeeded in finding new significant genes that were not marked in an earlier work on the same dataset.

When trying to learn a model for the prediction of an outcome given a set of covariates, a statistician has many estimation procedures in their toolbox. A few examples of these candidate learners are: least squares, least angle regression, random forests, and spline regression. Previous articles (van der Laan and Dudoit (2003); van der Laan et al. (2006); Sinisi et al. (2007)) theoretically validated the use of cross validation to select an optimal learner among many candidate learners. Motivated by this use of cross validation, we propose a new prediction method for creating a weighted combination of many candidate learners to build the super learner. This article proposes a fast algorithm for constructing a super learner in prediction which uses V-fold cross-validation to select weights to combine an initial set of candidate learners. In addition, this paper contains a practical demonstration of the adaptivity of this so called super learner to various true data generating distributions. This approach for construction of a super learner generalizes to any parameter which can be defined as a minimizer of a loss function.

Dose-response studies are commonly used in experiments in pharmaceutical research in order to investigate the dependence of the response on dose, i.e., a trend of the response level toxicity with respect to dose. In this paper we focus on dose-response experiments within a microarray setting in which several microarrays are available for a sequence of increasing dose levels. A gene is called differentially expressed if there is a monotonic trend (with respect to dose) in the gene expression. We review several testing procedures which can be used in order to test equality among the gene expression means against ordered alternatives with respect to dose, namely Williams' (Williams 1971 and 1972), Marcus' (Marcus 1976), global likelihood ratio test (Bartholomew 1961, Barlow et al. 1972, and Robertson et al. 1988), and M (Hu et al. 2005) statistics. Additionally we introduce a modification to the standard error of the M statistic. We compare the performance of these five test statistics. Moreover, we discuss the issue of one-sided versus two-sided testing procedures. False Discovery Rate (Benjamni and Hochberg 1995, Ge et al. 2003), and resampling-based Familywise Error Rate (Westfall and Young 1993) are used to handle the multiple testing issue. The methods above are applied to a data set with 4 doses (3 arrays per dose) and 16,998 genes. Results on the number of significant genes from each statistic are discussed. A simulation study is conducted to investigate the power of each statistic. A R library IsoGene implementing the methods is available from the first author.

Multiple testing procedures are commonly used in gene expression studies for the detection of differential expression, where typically thousands of genes are measured over at least two experimental conditions. Given the need for powerful testing procedures, and the attendant danger of false positives in multiple testing, the False Discovery Rate (FDR) controlling procedure of Benjamini and Hochberg (1995) has become a popular tool. When simultaneously testing hypotheses, suppose that R rejections are made, of which Fp are false positives. The Benjamini and Hochberg procedure ensures that the expectation of Fp/R is bounded above by some pre-specified proportion. In practice, the procedure is applied to a single experiment. In this paper we investigate the across-experiment variability of the proportion Fp/R as a function of three experimental parameters. The operational characteristics of the procedure when applied to dependent hypotheses are also considered.

Microarray studies often need to simultaneously examine thousands of genes to determine which are differentially expressed. One main challenge in those studies is to find suitable multiple testing procedures that provide accurate control of the error rates of interest and meanwhile are most powerful, that is, they return the longest list of truly interesting genes among competitors. Many multiple testing methods have been developed recently for microarray data analysis, especially resampling based methods, such as permutation methods, the null-centered and scaled bootstrap (NCSB) method, and the quantile-transformed-bootstrap-distribution (QTBD) method. Each of these methods has its own merits and limitations. Theoretically permutation methods can fail to provide accurate control of Type I errors when the so-called subset pivotality condition is violated. The NCSB method does not suffer from that limitation, but an impractical number of bootstrap samples are often needed to get proper control of Type I errors. The newly developed QTBD method has the virtues of providing accurate control of Type I errors under few restrictions. However, the relative practical performance of the above three types of multiple testing methods remains unresolved. This paper compares the above three resampling based methods according to the control of family wise error rates (FWER) through data simulations. Results show that among the three resampling based methods, the QTBD method provides relatively accurate and powerful control in more general circumstances.

We investigate an important issue of a meta-algorithm for selecting variables in the framework of microarray data. This wrapper method starts from any classification algorithm and weights each variable (i.e. gene) relative to its efficiency for classification. An optimization procedure is then inferred which exhibits important genes for the studied biological process.Theory and application with the SVM classifier were presented in Gadat and Younes, 2007 and we extend this method with CART. The classification error rates are computed on three famous public databases (Leukemia, Colon and Prostate) and compared with those from other wrapper methods (RFE, lo norm SVM, Random Forests). This allows the assessment of the statistical relevance of the proposed algorithm. Furthermore, a biological interpretation with the Ingenuity Pathway Analysis software outputs clearly shows that the gene selections from the different wrapper methods raise very relevant biological information, compared to a classical filter gene selection with T-test.

The investigation of microbial communities is an essential part of the study of the biosphere. Flexible molecular fingerprinting tools such as terminal-restriction fragment length polymorphism (T-RFLP) analysis are often applied in the studies to enable the characterization of the microbial population. However, such data have so far been primarily analyzed using conventional clustering methods. Here we introduce a Bayesian model-based method for the purpose of comparing microbial communities using T-RFLP data. Such datasets have in general several challenging features, e.g. sparseness, missing values and structurally zero-valued observations. These features are taken into account by developing a Bayesian latent class mixture model for the observations in our framework. To make inferences under the model we use a recent Markov chain Monte Carlo (MCMC) -based method for the Bayesian model selection. To assess the introduced method we analyze both simulated and real datasets. The simulations show that our approach compares preferably to standard statistical clustering tools, such as k-means, hierarchical clustering, and Autoclass. The developed tool is freely available as a software package T-BAPS at http://www.abo.fi/fak/mnf/mate/jc/software/t-baps.html.

We describe a general model for pairwise microsatellite allele matching probabilities. The model can be used for analysis of population substructure, and is particularly focused on relating genetic correlation to measurable covariates. The approach is intended for cases when the existence of subpopulations is uncertain and a priori assignment of samples to hypothesized subpopulations is difficult. Such a situation arises, for example, with western Arctic bowhead whales, where genetic samples are available only from a possibly mixed migratory assemblage. We estimate genetic structure associated with spatial, temporal, or other variables that may confound the detection of population structure. In the bowhead case, the model permits detection of genetic patterns associated with a temporally pulsed multi-population assemblage in the annual migration. Hypothesis tests for population substructure and for covariate effects can be carried out using permutation methods. Simulated and real examples illustrate the effectiveness and reliability of the approach and enable comparisons with other familiar approaches. Analysis of the bowhead data finds no evidence for two temporally pulsed subpopulations using the best available data, although a significant pattern found by other researchers using preliminary data is also confirmed here. Code in the R language is available from www.stat.colostate.edu/~geof/gammmp.html.

Array Comparative Genomic Hybridization (aCGH) is an array-based technology which provides simultaneous spot assays of relative genetic abundance (RGA) levels at multiple sites across the genome. These spot assays are spatially correlated with respect to genomic location and, as a result, the univariate tests conducted using data generated from these spot assays are also spatially correlated. In the context of multiple hypothesis testing, this spatial correlation complicates the question of how best to define a `discovery' and consequently, how best to estimate the false discovery rate (FDR) corresponding to a given rejection region.One can quantify the number of discoveries as the total number of spots for which the spot-based univariate test statistic falls within a given rejection region. Under this spot-based method, separate but correlated discoveries are identified. We show via a simulation study that the method of Benjamini and Hochberg (1995) can provide a reasonable estimate of the spot-wise FDR, but these results require that the simulated spot assays are categorized as true or false discoveries in a particular way. However, laboratory researchers may actually be interested in estimating a `regional' FDR, rather than a `local' spot-wise FDR. We describe an example of such circumstances, and present a method for estimating the (chromosome) arm-wise False Discovery Rate. In this framework, one can quantify the number of discoveries as the total number of chromosome arms for which at least one spot-based test statistic falls into a given rejection region. Defining the discoveries in this way, both the biological and testing objectives coincide. We provide results from a series of simulations which involved the analysis of preferentially re-sampled spot assay values from a real aCGH dataset.

Characterizing genetic variability in the human pathogenic Plasmodium species, the group of parasites that cause Malaria, may have broad global health implications. Specifically, discerning the combinations of mutations that lead to viable parasites and the population level frequencies of these clonal sequences will allow for targeted vaccine development and individualized treatment choices. This presents an analytical challenge, however, since haplotypic phase (i.e. the alignment of bases on a single DNA strand) is generally unobservable in multiply infected individuals. This manuscript describes an expectation maximization (EM) approach to maximum likelihood estimation of haplotype frequencies in this missing data setting. The approach is applied to a cohort of N=341 malaria infected children in Uganda, Cameroon and Sudan to characterize regional differences. A simulation study is also presented to characterize method performance and assess sensitivity to distributional assumptions.

Digital images obtained by the laser scanning of spotted microarrays often include saturated pixel values. These arise when the scan settings are sufficiently high and some pixels exceed the limit L=65535 and are instead set to L. Failure to adjust for this censoring leads to biased estimates of gene expression levels. To impute censored values, we propose a linear model based on the principal components of uncensored spots on the same array. This is computationally fast, flexible to adapt to distinctive spot shapes and profiles on different arrays, and is shown to be more effective than the polynomial-hyperbolic model in correcting for the bias. The application to biological data demonstrates the potential for enhancing the dynamic range of detection. Fortran90 subroutines implementing these methods are available at http://www.bioss.ac.uk/~chris.

Maximum Likelihood (ML) is used as a standard method for estimating divergence times in phylogenetic trees. The method is consistent and hence the precision can be improved by analyzing longer sequences. In this paper we show that the precision can be improved also by including more taxa to the existing tree. It is a theoretical study, complemented with simulations, showing that the gain in precision is faster with increasing sequence length than with increasing number of taxa.We further compare the results of estimating divergence times using Maximum Likelihood with the much faster and less complex estimation method of Mean Path Length (MPL), which works with the evolution model of Jukes-Cantor (1969). It is shown that MPL is as good as ML in estimating divergence times of nodes that are located near the root in the tree, but ML is better in estimating the divergence times of nodes lower down.

We present a Bayesian hierarchical model for detecting differentially expressed genes using a mixture prior on the parameters representing differential effects. We formulate an easily interpretable 3-component mixture to classify genes as over-expressed, under-expressed and non-differentially expressed, and model gene variances as exchangeable to allow for variability between genes. We show how the proportion of differentially expressed genes, and the mixture parameters, can be estimated in a fully Bayesian way, extending previous approaches where this proportion was fixed and empirically estimated. Good estimates of the false discovery rates are also obtained.Different parametric families for the mixture components can lead to quite different classifications of genes for a given data set. Using Affymetrix data from a knock out and wildtype mice experiment, we show how predictive model checks can be used to guide the choice between possible mixture priors. These checks show that extending the mixture model to allow extra variability around zero instead of the usual point mass null fits the data better.A software package for R is available.

Most genetic diseases are complex, i.e. associated to combinations of SNPs rather than individual SNPs. In the last few years, this topic has often been addressed in terms of SNP-SNP interaction patterns given as expressions linked by logical operators. Methods for multiple testing in high-dimensional settings can be applied when many SNPs are considered simultaneously. However, another less well-known multiple testing problem arises within a fixed subset of SNPs when the logic expression is chosen optimally. In this article, we propose a general asymptotic approach for deriving the distribution of the maximally selected chi-square statistic in various situations. We show how this result can be used for testing logic expressions - in particular SNP-SNP interaction patterns - while controlling for multiple comparisons. Simulations show that our method provides multiple testing adjustments when the logic expression is chosen such as to maximize the statistic. Its benefit is demonstrated through an application to a real dataset from a large population-based study considering allergy and asthma in KORA. An implementation of our method is available from the Comprehensive R Archive Network (CRAN) as R package 'SNPmaxsel'.