Systems of Preventive Cardiological Monitoring: Models, Algorithms, First Results, and Perspectives

Abstract

The results of work on creating methods, models, and computational algorithms for remote preventive health-monitoring systems are presented, in particular, cardiac preventive monitoring. The main attention is paid to the models and computational algorithms of preventive monitoring, the interaction of the computing kernels of a remote cluster with portable ECG recorders, implantable devices, and sensors. Computational kernels of preventive monitoring are a set of several thousand interacting automata of analog of Turing machines, recognizing the characteristic features and evolution of the hidden predictors of atrial fibrillation(AF), ventricular tachycardia or fibrillation (VT-VF), sudden cardiac death, and heart failure (HF) revealed by them. The estimation of the time for reaching the heart events boundaries is calculated on the basis of the evolution equations for the ECG multi-trajectories determined by recognizing automata. Evaluation time of heart event (HE) boundaries to achieve is calculated on the basis of the evolution equations for ECG multi-paths defined by recognizing machines. Ultimately, the computational cores reconstruct the ECG of the forecast and give temporary estimates of its achievement. Cloud computing cluster supports low-cost ECG ultra-portable recorders and does not limit the possibilities of using a more complex patient telemetry containing wearable and implantable devices: CRT and ICD, CardioMEMS HF System, and so on.

1. Introduction

The growing interest in remote cardiac-monitoring systems is associated primarily with the need for an early prognosis in the development of heart failure (HF) disease [1] and an early prognosis of developing against the background HF of such heart events as atrial fibrillation (AF), ventricular tachycardia or fibrillation (VT-VF), and sudden cardiac death [2]. In addition, on the assumption that remote monitoring systems can effectively cope with prognostic tasks, the next call to remote monitoring systems occurs namely the development and optimization of risk reduction strategies and the strategy of drug or device cardioversions. It is about the management of the patient’s cardiac condition through a drug or a device therapy, determining the degree of effectiveness of therapy and predicting the results of treatment [3]. However, despite the significant advances and new promising results [4] in the field of management and effective prevention, there are a significant number of problems associated with the lack of a physiologically justified mathematical model of management. There is some analog of the problems described earlier in the problems of engineering asset management (EAM). The basis of ЕАМ is Prognostics and Health Management (PHM) systems, including also algorithms of predictive analytics, big data, deep learning, and so on. Another component of ЕАМ is represented by Intelligent Maintenance Systems (IMS); here the main goal is to develop systems of preventive maintenance, self-maintenance, and self-recovery systems. Technical systems are somewhat simpler than biomechanical systems, so the mathematical PHM models, developed for the prognosis, are more formalized, including in the field of physics of failure. A simple mechanical transfer of PHM models is hardly possible, but some useful analogies are appropriate. In particular, in setting the tasks of calculating optimizing strategies for preventive maintenance, estimating the time to reach the boundaries of the mechanisms dysfunction and the failure boundary, and so on. In this chapter, some useful formulations and models of the PHM will be used as applied. The noted analogy can be supplemented by fact that used in cardioversion the devices are both wearable and implanted also need a prognosis of their failures and dysfunctions. As a result, the complex task of constructing PHM models for the system “biomechanical object plus implantable device or sensor” is relevant for the further development of preventive cardiac-monitoring systems.

Finally, the goal of PCM applications in medicine, using mHealth and eHealth platforms with telemetry transfer capabilities to a remote server and a large computational resource in the form of distributed computing systems, the cloud computing service, is reduced to the following:

1. based on the chronological database of each patient and current data, to calculate the parameters or predictors that characterize the evolution of HF;

2. on the same database, predict the appearance of hidden HF predictors (i.e., the appearance and evolution of those ECG signal characteristics that precede the appearance of HF);

3. to reveal the hidden predictors of cardiac events, namely AF, VT-VF, and sudden cardiac death;

4. to determine evolutionary equations for hidden predictors and estimate the time to reach the predicted characteristics of the ECG.

Similarly, remote systems of preventive monitoring should be aimed at solving management tasks, that is, constructing a strategy for variable parameters with drug therapy or programming CRT, IDC devices.

If to consider the problems noted earlier in the context of preventive monitoring, then the question of developing the systems of preventive maintenance strategies must be raised.

Thus, the prognosis of heart conditions is now reduced to the identification of predictors (ECG characteristics) cardiac events, which, in fact, have to answer the question of how likely the presence and numerical characteristics of ECG predictors provide the appearance of cardiac events. Such predictors include the length of the CT of the interval, the CT dispersion, the P-wave index, and so on. A new trend, known as nonlinear dynamics methods, adds to the existing predictors of the entropy characteristics commonly used in nonlinear dynamics for the classification of trajectories of dynamic systems. These include fractal, dimensional, and entropy characteristics, for example, information dimension, approximation entropy, and so on. The methods of symbolic dynamics are also used, the essence of which is reduced to the study of dynamical systems on the basis of the analysis of symbolic cascades. There are some results in the cases where the cascade is a topological Markov chain. At the same time, transferring results from cascades to a dynamic system requires that the dynamic symbolic system is the factor system.

In addition, there is no model of cardiac activity tied to a specific system of dynamic equations. Ultimately, the prediction achieved by measuring or calculating the predictors mentioned requires a number of conditional transitions, that is, the fulfillment of a multitude of conditions involving various facts from an anamnesis, the general state of the organism, the presence or absence of a range of diseases, which reduces the effectiveness of the prediction up to the phrase “positive prognosis“ and “negative prognosis.“.

