Unipolar MR elastography: Theory, numerical analysis and implementation

<jats:p>In MR elastography (MRE), zeroth moment balanced motion‐encoding gradients (MEGs) are incorporated into MRI sequences to induce a phase shift proportional to the local displacement caused by external actuation. To maximize the signal‐to‐noise ratio (SNR), fractional encoding is employed, i.e., the MEG duration is reduced below the wave period. Here, gradients encode primarily the velocity of the motion‐reducing encoding efficiency. Thus, in GRE‐MRE, T<jats:sub>2</jats:sub>* decay and motion sensitivity have to be balanced, imposing a lower limit on repetition times (TRs).</jats:p><jats:p>We propose to use a single trapezoidal gradient, a “unipolar gradient”, to directly encode spin displacement. Such gradients cannot be used in conventional sequences as they exhibit a large zeroth moment and dephase magnetization. By time‐reversing a spoiled SSFP sequence, the spoiling gradient becomes an efficient unipolar MEG. The proposed “unipolar MRE” technique benefits from this approach in three ways: first, displacement encoding is split over multiple TRs increasing motion sensitivity; second, spoiler and MEG coincide, allowing a reduction in TR; third, motion sensitivity of a typical unipolar lobe is of an order of magnitude higher than a bipolar MEG of equal duration.</jats:p><jats:p>In this work, motion encoding using unipolar MRE is analyzed using the extended phase graph (EPG) formalism with a periodic motion propagator. As an approximation, the two‐transverse TR approximation for diffusion‐weighted SSFP is extended to incorporate cyclic motion. A complex encoding efficiency metric is introduced to compare the displacement fields of unipolar and conventional GRE‐MRE sequences in both magnitude and phase. The derived theoretical encoding equations are used to characterize the proposed sequence using an extensive parameter study. Unipolar MRE is validated against conventional GRE‐MRE in a phantom study showing excellent agreement between measured displacement fields. In addition, unipolar MRE yields significantly increased octahedral shear strain‐SNR relative to conventional GRE‐MRE and allows for the recovery of high stiffness inclusions, where conventional GRE‐MRE fails.</jats:p>

section_content: In MR elastography (MRE), zeroth moment balanced motion-encoding gradients (MEGs) are incorporated into MRI sequences to induce a phase shift proportional to the local displacement caused by external actuation. To maximize the signal-tonoise ratio (SNR), fractional encoding is employed, i. e., the MEG duration is reduced below the wave period. Here, gradients encode primarily the velocity of the motion-reducing encoding efficiency. Thus, in GRE-MRE, T 2 * decay and motion sensitivity have to be balanced, imposing a lower limit on repetition times (TRs). 

We propose to use a single trapezoidal gradient, a "unipolar gradient", to directly encode spin displacement. Such gradients cannot be used in conventional sequences as they exhibit a large zeroth moment and dephase magnetization. By time-reversing a spoiled SSFP sequence, the spoiling gradient becomes an efficient unipolar MEG. 

The proposed "unipolar MRE" technique benefits from this approach in three ways: first, displacement encoding is split over multiple TRs increasing motion sensitivity; second, spoiler and MEG coincide, allowing a reduction in TR; third, motion sensitivity of a typical unipolar lobe is of an order of magnitude higher than a bipolar MEG of equal duration. 

In this work, motion encoding using unipolar MRE is analyzed using the extended phase graph (EPG) formalism with a periodic motion propagator. As an approximation, the two-transverse TR approximation for diffusion-weighted SSFP is extended to incorporate cyclic motion. A complex encoding efficiency metric is introduced to -------------------------------------------------------------------------------------------------------------------------------- 1 | INTRODUCTION Magnetic resonance elastography (MRE) allows the in vivo determination of the local tissue shear stiffness by time-resolved measurement of the 3D displacement field induced by an external actuator. 1 MRE utilizes the principle of phase contrast to encode the cyclic spin motion in the image phase. 2 Here, zeroth moment balanced gradients, ie, motion-encoding gradients (MEGs), are introduced before the readout to sensitize the image phase to displacements of the magnetization occurring while the gradient is switched on. Most MRE experiments currently rely on either spin echo (SE) or gradient-recalled echo (GRE) sequences. 1, 3, 4 Here, SE-based MRE benefits from the long echo time (TE) that allows maximization of the accumulated phase through matching of the MEG period to the wave period. Optimal motion sensitization with multiple gradient lobes and placement of MEGs before and after the echo pulse is possible. 4 However, fast readout strategies including echo-planar imaging (EPI) or spiral imaging need to be employed at the cost of T 2 * or frequency offset-related spatial blurring, eddy current-induced ghosting, or B0-inhomogeniety-induced geometric distortions. 5 In GRE-based MRE, a short Cartesian readout can be used due to the much shorter repetition time (TR). 7] [8] In addition, the loss of longitudinal magnetization due to the short TR as well as the use of RF-spoiling to reduce spin-history effects leads to additional reduction in SNR. 9, 10 [20] Instead of spoiling magnetization after the readout, refocusing gradients can be used to effectively reuse transverse magnetization in subsequent repetition periods. This balanced steady-state free precession (bSSFP) approach has been exploited for MRE, which greatly benefits from its high SNR. 6, 21, 22 However, the continuous refocusing of transverse magnetization in bSSFP is dependent on phase accrual during the TR interval, which leads to the typical banding artifacts in the presence of off-resonance. In the context of MRE, the interference of motion phase and inhomogeneity leads to nonlinearity in the encoding of the vibratory motion, rendering bSSFP-based MRE difficult for in vivo applications. 6, 21, 22 Nonetheless, bSSFP-based MRE has been reported for hepatic and cardiac MRE. 6, 23 5] [26] [27] These methods exploit additional refocusing gradients with integer multiples of the spoiling area before the readout of a conventional GRE sequence with consecutive balancing after the readout to refocus previously spoiled magnetization states while keeping the total spoiling area constant. These spoiling states, termed "configuration states", are carried over to the next TR period, where the subsequent RF pulse then partly flips the magnetization back to the longitudinal axis (storage), leaves it unaffected (dephasing configuration), or inverts it (rephasing configuration). 28 This leads to a plethora of discrete dephasing and rephasing configuration states, which can be refocused by a gradient that counteracts their accumulated spoiling. This work proposes a new MRE sequence, which is based on the readout of the SSFP SE using a time-reversed spoiled SSFP sequence. 29 Here, the spoiling gradient acts as a highly efficient MEG, which allows considerable reduction in TRs as well as encoding fractions. In addition, 3D slab encoding leads to increased SNR compared with conventional consecutive or interleaved multi-slice (MS) acquisitions. 30 In contrast to conventional GRE-or SE-based MRE sequences, which can be solely understood on the basis of phase contrast, the proposed technique depends on spin history effects and necessitates the analysis of echo pathways and their motion sensitivities. The objective of this work is 5-fold: first, a qualitative theoretical understanding of echo formation and motion encoding in unipolar MRE is provided; second, a periodic motion propagator in the presence of gradient fields is derived for the extended phase graph (EPG) formalism, which allows for the numerical, quantitative analysis of the motion encoding in unipolar MRE using EPGs; third, a complex encoding efficiency is proposed, which takes phase differences in encoding into account and hence allows direct comparison of the complex displacement field in both magnitude and phase of conventional and unipolar MRE despite their vast difference in encoding concept and sequence timing; fourth, a simplification of the full EPG-based simulation is attempted through adaptation of the two-transverse TR approximation, which was originally derived for DW-SSFP, allowing for fast encoding efficiency estimation 27, 31 ; and fifth, a comprehensive phantom study comparing unipolar MRE and conventional GRE-MRE is provided in order to demonstrate the feasibility of unipolar MRE and to validate theoretical findings.

