# Anatomical Phase Extraction (APE) Method: A Novel Method to Correct Detrimental Effects of Tissue-Inhomogeneity in Referenceless MR Thermometry—Preliminary Ex Vivo Investigation

Authors: Chien-Feng Judith Huang; Win-Li Lin; San-Chao Hwang; Ching Yao; Hsu Chang; Li-Wei Kuo
Journal: Computational and Mathematical Methods in Medicine
 (2021)

---

## Abstract

Purpose. We present a novel background tissue phase removing method, called anatomical phase extraction (APE), and to investigate the accuracy of temperature estimation and capability of reducing background artifacts compared with the conventional referenceless methods. Methods. Susceptibility variance was acquired by subtracting pretreatment baseline images taken at different locations (nine pretreatment baselines are acquired and called φ1 to φ9). The susceptibility phase data φS was obtained using the Wiener deconvolution algorithm. The background phase data φT was isolated by subtracting φS from the whole phase data. Finally, φT was subtracted from the whole phase data before applying the referenceless method. As a proof of concept, the proposed APE method was performed on ex vivo pork tenderloin and compared with other two referenceless temperature estimation approaches, including reweighted ℓ1 referenceless (RW-ℓ1) and ℓ2 referenceless methods. The proposed APE method was performed with four different baselines combination, namely, (φ1, φ5, φ2, φ4), (φ3, φ5, φ2, φ6), (φ7, φ5, φ8, φ4), and (φ9, φ5, φ8, φ6), and called APE experiment 1 to 4, respectively. The multibaseline method was used as a standard reference. The mean absolute error (MAE) and two-sample t-test analysis in temperature estimation of three regions of interest (ROI) between the multibaseline method and the other three methods, i.e., APE, RW-ℓ1, and ℓ2, were calculated and compared. Results. Our results show that the mean temperature errors of the APE method-experiment 1, APE method-experiment 2, APE method-experiment 3, APE method-experiment 4, and RW-ℓ1 and ℓ2 referenceless method are 1.02°C, 1.04°C, 1.00°C, 1.00°C, 4.75°C, and 13.65°C, respectively. The MAEs of the RW-ℓ1 and ℓ2 referenceless methods were higher than that of APE method. The APE method showed no significant difference (p>0.05), compared with the multibaseline method. Conclusion. The present work demonstrates the use of the APE method on referenceless MR thermometry to improve the accuracy of temperature estimation during MRI guided high-intensity focused ultrasound for ablation treatment.


---

## Body

## 1. Introduction

Recently, high-intensity focused ultrasound (HIFU) has received a surge of attention due to its ability to perform ablation or hyperthermia therapy without an incision [1]. HIFU treatment is superior to conventional resection surgery for several reasons, including reduced cost, shorter recovery time, and greater patient treatment acceptance [2]. During HIFU treatment, monitoring the temperature response is important to ensure that the required thermal dose is delivered to the target region while sparing the surrounding healthy tissues [3]. Magnetic resonance imaging (MRI) is well suited for this purpose because of its great sensitivity of temperature measurement and attractive soft-tissue contrast to clearly reveal protein denaturation and lesion formation [4]. Thus, MRI guided HIFU (MRgHIFU) treatment is widely used in several kinds of diseases, including uterine fibroids [5, 6], brain tumors [7], and essential tremor [8].MR thermometry can be carried out by several kinds of tissue contrasts, such as proton density [9], T1 relaxation time [10], diffusion coefficient of water molecules [11], and proton resonance frequency shift (PRFS) [12, 13]. Of these, the PRFS method is the most widely used approach because of its linear behavior, ease of measurement, and near tissue-type independence [3, 14]. The PRFS method is based on the chemical shift of water protons [12, 15] and utilizes gradient-echo pulse sequences to acquire MR phase images. With this method, one pretreatment phase image is acquired as a baseline and subsequently subtracted from sequentially acquired phase images during heating to obtain a phase difference and compute the thermal mapping image.In spite of its success, the PRFS method operates under the assumption that temperature variation is the sole contributor to phase variation. However, in practice, there are multiple sources of phase variations beside temperature variation, such as field drift [16], patient motion [17], heterogeneous fat/water distribution [18, 19], cavitation [20], oxygen concentration changes [21], and ultrasound transducer movement [22]. Among these potential error sources, ultrasound transducer movement is the dominant source of bias during MRgHIFU treatment. These error sources can cause misinterpretation of phase variation as temperature elevation, creating a temperature bias. Correction of this temperature bias is essential to avoid ineffective treatment or unexpected denaturation of adjacent healthy tissues.Several strategies have been proposed to overcome this issue. A typical method is the multibaseline method, which acquires multiple baseline images to form a lookup table and obtains the corresponding image for baseline subtraction during heating, at the expense of a prolonged scan time [23, 24]. Instead of acquiring multiple baseline images, Rieke et al. proposed a ℓ2 referenceless method which acquires baseline images by applying a two-dimensional polynomial fitting of the nonheated region [25], and Grissom et al. introduced a reweighted ℓ1 regression to this approach to enhance its convenience (hereafter referred to as the RW-ℓ1 referenceless method) [26]. However, tissue-inhomogeneity may cause background artifacts in the temperature maps of both referenceless methods, leading to misjudgment of the focal zone position. These artifacts limit the use of the referenceless methods in clinical applications.In this study, we introduced an anatomical phase extraction (APE) method to the referenceless MR thermometry whereby the background tissue phase data was isolated and removed before applying polynomial fitting. The APE method was designed to eliminate the detrimental effects of tissue-inhomogeneity. The performance of APE method was evaluated in ex vivo pork tenderloin by using the multibaseline method as standard reference and comparing with the conventional referenceless methods, includingℓ2 and RW-ℓ1.


## 2. Materials and Methods

### 2.1. MRgHIFU System Setup

The MRgHIFU system design and its experimental setup are shown in Figures1(a) and 1(b). A spherical HIFU-transducer (focal length: 12 cm, aperture radius: 8 cm, and operating frequency: 1.2 MHz) was mounted on an arc structure attached to the MRI patient table. A water bag filled the space between the HIFU transducer and the object to be ablated. The ultrasound power attenuation through the water bag to the ablated object was less than 1.5%. The MR-compatible arc structure can house mechanical mechanisms and associated devices, providing positional control along two degrees of freedom in an in-plane coordinate system. System software was developed with C and JAVA programming language and supported the necessary functions to perform MRgHIFU ablation procedures such as treatment plan, treatment execution with power and positioning controls, workstation-to-MRI scanner communication, target localization, temperature measurement, and ablated tissue necrosis monitoring. Although the fiber optic thermocouple is a common ground truth used to validate MR thermometry, the HIFU transmission would lead to temperature bias of fiber optic thermocouple measurement surrounding the focal region. Therefore, in this study, we did not measure exact temperature changes using the thermocouple device. Since our purpose was to develop a referenceless method which can provide the similar performance as the conventional method, we chose the multibaseline method as the reference standard.Figure 1
(a) Schematics and (b) photograph of experimental setup. Spherical HIFU-transducer was mounted on an arc structure attached to the MRI patient table. A water bag filled the space between the HIFU transducer and the object to be ablated. (c)Ex vivo experiments with pork tenderloin: ultrasound transducer performed ablation in nine different positions. Upper row: the time series of HIFU ablation. Lower row: the representative temperature maps of the nine ablation positions with the corresponding time. Note that HIFU ablation power transmitted in nine different positions was varying and increasing (detailed in Table 1), and thus, the temperature of successive points becomes hotter and hotter.
(a)(b)(c)


