The myokinetic stimulation interface: activation of proprioceptive neural responses with remotely actuated magnets implanted in rodent forelimb muscle

2.1.1. One-DoF force control of a single magnet—linear permanent magnet system

To generate a linear vibration in a magnet suspended in an elastic medium, a remote vibrating magnet (e.g. connected to a linear actuator) can be used (figure 2(A)). If the magnetization axes of the two magnets are aligned, the oscillatory movement of the latter generates a MF gradient that induces a force (and, in an elastic medium, a movement) on the target magnet. Quantitatively, the force generated on the target magnet by a remote magnet, with aligned magnetization axes but opposite directions, reads:

$$\begin{equation}F\left( h \right) = - \frac{{3{\mu _0}{m_{\text{r}}}m}}{{2\pi {h^4}}}\end{equation} \tag{ 1 }$$

where ${m_r}$ and $m$ correspond to the magnitude of the magnetic moment of the remote and the target magnets, respectively, and $h$ is their relative distance. Although equation (1) is nonlinear in $h$, for a small distance variation between the magnets (i.e. $\Delta h$) a linear approximation can be used:

$$\begin{align}F\left( {{h_0} + \Delta h} \right) &\approx - \frac{{3{\mu _0}{m_{\text{r}}}m}}{{2\pi h_0^4}} + \left( {\frac{{6{\mu _0}{m_{\text{r}}}m}}{{\pi h_0^5}}} \right)\nonumber\\ \Delta h &= {F_0} + \Delta F\end{align} \tag{ 2 }$$

where ${F_0}$ is the force applied in the stationary condition, and $\Delta F$ is the change in the applied force. Hence, a linear vibration, proportional to $\Delta F$, can be obtained by varying the distance between the magnets by $\Delta h$ (figure 2(A)). If only the vibratory stimulus is of interest, as the dynamic component of the force is entirely described by $\Delta F$, the term ${F_0}$ in equation (2) can be neglected.

Following this basic design, the linear permanent magnet system was implemented by means of a miniature permanent magnet (square cuboid 3.17 mm side × 9.5 mm length, N50, magnetic moment ${m_{\text{r}}} =$ 0.0746 Am²; CMS Magnetics, Garland, TX, USA) attached to the axis of a position-controlled voice coil motor (VCM, Equipment Solutions, Sunnyvale, CA, USA). The latter could generate vibrations with amplitude larger than 50 µm (peak to peak) at 90 Hz, and included a 150 nm resolution optical linear displacement sensor (Equipment Solutions, VCS-10, SCA-824, Sunnyvale, CA). The VCM was controlled through a dedicated driver on a PC.

2.1.2. One-DoF torque control of a single magnet—single-coil system

To implement a torsional vibration in a magnet suspended in an elastic medium, a single coil can be used (figure 2(B)). Quantitatively, the torque generated on the target magnet by a coil h distance from it, and with perpendicular magnetization axes, reads:

$$\begin{equation}T\left( {\,i\left( t \right)} \right) = {m_{\text{r}}}\left( {\frac{{{\mu _0}{m_{\text{c}}}}}{{2\pi {h^3}}}} \right) = {m_{\text{r}}}\left( {\frac{{{\mu _0}{k_{\text{c}}}i\left( t \right)}}{{2\pi {h^3}}}} \right)\end{equation} \tag{ 3 }$$

where ${m_{\text{r}}}$ is the target magnet magnetic moment, ${k_{\text{c}}}$ is a coil-dependent scale factor and $i\left( t \right)$ is the current supplied to the coil (with the equivalent dipole fitted on the coil surface). By applying a time-varying current, it is possible to induce a time-varying torque on the target magnet, thus generating a torsional vibration on the plane xz of the magnet (figure 2(B)).

The design of the single-coil system was based on equation (3), and used an electromagnetic coil (ITS-MS 3025, Intertec Components GmbH, Germany; characterized as in [18], table 1), a dedicated current driver (MC5004 P RS/CO, Faulhaber Minimotor SA, Switzerland), and a PC. The latter hosted the controller, implemented as a C++ application running on a real-time OS (Ubuntu Linux with the RT-PREEMPT patch, 500 Hz control frequency).

Table 1. Electromagnetic coils characteristics.