section_content: In GRE-based MRE, a spoiling gradient is played after the readout of the free induction decay (FID) to dephase transverse magnetization (Fig- ure 1A, C ). We propose to read out the F -1 configuration state (SSFP-echo) instead, which is refocused by time-reversing the spoiled GRE sequence (Figure 1B ), typically referred to as contrast-enhanced-FAST (CE-FAST), PSIF (Siemens; reversed FISP) or T2-CE-FFE (Philips). 34, 35 GURE 1 Depiction of echo formation in (A) spoiled (GRE) and (B) time-reversed spoiled SSFP. While the readout in spoiled GRE sequences is dominated by the FID, the first readout in time-reversed SSFP does not acquire any signal as the gradient spoils the signal of the FID. The second RF pulse splits the magnetization into a rephasing and a dephasing contribution. The former is refocused by the spoiling gradient and leads to the formation of the SSFP-echo (figure adapted from Scheffler, Ref. 32). Comparison of (C) conventional GRE-MRE with fractional bipolar motionencoding gradients (MEGs) and (D) the proposed unipolar MRE sequence that uses the spoiling gradient as the MEG and reads out the SSFP-echo instead. The gray-shaded gradients are zeroth-moment balanced over one TR period. In (E), the global sequence timing of unipolar MRE is shown, which follows the 3D Ristretto scheme. 12, 33 The innermost squares denote acquisitions according to (D), which are wave-phase interleaved (N P = 5, N D = 1). After a set of dummy shots to reach steady-state, successive k-lines are acquired in both slice-and phase-encoding direction (3D, slab excitation). After the full volume is acquired, the encoding direction is changed according to the Hadamard encoding scheme. A combination of TTL triggers and a correct choice of TR is used to ensure synchronization of the external wave generator and the unipolar MRE sequence transverse magnetization. The second RF pulse partly leaves the transverse magnetization unaffected (0°component, dephasing configuration) and partly inverts it (180°component, rephasing configuration). The latter is refocused by the spoiling gradient and forms the SSFP-echo. 32 7] [38] [39] [40] Here, this increased motion sensitivity is used for highly efficient elastography acquisitions. Since the bipolar MEG is replaced by a single gradient lobe, we will refer to the sequence as unipolar MRE. 

In order to acquire N P wave-phase offsets (with N P ≥2) to extract amplitude and phase of the displacement field, the sequence timing follows the Ristretto MRE scheme 33 for 3D encoding (Figure 1E ). Contrary to conventional 3D GRE-MRE schemes, 41 Ristretto interleaves the wave-phase acquisition allowing minimization of the TR. The wave phases are acquired in an interleaved fashion in the innermost sequence loop. The repetition time is constant and given by

with N D being a co-prime of N P ( gcd N D ; N P ðÞ ¼ 1). Each successive repetition is thus offset relative to the transistor-transistor logic (TTL) trigger played every N P repetitions, where the wave phase offsets amount to

Since RF pulse and net dephasing from imaging gradients is equal in each TR and the underlying motion is periodic in N P T R, a steady-state of N P distinct echoes forms. 42 the proposed sequence, the unipolar gradients for spoiling and motion sensitization coincide, necessitating the switching of the spoiling direction during sequence acquisition. To correct for phase contributions from imaging gradients, which are independent of the encoding direction, a reference phase image without motion sensitization is typically acquired in GRE-MRE. 7, 8 However, since motion sensitization and spoiling is facilitated with the same gradient, an encoding scheme that encodes with equal total gradient moment for all encoding directions needs to be employed to ensure equal image contrast as well as sufficient spoiling. To this end, Hadamard encoding is incorporated, which encodes along the four diagonals of a regular cube. 8, 43 The gradient direction matrix is given by

where the rows correspond to the three orthogonal gradient axes, while the columns denote successive encoding directions. Since the Hadamard scheme exhibits identical gradient moments in each of the three orthogonal gradient axes, isotropic image resolution is required to fulfill the gradient spoiling condition 30 (p. 353 )

for each spatial dimension separately. To ensure proper spoiling of the SSFP signal, the spoiling direction is changed only after the full k-space for all wave phase offsets has been acquired and additional dummy shots are played for each encoding direction separately to reestablish the steadystate (see Figure 1E ).
section_name: | THEORY

section_content: In phase contrast-based MRE, the phase accrual φ in each voxel due to the periodic motion is proportional to the time-varying displacement field

Here, d

→ 0 e iϕ is the complex displacement amplitude, Im ⋅ fgthe imaginary part, and ω ¼ 2πf the angular wave frequency. The proportionality factor is denoted by the encoding efficiency ϵ of the gradient waveform, which is given by 3

where γ is the gyromagnetic ratio and G τ ðÞdenotes the gradient waveform. While the phase of the complex encoding efficiency can often be neglected in conventional MRE, we specifically account for it to allow for absolute displacement field comparisons of conventional and unipolar MRE. 

The sampling of the N P equidistant phase offsets along the m-th direction can be stated as

where φ 0;n accounts for the phase accumulated by imaging gradients and the notation Ê •;m denotes the mth gradient direction. Thus the projection of the estimated displacement d → including imaging gradient bias along the mth encoding direction is obtained by inverse discrete Fourier transformation

which can then be solved for the unbiased d → using the pseudo-inverse of the encoding matrix. 8
section_name: | Lowest order encoding efficiency of unipolar MRE: SSFP spin-echo contribution

section_content: The complex encoding efficiency for a bipolar MEG (Figure 2A ), with duration T ¼ f -1 MEG, gradient strength G, and slope time t d played at T 0 relative to the TTL-trigger is given by The encoding efficiency of a single unipolar lobe (Figure 2D ) with constant gradient strength G, infinite slew, and finite duration T can be obtained by solving Equation 6 and can be written as ϵ single lobe ¼ Δk sinc q ðÞ e iπqþiωT0 ; (10) where q ¼ Tf denotes the encoding fraction, 6 Δk ¼ γGT is the gradient area of the unipolar lobe, and sinc x ðÞ¼ sin πx ðÞ =πx is the sinc-function. 

The primary contribution to the SSFP-echo is given by the spin echo (SSFP-SE) of two successive TR periods, which is sensitized to motion twice (Figure 2B ). The part of transverse magnetization, which is refocused, is complex conjugated at the time of the second RF pulse (Figure 2B ). 

Thus, the motion sensitization of the first unipolar gradient is effectively inverted (Figure 2C ). To obtain the correct encoding phase, the reference time with respect to the TTL trigger is set to be equal to the second RF pulse, thus the initial magnetization is excited at time -T R. The lowest order approximation to the encoding efficiency of unipolar MRE assuming infinite slew is then given by

where the timing of the Ristretto scheme as specified in Equation 1 is implied. In the limit of an infinitesimal encoding gradient with finite moment Δk (Figure 2E ), the encoding efficiency does not vanish but approaches
section_name: | Bipolar MEG

section_content: Generalization of the previous infinite-slew model to the real situation of a trapezoidal gradient lobe with finite slope time t d (Figure 2C ) leads to ϵ Unipolar ¼ 2iΔk⋅ sinc q d ðÞ sinc q-q d ðÞ sin π N D N P e iπ q-NDN -1 P ðÞ þiωT0 ;

where q d ¼ t d f is the fraction of the wave period that is lost to the gradient slewing. Taking the limit t d →0, the infinite-slew encoding efficiency Equation 11 is directly recovered. By replacing T→T=2 and T R →T=2 as well as temporal shifting by T 0 →T=2, the encoding efficiency of the bipolar MEG with finite slew is obtained. 

The above encoding efficiency Equation 13 is only valid under the assumption that each received signal is solely given by the SSFP-SE contribution, which is excited and refocused by the current and previous RF pulses. However, since no recovery period is included in time-reversed spoiled SSFP sequences compared with conventional SE sequences and the sequence TR is usually much shorter than T 2, spin history effects need to be taken into account. These will affect the encoding efficiency by constructive and destructive interference.
section_name: | Finite-slew unipolar MRE

section_content: In order to describe dephased magnetization states, EPG formalism is used. A comprehensive introduction to EPGs is given in the review article by Weigel, 28 which we follow here in terms of notation and conventions; and which should be consulted for explicit mathematical expressions of the operators. EPGs describe spoiled magnetization as a superposition of complex configuration states F AE and Z. There, the central assumption is uniformity of space so that spatial dependencies in the magnetization distribution are solely arising due to gradient fields. Each configuration state corresponds to a Fourier component k

RF pulses, gradient fields, motion and relaxation are described in terms of operators acting on these configuration states. RF pulses T ϕ α ðÞdistrib- In Figure 3, EPGs for unipolar MRE in the presence and absence of periodic motion are shown. Since the spoiling area is constant for all imaging shots, the y-axis counts the spatial frequencies in units of the spoiling area Δk in the direction of the spoiling gradient. The unipolar spoiling gradient, which is played before the readout in the time-reversed spoiled SSFP sequence, leads to dephasing of configuration states with k > 0 and to respective rephasing of configurations with k < 0. Due to the presence of periodic motion, the configuration states accumulate a phase depending on their k-value. This is depicted by the bending of the configuration states, which is proportional to the relative phase accrual with respect to the previous RF pulse. This should not be confused by an actual change in spatial frequency of the respective states. 

In addition to periodic motion, diffusion attenuates configuration states with high spatial frequency. 28 This is illustrated 2-fold: first, by reduced line thickness and opacity in the EPG with increasing |k| and, second, by the diffusion attenuation plot to the right. At the bottom of the graph, the unipolar MRE sequence is denoted together with its simplified operator notation. 

In the present work, two EPG-based analyses are proposed: (i) the EPG representation of unipolar MRE can be used for full numerical encoding efficiency calculations by driving a forward calculation according to the operator sequence denoted in Figure 3 to steady-state; and (ii) a reverse analysis of the echo signal contributions can be conducted, describing the lowest order perturbation theory solution of the EPG in terms of T 2decay. This approach is commonly known as the two-transverse TR-approximation, which has been proposed before for the understanding of diffusion contrast in DW-SSFP sequences, and which is extended here to incorporate periodic motion. 

Here, k → 0 is the initial wave vector of the configuration state. Assuming periodic motion (Equation 5 ), the general motion propagator for a single constant gradient lobe with infinite slew and finite duration T is given by

and

where

T is again the area of the gradient lobe. The propagator for longitudinal magnetization as well as transverse magnetization in the absence of a gradient field can be found by setting Δk → to zero and is thus given by J 0 only. Since a unipolar lobe causes a change in configuration state population, the full unipolar propagator is given by a concatenation of motion sensitization and subsequent shifting of the configuration state population according to the gradient operator S. The motion propagator for the general case of a trapezoidal lobe with finite slew is stated in Appendix A (see the supporting information). 

In addition to motion sensitization, the dephased magnetization states are subject to diffusion attenuation. 28 The 1D diffusion attenuation operator is given by

Thus, for transverse magnetization in the case of an infinite-slew gradient lobe, the operator becomes

with

Again, D 0 is the propagator in the absence of the gradient as well as the corresponding propagator for longitudinal magnetization states. The diffusion attenuation can be incorporated in the general motion propagator for unipolar MRE by applying the replacement J→ JD and J 0 →J 0 D 0, respectively. 

Diffusion attenuation is an important factor in unipolar MRE, as it leads to attenuation of high spatial frequency states since the exponential scales proportional to k 2. This ultimately reduces the sensitivity of the sequence to T 1 and T 2 relaxation, flip angle, as well as leading to the reduction of the nonlinearity of motion encoding.
section_name: | EPG formalism for unipolar MRE

