Ultrasound-guided navigation of a magnetic microrobot using acoustic phase analysis

Microrobots (MRs) have attracted significant interest for their potentialities in diagnosis and non-invasive intervention in hard-to-reach body areas. Fine control of biomedical MRs requires real-time feedback on their position and configuration. Ultrasound (US) imaging stands as a mature and advantageous technology for MRs tracking, but it suffers from disturbances due to low contrast resolution. To overcome these limitations and make US imaging suitable for closed-loop MR control, we propose a US contrast enhancement mechanism for MR visualization in heterogeneous and dynamic backgrounds (e.g., tissue). Our technique exploits the specific acoustic phase modulation produced by the MR characteristic motions. By applying this principle, we performed real-time visualization and position tracking of a magnetic MR rolling on a lumen boundary, both in static flow and opposing flow conditions, with an average error of 0.25 body-lengths. Overall, the reported results unveil countless possibilities to exploit the proposed approach as a robust feedback strategy for closed-loop control of medical MRs in-vivo. Phase-based a (a) In the dened the can assume two different states: rolling and idle. A rotating magnetic eld activates rolling while vibrations in the idle state are induced by a harmonic oscillation of the magnetic eld. Two different phase analysis strategies are put in place according to the MR status, namely rotations detection and vibrations detections. (b) Results of the trajectory tracking experiments. The hybrid approach allowed for robust tracking of the MR during the entire task. The estimated cartesian velocity faithfully reports MR behavior during the different navigation task states (1-5). The scale bar is 1mm.


Introduction
Recent progress in microrobotics provided relevant medical perspectives for non-invasive therapy and diagnosis [1][2][3] . However, the translation of microrobotics technologies to the clinics is hampered by the lack of suitable medical imaging strategies to provide robust and precise feedback for monitoring and control purposes 4,5 . At present, most of the advancements made in the field have been enabled by direct optical feedback through microscopy [6][7][8][9] . Accomplishing a comparable resolution and image quality in an in-vivo setting by medical imaging techniques is a major challenge, still to be addressed. Several medical imaging techniques have been considered for this purpose 10 , ranging from traditional techniques (e.g., MRI 11 or radiation-based 12 ) to innovative ones, such as photoacoustic 13 or magnetic particle imaging 14 . In this framework, medical ultrasound (US) stands as a mature technology that combines real-time capabilities with a good spatial resolution (100µ − 500µ ), deep tissue imaging (up to 25 far from the probe), no adverse health effects, and low equipment cost 15 , thus being a good candidate also for tracking microrobots (MRs) inside the body. The golden standard for US imaging is Brightness(B)-mode, which has been employed to monitor single agents and collective MR swarms during navigation and cargo delivery tasks [16][17][18][19][20][21][22] . In B-mode, image contrast is associated with objects echogenicity, defined as the ability to scatter US waves back to the source. The signal produced by the MR depends on the intensity of the backscattered waves and its ratio to the signal produced by the surrounding medium. Because of their small size, MRs cause poor scattering, resulting in weak US echoes. Most of the works reported in the state of the art have overcome the problem by visualizing MRs in controlled experimental conditions, where homogeneous and scarcely echogenic media (e.g., eye vitreous humor) were used to improve image contrast 23,24 . On the contrary, most biological tissues are highly heterogeneous and echogenic and produce high contrast imaging artifacts that hinder MR detection and visualization. In an in-vivo scenario, imaging artifacts can introduce large position tracking errors, produce dangerous instabilities, and compromise MRs control in closed-loop. To address this problem, US contrast enhancement could be pursued either with the inclusion of contrast agents 25 , such as microbubbles, or by exploiting the Doppler effect 15 . The first strategy suffers from the gas bubbles short half-life, whereas the latter looks more promising. In fact, moving objects produce acoustic phase lags in the backscattered waves, known as Doppler shifts, which are proportional to the objects' displacements 26 . According to this principle, researchers have recently investigated the possibility of exploiting ultrasound color Doppler imaging to visualize a swarm of microrobots in motion 27 . However, due to the weak Doppler signal produced by the swarm relative to the environment, they could not visualize the swarm directly. Instead, they exploited the local disturbance in the flow stream produced by the actuated swarm as an indirect localization method. In fact, traditional color Doppler is sensitive to all motions occurring in the imaging plane, hampering specific MR motions detection in dynamic backgrounds (e.g., biological tissues subject to physiological motions). Exploiting the specific acoustic phase modulation produced by characteristic MR motions, e.g., vibrations, can help to solve this limitation by enabling to distinguish MR movements from background motions, enabling MR direct visualization 28 . This strategy provides a contrastenhancing mechanism for improved MR imaging in echogenic and dynamic backgrounds, where traditional US imaging modalities fail. However, to make this contrast enhancement mechanism viable for MRs tracking, it should be extended to more common and widespread MR motion patterns, such as rolling, which is particularly efficient for MR navigation in different body districts (e.g., vascular system, urinary system, or gastrointestinal tract [29][30][31] ). To fill the aforementioned gap in US-based MR imaging, in this work the authors combine magnetically activated rolling locomotion of a cylindrical MR with realtime US tracking by detecting specific MR motion patterns through acoustic phase analysis. To simulate MR tracking in a realistic therapeutic task, a flexible detection algorithm was developed to detect the MR both during rolling locomotion (for target reaching) and during the idle state, i.e., when the MR is in the target position and target reaching has to be verified. During the locomotion state, acoustic phase analysis was used to detect the special motion signature produced by MR rotations. On the other hand, in the idle state, it was used to detect magnetically induced in place micro-vibrations. The approach was validated in a lumen of a tissue-mimicking phantom, both in static and counter flow conditions. The combination of the two detection strategies allowed to perform continuous US-guided navigation of the MR even in the presence of high contrast imaging artifacts by robustly tracking its centroid position over time and deriving features such as size and rotation frequency.

