Ivo Sbalzarini Group

We present a consistent mesh-free numerical scheme for solving the incompressible Navier–Stokes equations. Our method is based on entropically damped artificial compressibility for imposing the incompressibility constraint explicitly, and the Discretization-Corrected Particle Strength Exchange (DC-PSE) method to consistently discretize the differential operators on mesh-free particles. We further couple our scheme with Brinkman penalization to solve the Navier–Stokes equations in complex geometries. The method is validated using the 3D Taylor–Green vortex flow and the lid-driven cavity flow problem in 2D and 3D, where we also compare our method with hr-SPH and report better accuracy for DC-PSE. In order to validate DC-PSE Brinkman penalization, we study flow past obstacles, such as a cylinder, and report excellent agreement with previous studies.

We present an artificial neural network architecture, termed STENCIL-NET, for equation-free forecasting of spatiotemporal dynamics from data. STENCIL-NET works by learning a discrete propagator that is able to reproduce the spatiotemporal dynamics of the training data. This data-driven propagator can then be used to forecast or extrapolate dynamics without needing to know a governing equation. STENCIL-NET does not learn a governing equation, nor an approximation to the data themselves. It instead learns a discrete propagator that reproduces the data. It therefore generalizes well to different dynamics and different grid resolutions. By analogy with classic numerical methods, we show that the discrete forecasting operators learned by STENCIL-NET are numerically stable and accurate for data represented on regular Cartesian grids. A once-trained STENCIL-NET model can be used for equation-free forecasting on larger spatial domains and for longer times than it was trained for, as an autonomous predictor of chaotic dynamics, as a coarse-graining method, and as a data-adaptive de-noising method, as we illustrate in numerical experiments. In all tests, STENCIL-NET generalizes better and is computationally more efficient, both in training and inference, than neural network architectures based on local (CNN) or global (FNO) nonlinear convolutions.

We present a versatile open-source pipeline for simulating inhomogeneous reaction–diffusion processes in highly resolved, image-based geometries of porous media with reactive boundaries. Resolving realistic pore-scale geometries in numerical models is challenging and computationally demanding, as the scale differences between the sizes of the interstitia and the whole system can lead to prohibitive memory requirements. The present pipeline combines a level-set method with geometry-adapted sparse block grids on GPUs to efficiently simulate reaction–diffusion processes in image-based geometries. We showcase the method by applying it to fertilizer diffusion in soil, heat transfer in porous ceramics, and determining effective diffusion coefficients and tortuosity. The present approach enables solving reaction–diffusion partial differential equations in real-world geometries applicable to porous media across fields such as engineering, environmental science, and biology.

We present a meshfree collocation scheme to discretize intrinsic surface differential operators over scalar fields on smooth curved surfaces with given normal vectors and a non-intersecting tubular neighborhood. The method is based on discretization-corrected particle strength exchange (DC-PSE), which generalizes finite difference methods to meshfree point clouds. The proposed Surface DC-PSE method is derived from an embedding theorem, but we analytically reduce the operator kernels along surface normals to obtain a purely intrinsic computational scheme over surface point clouds. We benchmark Surface DC-PSE by discretizing the Laplace–Beltrami operator on a circle and a sphere, and we present convergence results for both explicit and implicit solvers. We then showcase the algorithm on the problem of computing Gauss and mean curvature of an ellipsoid and of the Stanford Bunny by approximating the intrinsic divergence of the normal vector field. Finally, we compare Surface DC-PSE with surface finite elements (SFEM) and diffuse-interface finite elements (DI FEM) in a validation case.

In this paper, we present a novel methodology for automatic adaptive weighting of Bayesian Physics-Informed Neural Networks (BPINNs), and we demonstrate that this makes it possible to robustly address multi-objective and multiscale problems. BPINNs are a popular framework for data assimilation, combining the constraints of Uncertainty Quantification (UQ) and Partial Differential Equation (PDE). The relative weights of the BPINN target distribution terms are directly related to the inherent uncertainty in the respective learning tasks. Yet, they are usually manually set a-priori, that can lead to pathological behavior, stability concerns, and to conflicts between tasks which are obstacles that have deterred the use of BPINNs for inverse problems with multiscale dynamics. The present weighting strategy automatically tunes the weights by considering the multitask nature of target posterior distribution. We show that this remedies the failure modes of BPINNs and provides efficient exploration of the optimal Pareto front. This leads to better convergence and stability of BPINN training while reducing sampling bias. The determined weights moreover carry information about task uncertainties, reflecting noise levels in the data and adequacy of the PDE model. We demonstrate this in numerical experiments in Sobolev training, and compare them to analytically ε-optimal baseline, and in a multiscale Lotka-Volterra inverse problem. We eventually apply this framework to an inpainting task and an inverse problem, involving latent field recovery for incompressible flow in complex geometries.

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.

Active polar fluids exhibit spontaneous flow when sufficient active stress is generated by internal molecular mechanisms. This is also referred to as an active Fréedericksz transition. Experiments have revealed the existence of competing in-plane and out-of-plane instabilities in three-dimensional active matter. So far, however, a theoretical model reconciling all observations is missing. In particular, the role of boundary conditions in these instabilities still needs to be explained. Here, we characterize the spontaneous flow transition in a symmetry-preserving three-dimensional active Ericksen-Leslie model, showing that the boundary conditions select the emergent behavior. Using nonlinear numerical solutions and linear perturbation analysis, we explain the mechanism for both in-plane and out-of-plane instabilities under extensile active stress for perpendicular polarity anchoring at the boundary, whereas parallel anchoring only permits in-plane flows under contractile stress or out-of-plane wrinkling under extensile stress.

We investigate a class of DNA mixture deconvolution algorithms based on variational inference, and we show that this can significantly reduce computational runtimes with little or no effect on the accuracy and precision of the result. In particular, we consider Stein Variational Gradient Descent (SVGD) and Variational Inference (VI) with an evidence lower-bound objective. Both provide alternatives to the commonly used Markov-Chain Monte-Carlo methods for estimating the model posterior in Bayesian probabilistic genotyping. We demonstrate that both SVGD and VI significantly reduce computational costs over the current state of the art. Importantly, VI does so without sacrificing precision or accuracy, presenting an overall improvement over previously published methods.

Mathematical definitions provide a precise, unambiguous way to formulate concepts. They also provide a common language between disciplines. Thus, they are the basis for a well-founded scientific discussion. In addition, mathematical definitions allow for deeper insights into the defined subject based on mathematical theorems that are incontrovertible under the given definition. Besides their value in mathematics, mathematical definitions are indispensable in other sciences like physics, chemistry, and computer science. In computer science, they help to derive the expected behavior of a computer program and provide guidance for the design and testing of software. Therefore, mathematical definitions can be used to design and implement advanced algorithms. One class of widely used algorithms in computer science is the class of particle-based algorithms, also known as particle methods. Particle methods can solve complex problems in various fields, such as fluid dynamics, plasma physics, or granular flows, using diverse simulation methods, including Discrete Element Methods (DEM), Molecular Dynamics (MD), Reproducing Kernel Particle Methods (RKPM), Particle Strength Exchange (PSE), and Smoothed Particle Hydrodynamics (SPH). Despite the increasing use of particle methods driven by improved computing performance, the relation between these algorithms remains formally unclear. In particular, particle methods lack a unifying mathematical definition and precisely defined terminology. This prevents the determination of whether an algorithm belongs to the class and what distinguishes the class. Here we present a rigorous mathematical definition for determining particle methods and demonstrate its importance by applying it to several canonical algorithms and those not previously recognized as particle methods. Furthermore, we base proofs of theorems about parallelizability and computational power on it and use it to develop scientific computing software. Our definition unified, for the first time, the so far loosely connected notion of particle methods. Thus, it marks the necessary starting point for a broad range of joint formal investigations and applications across fields.

We present a C++14 library for performance portability of scientific computing codes across CPU and GPU architectures. Our library combines generic data structures like vectors, multi-dimensional arrays, maps, graphs, and sparse grids with basic, reusable algorithms like convolutions, sorting, prefix sum, reductions, and scan. The memory layout of the data structures is adapted at compile-time using tuples with optional memory mirroring between CPU and GPU. We combine this transparent memory mapping with generic algorithms under two alternative programming interfaces: a CUDA-like kernel interface for multi-core CPUs, Nvidia GPUs, and AMD GPUs, as well as a lambda interface. We validate and benchmark the presented library using micro-benchmarks, showing that the abstractions introduce negligible performance overhead, and we compare performance against the current state of the art.

We provide an internal validation study of a recently published precise DNA mixture algorithm based on Hamiltonian Monte Carlo sampling (Susik et al., 2022). We provide results for all 428 mixtures analysed by Riman et al. (2021) and compare the results with two state-of-the-art software products: STRmix™  v2.6 and Euroformix v3.4.0. The comparison shows that the Hamiltonian Monte Carlo method provides reliable values of likelihood ratios (LRs) close to the other methods. We further propose a novel large-scale precision benchmark and quantify the precision of the Hamiltonian Monte Carlo method, indicating its improvements over existing solutions. Finally, we analyse the influence of the factors discussed by Buckleton et al. (2022).

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.

Particle methods are a widely used class of algorithms for computer simulation of complex phenomena in various fields, such as fluid dynamics, plasma physics, molecular chemistry, and granular flows, using diverse simulation methods, including Smoothed Particle Hydrodynamics (SPH), Particle-in-Cell (PIC) methods, Molecular Dynamics (MD), and Discrete Element Methods (DEM). Despite the increasing use of particle methods driven by improved computing performance, the relation between these algorithms remains formally unclear, and a unifying formal definition of particle methods is lacking. Here, we present a rigorous mathematical definition of particle methods and demonstrate its importance by applying it to various canonical and non-canonical algorithms, using it to prove a theorem about multi-core parallelizability, and designing a principled scientific computing software based on it. We anticipate that our formal definition will facilitate the solution of complex computational problems and the implementation of understandable and maintainable software frameworks for computer simulation.

