Display content:


In [1]:
% **Blog post code introduction**
% 
% Congrats on activating the "All cells" option in this interactive blog post =D
%
% Below, several new HTML blocks have appears prior to the figures, displaying the Octave/MATLAB code that was used to generate the figures in this blog post.
%
% If you want to reproduce the data on your own local computer, you simply need to have qMRLab installed in your Octave/MATLAB path and run the "startup.m" file, as is shown below.
%
% If you want to get under the hood and modify the code right now, you can do so in the Jupyter Notebook of this blog post hosted on MyBinder. The link to it is in the introduction above.
In [2]:

Variable Flip Angle T1 Mapping

Variable flip angle (VFA) T1 mapping (Christensen et al. 1974; Gupta 1977; Fram et al. 1987), also known as Driven Equilibrium Single Pulse Observation of T1 (DESPOT1) (Homer & Beevers 1985; Deoni et al. 2003), is a rapid quantitative T1 measurement technique that is widely used to acquire 3D T1 maps (e.g. whole-brain) in a clinically feasible time. VFA estimates T1 values by acquiring multiple spoiled gradient echo acquisitions, each with different excitation flip angles (θn for n = 1, 2, .., N and θiθj). The steady-state signal of this pulse sequence (Figure 1) uses very short TRs (on the order of magnitude of 10 ms) and is very sensitive to T1 for a wide range of flip angles.

VFA is a technique that originates from the NMR field, and was adopted because of its time efficiency and the ability to acquire accurate T1 values simultaneously for a wide range of values (Christensen et al. 1974; Gupta 1977). For imaging applications, VFA also benefits from an increase in SNR because it can be acquired using a 3D acquisition instead of multislice, which also helps to reduce slice profile effects. One important drawback of VFA for T1 mapping is that the signal is very sensitive to inaccuracies in the flip angle value, thus impacting the T1 estimates. In practice, the nominal flip angle (i.e. the value set at the scanner) is different than the actual flip angle experienced by the spins (e.g. at 3.0 T, variations of up to ±30%), an issue that increases with field strength. VFA typically requires the acquisition of another quantitative map, the transmit RF amplitude (B1+, or B1 for short), to calibrate the nominal flip angle to its actual value because of B1 inhomogeneities that occur in most loaded MRI coils (Sled & Pike 1998). The need to acquire an additional B1 map reduces the time savings offered by VFA over saturation-recovery techniques, and inaccuracies/imprecisions of the B1 map are also propagated into the VFA T1 map (Boudreau et al. 2017; Lee et al. 2017).

Figure 1. Simplified pulse sequence diagram of a variable flip angle (VFA) pulse sequence with a gradient echo readout. TR: repetition time, θn: excitation flip angle for the nth measurement, IMG: image acquisition (k-space readout), SPOIL: spoiler gradient.

Signal Modelling

The steady-state longitudinal magnetization of an ideal variable flip angle experiment can be analytically solved from the Bloch equations for the spoiled gradient echo pulse sequence {θn–TR}:

where Mz is the longitudinal magnetization, M0 is the magnetization at thermal equilibrium, TR is the pulse sequence repetition time (Figure 1), and θn is the excitation flip angle. The Mz curves of different T1 values for a range of θn and TR values are shown in Figure 2.

Figure 2. Variable flip angle technique signal curves (Eq. 1) for three different T1 values, approximating the main types of tissue in the brain at 3T.

In [3]:
%% MATLAB/OCTAVE CODE
% Adds qMRLab to the path of the environment

cd ../qMRLab
startup
warning: function /home/jovyan/work/t1_notebooks/qMRLab/src/Common/pulse/sinc.m shadows a core library function
warning: called from
    startup at line 1 column 1
warning: function /home/jovyan/work/t1_notebooks/qMRLab/External/NODDI_toolbox_v1.0/ParforProgMonv2/example.m shadows a core library function
warning: called from
    startup at line 1 column 1
loading struct
loading io
loading statistics
loading optim
loading image
In [4]:
%% MATLAB/OCTAVE CODE
% Code used to generate the data required for Figure 4 of the blog post

clear all

%% Setup parameters
% All times are in milliseconds
% All flip angles are in degrees

TR_range = 5:5:200;

params.EXC_FA = 1:90;

%% Calculate signals
%
% To see all the options available, run `help vfa_t1.analytical_solution`

for ii = 1:length(TR_range)
    params.TR = TR_range(ii);
    
    % White matter
    params.T1 = 900; % in milliseconds

    signal_WM(ii,:) = vfa_t1.analytical_solution(params);

    % Grey matter
    params.T1 = 1500;  % in milliseconds
    signal_GM(ii,:) = vfa_t1.analytical_solution(params);

    % CSF
    params.T1 = 4000;  % in milliseconds
    signal_CSF(ii,:) = vfa_t1.analytical_solution(params);
end
In [5]:
In [6]:

From Figure 2, it is clearly seen that the flip angle at which the steady-state signal is maximized is dependent on the T1 and TR values. This flip angle is a well known quantity, called the Ernst angle (Ernst & Anderson 1966), which can be solved analytically from Equation 1 using properties of calculus:

The closed-form solution (Equation 1) makes several assumptions which in practice may not always hold true if care is not taken. Mainly, it is assumed that the longitudinal magnetization has reached a steady state after a large number of TRs, and that the transverse magnetization is perfectly spoiled at the end of each TR. Bloch simulations – a numerical approach at solving the Bloch equations for a set of spins at each time point – provide a more realistic estimate of the signal if the number of repetition times is small (i.e. a steady-state is not achieved). As can be seen from Figure 3, the number of repetitions required to reach a steady state not only depends on T1, but also on the flip angle; flip angles near the Ernst angle need more TRs to reach a steady state. Preparation pulses or an outward-in k-space acquisition pattern are typically sufficient to reach a steady state by the time that the center of k-space is acquired, which is where most of the image contrast resides.

Figure 3. Signal curves simulated using Bloch simulations (orange) for a number of repetitions ranging from 1 to 150, plotted against the ideal case (Equation 1 – blue). Simulation details: TR = 25 ms, T1 = 900 ms, 100 spins. Ideal spoiling was used for this set of Bloch simulations (transverse magnetization was set to 0 at the end of each TR).

In [7]:
%% MATLAB/OCTAVE CODE
% Code used to generate the data required for Figure 4 of the blog post

clear all

%% Setup parameters
% All times are in milliseconds
% All flip angles are in degrees

% White matter
params.T1 = 900; % in milliseconds
params.T2 = 10000;
params.TR = 25;
params.TE = 5;
params.EXC_FA = 1:90;
Nex_range = 1:1:150;

%% Calculate signals
%
% To see all the options available, run `help vfa_t1.analytical_solution`

for ii = 1:length(Nex_range)
    params.Nex = Nex_range(ii);
    
    signal_analytical(ii,:) = vfa_t1.analytical_solution(params);

    [~, complex_signal] = vfa_t1.bloch_sim(params);
    signal_blochsim(ii,:) = abs(complex(complex_signal));
end
In [8]:
In [9]: