Stability in Conductivity Imaging from Partial Measurements of One Interior Current

We prove a stability result in the hybrid inverse problem of recovering the electrical conductivity from partial knowledge of one current density field generated inside a body by an imposed boundary voltage. The region where interior data stably reconstructs the conductivity is well defined by a combination of the exact and perturbed data.


Introduction
Since the 1980s, there have been numerous attempts in Electrial Impedance Tomography (EIT) toward a robust reconstruction of cross-sectional images of the conductivity distribution inside the human body [14]. Cross-sectional imaging of electrical conductivity distribution reveal important physiological and pathological information [14]. Therefore, high quality images could provide better differentiation of tissue or organs, resulting in enhanced diagnosis and treatment of numerous diseases [10,11,35]. In EIT, an electric potential u is generated inside a body Ω while maintaining a voltage f at the boundary of Ω. In the absence of sources the potential solves the following Dirichlet problem ∇ · σ∇u = 0 in Ω, The Dirichlet to Neumann map, or voltage to current map, is given by Λ σ : f → (σ∂u/∂ν)| ∂Ω , where ν denotes the unit outer normal to ∂Ω. The inverse EIT problem is to recover σ from knowledge of Λ σ . This problem is mathematically called Calderón's Problem [8] which is known to be severely ill-posed away from the boundary [1,3,21,28,29]. Such instability is reflected in the fact that EIT imaging suffers of poor spatial resolution and has not reached the stage of clinical applications. In recent years, there has been a grown interest in medical imaging modalities that combine the best imaging properties of different physics (e.g., electrical potential, magnetic resonance , optical waves, pressure waves, etc.) to generate high contrast and high spatial resolution images. These emerging techniques are known as coupled-physics imaging, sometimes called either hybrid imaging or multi-wave imaging. In some applications of non-invasive medical imaging there is need for high contrast and high resolution images. For example, in early cancer detection, high contrast discriminates between healthy and non-healthy tissue whereas high resolution is important to identify anomalies at an early stage [35]. In some situation current methodologies (e.g., EIT, optical tomography, ultrasound, magnetic resonance imaging) focus only in a particular type of wave that can either recover high resolution or high contrast images, but not both with the required accuracy [5]. For instance, using only electrical measurements we can potentially obtained high contrast but very low resolution.
Current density impedance imaging (CDII) is a multi-wave imaging modality that overcomes the poor spatial resolution of EIT [27]. In CDII, electrical measurements are combined with magnetic resonance measurements. In a CDII experiment, electric current is injected into the body and an MRI machine is used to measure the z-component of the magnetic field B generated from the electric current via the Biot-Savart law. Using two rotations of the body we can obtain the additional x and y components of B. Once B is known, using Amepère's Law we get the internal current density of the body, J = −σ∇u. This additional information enables improvement in the spatial resolution and mitigates for the instability of the EIT problem. For a review on the physical principles we refer the reader to the nice review paper [33,Scherzer and Widlak].
The use of current density to image electrical conductivity goes back to 1991 [30,34]. Since then, there has been extensive mathematical study of this problem from different points of view [7, 13, 16-19, 22-27, 30, 31]. Several recent works [4,9,22] used multiple boundary voltages to obtain redundant internal information that enables better stability results. In [22] Bal and Monard were able to take advantage of these additional information to prove Lipschitz stability under the assumption that the gradients of the solutions have maximal rank in R n at every point in the domain.
We focus on the problem of recovering the conductivity from a single boundary voltage potential and only one-dimensional information of the current density J (e.g. one component of J). Notice that since we are interested in obtaining an isotropic conductivity σ (i.e., a one variable function), access to the whole J would make the inverse problem mathematically over determined in dimensions n ≥ 2. We focus on recovering the conductivity from minimal information about the current density, while doing so we reveal geometric characteristics intrinsic of the CDII inverse problem. We show that we only need a well chosen component of J to recover σ in an stable way in any dimension. We study the full and partial data case an we describe the region where the recovery of σ is stable. As mentioned before, to obtain the whole current density J two rotations of the body are necessary. This is an obstacle in practice because patients can only be rotate once inside an MRI machine. A consequence of our result -that only one component of the current density is necessary to recover the conductivityis a better understating of geometric condition that imply that σ can be recover from two well chosen components of B. In particular, this implies that under some geometric assumptions only one rotation of the body is required to stably recover σ. We hope that the ideas presented in this work could lead to alternative approaches to performed the measurements in CDII.
In [24,25,Nachman,Tamasan and Timonov] proved that if the potential has no critical points in the domain, then information about the magnitude of the current density |J| uniquely determines the conductivity inside the body. A general approach presented by Kuchment and Steinhauer [19], dealt with the linearization of the problem. They proved that the linearization is stable, up to a finite dimensional kernel, which does not necessarily implies injectivity. For dimension n = 2, conditional stability to recover the voltage potential u was established in [31] by Tamasan, Timonov and Veras. In a joint work with Plamen Stefanov [23], we proved local Hölder stability for the linear and non-linear problem in any dimensions, under the assumption of no critical points. Our approach relies on a factorization of the linearization to prove stability for the linear problem by composition of stable operators. We would like to point out a small omission in [23, Theorem 1.2], were we need to add the assumption that the voltage potentials must be equal at the boundary. This is a natural assumption and it was clearly used in Equation (21) in there, but unfortunately omitted in the final version of the paper. In [18,Kim and Lee] used a curl-based approach of the CDII problem to proved Hölder stability for the planar case. They assume knowledge of all the current density, not only the magnitude, which makes the inverse problem linear. Their inversion is based on the method of characteristics (see also, [20]). In a similar spirit (i.e., using the whole current density), but for the anisotropic case, Bal, Guo and Monard [6] proved logarithmic stability in dimension n = 2.
In the partial data case to the author's knowledge there is only the result of uniqueness in dimension n = 2 by Nachman, Tamasan and Timonov [26]. They proved uniqueness in what we called the injectivity region. We show uniqueness and stability in what we called the stability region for any dimensions. When the part of the boundary Γ (where we have information) is simply connected the injectivity region and the stability region coincide under reasonable assumptions. However, when Γ is not simply connected they might differ as illustrated in Figure 2. It is unknown to the author if there is uniqueness for dimensions n ≥ 3, in injectivity region.
We prove global injectivity and stability of recovering the conductivity from the information about the voltage potential in a finite union of open connected sets and one component of the current density in a subset of the domain. We analyze the full data and partial data case. We use a non-linear version of the decomposition obtained in [23] and discover a "natural" component needed from the current density to obtain global stability.To the author's knowledge this is the first result where either injectivity or stability is establish from knowledge of only one component of the current density. As expected, the required component depends on the voltage potential inside the domain. However, if we assume an apriori condition that σ andσ are close in C 2,α (0 < α < 1), then we the stability region and the required component the current density could be "better controlled" from the boundary, as described in Theorem 1.2.
We work under the assumption of no critical points of the potential in the domain. In dimension two, the condition of no critical points can be satisfied under the assumption that the boundary illumination is almost two-to-one [24]. In dimensions n ≥ 3, it is still unclear how to chose boundary voltage potential to obtain no critical points in the interior. There are some recent works (e.g., [15], [2]) that are able to show global stability recoveries besides having a similar critical condition. The main obstacle in our analysis to include critical points is a better geometric understanding of the equipotential surfaces in dimensions n ≥ 3.
If σ and u describe, respectively, the electrical conductivity and electric potential of the body Ω, then the electrostatic approximation of Maxwell's equations are given by (2) and (3) were f is the voltage at the boundary. In this case, the current density vector field J : Knowledge of the current density is obtained from the magnetic field B (induced by the current injection). The current density and the magnetic field are related by the Ampère's Law where µ 0 ≈ 1, 3 × 10 −6 H/m is the vacuum magnetic permeability. For more details we refer to Scherzer and Widlak [33]. We address the stability of the inverse problem of recovering the conductivity σ from knowledge of f in an open subset Γ of ∂Ω, and from one component of J. In the case when Γ = ∂Ω and the internal information is given in all Ω we called it full data case. If Γ ∂Ω and the internal functional is known in some part of Ω (that in principle might depend on Γ) we called if partial data case. We give some preliminary notation that is useful to enunciate the results. Let u ∈ C 1 (Ω) be a function with non-vanishing gradient in Ω. For each p ∈ Ω we denote by x p the integral curve of u starting at p in the direction of ∇u. Assume |∇u| > 0 in Ω, then there exist a first time, denoted by t + p (resp. t − p ), such that the integral curve, starting from p and traveling in the direction (resp. opposite direction) to ∇u, hits the boundary transversely. We denote by the segment of the integral curve from t = 0 to t = t ± p . Also, denote by be the level curve of u = u(p) in Ω. We illustrated the previous notation in Figure 1. Figure 1: This picture illustrates the segment lp and the surface Σp.
The following definition describes the interior set in Ω were we can stably recover the conductivity σ from partial boundary information in Γ ⊂ ∂Ω and the set were we need information about the current density. Injectivity region: We say that p is visible from Γ, in transversal direction to ∇u, if The set of points in Ω that are visible from Γ in transversal direction to ∇u is called the injectivity region and denoted by I(Γ, u), i.e., Stability region: We say that the trajectory of p is visible from Γ, in parallel direction The set of point in Ω whose trayectory is visible from Γ in parallel direction to ∇u is called the stability region and denoted by S(Γ, u), i.e.,  Figure 2 illustrated the injectivity and stability regions from Γ, in the case that Γ is non connected. Clearly from definitions S(Γ, u) ⊂ I(Γ, u), when Γ is connected they can be equal I(Γ, u) = S(Γ, u) as shown in Figure 3.
The injectivity region I(Γ, u) from Γ is the set of all points that are visible from Γ in the transversal direction to ∇u; and the the stability region S(Γ, u) from Γ is the set of all points that are visible from Γ in transversal and parrallel direction to ∇u and for which one of the segments l ± p is contained in the injectivity region I(Γ, u +ũ). Here we called I(Γ, u +ũ) "injectivity region" because in dimension n = 2 is the region were σ can be uniquely determined from Γ and the magnitude of the current density as shown in [26]. To the author's knowledge there is still an open question whether in dimension n > 2 there is uniqueness in the region I(Γ, u +ũ) with some apriori information of the current density. We introduce the following notation to simplify the presentation.
Notation 1. For two vector fields w 0 and w in Ω, we denote by Π w 0 the orthogonal projection of w onto w 0 and Π ⊥ w 0 = Id − Π w 0 the orthogonal projection onto the orthogonal complement of w 0 , both in the Euclidean metric, i.e., In [26, Nachman, Tamasan and Timonov] proved in dimension n = 2, that if |J(σ)|− |J(σ)| = 0 in the injectivity region I(Γ, u +ũ) (in this case, I(Γ, u) = I(Γ, u +ũ)) then δσ = σ −σ = 0. The main theorem states is that in the stability region we can recover stably the conductivity from a "well chosen" component of the current density. Specifically, we can recover δσ from the gradient of the orthogonal projection of δJ onto the gradient field of ∇(u +ũ), i.e., This in particularly implies that by loss on one derivative, we can stably recover δσ from the scalar component of δJ in the direction of ∇(u +ũ), given by δJ · ∇(u +ũ) |∇(u +ũ)| as shown in Figure 4. To the author's knowledge this is the first stability result for the problem of CDII using only partial data, as well as the first time the recovery is done from one component of the current density. Theorem 1.1 (Global Stability with Partial Data). Let σ,σ ∈ C 2,α (Ω) be positive functions in Ω and let u ∈ C 2 (Ω) be σ-harmonic andũ ∈ C 2 (Ω) beσ-harmonic with u| ∂Ω ,ũ| ∂Ω ∈ C 2,α (∂Ω), 0 < α < 1. Let Γ and Γ be open sets of ∂Ω such that Γ is the union of finitely many connected component and Γ is compactly contained in Γ. Assume σ| Γ =σ| Γ , and u| Γ =ũ| Γ , and |∇(u +ũ)| = 0 in I(Γ , u +ũ). Then there exist C > 0, such that In particular, Remarks.

