We investigate the computational power of particle methods, a well-established class of algorit hms with applications in scientific computing and computer simulation. The computational power of a compute model determines the class of problems it can solve. Automata theory allows describing the computational power of abstract machines (automata) and the problems they can solve. At the top of the Chomsky hierarchy of formal languages and grammars are Turing machines, which resemble the concept on which most modern computers are built. Although particle methods can be interpreted as automata based on their formal definition, their computational power has so far not been studied. We address this by analyzing Turing completeness of particle methods. In particular, we prove two sets of restrictions under which a particle method is still Turing powerful, and we show when it loses Turing powerfulness. This contributes to understanding the theoretical foundations of particle methods and provides insight into the powerfulness of computer simulations. 17 pages, 23 appendix pages

We present an efficient raycasting algorithm for rendering Volumetric Depth Images (VDIs), and we show how it can be used in a remote visualization setting with VDIs generated and streamed from a remote server. VDIs are compact view-dependent volume representations that enable interactive visualization of large volumes at high frame rates by decoupling viewpoint changes from expensive rendering calculations. However, current rendering approaches for VDIs struggle with achieving interactive frame rates at high image resolutions. Here, we exploit the properties of perspective projection to simplify intersections of rays with the view-dependent frustums in a VDI and leverage spatial smoothness in the volume data to minimize memory accesses. Benchmarks show that responsive frame rates can be achieved close to the viewpoint of generation for HD display resolutions, providing high-fidelity approximate renderings of Gigabyte-sized volumes. We also propose a method to subsample the VDI for preview rendering, maintaining high frame rates even for large viewpoint deviations. We provide our implementation as an extension of an established open-source visualization library.

Abstract We present a C++ library for transparent memory and compute abstraction across CPU and GPU architectures. Our library combines generic data structures like vectors, multi‐dimensional arrays, maps, graphs, and sparse grids with basic generic algorithms like arbitrary‐dimensional convolutions, copying, merging, sorting, prefix sum, reductions, neighbor search, and filtering. The memory layout of the data structures is adapted at compile time using C++ tuples with optional memory double‐mapping between host and device and the capability of using memory managed by external libraries with no data copying. We combine this transparent memory layout with generic thread‐parallel algorithms under two alternative common interfaces: a CUDA‐like kernel interface and a lambda‐function interface. We quantify the memory and compute performance and portability of our implementation using micro‐benchmarks, showing that the abstractions introduce negligible performance overhead, and we compare performance against the current state of the art in a real‐world scientific application from computational fluid mechanics.

MOTIVATION: Access to unprecedented amounts of quantitative biological data allows us to build and test biochemically accurate reaction-diffusion models of intracellular processes. However, any increase in model complexity increases the number of unknown parameters and, thus, the computational cost of model analysis. To efficiently characterize the behavior and robustness of models with many unknown parameters remains, therefore, a key challenge in systems biology. RESULTS: We propose a novel computational framework for efficient high-dimensional parameter space characterization of reaction-diffusion models in systems biology. The method leverages the Lp-Adaptation algorithm, an adaptive-proposal statistical method for approximate design centering and robustness estimation. Our approach is based on an oracle function, which predicts for any given point in parameter space whether the model fulfills given specifications. We propose specific oracles to efficiently predict four characteristics of Turing-type reaction-diffusion models: bistability, instability, capability of spontaneous pattern formation and capability of pattern maintenance. We benchmark the method and demonstrate that it enables global exploration of a model's ability to undergo pattern-forming instabilities and to quantify robustness for model selection in polynomial time with dimensionality. We present an application of the framework to pattern formation on the endosomal membrane by the small GTPase Rab5 and its effectors, and we propose molecular mechanisms underlying this system. AVAILABILITY AND IMPLEMENTATION: Our code is implemented in MATLAB and is available as open source under https://git.mpi-cbg.de/mosaic/software/black-box-optimization/rd-parameter-space-screening. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

We present data structures and algorithms for native implementations of discrete convolution operators over Adaptive Particle Representations (APR) of images on parallel computer architectures. The APR is a content-adaptive image representation that locally adapts the sampling resolution to the image signal. It has been developed as an alternative to pixel representations for large, sparse images as they typically occur in fluorescence microscopy. It has been shown to reduce the memory and runtime costs of storing, visualizing, and processing such images. This, however, requires that image processing natively operates on APRs, without intermediately reverting to pixels. Designing efficient and scalable APR-native image processing primitives, however, is complicated by the APR's irregular memory structure. Here, we provide the algorithmic building blocks required to efficiently and natively process APR images using a wide range of algorithms that can be formulated in terms of discrete convolutions. We show that APR convolution naturally leads to scale-adaptive algorithms that efficiently parallelize on multi-core CPU and GPU architectures. We quantify the speedups in comparison to pixel-based algorithms and convolutions on evenly sampled data. We achieve pixel-equivalent throughputs of up to 1TB/s on a single Nvidia GeForce RTX 2080 gaming GPU, requiring up to two orders of magnitude less memory than a pixel-based implementation.

