A computer-implemented method for determining the electric propagation at the heart of a subject by external measurements

Abstract

The present invention relates to a computer-implemented method for determining the electrical propagation of a subject's heart by external measurements, wherein the external measurements are taken from a subject's torso by locating a plurality of electrodes that provide an electrical signal responsive to heart activity.

Description

[0001] A COMPUTER-IMPLEMENTED METHOD FOR DETERMINING THE ELECTRIC PROPAGATION AT THE HEART OF A SUBJECT BY EXTERNAL MEASUREMENTS DESCRIPTION

[0002] FIELD OF THE INVENTION

[0003] The present invention relates to a computer-implemented method for determining the electrical propagation of a subject's heart by external measurements, wherein the external measurements are taken from a subject's torso by locating a plurality of electrodes that provide the electrical signal providing information responsive to heart activity.

[0004] BACKGROUND

[0005] One of the most intense activities in medicine is to act efficiently in the application of therapies appropriate to the type of disease. In view of the great challenges, there are scientific fields where it is not possible to adequately establish the biological, electrical and mechanical behavior of organs as complex as the heart, since a direct access to the functioning organ is invasive and any surgical action can modify the operating conditions. Therefore, it is not possible to make direct measurements of potential or electrical activity at least across the tissues.

[0006] Importantly, current mapping technologies consider only the outer surface of the cardiac tissue: whereas invasive mapping allows to characterize the electrical activity on the endocardium surface, non-invasive mapping techniques allow to assess cardiac electrical activity on the epicardium surface. However, the cardiac tissue is three-dimensional, and counts with a midmyocardial substrate that can, in some cases, have a thickness of several millimeters. As an example, the importance of characterizing midmyocardial activity is paramount for an accurate evaluation of some cardiac diseases, such as ventricular arrhythmias related with transmural scars, in which the measurement of cardiac tissue surface electrical activity is not representative of the mechanisms sustaining the arrhythmia.Throughout this description, the term "source" will be used to refer to those regions of cardiac tissue that are capable of generating an electrical charge that allows electrical stimulation. These electrical charges propagate along a predetermined pathway along the cardiac tissue that trigger coordinated contraction of different parts of the heart.

[0007] Each person is characterized by a different distribution of sources. The characterization of different functional regions according to electrophysiological criteria is very important in the cardiological field.

[0008] The present invention overcomes the described drawbacks and makes it possible to carry out indirect measurements of the dynamic electrical activity and its propagation through the heart, using a plurality of measurements of the electrical potential by means of electrodes placed on the surface of the torso and adapted to detect the potential at such points.

[0009] DESCRIPTION OF THE INVENTION

[0010] The present invention solves the problem of how to determine electrical activity in the heart of a subject by indirect measurements that do not require invasive methods. This indirect determination makes it possible to determine the electrical pathways within the heart and, according to various application examples, also to distinguish between tissues with normal / healthy and abnormal / diseased electrophysiological behavior.

[0011] The impossibility of directly accessing each point of the tissue to measure the electrical potential is even greater when it is necessary to obtain measurements at points on the tissue that do not correspond to the surface. A surface measurement could be acquired with an electrode inserted surgically and / or implanted on the surface of an organ such as the heart, but measurement at a point inside the volume occupied by the organ is not possible since an intrusive measurement modifies the propagation conditions and the measurement is not representative of what occurs without the existence of the electrode. The method therefore allows at least non-intrusive measurements to be made at points in the volume delimiting the torso.

[0012] The solution to the problem requires a biological model in which the tissues and media of the subject's body are hypothetically considered as continuous media. In particular,when talking about electrical potential, this will refer to the extracellular potential. The considered biological model, is based on a coupling between intracellular and extracellular domains such that the potential difference across the cell membrane is considered the action potential (also called transmembrane voltage).

[0013] Another fact considered in the modelling of the problem is that there is no current flowing outside the surface bounding the torso of the subject.

[0014] When we talk about the torso of the individual, we mean a part of the body that includes the heart or at least a part of the heart.

[0015] Based on these assumptions, the method is a computer-implemented method for determining the electrical activity at the heart of a subject by external measurements, is carried out in a computer system and comprises the following steps:

[0016] defining a numerical model of the torso of a subject, the numerical model at least comprising a domain (Ω) representing a volume of the torso bounded by the torso surface (∂Ω) and by one or more closure surfaces, the domain being divided into a first sub-domain (Ω₁) representing at least one part of the heart and a second subdomain (Ω₂) representing the remaining volume, the numerical model at least comprising the electrical properties of the medium, preferably the conductivity (σ(x)), being x a point of the domain (Ω); obtaining as input data body surface potentials gᵢ measured on a plurality of N preestablished locations (yᵢ; i = 1,..., N) of the torso surface (∂Ω) of the subject, at

[0017]

[0018] least at nₜ ≥ 1 instant times.

[0019] The Ω domain is a volume bounded by a surface usually designated by ∂Ω. The surface bounding a volume will also be identified by the letter S as an alternative to d£l to simplify notation. The same applies to the volume, for instance, integrals extended to a volume V may use the notation dV or dΩ for the differential part.

[0020] The numerical model of the subject is a numerical model that at least represents the shape of the heart and at least part of the surrounding organs and / or tissues of the heart. The numerical model may be restricted to a part of the subject so that whenever reference is made to "torso", it includes a surrounding volume and need not include the head, arms or part of the body below the waist; however, the method isstill valid if these parts are also included in the numerical model with the only penalty being that the computational cost will be higher.

[0021] Any acquisition method allowing to estimate the shape of the torso of a subject, the position of the electrodes, and an estimation of the shape of the at least one part of the heart is valid for defining the numerical model.

[0022] The shape taken by the numerical model, in particular the first sub-domain (Ω₁) and the second sub-domain (Ω₂), according to a first example corresponds to a generic model having the usual shape and sizes. These usual shapes can be scaled to acquire dimensions of the subject.

[0023] According to another more preferred example, both the shape and the size of the first sub-domain (Ω₁) and the second sub-domain (Ω₂) are acquired by means of a 3D imaging process, for example by means of magnetic resonance imaging acquisition. In preferred cases, the shape of the domain and the electrical properties associated with each point in the domain are values acquired or stored in a database previously populated with subject data or with model data generated by other numerical models.

[0024] A specific way to identify different organs, such as for example the heart, is by means of image segmentation techniques which allows to distinguish the heart from other organs separated by fluids, fat or different tissues defining interfaces between limiting surfaces of the organs.

[0025] The segmentation process, according to examples, are preferably by means of automated, semiautomated, or manual techniques. Among automated and semiautomated methods are methods based on artificial intelligence where, after a learning phase, a segmentation module identifies each of the parts and organs establishing distinct segments or regions.

