A Monte Carlo method to account for sampling error in multi-species indicators

The usefulness of biodiversity indicators strongly increases if accompanied by measures of uncertainty. In the case of indicators that combine population indices of species, however, the inclusion of the uncertainty of the species indices has shown to be hard to realize, usually due to imperfections in monitoring programmes. Missing values and time series of different lengths preclude the use of analytical approaches, whereas bootstrapping across sites requires the raw abundance data on the site level, which may not always be available. Sometimes bootstrapping across species rather than sites is opted for, but this approach ignores the uncertainty attached to species indices. We developed a method to account for sampling error of species indices in the calculation of multi-species indicators based on Monte Carlo simulation of annual species indices. The construction of confidence intervals enables various trend assessments, like testing for linear or smooth trends, testing for changes between two time points, testing the significance of a suspected change-point and testing for differences between two multi-species indicators. Here, we compare our method with conventional methods and illustrate the benefits of our approach using Dutch breeding bird indicators.     

Bayesian phylogenetic inference via Markov chain Monte Carlo methods

We derive a Markov chain to sample from the posterior distribution for a phylogenetic tree given sequence information from the corresponding set of organisms, a stochastic model for these data, and a prior distribution on the space of trees. A transformation of the tree into a canonical cophenetic matrix form suggests a simple and effective proposal distribution for selecting candidate trees close to the current tree in the chain. We illustrate the algorithm with restriction site data on 9 plant species, then extend to DNA sequences from 32 species of fish. The algorithm mixes well in both examples from random starting trees, generating reproducible estimates and credible sets for the path of evolution.     

Estimation of sampling variance of molecular marker data using the bootstrap procedure

Knowledge of genetic relationships among genotypes is useful in a plant breeding program because it permits the organization of germplasm and provides for more efficient sampling. The genetic distance (GD) among genotypes can be estimated using random restriction fragment length polymorphisms (RFLPs) as molecular markers. Knowledge of the sampling variance associated with RFLP markers is needed to determine how many markers are required for a given level of precision in the estimate of GD. The sampling variance for GD among all pairs of 37 maize (Z. mays L.) inbred lines was estimated from 1202 RFLPs. The 1202 polymorphisms were generated from 251 enzyme-probe combinations (EPC). The sampling variance was used to determine how large a sample of RFLPs was required to provide a given level of precision. The coefficient of variation (CV) associated with GD has a nearly linear relationship between its expected standard deviation and mean. The magnitude of the decrease in the mean CV for GD with increasing numbers of bands was dependent upon the sampling unit; e.g., individual polymorphic bands vs EPC, and the degree of relatedness among the inbreds compared. The rate of reduction in mean CV with increasing sample size was the same regardless of the restriction enzyme used, BamHI, EcoRI or HindIII, when the bootstrap sampling units were individual polymorphic bands. In constrast, although the rate of reduction (slopes) was the same, the intercepts of the mean CVs were different when EPCs were used as the bootstrap sampling unit. This difference was due to the higher number of bands per EPC in BamHI (4.94) compared with EcoRI (4.83) and HindIII (4.63).     

Pulse pressure variation tracking using sequential Monte Carlo methods

The pulse pressure variation (PPV) is a measure of the respiratory effect on the variation of systemic arterial blood pressure (ABP) in patients receiving full mechanical ventilation. It is a promising predictor of increases in cardiac output due to an infusion of fluid. We describe a novel automatic algorithm to estimate the PPV of ABP signals recorded under full respiratory support. The algorithm utilizes our recently developed sequential Monte Carlo method (SMCM), which is called a maximum a-posteriori adaptive marginalized particle filter (MAM-PF). MAM-PF estimates the state-space model parameters of the ABP signal continuously and its upper and lower envelopes are derived as a combination of those parameter estimates. Then, the continuous PPV values can be easily obtained based on those estimated envelopes. We report the assessment results of the proposed algorithm on real ABP signals.     

Seriation in paleontological data using markov chain Monte Carlo methods

Given a collection of fossil sites with data about the taxa that occur in each site, the task in biochronology is to find good estimates for the ages or ordering of sites. We describe a full probabilistic model for fossil data. The parameters of the model are natural: the ordering of the sites, the origination and extinction times for each taxon, and the probabilities of different types of errors. We show that the posterior distributions of these parameters can be estimated reliably by using Markov chain Monte Carlo techniques. The posterior distributions of the model parameters can be used to answer many different questions about the data, including seriation (finding the best ordering of the sites) and outlier detection. We demonstrate the usefulness of the model and estimation method on synthetic data and on real data on large late Cenozoic mammals. As an example, for the sites with large number of occurrences of common genera, our methods give orderings, whose correlation with geochronologic ages is 0.95.     

Lower confidence bounds for prediction accuracy in high dimensions via AROHIL Monte Carlo

