diff --git a/advanced_topics/08_fslpy.ipynb b/advanced_topics/08_fslpy.ipynb
index 94680b2bb1909b7474ac69807fb0b355bfac4e07..546c1f7b39cafe34a67a6c4630595744e0c14a79 100644
--- a/advanced_topics/08_fslpy.ipynb
+++ b/advanced_topics/08_fslpy.ipynb
@@ -21,17 +21,12 @@
"perform analyses and image processing in conjunction with FSL.\n",
"\n",
"\n",
- "> **Note**: `fslpy` is distinct from `fslpython` - `fslpython` is the Python\n",
- "> environment that is baked into FSL. `fslpy` is a Python library which is\n",
- "> installed into the `fslpython` environment.\n",
- "\n",
- "\n",
"* [The `Image` class, and other data types](#the-image-class-and-other-data-types)\n",
" * [Creating images](#creating-images)\n",
" * [Working with image data](#working-with-image-data)\n",
" * [Loading other file types](#loading-other-file-types)\n",
" * [NIfTI coordinate systems](#nifti-coordinate-systems)\n",
- " * [Image processing](#image-processing)\n",
+ " * [Transformations and resampling](#transformations-and-resampling)\n",
"* [FSL wrapper functions](#fsl-wrapper-functions)\n",
" * [In-memory images](#in-memory-images)\n",
" * [Loading outputs into Python](#loading-outputs-into-python)\n",
@@ -51,6 +46,11 @@
" * [Working with atlases](#working-with-atlases)\n",
"\n",
"\n",
+ "> **Note**: `fslpy` is distinct from `fslpython` - `fslpython` is the Python\n",
+ "> environment that is baked into FSL. `fslpy` is a Python library which is\n",
+ "> installed into the `fslpython` environment.\n",
+ "\n",
+ "\n",
"Let's start with some standard imports and environment set-up:"
]
},
@@ -67,7 +67,8 @@
"import nibabel as nib\n",
"import numpy as np\n",
"import warnings\n",
- "warnings.filterwarnings(\"ignore\")"
+ "warnings.filterwarnings(\"ignore\")\n",
+ "np.set_printoptions(suppress=True, precision=4)"
]
},
{
@@ -97,6 +98,8 @@
" `fig` argument to plot overlays).\n",
" \"\"\"\n",
"\n",
+ " voxel = [int(round(v)) for v in voxel]\n",
+ "\n",
" data = np.asanyarray(data, dtype=np.float)\n",
" data[data <= 0] = np.nan\n",
"\n",
@@ -189,6 +192,14 @@
"- Ability to load metadata from BIDS sidecar files\n",
"\n",
"\n",
+ "> The `Image` class behaves differently to the `nibabel.Nifti1Image`. For\n",
+ "> example, when you create an `Image` object, the default behaviour is to load\n",
+ "> the image data into memory. This is configurable however; take a look at\n",
+ "> [the\n",
+ "> documentation](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.data.image.html#fsl.data.image.Image)\n",
+ "> to explore all of the options.\n",
+ "\n",
+ "\n",
"Some simple image processing routines are also provided - these are covered\n",
"[below](#image-processing).\n",
"\n",
@@ -213,8 +224,8 @@
"# load a FSL image - the file\n",
"# suffix is optional, just like\n",
"# in real FSL-land!\n",
- "img = Image(op.join(stddir, 'MNI152_T1_1mm'))\n",
- "print(img)"
+ "std1mm = Image(op.join(stddir, 'MNI152_T1_1mm'))\n",
+ "print(std1mm)"
]
},
{
@@ -233,7 +244,7 @@
"# load a nibabel image, and\n",
"# convert it into an FSL image\n",
"nibimg = nib.load(op.join(stddir, 'MNI152_T1_1mm.nii.gz'))\n",
- "img = Image(nibimg)"
+ "std1mm = Image(nibimg)"
]
},
{
@@ -249,10 +260,28 @@
"metadata": {},
"outputs": [],
"source": [
- "data = np.zeros((100, 100, 100))\n",
+ "data = np.zeros((182, 218, 182))\n",
"img = Image(data, xform=np.eye(4))"
]
},
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "If you have generated some data from another `Image` (or from a\n",
+ "`nibabel.Nifti1Image`) you can use the `header` option to set\n",
+ "the header information on the new image:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "img = Image(data, header=std1mm.header)"
+ ]
+ },
{
"cell_type": "markdown",
"metadata": {},
@@ -266,7 +295,8 @@
"metadata": {},
"outputs": [],
"source": [
- "img.save('empty.nii.gz')"
+ "img.save('empty')\n",
+ "!ls"
]
},
{
@@ -403,7 +433,7 @@
"\n",
"* The\n",
" [`fsl.data.bitmap.Bitmap`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.data.bitmap.html)\n",
- " class can be used to load a bitmap image (e.g. `jpg, `png`, etc) and\n",
+ " class can be used to load a bitmap image (e.g. `jpg`, `png`, etc) and\n",
" convert it to a NIfTI image.\n",
"* The\n",
" [`fsl.data.dicom.DicomImage`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.data.dicom.html)\n",
@@ -471,12 +501,27 @@
"\n",
"# Apply the world->voxel\n",
"# affine to the coordinates\n",
- "voxcoords = (np.dot(world2vox[:3, :3], mnicoords.T)).T + world2vox[:3, 3]\n",
- "\n",
- "# The code above is a bit fiddly, so\n",
- "# instead of figuring it out, you can\n",
- "# just use the transform() function:\n",
+ "voxcoords = (np.dot(world2vox[:3, :3], mnicoords.T)).T + world2vox[:3, 3]"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "The code above is a bit fiddly, so instead of figuring it out, you can just\n",
+ "use the\n",
+ "[`affine.transform`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.affine.html#fsl.transform.affine.transform)\n",
+ "function:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
"from fsl.transform.affine import transform\n",
+ "\n",
"voxcoords = transform(mnicoords, world2vox)\n",
"\n",
"# just to double check, let's transform\n",
@@ -615,14 +660,15 @@
"# To combine affines together, we\n",
"# have to list them in reverse -\n",
"# linear algebra is *weird*.\n",
- "funcvox2mni = concat(fsl2mni, func2std, vox2fsl)"
+ "funcvox2mni = concat(fsl2mni, func2std, vox2fsl)\n",
+ "print(funcvox2mni)"
]
},
{
"cell_type": "markdown",
"metadata": {},
"source": [
- "> Below we will use the\n",
+ "> In the next section we will use the\n",
"> [`fsl.transform.flirt.fromFlirt`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.flirt.html#fsl.transform.flirt.fromFlirt)\n",
"> function, which does all of the above for us.\n",
"\n",
@@ -648,15 +694,14 @@
"cell_type": "markdown",
"metadata": {},
"source": [
- "> Note that in the above example we are only applying a linear transformation\n",
- "> into MNI space - in reality you would also want to apply your non-linear\n",
- "> structural-to-standard transformation too. But this is left as [an exercise\n",
- "> for the\n",
- "> reader](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.fnirt.html).\n",
+ "Note that in the above example we are only applying a linear transformation\n",
+ "into MNI space - in reality you would also want to apply your non-linear\n",
+ "structural-to-standard transformation too. This is covered in the next\n",
+ "section.\n",
"\n",
"\n",
- "<a class=\"anchor\" id=\"image-processing\"></a>\n",
- "### Image processing\n",
+ "<a class=\"anchor\" id=\"transformations-and-resampling\"></a>\n",
+ "### Transformations and resampling\n",
"\n",
"\n",
"Now, it's all well and good to look at t-statistic values and voxel\n",
@@ -665,8 +710,13 @@
"this, we're going to resample our functional image into MNI space, so we can\n",
"overlay it on the MNI template. This can be done using some handy functions\n",
"from the\n",
+ "[`fsl.transform.flirt`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.flirt.html)\n",
+ "and\n",
"[`fsl.utils.image.resample`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.utils.image.resample.html)\n",
- "module:"
+ "modules.\n",
+ "\n",
+ "\n",
+ "Let's make sure we've got our source and reference images loaded:"
]
},
{
@@ -675,28 +725,51 @@
"metadata": {},
"outputs": [],
"source": [
- "from fsl.transform.flirt import fromFlirt\n",
- "from fsl.utils.image.resample import resampleToReference\n",
- "\n",
"featdir = op.join(op.join('08_fslpy', 'fmri.feat'))\n",
"tstat1 = Image(op.join(featdir, 'stats', 'tstat1'))\n",
- "std = Image(op.expandvars(op.join('$FSLDIR', 'data', 'standard', 'MNI152_T1_2mm')))\n",
+ "std = Image(op.expandvars(op.join('$FSLDIR', 'data', 'standard', 'MNI152_T1_2mm')))"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Now we'll load the `example_func2standard` FLIRT matrix, and adjust it so that\n",
+ "it transforms from functional *world* coordinates into standard *world*\n",
+ "coordinates - this is what is expected by the `resampleToReference` function,\n",
+ "used below:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from fsl.transform.flirt import fromFlirt\n",
"\n",
- "# Load the func2standard FLIRT matrix, and adjust it\n",
- "# so that it transforms from functional *world*\n",
- "# coordinates into standard *world* coordinates -\n",
- "# this is what is expected by the resampleToReference\n",
- "# function, used below\n",
"func2std = np.loadtxt(op.join(featdir, 'reg', 'example_func2standard.mat'))\n",
- "func2std = fromFlirt(func2std, tstat1, std, 'world', 'world')\n",
- "\n",
- "# All of the functions in the resample module\n",
- "# return a numpy array containing the resampled\n",
- "# data, and an adjusted voxel-to-world affine\n",
- "# transformation. But when using the\n",
- "# resampleToReference function, the affine will\n",
- "# be the same as the MNI152 2mm affine, so we\n",
- "# can ignore it.\n",
+ "func2std = fromFlirt(func2std, tstat1, std, 'world', 'world')"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Now we can use `resampleToReference` to resample our functional data into\n",
+ "MNI152 space. This function returns a `numpy` array containing the resampled\n",
+ "data, and an adjusted voxel-to-world affine transformation. But in this case,\n",
+ "we know that the data will be aligned to MNI152, so we can ignore the affine:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from fsl.utils.image.resample import resampleToReference\n",
+ "\n",
"std_tstat1 = resampleToReference(tstat1, std, func2std)[0]\n",
"std_tstat1 = Image(std_tstat1, header=std.header)"
]
@@ -725,15 +798,198 @@
"fig = ortho(std_tstat1, mnivoxels, cmap=plt.cm.inferno, fig=fig, cursor=True)"
]
},
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "In the example above, we resampled some data from functional space into\n",
+ "standard space using a linear transformation. But we all know that this is not\n",
+ "how things work in the real world - linear transformations are for kids. The\n",
+ "real world is full of lions and tigers and bears and warp fields.\n",
+ "\n",
+ "\n",
+ "The\n",
+ "[`fsl.transform.fnirt`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.fnirt.html#fsl.transform.fnirt.fromFnirt)\n",
+ "and\n",
+ "[`fsl.transform.nonlinear`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.nonlinear.html)\n",
+ "modules contain classes and functions for working with FNIRT-style warp fields\n",
+ "(modules for working with lions, tigers, and bears are still under\n",
+ "development).\n",
+ "\n",
+ "\n",
+ "Let's imagine that we have defined an ROI in MNI152 space, and we want to\n",
+ "project it into the space of our functional data. We can do this by combining\n",
+ "the nonlinear structural to standard registration produced by FNIRT with the\n",
+ "linear functional to structural registration generated by FLIRT. First of\n",
+ "all, we'll load images from each of the functional, structural, and standard\n",
+ "spaces:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "featdir = op.join('08_fslpy', 'fmri.feat')\n",
+ "func = Image(op.join(featdir, 'reg', 'example_func'))\n",
+ "struc = Image(op.join(featdir, 'reg', 'highres'))\n",
+ "std = Image(op.expandvars(op.join('$FSLDIR', 'data', 'standard', 'MNI152_T1_2mm')))"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Now, let's say we have obtained our seed location in MNI152 coordinates. Let's\n",
+ "convert them to MNI152 voxels just to double check:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "seedmni = [-48, -74, -9]\n",
+ "seedmnivox = transform(seedmni, std.getAffine('world', 'voxel'))\n",
+ "ortho(std.data, seedmnivox, cursor=True)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Now we'll load the FNIRT warp field, which encodes a nonlinear transformation\n",
+ "from structural space to standard space. FNIRT warp fields are often stored as\n",
+ "*coefficient* fields to reduce the file size, but in order to use it, we must\n",
+ "convert the coefficient field into a *deformation* (a.k.a. *displacement*)\n",
+ "field. This takes a few seconds:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from fsl.transform.fnirt import readFnirt\n",
+ "from fsl.transform.nonlinear import coefficientFieldToDeformationField\n",
+ "\n",
+ "struc2std = readFnirt(op.join(featdir, 'reg', 'highres2standard_warp'), struc, std)\n",
+ "struc2std = coefficientFieldToDeformationField(struc2std)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "We'll also load our FLIRT functional to structural transformation, adjust it\n",
+ "so that it transforms between voxel coordinate systems instead of the FSL\n",
+ "coordinate system, and invert so it can transform from structural voxels to\n",
+ "functional voxels:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from fsl.transform.affine import invert\n",
+ "func2struc = np.loadtxt(op.join(featdir, 'reg', 'example_func2highres.mat'))\n",
+ "func2struc = fromFlirt(func2struc, func, struc, 'voxel', 'voxel')\n",
+ "struc2func = invert(func2struc)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "Now we can transform our seed coordinates from MNI152 space into functional\n",
+ "space in two stages. First, we'll use our deformation field to transform from\n",
+ "MNI152 space into structural space:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "seedstruc = struc2std.transform(seedmni, 'world', 'voxel')\n",
+ "seedfunc = transform(seedstruc, struc2func)\n",
+ "\n",
+ "print('Seed location in MNI coordinates: ', seedmni)\n",
+ "print('Seed location in functional voxels:', seedfunc)\n",
+ "ortho(func.data, seedfunc, cursor=True)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "> FNIRT warp fields kind of work backwards - we can use them to transform\n",
+ "> reference coordinates into source coordinates, but would need to invert the\n",
+ "> warp field using `invwarp` if we wanted to transform from source coordinates\n",
+ "> into referemce coordinates.\n",
+ "\n",
+ "\n",
+ "Of course, we can also use our deformation field to resample an image from\n",
+ "structural space into MNI152 space. The `applyDeformation` function takes an\n",
+ "`Image` and a `DeformationField`, and returns a `numpy` array containing the\n",
+ "resampled data."
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "from fsl.transform.nonlinear import applyDeformation\n",
+ "\n",
+ "strucmni = applyDeformation(struc, struc2std)\n",
+ "\n",
+ "# remove low-valued voxels,\n",
+ "# just for visualisation below\n",
+ "strucmni[strucmni < 1] = 0\n",
+ "\n",
+ "fig = ortho(std.data, [45, 54, 45], cmap=plt.cm.gray)\n",
+ "fig = ortho(strucmni, [45, 54, 45], fig=fig)"
+ ]
+ },
+ {
+ "cell_type": "markdown",
+ "metadata": {},
+ "source": [
+ "The `premat` option to `applyDeformation` can be used to specify our linear\n",
+ "functional to structural transformation, and hence resample a functional image\n",
+ "into MNI152 space:"
+ ]
+ },
+ {
+ "cell_type": "code",
+ "execution_count": null,
+ "metadata": {},
+ "outputs": [],
+ "source": [
+ "tstatmni = applyDeformation(tstat1, struc2std, premat=func2struc)\n",
+ "tstatmni[tstatmni < 3] = 0\n",
+ "\n",
+ "fig = ortho(std.data, [45, 54, 45], cmap=plt.cm.gray)\n",
+ "fig = ortho(tstatmni, [45, 54, 45], fig=fig)"
+ ]
+ },
{
"cell_type": "markdown",
"metadata": {},
"source": [
"There are a few other useful functions tucked away in the\n",
- "[fsl.utils.image](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.utils.image.html)\n",
- "package, with more to be added in the future. The [`fsl.transform`]() package\n",
- "also contains a wealth of functionality for working with linear (FLIRT) and\n",
- "non-linear (FNIRT) transformations.\n",
+ "[`fsl.utils.image`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.utils.image.html)\n",
+ "and\n",
+ "[`fsl.transform`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.html)\n",
+ "packages, with more to be added in the future.\n",
"\n",
"\n",
"<a class=\"anchor\" id=\"fsl-wrapper-functions\"></a>\n",
@@ -761,7 +1017,7 @@
"metadata": {},
"outputs": [],
"source": [
- "from fsl.wrappers import bet, robustfov, LOAD\n",
+ "from fsl.wrappers import robustfov\n",
"\n",
"robustfov('08_fslpy/bighead', 'bighead_cropped')\n",
"\n",
@@ -805,6 +1061,8 @@
"metadata": {},
"outputs": [],
"source": [
+ "from fsl.wrappers import bet\n",
+ "\n",
"bet('bighead_cropped', 'bighead_cropped_brain', f=0.3, m=True, s=True)\n",
"\n",
"render('bighead_cropped -b 40 '\n",
@@ -878,6 +1136,8 @@
"metadata": {},
"outputs": [],
"source": [
+ "from fsl.wrappers import LOAD\n",
+ "\n",
"cropped = Image('bighead_cropped')\n",
"\n",
"# The loaded result is called \"output\",\n",
diff --git a/advanced_topics/08_fslpy.md b/advanced_topics/08_fslpy.md
index 1db751764b7892b694af22d8dcb463314fa0ea2f..043dfabaf0380154302c9d51374f81d20607f864 100644
--- a/advanced_topics/08_fslpy.md
+++ b/advanced_topics/08_fslpy.md
@@ -15,17 +15,12 @@ This practical highlights some of the most useful features provided by
perform analyses and image processing in conjunction with FSL.
-> **Note**: `fslpy` is distinct from `fslpython` - `fslpython` is the Python
-> environment that is baked into FSL. `fslpy` is a Python library which is
-> installed into the `fslpython` environment.
-
-
* [The `Image` class, and other data types](#the-image-class-and-other-data-types)
* [Creating images](#creating-images)
* [Working with image data](#working-with-image-data)
* [Loading other file types](#loading-other-file-types)
* [NIfTI coordinate systems](#nifti-coordinate-systems)
- * [Image processing](#image-processing)
+ * [Transformations and resampling](#transformations-and-resampling)
* [FSL wrapper functions](#fsl-wrapper-functions)
* [In-memory images](#in-memory-images)
* [Loading outputs into Python](#loading-outputs-into-python)
@@ -45,6 +40,11 @@ perform analyses and image processing in conjunction with FSL.
* [Working with atlases](#working-with-atlases)
+> **Note**: `fslpy` is distinct from `fslpython` - `fslpython` is the Python
+> environment that is baked into FSL. `fslpy` is a Python library which is
+> installed into the `fslpython` environment.
+
+
Let's start with some standard imports and environment set-up:
@@ -57,6 +57,7 @@ import nibabel as nib
import numpy as np
import warnings
warnings.filterwarnings("ignore")
+np.set_printoptions(suppress=True, precision=4)
```
@@ -78,6 +79,8 @@ def ortho(data, voxel, fig=None, cursor=False, **kwargs):
`fig` argument to plot overlays).
"""
+ voxel = [int(round(v)) for v in voxel]
+
data = np.asanyarray(data, dtype=np.float)
data[data <= 0] = np.nan
@@ -159,6 +162,14 @@ as:
- Ability to load metadata from BIDS sidecar files
+> The `Image` class behaves differently to the `nibabel.Nifti1Image`. For
+> example, when you create an `Image` object, the default behaviour is to load
+> the image data into memory. This is configurable however; take a look at
+> [the
+> documentation](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.data.image.html#fsl.data.image.Image)
+> to explore all of the options.
+
+
Some simple image processing routines are also provided - these are covered
[below](#image-processing).
@@ -178,8 +189,8 @@ stddir = op.expandvars('${FSLDIR}/data/standard/')
# load a FSL image - the file
# suffix is optional, just like
# in real FSL-land!
-img = Image(op.join(stddir, 'MNI152_T1_1mm'))
-print(img)
+std1mm = Image(op.join(stddir, 'MNI152_T1_1mm'))
+print(std1mm)
```
@@ -190,7 +201,7 @@ You can create an `Image` from an existing `nibabel` image:
# load a nibabel image, and
# convert it into an FSL image
nibimg = nib.load(op.join(stddir, 'MNI152_T1_1mm.nii.gz'))
-img = Image(nibimg)
+std1mm = Image(nibimg)
```
@@ -198,15 +209,26 @@ Or you can create an `Image` from a `numpy` array:
```
-data = np.zeros((100, 100, 100))
+data = np.zeros((182, 218, 182))
img = Image(data, xform=np.eye(4))
```
+If you have generated some data from another `Image` (or from a
+`nibabel.Nifti1Image`) you can use the `header` option to set
+the header information on the new image:
+
+
+```
+img = Image(data, header=std1mm.header)
+```
+
+
You can save an image to file via the `save` method:
```
-img.save('empty.nii.gz')
+img.save('empty')
+!ls
```
@@ -300,7 +322,7 @@ corresponding `nibabel` types:
* The
[`fsl.data.bitmap.Bitmap`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.data.bitmap.html)
- class can be used to load a bitmap image (e.g. `jpg, `png`, etc) and
+ class can be used to load a bitmap image (e.g. `jpg`, `png`, etc) and
convert it to a NIfTI image.
* The
[`fsl.data.dicom.DicomImage`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.data.dicom.html)
@@ -364,11 +386,18 @@ vox2world = std2mm.getAffine('voxel', 'world')
# Apply the world->voxel
# affine to the coordinates
voxcoords = (np.dot(world2vox[:3, :3], mnicoords.T)).T + world2vox[:3, 3]
+```
+
+
+The code above is a bit fiddly, so instead of figuring it out, you can just
+use the
+[`affine.transform`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.affine.html#fsl.transform.affine.transform)
+function:
+
-# The code above is a bit fiddly, so
-# instead of figuring it out, you can
-# just use the transform() function:
+```
from fsl.transform.affine import transform
+
voxcoords = transform(mnicoords, world2vox)
# just to double check, let's transform
@@ -468,9 +497,10 @@ from fsl.transform.affine import concat
# have to list them in reverse -
# linear algebra is *weird*.
funcvox2mni = concat(fsl2mni, func2std, vox2fsl)
+print(funcvox2mni)
```
-> Below we will use the
+> In the next section we will use the
> [`fsl.transform.flirt.fromFlirt`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.flirt.html#fsl.transform.flirt.fromFlirt)
> function, which does all of the above for us.
@@ -488,15 +518,14 @@ print('Peak activation (MNI voxels): ', mnivoxels)
```
-> Note that in the above example we are only applying a linear transformation
-> into MNI space - in reality you would also want to apply your non-linear
-> structural-to-standard transformation too. But this is left as [an exercise
-> for the
-> reader](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.fnirt.html).
+Note that in the above example we are only applying a linear transformation
+into MNI space - in reality you would also want to apply your non-linear
+structural-to-standard transformation too. This is covered in the next
+section.
-<a class="anchor" id="image-processing"></a>
-### Image processing
+<a class="anchor" id="transformations-and-resampling"></a>
+### Transformations and resampling
Now, it's all well and good to look at t-statistic values and voxel
@@ -505,33 +534,44 @@ at some images. Let's display our peak activation location in MNI space. To do
this, we're going to resample our functional image into MNI space, so we can
overlay it on the MNI template. This can be done using some handy functions
from the
+[`fsl.transform.flirt`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.flirt.html)
+and
[`fsl.utils.image.resample`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.utils.image.resample.html)
-module:
+modules.
-```
-from fsl.transform.flirt import fromFlirt
-from fsl.utils.image.resample import resampleToReference
+Let's make sure we've got our source and reference images loaded:
+
+```
featdir = op.join(op.join('08_fslpy', 'fmri.feat'))
tstat1 = Image(op.join(featdir, 'stats', 'tstat1'))
std = Image(op.expandvars(op.join('$FSLDIR', 'data', 'standard', 'MNI152_T1_2mm')))
+```
+
+
+Now we'll load the `example_func2standard` FLIRT matrix, and adjust it so that
+it transforms from functional *world* coordinates into standard *world*
+coordinates - this is what is expected by the `resampleToReference` function,
+used below:
+
+
+```
+from fsl.transform.flirt import fromFlirt
-# Load the func2standard FLIRT matrix, and adjust it
-# so that it transforms from functional *world*
-# coordinates into standard *world* coordinates -
-# this is what is expected by the resampleToReference
-# function, used below
func2std = np.loadtxt(op.join(featdir, 'reg', 'example_func2standard.mat'))
func2std = fromFlirt(func2std, tstat1, std, 'world', 'world')
+```
+
+
+Now we can use `resampleToReference` to resample our functional data into
+MNI152 space. This function returns a `numpy` array containing the resampled
+data, and an adjusted voxel-to-world affine transformation. But in this case,
+we know that the data will be aligned to MNI152, so we can ignore the affine:
+
+```
+from fsl.utils.image.resample import resampleToReference
-# All of the functions in the resample module
-# return a numpy array containing the resampled
-# data, and an adjusted voxel-to-world affine
-# transformation. But when using the
-# resampleToReference function, the affine will
-# be the same as the MNI152 2mm affine, so we
-# can ignore it.
std_tstat1 = resampleToReference(tstat1, std, func2std)[0]
std_tstat1 = Image(std_tstat1, header=std.header)
```
@@ -553,11 +593,137 @@ fig = ortho(std_tstat1, mnivoxels, cmap=plt.cm.inferno, fig=fig, cursor=True)
```
+In the example above, we resampled some data from functional space into
+standard space using a linear transformation. But we all know that this is not
+how things work in the real world - linear transformations are for kids. The
+real world is full of lions and tigers and bears and warp fields.
+
+
+The
+[`fsl.transform.fnirt`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.fnirt.html#fsl.transform.fnirt.fromFnirt)
+and
+[`fsl.transform.nonlinear`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.nonlinear.html)
+modules contain classes and functions for working with FNIRT-style warp fields
+(modules for working with lions, tigers, and bears are still under
+development).
+
+
+Let's imagine that we have defined an ROI in MNI152 space, and we want to
+project it into the space of our functional data. We can do this by combining
+the nonlinear structural to standard registration produced by FNIRT with the
+linear functional to structural registration generated by FLIRT. First of
+all, we'll load images from each of the functional, structural, and standard
+spaces:
+
+
+```
+featdir = op.join('08_fslpy', 'fmri.feat')
+func = Image(op.join(featdir, 'reg', 'example_func'))
+struc = Image(op.join(featdir, 'reg', 'highres'))
+std = Image(op.expandvars(op.join('$FSLDIR', 'data', 'standard', 'MNI152_T1_2mm')))
+```
+
+
+Now, let's say we have obtained our seed location in MNI152 coordinates. Let's
+convert them to MNI152 voxels just to double check:
+
+
+```
+seedmni = [-48, -74, -9]
+seedmnivox = transform(seedmni, std.getAffine('world', 'voxel'))
+ortho(std.data, seedmnivox, cursor=True)
+```
+
+
+Now we'll load the FNIRT warp field, which encodes a nonlinear transformation
+from structural space to standard space. FNIRT warp fields are often stored as
+*coefficient* fields to reduce the file size, but in order to use it, we must
+convert the coefficient field into a *deformation* (a.k.a. *displacement*)
+field. This takes a few seconds:
+
+
+```
+from fsl.transform.fnirt import readFnirt
+from fsl.transform.nonlinear import coefficientFieldToDeformationField
+
+struc2std = readFnirt(op.join(featdir, 'reg', 'highres2standard_warp'), struc, std)
+struc2std = coefficientFieldToDeformationField(struc2std)
+```
+
+We'll also load our FLIRT functional to structural transformation, adjust it
+so that it transforms between voxel coordinate systems instead of the FSL
+coordinate system, and invert so it can transform from structural voxels to
+functional voxels:
+
+
+```
+from fsl.transform.affine import invert
+func2struc = np.loadtxt(op.join(featdir, 'reg', 'example_func2highres.mat'))
+func2struc = fromFlirt(func2struc, func, struc, 'voxel', 'voxel')
+struc2func = invert(func2struc)
+```
+
+
+Now we can transform our seed coordinates from MNI152 space into functional
+space in two stages. First, we'll use our deformation field to transform from
+MNI152 space into structural space:
+
+
+```
+seedstruc = struc2std.transform(seedmni, 'world', 'voxel')
+seedfunc = transform(seedstruc, struc2func)
+
+print('Seed location in MNI coordinates: ', seedmni)
+print('Seed location in functional voxels:', seedfunc)
+ortho(func.data, seedfunc, cursor=True)
+```
+
+
+> FNIRT warp fields kind of work backwards - we can use them to transform
+> reference coordinates into source coordinates, but would need to invert the
+> warp field using `invwarp` if we wanted to transform from source coordinates
+> into referemce coordinates.
+
+
+Of course, we can also use our deformation field to resample an image from
+structural space into MNI152 space. The `applyDeformation` function takes an
+`Image` and a `DeformationField`, and returns a `numpy` array containing the
+resampled data.
+
+
+```
+from fsl.transform.nonlinear import applyDeformation
+
+strucmni = applyDeformation(struc, struc2std)
+
+# remove low-valued voxels,
+# just for visualisation below
+strucmni[strucmni < 1] = 0
+
+fig = ortho(std.data, [45, 54, 45], cmap=plt.cm.gray)
+fig = ortho(strucmni, [45, 54, 45], fig=fig)
+```
+
+
+The `premat` option to `applyDeformation` can be used to specify our linear
+functional to structural transformation, and hence resample a functional image
+into MNI152 space:
+
+
+```
+tstatmni = applyDeformation(tstat1, struc2std, premat=func2struc)
+tstatmni[tstatmni < 3] = 0
+
+fig = ortho(std.data, [45, 54, 45], cmap=plt.cm.gray)
+fig = ortho(tstatmni, [45, 54, 45], fig=fig)
+```
+
+
There are a few other useful functions tucked away in the
-[fsl.utils.image](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.utils.image.html)
-package, with more to be added in the future. The [`fsl.transform`]() package
-also contains a wealth of functionality for working with linear (FLIRT) and
-non-linear (FNIRT) transformations.
+[`fsl.utils.image`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.utils.image.html)
+and
+[`fsl.transform`](https://users.fmrib.ox.ac.uk/~paulmc/fsleyes/fslpy/latest/fsl.transform.html)
+packages, with more to be added in the future.
<a class="anchor" id="fsl-wrapper-functions"></a>
@@ -580,7 +746,7 @@ corresponding tool via the command-line:
```
-from fsl.wrappers import bet, robustfov, LOAD
+from fsl.wrappers import robustfov
robustfov('08_fslpy/bighead', 'bighead_cropped')
@@ -616,6 +782,8 @@ positional arguments, and pass the additional options as keyword arguments:
```
+from fsl.wrappers import bet
+
bet('bighead_cropped', 'bighead_cropped_brain', f=0.3, m=True, s=True)
render('bighead_cropped -b 40 '
@@ -673,6 +841,8 @@ files produced by the tool automatically loaded into memory for you:
```
+from fsl.wrappers import LOAD
+
cropped = Image('bighead_cropped')
# The loaded result is called "output",
diff --git a/advanced_topics/08_fslpy/fmri.feat/reg/example_func2highres.mat b/advanced_topics/08_fslpy/fmri.feat/reg/example_func2highres.mat
new file mode 100644
index 0000000000000000000000000000000000000000..dd61f138c86baeb47b3200251325c6fb2f661d69
--- /dev/null
+++ b/advanced_topics/08_fslpy/fmri.feat/reg/example_func2highres.mat
@@ -0,0 +1,4 @@
+0.99895165 0.0442250787 -0.01181694611 6.534548061
+-0.0439203422 0.9987243849 0.02491016319 9.692178016
+0.01290352651 -0.02436503913 0.9996197946 21.90296924
+0 0 0 1
diff --git a/advanced_topics/08_fslpy/fmri.feat/reg/highres.nii.gz b/advanced_topics/08_fslpy/fmri.feat/reg/highres.nii.gz
new file mode 100644
index 0000000000000000000000000000000000000000..0a4f8eb223ecb7605653ed088114dc622e0a810f
Binary files /dev/null and b/advanced_topics/08_fslpy/fmri.feat/reg/highres.nii.gz differ
diff --git a/advanced_topics/08_fslpy/fmri.feat/reg/highres2standard.mat b/advanced_topics/08_fslpy/fmri.feat/reg/highres2standard.mat
new file mode 100644
index 0000000000000000000000000000000000000000..033a907b8a7d549404141c72c2fa348f573da259
--- /dev/null
+++ b/advanced_topics/08_fslpy/fmri.feat/reg/highres2standard.mat
@@ -0,0 +1,4 @@
+1.056802026 -0.01924547726 0.02614687181 -36.51723948
+0.009055463297 0.9745460053 0.09056277052 -8.771603455
+-0.04315832679 -0.1680837227 1.136420957 -1.399839791
+0 0 0 1
diff --git a/advanced_topics/08_fslpy/fmri.feat/reg/highres2standard_warp.nii.gz b/advanced_topics/08_fslpy/fmri.feat/reg/highres2standard_warp.nii.gz
new file mode 100644
index 0000000000000000000000000000000000000000..df6f8601330538ded1fc12b1d93298b5b327626a
Binary files /dev/null and b/advanced_topics/08_fslpy/fmri.feat/reg/highres2standard_warp.nii.gz differ