Takumi Washio

Screening system for drug-induced arrhythmogenic risk combining a patch clamp and heart simulator

Finding the silent skipped beat: Predicting arrhythmia-causing drugs using a high-throughput hybrid heart simulator. To save time and cost for drug discovery, a paradigm shift in cardiotoxicity testing is required. We introduce a novel screening system for drug-induced arrhythmogenic risk that combines in vitro pharmacological assays and a multiscale heart simulator. For 12 drugs reported to have varying cardiotoxicity risks, dose-inhibition curves were determined for six ion channels using automated patch clamp systems. By manipulating the channel models implemented in a heart simulator consisting of more than 20 million myocyte models, we simulated a standard electrocardiogram (ECG) under various doses of drugs. When the drug concentrations were increased from therapeutic levels, each drug induced a concentration-dependent characteristic type of ventricular arrhythmia, whereas no arrhythmias were observed at any dose with drugs known to be safe. We have shown that our system combining in vitro and in silico technologies can predict drug-induced arrhythmogenic risk reliably and efficiently.

A three-dimensional simulation model of cardiomyocyte integrating excitation-contraction coupling and metabolism.

Recent studies have revealed that Ca(2+) not only regulates the contraction of cardiomyocytes, but can also function as a signaling agent to stimulate ATP production by the mitochondria. However, the spatiotemporal resolution of current experimental techniques limits our investigative capacity to understand this phenomenon. Here, we created a detailed three-dimensional (3D) cardiomyocyte model to study the subcellular regulatory mechanisms of myocardial energetics. The 3D cardiomyocyte model was based on the finite-element method, with detailed subcellular structures reproduced, and it included all elementary processes involved in cardiomyocyte electrophysiology, contraction, and ATP metabolism localized to specific loci. The simulation results were found to be reproducible and consistent with experimental data regarding the spatiotemporal pattern of cytosolic, intrasarcoplasmic-reticulum, and mitochondrial changes in Ca(2+); as well as changes in metabolite levels. Detailed analysis suggested that although the observed large cytosolic Ca(2+) gradient facilitated uptake by the mitochondrial Ca(2+) uniporter to produce cyclic changes in mitochondrial Ca(2+) near the Z-line region, the average mitochondrial Ca(2+) changes slowly. We also confirmed the importance of the creatine phosphate shuttle in cardiac energy regulation. In summary, our 3D model provides a powerful tool for the study of cardiac function by overcoming some of the spatiotemporal limitations of current experimental approaches.

Ordering strategies and related techniques to overcome the trade-off between parallelism and convergence in incomplete factorizations

This paper is concerned with the parallel implementation of the incomplete factorization preconditioned iterative method. Although the use of such parallel ordering as multicolor ordering may increase parallelism in factorization, it often slows convergence when used in the preconditioned method, and thus may offset the gain in speed obtained with parallelization. Further, the higher the parallelism of an ordering, the slower the convergence; the lower the parallelism, the faster the convergence. This well-known trade-off between parallelism and convergence is well explained by the property of compatibility, the level of which can be clearly seen when ordering is presented in graph form (S. Doi, A. Lichnewsky, A graph-theory approach for analyzing the effects of ordering on ILU preconditioning, INRIA report 1452, 1991). In any given method, the fewer the incompatible local graphs in an ordering (i.e., the lower the parallelism), the faster the convergence (S. Doi, Appl. Numer. Math. 7 (1991) 417–436; S. Doi, in: T. Nodera (Ed.), Advances in Numerical Methods for Large Sparse Sets of Linear Systems, 7 Keio University, 1991). An ordering with no incompatible local graphs, for example, such as that implemented on vector multiprocessors by using the nested dissection technique, will have excellent convergence, but its parallelism will be limited (S. Doi, A. Lichnewsky, Int. J. High Speed Comput. 2 (1990) 143–179). To attain a better balance, a certain degree of incompatibility is necessary. In this regard, increasing the number of colors in multicolor ordering can be a useful approach (S. Fujino, S. Doi, in: R. Beauwens (Ed.), Proceeding of the IMACS Internation Symposium on Iterative Methods in Linear Algebra, March 1991; S. Doi, A. Hoshi, Int. J. Comput. Meth. 44 (1992) 143–152). Two related techniques also presented here are the overlapped multicolor ordering (T. Washio, K. Hayami, SIAM J. Sci. Comput. 16 (1995) 631–650), and a fill-in strategy selectively applied to incompatible local graphs. Experiments conducted with an SX-5/16A vector parallel supercomputer show the relative effectiveness of increasing the number of colors and also of using this approach in combination with overlapping and with fill-ins.