From estimate (11) it follows
Notice that the r.h.s. of estimate (12) is independent of u +ũ, but it requires knowledge of the whole current density, i.e., n-dimensional information.
The stability estimate requires knowledge of J in a direction that cannot be controlled from the boundary illumination. However, if we assume an apriori condition that σ is close to σ 0 in C 2 (Ω). The level curves of u +ũ =const. will be close to the level curves u 0 =const. We could project δJ 0 = J(σ) − J(σ 0 ) onto the gradient field ∇u 0 instead of ∇(u + u 0 ) with a penalty term that will depend on the apriori closeness assumption. Hence we can work with an stability region and a projection of J that could be better controlled from the boundary. What we mean with that, is that we can choose voltages f as to induced level set of u 0 , by doing so, we would have a better understanding of the required direction of the current density J, as illustrated in Figure  5. We have the following local stability result that proves what we described above. Theorem 1.2. Let σ, σ 0 ∈ C 2,α (Ω) be positive functions in Ω. Let u ∈ C 2 (Ω) (resp. u 0 ∈ C 2 (Ω)) be σ-harmonic (resp. σ 0 -harmonic), with u| ∂Ω , u 0 | ∂Ω ∈ C 2,α (∂Ω) for some and |∇(u + u 0 )|, |∇(u 0 )| are non zero in Ω. There exist and > 0 such that if then for some C > 0, we have In particular,

