Volume 3, Issue 1

The increase in the number of available markers for many experimental populations has led to QTL studies with ever increasing marker numbers and densities. The resulting conundrum is that as marker density increases, so does the multiple testing problem. It is important to re-examine the detection of multiple QTL in light of increasing marker density. We explore through simulation whether existing methods have achieved the maximum possible power for detecting multiple QTL and whether increasing the marker density is an effective strategy for locating multiple QTL. In addition to existing methods, such as the maximum, the CET, and the Benjamini-Hochberg and Benjamini-Yekutieli procedures, we propose and evaluate the complete set of order statistics with their corresponding empirical joint distribution. We examine these statistics in conjunction with a novel application of the alpha-spending approach, providing a less conservative solution to the problem of controlling the false discovery rate (FDR) in multiple tests. We conducted a simulation study to assess the relative power of these approaches as well as their ability to control FDR. We find that several of the new approaches have a reasonable FDR, and can substantially improve the experimenter's ability to detect multiple QTL compared to existing approaches in many cases; however, the Benjamini-Hochberg procedure remains a very reasonable choice. The methods are applied to a nine-trait Oat vernalization dataset.

The problem of relating genotype (as represented by amino acid sequence) to phenotypes is distinguished from standard regression problems by the nature of sequence data. Here we investigate an instance of such a problem where the phenotype of interest is HIV-1 replication capacity and contiguous segments of protease and reverse transcriptase sequence constitutes genotype. A variety of data analytic methods have been proposed in this context. Shortcomings of select techniques are contrasted with the advantages afforded by tree-structured methods. However, tree-structured methods, in turn, have been criticized on grounds of only enjoying modest predictive performance. A number of ensemble approaches (bagging, boosting, random forests) have recently emerged, devised to overcome this deficiency. We evaluate random forests as applied in this setting, and detail why prediction gains obtained in other situations are not realized. Other approaches including logic regression, support vector machines and neural networks are also applied. We interpret results in terms of HIV-1 reverse transcriptase structure and function.

The problem of identifying differentially expressed genes in designed microarray experiments is considered. Lonnstedt and Speed (2002) derived an expression for the posterior odds of differential expression in a replicated two-color experiment using a simple hierarchical parametric model. The purpose of this paper is to develop the hierarchical model of Lonnstedt and Speed (2002) into a practical approach for general microarray experiments with arbitrary numbers of treatments and RNA samples. The model is reset in the context of general linear models with arbitrary coefficients and contrasts of interest. The approach applies equally well to both single channel and two color microarray experiments. Consistent, closed form estimators are derived for the hyperparameters in the model. The estimators proposed have robust behavior even for small numbers of arrays and allow for incomplete data arising from spot filtering or spot quality weights. The posterior odds statistic is reformulated in terms of a moderated t-statistic in which posterior residual standard deviations are used in place of ordinary standard deviations. The empirical Bayes approach is equivalent to shrinkage of the estimated sample variances towards a pooled estimate, resulting in far more stable inference when the number of arrays is small. The use of moderated t-statistics has the advantage over the posterior odds that the number of hyperparameters which need to estimated is reduced; in particular, knowledge of the non-null prior for the fold changes are not required. The moderated t-statistic is shown to follow a t-distribution with augmented degrees of freedom. The moderated t inferential approach extends to accommodate tests of composite null hypotheses through the use of moderated F-statistics. The performance of the methods is demonstrated in a simulation study. Results are presented for two publicly available data sets.

Likelihood-based cross-validation is a statistical tool for selecting a density estimate based on n i.i.d. observations from the true density among a collection of candidate density estimators. General examples are the selection of a model indexing a maximum likelihood estimator, and the selection of a bandwidth indexing a nonparametric (e.g. kernel) density estimator. In this article, we establish a finite sample result for a general class of likelihood-based cross-validation procedures (as indexed by the type of sample splitting used, e.g. V-fold cross-validation). This result implies that the cross-validation selector performs asymptotically as well (w.r.t. to the Kullback-Leibler distance to the true density) as a benchmark model selector which is optimal for each given dataset and depends on the true density. Crucial conditions of our theorem are that the size of the validation sample converges to infinity, which excludes leave-one-out cross-validation, and that the candidate density estimates are bounded away from zero and infinity. We illustrate these asymptotic results and the practical performance of likelihood-based cross-validation for the purpose of bandwidth selection with a simulation study. Moreover, we use likelihood-based cross-validation in the context of regulatory motif detection in DNA sequences.