Multi-scale simulations of cardiac electrophysiology and mechanics using the University of Tokyo heart simulator

The importance and need for an integrative mathematical modeling approach in the biological and medical fields is currently well recognized. Such an approach is crucial in understanding the complexity of hierarchical biological systems increasingly revealed by active researches in molecular and cellular biology. Particularly in cardiac functioning, modeling must cover such diverse phenomena as solid mechanics, fluid dynamics, electricity and biochemistry. Recent advancements in computational science and the development of high-performance computers have enabled the creation of multi-scale, multi-physics simulation heart models using the finite element method. Although whole heart or ventricular models of electrophysiology involving electro-mechanics with or without blood flow dynamics have been reported, to our knowledge no single model has yet succeeded in completely reproducing the behavior of the heart from the subcellular to whole organ levels. In this article, we present a brief methodology-focused review on some of the essential components for multi-scale, multi-physics heart modeling. A perspective of heart modeling in the era of high performance computing is also presented.

Fourier Analysis of GMRES( m ) Preconditioned by Multigrid

This paper deals with convergence estimates of GMRES(m) [Saad and Schultz, { SIAM J. Sci. Statist. Comput., 7 (1986), pp. 856--869] preconditioned by multigrid [Brandt, Math. Comp., 31 (1977), pp. 333--390], [Hackbusch, Multi-Grid Methods and Applications, Springer, Berlin, 1985]. Fourier analysis is a well-known and useful tool in the multigrid community for the prediction of two-grid convergence rates [Brandt, Math. Comp., 31 (1977), pp. 333--390], [Stuben and Trottenberg, in Multigrid Methods, Lecture Notes in Math. 960, K. Stuben and U. Trottenberg, eds., Springer, Berlin, pp. 1--176]. This analysis is generalized here to the situation in which multigrid is a preconditioner, since it is possible to obtain the whole spectrum of the two-grid iteration matrix. A preconditioned Krylov subspace acceleration method like GMRES(m) implicitly builds up a minimal residual polynomial. The determination of the polynomial coefficients is easily possible and can be done explicitly since, from Fourier analysis, a simple block-diagonal two-grid iteration matrix results. Based on the GMRES(m) polynomial, sharp theoretical convergence estimates can be obtained which are compared with estimates based on the spectrum of the iteration matrix. Several numerical scalar test problems are computed in order to validate the theoretical predictions.

Transmural and apicobasal gradients in repolarization contribute to T wave genesis in human surface ECG

The cellular basis of the T-wave morphology of surface ECG remains controversial in clinical cardiology. We examined the effect of action potential duration (APD) distribution on T-wave morphology using a realistic model of the human ventricle and torso. We developed a finite-element model of the ventricle consisting of ∼26 million elements, including the conduction system, each implemented with the ion current model of cardiomyocytes. This model was embedded in a torso model with distinct organ structures to obtain the standard ECG leads. The APD distribution was changed in the transmural direction by locating the M cells in either the endocardial or epicardial region. We also introduced apicobasal gradients by modifying the ion channel parameters. Both the transmural gradient (with M cells on the endocardial side) and the apicobasal gradient produced positive T waves, although a very large gradient was required for the apicobasal gradient. By contrast, T waves obtained with the transmural gradient were highly symmetric and, therefore, did not represent the true physiological state. Only combination of the transmural and the moderate apicobasal gradients produced physiological T waves in surface ECG. Positive T waves in surface ECG mainly originated from the transmural distribution of APD with M cells on the endocardial side, although the apicobasal gradient was also required to attain the physiological waveform.

A Parallel Multilevel Technique for Solving the Bidomain Equation on a Human Heart with Purkinje Fibers and a Torso Model

In this paper, we present a multigrid method and its implementation on parallel computers to solve the bidomain equation that appears in excitation propagation analysis of the human heart with the torso. The bidomain equation is discretized with the finite element method on a composite mesh composed of a fine voxel mesh around the heart and a coarse voxel mesh covering the torso. The extracellular potential problem on the torso is formulated as a variational problem with a constraint at the interface of the fine and coarse meshes. We show that this formulation naturally satisfies the conservation property of the electric currents and fits into the multilevel adaptive solution technique framework. We also present our special treatment of the Purkinje fiber network in the multigrid algorithm where it is modeled as multiway branching lines connected to the nodes in the voxel mesh of the heart. A parallel implementation of the proposed multigrid algorithm on distributed memory computers is presented and its performance is evaluated using real-life applications.

Parallel block preconditioning based on SSOR and MILU

Two kinds of parallel preconditioners for the solution of large sparse linear systems which arise from the 2-D 5-point finite difference discretization of a convection-diffusion equation are introduced. The preconditioners are based on the SSOR or MILU preconditioners and can be implemented on parallel computers with distributed memories. One is the block preconditioner, in which the interface components of the coefficient matrix between blocks are ignored to attain parallelism in the forward-backward substitutions. The other is the modified block preconditioner, in which the block preconditioner is modified by taking the interface components into account. The effect of these preconditioners on the convergence of preconditioned iterative methods and timing results on the parallel computer (Cenju) are presented.

Multiscale Heart Simulation with Cooperative Stochastic Cross-Bridge Dynamics and Cellular Structures

In this paper, we introduce a multiscale heart simulation method that integrates the three-scale phenomena from the microscopic stochastic sarcomere kinetics to the macroscopic heartbeat via the mesoscopic myocardial cell assembly. In sarcomere kinetics, the stochastic behavior of the myosin heads in the cooperative cross-bridge formation and the strain-dependent head rotations are directly simulated using a Monte Carlo (MC) method. This stochastic sarcomere model is coupled with contractile myofibril elements in a finite element (FE) cell assembly model. The cell assembly model is further coupled with a macroscopic FE heart model by means of the homogenization method. We bridge the large gap in the time step sizes between the MC and FE models with an idea based on impulse equilibrium between the sum of the stretches over the whole myosin arms and the contractile stress of the myofibril element. In the application of the homogenization method to the meso-macro coupling, a novel approach is introduced to d...

Approximation for Cooperative Interactions of a Spatially-Detailed Cardiac Sarcomere Model

We developed a novel ordinary differential equation (ODE) model, which produced results that correlated well with the Monte Carlo (MC) simulation when applied to a spatially-detailed model of the cardiac sarcomere. Configuration of the novel ODE model was based on the Ising model of myofilaments, with the “co-operative activation” effect introduced to incorporate nearest-neighbor interactions. First, a set of parameters was estimated using arbitrary Ca transient data to reproduce the combinational probability for the states of three consecutive regulatory units, using single unit probabilities for central and neighboring units in the MC simulation. The parameter set thus obtained enabled the calculation of the state transition of each unit using the ODE model with reference to the neighboring states. The present ODE model not only provided good agreement with the MC simulation results but was also capable of reproducing a wide range of experimental results under both steady-state and dynamic conditions including shortening twitch. The simulation results suggested that the nearest-neighbor interaction is a reasonable approximation of the cooperativity based on end-to-end interactions. Utilizing the modified ODE model resulted in a reduction in computational costs but maintained spatial integrity and co-operative effects, making it a powerful tool in cardiac modeling.