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}:

(5)\[\begin{equation}\tag{1} M_{z}(\theta_n) = M_0 \frac{1-e^{- \frac{TR}{T_1}}}{1-\text{cos}(\theta_n) e^{- \frac{TR}{T_1}}} \text{sin}(\theta_n) \end{equation}\]

where M_z is the longitudinal magnetization, M₀ is the magnetization at thermal equilibrium, TR is the pulse sequence repetition time (Figure 1), and θ_n is the excitation flip angle. The M_z curves of different T₁ values for a range of θ_n and TR values are shown in Figure 2.

Show code cell source Hide code cell source

from repo2data.repo2data import Repo2Data
import os 
import pickle
import matplotlib.pyplot as plt
import chart_studio.plotly as py
import plotly.graph_objs as go
import numpy as np
from plotly import __version__
from plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot
from IPython.display import display, HTML

data_req_path = os.path.join("..","..", "binder", "data_requirement.json")
repo2data = Repo2Data(data_req_path)
DATA_ROOT = os.path.join(repo2data.install()[0],"t1-book-neurolibre")

Show code cell output Hide code cell output

---- repo2data starting ----
/srv/conda/envs/notebook/lib/python3.8/site-packages/repo2data
Config from file :
../../binder/data_requirement.json
Destination:
../../data/qmrlab-t1-book

Info : ../../data/qmrlab-t1-book already downloaded

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

Show code cell source Hide code cell source

filename = os.path.join(DATA_ROOT,"02",'figure_2.pkl')

with open(filename, 'rb') as f:
    params, TR_range, signal_WM, signal_GM, signal_CSF = pickle.load(f)

config={'showLink': False, 'displayModeBar': False}

init_notebook_mode(connected=True)

data1 = [dict(
        visible = False,
        mode = 'lines',
        x = params["EXC_FA"],
        y = abs(np.squeeze(np.asarray(signal_WM[ii]))),
        name = 'T<sub>1</sub> = 0.9 s (White Matter)',
        text = 'T<sub>1</sub> = 0.9 s (White Matter)',
        hoverinfo = 'x+y+text') for ii in range(len(TR_range))]

data1[4]['visible'] = True

data2 = [dict(
        visible = False,
        mode = 'lines',
        x = params["EXC_FA"],
        y = abs(np.squeeze(np.asarray(signal_GM[ii]))),
        name = 'T<sub>1</sub> = 1.5 s (Grey Matter)',
        text = 'T<sub>1</sub> = 1.5 s (Grey Matter)',
        hoverinfo = 'x+y+text') for ii in range(len(TR_range))]

data2[4]['visible'] = True

data3 = [dict(
        visible = False,
        mode = 'lines',
        x = params["EXC_FA"],
        y = abs(np.squeeze(np.asarray(signal_CSF[ii]))),
        name = 'T<sub>1</sub> = 4.0 s (Cerebrospinal Fluid)',
        text = 'T<sub>1</sub> = 4.0 s (Cerebrospinal Fluid)',
        hoverinfo = 'x+y+text') for ii in range(len(TR_range))]

data3[4]['visible'] = True

data = data1 + data2 + data3

steps = []
for i in range(len(TR_range)):
    step = dict(
        method = 'restyle',  
        args = ['visible', [False] * len(data1)],
        label = str(TR_range[i])
        )
    step['args'][1][i] = True # Toggle i'th trace to "visible"
    steps.append(step)

sliders = [dict(
    x = 0,
    y = -0.02,
    active = 2,
    currentvalue = {"prefix": "TR value (ms): <b>"},
    pad = {"t": 50, "b": 10},
    steps = steps
)]