Results and Discussions
The proposed US contrast enhancement strategy employed for imaging and tracking is based on the specific acoustic phase modulation produced by MR motions. As a result of the Doppler effect, moving objects scatter US waves which are shifted in phase with respect to the incident wave. If a wave is propagating with wavelength and encounters a moving object, the acoustic phase shift in the backscattered echo is proportional to object displacement along the direction of wave propagation .
The acoustic frequency , given by the time derivative of the acoustic phase, is thus shifted proportionally to the object velocity .
Therefore, acoustic phase signal analysis allows accessing information on object displacement and velocity along the acoustic axis (i.e., the direction of wave propagation), independently on the amplitude of the backscattered echoes. This phenomenon can improve the detectability of moving objects even when they produce weak echoes with respect to the background. To exploit this principle as a contrast enhancement mechanism, we analyzed the rolling motion features. A cylindrical MR rolling on a boundary will appear in the imaging plane of the US probe as a rotating circle ( Figure 1). Given the MR with radius and angular velocity , the linear velocity of a generic point ( , ) with respect to a reference system placed in the center of rotation is given by According to eq (3), is null when = 0 (e.g., point 1 and any point along the acoustic axis y in Figure  1), and maximum for | | = (points 2 and 3). In fact, the linear velocity of points along the y axis is orthogonal to the US wave propagation direction and does not produce any phase shift in the backscattered echoes (curve 1 in Figure 1). On the other hand, the linear velocities of points with non-null x coordinate have a component parallel to the US wave propagation direction and produces a continuous phase drift in the backscattered echoes (curves 2 and 3). This characteristic velocity distribution can be exploited as a special signature to detect and visualize the rotating MR in the imaging plane, by properly processing the acoustic phase signal of backscattered waves.