### 2.2.Ex Vivo Experiments

Ex vivo pork tenderloin was used in this experiment. Nine HIFU ablations with varying power were conducted in different areas of the pork tenderloin (detailed spatial locations and power are given in Table 1) to simulate susceptibility variances caused by device repositioning in the clinical setting. The porcine specimen was obtained 2 hours before the ex vivo MR experiment and kept under the room temperature in the scan room. The positioning of experimental setup and anatomical localization scan took approximately 20 minutes, and the scan time for temperature mapping was approximately 18.5 minutes (1112 seconds). Thus, the total scan time for the whole experiment was approximately 38.5 minutes, and the porcine specimen was kept fresh during the experiment. HIFU energy was transmitted at positions one through nine at t=79,169,259,349,439,529,619,709,799 seconds, respectively, for a 30-second ablation and 60-second cooling time. The time series of HIFU ablation is also shown in Figure 1(c). For PRFS MR thermometry, the spoiled gradient echo sequence was implemented in a 1.5T MRI scanner (Symphony, Siemens, Erlangen, Germany) with the following parameters: repetitiontimeTR=37ms, echotimeTE=17.3ms, flipangle=18°, acquisition matrix size=128×77, fieldofview=256×256 mm, slicethickness=8 mm. The ablation points were preselected, and nine pretreatment baseline images were obtained for the conventional multibaseline thermometry method to acquire standard reference temperature images for comparison (detailed positions are described in Table 1 and relative positions shown in Figure 2(a)). The APE method has to choose four of these baselines to isolate the background tissue phase data. All images were processed and analyzed off-line using MATLAB (R2016b, The MathWorks, Inc, Natick, MA, USA).Table 1
Ex vivo HIFU ablation parameters.
LocationCoordinate (x, z)ShiftX (mm)ShiftZ (mm)HIFU powera (W)1(10, 10)1010402(0, 10)010603(-10, 10)1010804(10, 0)1001005(0, 0)001206(-10, 0)1001407(10, -10)10101608(0, -10)0101809(-10, -10)1010200aAblation duration: 30 seconds; cooling interval: 60 seconds.Figure 2
(a) Schematic diagram of relative positions one through nine. The distance between adjacent points is 10 mm. (b) Susceptibility phase signal images at the nine different positions acquired by the APE method. The contour of the central, or reference, positionφS,5 (in red) is superimposed on each susceptibility phase signal image to visualize the susceptibility variance between different positions. (c) Representative whole phase data (φ2), background tissue phase data (φT,2), and susceptibility phase data (φS,2). The tissue-inhomogeneity of φT,2 is apparent. When φT,2 is removed from φ2, the remaining φS,2 is smooth.
(a)(b)(c)


### 2.3. Principle

Letφj refer to the baseline phase image data when the moving objects (including the HIFU transducer, MRI RF coil, and water bag in our case) are at position j, we adopted three assumption for the APE method.Assumption 1.
During point-by-point ablation in MRgHIFU, only two kinds of phase signals exist in the baseline phase data: the background tissue phase data from the fixed object, hereafter calledφT,j, and the susceptibility phase data from moving objects (including the HIFU transducer, MRI RF coil, and water bag in our case), hereafter called φS,j. Thus, the baseline phase image data,φj, can be expressed as follows:(1)φj=φT,j+φS,j.Assumption 2.
When the moving objects are repositioned to different locations,φS is merely translated in the x and z directions, and the magnitude of φS does not change. As for φT, both magnitude and location of φT would not be affected by the moving objects repositioning (hence, φT is also referred to as the fixed phase data).Assumption 3.
The deleterious effects of tissue-inhomogeneity which would lead to background artifacts are contained withinφT.According to Assumption3, the background artifacts could be circumvented if φT can be removed before utilizing the referenceless methods. Thus, the core of the APE method is to isolate and remove φT,j from whole baseline phase data φj. To achieve this aim, we have to acquire φS,j first then φT,j can be derived from given φS,j and φj in Equation (1).


### 2.4. The Mathematical Model for the Image Degradation Process

To deriveφS,j from baseline images, a mathematical model for the image degradation process was used in the APE procedures. In our cases, susceptibility-shift process results in degradation of original (thermal) image. Image degradation process can be expressed as Equation (2) and detailed derivate in the appendix.(2)φS,ix,z∗Hi−cx,z=φS,ix,z−φS,cx,z,whereφS,ix,z and φS,cx,z denote the susceptibility phase data from moving objects when the moving objects are at position i and c, respectively. The distance between position i and c can be expressed as ∆x and ∆z (unit: pixel) along x and z directions, respectively. Hi−cx,z is point spread function (PSF) recording the devices shift effect between φi and φc(the matrix Hi−cx,z can also be considered as degrading operation), and the asterisk symbol (∗) indicates a two-dimensional spatial convolution. As Hi−cx,z and “φS,ix,z−φS,cx,z” are known (detailed in next section), Equation (2) could be cast as a deconvolution problem, and the Wiener deconvolution method can be used to acquire φS,ix,z.


### 2.5. Acquisition ofφS,ix,z

#### 2.5.1. PSF (Hi−c) Established

Hi−c is a Nx×Nz matrix recording the device shift effect between φi and φc. In the case of the transducer translated ShiftXi−cmm and ShiftZi−c (mm) in x direction and z direction, respectively, definition of Hi−c is(3)Hi−c=hkl,(4)hkl1ifk,l=0,0,−1ifk,l=∆x,∆z,0otherwise,(5)∆x=ShiftXi−cFOVx/Nx,∆z=ShiftZi−cFOVz/Nz,where the unit of∆x and ∆z is pixel.


#### 2.5.2. Acquisition of “φS,ix,z−φS,cx,z” Term

Recall thatφi and φc refer to the baseline phase image data when the moving objects are at position i and c, respectively. φi and φc can be expressed as(6)φix,z=φT,ix,z+φS,ix,z,φcx,z=φT,cx,z+φS,cx,z.As mentioned previously,φT is the fixed phase signal, so theoretically φT,i should equal φT,c and φT,ix,z−φT,cx,z=0. Thus, “φS,ix,z−φS,cx,z” term can be obtained by(7)φix,z−φcx,z=φT,ix,z+φS,ix,z−φT,cx,z+φS,cx,z=φT,ix,z−φT,cx,z+φS,ix,z−φS,cx,z=0+φS,ix,z−φS,cx,z=φS,ix,z−φS,cx,z.In other words, “φS,ix,z−φS,cx,z” term can be acquired by just subtracting two baseline images.


#### 2.5.3. Using Wiener Deconvolution Method to AcquireφS,i