In this article we describe and discuss implementation of a weighted simulation procedure, importance sampling, in the context of nonparametric linkage analysis. The objective is to estimate genome-wide p-values, i.e. the probability that the maximal linkage score exceeds given thresholds under the null hypothesis of no linkage. In order to reduce variance of the estimate for large thresholds, we simulate linkage scores under a distribution different from the null with an artificial disease locus positioned somewhere along the genome. To compensate for the fact that we simulate under the wrong distribution, the simulated scores are reweighted using a certain likelihood ratio. If the sampling distribution are properly chosen the variance of the corresponding estimate is reduced. This results in accurate genome-wide p-value estimates for a wide range of large thresholds with a substantially smaller cost adjusted relative efficiency with respect to standard unweighted simulation.We illustrate the performance of the method for several pedigree examples, discuss implementation including the amount of variance reduction and describe some possible generalizations.

Nuclear Magnetic Resonance (NMR) spectroscopy is a key experimental technique used to study protein structure, dynamics, and interactions. NMR methods face the bottleneck of spectral analysis, in particular determining the resonance assignments, which help define the mapping between atoms in the protein and peaks in the spectra. A substantial amount of noise in spectral data, along with ambiguities in interpretation, make this analysis a daunting task, and there exists no generally accepted measure of uncertainty associated with the resulting solutions. This paper develops a model-based inference approach that addresses the problem of characterizing uncertainty in backbone resonance assignment. We argue that NMR spectra are subject to random variation, and ignoring this stochasticity can lead to false optimism and erroneous conclusions. We propose a Bayesian statistical model that accounts for various sources of uncertainty and provides an automatable framework for inference. While assignment has previously been viewed as a deterministic optimization problem, we demonstrate the importance of considering all solutions consistent with the data, and develop an algorithm to search this space within our statistical framework. Our approach is able to characterize the uncertainty associated with backbone resonance assignment in several ways: 1) it quantifies of uncertainty in the individually assigned resonances in terms of their posterior standard deviations; 2) it assesses the information content in the data with a posterior distribution of plausible assignments; and 3) it provides a measure of the overall plausibility of assignments. We demonstrate the value of our approach in a study of experimental data from two proteins, Human Ubiquitin and Cold-shock protein A from E. coli. In addition, we provide simulations showing the impact of experimental conditions on uncertainty in the assignments.

Given a multiple testing situation, the null hypotheses that appear to have sufficiently low probabilities of truth may be rejected using a simple, nonparametric method based on decision theory. This applies not only to posterior levels of belief, but also to conditional probabilities in the sense of relative frequencies, as seen from their equality to local false discovery rates (dFDRs). This approach neither requires the estimation of probability densities, nor of their ratios. Decision theory can also inform the selection of false discovery rate weights. An application to gene expression microarrays is presented with a discussion of the applicability of the assumption of "clumpy dependence."

Motivation: Standard laboratory classification of the plasma cell dyscrasia monoclonal gammopathy of undetermined significance (MGUS) and the overt plasma cell neoplasm multiple myeloma (MM) is quite accurate, yet, for the most part, biologically uninformative. Most, if not all, cancers are caused by inherited or acquired genetic mutations that manifest themselves in altered gene expression patterns in the clonally related cancer cells. Microarray technology allows for qualitative and quantitative measurements of the expression levels of thousands of genes simultaneously, and it has now been used both to classify cancers that are morphologically indistinguishable and to predict response to therapy. It is anticipated that this information can also be used to develop molecular diagnostic models and to provide insight into mechanisms of disease progression, e.g., transition from healthy to benign hyperplasia or conversion of a benign hyperplasia to overt malignancy. However, standard data analysis techniques are not trivial to employ on these large data sets. Methodology designed to handle large data sets (or modified to do so) is needed to access the vital information contained in the genetic samples, which in turn can be used to develop more robust and accurate methods of clinical diagnostics and prognostics.Results: Here we report on the application of a panel of statistical and data mining methodologies to classify groups of samples based on expression of 12,000 genes derived from a high density oligonucleotide microarray analysis of highly purified plasma cells from newly diagnosed MM, MGUS, and normal healthy donors. The three groups of samples are each tested against each other. The methods are found to be similar in their ability to predict group membership; all do quite well at predicting MM vs. normal and MGUS vs. normal. However, no method appears to be able to distinguish explicitly the genetic mechanisms between MM and MGUS. We believe this might be due to the lack of genetic differences between these two conditions, and may not be due to the failure of the models. We report the prediction errors for each of the models and each of the methods. Additionally, we report ROC curves for the results on group prediction.Availability: Logistic regression: standard software, available, for example in SAS. Decision trees and boosted trees: C5.0 from www.rulequest.com. SVM: SVM-light is publicly available from svmlight.joachims.org. Naïve Bayes and ensemble of voters are publicly available from www.biostat.wisc.edu/~mwaddell/eov.html. Nearest Shrunken Centroids is publicly available from http://www-stat.stanford.edu/~tibs/PAM.