[0026] In a preferred embodiment, the first sub-domain (Ω₁) representing at least one part of the heart is a part or the whole heart volume, not only the surface of the heart. This option allows the method to provide at least indirect measurements at points on the heart that correspond to non-surface tissue points.According to an embodiment, the domain (Ω) is a closed domain wherein the torso surface (∂Ω) is within the domain (Ω).

[0027] Mathematical expressions will be used throughout the description. Mathematical expressions and equations are ways of linking variables. The same link can be expressed in different ways. For example, a circle is the geometric locus of points that verify that they have the same distance to a point, the center. This locus can be represented by the equation of a circle or also by expressing each coordinate by a parametric equation where both expressions define the same locus. For this reason, once a relation between variables or conditions imposed on such variables has been defined, these same relations or conditions, although expressed in a different way, will be considered equivalent.

[0028] The input data are N measurements of the torso surface potential. According to the preferred embodiment the measurements are acquired by means of electrodes placed on the skin surface of the subject at different preestablished locations yᵢ; i = 1,..., N. These locations yᵢ are located at the surface of the patient since the electrodes are just touching the skin of the patient registering the potential at said surface / skin. The preestablished locations where the subject's electrodes have been placed are positions that must also be located in the numerical model at the corresponding position in order to carry out the next step.

[0029] That is, the electric potential is directly measured at locations yᵢ at the torso, 1 ≤ i ≤ N, where N corresponds to the number of electrodes of the system. However, further indirect measurements can be obtained numerically from direct measurements.

[0030] The electric potential, according to a specific embodiment, can be interpolated, projected or estimated with data-based or Al-based methods to the hole torso surface. If the torso surface is discretized, for instance by means of a triangulation, some nodes are at the locations with measurements are direct measurements obtained by means of the electrodes and, other nodes has indirect measured values estimated from direct measurements as previously indicated. In any case, when referring to electric potential measurements of the torso surface will refer to any of the previous measurements, direct measurements with electrodes or indirect measurements that would involve an intermediate computational step departing from the direct measurement from the electrodes.determining functions G_{yᵢ}(x), i = 1,..., N solving numerically, for each G_{yᵢ}(x), the set of differential equations:

[0031] −∇·(σ(x) ∇G_{yᵢ}(x)) = δ_{yᵢ}(x), x ∈ ΩaG_{yᵢ}(x) + b∂ₙG_{yᵢ}(x) = c(x), x ∈ ∂Ωwherein G_{yᵢ}(x) is a function defined in the domain (Ω), wherein a, b ∈ R are two preestablished values, σ(x) is the conductivity function and c(x) is a preestablished function over the torso surface ∂Ω, wherein ∇· is the divergence operator, ∇ the gradient, δ_{yᵢ}(x) is the Dirac's delta distribution with center on the location yᵢ, and ∂ₙ denotes the outward normal derivative with respect to the torso surface ∂Ω; defining a first operator A defined such that, when A acts on a function f(x), represented by Af(x), f(x) being a scalar function defined on domain (Q) and representing the cardiac source, the first operator A at least determines a first contribution given by a numerical computation of ∫_Ω G_{yᵢ}(x)f(x)dV being f(x) = 0 in the second sub-domain (Ω₂) and, the integral is a volume integral over the domain (Ω); defining a second operator R defined such that, when R acts on a function (p(x), represented by R(p(x), (p(x) a scalar function defined on domain (Q) and representing the electric potential, the second operator R at least determines a second contribution given by a numerical computation of φ(x) + ∫_{∂Ω} ∂ₙ G_{yᵢ}(x) σ(x)φ(x)dS; the integral being a surface integral restricted to the torso surface (∂Ω), wherein φ(x) satisfies φ(yᵢ) = gᵢ, i = 1,..., N; obtaining the estimation of f(x) and φ(x) of the subject as the result of minimizing the contribution of at least a first quantity Af(x) and a second quantity −Rφ(x), measured under a predetermined norm ||·||, namely ||Af(x) − Rφ(x)||, satisfying that ∫_{Ω₁} f(x)dV = 0; determining the electric activity of the heart of a subject by at least one value of f

[0032]

