diff --git a/content/intro/intro.md b/Introduction to qMRI/01-intro.md similarity index 97% rename from content/intro/intro.md rename to Introduction to qMRI/01-intro.md index f1b4c06..c2d8196 100644 --- a/content/intro/intro.md +++ b/Introduction to qMRI/01-intro.md @@ -1,5 +1,5 @@ --- -title: Foreword +title: Basics of MRI and qMRI date: 2024-07-25 authors: - name: Agah Karakuzu diff --git a/content/intro/million_dollar.md b/Introduction to qMRI/02-million_dollar.md similarity index 100% rename from content/intro/million_dollar.md rename to Introduction to qMRI/02-million_dollar.md diff --git a/content/intro/journey.md b/Introduction to qMRI/03-journey.md similarity index 100% rename from content/intro/journey.md rename to Introduction to qMRI/03-journey.md diff --git a/content/intro/encode.md b/Introduction to qMRI/04-encode.md similarity index 100% rename from content/intro/encode.md rename to Introduction to qMRI/04-encode.md diff --git a/content/intro/sequences.md b/Introduction to qMRI/05-sequences.md similarity index 100% rename from content/intro/sequences.md rename to Introduction to qMRI/05-sequences.md diff --git a/content/intro/qmri_practice.md b/Introduction to qMRI/06-qmri_practice.md similarity index 100% rename from content/intro/qmri_practice.md rename to Introduction to qMRI/06-qmri_practice.md diff --git a/content/intro/img/int_fig1.jpg b/Introduction to qMRI/img/int_fig1.jpg similarity index 100% rename from content/intro/img/int_fig1.jpg rename to Introduction to qMRI/img/int_fig1.jpg diff --git a/content/intro/img/int_fig10.jpg b/Introduction to qMRI/img/int_fig10.jpg similarity index 100% rename from content/intro/img/int_fig10.jpg rename to Introduction to qMRI/img/int_fig10.jpg diff --git a/content/intro/img/int_fig11.jpg b/Introduction to qMRI/img/int_fig11.jpg similarity index 100% rename from content/intro/img/int_fig11.jpg rename to Introduction to qMRI/img/int_fig11.jpg diff --git a/content/intro/img/int_fig12.jpg b/Introduction to qMRI/img/int_fig12.jpg similarity index 100% rename from content/intro/img/int_fig12.jpg rename to Introduction to qMRI/img/int_fig12.jpg diff --git a/content/intro/img/int_fig13.jpg b/Introduction to qMRI/img/int_fig13.jpg similarity index 100% rename from content/intro/img/int_fig13.jpg rename to Introduction to qMRI/img/int_fig13.jpg diff --git a/content/intro/img/int_fig14.jpg b/Introduction to qMRI/img/int_fig14.jpg similarity index 100% rename from content/intro/img/int_fig14.jpg rename to Introduction to qMRI/img/int_fig14.jpg diff --git a/content/intro/img/int_fig15.jpg b/Introduction to qMRI/img/int_fig15.jpg similarity index 100% rename from content/intro/img/int_fig15.jpg rename to Introduction to qMRI/img/int_fig15.jpg diff --git a/content/intro/img/int_fig16.jpg b/Introduction to qMRI/img/int_fig16.jpg similarity index 100% rename from content/intro/img/int_fig16.jpg rename to Introduction to qMRI/img/int_fig16.jpg diff --git a/content/intro/img/int_fig17.jpg b/Introduction to qMRI/img/int_fig17.jpg similarity index 100% rename from content/intro/img/int_fig17.jpg rename to Introduction to qMRI/img/int_fig17.jpg diff --git a/content/intro/img/int_fig18.png b/Introduction to qMRI/img/int_fig18.png similarity index 100% rename from content/intro/img/int_fig18.png rename to Introduction to qMRI/img/int_fig18.png diff --git a/content/intro/img/int_fig19.png b/Introduction to qMRI/img/int_fig19.png similarity index 100% rename from content/intro/img/int_fig19.png rename to Introduction to qMRI/img/int_fig19.png diff --git a/content/intro/img/int_fig2.jpg b/Introduction to qMRI/img/int_fig2.jpg similarity index 100% rename from content/intro/img/int_fig2.jpg rename to Introduction to qMRI/img/int_fig2.jpg diff --git a/content/intro/img/int_fig20.jpg b/Introduction to qMRI/img/int_fig20.jpg similarity index 100% rename from content/intro/img/int_fig20.jpg rename to Introduction to qMRI/img/int_fig20.jpg diff --git a/content/intro/img/int_fig21.jpg b/Introduction to qMRI/img/int_fig21.jpg similarity index 100% rename from content/intro/img/int_fig21.jpg rename to Introduction to qMRI/img/int_fig21.jpg diff --git a/content/intro/img/int_fig22.jpg b/Introduction to qMRI/img/int_fig22.jpg similarity index 100% rename from content/intro/img/int_fig22.jpg rename to Introduction to qMRI/img/int_fig22.jpg diff --git a/content/intro/img/int_fig23.jpg b/Introduction to qMRI/img/int_fig23.jpg similarity index 100% rename from content/intro/img/int_fig23.jpg rename to Introduction to qMRI/img/int_fig23.jpg diff --git a/content/intro/img/int_fig24.jpg b/Introduction to qMRI/img/int_fig24.jpg similarity index 100% rename from content/intro/img/int_fig24.jpg rename to Introduction to qMRI/img/int_fig24.jpg diff --git a/content/intro/img/int_fig25.jpg b/Introduction to qMRI/img/int_fig25.jpg similarity index 100% rename from content/intro/img/int_fig25.jpg rename to Introduction to qMRI/img/int_fig25.jpg diff --git a/content/intro/img/int_fig3.jpg b/Introduction to qMRI/img/int_fig3.jpg similarity index 100% rename from content/intro/img/int_fig3.jpg rename to Introduction to qMRI/img/int_fig3.jpg diff --git a/content/intro/img/int_fig4.jpg b/Introduction to qMRI/img/int_fig4.jpg similarity index 100% rename from content/intro/img/int_fig4.jpg rename to Introduction to qMRI/img/int_fig4.jpg diff --git a/content/intro/img/int_fig5.jpg b/Introduction to qMRI/img/int_fig5.jpg similarity index 100% rename from content/intro/img/int_fig5.jpg rename to Introduction to qMRI/img/int_fig5.jpg diff --git a/content/intro/img/int_fig6.jpg b/Introduction to qMRI/img/int_fig6.jpg similarity index 100% rename from content/intro/img/int_fig6.jpg rename to Introduction to qMRI/img/int_fig6.jpg diff --git a/content/intro/img/int_fig7.jpg b/Introduction to qMRI/img/int_fig7.jpg similarity index 100% rename from content/intro/img/int_fig7.jpg rename to Introduction to qMRI/img/int_fig7.jpg diff --git a/content/intro/img/int_fig8.jpg b/Introduction to qMRI/img/int_fig8.jpg similarity index 100% rename from content/intro/img/int_fig8.jpg rename to Introduction to qMRI/img/int_fig8.jpg diff --git a/content/intro/img/int_fig9.jpg b/Introduction to qMRI/img/int_fig9.jpg similarity index 100% rename from content/intro/img/int_fig9.jpg rename to Introduction to qMRI/img/int_fig9.jpg diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/01-introduction.md b/T1 Mapping/Inversion Recovery T1 Mapping/01-introduction.md new file mode 100644 index 0000000..14eeff7 --- /dev/null +++ b/T1 Mapping/Inversion Recovery T1 Mapping/01-introduction.md @@ -0,0 +1,25 @@ +--- +title: Abstract +subtitle: Inversion Recovery +date: 2024-07-25 +authors: + - name: Mathieu Boudreau + affiliations: + - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada +numbering: + heading_2: false + figure: + template: Fig. %s +--- + +## Inversion Recovery T1 Mapping + +Widely considered the gold standard for T1 mapping, the inversion recovery technique estimates T1 values by fitting the signal recovery curve acquired at different delays after an inversion pulse (180°). In a typical inversion recovery experiment ([](#irFig1)), the magnetization at thermal equilibrium is inverted using a 180° RF pulse. After the longitudinal magnetization recovers through spin-lattice relaxation for predetermined delay (`inversion time`, TI), a 90° excitation pulse is applied, followed by a readout imaging sequence (typically a spin-echo or gradient-echo readout) to create a snapshot of the longitudinal magnetization state at that TI. + +Inversion recovery was first developed for NMR in the 1940s [@Hahn1949;@Drain1949], and the first T1 map was acquired using a saturation-recovery technique (90° as a preparation pulse instead of 180°) by [@Pykett1978]. Some distinct advantages of inversion recovery are its large dynamic range of signal change and an insensitivity to pulse sequence parameter imperfections [@Stikov2015]. Despite its proven robustness at measuring T1, inversion recovery is scarcely used in practice, because conventional implementations require repetition times (TRs) on the order of 2 to 5 T1 [@Steen1994], making it challenging to acquire whole-organ T1 maps in a clinically feasible time. Nonetheless, it is continuously used as a reference measurement during the development of new techniques, or when comparing different T1 mapping techniques, and several variations of the inversion recovery technique have been developed, making it practical for some applications [@Messroghli2004;@Piechnik2010]. + +```{figure} img/ir_pulsesequences.svg +:label: irFig1 + +Figure 1. Pulse sequence of an inversion recovery experiment. +``` \ No newline at end of file diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/02-IR_SignalModelling.ipynb b/T1 Mapping/Inversion Recovery T1 Mapping/02-IR_SignalModelling.ipynb new file mode 100644 index 0000000..f6ebffe --- /dev/null +++ b/T1 Mapping/Inversion Recovery T1 Mapping/02-IR_SignalModelling.ipynb @@ -0,0 +1,562 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "4c1c4b7f-505f-4112-8b4e-59fdd766eac2", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "---\n", + "title: Signal Modelling\n", + "subtitle: Inversion Recovery\n", + "date: 2024-07-25\n", + "authors:\n", + " - name: Mathieu Boudreau\n", + " affiliations:\n", + " - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada\n", + "numbering:\n", + " figure:\n", + " template: Fig. %s\n", + "---\n", + "\n", + "The steady-state longitudinal magnetization of an inversion recovery experiment can be derived from the Bloch equations for the pulse sequence {θ180 – TI – θ90 – (TR-TI)}, and is given by:\n", + "\n", + "\\begin{equation}\\label{eq:1}\n", + "M_{z}(TI) = M_0 \\frac{1-\\text{cos}(\\theta_{180})e^{- \\frac{TR}{T_1}} -[1-\\text{cos}(\\theta_{180})]e^{- \\frac{TI}{T_1}}}{1 - \\text{cos}(\\theta_{180}) \\text{cos}(\\theta_{90}) e^{- \\frac{TR}{T_1}}}\n", + "\\end{equation}\n", + "\n", + "where Mz is the longitudinal magnetization prior to the θ90 pulse. If the in-phase real signal is desired, it can be calculated by multiplying Eq. 1 by ksin(θ90)e-TE/T2, where k is a constant. This general equation can be simplified by grouping together the constants for each measurements regardless of their values (i.e. at each TI, same TE and θ90 are used) and assuming an ideal inversion pulse:\n", + "\n", + "\\begin{equation}\\label{eq:2}\n", + "M_z(TI) = C(1-2e^{- \\frac{TI}{T_1}} + e^{- \\frac{TR}{T_1}})\n", + "\\end{equation}\n", + "\n", + "where the first three terms and the denominator of Eq. 1 have been grouped together into the constant C. If the experiment is designed such that TR is long enough to allow for full relaxation of the magnetization (TR > 5T1), we can do an additional approximation by dropping the last term in Eq. 2:\n", + "\n", + "\\begin{equation}\\label{eq:3}\n", + "M_z(TI) = C(1-2e^{- \\frac{TI}{T_1}})\n", + "\\end{equation}\n", + "\n", + "The simplicity of the signal model described by Eq. 3, both in its equation and experimental implementation, has made it the most widely used equation to describe the signal evolution in an inversion recovery T1 mapping experiment. The magnetization curves are plotted in Figure 2 for approximate T1 values of three different tissues in the brain. Note that in many practical implementations, magnitude-only images are acquired, so the signal measured would be proportional to the absolute value of Eq. 3.\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "c2475527-acd3-47cf-a65b-10902e7d5866", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-output", + "hide-input" + ] + }, + "outputs": [], + "source": [ + "from repo2data.repo2data import Repo2Data\n", + "import os \n", + "import pickle\n", + "import matplotlib.pyplot as plt\n", + "import chart_studio.plotly as py\n", + "import plotly.graph_objs as go\n", + "import numpy as np\n", + "from plotly import __version__\n", + "from plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot\n", + "from IPython.display import display, HTML\n", + "from plotly import tools\n", + "\n", + "data_req_path = os.path.join(\"..\",\"..\", \"binder\", \"data_requirement.json\")\n", + "repo2data = Repo2Data(data_req_path)\n", + "DATA_ROOT = os.path.join(repo2data.install()[0],\"t1-book-neurolibre\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "64331b70-ae8a-42a6-af7f-b275c36c4580", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"01\",'figure_2.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " params, signal_WM, signal_GM, signal_CSF = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "\n", + "init_notebook_mode(connected=True)\n", + "\n", + "wm = go.Scatter(\n", + " x = params[\"TI\"],\n", + " y = signal_WM,\n", + " name = 'T1 = 0.9 s (White Matter)',\n", + " text = 'T1 = 0.9 s (White Matter)',\n", + " hoverinfo = 'x+y+text'\n", + ")\n", + "\n", + "gm = go.Scatter(\n", + " x = params[\"TI\"],\n", + " y = signal_GM,\n", + " name = 'T1 = 1.5 s (Grey Matter)',\n", + " text = 'T1 = 1.5 s (Grey Matter)',\n", + " hoverinfo = 'x+y+text'\n", + ")\n", + "\n", + "csf = go.Scatter(\n", + " x = params[\"TI\"],\n", + " y = signal_CSF,\n", + " name = 'T1 = 4.0 s (Cerebrospinal Fluid)',\n", + " text = 'T1 = 4.0 s (Cerebrospinal Fluid)',\n", + " hoverinfo = 'x+y+text'\n", + ")\n", + "\n", + "data = [wm, gm, csf]\n", + "\n", + "layout = go.Layout(\n", + " width=600,\n", + " height=350,\n", + " margin=go.layout.Margin(\n", + " l=100,\n", + " r=50,\n", + " b=60,\n", + " t=0,\n", + " ),\n", + " annotations=[\n", + " dict(\n", + " x=0.5004254919715793,\n", + " y=-0.175,\n", + " showarrow=False,\n", + " text='Inversion Time – TI (ms)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=-0.15,\n", + " y=0.50,\n", + " showarrow=False,\n", + " text='Long. Magnetization (Mz)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " ],\n", + " xaxis=dict(\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " yaxis=dict(\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " legend=dict(\n", + " x=0.55,\n", + " y=0.15,\n", + " traceorder='normal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=12,\n", + " color='#000'\n", + " ),\n", + " bordercolor='#000000',\n", + " borderwidth=2\n", + " ),\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'ir_fig_2.html', config = config)\n", + "\n", + "display(HTML('ir_fig_2.html'))" + ] + }, + { + "cell_type": "markdown", + "id": "ab705653-9cf7-4495-865b-90e95951e3d2", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "```{figure} img/plot1.png\n", + ":label: irPlot1\n", + "\n", + "Inversion recovery curves (Eq. 2) for three different T1 values, approximating the main types of tissue in the brain.\n", + "```\n", + "\n", + "Practically, Eq. 1 is the better choice for simulating the signal of an inversion recovery experiment, as the TRs are often chosen to be greater than 5T1 of the tissue-of-interest, which rarely coincides with the longest T1 present (e.g. TR may be sufficiently long for white matter, but not for CSF which could also be present in the volume). Equation 3 also assumes ideal inversion pulses, which is rarely the case due to slice profile effects. Figure 3 displays the inversion recovery signal magnitude (complete relaxation normalized to 1) of an experiment with TR = 5 s and T1 values ranging between 250 ms to 5 s, calculated using both equations." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "2f86f9e3-8855-4312-be1f-9637120665e1", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"01\",'figure_3.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " T1_range, signal_T1_Eq1, signal_T1_Eq3 = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "\n", + "init_notebook_mode(connected=True)\n", + "\n", + "data1 = [dict(\n", + " visible = False,\n", + " x = params[\"TI\"],\n", + " y = abs(signal_T1_Eq3[ii]),\n", + " name = '[Eq. 3] – Long TR approximation',\n", + " text = '[Eq. 3] – Long TR approximation',\n", + " hoverinfo = 'x+y+text') for ii in range(len(T1_range))]\n", + "\n", + "data1[3]['visible'] = True\n", + "\n", + "data2 = [dict(\n", + " visible = False,\n", + " x = params[\"TI\"],\n", + " y = abs(signal_T1_Eq1[ii]),\n", + " line = dict(\n", + " color = ('rgb(22, 96, 167)'),\n", + " dash = 'dash'),\n", + " name = '[Eq. 1] – General Equation',\n", + " text = '[Eq. 1] – General Equation',\n", + " hoverinfo = 'x+y+text') for ii in range(len(T1_range))]\n", + "\n", + "data2[3]['visible'] = True\n", + "\n", + "data = data1 + data2\n", + "\n", + "steps = []\n", + "for i in range(len(T1_range)):\n", + " step = dict(\n", + " method = 'restyle', \n", + " args = ['visible', [False] * len(data1)],\n", + " label = str(T1_range[i])\n", + " )\n", + " step['args'][1][i] = True # Toggle i'th trace to \"visible\"\n", + " steps.append(step)\n", + "\n", + "sliders = [dict(\n", + " x = 0,\n", + " y = -0.0,\n", + " active = 3,\n", + " currentvalue = {\"prefix\": \"T1 value (s): \"},\n", + " pad = {\"t\": 50, \"b\": 10},\n", + " steps = steps\n", + ")]\n", + "\n", + "layout = go.Layout(\n", + " width=580,\n", + " height=400,\n", + " margin=go.layout.Margin(\n", + " l=80,\n", + " r=40,\n", + " b=60,\n", + " t=10,\n", + " ),\n", + " annotations=[\n", + " dict(\n", + " x=0.5004254919715793,\n", + " y=-0.2,\n", + " showarrow=False,\n", + " text='Inversion Time – TI (ms)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=-0.14,\n", + " y=0.5,\n", + " showarrow=False,\n", + " text='Signal (magnitude)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " ],\n", + " xaxis=dict(\n", + " autorange=False,\n", + " range=[0, 5],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " yaxis=dict(\n", + " autorange=False,\n", + " range=[0, 1],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " legend=dict(\n", + " x=0.5,\n", + " y=0.5,\n", + " traceorder='normal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=12,\n", + " color='#000'\n", + " ),\n", + " bordercolor='#000000',\n", + " borderwidth=2\n", + " ), \n", + " sliders=sliders,\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'ir_fig_3.html', config = config)\n", + "display(HTML('ir_fig_3.html'))" + ] + }, + { + "cell_type": "markdown", + "id": "00b03f12", + "metadata": {}, + "source": [ + "```{figure} img/plot2.png\n", + ":label: irPlot2\n", + "\n", + "Signal recovery curves simulated using Eq. 3 (solid) and Eq. 1 (dotted) with a TR = 5 s for T1 values ranging between 0.25 to 5 s.\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "3405434d", + "metadata": { + "editable": false, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "\n", + "````{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 2.\n", + ":class: tip, dropdown\n", + "\n", + "```octave\n", + "\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "clear all\n", + "\n", + "%% Setup parameters\n", + "% All times are in seconds\n", + "% All flip angles are in degrees\n", + "\n", + "params.TR = 5.0;\n", + "params.TI = linspace(0.001, params.TR, 1000);\n", + " \n", + "params.TE = 0.004;\n", + "params.T2 = 0.040;\n", + " \n", + "params.EXC_FA = 90; % Excitation flip angle\n", + "params.INV_FA = 180; % Inversion flip angle\n", + "\n", + "params.signalConstant = 1;\n", + "\n", + "%% Calculate signals\n", + "%\n", + "% The option 'GRE-IR' selects the analytical equations for the\n", + "% gradient echo readout inversion recovery experiment The option\n", + "% '4' is a flag that selects the long TR approximation of the \n", + "% analytical solution (TR>>T1), Eq. 3 of the blog post.\n", + "%\n", + "% To see all the options available, run:\n", + "% `help inversion_recovery.analytical_solution`\n", + "\n", + "% White matter\n", + "params.T1 = 0.900; % in seconds\n", + "\n", + "signal_WM = inversion_recovery.analytical_solution(params, 'GRE-IR', 4);\n", + "\n", + "% Grey matter\n", + "params.T1 = 1.500; % in seconds\n", + "signal_GM = inversion_recovery.analytical_solution(params, 'GRE-IR', 4);\n", + "\n", + "% CSF\n", + "params.T1 = 4.000; % in seconds\n", + "signal_CSF = inversion_recovery.analytical_solution(params, 'GRE-IR', 4);\n", + "```\n", + "\n", + "````" + ] + }, + { + "cell_type": "markdown", + "id": "0a19cc2f", + "metadata": { + "editable": false, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "\n", + "```{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 3.\n", + ":class: tip, dropdown\n", + "\n", + "```octave\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "clear all\n", + "\n", + "%% Setup parameters\n", + "% All times are in seconds\n", + "% All flip angles are in degrees\n", + "\n", + "params.TR = 5.0;\n", + "params.TI = linspace(0.001, params.TR, 1000);\n", + " \n", + "params.TE = 0.004;\n", + "params.T2 = 0.040;\n", + " \n", + "params.EXC_FA = 90; % Excitation flip angle\n", + "params.INV_FA = 180; % Inversion flip angle\n", + "\n", + "params.signalConstant = 1;\n", + "\n", + "T1_range = 0.25:0.25:5; % in seconds\n", + "\n", + "%% Calculate signals\n", + "%\n", + "% The option 'GRE-IR' selects the analytical equations for the\n", + "% gradient echo readout inversion recovery experiment. The option\n", + "% '1' is a flag that selects full analytical solution equation \n", + "% (no approximation), Eq. 1 of the blog post. The option '4' is a\n", + "% flag that selects the long TR approximation of the analytical \n", + "% solution (TR>>T1), Eq. 3 of the blog post.\n", + "%\n", + "% To see all the options available, run:\n", + "% `help inversion_recovery.analytical_solution`\n", + "\n", + "for ii = 1:length(T1_range)\n", + " params.T1 = T1_range(ii);\n", + " \n", + " signal_T1_Eq1{ii} = inversion_recovery.analytical_solution(params, 'GRE-IR', 1);\n", + "\n", + " signal_T1_Eq3{ii} = inversion_recovery.analytical_solution(params, 'GRE-IR', 4);\n", + "end\n", + "```\n", + "\n", + "```" + ] + } + ], + "metadata": { + "celltoolbar": "Tags", + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.18" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/03-IR_DataFitting.ipynb b/T1 Mapping/Inversion Recovery T1 Mapping/03-IR_DataFitting.ipynb new file mode 100644 index 0000000..ed0f66f --- /dev/null +++ b/T1 Mapping/Inversion Recovery T1 Mapping/03-IR_DataFitting.ipynb @@ -0,0 +1,641 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "6644351b-4758-40f7-95e1-451cf676698f", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "---\n", + "title: Data Fitting\n", + "subtitle: Inversion Recovery\n", + "date: 2024-07-25\n", + "authors:\n", + " - name: Mathieu Boudreau\n", + " affiliations:\n", + " - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada\n", + "numbering:\n", + " figure:\n", + " template: Fig. %s\n", + "---" + ] + }, + { + "cell_type": "markdown", + "id": "d1efeaa4-45fa-409e-8c04-c93fdb7b01b3", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "Several factors impact the choice of the inversion recovery fitting algorithm. If only magnitude images are available, then a polarity-inversion is often implemented to restore the non-exponential magnitude curves (Figure 3) into the exponential form (Figure 2). This process is sensitive to noise due to the Rician noise creating a non-zero level at the signal null. If phase data is also available, then a phase term must be added to the fitting equation {cite:p}`Barral2010-qm`. Equation 3 must only be used to fit data for the long TR regime (TR > 5T1), which in practice is rarely satisfied for all tissues in subjects.\n", + "\n", + "Early implementations of inversion recovery fitting algorithms were designed around the computational power available at the time. These included the “null method” {cite:p}`Pykett1983`, assuming that each T1 value has unique zero-crossings (see Figure 2), and linear fitting of a rearranged version of Eq. 3 on a semi-log plot {cite:p}`Fukushima1981`. Nowadays, a non-linear least-squares fitting algorithm (e.g. Levenberg-Marquardt) is more appropriate, and can be applied to either approximate or general forms of the signal model (Eq. 3 or Eq. 1). More recent work {cite:p}`Barral2010-qm` demonstrated that T1 maps can also be fitted much faster (up to 75 times compared to Levenberg-Marquardt) to fit Eq. 1 – without a precision penalty – by using a reduced-dimension non-linear least squares (RD-NLS) algorithm. It was demonstrated that the following simplified 5-parameter equation can be sufficient for accurate T1 mapping:\n", + "\n", + "\\begin{equation}\\label{eq:4}\n", + "S(TI) = a + be^{- \\frac{TI}{T_1}}\n", + "\\end{equation}\n", + "\n", + "where a and b are complex values. If magnitude-only data is available, a 3-parameter model can be sufficient by taking the absolute value of Eq. 4. While the RD-NLS algorithms are too complex to be presented here (the reader is referred to the paper, (Barral et al. 2010)), the code for these algorithms [was released open-source](http://www-mrsrl.stanford.edu/~jbarral/t1map.html) along with the original publication, and is also available as a [qMRLab](https://github.com/qMRLab/qMRLab) T1 mapping model. One important thing to note about Eq. 4 is that it is general – no assumption is made about TR – and is thus as robust as Eq. 1 as long as all pulse sequence parameters other than TI are kept constant between each measurement. Figure 4 compares simulated data (Eq. 1) using a range of TRs (1.5T1 to 5T1) fitted using either RD-NLS & Eq. 4 or a Levenberg-Marquardt fit of Eq. 2." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "0a8f7c89-2950-466d-9615-8b0921d44a51", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-output", + "hide-input" + ] + }, + "outputs": [], + "source": [ + "from repo2data.repo2data import Repo2Data\n", + "import os \n", + "import pickle\n", + "import matplotlib.pyplot as plt\n", + "import chart_studio.plotly as py\n", + "import plotly.graph_objs as go\n", + "import numpy as np\n", + "from plotly import __version__\n", + "from plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot\n", + "from IPython.display import display, HTML\n", + "from plotly import tools\n", + "\n", + "data_req_path = os.path.join(\"..\",\"..\", \"binder\", \"data_requirement.json\")\n", + "repo2data = Repo2Data(data_req_path)\n", + "DATA_ROOT = os.path.join(repo2data.install()[0],\"t1-book-neurolibre\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "d0d95439-8f28-4d1d-98e6-1dde04cf5e11", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"01\",'figure_4.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " params, Mz_analytical, fitOutput_lm, fitOutput_barral, TR_range = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "\n", + "init_notebook_mode(connected=True)\n", + "\n", + "data1 = [dict(\n", + " visible = False,\n", + " mode = 'markers',\n", + " x = params[\"TI\"],\n", + " y = abs(np.squeeze(np.asarray(Mz_analytical[ii]))),\n", + " name = 'Simulated data',\n", + " text = 'Simulated data',\n", + " hoverinfo = 'x+y+text') for ii in range(len(TR_range))]\n", + "\n", + "data1[10]['visible'] = True\n", + "\n", + "data2 = [dict(\n", + " visible = False,\n", + " mode = 'lines',\n", + " x = params[\"TI\"],\n", + " y = abs(fitOutput_lm[ii]['c'] * (1 - 2*np.exp(-params['TI']/fitOutput_lm[ii]['T1']))),\n", + " name = '[C(1-2e-TI/T1)] Fitted T1: ' + str(round(fitOutput_lm[ii]['T1'])) + ' ms',\n", + " text = '[C(1-2e-TI/T1)]',\n", + " hoverinfo = 'x+y+text') for ii in range(len(TR_range))]\n", + "\n", + "data2[10]['visible'] = True\n", + "\n", + "data3 = [dict(\n", + " visible = False,\n", + " mode = 'lines',\n", + " x = params[\"TI\"],\n", + " y = abs((fitOutput_barral[ii]['ra']+fitOutput_barral[ii]['rb']*np.exp(-params['TI']/fitOutput_barral[ii]['T1']))),\n", + " name = '[a+be-TI/T1] Fitted T1: ' + str(round(fitOutput_barral[ii]['T1'])) + ' ms',\n", + " text = '[a+be-TI/T1]',\n", + " hoverinfo = 'x+y+text') for ii in range(len(TR_range))]\n", + "\n", + "data3[10]['visible'] = True\n", + "\n", + "\n", + "\n", + "data = data1 + data2 + data3\n", + "\n", + "steps = []\n", + "for i in range(len(TR_range)):\n", + " step = dict(\n", + " method = 'restyle', \n", + " args = ['visible', [False] * len(data1)],\n", + " label = str(TR_range[i])\n", + " )\n", + " step['args'][1][i] = True # Toggle i'th trace to \"visible\"\n", + " steps.append(step)\n", + "\n", + "sliders = [dict(\n", + " x = 0,\n", + " y = -0.02,\n", + " active = 10,\n", + " currentvalue = {\"prefix\": \"TR value (ms): \"},\n", + " pad = {\"t\": 50, \"b\": 10},\n", + " steps = steps\n", + ")]\n", + "\n", + "layout = go.Layout(\n", + " width=580,\n", + " height=450,\n", + " margin=go.layout.Margin(\n", + " l=80,\n", + " r=40,\n", + " b=60,\n", + " t=10,\n", + " ),\n", + " annotations=[\n", + " dict(\n", + " x=0.5004254919715793,\n", + " y=-0.18,\n", + " showarrow=False,\n", + " text='Inversion Time – TI (ms)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=-0.14,\n", + " y=0.5,\n", + " showarrow=False,\n", + " text='Signal (magnitude)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " ],\n", + " xaxis=dict(\n", + " autorange=False,\n", + " range=[0, params['TI'][-1]],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " yaxis=dict(\n", + " autorange=False,\n", + " range=[0, 1],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " legend=dict(\n", + " x=0.2,\n", + " y=0.9,\n", + " traceorder='normal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=12,\n", + " color='#000'\n", + " ),\n", + " bordercolor='#000000',\n", + " borderwidth=2\n", + " ), \n", + " sliders=sliders,\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'ir_fig_4.html', config = config)\n", + "display(HTML('ir_fig_4.html'))" + ] + }, + { + "cell_type": "markdown", + "id": "d2c9857d", + "metadata": {}, + "source": [ + "```{figure} img/plot3.png\n", + ":label: irPlot3\n", + "\n", + "Fitting comparison of simulated data (blue markers) with T1 = 1 s and TR = 1.5 to 5 s, using fitted using RD-NLS & Eq. 4 (green) and Levenberg-Marquardt & Eq. 2 (orange, long TR approximation).\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "d3accefe-ab3e-490f-9d97-db2f1b2b70b7", + "metadata": {}, + "source": [ + "

