B1 Mapping

1Double Angle techniques¶

Transmit and Receive RF field amplitudes¶

In an MRI experiment, magnetic field amplitude of the radiofrequency field (B₁) is an an important physical parameter that is a product of the interaction between RF coil design and the subject’s position (volume and relative position to the coil) and physical properties (electromagnetic permittivity ε and permeability μ). In any MRI experiments, two B₁ fields appear: the transmit RF field amplitude B₁⁺ and the receive RF field amplitude B₁^-. The latter, B₁^-, is often referred to in terms of the receive RF coil sensitivities and this is due to the principle of reciprocity Hoult & Richards, 1976Hoult, 2000 property of RF antennas. In terms of an MRI image, B₁^- is simply a multiplication factor that varies spatially but doesn’t change between pulse sequences if the subject remains motionless, and techniques to “flatten” images by estimating this field numerically have been developed Sled et al., 1998. Additionally, many quantitative MRI techniques calculate the ratio of images, which eliminates this B₁^- component from the resulting image.

B₁⁺, however, impacts the resulting MRI images in a much more complex way and is not a simple multiplication factor. B₁⁺ directly perturbs the system of spins by introducing energy in the system, which practically we quantify as the flip angle of an RF pulse. Two different B₁⁺ values will not have the same impact on voxel for different pulse sequences, as spin dynamics and steady-state conditions will vary. For example, let’s say you acquire a saturation recovery image and also a short TR steady-state gradient echo. For an actual flip angle = $\alpha \cdot B_{1}$, the magnetization after TE will be $M_{z}\text{sin}\left( \alpha \cdot B_{1} \right)\text{exp}\left(-\text{TE}/T_{2} \right)$ and $M_z$ prior to the RF pulse is given by Eq. 2.5. Thus, not only does B₁⁺ not appear as a simple multiplication factor, a change of B₁ will not impact this voxel for both sequences by the same ratio. Thus, knowledge of B₁⁺ through B₁ mapping can help us retroactively understand the dynamics of the spins accurately, and can play the role of a calibration factor for many quantitative MRI techniques (e.g. VFA, T₂/T₂^*, qMT, etc). In addition, B₁ mapping can also map the electromagnetic properties, but this won’t be discussed in this chapter.

In this chapter, we’ll be discussing a simple but widely used class of methods for B₁ mapping, the double angle method. For the sake of simplicity, and for consistency with the quantitative MRI literature, we’ll define B₁ = B₁⁺ and will explicitly state B₁^- when referring to the receive field.

Double Angle method(s)¶

The Double Angle method (DA or DAM) is a class of B₁ mapping techniques where two RF pulses at flip angles α and 2α are applied to a pulse sequence, and the ratio of the images are compared to expected output to produce a B₁ map. Several main pulse sequences (Figure 4.1) have been called the double angle method in the literature, and both have their own equations describing the relationship between the expected images and B₁. In this chapter, we’ll mostly explore there the α-180 method (Figure 4.1]), and then briefly explain the other method and its similarities/differences.

Gradient echo and 180 degree spin echo methods¶

This pulse sequence uses a 180 degree spin-echo refocusing pulse and acquires two images using an excitation pulse α and $2\alpha$. It assumes that there is full signal recovery (long TR), and because it refocuses T₂^*, it eliminates signal variability caused by B₀ in the resulting B₁ map Insko & Bolinger, 1993. Alternatively, a gradient echo could be used?

Assuming a refocusing pulse is used (i.e. isn’t dependent on B₁), we can develop the equation for a gradient echo and spin echo case.

Thus

and

Using a well known trigonometry identity (see Appendix A for derivation),

We can simplify Eq. 4.5,

And the true flip angle can be calculated from the ratio of these two magnetizations / signals / images:

Knowing that $\alpha = B_{1}$ alpha_nominal, B₁ is thus:

This equation is also used for α-180 spin echo pulses, however it assumes no dependency on of the refocusing pulse on B₁. Figure 4.3 explores this using Bloch simulations