Parameter	M. U.	Value
Size (diameter, length)	mm	30, 25
Magnetic moment scale factor	m2	753 ± 7
Dipole position (along z a )	mm	−2.67 ± 0.19
Resistance	Ω	40.52 ± 0.46
Inductance	mH	189 ± 3

^a axis perpendicular to the face of the coil.

2.1.3. Five-DoF control of a single permanent magnet—coil-array system

By generalizing the concept and the equations of the single-coil system it is possible to devise a coil-array architecture that provides full controllability over a target magnet. The force and torque generated on a target magnet with magnetic moment $\vec m$, due to a MF produced by a coil/electromagnet ${\vec B_{{\text{coil}}}}$, reads:

$$\begin{equation}\vec F = \nabla \left( {\vec m \cdot {{\vec B}_{{\text{coil}}}}\left( {\vec p} \right)} \right),{ }\vec T = \vec m \times {\vec B_{{\text{coil}}}}\left( {\vec p} \right)\end{equation} \tag{ 4 }$$

where $\vec p$ is the magnet position vector [19]. ${\vec B_{{\text{coil}}}}$ depends linearly on the applied current $i\left( t \right)$, and by adopting the magnetic dipole model [18], it reads:

$$\begin{align}{\vec B_{{\text{coil}}}} &= \frac{{{\mu _0}{k_{\text{c}}}}}{{4\pi }}\left( {\frac{{3\vec r\left( {\vec l \cdot \vec r} \right)}}{{{{\left| {\vec r} \right|}^5}}} - \frac{{\vec l}}{{{{\left| {\vec r} \right|}^3}}}} \right)i\left( t \right) = K\left( {\vec p} \right)i\left( t \right)\nonumber\\ &= \left[ {\begin{array}{*{20}{c}} {{k_x}\left( {\vec p} \right)} \\ {{k_y}\left( {\vec p} \right)} \\ {{k_z}\left( {\vec p} \right)} \end{array}} \right]i\left( t \right)\end{align} \tag{ 5 }$$

where $\vec r = \vec p - {\vec r_{{\text{coil}}}}$, ${\vec r_{{\text{coil}}}}$ is the position vector of the coil, ${k_{\text{c}}}$ is a scale factor and $\vec l$ is the unit vector representing the coil orientation.

From equations (4) and (5), it follows that the force $\vec F$ (and the torque $\vec T$) can be modulated through the coil current. In the case of C coils operating under certain assumptions [20], the compound MF generated in a point in space is simply defined as the sum of the fields generated by each coil:

$$\begin{align}{\vec B_{{\text{coils}}}} &= \sum\limits_{j = 1}^C {{{\vec B}_{{\text{coil}},j}}} = \sum\limits_{j = 1}^C {{K_j}} \left( {\vec p} \right){i_j}\nonumber\\& = \left[ {\begin{array}{*{20}{c}} {{k_{x,1}}\left( {\vec p} \right)}& \cdots &{{k_{x,C}}\left( {\vec p} \right)} \\ {{k_{y,1}}\left( {\vec p} \right)}& \cdots &{{k_{y,C}}\left( {\vec p} \right)} \\ {{k_{z,1}}\left( {\vec p} \right)}& \cdots &{{k_{z,C}}\left( {\vec p} \right)} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {{i_1}} \\ \vdots \\ {{i_C}} \end{array}} \right] = {\overline {\overline K}}\left( {\vec p} \right)\vec i\end{align} \tag{ 6 }$$

with ${\overline {\overline K}}\left( {\vec p} \right)$ a $3 \times C\,$matrix and $\vec i$ a C-dimensional vector. By substituting equation (6) in equation (4) we get:

$$\begin{equation*}\vec F = \nabla \left( {\vec m \cdot {{\vec B}_{{\text{coils}}}}} \right) = \nabla \left( {\vec m \cdot \overline {\overline K} \left( {\vec p} \right)\vec i{ }} \right) = {\overline {\overline M} _{\text{f}}}\vec i\end{equation*}$$

$$\begin{equation*}\vec T = \vec m \times {\vec B_{{\text{coils}}}} = \vec m \times \overline {\overline K} \left( {\vec p} \right)\vec i = {\overline {\overline M} _\tau }\vec i.\end{equation*}$$

These two systems, of three equations each, can be arranged in a unique system of six equations:

$$\begin{equation}\left[ {\begin{array}{*{20}{c}} {\vec F} \\ {\vec T} \end{array}} \right] = \left[ {\begin{array}{*{20}{c}} {{{\overline {\overline M} }_{\text{f}}}} \\ {{{\overline {\overline M} }_\tau }} \end{array}} \right]\vec i = \overline {\overline M} \vec i{ }\end{equation} \tag{ 7 }$$

with ${\overline {\overline M}}$ being a $6 \times C$ matrix. Thus, to generate specific forces and torques, e.g. vibration profiles, at a desired location within the workspace, it is sufficient to know (or identify) the entries of ${\overline {\overline M}}$ and to invert equation (7) to find $\vec i$. Notably, this requires $C$ to range from $6$ (when passive restoring forces are present, e.g. due to elasticity in the medium) to 8 (full independent control) [20, 21].

The design of the coil array system was based on these general equations. It consisted of $C = 8$ electromagnetic coils and drivers, two sensor boards each containing three 3D MF sensors (HMC5983, Honeywell; ±8 G max. range, 0.73 mG maximum resolution) (figure 3), and a PC with real-time OS. The electromagnetic coils were arranged around a hollow cylinder (d = 30 mm), representing the workspace of the device, and pointed towards its center (figure 3). The current of each coil was controlled by dedicated driver connected to the PC. As with the single-coil system, the latter hosted the controller, which both localized the magnet, based on the MF sampled by the sensors, and generated the current profiles to induce the desired forces/torques (periodic sinusoids, square waves, and generic polynomial signals). In this implementation, the desired signal shape, frequency, and amplitude could be chosen manually. The localization of the magnet was done once, to calibrate the system (or to identify ${\overline {\overline M}}$) (equation (7)) [22–24]; during normal operation, the position of the magnet was assumed to be time-invariant.

Zoom In Zoom Out Reset image size

2.2.1. Neural recordings

2.2.2. Animal preparation

2.2.3. Magnet implantation

After collecting the responses to these mechanical stimuli, the clip was removed from the distal tendon, and a cylindrical NdFeB permanent magnet (N50 magnetization, 1 mm diameter, 2 mm height, magnetic moment $m = \,$0.0016 Am²) was inserted into a pocket created in the muscle belly and covered by filter paper affixed with Vetbond tissue adhesive (3 m, St. Paul, MN) (figure 4(A)). While the adhesive was drying, the tendon clip was substituted by the linear permanent magnet system (figure 4(B)). For the linear permanent magnet system, the VCM was positioned so that the movement of the remote magnet moved the implanted magnet within the muscle belly. We used a micro-positioner to adjust the initial distance between the remote and the implanted magnet (${h_0}$) until successful magnetic coupling was achieved, a distance of ≈2–3 mm. The same static and small/large displacement frequency stimulus presentations were run as described above. It is worth noting that the displacements presented as input stimuli did not necessarily match with those at the (remote) muscle/magnet site, as the latter depend on equations (1) and (2) and on the dynamics of the system.

Zoom In Zoom Out Reset image size

2.2.4. Static displacement and vibrational response characterization

In-vitro tests were performed to verify the ability of the coil-array system to induce torsional and linear vibrations with different waveforms, on a remote target magnet. These tests were carried out with visualization by a high-speed video camera (Sony DSC-RX10M4; 1000 fps maximum frame rate, with a 1080 p maximum resolution), which allowed to direct observation of the displacement induced on the remote magnet in a single plane (xy in figure 5). A LED was used to synchronize the video with the actuation. The magnet, a commercial cylindrical NdFeB permanent magnet (N48 magnetization, 3 mm diameter, 8 mm height, magnetic moment $m = \,$0.0584 Am²), was rigidly attached on top of a cylinder ($d = 5\;{\text{mm}}$, $h = 10\,{\text{mm}}$) of viscoelastic material (Ecoflex 00–30, Smooth-On Inc., Macungie, Pennsylvania, USA), resembling the stiffness of muscular tissue [25] (cylinders mounted to each other on their curved side, figure 5, supplementary video 1 available online at stacks.iop.org/JNE/19/026048/mmedia).

Zoom In Zoom Out Reset image size