layout = go.Layout(
    width=580,
    height=450,
    margin=go.layout.Margin(
        l=80,
        r=40,
        b=60,
        t=10,
    ),
    annotations=[
        dict(
            x=0.5004254919715793,
            y=-0.18,
            showarrow=False,
            text='Excitation Flip Angle (°)',
            font=dict(
                family='Times New Roman',
                size=22
            ),
            xref='paper',
            yref='paper'
        ),
        dict(
            x=-0.15,
            y=0.5,
            showarrow=False,
            text='Long. Magnetization (M<sub>z</sub>)',
            font=dict(
                family='Times New Roman',
                size=22
            ),
            textangle=-90,
            xref='paper',
            yref='paper'
        ),
    ],
    xaxis=dict(
        autorange=False,
        range=[0, params['EXC_FA'][-1]],
        showgrid=False,
        linecolor='black',
        linewidth=2
    ),
    yaxis=dict(
        autorange=True,
        showgrid=False,
        linecolor='black',
        linewidth=2
    ),
    legend=dict(
        x=0.5,
        y=0.9,
        traceorder='normal',
        font=dict(
            family='Times New Roman',
            size=12,
            color='#000'
        ),
        bordercolor='#000000',
        borderwidth=2
    ), 
    sliders=sliders,
    plot_bgcolor='white'
)

fig = dict(data=data, layout=layout)

plot(fig, filename = 'vfa_fig_2.html', config = config)
display(HTML('vfa_fig_2.html'))

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

(6)\[\begin{equation}\tag{2} \theta_{Ernst} = \text{acos}(e^{- \frac{TR}{T_1}}) \end{equation}\]

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 T₁, 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, T₁ = 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).

Show code cell source Hide code cell source

filename = os.path.join(DATA_ROOT,"02",'figure_3.pkl')

with open(filename, 'rb') as f:
    params, Nex_range, signal_analytical, signal_blochsim = pickle.load(f)

config={'showLink': False, 'displayModeBar': False}

init_notebook_mode(connected=True)

data1 = [dict(
        visible = False,
        mode = 'lines',
        x = params["EXC_FA"],
        y = abs(np.squeeze(np.asarray(signal_analytical[ii]))),
        name = 'Analytical Solution',
        text = 'Analytical Solution',
        hoverinfo = 'x+y+text') for ii in range(len(Nex_range))]

data1[9]['visible'] = True

data2 = [dict(
        visible = False,
        mode = 'lines',
        x = params["EXC_FA"],
        y = abs(np.squeeze(np.asarray(signal_blochsim[ii]))),
        name = 'Bloch Simulation',
        text = 'Bloch Simulation',
        hoverinfo = 'x+y+text') for ii in range(len(Nex_range))]

data2[9]['visible'] = True

data = data1 + data2

steps = []
for i in range(len(Nex_range)):
    step = dict(
        method = 'restyle',  
        args = ['visible', [False] * len(data1)],
        label = str(Nex_range[i])
        )
    step['args'][1][i] = True # Toggle i'th trace to "visible"
    steps.append(step)

sliders = [dict(
    x = 0,
    y = -0.02,
    active = 9,
    currentvalue = {"prefix": "n<sup>th</sup> TR: <b>"},
    pad = {"t": 50, "b": 10},
    steps = steps
)]

layout = go.Layout(
    width=580,
    height=450,
    margin=go.layout.Margin(
        l=80,
        r=40,
        b=60,
        t=10,
    ),
    annotations=[
        dict(
            x=0.5004254919715793,
            y=-0.18,
            showarrow=False,
            text='Excitation Flip Angle (°)',
            font=dict(
                family='Times New Roman',
                size=22
            ),
            xref='paper',
            yref='paper'
        ),
        dict(
            x=-0.15,
            y=0.5,
            showarrow=False,
            text='Signal',
            font=dict(
                family='Times New Roman',
                size=22
            ),
            textangle=-90,
            xref='paper',
            yref='paper'
        ),
    ],
    xaxis=dict(
        autorange=False,
        range=[0, params['EXC_FA'][-1]],
        showgrid=False,
        linecolor='black',
        linewidth=2
    ),
    yaxis=dict(
        autorange=True,
        showgrid=False,
        linecolor='black',
        linewidth=2
    ),
    legend=dict(
        x=0.5,
        y=0.9,
        traceorder='normal',
        font=dict(
            family='Times New Roman',
            size=12,
            color='#000'
        ),
        bordercolor='#000000',
        borderwidth=2
    ), 
    sliders=sliders,
    plot_bgcolor='white'
)

fig = dict(data=data, layout=layout)

plot(fig, filename = 'vfa_fig_3.html', config = config)
display(HTML('vfa_fig_3.html'))