Figure 1. Acoustic phase shifts in the backscattered echoes induced by MR rotation. The dashed line represents the direction of wave propagation, namely the acoustic axis. Points in positions 2 and 3 continuously move respectively towards and out from the US probe, producing a continuous shift (with opposite sign) in the acoustic phase of received echoes. On the other hand, points in position 1 and all points along the acoustic axis y are not moving with respect to the probe. These rotating motion features produce a specific velocity distribution, which is exploited for MR detection and visualization.
To validate the proposed tracking strategy in an in-vivo like environment, we navigated a cylindrical magnetic MR (550 µ diameter, 990 µ length) through a tissue-mimicking phantom with a lumen (3 diameter) (Figure 2a). For more details regarding MR fabrication, the reader may refer to Section A of the Supplementary Information (Supplementary Figure S1). The magnetic MR was placed in the lumen and actuated through external inhomogeneous magnetic fields, generated by a cylindrical permanent magnet rotated through a robotic arm. The pulling force generated by the magnetic field gradients was exploited to migrate the MR at the lumen boundary and grant adherence to its surface, a precondition for rolling locomotion (Figure 2b). The rotation of the permanent magnet around its axis produced a magnetic torque aligning the magnetic moment of the MR with the external field, thus activating rolling locomotion. To demonstrate the robustness of the tracking approach in heterogeneous and dynamic backgrounds, a flow of blood-mimicking fluid was generated inside the lumen (3 −1 ) in the opposite direction of MR locomotion. For tracking the MR, the phantom was imaged with standard clinical US equipment while the raw Radio-Frequency (RF) data were processed by custom algorithms implemented on a desktop computer.

Figure 2. US-guided magnetic navigation platform. (a) The cylindrical magnetic microrobot is positioned inside the lumen of a tissue-mimicking phantom. The scale bar is 500 . (b) The MR is actuated by means of external magnetic fields produced by a permanent magnet mounted to the tip of a robot arm. The pulling force generated by the magnetic field gradients was exploited to migrate the MR at the lumen boundary and grant contact to its surface. The flow in the lumen is produced by a fluidic circuit supplied by a micropump. A standard US probe is used for continuous imaging of the phantom, and a Desktop computer is used for RF data processing, tracking, and visualization.
For the benefit of time efficiency, which is a requirement for real-time imaging, we developed a tracking algorithm working at two different levels: i) a coarse analysis, performed over the whole imaging plane, allows for fast MR localization by identifying the centroid position; ii) a finer and more computationally demanding analysis is then performed only on a reduced region of interest (ROI), defined in the identified centroid neighborhood, to derive MR features. The coarse analysis (localization) starts with raw RF data acquisition from the US transducer ( Figure 3). The RF data associated with consecutive frames constitute a 3D matrix called cineloop. The first two dimensions of the cineloop, ( ) and ( ), represent the lateral and axial dimensions for each frame, respectively. The third dimension ( ) represents the acquired frames index in time. The analytic acoustic signal * is obtained from the RF data cineloop through the Hilbert transform. For every pixel ( , ) in the cineloop, * along the time dimension is expressed as: ( ) is the instantaneous amplitude, converted into greyscale levels to obtain the pixel intensity for Bmode images 32 . The instantaneous phase ( ), which carries information on objects' displacements, is time derived to obtain the instantaneous frequency ( ). The average of ( ) over the acquired frames provides a 2D image representing the mean velocity for each pixel, which is defined as the Motion Image (MI) (Figure 3, in the center). Rotating objects in the MI produce the velocity distribution described in eq. (3). Therefore, to detect a rotating MR, the MI is cross-correlated with a velocity template, defined by eq.
(3). The template considers unitary rotation frequency since the maximum cross-correlation is expected in the MR centroid regardless of its rotation frequency. This condition stands true as long as no other circular object with comparable size to the MR is continuously rotating in the imaging plane. Once identified the MR centroid, a ROI centered in its position is defined, and the fine analysis is started. A Filtered Motion Image (FMI) is obtained from the ROI and is analyzed to estimate MR features. The MR diameter is estimated as the distance between the pixels in the FMI featuring the maximum and minimum intensity (velocity), respectively (as in points 2 and 3 of Figure 1). The MR rotation frequency is given by the maximum pixel intensity in the FMI, which corresponds to the tangential velocity on the MR surface (as again in points 2 or 3 of Figure 1), divided by the radius (half the estimated diameter). Finally, the FMI is overlapped with the B-mode image to perform enhanced imaging and tracking. This allows integrating morphological and anatomical information provided by B-mode with specific MR information provided by FMIs.

