Some inverse problems in electrocardiography Jocelyne Erhel and Edouard Canot Projects from M. Hayek and E. Deriaz ALADIN team - INRIA-Rennes • Models in electrocardiography and discretisation • Discrete ill-posed general least-squares problem • Time regularisation ERCIM workshop - February 2002 J. Erhel - 02/02 1
Electrocardiography 0 −20 −40 −60 −80 −100 −120 −140 −160 0 50 100 150 200 source of electrical current at the surface of the heart : Γ E propagation through the chest, bioelectric volume conductor : Ω measure of potential at the surface of the torso : Γ T J. Erhel - 02/02 2
Direct and inverse models in ECG Φ is the electrostatic potentiel σ is the electrical conductivity tensor Direct model ∇ . ( σ ∇ Φ) = 0 in Ω ∂ Φ (1) ∂n = 0 on Γ T Φ = Φ E on Γ E Well-posed problem Inverse model ∇ . ( σ ∇ Φ) = 0 in Ω ∂ Φ (2) ∂n = 0 on Γ T Φ = Φ T on Σ ⊂ Γ T Unique solution but not continuous : ill-posed problem J. Erhel - 02/02 3
Discretisation by Boundary Element Methods Well-suited for homogeneous and isotropic case ( σ scalar and constant) n nodes on Γ E and m ≥ n nodes on Γ T x 1 : discretisation of Φ E x 2 : discretisation of ∂ Φ E ∂n y : discretisation of Φ T � x 1 � A 11 A 12 A 13 = 0 (3) x 2 A 21 A 22 A 23 y Dense small linear system J. Erhel - 02/02 4
Discretisation by Finite Element Methods Well-suited for heterogeneous or anisotropic case ( σ tensor) n nodes on Γ E , N nodes in Ω and m ≥ n nodes on Γ T x 1 : discretisation of Φ E x 2 : discretisation of Φ in Ω y : discretisation of Φ T � x 1 � A 11 A 12 A 13 = 0 (4) x 2 A 21 A 22 A 23 y Sparse large linear system ( N large) J. Erhel - 02/02 5
Direct problem In BEM, we set N = n Square linear system of order N + m � � � � � � − A 11 A 12 A 13 x 2 = (5) x 1 − A 21 A 22 A 23 y which can be solved by ( A 23 − A 22 A − 1 12 A 13 ) y = ( A 22 A − 1 12 A 11 − A 21 ) x 1 (6) x 2 = − A − 1 12 ( A 11 x 1 + A 13 y ) Transfer matrix : y = Tx 1 , T ∈ R m × n J. Erhel - 02/02 6
Inverse problem N + n unknowns and N + m equations : least-squares problem y : measure of potential on the torso : errors y + e with e blank noise (0 ,µ 2 I ) Approach used in most papers and software min x 1 � y − Tx 1 � Our proposal min x,e � e � subject to Ax + Be = By x = ( x 1 ,x 2 ) ∈ R n + N , y ∈ R m A ∈ R ( N + m ) × ( n + N ) , B ∈ R ( N + m ) × m General Gauss-Markov linear model J. Erhel - 02/02 7
Regularisation Algorithms for general linear models Generalised Singular Value Decomposition GSVD Generalised QR factorisation GQR (Paige ’s algorithm) Iterative methods of type LSQR? Regularisation for discrete ill-posed problem x,e ( � e � 2 + λ 2 � Cx � 2 ) subject to Ax + Be = By min Ref : H. Zua and P.C. Hansen, 1990 Restricted SVD Algorithm similar to GQR Iterative methods of type LSQR? J. Erhel - 02/02 8
Some very preliminary results Discretisation by BEM - homogeneous and isotropic case - 2D Linear model min x � z − Bx � with z = By Tychonov regularisation and parameter selection using an L-curve Régularisation du problème inverse perturbé par un bruit d’amplitude 1e−3 1.5 solution régularisée solution analytique 1 0.5 ( φ c , q c ) 0 −0.5 Le potentiel électrique à la La dérivée normale du surface du potentiel électrique à la coeur ( φ c ) surface du coeur (q c ) −1 −1.5 0 20 40 60 80 100 120 140 nombre total de noeuds = 128 J. Erhel - 02/02 9
Time dependent problem Time dependent model Measures during a time interval Time discretisation ( t 1 , . . . ,t p ) Measure Y = ( y ( t 1 ) , . . . ,y ( t p )) Unknown X = ( x ( t 1 ) , . . . ,x ( t p )) Matrices A and B independent of time AX = B ( Y + E ) where A ∈ R ( N + m ) × n + N , B ∈ R ( N + m ) × m , X ∈ R ( n + N ) × p , Y ∈ R m × p ( I p ⊗ A ) vec ( X ) = ( I p ⊗ B )( vec ( Y ) + vec ( E )) Case of total blank noise p independent problems Ax ( t k ) = B ( y ( t k ) + e ( t k )) , k = 1 , . . . ,p J. Erhel - 02/02 10
Some very preliminary results 100 80 60 40 20 0 −20 −40 −60 −80 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Same algorithm as before for each time step J. Erhel - 02/02 11
Time dependent regularisation Approach based on the SVD of Y with m ≤ t Ref : F. Greensite and G. Huiskamp, 1998 Y = P ( S 0) Q T and Q = ( Q 1 Q 2 ) AXQ 1 = PS + EQ 1 and AXQ 2 = EQ 2 Preliminary method Solve min Z � AZ − PS � and take X = ZQ T 1 Solve p independent problems z ( t k ) � Az ( t k ) − ( PS )( t k ) � , k = 1 , . . . ,p min Regularisation for each ( PS )( t k ) Solution discarded if the Discrete Picard condition is not satisfied J. Erhel - 02/02 12
Some very preliminary results 10 10 0 0 −10 −10 −20 −20 −30 −30 −40 −40 −50 −50 −60 −60 −70 −70 0 10 20 30 40 50 60 70 80 0 10 20 30 40 50 60 70 80 Discrete Picard condition satisfied Discrete Picard condition not satisfied J. Erhel - 02/02 13
Some very preliminary results 100 80 60 40 20 0 −20 −40 −60 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 Solution using time regularisation J. Erhel - 02/02 14
Some perspectives • Regularisation of general linear model • Analysis of time regularisation • General time regularisation • Comparison with the classical approach J. Erhel - 02/02 15
Recommend
More recommend