Domain-specific languages (DSLs) are of increasing importance in scientific high-performance computing to reduce development costs, raise the level of abstraction, and, thus, ease scientific programming. However, designing DSLs is not easy, as it requires knowledge of the application domain and experience in language engineering and compilers. Consequently, many DSLs follow a weak approach using macros or text generators, which lack many of the features that make a DSL comfortable for programmers. Some of these features—e.g., syntax highlighting, type inference, error reporting—are easily provided by language workbenches, which combine language engineering techniques and tools in a common ecosystem. In this article, we present the Parallel Particle-Mesh Environment (PPME), a DSL and development environment for numerical simulations based on particle methods and hybrid particle-mesh methods. PPME uses the Meta Programming System, a projectional language workbench. PPME is the successor of the Parallel Particle-Mesh Language, a Fortran-based DSL that uses conventional implementation strategies. We analyze and compare both languages and demonstrate how the programmer’s experience is improved using static analyses and projectional editing, i.e., code-structure editing, constrained by syntax, as opposed to free-text editing. We present an explicit domain model for particle abstractions and the first formal type system for particle methods.

Embedded systems are often designed as complex architectures with numerous processing elements. Effectively programming such systems requires parallel programming models e.g. task-based or dataflow-based models. With these types of models, the mapping of the abstract application model to the existing hardware architecture plays a decisive role and is usually optimized to achieve an ideal resource footprint or a near-minimal execution time. However, when mapping several independent programs to the same platform, resource conflicts can arise. This can be circumvented by remapping some of the tasks of an application, which in turn affect its timing behavior, possibly leading to constraint violations. In this work we present a novel method to compute mappings that are robust against local task remapping. The underlying method is based on the bio-inspired design centering algorithm of Lp-Adaptation. We evaluate this with several benchmarks on different platforms and show that mappings obtained with our algorithm are indeed robust. In all experiments, our robust mappings tolerated significantly more run-time perturbations without violating constraints than mappings devised with optimization heuristics

Chemical reaction networks are ubiquitous in biology, and their dynamics is fundamentally stochastic. Here, we present the software library pSSAlib, which provides a complete and concise implementation of the most efficient partial-propensity methods for simulating exact stochastic chemical kinetics. pSSAlib can import models encoded in Systems Biology Markup Language, supports time delays in chemical reactions, and stochastic spatiotemporal reaction-diffusion systems. It also provides tools for statistical analysis of simulation results and supports multiple output formats. It has previously been used for studies of biochemical reaction pathways and to benchmark other stochastic simulation methods. Here, we describe pSSAlib in detail and apply it to a new model of the endocytic pathway in eukaryotic cells, leading to the discovery of a stochastic counterpart of the cut-out switch motif underlying early-to-late endosome conversion. pSSAlib is provided as a stand-alone command-line tool and as a developer API. We also provide a plug-in for the SBMLToolbox. The open-source code and pre-packaged installers are freely available from http://mosaic.mpi-cbg.de.

In image processing, the rapid approximate solution of variational problems involving generic data-fitting terms is often of practical relevance, for example in real-time applications. Variational solvers based on diffusion schemes or the Euler-Lagrange equations are too slow and restricted in the types of data-fitting terms. Here, we present a filter-based approach to reduce variational energies that contain generic data-fitting terms, but are restricted to specific regularizations. Our approach is based on reducing the regularization part of the variational energy, while guaranteeing non-increasing total energy. This is applicable to regularization-dominated models, where the data-fitting energy initially increases, while the regularization energy initially decreases. We present fast discrete filters for regularizers based on Gaussian curvature, mean curvature, and total variation. These pixel-local filters can be used to rapidly reduce the energy of the full model. We prove the convergence of the resulting iterative scheme in a greedy sense, and we show several experiments to demonstrate applications in image-processing problems involving regularization-dominated variational models.