One-DoF sinusoidal, square, and sawtooth, torsional ($\vec T$) and linear ($\vec F$) vibrations were generated at 20 and 90 Hz on the plane observed by the camera. We tested (a) torsional vibrations with a peak to peak torque amplitude of 900 $\mu$Nm, (b) linear vibrations with a peak to peak force amplitude of 64 mN, along the axis orthogonal to the magnetization axis of the magnet, and (c) of 40 mN along the axis parallel to the magnetization axis of the magnet. All vibrations had a zero mean (no DC component). For each experimental condition (three waveforms, two frequencies, and one amplitude) we ran five repetitions lasting for 3.5 s each (supplementary video 1). The videos, recorded at 500 fps, were processed frame by frame, using custom image processing algorithms developed in MATLAB (MathWorks, Natick, MA, USA) to reconstruct the magnet orientation and the position of its centroid, namely the planar absolute pose (blob analysis, figure 5). The position ($\rho$) and angular displacements ($\theta$) relative to the initial planar pose (when no actuation was delivered) were computed. The correlation coefficient, r, between the normalized power spectral densities (PSD) of the vibration signal (actuation signal) and of the magnet displacement (input and output waveforms, respectively), along the desired axes, was computed to assess frequency distortions in the generated vibrations.

For sinusoidal-shaped vibrations, the PSD of the magnet displacement signals were also used to extract some parameters, to quantify the coherence of the induced vibrations [18]. These were: (a) ${P_{{\text{pos,d}}}}$ and ${P_{{\text{pos,d|f}}}}$, namely the PSDs of the displacement signal along the desired actuation direction, in the whole spectrum and at the frequencies of interest (respectively), when a force was applied, (b) ${P_{{\text{pos,u}}}}$, the PSD of the displacement signal along the undesired (orthogonal) direction, in the whole frequency spectrum, when a force was applied, (c) ${P_{{\text{ang}}}}$ and ${P_{{\text{ang|f}}}}$, the PSD of the angular rotation, in the whole spectrum and at the frequencies of interest, respectively, when a torque was applied. In particular, the values of ${P_{{\text{pos,d|f}}}}$ and ${P_{{\text{ang|f}}}}$ were computed as the area of the main lobe of the PSD, centered at the frequencies of interest (i.e. 20 or 90 Hz, figure 5). These were used to assess the efficiency of the system according to the following equations:

$$\begin{equation}{\text{Torsional efficiency}}:{\text{ }}\gamma = \frac{{{P_{{\text{ang|f}}}}}}{{{P_{{\text{ang}}}}}}\end{equation} \tag{ 8 }$$

$$\begin{equation}{\text{Directional efficiency}}:{\text{ }}\eta = \frac{{{P_{{\text{pos,d|f}}}}}}{{{P_{{\text{pos,d}}}} + {P_{{\text{pos,u}}}}}}.\end{equation} \tag{ 9 }$$

The maximum efficiency was expected when the vibration was generated at the selected frequency (unique frequency in a sinusoid) and along the specified direction only.

Using the linear permanent magnet system, we recorded activity from a muscle sensory afferent neuron in the median nerve of the rat forelimb. By surgically degloving the forelimb and denervating the forepaw of rats, we silenced any potential confounding sensory input from the cutaneous receptors. We found that the linear permanent magnet system was able to elicit neural responses that were indistinguishable from responses that were elicited with matched mechanical actuation (figure 6(A)). In this neural recording, the magnetic actuation was able to reveal slowly adapting responsiveness to static displacement over one second (figure 6(B)). Similarly, the magnetic actuation produced neural responses that tracked with sinusoidal displacements across all tested frequencies from 20 to 150 Hz (figures 6(C)–(H)).

Zoom In Zoom Out Reset image size

At the 20 and 50 Hz displacements, the neurons fired at a rate close to one action potential per cycle (average of 104.0 spikes per 100 cycles at 20 Hz; and 100.0 spikes per 100 cycles at 50 Hz). This response activity reflected the spike-per-cycle activity of the large-amplitude mechanical displacements at the same frequencies (figure 7). As the frequency of the linear permanent magnet system actuation increased (70, 90, and 100 Hz), the number of average action potentials per cycle decreased steadily (average of 50.8 spikes per 100 cycles at 70 Hz; 61.9 spikes per 100 cycles at 90 Hz; and 35.8 spikes per 100 cycles at 100 Hz). At 150 Hz the spike-per-cycle activity of the linear permanent magnet actuated neural responses reflected the values obtained at 20 µm displacement with the direct mechanical displacements (15.6 spikes per 100 cycles at 150 Hz).