Figure 5:
We illustrate how we can controlled the visible region and the projection of the current density. The underlying idea is to chose voltages potential f at the boundary to induced specific level curves u0 =const. By doing so, we get a better understanding of the required direction of the current density to obtain an stable reconstruction. Here, δJ0 denote J(σ) − J(σ0) and w0 = Π ∇(u 0 ) (δJ0).
A consequence of Theorem 1.2 in dimension n = 3 is when the integral curves of u 0 are lines (i.e., equipotential surfaces are planes) we only need two components of the magnetic field B = (B x , B y , B z ) to stably recover the conductivity. Consider the following example.
Example 1. Assume that σ 0 = 1 and f = z. Then u 0 = z and by from Ampère's law Hence from Theorem 1.2 we can recover δσ 0 from δB x and δB y .
In particular, Theorem 1.2 gives some insight about necessary geometric conditions for the background equipotential surfaces needed to recover σ from only two components of B in an stable way.

Acknowledgments
I would like to thank Prof. Alexandru Tamasan for inviting me to talk at University of Central Florida, and for the fruitful conversations about the CDII inverse problem following that visit. In particular, for his important remarks on the assumption of no critical points for the voltage potential in dimensions n ≥ 3.

Preliminaries
The following decomposition of J(σ) − J(σ) reduced the problem of stability in the reconstruction of the conductivity σ to stability of recovering the voltage potential u and stability of the differential operator L given by Equation (16).
where L is a differential operator given by Proof. Consider the difference of two internal measurements If we dot product with ∇(u +ũ) we obtain 2(J(σ) − J(σ)) · ∇(u +ũ) = (σ −σ)|∇(u +ũ)| 2 + (σ +σ)∇(u −ũ) · ∇(u +ũ). (18) Solving for σ −σ we have On the other hand it follows easily from (2) Using (19) we substitute σ −σ in (20) and we get We notice that this decomposition is a non-linear version of the one obtained in [23] for the full data case. From there, we have the following representation of the operator L.
Proof. Clearly u +ũ trivially solves the eikonal equation c 2 |∇φ| 2 = 1 for the speed c = |∇(u +ũ)| −1 . Near x 0 , we can assume that u +ũ(x 0 ) = a; then u +ũ(x) is the signed distance from x to the level surface u +ũ = a. Choose local coordinates y on this level curve, and set y n = u 0 (x). Then y = (y , y n ) are boundary local coordinates to u +ũ = a and in those coordinates, the metric c −2 dx 2 takes the form and dx 2 = c 2 (dy n ) 2 + g αβ dy α dy β .