MOTIVATION: Implementation and development of statistical methods for high-dimensional data often require high-dimensional Monte Carlo simulations. Simulations are used to assess performance, evaluate robustness, and in some cases for implementation of algorithms. But simulation in high dimensions is often very complex, cumbersome and slow. As a result, performance evaluations are often limited, robustness minimally investigated and dissemination impeded by implementation challenges. This article presents a method for converting complex, slow high-dimensional Monte Carlo simulations into simpler, faster lower dimensional simulations. RESULTS: We implement the method by converting a previous Monte Carlo algorithm into this novel Monte Carlo, which we call AROHIL Monte Carlo. AROHIL Monte Carlo is shown to exactly or closely match pure Monte Carlo results in a number of examples. It is shown that computing time can be reduced by several orders of magnitude. The confidence bound method implemented using AROHIL outperforms the pure Monte Carlo method. Finally, the utility of the method is shown by application to a number of real microarray datasets.     

Enhancing backbone sampling in Monte Carlo simulations using internal coordinates normal mode analysis

Normal mode methods are becoming a popular alternative to sample the conformational landscape of proteins. In this study, we describe the implementation of an internal coordinate normal mode analysis method and its application in exploring protein flexibility by using the Monte Carlo method PELE. This new method alternates two different stages, a perturbation of the backbone through the application of torsional normal modes, and a resampling of the side chains. We have evaluated the new approach using two test systems, ubiquitin and c-Src kinase, and the differences to the original ANM method are assessed by comparing both results to reference molecular dynamics simulations. The results suggest that the sampled phase space in the internal coordinate approach is closer to the molecular dynamics phase space than the one coming from a Cartesian coordinate anisotropic network model. In addition, the new method shows a great speedup (∼5–7×), making it a good candidate for future normal mode implementations in Monte Carlo methods.     

Inferring synaptic inputs given a noisy voltage trace via sequential Monte Carlo methods

We discuss methods for optimally inferring the synaptic inputs to an electrotonically compact neuron, given intracellular voltage-clamp or current-clamp recordings from the postsynaptic cell. These methods are based on sequential Monte Carlo techniques ("particle filtering"). We demonstrate, on model data, that these methods can recover the time course of excitatory and inhibitory synaptic inputs accurately on a single trial. Depending on the observation noise level, no averaging over multiple trials may be required. However, excitatory inputs are consistently inferred more accurately than inhibitory inputs at physiological resting potentials, due to the stronger driving force associated with excitatory conductances. Once these synaptic input time courses are recovered, it becomes possible to fit (via tractable convex optimization techniques) models describing the relationship between the sensory stimulus and the observed synaptic input. We develop both parametric and nonparametric expectation-maximization (EM) algorithms that consist of alternating iterations between these synaptic recovery and model estimation steps. We employ a fast, robust convex optimization-based method to effectively initialize the filter; these fast methods may be of independent interest. The proposed methods could be applied to better understand the balance between excitation and inhibition in sensory processing in vivo.     

Estimation,variance and optimal sampling of gene diversity

An extension of Nei's analysis of diversity in a subdivided population is proposed for a haploid locus. The differentiation G _STbecomes a natural extension of Wright's F _STand generalizes Weir and Cockerham's parameter of co-ancestry by relaxing the assumption of identical correlation for all the alleles. Inter- and intrapopulation variances of the estimated diversities and differentiation are derived. Finally, the optimal sampling strategy for measuring G _STwhen a fixed number of individuals can be analysed is considered. It is shown that, at a given locus, there is a unique sample size per population which yields the smallest variance of G _ST,regardless of the number of populations studied. These theoretical developments are illustrated with an analysis of chloroplast DNA diversity in a forest tree. The results emphasize the necessity of sampling many populations, rather than many individuals per population, for an accurate measurement of the subdivision of gene diversity at a single locus.     

Monte Carlo sampling of near-native structures of proteins with applications

Since a protein's dynamic fluctuation inside cells affects the protein's biological properties, we present a novel method to study the ensemble of near-native structures (NNS) of proteins, namely, the conformations that are very similar to the experimentally determined native structure. We show that this method enables us to (i) quantify the difficulty of predicting a protein's structure, (ii) choose appropriate simplified representations of protein structures, and (iii) assess the effectiveness of knowledge-based potential functions. We found that well-designed simple representations of protein structures are likely as accurate as those more complex ones for certain potential functions. We also found that the widely used contact potential functions stabilize NNS poorly, whereas potential functions incorporating local structure information significantly increase the stability of NNS.     

De novo prediction of polypeptide conformations using dihedral probability grid Monte Carlo methodology.

We tested the dihedral probability grid Monte Carlo (DPG-MC) methodology to determine optimal conformations of polypeptides by applying it to predict the low energy ensemble for two peptides whose solution NMR structures are known: integrin receptor peptide (YGRGDSP, Type II beta-turn) and S3 alpha-helical peptide (YMSEDEL KAAEAAFKRHGPT). DPG-MC involves importance sampling, local random stepping in the vicinity of a current local minima, and Metropolis sampling criteria for acceptance or rejection of new structures. Internal coordinate values are based on side-chain-specific dihedral angle probability distributions (from analysis of high-resolution protein crystal structures). Important features of DPG-MC are: (1) Each DPG-MC step selects the torsion angles (phi, psi, chi) from a discrete grid that are then applied directly to the structure. The torsion angle increments can be taken as S = 60, 30, 15, 10, or 5 degrees, depending on the application. (2) DPG-MC utilizes a temperature-dependent probability function (P) in conjunction with Metropolis sampling to accept or reject new structures. For each peptide, we found close agreement with the known structure for the low energy conformational ensemble located with DPG-MC. This suggests that DPG-MC will be useful for predicting conformations of other polypeptides.     

