Stimulus-induced Rotary Saturation (SIRS): A potential method for the detection of neuronal currents with MRI

Technical Note

Stimulus-induced Rotary Saturation (SIRS): A potential method for the detection of

neuronal currents with MRI

Thomas Witzel

⁎

,

Fa-Hsuan Lin

, Bruce R. Rosen

, Lawrence L. Wald

a_{Speech and Hearing Bioscience and Technology Program, Harvard–MIT Division of Health Sciences and Technology, USA} b_{Department of Radiology, Athinoula A. Martinos Center for Biomedical Imaging, Massachusetts General Hospital, USA} c

Institute of Biomedical Engineering, National Taiwan University, Taipei, Taiwan

a b s t r a c t a r t i c l e i n f o Article history: Received 23 April 2007 Revised 15 April 2008 Accepted 1 May 2008 Available online 20 May 2008 Keywords:

fMRI

Neuronal current imaging Rotary saturation T1rho

Neuronal currents produce local transient and oscillatory magneticfields that can be readily detected by MEG. Previous work attempting to detect these magneticfields with MR focused on detecting local phase shifts and dephasing in T2or T2⁎-weighted images. For temporally biphasic and multi-phasic local currents the sensitivity of

these methods can be reduced through the cancellation of the accrued phase induced by positive and negative episodes of the neuronal current. The magnitude of the phase shift is also dependent on the distribution of the current within the voxel. Since spins on one side of a current source develop an opposite phase shift relative to those on the other side, there is likely to be significant cancellation within the voxel.

We introduce a potential method for detecting neuronal currents though their resonant T1ρsaturation during

a spin-lock preparation period. The method is insensitive to the temporal and spatial cancellation effects since it utilizes the multi-phasic nature of the neuronal currents and thus is not sensitive to the sign of the localfield. To produce a T1ρreduction, the Larmor frequency in the rotating frame, which is set byγB1lock(typically 20 Hz–

5 kHz), must match the major frequency components of the stimulus-induced neuronal currents. We validate the method in MRI phantom studies. The rotary saturation spectra showed a sharp resonance when a current dipole within the phantom was driven at the Larmor frequency in the rotating frame. A 7 min block-design experiment was found to be sensitive to a current dipole strength of 56 nAm, an approximate magneticfield of 1 nTat 1.5 mm from the dipole. This dipole moment is similar to that seen using the phase shift method in a similar experimental setup by Konn et al. [Konn, D., Gowland, P., Bowtell, R., 2003. MRI detection of weak magneticfields due to an extended current dipole in a conducting sphere: a model for direct detection of neuronal currents in the brain. Magn. Reson. Med. 50, 40–49], but is potentially less encumbered by temporal and spatial cancellation effects. Published by Elsevier Inc.

Introduction

Magnetoencephalography measures magnetic fields near the scalp and infers a distribution of dipole current sources within the cortex. The rationale is that postsynaptic currents occurring synchronously in a population of tens of thousands of cortical pyramidal cells can be modeled as a single current dipole. Because the cells are aligned almost perfectly parallel to each other, the dipolarfields of all the neurons add construc-tively. Unfortunately it is difficult to learn the microscopic distribution or even amplitude of the currents based on the detected MEGfield pattern since MEG measures fields at the

scalp (far from the dipole) where the structure of the dipole is not easily discernable. For example, it can be difficult to distinguish between a single larger dipole or a group of parallel small dipoles (Hämäläinen et al., 1993). Also, there could be cancellation among local dipoles in the folded cortex ( Hämä-läinen et al., 1993). Nevertheless, an estimate of the equivalent current dipole strengths can be made within the assumption of the MEG reconstruction and typical numbers range from 5 to 120 nAm with many robust stimulus paradigms giving dipole sources in the tens of nAm range (Hari, 1991). The temporal characteristics of the MEG-detected dipoles show a range of responses, most typically temporally asymmetric bipolar and multi-phasic responses with significant spectral power at frequencies below 100 Hz (Niedermeyer and Lopes da Silva, 1999). In addition to the transient evoked response, oscillatory fields can be evoked, for example using auditory click-trains (Lin et al., 2004). Thus in this case, a narrow-band peak at the stimulus repetition frequency is clearly visible in the MEG-⁎ Corresponding author. MGH, A.A. Martinos Center, Bldg. 149, 13th Street, Rm. 2301,

Charlestown MA 02129, USA.

E-mail address:twitzel@nmr.mgh.harvard.edu(T. Witzel). 1053-8119/$– see front matter. Published by Elsevier Inc. doi:10.1016/j.neuroimage.2008.05.010

Contents lists available atScienceDirect

NeuroImage

detectedfields. Similar narrow-band fields have been detected by MEG during visual paradigms with alternating checkerboard stimuli (Tallon-Baudry et al., 2001).

The direct detection of neuronal activity using MR imaging techniques has been pursued for over a decade (Bandettini et al., 2005; Bodurka and Bandettini, 2002; Bodurka et al., 1999; Konn et al., 2003, 2004; Parkes et al., 2007; Petridou et al., 2006; Singh, 1994; Xiong et al., 2003) and the progress has recently been reviewed (Bandettini et al., 2005). Most of the proposed techniques are based on the detection of phase shifts or dephasing caused by the extremely small changes in the local B0field induced by local neuronal currents. The phase

shift or dephasing of the detected MR magnetization accumu-lates during the TE period of either a gradient echo or spin echo sequence. The phase shift experienced by a given spin within the voxel,γΔBneuralTE, can be either positive or negative

de-pending on the sign of the z-component of thefield produced by the local current. The total dephasing measured in the image is proportional to the integral of this phase shift over the voxel volume and its temporal integral over the sensitive pe-riod of the sequence. The phantomfindings are generally en-couraging, with detectable local magneticfields as small as 0.2 nT (Bodurka and Bandettini, 2002) and current dipoles as small as 6.3 nAm (Konn et al., 2003), but the in vivo literature contains some positivefindings (Chow et al., 2006; Xiong et al., 2003), but mostly negativefindings (Chu et al., 2004), includ-ing a failure to reproduce the paradigm of Xiong et al. (Parkes et al., 2007).