Finally, we perform Wiener deconvolution [27] on Equation (2) to obtain φS,ix,z as follows:(8)φS,ix,z≈FT−1FTφS,ix,z−φS,cx,zFHi−c∙FTHi−c2FTHi−c2+K,whereFT and FT−1 mean the two-dimensional Fourier transform and inverse Fourier transform, respectively. K is a specified constant acting as a free parameter for optimizing the result.


### 2.6. Phase Unwrapping

All images underwent phase unwrapping prior to the APE method. We used Bruce Spottiswoode’s code from the MATLAB Central File Exchange [28], based on the algorithm proposed by Ghiglia et al. [29].


### 2.7. APE Procedures

The procedures of the APE method are displayed in a flowchart of Figure3 and are described below.(1)
Pretreatment phase images,φix,z,φcx,z, φmx,z, and φnx,z, were acquired where the moving objects were positioned in four different locations(2)
Withφix,z and φcx,z, φS,ix,z is derived according to Acquisition of φS,ix,z described above(3)
According to Assumption1, the background tissue phase φT,ix,z caused by the fixed object was isolated by subtracting φS,ix,z from φix,z(4)
To address the noises,φT,ix,z was subtracted from φix,z, φcx,z, φmx,z, and φn(x, z), respectively, to acquire φS,ix,z, φS,cx,z, φS,mx,z, and φS,nx,z. According to Assumption 2, we know that the magnitude of these four susceptibility phase data are equivalent, and there are merely translated in x and z directions. Thus, these four susceptibility phase data were all shifted to the location i and taking the mean as following equation to alleviate the noises(9)φS,ifinalx,z=14φS,ix,z+φS,cx−∆x,z−∆z+φS,mx−∆x,z+φS,nx,z−∆z.(5)
Finally, we subtractφS,ifinal from φi to isolate the fixed phase data, called φT,ifinalFigure 3
Flowchart of APE procedures.φi: baseline phase image when the moving objects are at position i (note that φi = φS,i + φT,i); matrix Hi−c: records the effect of φS,i−φS,c; φS,i−φS,c: equals the convolution of (φS,i, Hi−c); φS,i: susceptibility phase signal when the moving objects are at position i; φT,i: background tissue phase signal when the moving objects are at position i.


### 2.8. Conversion to Temperature

By introducing Assumption1 that MR phase data acquired during ablation, φjab, is comprised of two kinds of phase signals, φT,jab and φS,jab, into the PRFS equation proposed by Ishihara et al. [7], we calculated temperature change as follows:(10)φjab=φT,jab+φS,jab,ΔT=φjab−φjγαB0TE=φT,jab+φS,jab−φT,j+φS,jγαB0TE,whereγ is the gyromagnetic ratio,α(−0.01ppm/°C) is the temperature-dependent coefficient for tissue,B0 is the amplitude of the static magnetic field, φjab is the phase data acquired during ablation (comprised of φT,jab and φS,jab), and φj is the baseline phase data corresponding to φjab. As mentioned previously, φT is the fixed phase signal, so theoretically φT,jab should equal φT,j. Thus, the equation can be simplified as follows:(11)ΔT=φS,jab−φS,jγαB0TE.To acquire the baselineφS,j, we applied a fourth-order RW-ℓ1 referenceless method to φS,jab to remove phase variance in φS,jab due to temperature elevation. Then, temperature map can be obtained from Equation (11) with φS,jab and φS,j. Note that we applied a fourth-order RW-ℓ1 referenceless method to φS,jab to remove phase variance. This step would substantially reduce noise and computational errors generated by the Wiener deconvolution.


### 2.9. Validation of APE Method

To compare the performance of the APE method with other referenceless methods, we derived the temperature maps by using the APE method, multibaseline method, RW-ℓ1 [26], and ℓ2 referenceless method [25]. We used the code provided by Grissom et al. [22] for the RW-ℓ1 referenceless method and the APE method with fourth-order fitting in our data. Based on Rieke et al.’s recommendation [21], we also used fourth-order polynomial fitting in the ℓ2 referenceless method. Our APE method modified the RW-ℓ1 referenceless method by removing φT prior to polynomial fitting. We applied four different baseline combination sets to perform APE method. Let baselines which used in APE method as φi, φc, φm, and φn. The four baseline combinations sets are i,m,n,c=1,5,2,4, (3, 5, 2, 6), (7, 5, 8, 4), and (9, 5, 8, 6) used in APE experiment 1, APE experiment 2, APE experiment 3, and APE experiment 4, respectively. Since our purpose was to develop a referenceless method which can provide the similar performance as the conventional method, we chose the multibaseline method as the reference standard. The accuracy of each MR thermometry map was estimated using the temperature map obtained via the multibaseline method.To visualize the errors in temperature estimation, we subtracted the multibaseline method temperature map from the temperature maps obtained from the other three methods (APE and RW-ℓ1 and ℓ2 referenceless methods). Furthermore, to evaluate the temperature measurement errors of the three methods compared to the multibaseline method, we drew three rectangular regions-of-interest (ROI) and calculatedthe mean absolute temperature difference (or mean absolute error (MAE)) of each ROI. Two-sample t-test was performed to investigate the difference of the mean temperature values of ROIs between standard reference and each method. To examine the reproducibility of APE method, we also calculated the MAEs of each ROI of four APE experiments, respectively. Furthermore, we showed the representative temperature maps derived from these experiments and subtracted the temperature map obtained via experiment 1, experiment 2, experiment 3, experiment 4 from experiment 2, experiment 4, experiment 1, and experiment 3, respectively, to visualize the difference between the temperature maps acquired from these experiments.


## 2.1. MRgHIFU System Setup

The MRgHIFU system design and its experimental setup are shown in Figures1(a) and 1(b). A spherical HIFU-transducer (focal length: 12 cm, aperture radius: 8 cm, and operating frequency: 1.2 MHz) was mounted on an arc structure attached to the MRI patient table. A water bag filled the space between the HIFU transducer and the object to be ablated. The ultrasound power attenuation through the water bag to the ablated object was less than 1.5%. The MR-compatible arc structure can house mechanical mechanisms and associated devices, providing positional control along two degrees of freedom in an in-plane coordinate system. System software was developed with C and JAVA programming language and supported the necessary functions to perform MRgHIFU ablation procedures such as treatment plan, treatment execution with power and positioning controls, workstation-to-MRI scanner communication, target localization, temperature measurement, and ablated tissue necrosis monitoring. Although the fiber optic thermocouple is a common ground truth used to validate MR thermometry, the HIFU transmission would lead to temperature bias of fiber optic thermocouple measurement surrounding the focal region. Therefore, in this study, we did not measure exact temperature changes using the thermocouple device. Since our purpose was to develop a referenceless method which can provide the similar performance as the conventional method, we chose the multibaseline method as the reference standard.Figure 1
(a) Schematics and (b) photograph of experimental setup. Spherical HIFU-transducer was mounted on an arc structure attached to the MRI patient table. A water bag filled the space between the HIFU transducer and the object to be ablated. (c)Ex vivo experiments with pork tenderloin: ultrasound transducer performed ablation in nine different positions. Upper row: the time series of HIFU ablation. Lower row: the representative temperature maps of the nine ablation positions with the corresponding time. Note that HIFU ablation power transmitted in nine different positions was varying and increasing (detailed in Table 1), and thus, the temperature of successive points becomes hotter and hotter.
(a)(b)(c)


## 2.2.Ex Vivo Experiments

Ex vivo pork tenderloin was used in this experiment. Nine HIFU ablations with varying power were conducted in different areas of the pork tenderloin (detailed spatial locations and power are given in Table 1) to simulate susceptibility variances caused by device repositioning in the clinical setting. The porcine specimen was obtained 2 hours before the ex vivo MR experiment and kept under the room temperature in the scan room. The positioning of experimental setup and anatomical localization scan took approximately 20 minutes, and the scan time for temperature mapping was approximately 18.5 minutes (1112 seconds). Thus, the total scan time for the whole experiment was approximately 38.5 minutes, and the porcine specimen was kept fresh during the experiment. HIFU energy was transmitted at positions one through nine at t=79,169,259,349,439,529,619,709,799 seconds, respectively, for a 30-second ablation and 60-second cooling time. The time series of HIFU ablation is also shown in Figure 1(c). For PRFS MR thermometry, the spoiled gradient echo sequence was implemented in a 1.5T MRI scanner (Symphony, Siemens, Erlangen, Germany) with the following parameters: repetitiontimeTR=37ms, echotimeTE=17.3ms, flipangle=18°, acquisition matrix size=128×77, fieldofview=256×256 mm, slicethickness=8 mm. The ablation points were preselected, and nine pretreatment baseline images were obtained for the conventional multibaseline thermometry method to acquire standard reference temperature images for comparison (detailed positions are described in Table 1 and relative positions shown in Figure 2(a)). The APE method has to choose four of these baselines to isolate the background tissue phase data. All images were processed and analyzed off-line using MATLAB (R2016b, The MathWorks, Inc, Natick, MA, USA).Table 1
Ex vivo HIFU ablation parameters.
LocationCoordinate (x, z)ShiftX (mm)ShiftZ (mm)HIFU powera (W)1(10, 10)1010402(0, 10)010603(-10, 10)1010804(10, 0)1001005(0, 0)001206(-10, 0)1001407(10, -10)10101608(0, -10)0101809(-10, -10)1010200aAblation duration: 30 seconds; cooling interval: 60 seconds.Figure 2
(a) Schematic diagram of relative positions one through nine. The distance between adjacent points is 10 mm. (b) Susceptibility phase signal images at the nine different positions acquired by the APE method. The contour of the central, or reference, positionφS,5 (in red) is superimposed on each susceptibility phase signal image to visualize the susceptibility variance between different positions. (c) Representative whole phase data (φ2), background tissue phase data (φT,2), and susceptibility phase data (φS,2). The tissue-inhomogeneity of φT,2 is apparent. When φT,2 is removed from φ2, the remaining φS,2 is smooth.
(a)(b)(c)


## 2.3. Principle

Letφj refer to the baseline phase image data when the moving objects (including the HIFU transducer, MRI RF coil, and water bag in our case) are at position j, we adopted three assumption for the APE method.Assumption 1.
During point-by-point ablation in MRgHIFU, only two kinds of phase signals exist in the baseline phase data: the background tissue phase data from the fixed object, hereafter calledφT,j, and the susceptibility phase data from moving objects (including the HIFU transducer, MRI RF coil, and water bag in our case), hereafter called φS,j. Thus, the baseline phase image data,φj, can be expressed as follows:(1)φj=φT,j+φS,j.Assumption 2.
When the moving objects are repositioned to different locations,φS is merely translated in the x and z directions, and the magnitude of φS does not change. As for φT, both magnitude and location of φT would not be affected by the moving objects repositioning (hence, φT is also referred to as the fixed phase data).Assumption 3.
The deleterious effects of tissue-inhomogeneity which would lead to background artifacts are contained withinφT.According to Assumption3, the background artifacts could be circumvented if φT can be removed before utilizing the referenceless methods. Thus, the core of the APE method is to isolate and remove φT,j from whole baseline phase data φj. To achieve this aim, we have to acquire φS,j first then φT,j can be derived from given φS,j and φj in Equation (1).


## 2.4. The Mathematical Model for the Image Degradation Process

To deriveφS,j from baseline images, a mathematical model for the image degradation process was used in the APE procedures. In our cases, susceptibility-shift process results in degradation of original (thermal) image. Image degradation process can be expressed as Equation (2) and detailed derivate in the appendix.(2)φS,ix,z∗Hi−cx,z=φS,ix,z−φS,cx,z,whereφS,ix,z and φS,cx,z denote the susceptibility phase data from moving objects when the moving objects are at position i and c, respectively. The distance between position i and c can be expressed as ∆x and ∆z (unit: pixel) along x and z directions, respectively. Hi−cx,z is point spread function (PSF) recording the devices shift effect between φi and φc(the matrix Hi−cx,z can also be considered as degrading operation), and the asterisk symbol (∗) indicates a two-dimensional spatial convolution. As Hi−cx,z and “φS,ix,z−φS,cx,z” are known (detailed in next section), Equation (2) could be cast as a deconvolution problem, and the Wiener deconvolution method can be used to acquire φS,ix,z.


## 2.5. Acquisition ofφS,ix,z

### 2.5.1. PSF (Hi−c) Established

Hi−c is a Nx×Nz matrix recording the device shift effect between φi and φc. In the case of the transducer translated ShiftXi−cmm and ShiftZi−c (mm) in x direction and z direction, respectively, definition of Hi−c is(3)Hi−c=hkl,(4)hkl1ifk,l=0,0,−1ifk,l=∆x,∆z,0otherwise,(5)∆x=ShiftXi−cFOVx/Nx,∆z=ShiftZi−cFOVz/Nz,where the unit of∆x and ∆z is pixel.


### 2.5.2. Acquisition of “φS,ix,z−φS,cx,z” Term

Recall thatφi and φc refer to the baseline phase image data when the moving objects are at position i and c, respectively. φi and φc can be expressed as(6)φix,z=φT,ix,z+φS,ix,z,φcx,z=φT,cx,z+φS,cx,z.As mentioned previously,φT is the fixed phase signal, so theoretically φT,i should equal φT,c and φT,ix,z−φT,cx,z=0. Thus, “φS,ix,z−φS,cx,z” term can be obtained by(7)φix,z−φcx,z=φT,ix,z+φS,ix,z−φT,cx,z+φS,cx,z=φT,ix,z−φT,cx,z+φS,ix,z−φS,cx,z=0+φS,ix,z−φS,cx,z=φS,ix,z−φS,cx,z.In other words, “φS,ix,z−φS,cx,z” term can be acquired by just subtracting two baseline images.


### 2.5.3. Using Wiener Deconvolution Method to AcquireφS,i

Finally, we perform Wiener deconvolution [27] on Equation (2) to obtain φS,ix,z as follows:(8)φS,ix,z≈FT−1FTφS,ix,z−φS,cx,zFHi−c∙FTHi−c2FTHi−c2+K,whereFT and FT−1 mean the two-dimensional Fourier transform and inverse Fourier transform, respectively. K is a specified constant acting as a free parameter for optimizing the result.


## 2.5.1. PSF (Hi−c) Established