2Actual Flip Angle Imaging¶

Transmit radiofrequency field maps (B₁⁺, or B₁ for short) are used in diverse applications in MRI including: the study of electrical properties in tissues in vivo Sled & Pike, 1998Katscher et al., 2009, specific absorption rate (SAR) calculations Ibrahim et al., 2001, the calibration of quantitative T₁ Deoni, 2007Boudreau et al., 2017 and T₂ Sled & Pike, 2000 maps, better parameter estimation from magnetization transfer measurements Ropele et al., 2005Boudreau et al., 2018, B₁ shimming to improve image quality at whole-body ultra high fields Bergen et al., 2007, or quality control of RF coils Yarnykh, 2007. Several B₁ mapping techniques have been developed, and they can be broadly divided as magnitude-based and phase-based methods. The double angle method (DAM) is a saturation-recovery magnitude-based method that takes the ratio of the signal intensity of two magnitude images measured with different excitation flip angles Insko & Bolinger, 1993Stollberger & Wach, 1996. The Bloch-Siegert shift technique is a rapid phase-based method that encodes the B₁ information into phase signal Sacolick et al., 2010. The actual flip-angle imaging (AFI) is a magnitude-based B₁ mapping method that consists of a 3D acquisition that benefits from good anatomical coverage. In addition, this technique allows the acquisitions of whole-body (~7 min) and brain (~3 min) B₁ maps leading to a feasible implementation in clinics Yarnykh, 2007. On the other hand, the AFI pulse sequence has certain constraints that need to be considered for this B₁ mapping method to be widely deployed. Some of the limitations include the use of spoiler gradients that can give rise to prohibitive SAR values Sacolick et al., 2010, and the pulse sequence modifications on the MRI machine to implement the AFI method.

In this section, we will focus on presenting details about the AFI B₁ mapping method. We will cover signal modeling, data fitting, the benefits and the pitfalls of the technique. The figures are generated using the qMRLab module for this method.

Signal Modelling¶

The pulse sequence of the AFI method (Figure 4.5) is composed of two identical RF pulses and two different delays (TR₁ < TR₂). After each RF pulse, the signal intensity is acquired followed by a spoiler to destroy the residual transverse magnetization next to the following RF pulse. This method implements a pulsed steady-state signal with a gradient-echo acquisition, thus preventing the use of long repetition times Yarnykh, 2007. It has been demonstrated that if the delays TR₁ and TR₂ are sufficiently short (e.g. TR₁/TR₂ = 20 ms/100 ms), and the transverse magnetization is completely spoiled, the ratio of signal intensities (r = S₂/S₁) depends on the flip angle of applied pulses and is highly insensitive to T₁ Yarnykh, 2007.

The magnetization of an AFI experiment can be modeled under steady-state conditions by the implementation of a fast repetition of the sequence (TR₁ < TR₂ < T₁). The solution of the Bloch equations for the AFI method is given by Equations 1 and 2 that represent the longitudinal magnetization before the application of the RF pulses:

Mz_1,2 is the longitudinal magnetization of both pulses, M0 is the magnetization at thermal equilibrium, TR₁ is the delay time after the first pulse, TR₂ is the delay time after the second identical pulse (Figure 4.6), and θ is the excitation flip angle. The steady-state longitudinal magnetization Mz curves for different T₁ values for a range of $\theta_{n}$ and TR values are shown in Figure 4.6.

The analytical solution of the Bloch equations in a steady-state experiment (Eq. 4.10 and Eq. 4.11) makes several assumptions leading to practical challenges. First, it is assumed that the longitudinal magnetization has reached a steady state after a sufficiently large number of repetition times (TR), and that the transverse magnetization is perfectly spoiled prior to each pulse. To explore these properties, a numerical approach known as Bloch simulations is used to estimate the signal from an MRI experiment given a set of sequence parameters. Here, the Bloch simulations allow us to estimate the magnetization using a different number of sequence repetitions, and look at a special case when the steady-state is not achieved (due to a small number of sequence repetitions). As can be seen in Figure 4.7, the number of repetitions required to reach a steady-state depends on T₁ and the flip angle.