We present a mathematical model to extend the dynamic range of gene expression data measured by laser scanners. The strategy is based on the rather simple but novel idea of producing two images with different scanner sensitivities, obtaining two different sets of expression values: the first is a low-sensitivity measure to obtain high expression values which would be saturated in a high-sensitivity measure; the second, by the converse strategy, obtains additional information about the low-expression levels. Two mathematical models based on linear and gamma curves are presented for relating the two measurements to each other and producing a coherent and extended range of values. The procedure minimizes the quantization relative error and avoids the collateral effects of saturation. Since most of the current scanner devices are able to adjust the saturation level, the strategy can be considered as a universal solution, and not dependent on the image processing software used for reading the DNA chip. Various tests have been performed, on both proprietary and public domain data sets, showing a reduction of the saturation and quantization effects, not achievable by other methods, with a more complete description of gene-expression data and with a reasonable computational complexity.

Various discriminant methods have been applied for classification of tumors based on gene expression profiles, among which the nearest neighbor (NN) method has been reported to perform relatively well. Usually cross-validation (CV) is used to select the neighbor size as well as the number of variables for the NN method. However, CV can perform poorly when there is considerable uncertainty in choosing the best candidate classifier. As an alternative to selecting a single ``winner,'' we propose a weighting method to combine the multiple NN rules. Four gene expression data sets are used to compare its performance with CV methods. The results show that when the CV selection is unstable, the combined classifier performs much better.

The present article proposes general single-step multiple testing procedures for controlling Type I error rates defined as arbitrary parameters of the distribution of the number of Type I errors, such as the generalized family-wise error rate. A key feature of our approach is the test statistics null distribution (rather than data generating null distribution) used to derive cut-offs (i.e., rejection regions) for these test statistics and the resulting adjusted p-values. For general null hypotheses, corresponding to submodels for the data generating distribution, we identify an asymptotic domination condition for a null distribution under which single-step common-quantile and common-cut-off procedures asymptotically control the Type I error rate, for arbitrary data generating distributions, without the need for conditions such as subset pivotality. Inspired by this general characterization of a null distribution, we then propose as an explicit null distribution the asymptotic distribution of the vector of null value shifted and scaled test statistics. In the special case of family-wise error rate (FWER) control, our method yields the single-step minP and maxT procedures, based on minima of unadjusted p-values and maxima of test statistics, respectively, with the important distinction in the choice of null distribution. Single-step procedures based on consistent estimators of the null distribution are shown to also provide asymptotic control of the Type I error rate. A general bootstrap algorithm is supplied to conveniently obtain consistent estimators of the null distribution. The special cases of t- and F-statistics are discussed in detail. The companion articles focus on step-down multiple testing procedures for control of the FWER (van der Laan et al., 2004b) and on augmentations of FWER-controlling methods to control error rates such as tail probabilities for the number of false positives and for the proportion of false positives among the rejected hypotheses (van der Laan et al., 2004a). The proposed bootstrap multiple testing procedures are evaluated by a simulation study and applied to genomic data in the fourth article of the series (Pollard et al., 2004).

The present article proposes two step-down multiple testing procedures for asymptotic control of the family-wise error rate (FWER): the first procedure is based on maxima of test statistics (step-down maxT), while the second relies on minima of unadjusted p-values (step-down minP). A key feature of our approach is the characterization and construction of a test statistics null distribution (rather than data generating null distribution) for deriving cut-offs for these test statistics (i.e., rejection regions) and the resulting adjusted p-values. For general null hypotheses, corresponding to submodels for the data generating distribution, we identify an asymptotic domination condition for a null distribution under which the step-down maxT and minP procedures asymptotically control the Type I error rate, for arbitrary data generating distributions, without the need for conditions such as subset pivotality. Inspired by this general characterization, we then propose as an explicit null distribution the asymptotic distribution of the vector of null value shifted and scaled test statistics. Step-down procedures based on consistent estimators of the null distribution are shown to also provide asymptotic control of the Type I error rate. A general bootstrap algorithm is supplied to conveniently obtain consistent estimators of the null distribution.