Hi−c is a Nx×Nz matrix recording the device shift effect between φi and φc. In the case of the transducer translated ShiftXi−cmm and ShiftZi−c (mm) in x direction and z direction, respectively, definition of Hi−c is(3)Hi−c=hkl,(4)hkl1ifk,l=0,0,−1ifk,l=∆x,∆z,0otherwise,(5)∆x=ShiftXi−cFOVx/Nx,∆z=ShiftZi−cFOVz/Nz,where the unit of∆x and ∆z is pixel.


## 2.5.2. Acquisition of “φS,ix,z−φS,cx,z” Term

Recall thatφi and φc refer to the baseline phase image data when the moving objects are at position i and c, respectively. φi and φc can be expressed as(6)φix,z=φT,ix,z+φS,ix,z,φcx,z=φT,cx,z+φS,cx,z.As mentioned previously,φT is the fixed phase signal, so theoretically φT,i should equal φT,c and φT,ix,z−φT,cx,z=0. Thus, “φS,ix,z−φS,cx,z” term can be obtained by(7)φix,z−φcx,z=φT,ix,z+φS,ix,z−φT,cx,z+φS,cx,z=φT,ix,z−φT,cx,z+φS,ix,z−φS,cx,z=0+φS,ix,z−φS,cx,z=φS,ix,z−φS,cx,z.In other words, “φS,ix,z−φS,cx,z” term can be acquired by just subtracting two baseline images.


## 2.5.3. Using Wiener Deconvolution Method to AcquireφS,i

Finally, we perform Wiener deconvolution [27] on Equation (2) to obtain φS,ix,z as follows:(8)φS,ix,z≈FT−1FTφS,ix,z−φS,cx,zFHi−c∙FTHi−c2FTHi−c2+K,whereFT and FT−1 mean the two-dimensional Fourier transform and inverse Fourier transform, respectively. K is a specified constant acting as a free parameter for optimizing the result.


## 2.6. Phase Unwrapping

All images underwent phase unwrapping prior to the APE method. We used Bruce Spottiswoode’s code from the MATLAB Central File Exchange [28], based on the algorithm proposed by Ghiglia et al. [29].


## 2.7. APE Procedures

The procedures of the APE method are displayed in a flowchart of Figure3 and are described below.(1)
Pretreatment phase images,φix,z,φcx,z, φmx,z, and φnx,z, were acquired where the moving objects were positioned in four different locations(2)
Withφix,z and φcx,z, φS,ix,z is derived according to Acquisition of φS,ix,z described above(3)
According to Assumption1, the background tissue phase φT,ix,z caused by the fixed object was isolated by subtracting φS,ix,z from φix,z(4)
To address the noises,φT,ix,z was subtracted from φix,z, φcx,z, φmx,z, and φn(x, z), respectively, to acquire φS,ix,z, φS,cx,z, φS,mx,z, and φS,nx,z. According to Assumption 2, we know that the magnitude of these four susceptibility phase data are equivalent, and there are merely translated in x and z directions. Thus, these four susceptibility phase data were all shifted to the location i and taking the mean as following equation to alleviate the noises(9)φS,ifinalx,z=14φS,ix,z+φS,cx−∆x,z−∆z+φS,mx−∆x,z+φS,nx,z−∆z.(5)
Finally, we subtractφS,ifinal from φi to isolate the fixed phase data, called φT,ifinalFigure 3
Flowchart of APE procedures.φi: baseline phase image when the moving objects are at position i (note that φi = φS,i + φT,i); matrix Hi−c: records the effect of φS,i−φS,c; φS,i−φS,c: equals the convolution of (φS,i, Hi−c); φS,i: susceptibility phase signal when the moving objects are at position i; φT,i: background tissue phase signal when the moving objects are at position i.


## 2.8. Conversion to Temperature

By introducing Assumption1 that MR phase data acquired during ablation, φjab, is comprised of two kinds of phase signals, φT,jab and φS,jab, into the PRFS equation proposed by Ishihara et al. [7], we calculated temperature change as follows:(10)φjab=φT,jab+φS,jab,ΔT=φjab−φjγαB0TE=φT,jab+φS,jab−φT,j+φS,jγαB0TE,whereγ is the gyromagnetic ratio,α(−0.01ppm/°C) is the temperature-dependent coefficient for tissue,B0 is the amplitude of the static magnetic field, φjab is the phase data acquired during ablation (comprised of φT,jab and φS,jab), and φj is the baseline phase data corresponding to φjab. As mentioned previously, φT is the fixed phase signal, so theoretically φT,jab should equal φT,j. Thus, the equation can be simplified as follows:(11)ΔT=φS,jab−φS,jγαB0TE.To acquire the baselineφS,j, we applied a fourth-order RW-ℓ1 referenceless method to φS,jab to remove phase variance in φS,jab due to temperature elevation. Then, temperature map can be obtained from Equation (11) with φS,jab and φS,j. Note that we applied a fourth-order RW-ℓ1 referenceless method to φS,jab to remove phase variance. This step would substantially reduce noise and computational errors generated by the Wiener deconvolution.


## 2.9. Validation of APE Method

To compare the performance of the APE method with other referenceless methods, we derived the temperature maps by using the APE method, multibaseline method, RW-ℓ1 [26], and ℓ2 referenceless method [25]. We used the code provided by Grissom et al. [22] for the RW-ℓ1 referenceless method and the APE method with fourth-order fitting in our data. Based on Rieke et al.’s recommendation [21], we also used fourth-order polynomial fitting in the ℓ2 referenceless method. Our APE method modified the RW-ℓ1 referenceless method by removing φT prior to polynomial fitting. We applied four different baseline combination sets to perform APE method. Let baselines which used in APE method as φi, φc, φm, and φn. The four baseline combinations sets are i,m,n,c=1,5,2,4, (3, 5, 2, 6), (7, 5, 8, 4), and (9, 5, 8, 6) used in APE experiment 1, APE experiment 2, APE experiment 3, and APE experiment 4, respectively. Since our purpose was to develop a referenceless method which can provide the similar performance as the conventional method, we chose the multibaseline method as the reference standard. The accuracy of each MR thermometry map was estimated using the temperature map obtained via the multibaseline method.To visualize the errors in temperature estimation, we subtracted the multibaseline method temperature map from the temperature maps obtained from the other three methods (APE and RW-ℓ1 and ℓ2 referenceless methods). Furthermore, to evaluate the temperature measurement errors of the three methods compared to the multibaseline method, we drew three rectangular regions-of-interest (ROI) and calculatedthe mean absolute temperature difference (or mean absolute error (MAE)) of each ROI. Two-sample t-test was performed to investigate the difference of the mean temperature values of ROIs between standard reference and each method. To examine the reproducibility of APE method, we also calculated the MAEs of each ROI of four APE experiments, respectively. Furthermore, we showed the representative temperature maps derived from these experiments and subtracted the temperature map obtained via experiment 1, experiment 2, experiment 3, experiment 4 from experiment 2, experiment 4, experiment 1, and experiment 3, respectively, to visualize the difference between the temperature maps acquired from these experiments.


## 3. Results and Discussion

### 3.1. Results