Section 1 briefly describes the basic model of the propagation of the action potential (AP). The presented equation for the AP is purely a model, more realistic models; for example, Microdomain model [5] is not considered since the main purpose of this section is to demonstrate the result of numerical simulation of the propagation of AP on the basis of bundle cellular automata. Equations for ion currents are placed in the bundle of automata, which makes it possible to locally vary the parameters of the ion current models, the degree of anisotropy, and the local geometry of the cardio tissue [6]. In the context of this chapter, this section allows to refer to the microscopic theory of ion currents and action potential, without which the formal operation with this series of ECG wavelet coefficients is enriched by an innumerable number of possible scenarios for constructing recognizing automata.

Section 2 describes the evolution of wavelet coefficients in models of random walking on a multidimensional lattice and on a multidimensional continuum. With some realistic assumptions, it is possible to estimate the time to reach the boundaries of cardiac events or HF. The transition to high dimensions is also due to attempts of detection predictors of trajectories for the prognosis of rotors in AF.

Approximation entropy is given considerable attention in [7]. However, approximation entropy in principle is conceived as some estimate of Kolmogorov complexity. Moreover, K-complexity is an algorithmically unsolvable quantity, which encourages researchers to seek estimates for K-complexity in the form of various kinds of entropies, for example, Shannon entropy.

It should be noted that in reality, an estimate of not K-complexity but of conditional K-complexity is required, which induces the transitions from the time series of wavelet coefficients to vector processes.

There are other good reasons to go to the vector process in connection with the prognostic tasks and estimation of the time to reach. The automatic transfer of the random walk theory over multidimensional lattices and the multidimensional continuum [8, 9] is hampered by the fact that these models were created for problems in the theory of polymers, where the length of the jump is constant and equal to the length of the monomer. The situation is partly saved by the introduction of the distribution function of jumps along the lengths. For the visibility of finite formulas, it was necessary to introduce the averaged length of the jump. However, the main goal of this section is to construct multidimensional state spaces and then sets of trajectories on them. In such spaces of HF regions, the cardiac events are clearly distinguished by wavelet coefficient values. However, the Euclidean metric is not quite suitable for estimating the closeness of trajectories. Therefore, in Section 2, we give an example of the bifurcations of the distribution of density function on the basis of the elementary catastrophe theory.

Section 3 describes the construction of family of recognition automata and gives an example of an evolution equation for the internal states of automata. Here, the models described in Sections. 1 and 2 contribute to the understanding and further formalization of the principles of constructing recognizing automata and their operation. Since in algorithmic realizations of automata, only principles are laid down. All their further activities including reproduction, an increase in complexity, self-learning, and adaptation (in the context of personalized cardiology) should occur as a result of self-organization based on the principles laid down in the algorithms.

Section 5 describes the interaction of recognizing automata with peripheral devices, wearable or implanted, and also possibly other body sensors. It shows how the transmission and processing of ECG signals occurs and as a family of automata forms signaling automata and places them on portable devices. In turn, the automata, delivered to the device, control the calculations on a remote family of automata.

Section 6 is devoted to verification and discussion of the results and conditions for a correct prognosis.

Section 7 is devoted to the formalization of trajectories management tasks. The management model is outlined in the language of homotopy theory and the theory of infinite loop spaces. It considers the set of all admissible trajectories and the possibility of continuous deformation of some trajectories into others, as well as possible obstacles to such deformation. It is intuitively clear that the singularities that arouse in the model must be associated with the birth of filaments in models AF-VТ. However, this fact is not rigorously proved in this chapter. Hope is encouraged by the fact that filaments, like singularities in the base space, are of a homotopic nature in both cases and, therefore, are invariant with respect to any deformations, which means that they can be recorded by ECG analysis on the body.

The traditional approach to the diagnosis and prognosis of heart event (HE) based on an ECG can be found in [10, 11] and references cited therein. Summarizing the goal of PHM applications in preventive cardiac monitoring is reduced not only to accurate estimates of the time to reach the boundaries of cardiac events. One of PHM tasks is managing the state of complex systems and determining optimal management strategies to extend the lifetime of complex technical objects. Consequently, in the PHM cardiac applications, the PHM task is also the creation of models and management algorithms, that is, the retention of multitrajectories of cardiac states in classes of trajectories without cardiac events. How much such task is possible is discussed in the last section of this chapter in the framework of the topological model of trajectory management.

2. Basis

The construction of automata recognizing hidden, early predictors of cardiac events is determined by a set of requirements that follow from the basic models and representations. The basic principles should include the basic models describing the propagation of the action potential in cardio tissue [12, 13] and the models of ionic membrane currents [14]. Formula (1) represents a simplified equation for the action potential

$$\frac{\partial V\left(r,t\right)}{\partial t}=-D∆V\left(r,t\right)+I_{\mathit{ion}}\left(r\right)$$E1

V is the transmembrane potential, D is the homogeneous pseudo-diffusion constant of the intracellular gap junctional coupling, I_ion (r) is the total transmembrane ionic current, and $\mathit{∆}$ is the Laplace operator.

There are also more general models that take into account inhomogeneous cellular structure, mentioned in Section “Introduction.”

Using the models of ion transport through biological membranes based on models of statistical thermodynamics of nonequilibrium processes [13], a digital model of bundle cellular automata is created. The bundle cellular automata in their continuum limit represent a smooth fiber bundle [15]. In this case, the diffusion part of the “reaction–diffusion” equation is determined in the base manifold of the bundle. In the fibers, the anisotropic architecture of myocardial regions was modeled. Different fibers of bundle cellular automata corresponded to different degrees of anisotropy and different geometric configurations. The results of numerical simulation are shown in Figure 1.