section_content: The two-transverse TR approximation is an approximate solution to the problem of diffusion weighting in SSFP-based sequences. 31, 44 It determines the signal based on the assumption that only echoes with no more than two transverse TR-periods contribute to the signal amplitude of the SSFP-echo. This corresponds to the lowest order Taylor series expansion of the analytical solution in

In Figure 4, the EPG is drawn in reverse, ie, starting from an echo and tracing back the signal contributions in lowest order perturbation theory in E 2. This leads to the graphical depiction of the three distinct echo pathways of the two-transverse TR approximation. Solid paths are transverse magnetization states, whereas dotted lines denote longitudinal magnetization. Each path is preceded by infinitely many longitudinal magnetization-only paths to account for the equilibrium magnetization. Thus, every pathway can be characterized by the number m of longitudinal storage periods. Transverse magnetization is subject to motion sensitization according to the transverse motion propagators J T m þ 1 ðÞ T R ðÞ and J T0. J T0 corresponds to the motion sensitization of the refocusing pathway before readout, while Under consideration of the periodicity of the motion propagators

the derivation of the two-transverse TR approximation can be carried out by following the derivation of McNab and Miller, 44 which is laid out explicitly in Appendix C (see the supporting information). In conclusion, the echo signal is given by

where the SE and STE contributions are found to be

with

The final step, determination of the motion phase accumulated in steady-state, is performed by Fourier transformation of the argument of the 24 ) according to

Here, the integration is normalized such that the Fourier transformation reproduces amplitude and phase of a sinusoidal oscillation according to

To the knowledge of the authors, no closed-form solution of Equation 29 exists. Hence, the final solution is obtained numerically using MATLAB. 

3 | METHODS
section_name: | Extending the two-transverse TR approximation in the presence of periodic motion

section_content: The EPG simulations were based on the Stanford EDU MATLAB implementation 45 and run in MATLAB 2018a (MathWorks, Natick, Massachusetts). The unipolar motion propagator was implemented according to Equation 16 including diffusion attenuation and using the full trapezoidal gradient waveform. The steady-state signal behavior was simulated by repeated application of the EPG formalism to drive an initial equilibrium magnetization configuration to steady-state. The propagator ordering is given in abbreviated form in Figure 3 and was evaluated using the appropriate timing constraints given by the Ristretto MRE scheme (Equation 1 ). The analytical two-transverse TR approximation model was numerically solved using MATLAB as well. 

In order to determine the encoding efficiency in the steady-state, the unwrapped signal phase was obtained for the last N P elements, which
section_name: | EPG simulations

section_content: Unipolar MRE as well as conventional GRE-MRE were performed on a gel phantom with four cylindrical inclusions (CIRS Inc., Norfolk, Virginia) using a 3 T Philips Ingenia system (Philips Healthcare, Best, the Netherlands). The signal was received using a 15-channel head coil, and 35 Hz electromagnetic actuation was employed. 4 The phantom was positioned upright in the head coil with the actuator placed on top with primary actuation being in the top-down direction. Slice selection direction was chosen anterior-posterior (see also Figure 8 ). In the following, scan parameters are stated in pairs of "conventional/unipolar MRE"; 20/10 repeat measurements were acquired to investigate the dependence of the shear stiffness estimates with increasing SNR. All measurements were performed using fractional encoding 6 with encoding fractions of 17. 5%/7%, 11. 9/ 20. 9 mT/m bipolar/unipolar MEGs (duty cycle limitation in conventional MRE) resulting in an encoding efficiency of 4. 01/10. 8 rad/mm, 3 mm 3 isotropic resolution and 64 × 48 matrix size acquired using a 30°/20°flip angle, Cartesian readout, and third/first in-phase water-fat echo-time (TE 6. 9/2. 3 ms). Conventional GRE-MRE was performed as a MS acquisition of 12 interleaved slices resulting in a TR per slice of 110. 7 ms, while unipolar MRE was performed using slab-excitation and phase encoding of 36 slices with a TR of 5. 71 ms. The increased FOV in slice direction was chosen to eliminate fold-over artifacts due to production impurities with strong susceptibility mismatch in the bottom of the phantom. The conventional GRE-MRE sequence was accelerated using the Ristretto MRE 12, 33 scheme acquiring eight wave-phase offsets in reverse order (N D =7)

and excitation and readout of the 12 slices for one wave-phase offset within N W = 3 wave periods. Only five wave-phase offsets were acquired with unipolar MRE to conform with sequence-timing constraints. The spoiling gradient was chosen to result in 5 × 2π/Δx phase accumulation along each spatial direction, separately. The acquisition of all wave-phase offsets and encoding directions for a single slice took 14. 3/7. 0 seconds, respectively. Synchronization to the external wave generator was achieved using a TTL trigger signal. In order to achieve a magnetization and wave steady-state, 62 dummy shots were used for unipolar MRE, totaling to 1. 8 seconds per encoding direction (conventional MRE: 10 dummies, 1. 1 seconds/encoding direction). 

A second set of phantom scans was performed at 70 Hz actuation frequency. Unipolar MRE with five wave-phase offsets was acquired with
section_name: | Phantom comparison

section_content: In Figure 7, the results of the quantitative parameter mapping for (A) T 1 and T 2, and (B) ADC and FA are displayed. T 1 is found to be approximately constant over the whole phantom (2425 ± 33 ms), while T 2 differs significantly between background (193 ± 25 ms) and inclusions (632 ± 33 ms). The ADC is found to be approximately constant over the whole phantom ([1. 68 ± 0. 31] x 10 -3 mm 2 /s). FA is close to zero, consistent with the assumption of isotropic diffusion. ROI average and standard deviation for the inclusions, background and the whole phantom are given in Table S2. In Figure 7C, the magnitude and phase of the local encoding efficiency are shown based on the determined T 1 and T 2 map; however, diffusion effects are neglected. Substantial encoding efficiency differences between background (|ε| = 8. 77 ± 0. 30 rad/mm, ∠ε = 0. 192 ± 0. 008 rad) and inclusions (|ε| = 5. 76 ± 0. 13 rad/mm, ∠ε = 0. 276 ± 0. 009 rad) of ~52% can be observed. Figure 7D shows the encoding efficiency taking a constant diffusion coefficient of 1. 7 x 10 -3 mm 2 /s into account. Here, the encoding efficiency of the background S2. 