Zoom In Zoom Out Reset image size

Using the single-coil system we were able record solenoid-induced neural signals from the same sensory receptor represented in the direct mechanical and linear permanent magnet displacement recordings. Although the neural activation signals obtained across multiple vibratory frequencies (50, 70, 100, and 150 Hz) were complicated by electromagnetic interference due to the solenoid's proximity to the recording electrode we were able to reliably separate the neural signal from the background noise for the 50 and 70 Hz vibrational displacements. We found that the solenoid was able to activate neural activity at an average of 1.2 action potentials per vibratory cycle (five banks of 100 cycles each) at 50 Hz, and an average of 1.6 action potentials per vibratory cycle (five banks of 100 cycles each) at 70 Hz.

The several tests demonstrated that it was possible to produce directional and frequency selective vibrations along different magnet axes, using a five-DoF control, at both 20 and 90 Hz. When providing sinusoidal torsional vibration at 20 Hz as the control signal, the actual orientation $\theta$ followed the sinusoidal shape almost perfectly with amplitude ${\theta _{{\text{MAX}}}} \cong 7^\circ$ and oscillating frequency centered at 20 Hz (figure 8(A)). The torsional efficiency proved $\gamma \,$ = 0.898, meaning that approximately 90% of the power was concentrated around the desired frequency. The (undesired) displacements of the centroid of the magnet were confined to less than ±50 $\mu {\text{m}}$ (∼10% the arc subtending ${\theta _{{\text{MAX}}}}$), with a larger intensity on the direction orthogonal to the magnet axis.

Zoom In Zoom Out Reset image size

When compared to the 90 Hz signals, the 20 Hz sinusoids generally led to a larger efficiency together with a smaller dispersion of values around the mean (figure 8(B)). For the torsional vibrations, γ was equal to 0.897 ± 0.002 (average ± std. dev.) along the z rotational direction. Its value slightly decreased to 0.872 ± 0.017 at 90 Hz. Regarding the linear vibrations, $\eta$ proved equal to 0.897 ± 0.001 and 0.878 ± 0.006, respectively for the sinusoids at 20 and 90 Hz, for displacements orthogonal to the magnetization axis (${\eta _{\text{O}}}$), while it proved equal to 0.797 ± 0.03 and 0.58 ± 0.071, respectively for sinusoids at 20 and 90 Hz, for displacements parallel to the magnetization axis (${\eta _{\text{M}}}$).

With square wave and sawtooth waveforms, the comparison of the normalized theoretical and experimental PSD of the vibrations generated by the coil array system exhibited only minimal differences (representative examples are shown in figure 9). Irrespective of the shape (i.e. square wave, sawtooth) and frequency (i.e. 20 and 90 Hz), we retrieved a good overlap of the PSDs in the locations of the principal harmonics, for both torsional and linear vibrations. Indeed, the two signals appeared to share the same frequency peaks, although the output PSD exhibited a lower amplitude in the proximity of the harmonics at higher frequencies (figure 9).

Zoom In Zoom Out Reset image size

For square waves, irrespectively on the vibration frequency, the average correlation coefficient r proved 0.88 ($p$ < 0.05) and 0.94 ($p$ < 0.05), for linear and torsional vibrations, respectively.

For sawtooth waveforms, r was dependent on the vibration frequency and the actuation direction (for linear vibrations). In particular, in the case of 20 Hz vibrations, we found that r = 0.94 ($p$ < 0.05) and 0.92 ($p$ < 0.05), for linear vibrations orthogonal to the magnetization axis and torsional ones, respectively. These values dropped to r = 0.8 ($p$ < 0.05) and 0.65 ($p$ < 0.05) in the case of 90 Hz vibrations. Conversely, the correlation coefficient for linear vibrations along the magnetization axis increased from r = 0.26 ($p$ < 0.05) to r = 0.75 ($p$ < 0.05), respectively for 20 and 90 Hz.