The control of cell shape during cytokinesis requires a precise regulation of mechanical properties of the cell cortex. Only few studies have addressed the mechanisms underlying the robust production of unequal-sized daughters during asymmetric cell division. Here we report that unequal daughter-cell sizes resulting from asymmetric sensory organ precursor divisions in Drosophila are controlled by the relative amount of cortical branched Actin between the two cell poles. We demonstrate this by mistargeting the machinery for branched Actin dynamics using nanobodies and optogenetics. We can thereby engineer the cell shape with temporal precision and thus the daughter-cell size at different stages of cytokinesis. Most strikingly, inverting cortical Actin asymmetry causes an inversion of daughter-cell sizes. Our findings uncover the physical mechanism by which the sensory organ precursor mother cell controls relative daughter-cell size: polarized cortical Actin modulates the cortical bending rigidity to set the cell surface curvature, stabilize the division and ultimately lead to unequal daughter-cell size.

We present a parallel compositing algorithm for Volumetric Depth Images (VDIs) of large three-dimensional volume data. Large distributed volume data are routinely produced in both numerical simulations and experiments, yet it remains challenging to visualize them at smooth, interactive frame rates. VDIs are view-dependent piecewise constant representations of volume data that offer a potential solution. They are more compact and less expensive to render than the original data. So far, however, there is no method for generating VDIs from distributed data. We propose an algorithm that enables this by sort-last parallel generation and compositing of VDIs with automatically chosen content-adaptive parameters. The resulting composited VDI can then be streamed for remote display, providing responsive visualization of large, distributed volume data.

We present a computational scheme that derives a global polynomial level set parameterization for smooth closed surfaces from a regular surface-point set and prove its uniqueness. This enables us to approximate a broad class of smooth surfaces by affine algebraic varieties. From such a global polynomial level set parameterization, differential-geometric quantities like mean and Gauss curvature can be efficiently and accurately computed. Even fourth-order terms such as the Laplacian of mean curvature are approximated with high precision. The accuracy performance results in a gain of computational efficiency, significantly reducing the number of surface points required compared to classic alternatives that rely on surface meshes or embedding grids. We mathematically derive and empirically demonstrate the strengths and the limitations of the present approach, suggesting it to be applicable to a large number of computational tasks in numerical differential geometry.

We present a higher-order convergent numerical solver for active polar hydrodynamics in three-dimensional domains of arbitrary shape, along with a scalable open-source software implementation for shared- and distributed-memory parallel computers. This enables the computational study of the nonlinear dynamics of out-of-equilibrium materials from first principles. We numerically solve the nonlinear active Ericksen–Leslie hydrodynamic equations of three-dimensional (3D) active nematics using both a meshfree and a hybrid particle-mesh method in either the Eulerian or Lagrangian frame of reference. The solver is validated against a newly derived analytical solution in 3D and implemented using the OpenFPM software library for scalable scientific computing. We then apply the presented method to studying the transition of 3D active polar fluids to spatiotemporal chaos, the emergence of coherent angular motion in a 3D annulus, and chiral vortices in symmetric and asymmetric 3D shapes resembling dividing cells. Overall, this provides a robust and efficient open-source simulation framework for 3D active matter with verified numerical convergence and scalability on parallel computers.

Topological defects are distinguished objects in vector fields that occur in a wide range of applications, ranging from material sciences to cosmology to bio-medical imaging and fingerprint recognition. This thesis considers topological point defects, also known as singular points, of two-dimensional vector fields. Besides Euclidean vectors as representation of modulus and direction, this also includes nematic vectors that equally have a modulus but direction is replaced with a head-to-tail symmetric orientation. In both case, a singular point or topological defect is an isolated discontinuity in an otherwise continuous vector field. It is characterized by its index or topological charge, which attains integer values for polar and half-integer values for nematic vector fields. There are different yet equivalent approaches to define the index. They either base on homology groups and the Brouwer degree, or on the first fundamental group and the mapping degree, or relatedly on lifting of a loop path enclosing the singular point. The definition by lift used here translates changes in the vector field along a path into a summed change in orientation angle. This translates to topological defects in discretized vector fields, where topological charge is calculated as sum of finite angle differences along a loop path between discretization points. On closer inspection, this calculation is an estimation, and is guaranteed to yield the correct estimate only with additional assumptions, for example when the underlying continuous-domain vector field is smooth and sampled at sufficiently high spatial resolution. Otherwise, arbitrary locations and charges of topological defects are possible, which yield exactly the same discretized vectors by the periodicity of representative orientation angles. Besides, the estimated topological charge depends discontinuously on each of the discrete input vectors and exhibits discrete jumps. As application data typically is subject to noise and uncertainty, this raises the question how reliable are topological defects identified in it. The present thesis quantifies, how large perturbations of a vector field are admissible without alteration of topological defects and charges. To that end, it introduces a robustness measure for each edge in a discretization grid that are combined along loop paths. Replacing critical edges of minimal robustness within a loop path by other path segments around a minimally larger area allows targeted increase of robustness. This data-dependent method called expansion over the critical edge is iterated until a user-set robustness is satisfied. The final areas of this algorithm are shown to have minimal size and therefore maximal spatial resolution, which also adapts to the local quality of data. The areas are also given as the faces in the graph of sufficiently robust edges after deleting all vertices of degree 1 (leaves) and all their connected edges. The minimal robust areas turn out to be nested by inclusion according to their robustness threshold. This allows to tradeoff detection robustness of topological charges versus their localization accuracy, both within a selection of pre-defined loop path shapes, and for free data-dependent expansion over the critical edge. Differently from defect identification by pattern matching, there is no restriction on the charge detectable. Besides, the robustness is shown to detect the size of unordered cores of defects. Robust defect areas indicate possible defect dynamics comprising motion, defect pair generation and annihilation already from single time point data. The robustness is also applicable to irregular discretization grids thanks to its graph theoretic characterization, and an extension to curved surfaces is foreseeable. The robust data-dependent defect identification is exemplified on microscopy images of the fruit fly Drosophila melanogaster. During Dorsal Closure, a developmental process, a cell sheet called amnioserosa contracts in highly regulated manner, whereby forces are actively generated and propagated along filamentous proteins like actin. Thereby, activity level and visco-elastic properties of the tissue are linked to the topological defects in the actin orientation field. Robust detection of these reveals that the sum over robust charges is clearly positive in the hundreds, whereas the overall sum of charge without robustness consideration fluctuates around zero. Numerous charges are observed, but $\pm 1/2$ dominate and confirm the amnioserosa as nematic material despite polar molecular constituents like actin. The sizes of robust defects span three orders of magnitude, and the largest defects follow the shapes of biological cells. The size distribution decays by a power law with the power for positive defects being more negative. Time courses show slightly higher speed of motion for +1/2 defects than for -1/2 defects, an order of magnitude above material flow velocity. Experiments with a genetic modification in the protein Crumbs had shown excess contraction of the amnioserosa cell layer during development. Comparing defect velocity of these embryos to wildtype suggests that viscosity and rotational viscosity increase stronger than activity level. This hypothesis remains to be tested in a combination of experiments and simulations, yet it would not have been generated in the first place without consideration of robust defects. More generally, the presented robustness measure and optimal data-dependent identification of topological defects could benefit the analysis of defects in discretized vector fields in a variety of disciplines. The optimal data-dependent identification allows for example to calculate error distributions for charge and localization of defects. The size, shape, and nested inclusion of robust defects constitute new observables, that generate numerous follow-up questions already for the fruit fly and enable novel analyses.

Cell migration is crucial for organismal development and shapes organisms in health and disease. Although a lot of research has revealed the role of intracellular components and extracellular signaling in driving single and collective cell migration, the influence of physical properties of the tissue and the environment on migration phenomena in vivo remains less explored. In particular, the role of the extracellular matrix (ECM), which many cells move upon, is currently unclear. To overcome this gap, we use zebrafish optic cup formation, and by combining novel transgenic lines and image analysis pipelines, we study how ECM properties influence cell migration in vivo. We show that collectively migrating rim cells actively move over an immobile extracellular matrix. These cell movements require cryptic lamellipodia that are extended in the direction of migration. Quantitative analysis of matrix properties revealed that the topology of the matrix changes along the migration path. These changes in matrix topologies are accompanied by changes in the dynamics of cell-matrix interactions. Experiments and theoretical modeling suggest that matrix porosity could be linked to efficient migration. Indeed, interfering with matrix topology by increasing its porosity results in a loss of cryptic lamellipodia, less-directed cell-matrix interactions, and overall inefficient migration. Thus, matrix topology is linked to the dynamics of cell-matrix interactions and the efficiency of directed collective rim cell migration during vertebrate optic cup morphogenesis.

We propose a Gaussian jump process model on a regular Cartesian lattice for the diffusion part of the Reaction-Diffusion Master Equation (RDME). We derive the resulting Gaussian RDME (GRDME) formulation from analogy with a kernel-based discretization scheme for continuous diffusion processes and quantify the limits of its validity relative to the classic RDME. We then present an exact stochastic simulation algorithm for the GRDME, showing that the accuracies of GRDME and RDME are comparable, but exact simulations of the GRDME require only a fraction of the computational cost of exact RDME simulations. We analyze the origin of this speedup and its scaling with problem dimension. The benchmarks suggest that the GRDME is a particularly beneficial model for diffusion-dominated systems in three dimensional spaces, often occurring in systems biology and cell biology.

Adversarial attacks rely on the instability phenomenon appearing in general for all inverse problems, e.g., image classification and reconstruction, independently of the computational scheme or method used to solve the problem. We mathematically prove and empirically show that machine learning denoisers (MLD) are not excluded. That is to prove the existence of adversarial attacks given by noise patterns making the MLD run into instability, i.e., the MLD increases the noise instead of decreasing it. We further demonstrate that adversarial retraining or classic filtering do not provide an exit strategy for this dilemma. Instead, we show that adversarial attacks can be inferred by polynomial regression. Removing the underlying inferred polynomial distribution from the total noise distribution delivers an efficient technique yielding robust MLDs that make consistent computer vision tasks such as image segmentation or classification more reliable.

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.

Analysing mixed DNA profiles is a common task in forensic genetics. Due to the complexity of the data, such analysis is often performed using Markov Chain Monte Carlo (MCMC)-based genotyping algorithms. These trade off precision against execution time. When default settings (including default chain lengths) are used, as large as a 10-fold changes in inferred log-likelihood ratios (LR) are observed when the software is run twice on the same case. So far, this uncertainty has been attributed to the stochasticity of MCMC algorithms. Since LRs translate directly to strength of the evidence in a criminal trial, forensic laboratories desire LR with small run-to-run variability.