In Figure 8 S3. For 35 Hz, the OSS-SNR of unipolar MRE without averaging is higher than the OSS-SNR of the 20-fold averaged conventional scan for all regions. At 70 Hz, the single average unipolar MRE acquisition is comparable with an 8-to 10-fold averaged conventional MRE scan. 

Stiffness of inclusions 1 and 2 as well as the background are well recovered independent of the number of averages and technique. Conventional MRE qualitatively fails to predict the stiffness of inclusions 3 and 4, showing higher stiffness for the softer, third inclusion. Through averaging, the stiffness of inclusions 3 and 4 determined by conventional MRE increases significantly by 54% (inclusion 3) and 179% (inclusion 4) at 35 Hz and by 3. 8% and 10. 3% at 70 Hz, respectively. In conventional MRE and at 35 Hz, the stiffness of inclusion 3 reaches a plateau after 11 averages, while the stiffness of inclusion 4 increases even after 20 averages. With unipolar MRE, the change in stiffness is below 0. 3% for the first three inclusions as well as the background. Change in stiffness for inclusion 4 is 6. 2% at 35 Hz and 1. 9% at 70 Hz, and reaches a plateau after five averages in both cases.
section_name: | Phantom results

section_content: In this work, we have proposed a novel MR elastography sequence based on time-reversed spoiled SSFP, which utilizes the spoiling gradients as highly efficient unipolar MEGs. We have motivated the fundamental encoding principle of unipolar MRE through analogies to phase-contrast MRI and have described an extension of the two-transverse TR approximation developed for DW-SSFP to account for lowest order spin-history effects. 

The theoretical predictions of the proposed approximations have been compared with full EPG simulations in a comprehensive parameter study. It was found that longer T 1 and T 2 values as well as lower diffusion generally leads to a reduction in encoding efficiency. We attribute this finding to the destructive interference of an increasing number of refocused echo pathways. Contrary to conventional MRE, unipolar MRE encodes displacements d nonlinearly into the signal phase, as echo pathways with different motion sensitivities are mixed. However, we have identified a linear encoding domain, which is approximately given by Δk|d|≲10%. Since the contributions of the different echo pathways are dependent on relaxation and sequence parameters, it is safe to assume that the displacement range that is mapped linearly to the signal phase will be dependent on these parameters as well. This can also be appreciated in Figure S2, which is the short T 1 /T 2 equivalent of Figure 6 (T 1 = 300 ms, T 2 = 30 ms). There, the linear encoding domain is approximately twice as large and given by Δk|d|≲20%. Since encoding efficiency was found to be monotonically decreasing, unipolar MRE shows heightened sensitivity to low displacement amplitudes. In the case of conventional GRE-MRE, the encoding efficiency decreases as the duration of the encoding gradient is reduced, even if a finite moment for each of the two lobes is assumed. By contrast, unipolar MRE shows a finite encoding efficiency in the case of the infinitesimal model given by Equation 12. This demonstrates the uniqueness of unipolar MRE, which, unlike conventional GRE-MRE, splits encoding over two (and more) TR periods, thus effectively storing spin positions in the phase of the MR signal during the first TR and subtracting the updated position with the second TR. Analogous to SE sequences with single gradient lobes 47 or displacement encoding with stimulated echoes (DENSE), 55 it is precisely the mixing time between the two encoding steps that leads to the high motion sensitivity of the unipolar MRE sequence, and timereversed SSFP sequences in general. This relationship also manifests itself in the comparison of encoding efficiencies as a function of encoding fractions q in Figure 6A, where unipolar MRE shows a linear increase in encoding efficiency as a function of the MEG duration, whereas conventional bipolar MEGs follow a quadratic behavior in the small q domain. 

Note that the quadratic behavior is reminiscent of the encoding of displacement velocities, while the linear relationship arises from the direct encoding of displacement. 3 Since the gradient strength G is kept constant in the plot, the moment Δk vanishes for unipolar MRE with q → 0 and thus also its encoding efficiency. 

A direct consequence of encoding displacements in unipolar MRE (as opposed to encoding displacement velocities with conventional MRE) is the dependency of the encoding efficiency on sequence timing rather than encoding fraction (assuming Δk is kept constant). The SSFP-SE contribution to unipolar MRE is shown to decrease with increasing N P =N D ratio, as displacements between successive encoding/decoding steps become more and more similar. The full EPG simulations, on the other hand, show increasing encoding efficiency with increasing N P =N D ratio. 

We speculate that the interference of echo pathways is reduced as successive TR periods accumulate more and more similar motion phases. 

The parameter study also revealed an increase in encoding efficiency with increasing flip angle. Since an increase in flip angle leads to a reduction of the population of high-spatial frequency states, this finding is in general accordance with the dependency on T 1,T 2 and diffusion. Furthermore, it supports the hypothesis of destructive interference of motion encoding due to differences in motion sensitivity of echo pathways. Contrary to bSSFP-based MRE, 6, 21, 22 the encoding efficiency of unipolar MRE is independent of off-resonance, as dephasing and refocusing configuration states experience opposite phase accrual. Thus, each echo pathway shows the same dependency on off-resonance Δω, which is given by The x-axis denotes the average octahedral shear strain-SNR (OSS-SNR) for each ROI. Sample averaging was performed to increase the SNR in both unipolar and conventional MRE to demonstrate the influence of SNR on the reconstructed shear stiffness. Here, the left-most data-points in each series correspond to the single-scan case, where each successive point increases the number of sample averages up to N Avg = 20. While the stiffness of soft inclusions (numbers 1 and 2) and the background are found to be independent of the OSS-SNR, stiff inclusions (numbers 3 and 4) depend on the number of sample averages for conventional MRE at 35 Hz. Due to the increased OSS-SNR of unipolar MRE, even without averaging, unipolar MRE can directly capture the stiffness of these inclusions exp iΔω T E -T R ðÞ fg. Thus, in equivalence to conventional spoiled GRE-MRE, off-resonance leads to an image phase offset, which is independent of the wave-phase offset and is removed once the Fourier transformation over wave phases is performed. 

Since unipolar MRE uses the spoiling gradients as MEGs, the sequence TR can be considerably reduced compared with conventional GRE-MRE. 