In practice, gradient and RF spoiling are important parameters to consider in an AFI experiment. A combination of both Zur et al., 1991Bernstein et al., 2004 is typically recommended, and Figure 4.8 shows how this better approximates the ideal spoiling case.

Data Fitting¶

The ratio of Eq. 4.10 and Eq. 4.11, gives rise to Eq. 4.12 that depends on the parameters T₁, TR₁, TR₂ and the excitation flip angle (θ).

Eq. 4.12 can be simplified if the Taylor series expansion of the exponential function is used, followed by a first-order approximation to its terms. For this expansion, short TR₁ and TR₂ (TR₁ < T₁ and TR₂ < T₁) are assumed to approximate the signal intensities ratio (Eq. 4.13) where n = TR₂/TR₁.

Finally, a measure of the actual flip-angle (θ) can be achieved by solving Eq. 4.13 to obtain Eq. 4.14, which only depends on the signal intensities ratio (r = S₂/S₁) and the parameters TR₁ and TR₂.

The actual flip-angle is estimated using an approximation (Eq. 4.13) of a complete analytical solution (Eq. 4.12), and the nature of this approximation makes it worthwhile to assess the accuracy of the signal intensities ratio between both equations. Next, a set of simulations are displayed to analyze how the choice of r is affected by T₁, TR₁ and TR₂. First, the effect of the relaxation time T₁ is simulated in Figure 4.9 for both the approximation and the complete analytical solution.

The signal ratio r is highly insensitive to the relaxation time T₁, except for the low T₁ values and large flip angles (>70°). This shows that the Taylor expansion is a good approximation to the signal ratio r because it is possible to get rid of the inverse quadratic T₁ dependance by taking the first-order terms of the expansion.

The effect of the TR₁ parameter on the signal ratio is shown in Figure 4.10. To assess the influence of the repetition time, we fix n=5 and vary the parameter TR₁ in accordance to the relation n = TR₂/TR₁. As TR₁ increases (> 50 ms), the approximated ratio r slightly deviates from the analytical approach. Although the deviation is slight only at high flip angles, a good signal ratio approximation can be achieved for a wide range of flip angles and repetition times.

Finally, the effect of the parameter n on the signal ratio r (Figure 4.11) does not seem to significantly affect the signal ratio between the approximated equation and the analytical approach. However, the parameter n has a major impact on the sensitivity of the AFI method to variations in the flip angle. Figure 4.11 shows that the increase of the parameter n (= TR₂/TR₁) allows for improvement of the dynamic range of flip angles measurements. These simulations have shown that an optimal implementation of the AFI method requires a careful selection of sequence parameters.

Figure 4.12 displays an example AFI dataset and its corresponding field B₁ map in a healthy human brain. Although not clearly visible, both AFI images present a small Gibbs ringing artifact that is propagated and amplified due to the AFI calculation consisting of the division of both images Boudreau et al., 2017. The ringing artifact is clearly seen in the unfiltered/raw B₁ field map shown in Figure 4.12 (right).

The ringing artifact shown in Figure 4.12 can be attenuated by implementing a smoothing process. Figure 4.13 shows the raw (left) and the filtered (right) B₁ map where a median filter was used to smooth the field map.

Benefits and Pitfalls¶