[0033] (x), (p(x) or both at least at one instant time.Functions G_{y₁}(x) are the result of solving a system of differential equations containing Dirac's delta distributions centered at the positions of the model corresponding to positions of the subject where the electrodes are located and, allow to establish a plurality of known functions that enter as arguments of the operators A and R as defined.

[0034] The resolution of the system of equations that allows the calculation of the functions G_{yᵢ}(x) requires a discretization of the numerical model in the whole domain Ω. The discretization of the domain depends on the numerical technique used to solve the system of differential equations.

[0035] Among the methods of solving partial differential equations are: the finite element method, the finite volume method, the finite difference method, spectral methods, etc. The preferred method, which will be described below, is the finite element method, since it allows to represent complex shapes in a simple way, such as the shapes adopted by the different organs of the body and in particular the heart.

[0036] Once the functions G_{yᵢ}(x) have been calculated, they are considered known in the subsequent calculations involving the iterative method of searching for the minimum of an objective function.

[0037] Let a, b ∈ R be two preestablished values and c(x) is a preestablished function over the torso surface ∂Ω. As it will be disclosed, according to an embodiment, the chosen values are a = 0, b = 1 and c = −1 / Sr where Sr is the area of the surface of the boundary. Any other values may be preestablished since resulting functions G_{yᵢ}(x) will be different but, at the end of the method the solution for f(x) and φ(x) would not depend on such preestablished values.

[0038] The implementation of the numerical computation of the operators A and R requires a discretization of the domain Ω. This discretization does not need to be the same discretization as the one used when determining the G_{yᵢ}(x) functions. However, according to a preferred example the discretization is the same. This avoids the computational cost of determining a second discretization and also the cost of expressing the G_{yᵢ}(x) functions, already determined in the first discretization, in the second discretization. The latter transformation would require an interpolationprocess which is avoided by using the same discretization in the two steps.

[0039] In turn, the operator A has a first contribution by evaluating an integral over the whole volume where the function f(x) representing the sources, i.e. the charge causing an electrical potential, must be known.

[0040] On the other hand, the operator R has a second contribution that only takes into account points of the surface and, for its evaluation, it is necessary to know the value of the potential φ(x).

[0041] Although the potential φ(x) and the value of the sources f(x) must be known to carry out the calculation of the integrals that define one and the other operator, such operators are found within a quantity to be minimized by a numerical method, the quantity measured under the predetermined norm ||·||. That is, when the term "quantity" is used, it does not necessarily have to be a scalar, it can be mainly vectorial or even a matrix. The latter case is obtained when the norm is minimized over a set of time instants, the time instants corresponding to the time discretization. Of course, the norm applied to any quantity is a scalar value.

[0042] This implies that, according to a preferred embodiment, the numerical method is an iterative method in which initially a value of / °(x) is proposed as the first approximation of the solution, and it is the iterative process that establishes new convergent approximations to the solution at least for f(x).

[0043] Initial values of the iterative process can be, for example, the null value for the sources in the whole domain. The iterative process will apply successive corrections for one or more variables until the solution is reached after convergence.

[0044] An iterative method includes a stopping criterion. This stopping criterion may be based on a maximum number of iterations or a maximum computation time. Additionally, the convergence of the method is assessed when the norm of an error estimate or a residual value is less than a predetermined value. If both stopping criteria are imposed, then the iterative method stops when either of the stopping criteria is verified. In another embodiment, the iterative method stops when more than one of the stopping criteria are verified.The quantity to be minimized by a numerical method, in an alternative embodiment, is solved using a solver, in particular a direct solver for example based on matrix factorization.

[0045] During the calculation of the operator R, it is necessary to know the values of φ(x) at least in an estimated form. According to the preferred example, the values of φ(x) are estimated from the measurements at φ(yᵢ) which are extrapolated to those points of the discretization involved in the calculation of the operator R. A specific way of interpolation is to use the Laplacian operator to propagate the values at the positions yᵢ to the rest of positions where the value of φ(x) is needed.

[0046] The advantage of this example realization is that the estimated values of φ(x) obtained by interpolation do not change from one iteration to another, they are like if they were values coming from a measurement as well making the calculation more efficient.

[0047] According to another example of realization, the values of φ(x) in positions other than the positions yᵢ are part of the variables to be solved in the minimization problem so that these, in an iterative process, will be modified at each iteration converging towards the solution and, the operator R must be updated in each iteration. In this case the approximation of φ(x) is not by interpolation.

[0048] After convergence the method provides the electrical activity of the heart since the solution to the minimizing method is the function f(x) providing the sources at any location of the domain (Ω), the electric potential φ(x) also at any location of the domain (Ω), or both.

[0049] The values of f(x) and φ(x) are obtained directly at each node of the discretization although the discretization method also provides the approximate values in the rest of the domain. For example, a finite element discretization also defines the basis on which any value in the domain is expressed as a function of the nodal values in the element.

[0050] An efficient way to operate with the nodal values is to make use of vector and matrixtechniques as it allows an efficient implementation in the computation by a computer system. The resolution of an embodiment using finite element and vector techniques will be the example of choice when describing a detailed embodiment later.

[0051] Once the objective function to be minimized has been defined, it is possible to use available solvers, for instance to carry out an iterative process that converges to the solution or an alternative direct method for solving the minimization problem.

[0052] According to an embodiment of any of the previously disclosed embodiments, in the minimizing step the norm expression to be minimized further comprises the contribution of a third quantity λh(f(x)), being λ a real regularization parameter, and h(·) being a preestablished real regularization function.

[0053] The use of a third quantity contributing to the minimization step allows to adjust the solution when there is mismatch model, electrical noise or numerical errors to stabilize problems that would otherwise end up being degenerate problems.

[0054] Values used as the parameter A are in a predetermined range. In a preferred embodiment, parameter A is the regularization parameter used in invers problems controlling the equilibrium between the combination of the first and the second quantities and, the third quantity. In an embodiment parameter A is determined as the value of maximum curvature of the L-curve method. It has been found as adequate values of A in the range [10 “7, 10-1].

[0055] Examples of preestablishing real regularization function h(·) are choosing any of the classical functional norms such as the norm L₂, L₁, etc.

[0056] Another type of regularization consists in taking a set of linearly independent functions capable of describing the cardiac sources through a linear combination. Then, an optimization problem is formulated aimed at estimating the coefficients of the linear combination. The set of linearly independent functions, according to some examples, is determined by the geometric properties of the heart, the characteristics of the cardiac sources, learned a priori through data-driven algorithms or artificial intelligence-based methods or taken from a dictionary of multiple functions.According to an embodiment of the previously disclosed embodiments introducing the contribution of a third quantity, the real regularization function is a function responsive to:

[0057] a) the activation time at point x, b) the conductivity of the medium at point x, c) a weighted measure w(x) / (x) wherein w(x) is a real function representing a dysfunctional degree of the tissue at location x, the value ofw(x) ranging from 1 to a predetermined f> real value, f> > 1, d

[0058]

[0059] ) a combination of any of them and, in particular values under a norm.

[0060] In the most general case, the minimization process takes into account the functions / (%) and < >(%) over a period of time. That is, the two functions can be expressed as (x, t) and < >(%, t) where (%) and < >(%) are sampled at a given set of time instants. When the spatial discretization does not evolve with time, the nodal values for each time instant allow to assess in a simple way the value of an indirectly measured signal, the signal corresponding to the value of the source (%) at the nodal position and the signal corresponding to the value of the potential < >(%) also at the nodal position. The resulting signal is a sampled function at the given set of time instants where (%) and < >(%) are sampled. Throughout the description, to simplify the notation, when (%) and < >(%) are used, it will be understood that both functions can depend on time even though it is not explicitly stated in the dependence, unless otherwise stated.

[0061] The processing of these signals is the one that allows to establish for example the activation time at a point x.

[0062] The conductivity values are assessed e.g. by using databases where known values are stored depending on the tissue. For example, after a segmentation process on an image obtained by magnetic resonance imaging, each volume identified as a certain organ, that organ or part of that organ has a type of tissue for which conductivity values are available. At the time of defining the numerical model, these values are assigned depending on the tissue identified. In case the conductivity is not given to the system, the system will assume a constant conductivity value, typically a conductivity value of 1, as it does not affect the mathematical properties of the problem.

[0063] According to the third option, information on the degree of tissue dysfunction isavailable. This degree of dysfunction is expressed by a parameter w(x) that varies between 1 and f> where f> is a value greater than 1 and therefore positive. The degree of dysfunction is the function w(x) which takes these values and weights the function / (%) by giving it a higher weight.

[0064] According to an embodiment of the previously disclosed embodiment, the real regularization function is a function responsive to the weighted measure, wherein: - w(x) = 1 when there is no dysfunction of the tissue at location x,

[0065] - w(x) takes a value greater than 1 and less than f> when the tissue at location x has a condition that negatively affects the magnitude of the allowed electrical activity and, - w(x) = f> when the tissue at location x has a condition that prevents the electrical activity.

[0066] According to an embodiment of the previously disclosed embodiment, the condition of the tissue affecting the magnitude of the allowed electrical activity is characterized by the degree of fibrosis.

