 # Regularization of EIT reconstruction based on multi-scales wavelet transforms

• , Benjamin Schullcke , Sabine Krueger-Ziolek and Knut Moeller

## Abstract

Electrical Impedance Tomography (EIT) intends to obtain the conductivity distribution of a domain from the electrical boundary conditions. This is an ill-posed inverse problem usually solved on finite element meshes. Wavelet transforms are widely used for medical image reconstruction. However, because of the irregular form of the finite element meshes, the canonical wavelet transforms is impossible to perform on meshes. In this article, we present a framework that combines multi-scales wavelet transforms and finite element meshes by viewing meshes as undirected graphs and applying spectral graph wavelet transform on the meshes.

## 1 Introduction

Electrical Impedance Tomography (EIT) is a radiation-free imaging method . A commonly used imaging strategy is called time-difference reconstruction. It attempts to reveal the conductivity distribution changes inside the human body between two time points through the electrical data obtained via electrodes attached on skin. In time-difference lung EIT typically 16 electrodes are placed equidistantly on a horizontal chest plane. For each EIT frame, currents are successively injected into the human body through adjacent electrodes. A conductivity image is then reconstructed by collecting the voltage measurements recorded from the remaining electrodes.

We denote the conductivity changes of the domain between two time steps by Δs and the measured voltage changes on the electrodes by a vector ΔV. Employing the finite element model (FEM) with M elements, the conductivity change Δs is represented by a M × 1 vector. For time-difference imaging, the reconstruction problem can be linearized . Approximately, there exists the following relation:

(1)JEΔsΔV

where JE denotes the Jacobian matrix of elements calculated at the constant conductivity 1:

(2)JijE=ΔViΔsj|1

This Jacobian matrix of elements could be calculated by EIDORS toolbox . The Jacobian matrix JE is ill-conditioned because the degree of freedom of this inverse problem is too large. Solving the conductivity changes from the above relation leads to an ill-posed inverse problem. Especially, the solution Δs from equation (eq 1) is unstable. To circumvent this difficulty, additional prior information is introduced to restrict the flexibility of the solution . In this study, we employ a prior illustrating the sparsity hypothesis. Explicitly, we solve Δs from the following regularized optimization problem instead of (eq 1):

(3)Δs^=argminΔs12ΔVJEΔs22+αRΔs1

where α is a regularization parameter that controls the trade-off between the regularization term and the fidelity term. According to sparse regularization theory, the regularization term RΔs1 promotes a solution Δs^ with RΔs^ being sparse.

Triangular finite element meshes are commonly used in EIT reconstructions. Triangular meshes can simulate irregular domains as well as the electrical properties around electrodes without approximation. Wavelet transforms have been widely used in medical imaging such as Computed Tomography or Magnetic Resonance Imaging and shown a great success in these imaging techniques. However, canonical wavelet transforms are not appropriate to apply on general finite element meshes. Instead, in this article we used spectral graph wavelet transforms and view the triangular meshes as graphs.

## 2 Method

### 2.1 Nodal Jacobian matrix

To make triangular finite element meshes compatible with graph wavelet transforms, we need to view the meshes as undirected graphs. Instead of reconstructing the conductivity changes on triangular elements we pursue a solution of (eq 3) on the nodes. This relies on reforming the element-based Jacobian JE to a node-based Jacobian matrix, i.e. nodal Jacobian. We denote this nodal Jacobian matrix as JN. The entry JijN of Nodal Jacobian records the i-th voltage measurement changes with respect to small conductivity changes on the j-th node. For 2D meshes, we identified the conductivity on each triangular element as the linear interpolation of the conductivities on its three vertices. According to this assumption, the nodal Jacobian can be calculated from the element based Jacobian JE (see Algorithm 1). Using this nodal Jacobian, conductivity changes on nodes can be reconstructed. To complete the reconstruction process, we need additionally interpolate the conductivity information on nodes back to elements.

Algorithm 1

Nodal Jacobian .

 Input: the Jacobian with respect to elements: JE. Output: the Jacobian with respect to nodes: JN. For each node n in the mesh: elems = list of elements using node n, J:,nN=∑i∈elems13⁢J :,iE where J :, i denotes the ith column of matrix J End for