\n", + "Figure 5 displays an example brain dataset from an inversion recovery experiment, along with the T1 map fitted using the RD-NLS technique.\n", + "

" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "f2d2921d-e5a3-4254-a3e3-3fa9c8a70496", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"01\",'figure_5.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " T1_map, TI_0030, TI_0530, TI_1030, TI_1530, xAxis, yAxis = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "\n", + "init_notebook_mode(connected=True)\n", + "\n", + "trace1 = go.Heatmap(x = xAxis,\n", + " y = yAxis,\n", + " z=TI_0030,\n", + " colorscale='gray',\n", + " showscale = False,\n", + " visible=False,\n", + " name = 'Signal')\n", + "trace2 = go.Heatmap(x = xAxis,\n", + " y = yAxis,\n", + " z=TI_0530,\n", + " colorscale='gray',\n", + " showscale = False,\n", + " visible=False,\n", + " name = 'Signal')\n", + "trace3 = go.Heatmap(x = xAxis,\n", + " y = yAxis,\n", + " z=TI_1030,\n", + " colorscale='gray',\n", + " showscale = False,\n", + " visible=True,\n", + " name = 'Signal')\n", + "trace4 = go.Heatmap(x = xAxis,\n", + " y = yAxis,\n", + " z=TI_1530,\n", + " colorscale='gray',\n", + " visible=False,\n", + " showscale = False,\n", + " name = 'Signal')\n", + "trace5 = go.Heatmap(x = xAxis,\n", + " y = yAxis,\n", + " z=T1_map,\n", + " colorscale='Portland',\n", + " xaxis='x2',\n", + " yaxis='y2',\n", + " visible=True,\n", + " name = 'T1 values (ms)')\n", + "\n", + "data=[trace1, trace2, trace3, trace4, trace5]\n", + "\n", + "\n", + "updatemenus = list([\n", + " dict(active=2,\n", + " x = 0.12,\n", + " xanchor = 'left',\n", + " y = -0.15,\n", + " yanchor = 'bottom',\n", + " direction = 'up',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=16\n", + " ),\n", + " buttons=list([ \n", + " dict(label = '30 ms',\n", + " method = 'update',\n", + " args = [{'visible': [True, False, False, False, True]},\n", + " ]),\n", + " dict(label = '530 ms',\n", + " method = 'update',\n", + " args = [{'visible': [False, True, False, False, True]},\n", + " ]),\n", + " dict(label = '1030 ms',\n", + " method = 'update',\n", + " args = [{'visible': [False, False, True, False, True]},\n", + " ]),\n", + " dict(label = '1530 ms',\n", + " method = 'update',\n", + " args = [{'visible': [False,False, False, True, True]},\n", + " ])\n", + " ]),\n", + " )\n", + "])\n", + "\n", + "layout = dict(\n", + " width=560,\n", + " height=345,\n", + " margin = dict(\n", + " t=40,\n", + " r=50,\n", + " b=10,\n", + " l=50),\n", + " annotations=[\n", + " dict(\n", + " x=0.06,\n", + " y=1.15,\n", + " showarrow=False,\n", + " text='MR Image',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=26\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=0.6,\n", + " y=1.15,\n", + " showarrow=False,\n", + " text='T1 map',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=26\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=1.22,\n", + " y=1.15,\n", + " showarrow=False,\n", + " text='T1 (ms)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=26\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=0.02,\n", + " y=-0.15,\n", + " showarrow=False,\n", + " text='TI:',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " ],\n", + " xaxis = dict(range = [0,127], autorange = False,\n", + " showgrid = False, zeroline = False, showticklabels = False,\n", + " ticks = '', domain=[0, 0.6]),\n", + " yaxis = dict(range = [0,127], autorange = False,\n", + " showgrid = False, zeroline = False, showticklabels = False,\n", + " ticks = '', domain=[0, 1]),\n", + " xaxis2 = dict(range = [0,127], autorange = False,\n", + " showgrid = False, zeroline = False, showticklabels = False,\n", + " ticks = '', domain=[0.4, 1]),\n", + " yaxis2 = dict(range = [0,127], autorange = False,\n", + " showgrid = False, zeroline = False, showticklabels = False,\n", + " ticks = '', domain=[0, 1], anchor='x2'),\n", + " showlegend = False,\n", + " autosize = False,\n", + " updatemenus=updatemenus,\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'ir_fig_5.html', config = config)\n", + "display(HTML('ir_fig_5.html'))" + ] + }, + { + "cell_type": "markdown", + "id": "a87ab22c", + "metadata": {}, + "source": [ + "```{figure} img/plot4.png\n", + ":label: irPlot4\n", + "\n", + "Example inversion recovery dataset of a healthy adult brain (left). Inversion times used to acquire this magnitude image dataset were 30 ms, 530 ms, 1030 ms, and 1530 ms, and the TR used was 1550 ms. The T1 map (right) was fitted using a RD-NLS algorithm.\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "023c58aa", + "metadata": {}, + "source": [ + "\n", + "```{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 4.\n", + ":class: tip, dropdown\n", + "\n", + "\n", + "```octave\n", + "\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "clear all\n", + "\n", + "%% Setup parameters\n", + "% All times are in milliseconds\n", + "% All flip angles are in degrees\n", + "\n", + "params.TI = 50:50:1500;\n", + "TR_range = 1500:50:5000;\n", + "\n", + "params.EXC_FA = 90;\n", + "params.INV_FA = 180;\n", + "\n", + "params.T1 = 1000;\n", + "\n", + "%% Calculate signals\n", + "%\n", + "% The option 'GRE-IR' selects the analytical equations for the gradient echo readout inversion recovery experiment\n", + "% The option '1' is a flag that selects full analytical solution equation (no approximation), Eq. 1 of the blog post.\n", + "%\n", + "% To see all the options available, run `help inversion_recovery.analytical_solution`\n", + "\n", + "for ii = 1:length(TR_range)\n", + " params.TR = TR_range(ii);\n", + " Mz_analytical(ii,:) = inversion_recovery.analytical_solution(params, 'GRE-IR', 1);\n", + "end\n", + "\n", + "%% Fit data using Levenberg-Marquardt with the long TR approximation equation\n", + "%\n", + "% The option '4' is a flag that selects the long TR approximation of the analytical solution (TR>>T1), Eq. 3 of the blog post.\n", + "%\n", + "% To see all the options available, run `help inversion_recovery.fit_lm`\n", + "\n", + "\n", + "for ii=1:length(TR_range)\n", + " fitOutput_lm{ii} = inversion_recovery.fit_lm(Mz_analytical(ii,:), params, 4);\n", + " T1_lm(ii) = fitOutput_lm{ii}.T1;\n", + "end\n", + "\n", + "%% Fit data using the RDLS method (Barral), Eq. 4 of the blog post.\n", + "%\n", + "\n", + "% Create a qMRLab inversion recovery model object and load protocol values\n", + "irObj = inversion_recovery();\n", + "irObj.Prot.IRData.Mat = params.TI';\n", + "\n", + "for ii=1:length(TR_range)\n", + "\n", + " data.IRData = Mz_analytical(ii,:);\n", + "\n", + " fitOutput_barral{ii} = irObj.fit(data);\n", + "\n", + " T1_barral(ii) = fitOutput_barral{ii}.T1;\n", + "\n", + "end\n", + "```\n", + "\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "984c4b97", + "metadata": {}, + "source": [ + "\n", + "```{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 5.\n", + ":class: tip, dropdown\n", + "\n", + "```octave\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "\n", + "clear all\n", + "\n", + "%% MATLAB/OCTAVE CODE\n", + "\n", + "% Load data into environment, and rotate mask to be aligned with IR data\n", + "load('data/ir_dataset/IRData.mat');\n", + "load('data/ir_dataset/IRMask.mat');\n", + "\n", + "IRData = data;\n", + "Mask = imrotate(Mask,180);\n", + "clear data\n", + "\n", + "% Format qMRLab inversion_recovery model parameters, and load them into the Model object\n", + "Model = inversion_recovery; \n", + "TI = [30; 530; 1030; 1530];\n", + "Model.Prot.IRData.Mat = [TI];\n", + "\n", + "% Format data structure so that they may be fit by the model\n", + "data = struct();\n", + "data.IRData= double(IRData);\n", + "data.Mask= double(Mask);\n", + "\n", + "FitResults = FitData(data,Model,0); % The '0' flag is so that no wait bar is shown.\n", + "\n", + "% Code used to re-orient the images to make pretty figures, and to assign variables with the axis lengths.\n", + "\n", + "T1_map = imrotate(FitResults.T1.*Mask,-90);\n", + "xAxis = [0:size(T1_map,2)-1];\n", + "yAxis = [0:size(T1_map,1)-1];\n", + "\n", + "% Raw MRI data at different TI values\n", + "TI_0030 = imrotate(squeeze(IRData(:,:,:,1).*Mask),-90);\n", + "TI_0530 = imrotate(squeeze(IRData(:,:,:,2).*Mask),-90);\n", + "TI_1030 = imrotate(squeeze(IRData(:,:,:,3).*Mask),-90);\n", + "TI_1530 = imrotate(squeeze(IRData(:,:,:,4).*Mask),-90);\n", + "```\n", + "\n", + "```" + ] + } + ], + "metadata": { + "celltoolbar": "Tags", + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.18" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/04-IR_BenefitsAndPitfalls.ipynb b/T1 Mapping/Inversion Recovery T1 Mapping/04-IR_BenefitsAndPitfalls.ipynb new file mode 100644 index 0000000..e194d2e --- /dev/null +++ b/T1 Mapping/Inversion Recovery T1 Mapping/04-IR_BenefitsAndPitfalls.ipynb @@ -0,0 +1,472 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "537b8c4a-fbcc-45cb-804c-4a0d299750be", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "---\n", + "title: Benefits and Pitfalls\n", + "subtitle: Inversion Recovery\n", + "date: 2024-07-25\n", + "authors:\n", + " - name: Mathieu Boudreau\n", + " affiliations:\n", + " - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada\n", + "numbering:\n", + " figure:\n", + " template: Fig. %s\n", + "---" + ] + }, + { + "cell_type": "markdown", + "id": "f8d1db4b-95c2-4a61-9217-1be991784701", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "The conventional inversion recovery experiment is considered the gold standard T1 mapping technique for several reasons: \n", + "* A typical protocol has a long TR value and a sufficient number of inversion times for stable fitting (typically 5 or more) covering the range [0, TR]. \n", + "* It offers a wide dynamic range of signals ([up to `-kM0`, `kM0`]), allowing a number of inversion times where high SNR is available to sample the signal recovery curve {cite:p}`Fukushima1981`. \n", + "* T1 maps produced by inversion recovery are largely insensitive to inaccuracies in excitation flip angles and imperfect spoiling {cite:p}`Stikov2015`, as all parameters except TI are constant for each measurement and only a single acquisition is performed (at TI) during each TR. \n", + "\n", + "One important protocol design consideration is to avoid acquiring at inversion times where the signal for T1 values of the tissue-of-interest is nulled, as the magnitude images at this TI time will be dominated by Rician noise which can negatively impact the fit under low SNR circumstances (Figure 6). Inversion recovery can also often be acquired using commonly available standard pulse sequences available on most MRI scanners by setting up a customized acquisition protocol, and does not require any additional calibration measurements. For an example, please visit the interactive preprint of the ISMRM Reproducible Research Group 2020 Challenge on inversion recovery T1 mapping {cite:p}`Boudreau2023`. " + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "4130dad2-b8dd-43b6-96da-5a457cf29ea8", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-output", + "hide-input" + ] + }, + "outputs": [], + "source": [ + "from repo2data.repo2data import Repo2Data\n", + "import os \n", + "import pickle\n", + "import matplotlib.pyplot as plt\n", + "import chart_studio.plotly as py\n", + "import plotly.graph_objs as go\n", + "import numpy as np\n", + "from plotly import __version__\n", + "from plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot\n", + "from IPython.display import display, HTML\n", + "\n", + "data_req_path = os.path.join(\"..\",\"..\", \"binder\", \"data_requirement.json\")\n", + "repo2data = Repo2Data(data_req_path)\n", + "DATA_ROOT = os.path.join(repo2data.install()[0],\"t1-book-neurolibre\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "369dd138-a820-43e3-b0fe-fdc968ca41e6", + "metadata": { + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"01\",'figure_6.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " TR_range, TI_lowres, TI_highres, T1_mean, T1_std, data_mean, data_std, data_noiseless = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "\n", + "init_notebook_mode(connected=True)\n", + "\n", + "data1 = [dict(\n", + " visible = False,\n", + " x = np.squeeze(np.asarray(TI_lowres[ii,:])),\n", + " y = np.squeeze(np.asarray(data_mean[ii,:])),\n", + " error_y=dict(\n", + " type='data',\n", + " color = ('rgb(22, 96, 167)'),\n", + " array=np.squeeze(np.asarray(data_std[ii,:])),\n", + " visible=True\n", + " ),\n", + " line = dict(\n", + " color = ('rgb(22, 96, 167)'),\n", + " dash = 'dot'),\n", + " mode = 'markers',\n", + " name = 'Monte Carlo simulated signal',\n", + " text = 'Monte Carlo simulated signal',\n", + " hoverinfo = 'x+y+text') for ii in range(len(TR_range))]\n", + "\n", + "data1[28]['visible'] = True\n", + "\n", + "data2 = [dict(\n", + " visible = False,\n", + " x = np.squeeze(np.asarray(TI_highres[ii,:])),\n", + " y = np.squeeze(np.asarray(data_noiseless[ii,:])),\n", + " line = dict(\n", + " color = ('rgb(247, 152, 19)'),\n", + " ),\n", + " name = 'Noiseless signal',\n", + " text = 'Noiseless signal',\n", + " hoverinfo = 'x+y+text') for ii in range(len(TR_range))]\n", + "\n", + "data2[28]['visible'] = True\n", + "\n", + "data_meanT1 = [dict(\n", + " visible = False,\n", + " x = TR_range,\n", + " y = T1_mean,\n", + " name = 'Mean T1 (s)',\n", + " text = 'Mean T1 (s)',\n", + " hoverinfo = 'x+y+text',\n", + " xaxis='x2',\n", + " yaxis='y2') for ii in range(len(TR_range))]\n", + "\n", + "data_meanT1[15]['visible'] = True\n", + "\n", + "data_stdT1 = [dict(\n", + " visible = False,\n", + " x = TR_range,\n", + " y = T1_std,\n", + " line = dict(\n", + " color = ('rgb(222, 22, 22)'),\n", + " ),\n", + " name = 'STD T1 (s)',\n", + " text = 'STD T1 (s)',\n", + " hoverinfo = 'x+y+text',\n", + " xaxis='x2',\n", + " yaxis='y3') for ii in range(len(TR_range))]\n", + "\n", + "data_stdT1[28]['visible'] = True\n", + "\n", + "data = data2 + data1 + data_meanT1 + data_stdT1\n", + "\n", + "steps = []\n", + "for i in range(len(TR_range)):\n", + " step = dict(\n", + " method = 'restyle', \n", + " args = ['visible', [False] * len(data1)],\n", + " label = str(TR_range[i])\n", + " )\n", + " step['args'][1][i] = True # Toggle i'th trace to \"visible\"\n", + " steps.append(step)\n", + "\n", + "sliders = [dict(\n", + " x = 0,\n", + " y = -0.02,\n", + " active = 28,\n", + " currentvalue = {\"prefix\": \"TR value (ms): \"},\n", + " pad = {\"t\": 50, \"b\": 10},\n", + " steps = steps\n", + ")]\n", + "\n", + "layout = go.Layout(\n", + " width=540,\n", + " height=540,\n", + " margin = dict(\n", + " t=0,\n", + " r=25,\n", + " b=100,\n", + " l=75),\n", + " annotations=[\n", + " dict(\n", + " x=0.5004254919715793,\n", + " y=-0.17,\n", + " showarrow=False,\n", + " text='Inversion Time – TI (ms)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=26\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=-0.15,\n", + " y=0.5,\n", + " showarrow=False,\n", + " text='Signal (magnitude)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=26\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=0.76,\n", + " y=0.77,\n", + " showarrow=False,\n", + " text='TR (ms)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=14\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=0.40,\n", + " y=0.35,\n", + " showarrow=False,\n", + " text='Mean T1 (ms)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=14\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=1.00,\n", + " y=0.35,\n", + " showarrow=False,\n", + " text='STD T1 (ms)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=14\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " )\n", + " ],\n", + " xaxis=dict(\n", + " autorange=False,\n", + " range=[0, 5000],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " yaxis=dict(\n", + " autorange=False,\n", + " range=[0, 1],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " xaxis2=dict(\n", + " domain=[0.5, 0.90],\n", + " anchor='y2',\n", + " mirror = True,\n", + " side='top',\n", + " ticks='inside',\n", + " showline=True,\n", + " ),\n", + " yaxis2=dict(\n", + " autorange=False,\n", + " range=[500, 1300],\n", + " domain=[0.05, 0.65],\n", + " anchor='x2',\n", + " mirror = True,\n", + " ticks='inside',\n", + " showline=True,\n", + " ),\n", + " yaxis3=dict(\n", + " autorange=False,\n", + " range=[0, 190],\n", + " domain=[0.05, 0.65],\n", + " anchor='x2',\n", + " overlaying='y2',\n", + " side='right',\n", + " ticks='inside',\n", + " ),\n", + " legend=dict(\n", + " x=0.3,\n", + " y=1.35,\n", + " traceorder='normal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=12,\n", + " color='#000'\n", + " ),\n", + " bordercolor='#000000',\n", + " borderwidth=2\n", + " ), \n", + " sliders=sliders,\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'ir_fig_6.html', config = config)\n", + "display(HTML('ir_fig_6.html'))" + ] + }, + { + "cell_type": "markdown", + "id": "44739c40-81bd-412b-ab14-2cf46ba03d1d", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [] + }, + "source": [ + "\n", + "```{figure} img/plot5.png\n", + ":label: irPlot5\n", + "\n", + "Figure 6. Monte Carlo simulations (mean and standard deviation (STD), blue markers) and fitted T1 values (mean and STD, red and green respectively) generated for a T1 value of 900 ms and 5 TI values linearly spaced across the TR (ranging from 1 to 5 s). A bump in T1 STD occurs near TR = 3000 ms, which coincides with the TR where the second TI is located near a null point for this T1 value.\n", + "```\n", + "\n", + "Despite a widely acknowledged robustness for measuring accurate T1 maps, inversion recovery is not often used in studies. An important drawback of this technique is the need for long TR values, generally on the order of a few T1 for general models (e.g. Eqs. 1 and 4), and up to 5T1 for long TR approximated models (Eq. 3). It takes about to 10-25 minutes to acquire a single-slice T1 map using the inversion recovery technique, as only one TI is acquired per TR (2-5 s) and conventional cartesian gradient readout imaging acquires one phase encode line per excitation (for a total of ~100-200 phase encode lines). The long acquisition time makes it challenging to acquire whole-organ T1 maps in clinically feasible protocol times. Nonetheless, it is useful as a reference measurement for comparisons against other T1 mapping methods, or to acquire a single-slice T1 map of a tissue to get T1 estimates for optimization of other pulse sequences." + ] + }, + { + "cell_type": "markdown", + "id": "7bf30faf-3f00-4c6e-9eed-c6b880307009", + "metadata": {}, + "source": [ + "\n", + "```{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 6.\n", + ":class: tip, dropdown\n", + "\n", + "```octave\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "clear all\n", + "\n", + "%% Setup parameters\n", + "% All times are in seconds\n", + "% All flip angles are in degrees\n", + "\n", + "TR_range = 1000:100:5000; % in seconds\n", + "\n", + "x = struct;\n", + "x.T1 = 900; % in seconds\n", + "\n", + "Opt.SNR = 25;\n", + "Opt.M0 = 1;\n", + "Opt.FAexcite = 90; % Excitation flip angle\n", + "Opt.FAinv = 180; % Inversion flip angle\n", + "\n", + "%% Monte Carlo data simulation\n", + "% Simulate noisy signal data 1,000 time, fit the data, then calculate the means and standard deviations of the data and fitted T1\n", + "% Data is calculated by calculating the a and b values of Eq. 4 from the full analytical equations (Eq. 1)\n", + "\n", + "Model = inversion_recovery; \n", + "\n", + "for ii = 1:length(TR_range)\n", + " Opt.TR = TR_range(ii);\n", + " Opt.T1 = x.T1;\n", + " TI_lowres(ii,:) = linspace(0.05, Opt.TR, 6)';\n", + " Model.Prot.IRData.Mat = [TI_lowres(ii,:)];\n", + " [ra,rb] = Model.ComputeRaRb(x,Opt);\n", + " x.rb = rb;\n", + " x.ra = ra;\n", + " for jj = 1:1000\n", + " [FitResult{ii,jj}, noisyData{ii,jj}] = Model.Sim_Single_Voxel_Curve(x,Opt,0); \n", + " fittedT1(ii,jj) = FitResult{ii,jj}.T1;\n", + " noisyData_array(ii,jj,:) = noisyData{ii,jj}.IRData;\n", + " end\n", + " \n", + " for kk=1:length(TI_lowres(ii,:))\n", + " data_mean(ii,kk) = mean(noisyData_array(ii,:,kk));\n", + " data_std(ii,kk) = std(noisyData_array(ii,:,kk));\n", + " end\n", + " \n", + " T1_mean(ii) = mean(fittedT1(ii,:));\n", + " T1_std(ii) = std(fittedT1(ii,:));\n", + "end\n", + "\n", + "%% Calculate the noiseless data at a higher TI resolution to plot the ideal signal curve.\n", + "%\n", + "\n", + "for ii = 1:length(TR_range)\n", + " TI_highres(ii,:) = linspace(0.05, TR_range(ii), 500);\n", + " Model.Prot.IRData.Mat = [TI_highres(ii,:)];\n", + " Opt.TR = TR_range(ii);\n", + " Opt.T1 = x.T1;\n", + " [ra,rb] = Model.ComputeRaRb(x,Opt);\n", + " x.rb = rb;\n", + " x.ra = ra;\n", + "\n", + " data_noiseless(ii,:) = Model.equation(x);\n", + "end\n", + "```\n", + "\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "7c77eab6-945c-4086-a932-1778592a606c", + "metadata": {}, + "source": [ + "```{admonition} References\n", + ":class: seealso\n", + "\n", + "```{bibliography}\n", + ":filter: docname in docnames\n", + "```\n", + "\n", + "```" + ] + } + ], + "metadata": { + "celltoolbar": "Tags", + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.18" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/IR_OtherMethods.md b/T1 Mapping/Inversion Recovery T1 Mapping/IR_OtherMethods.md new file mode 100644 index 0000000..8a0e1c5 --- /dev/null +++ b/T1 Mapping/Inversion Recovery T1 Mapping/IR_OtherMethods.md @@ -0,0 +1,16 @@ +--- +title: Other Saturation-Recovery T1 Mapping techniques +subtitle: Inversion Recovery +date: 2024-07-25 +authors: + - name: Mathieu Boudreau + affiliations: + - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada +numbering: + figure: + template: Fig. %s +--- + +Several variations of the inversion recovery pulse sequence were developed to overcome challenges like those specified above. Amongst them, the Look-Locker technique {cite:p}`Look1970` stands out as one of the most widely used in practice. Instead of a single 90° acquisition per TR, a periodic train of small excitation pulses θ are applied after the inversion pulse, {θ180 – 𝛕 – θ – 𝛕 – θ – ...}, where 𝛕 = TR/n and n is the number of sampling acquisitions. This pulse sequence samples the inversion time relaxation curve much more efficiently than conventional inversion recovery, but at a cost of lower SNR. However, because the magnetization state of each TI measurement depends on the previous series of θ excitation, it has higher sensitivity to B1-inhomogeneities and imperfect spoiling compared to inversion recovery {cite:p}`Gai2013,Stikov2015`. Nonetheless, Look-Locker is widely used for rapid T1 mapping applications, and variants like MOLLI (Modified Look-Locker Inversion recovery) and ShMOLLI (Shortened MOLLI) are widely used for cardiac T1 mapping {cite:p}`Messroghli2004,Piechnik2010`. + +Another inversion recovery variant that’s worth mentioning is saturation recovery, in which the inversion pulse is replaced with a saturation pulse: {θ90 – TI – θ90}. This technique was used to acquire the very first T1 map {cite:p}`Pykett1978`. Unlike inversion recovery, this pulse sequence does not need a long TR to recover to its initial condition; every θ90 pulse resets the longitudinal magnetization to the same initial state. However, to properly sample the recovery curve, TIs still need to reach the order of ~T1, the dynamic range of signal potential is cut in half ([0, M0]), and the short TIs (which have the fastest acquisition times) have the lowest SNRs. \ No newline at end of file diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/equation1.png b/T1 Mapping/Inversion Recovery T1 Mapping/img/equation1.png new file mode 100644 index 0000000..15b1f42 Binary files /dev/null and b/T1 Mapping/Inversion Recovery T1 Mapping/img/equation1.png differ diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/equation2.png b/T1 Mapping/Inversion Recovery T1 Mapping/img/equation2.png new file mode 100644 index 0000000..85b107c Binary files /dev/null and b/T1 Mapping/Inversion Recovery T1 Mapping/img/equation2.png differ diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/equation3.png b/T1 Mapping/Inversion Recovery T1 Mapping/img/equation3.png new file mode 100644 index 0000000..d421ce1 Binary files /dev/null and b/T1 Mapping/Inversion Recovery T1 Mapping/img/equation3.png differ diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/equation4.png b/T1 Mapping/Inversion Recovery T1 Mapping/img/equation4.png new file mode 100644 index 0000000..4e21645 Binary files /dev/null and b/T1 Mapping/Inversion Recovery T1 Mapping/img/equation4.png differ diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/equations_all.svg b/T1 Mapping/Inversion Recovery T1 Mapping/img/equations_all.svg new file mode 100644 index 0000000..9361a87 --- /dev/null +++ b/T1 Mapping/Inversion Recovery T1 Mapping/img/equations_all.svg @@ -0,0 +1,700 @@ + + + + + + + + + + image/svg+xml + + + + + + + (1) + + + + + + + + + diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/ir_pulsesequences.png b/T1 Mapping/Inversion Recovery T1 Mapping/img/ir_pulsesequences.png new file mode 100644 index 0000000..718c54b Binary files /dev/null and b/T1 Mapping/Inversion Recovery T1 Mapping/img/ir_pulsesequences.png differ diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/ir_pulsesequences.svg b/T1 Mapping/Inversion Recovery T1 Mapping/img/ir_pulsesequences.svg new file mode 100644 index 0000000..557c8c3 --- /dev/null +++ b/T1 Mapping/Inversion Recovery T1 Mapping/img/ir_pulsesequences.svg @@ -0,0 +1,463 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + image/svg+xml + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + IMG + + + TIn + TR + + 180° + 90° + + 180° + + + diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/plot1.png b/T1 Mapping/Inversion Recovery T1 Mapping/img/plot1.png new file mode 100644 index 0000000..d17443e Binary files /dev/null and b/T1 Mapping/Inversion Recovery T1 Mapping/img/plot1.png differ diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/plot2.png b/T1 Mapping/Inversion Recovery T1 Mapping/img/plot2.png new file mode 100644 index 0000000..e482205 Binary files /dev/null and b/T1 Mapping/Inversion Recovery T1 Mapping/img/plot2.png differ diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/plot3.png b/T1 Mapping/Inversion Recovery T1 Mapping/img/plot3.png new file mode 100644 index 0000000..1e799d0 Binary files /dev/null and b/T1 Mapping/Inversion Recovery T1 Mapping/img/plot3.png differ diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/plot4.png b/T1 Mapping/Inversion Recovery T1 Mapping/img/plot4.png new file mode 100644 index 0000000..b7de956 Binary files /dev/null and b/T1 Mapping/Inversion Recovery T1 Mapping/img/plot4.png differ diff --git a/T1 Mapping/Inversion Recovery T1 Mapping/img/plot5.png b/T1 Mapping/Inversion Recovery T1 Mapping/img/plot5.png new file mode 100644 index 0000000..7e1ff9f Binary files /dev/null and b/T1 Mapping/Inversion Recovery T1 Mapping/img/plot5.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/01-Introduction.md b/T1 Mapping/Variable Flip Angle T1 Mapping/01-Introduction.md new file mode 100644 index 0000000..026866f --- /dev/null +++ b/T1 Mapping/Variable Flip Angle T1 Mapping/01-Introduction.md @@ -0,0 +1,25 @@ +--- +title: Abstract +subtitle: Variable Flip Angle +date: 2024-07-25 +authors: + - name: Mathieu Boudreau + affiliations: + - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada +numbering: + heading_2: false + figure: + template: Fig. %s +--- + +## Variable Flip Angle T1 Mapping + +Variable flip angle (VFA) T1 mapping {cite:p}`Christensen1974,Fram1987,Gupta1977`, also known as Driven Equilibrium Single Pulse Observation of T1 (DESPOT1) {cite:p}`Homer1985,Deoni2003`, 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 {cite:p}`Christensen1974,Gupta1977`. 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 {cite:p}`Sled1998`. 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 {cite:p}`Boudreau2017,Lee2017`. + +```{figure} img/vfa_pulsesequence.png +:label: vfaFig1 + +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. +``` \ No newline at end of file diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/02-VFA_SignalModelling.ipynb b/T1 Mapping/Variable Flip Angle T1 Mapping/02-VFA_SignalModelling.ipynb new file mode 100644 index 0000000..a5c3c14 --- /dev/null +++ b/T1 Mapping/Variable Flip Angle T1 Mapping/02-VFA_SignalModelling.ipynb @@ -0,0 +1,789 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "a1532bfe-fb4c-4258-a2c6-869b645579e8", + "metadata": {}, + "source": [ + "---\n", + "title: Signal Modelling\n", + "subtitle: Variable Flip Angle\n", + "date: 2024-07-25\n", + "authors:\n", + " - name: Mathieu Boudreau\n", + " affiliations:\n", + " - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada\n", + "numbering:\n", + " figure:\n", + " template: Fig. %s\n", + "---\n", + "\n", + "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}:\n", + "\n", + "\\begin{equation}\\tag{1}\n", + "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)\n", + "\\end{equation}\n", + "\n", + "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.\n" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "c22690d8-66dd-401f-86be-e3f3968f828d", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-output", + "hide-input" + ] + }, + "outputs": [], + "source": [ + "from repo2data.repo2data import Repo2Data\n", + "import os \n", + "import pickle\n", + "import matplotlib.pyplot as plt\n", + "import chart_studio.plotly as py\n", + "import plotly.graph_objs as go\n", + "import numpy as np\n", + "from plotly import __version__\n", + "from plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot\n", + "from IPython.display import display, HTML\n", + "\n", + "data_req_path = os.path.join(\"..\",\"..\", \"binder\", \"data_requirement.json\")\n", + "repo2data = Repo2Data(data_req_path)\n", + "DATA_ROOT = os.path.join(repo2data.install()[0],\"t1-book-neurolibre\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "6ad4c416-1995-48d5-a9dd-a248b7172585", + "metadata": { + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"02\",'figure_2.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " params, TR_range, signal_WM, signal_GM, signal_CSF = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "\n", + "init_notebook_mode(connected=True)\n", + "\n", + "data1 = [dict(\n", + " visible = False,\n", + " mode = 'lines',\n", + " x = params[\"EXC_FA\"],\n", + " y = abs(np.squeeze(np.asarray(signal_WM[ii]))),\n", + " name = 'T1 = 0.9 s (White Matter)',\n", + " text = 'T1 = 0.9 s (White Matter)',\n", + " hoverinfo = 'x+y+text') for ii in range(len(TR_range))]\n", + "\n", + "data1[4]['visible'] = True\n", + "\n", + "data2 = [dict(\n", + " visible = False,\n", + " mode = 'lines',\n", + " x = params[\"EXC_FA\"],\n", + " y = abs(np.squeeze(np.asarray(signal_GM[ii]))),\n", + " name = 'T1 = 1.5 s (Grey Matter)',\n", + " text = 'T1 = 1.5 s (Grey Matter)',\n", + " hoverinfo = 'x+y+text') for ii in range(len(TR_range))]\n", + "\n", + "data2[4]['visible'] = True\n", + "\n", + "data3 = [dict(\n", + " visible = False,\n", + " mode = 'lines',\n", + " x = params[\"EXC_FA\"],\n", + " y = abs(np.squeeze(np.asarray(signal_CSF[ii]))),\n", + " name = 'T1 = 4.0 s (Cerebrospinal Fluid)',\n", + " text = 'T1 = 4.0 s (Cerebrospinal Fluid)',\n", + " hoverinfo = 'x+y+text') for ii in range(len(TR_range))]\n", + "\n", + "data3[4]['visible'] = True\n", + "\n", + "data = data1 + data2 + data3\n", + "\n", + "steps = []\n", + "for i in range(len(TR_range)):\n", + " step = dict(\n", + " method = 'restyle', \n", + " args = ['visible', [False] * len(data1)],\n", + " label = str(TR_range[i])\n", + " )\n", + " step['args'][1][i] = True # Toggle i'th trace to \"visible\"\n", + " steps.append(step)\n", + "\n", + "sliders = [dict(\n", + " x = 0,\n", + " y = -0.02,\n", + " active = 2,\n", + " currentvalue = {\"prefix\": \"TR value (ms): \"},\n", + " pad = {\"t\": 50, \"b\": 10},\n", + " steps = steps\n", + ")]\n", + "\n", + "layout = go.Layout(\n", + " width=580,\n", + " height=450,\n", + " margin=go.layout.Margin(\n", + " l=80,\n", + " r=40,\n", + " b=60,\n", + " t=10,\n", + " ),\n", + " annotations=[\n", + " dict(\n", + " x=0.5004254919715793,\n", + " y=-0.18,\n", + " showarrow=False,\n", + " text='Excitation Flip Angle (°)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=-0.15,\n", + " y=0.5,\n", + " showarrow=False,\n", + " text='Long. Magnetization (Mz)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " ],\n", + " xaxis=dict(\n", + " autorange=False,\n", + " range=[0, params['EXC_FA'][-1]],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " yaxis=dict(\n", + " autorange=True,\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " legend=dict(\n", + " x=0.5,\n", + " y=0.9,\n", + " traceorder='normal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=12,\n", + " color='#000'\n", + " ),\n", + " bordercolor='#000000',\n", + " borderwidth=2\n", + " ), \n", + " sliders=sliders,\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'vfa_fig_2.html', config = config)\n", + "display(HTML('vfa_fig_2.html'))" + ] + }, + { + "cell_type": "markdown", + "id": "5f6de3d6-fbf2-4e50-904f-693cbbed6dcb", + "metadata": {}, + "source": [ + "```{figure} img/plot1.png\n", + ":label: vfaPlot1\n", + "\n", + "Variable flip angle technique signal curves (Eq. 1) for three different T1 values, approximating the main types of tissue in the brain at 3T.\n", + "```\n", + "\n", + "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 {cite:p}`Ernst1966`, which can be solved analytically from Equation 1 using properties of calculus:\n", + "\n", + "\n", + "\\begin{equation}\\tag{2}\n", + "\\theta_{Ernst} = \\text{acos}(e^{- \\frac{TR}{T_1}})\n", + "\\end{equation}\n", + "\n", + "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." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "91fc5cf8-1e76-4629-95b9-63723e2f80e4", + "metadata": { + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"02\",'figure_3.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " params, Nex_range, signal_analytical, signal_blochsim = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "\n", + "init_notebook_mode(connected=True)\n", + "\n", + "data1 = [dict(\n", + " visible = False,\n", + " mode = 'lines',\n", + " x = params[\"EXC_FA\"],\n", + " y = abs(np.squeeze(np.asarray(signal_analytical[ii]))),\n", + " name = 'Analytical Solution',\n", + " text = 'Analytical Solution',\n", + " hoverinfo = 'x+y+text') for ii in range(len(Nex_range))]\n", + "\n", + "data1[9]['visible'] = True\n", + "\n", + "data2 = [dict(\n", + " visible = False,\n", + " mode = 'lines',\n", + " x = params[\"EXC_FA\"],\n", + " y = abs(np.squeeze(np.asarray(signal_blochsim[ii]))),\n", + " name = 'Bloch Simulation',\n", + " text = 'Bloch Simulation',\n", + " hoverinfo = 'x+y+text') for ii in range(len(Nex_range))]\n", + "\n", + "data2[9]['visible'] = True\n", + "\n", + "data = data1 + data2\n", + "\n", + "steps = []\n", + "for i in range(len(Nex_range)):\n", + " step = dict(\n", + " method = 'restyle', \n", + " args = ['visible', [False] * len(data1)],\n", + " label = str(Nex_range[i])\n", + " )\n", + " step['args'][1][i] = True # Toggle i'th trace to \"visible\"\n", + " steps.append(step)\n", + "\n", + "sliders = [dict(\n", + " x = 0,\n", + " y = -0.02,\n", + " active = 9,\n", + " currentvalue = {\"prefix\": \"nth TR: \"},\n", + " pad = {\"t\": 50, \"b\": 10},\n", + " steps = steps\n", + ")]\n", + "\n", + "layout = go.Layout(\n", + " width=580,\n", + " height=450,\n", + " margin=go.layout.Margin(\n", + " l=80,\n", + " r=40,\n", + " b=60,\n", + " t=10,\n", + " ),\n", + " annotations=[\n", + " dict(\n", + " x=0.5004254919715793,\n", + " y=-0.18,\n", + " showarrow=False,\n", + " text='Excitation Flip Angle (°)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=-0.15,\n", + " y=0.5,\n", + " showarrow=False,\n", + " text='Signal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " ],\n", + " xaxis=dict(\n", + " autorange=False,\n", + " range=[0, params['EXC_FA'][-1]],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " yaxis=dict(\n", + " autorange=True,\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " legend=dict(\n", + " x=0.5,\n", + " y=0.9,\n", + " traceorder='normal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=12,\n", + " color='#000'\n", + " ),\n", + " bordercolor='#000000',\n", + " borderwidth=2\n", + " ), \n", + " sliders=sliders,\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'vfa_fig_3.html', config = config)\n", + "display(HTML('vfa_fig_3.html'))" + ] + }, + { + "cell_type": "markdown", + "id": "fe879588-3b24-40a5-8160-eee5a5dff39f", + "metadata": {}, + "source": [ + "```{figure} img/plot2.png\n", + ":label: vfaPlot2\n", + "\n", + "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).\n", + "```\n", + "\n", + "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 {cite:p}`Handbook2004,Zur1991` 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 {cite:p}`Yarnykh2010`." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "55252e83-7b0c-407d-930d-be0483be3c1e", + "metadata": { + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"02\",'figure_4.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " params, Nex_range, signal_ideal_spoil, signal_optimal_crush_and_rf_spoil, signal_no_gradient_and_rf_spoil = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "init_notebook_mode(connected=True)\n", + "\n", + "data1 = [dict(\n", + " visible = False,\n", + " mode = 'lines',\n", + " x = params[\"EXC_FA\"],\n", + " y = abs(np.squeeze(np.asarray(signal_ideal_spoil[ii]))),\n", + " name = 'Ideal Spoiling',\n", + " text = 'Ideal Spoiling',\n", + " hoverinfo = 'x+y+text') for ii in range(len(Nex_range))]\n", + "\n", + "data1[10]['visible'] = True\n", + "\n", + "data2 = [dict(\n", + " visible = False,\n", + " mode = 'lines',\n", + " x = params[\"EXC_FA\"],\n", + " y = abs(np.squeeze(np.asarray(signal_optimal_crush_and_rf_spoil[ii]))),\n", + " name = 'Gradient & RF Spoiling',\n", + " text = 'Gradient & RF Spoiling',\n", + " hoverinfo = 'x+y+text') for ii in range(len(Nex_range))]\n", + "\n", + "data2[10]['visible'] = True\n", + "\n", + "data3 = [dict(\n", + " visible = False,\n", + " mode = 'lines',\n", + " x = params[\"EXC_FA\"],\n", + " y = abs(np.squeeze(np.asarray(signal_no_gradient_and_rf_spoil[ii]))),\n", + " name = 'No Spoiling',\n", + " text = 'No Spoiling',\n", + " hoverinfo = 'x+y+text') for ii in range(len(Nex_range))]\n", + "\n", + "data3[10]['visible'] = True\n", + "\n", + "data = data1 + data2+ data3\n", + "\n", + "steps = []\n", + "for i in range(len(Nex_range)):\n", + " step = dict(\n", + " method = 'restyle', \n", + " args = ['visible', [False] * len(data1)],\n", + " label = str(Nex_range[i])\n", + " )\n", + " step['args'][1][i] = True # Toggle i'th trace to \"visible\"\n", + " steps.append(step)\n", + "\n", + "sliders = [dict(\n", + " x = 0,\n", + " y = -0.02,\n", + " active = 10,\n", + " currentvalue = {\"prefix\": \"nth TR: \"},\n", + " pad = {\"t\": 50, \"b\": 10},\n", + " steps = steps\n", + ")]\n", + "\n", + "layout = go.Layout(\n", + " width=580,\n", + " height=450,\n", + " margin=go.layout.Margin(\n", + " l=80,\n", + " r=40,\n", + " b=60,\n", + " t=10,\n", + " ),\n", + " annotations=[\n", + " dict(\n", + " x=0.5004254919715793,\n", + " y=-0.18,\n", + " showarrow=False,\n", + " text='Excitation Flip Angle (°)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=-0.15,\n", + " y=0.5,\n", + " showarrow=False,\n", + " text='Signal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " ],\n", + " xaxis=dict(\n", + " autorange=False,\n", + " range=[0, params['EXC_FA'][-1]],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " yaxis=dict(\n", + " autorange=True,\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " legend=dict(\n", + " x=0.5,\n", + " y=0.9,\n", + " traceorder='normal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=12,\n", + " color='#000'\n", + " ),\n", + " bordercolor='#000000',\n", + " borderwidth=2\n", + " ), \n", + " sliders=sliders,\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'vfa_fig_4.html', config = config)\n", + "display(HTML('vfa_fig_4.html'))\n", + "\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "id": "ff310736", + "metadata": {}, + "source": [ + "```{figure} img/plot3.png\n", + ":label: vfaPlot3\n", + "\n", + "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, T1 = 900 ms, Te = 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ɸ0 = ½ ɸ0(n2 + n + 2) {cite:p}`Handbook2004` with ɸ0 = 117° {cite:p}`Zur1991` after each TR.\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "42aa7bdb", + "metadata": {}, + "source": [ + "\n", + "```{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 1.\n", + ":class: tip, dropdown\n", + "\n", + "```octave\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "clear all\n", + "\n", + "%% Setup parameters\n", + "% All times are in milliseconds\n", + "% All flip angles are in degrees\n", + "\n", + "TR_range = 5:5:200;\n", + "\n", + "params.EXC_FA = 1:90;\n", + "\n", + "%% Calculate signals\n", + "%\n", + "% To see all the options available, run `help vfa_t1.analytical_solution`\n", + "\n", + "for ii = 1:length(TR_range)\n", + " params.TR = TR_range(ii);\n", + " \n", + " % White matter\n", + " params.T1 = 900; % in milliseconds\n", + "\n", + " signal_WM(ii,:) = vfa_t1.analytical_solution(params);\n", + "\n", + " % Grey matter\n", + " params.T1 = 1500; % in milliseconds\n", + " signal_GM(ii,:) = vfa_t1.analytical_solution(params);\n", + "\n", + " % CSF\n", + " params.T1 = 4000; % in milliseconds\n", + " signal_CSF(ii,:) = vfa_t1.analytical_solution(params);\n", + "end\n", + "\n", + "```\n", + "\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "3d1b8cb3", + "metadata": {}, + "source": [ + "\n", + "```{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 3.\n", + ":class: tip, dropdown\n", + "\n", + "```octave\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "clear all\n", + "\n", + "%% Setup parameters\n", + "% All times are in milliseconds\n", + "% All flip angles are in degrees\n", + "\n", + "% White matter\n", + "params.T1 = 900; % in milliseconds\n", + "params.T2 = 10000;\n", + "params.TR = 25;\n", + "params.TE = 5;\n", + "params.EXC_FA = 1:90;\n", + "Nex_range = 1:1:150;\n", + "\n", + "%% Calculate signals\n", + "%\n", + "% To see all the options available, run `help vfa_t1.analytical_solution`\n", + "\n", + "for ii = 1:length(Nex_range)\n", + " params.Nex = Nex_range(ii);\n", + " \n", + " signal_analytical(ii,:) = vfa_t1.analytical_solution(params);\n", + "\n", + " [~, complex_signal] = vfa_t1.bloch_sim(params);\n", + " signal_blochsim(ii,:) = abs(complex(complex_signal));\n", + "end\n", + "\n", + "```\n", + "\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "17fcbae5", + "metadata": {}, + "source": [ + "\n", + "```{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 4.\n", + ":class: tip, dropdown\n", + "\n", + "```octave\n", + "\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "clear all\n", + "\n", + "%% Setup parameters\n", + "% All times are in milliseconds\n", + "% All flip angles are in degrees\n", + "\n", + "% White matter\n", + "params.T1 = 900; % in milliseconds\n", + "params.T2 = 100;\n", + "params.TR = 25;\n", + "params.TE = 5;\n", + "params.EXC_FA = 1:90;\n", + "Nex_range = [1:9, 10:10:100];\n", + "\n", + "%% Calculate signals\n", + "%\n", + "% To see all the options available, run `help vfa_t1.analytical_solution`\n", + "\n", + "for ii = 1:length(Nex_range)\n", + " params.Nex = Nex_range(ii);\n", + " \n", + " params.crushFlag = 1;\n", + " \n", + " [~, complex_signal] = vfa_t1.bloch_sim(params);\n", + " signal_ideal_spoil(ii,:) = abs(complex_signal);\n", + " \n", + " \n", + " params.inc = 117;\n", + " params.partialDephasing = 1;\n", + " params.partialDephasingFlag = 1;\n", + " params.crushFlag = 0;\n", + " \n", + " [~, complex_signal] = vfa_t1.bloch_sim(params);\n", + " signal_optimal_crush_and_rf_spoil(ii,:) = abs(complex_signal);\n", + " \n", + " params.inc = 0;\n", + " params.partialDephasing = 0;\n", + "\n", + " [~, complex_signal] = vfa_t1.bloch_sim(params);\n", + " signal_no_gradient_and_rf_spoil(ii,:) = abs(complex_signal);\n", + "end\n", + "```\n", + "\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "da5d20eb-81ba-4b4b-8345-dbf807f239e9", + "metadata": {}, + "source": [ + "```{admonition} References\n", + ":class: seealso\n", + "\n", + "```{bibliography}\n", + ":filter: docname in docnames\n", + "```\n", + "\n", + "```" + ] + } + ], + "metadata": { + "celltoolbar": "Tags", + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.18" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/03-VFA_DataFitting.ipynb b/T1 Mapping/Variable Flip Angle T1 Mapping/03-VFA_DataFitting.ipynb new file mode 100644 index 0000000..c83a139 --- /dev/null +++ b/T1 Mapping/Variable Flip Angle T1 Mapping/03-VFA_DataFitting.ipynb @@ -0,0 +1,958 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "id": "cbccc8bd-7bfe-46ba-b37d-49c518746271", + "metadata": {}, + "source": [ + "---\n", + "title: Data Fitting\n", + "subtitle: Variable Flip Angle\n", + "date: 2024-07-25\n", + "authors:\n", + " - name: Mathieu Boudreau\n", + " affiliations:\n", + " - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada\n", + "numbering:\n", + " figure:\n", + " template: Fig. %s\n", + "---\n", + "\n", + "At first glance, one could be tempted to fit VFA data using a non-linear least squares fitting algorithm such as Levenberg-Marquardt with Eq. 1, which typically only has two free fitting variables (T1 and M0). Although this is a valid way of estimating T1 from VFA data, it is rarely done in practice because a simple refactoring of Equation 1 allows T1 values to be estimated with a linear least square fitting algorithm, which substantially reduces the processing time. Without any approximations, Equation 1 can be rearranged into the form y = mx+b {cite:p}`Gupta1977`:\n", + "\n", + "\\begin{equation}\\tag{3}\n", + "\\frac{S_n}{ \\text{sin}(\\theta_n)} = e^{- \\frac{TR}{T_1}} \\frac{S_n}{ \\text{tan}(\\theta_n)} + C (1-e^{- \\frac{TR}{T_1}})\n", + "\\end{equation}\n", + "\n", + "As the third term does not change between measurements (it is constant for each θn), it can be grouped into the constant for a simpler representation:\n", + "\n", + "\\begin{equation}\\tag{4}\n", + "\\frac{S_n}{ \\text{sin}(\\theta_n)} = e^{- \\frac{TR}{T_1}} \\frac{S_n}{ \\text{tan}(\\theta_n)} + C\n", + "\\end{equation}\n", + "\n", + "With this rearranged form of Equation 1, T1 can be simply estimated from the slope of a linear regression calculated from Sn/sin(θn) and Sn/tan(θn) values:\n", + "\n", + "\\begin{equation}\\tag{5}\n", + "T_1 = - \\frac{TR}{ \\text{ln}(slope)}\n", + "\\end{equation}\n", + "\n", + "If data were acquired using only two flip angles – a very common VFA acquisition protocol – then the slope can be calculated using the elementary slope equation. Figure 5 displays both Equation 1 and 4 plotted for a noisy dataset." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "95947bb3-31f4-4f78-891e-1da8c2e16882", + "metadata": { + "editable": true, + "slideshow": { + "slide_type": "" + }, + "tags": [ + "hide-output", + "hide-input" + ] + }, + "outputs": [], + "source": [ + "from repo2data.repo2data import Repo2Data\n", + "import os \n", + "import pickle\n", + "import matplotlib.pyplot as plt\n", + "import chart_studio.plotly as py\n", + "import plotly.graph_objs as go\n", + "import numpy as np\n", + "from plotly import __version__\n", + "from plotly.offline import download_plotlyjs, init_notebook_mode, plot, iplot\n", + "from IPython.display import display, HTML\n", + "from plotly import tools\n", + "\n", + "data_req_path = os.path.join(\"..\",\"..\", \"binder\", \"data_requirement.json\")\n", + "repo2data = Repo2Data(data_req_path)\n", + "DATA_ROOT = os.path.join(repo2data.install()[0],\"t1-book-neurolibre\")" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "58c59812-8714-434e-8b7a-8e8177802170", + "metadata": { + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"02\",'figure_5.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " params, data_mean, data_mean_div_sin, data_mean_div_tan, data_std, data_std_div_sin, data_std_div_tan, params_highres, signal_WM, signal_WM_div_sin, signal_WM_div_tan = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "\n", + "init_notebook_mode(connected=True)\n", + "\n", + "data1 = dict(\n", + " visible = True,\n", + " x = params_highres[\"EXC_FA\"],\n", + " y = signal_WM,\n", + " name = 'Analytical Solutions',\n", + " text = params[\"EXC_FA\"],\n", + " mode = 'lines', \n", + " line = dict(\n", + " color = ('rgb(0, 0, 0)'),\n", + " dash = 'dot'),\n", + " hoverinfo='none')\n", + "\n", + "data2 = dict(\n", + " visible = True,\n", + " x = signal_WM_div_tan,\n", + " y = signal_WM_div_sin,\n", + " name = 'Analytical Solutions',\n", + " text = params_highres[\"EXC_FA\"],\n", + " mode = 'lines',\n", + " xaxis='x2',\n", + " yaxis='y2',\n", + " line = dict(\n", + " color = ('rgb(0, 0, 0)'),\n", + " dash = 'dot'\n", + " ),\n", + " hoverinfo='none',\n", + " showlegend=False)\n", + "\n", + "data3 = dict(\n", + " visible = True,\n", + " x = params[\"EXC_FA\"],\n", + " y = data_mean,\n", + " name = 'Nonlinear Form - Noisy',\n", + " text = [\"Flip angle: \" + str(x) + \"°\" for x in params[\"EXC_FA\"]],\n", + " mode = 'markers',\n", + " hoverinfo = 'y+text',\n", + " line = dict(\n", + " color = ('rgb(22, 96, 167)'),\n", + " ),\n", + " error_y=dict(\n", + " type='data',\n", + " array=data_std,\n", + " visible=True,\n", + " color = ('rgb(142, 192, 240)')\n", + " ))\n", + "\n", + "data4 = dict(\n", + " visible = True,\n", + " x = data_mean_div_tan,\n", + " y = data_mean_div_sin,\n", + " name = 'Linear Form - Noisy',\n", + " text = [\"Flip angle: \" + str(x) + \"°\" for x in params[\"EXC_FA\"]],\n", + " mode = 'markers',\n", + " xaxis='x2',\n", + " yaxis='y2',\n", + " hoverinfo = 'x+y+text',\n", + " line = dict(\n", + " color = ('rgb(205, 12, 24)'),\n", + " ),\n", + " error_x=dict(\n", + " type='data',\n", + " array=data_std_div_tan,\n", + " visible=True,\n", + " color = ('rgb(248, 135, 142)')\n", + " ),\n", + " error_y=dict(\n", + " type='data',\n", + " array=data_std_div_sin,\n", + " visible=True,\n", + " color = ('rgb(248, 135, 142)')\n", + " ))\n", + "\n", + "data = [data1, data2, data3, data4]\n", + "\n", + "layout = go.Layout(\n", + " width=580,\n", + " height=450,\n", + " margin=go.layout.Margin(\n", + " l=80,\n", + " r=80,\n", + " b=60,\n", + " t=60,\n", + " ),\n", + " annotations=[\n", + " dict(\n", + " x=0.5004254919715793,\n", + " y=-0.14,\n", + " showarrow=False,\n", + " text='Excitation Flip Angle (θn)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22,\n", + " color=('rgb(21, 91, 158)')\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=-0.17,\n", + " y=0.5,\n", + " showarrow=False,\n", + " text='Signal (Sn)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22,\n", + " color=('rgb(21, 91, 158)')\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=0.5004254919715793,\n", + " y=1.15,\n", + " showarrow=False,\n", + " text='Sn / tan(θn)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22,\n", + " color=('rgb(169, 10, 20)') \n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=1.16,\n", + " y=0.5,\n", + " showarrow=False,\n", + " text='Sn / sin(θn)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22,\n", + " color=('rgb(169, 10, 20)') \n", + " ),\n", + " xref='paper',\n", + " yref='paper',\n", + " textangle=-90,\n", + " ),\n", + " ],\n", + " xaxis=dict(\n", + " autorange=False,\n", + " range=[params['EXC_FA'][0], params['EXC_FA'][-1]],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " yaxis=dict(\n", + " autorange=True,\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " xaxis2=dict(\n", + " autorange=False,\n", + " range=[0, 1],\n", + " showgrid=False,\n", + " mirror=True,\n", + " overlaying= 'x',\n", + " anchor= 'y2',\n", + " side= 'top',\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " yaxis2=dict(\n", + " autorange=False,\n", + " range=[0, 1],\n", + " showgrid=False,\n", + " overlaying= 'y',\n", + " anchor= 'x',\n", + " side= 'right',\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " legend=dict(\n", + " x=0.32,\n", + " y=0.98,\n", + " traceorder='normal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=12,\n", + " color='#000'\n", + " ),\n", + " bordercolor='#000000',\n", + " borderwidth=2\n", + " ),\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'vfa_fig_5.html', config = config)\n", + "display(HTML('vfa_fig_5.html'))" + ] + }, + { + "cell_type": "markdown", + "id": "9cf1205d-f8f0-4d7a-a13c-78db3c3093d1", + "metadata": {}, + "source": [ + "```{figure} img/plot4.png\n", + ":label: vfaPlot4\n", + "\n", + "Mean and standard deviation of the VFA signal plotted using the nonlinear form (Equation 1 – blue) and linear form (Equation 4 – red). Monte Carlo simulation details: SNR = 25, N = 1000. VFA simulation details: TR = 25 ms, T1 = 900 ms.\n", + "```\n", + "\n", + "There are two important imaging protocol design considerations that should be taken into account when planning to use VFA: (1) how many and which flip angles to use to acquire VFA data, and (2) correcting inaccurate flip angles due to transmit RF field inhomogeneity. Most VFA experiments use the minimum number of required flip angles (two) to minimize acquisition time. For this case, it has been shown that the flip angle choice resulting in the best precision for VFA T1 estimates for a sample with a single T1 value (i.e. single tissue) are the two flip angles that result in 71% of the maximum possible steady-state signal (i.e. at the Ernst angle) {cite:p}`Deoni2003,Schabel2008`.\n", + "\n", + "Time allowing, additional flip angles are often acquired at higher values and in between the two above, because greater signal differences between tissue T1 values are present there (e.g. Figure 2). Also, for more than two flip angles, Equations 1 and 4 do not have the same noise weighting for each fitting point, which may bias linear least-square T1 estimates at lower SNRs. Thus, it has been recommended that low SNR data should be fitted with either Equation 1 using non-linear least-squares (slower fitting) or with a weighted linear least-squares form of Equation 4 {cite:p}`Chang2008`.\n", + "\n", + "Accurate knowledge of the flip angle values is very important to produce accurate T1 maps. Because of how the RF field interacts with matter {cite:p}`Sled1998`, the excitation RF field (B1+, or B1 for short) of a loaded RF coil results in spatial variations in intensity/amplitude, unless RF shimming is available to counteract this effect (not common at clinical field strengths). For quantitative measurements like VFA which are sensitive to this parameter, the flip angle can be corrected (voxelwise) relative to the nominal value by multiplying it with a scaling factor (B1) from a B1 map that is acquired during the same session:\n", + "\n", + "\\begin{equation}\\tag{6}\n", + "\\theta_{corrected} = B_1 \\theta_{nominal}\n", + "\\end{equation}\n", + "\n", + "B1 in this context is normalized, meaning that it is unitless and has a value of 1 in voxels where the RF field has the expected amplitude (i.e. where the nominal flip angle is the actual flip angle). Figure 6 displays fitted VFA T1 values from a Monte Carlo dataset simulated using biased flip angle values, and fitted without/with B1 correction." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "7a542739-7c65-4490-bb77-bc9ff0ce08ff", + "metadata": { + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"02\",'figure_6.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " B1Range, mean_T1_noB1Correction, mean_T1_withB1Correction, std_T1_noB1Correction, std_T1_withB1Correction = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "\n", + "init_notebook_mode(connected=True)\n", + "\n", + "data1 = dict(\n", + " visible = True,\n", + " x = B1Range,\n", + " y = mean_T1_noB1Correction,\n", + " name = 'Nominal flip angles',\n", + " text = 'Nominal flip angles',\n", + " mode = 'lines+markers',\n", + " hoverinfo = 'x+y+text',\n", + " line = dict(\n", + " color = ('rgb(22, 96, 167)'),\n", + " ),\n", + " error_y=dict(\n", + " type='data',\n", + " array=std_T1_noB1Correction,\n", + " visible=True,\n", + " color = ('rgb(142, 192, 240)')\n", + " ))\n", + "\n", + "data2 = dict(\n", + " visible = True,\n", + " x = B1Range,\n", + " y = mean_T1_withB1Correction,\n", + " name = 'B1-corrected flip angles',\n", + " text = 'B1-corrected flip angles',\n", + " mode = 'lines+markers',\n", + " hoverinfo = 'x+y+text',\n", + " line = dict(\n", + " color = ('rgb(205, 12, 24)'),\n", + " ),\n", + " error_y=dict(\n", + " type='data',\n", + " array=std_T1_withB1Correction,\n", + " visible=True,\n", + " color = ('rgb(248, 135, 142)')\n", + " ))\n", + "\n", + "data = [data1, data2]\n", + "\n", + "layout = go.Layout(\n", + " width=580,\n", + " height=450,\n", + " margin=go.layout.Margin(\n", + " l=80,\n", + " r=80,\n", + " b=60,\n", + " t=60,\n", + " ),\n", + " annotations=[\n", + " dict(\n", + " x=0.5004254919715793,\n", + " y=-0.14,\n", + " showarrow=False,\n", + " text='B1 (n.u.)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=-0.17,\n", + " y=0.5,\n", + " showarrow=False,\n", + " text='T1 (s)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=22\n", + " ),\n", + " textangle=-90,\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " ],\n", + " xaxis=dict(\n", + " autorange=False,\n", + " range=[B1Range[0], B1Range[-1]],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " yaxis=dict(\n", + " autorange=False,\n", + " range=[0, max(mean_T1_noB1Correction)],\n", + " showgrid=False,\n", + " linecolor='black',\n", + " linewidth=2\n", + " ),\n", + " legend=dict(\n", + " x=0.32,\n", + " y=0.98,\n", + " traceorder='normal',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=12,\n", + " color='#000'\n", + " ),\n", + " bordercolor='#000000',\n", + " borderwidth=2\n", + " ),\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'vfa_fig_6.html', config = config)\n", + "display(HTML('vfa_fig_6.html'))" + ] + }, + { + "cell_type": "markdown", + "id": "94a13435-26d7-4e8c-9d1d-e198b224746c", + "metadata": {}, + "source": [ + "```{figure} img/plot5.png\n", + ":label: vfaPlot5\n", + "\n", + "Mean and standard deviations of fitted VFA T1 values for a set of Monte Carlo simulations (SNR = 100, N = 1000), simulated using a wide range of biased flip angles and fitted without (blue) or with (red) B1 correction. Simulation parameters: TR = 25 ms, T1 = 900 ms, θnominal = 6° and 32° (optimized values for this TR/T1 combination). Notice how even after B1 correction, fitted T1 values at B1 values far from the nominal case (B1 = 1) exhibit larger variance, as the actual flip angles of the simulated signal deviate from the optimal values for this TR/T1 (Deoni et al. 2003).\n", + "```\n", + "\n", + "

\n", + "Figure 7 displays an example VFA dataset and a B1 map in a healthy brain, along with the T1 map estimated using a linear fit (Equations 4 and 5).\n", + "

" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "id": "e400a24f-4b21-4200-9f7a-ba0646c2895d", + "metadata": { + "tags": [ + "hide-input" + ] + }, + "outputs": [], + "source": [ + "filename = os.path.join(DATA_ROOT,\"02\",'figure_7.pkl')\n", + "\n", + "with open(filename, 'rb') as f:\n", + " T1_map, FA_03, FA_20, B1map, xAxis, yAxis = pickle.load(f)\n", + "\n", + "config={'showLink': False, 'displayModeBar': False}\n", + "\n", + "init_notebook_mode(connected=True)\n", + "\n", + "trace1 = go.Heatmap(x = xAxis,\n", + " y = yAxis,\n", + " z=FA_03,\n", + " colorscale='gray',\n", + " showscale = False,\n", + " visible=False,\n", + " name = 'Signal')\n", + "trace2 = go.Heatmap(x = xAxis,\n", + " y = yAxis,\n", + " z=FA_20,\n", + " colorscale='gray',\n", + " showscale = False,\n", + " visible=True,\n", + " name = 'Signal')\n", + "trace3 = go.Heatmap(x = xAxis,\n", + " y = yAxis,\n", + " z=B1map,\n", + " zmin=0.7,\n", + " zmax=1.3,\n", + " colorscale='balance',\n", + " showscale = False,\n", + " visible=False,\n", + " name = 'B1 values')\n", + "trace5 = go.Heatmap(x = xAxis,\n", + " y = yAxis,\n", + " z=T1_map,\n", + " zmin=0.0,\n", + " zmax=5000,\n", + " colorscale='Portland',\n", + " xaxis='x2',\n", + " yaxis='y2',\n", + " visible=True,\n", + " name = 'T1 values (ms)')\n", + "\n", + "data=[trace1, trace2, trace3, trace5]\n", + "\n", + "\n", + "updatemenus = list([\n", + " dict(active=1,\n", + " x = 0.09,\n", + " xanchor = 'left',\n", + " y = -0.15,\n", + " yanchor = 'bottom',\n", + " direction = 'up',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=16\n", + " ),\n", + " buttons=list([ \n", + " dict(label = '3 deg',\n", + " method = 'update',\n", + " args = [{'visible': [True, False, False, True]},\n", + " ]),\n", + " dict(label = '20 deg',\n", + " method = 'update',\n", + " args = [{'visible': [False, True, False, True]},\n", + " ]),\n", + " dict(label = 'B1 map',\n", + " method = 'update',\n", + " args = [{'visible': [False, False, True, True]},\n", + " ])\n", + " ])\n", + " )\n", + "])\n", + "\n", + "layout = dict(\n", + " width=560,\n", + " height=345,\n", + " margin = dict(\n", + " t=40,\n", + " r=50,\n", + " b=10,\n", + " l=50),\n", + " annotations=[\n", + " dict(\n", + " x=0.055,\n", + " y=1.15,\n", + " showarrow=False,\n", + " text='Input Data',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=26\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=0.6,\n", + " y=1.15,\n", + " showarrow=False,\n", + " text='T1 map',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=26\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " dict(\n", + " x=1.22,\n", + " y=1.15,\n", + " showarrow=False,\n", + " text='T1 (ms)',\n", + " font=dict(\n", + " family='Times New Roman',\n", + " size=26\n", + " ),\n", + " xref='paper',\n", + " yref='paper'\n", + " ),\n", + " ],\n", + " xaxis = dict(range = [0,127], autorange = False,\n", + " showgrid = False, zeroline = False, showticklabels = False,\n", + " ticks = '', domain=[0, 0.58]),\n", + " yaxis = dict(range = [0,127], autorange = False,\n", + " showgrid = False, zeroline = False, showticklabels = False,\n", + " ticks = '', domain=[0, 1]),\n", + " xaxis2 = dict(range = [0,127], autorange = False,\n", + " showgrid = False, zeroline = False, showticklabels = False,\n", + " ticks = '', domain=[0.40, 0.98]),\n", + " yaxis2 = dict(range = [0,127], autorange = False,\n", + " showgrid = False, zeroline = False, showticklabels = False,\n", + " ticks = '', domain=[0, 1], anchor='x2'),\n", + " showlegend = False,\n", + " autosize = False,\n", + " updatemenus=updatemenus,\n", + " plot_bgcolor='white'\n", + ")\n", + "\n", + "\n", + "fig = dict(data=data, layout=layout)\n", + "\n", + "plot(fig, filename = 'vfa_fig_7.html', config = config)\n", + "display(HTML('vfa_fig_7.html'))" + ] + }, + { + "cell_type": "markdown", + "id": "6b0850e8", + "metadata": {}, + "source": [ + "```{figure} img/plot6.png\n", + ":label: vfaPlot6\n", + "\n", + "Example variable flip angle dataset and B1 map of a healthy adult brain (left). The relevant VFA protocol parameters used were: TR = 15 ms, θnominal = 3° and 20°. The T1 map (right) was fitted using a linear regression (Equations 4 and 5).\n", + "```\n", + "\n" + ] + }, + { + "cell_type": "markdown", + "id": "fe9a024e", + "metadata": {}, + "source": [ + "\n", + "```{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 5.\n", + ":class: tip, dropdown\n", + "\n", + "```octave\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "clear all\n", + "\n", + "%% Setup parameters\n", + "% All times are in milliseconds\n", + "% All flip angles are in degrees\n", + "\n", + "params.EXC_FA = [1:4,5:5:90];\n", + "\n", + "%% Calculate signals\n", + "%\n", + "% To see all the options available, run `help vfa_t1.analytical_solution`\n", + "\n", + "params.TR = 0.025;\n", + "params.EXC_FA = [2:9,10:5:90];\n", + "\n", + "% White matter\n", + "x.M0 = 1;\n", + "x.T1 = 0.900; % in milliseconds\n", + "\n", + "Model = vfa_t1; \n", + "\n", + "Opt.SNR = 25;\n", + "Opt.TR = params.TR;\n", + "Opt.T1 = x.T1;\n", + "\n", + "clear Model.Prot.VFAData.Mat(:,1) \n", + "Model.Prot.VFAData.Mat = zeros(length(params.EXC_FA),2);\n", + "Model.Prot.VFAData.Mat(:,1) = params.EXC_FA';\n", + "Model.Prot.VFAData.Mat(:,2) = Opt.TR;\n", + "\n", + "for jj = 1:1000\n", + " [FitResult{jj}, noisyData{jj}] = Model.Sim_Single_Voxel_Curve(x,Opt,0); \n", + " fittedT1(jj) = FitResult{jj}.T1;\n", + " noisyData_array(jj,:) = noisyData{jj}.VFAData;\n", + " noisyData_array_div_sin(jj,:) = noisyData_array(jj,:) ./ sind(Model.Prot.VFAData.Mat(:,1))';\n", + " noisyData_array_div_tan(jj,:) = noisyData_array(jj,:) ./ tand(Model.Prot.VFAData.Mat(:,1))';\n", + "end\n", + " \n", + "for kk=1:length(params.EXC_FA)\n", + " data_mean(kk) = mean(noisyData_array(:,kk));\n", + " data_std(kk) = std(noisyData_array(:,kk));\n", + " \n", + " data_mean_div_sin(kk) = mean(noisyData_array_div_sin(:,kk));\n", + " data_std_div_sin(kk) = std(noisyData_array_div_sin(:,kk));\n", + " \n", + " data_mean_div_tan(kk) = mean(noisyData_array_div_tan(:,kk));\n", + " data_std_div_tan(kk) = std(noisyData_array_div_tan(:,kk));\n", + "end\n", + "\n", + "\n", + "%% Setup parameters\n", + "% All times are in milliseconds\n", + "% All flip angles are in degrees\n", + "\n", + "params_highres.EXC_FA = [2:1:90];\n", + "\n", + "%% Calculate signals\n", + "%\n", + "% To see all the options available, run `help vfa_t1.analytical_solution`\n", + "\n", + "params_highres.TR = params.TR * 1000; % in milliseconds\n", + " \n", + "% White matter\n", + "params_highres.T1 = x.T1*1000; % in milliseconds\n", + "\n", + "signal_WM = vfa_t1.analytical_solution(params_highres);\n", + "signal_WM_div_sin = signal_WM ./ sind(params_highres.EXC_FA);\n", + "signal_WM_div_tan = signal_WM ./ tand(params_highres.EXC_FA);\n", + "```\n", + "\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "cf6243e6", + "metadata": {}, + "source": [ + "\n", + "```{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 6.\n", + ":class: tip, dropdown\n", + "\n", + "```octave\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "clear all\n", + "\n", + "%% Setup parameters\n", + "% All times are in seconds\n", + "% All flip angles are in degrees\n", + "\n", + "params.TR = 0.025; % in seconds\n", + "\n", + "% White matter\n", + "params.T1 = 0.900; % in seconds\n", + "\n", + "% Calculate optimal flip angles for a two flip angle VFA experiment (for this T1 and TR)\n", + "% The will be the nominal flip angles (the flip angles assumed by the \"user\", before a \n", + "% \"realistic\"B1 bias is applied)\n", + "\n", + "nominal_EXC_FA = vfa_t1.find_two_optimal_flip_angles(params); % in degrees\n", + "disp('Nominal flip angles:')\n", + "disp(nominal_EXC_FA)\n", + "\n", + "% Range of B1 values biasing the excitation flip angle away from their nominal values\n", + "B1Range = 0.1:0.1:2;\n", + "\n", + "x.M0 = 1;\n", + "x.T1 = params.T1; % in seconds\n", + "\n", + "Model = vfa_t1; \n", + "\n", + "Opt.SNR = 100;\n", + "Opt.TR = params.TR;\n", + "Opt.T1 = x.T1;\n", + "\n", + "% Monte Carlo signal simulations\n", + "for ii = 1:1000\n", + " for jj = 1:length(B1Range)\n", + " B1 = B1Range(jj);\n", + " actual_EXC_FA = B1 * nominal_EXC_FA;\n", + " \n", + " params.EXC_FA = actual_EXC_FA;\n", + "\n", + " clear Model.Prot.VFAData.Mat(:,1)\n", + " Model.Prot.VFAData.Mat = zeros(length(params.EXC_FA),2);\n", + " Model.Prot.VFAData.Mat(:,1) = params.EXC_FA';\n", + " Model.Prot.VFAData.Mat(:,2) = Opt.TR;\n", + "\n", + " [FitResult{ii,jj}, noisyData{ii,jj}] = Model.Sim_Single_Voxel_Curve(x,Opt,0); \n", + " noisyData_array(ii,jj,:) = noisyData{ii,jj}.VFAData;\n", + " end\n", + "end\n", + "%\n", + "Model = vfa_t1; \n", + " \n", + "FlipAngle = nominal_EXC_FA';\n", + "TR = params.TR .* ones(size(FlipAngle));\n", + "\n", + "Model.Prot.VFAData.Mat = [FlipAngle TR];\n", + "\n", + "data.VFAData(:,:,1,1) = noisyData_array(:,:,1);\n", + "data.VFAData(:,:,1,2) = noisyData_array(:,:,2);\n", + "data.Mask = repmat(ones(size(B1Range)),[size(noisyData_array,1),1]);\n", + "\n", + "data.B1map = repmat(ones(size(B1Range)),[size(noisyData_array,1),1]);\n", + "FitResults_noB1Correction = FitData(data,Model,0);\n", + "\n", + "data.B1map = repmat(B1Range,[size(noisyData_array,1),1]);\n", + "FitResults_withB1Correction = FitData(data,Model,0);\n", + "\n", + "%%\n", + "%\n", + "\n", + "mean_T1_noB1Correction = mean(FitResults_noB1Correction.T1);\n", + "mean_T1_withB1Correction = mean(FitResults_withB1Correction.T1);\n", + "std_T1_noB1Correction = std(FitResults_noB1Correction.T1);\n", + "std_T1_withB1Correction = std(FitResults_withB1Correction.T1);\n", + "```\n", + "\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "ea032671", + "metadata": {}, + "source": [ + "```{admonition} Click here to view the qMRLab (MATLAB/Octave) code that generated Figure 7.\n", + ":class: tip, dropdown\n", + "\n", + "```octave\n", + "\n", + "% Verbosity level 0 overrides the disp function and supresses warnings.\n", + "% Once executed, they cannot be restored in this session\n", + "% (kernel needs to be restarted or a new notebook opened.)\n", + "VERBOSITY_LEVEL = 0;\n", + "\n", + "if VERBOSITY_LEVEL==0\n", + " % This hack was used to supress outputs from external tools\n", + " % in the Jupyter Book.\n", + " function disp(x)\n", + " end\n", + " warning('off','all')\n", + "end\n", + "\n", + "try\n", + " cd qMRLab\n", + "catch\n", + " try\n", + " cd ../../../qMRLab\n", + " catch\n", + " cd ../qMRLab\n", + " end\n", + "end\n", + "\n", + "startup\n", + "clear all\n", + "\n", + "%% MATLAB/OCTAVE CODE\n", + "\n", + "% Load data into environment, and rotate mask to be aligned with IR data\n", + "load('data/vfa_dataset/VFAData.mat');\n", + "load('data/vfa_dataset/B1map.mat');\n", + "load('data/vfa_dataset/Mask.mat');\n", + "\n", + "% Format qMRLab vfa_t1 model parameters, and load them into the Model object\n", + "Model = vfa_t1; \n", + "FlipAngle = [ 3; 20];\n", + "TR = [0.015; 0.0150];\n", + "\n", + "Model.Prot.VFAData.Mat = [FlipAngle, TR];\n", + "\n", + "% Format data structure so that they may be fit by the model\n", + "data = struct();\n", + "data.VFAData= double(VFAData);\n", + "data.B1map= double(B1map);\n", + "data.Mask= double(Mask);\n", + "\n", + "FitResults = FitData(data,Model,0); % The '0' flag is so that no wait bar is shown.\n", + "\n", + "T1_map = imrotate(FitResults.T1.*Mask,-90);\n", + "T1_map(T1_map>5)=0;\n", + "T1_map = T1_map*1000; % Convert to ms\n", + "\n", + "xAxis = [0:size(T1_map,2)-1];\n", + "yAxis = [0:size(T1_map,1)-1];\n", + "\n", + "% Raw MRI data at different TI values\n", + "FA_03 = imrotate(squeeze(VFAData(:,:,:,1).*Mask),-90);\n", + "FA_20 = imrotate(squeeze(VFAData(:,:,:,2).*Mask),-90);\n", + "B1map = imrotate(squeeze(B1map.*Mask),-90);\n", + "```\n", + "\n", + "```" + ] + }, + { + "cell_type": "markdown", + "id": "3b6ac720-ef06-4125-9c32-011a55e7ba2b", + "metadata": {}, + "source": [ + "```{admonition} References\n", + ":class: seealso\n", + "\n", + "```{bibliography}\n", + ":filter: docname in docnames\n", + "```\n", + "\n", + "```" + ] + } + ], + "metadata": { + "celltoolbar": "Tags", + "kernelspec": { + "display_name": "Python 3 (ipykernel)", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.8.18" + } + }, + "nbformat": 4, + "nbformat_minor": 5 +} diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/04-VFA_BenefitsAndPitfalls.md b/T1 Mapping/Variable Flip Angle T1 Mapping/04-VFA_BenefitsAndPitfalls.md new file mode 100644 index 0000000..f3734a1 --- /dev/null +++ b/T1 Mapping/Variable Flip Angle T1 Mapping/04-VFA_BenefitsAndPitfalls.md @@ -0,0 +1,16 @@ +--- +title: Benefits and Pitfalls +subtitle: Variable Flip Angle +date: 2024-07-25 +authors: + - name: Mathieu Boudreau + affiliations: + - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada +numbering: + figure: + template: Fig. %s +--- + +It has been well reported in recent years that the accuracy of VFA T1 estimates is very sensitive to pulse sequence implementations {cite:p}`Baudrexel2017,Lutti2013,Stikov2015`, and as such is less robust than the gold standard inversion recovery technique. In particular, the signal bias resulting from insufficient spoiling can result in inaccurate T1 estimates of up to 30% relative to inversion recovery estimated values {cite:p}`Stikov2015`. VFA T1 map accuracy and precision is also strongly dependent on the quality of the measured B1 map {cite:p}`Lee2017`, which can vary substantially between implementations {cite:p}`Boudreau2017`. Modern rapid B1 mapping pulse sequences are not as widely available as VFA, resulting in some groups attempting alternative ways of removing the bias from the T1 maps like generating an artificial B1 map through the use of image processing techniques {cite:p}`Liberman2013` or omitting B1 correction altogether {cite:p}`Yuan2012`. The latter is not recommended, because most MRI scanners have default pulse sequences that, with careful protocol settings, can provide B1 maps of sufficient quality very rapidly {cite:p}`Boudreau2017,Samson2006,Wang2005`. + +Despite some drawbacks, VFA is still one of the most widely used T1 mapping methods in research. Its rapid acquisition time, rapid image processing time, and widespread availability makes it a great candidate for use in other quantitative imaging acquisition protocols like quantitative magnetization transfer imaging {cite:p}`Cercignani2005,Yarnykh2002` and dynamic contrast enhanced imaging {cite:p}`Li2018,Sung2013`. \ No newline at end of file diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation1.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation1.png new file mode 100644 index 0000000..01e594d Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation1.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation2.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation2.png new file mode 100644 index 0000000..d5628ab Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation2.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation3.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation3.png new file mode 100644 index 0000000..5ab882d Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation3.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation4.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation4.png new file mode 100644 index 0000000..cdb6fac Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation4.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation5.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation5.png new file mode 100644 index 0000000..3903466 Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation5.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation6.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation6.png new file mode 100644 index 0000000..c7bd1f4 Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equation6.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/equations_all.svg b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equations_all.svg new file mode 100644 index 0000000..b637c34 --- /dev/null +++ b/T1 Mapping/Variable Flip Angle T1 Mapping/img/equations_all.svg @@ -0,0 +1,1916 @@ + + + + + + + + + + image/svg+xml + + + + + + + (1) + + + + + + + + + + + + + + + + + + + + + + + + + + + diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot1.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot1.png new file mode 100644 index 0000000..cd62eaa Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot1.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot2.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot2.png new file mode 100644 index 0000000..fb115ba Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot2.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot3.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot3.png new file mode 100644 index 0000000..4967519 Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot3.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot4.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot4.png new file mode 100644 index 0000000..5452394 Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot4.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot5.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot5.png new file mode 100644 index 0000000..1c23719 Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot5.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot6.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot6.png new file mode 100644 index 0000000..aabab04 Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/plot6.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/vfa_pulsesequence.png b/T1 Mapping/Variable Flip Angle T1 Mapping/img/vfa_pulsesequence.png new file mode 100644 index 0000000..15554fa Binary files /dev/null and b/T1 Mapping/Variable Flip Angle T1 Mapping/img/vfa_pulsesequence.png differ diff --git a/T1 Mapping/Variable Flip Angle T1 Mapping/img/vfa_pulsesequence.svg b/T1 Mapping/Variable Flip Angle T1 Mapping/img/vfa_pulsesequence.svg new file mode 100644 index 0000000..452dd11 --- /dev/null +++ b/T1 Mapping/Variable Flip Angle T1 Mapping/img/vfa_pulsesequence.svg @@ -0,0 +1,365 @@ + + + + + + + + + + + + + + + + + + + + image/svg+xml + + + + + + + + + + + + + + + + + + + + + + + + + IMG + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + SPOIL + + TR + + θn + θn + + + diff --git a/T1 Mapping/foreword.md b/T1 Mapping/foreword.md new file mode 100644 index 0000000..59b0636 --- /dev/null +++ b/T1 Mapping/foreword.md @@ -0,0 +1,13 @@ +--- +title: About this chapter +date: 2024-07-25 +authors: + - name: Mathieu Boudreau + affiliations: + - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada +numbering: + figure: + template: Fig. %s +--- + +Lorem ipsum dolor sit amet, consectetur adipiscing elit. Morbi aliquam, augue nec varius tincidunt, diam turpis pharetra ex, quis egestas tellus libero in erat. Etiam libero purus, volutpat ac tempus vel, rutrum ut orci. Mauris et ipsum gravida, porta augue sed, rutrum nibh. Vestibulum mollis nibh eget ipsum ullamcorper, vel luctus nisl luctus. Mauris dictum blandit ullamcorper. Nulla scelerisque a mauris ac accumsan. Fusce vel placerat arcu. Donec consectetur aliquet purus, sed pulvinar diam luctus cursus. Mauris blandit ex vel est eleifend accumsan. Nulla ligula odio, iaculis id justo eget, facilisis mollis ante. Etiam sollicitudin scelerisque ante vel suscipit. Nullam condimentum nunc at augue interdum, ac mollis orci blandit. In convallis erat viverra tellus interdum efficitur. Maecenas commodo facilisis odio quis elementum. \ No newline at end of file diff --git a/T1 Mapping/mp2rage T1 Mapping/01-Abstract.md b/T1 Mapping/mp2rage T1 Mapping/01-Abstract.md new file mode 100644 index 0000000..0951e0d --- /dev/null +++ b/T1 Mapping/mp2rage T1 Mapping/01-Abstract.md @@ -0,0 +1,17 @@ +--- +title: Abstract +subtitle: MP2RAGE +date: 2024-07-25 +authors: + - name: Mathieu Boudreau + affiliations: + - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada +numbering: + heading_2: false + figure: + template: Fig. %s +--- + +Dictionary-based MRI techniques capable of generating T1 maps are increasing in popularity, due to their growing availability on clinical scanners, rapid scan times, and fast post-processing computation time, thus making quantitative T1 mapping accessible for clinical applications. Generally speaking, dictionary-based quantitative MRI techniques use numerical dictionaries—databases of pre-calculated signal values simulated for a wide range of tissue and protocol combinations—during the image reconstruction or post-processing stages. Popular examples of dictionary-based techniques that have been applied to T1 mapping are MR Fingerprinting (MRF) (Ma et al. 2013), certain flavours of compressed sensing (CS) (Doneva et al. 2010; Li et al. 2012), and Magnetization Prepared 2 Rapid Acquisition Gradient Echoes (MP2RAGE) (Marques et al. 2010). Dictionary-based techniques can usually be classified into one of two categories: techniques that use information redundancy from parametric data to assist in accelerated imaging (e.g. CS, MRF), or those that use dictionaries to estimate quantitative maps using the MR images after reconstruction. Because MP2RAGE is a technique implemented primarily for T1 mapping, and it is becoming increasingly available as a standard pulse sequence on many MRI systems, the remainder of this section will focus solely on this technique. However, many concepts discussed are shared by other dictionary-based techniques. + +MP2RAGE is an extension of the conventional MPRAGE pulse sequence widely used in clinical studies (Haase et al. 1989; Mugler & Brookeman 1990). A simplified version of the MP2RAGE pulse sequence is shown in Figure 1. MP2RAGE can be seen as a hybrid between the inversion recovery and VFA pulse sequences: a 180° inversion pulse is used to prepare the magnetization for T1 sensitivity at the beginning of each TRMP2RAGE, and then two images are acquired at different inversion times using gradient recalled echo (GRE) imaging blocks with low flip angles and short repetition times (TR). During a given GRE imaging block, each excitation pulse is followed by a constant in-plane (“y”) phase encode weighting (varied for each TRMP2RAGE), but with different 3D (“z”) phase encoding gradients (varied at each TR). The center of k-space for the 3D phase encoding direction is acquired at the TI for each GRE imaging block. The main motivation for developing the MP2RAGE pulse sequence was to provide a metric similar to MPRAGE, but with self-bias correction of the static (B0) and receive (B1-) magnetic fields, and a first order correction of the transmit magnetic field (B1+). However, because two images at different TIs are acquired (unlike MPRAGE, which only acquires data at a single TI), information about the T1 values can also be inferred, thus making it possible to generate quantitative T1 maps using this data. \ No newline at end of file diff --git a/T1 Mapping/mp2rage T1 Mapping/02-SignalModelling.md b/T1 Mapping/mp2rage T1 Mapping/02-SignalModelling.md new file mode 100644 index 0000000..60f0542 --- /dev/null +++ b/T1 Mapping/mp2rage T1 Mapping/02-SignalModelling.md @@ -0,0 +1,53 @@ +--- +title: Signal Modelling +subtitle: MP2RAGE +date: 2024-07-25 +authors: + - name: Mathieu Boudreau + affiliations: + - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada +numbering: + heading_2: false + figure: + template: Fig. %s +--- + + + +Prior to considering the full signal equations, we will first introduce the equation for the MP2RAGE parameter (SMP2RAGE) that is calculated in addition to the T1 map. For complex data (magnitude and phase, or real and imaginary), the MP2RAGE signal (SMP2RAGE) is calculated from the images acquired at two TIs (SGRE,TI1 and SGRE,TI2) using the following expression (Marques et al. 2010): + +```{figure} img/equation1.png +:label: mp2ragefig1 +``` + +This value is bounded between [-0.5, 0.5], and helps reduce some B0 inhomogeneity effects using the phase data. For real data, or magnitude data with polarity restoration, this metric is instead calculated as: + + +```{figure} img/equation2.png +:label: mp2ragefig2 +``` + +Because MP2RAGE is a hybrid of pulse sequences used for inversion recovery and VFA, the resulting signal equations are more complex. Typically, a steady state is not achieved during the short train of GRE imaging blocks, so the signal at the center of k-space for each readout (which defines the contrast weighting) will depend on the number of phase-encoding steps. For simplicity, the equations presented here assume that the 3D phase-encoding dimension is fully sampled (no partial Fourier or parallel imaging acceleration). For this case (see appendix of (Marques et al. 2010) for derivation details), the signal equations are: + + +```{figure} img/equation3.png +:label: mp2ragefig3 +``` + +```{figure} img/equation4.png +:label: mp2ragefig14 +``` + +where B1- is the receive field sensitivity, “eff” is the adiabatic inversion pulse efficiency, ER = exp(-TR/T1), EA = exp(-TA/T1), EB = exp(-TB/T1), EC = exp(-TC/T1). The variables TA, TB, and TC are the three different delay times (TA: time between inversion pulse and beginning of the GRE1 block, TB: time between the end of GRE1 and beginning of GRE2, TC: time between the end of GRE2 and the end of the TR). If no k-space acceleration is used (e.g. no partial Fourier or parallel imaging acceleration), then these values are TA = TI1 - (n/2)TR, TB = TI2 - (TI1 + nTR), and TC = TRMP2RAGE - (TI2 + (n/2)TR), where n is the number of voxels acquired in the 3D phase encode direction varied within each GRE block. The value mz,ss is the steady-state longitudinal magnetization prior to the inversion pulse, and is given by: + + +```{figure} img/equation5.png +:label: mp2ragefig5 +``` + + +```{figure} img/equation6.png +:label: mp2ragefig6 +``` + +From Eqs. (3–6), it is evident that the MP2RAGE parameter SMP2RAGE (Eqs. (1.11) and (1.12)) cancels out the effects of receive field sensitivity, T2*, and M0. The signal sensitivity related to the transmit field (B1+), hidden in Eqs. (3–6) within the flip angle values θ1 and θ2, can also be reduced by careful pulse sequence protocol design (Marques et al. 2010), but not entirely eliminated (Marques & Gruetter 2013). \ No newline at end of file diff --git a/T1 Mapping/mp2rage T1 Mapping/03-DataFitting.md b/T1 Mapping/mp2rage T1 Mapping/03-DataFitting.md new file mode 100644 index 0000000..ee985d2 --- /dev/null +++ b/T1 Mapping/mp2rage T1 Mapping/03-DataFitting.md @@ -0,0 +1,33 @@ +--- +title: Data Fitting +subtitle: MP2RAGE +date: 2024-07-25 +authors: + - name: Mathieu Boudreau + affiliations: + - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada +numbering: + figure: + template: Fig. %s +--- + +Dictionary-based techniques such as MP2RAGE do not typically use conventional minimization algorithms (e.g. Levenberg-Marquardt) to fit signal equations to observed data. Instead, the MP2RAGE technique uses pre-calculated signal values for a wide range of parameter values (e.g. T1), and then interpolation is done within this dictionary of values to estimate the T1 value that matches the observed signal. This approach results in rapid post-processing times because the dictionaries can be simulated/generated prior to scanning and interpolating between these values is much faster than most fitting algorithms. This means that the quantitative image can be produced and displayed directly on the MRI scanner console rather than needing to be fitted offline. + +MP2RAGE is an extension of the conventional MPRAGE pulse sequence widely used in clinical studies (Haase et al. 1989; Mugler & Brookeman 1990). A simplified version of the MP2RAGE pulse sequence is shown in Figure 1. MP2RAGE can be seen as a hybrid between the inversion recovery and VFA pulse sequences: a 180° inversion pulse is used to prepare the magnetization for T1 sensitivity at the beginning of each TRMP2RAGE, and then two images are acquired at different inversion times using gradient recalled echo (GRE) imaging blocks with low flip angles and short repetition times (TR). During a given GRE imaging block, each excitation pulse is followed by a constant in-plane (“y”) phase encode weighting (varied for each TRMP2RAGE), but with different 3D (“z”) phase encoding gradients (varied at each TR). The center of k-space for the 3D phase encoding direction is acquired at the TI for each GRE imaging block. The main motivation for developing the MP2RAGE pulse sequence was to provide a metric similar to MPRAGE, but with self-bias correction of the static (B0) and receive (B1-) magnetic fields, and a first order correction of the transmit magnetic field (B1+). However, because two images at different TIs are acquired (unlike MPRAGE, which only acquires data at a single TI), information about the T1 values can also be inferred, thus making it possible to generate quantitative T1 maps using this data. + +```{figure} img/plot1.png +:label: mp2rageplot1 + +T1 lookup table as a function of B1 and SMP2RAGE value. Inversion times used to acquire this magnitude image dataset were 800 ms and 2700 ms, the flip angles were 4° and 5° (respectively), TRMP2RAGE = 6000 ms, and TR = 6.7 ms. The code that was used were shared open sourced by the authors of the original MP2RAGE paper (https://github.com/JosePMarques/MP2RAGE-related-scripts). +``` + +To produce T1 maps with good accuracy and precision using dictionary-based interpolation methods, it is important that the signal curves are unique for each parameter value. MP2RAGE can produce good T1 maps by using a dictionary with only dimensions (T1, SMP2RAGE), since SMP2RAGE is unique for each T1 value for a given protocol (Marques et al. 2010). However, as was noted above, SMP2RAGE is also sensitive to B1 because of θ1 and θ2 in Eqs. (1.13–1.16). The B1-sensitivity can be reduced substantially with careful MP2RAGE protocol optimization (Marques et al. 2010), and further improved by including B1 as one of the dictionary dimensions [T1, B1, SMP2RAGE] (Figure 1.15). This requires an additional acquisition of a B1 map (Marques & Gruetter 2013), which lengthens the scan time. + + +```{figure} img/plot2.png +:label: mp2rageplot2 + +Example MP2RAGE dataset of a healthy adult brain at 7T and T1 map. Inversion times used to acquire this magnitude image dataset were 800 ms and 2700 ms, the flip angles were 4° and 5° (respectively), TRMP2RAGE = 6000 ms, and TR = 6.7 ms. The dataset and code that was used were shared open sourced by the authors of the original MP2RAGE paper (https://github.com/JosePMarques/MP2RAGE-related-scripts). +``` + +The MP2RAGE pulse sequence is increasingly being distributed by MRI vendors, thus typically a data fitting package is also available to reconstruct the T1 maps online. Alternatively, several open source packages to create T1 maps from MP2RAGE data are available online (Marques 2017; de Hollander 2017), and for new users these are recommended—as opposed to programming one from scratch—as there are many potential pitfalls (e.g. adjusting the equations to handle partial Fourier or parallel imaging acceleration). \ No newline at end of file diff --git a/T1 Mapping/mp2rage T1 Mapping/04-BenefitsAndPitfalls.md b/T1 Mapping/mp2rage T1 Mapping/04-BenefitsAndPitfalls.md new file mode 100644 index 0000000..1078fc0 --- /dev/null +++ b/T1 Mapping/mp2rage T1 Mapping/04-BenefitsAndPitfalls.md @@ -0,0 +1,18 @@ +--- +title: Benefits and Pitfalls +subtitle: MP2RAGE +date: 2024-07-25 +authors: + - name: Mathieu Boudreau + affiliations: + - NeuroPoly Lab, Polytechnique Montreal, Quebec, Canada +numbering: + figure: + template: Fig. %s +--- + +This widespread availability and its turnkey acquisition/fitting procedures are a main contributing factor to the growing interest for including quantitative T1 maps in clinical and neuroscience studies. T1 values measured using MP2RAGE showed high levels of reproducibility for the brain in an inter- and intra-site study at eight sites (same MRI hardware/software and at 7 T) of two subjects (Voelker et al. 2016). Not only does MP2RAGE have one of the fastest acquisition and post-processing times among quantitative T1 mapping techniques, but it can accomplish this while acquiring very high resolution T1 maps (1 mm isotropic at 3T and submillimeter at 7T, both in under 10 min (Fujimoto et al. 2014)), opening the doors to cortical studies which greatly benefit from the smaller voxel size (Waehnert et al. 2014; Beck et al. 2018; Haast et al. 2018). +