We present a statistical learning framework for robust identification of differential equations from noisy spatio-temporal data. We address two issues that have so far limited the application of such methods, namely their robustness against noise and the need for manual parameter tuning, by proposing stability- based model selection to determine the level of regularization required for reproducible inference. This avoids manual parameter tuning and improves robustness against noise in the data. Our stability selection approach, termed PDE-STRIDE, can be combined with any sparsity-promoting regression method and provides an interpretable criterion for model component importance. We show that the particular combination of stability selection with the iterative hard-thresholding algorithm from compressed sensing provides a fast and robust framework for equation inference that outperforms previous approaches with respect to accuracy, amount of data required, and robustness. We illustrate the performance of PDE-STRIDE on a range of simulated benchmark problems, and we 2 demonstrate the applicability of PDE-STRIDE on real-world data by considering purely data-driven inference of the protein interaction network for embryonic polarization in Caenorhabditis elegans. Using fluorescence microscopy images of C. elegans zygotes as input data, PDE-STRIDE is able to learn the molecular interactions of the proteins.

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.

Quantifying molecular dynamics within the context of complex cellular morphologies is essential toward understanding the inner workings and function of cells. Fluorescence recovery after photobleaching (FRAP) is one of the most broadly applied techniques to measure the reaction diffusion dynamics of molecules in living cells. FRAP measurements typically restrict themselves to single-plane image acquisition within a subcellular-sized region of interest due to the limited temporal resolution and undesirable photobleaching induced by 3D fluorescence confocal or widefield microscopy. Here, an experimental and computational pipeline combining lattice light sheet microscopy, FRAP, and numerical simulations, offering rapid and minimally invasive quantification of molecular dynamics with respect to 3D cell morphology is presented. Having the opportunity to accurately measure and interpret the dynamics of molecules in 3D with respect to cell morphology has the potential to reveal unprecedented insights into the function of living cells.

We characterize and remedy a failure mode that may arise from multi-scale dynamics with scale imbalances during training of deep neural networks, such as physics informed neural networks (PINNs). PINNs are popular machine-learning templates that allow for seamless integration of physical equation models with data. Their training amounts to solving an optimization problem over a weighted sum of data-fidelity and equation-fidelity objectives. Conflicts between objectives can arise from scale imbalances, heteroscedasticity in the data, stiffness of the physical equation, or from catastrophic interference during sequential training. We explain the training pathology arising from this and propose a simple yet effective inverse Dirichlet weighting strategy to alleviate the issue. We compare with Sobolev training of neural networks, providing the baseline of analytically ε-optimal training. We demonstrate the effectiveness of inverse Dirichlet weighting in various applications, including a multi-scale model of active turbulence, where we show orders of magnitude improvement in accuracy and convergence over conventional PINN training. For inverse modeling using sequential training, we find that inverse Dirichlet weighting protects a PINN against catastrophic forgetting.

We present an autotuning approach for compile-time optimization of numerical discretization methods insimulations of partial differential equations. Our approach is based on data-driven regression of performancemodels for numerical methods. We use these models at compile time to automatically determine the parameters(e.g., resolution, time step size, etc.) of numerical simulations of continuum spatio-temporal models in order tooptimize the tradeoff between simulation accuracy and runtime. The resulting autotuner is developed for thecompiler of a Domain-Specific Language (DSL) for numerical simulations. The abstractions in the DSL enablethe compiler to automatically determine the performance models and know which discretization parametersto tune. We demonstrate that this high-level approach can explore a large space of possible simulations, withsimulation runtimes spanning multiple orders of magnitude. We evaluate our approach in two test cases:the linear diffusion equation and the nonlinear Gray-Scott reaction–diffusion equation. The results show thatour model-based autotuner consistently finds configurations that outperform those found by state-of-the-artgeneral-purpose autotuners. Specifically, our autotuner yields simulations that are on average 4.2x faster thanthose found by the best generic exploration algorithms, while using 16x less tuning time. Compared to manualtuning by a group of researchers with varying levels of expertise, the autotuner was slower than the best usersby not more than a factor of 2, whereas it was able to significantly outperform half of them.

We show how the finite sizes of unordered defect cores in discretized orientation and vector fields can reliably be estimated using a robustness measure for topological defects. Topological defects, or singular points, in vector and orientation fields are considered in applications from material science to life sciences to fingerprint recognition. Their identification from dis- cretized two-dimensional fields must deal with discontinuities, since the estimated topological charge jumps in (half-)integer steps upon orientation changes above a certain threshold. We use a recently proposed robustness measure [Hoffmann & Sbalzarini, Phys. Rev. E 103(1), 012602 (2021)] that exploits this effect to quantify the influence of noise in a vector field, and of the path chosen for defect estimation, on the detection reliability in two-dimensional discrete domains. Here, we show how this robustness measure can be used to quantify the sizes of unordered regions surrounding a defect, which are known as unordered cores. We suggest that the size of an unordered core can be identified as the smallest path radius of sufficient robustness. The resulting robust core-size estimation complements singular point and index estimation and may serve as uncertainty quantification of defect localization, or as an additional feature for defect characterization.

We present a user-friendly and intuitive C++ expression system to implement numerical simulations of continuum biological hydrodynamics. The expression system allows writing simulation programs in near-mathematical notation and makes codes more readable, more compact, and less error-prone. It also cleanly separates the implementation of the partial differential equation model from the implementation of the numerical methods used to discretize it. This allows changing either of them with minimal changes to the source code. The presented expression system is implemented in the high-performance computing platform OpenFPM, supporting simulations that transparently parallelize on multi-processor computer systems. We demonstrate that our expression system makes it easier to write scalable codes for simulating biological hydrodynamics in space and time. We showcase the present framework in numerical simulations of active polar fluids, as well as in classic simulations of fluid dynamics from the incompressible Navier-Stokes equations to Stokes flow in a ball. The presented expression system accelerates scalable simulations of spatio-temporal models that encode the physics and material properties of tissues in order to algorithmically study morphogenesis.

Mathematical modeling and simulation has emerged as a fundamental means to understand physical process around us with countless real-world applications in applied science and engineering problems. However, heavy reliance on first principles, symmetry relations, and conservation laws has limited its applicability to a few scientific domains and even few real-world scenarios. Especially in disciplines like biology the underlying living constituents exhibit a myriad of complexities like non-linearities, non-equilibrium physics, self-organization and plasticity that routinely escape mathematical treatment based on governing laws. Meanwhile, recent decades have witnessed rapid advancement in computing hardware, sensing technologies, and algorithmic innovations in machine learning. This progress has helped propel data-driven paradigms to achieve unprecedented practical success in the fields of image processing and computer vision, natural language processing, autonomous transport, and etc. In the current thesis, we explore, apply, and advance statistical and machine learning strategies that help bridge the gap between data and mathematical models, with a view toward modeling and simulation of spatiotemporal processes in biology. As first, we address the problem of learning interpretable mathematical models of biologial process from limited and noisy data. For this, we propose a statistical learning framework called PDE-STRIDE based on the theory of stability selection and ℓ0-based sparse regularization for parsimonious model selection. The PDE-STRIDE framework enables model learning with relaxed dependencies on tuning parameters, sample-size and noise-levels. We demonstrate the practical applicability of our method on real-world data by considering a purely data-driven re-evaluation of the advective triggering hypothesis explaining the embryonic patterning event in the C. elegans zygote. As a next natural step, we extend our PDE-STRIDE framework to leverage prior knowledge from physical principles to learn biologically plausible and physically consistent models rather than models that simply fit the data best. For this, we modify the PDE-STRIDE framework to handle structured sparsity constraints for grouping features which enables us to: 1) enforce conservation laws, 2) extract spatially varying non-observables, 3) encode symmetry relations associated with the underlying biological process. We show several applications from systems biology demonstrating the claim that enforcing priors dramatically enhances the robustness and consistency of the data-driven approaches. In the following part, we apply our statistical learning framework for learning mean-field deterministic equations of active matter systems directly from stochastic self-propelled active particle simulations. We investigate two examples of particle models which differs in the microscopic interaction rules being used. First, we consider a self-propelled particle model endowed with density-dependent motility character. For the chosen hydrodynamic variables, our data-driven framework learns continuum partial differential equations that are in excellent agreement with analytical derived coarse-grain equations from Boltzmann approach. In addition, our structured sparsity framework is able to decode the hidden dependency between particle speed and the local density intrinsic to the self-propelled particle model. As a second example, the learning framework is applied for coarse-graining a popular stochastic particle model employed for studying the collective cell motion in epithelial sheets. The PDE-STRIDE framework is able to infer novel PDE model that quantitatively captures the flow statistics of the particle model in the regime of low density fluctuations. Modern microscopy techniques produce GigaBytes (GB) and TeraBytes (TB) of data while imaging spatiotemporal developmental dynamics of living organisms. However, classical statistical learning based on penalized linear regression models struggle with issues like accurate computation of derivatives in the candidate library and problems with computational scalability for application to “big” and noisy data-sets. For this reason we exploit the rich parameterization of neural networks that can efficiently learn from large data-sets. Specifically, we explore the framework of Physics-Informed Neural Networks (PINN) that allow for seamless integration of physics priors with measurement data. We propose novel strategies for multi-objective optimization that allow for adapting PINN architecture to multi-scale modeling problems arising in biology. We showcase application examples for both forward and inverse modeling of mesoscale active turbulence phenomenon observed in dense bacterial suspensions. Employing our strategies, we demonstrate orders of magnitude gain in accuracy and convergence in comparison with conventional formulation for solving multi-objective optimization in PINNs. In the concluding chapter of the thesis, we skip model interpretability and focus on learning computable models directly from noisy data for the purpose of pure dynamics forecasting. We propose STENCIL-NET, an artificial neural network architecture that learns solution adaptive spatial discretization of an unknown PDE model that can be stably integrated in time with negligible loss in accuracy. To support this claim, we present numerical experiments on long-term forecasting of chaotic PDE solutions on coarse spatio-temporal grids, and also showcase de-noising application that help decompose spatiotemporal dynamics from the noise in an equation-free manner.