Understanding the microscopic spatial distribution and time course of the neuronal currents is critical to assessing the potential ability of these MR-based methods to detect neu-ronal activation. The phasensitive portion of the MR se-quence must be properly timed relative to the transient or oscillating neuronal current to ensure that cancellation of the induced phase during a biphasic or multi-phasic current is minimized. For example, the phase-sensitive portion of the sequence (the TE period of a gradient-recalled echo sequence) might be limited to a single lobe of a biphasic response. Or, for a spin echo sequence, the timing of the 180° pulse can be chosen to occur at the transition between the two phases of the re-sponse allowing the biphasic signal to be detected without cancellation (Singh, 1994). Finally, the conventional methods using T2and T2⁎ are only sensitive to the component of the local

field along the applied external B0field.

The spatial distribution of the current within the voxel can also lead to cancellation of the phase shift induced in the MR signal. Even a DC current in a small, straight wire running through the center of the voxel will not induce a net phase shift in the MR signal in that voxel; the current must be offset from the voxel center to produce a net phase accrual (Bandettini et al., 2005). This arises since the phase is sensitive to the sign of the magneticfield and thus the integral of the phase shift over the voxel is zero for this symmetric case. In other cases, the magneticfield is expected to cancel for points within the voxel. For example, for a uniform distribution of discrete par-allel wires there are some areas where the magneticfield itself cancels, but the integral of the phase shift over the region cancels everywhere but at the edge of the distribution (Park and Lee, 2007). In this case the signal will experience mag-nitude dephasing but not a net phase shift. In a recent simu-lation study, Park et al. (2007) numerically calculated the effects of neuronal magneticfields on voxel magnitude and phase, simulating both dendritic as well as axonal currents.

They estimate average neuronalfields of 0.3 nT for simulta-neous dendritic currents.

Recently, the use of ultra-low-field MRI has been proposed to detect neuronal currents (Kraus, 2006; Kraus et al., 2008). An extremely low B0measurementfield allows resonant

inter-actions between the neuronal currents and spin magnetization ifγB0is brought into the frequency range of the neuronal

cur-rents. Although theirfirst measurements were at significantly higherfields, Kraus and colleagues proposed static measure-mentfields weaker than 0.025 G using a SQUID magnetome-ter to detect the MR signals at a Larmor frequency of less than 100 Hz (Matlachov et al., 2004; Volegov et al., 2004). In this scenario, the transient neuronalfields can possibly act as re-sonant excitation pulse providing the initial excitation of the proton magnetization. Thus, the neuronal current supplies the local B1field for excitation and no transverse magnetization

should be detected in locations without neuronal currents. Contrary to the high-field MRI-based phase shift approaches mentioned above, multi-phasic or oscillatory currents are be-neficial for increasing the tip angle and thus the detected response— provided that they have significant spectral power at the Larmor frequency. However, like the high-field phase shift methods, the technique is still susceptible to cancellation effects from the spatial symmetry of the currents within a voxel. Also, the method is sensitive only to the component of the neuronal field perpendicular to the z axis and the effect is still dependent on the sign of the neuronal current. Thus, a small, straight wire running through the center of the voxel along the static mag-netizingfield direction produces an excitation giving transverse magnetization along opposite directions in the x–y plane on opposite sides of the wire. Thus, no net signal from the voxel is expected for this situation. The low-field method is also ad-versely affected by the small Boltzmann population and other technical confounds of the low-frequency MRI apparatus.

In this work, we introduce a new approach to effectively lower the Larmor frequency below 100 Hz in high-field MR systems based on the spin-locking mechanism (Abragam, 1961). During a spin-lock, the magnetization is locked in the rotating frame by an RFfield, the so-called spin-lock field B1lock. The effective resonance frequency in the rotating frame

during the spin-lock condition is very close toω1lock=γB1lock,

deviating only by a potential Bloch–Siegert shift (Bloch and Siegert, 1940). This latter small shift in resonance frequency results from the counter-rotatingfield component of a linearly polarized time varying field. The Larmor frequency in the rotating frame can be chosen by setting the RF power of the spin-lock pulse. In this work, we acquire phantom data with γB1lockranging from 20 Hz to 200 Hz. T1ρrelaxation occurs in

the spin-lock state from external magneticfield fluctuations with a z component at the resonant frequencyγB1lock and

harmonics.

Redfield described the saturation of the spin-locked mag-netization with an external coil and audio-frequency source as “rotary saturation” (Redfield, 1955). Redfield analyzed rotary saturation in solids using B1lockfields large compared to the

localfields. They derived a steady state magnetization change in the detected spin-locked magnetization M_ρcompared to its initial value M0for the solid state which is given by:

M0−Mρ



=M0¼ γ2B2rotarysatπT1ρfðω1lockÞ ð1Þ

