VoroCrust: Voronoi Meshing Without Clipping.

Polyhedral meshes are increasingly becoming an attractive option with particular advantages over traditional meshes for certain applications. What has been missing is a robust polyhedral meshing algorithm that can handle broad classes of domains exhibiting arbitrary curved boundaries and sharp features. In addition, the power of primal-dual mesh pairs, exemplified by Voronoi-Delaunay meshes, has been recognized as an important ingredient in numerous formulations. The VoroCrust algorithm is the first provably correct algorithm for conforming Voronoi meshing for non-convex and possibly non-manifold domains with guarantees on the quality of both surface and volume elements. A robust refinement process estimates a suitable sizing field that enables the careful placement of Voronoi seeds across the surface circumventing the need for clipping and avoiding its many drawbacks. The algorithm has the flexibility of filling the interior by either structured or random samples, while all sharp features are preserved in the output mesh. We demonstrate the capabilities of the algorithm on a variety of models and compare against state-of-the-art polyhedral meshing methods based on clipped Voronoi cells establishing the clear advantage of VoroCrust output.

Statistical Framework for Uncertainty Quantification in Computational Molecular Modeling.

As computational modeling, simulation, and predictions are becoming integral parts of biomedical pipelines, it behooves us to emphasize the reliability of the computational protocol. For any reported quantity of interest (QOI), one must also compute and report a measure of the uncertainty or error associated with the QOI. This is especially important in molecular modeling, since in most practical applications the inputs to the computational protocol are often noisy, incomplete, or low-resolution. Unfortunately, currently available modeling tools do not account for uncertainties and their effect on the final QOIs with sufficient rigor. We have developed a statistical framework that expresses the uncertainty of the QOI as the probability that the reported value deviates from the true value by more than some user-defined threshold. First, we provide a theoretical approach where this probability can be bounded using Azuma-Hoeffding like inequalities. Second, we approximate this probability empirically by sampling the space of uncertainties of the input and provide applications of our framework to bound uncertainties of several QOIs commonly used in molecular modeling. Finally, we also present several visualization techniques to effectively and quantitavely visualize the uncertainties: in the input, final QOIs, and also intermediate states.

Sampling Conditions for Conforming Voronoi Meshing by the VoroCrust Algorithm.

We study the problem of decomposing a volume with a smooth boundary into a collection of Voronoi cells. Unlike the dual problem of conforming Delaunay meshing, a principled solution to this problem for generic smooth surfaces remained elusive. VoroCrust leverages ideas from weighted α-shapes and the power crust algorithm to produce unweighted Voronoi cells conforming to the surface, yielding the first provably-correct algorithm for this problem. Given a κ-sparse ε-sample, we work with the balls of radius δ times the local feature size centered at each sample. The corners of the union of these balls on both sides of the surface are the Voronoi sites and the interface of their cells is a watertight surface reconstruction embedded in the dual shape of the union of balls. With the surface protected, the enclosed volume can be further decomposed by generating more sites inside it. Compared to clipping-based algorithms, VoroCrust cells are full Voronoi cells, with convexity and fatness guarantees. Compared to the power crust algorithm, VoroCrust cells are not filtered, are unweighted, and offer greater flexibility in meshing the enclosed volume by either structured or randomly genenerated samples.

Disk Density Tuning of a Maximal Random Packing.

We introduce an algorithmic framework for tuning the spatial density of disks in a maximal random packing, without changing the sizing function or radii of disks. Starting from any maximal random packing such as a Maximal Poisson-disk Sampling (MPS), we iteratively relocate, inject (add), or eject (remove) disks, using a set of three successively more-aggressive local operations. We may achieve a user-defined density, either more dense or more sparse, almost up to the theoretical structured limits. The tuned samples are conflict-free, retain coverage maximality, and, except in the extremes, retain the blue noise randomness properties of the input. We change the density of the packing one disk at a time, maintaining the minimum disk separation distance and the maximum domain coverage distance required of any maximal packing. These properties are local, and we can handle spatially-varying sizing functions. Using fewer points to satisfy a sizing function improves the efficiency of some applications. We apply the framework to improve the quality of meshes, removing non-obtuse angles; and to more accurately model fiber reinforced polymers for elastic and failure simulations.

All-Hex Meshing of Multiple-Region Domains without Cleanup.

In this paper, we present a new algorithm for all-hex meshing of domains with multiple regions without post-processing cleanup. Our method starts with a strongly balanced octree. In contrast to snapping the grid points onto the geometric boundaries, we move points a slight distance away from the common boundaries. Then we intersect the moved grid with the geometry. This allows us to avoid creating any flat angles, and we are able to handle two-sided regions and more complex topologies than prior methods. The algorithm is robust and cleanup-free; without the use of any pillowing, swapping, or smoothing. Thus, our simple algorithm is also more predictable than prior art.