[0067] One way to determine the value of w(x) is to establish a relationship between the condition that negatively affects electrical activity. For example, if there is tissue damage such as a certain degree of fibrosis, the quatification of fibrosis is functionally linked to w(x), for example by making it proportional or by a predetermined function or, more specifically by making it linearly proportional or by a predetermined monotonically increasing or monotonically decreasing function.

[0068] Quantification of the degree of fibrosis and its relationship to w(x) can be performed on biological samples that are not from the patient but allow a relationship to be established between the condition that adversely affects the magnitude that allows electrical activity and w(x). The negatively affecting condition can for example be assessed through magnetic resonance imaging of the subject.

[0069] According to an embodiment of any of the previously disclosed embodiments, the numerical model is generated from a 3D image of the anatomy of the torso and the heart, preferably from a computer tomography or a magnetic resonance image acguisition from the subject.This example has already been cited above. The acquisition of a 3D image of the subject's torso makes it possible to establish the anatomy or shape of the organs, in particular the heart, and therefore the geometry and properties of the computational domain. Typically, the3D image from CT (Computed Tomography) or MRI (Magnetic Resonance Imaging) is a volume or collection of planes representing scalar values associated to some properties of the tissues. On these grayscale values, segmentation algorithms can be applied aimed at defining the geometries or shapes of the different tissues and organs. An example of an MRI acquisition procedure can determine the torso geometry from scouting sequences, the heart shape from CINE sequences typically based on bSSFP (balanced Steady-State Free Precession), and the determination and characterization of fibrotic areas by mean of LGE-RMI (Late Gadolinium Enhancement) acquisition.

[0070] According to another embodiment, the 3D geometry is obtained by photogrammetry. According to an embodiment applicable to the previous embodiment, a database with a plurality of previously acquired hearts and corresponding torsos are available. With the torso shape, a selection of the heart shape is carried out, among the heart shapes stored in the database, and this is adjusted to the dimensions corresponding to the torso dimensions.

[0071] The geometry of the 3D geometrical model is required in the system. A preferred example of the geometry representation makes use of elements with nodes connected by triangles for surfaces descriptions or by tetrahedra for filled volumes.

[0072] According to an embodiment of any of the previously disclosed embodiments, the first sub-domain (Ω₁) is a domain of the numerical model comprising the right atria, the left atria, the right ventricle, the left ventricle, or any combination of them.

[0073] Setting the first subdomain

[0074]

[0075] as at least a part of the heart allows to limit the computer domain where / (%) is nonzero. This allows the evaluation of the integrals in operators A and R to be extended to a smaller domain and thus increases the computational efficiency of the computation even when the entire domain is very large.

[0076] According to an embodiment of any of the previously disclosed embodiments, the firstsub-domain (Ω₁) and the second sub-domain (Ω₂) are discretized according to a first discretization for numerically solving the set of differential eguations, and according to a second discretization for determining the first and the second contribution, wherein preferably the first discretization and the second discretization are the same discretization.

[0077] As described above, the functions Gy(x) and the operators A and R can be expressed in different discretizations but, nevertheless, the use of the same discretization makes the calculation more efficient, mainly because the functions Gy(x), the ones initially calculated, do not have to be interpolated and are already directly expressed in the second discretization. This not only avoids additional computer cost in the interpolation process but also reduces rounding and approximation errors due to the use of different discretizations.

[0078] According to an embodiment of any of the previously disclosed embodiments, functions Gy(x), i = 1,..., N are determined solving the set of differential equations by a finite element method, where N correspond to the torso location where measurements are available

[0079] The finite element method is the preferred method of discretization since it allows the subdomains to be adapted to complex shapes with known computer tools, e.g. for mesh generation.

[0080] The finite element method provides a discretization technique that also defines the way to approximate the solution at locations of an element other than nodes.

[0081] According to an embodiment of any of the previously disclosed embodiments, the second sub-domain (Ω₂) is segmented into a plurality of sub-volumes, the second sub-domain (Ω₂) being a direct sum of the k sub-volumes, Γⱼ; j = 1,..., k + 1, k being a predetermined positive integer, is the set of boundaries between the plurality of sub-volumes, further including the torso surface (∂Ω), the conductivity σ(x) is taken as a constant value at each sub-volume, being σⱼ⁺ the conductivity inside the sub-volume Γⱼ and σⱼ⁻ the conductivity outside the sub-

[0082]

[0083] volume Γⱼ;the first operator A further comprises a third contribution given by a numerical computation of ∑ⱼ₌₁ᵏ Bⱼ = ∑ⱼ₌₁ᵏ (1 / σⱼ⁻ − 1 / σⱼ⁺) ∫_Γⱼ G_{yᵢ}(x) ∂ₙ(σ(x) ∇φ(x)) dS.

[0084]

[0085] j=i j=i ' J J '

[0086] The general expressions of operators A and R include the value of the conductivity <j(x) as a function of x as a parameter within the integral. However, a simple way to identify distinct tissue regions is to define sub-volumes and for each sub-volume populate a constant value of the conductivity. This strategy results in discontinuities at the boundary between sub-volumes of the conductivity value. The way to address the problem of the existence of discontinuities in the integrand, according to this embodiment, is to take into account an additional term for each boundary that takes into account the effect of the discontinuity on the conductivity value, the value Bj with index j identifying an integral value over one of the boundaries between sub-volumes identifying an inhomogeneous medium while the integrals defining the operators A and R are calculated in each of the sub-volumes where the value of the σ(x) is now constant and may even go outside the integral. Operator A further comprises the third contribution value of: ∑ⱼ₌₁ᵏ Bⱼ being the sum of all integral over the inhomogeneous regions identifying the heterogeneous medium.

[0087] According to an embodiment of any of the previously disclosed embodiments, in the minimizing step the norm expression to be minimized further comprises the contribution of a fourth guantity Me

[0088] being j a mute variable of the summation identifying an specific location zⱼ of measurement of the intracardiac tissue among a total number of measurements Me, <

[0089]

[0090] p(zj) the potential estimated at Zj and, g(zj) being the signal value measured at Zj.

[0091] According to this embodiment, a set of intracardiac point zⱼ potential measurements ĝ(zⱼ) are available, which are either acquired by means of an electrode in contact with the intracardiac tissue or are measurements obtained from databases storing such measurements previously acquired or estimated.

[0092] These measurements are potential values known from the acquired data; however,the method does not impose such values on the solution of the numerical problem but defines them as objectives to be achieved in the minimization process. That is, the solution will be the result of minimizing an objective function defined by the terms included within the norm, among which there is now a new term, the term that measures a distance between the values of the measurement and those of the estimated solution to be solved.

