Neuroimaging in Python-Empty.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Import standard modules"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "from matplotlib import pyplot as plt\n",
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Import NiPy modules\n",
    "\n",
    "If you've installed python via miniconda, most packages are available by running `conda install package` in the terminal. In general, if a package is available through conda that's the easiest way to install it. Sometimes packages aren't available through conda, in which case you can use the python package installer, `pip`. If you're running this tutorial on `syzygy.ca`, conda is not installed so all packages are installed through pip. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "#!pip install --user nipy nibabel dipy nilearn nipype"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import nipy\n",
    "import nibabel as nib\n",
    "import dipy\n",
    "import nilearn\n",
    "import nipype"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Load up some sample data\n",
    "\n",
    "First, we'll load up the MNI T1 template as \"t1img\". `nib.load()` returns an image object. The exact structure of the image object will depend on what type of file you load. Mostly you're\n",
    "\n",
    "Read more about nibabel image objects here: http://nipy.org/nibabel/nibabel_images.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Loading the image object is fast because it doesn't load the image data automatically. This lets you quickly access header information if that's all you want. To load the image data as a numpy array just call `get_data()`"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Visualizing a slice is easy"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Note: Choose colormaps carefully! \n",
    "\n",
    "https://matplotlib.org/users/colormaps.html"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### t1data and maskdata are both just numpy 3D arrays, so we perform normal array operations"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise (optional): Parse an atlas xml file\n",
    "\n",
    "* Parse the xml file using xml.etree.ElementTree\n",
    "* Determine how to access the index and text of each element\n",
    "* Loop over all the elements and place the data into a dictionary, with the label text as key and the index as value\n",
    "* Plot an ROI by calling it by it's string, ie 'Body of corpus callosum'\n",
    "* Evaluate the mean T1 in that ROI, again calling it by name"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Saving modified images as Nifti files\n",
    "Let's make a copy of t1data and zero out voxels where maskdata==48"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Save the modified image\n",
    "\n",
    "Create a new nifti image object with `nib.Nifti1Image()`. It needs a image array and an affine matrix"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Nibabel has tools for opening multiple formats including PAR/REC\n",
    "Loading a PAR/REC image with nibabel and then saving to nifti is, in my experience, more robust than dcm2nii"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise (optional): Make a  pandas dataframe with the mean values for each ROI\n",
    "* Calculate the mean of each roi in the JHU atlas. Put in in a list or a dict\n",
    "* Create a new `pd.DataFrame()` with the list/dict as the data, the mask labels as the index, and a relevant string as the column header "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Not very \"pythonic\". Avoid loops when you can do it in one line!"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Do this for multiple subjects and you've started building a database that you can run stats on, save to Excel, etc."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Manipulating the affine\n",
    "The affine matrix describes the relationship between the image matrix and the \"real world\" coordinate system. Image orientation manipulation should be approached with caution, but it's often necessary."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Notice the difference?\n",
    "\n",
    "Hard to spot, but the canonical oriented image is flipped left/right from the standard image. If you have a sagitally acquired image, the image will be reoriented in more dimensions."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "## Exercise: Rotate the image and the affine matrix by 90 degrees. Save and open in FSL. Are the orientation labels correct?\n",
    "\n",
    "http://nipy.org/nibabel/coordinate_systems.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise: Tissue segmentation\n",
    "\n",
    "http://nipy.org/dipy/examples_built/tissue_classification.html#example-tissue-classification\n",
    "\n",
    "Let's run a tissue classification algorithm on our T1 image"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Nipype: Integrate imaging software into your Python code\n",
    "\n",
    "https://miykael.github.io/nipype_tutorial/\n",
    "\n",
    "Nipype is an amazing utility that allows you to link together different non-Python imaging tools within your Python analysis. It can be used simply to (for example) add FSL's eddy correction processing to you Python DTI processing; or it can be used to organize your entire processing workflow.\n",
    "\n",
    "Personally I don't use the full workflow functionality. Instead I use it to call FSL and SPM functions from my Python scripts. It's much easier to write everything in pure Python instead of switching back and forth between FSL command line tools and SPM in Matlab. Nipype has an amazing number of interfaces, including:\n",
    "* FSL\n",
    "* SPM\n",
    "* FreeSurfer\n",
    "* niftireg\n",
    "* ANTS\n",
    "* AFNI\n",
    "* Camino\n",
    "* MIPAV\n",
    "* DTK\n",
    "\n",
    "Full list here: http://nipype.readthedocs.io/en/latest/documentation.html\n",
    "\n",
    "For the Nipype interfaces to work, you just have to have the software you want to use install and certain environmental variables set so that Python knows where to look for the programs. For example,"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "None of these software tools are on the syzygy server. This example will only work on a computer with FSL installed."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`my_reg_file` is a string describing a path to a file. We can now use this variable as the input to another nipype interface, or we can load the nifti file into Python with the tools we learned previously."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise: DTI processing\n",
    "http://nipy.org/dipy/examples_built/reconst_dti.html#example-reconst-dti"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "import dipy.reconst.dti as dti\n",
    "from dipy.data import read_stanford_hardi"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The gtab (**g**radient **tab**le) is the object that holds the gradient information. This example gives it to us automatically. If you created your dti nifti file with dcm2nii, it probably also gave you  .bvec and .bval files. In that case, you could create the gtab with\n",
    "```\n",
    "from dipy.io import read_bvals_bvecs\n",
    "bvals, bvecs = read_bvals_bvecs(fbval, fbvec)\n",
    "from dipy.core.gradients import gradient_table\n",
    "gtab = gradient_table(bvals, bvecs)\n",
    "```\n",
    "\n",
    "If you're reading data straight from a PAR/REC file, \n",
    "```\n",
    "dtipar = nib.parrec.load(dti_par)\n",
    "bvals, bvecs = dtipar.header.get_bvals_bvecs()\n",
    "gtab = gradient_table(bvals, bvecs)\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create a TensorModel object with our gradient information"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That's it! We can save the FA as a nifti, or do some numerical analysis on it with numpy, or make publication images, or whatever. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise: Add an eddy correction step to the above pipeline using a Nipype interface"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise: fMRI analysis with ICA\n",
    "\n",
    "http://nilearn.github.io/connectivity/resting_state_networks.html"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "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.6.0"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}