Figure 1.

The results of numerical simulation.

The main result obtained by numerical simulation is reduced to demonstrating the growth of wavefront fluctuations and fluctuations in the curvature of the wavefront after repeated passage of the front along the same macroscopic myocardial regions, that is, the model of bundle cellular automata, and in the continuum limit of smooth fiber bundle is an acceptable approximation in the problem of macroscopic description of the evolution of the wavefront of the AP, taking into account the changes taking place on the cellular, microscopic scale, when changing the microgeometry of cell bonds and the variation in the parameters of the microscopic model of intercellular conductivity by multiplicative noise, as was done in the work [16] in the model of Khodzhin-Huxley axon with Markov dichotomous voltage noise.

Models of the AP by individual cells and its conduction from cell to cell through intercellular gap junctions are discussed in the work in [12].

3. Model

Another basis for constructing prediction algorithms and recognizing automata is the random walk model on a multidimensional lattice or in a multidimensional continuum. Thus, at each set of R-R intervals of a fixed ECG length, the signal is represented as a set of its wavelet coefficients

$$\left\{{}_{\mathit{\text{Hist}}}{}^k{W}_{i,j}^N\right\},N=1,2,3,\dots ,N^{\ast }$$E2

$N$ is the number of R-R cycle, $i,j$ are indexes of wavelet decomposition, $\mathit{\text{Hist}}$ is the heart rate histogram column index, and $k$ is the numbering of vectors from wavelet coefficients of dimension $N^{\ast }$.

Moreover, for all fixed indices except $N$, stochastic process with discrete time$,N$ cascade, is determined. Further, fixing the limiting value $N$ as $N^{\ast }$ is determined by a set of vectors of dimension $N^{\ast }$, chosen from the consistent values of the process under consideration with discrete time. As a result, the space ${\mathit{R}}^{{\mathit{N}}^{\mathit{\ast }}}$ is determined, consisting of all possible finite segments of dimension размерности $N^{\ast }$

$$\left\{{\mathit{R}}_{\mathit{k}}\right\}\stackrel{\scriptscriptstyle\mathrm{def}}{=}\left\{{}_{\mathit{\text{Hist}}}{}^k{W}_{i,j}^N:k{N}^{\ast }\le N\le \left(k+1\right)N^{\ast },k=0,1,2,3\dots \right\},{\mathit{R}}_{\mathit{k}}\mathbf{\in }{\mathit{R}}^{{\mathit{N}}^{\mathbf{\ast }}}$$E3

The total number of such spaces ${\mathit{R}}^{{\mathit{N}}^{\mathbf{\ast }}}$ corresponds to the product of the number of elements of sets $\left\{\mathit{\text{Hist}}\right\},\left\{i\right\},\left\{j\right\}.$

Thus, the formed product of spaces contains the set of all admissible values of any ECG signals (set of admissible states). In this set, there is a subset containing sequences of cardiac events (CE) и HF predictors. Such a subset is formed by wavelet coefficients from the ECG database of cardiac events и HF predictors. Thus, each finite segment of any process is typed by elementary transitions in each coordinate ${\mathit{R}}_{\mathit{k}}\mathbf{\in }{\mathit{R}}^{{\mathit{N}}^{\mathbf{\ast }}}$, if the probability of such transition is known. For the transition, one vector of dimension $N^{\ast }$ to another is required,$N^{\ast }$ transitions. As a result, the initial vector goes into another vector. If the probabilities of elementary transitions are known, the probability from vector to vector is determined by the product of elementary coordinate-wise probabilities. Successive transitions from one vector to another form of the trajectory in space of the product. In each product space, the probability density of the transition is determined, expressed as a Feynman path integral. A vector in each fixed product space is defined as a state, and a sequence of vectors determines the state trajectory. The set of states in the product of spaces defines a multistate, a multitrajectory, respectively.

The introduction of the probability density of transition from one state to another is interpreted as a random walk of a trajectory in state spaces. Thus, the problem of predicting cardiac events is reduced to the calculation of the probability of transition to the boundary of the region obtained by mapping the ECG from the base of cardiac events into a product of spaces. In addition to a subset of cardiac events clearly expressed by their ECG (AF-VF) in the space of events, regions corresponding to various changes in the point-wise Holder regularity are also separated relatively to the regularity characteristics for the ECG norm. The procedure of allocation of such regions is performed during the monitoring and based on two microlocal spaces [17]. Therefore, the prognosis is reduced to calculating the probability density function (PDF) for transition probabilities from one state to another (${\mathit{R}}_{\mathit{k}}\mathit{\to }{\mathit{R}}_{\mathit{k}\mathit{+}\mathbf{1}})$. Since the regions of cardiac events and the irregular behavior of the ECG signal are separated, the probability of transition of the state to the boundary of the “critical” regions in L (L is the continuous analog of the number of cardiac cycles) steps is calculated. At the difficulty of calculating the probability density function, calculations are carried out for the moments of the PDF or the cumulants. The time required to reach the “critical region” is determined from explicit analytical expressions for the second moment of PDF.

Evolutionary equations for the probability density are derived from the Feynman representation of the PDF under certain limiting assumptions, in particular, the Fokker-Planck equation, and the equations for single-step processes. The calculation of the Feynman integral also reduces to solutions of the Hamilton-Jacobi equation, the Schrödinger equation, and the WKB method [18].

However, in the multidimensional case, the solution of the equations is possible only numerically and generates some computational problems.

In [19], on the basis of modification of random walk models [8, 9], the PDF of transition probabilities (probability function of end-to-end distribution) представлено as the Feynman path integral