Biological filaments, such as actin filaments, microtubules, and cilia, are often imaged using different light-microscopy techniques. Reconstructing the filament curve from the acquired images constitutes the filament segmentation problem. Since filaments have lower dimensionality than the image itself, there is an inherent trade-off between tracing the filament with sub-pixel accuracy and avoiding noise artifacts. Here, we present a globally optimal filament segmentation method based on B-spline vector level-sets and a generalized linear model for the pixel intensity statistics. We show that the resulting optimization problem is convex and can hence be solved with global optimality. We introduce a simple and efficient algorithm to compute such optimal filament segmentations, and provide an open-source implementation as an ImageJ/Fiji plugin. We further derive an information-theoretic lower bound on the filament segmentation error, quantifying how well an algorithm could possibly do given the information in the image. We show that our algorithm asymptotically reaches this bound in the spline coefficients. We validate our method in comprehensive benchmarks, compare with other methods, and show applications from fluorescence, phase-contrast, and dark-field microscopy.

The distributed resampling algorithm with non-proportional allocation (RNA) [1] is key to implementing particle filtering applications on parallel computer systems. We extend the original work by Bolić et al. by introducing an adaptive RNA (ARNA) algorithm, improving RNA by dynamically adjusting the particle-exchange ratio and randomizing the process ring topology. This improves the runtime performance of ARNA by about 9% over RNA with 10% particle exchange. ARNA also significantly improves the speed at which information is shared between processing elements, leading to about 20-fold faster convergence. The ARNA algorithm requires only a few modifications to the original RNA, and is hence easy to implement.

With the advent of multicore processors, numerical and mathematical software relies on parallelism in order to benefit from hardware performance increases. We present the design and use of a Fortran 2003 wrapper for POSIX threads, called forthreads. Forthreads is complete in the sense that is provides native Fortran 2003 interfaces to all pthreads routines where possible. We demonstrate the use and efficiency of forthreads for SIMD parallelism and task parallelism. We present forthreads/MPI implementations that enable hybrid shared-/distributed-memory parallelism in Fortran 2003. Our benchmarks show that forthreads offers performance comparable to that of OpenMP, but better thread control and more freedom. We demonstrate the latter by presenting a multithreaded Fortran 2003 library for POSIX Internet sockets, enabling interactive numerical simulations with runtime control.

We present a variational model with local weighted Gaussian curvature as regularizer. We show its convexity for an area-weight function and provide a closed-form solution for this case. The corresponding regularization coefficient has a theoretical bound. Moreover, we prove that the model is convex for a wide range of weight functions and show that it can be efficiently solved using splitting techniques. Finally, we demonstrate several applications of the model in image de-noising, smoothing, texture decomposition, image sharpening, and regularization-coefficient optimization.

BACKGROUND: Analyzing spatial distributions of objects in images is a fundamental task in many biological studies. The relative arrangement of a set of objects with respect to another set of objects contains information about potential interactions between the two sets of objects. If they do not "feel" each other's presence, their spatial distributions are expected to be independent of one another. Spatial correlations in their distributions are indicative of interactions and can be modeled by an effective interaction potential acting between the points of the two sets. This can be used to generalize co-localization analysis to spatial interaction analysis. However, no user-friendly software for this type of analysis was available so far. RESULTS: We present an ImageJ/Fiji plugin that implements the complete workflow of spatial pattern and interaction analysis for spot-like objects. The plugin detects objects in images, infers the interaction potential that is most likely to explain the observed pattern, and provides statistical tests for whether an inferred interaction is significant given the number of objects detected in the images and the size of the space within which they can distribute. We benchmark and demonstrate the present software using examples from confocal and PALM single-molecule microscopy. CONCLUSIONS: The present software greatly simplifies spatial interaction analysis for point patterns, and makes it available to the large user community of ImageJ and Fiji. The presented showcases illustrate the usage of the software.

Abstract We introduce a new class of data-fitting energies that couple image segmentation with image restoration. These functionals model the image intensity using the statistical framework of generalized linear models. By duality, we establish an information-theoretic interpretation using Bregman divergences. We demonstrate how this formulation couples in a principled way image restoration tasks such as denoising, deblurring (deconvolution), and inpainting with segmentation. We present an alternating minimization algorithm to solve the resulting composite photometric/geometric inverse problem. We use Fisher scoring to solve the photometric problem and to provide asymptotic uncertainty estimates. We derive the shape gradient of our data-fitting energy and investigate convex relaxation for the geometric problem. We introduce a new alternating split-Bregman strategy to solve the resulting convex problem and present experiments and comparisons on both synthetic and real-world images.

