Cactus Tools for Grid Applications

Cactus is an open source problem solving environment designed for scientists and engineers. Its modular structure facilitates parallel computation across different architectures and collaborative code development between different groups. The Cactus Code originated in the academic research community, where it has been developed and used over many years by a large international collaboration of physicists and computational scientists. We discuss here how the intensive computing requirements of physics applications now using the Cactus Code encourage the use of distributed and metacomputing, and detail how its design makes it an ideal application test-bed for Grid computing. We describe the development of tools, and the experiments which have already been performed in a Grid environment with Cactus, including distributed simulations, remote monitoring and steering, and data handling and visualization. Finally, we discuss how Grid portals, such as those already developed for Cactus, will open the door to global computing resources for scientific users.     

An application of fuzzy logic to build ecological sympatry networks

In recent years sympatry networks have been proposed as a mean to perform biogeographic analysis, but their computation posed practical difficulties that limited their use. We propose a novel approach, bringing closer the application of well-established network analysis tools to the study of sympatry patterns using both geographic and environmental data associated with the occurrence of species. Our proposed algorithm, SGraFuLo, combines the use of fuzzy logic and numerical methods to directly compute the network of interest from point locality records, without the need of specialized tools, such as geographic information systems, thereby simplifying the process for end users. By posing the problem in matrix terms, SGraFuLo is able to achieve remarkable efficiency even for large datasets, taking advantage of well established scientific computing algorithms. We present sympatry networks constructed using real-world data collected in Mexico and Central America and highlight the potential of our approach in the analysis of overlapping niches of species that could have important applications even in evolutionary studies. We also present details on the design and implementation of the algorithm, as well as experiments that show its efficiency. The source code is freely released and datasets are also available to support the reproducibility of our results.     

Synthesis of application-level utility functions for autonomic self-assessment

We present an approach to self-assessment for Autonomic Computing, based on the synthesis of utility functions, at the level of an autonomic application, or even a single task or feature performed by that application. This paper describes the fundamental steps of our approach: instrumentation of the application; collection of exhaustive samples of runtime data about relevant quality attributes of the application, as well as characteristics of its runtime environment; synthesis of a utility function through statistical correlation over the collected data points; and embedding of code corresponding to the equation of the synthesized utility function within the application, which enables the computation of utility values at run time. We describe a number of case studies, with their results and implications, to motivate and discuss the significance of application-level utility, illustrate our statistical synthesis method, and present our framework for instrumentation, monitoring, and utility function embedding/evaluation.     

Faster sequential genetic linkage computations.

Linkage analysis using maximum-likelihood estimation is a powerful tool for locating genes. As available data sets have grown, the computation required for analysis has grown exponentially and become a significant impediment. Others have previously shown that parallel computation is applicable to linkage analysis and can yield order-of-magnitude improvements in speed. In this paper, we demonstrate that algorithmic modifications can also yield order-of-magnitude improvements, and sometimes much more. Using the software package LINKAGE, we describe a variety of algorithmic improvements that we have implemented, demonstrating both how these techniques are applied and their power. Experiments show that these improvements speed up the programs by an order of magnitude, on problems of moderate and large size. All improvements were made only in the combinatorial part of the code, without restoring to parallel computers. These improvements synthesize biological principles with computer science techniques, to effectively restructure the time-consuming computations in genetic linkage analysis.     

遗传算法在转录因子结合位点识别中的应用

遗传算法是模拟生物进化过程的计算模型,是一种全局优化搜索算法。将遗传算法与转录因子结合位点识别问题相结合的新方法,以一致性序列模型作为保守motif的描述模型,通过对motif序列与待测序列的比对问题进行编码,将其转化成搜索空间中的优化问题,利用遗传算法来搜索最优解,预测转录因子的结合位点。实验结果表明,这种新的方法是有效的,它在占用少量内存的情况下能够准确地识别出待测转录因子结合位点。     