where, Brotarysatis the applied rotary saturationfield (a

T1ρis the spin-lattice relaxation in the spin-locked state and

f(ω1lock) is the Lorentzian shape function which has a value of

1 s on resonance, whenγB1lock=γBrotarysat(Abragam, 1961).

Since Eq. (1) was derived for solids and describes a steady state effect, it cannot be used to accurately model the SIRS mechanism. In order to model the effect in liquids and for transitory applications of a rotary saturationfield, we model the rotary saturation effect as a coherent rotation away from the B1lockdirection by the rotary saturationfield (i.e., the neuronal

field). The equations of motion for the magnetization in the spin-lock state are formulated in the doubly rotating frame as

dM dt ¼ γM Heff 1− ωrotarysat γjHeffj   þ Brotarysat   −R M−Mρ
; with R¼ 1 T4 2 0 0 0 1 T₂4 0 0 0 1 T1ρ 0 B B B B B B B @ 1 C C C C C C C A ; ð2Þ

where Heffis the effectivefield formed from the applied B1lock

and the off-resonance pseudo-field Δω/γ, Mρ denotes the

equilibrium magnetization in the rotating frame, which is very small since B1lockis small. Here Heffis the relevant axis toward

which the magnetization is at equilibrium, and is not the static magnetic field direction. If both the spin-lock RF and the Rotary Saturationfield Brotarysatis on resonance, then Heffis

B1lock. Brotarysatis the audio-frequencyfield applied orthogonal

to the spin-lock axis (here B1lockis taken to be along the x axis

in the singly rotating frame). To preserve the analogy to the singly rotating from, the axes of the doubly rotating frame are renamed so that the z axis is the equilibrium direction. Thus in the doubly rotating frame, the z axis is the direction B1lock.

Then Brotarysatmust be applied orthogonal to B1lockto produce

a rotation. The magnetization relaxes along Heffwith the time

constant T_1ρ, and in the other two directions with the time constant T2⁎. We propose to utilize spin-locking to sensitize

the spins to neuronal magneticfields oscillating at the Larmor frequency in the rotating frame, whereγB1lockranges between

20 and 500 Hz, as set by the applied B1lock. Oscillating

neuronal currents with spectral power at γB1lock are then

capable of producing rotary saturation of the spin-locked magnetization. Thus, the neuronal currents produce a reso-nant saturation effect on the MR signal in the rotating frame during the spin-lock state, and the spin-lock state can be tuned to be sensitive to a frequency of interest by adjusting γB1lock. In practice, we anticipate choosing γB1lockto the

ex-pected frequency of the neuronal oscillations to sensitize the sequence to these frequencies. Thus the neuronal currents induced by the stimulus create the rotary saturation effect. In analogy with the classic rotary saturation experiment of

Redfield (1955), we refer to this potential method as Stimulus-induced Rotary Saturation (SIRS).

To validate the potential of SIRS to detect local current dipoles associated with neuronal activation, we performed two different experiments. Thefirst is the measurement of rotary saturation spectra to characterize the rotary saturation effect from a current dipole. Secondly, we performed a block-design experiment in order to establish detection limits of the method compared to conventional phase-sensitive methods. Methods

The measurements were performed in a 17-cm diameter spherical phantomfilled with a solution of saline (0.9% NaCl) and 0.25 ml of Gadopentetate dimeglumine (Magnevist, Berlex Laboratories) yielding a T1measured by inversion recovery of

1100 ms. The current dipole had a length of 1 cm was formed from 0.89 mm diameter insulated copper wire exposed only at the tips. The current dipole produced by this phantom is therefore the product of the current and the dipole length. A photograph of the dipole and a schematic of the driving circuit are shown in Fig. 1. The return path for the current was therefore an ionic current distributed within the saline. The dipole was placed 3.9 cm from the edge (4.6 cm from center) of the phantom and the phantom was aligned so that the current dipole was perpendicular to the B0field. The dipole was at the

end of an insulated rigid copper coaxial cable that extended approximately 2 cm out of the phantom. Thus, with the

Fig. 1. Circuit diagram of the apparatus for applying oscillating currents to the dipole phantom during the spin-locking. The function generator is gated by a logic pulse generated by the pulse sequence so that its output is non-zero only during TSL. The output of the function generator then goes through an isolation transformer to eliminate DC and common-mode currents and drives the dipole phantom through a series resistor (10 kΩ) and series RF chokes to block RF pickup on the lines. The current in the phantom is measured by monitoring the voltage drop across the 10 kΩ series resistor using an optically isolated operational amplifier with a gain of 10 and 500 kΩ input impedance. The RMS amplitude of the output of the operational amplifier is measured with a digital oscilloscope.

exception of the two 0.5-cm portions of the dipole, the current was confined to a semi-rigid coaxial cable. Outside of the phantom, the coax was wound into a four-turn solenoid and tuned with a parallel capacitor to form a trap circuit resonant at the Larmor frequency (123 MHz) to reject RF common-mode currents that might be induced by the RF onto the dipole supply leads.

The magnetic field strengths generated by the current dipole were calculated from the current dipole moment using the equations published by Heller et al. (2004). The exact dipole location in the sphere was determined based on a high-resolution, 3D Turbo Spin Echo image of the phantom. The calculation was performed for a 2 cm × 2 cm region of interest centered at the current dipole. Since the equations ofHeller et al. (2004) assume an infinitely small dipole, we modeled our extended dipole by superposition of 1000 dipoles distributed with 10-μm spacing along the extent of our dipole. The dipole is connected to a digital frequency synthesizer (Global Specialties model 2003) via an audio transformer to eliminate DC offsets and a twisted pair cable with 1μH series RF chokes to prevent the conduction of signals at RF frequencies. The frequency syn-thesizer uses a gating pulse controlled by the MR pulse se-quence in order to supply current only during the spin-lock portion of the sequence. The current through the phantom was determined by measuring the voltage across a 10-kW series resistor using an optically isolated differential probe with ap-proximately 500 kW input impedance.

A 3 T whole-body MR scanner (Siemens Medical Solutions, Erlangen, Germany) was used with a transmit and receive birdcage head coil. The spin-lock prepared EPI pulse sequence (Borthakur et al., 2006) used is shown inFig. 2. It consists of a 90° hard pulse in y direction, immediately followed by a rectangular spin-lock pulse in x direction, which remains on for the duration TSL. The spin-lock period is followed by a hard 90° pulse in−y direction which returns the spin-locked mag-netization to the longitudinal axis for detection by the stan-dard gradient echo EPI sequence. All remaining transverse magnetization is spoiled by immediately following spoiler gradients. The entire spin-lock preparation is non-selective. Following the preparation a single slice EPI sequence is used to

acquire the T1ρ-weighted image (64 × 64 matrix, 20-cm FOV,

3.14 × 3.1 mm, 4-mm thickness with a TE of 19 ms, TR of 1 s, and a BW of 2298 Hz). For the studies addressing the sensitivity of the technique to dipolarfields, the spin-lock duration, TSL, was 100 ms. In order to determine the T1ρof the phantom, TSL was

set to 0, 15, 30, 50, and 70 ms, and the resulting image intensity wasfit to the relation I=M0exp(− TSL/T1ρ).

Bloch equation simulation

The simulation starts with the magnetization aligned along the spin-lock field. During the spin-lock phase the audio-frequency saturation field acts to rotate the magnetization away from the spin-lockfield. The Bloch simulator tracks all three components of the magnetization in Eq. (2), however only the component along the spin-lockfield is of interest since this component will be returned to the z axis for subsequent imaging. The other components of the magnetiza-tion are dephased by the spoiler gradients following the 90° pulse. We simulated the Bloch equations in the doubly rotating frame using an MATLAB code solver with parameters matching the phantom (T1= 1.1 s, T2= 0.95 s, T1ρ= 0.98 s, T2⁎=50 ms) and

human brain gray matter at 1.5 T (T1= 0.9 s, T2= 80 ms, T1ρ=

147 ms, T2⁎=50 ms). Here T1ρwas taken from the R2=1/T2and

R1= 1 / T1values using the equation T1ρ= 2 / (R2+ R1) (Abragam,

1961). The rotary saturation spectra simulate the case where a Brotarysatfield of 1 nT and 5 nT is swept in frequency past the

γB1lockresonance value. Spectra were obtained for B1lock= 10,

20 and 40 Hz for both the phantom and human relaxation parameters with a spin-lock duration of 100 ms.

Rotary saturation spectra

In order to measure the rotary saturation spectra, the frequency of the synthesizer was swept from below resonance to above the expected rotary saturation resonance in 1 Hz steps. Each frequency was repeated for 20 images. In order to exclude potential systematic errors, the frequencies were alternately swept up from 20 Hz to 60 Hz or down from 60 Hz to 20 Hz. The spectra taken withω1lockset to 10 Hz and 20 Hz

were swept from 2 Hz to 60 Hz. The entire rotary saturation spectrum measurement was repeated for each B1lock

ampli-tude in order to experimentally confirm the resonance frequency in the rotating frame. Due to capacitive effects in the phantom, the impedance of the phantom, and therefore the current amplitude, depends on the frequency of the applied current. The currents reported for the rotary saturation spectra refer to the current measured when the synthesizer was set to the resonance frequencyγB1lock. Since the other

measurements were performed at afixed frequency this effect was a limited disruption. For a 50-Hz driving source, the impedance of the phantom was measured to be 950Ω. Block-design experiment

In the block detection experiment the gating of the function generator was activated to produce the oscillating current for “on” blocks and deactivated for “off” blocks. Thus, a driving current was present only during the spin-locking period of the “on” blocks. The function generator was set to the resonance frequency previously determined with a rotary saturation spectrum measurement. AlthoughγB1lockwas set to be 40 Hz

based on the scanner's transmit adjustment procedure, the Fig. 2. Pulse sequence diagram for the spin-lock prepared EPI sequence. After the initial

90° excitation, the magnetization is spin-locked along the y axis by a continuous wave RF pulse applied for a time TSL. Immediately following the spin-lock period, TSL, the magnetization is rotated back to the z axis by the 90°− xpulse. Then the traditional EPI

sequence images the magnetization which experienced the spin-lock state. T1ρ

relaxation processes during TSL are therefore expected to diminish the signal intensity in the EPI image. The sinusoidal current input to the phantom is applied during the spin-lock period and is shown schematically in the bottom trace.

minimum was detected at 37 Hz. All block-design measurements were thus performed with the function generator set to 37 Hz.

The on-resonance and off-resonance blocks lasted for 20 TR periods of 1 s each and were repeated 10 times during each run. Total run time was 6 min, 40 s for the blocks plus 20 dummy scans to ensure that the magnetization was in steady state. For high current levels two runs were collected, for the low current levels six runs were collected and analyzed together. Statistical significance was established using a two-sided t-test between the“on” and “off” blocks for each voxel in the image.

Block detection experiments were performed at function generator voltage settings of 4 V, 0.4 V and 0.2 V, yielding a

measured current (from the voltage drop across the series 10 kΩ resistor) of 118 μA, 11.5 μA, and 5.6 μA, respectively, corresponding to current dipoles of 1180 nAm, 115 nAm, and 56 nAm respectively.

Results

Spin-lock images

Fig. 3 shows representative spin-lock images for both the phantom withω1lock= 40 Hz as well as a human image obtained

withω1lock= 40 Hz and 10 Hz.

Fig. 3. Representative spin-locked prepared images from the sequence ofFig. 2. (A) phantom image spin-locked withγB1lock= 40 Hz shown with the location of the current dipole

overlaid in yellow. (B) Time series standard deviation from a time series of 420 of the spin-lock prepped images. (C) Brain image obtained with the spin-lock prepared sequence with γB1lock= 40 Hz and D)γB1lock= 10 Hz.