We present a design and implementation of distributed sparse block grids that transparently scale from a single CPU to multi-GPU clusters. We support dynamic sparse grids as, e.g., occur in computer graphics with complex deforming geometries and in multi-resolution numerical simulations. We present the data structures and algorithms of our approach, focusing on the optimizations required to render them computationally efficient on CPUs and GPUs alike. We provide a scalable implementation in the OpenFPM software library for HPC. We benchmark our implementation on up to 16 Nvidia GTX 1080 GPUs and up to 64 Nvidia A100 GPUs showing state-of-the-art scalability (68% to 96% parallel efficiency) on three benchmark problems. On a single GPU, our implementation is 14 to 140-fold faster than on a multi-core CPU.

We generalize Biggs Theorem to the case of directed cycles of multi-digraphs allowing to compute the dimension of the directed cycle space independently of the graph representation with linear runtime complexity. By considering two-dimensional CW complex of elementary cycles and deriving formulas for the Betti numbers of the associated cellular homology groups, we extend the list of representation independent topological inavariants measuring the graph structure. We prove the computation of the 2nd Betti number to be sharp #Phard in general and present specific representation invariant sub-fillings yielding efficiently computable homology groups. Finally, we suggest howto use the provided structural measures to shed new light on graph theoretical problems as graph embeddings, discrete Morse theoryand graph clustering.

We present a novel method for model-based segmentation of extended, blob-like objects on curved surfaces. Our method addresses several challenges arising when imaging curved biological membrane, such as out-of-membrane signal and geometry-induced background variations. We use a particle-based reconstruction of the membrane geometry, moment-conserving intensity interpolation from pixels to surface particles, and model-based in-surface segmentation. Our method denoises and deconvolves images, corrects for background variations, and quantifies the number, size, and intensity of segmented objects. We benchmark the accuracy of the method and present two applications to (1) neuroepithelial focal adhesion sites during optic cup morphogenesis in zebrafish and (2) reconstituted membrane domains bearing the small GTPase Rab5 on spherical beads.

We propose a statistical learning framework based on group-sparse regression that can be used to (i) enforce conservation laws, (ii) ensure model equivalence, and (iii) guarantee symmetries when learning or inferring differential-equation models from data. Directly learning interpretable mathematical models from data has emerged as a valuable modeling approach. However, in areas such as biology, high noise levels, sensor-induced correlations, and strong intersystem variability can render data-driven models nonsensical or physically inconsistent without additional constraints on the model structure. Hence, it is important to leverage prior knowledge from physical principles to learn biologically plausible and physically consistent models rather than models that simply fit the data best. We present the group iterative hard thresholding algorithm and use stability selection to infer physically consistent models with minimal parameter tuning. We show several applications from systems biology that demonstrate the benefits of enforcing priors in data-driven modeling.

SNT is an end-to-end framework for neuronal morphometry and whole-brain connectomics that supports tracing, proof-editing, visualization, quantification and modeling of neuroanatomy. With an open architecture, a large user base, community-based documentation, support for complex imagery and several model organisms, SNT is a flexible resource for the broad neuroscience community. SNT is both a desktop application and multi-language scripting library, and it is available through the Fiji distribution of ImageJ.

We develop a theoretical foundation for a time-series analysis method suitable for revealing the spectrum of diffusion coefficients in mixed Brownian systems, where no prior knowledge of particle distinction is required. This method is directly relevant for particle tracking in biological systems, where diffusion processes are often non-uniform. We transform Brownian data onto the logarithmic domain, where the coefficients for individual modes of diffusion appear as distinct spectral peaks in the probability density. We refer to the method as the logarithmic measure of diffusion, or simply as the logarithmic measure. We provide a general protocol for deriving analytical expressions for the probability densities on the logarithmic domain. The protocol is applicable for any number of spatial dimensions with any number of diffusive states. The analytical form can be fitted to data to reveal multiple diffusive modes. We validate the theoretical distributions and benchmark the accuracy and sensitivity of the method by extracting multi-modal diffusion coefficients from 2D Brownian simulations of poly-disperse filament bundles. Bundling the filaments allows us to control the system non-uniformity and hence quantify the sensitivity of the method. By exploiting the anisotropy of the simulated filaments, we generalize the logarithmic measure to rotational diffusion. By fitting the analytical forms to simulation data, we confirm the method's theoretical foundation. An error analysis in the single-mode regime shows that the proposed method is comparable in accuracy to the standard mean squared displacement approach for evaluating diffusion coefficients. For the case of multi-modal diffusion, we compare the logarithmic measure against other more sophisticated methods, showing that both model selectivity and extraction accuracy are comparable for small data sets. Therefore we suggest that the logarithmic measure, as a method for multi-modal diffusion coefficient extraction, is ideally suited for small data sets, a condition often confronted in the experimental context. Finally, we critically discuss the proposed benefits of the method and its information content.

The classical NP–hard feedback arc set problem (FASP) and feedback vertex set problem (FVSP) ask for a minimum set of arcs ε ⊆ E or vertices ν ⊆ V whose removal G ∖ ε, G ∖ ν makes a given multi–digraph G=(V, E) acyclic, respectively. Though both problems are known to be APX–hard, constant ratio approximations or proofs of inapproximability are unknown. We propose a new universal O(|V||E|4)–heuristic for the directed FASP. While a ratio of r ≈ 1.3606 is known to be a lower bound for the APX–hardness, at least by empirical validation we achieve an approximation of r ≤ 2. Most of the relevant applications, such as circuit testing, ask for solving the FASP on large sparse graphs, which can be done efficiently within tight error bounds with our approach.

Live in situ visualization of numerical simulations – interactive visualization while the simulation is running – can enable new modes of interaction, including computational steering. Designing easy-to-use distributed in situ architectures, with viewing latency low enough, and frame rate high enough, for interactive use, is challenging. Here, we propose a fully asynchronous, hybrid CPU–GPU in situ architecture that emphasizes interactivity. We also present a transparent implementation of this architecture embedded into the OpenFPM simulation framework. The benchmarks show that our architecture minimizes visual latencies, and achieves frame rates between 6 and 60 frames/second – depending on simulation data size and degree of parallelism – by changing only a few lines of an existing simulation code.

We present Bionic Tracking, a novel method for solving biological cell tracking problems with eye tracking in virtual reality using commodity hardware. Using gaze data, and especially smooth pursuit eye movements, we are able to track cells in time series of 3D volumetric datasets. The problem of tracking cells is ubiquitous in developmental biology, where large volumetric microscopy datasets are acquired on a daily basis, often comprising hundreds or thousands of time points that span hours or days. The image data, however, is only a means to an end, and scientists are often interested in the reconstruction of cell trajectories and cell lineage trees. Reliably tracking cells in crowded three-dimensional space over many time points remains an open problem, and many current approaches rely on tedious manual annotation or curation. In the Bionic Tracking approach, we substitute the usual 2D point-and-click interface for annotation or curation with eye tracking in a virtual reality headset, where users follow cells with their eyes in 3D space in order to track them. We detail the interaction design of our approach and explain the graph-based algorithm used to connect different time points, also taking occlusion and user distraction into account. We demonstrate Bionic Tracking using examples from two different biological datasets. Finally, we report on a user study with seven cell tracking experts, highlighting the benefits and limitations of Bionic Tracking compared to point-and-click interfaces.

The design principles dictating the spatio-temporal organisation of eukaryotic cells, and in particular the mechanisms controlling the self-organisation and dynamics of membrane-bound organelles such as the Golgi apparatus, remain elusive. Although this organelle was discovered 120 years ago, such basic questions as whether vesicular transport through the Golgi occurs in an anterograde (from entry to exit) or retrograde fashion are still strongly debated. Here, we address these issues by studying a quantitative model of organelle dynamics that includes: de-novo compartment generation, inter-compartment vesicular exchange, and biochemical conversion of membrane components. We show that anterograde or retrograde vesicular transports are asymptotic behaviors of a much richer dynamical system. Indeed, the structure and composition of cellular compartments and the directionality of vesicular exchange are intimately linked. They are emergent properties that can be tuned by varying the relative rates of vesicle budding, fusion and biochemical conversion.

Proteins can self-organize into spatial patterns via non-linear dynamic interactions on cellular membranes. Modelling and simulations have shown that small GTPases can generate patterns by coupling guanine nucleotide exchange factors (GEF) to effectors, generating a positive feedback of GTPase activation and membrane recruitment. Here, we reconstituted the patterning of the small GTPase Rab5 and its GEF/effector complex Rabex5/Rabaptin5 on supported lipid bilayers. We demonstrate a 'handover' of Rab5 from Rabex5 to Rabaptin5 upon nucleotide exchange. A minimal system consisting of Rab5, RabGDI and a complex of full length Rabex5/Rabaptin5 was necessary to pattern Rab5 into membrane domains. Rab5 patterning required a lipid membrane composition mimicking that of early endosomes, with PI(3)P enhancing membrane recruitment of Rab5 and acyl chain packing being critical for domain formation. The prevalence of GEF/effector coupling in nature suggests a possible universal system for small GTPase patterning involving both protein and lipid interactions.

The cell cortex, a thin film of active material assembled below the cell membrane, plays a key role in cellular symmetry-breaking processes such as cell polarity establishment and cell division. Here, we present a minimal model of the self-organization of the cell cortex that is based on a hydrodynamic theory of curved active surfaces. Active stresses on this surface are regulated by a diffusing molecular species. We show that coupling of the active surface to a passive bulk fluid enables spontaneous polarization and the formation of a contractile ring on the surface via mechanochemical instabilities. We discuss the role of external fields in guiding such pattern formation. Our work reveals that key features of cellular symmetry breaking and cell division can emerge in a minimal model via general dynamic instabilities.

We propose an approach to the numerical simulation of thin-film flows based on the lattice Boltzmann method. We outline the basic features of the method, show in which limits the expected thin-film equations are recovered, and perform validation tests. The numerical scheme is applied to the viscous Rayleigh-Taylor instability of a thin film and to the spreading of a sessile drop toward its equilibrium contact angle configuration. We show that the Cox-Voinov law is satisfied and that the effect of a tunable slip length on the substrate is correctly captured. We address, then, the problem of a droplet sliding on an inclined plane, finding that the Capillary number scales linearly with the Bond number, in agreement with experimental results. At last, we demonstrate the ability of the method to handle heterogenous and complex systems by showcasing the controlled dewetting of a thin film on a chemically structured substrate.