Computational reconstitution of spine calcium transients from individual proteins.

We have built a stochastic model in the program MCell that simulates Ca(2+) transients in spines from the principal molecular components believed to control Ca(2+) entry and exit. Proteins, with their kinetic models, are located within two segments of dendrites containing 88 intact spines, centered in a fully reconstructed 6 × 6 × 5 µm(3) cube of hippocampal neuropil. Protein components include AMPA- and NMDA-type glutamate receptors, L- and R-type voltage-dependent Ca(2+) channels, Na(+)/Ca(2+) exchangers, plasma membrane Ca(2+) ATPases, smooth endoplasmic reticulum Ca(2+) ATPases, immobile Ca(2+) buffers, and calbindin. Kinetic models for each protein were taken from published studies of the isolated proteins in vitro. For simulation of electrical stimuli, the time course of voltage changes in the dendritic spine was generated with the desired stimulus in the program NEURON. Voltage-dependent parameters were then continuously re-adjusted during simulations in MCell to reproduce the effects of the stimulus. Nine parameters of the model were optimized within realistic experimental limits by a process that compared results of simulations to published data. We find that simulations in the optimized model reproduce the timing and amplitude of Ca(2+) transients measured experimentally in intact neurons. Thus, we demonstrate that the characteristics of individual isolated proteins determined in vitro can accurately reproduce the dynamics of experimentally measured Ca(2+) transients in spines. The model will provide a test bed for exploring the roles of additional proteins that regulate Ca(2+) influx into spines and for studying the behavior of protein targets in the spine that are regulated by Ca(2+) influx.

Robust All-quad Meshing of Domains with Connected Regions.

In this paper, we present a new algorithm for all-quad meshing of non-convex domains, with connected regions. Our method starts with a strongly balanced quadtree. In contrast to snapping the grid points onto the geometric boundaries, we move points a slight distance away from the common boundaries. Then we intersect the moved grid with the geometry. This allows us to avoid creating any flat quads, and we are able to handle two-sided regions and more complex topologies than prior methods. The algorithm is provably correct, robust, and cleanup-free; meaning we have angle and edge length bounds, without the use of any pillowing, swapping, or smoothing. Thus, our simple algorithm is also more predictable than prior art.

Extracellular sheets and tunnels modulate glutamate diffusion in hippocampal neuropil.

Although the extracellular space in the neuropil of the brain is an important channel for volume communication between cells and has other important functions, its morphology on the micron scale has not been analyzed quantitatively owing to experimental limitations. We used manual and computational techniques to reconstruct the 3D geometry of 180 µm(3) of rat CA1 hippocampal neuropil from serial electron microscopy and corrected for tissue shrinkage to reflect the in vivo state. The reconstruction revealed an interconnected network of 40-80 nm diameter tunnels, formed at the junction of three or more cellular processes, spanned by sheets between pairs of cell surfaces with 10-40 nm width. The tunnels tended to occur around synapses and axons, and the sheets were enriched around astrocytes. Monte Carlo simulations of diffusion within the reconstructed neuropil demonstrate that the rate of diffusion of neurotransmitter and other small molecules was slower in sheets than in tunnels. Thus, the non-uniformity found in the extracellular space may have specialized functions for signaling (sheets) and volume transmission (tunnels).

Single-particle reconstruction using L(2)-gradient flow.

In this paper, we present an iterative algorithm for reconstructing a three-dimensional density function from a set of two dimensional electron microscopy images. By minimizing an energy functional consisting of a fidelity term and a regularization term, an L(2)-gradient flow is derived. The flow is integrated by a finite element method in the spatial direction and an explicit Euler scheme in the temporal direction. Our method compares favorably with those of the weighted back projection, Fourier method, algebraic reconstruction technique and simultaneous iterative reconstruction technique.

An Automatic 3D Mesh Generation Method for Domains with Multiple Materials.

This paper describes an automatic and efficient approach to construct unstructured tetrahedral and hexahedral meshes for a composite domain made up of heterogeneous materials. The boundaries of these material regions form non-manifold surfaces. In earlier papers, we developed an octree-based isocontouring method to construct unstructured 3D meshes for a single-material (homogeneous) domain with manifold boundary. In this paper, we introduce the notion of a material change edge and use it to identify the interface between two or several different materials. A novel method to calculate the minimizer point for a cell shared by more than two materials is provided, which forms a non-manifold node on the boundary. We then mesh all the material regions simultaneously and automatically while conforming to their boundaries directly from volumetric data. Both material change edges and interior edges are analyzed to construct tetrahedral meshes, and interior grid points are analyzed for proper hexahedral mesh construction. Finally, edge-contraction and smoothing methods are used to improve the quality of tetrahedral meshes, and a combination of pillowing, geometric flow and optimization techniques is used for hexahedral mesh quality improvement. The shrink set of pillowing schemes is defined automatically as the boundary of each material region. Several application results of our multi-material mesh generation method are also provided.