This article shows that any single-step or stepwise multiple testing procedure (asymptotically) controlling the family-wise error rate (FWER) can be augmented into procedures that (asymptotically) control tail probabilities for the number of false positives and the proportion of false positives among the rejected hypotheses. Specifically, given any procedure that (asymptotically) controls the FWER at level alpha, we propose simple augmentation procedures that provide (asymptotic) level-alpha control of: (i) the generalized family-wise error rate, i.e., the tail probability, gFWER(k), that the number of Type I errors exceeds a user-supplied integer k, and (ii) the tail probability, TPPFP(q), that the proportion of Type I errors among the rejected hypotheses exceeds a user-supplied value 0

We present a statistical approach to scoring changes in activity of metabolic pathways from gene expression data. The method identifies the biologically relevant pathways with corresponding statistical significance. Based on gene expression data alone, only local structures of genetic networks can be recovered. Instead of inferring such a network, we propose a hypothesis-based approach. We use given knowledge about biological networks to improve sensitivity and interpretability of findings from microarray experiments.Recently introduced methods test if members of predefined gene sets are enriched in a list of top-ranked genes in a microarray study. We improve this approach by defining scores that depend on all members of the gene set and that also take pairwise co-regulation of these genes into account. We calculate the significance of co-regulation of gene sets with a nonparametric permutation test. On two data sets the method is validated and its biological relevance is discussed. It turns out that useful measures for co-regulation of genes in a pathway can be identified adaptively.We refine our method in two aspects specific to pathways. First, to overcome the ambiguity of enzyme-to-gene mappings for a fixed pathway, we introduce algorithms for selecting the best fitting gene for a specific enzyme in a specific condition. In selected cases, functional assignment of genes to pathways is feasible. Second, the sensitivity of detecting relevant pathways is improved by integrating information about pathway topology. The distance of two enzymes is measured by the number of reactions needed to connect them, and enzyme pairs with a smaller distance receive a higher weight in the score calculation.

We propose a family-based association test, FBAT-PC, for studies with quantitative traits that are measured repeatedly. The traits may be influenced by partially or completely unknown factors that may vary for each measurement. Using generalized principal component analysis, FBAT-PC amplifies the genetic effects of each measurement by constructing an overall phenotype with maximal heritability. Analytically, and in the simulation studies, we compare FBAT-PC with standard methodology and assess both the heritability of the overall phenotype and the power of FBAT-PC. Compared to univariate analysis, FBAT-PC achieves power gains of up to 200%. Applications of FBAT-PC to an osteoporosis study and to an asthma study show the practical relevance of FBAT-PC. FBAT-PC has been implemented in the software package PBAT and is freely available at http://www.biostat.harvard.edu/~clange/default.htm.

van der Laan and Dudoit (2003) provide a road map for estimation and performance assessment where a parameter of interest is defined as the risk minimizer for a suitable loss function and candidate estimators are generated using a loss function. After briefly reviewing this approach, this article proposes a general deletion/substitution/addition algorithm for minimizing, over subsets of variables (e.g., basis functions), the empirical risk of subset-specific estimators of the parameter of interest. This algorithm provides us with a new class of loss-based cross-validated algorithms in prediction of univariate outcomes, which can be extended to handle multivariate outcomes, conditional density and hazard estimation, and censored outcomes such as survival. In the context of regression, using polynomial basis functions, we study the properties of the deletion/substitution/addition algorithm in simulations and apply the method to detect transcription factor binding sites in yeast gene expression experiments.

We present a new approach to molecular classification based on mRNA comparisons. Our method, referred to as the top-scoring pair(s) (TSP) classifier, is motivated by current technical and practical limitations in using gene expression microarray data for class prediction, for example to detect disease, identify tumors or predict treatment response. Accurate statistical inference from such data is difficult due to the small number of observations, typically tens, relative to the large number of genes, typically thousands. Moreover, conventional methods from machine learning lead to decisions which are usually very difficult to interpret in simple or biologically meaningful terms. In contrast, the TSP classifier provides decision rules which i) involve very few genes and only relative expression values (e.g., comparing the mRNA counts within a single pair of genes); ii) are both accurate and transparent; and iii) provide specific hypotheses for follow-up studies. In particular, the TSP classifier achieves prediction rates with standard cancer data that are as high as those of previous studies which use considerably more genes and complex procedures. Finally, the TSP classifier is parameter-free, thus avoiding the type of over-fitting and inflated estimates of performance that result when all aspects of learning a predictor are not properly cross-validated.