Given that time-reversed spoiled SSFP sequences cannot be used with interleaved MS acquisitions, the sequence TR is generally much shorter than typical T 2 values. A Taylor series expansion of the echo signal in exp -T R =T 2 fg is at first sight not feasible. Nonetheless, the two-transverse TR approximation proves to be valuable in describing the encoding efficiency of the unipolar MRE acquisition, especially within the linear encoding regime. This can be appreciated from the closeness of predicted encoding efficiency of the two-transverse TR approximation and the full EPG simulations in Figures 5 and 6, as well as the close coincidence of the displacement fields in the phantom experiments depicted in Figure 8. While differences between the full EPG and the two-transverse TR approximation certainly persist in theory, we believe that the attenuation of configuration states is not only due to diffusion and transverse relaxation but is complemented by intra-voxel phase dispersion (see below) additionally driving the sequence closer to the two-transverse TR solution. 

The propagator definition in Equation 16 assumes periodic bulk motion. From conventional MRE it is known that the variation of the displacement field within one voxel leads to intra-voxel phase dispersion (IVPD) and thus signal attenuation. 56, 57 In the 1D case of a propagating planewave with wavelength λ, signal attenuation by IVPD in the conventional MRE case can be estimated to be

where S 0 denotes the expected signal without IVPD and S the attenuated signal. 8, 57 In the case of unipolar MRE, intra-voxel variations of the displacement field can be treated using the spatially resolved EPG, which assumes each spatial position as independent EPG simulations. 58 The echo signal can then be approximately expressed as the average of sub-voxel signals, each being a superposition of refocused echo pathways. 

Since summation over echo pathways and spatial averaging can be exchanged, IVPD in unipolar MRE can not only be discussed as a total effect of signal attenuation, it can also be studied for each echo pathway separately. As the motion sensitivity of configuration states increases proportionally to their spatial frequency, each echo pathway accumulates a different motion-dependent phase. From Equation 31 it becomes clear that echo pathways involving higher order configuration states will exhibit stronger attenuation due to their increased phase dispersion. Thus, similarto diffusion, the presence of IVPD is expected to reduce the contributions of higher order configuration states to the steady-state signal and therefore reduces the influence of tissue, experiment and sequence properties on the encoding efficiency. 

In the derivation and analysis of unipolar MRE, the assumption of isotropic diffusion was made. While the DTI measurements support this assumption in the phantom used, a discussion of anisotropic diffusion is necessary when unipolar MRE is employed in vivo. McNab and Miller extended the two-transverse TR approximation for anisotropic DW-SSFP by replacing the exponent of the diffusion propagator by its respective anisotropic form. 44 Weigel et al later introduced anisotropic diffusion in the context of the EPG formalism. 59 Both approaches can be equally used here. Since the encoding direction and thus diffusion sensitization is held constant in successive TRs, the diffusion coefficient D that enters into the 1D diffusion propagator (Equation 20 ) can be simply replaced by the projection of the diffusion tensor D on the encoding direction e

Thus, in the presence of anisotropic diffusion, encoding direction-dependent attenuation of configuration states, and with these direction-dependent encoding efficiencies, are expected to occur. 

In the nonlinear regime (here Δkd jj>10%), full EPG simulations show a monotonous decline in encoding efficiency, which leads to systematic underestimation of large displacements (Figure 6C ). Since the dependency on the displacement amplitude is found to be monotonous within the relevant displacement domain, a bidirectional relationship is likely to exist. With exact knowledge of the B þ 1 field as well as tissue parameters within the ROI, full EPG simulations with different displacement amplitudes can be performed to correct the displacement fields, eg, via a lookup-table approach. This also holds for the variation in phase angle (as observed in Figure S4 ). 

The sequence has been validated in a phantom study against a conventional GRE-MRE sequence at two driving frequencies, where displacement fields were found to be in excellent agreement. We attribute the remaining differences to differences in the slice profile of MS and 3D acquisitions as well as to residual differences in encoding efficiency due to local relaxation parameter variations. Especially at low driving frequency, the displacement noise is visually higher in the conventional MRE scan than in unipolar MRE, which is also confirmed by the OSS-SNR calculations. 

The low SNR of the conventional MRE scan leads to underestimation of the stiffness of the two stiffest inclusions. Here, conventional MRE without averaging even fails to qualitatively predict the elevated stiffness of inclusion 4 compared with inclusion 3 (Figure 9 ). The stiffness of inclusions 1 and 2 as well as the phantom background were found to be in excellent agreement between unipolar and conventional MRE. 

Through averaging, the SNR of the conventional MRE scan could be successfully increased to the point of close agreement of recovered stiffness values with unipolar MRE. At 35 Hz, however, even 20-fold averaging of the conventional MRE scan did not provide sufficient SNR to allow estimation of the stiffness of the two stiffest inclusions. Minor differences in the order of 0. 2 kPa persist, however, for inclusion 3 at 35 Hz, and also the two stiffest inclusions at 70 Hz. We believe these differences to be again attributable to slice profile differences, since the high stiffnessof the inclusions leads to very low spatial derivatives and thus renders the inversion especially sensitive to geometrical mismatches. 

In this study, conventional MRE was performed with fewer slices compared with unipolar MRE (12 vs. 36 slices) to keep scan durations within acceptable ranges. Since unipolar MRE uses slab excitation of the whole phantom and phase encoding in the slice-encoding direction, SNR of 3D unipolar MRE is ffiffiffiffiffiffi N S p -fold higher than a respective 2D scan with equal TR, where N S denotes the number of phase-encoding steps in the slice direction. Accordingly, the SNR of the unipolar MRE scan was √3 higher, equaling 3-fold averaging. An OSS-SNR efficiency adjusted for the slice acquisition technique can be defined, which relates the OSS-SNR to the total acquisition duration T Tot and accounts for the number of averages N Avg and slices N S. For conventional MRE, where MS acquisition is performed, the OSS-SNR efficiency can be defined as

The use of slab excitation (3D) in unipolar MRE leads to In the phantom study, the actuation strength was specifically chosen to remain within the linear encoding regime. While in conventional MRE the encoding efficiency can be reduced by gradient scaling or lowering the encoding fraction, 6 in unipolar MRE spoiling and motion sensitivity are coupled. Thus, a reduction of encoding efficiency to remain within the linear encoding regime is not feasible. Changes to the duration of the encoding gradient at constant zeroth moment only leads to minor reductions in encoding efficiency. Thus, if linear encoding is intended, displacement amplitudes must be chosen accordingly. 

The motion propagator in Equation 18 cannot only be used for the simulation of unipolar MRE, but for sequences with periodic motion in general. Since the propagator describes motion sensitivity of configuration states subject to a generalized trapezoidal gradient, concatenation of the propagator allows performing EPG simulations with bipolar or flow-compensated gradients. Hence, it can be employed to study spin-history effects in conventional 3D GRE-MRE with short TRs or the influence of RF spoiling and its efficiency in canceling higher order echo pathways. 