[0093] According to an embodiment of any of the previously disclosed embodiments, the method further comprises determining regions of the heart with biological tissue at which the initiation of the pulse occurs, wherein the method further comprises: - for each node at location yᵢ; i = 1, ..., Iₜ ≤ N_H, being Iₜ predetermined and among the plurality of nodes of the second discretization NH, determining the Local Activation Time (LAT); selecting the at least one node location yshaving the minimum Local Activation Time (LAT). determining the regions of the heart with biological tissue at which the initiation of

[0094]

[0095] the pulse occurs as the region defined by the least one node location ys.

[0096] According to an embodiment of the previous embodiment, the Local Activation Time (LAT) is determined according to the following steps:

[0097] selecting a function of time fs(t) either the function φ(yᵢ, t) or the function f(yᵢ, t) for node at location yᵢ; computing the instant time of the maximum negative slope of the temporal signal fs(t), being the Local Activation Time (LAT) the time fs(t) that has lapsed until the

[0098]

[0099] computed instant time.

[0100] There is a plurality of methods for determining the Local Activation Time (LAT).

[0101] According to any of the last two embodiments, knowledge of the electrical potential and its evolution over time are used to establish the earliest activation site, i.e., to determine the region formed by the few cells that cause the cardiac electrical activity to start and subsequently propagate throughout the rest of the cardiac tissue.

[0102] In cardiac tissue the pulse is initiated from a few cells that generate electrical charge and this causes an elevation of the electrical potential. After the elevation of the electrical potential there is a sharp drop which then recovers, stabilizing the electricalpotential value at that point. The signal shows a strong variation first with a positive peak and then with a negative peak, in a short period of time, after which the electrical potential stabilizes.

[0103] In many cases, the onset of the signal is located in the wavefront and, the wave-front position is approximated by the point where the maximum decrease occurs, i.e. the point where the derivative is negative and in absolute value reaches the maximum.

[0104] The method looks for this behavior in a plurality of points, preferably in all of them, and also mainly in the domain where the heart is, the first subdomain, since this is where the pulse is produced.

[0105] The time value where the condition of the maximum decrease of, the electrical potential or the value of the source, is obtained at different time instants in each place of the heart. The points where this signal is produced first are selected, interpreting that it is at these points where it starts and that the rest of the points show the propagation of this signal causing the heart to contract following a certain propagation pattern.

[0106] It has been observed that the strong variation in the electrical potential signal that occurs in a subject is not like the one obtained in the solutions of f

[0107]

[0108] f(x) and φ(x) according to the disclosed method since, when the response obtained in the simulation is observed, the signal is smoother. It is interpreted that the solution to the numerical problem is filtered by a low-pass filter eliminating strong local variations and delaying the estimated time instants when compared with values measured on the patient.

[0109] Surprisingly, however, although the simulation values are delayed in time, what is determined is the place where the first changes occur, so the accuracy in determining the place where the pulse onset is established and therefore the identification of the tissue causing the pulse is very accurate.

[0110] According to a further embodiment applied to the previous embodiment, Iₜ is equal to the total number of nodes of the first sub-domain (Ω₁).According to this embodiment, the number of nodes where the computation is executed is restricted since the origin of the pulse is expected to occur in a part of the heart.

[0111] According to an embodiment of any of the previously disclosed embodiments, the method further comprises:

[0112] determining at least one wavefront of propagation of an electric impulse, the location of the wave being those regions with a Local Activation Time (LAT) that are close, being close if they have a value of the Local Activation Time (LAT) that does not differ more than a certain predetermined threshold value.

[0113] According to this embodiment a map of LATs may be generated determining the Local Activation Time at each point.

[0114] In particular, after the previous step, according to a further embodiment, the method also comprises:

[0115] determining whether a set of wavefronts of propagations follows a path of propagation converging in a loop.

[0116] According to this specific method, reentries of the electrical impulse are detected allowing to identify arrhythmias, in particular regular arrhythmias.

[0117] Once the Local Activation Time (LAT) has been computed for each node or a plurality of nodes, wave fronts can be determined traveling through the heart tissue. The evolution of the wave fronts are useful for identifying the origin and direction of electrical impulses and detecting abnormalities like reentrant circuits, which are loops of electrical activity that can cause arrhythmias. The evaluation of the arrhythmias are carried out by a medical doctor under his / her criterion in view of the evolution provided by the method.

[0118] According to an embodiment of any of the previously disclosed embodiments, the method further comprises for determining the conduction velocity of the electrical impulse:

[0119] - for each node, computing at least one distance between at least one neighboring node;determining the speed of propagation of the electrical impulse between two nodes as the rate between the distance and the Local Activation Time (LAT); determining regions with a conduction velocity having a value measured according to a predetermined norm as being below a predetermined positive threshold value.

[0120] This approach results in a conduction velocity map that allow to represent the conduction velocities across the heart's surface.

[0121] Conduction velocity maps are crucial for identifying and analyzing slow conducting regions and / or reentrant circuits, which cause arrhythmias like atrial fibrillation and ventricular tachycardia. These maps highlight areas with slow conduction velocities, indicating potential reentrant circuit formation, aiding in targeted interventions such as ablation therapy. They also assess treatment effectiveness by comparing pre- and post-procedure conduction velocities.

[0122] According to an embodiment of any of the previously disclosed embodiments, wherein the method further comprises:

[0123] - transforming the electrical signal at each node in a period of time to a phase function, preferably by applying a Hilbert transform and, in the interval [−π, π];

[0124] - determining if the domain comprises a spiral or a rotors with a central core.

[0125] Other non-limiting examples of electroanatomic maps that can be constructed from the estimated electrical signals are voltage, rotor histogram, dominant frequency, singularity points, slew-rate, fractionation, entropy or propagation maps.

[0126] A phase map of the domain allows to identify reentrant electrical circuits. The transform converts the time-domain signals into a phase representation that cycles continuously from −π to π. The transformation can be done using methods like the Hilbert transform, which extracts the instantaneous phase of the signal.

[0127] This transformation makes it possible to identify patterns in the electrical activity, such as spirals or rotors, which provides valuable information to the medical doctor.A second aspect of the invention is a data processing system comprising means for carrying out the steps of any of the methods previously disclosed. In particular, the system comprises:

[0128] at least one first electrode configured to be located over the skin of a patient and adapted to acquire a signal responsive to the electric potential of the skin at its location;

[0129] a processor at least configured to receive the signal of the at least one electrode; wherein the processor is adapted for carrying out the steps of the method of any of the previously disclosed embodiments.

[0130] An embodiment of the invention according to the second aspect of the invention is a system further comprising a set of electrodes wherein such electrodes are in electrical contact with the surface of the patient's skin on which measurements are to be taken. For instance, the electrodes comprise a surface adapted to be in contact with the surface of the patient's skin and are held in contact by an adhesive.