Proteins wriggle

We propose an algorithmic strategy for improving the efficiency of Monte Carlo searches for the low-energy states of proteins. Our strategy is motivated by a model of how proteins alter their shapes. In our model, when proteins fold under physiological conditions, their backbone dihedral angles change synchronously in groups of four or more to avoid steric clashes and respect the kinematic conservation laws. They wriggle; they do not thrash. We describe a simple algorithm that can be used to incorporate wriggling in Monte Carlo simulations of protein folding. We have tested this wriggling algorithm against a code in which the dihedral angles are varied independently (thrashing). Our standard of success is the average root-mean-square distance (rmsd) between the alpha-carbons of the folding protein and those of its native structure. After 100,000 Monte Carlo sweeps, the relative decrease in the mean rmsd, as one switches from thrashing to wriggling, rises from 11% for the protein 3LZM with 164 amino acids (aa) to 40% for the protein 1A1S with 313 aa and 47% for the protein 16PK with 415 aa. These results suggest that wriggling is useful and that its utility increases with the size of the protein. One may implement wriggling on a parallel computer or a computer farm.     

Fast Acceleration of 2D Wave Propagation Simulations Using Modern Computational Accelerators

Recent developments in modern computational accelerators like Graphics Processing Units (GPUs) and coprocessors provide great opportunities for making scientific applications run faster than ever before. However, efficient parallelization of scientific code using new programming tools like CUDA requires a high level of expertise that is not available to many scientists. This, plus the fact that parallelized code is usually not portable to different architectures, creates major challenges for exploiting the full capabilities of modern computational accelerators. In this work, we sought to overcome these challenges by studying how to achieve both automated parallelization using OpenACC and enhanced portability using OpenCL. We applied our parallelization schemes using GPUs as well as Intel Many Integrated Core (MIC) coprocessor to reduce the run time of wave propagation simulations. We used a well-established 2D cardiac action potential model as a specific case-study. To the best of our knowledge, we are the first to study auto-parallelization of 2D cardiac wave propagation simulations using OpenACC. Our results identify several approaches that provide substantial speedups. The OpenACC-generated GPU code achieved more than speedup above the sequential implementation and required the addition of only a few OpenACC pragmas to the code. An OpenCL implementation provided speedups on GPUs of at least faster than the sequential implementation and faster than a parallelized OpenMP implementation. An implementation of OpenMP on Intel MIC coprocessor provided speedups of with only a few code changes to the sequential implementation. We highlight that OpenACC provides an automatic, efficient, and portable approach to achieve parallelization of 2D cardiac wave simulations on GPUs. Our approach of using OpenACC, OpenCL, and OpenMP to parallelize this particular model on modern computational accelerators should be applicable to other computational models of wave propagation in multi-dimensional media.     

VDA, a method of choosing a better algorithm with fewer validations

The multitude of bioinformatics algorithms designed for performing a particular computational task presents end-users with the problem of selecting the most appropriate computational tool for analyzing their biological data. The choice of the best available method is often based on expensive experimental validation of the results. We propose an approach to design validation sets for method comparison and performance assessment that are effective in terms of cost and discrimination power.Validation Discriminant Analysis (VDA) is a method for designing a minimal validation dataset to allow reliable comparisons between the performances of different algorithms. Implementation of our VDA approach achieves this reduction by selecting predictions that maximize the minimum Hamming distance between algorithmic predictions in the validation set. We show that VDA can be used to correctly rank algorithms according to their performances. These results are further supported by simulations and by realistic algorithmic comparisons in silico.VDA is a novel, cost-efficient method for minimizing the number of validation experiments necessary for reliable performance estimation and fair comparison between algorithms.Our VDA software is available at http://sourceforge.net/projects/klugerlab/files/VDA/     

Hybrid framework for feasible modeling of an edge illumination X-ray phase-contrast imaging system at a human scale