There are several important issues to be addressed for gene expression temporal patterns' analysis: first, the correlation structure of multidimensional temporal data; second, the numerous sources of variations with existing high level noise; and last, gene expression mostly involves heterogeneous multiple dynamic patterns. We propose a Hierarchical Bayesian Neural Network model to account for the input correlations of time course gene array data. The variations in absolute gene expression levels and the noise can be estimated with the hierarchical Bayesian setting. The network parameters and the hyperparameters were simultaneously optimized with Monte Carlo Markov Chain simulation. Results show that the proposed model and algorithm can well capture the dynamic feature of gene expression temporal patterns despite the high noise levels, the highly correlated inputs, the overwhelming interactions, and other complex features typically present in microarray data. We test and demonstrate the proposed models with yeast cell cycle temporal data sets. The model performance of Hierarchical Bayesian Neural Network was compared to other popular machine learning methods such as Nearest Neighbor, Support Vector Machine, and Self Organized Map.

A genome-wide location analysis method has been introduced as a means to simultaneously study protein-DNA binding interactions for a large number of genes on a microarray platform. Identification of interactions between transcription factors (TF) and genes provide insight into the mechanisms that regulate a variety of cellular responses. Drawing proper inferences from the experimental data is key to finding statistically significant TF-gene binding interactions. We describe how the analysis and interpretation of genome-wide location data can be fit into a traditional statistical modeling framework that considers the data across all arrays and formulizes appropriate hypothesis tests. The approach is illustrated with data from a yeast transcription factor binding experiment that illustrates how identified TF-gene interactions can enhance initial exploration of transcriptional regulatory networks. Examples of five kinds of transcriptional regulatory structure are also demonstrated. Some stark differences with previously published results are explored.

Under the assumption that a significant motivation for sequencing the genomes of mammals is the resulting ability to help us locate and characterize functional DNA segments shared with humans, we have developed a statistical analysis to quantify the expected advantage. Examining uncertainty in terms of the width of a confidence interval, we show that uncertainty in the rate of nucleotide mutation can be shrunk by a factor of nearly four when nine mammals; human, chimpanzee, baboon, cat, dog, cow, pig, rat, mouse; are used instead of just two; human and mouse. Contrastingly, we show confidence interval shrinkage by a factor of only 1.5 for measurements of the distribution of nucleotides at an aligned sequence site. These additional genomes should greatly help in identifying conserved DNA sites, but would be much less effective at precisely describing the expected pattern of nucleotides at those sites.

A common practice in protein sequence alignment is to try several scoring matrices until ``something interesting'' is found. This leads to a multiple testing problem making p- and E-values hard to interpret. We focus on local alignment and propose to use logistic copula functions to model explicitly the dependence structure of scores obtained using different scoring matrices. By doing this, we obtain p-value correction factors when using more than one scoring matrix on the same sequences. Furthermore the parameter of the logistic copula can be interpreted as measure of dependence, providing insight concerning the relatedness of the scores from different matrices.

Phenotype and/or genotype misclassification can: significantly increase type II error probabilities for genetic case/control association, causing decrease in statistical power; and produce inaccurate estimates of population frequency parameters. We present a method, the likelihood ratio test allowing for errors (LRTae) that incorporates double-sample information for phenotypes and/or genotypes on a sub-sample of cases/controls. Population frequency parameters and misclassification probabilities are determined using a double-sample procedure as implemented in the Expectation-Maximization (EM) method. We perform null simulations assuming a SNP marker or a 4-allele (multi-allele) marker locus. To compare our method with the standard method that makes no adjustment for errors (LRTstd), we perform power simulations using a 2^k factorial design with high and low settings of: case/control samples, phenotype/genotype costs, double-sampled phenotypes/genotypes costs, phenotype/genotype error, and proportions of double-sampled individuals. All power simulations are performed fixing equal costs for the LRTstd and LRTae methods. We also consider case/control ApoE genotype data for an actual Alzheimer's study.The LRTae method maintains correct type I error proportions for all null simulations and all significance level thresholds (10%, 5%, 1%). LRTae average estimates of population frequencies and misclassification probabilities are equal to the true values, with variances of 10e-7 to 10e-8. For power simulations, the median power difference LRTae-LRTstd at the 5% significance level is 0.06 for multi-allele data and 0.01 for SNP data. For the ApoE data example, the LRTae and LRTstd p-values are 5.8 x 10e-5 and 1.6 x 10e-3, respectively. The increase in significance is due to adjustment in the LRTae for misclassification of the most commonly reported risk allele. We have developed freely available software that performs our LRTae statistic.