[0131] When at least one first electrode is a plurality of electrodes, this plurality of electrodes is a first set of electrodes intended to acquire electrical potential measurements at distinct points on the patient's skin.

[0132] According to another embodiment of the system, it further comprises at least one second electrode configured to be located in an intracardiac tissue location of the patient and adapted to acquire a signal responsive to the electric potential of the tissue at its location and, the processor is further adapted to carry out a method according to the embodiment using, in the step minimizing the norm expression to be minimized, the

[0133]

[0134] contribution of a fourth quantity ∑ⱼ₌₁^Mₑ (φ(zⱼ) − ĝ(zⱼ))²

[0135] According to another embodiment of the system, the measurements of at least one second electrode or from a plurality of second electrodes, are retrieved from a data base storing such measurements.

[0136] A third aspect of the invention is a computer program product comprising instructions which, when the program is executed by a computer, cause the computer to carry out steps of any of the disclosed methods.DESCRIPTION OF THE DRAWINGS

[0137] These and other features and advantages of the invention will be seen more clearly from the following detailed description of a preferred embodiment provided only by way of illustrative and non-limiting example in reference to the attached drawings.

[0138] Figure 1 This figure shows the shape of a heart over-imposing the electrical potential simulated according to a first embodiment of the method. The heart is identified with roman letter (I). A small region is zoomed in and identified with roman letters (II). The normalized signals are presented in the right lower part of the figure and identified with roman letters (III).

[0139] Figure 2. Figure 2 shows a sectional view of the simulated heart with three locations (I, II and III) at a certain instant time.

[0140] Figure 3. Figure 3 shows three signals, one graphic per location (I, II and III) identified in figure 2.

[0141] Figure 4. Figure 4 shows the shape and values of functions Gy.(x) according to the measurement of the electric potential of an electrode at location y. Functions Gy(x) are shown on domain Ω wherein the first subdomain (Ω1) is shown as ΩH, the heart domain.

[0142] Figure 5. Figure 5 shows three different views for highlighting the importance of the method estimating the electrical activity in the mid-myocardial tissue.

[0143] DETAILED DESCRIPTION OF THE INVENTION

[0144] As will be appreciated by one skilled in the art, aspects of the present invention may be embodied as a system, method or computer program product.

[0145] The specific embodiment applies to a subject on which electrical propagation is to be determined, specifically determining the origin of the signal causing the pulse actingon his heart and the mode of propagation of said signal in the tissue of the heart.

[0146] As a first step, according to this embodiment, the subject is subjected to a 3D imaging test, an image obtained by magnetic resonance techniques, resulting in a 3D image formed by grayscale voxels. The image is acquired over the torso and includes the heart.

[0147] A segmentation algorithm is used to differentiate the organs and tissues. This segmentation makes it possible to generate a numerical model in which a first subdomain Ω1contains the heart and a second sub-domain Ω2contains the rest of the domain including for instance other organs, with different tissues.

[0148] Figure 4 shows a torso (T) being modeled as domain (Ω) wherein the heart is a subdomain, represented as ΩHin figure 4 and identified as a first sub-domain (Ω1) along the description. The second subdomain (Ω2) is the remaining region of Ω, that is

[0149] The numerical model makes use of a finite element discretization of one (Ω1) and the other sub-domain (Ω2) such that the shape of the tissue structures of the torso region are defined by elements and the discretization also comprises nodes in which the values of various functions are determined, in particular the solutions that determine values of interest.

[0150] In a further stage a plurality of electrodes are incorporated on the surface of the torso, the plurality of electrodes located at N different locations yi; i = 1,..., N for the measurement of the potential over a certain period of time. These electrodes are for example attached on the skin of the subject in such a way that they are in electrical contact with the skin.

[0151] The measurement of the potentials giat the predetermined positions (yi; i = 1,..., N) generates N readout signals which are further processed.

[0152] If the signals are captured at a single instant, the method provides the distribution of potential sources and the distribution of the potential at the instant of measurement. If the measurement is carried out over time, the values of the sources and of thepotential are obtained for each time instant at which a sample of the measurement has been carried out.

[0153] In this embodiment, the discretization of the surface is a closed triangulation comprising nodes and edges of the torso surface. The volume is discretized using a tetrahedralization matching at the surface with the triangles of the discretization of such surface. It can include the neck, arms, and other body parts. In this embodiment, the locations of the electrodes correspond to the location of nodes of the discretization. The heart is a tetrahedralization, also comprising nodes and edges forming a closed triangulation at each face of the elements, representing the endoepicardium, being the first sub-domain (Ω1). In this case, the heart geometry is the right atria, left atria, right ventricle and left ventricle. According to other embodiments, the elements of the first sub-domain (Ω1) may be any combination of the cited four chambers.

[0154] According to another embodiment wherein a plurality of regions has different conductivity, the discretization may also use different elements.

[0155] Once the set of electrodes are already installed on the skin surface and the numerical model generated, the next step is the acquisition of the body surface potentials measured at the electrode locations with a measurement system.

[0156] The signals, according to an embodiment, are filtered with a proper algorithm, in order to clean up the signal from noise. The signals are assumed to be acquired in a time window of clinical interest.

[0157] According to this embodiment, the functions Gy.(x); i = 1,..., N and the operators A and R are expressed in matrix form, resulting in a computationally efficient method when running in a computer system.

[0158] Body surface potentials measured by the electrodes are defined as g ∈ E, being E the set of electrodes. These measurements along nttime samples are stored into a matrix G of nErows and ntcolumns where nEis the number of electrodes (equal to N).

[0159] According to this embodiment, the torso (T) of the subject being modeled by thenumerical model is discretized using a tetrahedralization, each face of the tetrahedron being a triangle, defining a discretization of finite elements. The discretization comprises the location of the electrodes and the discretization is extended at least to the first sub-domain (Ω1) comprising at least one part of the heart and to the second sub-domain (Ω2).

[0160] Once the tetrahedral mesh of the discretization has been generated according to the shape of the tissue structures (according to the physiology), functions Gyi(x) are calculated for certain boundary conditions. In particular, in this embodiment a = 0, b = 1 and c = —1 / Sr where Sr denotes the area of the surface of the boundary. The resulting boundary condition is expressed as:

[0161] ∂nGy(x) = −1 / Sr

[0162] The gradient potential in the normal direction of the surface is zero since, by hypothesis the body of a subject is electrically isolated and therefore the electrical field does not exit from it.

[0163] The solution of the differential equation −∇·(σ(x) ∇Gy(x)) = δy(x), results in a matrix Gyof one row and nΩcolumns for i = 1,...,nE, where nΩis the total number of nodes of the body tetrahedralization.