Figure 3. Acoustic phase processing for enhanced contrast US imaging and tracking. Algorithm flow chart: Ultrasound RF data are acquired in the form of a cineloop from the probe. For each pixel in the cineloop, the analytic acoustic signal in the time dimension (n) is obtained through the Hilbert transform from the raw RF data. The instantaneous amplitude A is used for B-mode imaging, while the instantaneous phase is analyzed for MR tracking. Fast MR localization is achieved through a coarse phase analysis over the whole image plane. Then a fine analysis is performed in a ROI defined around the centroid, and features as MR diameter and rotating frequency are estimated by analyzing the velocity distribution. Finally, the selected ROI (Filtered Motion Image) is overlapped on a B-mode image for enhanced imaging and tracking. Illustrative outcomes of the process along the pipeline are displayed by the brighter arrows in the flow chart.
First, we validated the performances of the algorithm in estimating MR features for variable diameters and rotation frequencies. The algorithm resulted accurate with a maximum estimation error of 3% and 13%, for diameter and frequency, respectively (Section B of Supplementary Information and Supplementary Figure S2). Then, we conducted navigation experiments to demonstrate the method potential for dynamic tracking of MR trajectories in the presence of high contrast imaging artifacts. To this aim, we rolled the MR on the boundary of the fluid-filled lumen of the tissue-mimicking phantom under continuous US imaging (Figure 4a), first in static fluid conditions and then with an induced flow, opposing the rolling direction. In all the experiments, the MR was actuated with a magnetic field rotating at the constant frequency of 1.5Hz and monitored over a 10s acquisition time with a tracking frame rate of 3 . During navigation, we compared the performances of the proposed algorithm with those of a Bmode algorithm based on image differentiation, particularly suitable for tracking rotating objects 22 (Figure  4b). By analyzing the tracked trajectory points in the static fluid experiments, we observed that the Bmode tracking algorithm intermittently failed in localizing the MR. This was caused by background elements, e.g., the lumen boundaries, which showed considerably higher contrast than the MR (image inset in Figure 4b). In fact, most of the tracked points were concentrated in few points along the lumen boundary, corresponding to the highest contrast pixels in the image. On the other hand, the phase analysis-based algorithm allowed to robustly track the MR positions along the entire trajectory. By overlapping the FMIs acquired during different time instants of the MR trajectory with a B-mode image, we could verify two aspects: (i) despite the low contrast of the MR relatively to the phantom (image inset in Figure 4b), the specific phase signature produced by its motion allowed for enhanced visualization inside the lumen; (ii) throughout the entire trajectory, the MR centroid was precisely localized. In particular, assuming pure rolling and rigid bodies, the average centroid tracking error was 135 µ (about 0.

Results of trajectory tracking in the presence of an induced flow. The B-mode tracker is disturbed by high contrast particles in the flow (bubbles), while phase-based tracker can follow the MR robustly. The presence of an opposite flow introduces slip in the rolling motion, reducing the MR cartesian velocity until it eventually approaches zero when the MR runs across the inclined section of the lumen (at 10s). The scale bar is 1mm.
Lastly, we conducted additional experiments to validate the proposed tracking method in a multi-modal MR navigation task ( Figure 5 and Supplementary Video). For this purpose, we developed a flexible tracking algorithm to detect the MR not only during rolling locomotion but also when idle in-place and not rotating. Specifically, a different phase analysis strategy was put in place depending on the MR state (either rolling or idle state): during the rolling state, the MR was actuated with a rotating magnetic field and detected through its characteristic rotations, with the algorithm proposed in this paper; during the idle state, the magnet was rotated to produce a harmonic motion over a circular sector to induce in-place micro-vibrations of the MR, at a frequency of 5 ( Figure 5a). These vibrations allowed to track the MR through frequency analysis of the acoustic phase signal, as previously reported by the authors 28 . When passing from the rolling state (i.e., employed to reach the target location) to the idle state (i.e., once reached the target), the robot arm controller communicates the MR status to the tracker, which can switch detection modality accordingly. We designed a simple task for experimental validation: the MR is navigated from the start point (1) to the target point a by rolling motion (2). When the target a is reached (3), the MR waits idly in place (vibration state) for 2 . After this stage, in which the reaching of the target location can be verified (i.e., in a real in vivo environment by identifying potential anatomical markers), the MR is actuated again in rolling mode (4) towards target b. Once target point b is reached, the MR waits idly again for 2 and finally completes the task (5). The FMIs overlapped with B-mode show the MR positions in the lumen at different time instants, demonstrating continuous tracking over the entire trajectory (Figure 5b). The position tracking error was in line with previous experiments (less than 0.25 body-lengths). The time evolution of cartesian velocity faithfully reflected the MR behavior during the different navigation task states: velocity was almost constant during rolling through task states 2 and 4. During task states 3 and 5, the velocity was close to zero as the MR was waiting in place idle, and the tracked trajectory points concentrated in the targets a and b.

Conclusion
In conclusion, we proposed a simple yet powerful strategy to improve MR visualization with US in heterogeneous and dynamic media, such as biological tissues. Particularly, we exploited the specific acoustic phase modulation produced by characteristic MR motions as a contrast-enhancing mechanism for robust real-time imaging in echogenic backgrounds, where standard modalities fail. Position tracking performances were validated through navigation experiments conducted in a tissue-mimicking phantom, first in static fluid conditions and then with an opposing flow. Furthermore, we proposed and validated a flexible phase analysis approach to adapt the MR detection strategy according to its different states during a navigation task (e.g., rotations or vibrations). In all experiments, the proposed approach provided direct MR visualization and robust position estimation over the entire trajectories, with an average tracking error of about 0.25 body-lengths, outperforming B-mode, which provided poor MR visualization and produced substantial tracking errors (up to 40 body-lengths) due to background disturbances. Furthermore, we could estimate MR diameter and rotation frequency with a maximum estimation error of 3% and 13%, respectively. Although the MR presented in this study featured fixed size and constant rotation frequency, real-time estimation of these features can be a powerful tool when applied to shape morphing MR 33 or nanoparticles swarms 34 . Overall, the experiments validated acoustic phase analysis as a successful approach for visualizing and tracking MRs in high echogenic and dynamic backgrounds and during different task states, offering new insights into smart combinations of actuation and imaging at the microscale. In fact, the proposed contrast-enhancing mechanism is not only limited to detecting rotating or vibrating motion patterns. Instead, it is more generally based on detecting any displacement along the direction of wave propagation, a typical condition of several microrobotics tasks, including different navigation and function activation strategies 1,35,36 . These perspectives unveil countless possibilities to exploit acoustic phase analysis as a robust feedback strategy for closed-loop control of medical MRs invivo.

Methods
Tissue-mimicking phantom: The phantom was designed to mimic the size and acoustic properties of a tract of the human medium artery (3 − 4 in diameter) with the surrounding soft tissue. The tissue-mimicking phantom was made from agarose and soy milk used as a scatter-enhancing agent. Agarose powder (Sigma-Aldrich) was dissolved in deionized and degassed water (dd-H2O) -soy milk (5% v/v) solution and kept at 90 ° for 1 ℎ under continuous stirring. Selecting the proper agarose concentration (2% v/v) produces mechanical and acoustic properties that mimic human tissues 37 . Physical reticulation occurred at room temperature in the target mold (standard Petri dish). A 3mm diameter rubber tube was embedded in the phantom before reticulation to obtain a lumen for navigation experiments. After reticulation the tube was removed to generate the desired lumen in the phantom.
Blood-mimicking fluid: A fluid that mimics the viscosity and acoustic properties of blood was obtained from an aqueous glycerol solution (60% v/v) 38 .
Rolling Microrobot: To perform controlled rolling over the lumen internal surface through external magnetic fields, we required a sub-millimeter cylindrical MR with remanent magnetization along the radial direction. For this purpose, we employed extrusion-based printing of a UV curable magnetic ink (Supplementary Figure S1). To achieve uniform radial magnetization, the printed cylindrical string was magnetized radially by an impulse magnetizer with a peak field intensity of 1.8 (T-Series, Magnet-Physik Dr. Steingroever GmbH, Germany). The final length of the cylindrical MR was defined by cutting the magnetized string into smaller segments with required aspect ratio (2:1) by using a scalpel under microscope guidance. The resulting MR length and diameter were measured respectively around 990 µ and 550 µ , using a digital microscope (KH-7700, Hirox Co., Ltd, Japan), and the remanent magnetization was measured 27.5 ± 1.7 −1 , using a Vibrating Sample Magnetometer (Model10 VSM, MicroSense, USA).
Ultrasound-guided magnetic navigation in phantom: To validate the proposed tracking strategy in an invivo like environment, we built an experimental platform to navigate the magnetic MR through the lumen of the tissue-mimicking phantom. The magnetic MR was placed in the lumen and actuated through external inhomogeneous magnetic fields, generated by a cylindrical permanent magnet (6 in diameter, 7 in height, NdFeB, diametral magnetization, grade N35) positioned 10 far away from the phantom and rotated through a robotic arm. The pulling force generated by the magnetic field gradients was exploited to migrate the MR at the lumen boundary and grant adherence to its surface, a precondition for rolling locomotion. The rotation of the permanent magnet around its axis produced a magnetic torque aligning the magnetic moment of the MR with the external field, thus activating rolling locomotion.    Figure 1 Acoustic phase shifts in the backscattered echoes induced by MR rotation. The dashed line represents the direction of wave propagation, namely the acoustic axis. Points in positions 2 and 3 continuously move respectively towards and out from the US probe, producing a continuous shift (with opposite sign) in the acoustic phase of received echoes. On the other hand, points in position 1 and all points along the acoustic axis y are not moving with respect to the probe. These rotating motion features produce a speci c velocity distribution, which is exploited for MR detection and visualization. US-guided magnetic navigation platform. (a) The cylindrical magnetic microrobot is positioned inside the lumen of a tissue-mimicking phantom. The scale bar is 500 μm (b) The MR is actuated by means of external magnetic elds produced by a permanent magnet mounted to the tip of a robot arm. The pulling force generated by the magnetic eld gradients was exploited to migrate the MR at the lumen boundary and grant contact to its surface. The ow in the lumen is produced by a uidic circuit supplied by a micropump. A standard US probe is used for continuous imaging of the phantom, and a Desktop computer is used for RF data processing, tracking, and visualization. Acoustic phase processing for enhanced contrast US imaging and tracking. Algorithm owchart:

Figures
Ultrasound RF data are acquired in the form of a cineloop from the probe. For each pixel in the cineloop, the analytic acoustic signal in the time dimension (n) is obtained through the Hilbert transform from the raw RF data. The instantaneous amplitude A is used for B-mode imaging, while the instantaneous phase φ is analyzed for MR tracking. Fast MR localization is achieved through a coarse phase analysis over the whole image plane. Then a ne analysis is performed in a ROI de ned around the centroid, and features as MR diameter and rotating frequency are estimated by analyzing the velocity distribution. Finally, the selected ROI (Filtered Motion Image) is overlapped on a B-mode image for enhanced imaging and tracking. Illustrative outcomes of the process along the pipeline are displayed by the brighter arrows in the ow chart. tracking fails due to objects with higher contrast than the MR (e.g., lumen boundaries). The phase-based approach succeeds in tracking the MR along the entire trajectory. The almost constant cartesian velocity suggests that the MR performs a nearly pure rolling motion. (c) Results of trajectory tracking in the presence of an induced ow. The B-mode tracker is disturbed by high contrast particles in the ow (bubbles), while phase-based tracker can follow the MR robustly. The presence of an opposite ow introduces slip in the rolling motion, reducing the MR cartesian velocity until it eventually approaches zero when the MR runs across the inclined section of the lumen (at 10s). The scale bar is 1mm. Figure 5 Phase-based tracking in a multi-modal navigation task. (a) In the de ned navigation task, the MR can assume two different states: rolling and idle. A rotating magnetic eld activates rolling while vibrations in the idle state are induced by a harmonic oscillation of the magnetic eld. Two different phase analysis strategies are put in place according to the MR status, namely rotations detection and vibrations detections. (b) Results of the trajectory tracking experiments. The hybrid approach allowed for robust tracking of the MR during the entire task. The estimated cartesian velocity faithfully reports MR behavior during the different navigation task states (1)(2)(3)(4)(5). The scale bar is 1mm.

Supplementary Files
This is a list of supplementary les associated with this preprint. Click to download.