$$P\left({\mathit{R}}_{\mathbf{0}},{\mathit{R}}_{\mathit{L}},L\right)={\int }_{\mathit{r}\left(0\right)={\mathit{R}}_{\mathbf{0}}}^{\mathit{r}\left(L\right)={\mathit{R}}_{\mathit{L}}}\mathcal{D}\left[\mathit{r}\left(s\right)\right]\text{exp}\left[-{\int }_0^L\mathit{ds}{⟨\xi ⟩}_L{\mathit{r}}_{\mathit{L}}^{\mathbf{2}}\left(\mathit{s}\right)\right]$$E4

${\mathit{r}}_{\mathit{L}}\left(\mathit{s}\right)$ is the parameterization of polygonal $\left\{∆{\mathit{R}}_i;\mathrm{i}=1,2,3\dots \right\}$ and L is the continuous analog of the number of cardiac cycles.

On the basis of such representations, second moments of the PDF are analytically calculated, which allows to determine the average time to reach the HE boundaries, represented by formulas (6)–(8). Presented cases of analytical expressions take into account not all possible scenarios for the evolution of states, in particular, the fact of the appearance of obstacles to the realization of certain types of trajectories, as a result of degradation of the cardio tissue, is not taken into account.

To solve these problems, we introduce a one-step cascade operator, defined as follows:

$({\mathit{r}}_{\mathbf{1}}\mathit{,}\dots {\mathit{r}}_{\mathit{m}}\mathit{,}\dots {\mathit{r}}_{{\mathit{N}}^{\mathit{\ast }\mathit{,}}}{\mathit{r}}_{{\mathit{N}}^{\mathit{\ast }\mathit{,}}\mathit{+}\mathbf{1}}$,…) - cascade

$${\mathrm{\Omega }}^{\ast }\left(..m_i..\right)=m_i+1$$E5

$${\left({\mathrm{\Omega }}^{\ast }\right)}^{N^{\ast }}{\mathit{R}}_{\mathit{k}}={\mathit{R}}_{\mathit{k}+\mathbf{1}},\forall \mathit{k}$$

The introduced operator has an analog in the symbolic space on the basis of which the recognizing automata, which is presented in the next section, are constructed.

Within the framework of the presented model of the segment wander from the wavelet coefficients of the ECG signal and in the performance of the Markov property, the prediction and estimates of the time to reach the boundaries of the regions of dysfunctions or cardiac events of AF-HF type are obtained by inferring the probability density from the Feyman representation. In some cases, the Fokker-Planck equations are obtained. By solving these equations, or by calculating the end-to-end distribution function by the methods described in the paper [8, 9], the following estimates of the time to reach the boundaries (L_time) events are obtained:

1. Models of the free random walk.

$$L_{\mathit{\text{time}}}=\frac{<{‖{\mathit{R}}^{\mathit{C}}‖}^2{>}^2}{\left(N^{\ast }-1\right)⟨\xi ⟩},$$E6

1. Models of random walk with constraints. Under condition of confirmation of hypothesis of model of random walk with constraints, the estimation of L_time is determined as a solution of equation

$$<{‖{\mathit{R}}^{\mathit{c}}‖}^2>=2\left\{\Theta L_{\mathit{\text{time}}}-{\Theta }^2\left[1-e^{\raisebox{1ex}{$-L_{\mathit{\text{time}}}$}\!\left/ \!\raisebox{-1ex}{$\Theta $}\right.}\right]\right\},\Theta =\frac{\mathit{\text{const}}}{⟨\xi ⟩\left(N^{\ast }-1\right)}$$E7

1. Models of random walk in a non-simply-connected domain. Under condition of confirmation of a hypothesis of model of random walk in a non-simply-connected domain, estimates of L_time is determined as a solution of equation

$$<{‖{\mathit{R}}^{\mathit{c}}‖}^2>={⟨\xi ⟩}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$\left(N^{\ast }+2\right)$}\right.}{\left(\frac{N^{\ast }+2}{3}\sqrt{\frac{2⟨\xi ⟩}{N^{\ast }}}{L}_{\mathit{\text{time}}}\right)}^{\raisebox{1ex}{$6$}\!\left/ \!\raisebox{-1ex}{$N^{\ast }+2$}\right.}$$E8

where ${\mathit{R}}^{\mathit{C}}\mathbf{\in }\left\{\mathit{\text{border of}}\mathit{HE}\right.,\mathit{HF}\}$, $<{‖{\mathit{R}}^{\mathit{c}}‖}^2>$ 2-th moment probability function of end-to-end distribution; $⟨\xi ⟩$ is the average length of an elementary transition.

Similar to the previous ones, it is advisable to construct the space of events and the symbolic space for a continuous analog of the wavelet transform. The transition to the continual version requires much computing power, but it has some advantages. Figure 2 shows the value of the wavelet coefficients of the continuous wavelet transformation. Only in this case is the development successful of the multiscale fluctuations of the QRS complex by the scaling variable in a patient with persistent AF. The prognosis of the development of such fluctuations and the time to reach the value by the wavelet coefficients of fields of strong irregularity according to Hölder are also modeled in terms of a random walk on a multidimensional lattice or in a continuum.

Figure 2.

Continuous wavelet transformation of ECG with AF and ethalon.