+ +

+Despite these benefits, MP2RAGE and similar dictionary-based techniques have certain limitations that are important to consider before deciding to incorporate them in a study. Good reproducibility of the quantitative T1 map is dependent on using one pre-calculated dictionary. If two different dictionaries are used (e.g. cross-site with different MRI vendors), the differences in the dictionary interpolations will likely result in minor differences in T1 estimates for the same data. Also, although the B1-sensitivity of the MP2RAGE T1 maps can be reduced with proper protocol optimization, it can be substantial enough that further correction using a measured B1 map should be done (Marques & Gruetter 2013; Haast et al. 2018). However B1 mapping brings an additional potential source of error, so carefully selecting a B1 mapping technique and accompanying post-processing method (e.g. filtering) should be done before integrating it in a T1 mapping protocol (Boudreau et al. 2017). Lastly, the MP2RAGE equations (and thus, dictionaries) assume monoexponential longitudinal relaxation, and this has been shown to result in suboptimal estimates of the long T1 component for a biexponential relaxation model (Rioux et al. 2016), an effect that becomes more important at higher fields. \ No newline at end of file diff --git a/T1 Mapping/mp2rage T1 Mapping/img/equation1.png b/T1 Mapping/mp2rage T1 Mapping/img/equation1.png new file mode 100644 index 0000000..b82bd62 Binary files /dev/null and b/T1 Mapping/mp2rage T1 Mapping/img/equation1.png differ diff --git a/T1 Mapping/mp2rage T1 Mapping/img/equation2.png b/T1 Mapping/mp2rage T1 Mapping/img/equation2.png new file mode 100644 index 0000000..01f6525 Binary files /dev/null and b/T1 Mapping/mp2rage T1 Mapping/img/equation2.png differ diff --git a/T1 Mapping/mp2rage T1 Mapping/img/equation3.png b/T1 Mapping/mp2rage T1 Mapping/img/equation3.png new file mode 100644 index 0000000..a4c66b3 Binary files /dev/null and b/T1 Mapping/mp2rage T1 Mapping/img/equation3.png differ diff --git a/T1 Mapping/mp2rage T1 Mapping/img/equation4.png b/T1 Mapping/mp2rage T1 Mapping/img/equation4.png new file mode 100644 index 0000000..8e8e3a7 Binary files /dev/null and b/T1 Mapping/mp2rage T1 Mapping/img/equation4.png differ diff --git a/T1 Mapping/mp2rage T1 Mapping/img/equation5.png b/T1 Mapping/mp2rage T1 Mapping/img/equation5.png new file mode 100644 index 0000000..d1286e4 Binary files /dev/null and b/T1 Mapping/mp2rage T1 Mapping/img/equation5.png differ diff --git a/T1 Mapping/mp2rage T1 Mapping/img/equation6.png b/T1 Mapping/mp2rage T1 Mapping/img/equation6.png new file mode 100644 index 0000000..dc99697 Binary files /dev/null and b/T1 Mapping/mp2rage T1 Mapping/img/equation6.png differ diff --git a/T1 Mapping/mp2rage T1 Mapping/img/mp2rage_pulsesequence.png b/T1 Mapping/mp2rage T1 Mapping/img/mp2rage_pulsesequence.png new file mode 100644 index 0000000..7b84416 Binary files /dev/null and b/T1 Mapping/mp2rage T1 Mapping/img/mp2rage_pulsesequence.png differ diff --git a/T1 Mapping/mp2rage T1 Mapping/img/mp2rage_pulsesequence.svg b/T1 Mapping/mp2rage T1 Mapping/img/mp2rage_pulsesequence.svg new file mode 100644 index 0000000..d887016 --- /dev/null +++ b/T1 Mapping/mp2rage T1 Mapping/img/mp2rage_pulsesequence.svg @@ -0,0 +1,870 @@ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + image/svg+xml + + + + + + + + + + + + + + + + + + + + + + + + + + θ1 + + + + θ1 + + + + θ1 + + + + θ1 + + + + θ1 + + + + + + + + + + + + + + + + + + + + + 180° + 180° + + + + θ2 + + + + θ2 + + + + θ2 + + + + θ2 + + + + θ2 + + + + + + + + + + + + + + + + + + + + + + + + + TR + + TRMP2RAGE + nTR + + + TI1 + + TI2 + 1st GRE image + 2nd GRE image + + diff --git a/T1 Mapping/mp2rage T1 Mapping/img/plot1.png b/T1 Mapping/mp2rage T1 Mapping/img/plot1.png new file mode 100644 index 0000000..6885013 Binary files /dev/null and b/T1 Mapping/mp2rage T1 Mapping/img/plot1.png differ diff --git a/T1 Mapping/mp2rage T1 Mapping/img/plot2.png b/T1 Mapping/mp2rage T1 Mapping/img/plot2.png new file mode 100644 index 0000000..2b45b9c Binary files /dev/null and b/T1 Mapping/mp2rage T1 Mapping/img/plot2.png differ diff --git a/content/banner.jpg b/banner.jpg similarity index 100% rename from content/banner.jpg rename to banner.jpg diff --git a/book.bib b/bibliography/book.bib similarity index 100% rename from book.bib rename to bibliography/book.bib diff --git a/bibliography/ir.bib b/bibliography/ir.bib new file mode 100644 index 0000000..b4169c0 --- /dev/null +++ b/bibliography/ir.bib @@ -0,0 +1,531 @@ +@incollection{Boudreau2019, +title = {Chapter 2 - Quantitative T1 and T1ρ Mapping}, +editor = {Nicole Seiberlich and Vikas Gulani and Fernando Calamante and Adrienne Campbell-Washburn and Mariya Doneva and Houchun Harry Hu and Steven Sourbron}, +series = {Advances in Magnetic Resonance Technology and Applications}, +publisher = {Academic Press}, +volume = {1}, +pages = {19-45}, +year = {2020}, +booktitle = {Quantitative Magnetic Resonance Imaging}, +issn = {2666-9099}, +doi = {10.1016/B978-0-12-817057-1.00004-4}, +author = {Mathieu Boudreau and Kathryn E. Keenan and Nikola Stikov}, +keywords = {, , Relaxometry, Inversion recovery, Variable flip angle, VFA, DESPOT1, MP2RAGE}, +abstract = {The spin-lattice relaxation time (T1) is one of the fundamental parameters of Magnetic Resonance Imaging (MRI). It characterizes the rate at which the longitudinal magnetization component recovers to its equilibrium state. T1 has been established as a sensitive biomarker for tissue characterization, and is also an important consideration in pulse sequence development and optimization. The development of T1 mapping techniques has been an active field of research since the early days of NMR, and a wide range of strategies have emerged for measuring T1, each with unique benefits and pitfalls. This chapter covers several strategies for measuring the T1 parameter. First, we introduce the inversion recovery (IR) T1 mapping technique. Limited by its long acquisition time, IR is widely used as a reference measure due to its high accuracy and robustness. The second section is on the variable flip angle (VFA, a.k.a. DESPOT1) technique, which leverages the steady-state signal sensitivity with respect to the excitation flip angles and T1. Capable of rapid 3D imaging, VFA is dependent on another quantitative MRI map, the radio-frequency transmit field (B1+, or B1). The third section presents a widely used dictionary-based T1 mapping technique, MP2RAGE. Boasting rapid acquisition and reconstruction times, MP2RAGE is a promising technique to get T1 maps used in the clinic. The final section presents an overview of spin-lattice relaxation time in the rotating frame (T1ρ) and its acquisition techniques. T1ρ is more sensitive to molecular motion in the kHz range, and has been extensively used to study the musculoskeletal system.} +} + +@article{Karakuzu2020-ul, + author = {Karakuzu, Agah and Boudreau, Mathieu and Duval, Tanguy and +Boshkovski, Tommy and Leppert, Ilana and Cabana, Jean-Fran{\c +c}ois and Gagnon, Ian and Beliveau, Pascale and Pike, G and +Cohen-Adad, Julien and Stikov, Nikola}, + copyright = {http://creativecommons.org/licenses/by/4.0/}, + doi = {10.21105/joss.02343}, + journal = {J. Open Source Softw.}, + month = {September}, + number = {53}, + pages = {2343}, + publisher = {The Open Journal}, + title = {{qMRLab}: Quantitative {MRI} analysis, under one umbrella}, + volume = {5}, + year = {2020} +} + +@article{Boudreau2023, + doi = {10.55458/neurolibre.00014}, + year = {2023}, + publisher = {NeuroLibre}, + author = {Mathieu Boudreau and Agah Karakuzu and Julien Cohen-Adad and Ecem Bozkurt and Madeline Carr and Marco Castellaro and Luis Concha and Mariya Doneva and Seraina Dual and Alex Ensworth and Alexandru Foias and Véronique Fortier and Refaat E. Gabr and Guillaume Gilbert and Carri K. Glide-Hurst and Matthew Grech-Sollars and Siyuan Hu and Oscar Jalnefjord and Jorge Jovicich and Kübra Keskin and Peter Koken and Anastasia Kolokotronis and Simran Kukran and Nam. G. Lee and Ives R. Levesque and Bochao Li and Dan Ma and Burkhard Mädler and Nyasha Maforo and Jamie Near and Erick Pasaye and Alonso Ramirez-Manzanares and Ben Statton and Christian Stehning and Stefano Tambalo and Ye Tian and Chenyang Wang and Kilian Weis and Niloufar Zakariaei and Shuo Zhang and Ziwei Zhao and Nikola Stikov}, + title = {Results of the ISMRM 2020 joint Reproducible Research & Quantitative MR study groups reproducibility challenge on phantom and human brain T1 mapping}, + journal = {NeuroLibre Reproducible Preprints} } + +@article{Stikov2015, + abstract = {Purpose There are many T1 mapping methods available, each of them validated in phantoms and reporting excellent agreement with literature. However, values in literature vary greatly, with T1 in white matter ranging from 690 to 1100 ms at 3 Tesla. This brings into question the accuracy of one of the most fundamental measurements in quantitative MRI. Our goal was to explain these variations and look into ways of mitigating them. Theory and Methods We evaluated the three most common T1 mapping methods (inversion recovery, Look-Locker, and variable flip angle) through Bloch simulations, a white matter phantom and the brains of 10 healthy subjects (single-slice). We pooled the T1 histograms of the subjects to determine whether there is a sequence-dependent bias and whether it is reproducible across subjects. Results We found good agreement between the three methods in phantoms, but poor agreement in vivo, with the white matter T1 histogram peak in healthy subjects varying by more than 30\% depending on the method used. We also found that the pooled brain histograms displayed three distinct white matter peaks, with Look-Locker consistently underestimating, and variable flip angle overestimating the inversion recovery T1 values. The Bloch simulations indicated that incomplete spoiling and inaccurate B1 mapping could account for the observed differences. Conclusion We conclude that the three most common T1 mapping protocols produce stable T1 values in phantoms, but not in vivo. To improve the accuracy of T1 mapping, we recommend that sites perform in vivo validation of their T1 mapping method against the inversion recovery reference method, as the first step toward developing a robust calibration scheme. Magn Reson Med 73:514–522, 2015. © 2014 Wiley Periodicals, Inc.}, + author = {Stikov, Nikola and Boudreau, Mathieu and Levesque, Ives R. and Tardif, Christine L. and Barral, Joëlle K. and Pike, G. Bruce}, + doi = {10.1002/mrm.25135}, + eprint = {https://onlinelibrary.wiley.com/doi/pdf/10.1002/mrm.25135}, + journal = {Magnetic Resonance in Medicine}, + keywords = {relaxometry, T1 mapping, quantitative MRI, accuracy, precision, inversion recovery, Look-Locker, variable flip angle, B1 mapping}, + number = {2}, + pages = {514-522}, + title = {On the accuracy of T1 mapping: Searching for common ground}, + volume = {73}, + year = {2015} +} + +@book{Fukushima1981, + title={Experimental pulse NMR: a nuts and bolts approach}, + author={Fukushima, Eiichi}, + year={1981}, + publisher={CRC Press}, + doi = {10.1201/9780429493867} +} + +@article{Barral2010-qm, + author = {Barral, Jo{\"e}lle K and Gudmundson, Erik and Stikov, Nikola and +Etezadi-Amoli, Maryam and Stoica, Petre and Nishimura, Dwight G}, + doi = {10.1002/mrm.22497}, + journal = {Magn. Reson. Med.}, + language = {en}, + month = {October}, + number = {4}, + pages = {1057--1067}, + title = {A robust methodology for in vivo {T1} mapping}, + volume = {64}, + year = {2010} +} + +@article{Pykett1983, + title={Measurement of spin-lattice relaxation times in nuclear magnetic resonance imaging}, + author={Pykett, IL and Rosen, BR and Buonanno, FS and Brady, TJ}, + journal={Physics in Medicine \& Biology}, + volume={28}, + number={6}, + pages={723}, + year={1983}, + publisher={IOP Publishing}, + doi = {10.1088/0031-9155/28/6/012} +} + +@article{Steen1994, + title={Precise and accurate measurement of proton T1 in human brain in vivo: validation and preliminary clinical application}, + author={Steen, R Grant and Gronemeyer, Suzanne A and Kingsley, Peter B and Reddick, Wilbum E and Langston, James S and Taylor, June S}, + journal={Journal of Magnetic Resonance Imaging}, + volume={4}, + number={5}, + pages={681--691}, + year={1994}, + publisher={Wiley Online Library}, + doi = {10.1002/jmri.1880040511} +} + +@article{Pykett1978, + title={A line scan image study of a tumorous rat leg by NMR}, + author={Pykett, IL and Mansfield, Peter}, + journal={Physics in Medicine \& Biology}, + volume={23}, + number={5}, + pages={961}, + year={1978}, + publisher={IOP Publishing}, + doi = {10.1088/0031-9155/23/5/012} +} + +@article{Piechnik2010, + title={Shortened Modified Look-Locker Inversion recovery (ShMOLLI) for clinical myocardial T1-mapping at 1.5 and 3 T within a 9 heartbeat breathhold}, + author={Piechnik, Stefan K and Ferreira, Vanessa M and Dall'Armellina, Erica and Cochlin, Lowri E and Greiser, Andreas and Neubauer, Stefan and Robson, Matthew D}, + journal={Journal of cardiovascular magnetic resonance}, + volume={12}, + number={1}, + pages={1--11}, + year={2010}, + publisher={BioMed Central}, + doi = {10.1186/1532-429X-12-69} +} + +@article{Messroghli2004, + title={Modified Look-Locker inversion recovery (MOLLI) for high-resolution T1 mapping of the heart}, + author={Messroghli, Daniel R and Radjenovic, Aleksandra and Kozerke, Sebastian and Higgins, David M and Sivananthan, Mohan U and Ridgway, John P}, + journal={Magnetic Resonance in Medicine: An Official Journal of the International Society for Magnetic Resonance in Medicine}, + volume={52}, + number={1}, + pages={141--146}, + year={2004}, + publisher={Wiley Online Library}, + doi = {10.1002/mrm.20110} +} + +@article{Look1970, + title={Time saving in measurement of NMR and EPR relaxation times}, + author={Look, David C and Locker, Donald R}, + journal={Review of scientific instruments}, + volume={41}, + number={2}, + pages={250--251}, + year={1970}, + publisher={American Institute of Physics}, + doi = {10.1063/1.1684482} +} + +@article{Hahn1949, + title={An accurate nuclear magnetic resonance method for measuring spin-lattice relaxation times}, + author={Hahn, Erwin L}, + journal={Physical Review}, + volume={76}, + number={1}, + pages={145}, + year={1949}, + publisher={APS}, + doi = {10.1103/PhysRev.76.145} +} + +@article{Gai2013, + title={Modified Look-Locker T1 evaluation using Bloch simulations: human and phantom validation}, + author={Gai, Neville D and Stehning, Christian and Nacif, Marcelo and Bluemke, David A}, + journal={Magnetic resonance in medicine}, + volume={69}, + number={2}, + pages={329--336}, + year={2013}, + publisher={Wiley Online Library}, + doi = {10.1002/mrm.24251} +} + +@article{Drain1949, + title={A direct method of measuring nuclear spin-lattice relaxation times}, + author={Drain, LE}, + journal={Proceedings of the Physical Society. Section A}, + volume={62}, + number={5}, + pages={301}, + year={1949}, + publisher={IOP Publishing}, + doi = {10.1088/0370-1298/62/5/306} +} + +@article{Baudrexel2017, + doi = {10.1002/mrm.26979}, + year = {2017}, + month = oct, + publisher = {Wiley}, + volume = {79}, + number = {6}, + pages = {3082--3092}, + author = {Simon Baudrexel and Ulrike N\"{o}th and Jan-R\"{u}diger Sch\"{u}re and Ralf Deichmann}, + title = {T1 mapping with the variable flip angle technique: A simple correction for insufficient spoiling of transverse magnetization}, + journal = {Magnetic Resonance in Medicine} +} + +@book{Handbook2004, + doi = {10.1016/b978-0-12-092861-3.x5000-6}, + author={Bernstein, Matt A and King, Kevin F and Zhou, Xiaohong Joe}, + year = {2004}, + publisher = {Elsevier}, + title = {Handbook of {MRI} Pulse Sequences} +} + +@article{Boudreau2017, + doi = {10.1002/jmri.25692}, + year = {2017}, + month = mar, + publisher = {Wiley}, + volume = {46}, + number = {6}, + pages = {1673--1682}, + author = {Mathieu Boudreau and Christine L. Tardif and Nikola Stikov and John G. Sled and Wayne Lee and G. Bruce Pike}, + title = {B1 mapping for bias-correction in quantitative T1 imaging of the brain at 3T using standard pulse sequences}, + journal = {Journal of Magnetic Resonance Imaging} +} + +@article{Cercignani2005, + doi = {10.1016/j.neuroimage.2005.04.031}, + year = {2005}, + month = aug, + publisher = {Elsevier {BV}}, + volume = {27}, + number = {2}, + pages = {436--441}, + author = {Mara Cercignani and Mark R. Symms and Klaus Schmierer and Philip A. Boulby and Daniel J. Tozer and Maria Ron and Paul S. Tofts and Gareth J. Barker}, + title = {Three-dimensional quantitative magnetisation transfer imaging of the human brain}, + journal = {{NeuroImage}} +} + +@article{Chang2008, + doi = {10.1002/mrm.21669}, + year = {2008}, + month = jul, + publisher = {Wiley}, + volume = {60}, + number = {2}, + pages = {496--501}, + author = {Lin-Ching Chang and Cheng Guan Koay and Peter J. Basser and Carlo Pierpaoli}, + title = {Linear least-squares method for unbiased estimation of T1 from {SPGR} signals}, + journal = {Magnetic Resonance in Medicine} +} + +@article{Christensen1974, + title={Optimal determination of relaxation times of Fourier transform nuclear magnetic resonance. Determination of spin-lattice relaxation times in chemically polarized species}, + author={Christensen, Kenner A and Grant, David M and Schulman, Edward M and Walling, Cheves}, + journal={The Journal of Physical Chemistry}, + volume={78}, + number={19}, + pages={1971--1977}, + year={1974}, + doi = {10.1021/j100612a022}, + publisher={ACS Publications} +} + +@article{Feinberg2010, + doi = {10.1371/journal.pone.0015710}, + year = {2010}, + month = dec, + publisher = {Public Library of Science ({PLoS})}, + volume = {5}, + number = {12}, + pages = {e15710}, + author = {David A. Feinberg and Steen Moeller and Stephen M. Smith and Edward Auerbach and Sudhir Ramanna and Matt F. Glasser and Karla L. Miller and Kamil Ugurbil and Essa Yacoub}, + editor = {Pedro Antonio Valdes-Sosa}, + title = {Multiplexed Echo Planar Imaging for Sub-Second Whole Brain {FMRI} and Fast Diffusion Imaging}, + journal = {{PLoS} {ONE}} +} + +@article{Ernst1966, + doi = {10.1063/1.1719961}, + year = {1966}, + month = jan, + publisher = {{AIP} Publishing}, + volume = {37}, + number = {1}, + pages = {93--102}, + author = {R. R. Ernst and W. A. Anderson}, + title = {Application of Fourier Transform Spectroscopy to Magnetic Resonance}, + journal = {Review of Scientific Instruments} +} + +@article{Fram1987, + doi = {10.1016/0730-725x(87)90021-x}, + year = {1987}, + month = jan, + publisher = {Elsevier {BV}}, + volume = {5}, + number = {3}, + pages = {201--208}, + author = {Evan K. Fram and Robert J. Herfkens and G.Allan Johnson and Gary H. Glover and John P. Karis and Ann Shimakawa and Tom G. Perkins and Norbert J. Pelc}, + title = {Rapid calculation of T1 using variable flip angle gradient refocused imaging}, + journal = {Magnetic Resonance Imaging} +} + +@article{Gupta1977, + doi = {10.1016/0022-2364(77)90138-x}, + year = {1977}, + month = jan, + publisher = {Elsevier {BV}}, + volume = {25}, + number = {1}, + pages = {231--235}, + author = {Raj K Gupta}, + title = {A new look at the method of variable nutation angle for the measurement of spin-lattice relaxation times using fourier transform {NMR}}, + journal = {Journal of Magnetic Resonance (1969)} +} + +@article{Homer1985, + doi = {10.1016/0022-2364(85)90318-x}, + year = {1985}, + month = jun, + publisher = {Elsevier {BV}}, + volume = {63}, + number = {2}, + pages = {287--297}, + author = {John Homer and Martin S Beevers}, + title = {Driven-equilibrium single-pulse observation of T1 relaxation. A reevaluation of a rapid {\textquotedblleft}new{\textquotedblright} method for determining {NMR} spin-lattice relaxation times}, + journal = {Journal of Magnetic Resonance (1969)} +} + +@article{Lee2017, + doi = {10.3389/fnins.2017.00106}, + year = {2017}, + month = mar, + publisher = {Frontiers Media {SA}}, + volume = {11}, + author = {Yoojin Lee and Martina F. Callaghan and Zoltan Nagy}, + title = {Analysis of the Precision of Variable Flip Angle T1 Mapping with Emphasis on the Noise Propagated from {RF} Transmit Field Maps}, + journal = {Frontiers in Neuroscience} +} + +@article{Liberman2013, + doi = {10.1002/jmri.24373}, + year = {2013}, + month = nov, + publisher = {Wiley}, + volume = {40}, + number = {1}, + pages = {171--180}, + author = {Gilad Liberman and Yoram Louzoun and Dafna Ben Bashat}, + title = {T1 Mapping using variable flip angle {SPGR} data with flip angle correction}, + journal = {Journal of Magnetic Resonance Imaging} +} + +@article{Li2018, + doi = {10.1088/1361-6560/aad519}, + year = {2018}, + month = aug, + publisher = {{IOP} Publishing}, + volume = {63}, + number = {16}, + pages = {16NT01}, + author = {Z F Li and W Zhao and T F Qi and C Gao and Q Gu and J S Zhao and T S Koh}, + title = {A simple B1correction method for dynamic contrast-enhanced {MRI}}, + journal = {Physics in Medicine and Biology} +} + + +@inproceedings{Lutti2013, + title={Optimizing the accuracy of T1 mapping accounting for RF non-linearities and spoiling characteristics in FLASH imaging}, + author={Lutti, Antoine and Weiskopf, Nikolaus}, + booktitle={Proceedings of the 21st Annual Meeting of ISMRM, Salt Lake City, Utah, USA}, + pages={2478}, + year={2013} +} + +@article{Samson2006, + doi = {10.1016/j.mri.2005.10.025}, + year = {2006}, + month = apr, + publisher = {Elsevier {BV}}, + volume = {24}, + number = {3}, + pages = {255--263}, + author = {Rebecca S. Samson and Claudia A.M. Wheeler-Kingshott and Mark R. Symms and Daniel J. Tozer and Paul S. Tofts}, + title = {A simple correction for B1 field errors in magnetization transfer ratio measurements}, + journal = {Magnetic Resonance Imaging} +} + +@article{Schabel2008, + doi = {10.1088/0031-9155/54/1/n01}, + year = {2008}, + month = dec, + publisher = {{IOP} Publishing}, + volume = {54}, + number = {1}, + pages = {N1--N8}, + author = {Matthias C Schabel and Glen R Morrell}, + title = {Uncertainty in T1 mapping using the variable flip angle method with two flip angles}, + journal = {Physics in Medicine and Biology} +} + +@article{Sled1998, + doi = {10.1109/42.730409}, + year = {1998}, + publisher = {Institute of Electrical and Electronics Engineers ({IEEE})}, + volume = {17}, + number = {4}, + pages = {653--662}, + author = {J.G. Sled and G.B. Pike}, + title = {Standing-wave and {RF} penetration artifacts caused by elliptic geometry: an electrodynamic analysis of {MRI}}, + journal = {{IEEE} Transactions on Medical Imaging} +} + +@article{Sung2013, + doi = {10.1002/jmri.23996}, + year = {2013}, + month = jan, + publisher = {Wiley}, + volume = {38}, + number = {2}, + pages = {454--459}, + author = {Kyunghyun Sung and Bruce L. Daniel and Brian A. Hargreaves}, + title = {Transmit field inhomogeneity and T1 estimation errors in breast {DCE}-{MRI} at 3 tesla}, + journal = {Journal of Magnetic Resonance Imaging} +} + +@article{Wang2005, + doi = {10.1002/mrm.20377}, + year = {2005}, + month = feb, + publisher = {Wiley}, + volume = {53}, + number = {3}, + pages = {666--674}, + author = {Jinghua Wang and Maolin Qiu and R. Todd Constable}, + title = {In vivo method for correcting transmit/receive nonuniformities with phased array coils}, + journal = {Magnetic Resonance in Medicine} +} + +@article{Yarnykh2010, + doi = {10.1002/mrm.22394}, + year = {2010}, + month = may, + publisher = {Wiley}, + volume = {63}, + number = {6}, + pages = {1610--1626}, + author = {Vasily L. Yarnykh}, + title = {Optimal radiofrequency and gradient spoiling for improved accuracy of T1 and B1 measurements using fast steady-state techniques}, + journal = {Magnetic Resonance in Medicine} +} + +@article{Yarnykh2002, + doi = {10.1002/mrm.10120}, + year = {2002}, + month = apr, + publisher = {Wiley}, + volume = {47}, + number = {5}, + pages = {929--939}, + author = {Vasily L. Yarnykh}, + title = {Pulsed Z-spectroscopic imaging of cross-relaxation parameters in tissues for human {MRI}: Theory and clinical applications}, + journal = {Magnetic Resonance in Medicine} +} + +@article{Yuan2012, + author = {Jing Yuan and Steven Kwok Keung Chow and David Ka Wai Yeung and Anil T Ahuja and Ann D King}, + title = {Quantitative evaluation of dual-flip-angle T 1 mapping on DCE-MRI kinetic parameter estimation in head and neck}, + journal = {Quantitative Imaging in Medicine and Surgery}, + volume = {2}, + number = {4}, + year = {2012}, + issn = {2223-4306}, + doi = {10.3978/j.issn.2223-4292.2012.11.04} +} + +@article{Zur1991, + doi = {10.1002/mrm.1910210210}, + year = {1991}, + month = oct, + publisher = {Wiley}, + volume = {21}, + number = {2}, + pages = {251--263}, + author = {Y. Zur and M. L. Wood and L. J. Neuringer}, + title = {Spoiling of transverse magnetization in steady-state sequences}, + journal = {Magnetic Resonance in Medicine} +} + +@article{Deoni2003, + doi = {10.1002/mrm.10407}, + year = {2003}, + month = mar, + publisher = {Wiley}, + volume = {49}, + number = {3}, + pages = {515--526}, + author = {Sean C.L. Deoni and Brian K. Rutt and Terry M. Peters}, + title = {Rapid combined T1 and T2 mapping using gradient recalled acquisition in the steady state}, + journal = {Magnetic Resonance in Medicine} +} +@misc{Karakuzu2022-nlwf, + title={NeuroLibre : A preprint server for full-fledged reproducible neuroscience}, + url={osf.io/h89js}, + DOI={10.31219/osf.io/h89js}, + publisher={OSF Preprints}, + author={Karakuzu, Agah and DuPre, Elizabeth and Tetrel, Loic and Bermudez, Patrick and Boudreau, Mathieu and Chin, Mary and Poline, Jean-Baptiste and Das, Samir and Bellec, Pierre and Stikov, Nikola}, + year={2022}, + month={Apr} +} + +@article{Dupre2022-iro, + title={Beyond advertising: New infrastructures for publishing integrated research objects}, + author={DuPre, Elizabeth and Holdgraf, Chris and Karakuzu, Agah and Tetrel, Lo{\"\i}c and Bellec, Pierre and Stikov, Nikola and Poline, Jean-Baptiste}, + journal={PLOS Computational Biology}, + doi = {10.1371/journal.pcbi.1009651}, + volume={18}, + number={1}, + pages={e1009651}, + year={2022}, + publisher={Public Library of Science San Francisco, CA USA} +} + +@article{Harding2023-conp, + title = {The {Canadian} {Open} {Neuroscience} {Platform}—{An} open science framework for the neuroscience community}, + volume = {19}, + url = {10.1371/journal.pcbi.1011230}, + doi = {10.1371/journal.pcbi.1011230}, + abstract = {The Canadian Open Neuroscience Platform (CONP) takes a multifaceted approach to enabling open neuroscience, aiming to make research, data, and tools accessible to everyone, with the ultimate objective of accelerating discovery. Its core infrastructure is the CONP Portal, a repository with a decentralized design, where datasets and analysis tools across disparate platforms can be browsed, searched, accessed, and shared in accordance with FAIR principles. Another key piece of CONP infrastructure is NeuroLibre, a preprint server capable of creating and hosting executable and fully reproducible scientific publications that embed text, figures, and code. As part of its holistic approach, the CONP has also constructed frameworks and guidance for ethics and data governance, provided support and developed resources to help train the next generation of neuroscientists, and has fostered and grown an engaged community through outreach and communications. In this manuscript, we provide a high-level overview of this multipronged platform and its vision of lowering the barriers to the practice of open neuroscience and yielding the associated benefits for both individual researchers and the wider community.}, + number = {7}, + journal = {PLOS Computational Biology}, + author = {Harding, Rachel J. and Bermudez, Patrick and Bernier, Alexander and Beauvais, Michael and Bellec, Pierre and Hill, Sean and Karakuzu, Agah and Knoppers, Bartha M. and Pavlidis, Paul and Poline, Jean-Baptiste and Roskams, Jane and Stikov, Nikola and Stone, Jessica and Strother, Stephen and Consortium, CONP and Evans, Alan C.}, + month = jul, + year = {2023}, + note = {Publisher: Public Library of Science}, + pages = {1--14}, +} \ No newline at end of file diff --git a/_toc.yml b/extra.yml similarity index 55% rename from _toc.yml rename to extra.yml index 21a5dfa..04f898d 100644 --- a/_toc.yml +++ b/extra.yml @@ -1,20 +1,19 @@ format: jb-book root: content/foreword.md parts: - - caption: Basics of MRI and qMRI - chapters: - - file: content/intro/intro - - file: content/intro/million_dollar - - file: content/intro/journey - - file: content/intro/encode - - file: content/intro/sequences - - file: content/intro/qmri_practice - - caption: T1 Mapping - chapters: - - file: content/t1/into - - file: content/t1/ir - - file: content/t1/vfa - - file: content/t1/mp2rage + caption: Basics of MRI and qMRI + chapters: + - file: content/intro/intro + - file: content/intro/million_dollar + - file: content/intro/journey + - file: content/intro/encode + - file: content/intro/sequences + - file: content/intro/qmri_practice + caption: T1 Mapping + chapters: + - content/t1/t1_main.md + sections: + - file: content/t1/ir/introduction # - caption: T2 Mapping # chapters: # - file: path/to/part2/chapter1 diff --git a/content/foreword.md b/index.md similarity index 93% rename from content/foreword.md rename to index.md index b1aab72..df33855 100644 --- a/content/foreword.md +++ b/index.md @@ -1,13 +1,7 @@ --- title: qMRI MOOC -subtitle: A massive open online course powered by qMRLab +description: A massive open online course powered by qMRLab date: 2024-07-25 -editors: - - name: Agah Karakuzu - orcid: 0000-0001-7283-271X - affiliations: - - NeuroPoly, Polytechnique Montreal, Quebec, Canada - - Montreal Heart Institute, Montreal, Quebec, Canada --- ```{figure} banner.jpg diff --git a/myst.yml b/myst.yml index b1a65be..d4dad34 100644 --- a/myst.yml +++ b/myst.yml @@ -14,7 +14,7 @@ project: - quantitative MRI - education github: https://github.com/qmrlab/moooc - banner: content/banner.jpg + banner: banner.jpg venue: qMRLab subject: MOOC numbering: @@ -47,8 +47,11 @@ project: RF: Radiofrequency FA: Flip Angle VFA: Variable Flip Angle + exclude: + - README.md bibliography: - - book.bib + - bibliography/book.bib + - bibliography/ir.bib site: template: book-theme logo: public/icon.png