Estimation of admixture proportions: a likelihood-based approach using Markov chain Monte Carlo.

When populations are separated for long periods and then brought into contact for a brief episode in part of their range, this can result in genetic admixture. To analyze this type of event we considered a simple model under which two parental populations (P1 and P2) mix and create a hybrid population (H). After that event, the three populations evolve under pure drift without exchange during T generations. We developed a new method, which allows the simultaneous estimation of the time since the admixture event (scaled by the population size t(i) = T/N(i), where N(i) is the effective population size of population i) and the contribution of one of two parental populations (which we call p1). This method takes into account drift since the admixture event, variation caused by sampling, and uncertainty in the estimation of the ancestral allele frequencies. The method is tested on simulated data sets and then applied to a human data set. We find that (i) for single-locus data, point estimates are poor indicators of the real admixture proportions even when there are many alleles; (ii) biallelic loci provide little information about the admixture proportion and the time since admixture, even for very small amounts of drift, but can be powerful when many loci are used; (iii) the precision of the parameters' estimates increases with sample size n = 50 vs. n = 200 but this effect is larger for the t(i)'s than for p1; and (iv) the increase in precision provided by multiple loci is quite large, even when there is substantial drift (we found, for instance, that it is preferable to use five loci than one locus, even when drift is 100 times larger for the five loci). Our analysis of a previously studied human data set illustrates that the joint estimation of drift and p1 can provide additional insights into the data.     

Signal detection in averaged evoked potentials: Monte Carlo comparison of the sensitivity of different methods

Many clinical and research applications rely on detecting evoked potential (EP) signal or EP differences between conditions. Statistical methods for objective signal detection should be sensitive to the presence of signal, but must provide the user strict control on tolerated false alarm rate. The respective sensitivities of 6 signal detection methods were compared through several Monte Carlo simulations involving 2 autocorrelation structures, 5% and 1% significance levels, 8, 10 or 12 replications per study, and increasing signal to noise ratio. The signal detection methods compared were: (1) the Record Orthogonality Test by Permutations (ROT-p), a variant of the Residual Orthogonality Test (Achim et al., 1988), that provides an unbiased estimate of the energy of the signal present in the averaged data, (2) the Tsum² permutation test of Karniski et al. (1994), (3) a Principal Component Analysis method (PC1) consisting of a t test on the weights of the first principal component, (4) multiple t tests on amplitudes with empirical adjustment for global false alarm rate, and (5–6) the test of Guthrie and Buchwald (1991) on length of consecutive t tests significant at P < 0.05 or 0.01 per-test. The first 3 methods did not exceed their nominal false alarm rate and clearly outperformed the last 3, with the ROT-p method being significantly more sensitive than all others under almost all conditions.     

Sequence specific resonance assignment via Multicanonical Monte Carlo search using an ABACUS approach

ABACUS [Grishaev et al. (2005) Proteins 61:36-43] is a novel protocol for automated protein structure determination via NMR. ABACUS starts from molecular fragments defined by unassigned J-coupled spin-systems and involves a Monte Carlo stochastic search in assignment space, probabilistic sequence selection, and assembly of fragments into structures that are used to guide the stochastic search. Here, we report further development of the two main algorithms that increase the flexibility and robustness of the method. Performance of the BACUS [Grishaev and Llinás (2004) J Biomol NMR 28:1-101] algorithm was significantly improved through use of sequential connectivities available from through-bond correlated 3D-NMR experiments, and a new set of likelihood probabilities derived from a database of 56 ultra high resolution X-ray structures. A Multicanonical Monte Carlo procedure, Fragment Monte Carlo (FMC), was developed for sequence-specific assignment of spin-systems. It relies on an enhanced assignment sampling and provides the uncertainty of assignments in a quantitative manner. The efficiency of the protocol was validated on data from four proteins of between 68-116 residues, yielding 100% accuracy in sequence specific assignment of backbone and side chain resonances.     

Joint oligogenic segregation and linkage analysis using bayesian Markov chain Monte Carlo methods

One of the most challenging areas in human genetics is the dissection of quantitative traits. In this context, the efficient use of available data is important, including, when possible, use of large pedigrees and many markers for gene mapping. In addition, methods that jointly perform linkage analysis and estimation of the trait model are appealing because they combine the advantages of a model-based analysis with the advantages of methods that do not require prespecification of model parameters for linkage analysis. Here we review a Markov chain Monte Carlo approach for such joint linkage and segregation analysis, which allows analysis of oligogenic traits in the context of multipoint linkage analysis of large pedigrees. We provide an outline for practitioners of the salient features of the method, interpretation of the results, effect of violation of assumptions, and an example analysis of a two-locus trait to illustrate the method.