{
 "metadata": {
  "name": ""
 },
 "nbformat": 3,
 "nbformat_minor": 0,
 "worksheets": [
  {
   "cells": [
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "# Issue 1: Infinite Values in the Final Z-Scores\n",
      "There are infinte values in the data when we don't separate the distribution into two:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "import matplotlib\n",
      "import matplotlib.pylab as plt\n",
      "import numpy as np\n",
      "%matplotlib inline\n",
      "import nibabel as nib\n",
      "from nilearn.plotting import plot_stat_map, plot_roi\n",
      "from scipy.spatial.distance import pdist\n",
      "from scipy.stats import norm, t\n",
      "import seaborn as sns\n",
      "\n",
      "all_copes_file = \"/home/vanessa/Desktop/tfMRI_LANGUAGE_STORY.nii_tstat1.nii.gz\"\n",
      "all_copes = nib.load(all_copes_file)\n",
      "\n",
      "data = all_copes.get_data()\n",
      "\n",
      "# Let's select just the nonzero voxels\n",
      "nonzero = data[data!=0]\n",
      "\n",
      "# Degrees of freedom is number of subjects -2\n",
      "dof = 486 -2\n",
      "Zdata = norm.ppf(t.cdf(data, dof))\n",
      "\n",
      "# Are there inf in the data?\n",
      "infs = np.isinf(Zdata).any()\n",
      "\n",
      "print \"Are there infs in the data?: %s\" %(infs)\n",
      "print \"There are %s infs in the data\" %(len(np.where(np.isinf(Zdata))[0]))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Are there infs in the data?: True\n",
        "There are 7982 infs in the data\n"
       ]
      }
     ],
     "prompt_number": 35
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "# Why are there infs in the data?\n",
      "This first issue can be traced to the case when the p_values are equal to 1.0"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "p_values = t.cdf(data, dof)\n",
      "Zdata = norm.ppf(p_values)\n",
      "index = np.where(np.isinf(Zdata))\n",
      "\n",
      "print \"The inf values occur where p-values are %s\" %(np.unique(p_values[index]))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "The inf values occur where p-values are [ 1.]\n"
       ]
      }
     ],
     "prompt_number": 36
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "# How can we fix this?\n",
      "This can be helped by fudging the number a bit: instead of 1.0 we can fill in 0.999999"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "p_values[p_values==1.0] = 0.999999999\n",
      "Zdata = norm.ppf(p_values)\n",
      "\n",
      "# Are there inf in the data?\n",
      "infs = np.isinf(Zdata).any()\n",
      "\n",
      "print \"Are there infs in the data?: %s\" %(infs)\n",
      "print \"There are %s infs in the data\" %(len(np.where(np.isinf(Zdata))[0]))"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [
      {
       "output_type": "stream",
       "stream": "stdout",
       "text": [
        "Are there infs in the data?: False\n",
        "There are 0 infs in the data\n"
       ]
      }
     ],
     "prompt_number": 37
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "# Issue 2: The Z-Score Distribution is Truncated\n",
      "Here I thought we were in the clear! But then I looked at the distribution of the map:"
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [
      "Z_nii = nib.nifti1.Nifti1Image(Zdata,affine=all_copes.get_affine(),header=all_copes.get_header())\n",
      "nib.save(Z_nii,\"/home/vanessa/Desktop/Zdata.nii\")"
     ],
     "language": "python",
     "metadata": {},
     "outputs": [],
     "prompt_number": 45
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "On the left is the Zdata that we just made - if it were a lizard, someone cut part of his right tail off! The map on the right is generated by first splitting at zero.\n",
      "\n",
      "\n",
      "![zmap_comparison](http://www.vbmis.com/bmi/share/chris/t_to_z/comparison.png)"
     ]
    },
    {
     "cell_type": "markdown",
     "metadata": {},
     "source": [
      "Apologies for the different scalings for the two distributions, I just opened them up in MriCron :)\n",
      "\n",
      "The concern is that we are losing potentially the strongest activations in the conversion."
     ]
    },
    {
     "cell_type": "code",
     "collapsed": false,
     "input": [],
     "language": "python",
     "metadata": {},
     "outputs": []
    }
   ],
   "metadata": {}
  }
 ]
}