Here, we present a hybrid approach for simulating an edge illumination X-ray phase-contrast imaging (EIXPCi) set-up using graphics processor units (GPU) with a high degree of accuracy. In this study, the applicability of pixel, mesh and non-uniform rational B-splines (NURBS) objects to carry out realistic maps of X-ray phase-contrast distribution at a human scale is accounted for by using numerical anthropomorphic phantoms and a very fast and robust simulation framework which integrates total interaction probabilities along selected X-ray paths. We exploit the mathematical and algorithmic properties of NURBS and describe how to represent human scale phantoms in an edge illumination X-ray phase-contrast model. The presented implementation allows the modeling of a variety of physical interactions of x-rays with different mathematically described objects and the recording of quantities, e.g. path integrals, interaction sites and deposited energies. Furthermore, our efficient, scalable and optimized hybrid Monte Carlo and ray-tracing projector can be used in iterative reconstruction algorithms on multi GPU heterogeneous systems. The preliminary results of our innovative approach show the fine performance of an edge illumination X-ray phase-contrast medical imaging system on various human-like soft tissues with noticeably reduced computation time. Our approach to the EIXPCi modeling confirms that building a true imaging system at a human scale should be possible and the simulations presented here aim at its future development.     

Secure outsourcing algorithms of modular exponentiations with optimal checkability based on a single untrusted cloud server

Modular exponentiation is an expensive discrete-logarithm operation, difficult for resource-constrained users to perform locally. Fortunately, thanks to burgeoning cloud computing, users are willing to securely outsourcing modular exponentiations to cloud servers to reduce computation overhead. In this paper, we contrive a fully verifiable secure outsourcing scheme for modular exponentiation with only a single server, named MExp. MExp not only prevents users’ private information leakage during outsourcing by our new logical division method, but also eliminates collusion attacks occurring in algorithms with two untrusted servers. Moreover, our MExp allows outsourcers to detect any misbehavior with a probability of 1, which shows significant improvement in checkability when compare to other single-server-based schemes. With a view to reducing computation overhead, MExp is extended to multiple modular exponentiations, named M2Exp. The algorithm significantly diminishes the local costs of multiple modular exponentiation calculations and the checkability is still 1. Compared with existing state-of-the-art schemes, MExp and M2Exp have outstanding performance in both efficiency and checkability. Finally, MExp and M2Exp are applied to Cramer–Shoup encryptions and Schnorr signatures.     

Bayesian analysis of an admixture model with mutations and arbitrarily linked markers

We introduce here a Bayesian analysis of a classical admixture model in which all parameters are simultaneously estimated. Our approach follows the approximate Bayesian computation (ABC) framework, relying on massive simulations and a rejection-regression algorithm. Although computationally intensive, this approach can easily deal with complex mutation models and partially linked loci, and it can be thoroughly validated without much additional computation cost. Compared to a recent maximum-likelihood (ML) method, the ABC approach leads to similarly accurate estimates of admixture proportions in the case of recent admixture events, but it is found superior when the admixture is more ancient. All other parameters of the admixture model such as the divergence time between parental populations, the admixture time, and the population sizes are also well estimated, unlike the ML method. The use of partially linked markers does not introduce any particular bias in the estimation of admixture, but ML confidence intervals are found too narrow if linkage is not specifically accounted for. The application of our method to an artificially admixed domestic bee population from northwest Italy suggests that the admixture occurred in the last 10-40 generations and that the parental Apis mellifera and A. ligustica populations were completely separated since the last glacial maximum.     

A communication-efficient linear system solver for large eddy simulation of jet engine noise

