mri_plot.ipynb

{
  "metadata": {
    "name": "Basic (f)MRI plotting"
  }, 
  "nbformat": 3, 
  "nbformat_minor": 0, 
  "worksheets": [
    {
      "cells": [
        {
          "cell_type": "heading", 
          "level": 1, 
          "metadata": {}, 
          "source": [
            "Basic (f)MRI plotting"
          ]
        }, 
        {
          "cell_type": "markdown", 
          "metadata": {}, 
          "source": [
            "\n\n", 
            "When running an fMRI data analysis it is often necessary to visualize results\nin their original dataspace, typically as an overlay on some anatomical brain\nimage. PyMVPA has the ability to export results into the NIfTI format, and via\nthis data format it is compatible with virtually any MRI visualization software.\n\n", 
            "However, sometimes having a scriptable plotting facility within Python is\ndesired. There are a number of candidate tools for this purpose (e.g. ", 
            "[Mayavi](http://code.enthought.com/projects/mayavi/)),\nbut also PyMVPA itself offers some basic MRI plotting.\n\n", 
            "In this example, we are showing a quick-and-dirty plot of a voxel-wise\nANOVA measure, overlaid on the respective brain anatomy. Note that the plotting\nis not specific to ANOVAs. Any feature-wise measure can be plotted this way.\n\n", 
            "We start with basic steps: loading PyMVPA and the example fMRI dataset, only\nselect voxels that correspond to some pre-computed gray matter mask, do basic\npreprocessing, and estimate ANOVA scores. This has already been described\nelsewhere, hence we only provide the code here for the sake of completeness."
          ]
        }, 
        {
          "cell_type": "code", 
          "collapsed": false, 
          "input": [
            "from mvpa2.suite import *\n\n# load PyMVPA example dataset\ndatapath = os.path.join(mvpa2.cfg.get('location', 'tutorial data'), 'haxby2001')\ndataset = load_tutorial_data(roi='gray')\n\n# do chunkswise linear detrending on dataset\npoly_detrend(dataset, chunks_attr='chunks')\n\n# exclude the rest conditions from the dataset, since that should be\n# quite different from the 'active' conditions, and make the computation\n# below pointless\ndataset = dataset[dataset.sa.targets != 'rest']\n\n# define sensitivity analyzer to compute ANOVA F-scores on the remaining\n# samples\nsensana = OneWayAnova()\nsens = sensana(dataset)"
          ], 
          "language": "python", 
          "metadata": {}, 
          "outputs": []
        }, 
        {
          "cell_type": "markdown", 
          "metadata": {}, 
          "source": [
            "\n\n", 
            "The measure is computed, and we can look at the actual plotting. Typically, it\nis useful to pre-define some common plotting arguments, for example to ensure\nconsistency throughout multiple figures. This following sets up which backround\nimage to use (`background`), which portions of the image to plot\n(`background_mask`), and which portions of the overlay images to plot\n(`overlay_mask`). All these arguments are actually NIfTI images of the same\ndimensions and orientation as the to be plotted F-scores image. the remaining\nsettings configure the colormaps to be used for plotting and trigger\ninteractive plotting."
          ]
        }, 
        {
          "cell_type": "code", 
          "collapsed": false, 
          "input": [
            "mri_args = {\n    'background' : os.path.join(datapath, 'sub001', 'anatomy', 'highres001.nii.gz'),\n    'background_mask' : os.path.join(datapath, 'sub001', 'masks', 'orig', 'brain.nii.gz'),\n    'overlay_mask' : os.path.join(datapath, 'sub001', 'masks', 'orig', 'gray.nii.gz'),\n    'cmap_bg' : 'gray',\n    'cmap_overlay' : 'autumn', # YlOrRd_r # pl.cm.autumn\n    'interactive' : cfg.getboolean('examples', 'interactive', True),\n    }"
          ], 
          "language": "python", 
          "metadata": {}, 
          "outputs": []
        }, 
        {
          "cell_type": "markdown", 
          "metadata": {}, 
          "source": [
            "\n\n", 
            "All that remains to do is a single call to ", 
            "`plot_lightbox()`. We pass it the\nF-score vector. ", 
            "`map2nifti` uses the mapper in our original dataset to project\nit back into the functional MRI volume space. We treshold the data with the\ninterval [0, +inf] (i.e. all possible values and F-Score can have), and select\na subset of slices to be plotted. That's it."
          ]
        }, 
        {
          "cell_type": "code", 
          "collapsed": false, 
          "input": [
            "fig = plot_lightbox(overlay=map2nifti(dataset, sens),\n              vlim=(0, None), slices=range(25,29), **mri_args)"
          ], 
          "language": "python", 
          "metadata": {}, 
          "outputs": []
        }, 
        {
          "cell_type": "markdown", 
          "metadata": {}, 
          "source": [
            "\n\n", 
            "The resulting figure would look like this:\n\n\\[Visit [http://pymvpa.org/examples/mri_plot.html](http://pymvpa.org/examples/mri_plot.html) to view this figure\\]\n\n", 
            "In interactive mode it is possible to click on the histogram to adjust the\nthresholding of the overlay volumes. Left-click sets the value corresponding\nto the lowest value in the colormap, and right-click set the value for the upper\nend of the colormap. Try right-clicking somewhere at the beginning of the x-axis\nand left on the end of the x-axis."
          ]
        }
      ], 
      "metadata": {}
    }
  ]
}