[0164] It is noted that, due to cardiac sources f(x) are set to zero outside the heart, Gyresults in a smaller matrix GyΩof one row and nHcolumns where nHis the number of nodes of the heart volume.

[0165] Once GyΩhas been computed operators A and R must be built.

[0166] Operator A acting on f is Af and defined by the integral value

[0167] f Gy.(x)f(x)dV

[0168]

[0169] wherein using the tetrahedralization (the discretization used in this embodiment), the integral is computed using quadrature formulae producing matrix A of nErows and nHcolumns wherein Af is approximated by matrix A multiplied by vector f of nHcomponents. It is noted that each row of A has a "weighted" function Gy.(x) restricted to the heart nodes.When the domain Ω comprises several regions of different piecewise constant conductivities, operator A is incremented with ΣjMBwherein MBis a matrix of nErows and nHcolumns built as the weights the quadrature schemes used for volume and surface integrals required to compute the term Bjused to build the matrix expression MB, given by

[0170] k k..

[0171] Σj=1kBj= Σj=1k(1 / σj−− 1 / σj+) ∫ΓGy(x) ∂n(σ(x)∇φ(x))d ∂Ω

[0172] j

[0173]

[0174] =i j=i ' J J 'Jri

[0175] for k the number of boundaries between regions. Regions may have different conductivity. If regions have the same conductivity these terms are zero and operator A does no need to be incremented.