Estimation of the association between haplotypes and disease from a case-control study is considered. Assuming a single ``disease haplotype'' leads to the increased risk, attention focusses on the relative risks associated with a single copy, or two copies of the disease haplotype, relative to individuals with no copies. In this setting, case frequencies of the haplotype pairs are in Hardy-Weinberg Equilibrium (HWE) only if the combined influence of the two copies of the disease haplotype on risk is multiplicative. Thus, imputation cannot rely on the assumption of HWE for cases. A method is presented for obtaining estimates of the relative risks, making use of the EM algorithm and the assumption of HWE only for controls. The method accounts for the additional variation in the estimates due to the imputation of expected frequencies of haplotype pairs from ambiguous genotypes. A simulation study shows that the resulting confidence intervals have nominal coverage, and that the methods based on the assumption of HWE for both cases and controls can lead to bias.

The assessment of residue conservation in a multiple sequence alignment is a central issue in bioinformatics. Conserved residues and regions are used to determine structural and functional motifs or evolutionary relationships between the sequences of a multiple sequence alignment. For this reason, residue conservation is a valuable measure for database and motif search or for estimating the quality of alignments. In this paper, we present statistical methods for identifying conserved residues in multiple sequence alignments. While most earlier studies examine the positional conservation of the alignment, we focus on the detection of individual conserved residues at a position. The major advantages of multiple comparison methods originate from their ability to select conserved residues simultaneously and to consider the variability of the residue estimates. Large-scale simulations were used for the comparative analysis of the methods. Practical performance was studied by comparing the structurally and functionally important residues of Src homology 2 (SH2) domains to the assignments of the conservation indices. The applicability of the indices was also compared in three additional protein families comprising different degrees of entropy and variability in alignment positions. The results indicate that statistical multiple comparison methods are sensitive and reliable in identifying conserved residues.

Cross-study validation of gene expression investigations is critical in genomic analysis. We developed an R package and associated object definitions to merge and visualize multiple gene expression datasets. Our merging functions use arbitrary character IDs and generate objects that can efficiently support a variety of joint analyses. Visualization tools support exploration and cross-study validation of the data, without requiring normalization across platforms. Tools include ``integrative correlation'' plots that is, scatterplots of all pairwise correlations in one study against the corresponding pairwise correlations of another, both for individual genes and all genes combined. Gene-specific plots can be used to identify genes whose changes are reliably measured across studies. Visualizations also include scatterplots of gene-specific statistics quantifying relationships between expression and phenotypes of interest, using linear, logistic and Cox regression.

Gametic models for fitting breeding values at QTL as random effects in outbred populations have become popular because they require few assumptions about the number and distribution of QTL alleles segregating. The covariance matrix of the gametic effects has an inverse that is sparse and can be constructed rapidly by a simple algorithm, provided that all individuals have marker data, but not otherwise. An equivalent model, in which the joint distribution of QTL breeding values and marker genotypes is considered, was shown to generate a covariance matrix with a sparse inverse that can be constructed rapidly with a simple algorithm. This result makes more feasible including QTL as random effects in analyses of large pedigrees for QTL detection and marker assisted selection. Such analyses often use algorithms that rely upon sparseness of the mixed model equations and require the inverse of the covariance matrix, but not the covariance matrix itself. With the proposed model, each individual has two random effects for each possible unordered marker genotype for that individual. Therefore, individuals with marker data have two random effects, just as with the gametic model. To keep the notation and the derivation simple, the method is derived under the assumptions of a single linked marker and that the pedigree does not contain loops. The algorithm could be applied, as an approximate method, to pedigrees that contain loops.

With the rapid growth of entire genome data, reconstructing the phylogenetic relationship among different genomes has become a hot topic in comparative genomics. Maximum likelihood approach is one of the various approaches, and has been very successful. However, there is no reported study for any applications in the genome tree-making mainly due to the lack of an analytical form of a probability model and/or the complicated calculation burden. In this paper we studied the mathematical structure of the stochastic model of genome evolution, and then developed a simplified likelihood function for observing a specific phylogenetic pattern under four genome situation using gene content information. We use the maximum likelihood approach to identify phylogenetic trees. Simulation results indicate that the proposed method works well and can identify trees with a high correction rate. Real data application provides satisfied results. The approach developed in this paper can serve as the basis for reconstructing phylogenies of more than four genomes.