Let us return to the one-dimensional processes of wavelet coefficients with a fixed quadruple of indices. The PDF of these processes determines the probability of an elementary jump under the action of the shift operator by one step, that is, eventually the probability of a transition for L steps, as in the case of the vector process considered earlier. The type of distribution density thus affects the construction of the trajectory of the vector process. Here, it is necessary to note the following: in those cases when the one-dimensional process under consideration is a Markov diffusion process, and on the part of the wavelet coefficients of the ECG signal it is. The distribution density function here is a stationary solution of the Fokker-Planck equation. Also assume, for example, that the stationary solution is expressed as an exponential function of some potential function [16].

$$P_s\left(x\right)=\mathit{\text{const}}\bullet \mathit{exp}\left(U\left(x\right)\right)$$E9

where x is some wavelet cascade in a continuum limit.

In this case, it is possible to calculate possible rapid transformations of the PDF of one type, another type using elementary constructions from the theory of catastrophes (Poston Stewart) For this, the potential function or the neighborhood of its maximum is approximated by a polynomial of the fourth order. As a result of standard transformations (Poston), the so-called bifurcation set on the plane of the approximating polynomial coefficients is constructed, as shown in Figure 3.

Figure 3.

The sets of bifurcations or catastrophes and the changing of type of the PDF.

4. Recognizing automata

For each cascade $\left\{{}_{\mathit{\text{Hist}}}^k,W_{i,j}^N\right\},N=1,2,3,\dots ,N$*, an additional symbolic affine space $S^{n+1}$ is defined, the dimension of which is determined as n + 1 number of columns of the cascade histogram, that is, each histogram of the wavelet coefficient values with a fixed triplet of indices (Hist, i, j) is a point in the symbolic space $S^n$. This point determines the internal state of the automata. An elementary shift operator, shifting the segment $\left\{{}_{\mathit{\text{Hist}}}^k,W_{i,j}^N\right\},N=1,2,3,\dots ,\infty$ by one in the direction of index $N$* growth, is also defined at preservation of a segment of length $N$*. In this case, the internal state of the automata changes, the point in space shifts. The state of the automata changes with an elementary shift by the transition of elementary value from one coordinate to another, that is, in the case of a shift, the transition of elementary value from one column of the histogram to another column is carried out. By construction, the elementary transition is equivalent to the shift operator of some initial segment in segment $\left\{{}_{\mathit{\text{Hist}}}^k,W_{i,j}^N\right\},N=1,2,3,\dots ,\infty$ by a unit step. In a symbolic space, the elementary shift corresponds to the multiplication of affine matrix of the form

$${\mathrm{\Omega }}_{i,k}\left(\begin{array}{c}{n}_1\\ n_2\\ \begin{array}{c}.\\ n_i\\ \begin{array}{c}{n}_k\\ 1\end{array}\end{array}\end{array}\right)\stackrel{\scriptscriptstyle\mathrm{def}}{=}\left(\begin{array}{ccc}\begin{array}{cc}1& 0\\ 0& 1\end{array}& \begin{array}{cc}0& 0\\ 0& 0\end{array}& \begin{array}{cc}0& 0\\ 0& 0\end{array}\\ \begin{array}{cc}.& .\\ 0& 0\end{array}& \begin{array}{cc}.& .\\ 1& 0\end{array}& \begin{array}{cc}.& .\\ 0& 1\end{array}\\ \begin{array}{cc}0& 0\\ 0& 0\end{array}& \begin{array}{cc}0& 0\\ 0& 0\end{array}& \begin{array}{cc}1& -1\\ 0& 1\end{array}\end{array}\right)\left(\begin{array}{c}{n}_1\\ n_2\\ \begin{array}{c}.\\ n_i\\ \begin{array}{c}{n}_k\\ 1\end{array}\end{array}\end{array}\right)=\left(\begin{array}{c}{n}_1\\ n_2\\ \begin{array}{c}.\\ n_i+1\\ \begin{array}{c}{n}_k-1\\ 1\end{array}\end{array}\end{array}\right)$$E10

on the set of coordinates of the state point. The consistent application of the shift operator of $N^{\ast }$ times ensures a transition from one segment of length $N^{\ast }$ to the next segment of the length N* in the state space. In a symbolic space, such a transition for $N^{\ast }$ steps will correspond to a sequence of matrices of the form.

$${\Omega }_{N^{\ast }}\stackrel{\scriptscriptstyle\mathrm{def}}{=}\left(\prod _1^{N\ast }{\Omega }_{i,k}\right) \mathit{for} \forall \left(i,k\right)$$E11

The $N^{\ast }$-fold multiplication of the matrices ${\Omega }_{N^{\ast }}$ determines the matrix of frequency for each elementary transition and, accordingly, the transition probabilities matrix ${\widehat{w}}_{m,k}$.

Assuming that the transition probabilities obey the stationary properties, the prognosis is carried out as follows:

First, in a symbolic space, the problem of the wandering of a point for $N^{\ast }$ elementary steps is formulated. Further constructions of model and prediction algorithms have many options, and this set is determined by the characteristics of the observed signal as sequences of indexed N and fixed triple of remaining indices. These properties include the Markov property, stationarity, and ergodicity of each cascade of wavelet coefficients.

The presence of these properties must be checked either by direct calculations based on manipulations with wavelet decompositions or by the calculation of the entropic, dimensional (information, capacitive, dimensionality, etc.) characteristics of the observed trajectories, their correlation radius. In addition, it is necessary to analyze the conditions of point-wise Hölder regularity of the observed and predicted trajectory on the basis of two microlocalizations.