Fig. 4. Magnitude of the z-component of the magneticfield associated with the current dipole phantom. Left, axial image of the z-component of the field for a 5.6 μA current in the dipole (56 nAm dipole) based on the distributed dipole model. Afield of 1 nT is estimated for a distance of 1.5 mm from the wire. Right, plot along the y axis showing both the field values derived from the distributed dipole model and the Biot–Savart law for a long thin wire.

Magneticfield calculation

Magneticfield calculationFig. 4shows the magnitude of the z component of the magneticfield associated with a 56 nAm current dipole shown in an axial plane centered on the dipole.

Fig. 3also shows a plot of magneticfield of the dipole along with a reference plot derived from the Biot–Savart law for a thin, long straight wire. Thefield is calculated to be 1 nT at a distance of 1.5 mm from the dipole. At this distance the 1 cm extended dipole looks approximately like an infinite wire. Bloch simulation

Fig. 5shows the simulation results for the phantom and the human gray matter relaxation parameters. Rotary saturation spectra are shown forγB1lockof 10, 20 and 40 Hz

and applied audio-field strengths of 1 nT and 5 nT. In each case, as the frequency of the simulated B1rotarysat field is

swept past the resonance condition atγB1lock, a reduction of

signal is observed which peaks when B1rotarysat=γB1lock. The