A common experimental strategy utilizing microarrays is to develop a signature of genes responding to some treatment in a model system, and then ask whether the same genes respond in an analogous way in a more natural and uncontrolled environment. In statistical terms, the question posed is whether genes score similarly on some statistical test in two independent data sets. Approaches to this problem ignoring gene/gene correlations common to all microarray data sets are known to give overstated statistical confidence levels. Permutation approaches have been proposed to give more accurate confidence levels, but can not be applied when sample sizes are small. Here we argue that the product moment correlation between test statistics in the two experiments is an ideal measure for summarizing concordance between the experiments, as confidence levels accounting for intergene correlations depend only on a single number -- the average squared correlation between gene pairs in the data set. The resulting null standard deviation is shown to vary by less than a factor of two over six distinct experimental data sets, suggesting that a universal constant may be used for this quantity. We show how a hidden assumption of the permutation approach may lead to incorrect p-values, while the analytic approach presented here is shown to be resistant to this assumption.

Partial Least Squares (PLS) dimension reduction is known to give good prediction accuracy in the context of classification with high-dimensional microarray data. In this paper, the classification procedure consisting of PLS dimension reduction and linear discriminant analysis on the new components is compared with some of the best state-of-the-art classification methods. Moreover, a boosting algorithm is applied to this classification method. In addition, a simple procedure to choose the number of PLS components is suggested. The connection between PLS dimension reduction and gene selection is examined and a property of the first PLS component for binary classification is proved. In addition, we show how PLS can be used for data visualization using real data. The whole study is based on 9 real microarray cancer data sets.

We present a series of statistical solutions to challenges that commonly arise in the production and analysis of genomic tag libraries. Tag libraries are collections of fragments of DNA or RNA, with each unique fragment often present in millions or billions of copies. Inferences can be made from data obtained by sequencing a subset of the library. The statistical approaches outlined in this paper are divided into three parts. First, we demonstrate the application of classical capture-recapture theory to the question of library complexity, i.e. the number of unique fragments in the library. Simulation studies verify the accuracy, for sample sizes of magnitudes typical in genomic studies, of the formulas we use to make our estimates. Second, we present a straightforward statistical cost analysis of tag experiments designed to uncover either disease-causing pathogens or new genes. Third, we develop a hidden Markov model approach to karyotyping a sample using a tag library derived from the sample's genomic DNA. While the resolution of the approach depends upon the number of tags sequenced from the library, we show via simulation that copy number alterations can be reliably detected for lengths as small as 1 Mb, even when a moderate number of tags are sequenced. Simulations predict very good specificity as well. Finally, all three of our approaches are applied to data from real tag library experiments. The hidden Markov model results are in line with what was expected from simulation, and genomic alterations found by applying the method to a cancer cell line library are confirmed using PCR.The methods and data described in this paper are contained in an R package, tagAnalysis, freely available at http://meyerson.dfci.harvard.edu/~tl974/tagAnalysis.

Recent analyses have shown that the relationship between intensity measurements from high density oligonucleotide microarrays and known concentration is non linear. Thus many measurements of so-called gene expression are neither measures of transcript nor mRNA concentration as might be expected.Intensity as measured in such microarrays is a measurement of fluorescent dye attached to probe-target duplexes formed during hybridization of a sample to the probes on the microarray. We develop several dynamic adsorption models relating fluorescent dye intensity to target RNA concentration, the simplest of which is the equilibrium Langmuir isotherm, or hyperbolic response function. Using data from the Affymerix HG-U95A Latin Square experiment, we evaluate various physical models, including equilibrium and non-equilibrium models, by applying maximum likelihood methods. We show that for these data, equilibrium Langmuir isotherms with probe dependent parameters are appropriate. We describe how probe sequence information may then be used to estimate the parameters of the Langmuir isotherm in order to provide an improved measure of absolute target concentration.