The muscle sensory receptors, joint capsule receptors, and the skin around the joints all play a role in the sense of the movement and position of our limbs and bodies in space [4, 6]. In particular, the type Ia muscle spindle fibers are considered to be the sensory receptors responsible for signalling limb movement and position and are responsive to static displacements and a wide range of vibratory frequencies [29]. With the use of the linear permanent magnet system we were able to elicit neural responses from a slowly adapting muscle sensory receptor. We found that we could activate action potentials in a single afferent receptor with static displacements and vibration induced by the coupling between a permanent magnet that was implanted in the muscle and a remote permanent magnet attached to a linear actuator. The ability to apply a static displacement using magnetic actuation suggests that different onset velocities can be provided to signal rate of change of movement to proprioceptive afferents. Furthermore, the frequencies that we were able to achieve spanned from 20 to 150 Hz; frequencies that are known to activate the type Ia muscle spindle receptors [30]. Specifically, illusions of limb movement can be triggered with vibration of the muscle at frequencies above 70 Hz [7, 8, 30, 31]. In our previous work, we have shown that 90 Hz vibration of deep muscles in human patients with amputation and a neural-machine interface (targeted reinnervation) improves functional prosthetic outcomes through kinaesthetic perception of hand grasping movements [14, 27]. The linear permanent magnet system provides evidence that muscle sensory action potential responses can be remotely activated at 90 Hz frequency. These results provide evidence in an in-vivo animal preparation that remote actuation of permanent magnets implanted in the muscle can be utilized for effective kinaesthetic prosthetic feedback.

A key feature of the vibrational input from the remotely actuated implanted magnets was that the displacements required to trigger muscle sensory afferent responses appeared to be very small. We used the VCM optical position sensor to physically measure the displacements that were applied directly to the muscle tendon in the control experiments. We could then compare that to the number of action potentials activated by the movement of the implanted magnet for each bank of 100 sinusoids to the number of action potentials elicited in the same neuron by the measured mechanical displacement. With that information, we could infer the displacement range at which the magnets were likely operating (figure 10). At the lower frequencies (20 and 50 Hz) we were able to elicit action potential responses that tracked one to one with each sinusoidal cycle (see: figure 7 red line and open circles). This was similar to measured mechanical displacements eliciting one to one firing that ranged from 221 µm down to 97 µm (see: figure 10). It appears that at 20 and 50 Hz, the magnitude of the implanted magnet displacement was sufficient to achieve a one-to-one response. However, because that critical lower bound was not expressly identified, we can only infer that the amplitude of the implanted magnet is somewhere in the 21–97 µm range for frequencies of 70 Hz and above. Conversely, if the receptor had a response plateau of one impulse per cycle, we are equally unable to determine the upper bound of displacement to achieve the one-to-one response. Moreover, as the frequency of remote magnetic activation increased, the relative rate of firing per 100 cycles decreased, suggesting that the distance the magnet moved likely attenuated (see: figure 7 red line and open circles). At the highest frequency the spikes-per-100-cycles was equivalent to those elicited by the mechanical vibration at 21 µm displacement (see: figure 7 red line and open circles). Moreover, with frequencies higher than 150 Hz, activating Pacinian corpuscles in the muscle is a possibility. However, with magnet displacement attenuation, further experiments would need to be conducted to determine if these receptors could be activated ([32]). Together, these results suggest that the remotely actuated implanted magnets were potentially vibrating within the range of 20–100 µm (figure 10 light red band), although the implanted magnet could be vibrating at a higher displacement if the receptor response plateaued at one impulse per cycle. Similarly, this approach shows the potential for remotely actuating the deep muscle sensory afferents without the necessity of pushing through the skin [14, 27].

Zoom In Zoom Out Reset image size

We were able to demonstrate that a remote solenoid could be used to actuate an implanted permanent magnet to induce single-unit neural responses within the muscle. Solenoid actuation alone is an important proof of concept result because of the increased flexibility of the stimulation parameters (discussed in further detail below). Most importantly, we were able to show neural activation at the 70 Hz frequency, which is the baseline frequency for activation of active kinaesthetic percepts [30]. Nonetheless, the proximity of the coil to the recording electrode made it difficult to effectively record the micro-volt level single unit responses from the axon hook electrode. Although the noise prevented similar recordings with the coil-array system, this work invites further studies in which electromagnetic interference is prevented by design, e.g. by using shields and/or appropriate materials, or robust localization algorithms. It is worth mentioning though, that electromagnetic emissions would not necessarily hamper the implementation of the myokinetic stimulation interface, since their (known) frequency content could be tracked and rejected. In addition, substituting the coils with air-core coils could largely reduce those emissions, favouring the integration of the two systems [21].