### 2.2 Spectral graph wavelets

In the proposed method, the finite element mesh was viewed as an undirected connected graph G = {V, E, A} consisting of a set of indexed vertices V with |V| = N vertices, a set of edges E, and the adjacency matrix A. With this graph representation, the spectral graph wavelets for finite element mesh can be constructed. We refer to the original article for the detailed theoretical background as well as the construction steps of the spectral graph wavelet . Indeed, the wavelet transforms are generated from two functions , namely, the scaling function φ0 and the mother wavelet transform φ. For any scale t > 0, a wavelet transform φt can be derived from the mother wavelet transform. Intuitively, in the frequency domain, the scaling function φ0 acts as a low-pass filter, while the wavelet transform φt performs like a band-pass filter with the frequency band determined by the scale t. Each of these transforms can be represented as a N × N matrix. The i-th column of a wavelet transform φt is the wavelet with scale t centered at the i-th node of the graph. Some selected graph wavelets locate on the same center are plotted in Figure 1. In this paper, we employ the wavelet transform φ0 and the wavelet transforms φt with scales t = [0.5, 0.25, 0.125]. Figure 1

Some selected graph wavelets. These plots represent the wavelets located at the 500-th node with different scales.

### 2.3 Multi-scales wavelet transform

Multi-scales wavelet transforms are commonly used in image processing. They provide a multi-resolution view of the underlying problem. In canonical 2D image processing, a multi-scales wavelet transform is constructed through down-sampling the original domain and the wavelet transforms scale-by-scale. We adapt this strategy to the graph wavelet transforms developed in Section 2.2. Indeed, such multi-scales wavelet transform can be considered as a filter bank. Given a set of scales t ∈ [2−1, ⋯ 2l], this filter bank consists of wavelets selected from the wavelet transforms φt and φ0. We select all those N wavelets with the scale 2l and select around 2jN wavelets with scale 2l + j. Here, N is the number of nodes in mesh. This process is realized by downsampling the mesh nodes and forming a chain of node sets VlVl1V0. In this article, we fix l = 3. Regarding the mesh in Figure 1, the downsampled node sets are plotted in Figure 2. For j = 0, 1, 2, the wavelets of scale 2l + j with centers in Vj are selected. For j = 3, the wavelets of φ0 with centers in V0 are selected. Indeed, since the wavelets are identified as columns of matrices φt, we need only select the corresponding columns in wavelet transform φt for each scale. A multi-scales wavelet transform W will be obtained by setting its column as the selected wavelets scale-by-scale. Figure 2

The down-sampled node sets from the mesh used in Figure 1. With the decreasing of wavelet scales t, the denser node sets are used to select the wavelets from the corresponding φt.

### 2.4 Regularization with multi-scales wavelet transform

Combining the above tools, we introduce a sparse solver using the multi-scales wavelet transform as regularization matrix. Given the multi-scales wavelet transform W constructed in Section 2.3, the regularized EIT inverse problem (eq 3) has been rewritten as:

(4)Δs^=argminΔs12ΔVJNΔs22+αWTΔs1

Here we abused Δs to denote the conductivity changes on nodes. This optimization problem can be solved by Split Bregman iterative solver .

We summarize the proposed sparse solver with wavelet-based regularization (SSWR) by flowcharts in Figure 3. Figure 3

Workflow of the proposed framework.

## 3 Simulation results

The performance of this wavelet based regularization framework has been evaluated through simulation. The simulations were executed on a circular finite element domain with radius 1 (see. Figure 4A). On the boundary, 16 electrodes were attached equidistantly. The background conductivity is fixed to s0 = 1.0 Sm−1. Alternative currents of 10 mA were injected successively and adjacently through the electrodes, the voltages between adjacent electrodes were recorded and denoted by Vh. Two contrasts with conductivities 1.6 Sm−1 (blue contrast) and 0.6 Sm−1 (red contrast) are embedded into the background (see Figure 4B). After embedding, another voltage measurement (denoted by Vih) was simulated. The voltage difference was calculated by ΔV = VihVh . In addition, 10% white noise was added to ΔV. Figure 4