The methodological advancement in microarray data analysis on the basis of false discovery rate (FDR) control, such as the q-value plots, allows the investigator to examine the FDR from several perspectives. However, when FDR control at the ``customary" levels 0.01, 0.05, or 0.1 does not provide fruitful findings, there is little guidance for making the trade off between the significance threshold and the FDR level by sound statistical or biological considerations. Thus, meaningful statistical significance criteria that complement the existing FDR methods for large-scale multiple tests are desirable. Three statistical significance criteria, the profile information criterion, the total error proportion, and the guide-gene driven selection, are developed in this research. The first two are general significance threshold criteria for large-scale multiple tests; the profile information criterion is related to the recent theoretical studies of the connection between FDR control and minimax estimation, and the total error proportion is closely related to the asymptotic properties of FDR control in terms of the total error risk. The guide-gene driven selection is an approach to combining statistical significance and the existing biological knowledge of the study at hand. Error properties of these criteria are investigated theoretically and by simulation. The proposed methods are illustrated and compared using an example of genomic screening for novel Arf gene targets. Operating characteristics of q-value and the proposed significance threshold criteria are investigated and compared in a simulation study that employs a model mimicking a gene regulatory pathway. A guideline for using these criteria is provided. Splus/R code is available from the corresponding author upon request.

We demonstrate a concept and implementation of a compendium for the classification of high-dimensional data from microarray gene expression profiles. A compendium is an interactive document that bundles primary data, statistical processing methods, figures, and derived data together with the textual documentation and conclusions. Interactivity allows the reader to modify and extend these components. We address the following questions: how much does the discriminatory power of a classifier depend on the choice of the algorithm that was used to identify it; what alternative classifiers could be used just as well; how robust is the result. The answers to these questions are essential prerequisites for validation and biological interpretation of the classifiers. We show how to use this approach by looking at these questions for a specific breast cancer microarray data set that first has been studied by Huang et al. (2003).

Markov models provide a natural framework for modeling cellular and molecular level changes over time. Kalbfleisch and Lawless propose using a Chi-squared statistic for assessing the appropriateness of assuming a first-order, homogeneous Markov process. While this statistic provides a global test of the Markov assumption, it does not permit identification of individual departures. We consider two approaches for discovering specific departures from the Markov assumption. First, we propose a diagnostic that tests whether the number of observed transitions out of a given state at a given time point is different than expected. Second, we construct statistics based on the number of observations in each state at each time point. In both cases, we construct multiple correlated statistics and testing is achieved through simulations. These approaches are applied to HIV genetics sequences measured over time.

Accurate systems biology modeling requires a complete catalog of protein complexes and their constituent proteins. We discuss a graph-theoretic/statistical algorithm for local dynamic modeling of protein complexes using data from affinity purification-mass spectrometry experiments. The algorithm readily accommodates multicomplex membership by individual proteins and dynamic complex composition, two biological realities not accounted for in existing topological descriptions of the overall protein network. A likelihood-based objective function guides the protein complex modeling algorithm. With an accurate complex membership catalog in place, systems biology can proceed with greater precision.

This note is a response to the article "relating HIV-1 sequence variation to replication capacity via trees and forests" that appeared in Statistical Applications in Genetics and Molecular Biology, 3(1): http://www.bepress.com/sagmb/vol3/iss1/art2,2004

This note is the authors' response to the Reader's Reaction provided by Foulkes and De Gruttola as published in Statistical Applications in Genetics and Molecular Biology

Several important syndromes are caused by deleterious germline mutations of individual genes. In both clinical and research applications it is useful to evaluate the probability that an individual carries an inherited genetic variant of these genes, and to predict the risk of disease for that individual, using information on his/her family history. Mendelian risk prediction models accomplish these goals by integrating Mendelian principles and state-of-the-art statistical models to describe phenotype/genotype relationships. Here we introduce an R library called BayesMendel that allows implementation of Mendelian models in research and counseling settings. BayesMendel is implemented in an object-oriented structure in the language R and distributed freely as an open source library. In its first release, it includes two major cancer syndromes: the breast-ovarian cancer syndrome and the hereditary non-polyposis colorectal cancer syndrome, along with up-to-date estimates of penetrance and prevalence for the corresponding genes. Input genetic parameters can be easily modified by users. BayesMendel can also serve as a generic tool for genetic epidemiologists to flexibly implement their own Mendelian models for novel syndromes and local subpopulations, without reprogramming complex statistical analyses and prediction tools.

We discuss our concerns regarding the reliability of data generated by spotted cDNA microarrays. Two types of error we highlight are cross-hybridization artifact due to sequence homologies and sequence errors in the cDNA used for spotting on microarrays. We feel that statisticians who analyze microarray data should be aware of these sources of unreliability intrinsic to cDNA microarray design and use.