Susceptibility phase dataφS, which is caused by the moving objects at different positions and derived from the APE method outlined above, is shown in Figure 2(b). The φS shift, or susceptibility variance, can be visualized by comparing the contours of the central, or reference, position φS,5 (in red) with the contours of other positions φS,i (in blue). Although the variance seems nonsignificant, this susceptibility variance would nevertheless result in noticeable and clinically relevant errors in temperature measurement. Figure 2(c) shows the susceptibility phase data, φS,2, and the background tissue phase data, φT,2, which is separated from whole phase data φ2. The tissue-inhomogeneity contained in φ2 and φT,2 is obviously observed, whereas the remaining φS,2 becomes much smoother.The ROI surrounding the ablation points (ROI1 and ROI2) and the background (ROI3) are shown in Figure 4(a). Each ROI contains 12 pixels. The absolute errors in temperature measurement of the 12 pixels within each ROI were averaged to obtain the MAE of each ROI. A total of 230 frames were captured for each ROI, and the MAEs of these frames for each ROI of each method are plotted in Figures 4(b)–4(g). The mean temperature error was derived by averaging all the MAEs (Table 2). The root mean square errors (RMSEs) between multibaseline method and the other methods of each ROI are also shown in Table 2 to compare the performance of each method. The mean temperature errors of the APE method-experiment 1, APE method-experiment 2, APE method-experiment 3, APE method-experiment 4, and RW-ℓ1 and ℓ2 referenceless method are 1.02°C, 1.04°C, 1.00°C, 1.00°C, 4.75°C, and 13.65°C, respectively. The MAEs and RMSEs of the RW-ℓ1 and ℓ2 referenceless method were higher than APE experiments, especially in ROI3 where the proposed APE method achieved the lowest MAEs and RMSEs of mostly less than 1°C. As compared with multibaseline method, both RW-ℓ1 and ℓ2 referenceless methods showed significant difference (p<0.05). On the contrary, there was no significant difference between multibaseline method and each APE experiment (p>0.05).Figure 4
(a) A representative temperature map with the three selected ROIs. Each ROI contains 12 pixels, and the absolute errors in temperature measurement of the 12 pixels within each ROI were averaged to obtain the MAE. (b–g) The MAEs of a total of 230 frames among these three ROIs for the proposed APE #experiment 1 (b), proposed APE #experiment 2 (c), proposed APE #experiment 3 (e), proposed APE #experiment 4 (f), and RW-ℓ1 (d) and ℓ2referenceless (g) methods, respectively.Table 2
Mean values and standard deviations for absolute temperature estimation errors, RMSEs, and two-samplet-test analysis between multibaseline and proposed APE, RW−ℓ1 referenceless, and ℓ2 referenceless method, respectively.
ROI1ROI2ROI3p valueProposed APE (°C)#experiment 11.14±0.69 (RMSE=1.33)1.08±1.04 (RMSE=1.50)0.83±0.37 (RMSE=0.90)0.35#experiment 21.15±0.68 (RMSE=1.33)1.11±1.07 (RMSE=1.53)0.88±0.37 (RMSE=0.95)0.40#experiment 31.17±0.69 (RMSE=1.36)0.96±0.91 (RMSE=1.32)0.86±0.38 (RMSE=0.94)0.41#experiment 41.14±0.69 (RMSE=1.34)1.06±1.02 (RMSE=1.47)0.79±0.35 (RMSE=0.87)0.32RW-ℓ1 referenceless (°C)1.37±0.84 (RMSE=1.61)3.00±0.37 (RMSE=3.03)9.89±0.98 (RMSE=9.93)< 0.05ℓ2 referenceless (°C)2.76±0.78 (RMSE=2.86)1.98±1.49 (RMSE=2.48)36.20±1.34 (RMSE=36.22)< 0.05Figures5(a)–5(d) show the temperature maps derived from the multibaseline, proposed APE, RW-ℓ1, and ℓ2 referenceless thermometry methods, respectively. Figure 5(e) shows spatial profiles through the heated area of the respective temperature maps. The multibaseline and proposed APE methods demonstrated considerable similarity in both the temperature map and the spatial profile. On the contrary, the ℓ2 method resulted in substantial background artifacts and noticeable discrepancies in both the temperature map and the spatial profile when compared to the multibaseline method. Although the RW-ℓ1 method performed better than the ℓ2 method, background artifacts were still remarkable.Figure 5
Temperature maps derived from the (a) multibaseline, (b) proposed APE, (c) RW-ℓ1 referenceless, and (d) ℓ2referenceless method, respectively. (e) Comparison of spatial profiles through the heated area (dotted-white line) derived from (a) to (d).To visualize the error in temperature measurement of each method, Figure6 shows the representative difference in temperature maps acquired by subtracting the temperature map obtained via each method (APE, RW-ℓ1, and ℓ2 referenceless) from the multibaseline temperature map. The difference map of APE method exhibited negligible errors. In contrast, both RW-ℓ1 and ℓ2 referenceless methods resulted in noticeable background artifacts in the difference map, with the ℓ2 method demonstrating the largest background artifacts. The background artifacts of the RW-ℓ1 method were still notable, despite being smaller compared to the ℓ2 method.Figure 6
(a) Temperature maps derived from the multibaseline (standard reference) method. (b) Temperature maps derived from the proposed APE (top row), RW-ℓ1 referenceless (middle row), and ℓ2referenceless (bottom row) methods. (c) The difference between temperature maps of the multibaseline method and the other three methods. Background artifacts in the RW-ℓ1 referenceless and ℓ2 referenceless methods are strikingly obvious. The proposed APE method markedly reduced these background artifacts.To validate the reproducibility of APE method, Figure7 shows the representative temperature map of APE experiments 1 to 4. It is obvious that these temperature maps demonstrate considerable similarity, and the noises in their difference maps are negligible. Note that the color bar range of difference maps in Figure 7 has been narrowed to 1°C to 4°C because the difference values are too small to identify.Figure 7
The representative temperature map of APE experiments 1 to 4 and the difference map acquired from temperature map difference of experiment 1-experiment 2, experiment 2-experiment 3, experiment 3-experiment 4, and experiment 1-experiment 4. It is obvious that these temperature maps demonstrate considerable similarity, and the noises in their different maps are negligible.


### 3.2. Discussion