High-fidelity computational fluid dynamics (CFD) tools, such as the large eddy simulation technique, have become feasible in aiding the field of computational aeroacoustics (CAA) to compute noise on petascale computing platforms. CAA poses significant challenges for researchers because the computational schemes used in the CFD tools should have high accuracy, good spectral resolution, and low dispersion and diffusion errors. A high-order compact finite difference scheme, which is implicit in space, can be used for such simulations because it fulfills the requirements for CAA. Usually, this method is parallelized using a transposition scheme; however, that approach has a high communication overhead. In this paper, we discuss the use of a parallel tridiagonal linear system solver based on the truncated SPIKE algorithm for reducing the communication overhead in our large eddy simulations. We present theoretical performance analysis and report experimental results collected on two parallel computing platforms.     

Computational modeling of multicellular constructs with the material point method

Computational modeling of the mechanics of cells and multicellular constructs with standard numerical discretization techniques such as the finite element (FE) method is complicated by the complex geometry, material properties and boundary conditions that are associated with such systems. The objectives of this research were to apply the material point method (MPM), a meshless method, to the modeling of vascularized constructs by adapting the algorithm to accurately handle quasi-static, large deformation mechanics, and to apply the modified MPM algorithm to large-scale simulations using a discretization that was obtained directly from volumetric confocal image data. The standard implicit time integration algorithm for MPM was modified to allow the background computational grid to remain fixed with respect to the spatial distribution of material points during the analysis. This algorithm was used to simulate the 3D mechanics of a vascularized scaffold under tension, consisting of growing microvascular fragments embedded in a collagen gel, by discretizing the construct with over 13.6 million material points. Baseline 3D simulations demonstrated that the modified MPM algorithm was both more accurate and more robust than the standard MPM algorithm. Scaling studies demonstrated the ability of the parallel code to scale to 200 processors. Optimal discretization was established for the simulations of the mechanics of vascularized scaffolds by examining stress distributions and reaction forces. Sensitivity studies demonstrated that the reaction force during simulated extension was highly sensitive to the modulus of the microvessels, despite the fact that they comprised only 10.4% of the volume of the total sample. In contrast, the reaction force was relatively insensitive to the effective Poisson's ratio of the entire sample. These results suggest that the MPM simulations could form the basis for estimating the modulus of the embedded microvessels through a parameter estimation scheme. Because of the generality and robustness of the modified MPM algorithm, the relative ease of generating spatial discretizations from volumetric image data, and the ability of the parallel computational implementation to scale to large processor counts, it is anticipated that this modeling approach may be extended to many other applications, including the analysis of other multicellular constructs and investigations of cell mechanics.     

Super pairwise alignment (SPA): an efficient approach to global alignment for homologous sequences.

Sequence analysis is the basis of bioinformatics, while sequence alignment is a fundamental task for sequence analysis. The widely used alignment algorithm, Dynamic Programming, though generating optimal alignment, takes too much time due to its high computation complexity O(N(2)). In order to reduce computation complexity without sacrificing too much accuracy, we have developed a new approach to align two homologous sequences. The new approach presented here, adopting our novel algorithm which combines the methods of probabilistic and combinatorial analysis, reduces the computation complexity to as low as O(N). The computation speed by our program is at least 15 times faster than traditional pairwise alignment algorithms without a loss of much accuracy. We hence named the algorithm Super Pairwise Alignment (SPA). The pairwise alignment execution program based on SPA and the detailed results of the aligned sequences discussed in this article are available upon request.     

Automatic classification and robust identification of vestibulo-ocular reflex responses: from theory to practice