Thus, the vector $\left(n_1,\dots n_k\right)$ is defined as the internal state of the automaton, and this vector corresponds to the frequency histogram of the vector $\left\{{}_{\mathit{\text{Hist}}}{}^k{W}_{i,j}^N;N=1,2,3,\dots ,N^{\ast }\right\}\subset R^{N^{\ast }}$, Matrix operator ${\mathit{\Omega }}_{{\mathit{N}}^{\mathit{\ast }}}$ determines the transition frequencies. If to move from frequency representations to probabilistic standard renormalization $\left(n_i\to p_i\right)$, then the evolution equation for $p_i$ looks as follows in the continuous representation in time:

$$\frac{\partial p_m}{\partial t}=\sum _k\left[{\widehat{w}}_{k,m}{p}_k-{\widehat{w}}_{m,k}{p}_m\right]$$E12

where ${\widehat{w}}_{k,m}$ is the transition probabilities defined by the operator ${\mathit{\Omega }}_{{\mathit{N}}^{\mathit{\ast }}}$. Eq. (12) is the so-called kinetic equation, the master equation or, by its nature, the balance Equation [18]. The change in the internal state of the automaton is described by the solution of Eq. (12). On the other hand, the change in the internal state of automaton is defined as the walk of a point along $N^{\ast }-1$ dimensional simplex ${\Sigma }^{N^{\ast }-1}$, defined by condition $\sum _i{n}_i=N^{\ast }$ and transitions ${\Omega }_{i,k}$ are defined on the one-dimensional faces of the simplex.

The solution of Eq. (12) gives sufficient information for solving the problem of random walk in the space $\mathit{W}\stackrel{\scriptscriptstyle\mathrm{def}}{=}\left\{{}_{\mathit{\text{Hist}}}^k,W_{i,j}^N\right\}$ in the representation of the probabilities of the transition from one vector of space to another in the form of a Feynman path integral. In many cases, to implement the tasks of the prognosis и and time estimates of reaching the boundary of the region of cardiac events, it is only sufficient to solve Eq. (12), that is, the solution of the problem of a walk on the simplex ${\Sigma }^{N^{\ast }-1}$. This is possible in cases where the early predictors are expressed in the singularities of walking on the simplex ${\Sigma }^{N^{\ast }-1}$. As in the space of trajectories, such predictors distinguish regions of the simplex for such predictors.

However, in the space formed by the cascade wavelet coefficients of the ECG, from the construction itself follows that in the normal operation of the heart in the transition from one state vector to another, the solution of equation must be stationary, if not for all cascades, then for some subset of them. Taking into account biological rhythms and other cyclic processes in the body from experimental tests, it follows that in fact the solutions of Eq. (12) satisfy the stationarity condition (13) in the “average”

$${\widehat{w}}_{k,m}=\sum _m{\widehat{w}}_{m,k}$$E13

That is, fluctuations of the stationary solution whose amplitude is determined from the chronological database are allowed. Thus, the loss of stability of the stationary solution with subsequent transitions to another stationary solution changes the nature of a walk in space $\mathit{W}$ and leads the trajectory to the boundary of the region of cardiac events. In this case, the nature of the fluctuations and their amplitude varies, and the analogy of such changes is the dynamics of the fluctuations during phase transitions. In addition, this is only one of the prognosis scenarios. Another scenario corresponds to the violation of the stationarity conditions and the time-dependent solutions of Eq. (12). This scenario will be discussed in the next section.

In conclusion, a few words about early predictors are discussed. In order for this class of predictors not to be empty, it suffices to point out violations of stationarity conditions and conditions of point-wise Holder regularity, in particular, on the change of Holder exponent of the ECG signal at the time $t_0$. The conditions for such events, the so-called two-microlocal conditions, are also defined on the wavelet coefficients [17] and, consequently, are present in space in the form of regions determined by various kinds of inequalities that limit the set of admissible values of wavelet coefficients.

Let us give some examples. When the stationarity condition is satisfied, the prediction problem reduces to the problem of the wandering of point on the multidimensional lattice $Z^n$. Here, the method of constructing a symbolic space allows to correctly move to the continuous model and uses the technique Feynman path integral. In this case, the probability of reaching a pre-assigned point of symbolic or basic state space is estimated. The time estimates in terms of the number of elementary steps to achieve a pre-determined state are given subsequently. These states also include the ranges of values of wavelet coefficients that do not satisfy the conditions of homogeneous or point-wise Hölder regularity. As can be seen from the formulas, the estimates vary depending on the characteristics of the observed ECG signal or its wavelet coefficients. For example, if the effect of excluded volume is taken into account, when forbidden trajectories appear in a set of admissible trajectories. Such a phenomenon is possible when the fragments of the cardiac myocyte sequences are turned off as a result of local myocardial degradation, when groups of conducting cardiac myocytes no longer can fully or partially perform their conductive functions. The effect of excluded volume significantly changes the properties of the considered processes.

For example, the Chapman-Kolmogorov equation becomes unjust, the system is no longer a Markovian system, and so on. Such phenomena on the one hand are themselves predictors; on the other hand, for correct estimates, it is necessary to introduce the transition probabilities and PDF containing more variables.

In this situation, the change in the properties of a myocardial tissue is associated with the discrepancy between the observed and predicted ECG parameters of the wavelet coefficients, which is a signal for automatic complication of the model by the birth of new automata. In concrete example, the appearance of forbidden trajectories (the effect of excluded volume) is a command to construct three-dimensional histograms and multidimensional transition functions. However, the further principles of the operation of automata remain unchanged. The same multiplication of automata can occur at searching for hidden predictors. For example, in estimating the change in the properties of regularity, smoothness, and so on, trajectories as a criterion are often inequalities that contain sums over the time index j.