CRYO-ELECTRON MICROSCOPY DATA DENOISING BASED ON THE GENERALIZED DIGITIZED TOTAL VARIATION METHOD.

The energy functional used in digitalized total variation method is expanded to a general form and a generalized digitized total variation (GDTV) denoising method is obtained. We further expand this method from 2-dimensional (2D) image to 3-dimensional (3D) image processing field. Cryo-electron microscopy (cryo EM) and single particle reconstruction are becoming part of standard collection of structural techniques used for studying macromolecular assemblies. Howerver, the 3D data obtained suffers greatly from noise and degradation for the low dose electron radiation. Thus, image enhancement and noise reduction are theoretically necessary to improve the data for the subsequent segmentation and/or structure skeletonization. Although there are several methods to tackle this problem, we are pleased to find that GDTV method is very efficient and can achieve better performance. Comparative experiments are carried out to verify the enhancement achieved by the GDTV method.

Quality Multi-domain Meshing for Volumetric Data.

Multi-domain meshing from volumetric data is of great importance in many fields like medicine, biology and geology. This paper proposes a new approach to produce a high quality mesh with separated multiple domains. A point cloud is generated from a preliminary mesh representing the boundary between different domains from the discrete volumetric representation used as input. A higher-order level-set method is employed to produce a quality sub-mesh from this point cloud and geometric flow is used as smoothing mechanism. A new approach to detect and curate intersections within an assembly of these 2-manifold sub-meshes by utilizing the intermediate volumetric representation is developed. The separation between sub-meshes can be controlled by the user using a gap threshold parameter. The resulting high quality multi-domain mesh is free from self- and inter-domain intersections and can be further utilized in finite element and boundary element computations. The proposed pipeline has been efficiently implemented and sample meshes have been provided for visualization.

A Fast Variational Method for the Construction of Resolution Adaptive C-Smooth Molecular Surfaces.

We present a variational approach to smooth molecular (proteins, nucleic acids) surface constructions, starting from atomic coordinates, as available from the protein and nucleic-acid data banks. Molecular dynamics (MD) simulations traditionally used in understanding protein and nucleic-acid folding processes, are based on molecular force fields, and require smooth models of these molecular surfaces. To accelerate MD simulations, a popular methodology is to employ coarse grained molecular models, which represent clusters of atoms with similar physical properties by psuedo- atoms, resulting in coarser resolution molecular surfaces. We consider generation of these mixed-resolution or adaptive molecular surfaces. Our approach starts from deriving a general form second order geometric partial differential equation in the level-set formulation, by minimizing a first order energy functional which additionally includes a regularization term to minimize the occurrence of chemically infeasible molecular surface pockets or tunnel-like artifacts. To achieve even higher computational efficiency, a fast cubic B-spline C(2) interpolation algorithm is also utilized. A narrow band, tri-cubic B-spline level-set method is then used to provide C(2) smooth and resolution adaptive molecular surfaces.

Three-dimensional geometric modeling of membrane-bound organelles in ventricular myocytes: bridging the gap between microscopic imaging and mathematical simulation.

A general framework of image-based geometric processing is presented to bridge the gap between three-dimensional (3D) imaging that provides structural details of a biological system and mathematical simulation where high-quality surface or volumetric meshes are required. A 3D density map is processed in the order of image pre-processing (contrast enhancement and anisotropic filtering), feature extraction (boundary segmentation and skeletonization), and high-quality and realistic surface (triangular) and volumetric (tetrahedral) mesh generation. While the tool-chain described is applicable to general types of 3D imaging data, the performance is demonstrated specifically on membrane-bound organelles in ventricular myocytes that are imaged and reconstructed with electron microscopic (EM) tomography and two-photon microscopy (T-PM). Of particular interest in this study are two types of membrane-bound Ca(2+)-handling organelles, namely, transverse tubules (T-tubules) and junctional sarcoplasmic reticulum (jSR), both of which play an important role in regulating the excitation-contraction (E-C) coupling through dynamic Ca(2+) mobilization in cardiomyocytes.

Bio-molecule Surfaces Construction via a Higher-Order Level-Set Method.

We present a general framework for a higher-order spline level-set (HLS) method and apply this to bio-molecule surfaces construction. Starting from a first order energy functional, we obtain a general level set formulation of geometric partial differential equation, and provide an efficient approach to solve this partial differential equation using a C(2) spline basis. We also present a fast cubic spline interpolation algorithm based on convolution and the Z-transform, which exploits the local relationship of interpolatory cubic spline coefficients with respect to given function data values. One example of our HLS method is demonstrated, which is the construction of bio-molecule surfaces (an implicit solvation interface) with their individual atomic coordinates and solvated radii as prerequisite.