line width and degree of saturation in the Bloch simulations were found to be independent of γB1lock. The effect is

roughly linear with Brotarysatfor small Brotarysat. For neuronal

oscillations at 10, 20 and 40 Hz, a Brotarysat field of 6.8 nT

reduced the magnetization by 0.1% for the 100 ms spin-lock duration experiment using gray matter relaxation param-eters at 1.5 T. The neuronalfield causing a 0.1% magnetization change decreased to 4.0 nT and 2.2 nT for spin-lock durations of 200 ms and 500 ms respectively.

Rotary saturation spectra

The phantom data (Fig. 6) show a decrease in image ampli-tude when the applied current in the phantom oscillates close to the Larmor frequency in the rotating frame, γB1lock. The

observed line width of the saturation spectra is about 6 Hz FWHM. The peak often appeared at a slightly (2 Hz) higher frequency than expected from the B1lock level determined

from the instrument's standard transmit calibration, possibly due to the Bloch–Siegert shift (Bloch and Siegert, 1940) or a mis-calibrated B1 level. For high dipole current strengths

(400μA), the rotary saturation effect achieved a 30%

ampli-tude decrease. For a small dipole current (16μA), the shape of the rotary saturation spectrum can still be observed in voxels immediately neighboring the dipole (Fig. 4). The amplitude minimum still occurs close to the selected Larmor frequency Fig. 5. Simulation of the rotary saturation effect based on the Bloch equations in the rotating frame. The simulations show the fraction of the spin-locked magnetization as a function of the frequency of the applied Brotarysatfield. When the frequency of Brotarysatreaches the resonance condition in the rotating frame (equal toγB1lock) the spin-locked magnetization

is rotated away from the spin-lockfield B1lock. In the spin-lock prepared experiment, only the spin-lock magnetization is returned to the z axis for subsequent imaging. Simulations

were performed for human gray matter at 1.5 T (T1= 900 ms, T2= 80 ms, T2⁎=50 ms) and for the phantom used in this study (T1= 1100 ms, T2= 950 ms, T2⁎=100 ms). The effect is

shown for rotary saturation atωrotarysatof 10 Hz, 20 Hz and 40 Hz (black, blue and red curves respectively) using Brotarysatof 1 nT (solid line) and 5 nT (dashed line).

Fig. 6. Image intensity as a percentage of the mean near the current dipole as a function of the frequency of the applied current to the dipole. (A) Intensity changes in a single voxel adjacent to the dipole using a 400μA applied current (4 μAm dipole strength),10 averages. The black curve shows the rotary saturation effect with the spin-lockfield set to γB1lock=

40 Hz and the red trace shows a similar measurement withγB1lock= 50 Hz. (B) Intensity

changes in a 5 voxel ROI immediately below dipole using a16μA applied current (160 nA dipole strength), 160 averages. The red curve shows the rotary saturation effect with the spin-lockfield set to γB1lock= 40 Hz and the black trace shows a 5 pixel ROI distant from the

in the rotating frame as expected. The observed line width of the rotary saturation spectra is also about 6 Hz.

Block detection experiment

The t-test maps from the block-design runs at each of the three current levels are shown inFig. 6. The t-test results are overlaid on a high-resolution Turbo Spin Echo image, but an example of an individual spin-lock-prepped echoplanar image is shown inFig. 3A. The two runs performed with 4 V function generator setting (118μA current through the dipole) resulted in a large number of significant (pb0.001) voxels above and below the dipole, in approximate agreement with the expected distribution of the z component of the dipolarfield. The six runs performed with 0.4 V function generator setting (11.5μA in the dipole) showed significant voxels only within 5 mm above and below the dipole. The lowest current experiment performed with a 0.2 V function generator setting (5.6μA in the dipole) yielded three significant voxels, all immediately below the dipole. The percentage signal change observed was 0.15% for the 5.6μA current in the dipole.

Discussion

The validation of this methodology is based on the equivalent current dipole model which is commonly utilized in MEG data analysis. The rationale is that postsynaptic currents occurring synchronously in a population of tens of thousands of cortical pyramidal cells can be modeled as a single current dipole because the cells themselves are almost perfectly par-allel to each other and thus the dipolarfields of all the neurons add constructively. From the farfield it is impossible to dis-tinguish between a single larger dipole or a group of parallel small dipoles. Thus, this phantom study is an appropriate benchmark for the MEG experiment. Here our phantom tries to match the magnetic fields observed in MEG using a single artificial dipole. For the MR experiment, however, the dipolar fields are detected locally by their effect on the MR signal, in our case due to their rotary saturation effect on neighboring water. In this case, the relevant magneticfields are not the far-fields of the dipole, as in MEG, but the localfields among the neurons. Nonetheless, a single current dipole provides an easy-to-construct approximation which at least matches the neural situation very closely in the farfield but likely only to an order of magnitude in the localfields of activated cortex.

While we show that the rotary saturation method appears to have adequate sensitivity to local magneticfields using a current dipole phantom, T2⁎-based phase and amplitude

me-thods are also adequate by these measures. For example, Konn et al. are able to see a minimum dipole moment of 6.3 nAm (current of 1.3μA and dipole size of 5 mm) using the phase shift in EPI in a phantom similar to ours using 62 min of averaging (Konn et al., 2003). Normalizing for averaging time, this corresponds to a dipole moment of 77 nAm (10μA current) in a minute of scan time.Pell et al. (2006) used the magnitude changes and observed a minimum current of 1.7μA in a 1.8 min run; a normalized sensitivity of 2.3μA in a minute of scan time. The block-design results shown in Fig. 7 show a non-Bonferroni- corrected p-value of 10− 7 for a 56 nA dipole (5.6μA current) in 40 min of imaging time. However, unlike the phase and magnitude EPI experiments, the signal change in the SIRS experiment does not scale linearly with current. The method detects the reduction in magnetization as the spin-locked magnetization is rotated away from the spin-lock axis, thus scales as 1−cos(φ)≈φ2_{/ 2 whereφ is the small tip angle}