The experimental conditions for the characterization of the coil-array system were chosen in order to gather overall insights on the capabilities of the system, and testing them represented an important comparison to the in-vivo tests. In particular, the 90 Hz vibration was motivated by the neurophysiology of the kinaesthetic illusion [14, 27], whereas the 20 Hz corresponded to the lowest frequency tested with the direct mechanical displacement in the in-vivo experiments. The chosen amplitudes (40 and 64 ${\text{mN}}$, for linear vibrations along the magnetization and the orthogonal axis, respectively, and 420 $\mu {\text{m}}$ for torsional vibrations) were the maximum allowed by the coil-array system in all involved experimental conditions, in order to allow for direct comparison. Furthermore, due to visibility limitations of the camera setup (figure 5), we placed the remote magnet on top of the tissue phantom, rather than placing it inside. Finally, the sawtooth and the square waveforms were tested to generalize the physical capabilities of the coil-array system.

The efficiency of the torsional (γ) and linear (orthogonal to the magnetization axis, ${\eta _{\text{O}}}$) sinusoidal vibrations, was large (∼85%–90%) and similar to our previous work [18] (figure 8(B)). This was also true with fewer coils (8 instead of 12, which implies a lower isotropy and controllability of the system), as well as with the different viscoelastic properties of the silicone gel (Ecoflex 00–30 vs. ballistic gelatine). The silicone gel (which was not modelled in the control equations) may have reduced the efficiency from altering the PSD distribution in correspondence to the main frequency components and also at low frequencies (<2 Hz, figure 9(A)). Such a measured spread may also have been an artifact of the video analysis, considering that the frequency resolution of the Fast Fourier Transform was inversely proportional to the size of the temporal window. In particular, having longer acquisitions of the vibrations, or faster acquisition rates, would allow for a better estimate of the PSD, reducing the size of the spectral side lobes. Thus, this limitation may be partially attributed to the available camera, which allowed high speed recordings (500 fps) for no more than 3.5 s. In addition, it is possible that other sources of noise in the video processing may have altered the measures. The main confound was the change in the lightning conditions from the different reflection of the light on the surface of the magnet while moving/vibrating. This was exacerbated by the uneven light across frames due to the 50 Hz power line switching, which led to additional distortion of the measured spectra.

Although the system exhibited relatively large efficiencies for torsional (γ) and linear (orthogonal to the magnetization axis, ${\eta _{\text{O}}}$) vibration, we failed to demonstrate similar efficiencies for linear vibrations parallel to the magnetization axis (${\eta _{\text{M}}}$), especially at 90 Hz (figure 8(B)). This outcome was counterintuitive. Considering the geometrical configuration of the coils in the coil-array system, the lower surface offered by the cylinder base compared to the lateral face, and most of all, standing to equation (5), we expected a larger current/displacement efficiency along the axis of magnetization with respect to the orthogonal axis, as in [18]. Here instead, it was difficult to produce significant displacements along the magnetization axis (i.e. unexpectedly large currents were required, with values close to the maximum produced by the coils). This result may be reflected in an unexpectedly low efficiency ${\eta _{\text{M}}}$. The difficulty in producing displacements along the magnetization axis was the reason choosing a lower amplitude (40 mN vs. 64 mN, peak to peak) to be tested for that axis. However, since this displacement was rather small, and the resolution of the video camera was limited, any environmental noise, undesired movement, or inaccuracies in the video analysis likely had an outsized impact on that particular measurement. The tests with the saw-tooth and square waves generalized the capabilities demonstrated with the sinusoidal vibrations. The imperfect match between the input and output PSDs of the waveforms (figure 9) may be attributed (again) to the viscoelasticity of the medium but also because the input signal was represented in its ideal form (pure tones) rather than with its actual (measured) spectrum.

Motivated by the success of magnetic actuation in generating relevant neural activations in the frequency range of the kinaesthetic illusion, future refinements of the proposed stimulation system would allow the development of bidirectional interface for a prosthetic limb. As a next step we foresee the complete integration of the localizer in the stimulation interface, recalling that in the present work it was not used while inducing vibrations. This step is necessary to induce precise vibrations during muscle deformations caused by voluntary contraction or elongation. In this regard, system performance in terms of spatial resolution (submillimetre) and bandwidth (around 50 Hz for the tracking), seems to not represent a bottleneck for the clinical translation. Final refinements, as well as studies with human participants, could contribute in the future to a better understanding of proprioception.