Sufficient spoiling is likely the most challenging parameter to control for in a VFA experiment. A combination of both gradient spoiling and RF phase spoiling [Bernstein et al., 2004, Zur et al., 1991] are typically recommended (Figure 4). It has also been shown that the use of very strong gradients, introduces diffusion effects (not considered in Figure 4), further improving the spoiling efficacy in the VFA pulse sequence [Yarnykh, 2010].

Figure 4. Signal curves estimated using Bloch simulations for three categories of signal spoiling: (1) ideal spoiling (blue), gradient & RF Spoiling (orange), and no spoiling (green). Simulations details: TR = 25 ms, T₁ = 900 ms, T_e = 100 ms, TE = 5 ms, 100 spins. For the ideal spoiling case, the transverse magnetization is set to zero at the end of each TR. For the gradient & RF spoiling case, each spin is rotated by different increments of phase (2𝜋 / # of spins) to simulate complete decoherence from gradient spoiling, and the RF phase of the excitation pulse is ɸ_n = ɸ_n-1 + nɸ₀ = ½ ɸ₀(n² + n + 2) [Bernstein et al., 2004] with ɸ₀ = 117° [Zur et al., 1991] after each TR.

Show code cell source Hide code cell source

filename = os.path.join(DATA_ROOT,"02",'figure_4.pkl')

with open(filename, 'rb') as f:
    params, Nex_range, signal_ideal_spoil, signal_optimal_crush_and_rf_spoil, signal_no_gradient_and_rf_spoil = pickle.load(f)

config={'showLink': False, 'displayModeBar': False}
init_notebook_mode(connected=True)

data1 = [dict(
        visible = False,
        mode = 'lines',
        x = params["EXC_FA"],
        y = abs(np.squeeze(np.asarray(signal_ideal_spoil[ii]))),
        name = 'Ideal Spoiling',
        text = 'Ideal Spoiling',
        hoverinfo = 'x+y+text') for ii in range(len(Nex_range))]

data1[10]['visible'] = True

data2 = [dict(
        visible = False,
        mode = 'lines',
        x = params["EXC_FA"],
        y = abs(np.squeeze(np.asarray(signal_optimal_crush_and_rf_spoil[ii]))),
        name = 'Gradient & RF Spoiling',
        text = 'Gradient & RF Spoiling',
        hoverinfo = 'x+y+text') for ii in range(len(Nex_range))]

data2[10]['visible'] = True

data3 = [dict(
        visible = False,
        mode = 'lines',
        x = params["EXC_FA"],
        y = abs(np.squeeze(np.asarray(signal_no_gradient_and_rf_spoil[ii]))),
        name = 'No Spoiling',
        text = 'No Spoiling',
        hoverinfo = 'x+y+text') for ii in range(len(Nex_range))]

data3[10]['visible'] = True

data = data1 + data2+ data3

steps = []
for i in range(len(Nex_range)):
    step = dict(
        method = 'restyle',  
        args = ['visible', [False] * len(data1)],
        label = str(Nex_range[i])
        )
    step['args'][1][i] = True # Toggle i'th trace to "visible"
    steps.append(step)

sliders = [dict(
    x = 0,
    y = -0.02,
    active = 10,
    currentvalue = {"prefix": "n<sup>th</sup> TR: <b>"},
    pad = {"t": 50, "b": 10},
    steps = steps
)]

layout = go.Layout(
    width=580,
    height=450,
    margin=go.layout.Margin(
        l=80,
        r=40,
        b=60,
        t=10,
    ),
    annotations=[
        dict(
            x=0.5004254919715793,
            y=-0.18,
            showarrow=False,
            text='Excitation Flip Angle (°)',
            font=dict(
                family='Times New Roman',
                size=22
            ),
            xref='paper',
            yref='paper'
        ),
        dict(
            x=-0.15,
            y=0.5,
            showarrow=False,
            text='Signal',
            font=dict(
                family='Times New Roman',
                size=22
            ),
            textangle=-90,
            xref='paper',
            yref='paper'
        ),
    ],
    xaxis=dict(
        autorange=False,
        range=[0, params['EXC_FA'][-1]],
        showgrid=False,
        linecolor='black',
        linewidth=2
    ),
    yaxis=dict(
        autorange=True,
        showgrid=False,
        linecolor='black',
        linewidth=2
    ),
    legend=dict(
        x=0.5,
        y=0.9,
        traceorder='normal',
        font=dict(
            family='Times New Roman',
            size=12,
            color='#000'
        ),
        bordercolor='#000000',
        borderwidth=2
    ), 
    sliders=sliders,
    plot_bgcolor='white'
)

fig = dict(data=data, layout=layout)

plot(fig, filename = 'vfa_fig_4.html', config = config)
display(HTML('vfa_fig_4.html'))

Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 1.

% Verbosity level 0 overrides the disp function and supresses warnings.
% Once executed, they cannot be restored in this session
% (kernel needs to be restarted or a new notebook opened.)
VERBOSITY_LEVEL = 0;

if VERBOSITY_LEVEL==0
    % This hack was used to supress outputs from external tools
    % in the Jupyter Book.
    function disp(x)
    end
    warning('off','all')
end

try
    cd qMRLab
catch
    try
        cd ../../../qMRLab
    catch
        cd ../qMRLab
    end
end

startup
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

Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 3.

% Verbosity level 0 overrides the disp function and supresses warnings.
% Once executed, they cannot be restored in this session
% (kernel needs to be restarted or a new notebook opened.)
VERBOSITY_LEVEL = 0;

if VERBOSITY_LEVEL==0
    % This hack was used to supress outputs from external tools
    % in the Jupyter Book.
    function disp(x)
    end
    warning('off','all')
end

try
    cd qMRLab
catch
    try
        cd ../../../qMRLab
    catch
        cd ../qMRLab
    end
end

startup
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

Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 4.

% Verbosity level 0 overrides the disp function and supresses warnings.
% Once executed, they cannot be restored in this session
% (kernel needs to be restarted or a new notebook opened.)
VERBOSITY_LEVEL = 0;

if VERBOSITY_LEVEL==0
    % This hack was used to supress outputs from external tools
    % in the Jupyter Book.
    function disp(x)
    end
    warning('off','all')
end

try
    cd qMRLab
catch
    try
        cd ../../../qMRLab
    catch
        cd ../qMRLab
    end
end

startup
clear all

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

% White matter
params.T1 = 900; % in milliseconds
params.T2 = 100;
params.TR = 25;
params.TE = 5;
params.EXC_FA = 1:90;
Nex_range = [1:9, 10:10:100];

%% 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);
    
    params.crushFlag = 1;
    
    [~, complex_signal] = vfa_t1.bloch_sim(params);
    signal_ideal_spoil(ii,:) = abs(complex_signal);
    
    
    params.inc = 117;
    params.partialDephasing = 1;
    params.partialDephasingFlag = 1;
    params.crushFlag = 0;
    
    [~, complex_signal] = vfa_t1.bloch_sim(params);
    signal_optimal_crush_and_rf_spoil(ii,:) = abs(complex_signal);
    
    params.inc = 0;
    params.partialDephasing = 0;

    [~, complex_signal] = vfa_t1.bloch_sim(params);
    signal_no_gradient_and_rf_spoil(ii,:) = abs(complex_signal);
end

References

1(1,2)

Matt A Bernstein, Kevin F King, and Xiaohong Joe Zhou. Handbook of MRI Pulse Sequences. Elsevier, 2004. doi:10.1016/b978-0-12-092861-3.x5000-6.

2

R. R. Ernst and W. A. Anderson. Application of fourier transform spectroscopy to magnetic resonance. Review of Scientific Instruments, 37(1):93–102, January 1966. doi:10.1063/1.1719961.

3

Vasily L. Yarnykh. Optimal radiofrequency and gradient spoiling for improved accuracy of t1 and b1 measurements using fast steady-state techniques. Magnetic Resonance in Medicine, 63(6):1610–1626, May 2010. doi:10.1002/mrm.22394.

4(1,2)

Y. Zur, M. L. Wood, and L. J. Neuringer. Spoiling of transverse magnetization in steady-state sequences. Magnetic Resonance in Medicine, 21(2):251–263, October 1991. doi:10.1002/mrm.1910210210.