Thus, when the stationarity conditions are fulfilled, the prognosis is reduced to estimates to reach the critical regions by a trajectory or a class of trajectories. Depending on the characteristics of the process, the estimate may vary, and a correction of the prognosis is necessary. It is this fact that determines the monitoring regimes, their frequency, and duration. In this case, the degree of deviation of the observed trajectory from the predicted trajectory, or its characteristics: moments, properties, conditions of Holder, and so on, is also estimated.

If the quasi-stationary conditions are violated at the first step, taking into account the trends of transition probabilities, a new PDF is determined on the basis of the kinetic equation; then, taking into account the changes in the PDF and the trends of transition probabilities, trajectories or their new characteristics are recalculated. However, under the conditions of fulfilling the quasi-stationary conditions, a change in the structure of the transition probabilities, an increase in the amplitude of their fluctuations, and a change in the structure of the set of transition functions are possible. In this case, the algorithms calculate the change in the positions of the values in the space of the PDF approximation parameters, their proximity to the bifurcation set. Taking into account the approach speed of the approximation parameters with the bifurcation set, a further prediction in the state space is corrected. Trajectories in this case means a set of trajectories in a set of spaces, and the earliest signs can appear for the coefficients of only one class of trajectories with a fixed triplet of indices. Subsequently, predictors can appear on the remaining classes of trajectories with other fixed indices.

5. Interaction with devices

The interaction of the remote monitoring system with implantable and portable devices is built according to the following scheme:

1. Initially, there is a set of statistical data in a chronological database. In view of the peculiarities of the set of statistics, it is necessary that the samples for all fixed indices be representative both at rest and in the state of motion.

2. At such set, the automata create a database of valid trajectories in small dimensions.

3. With the full set of chronological database, automata begin calculations according to the algorithms described earlier.

The further mode of preventive monitoring is determined by automata and is based on the requirement of sufficiency of statistics. After that, the automata determine the further monitoring strategy. During the operation of automata, numerous additional hidden predictors are identified, the earliest. If the development of hidden predictors leads to the emergence and development of existing predictors, a minimal subclass is allocated as signal and management automata from the whole class of automata. Their volume should correspond to small computing resources of implantable or portable devices. Further, this minimal subset is transferred to the device’s memory and further serves as a signal device that controls the calculations on the remote server. The described scenario allows to optimize costs and the preventive monitoring modes.

6. Verification

The experimental verification of the capabilities of the set of recognizing automata presented in this chapter is carried out continuously during the last 4 years [19, 20, 21]. The ECG signal is selected as an initial observed signal for the analysis and prediction of cardiac events. The standard scheme for measuring the ECG signal by a recorder in 12 leads with a variable sampling rate of 1–2 kHz and a 24-bit resolution is considered.

The wearable ECG set has DSP compute block on-board to partially offload cloud infrastructure and to monitor cardio events in real time. The particular set of automata computed locally depends on power demand-reaction time tradeoff, which in turn depends on particular patient’s case. In any case, all the collected data are compressed and transferred to the cloud. Currently, Wi-Fi is used for communication. Mobile connectivity in spaces where Wi-Fi infrastructure is not present is achieved via mobile Wi-Fi tethering with smartphone (mostly to get rid of multiple sim cards burden). Also, the device carries Bluetooth LE which is used for settings transfer and standard on-site real-time monitoring via tablet software or PC software if one carries BLE receiver. A schematic diagram of the cluster work is presented in Figure 4.

Figure 4.

A schematic diagram of work of remote preventive cardiac-monitoring cluster: ECG device and on-board recognizing automata for the realization of signal function and management by recognizing automata in the cloud (I); database, chronological database of patient, database of HE (II); smartphone for text-graphic messages of cluster (III); 1—space W; 2—the set of interacting automata parallel to the processing of W-cascades and defining ${\widehat{w}}_{k,m}$; 3—two-microlocal analysis; 4—check of the quasi-stationarity; 5—prognosis in the symbolic space ${\Sigma }^{N^{\ast }-1}$ using automata from block 2; 6—check of stability of the prognosis; 7—constructing the prognosis in the state space W.

The check consists of two stages. At the first stage, a chronological ECG database for patients is used to predict cardiac events. In doing so, the ECG is used both for cardiac events and for the ECG of the control group. Chronological databases of patients with a long history from the occurrence of cardiac events, and their subsequent treatment of drug or catheter ablation were also used. Monitoring was carried out both before ablation and after it during a long period of drug support. At the first stage of checking, the effectiveness of recognizing automata and the ability of automata to predict cardiac events were tested. In those cases when the automata did not predict a cardiac event, the automata returned to the beginning of the recording, they became more complicated, and the process was repeated. It should be noted that the main purpose of the described experiments and the basic principles laid down in the algorithms of automata are aimed at preventive monitoring, that is, on the detection of the earliest predictors of cardiac events with subsequent time estimate of the evolution of these predictors until the appearance of later predictors, already known, such as dispersion QT interval, P wave index, increased QT interval duration, and so on. The main problem, because of which there were gaps of cardiac events by automata, is as follows. If to analyze the results of each automaton individually, then during the evolution of their states, there was no approximation of the state trajectory to the boundaries of the regions of cardiac events in any of the selected metrics, but the event was happening. The analysis showed that the predictor of the event in these situations is not a violation of the quasi-stationary conditions, but a change or a mismatch in the structure of transition probabilities for some subset of automata. The revealed mismatch can be characterized in terms of conditional entropy. Ultimately, conditionally entropic characteristics were used in the approximation of conditional K-complexity, since the calculation of K-complexity is an algorithmically unsolvable problem. Another way to solve the problem is to reduce the complexity of automata in dimension, which has an analogy in the transition from single-particle PDF to multiparticle PDF.