The APE method improves upon the referenceless thermometry methods by eliminating background artifacts and temperature measurement errors in the heated region caused by tissue-inhomogeneity. By removing the background tissue phase dataφT before implementing the referenceless methods, we reduced the adverse effects of tissue-inhomogeneity and obtained a temperature map comparable to that of the multibaseline method. Moreover, by requiring only four pretreatment baseline images instead of an entire baseline library required by the multibaseline method, our proposed APE method could potentially reduce scan time and improve patient comfortability.Theoretically, baselines of different methods should be compared with that of the multibaseline method to determine the accuracy of each method. However, our APE method acquires a different baseline compared to other thermometry methods. Generally, baselines comprise of bothφS and φT regardless of whether they were acquired via scanning or polynomial fitting. Our strategy removes φT to resolve the undesirable consequences of tissue-inhomogeneity; thus, our baseline and floating images do not contain φT. Therefore, comparing baselines to determine the accuracy of each thermometry method is not feasible in our case. Instead, to assess the accuracy of each method, we qualitatively compared the temperature maps obtained by each method and quantified the temperature measurement errors by calculating the MAEs.Repositioning of the moving objects results in shifts inφS or susceptibility variation. Susceptibility variation in point-by-point ablation (also referred to as sequential sonications) causes errors in temperature measurement in the PRFS method. To our knowledge, φTand φS always coexist in phase data and are difficult to distinguish. In this study, we developed a method to separate these two phase sources and remove φT from whole phase data, thereby visualizing φS. Consistent with Assumption 2, our results indicate that when the moving objects are repositioned to different locations, the magnitude of φS does not change; rather, φS is merely translated in the x and z directions.Although both RW-ℓ1 and ℓ2 referenceless methods suffer from negative effects of tissue-inhomogeneity, tissue-inhomogeneity has a smaller influence on the RW-ℓ1 method with regard to the area of the background artifact. The RW-ℓ1 method optimizes the ℓ2 referenceless method by automatically excluding the heated pixels and noise from the polynomial fitting procedure. This may have partially alleviated the effects of tissue-inhomogeneity by unintentionally excluding some, but not all, pixels affected by tissue-inhomogeneity. Despite the improvements of the RW-ℓ1 method, background artifacts and temperature measurement errors are still nonnegligible, limiting its clinical applications. Thus, φT isolation and removal are still necessary.Theoretically, our APE algorithm can separateφS and φT using any two baselines. However, we found that φS would slightly deform along the shift direction, resulting in small inconsistencies in φT, which should be fixed. To address this problem and diminish noise, we take the mean of φS,ix,z, φS,cx−∆x,z−∆z, φS,mx−∆x,z, and φS,nx,z−∆z to acquire φS,ifinalx,z (as mentioned in APE Procedures’ step 4) and then subtract φS,ifinalx,z from φix,z to acquire the fixed phase data φT,ifinal.Our study has some limitations. First, the proposed APE method was developed for two-dimensional susceptibility movement, whereas it is not suitable for susceptibility shifts in three dimensions or cylindrical directions. Thus, it is not viable in three-dimensional or cylindrical positioning devices. Currently, the transducers are typically repositioned along two dimensions (e.g.,x-z directions) in most HIFU ablation treatments, so the proposed APE method is applicable for most of clinical use. In further study, it would be beneficial to improve the proposed APE method to a three-dimensional algorithm for providing more clinical feasibility. Second, the shift of moving objects should be detected when the APE method is performed, but there is not a well-established method to detect these shifts to the best of our knowledge. Thus, it is necessary to develop an approach to record the moving objects shifts in MRgHIFU treatment for strengthening its further clinical use.This study described an APE method which reduces the effect of tissue-inhomogeneity on referenceless methods and improve the accuracy of PRFS-based temperature measurement. This strategy expands the potential clinical applications of the referenceless methods to beyond homogenous tissues. With our APE method, the referenceless methods can potentially be used in clinical MRgHIFU treatment. Compared to the conventional multibaseline method, the APE method demonstrates comparable accuracy and not only avoids unnecessary risk of treatment interruptions arising from patient movements but also decreases preparation time spent obtaining baseline images. Further studies applying the APE method toin vivo experiments are necessary to validate its feasibility.


## 3.1. Results

Susceptibility phase dataφS, which is caused by the moving objects at different positions and derived from the APE method outlined above, is shown in Figure 2(b). The φS shift, or susceptibility variance, can be visualized by comparing the contours of the central, or reference, position φS,5 (in red) with the contours of other positions φS,i (in blue). Although the variance seems nonsignificant, this susceptibility variance would nevertheless result in noticeable and clinically relevant errors in temperature measurement. Figure 2(c) shows the susceptibility phase data, φS,2, and the background tissue phase data, φT,2, which is separated from whole phase data φ2. The tissue-inhomogeneity contained in φ2 and φT,2 is obviously observed, whereas the remaining φS,2 becomes much smoother.The ROI surrounding the ablation points (ROI1 and ROI2) and the background (ROI3) are shown in Figure 4(a). Each ROI contains 12 pixels. The absolute errors in temperature measurement of the 12 pixels within each ROI were averaged to obtain the MAE of each ROI. A total of 230 frames were captured for each ROI, and the MAEs of these frames for each ROI of each method are plotted in Figures 4(b)–4(g). The mean temperature error was derived by averaging all the MAEs (Table 2). The root mean square errors (RMSEs) between multibaseline method and the other methods of each ROI are also shown in Table 2 to compare the performance of each method. The mean temperature errors of the APE method-experiment 1, APE method-experiment 2, APE method-experiment 3, APE method-experiment 4, and RW-ℓ1 and ℓ2 referenceless method are 1.02°C, 1.04°C, 1.00°C, 1.00°C, 4.75°C, and 13.65°C, respectively. The MAEs and RMSEs of the RW-ℓ1 and ℓ2 referenceless method were higher than APE experiments, especially in ROI3 where the proposed APE method achieved the lowest MAEs and RMSEs of mostly less than 1°C. As compared with multibaseline method, both RW-ℓ1 and ℓ2 referenceless methods showed significant difference (p<0.05). On the contrary, there was no significant difference between multibaseline method and each APE experiment (p>0.05).Figure 4
(a) A representative temperature map with the three selected ROIs. Each ROI contains 12 pixels, and the absolute errors in temperature measurement of the 12 pixels within each ROI were averaged to obtain the MAE. (b–g) The MAEs of a total of 230 frames among these three ROIs for the proposed APE #experiment 1 (b), proposed APE #experiment 2 (c), proposed APE #experiment 3 (e), proposed APE #experiment 4 (f), and RW-ℓ1 (d) and ℓ2referenceless (g) methods, respectively.Table 2
Mean values and standard deviations for absolute temperature estimation errors, RMSEs, and two-samplet-test analysis between multibaseline and proposed APE, RW−ℓ1 referenceless, and ℓ2 referenceless method, respectively.
ROI1ROI2ROI3p valueProposed APE (°C)#experiment 11.14±0.69 (RMSE=1.33)1.08±1.04 (RMSE=1.50)0.83±0.37 (RMSE=0.90)0.35#experiment 21.15±0.68 (RMSE=1.33)1.11±1.07 (RMSE=1.53)0.88±0.37 (RMSE=0.95)0.40#experiment 31.17±0.69 (RMSE=1.36)0.96±0.91 (RMSE=1.32)0.86±0.38 (RMSE=0.94)0.41#experiment 41.14±0.69 (RMSE=1.34)1.06±1.02 (RMSE=1.47)0.79±0.35 (RMSE=0.87)0.32RW-ℓ1 referenceless (°C)1.37±0.84 (RMSE=1.61)3.00±0.37 (RMSE=3.03)9.89±0.98 (RMSE=9.93)< 0.05ℓ2 referenceless (°C)2.76±0.78 (RMSE=2.86)1.98±1.49 (RMSE=2.48)36.20±1.34 (RMSE=36.22)< 0.05Figures5(a)–5(d) show the temperature maps derived from the multibaseline, proposed APE, RW-ℓ1, and ℓ2 referenceless thermometry methods, respectively. Figure 5(e) shows spatial profiles through the heated area of the respective temperature maps. The multibaseline and proposed APE methods demonstrated considerable similarity in both the temperature map and the spatial profile. On the contrary, the ℓ2 method resulted in substantial background artifacts and noticeable discrepancies in both the temperature map and the spatial profile when compared to the multibaseline method. Although the RW-ℓ1 method performed better than the ℓ2 method, background artifacts were still remarkable.Figure 5
Temperature maps derived from the (a) multibaseline, (b) proposed APE, (c) RW-ℓ1 referenceless, and (d) ℓ2referenceless method, respectively. (e) Comparison of spatial profiles through the heated area (dotted-white line) derived from (a) to (d).To visualize the error in temperature measurement of each method, Figure6 shows the representative difference in temperature maps acquired by subtracting the temperature map obtained via each method (APE, RW-ℓ1, and ℓ2 referenceless) from the multibaseline temperature map. The difference map of APE method exhibited negligible errors. In contrast, both RW-ℓ1 and ℓ2 referenceless methods resulted in noticeable background artifacts in the difference map, with the ℓ2 method demonstrating the largest background artifacts. The background artifacts of the RW-ℓ1 method were still notable, despite being smaller compared to the ℓ2 method.Figure 6
(a) Temperature maps derived from the multibaseline (standard reference) method. (b) Temperature maps derived from the proposed APE (top row), RW-ℓ1 referenceless (middle row), and ℓ2referenceless (bottom row) methods. (c) The difference between temperature maps of the multibaseline method and the other three methods. Background artifacts in the RW-ℓ1 referenceless and ℓ2 referenceless methods are strikingly obvious. The proposed APE method markedly reduced these background artifacts.To validate the reproducibility of APE method, Figure7 shows the representative temperature map of APE experiments 1 to 4. It is obvious that these temperature maps demonstrate considerable similarity, and the noises in their difference maps are negligible. Note that the color bar range of difference maps in Figure 7 has been narrowed to 1°C to 4°C because the difference values are too small to identify.Figure 7
The representative temperature map of APE experiments 1 to 4 and the difference map acquired from temperature map difference of experiment 1-experiment 2, experiment 2-experiment 3, experiment 3-experiment 4, and experiment 1-experiment 4. It is obvious that these temperature maps demonstrate considerable similarity, and the noises in their different maps are negligible.