B₁ mapping is of interest for diverse MRI applications, and several mapping techniques have been developed. The DAM method consists of acquiring two scans at two different flip angles. To avoid the dependence of the signal on T₁, long repetition times are required to allow the recovery of the longitudinal magnetization between pulses Yarnykh, 2007Insko & Bolinger, 1993. The AFI method overcomes this practical limitation by repeating the pulse sequence at a fast rate to achieve a pulsed state of magnetization and shorter time delays between pulses. In addition, due to scan-time constraints, B₁ mapping methods are often implemented in 2D Chavez & Stanisz, 2012. However, the accuracy of the measurements of 2D B₁ mapping techniques is compromised by the slice profile effects due to the problem of nonuniform excitation across slices Yarnykh, 2007Chavez & Stanisz, 2012. The AFI method on the other hand, adresses this issue using a fast 3D implementation leading to scans with an excellent anatomical coverage in clinically feasible times, with an increase in signal-to-noise ratio compared to 2D multislice acquisitions.

The performance of the AFI method is based on the following assumptions. First, the two images acquired at different times should be registered to avoid motion effects. It is also assumed that the signal is insensitive to the main magnetic field non-uniformities and chemical shift effects that are canceled out when taking the signal ratio r Yarnykh, 2007. Despite some clear advantages over other B₁ mapping techniques, the application of spoiler gradients to mitigate the T₁ dependence can be a limitation due to significant levels of RF power depositions Sacolick et al., 2010. Furthermore, it is necessary to adapt the AFI pulse sequence to different scanner manufacturers, and programming experience is required to bring this technique into the clinic.

3Filtering¶

The behaviour of electromagnetic fields produced by RF antennas are bound by the laws of physics. The Maxwell equations impose many limitations on how these fields can not only vary spatially and temporally, but how the electric and magnetic fields are linked. While propagating magnetic fields interface of boundary between materials can be discontinuous (a result of Maxwell’s equations), it’s been shown in the context of MRI and tissues that the magnetic field amplitudes are expected to be smoothly varying when using clinical MRIs Sled & Pike, 1998Sled et al., 1998. At ultra-high fields, standing wave artifacts can lead to more B₁ variations and even signal nulls, however the field amplitude nonetheless varies continuously Uğurbil, 2018Vaughan et al., 2001Yang et al., 2002. Thus, for both B₁⁺ and B₁^-, their amplitude is expected to be a smoothly varying multiplicative field, and at clinical field strength it’s also expected to be a slowly or low frequency varying field.

In practice, measured B₁⁺ maps are rarely perfectly smooth over the anatomy-of-interest being imaged. Figure 4.14 shows a comparison of measured B₁ maps in the brain produced by three methods: double angle, actual flip angle imaging (AFI), and Bloch-Siegert shift.

The overall “shape” of the B₁ map is the same for all three maps, and this nonuniformity pattern is expected due to the elliptical shape of the brain and its electromagnetic properties Sled & Pike, 1998. We see in the B₁ maps of Figure 4.14 that there is some noise, some distinguishable anatomical structures (caused by T₁ sensitivity and/or k-space propagation susceptibility effects), and in one case (AFI) an artifact caused by Gibbs ringing in the acquired images. All of these variations are not present in the actual B₁⁺ field that the spins experience during a pulse sequence, and so using this “raw” B₁ map to calibration flip angles or RF power for other quantitative MRI techniques (eg. variable flip angle T₁ mapping, quantitative magnetization transfer) risks introducing errors during the correction.

Although not a perfect solution, researchers often smoothen their B₁ maps Yarnykh, 2007Lutti et al., 2010Boudreau et al., 2017 in an effort to mitigate the error propagation from the B₁⁺ map noise and artifacts prior to use for other techniques. This chapter will discuss some common ways this B₁ map smoothing is achieved, show some examples of their benefits and weaknesses, and discuss some best practices.

Filters and smoothing¶

There are two main ways that field maps are smoothened in practice: filters and fitting. The study of filters is typically presented in a signal processing context, however its basic principles (in particular, convolutions) are observed in many other fields of study, in particular physics.

We’ll begin by providing a very brief overview of some key filtering properties, then move on to some illustrative 1D examples related to MRI situations before finally returning to their applications in actual B₁ maps.