Excessive alcohol consumption is often linked to anxiety states and has a major relay center in the anterior part of bed nucleus of stria terminalis (BNST). We analyzed the impact of (i) genetic predisposition to high alcohol preference and consumption, and (ii) alcohol intake on anterior BNST, namely anterolateral (AL), anteromedial (AM), and anteroventral (lateral + medial subdivisions: AVl, AVm) subnuclei. We used two rat lines selectively bred for low- and high-alcohol preference and consumption, named Sardinian alcohol-non preferring (sNP) and -preferring (sP), respectively, the latter showing also inherent anxiety-related behaviors. We analyzed the modulation of calcitonin gene-related peptide (CGRP; exerting anxiogenic effects in BNST), neuropeptide Y (NPY; exerting mainly anxiolytic effects), and microglia activation (neuroinflammation marker, thought to increase anxiety). Calcitonin gene-related peptide-immunofluorescent fibers/terminals did not differ between alcohol-naive sP and sNP rats. Fiber/terminal NPY-immunofluorescent intensity was lower in BNST-AM and BNST-AVm of alcohol-naive sP rats. Activation of microglia (revealed by morphological analysis) was decreased in BNST-AM and increased in BNST-AVm of alcohol-naive sP rats. Prolonged (30 consecutive days), voluntary alcohol intake under the homecage 2-bottle "alcohol vs. water" regimen strongly increased CGRP intensity in BNST of sP rats in a subnucleus-specific manner: in BNST-AL, BNST-AVm, and BNST-AM. CGRP area sum, however, decreased in BNST-AM, without changes in other subnuclei. Alcohol consumption increased NPY expression, in a subnucleus-specific manner, in BNST-AL, BNST-AVl, and BNST-AVm. Alcohol consumption increased many size/shapes parameters in microglial cells, indicative of microglia de-activation. Finally, microglia density was increased in ventral anterior BNST (BNST-AVl, BNST-AVm) by alcohol consumption. In conclusion, genetic predisposition of sP rats to high alcohol intake could be in part mediated by anterior BNST subnuclei showing lower NPY expression and differential microglia activation. Alcohol intake in sP rats produced complex subnucleus-specific changes in BNST, affecting CGRP/NPY expression and microglia and leading to hypothesize that these changes might contribute to the anxiolytic effects of voluntarily consumed alcohol repeatedly observed in sP rats.

Mechanochemical processes in thin biological structures, such as the cellular cortex or epithelial sheets, play a key role during the morphogenesis of cells and tissues. In particular, they are responsible for the dynamical organization of active stresses that lead to flows and deformations of the material. Consequently, advective transport redistributes force-generating molecules and thereby contributes to a complex mechanochemical feedback loop. It has been shown in fixed geometries that this mechanism enables patterning, but the interplay of these processes with shape changes of the material remains to be explored. In this work, we study the fully self-organized shape dynamics using the theory of active fluids on deforming surfaces and develop a numerical approach to solve the corresponding force and torque balance equations. We describe the spontaneous generation of nontrivial surface shapes, shape oscillations, and directed surface flows that resemble peristaltic waves from self-organized, mechanochemical processes on the deforming surface. Our approach provides opportunities to explore the dynamics of self-organized active surfaces and can help to understand the role of shape as an integral element of the mechanochemical organization of morphogenetic processes.

Modern microscopes create a data deluge with gigabytes of data generated each second, and terabytes per day. Storing and processing this data is a severe bottleneck, not fully alleviated by data compression. We argue that this is because images are processed as grids of pixels. To address this, we propose a content-adaptive representation of fluorescence microscopy images, the Adaptive Particle Representation (APR). The APR replaces pixels with particles positioned according to image content. The APR overcomes storage bottlenecks, as data compression does, but additionally overcomes memory and processing bottlenecks. Using noisy 3D images, we show that the APR adaptively represents the content of an image while maintaining image quality and that it enables orders of magnitude benefits across a range of image processing tasks. The APR provides a simple and efficient content-aware representation of fluosrescence microscopy images.

We quantified cell population increase in the quail embryo enteric nervous system (ENS) from E2.5 (about 1500 cells) to E12 (about 8 million cells). We then probed ENS proliferative capacity by grafting to the chorio-allantoic membrane large (600 cells) and small (40 cells) populations of enteric neural crest (ENC) cells with aneural gut. This demonstrated that ENC cells show an extremely high capacity to regulate their proliferation while forming the ENS. Previous mathematical models and clonal label experiments revealed that a few dominant ENS "superstar" cell clones emerge but most clones are small. The model implied that "superstars" arise stochastically, but the same outcome could arise if "superstars" were pre-determined. We investigated these two modes mathematically and by grafting experiments with large and small numbers of ENCs, each including one EGFP-labelled ENC cell. The stochastic model predicts that the frequency of "superstar" detection increases as the ENC population decreases, the pre-determined model does not. Experimentally, as predicted by the stochastic model, the frequency of "superstar" detection increased with small ENC cell number. We conclude that ENS "superstar" clones achieve this status stochastically. Clonal dominance implies that clonal diversity is greatly reduced and in this case, somatic mutations may affect the phenotype. We suggest that somatic mutations coupled with loss of clonal diversity may contribute to variable penetrance and expressivity in individuals with genetically identical ENS pathologies.

An implicit aim in cellular infection biology is to understand the mechanisms how viruses, microbes, eukaryotic parasites, and fungi usurp the functions of host cells and cause disease. Mechanistic insight is a deep understanding of the biophysical and biochemical processes that give rise to an observable phenomenon. It is typically subject to falsification, that is, it is accessible to experimentation and empirical data acquisition. This is different from logic and mathematics, which are not empirical, but built on systems of inherently consistent axioms. Here, we argue that modeling and computer simulation, combined with mechanistic insights, yields unprecedented deep understanding of phenomena in biology and especially in virus infections by providing a way of showing sufficiency of a hypothetical mechanism. This ideally complements the necessity statements accessible to empirical falsification by additional positive evidence. We discuss how computational implementations of mathematical models can assist and enhance the quantitative measurements of infection dynamics of enveloped and non-enveloped viruses and thereby help generating causal insights into virus infection biology.

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.

Stochasticity in biological systems often referred to as gene expression noise is ubiquitous. The main sources of this noise come about in two ways. The implicit randomness of the biochemical reactions generates intrinsic noise inside the cell. Other cellular processes are themselves products that vary over time and from each cell to another, producing the so-called extrinsic noise. Controlling the mean expression level of a gene while reducing its noise is a challenge in many applications of Synthetic Biology. In previous works, we proposed a gene synthetic circuit to reduce gene expression noise while achieving a desired mean expression level. The circuit combines a negative feedback loop and a cell-to-cell communication mechanism based on quorum sensing. In this work, we use a multi-objective optimization design approach to find the best values for the tunable-in-the-lab parameters that, for a given desired mean expression value, achieve minimization of gene expression noise caused by intrinsic and extrinsic fluctuations. Our approach allows tuning the circuit parameters required to minimize noise effects, providing results which prove in accordance with genome-wide experimental data reported in the literature. Exploring different scenarios, either considering only intrinsic noise or considering both extrinsic and intrinsic ones, we find that the design strategies obtained for both cases are not transferable. Thus, designing the circuit parameters taking into account only intrinsic noise yields a sub-optimal design with decreased performance when evaluated in a scenario where both extrinsic and intrinsic noise are present.

Stochastic fluctuations in gene expression trigger both beneficial and harmful consequences for cell behavior. Therefore, achieving a desired mean protein expression level while minimizing noise is of interest in many applications, including robust protein production systems in industrial biotechnology. Here, we consider a synthetic gene circuit combining intracellular negative feedback and cell-to-cell communication based on quorum sensing. Accounting for both intrinsic and extrinsic noise, stochastic simulations allow us to analyze the capability of the circuit to reduce noise strength as a function of its parameters. We obtain mean expression levels and noise strengths for all species under different scenarios, showing good agreement with system-wide available experimental data of protein abundance and noise in Escherichia coli. Our in silico experiments, validated by preliminary in vivo results, reveal significant noise attenuation in gene expression through the interplay between quorum sensing and negative feedback and highlight the differential role that they play in regard to intrinsic and extrinsic noise.

Biologically functional DNA hairpins are found in archaea, prokaryotes and eukaryotes, playing essential roles in various DNA transactions. However, during DNA replication, hairpin formation can stall the polymerase and is therefore prevented by the single-stranded DNA binding protein (SSB). Here, we address the question how hairpins maintain their functional secondary structure despite SSB's presence. As a model hairpin, we used the recombinogenic form of the attC site, essential for capturing antibiotic-resistance genes in the integrons of bacteria. We found that attC hairpins have a conserved high GC-content near their apical loop that creates a dynamic equilibrium between attC fully opened by SSB and a partially structured attC-6-SSB complex. This complex is recognized by the recombinase IntI, which extrudes the hairpin upon binding while displacing SSB. We anticipate that this intriguing regulation mechanism using a base pair distribution to balance hairpin structure formation and genetic stability is key to the dissemination of antibiotic resistance genes among bacteria and might be conserved among other functional hairpins.

The polarization of cells and tissues is fundamental for tissue morphogenesis during biological development and regeneration. A deeper understanding of biological polarity pattern formation can be gained from the consideration of pattern reorganization in response to an opposing instructive cue, which we here consider using the example of experimentally inducible body axis inversions in planarian flatworms. We define a dynamically diluted alignment model linking three processes: entrainment of cell polarity by a global signal, local cell-cell coupling aligning polarity among neighbours, and cell turnover replacing polarized cells by initially unpolarized cells. We show that a persistent global orienting signal determines the final mean polarity orientation in this stochastic model. Combining numerical and analytical approaches, we find that neighbour coupling retards polarity pattern reorganization, whereas cell turnover accelerates it. We derive a formula for an effective neighbour coupling strength integrating both effects and find that the time of polarity reorganization depends linearly on this effective parameter and no abrupt transitions are observed. This allows us to determine neighbour coupling strengths from experimental observations. Our model is related to a dynamic 8-Potts model with annealed site-dilution and makes testable predictions regarding the polarization of dynamic systems, such as the planarian epithelium.

We present the Parallel Particle-Mesh Environment (PPME), a domainspecific language (DSL) and development environment for numerical simulations using particles and hybrid particle-mesh methods. PPME is the successor of the Parallel Particle-Mesh Language (PPML), a Fortran-based DSL that provides high-level abstractions for the development of distributed-memory particle-mesh simulations. On top of PPML, PPME provides a complete development environment for particle-based simulations usin state-of-the-art language engineering and compiler construction techniques. Relying on a novel domain metamodel and formal type system for particle methods, it enables advanced static code correctness checks at the level of particle abstractions, complementing the low-level analysis of the compiler. Furthermore, PPME adopts Herbie for improving the accuracy of floating-point expressions and supports a convenient high-level mathematical notation for equations and differential operators. For demonstration purposes, we discuss an example from Discrete Element Methods (DEM) using the classic Silbert model to simulate granular flows.