Simulation results. (A) shows the background of the simulated phantom. (B) shows the ground truth. Two contrasts are embedded into the background. (C) and (D) are the reconstructed images from the proposed SSWR solver and the one step GN solver with Tikhonov regularization respectively. The reconstructions are performed on a coarser mesh (C and D).

The reconstructions were applied on the coarser mesh shown in Figure 4C and D to avoid the so-called “inverse crime” problem. This coarser mesh is the same as the meshes used for Figures 1 and 2. Indeed, the multi-scale wavelet transform matrix was obtained as described in Section 2.3 using the down-sampled node sets depicted in Figure 2. In Split Bregman iterative solver, the number of iterations was fixed to be 5. The regularization parameter is set to be α = 10−5. In simulation, computation time for SSWR solver takes around 0.4 s to accomplish the reconstruction. Please insert your manuscript here.

The proposed SSWR solver has been compared with the standard one step Gauss-Newton method (GN) . GN was combined with standard Tikhonov regularizer which employs the identity matrix as regularization matrix. The regularization parameter used for Tikhonov regularization is set to be 5 ⋅ 10−3. The reconstructed images are plotted in Figure 4D. All the parameters are heuristically chosen for different solvers to get the best performance.

## 4 Discussion

The reconstructed images (see Figure 4C and D) demonstrated that the proposed wavelet-based regularization generates reconstruction with fewer artefacts. In particular, the image reconstructed by SSWR is more smooth and homogenous in contrast value. Moreover, the shape of the contrasts is better detected by the SSWR solver.

In this study, all parameters, such as the regularization parameters and the number of iteration for SSWR solver, are determined heuristically. A strategy to find the optimal parameter setting is preferred. Another shortcoming is the computational time. The proposed SSWR solver is an iterative solver. It is very efficient, but not enough for real-time imaging.

## 5 Conclusion

In this preliminary study, a multi-scales wavelet transform is constructed in the context of irregular finite element mesh. This is realized by applying the spectral graph wavelet transforms and down-sampling the nodes. Sparse regularization using this multi-scales wavelet transform was applied in simulation. Simulation results indicated that the proposed wavelet-based sparse solver obtains better images. In particular, the image reconstructed by the wavelet-based sparse solver avoids small turbulences and the shape of contrasts can be better recovered.

## Author’s Statement

Research funding: This work was supported by the BMBF grant no. 03FH038I3 (MOSES) and by the European Commission under grant no. FP7-PIRSES 318943 (eTime). Conflict of interest: Authors state no conflict of interest. Material and Methods: Informed consent: Informed consent is not applicable. Ethical approval: The conducted research is not related to either human or animal use.

## References

 Zhao Z, Fischer R, Frerichs I, Müller-Lisse U, Möller K. Regional ventilation in cystic fibrosis measured by electrical impedance tomography. J Cyst Fibros. 2012;11:412–8.10.1016/j.jcf.2012.03.011Search in Google Scholar PubMed

 Holder DS. Electrical impedance tomography: methods, history and applications. CRC Press; 2004.10.1201/9781420034462.ch4Search in Google Scholar

 Adler A, Lionheart WR. Uses and abuses of EIDORS: an extensible software base for EIT. Physiol Meas. 2006;27:S25.10.1088/0967-3334/27/5/S03Search in Google Scholar PubMed

 Brad Graham AA. A nodal jacobian algorithm for reduced complexity EIT reconstructions. Int J Inform Syst Sci. 2006;2:453–68.Search in Google Scholar

 Hammond DK, Vandergheynst P, Gribonval R. Wavelets on graphs via spectral graph theory. Appl Comput Harmon A. 2011;30:129–50.10.1016/j.acha.2010.04.005Search in Google Scholar

 Goldstein T, Osher S. The split Bregman method for L1-regularized problems. SIAM J Imaging Sci. 2009;2:323–43.10.1137/080725891Search in Google Scholar

Published Online: 2016-9-30
Published in Print: 2016-9-1 