Filtering is presented as a convolution process to produce an output that is smoother, meaning less sharp edges. A convolution is the multiplication of a kernel (a predetermined function or property, such as the mean, median, Gaussian function, etc) that is shifted at each point of the signal or image, and the summed value of this multiplication is assigned to the time or spatial point where it was applied. Figure 4.16 illustrates this for the mean using a three-position mean as a kernel:

In the context of MRI, the mean is not the best choice for a filter, as it is sensitive to high values relative to the base signal. The median is a better choice, which we’ll demonstrate in the next section. In terms of equations, the convolution is shown using the symbol $\otimes$, such that analytically it is represented as:

where $f\left( t \right)$ is the signal of interest and $g\left( t \right)$ is the kernel. Not every kernel will lead to smoothing (reduction of high frequencies) of the signal of interest when convolved, however the Gaussian distribution is one such smoothing function:

where $x_{0}$ is the center position of the distribution, and σ is a measure of the width. The convolution using this function with a 9-point sample for different widths is shown in Figure 4.16.

One property of the convolution is that the convolution of two functions is the multiplication of the Fourier transforms of each function following by an inverse Fourier transform:

Although convolutions can be computed this way and may me more conceptually clear, particularly the role of the kernel, practically this ends up being slower than the convolution method when using only a small number of samples for the kernel.

As for “smoothing” the signal using fitting, splines are simply a piecewise fit of your signal to some function with a continuity condition set at different points throughout the signal. Typically this is done using polynomials, such as a third-order polynomial: $\text{a}x^{3} + \text{b}x^{2} +\text{c}x^{1} + \text{d}x^{0}$. There are a lot of different algorithms and ways to do this, which is out of scope for this work.

B1 map examples¶

Let’s revisit our initial B₁ maps in Figure 4.14 and see how they respond to the filters we’ve explored in the previous section. The double angle B₁ maps were mostly impacted by noise and structural T₁ patterns, AFI had some artifacts that were caused by Gibbs ringing in the raw images, and the Bloch-Siegert B₁ map had an artifact caused by a phase pole at the end of a fringe line. Figure 4.22 shows each of the B₁ map and the filtered maps using the median, Gaussian, and spline filtering techniques.

All three methods worked well with the double angle B₁ map, and the outputs of the median and Gaussian are most similar. The top right corner of the spline-filtered double angle B₁ map has higher intensity, likely due to an edge effect as discussed in the 1D example section. For AFI, median and gaussian filters removed most of the repeated variations, whereas spline-filtering didn’t at the medium filter strength. Lastly, for Bloch-Siegert, the median filter performed well at removing the noise and smoothing out the artifact, though some still remains. For the Gaussian and spline cases, there was a single pixel in the left that had very high value in the unfiltered images and this led to a spreading of high B₁ values to nearby voxels, something that didn’t occur for the median filter case. If either of these filters were used in an automated pipeline without quality control, inaccurate B₁ values would have been spread, which is undesireable.

Recommendations, benefits, and pitfalls¶

Overall, median, Gaussian, and spline filters perform at smoothing noisy B₁ maps. If image artifacts exists in the B₁ maps, then the choice of filter could impact the output B₁ map. In all our B₁ map examples (Figure 4.22), a 5x5 voxel median filter would have performed sufficiently all while avoiding spreading errors. This may be a relatively safe filter to try first for clinical use in the brain at clinical field strengths. In other field strengths or anatomies, or if different artifacts exist, this may not always be the case. Good care should always be applied when selecting a filter; know why you are using it, what its potential drawbacks are, and look for error spreading in the output B₁ map.

If your filtered B₁ map is intended for use at boundary edges, such as grey matter, extra precautions should be taken when applying filters and doing quality control. Know how your filters handle edges, and if needed choose to mirror or extrapolate B₁ values beyond the masked region of interest. Quality control is important, as there can be substantial edge artifacts when using filters.

Finally, remember that using the unfiltered B₁ map is also a choice, and many researchers use these. It’s important to report if you filtered or not your B₁ maps when reporting them in your research.

4Appendix A¶