The design of systems or models that work robustly under uncertainty and environmental fluctuations is a key challenge in both engineering and science. This is formalized in the design-centering problem, which is defined as finding a design that fulfills given specifications and has a high probability of still doing so if the system parameters or the specifications fluctuate randomly. Design centering is often accompanied by the problem of quantifying the robustness of a system. Here we present a novel adaptive statistical method to simultaneously address both problems. Our method, L p -Adaptation, is inspired by the evolution of robustness in biological systems and by randomized schemes for convex volume computation. It is able to address both problems in the general, non-convex case and at low computational cost. We describe the concept and the algorithm, test it on known benchmarks, and demonstrate its real-world applicability in electronic and biological systems. In all cases, the present method outperforms the previous state of the art. This enables re-formulating optimization problems in engineering and biology as design centering problems, taking global system robustness into account.

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.

Bile, the central metabolic product of the liver, is transported by the bile canaliculi network. The impairment of bile flow in cholestatic liver diseases has urged a demand for insights into its regulation. Here, we developed a predictive 3D multi-scale model that simulates fluid dynamic properties successively from the subcellular to the tissue level. The model integrates the structure of the bile canalicular network in the mouse liver lobule, as determined by high-resolution confocal and serial block-face scanning electron microscopy, with measurements of bile transport by intravital microscopy. The combined experiment-theory approach revealed spatial heterogeneities of biliary geometry and hepatocyte transport activity. Based on this, our model predicts gradients of bile velocity and pressure in the liver lobule. Validation of the model predictions by pharmacological inhibition of Rho kinase demonstrated a requirement of canaliculi contractility for bile flow in vivo. Our model can be applied to functionally characterize liver diseases and quantitatively estimate biliary transport upon drug-induced liver injury.

In neoplastic cell growth, clones and subclones are variable both in size and mutational spectrum. The largest of these clones are believed to represent those cells with mutations that make them the most "fit," in a Darwinian sense, for expansion in their microenvironment. Thus, the degree of quantitative clonal expansion is regarded as being determined by innate qualitative differences between the cells that originate each clone. Here, using a combination of mathematical modelling and clonal labelling experiments applied to the developmental model system of the forming enteric nervous system, we describe how cells which are qualitatively identical may consistently produce clones of dramatically different sizes: most clones are very small while a few clones we term "superstars" contribute most of the cells to the final population. The basis of this is minor stochastic variations ("luck") in the timing and direction of movement and proliferation of individual cells, which builds a local advantage for daughter cells that is cumulative. This has potentially important consequences. In cancers, especially before strongly selective cytotoxic therapy, the assumption that the largest clones must be the cells with deterministic proliferative ability may not always hold true. In development, the gradual loss of clonal diversity as "superstars" take over the population may erode the resilience of the system to somatic mutations, which may have occurred early in clonal growth.

During metaphase, forces on kinetochores are exerted by k-fibres, bundles of microtubules that end at the kinetochore. Interestingly, non-kinetochore microtubules have been observed between sister kinetochores, but their function is unknown. Here we show by laser-cutting of a k-fibre in HeLa and PtK1 cells that a bundle of non-kinetochore microtubules, which we term 'bridging fibre', bridges sister k-fibres and balances the interkinetochore tension. We found PRC1 and EB3 in the bridging fibre, suggesting that it consists of antiparallel dynamic microtubules. By using a theoretical model that includes a bridging fibre, we show that the forces at the pole and at the kinetochore depend on the bridging fibre thickness. Moreover, our theory and experiments show larger relaxation of the interkinetochore distance for cuts closer to kinetochores. We conclude that the bridging fibre, by linking sister k-fibres, withstands the tension between sister kinetochores and enables the spindle to obtain a curved shape.

Viruses spread between cells, tissues, and organisms by cell-free and cell-cell mechanisms, depending on the cell type, the nature of the virus, or the phase of the infection cycle. The mode of viral transmission has a large impact on disease development, the outcome of antiviral therapies or the efficacy of gene therapy protocols. The transmission mode of viruses can be addressed in tissue culture systems using live-cell imaging. Yet even in relatively simple cell cultures, the mechanisms of viral transmission are difficult to distinguish. Here we present a cross-platform software framework called "Infectio," which is capable of simulating transmission phenotypes in tissue culture of virtually any virus. Infectio can estimate interdependent biological parameters, for example for vaccinia virus infection, and differentiate between cell-cell and cell-free virus spreading. Infectio assists in elucidating virus transmission mechanisms, a feature useful for designing strategies of perturbing or enhancing viral transmission. The complexity of the Infectio software is low compared to that of other software commonly used to quantitate features of cell biological images, which yields stable and relatively error-free output from Infectio. The software is open source (GPLv3 license), and operates on the major platforms (Windows, Mac, and Linux). The complete source code can be downloaded from http://infectio.github.io/index.html. IMPORTANCE Infectio presents a generalized platform to analyze virus infection spread between cells. It allows the simulation of plaque phenotypes from image-based assays. Viral plaques are the result of virus spreading from primary infected cells to neighboring cells. This is a complex process and involves neighborhood effects at cell-cell contact sites or fluid dynamics in the extracellular medium. Infectio differentiates between two major modes of virus transmission between cells, allowing in silico testing of hypotheses about spreading mechanisms of any virus which can be grown in cell cultures, based on experimentally measured parameters, such as infection intensity or cell killing. The results of these tests can be compared with experimental data and allow interpretations with regard to biophysical mechanisms. Infectio also facilitates characterizations of the mode of action of therapeutic agents, such as oncolytic viruses or other infectious or cytotoxic agents.

Detection and quantification of fluorescently labeled molecules in subcellular compartments is a key step in the analysis of many cell biological processes. Pixel-wise colocalization analyses, however, are not always suitable, because they do not provide object-specific information, and they are vulnerable to noise and background fluorescence. Here we present a versatile protocol for a method named 'Squassh' (segmentation and quantification of subcellular shapes), which is used for detecting, delineating and quantifying subcellular structures in fluorescence microscopy images. The workflow is implemented in freely available, user-friendly software. It works on both 2D and 3D images, accounts for the microscope optics and for uneven image background, computes cell masks and provides subpixel accuracy. The Squassh software enables both colocalization and shape analyses. The protocol can be applied in batch, on desktop computers or computer clusters, and it usually requires <1 min and <5 min for 2D and 3D images, respectively. Basic computer-user skills and some experience with fluorescence microscopy are recommended to successfully use the protocol.

Particle tracking is of key importance for quantitative analysis of intracellular dynamic processes from time-lapse microscopy image data. Because manually detecting and following large numbers of individual particles is not feasible, automated computational methods have been developed for these tasks by many groups. Aiming to perform an objective comparison of methods, we gathered the community and organized an open competition in which participating teams applied their own methods independently to a commonly defined data set including diverse scenarios. Performance was assessed using commonly defined measures. Although no single method performed best across all scenarios, the results revealed clear differences between the various approaches, leading to notable practical conclusions for users and developers.

Understanding biological phenomena requires a systemic approach that incorporates different mechanisms acting on different spatial and temporal scales, since in organisms the workings of all components, such as organelles, cells, and organs interrelate. This inherent interdependency between diverse biological mechanisms, both on the same and on different scales, provides the functioning of an organism capable of maintaining homeostasis and physiological stability through numerous feedback loops. Thus, developing models of organisms and their constituents should be done within the overall systemic context of the studied phenomena. We introduce such a method for modeling growth and regeneration of livers at the organ scale, considering it a part of the overall multi-scale biochemical and biophysical processes of an organism. Our method is based on the earlier discovered general growth law, postulating that any biological growth process comprises a uniquely defined distribution of nutritional resources between maintenance needs and biomass production. Based on this law, we introduce a liver growth model that allows to accurately predicting the growth of liver transplants in dogs and liver grafts in humans. Using this model, we find quantitative growth characteristics, such as the time point when the transition period after surgery is over and the liver resumes normal growth, rates at which hepatocytes are involved in proliferation, etc. We then use the model to determine and quantify otherwise unobservable metabolic properties of livers.

This essay provides an introduction to the terminology, concepts, methods, and challenges of image-based modeling in biology. Image-based modeling and simulation aims at using systematic, quantitative image data to build predictive models of biological systems that can be simulated with a computer. This allows one to disentangle molecular mechanisms from effects of shape and geometry. Questions like "what is the functional role of shape" or "how are biological shapes generated and regulated" can be addressed in the framework of image-based systems biology. The combination of image quantification, model building, and computer simulation is illustrated here using the example of diffusion in the endoplasmic reticulum.

Viruses spread between cells, tissues, and organisms by cell-free and cell-cell transmissions. Both mechanisms enhance disease development, but it is difficult to distinguish between them. Here, we analyzed the transmission mode of human adenovirus (HAdV) in monolayers of epithelial cells by wet laboratory experimentation and a computer simulation. Using live-cell fluorescence microscopy and replication-competent HAdV2 expressing green fluorescent protein, we found that the spread of infection invariably occurred after cell lysis. It was affected by convection and blocked by neutralizing antibodies but was independent of second-round infections. If cells were overlaid with agarose, convection was blocked and round plaques developed around lytic infected cells. Infected cells that did not lyse did not give rise to plaques, highlighting the importance of cell-free transmission. Key parameters for cell-free virus transmission were the time from infection to lysis, the dose of free viruses determining infection probability, and the diffusion of single HAdV particles in aqueous medium. With these parameters, we developed an in silico model using multiscale hybrid dynamics, cellular automata, and particle strength exchange. This so-called white box model is based on experimentally determined parameters and reproduces viral infection spreading as a function of the local concentration of free viruses. These analyses imply that the extent of lytic infections can be determined by either direct plaque assays or can be predicted by calculations of virus diffusion constants and modeling.

We present a discrete, unsupervised multi-regioncompetition 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.