Finite element analysis of the time-dependent Smoluchowski equation for acetylcholinesterase reaction rate calculations.

This article describes the numerical solution of the time-dependent Smoluchowski equation to study diffusion in biomolecular systems. Specifically, finite element methods have been developed to calculate ligand binding rate constants for large biomolecules. The resulting software has been validated and applied to the mouse acetylcholinesterase (mAChE) monomer and several tetramers. Rates for inhibitor binding to mAChE were calculated at various ionic strengths with several different time steps. Calculated rates show very good agreement with experimental and theoretical steady-state studies. Furthermore, these finite element methods require significantly fewer computational resources than existing particle-based Brownian dynamics methods and are robust for complicated geometries. The key finding of biological importance is that the rate accelerations of the monomeric and tetrameric mAChE that result from electrostatic steering are preserved under the non-steady-state conditions that are expected to occur in physiological circumstances.

Patient-Specific Vascular NURBS Modeling for Isogeometric Analysis of Blood Flow.

We describe an approach to construct hexahedral solid NURBS (Non-Uniform Rational B-Splines) meshes for patient-specific vascular geometric models from imaging data for use in isogeometric analysis. First, image processing techniques, such as contrast enhancement, filtering, classification, and segmentation, are used to improve the quality of the input imaging data. Then, lumenal surfaces are extracted by isocontouring the preprocessed data, followed by the extraction of vascular skeleton via Voronoi and Delaunay diagrams. Next, the skeleton-based sweeping method is used to construct hexahedral control meshes. Templates are designed for various branching configurations to decompose the geometry into mapped meshable patches. Each patch is then meshed using one-to-one sweeping techniques, and boundary vertices are projected to the lumenal surface. Finally, hexahedral solid NURBS are constructed and used in isogeometric analysis of blood flow. Piecewise linear hexahedral meshes can also be obtained using this approach. Examples of patient-specific arterial models are presented.

AUTOMATIC STRUCTURE INTERPRETATION OF SINGLE PARTICLE CRYO-ELECTRON MICROSCOPY: FROM IMAGES TO PSUEDO-ATOMIC MODELS.

Three dimensional Electron Microscopy (EM) and in particular single particle reconstruction using cryo-EM, has rapidly advanced over recent years, such that increasingly several macromolecular complexes can be resolved at subnanometer resolution (6-10 Å). This paper reviews some of the main volumetric image and geometric post-processing steps once a three dimensional EM map (henceforth a 3D map) has been reconstructed from single particle Cryo-EM, as essential steps in an enhanced and automated computational structure interpretation pipeline. In particular the paper addresses automated filtering, critical point calculations, symmetric and non-symmetric molecular domain segmentation, molecular surface selection, curation, and protein secondary structure (α- helices and ß-sheets) elucidation from 3D maps.

Application of new multi-resolution methods for the comparison of biomolecular electrostatic properties in the absence of global structural similarity.

In this paper we present a method for the multi-resolution comparison of biomolecular electrostatic potentials without the need for global structural alignment of the biomolecules. The underlying computational geometry algorithm uses multi-resolution attributed contour trees (MACTs) to compare the topological features of volumetric scalar fields. We apply the MACTs to compute electrostatic similarity metrics for a large set of protein chains with varying degrees of sequence, structure, and function similarity. For calibration, we also compute similarity metrics for these chains by a more traditional approach based upon 3D structural alignment and analysis of Carbo similarity indices. Moreover, because the MACT approach does not rely upon pairwise structural alignment, its accuracy and efficiency promises to perform well on future large-scale classification efforts across groups of structurally-diverse proteins. The MACT method discriminates between protein chains at a level comparable to the Carbo similarity index method; i.e., it is able to accurately cluster proteins into functionally-relevant groups which demonstrate strong dependence on ligand binding sites. The results of the analyses are available from the linked web databases http://ccvweb.cres.utexas.edu/MolSignature/ and http://agave.wustl.edu/similarity/. The MACT analysis tools are available as part of the public domain library of the Topological Analysis and Quantitative Tools (TAQT) from the Center of Computational Visualization, at the University of Texas at Austin (http://ccvweb.csres.utexas.edu/software). The Carbo software is available for download with the open-source APBS software package at http://apbs.sf.net/.

Discrete Surface Modelling Using Partial Differential Equations.

We use various nonlinear partial differential equations to efficiently solve several surface modelling problems, including surface blending, N-sided hole filling and free-form surface fitting. The nonlinear equations used include two second order flows, two fourth order flows and two sixth order flows. These nonlinear equations are discretized based on discrete differential geometry operators. The proposed approach is simple, efficient and gives very desirable results, for a range of surface models, possibly having sharp creases and corners.