The Vestibulo-Ocular Reflex (VOR) stabilizes images of the world on our retinae when our head moves. Basic daily activities are thus impaired if this reflex malfunctions. During the past few decades, scientists have modeled and identified this system mathematically to diagnose and treat VOR deficits. However, traditional methods do not analyze VOR data comprehensively because they disregard the switching nature of nystagmus; this can bias estimates of VOR dynamics. Here we propose, for the first time, an automated tool to analyze entire VOR responses (slow and fast phases), without a priori classification of nystagmus segments. We have developed GNL-HybELS (Generalized NonLinear Hybrid Extended Least Squares), an algorithmic tool to simultaneously classify and identify the responses of a multi-mode nonlinear system with delay, such as the horizontal VOR and its alternating slow and fast phases. This algorithm combines the procedures of Generalized Principle Component Analysis (GPCA) for classification, and Hybrid Extended Least Squares (HybELS) for identification, by minimizing a cost function in an optimization framework. It is validated here on clean and noisy VOR simulations and then applied to clinical VOR tests on controls and patients. Prediction errors were less than 1 deg for simulations and ranged from .69 deg to 2.1 deg for the clinical data. Nonlinearities, asymmetries, and dynamic parameters were detected in normal and patient data, in both fast and slow phases of the response. This objective approach to VOR analysis now allows the design of more complex protocols for the testing of oculomotor and other hybrid systems.     

Bacterial community reconstruction using compressed sensing

Bacteria are the unseen majority on our planet, with millions of species and comprising most of the living protoplasm. We propose a novel approach for reconstruction of the composition of an unknown mixture of bacteria using a single Sanger-sequencing reaction of the mixture. Our method is based on compressive sensing theory, which deals with reconstruction of a sparse signal using a small number of measurements. Utilizing the fact that in many cases each bacterial community is comprised of a small subset of all known bacterial species, we show the feasibility of this approach for determining the composition of a bacterial mixture. Using simulations, we show that sequencing a few hundred base-pairs of the 16S rRNA gene sequence may provide enough information for reconstruction of mixtures containing tens of species, out of tens of thousands, even in the presence of realistic measurement noise. Finally, we show initial promising results when applying our method for the reconstruction of a toy experimental mixture with five species. Our approach may have a potential for a simple and efficient way for identifying bacterial species compositions in biological samples. All supplementary data and the MATLAB code are available at www.broadinstitute.org/?orzuk/publications/BCS/.     

Flux modules in metabolic networks

The huge number of elementary flux modes in genome-scale metabolic networks makes analysis based on elementary flux modes intrinsically difficult. However, it has been shown that the elementary flux modes with optimal yield often contain highly redundant information. The set of optimal-yield elementary flux modes can be compressed using modules. Up to now, this compression was only possible by first enumerating the whole set of all optimal-yield elementary flux modes. We present a direct method for computing modules of the thermodynamically constrained optimal flux space of a metabolic network. This method can be used to decompose the set of optimal-yield elementary flux modes in a modular way and to speed up their computation. In addition, it provides a new form of coupling information that is not obtained by classical flux coupling analysis. We illustrate our approach on a set of model organisms.     

Effective and efficient data sampling using bitmap indices

With growing computational capabilities of parallel machines, scientific simulations are being performed at finer spatial and temporal scales, leading to a data explosion. The growing sizes are making it extremely hard to store, manage, disseminate, analyze, and visualize these datasets, especially as neither the memory capacity of parallel machines, memory access speeds, nor disk bandwidths are increasing at the same rate as the computing power. Sampling can be an effective technique to address the above challenges, but it is extremely important to ensure that dataset characteristics are preserved, and the loss of accuracy is within acceptable levels. In this paper, we address the data explosion problems by developing a novel sampling approach, and implementing it in a flexible system that supports server-side sampling and data subsetting. We observe that to allow subsetting over scientific datasets, data repositories are likely to use an indexing technique. Among these techniques, we see that bitmap indexing can not only effectively support subsetting over scientific datasets, but can also help create samples that preserve both value and spatial distributions over scientific datasets. We have developed algorithms for using bitmap indices to sample datasets. We have also shown how only a small amount of additional metadata stored with bitvectors can help assess loss of accuracy with a particular subsampling level. Some of the other properties of this novel approach include: (1) sampling can be flexibly applied to a subset of the original dataset, which may be specified using a value-based and/or a dimension-based subsetting predicate, and (2) no data reorganization is needed, once bitmap indices have been generated. We have extensively evaluated our method with different types of datasets and applications, and demonstrated the effectiveness of our approach.