## 3.2. Discussion

The APE method improves upon the referenceless thermometry methods by eliminating background artifacts and temperature measurement errors in the heated region caused by tissue-inhomogeneity. By removing the background tissue phase dataφT before implementing the referenceless methods, we reduced the adverse effects of tissue-inhomogeneity and obtained a temperature map comparable to that of the multibaseline method. Moreover, by requiring only four pretreatment baseline images instead of an entire baseline library required by the multibaseline method, our proposed APE method could potentially reduce scan time and improve patient comfortability.Theoretically, baselines of different methods should be compared with that of the multibaseline method to determine the accuracy of each method. However, our APE method acquires a different baseline compared to other thermometry methods. Generally, baselines comprise of bothφS and φT regardless of whether they were acquired via scanning or polynomial fitting. Our strategy removes φT to resolve the undesirable consequences of tissue-inhomogeneity; thus, our baseline and floating images do not contain φT. Therefore, comparing baselines to determine the accuracy of each thermometry method is not feasible in our case. Instead, to assess the accuracy of each method, we qualitatively compared the temperature maps obtained by each method and quantified the temperature measurement errors by calculating the MAEs.Repositioning of the moving objects results in shifts inφS or susceptibility variation. Susceptibility variation in point-by-point ablation (also referred to as sequential sonications) causes errors in temperature measurement in the PRFS method. To our knowledge, φTand φS always coexist in phase data and are difficult to distinguish. In this study, we developed a method to separate these two phase sources and remove φT from whole phase data, thereby visualizing φS. Consistent with Assumption 2, our results indicate that when the moving objects are repositioned to different locations, the magnitude of φS does not change; rather, φS is merely translated in the x and z directions.Although both RW-ℓ1 and ℓ2 referenceless methods suffer from negative effects of tissue-inhomogeneity, tissue-inhomogeneity has a smaller influence on the RW-ℓ1 method with regard to the area of the background artifact. The RW-ℓ1 method optimizes the ℓ2 referenceless method by automatically excluding the heated pixels and noise from the polynomial fitting procedure. This may have partially alleviated the effects of tissue-inhomogeneity by unintentionally excluding some, but not all, pixels affected by tissue-inhomogeneity. Despite the improvements of the RW-ℓ1 method, background artifacts and temperature measurement errors are still nonnegligible, limiting its clinical applications. Thus, φT isolation and removal are still necessary.Theoretically, our APE algorithm can separateφS and φT using any two baselines. However, we found that φS would slightly deform along the shift direction, resulting in small inconsistencies in φT, which should be fixed. To address this problem and diminish noise, we take the mean of φS,ix,z, φS,cx−∆x,z−∆z, φS,mx−∆x,z, and φS,nx,z−∆z to acquire φS,ifinalx,z (as mentioned in APE Procedures’ step 4) and then subtract φS,ifinalx,z from φix,z to acquire the fixed phase data φT,ifinal.Our study has some limitations. First, the proposed APE method was developed for two-dimensional susceptibility movement, whereas it is not suitable for susceptibility shifts in three dimensions or cylindrical directions. Thus, it is not viable in three-dimensional or cylindrical positioning devices. Currently, the transducers are typically repositioned along two dimensions (e.g.,x-z directions) in most HIFU ablation treatments, so the proposed APE method is applicable for most of clinical use. In further study, it would be beneficial to improve the proposed APE method to a three-dimensional algorithm for providing more clinical feasibility. Second, the shift of moving objects should be detected when the APE method is performed, but there is not a well-established method to detect these shifts to the best of our knowledge. Thus, it is necessary to develop an approach to record the moving objects shifts in MRgHIFU treatment for strengthening its further clinical use.This study described an APE method which reduces the effect of tissue-inhomogeneity on referenceless methods and improve the accuracy of PRFS-based temperature measurement. This strategy expands the potential clinical applications of the referenceless methods to beyond homogenous tissues. With our APE method, the referenceless methods can potentially be used in clinical MRgHIFU treatment. Compared to the conventional multibaseline method, the APE method demonstrates comparable accuracy and not only avoids unnecessary risk of treatment interruptions arising from patient movements but also decreases preparation time spent obtaining baseline images. Further studies applying the APE method toin vivo experiments are necessary to validate its feasibility.


## 4. Conclusion

In this work, we proposed the APE method, a strategy to eliminate the detrimental consequences of tissue-inhomogeneity in referenceless methods. The theoretical aspects of the strategy were described, and an ex vivo experiment was performed for qualitative and quantitative evaluation. Our results indicate that the accuracy of the APE method is comparable with the conventional multibaseline method, implying that the APE method could potentially reduce the total scan time, thereby improving patient comfortability and further reducing the risk of treatment interruption.



---

*Source: 5566775-2021-08-11.xml*