induced by the rotary saturationfield in the doubly rotating frame. Sinceφ is proportional to Brotarysat, and thus the current,

the effect is expected to scale with the square of the current. Thus, for the SIRS experiment a 1 min experiment would be expected to detect a 2.5 fold smaller current than the 40 min experiment. This suggests the ability to detect a 22.4 nAm dipole or 2.2μA current in the normalized 1 min of imaging time. While larger than that of the other two studies, our result occurred with a Bonferroni corrected p-value of approximately 2 × 10− 5compared to the threshold for statistical significance of p = 10− 3used in the Pell result.

The failure of the T2⁎-weighted imaging to detect actual

neuronal currents in human experiments suggests that the single current dipole phantom does not adequately mimic the local phenomena near the neurons. Our best insight into whether or not the SIRS method (or other MR methods) is sufficiently sensitive to allow detection of neuronal fields comes from microSQUID measurements of evoked fields near the cortex. Okada et al. has measured the magneticfields 1.2 mm above the cortex of juvenile swine (Okada et al., 1999). They measurefields in excess of 1 pT in response to snout stimulation. While the current dipole model used provides a diffuse current return path similar to that expected in vivo, practical implementation of the phantom causes problems which make it difficult to use DC current measurements. When a DC current

Fig. 7. Two-sided t-test results (not Bonferroni corrected) from the block-design paradigm shown with a lower threshold of pb10− 3_{and a color scale saturation level of pb10}− 8_{. All}

block designs were 20 s per period, andγB1lock= 40 Hz. T-test results are overlaid on a high-resolution Turbo Spin Echo (TSE) image which allows visualization of the current dipole's

wires. A) 118μA applied current, B) 11.5 μA applied current, C) 5.6 μA applied current. The statistical parametric map is overlaid over a high-resolution fast spin echo image of the phantom to visualize the location of the dipole wire.

is used, electrolysis at the electrode tips occurs causing gas bubbles and alteration of the current. Eventually the tips must be sanded. Although platinum exhibits reduced electrolysis, unlike copper it is severely diamagnetically different from water and altered the EPI images. Thus, we did not perform a direct comparison to the DC current measurements in the literature using T2⁎-weighted EPI phase and amplitude

mea-sures and rely on the extensive literature for the sensitivity of these methods.

A major difference between the macroscopic phantom's local fields and those within the microscopic distribution of neurons is likely caused by the spatial cancellation effects in bothfield and MR phase arising from the multiple parallel neuronal currents in neighborhood of activated cortex. We break this down into two effects. First is the cancellation of the neuronal magneticfields at some spatial locations within the distribution of synchronous neurons. For water at spatial locations where this cancellation occurs, it is unlikely that any MR method will be able to detect a signal change. For example, a region with uniformly spaced, current-carrying discrete wires with aligned currents, the magneticfield cancels at spatial locations half way between the wires and is reduced near these locations. The second concern is based on how the survivingfields affect the MR signal. The signal changes are determined by integrating the local effects over both time and space. Spatially, the magnitude and phase changes in a T2⁎-weighed experiment are determined

by integrating over the voxel. In this case, the effect can be reduced by the symmetry of the voxel with respect to the current distribution. For example, a single current-carrying wire centered on the voxel will not produce a net phase shift, and will exhibit a reduced magnitude reduction. Park et al. (2007)

estimate the signal and magnitude changes for a region of activated cortex larger than a voxel. In this case, the voxels near the edge of the activated region have the strongest effect.

Considering time-domain effects is important for neuronal currents which reverse polarity during the measurement period. In the case of a symmetric multi-phasic current time course, there can be complete elimination of the accrued phase shift in the MR signal. MEG experiments suggest that at least partial cancellation will occur for many stimuli. This effect can be reduced by carefully timing the MR pulse sequence relative to the evoked response so that the phase-sensitive part of the MR experiment only occurs during one phase of the response, or exactly straddles the refocusing pulse of a spin echo experi-ment (Singh, 1994).

An additional source of cancellation of the neuronal mag-neticfields can derive from a lack of phase coherence among the neuronal currents. This leads to incomplete superposition among the individual neuronsfields. This is a problem com-mon to all neuronal magneticfield measurements including MEG since the distal nature of the MEG detectors is unlikely to help a lack of phase superposition. Although phase cancella-tion among the localfields undoubtedly occurs, a 1 nT estimate for the strength of local neuronalfields within millimeters of the activation is supported by MicroSQUID measurements close to the cortex that indicate the presence of neuronal magneticfields above 1 nT 1 mm above active cortex (Okada et al., 1999; Ikeda et al., 2005).

The rotary saturation approach presented here cannot offer improved detection in areas where thefields from the local neuronal currents cancel, but likely has advantageous spatial and temporal averaging properties relative to MR phase mea-surements and T2⁎-weighted imaging. The saturation method

does not suffer from intra-voxel cancellation due to the phase of the saturationfield; local water molecules will be saturated by the neuronalfields regardless of the phase of these fields, and the averaged effect over the voxel or time period of the spin-locking will always further reduce the signal intensity.

A further limitation that is common to both our approach and the T2⁎-based methods is that the neuronal magnetic field

needs to be parallel to B0in order to have an effect. When

investigating the cerebral cortex and assuming large pyrami-dal neurons to be the main generator of the neuronal magnetic field, only certain regions of the cortex satisfy this condition for a given head position, thus, in an in vivo experiment the target area needs to be carefully selected and the subjects head ap-propriately positioned. A similar “orientation blindness” occurs in MEG which is insensitive to radial dipoles. For the Rotary Saturation experiment, the neuronal fields must be oriented orthogonal to the spin-lockfield. If needed the ex-periment could be repeated with the spin-lock along the other axis to eliminate the orientation blind spot.

Furthermore, there are additional limitations to an in vivo rotary saturation experiment. The spin-locking period is limited in time by both SAR restrictions and the T_1ρtime. The method, as implemented, was limited to single slice measure-ments, although multislice spin-locking sequences have been developed for other purposes (Wheaton et al., 2004). In a mo-tional narrowed regime, such as the bulk water that is usually imaged in MRI, the T_1ρtime is almost identical to T2when the