Incoming Simian Virus 40 particles bind to their cellular receptor, the glycolipid GM1, in the plasma membrane and thereby induce membrane deformation beneath the virion leading to endocytosis and infection. Efficient membrane deformation depends on receptor lipid structure and the organization of binding sites on the internalizing particle. To determine the role of receptor diffusion, concentration and the number of receptors required for stable binding in this interaction, we analyze the binding of SV40 to GM1 in supported membrane bilayers by computational modeling based on experimental data. We measure the diffusion rates of SV40 virions in solution by fluorescence correlation spectroscopy and of the receptor in bilayers by single molecule tracking. Quartz-crystal microbalance with dissipation (QCM-D) is used to measure binding of SV40 virus-like particles to bilayers containing the viral receptor GM1. We develop a phenomenological stochastic dynamics model calibrated against this data, and use it to investigate the early events of virus attachment to lipid membranes. Our results indicate that SV40 requires at least 4 attached receptors to achieve stable binding. We moreover find that receptor diffusion is essential for the establishment of stable binding over the physiological range of receptor concentrations and that receptor concentration controls the mode of viral motion on the target membrane. Our results provide quantitative insight into the initial events of virus-host interaction at the nanoscopic level.

We present the energy minimization of atomic clusters as a promising problem class for continuous black box optimization benchmarks. Finding the arrangement of atoms that minimizes a given potential energy is a specific instance of the more general class of geometry optimization or packing problems, which are generally NP-complete. Atomic clusters are a well-studied subject in physics and chemistry. From the large set of available cluster optimization problems, we propose two specific instances: Cohn-Kumar clusters and Lennard-Jones clusters. The potential energies of these clusters are governed by distance-dependent pairwise interaction potentials. The resulting collection of landscapes is composed of smooth and rugged single-funnel topologies, as well as tunable double-funnel topologies. In addition, all problems possess a feature that is not covered by the synthetic functions in current black box optimization test suites: isospectral symmetry. This property implies that any atomic arrangement is uniquely defined by the pairwise distance spectrum, rather than the absolute atomic positions. We hence suggest that the presented problem instances should be included in black box optimization benchmark suites.

We present a discrete, unsupervised multi-region competition algorithm for image segmentation over different energy functionals. The number of regions present in an image does not need to be known a priori, nor their photometric properties. The algorithm jointly estimates the number of regions, their photometries, and their contours. The required regularization is provided by defining a region as a connected set of pixels. The evolving contours in the image are represented by computational particles that move as driven by an energy-minimization algorithm. We present an efficient discrete algorithm that allows minimizing a range of well-known energy functionals under the topological constraint of regions being connected components. The presented framework and algorithms are implemented in the open-source Insight Toolkit (ITK) image-processing library.

We present a unifying framework for continuous optimization and sampling. This framework is based on Gaussian Adaptation (GaA), a search heuristic developed in the late 1960's. It is a maximum-entropy method that shares several features with the (1+1)-variant of the Covariance Matrix Adaptation Evolution Strategy (CMA-ES). The algorithm samples single candidate solutions from a multivariate normal distribution and continuously adapts the first and second moments. We present modifications that turn the algorithm into both a robust continuous black-box optimizer and, alternatively, an adaptive Random Walk Monte Carlo sampler. In black-box optimization, sample-point selection is controlled by a monotonically decreasing, fitness-dependent acceptance threshold. We provide general strategy parameter settings, stopping criteria, and restart mechanisms that render GaA quasi parameter free. We also introduce Metropolis GaA (M-GaA), where sample-point selection is based on the Metropolis acceptance criterion. This turns GaA into a Monte Carlo sampler that is conceptually similar to the seminal Adaptive Proposal (AP) algorithm. We evaluate the performance of Restart GaA on the CEC 2005 benchmark suite. Moreover, we compare the efficacy of M-GaA to that of the Metropolis-Hastings and AP algorithms on selected target distributions.

BACKGROUND: Sub-cellular structures interact in numerous direct and indirect ways in order to fulfill cellular functions. While direct molecular interactions crucially depend on spatial proximity, other interactions typically result in spatial correlations between the interacting structures. Such correlations are the target of microscopy-based co-localization analysis, which can provide hints of potential interactions. Two complementary approaches to co-localization analysis can be distinguished: intensity correlation methods capitalize on pattern discovery, whereas object-based methods emphasize detection power. RESULTS: We first reinvestigate the classical co-localization measure in the context of spatial point pattern analysis. This allows us to unravel the set of implicit assumptions inherent to this measure and to identify potential confounding factors commonly ignored. We generalize object-based co-localization analysis to a statistical framework involving spatial point processes. In this framework, interactions are understood as position co-dependencies in the observed localization patterns. The framework is based on a model of effective pairwise interaction potentials and the specification of a null hypothesis for the expected pattern in the absence of interaction. Inferred interaction potentials thus reflect all significant effects that are not explained by the null hypothesis. Our model enables the use of a wealth of well-known statistical methods for analyzing experimental data, as demonstrated on synthetic data and in a case study considering virus entry into live cells. We show that the classical co-localization measure typically under-exploits the information contained in our data. CONCLUSIONS: We establish a connection between co-localization and spatial interaction of sub-cellular structures by formulating the object-based interaction analysis problem in a spatial statistics framework based on nearest-neighbor distance distributions. We provide generic procedures for inferring interaction strengths and quantifying their relative statistical significance from sets of discrete objects as provided by image analysis methods. Within our framework, an interaction potential can either refer to a phenomenological or a mechanistic model of a physico-chemical interaction process. This increased flexibility in designing and testing different hypothetical interaction models can be used to quantify the parameters of a specific interaction model or may catalyze the discovery of functional relations.