Figure 5 shows changes in the states of automata under the action of the operator ${\mathit{\Omega }}_{{\mathit{N}}^{\mathit{\ast }}}=\left(\prod _{\mathbf{1}}^{\mathit{N}\mathit{\ast }}{\mathit{\Omega }}_{\mathit{i}\mathit{,}\mathit{k}}\right)$ considered in Section 3. The figure reflects the mixing nature of the actions of the operator ${\mathit{\Omega }}_{{\mathit{N}}^{\mathit{\ast }}}$ and the temporal evolution of transition functions, entropic, informational, and dimensional characteristics of which determine the earliest predictors.

Figure 5.

Changes in the states of automata.

By changing the consistency of the state trajectories of a certain subset of automata, it is implied in this case that early predictors are defined as differences in the structure of operators ${\mathit{\Omega }}_{{\mathit{N}}^{\mathit{\ast }}}$. However, here we are talking about the earliest predictors or automata, distinguishing trajectories, leading to cardiac events from trajectories without cardiac events.

If to consider traditional predictors associated with the estimation of the duration of intervals, the characteristics of the QRS complex of the alteration, P- wave, dispersion of P-wave index, then in this case, the automata predict rather successfully the evolution of the listed predictors.

This is important for drug treatment of persistent AF and optimization of drug therapy, as for many other cardiac applications.

7. Task of management of the trajectories

We now return to discuss the management problem mentioned in Section “Introduction.” The set of admissible trajectories is sufficiently variable. Depending on many factors, these changes are associated with cardiovascular degradation, changes in conductivity at cell scales, changes in the architecture of a set of conductive paths between sets of conducting cells, and so on. The very set of trajectories is so factorized into equivalence classes, regarding the action of groups of process symmetries, the type of the PDF process, and the set of transition probabilities. In view of the factors mentioned earlier, there are prohibitions on transitions from one state P to another. Thus, some trajectories in classes become forbidden when all these factors are taken into account, or the probability of some trajectories becomes small.

In fact, the task of controlling trajectories reduces to changing the class of trajectory or to the task of keeping a trajectory in a given class by means of variable management parameters. The management parameters include all parameters on the macro- and microlevel, which can be varied in various ways. Such methods include drug therapy with AF-AT events, and AV and VV programming of the CRT device [22, 23]. The same goals are pursued with ablation or defibrillation. Within the framework of the presented model, all possible ways to change the class of trajectories to the class of trajectories that do not terminate HE are formalized as shown in Figure 6, where the management loop is mapped into a state space or trajectories.

Figure 6.

Management loop (yellow color) is mapped into a state space of trajectories (blue color).

Passing to formal language of homotopy theory and theory of infinite loop spaces [24], the process of management defined by the mapping

$$S:\partial I^{k+1}\to \partial C^k$$E14

k + 1 is the dimension cube, k is the number of management parameters, k + 1 is the parameter—time.

At each fixed time t and with the variations of other к parameters, a mapping the boundary of k—dimension cube to the area homotopically equivalent k-dimensional sphere in K is defined

$$S^{-}:{\partial I}^k\to \partial C^k\subset K$$E15

Definition 1.

The management task is solvable if and only if.

1. The set of homotopy classes $[{\partial I}^k,$ $\partial C^k]$ is trivial or.

2. The mapping $S^{-}$ belongs to the trivial element of the set of homotopy classes.

In these examples, the set of homotopy classes has a group structure and is defined as the k-homotopy group of the k-dimensional sphere $\left[{\partial I}^k,\partial C^k\right]={\pi }_k\left(\partial C^k\right)$.

In other words, the management problem is solvable if, with the help of a variation of management parameters, the observed trajectory can be deformed into a predefined trajectory if and only if there is no forbidden trajectory or other topological obstacles between them. A little more optimism is given by the following statement: if the management problem is not solvable with this set of management parameters, then changing the number of management-led parameters can possibly translate the problem into a class of solvable management problems. In the version of recognizing automata, the management problem and the calculation of the management strategy are reduced to determining the influence of the management parameters on the set of admissible transition probabilities in the predefined class of trajectories.

8. Conclusion

Several years of experimenting with recognizing automata and the preliminary obtained results allow to make optimistic conclusions about preventive medicine, in particular, preventive monitoring. The development of mHealth platforms, portable and implantable devices, on-body sensors, and so on, available wireless data transmission systems and available computing power allow solving very complex prognosis in a very short time and in a number of cases in real time. Thus, the development of preventive monitoring systems creates a trend toward changing paradigms, at least in the field of cardiology. If early predictors exist, then a natural question arises, but is there any possibility of influencing the character and speed of development of early, preventive predictors via preventive drug therapy, different diets, and so on in the direction of reversibility of the current situation. That is, does the class of trajectories exist when early appeared predictors are eliminated by preventive maintenance, otherwise when the situation is physically reversible? In the PHM/IMS applications to technical objects, the term “self-maintenance, self-recovery” appears. Medicine is more conservative, and yet the above analogy is appropriate. In the context of preventive cardiomonitoring, this means that the number of age groups that preventive control is recommended increases markedly and begins on average from 30 to 40 years. This is also indicated by the statistics of the growth of heart diseases, which currently has the nature of a pandemic, as well as statistics on the rejuvenation of heart diseases [25]. The concrete ways of creating preventive monitoring systems are now quite realistic and are reduced to the realization of the fact that the chronological basis of the individual’s ECG data is needed to identify early predictors.

The database is updated periodically. The question of the refresh rate is solved by the system of recognizing automata.