spin-lock is performed on resonance (Moran and Hamilton, 1995). Practically this limits the spin-lock duration to about 100 ms. Furthermore, the selected Larmor frequency in the rotating frame has an upper limit due to the SAR restrictions and a lower limit imposed by the ability to create a spin-lock. However, our initial tests with in vivo imaging have found that Larmor frequencies in the rotating frame as low as 10 Hz are feasible (seeFig. 3). Additionally, the method is currently sen-sitive to only a single frequency of neuronal currents but could, in principle, be extended to at least several frequencies.

We designed our phantom experiment to demonstrate the SIRS effect in the frequency range of evoked responses com-monly observed in MEG/EEG measurements of the human brain. Virtually all magnetic fields originating from postsy-naptic currents have frequencies of less than 120 Hz, with the bulk offields being below 20 Hz. Examples for experiments that evoke narrow-band oscillations are auditory click-trains (Lin et al., 2004), somatosensory stimulation with rapid elec-trical or mechanical pulse trains (Prevec and Ribaric-Jankes, 1996) and visual stimulation with periodic rapidly changing intensity (Muller et al., 1997). It is also possible to evoke or induce so-called gamma oscillations (30 Hz–60 Hz) with spe-cific stimuli, which will give rise to narrow-band magnetic fields occurring often from the sensory cortex associated with the modality.

Contamination by the BOLD effect is a problem common to all MR neuronal current methods and several studies have attempted to address this (Bandettini et al., 2005). For exam-ple, eliminating BOLD contamination might be achieved by timing to make the experiment sensitive only to the early response and not the later BOLD signal. This requires a long Inter-Stimulus Interval (ISI) and renders the experiment inefficient. A second approach is to saturate the BOLD effect by attempting to have similar BOLD responses for all condi-tions. For the SIRS method this might involve altering the stimulus frequency by approximately 10 Hz so that the Rotary

saturation effect is on resonance for one condition and off resonance for the other. This gives full SIRS contrast but little BOLD changes if the two stimuli have similar BOLD response.

A third method is to attempt to reduce the BOLD effect by performing the study at low B0 fields. Unlike BOLD, most

neuronal current methods (including SIRS) are expected to have little dependence on B0. The contrast itself is not

altered for SIRS except through secondary effects such as relaxation time changes. Since T2and thus T1ρare expected

to lengthen at lower B0, SIRS contrast should slightly

im-prove with experiments at lower B0field. Similarly the

re-duced Specific Absorption Rate (SAR) constraints will allow longer TSL and more sensitive detection. The time series SNR is expected to have only modest decreases as the B0 is

lowered due to the dominance of physiological noise at high field for conventional voxel volumes (Triantafyllou et al., 2005). Finally, there is some hope that pharmacologic preparation might help for neuronal current studies in animal models. Here agents are being sought which act only on astrocytes and affect theflow response but not electro-physiology (Baslow et al., 2005).

The sensitivity of the method is limited by T1ρand the

spin-lock duration TSL. T1ρis influenced by many parameters

in-cluding the B0field, the spin-lock field amplitude (B1lock), the

temperature, and the tissue type.Borthakur et al. (2006)found T_1ρ for gray and white matter to be 99 ms and 85 ms respectively atω1lock= 500 Hz at 1.5 T. The corresponding T2

values were 75 ms and 66 ms. Thus, the T_1ρin the brain is similar in magnitude to T2, but not identical.

Conclusion

We introduce a new method by which neuronal currents might be construed to affect MR signals. We demonstrate that the method has the sensitivity to detect the magneticfields of dipole current sources of similar magnitude and frequencies to that observed in MEG experiments. The method is an analog to the excitation of MR signal by neuronalfields proposed by

Kraus (2006)for ultra-low-field MRI, but has the advantage of being applicable to standard high-field MR instrumentation. While we have not yet demonstrated its ability to detect neuronalfields in vivo, the method has many favorable spatial and temporal averaging properties which may be holding back conventional MR phase and T2⁎ magnitude-based studies.

Acknowledgments

We would like to acknowledge Dr. Joseph Mandeville and Dr. Jonathan Polimeni for the stimulating discussions and helpful assistance with the implementation of the experiment. This project wasfinancially supported by NIH/NIDCD 5T32DC000038-15, NCRR P41RR14075, NIH 1R21EB007298, NSC 96-2320-B-002-085, NHRI-EX97-9715EC and the MIND Institute.

References

Abragam, A., 1961. Principles of Nuclear Magnetism. Oxford University Press. Bandettini, P.A., Petridou, N., Bodurka, J., 2005. Direct detection of neuronal activity with

MRI: fantasy, possibility, or reality? Appl. Magn. Reson 29.

Baslow, M.H., Dyakin, V.V., Nowak, K.L., Hungund, B.L., Guilfoyle, D.N., 2005. 2-PMPA, a NAAG peptidase inhibitor, attenuates magnetic resonance BOLD signals in brain of

anesthetized mice: evidence of a link between neuron NAAG release and hyperemia. J. Mol. Neurosci. 26, 1–15.

Bloch, F., Siegert, A., 1940. Magnetic resonance for nonrotatingfields. Phys. Rev. 57, 522. Bodurka, J., Bandettini, P.A., 2002. Toward direct mapping of neuronal activity: MRI detection of ultraweak, transient magneticfield changes. Magn. Reson. Med. 47, 1052–1058.

Bodurka, J., Jesmanowicz, A., Hyde, J.S., Xu, H., Estkowski, L., Li, S.J., 1999. Current-induced magnetic resonance phase imaging. J. Magn. Reson. 137, 265–271. Borthakur, A., Hulvershorn, J., Gualtieri, E., Wheaton, A.J., Charagundla, S., Elliott, M.A.,

Reddy, R., 2006. A pulse sequence for rapid in vivo spin-locked MRI. J. Magn. Reson. Imaging 23, 591–596.