Previous work by Scheffler et al has described an encoding concept similar to unipolar MRE to induce an "oscillating steady-state" in a spoiled SSFP sequence with readout of the F -2,F -1,F 0 and F 1 configuration states. 60 and thus makes assumptions on the contributing echo pathways. On the other hand, Scheffler et al developed a signal model for small perturbations of the accumulated phase during a TR period without making assumptions on relaxation properties. Transferring their approach to unipolar MRE, it becomes a model to describe motion encoding for small displacements in the case of N P = 2. An extension of their work to more than two states is beyond the scope of this work; however, it would provide a valuable addition to the understanding of encoding in unipolar MRE in the absence of diffusion or IVPD. 

In the present work, Hadamard encoding has been used to allow for the extraction of the 3D volumetric displacement field. It should be noted that unipolar MRE is not restricted to Hadamard encoding alone-all encoding schemes lacking a separate reference scan without motion sensitization can be used. Clinically, MR elastography is often performed on a small set of discontinuous 2D slices with only through-plane motion sensitization. 63, 64 In unipolar MRE, slab excitation and phase encoding in the slice direction can be replaced by single slice excitation. Interleaved MS acquisitions cannot be paired with unipolar MRE as it is described here, as spoiling gradients in consecutive TR periods will accumulate and change the overall spoiling area as well as the motion sensitization scheme. 

While the present work has focused on the theoretical foundations and phantom validation of unipolar MRE, future work is required to assess the utility of the method in vivo. The high intended motion sensitivity of unipolar MRE requires careful consideration and provisions depending on the target anatomy, in particular if physiological and patient motion is present. Accordingly, unipolar MRE is expected to require breath-holding in the case of thoracic and abdominal applications. To address cardiac-related motion, artifact reduction techniques such as averaging or gating may need to be incorporated. 65 In brain and spine, the pulsatile CSF flow is expected to be encoded by unipolar MRE, potentially interfering with displacement encoding. 36, 42, 66 However, since SNR is much higher in unipolar MRE than in a comparable conventional GRE-MRE sequence, a combination of compressed sensing and self-gating techniques is envisaged. 13, 65, 67, 68 Overall, the influence of physiological and patient motion on the steady-state signal and resulting changes in encoding efficiency remains to be understood. In addition, in vivo parametric mapping of T 1,T 2 and ADC would be required for the tissue of interest to estimate their impact on the local encoding efficiency. While the phantom results and theoretical considerations suggest that relaxation effects are not critical, these results need to be confirmed in vivo. 

Given the considerations above, a first application of unipolar MRE could focus on brain, where mechanical attenuation by the skull results in low displacement amplitudes necessitating high encoding efficiency and SNR. 

In the future, unipolar MRE might prove valuable for multi-frequency MRE, when different actuation frequencies are acquired consecutively. This is especially the case for low driver frequencies as the SNR is significantly increased with unipolar MRE compared with conventional MRE with fractional encoding. At high frequencies, however, unipolar MRE loses its advantage over conventional MRE as encoding fractions typically become much larger and acceleration is limited by the combined duration of RF pulse and readout, not by the MEG duration.
section_name: | DISCUSSION

section_content: We have proposed a novel MR elastography sequence termed "unipolar MRE", which utilizes the spoiling gradients of a time-reversed SSFP sequence to acquire the full 3D shear displacement field. Compared with conventional GRE-MRE, unipolar MRE offers significantly reduced TRs with comparable or even increased displacement encoding efficiency. Accordingly, 4-to 8-fold faster acquisitions at equal or higher OSS-SNR can be achieved. 

The steady-state echo is given by a superposition of these terms, yielding

Insertion of the STE definition into the infinite sum leads to

An auxiliary term A can be identified by collecting all terms dependent on m. Shifting the summation index m→m-2, we arrive at

Splitting the infinite sum into two

allows to make use of the periodicity of the motion propagators implied by the Ristretto sequence timing, i. e. 

J T m þ N P ðÞ T R ðÞ ¼ J T mT R ðÞ :

The summation over m can again be performed, leading to

A ¼ E 1 cos α ðÞ ðÞ 2 1-E 1 cos α ðÞ ðÞ NP κ N P ðÞ :

Hence, the STE contributions with m > 1 are given by
section_name: | CONCLUSION

section_name: ACKNOWLEDGEMENTS

section_content: The authors would like to thank Ayse Sila Dokumaci at King's College London and Jurgen Henk Runge at King's College London and University Medical Center Amsterdam for their support with initial experiments.
section_name: ACKNOWLEDGEMENTS

section_name: | RESULTS

section_content: In addition to the acquisition of the MRE datasets, quantitative parameter mapping was performed to obtain local estimates for T 1,T 2, the apparent diffusion coefficient (ADC) and the fractional anisotropy (FA) of the phantom in one representative slice in order to calculate the local, complex encoding efficiency for the unipolar MRE acquisition. T 1 and T 2 estimates were obtained on the same 3 T Philips Ingenia system using a spoiled GRE-based MR fingerprinting (FISP-MRF) 46 sequence with a constant TR of 16 ms, fully sampled multi-shot spiral acquisition with 15 interleaves, a repetition delay of six seconds to ensure relaxation to equilibrium, and 1000 variable excitations after adiabatic inversion resulting in a total scan duration of 5. 5 minutes. The flip angle train is given in Figure S1. Dictionary-based matching was performed to determine T 1 and T 2 locally using MATLAB 2018a (MathWorks). The dictionary was created based on the Stanford EDU MATLAB EPG implementation. 45 The dictionary resolution is given in Table S1. SE-based diffusion tensor imaging (DTI) with unipolar encoding gradients, six encoding directions and five b-values (linearly interpolated, 0 to 1000 ms/mm 2 ) was performed to obtain an estimate of the ADC and FA in a single slice. 47, 48 The scan was performed on a 1. 5 T Philips Achieva system (Philips Healthcare). The signal was received using an eight-channel head coil. Sequence parameters were 3 x 3 mm 2 in-plane resolution, 8 mm slice thickness, 10 averages, single-shot EPI readout, matrix size 64 x 45, FOV 196 x 136, flip angle 90°/180°, 60 ms echo time, and a TR of 1000 ms, resulting in 4 minutes 11 seconds scan duration. 