[0176] The other operator to be discretized is operator R that is applied to the discrete vector (p, the electric potential at the nodes.

[0177] Operator R may be expressed from an intermediate operator Q acting over g, Qg, wherein

[0178] Qg ≈ ∫∂Ω∂nGy(x)σ(x)g̅(x)dS

[0179]

[0180] Jan

[0181] where g̅(x) denotes function g restricted to the surface S or also denoted as ∂Ω. Operator Q is a square matrix of nTrows and nTcolumns denoting the quadrature matrix computing ∂nGy(x), restricted to the surface, numerically. Therefore, R = I + Q, being I the identity matrix having 1's on the index that corresponds to the electrode locations and is the operator computing g̅(x) +

[0182]

[0183] ∂nGy(x)σ(x)g̅(x)dS when applied to the discretized version of g(x), that is, expressed in vector form.

[0184] One of the calculated values is the cardiac sources f at the nodes of the discretization and then the following condition is imposed:

[0185] ∫Ωf(x)dV = − ∫∂Ω∂n(σ(x)∇φ(x))dS

[0186] <

[0187]

[0188] a Jan

[0189] already identified as the Existence condition. According to an embodiment, additional restrictive conditions are taken where the right term of the integral is constant at each time instant and then can be discretized using a numerical scheme as a vector vtof 1 row and ntcolumns. Then, the left term can be discretized using a quadrature scheme such that

[0190] REF ≈ ∫Ωf(x)dV

[0191]

[0192] ■J 0.where REis a vector of one row and nHcolumns and F is a matrix of nHrows and ntcolumns stacking the sources along the time evolution and, as a result, the existence condition is transformed into a linear restriction of the form REF = vt.

[0193] Once all the cited matrices are built, the following discrete optimization problem is solved: arg minF||AF − RG|| such that REF = vt.

[0194] Now, returning to figure 1, it shows a spatial and temporal visualization of the cardiac sources determined according to the disclosed method.

[0195] At the left side and identified with roman number (I), there is a representation of the cardiac sources at a fixed time instant. The heart is shown according to the numerical model defining the shape and properties of the heart.

[0196] The same figure shows a wave front produced by the cardiac sources over the endo-epicardial surface. The representation is normalized for a better visualization.

[0197] In the center and at the upper part of figure 1, identified with roman number (II), there is a zoomed in view of the mid-myocardial region, represented with a point could. This image shows how the wave front is produced transmurally.

[0198] In the right side of the same figure, identified with roman number (III), a graphic shows the temporal behavior of the cardiac source for a node in the aforementioned cloud of points. Two separated variables are represented: the action potential (identified as TVMs) and the cardiac source. The horizontal axis represents time, and the vertical axis displays arbitrary units for the normalized variables.

[0199] Figure 2 shows a map of the Local Activation Times (LATs), that is, the time that elapses from the origin of time until the moment at which the steepest decreasing slope of the electric potential is measured. These instants, in the heart domain, occur between instants 0 to 110 milliseconds. Three separated locations are identified as (I, II, and III).

[0200] Figures 2 and 3 allow to visualize at least three variables defined on the heart (actionpotential, electric potential, and cardiac sources).

[0201] The three nodes identified in figure 2 are (I) a node near the origin of activation, (II) a node at an intermediate point of propagation, and (III) a node near the latest place of activation.

[0202] Figure 3 shows three graphics, the different variables visualized for the three nodes at the chosen location: the action potential, the cardiac source, and the electric potential, also called electrograms, or EGMs. For each plot, the horizontal axis represents time, and the vertical axis displays arbitrary units for the normalized variables obtained by the disclosed method.

[0203] Figure 4 shows a visualization of functions Gy.(x), the functions that relate the cardiac sources to the electric potential on an electrode. Different regions of interest can be observed. The full domain, denoted as 1, is the volume encompassed by the torso surface T = d£l. Furthermore, the heart region is included in the domain £lHc 1. The location of an electrode y G T is also high-lighted, and the contour plot of the Gy(x) function related to this electrode is shown. The closest regions have a higher contribution whereas the farthest regions have a lower contribution.

[0204] Figure 5 shows the relevance of the estimation of the electrical activity in the mid-myocardial tissue according to the disclosed method. Figure 5 shows the ectopic beat originating from the septum and is identified on various heart structures. Figure 5A, shows the activation times on the epicardial surface, corresponding to the classical ECGI model. Figure 5B, shows a new internal surface allowing for the representation of activation on an endo-epicardial surface. This corresponds to more advanced techniques like the Equivalent Double Layer, which aims to estimate action potentials. Figure 5C shows the entirety of the heart volume, i.e., the mid-myocardial tissue, and the activation time is represented over a point cloud. This corresponds to the proposed volumetric ECGI approach according to the disclosed embodiment.

Claims

CLAIMS1.- A computer-implemented method for determining the electric activity at the heart of a subject by external measurements, the method, carried out in a computer system, comprising the flowing steps:defining a numerical model of the torso of a subject, the numerical model at least comprising a domain (fl) representing a volume of the torso bounded by the torso surface (dfl) and by one or more closure surfaces, the domain being divided into a first sub-domain (fl representing at least one part of the heart and a second subdomain (0.2) representing the remaining volume, the numerical model at least comprising the electrical properties of the medium, preferably the conductivity (<j(x)), being x a point of the domain (fl);obtaining as input data body surface potentials gimeasured on a plurality of N preestablished locations (yi; i = 1,...,N) of the torso surface (∂Ω) of the subject, at least at nt≥ 1 instant times;determining functions Gy(x), i = 1,...,N solving numerically, for each Gy(x), the set of differential equations:-V-fa(x)VGy.(x)) = 5y.(x), x G fl aaGy(x) + b∂nGy(x) = c(x), x ∈ ∂Ω wherein Gy(x) is a function defined in the domain (Ω),wherein a, b ∈ R are two preestablished values, σ(x) is the conductivity function and c(x) is a preestablished function over the torso surface ∂Ω,wherein ∇ · is the divergence operator, ∇ is the gradient operator, δy(x) is the Dirac's delta distribution with center on the location yi, and ∂ndenotes the outward normal derivative with respect to the torso surface ∂Ω;defining a first operator A defined such that, when A acts on a function f(x), represented by Af(x), f(x) being a scalar function defined on domain (Ω) and representing the cardiac source, the first operator A at least determines a first contribution given by a numerical computation off Gy.(x) / (x)dVJo.being f(x) = 0 in the second sub-domain (Ω2) and, the integral is a volume integral over the domain (Ω);defining a second operator R defined such that, when R acts on a function φ(x), represented by Rφ(x), φ(x) a scalar function defined on domain (Ω) and representing the electric potential, the second operator R at least determines asecond contribution given by a numerical computation ofφ(x) + ∫∂Ω∂nGy(x)σ(x)φ(x)dS;the integral being a surface integral restricted to the torso surface (∂Ω), wherein φ(x) satisfies φ(yi) = gi, i = 1 ... N;obtaining the estimation of f(x) and φ(x) of the subject as the result of minimizing, the contribution of at least a first quantity Af(x) and a second quantity −Rφ(x), measured under a predetermined norm ||·||, namely ||Af(x) − Rφ(x)||,satisfying that∫Ωf(x)dV = 0;determining the electric activity of the heart of a subject by at least one value of f(x), φ(x) or both at least at one instant time.2.- A method according to claim 1, wherein the minimizing step further comprises the contribution of λh(f(x)), being λ a real regularization parameter, a real number, and h(·) being a preestablished real regularization function.3.- A method according to any of the previous claims, wherein the numerical model is generated from an acquisition measurement of the torso and the shape of the heart.4.- A method according to any of the previous claims, wherein first sub-domain (Ω1) is a domain of the numerical model comprising the right atria, the left atria, the right ventricle, the left ventricle, or any combination of them.5.- A method according to any of the previous claims, wherein the first subdomain (Ω1) and the second sub-domain (Ω2) are discretized according to a first discretization for numerically solving the set of differential equations, and according to a second discretization for determining the first and the second contribution, wherein preferably the first discretization and the second discretization are the same discretization.6.- A method according to the previous claim, wherein functions Gy.(x), i = 1,..., N are determined solving the set of differential equations by a finite element method.7.- A method according to any of the previous claims, whereinthe second sub-domain (Ω2) is segmented into a plurality of k sub-volumes, the second sub-domain (Ω2) being a direct sum of the k sub-volumes, k being a predetermined positive integer,Γj; j = 1,..., k + 1, is the set of boundaries between the plurality of sub-volumes, further including the torso surface (∂Ω),the conductivity σ(x) is taken as a constant value at each sub-volume, being σj+the conductivity inside the sub-volume Γjand σj−the conductivity outside the subvolume Γj;the first operator A further comprises a third contribution given by a numerical computation ofk k..Σj=1kBj= Σj=1k(1 / σj−− 1 / σj+) ∫ΓGy(x) ∂n(σ(x)∇φ(x))d ∂Ω.j=i j=i '7 7'Jri8.- A method according to any of the previous claims, wherein the method further comprises determining regions of the heart with biological tissue at which the initiation of the pulse occurs, wherein the method further comprises:for each node at location yi; i = 1,..., It≤ NH, being Itpredetermined and among the plurality of nodes of the second discretization NH, determining the Local Activation Time (LAT);selecting the at least one node location yshaving the minimum Local Activation Time (LAT).determining the regions of the heart with biological tissue at which the initiation of the pulse occurs as the region defined by the least one node location ys.9.- A method according to the previous claim, wherein the Local Activation Time (LAT) is determined according to the following steps:selecting a function of time fs(t) either the function φ(yi, t) or the function f(yi, t) for node at location yi;computing the instant time of the maximum negative slope of the temporal signal being the Local Activation Time (LAT) the time fs(t) that has lapsed until the computed instant time.10.- A method according to claim 8 or 9, wherein Itis equal to the total number of nodes of the first sub-domain (fli).11.- A method according to any of claims 8 to 10, wherein the method further comprises:determining at least one wavefront of propagation of an electric impulse, the location of the wavefront being those regions with a Local Activation Time (LAT) that are close, being close if they have a value of the Local Activation Time (LAT) that does not differ more than a certain predetermined threshold value.12.- A method according to the previous claim, wherein the method further comprises:determining whether a set of wavefronts of propagations follows a path of propagation converging in a loop.13.- A method according to any of claims 8 to 12, wherein the method further comprises for determining the conduction velocity of the electrical impulse:for each node, computing at least one distance between at least one neighboring node;determining the speed of propagation of the electrical impulse between two nodes as the rate between the distance and the Local Activation Time (LAT); determining regions with a conduction velocity having a value measured according to a predetermined norm as being below a predetermined positive threshold value.14.- A method according to any of claims 8 to 13, wherein the method further comprises:- transforming the electrical signal at each node in a period of time to a phase function, preferably by applying a Hilbert transform and, in the interval [−π, π];- determining if the domain comprises a spiral or a rotors with a central core.15.- A method according to any of previous claims, wherein in the minimizing step the norm expression to be minimized further comprises the contribution of a fourth quantityMej=1being j a mute variable of the summation identifying an specific location Zj of measurement of the intracardiac tissue among a total number of measurements Me, <p z ) the potential estimated at Zj and, g(z ) being the signal value measured at Zj.16.- A system comprisingat least one first electrode configured to be located over the skin of a patient and adapted to acquire a signal responsive to the electric potential of the skin at its location;a processor at least configured to receive the signal of the at least one electrode; wherein the processor is adapted for carrying out the steps of the method of any of claims 1 to 15.17.- A system according to claim 16 further comprising at least one second electrode configured to be located in an intracardiac tissue location of the patient and adapted to acquire a signal responsive to the electric potential of the tissue at its location, wherein the processor is further adapted to carry out a method according to claim 14.18.- A computer program product comprising instructions which, when the program is executed by a computer, cause the computer to carry out steps of any of claims 1 to

Images