The proteome of the radiation- and desiccation-resistant bacterium D. radiodurans features a group of proteins that contain significant intrinsically disordered regions that are not present in non-extremophile homologues. Interestingly, this group includes a number of housekeeping and repair proteins such as DNA polymerase III, nudix hydrolase and rotamase. Here, we focus on a member of the nudix hydrolase family from D. radiodurans possessing low-complexity N- and C-terminal tails, which exhibit sequence signatures of intrinsic disorder and have unknown function. The enzyme catalyzes the hydrolysis of oxidatively damaged and mutagenic nucleotides, and it is thought to play an important role in D. radiodurans during the recovery phase after exposure to ionizing radiation or desiccation. We use molecular dynamics simulations to study the dynamics of the protein, and study its hydration free energy using the GB/SA formalism. We show that the presence of disordered tails significantly decreases the hydration free energy of the whole protein. We hypothesize that the tails increase the chances of the protein to be located in the remaining water patches in the desiccated cell, where it is protected from the desiccation effects and can function normally. We extrapolate this to other intrinsically disordered regions in proteins, and propose a novel function for them: intrinsically disordered regions increase the "surface-properties" of the folded domains they are attached to, making them on the whole more hydrophilic and potentially influencing, in this way, their localization and cellular activity.

We extend the evolution strategy with covariance matrix adaptation (CMA-ES) by collaborative concepts from particle swarm optimization (PSO). The proposed particle swarm CMA-ES (PS-CMA-ES) algorithm is a hybrid real-parameter algorithm that combines the robust local search performance of CMA-ES with the global exploration power of PSO using multiple CMA-ES instances to explore different parts of the search space in parallel. Swarm intelligence is introduced by considering individual CMA-ES instances as lumped particles that communicate with each other. This includes non-local information in CMA-ES, which improves the search direction and the sampling distribution. We evaluate the performance of PS-CMA-ES on the IEEE CEC 2005 benchmark test suite. The new PS-CMA-ES algorithm shows superior performance on noisy problems and multi-funnel problems with non-convex underlying topology.

We present pCMALib, a parallel software library that implements the Evolution Strategy with Covariance Matrix Adaptation (CMA-ES). The library is written in Fortran 90/95 and uses the Message Passing Interface (MPI) for efficient parallelization on shared and distributed memory machines. It allows single CMA-ES optimization runs, embarrassingly parallel CMA-ES runs, and coupled parallel CMA-ES runs using a cooperative island model. As one instance of an island model CMA-ES, the recently presented Particle Swarm CMA-ES (PS-CMA-ES) is included using collaborative concepts from Swarm Intelligence for the migration model. Special attention has been given to an efficient design of the MPI communication protocol, a modular software architecture, and a user-friendly programming interface. The library includes a Matlab interface and is supplemented with an efficient Fortran implementation of the official CEC 2005 set of 25 real-valued benchmark functions. This is the first freely available Fortran implementation of this standard benchmark test suite. We present test runs and parallel scaling benchmarks on Linux clusters and multi-core desktop computers, showing good parallel efficiencies and superior computational performance compared to the reference implementation.

Single- and multicriteria evolution strategies are implemented to optimize micro-fluidic devices, namely the shape of a microchannel used for bioanalysls and the mixing rate in a micromixer used for medical applications. First, multimembered evolution strategies employing mutative step size adaptation are combined with the Strength Pareto approach. In order to support targetting, an extension of the Strength Pareto evolutionary algorithm is proposed. Applied on the optimization of the microchannel, these algorithms suggest a novel design with improved properties over traditional designs. A comparison with a gradient method is presented. Second, an evolution strategy with derandomized self-adaptation of the mutation distribution is used to optimize the micromixer. The results agree well with dynamical systems theory.