For comparability, the FOV of the unipolar MRE scan was cropped to the 12 central slices to match the conventional MRE scan prior to postprocessing. The image phase was spatially unwrapped using the statistical cost, network-flow algorithm for phase unwrapping (SNAPHU), [49] [50] [51] [52] temporally aligned, and decoded using the pseudo-inverse of the encoding matrix. 8 Cumulative averaging of the decoded wave phase images was performed prior to Fourier-transformation. The complex displacement fields were estimated from the first Fourier component and corrected for the encoding efficiency of the conventional and unipolar MRE scans. The phase of the complex displacement field of the conventional MRE scans was corrected for inter-slice phase shifts according to the Ristretto MRE scheme. 12, 33 The encoding efficiency of the unipolar MRE acquisitions was determined using the two-transverse TR approximation assuming homogenous relaxation and diffusion (T 1 = 2. 4 s, T 2 = 190 ms and D = 1. 7 x 10 -3 mm 2 /s). The complex displacement fields were then further processed without smoothing using a FEM-based reconstruction algorithm 53 to obtain the complex shear stiffness, where we report the averaged magnitude over the innermost six slices. In addition, the wave-phase images were used to estimate the displacement noise in each voxel, allowing for the determination of the octahedral shear strain-SNR (OSS-SNR). 54 OSS-SNR values are reported for the inclusions and for the background as ROI-averaged values.

section_content: In Figure 5, the encoding efficiency of unipolar MRE is shown as a function of tissue parameters T 1,T 2, and the diffusion coefficient D within the phase-contrast approximation (red, SSFP-SE contribution), the two-transverse TR approximation (yellow), and the full EPG simulation (violet). If applicable, full EPG simulations neglecting diffusion attenuation are shown in dotted violet. Sequence parameters were held constant at the values displayed by the vertical dashed lines. In Figure 5A, the encoding efficiency is shown as a function of the diffusion coefficient. The phase-contrast model is independent of diffusion by construction. For typical diffusion coefficients in vivo (D = ~2 x 10 -3 mm 2 /s), the encoding efficiency obtained from full EPG simulations is within 3% of the two-transverse TR approximation and converges for strong diffusion. Neglecting diffusion, the encoding efficiency obtained from full EPG simulations deviates by ~20%. 

In Figure 5B-D, the encoding efficiency is shown as a function of T 1 and T 2 under the assumption of finite diffusion (D = 1. 7 x 10 -3 mm 2 /s). 

Figure 5B shows the case of long T 2 relaxation (190 ms), where full EPG and the two-transverse TR approximation show little dependency on T 1 with variations below 3%. Figure 5D shows the case of short T 2 relaxation (30 ms). For long T 1 (>1000 ms), full EPG with and without diffusion as well as the two-transverse TR approximation become equivalent. At very short T 1 (<100 ms), maximal deviations of up to 5% are found. Contrary to the long T 2 case, diffusion only plays a minor role, as shown by the near equivalence of full EPG with and without diffusion. This can also be appreciated in Figure 5C, which shows the T 2 dependency for the case of long T 1 (2400 ms). Below 50 ms, the two-transverse TR approximation and the full EPG with and without diffusion are closely aligned (<1% deviation). For increasing T 2, full EPG simulations, especially without diffusion attenuation, show decreasing encoding efficiency compared with the two-transverse TR approximation. Again, the encoding efficiency of both phase-contrast and the two-transverse TR approximation are independent of T 2 by design. The dependency of the encoding efficiency phase angle on the same set of parameters is shown in Figure S3. 

In Figure 6, the encoding efficiency of unipolar MRE is compared for different sequence parameters. Figure 6A shows the encoding efficiency as a function of encoding fraction q and compares the different unipolar MRE approximations against conventional phase-contrast MRE with bipolar MEGs (blue). The red line indicates the encoding fraction limit imposed by the sequence timing (here N P =5,N D = 1), which limits the duration of the encoding gradient. In the small q regime (<10%), the encoding efficiency of unipolar MRE is proportional to q, whereas conventional bipolar MRE is proportional to q 2. In Figure 6B, the dependency on the sequence-timing parameters N P and N D is shown. The SSFP-SE contribution shows its maximal encoding efficiency at N P = 2, where the phase contributions from two consecutive TRs interfere
section_name: | EPG simulations

section_content: This project received funding from the European Union's Horizon 2020 research and innovation programme under grant agreement no. 668039.
section_name: FUNDING INFORMATION

section_content: Additional supporting information may be found online in the Supporting Information section at the end of the article. 

How to cite this article: Guenthner C, Sethi S, Troelstra M, et al. Unipolar MR elastography: Theory, numerical analysis and implementation. NMR in Biomedicine. 20 ;33:e4138. https://doi. org/10. 1002/nbm. 4138
section_name: SUPPORTING INFORMATION

section_content: The motion propagator for a trapezoidal gradient with finite slope time t d, duration T, and gradient area Δk is given by

The propagator in the absence of the gradient is equivalent to J 0 in Eq. 19.
section_name: General MRE Propagator for a Trapezoidal Gradient Lobe

section_content: To evaluate the two-transverse TR approximation, the motion propagators for storage and transverse magnetization periods need to be defined. 

The propagator directly preceeding readout is given by

Here, the notation •j rule should be interpreted as a replacement rule, which is applied to the expression on its left. 

The motion propagator directly following initial RF excitation is given by

The first term denotes initial sensitization through the unipolar MEG, while the second term accounts for phase accrual due to magnetization being spoiled. 

Configuration states in longitudinal storage experience motion sensitization over the duration of the storage period mT R. The motion propagator is thus given by

Derivation of the Two-Transverse TR Approximation in the Presence of Bulk Periodic Motion

For the derivation of the two-transverse TR approximation in the presence of periodic bulk motion, we follow McNab and Miller, Ref. 44 except for the respective additions for unipolar MRE. Each echo pathway is preceded by infinitely many longitudinal storage states, which have to be taken into account to capture T 1 recovery. These contributions are given by the factor

where M 0 is the equilibrium magnetization. The solution of the infinite series is given by the geometric series identity 69

Since k ¼ 0 holds for these states, they are not subject to motion sensitization. Taking the remaining transverse periods and the RF pulses into account, the spin-echo contribution can be stated as

The stimulated echo (STE) contributions are given by two distinct terms. The STE contributions are given by (STE 1 : one storage period)
section_name: APPENDIX B Definition of Motion Propagators for the Two-Transverse TR Approximation