Stability with Partial Data
The proof of stability with partial data is based on establishing stability estimates for the potential and the operator L in decomposition (15) as in [23]. The visibility of S(Γ , u +ũ) from Γ in parallel direction to ∇(u +ũ) allows to show stability for the potential, while the the visibility in the transversal direction is used to prove stability of the operator L. From Proposition 2.2 we have that the restriction of L to the level curves of u 1 +ũ 1 is elliptic, we use local boundary normal coordinates along the level curves of u 1 +ũ 1 and the parallel condition to extend the analysis to the boundary and obtain information from Γ .
To prove estimates (27) and (28), we will need to extend σ,σ, u,ũ in two different ways to a larger domain Ω 1 containing I ⊃ V . We need this extensions because we cannot guarantee that neither the level curves nor the integral curves of u +ũ are transversal to the boundary, hence is not clear how to use initial boundary condition when the surfaces are non-transversal. The first extension takes advantage of the boundary information while the second extension uses the regularity of the functions to extend the geometry and the operators in decomposition (15). This two extensions are similar in the proof of both estimates so we present them before the proof.

Proof of inequality
Denote by x(t + 0 ) = x + 0 the first point on that the integral curve, starting from x 0 and traveling in the same direction of ∇(u 2 +ũ 2 ), hits the boundary Ω transversely. Since x 0 belongs to the the visible from Γ, then x + 0 ∈ Γ. Also denote by x(t + 1 ) = x + 1 the first point on that the integral curve, starting from x 0 and traveling in the same direction of ∇(u 2 +ũ 2 ), hits the boundary Ω 1 . By possibly making Ω 1 closer to the boundary, we can assume that the the integral curve hits Ω 1 right after exiting Γ at x + 0 (i.e., the integral curve does not come back to Ω before hitting Ω 1 at x + 1 . We denote by b 0 a fixed number such that t + 0 < b 0 < t + 1 . Now, if we travel from x 0 in the direction opposite to ∇(u 2 +ũ 2 ) then either, the integral curve hits the boundary at x(t − 0 ) = x − 0 ∈ Γ and we can similarly defined x − 1 as the first point were the integral curve traveling in the opposite direction to the flow hits Ω 1 , or there exist t m such that x(t m ) = min V (u +ũ). In the first case we define denote a 0 such that t − 1 < a 0 < t − 0 , in the second case we let a 0 = t m .
Consider a tubular neighborhood of the integral curve x(t) by, where δ 0 > 0 is small enough so that T x 0 ∩ {y n = y n 0 + b 0 } is contained in Ω 1 \ Ω. Let h = (σ 1 −σ 1 ), notice that h satisfy the following first order equation Since in the extension of h to Ω 1 the differential commutes with the extension (30), we can extend (31) to get were recall that H is the characteristic function in Ω. In local coordinates in T x 0 , since h = 0 in Ω 1 \ ω, we can write (32) as were G = −(∇ · (σ +σ)∇(u −ũ)) and c = |∇(u +ũ)| −1 . Using an integral factor we can solve (33) and get h(y , y n ) = 1 µ(y , y n ) were µ(y , y n ) = exp y n b 0 c 2 ∆(u +ũ)(y , t)Hdt .
Notice that since ∇(u 2 +ũ 2 ) = 0 in V then 0 < C 1 < µ < C 2 for some positive constants C 1 and C 2 . Then using the Cauchy inequality we get that for δ 0 ≥ δ > 0, h(y) 2 L 2 (Tx 0 ) = |y −y 0 |<δ b a 1 µ(y , y n ) y n a c 2 G(y , t)Hµ(y , t)dt By the compactness of V , we can find T x 0 , T x 1 , . . . , T xm such that their union covers V except the n − 1 dimensional set {x ∈ V : (u +ũ)(x) = min V u +ũ}, which has measure zero. Since the right-hand side in (35) does not depends on this collection of tubular neighborhoods we use a partition of unity subordinated to this covering to conclude Proof of inequality (28): let x 0 ∈ I , and denote by Σ 0 the level surface of u 2 +ũ 2 in Ω 1 containing x 0 . Notice that Σ 0 is bounded and closed in Ω 1 , hence is compact in R 2 . Its restriction to the interior of Ω 1 is an open surface (locally given by u 2 +ũ 2 = const. with ∇(u 2 +ũ 2 ) = 0). Note that any such level surface may have points on Γ, where it is not transversal to Γ. Let y = (y , y n ) be local boundary normal coordinates for x 0 ∈ Σ 0 as in (21). By compactness and since x 0 belongs to the the visible part from Γ in transversal direction, we can define these coordinates to an open U 0 neighborhood of Σ 0 ∩ Ω contained in Ω 1 . In these coordinates we can write this open neighborhood as U 0 =Σ 0 × (a 0 − 0 , a 0 + 0 ), forΣ 0 = Σ 0 ∩Ω, where Ω Ω Ω 1 ; a 0 = (u +ũ)(x 0 ); and 0 < min{dist(∂Ω, ∂Ω), dist(∂Ω, ∂Ω 1 )}.
Denote by v = u 1 −ũ 1 . From (22), the local representation of L in U 0 we have (36) We now use Poincaré inequality onΣ 0 to obtain that for each x 0 ∈ I there exist 0 such that for all 0 < < 0 From (36) and (37) and since v = 0 in Ω 1 \ Ω we obtain using again compactness of I we have we implies (28).
Proof of Theorem 1.2. Denote by δJ 0 = J(σ) − J(σ 0 ). From Theorem 1.1 we have that there exist a C > independent of such that and notice from (18) the second term in the r.h.s. of (41) can be written as Using (20) and Cauchy-Schwartz inequality we obtain the following estimate for C independent of . Now since σ + σ 0 ∞ > δ > 0 in Ω we obtain Using Schauder theory for (20) we see easily that u − u 0 C 2 (Ω) = O( ). Combining this fact with estimates (42) and (43) Finally, using (45) for < 1/C we obtain with implies (14) and proves Theorem 1.2.

Conclusions
We have shown that one component of the current density is enough to stably recover the conductivity in the problem of CDII. We have the restriction of working under the assumption of no critical point of the voltage potential inside the domain. At the moment is unknown how to choose illumination to satisfy those condition for n ≥ 3. Moreover, is has been shown in similar coupled-physics inverse problems (e.g., coupling microwave and ultrasound [2], coupling elastrography and magnetic resonance [15]) that such degenerate condition can be avoided from the analysis. A better understanding of the equipotencial lines in dimension n ≥ 3 is necessary to carry over the analysis presented here to the case of critical points. For practical applications it would be interesting to explore the relation between the needed component of the current density and the magnetic field B. As shown in the introduction, under certain geometric restrictions we can recover σ from only two components of B. In such cases, only one rotation of the body is necessary to obtain the current density, making it feasible for human applications. The geometric properties of the equipotential surfaces could be approximated by means of proper boundary illumination. The stability of such approximation is also a matter to consider for practical purposes.
As mentioned in the introduction, Nachman, Tamasan and Timonov, proved in [26] for dimension n = 2 uniqueness in the in what we called the injectivity region. We show uniqueness and stability in what we called the stability region for any dimensions. When the part of the boundary Γ (where we have information) is simply connected the injectivity region and the stability region coincide (under reasonable geometric assumptions of the domain). However, when Γ is not simply connected they might defer as illustrated in Figure 2. It is unclear to the author if there is uniqueness in general in the injectivity region for dimensions n ≥ 3.
Finally, it might be possible to extend this analysis to the anisotropic case. For example in dimension n = 2, Bal, Guo and Monard, proved in [6] Lipschitz stability when the complete current density is known. Given that only an appropriate component of the current density is required to show Hölder stability in the isotropic case, minimal requirements about the internal functional might give Hölder stability for the anisotropic case as well.