The plasma membrane of cells is a dynamic mixture of different lipids, proteins, and sugars. In a recent issue of Cell, Gowrishankar et al. (2012) propose a model for how the actin cortex may generate and regulate lateral heterogeneity in the plasma membrane by actively clustering cell surface molecules.

Accurate extraction of cell outlines from microscopy images is essential for analysing the dynamics of migrating cells. Phase-contrast microscopy is one of the most common and convenient imaging modalities for observing cell motility because it does not require exogenous labelling and uses only moderate light levels with generally negligible phototoxicity effects. Automatic extraction and tracking of high-resolution cell outlines from phase-contrast images, however, is difficult due to complex and non-uniform edge intensity. We present a novel image-processing method based on refined level-set segmentation for accurate extraction of cell outlines from high-resolution phase-contrast images. The algorithm is validated on synthetic images of defined noise levels and applied to real image sequences of polarizing and persistently migrating keratocyte cells. We demonstrate that the algorithm is able to reliably reveal fine features in the cell edge dynamics.

We consider the problem of inferring the unknown parameters of a stochastic biochemical network model from a single measured time-course of the concentration of some of the involved species. Such measurements are available, e.g., from live-cell fluorescence microscopy in image-based systems biology. In addition, fluctuation time-courses from, e.g., fluorescence correlation spectroscopy (FCS) provide additional information about the system dynamics that can be used to more robustly infer parameters than when considering only mean concentrations. Estimating model parameters from a single experimental trajectory enables single-cell measurements and quantification of cell-cell variability. We propose a novel combination of an adaptive Monte Carlo sampler, called Gaussian Adaptation (GaA), and efficient exact stochastic simulation algorithms (SSA) that allows parameter identification from single stochastic trajectories. We benchmark the proposed method on a linear and a non-linear reaction network at steady state and during transient phases. In addition, we demonstrate that the present method also provides an ellipsoidal volume estimate of the viable part of parameter space and is able to estimate the physical volume of the compartment in which the observed reactions take place.

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.

Molecular discreteness is apparent in small-volume chemical systems, such as biological cells, leading to stochastic kinetics. Here we present a theoretical framework to understand the effects of discreteness on the steady state of a monostable chemical reaction network. We consider independent realizations of the same chemical system in compartments of different volumes. Rate equations ignore molecular discreteness and predict the same average steady-state concentrations in all compartments. However, our theory predicts that the average steady state of the system varies with volume: if a species is more abundant than another for large volumes, then the reverse occurs for volumes below a critical value, leading to a concentration inversion effect. The addition of extrinsic noise increases the size of the critical volume. We theoretically predict the critical volumes and verify, by exact stochastic simulations, that rate equations are qualitatively incorrect in sub-critical volumes.

Stochastic reaction-diffusion systems frequently exhibit behavior that is not predicted by deterministic simulation models. Stochastic simulation methods, however, are computationally expensive. We present a more efficient stochastic reaction-diffusion simulation algorithm that samples realizations from the exact solution of the reaction-diffusion master equation. The present algorithm, called partial-propensity stochastic reaction-diffusion (PSRD) method, uses an on-lattice discretization of the reaction-diffusion system and relies on partial-propensity methods for computational efficiency. We describe the algorithm in detail, provide a theoretical analysis of its computational cost, and demonstrate its computational performance in benchmarks. We then illustrate the application of PSRD to two- and three-dimensional pattern-forming Gray-Scott systems, highlighting the role of intrinsic noise in these systems.

We propose a new technique to measure the volume of adherent migrating cells. The method is based on a negative staining where a fluorescent, non-cell-permeant dye is added to the extracellular medium. The specimen is observed with a conventional fluorescence microscope in a chamber of uniform height. Given that the fluorescence signal depends on the thickness of the emitting layer, the objects excluding the fluorescent dye (i.e., cells) appear dark, and the decrease of the fluorescent signal with respect to the background is expected to give information about the height and the volume of the object. Using a glass microfabricated pattern with steps of defined heights, we show that the drop in fluorescence intensity is indeed proportional to the height of the step and obtain calibration curves relating fluorescence intensity to height. The technique, termed the fluorescence displacement method, is further validated by comparing our measurements with the ones obtained by atomic force microscopy (AFM). We apply our method to measure the real-time volume dynamics of migrating fish epidermal keratocytes subjected to osmotic stress. The fluorescence displacement technique allows fast and precise monitoring of cell height and volume, thus providing a valuable tool for characterizing the three-dimensional behaviour of migrating cells.

Influenza A virus (IAV) enters host cells by endocytosis followed by acid-activated penetration from late endosomes (LEs). Using siRNA silencing, we found that histone deacetylase 8 (HDAC8), a cytoplasmic enzyme, efficiently promoted productive entry of IAV into tissue culture cells, whereas HDAC1 suppressed it. HDAC8 enhanced endocytosis, acidification, and penetration of the incoming virus. In contrast, HDAC1 inhibited acidification and penetration. The effects were connected with dramatic alterations in the organization of the microtubule system, and, as a consequence, a change in the behavior of LEs and lysosomes (LYs). Depletion of HDAC8 caused loss of centrosome-associated microtubules and loss of directed centripetal movement of LEs, dispersing LE/LYs to the cell periphery. For HDAC1, the picture was the opposite. To explain these changes, centrosome cohesion emerged as the critical factor. Depletion of HDAC8 caused centrosome splitting, which could also be induced by depleting a centriole-linker protein, rootletin. In both cases, IAV infection was inhibited. HDAC1 depletion reduced the splitting of centrosomes, and enhanced infection. The longer the distance between centrosomes, the lower the level of infection. HDAC8 depletion was also found to inhibit infection of Uukuniemi virus (a bunyavirus) suggesting common requirements among late penetrating enveloped viruses. The results established class I HDACs as powerful regulators of microtubule organization, centrosome function, endosome maturation, and infection by IAV and other late penetrating viruses.

Prokaryotic environmental adaptations occur at different levels within cells to ensure the preservation of genome integrity, proper protein folding and function as well as membrane fluidity. Although specific composition and structure of cellular components suitable for the variety of extreme conditions has already been postulated, a systematic study describing such adaptations has not yet been performed. We therefore explored whether the environmental niche of a prokaryote could be deduced from the sequence of its proteome. Finally, we aimed at finding the precise differences between proteome sequences of prokaryotes from different environments.

Stochastic effects from correlated noise non-trivially modulate the kinetics of non-linear chemical reaction networks. This is especially important in systems where reactions are confined to small volumes and reactants are delivered in bursts. We characterise how the two noise sources confinement and burst modulate the relaxation kinetics of a non-linear reaction network around a non-equilibrium steady state. We find that the lifetimes of species change with burst input and confinement. Confinement increases the lifetimes of all species that are involved in any non-linear reaction as a reactant. Burst monotonically increases or decreases lifetimes. Competition between burst-induced and confinement-induced modulation may hence lead to a non-monotonic modulation. We quantify lifetime as the integral of the time autocorrelation function (ACF) of concentration fluctuations around a non-equilibrium steady state of the reaction network. Furthermore, we look at the first and second derivatives of the ACF, each of which is affected in opposite ways by burst and confinement. This allows discriminating between these two noise sources. We analytically derive the ACF from the linear Fokker-Planck approximation of the chemical master equation in order to establish a baseline for the burst-induced modulation at low confinement. Effects of higher confinement are then studied using a partial-propensity stochastic simulation algorithm. The results presented here may help understand the mechanisms that deviate stochastic kinetics from its deterministic counterpart. In addition, they may be instrumental when using fluorescence-lifetime imaging microscopy (FLIM) or fluorescence-correlation spectroscopy (FCS) to measure confinement and burst in systems with known reaction rates, or, alternatively, to correct for the effects of confinement and burst when experimentally measuring reaction rates.

Mesoscopic oscillatory reaction systems, for example in cell biology, can exhibit stochastic oscillations in the form of cyclic random walks even if the corresponding macroscopic system does not oscillate. We study how the intrinsic noise from molecular discreteness influences the frequency spectrum of mesoscopic oscillators using as a model system a cascade of coupled Brusselators away from the Hopf bifurcation. The results show that the spectrum of an oscillator depends on the level of noise. In particular, the peak frequency of the oscillator is reduced by increasing noise, and the bandwidth increased. Along a cascade of coupled oscillators, the peak frequency is further reduced with every stage and also the bandwidth is reduced. These effects can help understand the role of noise in chemical oscillators and provide fingerprints for more reliable parameter identification and volume measurement from experimental spectra.

Several real-world systems, such as gene expression networks in biological cells, contain coupled chemical reactions with a time delay between reaction initiation and completion. The non-Markovian kinetics of such reaction networks can be exactly simulated using the delay stochastic simulation algorithm (dSSA). The computational cost of dSSA scales with the total number of reactions in the network. We reduce this cost to scale at most with the smaller number of species by using the concept of partial reaction propensities. The resulting delay partial-propensity direct method (dPDM) is an exact dSSA formulation for well-stirred systems of coupled chemical reactions with delays. We detail dPDM and present a theoretical analysis of its computational cost. Furthermore, we demonstrate the implications of the theoretical cost analysis in two prototypical benchmark applications. The dPDM formulation is shown to be particularly efficient for strongly coupled reaction networks, where the number of reactions is much larger than the number of species.

Reaction-diffusion processes on complex deforming surfaces are fundamental to a number of biological processes ranging from embryonic development to cancer tumor growth and angiogenesis. The simulation of these processes using continuum reaction-diffusion models requires computational methods capable of accurately tracking the geometric deformations and discretizing on them the governing equations. We employ a Lagrangian level-set formulation to capture the deformation of the geometry and use an embedding formulation and an adaptive particle method to discretize both the level-set equations and the corresponding reaction-diffusion. We validate the proposed method and discuss its advantages and drawbacks through simulations of reaction-diffusion equations on complex and deforming geometries.

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.

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.

We present the partial-propensity stochastic simulation algorithm with composition-rejection sampling (PSSA-CR). It is an exact formulation of the stochastic simulation algorithm (SSA) for well-stirred systems of coupled chemical reactions. The new formulation is a partial-propensity variant [R. Ramaswamy, N. Gonzalez-Segredo, and I. F. Sbalzarini, J. Chem. Phys. 130, 244104 (2009)] of the composition- rejection SSA [A. Slepoy, A. P. Thompson, and S. J. Plimpton, J. Chem. Phys. 128, 205101 (2008)]. The computational cost of this new formulation is bounded by a constant for weakly coupled reaction networks, and it increases at most linearly with the number of chemical species for strongly coupled reaction networks. PSSA-CR thus combines the advantages of partial-propensity methods and the composition-rejection SSA, providing favorable scaling of the computational cost for all classes of reaction networks.