Chow, L.S., Cook, G.G., Whitby, E., Paley, M.N., 2006. Investigating direct detection of axonfiring in the adult human optic nerve using MRI. NeuroImage 30, 835–846.

Chu, R., de Zwart, J.A., van Gelderen, P., Fukunaga, M., Kellman, P., Holroyd, T., Duyn, J.H., 2004. Hunting for neuronal currents: absence of rapid MRI signal changes during visual-evoked response. NeuroImage 23, 1059–1067.

Hämäläinen, M., Hari, R., Ilmoniemi, R.J., Knuutila, J., Lounasmaa, O.V., 1993. Magnetoencephalography; theory, instrumentation, and applications to noninva-sive studies of the working human brain. Rev. Mod. Phys. 65, 413.

Hari, R., 1991. On brain's magnetic responses to sensory stimuli. J. Clin. Neurophysiol. 8, 157–169.

Heller, L., Ranken, D., Best, E., 2004. The magneticfield inside special conducting geometries due to internal current. IEEE Trans. Biomed. Eng. 51, 1310–1318. Ikeda, H., Wang, Y., Okada, Y.C., 2005. Origins of the somatic N20 and high-frequency

oscillations evoked by trigeminal stimulation in the piglets. Clin. Neurophysiol. 116, 827–841.

Konn, D., Gowland, P., Bowtell, R., 2003. MRI detection of weak magneticfields due to an extended current dipole in a conducting sphere: a model for direct detection of neuronal currents in the brain. Magn. Reson. Med. 50, 40–49.

Konn, D., Leach, S., Gowland, P., Bowtell, R., 2004. Initial attempts at directly detect-ing alpha wave activity in the brain usdetect-ing MRI. Magn. Reson. Imagdetect-ing 22, 1413–1427.

Kraus Jr., R.H., 2006. Lowfield SQUID MRI. Workshop on Multi-modal Functional Neuroimaging. Cortona, Italy.

Kraus Jr., R.H., Volegov, P., Matlachov, A., Espy, M., 2008. Toward direct neural current imaging by resonant mechanisms at ultra-lowfield. NeuroImage 39, 310–317. Lin, F.H., Witzel, T., Hämäläinen, M.S., Dale, A.M., Belliveau, J.W., Stufflebeam, S.M., 2004.

Spectral spatiotemporal imaging of cortical oscillations and interactions in the human brain. NeuroImage 23, 582–595.

Matlachov, A.N., Volegov, P.L., Espy, M.A., George, J.S., Kraus Jr., R.H., 2004. SQUID detected NMR in microtesla magneticfields. J. Magn. Reson. 170, 1–7.

Moran, P.R., Hamilton, C.A., 1995. Near-resonance spin-lock contrast. Magn. Reson. Imaging 13, 837–846.

Muller, M.M., Teder, W., Hillyard, S.A., 1997. Magnetoencephalographic recording of steady-state visual evoked cortical activity. Brain Topogr. 9, 163–168.

Niedermeyer, E., Lopes da Silva, F.H., 1999. Electroencephalography: Basic Principles, Clinical Applications, and Related Fields, 4th ed. Williams & Wilkins, Baltimore. Okada, Y., Lahteenmaki, A., Xu, C., 1999. Comparison of MEG and EEG on the basis of

somatic evoked responses elicited by stimulation of the snout in the juvenile swine. Clin. Neurophysiol. 110, 214–229.

Park, T.S., Lee, S.Y., 2007. Effects of neuronal magneticfields on MRI: numerical analysis with axon and dendrite models. NeuroImage.

Parkes, L.M., de Lange, F.P., Fries, P., Toni, I., Norris, D.G., 2007. Inability to directly detect magneticfield changes associated with neuronal activity. Magn. Reson. Med. 57, 411–416.

Pell, G.S., Abbot, D.F., Fleming, S.W., Prichard, J.W., Jackson, G.D., 2006. Further steps toward direct magnetic resonance (MR) imaging detection of neuronal action currents: optimization of MR sensitivity to transient and weak currents in a conductor. MRM 55 (5), 1038–1046.

Petridou, N., Plenz, D., Silva, A.C., Loew, M., Bodurka, J., Bandettini, P.A., 2006. Direct magnetic resonance detection of neuronal electrical activity. Proc. Natl. Acad. Sci. U. S. A. 103, 16015–16020.

Prevec, T.S., Ribaric-Jankes, K., 1996. Can somatosensory system generate frequency following response? Pflugers Arch. 431, R301–R302.

Redfield, A.G., 1955. Nuclear magnetic resonance saturation and rotary saturation in solids. Phys. Rev. 98, 1787.

Singh, M., 1994. Sensitivity of MR phase shift to detect evoked neuromagneticfields inside the head. IEEE Trans. Nucl. Sci. 41, 349–351.

Tallon-Baudry, C., Bertrand, O., Fischer, C., 2001. Oscillatory synchrony between human extrastriate areas during visual short-term memory maintenance. J. Neurosci. 21, RC177.

Triantafyllou, C., Hoge, R.D., Krueger, G., Wiggins, C.J., Potthast, A., Wiggins, G.C., Wald, L.L., 2005. Comparison of physiological noise at 1.5 T, 3 T and 7 T and optimization of fMRI acquisition parameters. NeuroImage 26, 243–250.

Volegov, P., Matlachov, A.N., Espy, M.A., George, J.S., Kraus Jr., R.H., 2004. Simultaneous magnetoencephalography and SQUID detected nuclear MR in microTesla magnetic fields. Magn. Reson. Med. 52, 467–470.

Wheaton, A.J., Borthakur, A., Charagundla, S.R., Reddy, R., 2004. Pulse sequence for multislice T1rho-weighted MRI. Magn. Reson. Med. 51, 362–369.

Xiong, J., Fox, P.T., Gao, J.H., 2003. Directly mapping magneticfield effects of neuronal activity by magnetic resonance imaging. Hum. Brain Mapp. 20, 41–49.