Live imaging of subcellular structures is indispensible to advance our understanding of cellular processes. The blurred digital images acquired in light microscopy are, however, complex to analyze, and identification and reconstruction of subcellular structures from such images remains a major challenge. We present a novel, model-based image analysis algorithm to reconstruct outlines of subcellular structures using a sub-pixel representation. The algorithm explicitly accounts for the optical properties of the microscope. We validate the reconstruction performance on synthetic data and apply the new method to fluorescence microscopy images of endosomes identified by the GTPase EGFP-Rab5. The benefits of the new algorithm are outlined by comparison to standard techniques. We demonstrate that the new algorithm leads to better discrimination between different endosomal virus entry pathways and to more robust, accurate, and self-consistent quantification of endosome shape features. This allows establishing a set of features that quantify endosome morphology and robustly capture the dynamics of endosome fusion.

The concept of high-resolution shapes (also referred to as folds or states, depending on the context) of a polymer chain plays a central role in polymer science, structural biology, bioinformatics, and biopolymer dynamics. However, although the idea of shape is intuitively very useful, there is no unambiguous mathematical definition for this concept. In the present work, the distributions of high-resolution shapes within the ideal random-walk ensembles with N=3,...,6 beads (or up to N=10 for some properties) are investigated using a systematic (grid-based) approach based on a simple working definition of shapes relying on the root-mean-square atomic positional deviation as a metric (i.e., to define the distance between pairs of structures) and a single cutoff criterion for the shape assignment. Although the random-walk ensemble appears to represent the paramount of homogeneity and randomness, this analysis reveals that the distribution of shapes within this ensemble, i.e., in the total absence of interatomic interactions characteristic of a specific polymer (beyond the generic connectivity constraint), is significantly inhomogeneous. In particular, a specific (densest) shape occurs with a local probability that is 1.28, 1.79, 2.94, and 10.05 times (N=3,...,6) higher than the corresponding average over all possible shapes (these results can tentatively be extrapolated to a factor as large as about 10(28) for N=100). The qualitative results of this analysis lead to a few rather counterintuitive suggestions, namely, that, e.g., (i) a fold classification analysis applied to the random-walk ensemble would lead to the identification of random-walk "folds;" (ii) a clustering analysis applied to the random-walk ensemble would also lead to the identification random-walk "states" and associated relative free energies; and (iii) a random-walk ensemble of polymer chains could lead to well-defined diffraction patterns in hypothetical fiber or crystal diffraction experiments. The inhomogeneous nature of the shape probability distribution identified here for random walks may represent a significant underlying baseline effect in the analysis of real polymer chain ensembles (i.e., in the presence of specific interatomic interactions). As a consequence, a part of what is called a polymer shape may actually reside just "in the eye of the beholder" rather than in the nature of the interactions between the constituting atoms, and the corresponding observation-related bias should be taken into account when drawing conclusions from shape analyses as applied to real structural ensembles.

We introduce an alternative formulation of the exact stochastic simulation algorithm (SSA) for sampling trajectories of the chemical master equation for a well-stirred system of coupled chemical reactions. Our formulation is based on factored-out, partial reaction propensities. This novel exact SSA, called the partial-propensity direct method (PDM), is highly efficient and has a computational cost that scales at most linearly with the number of chemical species, irrespective of the degree of coupling of the reaction network. In addition, we propose a sorting variant, SPDM, which is especially efficient for multiscale reaction networks.

We study the thermophoretic motion of water nanodroplets confined inside carbon nanotubes using molecular dynamics simulations. We find that the nanodroplets move in the direction opposite the imposed thermal gradient with a terminal velocity that is linearly proportional to the gradient. The translational motion is associated with a solid body rotation of the water nanodroplet coinciding with the helical symmetry of the carbon nanotube. The thermal diffusion displays a weak dependence on the wetting of the water-carbon nanotube interface. We introduce the use of the moment scaling spectrum (MSS) in order to determine the characteristics of the motion of the nanoparticles inside the carbon nanotube. The MSS indicates that affinity of the nanodroplet with the walls of the carbon nanotubes is important for the isothermal diffusion and hence for the Soret coefficient of the system.

Biological trajectories can be characterized by transient patterns that may provide insight into the interactions of the moving object with its immediate environment. The accurate and automated identification of trajectory motifs is important for the understanding of the underlying mechanisms. In this work, we develop a novel trajectory segmentation algorithm based on supervised support vector classification. The algorithm is validated on synthetic data and applied to the identification of trajectory fingerprints of fluorescently tagged human adenovirus particles in live cells. In virus trajectories on the cell surface, periods of confined motion, slow drift, and fast drift are efficiently detected. Additionally, directed motion is found for viruses in the cytoplasm. The algorithm enables the linking of microscopic observations to molecular phenomena that are critical in many biological processes, including infectious pathogen entry and signal transduction.

We present a computational particle method for the simulation of isotropic and anisotropic diffusion on curved biological surfaces that have been reconstructed from image data. The method is capable of handling surfaces of high curvature and complex shape, which are often encountered in biology. The method is validated on simple benchmark problems and is shown to be second-order accurate in space and time and of high parallel efficiency. It is applied to simulations of diffusion on the membrane of endoplasmic reticula (ER) in live cells. Diffusion simulations are conducted on geometries reconstructed from real ER samples and are compared to fluorescence recovery after photobleaching experiments in the same ER samples using the transmembrane protein tsO45-VSV-G, C-terminally tagged with green fluorescent protein. Such comparisons allow derivation of geometry-corrected molecular diffusion constants for membrane components from fluorescence recovery after photobleaching data. The results of the simulations indicate that the diffusion behavior of molecules in the ER membrane differs significantly from the volumetric diffusion of soluble molecules in the lumen of the same ER. The apparent speed of recovery differs by a factor of approximately 4, even when the molecular diffusion constants of the two molecules are identical. In addition, the specific shape of the membrane affects the recovery half-time, which is found to vary by a factor of approximately 2 in different ER samples.

The lateral mobility of individual murine polyoma virus-like particles (VLPs) bound to live cells and artificial lipid bilayers was studied by single fluorescent particle tracking using total internal reflection fluorescence microscopy. The particle trajectories were analyzed in terms of diffusion rates and modes of motion as described by the moment scaling spectrum. Although VLPs bound to their ganglioside receptor in lipid bilayers exhibited only free diffusion, analysis of trajectories on live 3T6 mouse fibroblasts revealed three distinct modes of mobility: rapid random motion, confined movement in small zones (30-60 nm in diameter), and confined movement in zones with a slow drift. After binding to the cell surface, particles typically underwent free diffusion for 5-10 s, and then they were confined in an actin filament-dependent manner without involvement of clathrin-coated pits or caveolae. Depletion of cholesterol dramatically reduced mobility of VLPs independently of actin, whereas inhibition of tyrosine kinases had no effect on confinement. The results suggested that clustering of ganglioside molecules by the multivalent VLPs induced transmembrane coupling that led to confinement of the virus/receptor complex by cortical actin filaments.

The determination of diffusion coefficients from fluorescence recovery data is often complicated by geometric constraints imposed by the complex shapes of intracellular compartments. To address this issue, diffusion of proteins in the lumen of the endoplasmic reticulum (ER) is studied using cell biological and computational methods. Fluorescence recovery after photobleaching (FRAP) experiments are performed in tissue culture cells expressing GFP-KDEL, a soluble, fluorescent protein, in the ER lumen. The three-dimensional (3D) shape of the ER is determined by confocal microscopy and computationally reconstructed. Within these ER geometries diffusion of solutes is simulated using the method of particle strength exchange. The simulations are compared to experimental FRAP curves of GFP-KDEL in the same ER region. Comparisons of simulations in the 3D ER shapes to simulations in open 3D space show that the constraints imposed by the spatial confinement result in two- to fourfold underestimation of the molecular diffusion constant in the ER if the geometry is not taken into account. Using the same molecular diffusion constant in different simulations, the observed speed of fluorescence recovery varies by a factor of 2.5, depending on the particular ER geometry and the location of the bleached area. Organelle shape considerably influences diffusive transport and must be taken into account when relating experimental photobleaching data to molecular diffusion coefficients. This novel methodology combines experimental FRAP curves with high accuracy computer simulations of diffusion in the same ER geometry to determine the molecular diffusion constant of the solute in the particular ER lumen.

This paper presents a computationally efficient, two-dimensional, feature point tracking algorithm for the automated detection and quantitative analysis of particle trajectories as recorded by video imaging in cell biology. The tracking process requires no a priori mathematical modeling of the motion, it is self-initializing, it discriminates spurious detections, and it can handle temporary occlusion as well as particle appearance and disappearance from the image region. The efficiency of the algorithm is validated on synthetic video data where it is compared to existing methods and its accuracy and precision are assessed for a wide range of signal-to-noise ratios. The algorithm is well suited for video imaging in cell biology relying on low-intensity fluorescence microscopy. Its applicability is demonstrated in three case studies involving transport of low-density lipoproteins in endosomes, motion of fluorescently labeled Adenovirus-2 particles along microtubules, and tracking of quantum dots on the plasma membrane of live cells. The present automated tracking process enables the quantification of dispersive processes in cell biology using techniques such as moment scaling spectra.

Polarized cells frequently use diffusion barriers to separate plasma membrane domains. It is unknown whether diffusion barriers also compartmentalize intracellular organelles. We used photobleaching techniques to characterize protein diffusion in the yeast endoplasmic reticulum (ER). Although a soluble protein diffused rapidly throughout the ER lumen, diffusion of ER membrane proteins was restricted at the bud neck. Ultrastructural studies and fluorescence microscopy revealed the presence of a ring of smooth ER at the bud neck. This ER domain and the restriction of diffusion for ER membrane proteins through the bud neck depended on septin function. The membrane-associated protein Bud6 localized to the bud neck in a septin-dependent manner and was required to restrict the diffusion of ER membrane proteins. Our results indicate that Bud6 acts downstream of septins to assemble a fence in the ER membrane at the bud neck. Thus, in polarized yeast cells, diffusion barriers compartmentalize the ER and the plasma membrane along parallel lines.