{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Multi-subject dictionary learning & CanICA"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Reference:**\n",
    "\n",
    "  * [Dictionary Learning and ICA for doing group analysis of resting-state fMRI](http://nilearn.github.io/auto_examples/03_connectivity/plot_compare_resting_state_decomposition.html)\n",
    "  * [nilearn.decomposition.DictLearning](http://nilearn.github.io/modules/generated/nilearn.decomposition.DictLearning.html)\n",
    "  * [nilearn.decomposition.CanICA](http://nilearn.github.io/modules/generated/nilearn.decomposition.CanICA.html)\n",
    "  * [CanICA](https://github.com/GaelVaroquaux/canica)\n",
    "\n",
    "This example applies dictionary learning and ICA to resting-state data, visualizing resulting components using atlas plotting tools.\n",
    "\n",
    "Dictionary learning is a sparsity based decomposition method for extracting spatial maps. It extracts maps that are naturally sparse and usually cleaner than ICA.\n",
    "CanICA is an ICA method for group-level analysis of fMRI data. Compared to other strategies, it brings a well-controlled group model, as well as a thresholding algorithm controlling for specificity and sensitivity with an explicit model of the signal. \n",
    "\n",
    "  * G. Varoquaux et al. “A group model for stable multi-subject ICA on fMRI datasets”, NeuroImage Vol 51 (2010), p. 288-299\n",
    "  * G. Varoquaux et al. “ICA-based sparse features recovery from fMRI datasets”, IEEE ISBI 2010, p. 1177\n",
    "  * Gael Varoquaux et al. Multi-subject dictionary learning to segment an atlas of brain spontaneous activity Information Processing in Medical Imaging, 2011, pp. 562-573, Lecture Notes in Computer Science"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/datagrid/personal/borovec/Applications/vEnv2/lib/python2.7/site-packages/h5py/__init__.py:36: FutureWarning: Conversion of the second argument of issubdtype from `float` to `np.floating` is deprecated. In future, it will be treated as `np.float64 == np.dtype(float).type`.\n",
      "  from ._conv import register_converters as _register_converters\n",
      "/mnt/datagrid/personal/borovec/Applications/vEnv2/lib/python2.7/site-packages/IPython/html.py:14: ShimWarning: The `IPython.html` package has been deprecated. You should import from `notebook` instead. `IPython.html.widgets` has moved to `ipywidgets`.\n",
      "  \"`IPython.html.widgets` has moved to `ipywidgets`.\", ShimWarning)\n"
     ]
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "%load_ext autoreload\n",
    "%autoreload 2\n",
    "import os, sys\n",
    "import numpy as np\n",
    "import nibabel as nib\n",
    "from skimage import io\n",
    "import matplotlib.pylab as plt\n",
    "sys.path += [os.path.abspath('.'), os.path.abspath('..')]  # Add path to root\n",
    "import notebooks.notebook_utils as uts"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## load datset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "loading dataset: (True) exists -> /datagrid/Medical/microscopy/drosophila/synthetic_data/atomicPatternDictionary_v0/datasetFuzzy_raw\n",
      "loading (True) <- /datagrid/Medical/microscopy/drosophila/synthetic_data/atomicPatternDictionary_v0/dictionary/atlas.png\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAEICAYAAAB/KknhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAEjFJREFUeJzt3X2QXXV9x/H3xxASnkII0BASxmBBMHUk2JWHkaISeSyV/MEwICNrJ87OtLaFaiugo1Mcx0E7I9rRod0CElsUEMVgxgdipKPTSmAj4TEiUYIkTVgQ0gCtIcFv/zi/rZftPpzdex528/u8Zu7sPeeee883ufu5v9/v3LO/o4jAzPLyurYLMLPmOfhmGXLwzTLk4JtlyME3y5CDb5YhB3+Kk7RZ0rvbrmM0kj4q6Ya267CJcfBbIOnfJL0gadaw9TdL+lRbdU1GRHw6Ij7Q9H4lfVfSSx23VyQ93HQd05WD3zBJi4E/AgJ4T6vFTGMRcW5EHDh0A/4D+HrbdU0XDn7zLgPuBW4GeodWSuoDLgU+klqwbw9/oqSTJP1E0g5J2yR9UdK+6TFJuk7SoKSdkh6W9Ob02HslPTRaQcOHE5L+TtK/pvuLJYWkXkm/kvScpI+NtG1afp+kpyT9WtLHOl97eI9G0jslbelYPlLSNyQ9K+lJSX9V5j+048P0K2W2Nwe/DZcBt6Tb2ZLmA0REf1r32dSK/ckIz30V+GvgMOBUYBnw5+mxs4DTgTcCBwMXAb9Or/3ViHhLl3WfBhyX9vkJSW8avoGkJcD1wPuAI4FDgUVlXlzS64BvAw8CC9N+rpB0domnXwb8OCI2l9mXOfiNknQa8Hrg9ohYD/wCeG/Z50fE+oi4NyL2pF/yfwLekR7eDRwEHA8oIjZGxLYKy78mIv4nIh6kCOcJI2xzIbA6In4UEbuAjwO/Lfn6bwMOj4hPRsQrEfFL4J+Bi0s89zKKHpSV5OA3qxe4OyKeS8tfpaO7Px5Jb5S0WtJ2STuBT1O0/kTED4EvAl8CBiX1S5pTYe3bO+7/N3DgCNscCTw9tBARL5N6HSW8HjgyDWN2SNoBfBSYP9aT0ofpEcAdJfdjOPiNkbQfRff7HSm42ym67SdIGmo9x/tTyeuBnwHHRsQcimBo6MGI+IeI+ENgCUWX/29LlvcysH/H8hElnzfcNuCooQVJ+1N098vs52ngyYiY23E7KCLOG2efvcA3I+KlSdacJQe/OcspxuhLgKXp9ibgxxRdVYBngDeM8RoHATuBlyQdD/zZ0AOS3ibpZEkzKQL2G1I3W9L7JW0e43U3ABdLmimph6LLPhl3AOdLOi0ddPwkr/0d2wCcJ2mepCOAKzoeuw94UdKVkvaTNEPSmyW9bbSddXyY3jzJerPl4DenF/hyRPwqIrYP3Si655dK2ge4EViSurrfGuE1/obimMCLFOPf2zoem5PWvQA8RdHF/vv02FHAv49R28eB30/PvYZiCDJhEfEo8MH0/G3p9bZ0bPIvFMcHNgN3d9YfEa8C51N8ID4JPAfcQHGgcjTLgR3APZOpN2fyRBx7P0l3A5dHxMYW9r0Z+EBE/KDpfdvo9mm7AKtfRJzVdg02tbirb5Yhd/XNMtRViy/pHEmPS9ok6aqqijKzek26xZc0A/g5cCbFkdv7gUsi4rHRnrOvZsVsDpjU/sxsfL/hZV6JXRpvu24O7p0EbEqnViLpVuACYNTgz+YATtayLnZpZmNZF2tLbddNV38hHadnUrT6C4dvJKlP0oCkgd3s6mJ3ZlaV2o/qR0R/RPRERM9MZo3/BDOrXTfB30rHedkUf365tbtyzKwJ3QT/fuBYSUen87IvBu6qpiwzq9OkD+5FxB5JfwF8H5gB3JTO1TazKa6rU3Yj4jvAdyqqxcwa4lN2zTLk4JtlyME3y5CDb5YhB98sQw6+WYYcfLMMOfhmGXLwzTKU5WSb//WdY0ptd/B5m2quxKwdbvHNMuTgm2Uom65+2e79WM9x19/2Fm7xzTLk4JtlyME3y5CDb5YhB98sQw6+WYYcfLMMOfhmGXLwzTLk4JtlKJtTdi0P5z664zXL3/2DuS1VMrW5xTfL0LjBl3STpEFJj3SsmydpjaQn0s9D6i3TzKqkiBh7A+l04CXgKxHx5rTus8DzEXGtpKuAQyLiyvF2Nkfz4mQtq6DsiZvMX+eNxX+pV7/h3fZu5dDtXxdr2RnPa7ztxm3xI+JHwPPDVl8ArEz3VwLLJ1yhmbVmsgf35kfEtnR/OzB/tA0l9QF9ALPZf5K7M7MqdX1UPyJC0qjjhYjoB/qh6Op3u782uXtfr6q79ja6yR7Vf0bSAoD0c7C6ksysbpMN/l1Ab7rfC6yqphwza0KZr/O+BvwEOE7SFkkrgGuBMyU9Abw7LZvZNDHuGD8iLhnloXa+l7O9SpPj+s595fDV3lh85p5Zhhx8swxl80c6nV/F+RJa7ZkqX9nl/sc8bvHNMuTgm2XIwTfLUDZj/E5Vj92f6zu10tebiMP6f9LavsuaKuN6+x23+GYZcvDNMjTuRBxVanMijiq02aWfjKk4DJgO3f7p/NVeZRNxmNnex8E3y1CWR/XLmm5d++HGqn8yw4Dv/+eG0tuefeTSCb++NcctvlmGHHyzDDn4ZhnyGH8Mw8fB03nMP9mv9iYyrh/teR7vTz1u8c0y5OCbZchd/TFM5679cMP/LVPxrD5rjlt8sww5+GYZcvDNMuQxfiY8ph/bdP6LvMlwi2+WoTKX0DpK0j2SHpP0qKTL0/p5ktZIeiL9PKT+cs2sCmW6+nuAD0fETyUdBKyXtAZ4P7A2Iq6VdBVwFXBlfaXaRLXZve88W286TL6Rm3Fb/IjYFhE/TfdfBDYCC4ELgJVps5XA8rqKNLNqTejgnqTFwInAOmB+RGxLD20H5o/ynD6gD2A2+0+2TjOrUOmDe5IOBL4BXBEROzsfi2LivhEn74uI/ojoiYiemczqqlgzq0apFl/STIrQ3xIR30yrn5G0ICK2SVoADNZVZFv2pr/Oa9pof53n8f7UUOaovoAbgY0R8bmOh+4CetP9XmBV9eWZWR3KtPhvB94HPCxp6GP8o8C1wO2SVgBPARfVU6KZVc3z6lcgl0to1THZZltd/731TD3Pq29mo3LwzTLkrn6Lqp73frqru9u/t3bvO7mrb2ajcvDNMuTgm2XIE3G0KMdx/FiGj8HLjvlzGLtXzS2+WYYcfLMMuatvU5a78PVxi2+WIQffLEMOvlmGHHyzDDn4Zhly8M0y5OCbZcjBN8uQg2+WIQffLEMOvlmGHHyzDDn4Zhly8M0y5OCbZajMtfNmS7pP0oOSHpV0TVp/tKR1kjZJuk3SvvWXa2ZVKNPi7wLOiIgTgKXAOZJOAT4DXBcRxwAvACvqK9PMqjRu8KPwUlqcmW4BnAHckdavBJbXUqGZVa7UGF/SjHSl3EFgDfALYEdE7EmbbAEWjvLcPkkDkgZ2s6uKms2sS6WCHxGvRsRSYBFwEnB82R1ERH9E9EREz0xmTbJMM6vShI7qR8QO4B7gVGCupKHJOhcBWyuuzcxqUuao/uGS5qb7+wFnAhspPgAuTJv1AqvqKtLMqlVmeu0FwEpJMyg+KG6PiNWSHgNulfQp4AHgxhrrNLMKjRv8iHgIOHGE9b+kGO+b2TTjM/fMMuTgm2XIwTfLkINvliEH3yxDDr5ZhnyZ7Dqc8pbf3b/3ofbqMBuFW3yzDDn4ZhnKs6vf2RUfrrNrPtZ2k92Xu/42BbjFN8uQg2+WIQffLEP5jPHLjterGNeXfX2P960lbvHNMuTgm2Vo7+3q191lr4K/6rOWuMU3y5CDb5YhB98sQw6+WYYcfLMMOfhmGXLwzTLk4JtlqHTw06WyH5C0Oi0fLWmdpE2SbpO0b31lmlmVJnLm3uUUF8uck5Y/A1wXEbdK+kdgBXB9xfVN3vCz4KbimXw+U89aUqrFl7QI+GPghrQs4AzgjrTJSmB5HQWaWfXKdvU/D3wE+G1aPhTYERF70vIWYOFIT5TUJ2lA0sBudnVVrJlVY9zgSzofGIyI9ZPZQUT0R0RPRPTMZNZkXsLMKlZmjP924D2SzgNmU4zxvwDMlbRPavUXAVvrK7MCVUyiWfVrNGjT504pve0xH7q3xkpsKhi3xY+IqyNiUUQsBi4GfhgRlwL3ABemzXqBVbVVaWaV6uZ7/CuBD0naRDHmv7GaksysboqIxnY2R/PiZC1rbH+5m0j3fjTu9k8v62ItO+N5jbedz9wzy5CDb5ahvXfOvUxV0b0f7fXc7d97uMU3y5CDb5YhB98sQw6+WYYcfLMMOfhmGXLwzTLk4JtlyME3y5CDb5YhB98sQw6+WYYcfLMMOfhmGXLwzTLk4JtlyME3y5CDb5YhB98sQw6+WYYcfLMMOfhmGSo1vbakzcCLwKvAnojokTQPuA1YDGwGLoqIF+op08yqNJEW/10RsTQietLyVcDaiDgWWJuWzWwa6KarfwGwMt1fCSzvvhwza0LZ4Adwt6T1kvrSuvkRsS3d3w7MH+mJkvokDUga2M2uLss1syqUvYTWaRGxVdLvAWsk/azzwYgISSNedjci+oF+KK6W21W1ZlaJUi1+RGxNPweBO4GTgGckLQBIPwfrKtLMqjVu8CUdIOmgofvAWcAjwF1Ab9qsF1hVV5FmVq0yXf35wJ2Shrb/akR8T9L9wO2SVgBPARfVV6aZVUkRzQ2752henKxlje0vd1VcMtuXxp5e1sVadsbzGm87n7lnliEH3yxDDr5Zhsp+j2/TUOf4fCLjfY/r935u8c0y5OCbZchd/Uy4+26d3OKbZcjBN8uQg2+WIQffLEMOvlmGHHyzDDn4Zhly8M0y5OCbZcjBN8uQg2+WIQffLEMOvlmGHHyzDDn4Zhly8M0y5OCbZcjBN8tQqeBLmivpDkk/k7RR0qmS5klaI+mJ9POQuos1s2qUbfG/AHwvIo4HTgA2AlcBayPiWGBtWjazaaDM1XIPBk4HbgSIiFciYgdwAbAybbYSWF5XkWZWrTIt/tHAs8CXJT0g6YZ0uez5EbEtbbOd4qq6/4+kPkkDkgZ2s6uaqs2sK2WCvw/wVuD6iDgReJlh3fooLrk74mV3I6I/Inoiomcms7qt18wqUCb4W4AtEbEuLd9B8UHwjKQFAOnnYD0lmlnVxg1+RGwHnpZ0XFq1DHgMuAvoTet6gVW1VGhmlSt7JZ2/BG6RtC/wS+BPKT40bpe0AngKuKieEs2saqWCHxEbgJ4RHlpWbTlm1gSfuWeWIQffLEMOvlmGHHyzDDn4Zhly8M0ypOJs24Z2Jj1L8Z3/YcBzje14ZFOhBnAdw7mO15poHa+PiMPH26jR4P/fTqWBiBjpvICsanAdrqOtOtzVN8uQg2+WobaC39/SfjtNhRrAdQznOl6rljpaGeObWbvc1TfLkINvlqFGgy/pHEmPS9okqbFZeSXdJGlQ0iMd6xqfHlzSUZLukfSYpEclXd5GLZJmS7pP0oOpjmvS+qMlrUvvz21p/oXaSZqR5nNc3VYdkjZLeljSBkkDaV0bvyONTGXfWPAlzQC+BJwLLAEukbSkod3fDJwzbF0b04PvAT4cEUuAU4APpv+DpmvZBZwREScAS4FzJJ0CfAa4LiKOAV4AVtRcx5DLKaZsH9JWHe+KiKUd35u38TvSzFT2EdHIDTgV+H7H8tXA1Q3ufzHwSMfy48CCdH8B8HhTtXTUsAo4s81agP2BnwInU5whts9I71eN+1+UfpnPAFYDaqmOzcBhw9Y1+r4ABwNPkg6611lHk139hcDTHctb0rq2lJoevC6SFgMnAuvaqCV1rzdQTJK6BvgFsCMi9qRNmnp/Pg98BPhtWj60pToCuFvSekl9aV3T70tXU9lPhA/uMfb04HWQdCDwDeCKiNjZRi0R8WpELKVocU8Cjq97n8NJOh8YjIj1Te97BKdFxFsphqIflHR654MNvS9dTWU/EU0GfytwVMfyorSuLa1MDy5pJkXob4mIb7ZZC0AUV0W6h6JLPVfS0DyMTbw/bwfeI2kzcCtFd/8LLdRBRGxNPweBOyk+DJt+Xxqbyr7J4N8PHJuO2O4LXEwxRXdbGp8eXJIoLkW2MSI+11Ytkg6XNDfd34/iOMNGig+AC5uqIyKujohFEbGY4vfhhxFxadN1SDpA0kFD94GzgEdo+H2JJqeyr/ugybCDFOcBP6cYT36swf1+DdgG7Kb4VF1BMZZcCzwB/ACY10Adp1F00x4CNqTbeU3XArwFeCDV8QjwibT+DcB9wCbg68CsBt+jdwKr26gj7e/BdHt06Hezpd+RpcBAem++BRxSRx0+ZdcsQz64Z5YhB98sQw6+WYYcfLMMOfhmGXLwzTLk4Jtl6H8BEGZDb3ifE5cAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f8634425e90>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# uts.DEFAULT_PATH = '/datagrid/Medical/microscopy/drosophila/synthetic_data/atomicPatternDictionary_v0'\n",
    "p_dataset = os.path.join(uts.DEFAULT_PATH, uts.SYNTH_DATASETS_FUZZY[0])\n",
    "print ('loading dataset: ({}) exists -> {}'.format(os.path.exists(p_dataset), p_dataset))\n",
    "\n",
    "p_atlas = os.path.join(uts.DEFAULT_PATH, 'dictionary/atlas.png')\n",
    "atlas_gt = io.imread(p_atlas)\n",
    "nb_patterns = len(np.unique(atlas_gt))\n",
    "print ('loading ({}) <- {}'.format(os.path.exists(p_atlas), p_atlas))\n",
    "plt.imshow(atlas_gt, interpolation='nearest')\n",
    "_ = plt.title('Atlas; unique %i' % nb_patterns)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "('loaded # images: ', 800)\n",
      "('image shape:', (64, 64))\n"
     ]
    }
   ],
   "source": [
    "list_imgs = uts.load_dataset(p_dataset)\n",
    "print ('loaded # images: ', len(list_imgs))\n",
    "img_shape = list_imgs[0].shape\n",
    "print ('image shape:', img_shape)\n",
    "\n",
    "nii_images = [nib.Nifti1Image(np.expand_dims(img, axis=0), affine=np.eye(4)) for img in list_imgs]\n",
    "mask_full = nib.Nifti1Image(np.expand_dims(np.ones(atlas_gt.shape, dtype=np.int8), axis=0), affine=np.eye(4))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "('shape:', (1, 64, 64))\n",
      "('values:', array([1], dtype=int8))\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAD8CAYAAABXXhlaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAADLVJREFUeJzt3X/oXfV9x/Hna/lZ7Y+Y1oVgZHEYKv4xY/niD5SymlmyrtT8IaKUEUYg/7hhWaHTDQaF/VH/qfWPMQjVNX+4qrN1ESltXWoZgxH9WrWNptbUKSZE021Ku8LSxL73xz0pX8M3fm9yz7m34fN8wJfvPeee63nj/T7vr1zOSVUhqS2/M+sBJE2f4UsNMnypQYYvNcjwpQYZvtQgw5caNFH4SbYmeSnJwSR39jWUpGHlbL/Ak2QZ8BPgRuAQ8DRwW1W92N94koawfILbXgUcrKpXAJI8CNwEnDb8lVlVqzl/gl1Kei//xy/5VR3LUttNEv5FwOsLlg8BV7/XDVZzPldnywS7lPRe9tXesbabJPyxJNkJ7ARYzXlD707SGCb5cO8wcPGC5Q3dunepql1VNVdVcytYNcHuJPVlkvCfBjYluSTJSuBW4LF+xpI0pLN+qV9VJ5L8OfAdYBlwf1W90NtkkgYz0Xv8qvoW8K2eZpE0JX5zT2qQ4UsNMnypQYYvNcjwpQYZvtQgw5caZPhSgwxfapDhSw0yfKlBhi81yPClBhm+1CDDlxpk+FKDDF9qkOFLDTJ8qUGGLzXI8KUGGb7UIMOXGmT4UoMMX2rQkuEnuT/J0ST7F6xbm+SJJC93vy8YdkxJfRrnGf9rwNZT1t0J7K2qTcDeblnSOWLJ8Kvq34D/OWX1TcDu7vJuYFvPc0ka0Nm+x19XVUe6y28A63qaR9IUTPzhXlUVUKe7PsnOJPNJ5o9zbNLdSerB2Yb/ZpL1AN3vo6fbsKp2VdVcVc2tYNVZ7k5Sn842/MeA7d3l7cCefsaRNA3j/HPe14H/AD6a5FCSHcCXgBuTvAz8Ubcs6RyxfKkNquq201y1pedZJE2J39yTGmT4UoMMX2qQ4UsNMnypQYYvNcjwpQYZvtQgw5caZPhSgwxfapDhSw0yfKlBhi81yPClBhm+1CDDlxpk+FKDDF9qkOFLDTJ8qUGGLzXI8KUGGb7UIMOXGjTOKbQuTvJkkheTvJDkjm792iRPJHm5+33B8ONK6sM4z/gngM9X1eXANcDtSS4H7gT2VtUmYG+3LOkcsGT4VXWkqn7QXf4FcAC4CLgJ2N1tthvYNtSQkvp1Ru/xk2wErgT2Aeuq6kh31RvAul4nkzSYscNP8n7gG8DnqurnC6+rqgLqNLfbmWQ+yfxxjk00rKR+jBV+khWMon+gqr7ZrX4zyfru+vXA0cVuW1W7qmququZWsKqPmSVNaJxP9QPcBxyoqi8vuOoxYHt3eTuwp//xJA1h+RjbXAf8KfCjJM916/4a+BLwcJIdwGvALcOMKKlvS4ZfVf8O5DRXb+l3HEnT4Df3pAYZvtQgw5caZPhSgwxfapDhSw0yfKlBhi81yPClBhm+1CDDlxpk+FKDDF9qkOFLDTJ8qUGGLzXI8KUGGb7UIMOXGmT4UoMMX2qQ4UsNMnypQYYvNcjwpQaNc+681UmeSvJ8kheSfLFbf0mSfUkOJnkoycrhx5XUh3Ge8Y8BN1TVFcBmYGuSa4C7gXuq6lLgLWDHcGNK6tOS4dfI/3aLK7qfAm4AHunW7wa2DTKhpN6N9R4/ybLuTLlHgSeAnwJvV9WJbpNDwEXDjCipb2OFX1XvVNVmYANwFXDZuDtIsjPJfJL54xw7yzEl9emMPtWvqreBJ4FrgTVJTp5mewNw+DS32VVVc1U1t4JVEw0rqR/jfKp/YZI13eX3ATcCBxg9ANzcbbYd2DPUkJL6tXzpTVgP7E6yjNEDxcNV9XiSF4EHk/wd8Cxw34BzSurRkuFX1Q+BKxdZ/wqj9/uSzjF+c09qkOFLDTJ8qUGGLzXI8KUGGb7UIMOXGmT4UoMMX2qQ4UsNMnypQYYvNcjwpQYZvtQgw5caZPhSgwxfapDhSw0yfKlBhi81yPClBhm+1CDDlxpk+FKDDF9q0Njhd6fKfjbJ493yJUn2JTmY5KEkK4cbU1KfzuQZ/w5GJ8s86W7gnqq6FHgL2NHnYJKGM1b4STYAfwJ8tVsOcAPwSLfJbmDbEANK6t+4z/hfAb4A/Lpb/jDwdlWd6JYPARf1PJukgSwZfpJPA0er6pmz2UGSnUnmk8wf59jZ/Cck9WzJ02QD1wGfSfIpYDXwQeBeYE2S5d2z/gbg8GI3rqpdwC6AD2Zt9TK1pIks+YxfVXdV1Yaq2gjcCnyvqj4LPAnc3G22Hdgz2JSSejXJv+P/FfCXSQ4yes9/Xz8jSRraOC/1f6Oqvg98v7v8CnBV/yNJGprf3JMaZPhSgwxfapDhSw0yfKlBhi81yPClBhm+1CDDlxpk+FKDDF9qkOFLDTJ8qUGGLzXI8KUGGb7UIMOXGmT4UoMMX2qQ4UsNMnypQYYvNcjwpQYZvtQgw5caNNaZdJK8CvwCeAc4UVVzSdYCDwEbgVeBW6rqrWHGlNSnM3nG/0RVba6quW75TmBvVW0C9nbLks4Bk7zUvwnY3V3eDWybfBxJ0zBu+AV8N8kzSXZ269ZV1ZHu8hvAut6nkzSIcc+We31VHU7yu8ATSX688MqqqiS12A27B4qdAKs5b6JhJfVjrGf8qjrc/T4KPMro9NhvJlkP0P0+eprb7qqquaqaW8GqfqaWNJElw09yfpIPnLwMfBLYDzwGbO822w7sGWpISf0a56X+OuDRJCe3/6eq+naSp4GHk+wAXgNuGW5MSX1aMvyqegW4YpH1/w1sGWIoScPym3tSgwxfapDhSw0yfKlBhi81yPClBhm+1CDDlxpk+FKDDF9qkOFLDTJ8qUGGLzXI8KUGGb7UIMOXGmT4UoMMX2qQ4UsNMnypQYYvNcjwpQYZvtQgw5caZPhSg8YKP8maJI8k+XGSA0muTbI2yRNJXu5+XzD0sJL6Me4z/r3At6vqMkan0zoA3AnsrapNwN5uWdI5YJyz5X4I+DhwH0BV/aqq3gZuAnZ3m+0Gtg01pKR+jfOMfwnwM+Afkzyb5Kvd6bLXVdWRbps3GJ1VV9I5YJzwlwMfA/6hqq4EfskpL+urqoBa7MZJdiaZTzJ/nGOTziupB+OEfwg4VFX7uuVHGD0QvJlkPUD3++hiN66qXVU1V1VzK1jVx8ySJrRk+FX1BvB6ko92q7YALwKPAdu7dduBPYNMKKl3y8fc7i+AB5KsBF4B/ozRg8bDSXYArwG3DDOipL6NFX5VPQfMLXLVln7HkTQNfnNPapDhSw0yfKlBhi81yPClBhm+1CDDlxqU0dfsp7Sz5GeMvuzzEeC/prbjxf02zADOcSrneLczneP3qurCpTaaavi/2WkyX1WLfSGoqRmcwzlmNYcv9aUGGb7UoFmFv2tG+13ot2EGcI5TOce7DTLHTN7jS5otX+pLDZpq+Em2JnkpycEkUzsqb5L7kxxNsn/BuqkfHjzJxUmeTPJikheS3DGLWZKsTvJUkue7Ob7Yrb8kyb7u/nmoO/7C4JIs647n+Pis5kjyapIfJXkuyXy3bhZ/I1M5lP3Uwk+yDPh74I+By4Hbklw+pd1/Ddh6yrpZHB78BPD5qrocuAa4vft/MO1ZjgE3VNUVwGZga5JrgLuBe6rqUuAtYMfAc5x0B6NDtp80qzk+UVWbF/zz2Sz+RqZzKPuqmsoPcC3wnQXLdwF3TXH/G4H9C5ZfAtZ3l9cDL01rlgUz7AFunOUswHnAD4CrGX1RZPli99eA+9/Q/THfADwOZEZzvAp85JR1U71fgA8B/0n32duQc0zzpf5FwOsLlg9162ZlpocHT7IRuBLYN4tZupfXzzE6SOoTwE+Bt6vqRLfJtO6frwBfAH7dLX94RnMU8N0kzyTZ2a2b9v0ytUPZ++Ee73148CEkeT/wDeBzVfXzWcxSVe9U1WZGz7hXAZcNvc9TJfk0cLSqnpn2vhdxfVV9jNFb0duTfHzhlVO6XyY6lP2ZmGb4h4GLFyxv6NbNyliHB+9bkhWMon+gqr45y1kAanRWpCcZvaRek+TkcRincf9cB3wmyavAg4xe7t87gzmoqsPd76PAo4weDKd9v0x0KPszMc3wnwY2dZ/YrgRuZXSI7lmZ+uHBk4TRqcgOVNWXZzVLkguTrOkuv4/R5wwHGD0A3DytOarqrqraUFUbGf09fK+qPjvtOZKcn+QDJy8DnwT2M+X7paZ5KPuhPzQ55UOKTwE/YfR+8m+muN+vA0eA44weVXcwei+5F3gZ+Fdg7RTmuJ7Ry7QfAs91P5+a9izAHwDPdnPsB/62W//7wFPAQeCfgVVTvI/+EHh8FnN0+3u++3nh5N/mjP5GNgPz3X3zL8AFQ8zhN/ekBvnhntQgw5caZPhSgwxfapDhSw0yfKlBhi81yPClBv0/IkY+Muem1BIAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f862dc956d0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "from nilearn.input_data import NiftiMasker\n",
    "masker = NiftiMasker(low_pass=0.5, verbose=0)\n",
    "masker.fit(nii_images)\n",
    "print ('shape:', masker.mask_img_.get_data().shape)\n",
    "print ('values:', np.unique(masker.mask_img_.get_data()))\n",
    "_= plt.imshow(masker.mask_img_.get_data()[0])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## CanICA"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "CanICA is an ICA method for group-level analysis of fMRI data. Compared to other strategies, it brings a well-controlled group model, as well as a thresholding algorithm controlling for specificity and sensitivity with an explicit model of the signal. \n",
    "\n",
    "[nilearn.decomposition.CanICA](http://nilearn.github.io/modules/generated/nilearn.decomposition.CanICA.html)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/mnt/datagrid/personal/borovec/Applications/vEnv2/lib/python2.7/site-packages/nilearn/signal.py:139: UserWarning: Detrending of 3D signal has been requested but would lead to zero values. Skipping.\n",
      "  warnings.warn('Detrending of 3D signal has been requested but '\n",
      "/mnt/datagrid/personal/borovec/Applications/vEnv2/lib/python2.7/site-packages/nilearn/signal.py:51: UserWarning: Standardization of 3D signal has been requested but would lead to zero values. Skipping.\n",
      "  warnings.warn('Standardization of 3D signal has been requested but '\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "CanICA(detrend=True, do_cca=True, high_pass=None, low_pass=None,\n",
       "    mask=<nibabel.nifti1.Nifti1Image object at 0x7f8634425510>,\n",
       "    mask_args=None, mask_strategy='background',\n",
       "    memory=Memory(cachedir=None), memory_level=0, n_components=7, n_init=5,\n",
       "    n_jobs=1, random_state=None, smoothing_fwhm=6, standardize=True,\n",
       "    t_r=None, target_affine=None, target_shape=None, threshold='auto',\n",
       "    verbose=0)"
      ]
     },
     "execution_count": 10,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from nilearn.decomposition import CanICA\n",
    "canica = CanICA(mask=mask_full, n_components=nb_patterns, \n",
    "                mask_strategy='background',\n",
    "                threshold='auto', n_init=5, verbose=0)\n",
    "canica.fit(nii_images)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x7f862c411ed0>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAD8CAYAAABXXhlaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAD5xJREFUeJzt3X+MFOd9x/H3p5gfNomNCS6lYBUikC0qYRxdDchWRUydUNeK+4dlxY0qVCGdVLmV06aKoZWqRGql+J8Q/1FFOhU3/OEGO04cEEVx6BVaRbKxzzV2+BHChWIZAj7TgG1FKgXn2z92Dpbt/Zi7nZnd2efzktDNzM7ufLm97z7f53lmZxQRmFlafq3TAZhZ9Zz4Zgly4pslyIlvliAnvlmCnPhmCXLimyWorcSXtFHScUnDkrYUFZSZlUvTPYFH0gzgp8ADwGngNeCxiDhaXHhmVoYb2njuPcBwRJwEkLQTeBgYN/FnaXbMYe7krzz3xutW5y39sI0wu8vFI+38ys0m9j/8kv+NS5psv3b+ChcD7zStnwbWTPSEOcxljTZM/sqrVl23+tD2f596dF1qz2/f2ukQrIcdjMFc+5Xe/EjqB/oB5nBT2YczsxzaGdw7A9zetL4k23adiBiIiL6I6JvJ7DYOZ2ZFaSfxXwNWSFomaRbweWB3MWGZWZmmXepHxBVJfwa8BMwAnomII4VFZmalaauPHxF7gb0FxWJmFemeuaW110bye2kUv9VDRy5cXfYIv3WKT9k1S5AT3yxBTnyzBDnxzRLkxDdLkBPfLEFdM53Xy1N41r7hbWvbfo3lf/FKAZH0Brf4Zgly4pslqGtKfbMiyvnpvn5q3QC3+GYJcuKbJcilfgc1f2EH0vjSTtnl/HQ1x5VC2e8W3yxBTnyzBDnxzRLkPr6Vrlv79eNpjbcX+/xu8c0S5MQ3S5BL/Q6q+/Td+3uX59vxeLlxlK0Xp/rc4pslyIlvliAnvlmC3MevWJX9+pd+fqjt11h76JECIukdvTLVN2mLL+kZSSOSDjdtmy9pn6QT2c96j1KZJSZPqf8tYGPLti3AYESsAAazdTOriUlL/Yj4D0lLWzY/DKzPlncAB4AnC4yrZ9R9yq4It91x/urye8cXdDASGzXdwb2FEXE2Wz4HLCwoHjOrQNuj+hERQIz3uKR+SUOShi5zqd3DmVkBpjuq/66kRRFxVtIiYGS8HSNiABgAuFnzx/2A6CUu79NR17P6ptvi7wY2ZcubgF3FhGNmVcgznfdt4GXgDkmnJW0GvgY8IOkE8HvZupnVRJ5R/cfGeWhDwbGYWUV85l4BerlP/8rqF64uF3EWX/PUHnh6r1N8rr5Zgpz4ZglyqT9NvVze2/TUaWrPLb5Zgpz4Zgly4pslqGv6+M195tZ7ynUL9+utV7jFN0uQE98sQV1T6ncLl/PV8kU6OsMtvlmCnPhmCerKUr+13C5zlN+lfX7NX9gBX3q7ztzimyXIiW+WICe+WYK6so/fajpn9T379u+M+9gtDw63HZNZnbnFN0uQE98sQbUo9Zs1l/3v712e+3ku782ucYtvliAnvlmCnPhmCapdH79Za7+9uc+fap/+pZ8f6nQI0+Zr7lcnzy20bpe0X9JRSUckPZFtny9pn6QT2U+f9G5WE3lK/SvAlyJiJbAWeFzSSmALMBgRK4DBbN3MaiDPvfPOAmez5Q8lHQMWAw8D67PddgAHgCdLiTKnVMv7qvjbeL1jSoN7kpYCdwMHgYXZhwLAOWBhoZGZWWlyJ76kjwHfBb4YER80PxYRAcQ4z+uXNCRp6DKX2grWzIqRK/ElzaSR9M9GxPeyze9KWpQ9vggYGeu5ETEQEX0R0TeT2UXEbGZtyjOqL2A7cCwivt700G5gU7a8CdhVfHhmVoY88/j3An8M/FjS6CTxXwNfA56XtBl4G3i0nBDNrGh5RvV/BGichzcUG46ZVaHWZ+6ZdZNuvzV2M5+rb5YgJ75Zgpz4Zgly4pslyIlvliAnvlmCPJ1n1oY6TeE1c4tvliAnvlmCXOrXXJ2vsVdHdS3tW7nFN0uQE98sQU58swQ58c0S5MQ3S5AT3yxBns4za9ErU3YTcYtvliAnvlmCXOqbkUZ538wtvlmCnPhmCXLimyUonT7+2lWdjmByr7w15ad89jdXX7feS9/Wa+53D29b2/Zr2DV57p03R9Krkt6UdETSV7PtyyQdlDQs6TlJs8oP18yKkKfUvwTcHxF3AauBjZLWAk8B2yJiOXAB2FxemGZWJDVubZ9zZ+km4EfAnwL/AvxGRFyRtA74SkR8dqLn36z5sUYl3m6vDuV8XtMo+6dioi5Ba/dhPO/vXd52HLc8ONz2a9g1B2OQD+IX493r8qpcg3uSZmR3yh0B9gE/Ay5GxJVsl9PA4ukGa2bVypX4EfFRRKwGlgD3AHfmPYCkfklDkoYuc2maYZpZkaY0nRcRF4H9wDpgnqTRWYElwJlxnjMQEX0R0TeT2W0Fa2bFmHQ6T9JtwOWIuCjpRuABGgN7+4FHgJ3AJmBXmYGOqZf69K0m+r8V0P/P24+fiPvn9ZVnHn8RsEPSDBoVwvMRsUfSUWCnpL8D3gC2lxinmRVo0sSPiLeAu8fYfpJGf9/MaqZ+Z+71cnlvVhGfq2+WICe+WYLqUeq7vL9e8++j5DP8rDe5xTdLkBPfLEFOfLMEdWcf3336/Fp/V+7zWw5u8c0S5MQ3S1B3lvpWuPP963Lvu2Dg5RIjsW7gFt8sQU58swQ58c0S1D19fE/hFWIqffl2X8NjAfXlFt8sQU58swR1T6lvtdPcJXDZXy9u8c0S5MQ3S5BL/Zo7v2pup0MA/v9MgEv/7uYW3yxBTnyzBDnxzRLkPr6VwlN93S13i5/dKvsNSXuy9WWSDkoalvScpFnlhWlmRZpKqf8EcKxp/SlgW0QsBy4Am4sMzMzKk6vUl7QE+APg74G/lCTgfuCPsl12AF8BvllCjNaiW6bw8nLZ333ytvjfAL4M/Cpb/wRwMSKuZOungcUFx2ZmJZk08SU9BIxExOvTOYCkfklDkoYuc2k6L2FmBctT6t8LfE7Sg8Ac4GbgaWCepBuyVn8JcGasJ0fEADAAcLPmRyFRm1lbJk38iNgKbAWQtB74q4j4gqTvAI8AO4FNwK62Imm+HrwvynGduvXprfu1cwLPkzQG+oZp9Pm3FxOSmZVtSifwRMQB4EC2fBK4p/iQzKxsPnPPKuVv8XUHn6tvliAnvlmCurPUb73ja4Kj/B7JtzK5xTdLkBPfLEFOfLMEdWcfv5XP6jMrlFt8swQ58c0SVI9Sv1kPl/2ewrOquMU3S5AT3yxBTnyzBNWvj9/Mp/aaTYtbfLMEOfHNElTvUr9Va+nfrAu7AZ6+s05xi2+WICe+WYJ6q9SfyETdgE5ZtW7yfcxK4BbfLEFOfLMEOfHNEuTEN0tQrsE9SaeAD4GPgCsR0SdpPvAcsBQ4BTwaERfKCdPMijSVFv/TEbE6Ivqy9S3AYESsAAazdTOrgXZK/YeBHdnyDuAP2w/HzKqQN/ED+KGk1yX1Z9sWRsTZbPkcsLDw6MysFHlP4LkvIs5I+nVgn6SfND8YESEpxnpi9kHRDzCHm9oK1syKkavFj4gz2c8R4EUat8d+V9IigOznyDjPHYiIvojom8nsYqI2s7ZMmviS5kr6+Ogy8BngMLAb2JTttgnYVVaQZlasPKX+QuBFSaP7/3NE/EDSa8DzkjYDbwOPlhemmRVp0sSPiJPAXWNs/29gQxlBmVm5fOaeWYKc+GYJcuKbJciJb5agdK7A04UWDLx83fr5/t6/Ik/r/9k6wy2+WYKc+GYJcuKbJciJb5YgJ75Zgjyq30WaR7x7aYTfI/ndxy2+WYKc+GYJcuKbJch9/C5V9/6++/XdzS2+WYKc+GYJcqlfA3X4Mo9L+3pxi2+WICe+WYKc+GYJch+/hrplqs/9+vpyi2+WICe+WYJc6tdc3nJ7Kl0Cl/C9L1eLL2mepBck/UTSMUnrJM2XtE/SieznrWUHa2bFyFvqPw38ICLupHE7rWPAFmAwIlYAg9m6mdWAIsa8rf21HaRbgEPAJ6NpZ0nHgfURcTa7TfaBiLhjote6WfNjjXy7PbOyHIxBPohfaLL98rT4y4D3gH+S9Iakf8xul70wIs5m+5yjcVddM6uBPIl/A/Ap4JsRcTfwS1rK+qwSGLN0kNQvaUjS0GUutRuvmRUgT+KfBk5HxMFs/QUaHwTvZiU+2c+RsZ4cEQMR0RcRfTOZXUTMZtamSRM/Is4B70ga7b9vAI4Cu4FN2bZNwK5SIjSzwuWdx/9z4FlJs4CTwJ/Q+NB4XtJm4G3g0XJCNLOi5Ur8iDgE9I3xkIfozWrIp+yaJciJb5YgJ75Zgpz4Zgly4pslyIlvliAnvlmCJv12XqEHk96jcbLPAuB8ZQceWzfEAI6jleO43lTj+K2IuG2ynSpN/KsHlYYiYqwTgpKKwXE4jk7F4VLfLEFOfLMEdSrxBzp03GbdEAM4jlaO43qlxNGRPr6ZdZZLfbMEVZr4kjZKOi5pWFJlV+WV9IykEUmHm7ZVfnlwSbdL2i/pqKQjkp7oRCyS5kh6VdKbWRxfzbYvk3Qwe3+ey66/UDpJM7LrOe7pVBySTkn6saRDkoaybZ34G6nkUvaVJb6kGcA/AL8PrAQek7SyosN/C9jYsq0Tlwe/AnwpIlYCa4HHs99B1bFcAu6PiLuA1cBGSWuBp4BtEbEcuABsLjmOUU/QuGT7qE7F8emIWN00fdaJv5FqLmUfEZX8A9YBLzWtbwW2Vnj8pcDhpvXjwKJseRFwvKpYmmLYBTzQyViAm4D/BNbQOFHkhrHerxKPvyT7Y74f2AOoQ3GcAha0bKv0fQFuAf6LbOytzDiqLPUXA+80rZ/OtnVKRy8PLmkpcDdwsBOxZOX1IRoXSd0H/Ay4GBFXsl2qen++AXwZ+FW2/okOxRHADyW9Lqk/21b1+1LZpew9uMfElwcvg6SPAd8FvhgRH3Qiloj4KCJW02hx7wHuLPuYrSQ9BIxExOtVH3sM90XEp2h0RR+X9LvND1b0vrR1KfupqDLxzwC3N60vybZ1Sq7LgxdN0kwaSf9sRHyvk7EARMRFYD+NknqepNHrMFbx/twLfE7SKWAnjXL/6Q7EQUScyX6OAC/S+DCs+n1p61L2U1Fl4r8GrMhGbGcBn6dxie5Oqfzy4JIEbAeORcTXOxWLpNskzcuWb6QxznCMxgfAI1XFERFbI2JJRCyl8ffwbxHxharjkDRX0sdHl4HPAIep+H2JKi9lX/agScsgxYPAT2n0J/+mwuN+GzgLXKbxqbqZRl9yEDgB/Cswv4I47qNRpr1F436Eh7LfSaWxAKuAN7I4DgN/m23/JPAqMAx8B5hd4Xu0HtjTiTiy472Z/Tsy+rfZob+R1cBQ9t58H7i1jDh85p5Zgjy4Z5YgJ75Zgpz4Zgly4pslyIlvliAnvlmCnPhmCXLimyXo/wC/FyMy1axLFAAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f862c664c50>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "components = np.argmax(canica.components_, axis=0) + 1\n",
    "atlas = components.reshape(atlas_gt.shape)\n",
    "plt.imshow(atlas)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x7f862c31c190>"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAD8CAYAAABXXhlaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAADxFJREFUeJzt3V2sHPV9xvHvg/FLIHWME2odbFS7wg3yRTHpETYFRY5diENRzEWEIFHlVpZ8QyuipgqmVaumLxLchHDRprIKxRe8BkJsWSTEPcWqqraG42KIXwA7rhF2bE7aYJlGqmPDrxc7NuvTs+eMz87Mnj2/5yNZOzM7u/ODPc/O/z8z+x9FBGaWyyW9LsDMmufgmyXk4Jsl5OCbJeTgmyXk4Jsl5OCbJdRV8CWtlfSmpEOSNlVVlJnVS5O9gEfSDOAt4BbgKPAKcHdE7K+uPDOrw6VdvPYG4FBEHAaQ9BSwDugY/FmaHXO4vItNmtl4/pef84s4rYnW6yb4C4F32uaPAivGe8EcLmeF1nSxSTMbz64YKrVeN8EvRdJGYCPAHC6re3NmVkI3B/eOAVe3zS8qll0gIjZHxGBEDM5kdhebM7OqdBP8V4ClkpZImgXcBWyrpiwzq9Okm/oRcVbS7wMvAjOARyNiX2WVmVltuurjR8QLwAsV1WJmDan94N5U9JU3jpZa7/FrF9VciVlv+JJds4QcfLOEpm1Tv2xzfrLv4W6A9TPv8c0ScvDNEnLwzRKaVn38Kvr1k9mW+/vWb7zHN0vIwTdLyME3S8jBN0vIwTdLqK+P6jd5FH88o+vwUX6b6rzHN0vIwTdLyME3S8jBN0vIwTdLyME3S8jBN0vIwTdLyME3S8jBN0uory/ZtTwOPvYbk3rd0t/dXXEl08OEe3xJj0oakbS3bdl8STskHSwer6i3TDOrUpmm/mPA2lHLNgFDEbEUGCrmzaxPTNjUj4h/lrR41OJ1wKpieguwE7ivwrrMJt28L/Me2bsAkz24tyAijhfTJ4AFFdVjZg3o+qh+RAQQnZ6XtFHSsKThM5zudnNmVoHJHtV/V9JARByXNACMdFoxIjYDmwHman7HL4h+5oE3qlFF036y28rW9J/sHn8bsL6YXg9sraYcM2tCmdN5TwL/Bnxa0lFJG4AHgFskHQR+q5g3sz5R5qj+3R2eWlNxLWbWkL6+cm9033qqDL6Z0cydAx2fO7PqeMfnmuzX20d8rb5ZQg6+WUJ93dTvpSfvuvX89CXLm932h3v2d/0eL/5kT6n1Pn9V5/+48Zr3pdc7Uuotatfe5chwas97fLOEHHyzhBx8s4SmVR+//fReHaf22vv1vXTJ8mUdn2vv/5ftx4+n/T1uf+sLXb+fTQ3e45sl5OCbJTStmvrtqvjF3HhN6qnqwpq7b+rXbdnin5yf3n/kqh5Wkov3+GYJOfhmCU3bpv5k9WPzvinbf+37F8xP16P8GQbp8B7fLCEH3ywhB98sIQffLCEH3ywhB98sIQffLCEH3ywhB98sIQffLCFfsjvNfP+FJ3pdwqS1/1IPevdrvel4ie5oZW6hdbWklyTtl7RP0r3F8vmSdkg6WDxeUX+5ZlaFMk39s8DXImIZsBK4R9IyYBMwFBFLgaFi3sz6QJl75x0HjhfT70s6ACwE1gGritW2ADuB+2qpskHtY9b5l3oXmq6/xsvoog7uSVoMXA/sAhYUXwoAJ4AFlVZmZrUpHXxJHweeA74aEafan4uIAKLD6zZKGpY0fIbTXRVrZtUoFXxJM2mF/vGI+G6x+F1JA8XzA8DIWK+NiM0RMRgRgzOZXUXNZtalMkf1BTwCHIiIb7Y9tQ1YX0yvB7ZWX56Z1aHMefybgN8BfiTp3LCtfww8ADwjaQPwNnBnPSWaWdXKHNX/F0Adnl5TbTlm1gRfuTeO0bej9um96SvD1XrtfK2+WUIOvllCbupfhNFN/3Oa7gJ0qgPg81ctPz9dxd1ym3Rm1fELFzxW3490sjXtR/Me3ywhB98sIQffLCH38SswXp+7l/q5vw8X9sNH389uMu9hH/Ee3ywhB98sIbV+UduMuZofK+SrfKea8boB7d2F0WbuHOh62//vFJ51ZVcMcSp+1ukS+/O8xzdLyME3S8jBN0vIp/Ns3H78eNw/71/e45sl5OCbJeTgmyXk4Jsl5OCbJeTgmyXk4Jsl5OCbJeTgmyXk4JslVObeeXMkvSzpNUn7JH2jWL5E0i5JhyQ9LWlW/eWaWRXK7PFPA6sj4jpgObBW0krgQeChiLgGeA/YUF+ZZlalCYMfLf9TzM4s/gWwGni2WL4FuKOWCs2scqX6+JJmFHfKHQF2AD8GTkbE2WKVo8DCeko0s6qVCn5EfBARy4FFwA3AtWU3IGmjpGFJw2c4PckyzaxKF3VUPyJOAi8BNwLzJJ37Pf8i4FiH12yOiMGIGJzJ7K6KNbNqlDmqf6WkecX0x4BbgAO0vgC+VKy2HthaV5FmVq0yI/AMAFskzaD1RfFMRGyXtB94StJfAa8Cj9RYp5lVaMLgR8TrwPVjLD9Mq79vZn3GV+6ZJeTgmyXk4JsllH547VNfXtno9uY+8e+Nbs9sLN7jmyXk4Jsl5OCbJZSmj990X76T9jrc37de8R7fLCEH3yyhadvUnypN+/GMrtFNf2uK9/hmCTn4Zgk5+GYJTas+fj/068fjU33WFO/xzRJy8M0ScvDNEnLwzRJy8M0ScvDNEnLwzRJy8M0ScvDNEnLwzRIqHfziVtmvStpezC+RtEvSIUlPS5pVX5lmVqWL2ePfS+tmmec8CDwUEdcA7wEbqizMzOpT6kc6khYBvw38NfCHkgSsBr5crLIF+HPg2zXUaJP0zp/+ZtfvcfVf/msFldhUU3aP/y3g68CHxfwngZMRcbaYPwosrLg2M6vJhMGXdDswEhG7J7MBSRslDUsaPsPpybyFmVWsTFP/JuCLkm4D5gBzgYeBeZIuLfb6i4BjY704IjYDmwHman5UUrWZdUUR5bMoaRXwRxFxu6TvAM9FxFOS/g54PSL+drzXz9X8WKE1XRVcVj8OylHF4BtV9Os7cX9/6tsVQ5yKn2mi9bo5j38frQN9h2j1+R/p4r3MrEEXNfRWROwEdhbTh4Ebqi/JzOp2UU39bjXZ1B9tKjb9p3rTfiJu+k89TTT1zaxPOfhmCU2r4bXH06lZXXcXwMNk21TkPb5ZQg6+WUIOvllCaU7nTSe9PIXXiU/tTQ0+nWdmHTn4Zgk5+GYJOfhmCTn4Zgk5+GYJOfhmCTn4Zgk5+GYJOfhmCTn4Zgk5+GYJOfhmCTn4Zgk5+GYJOfhmCTn4ZgmVGmVX0hHgfeAD4GxEDEqaDzwNLAaOAHdGxHv1lGlmVbqYPf7nImJ5RAwW85uAoYhYCgwV82bWB7pp6q8DthTTW4A7ui/HzJpQNvgB/FDSbkkbi2ULIuJ4MX0CWFB5dWZWi7J30rk5Io5J+mVgh6Q32p+MiJA05nC9xRfFRoA5XNZVsWZWjVJ7/Ig4VjyOAM/Tuj32u5IGAIrHkQ6v3RwRgxExOJPZ1VRtZl2ZcI8v6XLgkoh4v5i+FfgLYBuwHnigeNxaZ6H2kfYx7H2bbJuMMk39BcDzks6t/0RE/EDSK8AzkjYAbwN31lemmVVpwuBHxGHgujGW/zfg2+KY9SFfuWeWkINvlpCDb5aQg2+WUNkLeGyKGn1Krc7Tez59N314j2+WkINvlpCb+tNM1Vf1uXk/PXmPb5aQg2+WkJv605ib6daJ9/hmCTn4Zgk5+GYJOfhmCTn4Zgk5+GYJOfhmCTn4Zgk5+GYJOfhmCTn4Zgk5+GYJOfhmCTn4ZgmVCr6keZKelfSGpAOSbpQ0X9IOSQeLxyvqLtbMqlF2j/8w8IOIuJbW7bQOAJuAoYhYCgwV82bWByYMvqRPAJ8FHgGIiF9ExElgHbClWG0LcEddRZpZtcrs8ZcAPwX+QdKrkv6+uF32gog4XqxzgtZddc2sD5QJ/qXAZ4BvR8T1wM8Z1ayPiABirBdL2ihpWNLwGU53W6+ZVaBM8I8CRyNiVzH/LK0vgnclDQAUjyNjvTgiNkfEYEQMzmR2FTWbWZcmDH5EnADekfTpYtEaYD+wDVhfLFsPbK2lQjOrXNlRdv8AeFzSLOAw8Hu0vjSekbQBeBu4s54SzaxqpYIfEXuAwTGeWlNtOWbWBF+5Z5aQg2+WkINvlpCDb5aQg2+WkINvlpCDb5aQWpfZN7Qx6ae0Lvb5FPBfjW14bFOhBnAdo7mOC11sHb8SEVdOtFKjwT+/UWk4Isa6IChVDa7DdfSqDjf1zRJy8M0S6lXwN/dou+2mQg3gOkZzHReqpY6e9PHNrLfc1DdLqNHgS1or6U1JhyQ1NiqvpEcljUja27as8eHBJV0t6SVJ+yXtk3RvL2qRNEfSy5JeK+r4RrF8iaRdxefzdDH+Qu0kzSjGc9zeqzokHZH0I0l7JA0Xy3rxN9LIUPaNBV/SDOBvgC8Ay4C7JS1raPOPAWtHLevF8OBnga9FxDJgJXBP8f+g6VpOA6sj4jpgObBW0krgQeChiLgGeA/YUHMd59xLa8j2c3pVx+ciYnnb6bNe/I00M5R9RDTyD7gReLFt/n7g/ga3vxjY2zb/JjBQTA8AbzZVS1sNW4FbelkLcBnwH8AKWheKXDrW51Xj9hcVf8yrge2AelTHEeBTo5Y1+rkAnwD+k+LYW511NNnUXwi80zZ/tFjWKz0dHlzSYuB6YFcvaima13toDZK6A/gxcDIizharNPX5fAv4OvBhMf/JHtURwA8l7Za0sVjW9OfS2FD2PrjH+MOD10HSx4HngK9GxKle1BIRH0TEclp73BuAa+ve5miSbgdGImJ309sew80R8RlaXdF7JH22/cmGPpeuhrK/GE0G/xhwddv8omJZr5QaHrxqkmbSCv3jEfHdXtYCEK27Ir1Eq0k9T9K5cRib+HxuAr4o6QjwFK3m/sM9qIOIOFY8jgDP0/oybPpz6Woo+4vRZPBfAZYWR2xnAXfRGqK7VxofHlySaN2K7EBEfLNXtUi6UtK8YvpjtI4zHKD1BfClpuqIiPsjYlFELKb19/BPEfGVpuuQdLmkXzo3DdwK7KXhzyWaHMq+7oMmow5S3Aa8Ras/+ScNbvdJ4Dhwhta36gZafckh4CDwj8D8Buq4mVYz7XVgT/HvtqZrAX4deLWoYy/wZ8XyXwVeBg4B3wFmN/gZrQK296KOYnuvFf/2nfvb7NHfyHJguPhsvgdcUUcdvnLPLCEf3DNLyME3S8jBN0vIwTdLyME3S8jBN0vIwTdLyME3S+j/AExFFzMFoTZoAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f862c3db050>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "max_ptn = np.max(canica.components_, axis=0).reshape(atlas_gt.shape)\n",
    "atlas[max_ptn < np.mean(max_ptn[max_ptn > 0])] = 0\n",
    "plt.imshow(atlas)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Dictionary Learning"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Perform a map learning algorithm based on spatial component sparsity, over a CanICA initialization. This yields more stable maps than CanICA.\n",
    "\n",
    "[nilearn.decomposition.DictLearning](http://nilearn.github.io/modules/generated/nilearn.decomposition.DictLearning.html) "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "DictLearning(alpha=10, batch_size=20, detrend=True, dict_init=None,\n",
       "       high_pass=None, low_pass=None,\n",
       "       mask=<nibabel.nifti1.Nifti1Image object at 0x7f8634425510>,\n",
       "       mask_args=None, mask_strategy='background',\n",
       "       memory=Memory(cachedir=None), memory_level=0, method='cd',\n",
       "       n_components=7, n_epochs=10, n_jobs=1, random_state=None,\n",
       "       reduction_ratio='auto', smoothing_fwhm=4, standardize=True,\n",
       "       t_r=None, target_affine=None, target_shape=None, verbose=0)"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "from nilearn.decomposition import DictLearning\n",
    "dict_learn = DictLearning(mask=mask_full, n_components=nb_patterns,\n",
    "                          mask_strategy='background',\n",
    "                          verbose=0, n_epochs=10)\n",
    "dict_learn.fit(nii_images)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x7f862c25f3d0>"
      ]
     },
     "execution_count": 16,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAD8CAYAAABXXhlaAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAADYFJREFUeJzt3X/oXfV9x/HnaxqTalt/VBcyI9WhVBzUWL74A6VYnZ2zUv1DpK6MMAL5xw27dbS6wVhhg/pPrYxRCNOZP1zV2rqISK3LlDFoo1/nj6qpNXWKydS0m65dYamx7/1xT7av4Ru/N997zr3JPs8HfLnnfM45OW9y7uueH/fcz0lVIaktvzLrAiRNn8GXGmTwpQYZfKlBBl9qkMGXGmTwpQZNFPwklyd5IcmOJDf2VZSkYWW5N/AkOQL4IXAZsBN4HLiuqp7vrzxJQzhygmXPBXZU1UsASe4CrgIOGPyjsrJWccwEq5T0Xv6bn/OL2pOl5psk+CcDry4Y3wmc914LrOIYzsulE6xS0nvZVlvHmm+S4I8lyUZgI8Aqjh56dZLGMMnFvV3AKQvG13Zt71JVm6pqrqrmVrBygtVJ6sskwX8cOCPJaUmOAj4D3N9PWZKGtOxD/aram+T3gYeAI4Dbq+q53iqTNJiJzvGr6kHgwZ5qkTQl3rknNcjgSw0y+FKDDL7UIIMvNcjgSw0y+FKDDL7UIIMvNcjgSw0y+FKDDL7UIIMvNcjgSw0y+FKDDL7UIIMvNWjwXnYPSed/dLz5vvfMsHVIM+IeX2qQwZcaZPClBhl8qUEGX2qQwZcaZPClBi0Z/CS3J9md5NkFbSckeTjJi93r8cOWKalP4+zx7wAu36/tRmBrVZ0BbO3GJR0mlgx+Vf0T8B/7NV8FbO6GNwNX91yXpAEt9xx/dVW91g2/DqzuqR5JUzDxxb2qKqAOND3JxiTzSebfZs+kq5PUg+UG/40kawC6190HmrGqNlXVXFXNrWDlMlcnqU/LDf79wPpueD2wpZ9yJE3DOF/nfR34LvCRJDuTbAC+DFyW5EXgN7txSYeJJX+PX1XXHWDSpT3XImlK2uyIww421Dhv2ZUaZPClBhl8qUEGX2qQwZcaZPClBhl8qUEGX2qQwZca1OadewN76N+eWtZyf/Xmhxdtv+PWK941fuKm7y7r35f2cY8vNcjgSw3yUL8HfR3a739IfyA/2XjBAad5GqBxuMeXGmTwpQYZfKlBnuMvUx/n9eOe00t9c48vNcjgSw3yUH9gy/3KbrkWftXX91d7O245f+x5T//D7/W6bvXLPb7UIIMvNcjgSw3yHH8As/zKzlt2NY5xHqF1SpJHkjyf5LkkN3TtJyR5OMmL3evxw5crqQ/jHOrvBT5fVWcB5wPXJzkLuBHYWlVnAFu7cUmHgXGenfca8Fo3/LMk24GTgauAi7vZNgOPAl8cpMpD3IE60JgGD+21HAd1cS/JqcA5wDZgdfehAPA6sLrXyiQNZuzgJ3k/8E3gc1X104XTqqqAOsByG5PMJ5l/mz0TFSupH2MFP8kKRqG/s6q+1TW/kWRNN30NsHuxZatqU1XNVdXcClb2UbOkCS15jp8kwG3A9qr6yoJJ9wPrgS93r1sGqfAQ9Vu/tu6A06587s1B1z2r83pvw/3/Y5zv8S8Efhf4fpJ9v0X9E0aBvyfJBuAV4NphSpTUt3Gu6v8zkANMvrTfciRNg3fuDeCB3/i/e5lOxK/bdOjxXn2pQQZfapDBlxpk8KUGGXypQQZfapDBlxpk8KUGGXypQQZfapDBlxpk8KUGGXypQQZfapDBlxpk8KUGGXypQQZfapDBlxpk8KUGGXypQQZfapDBlxpk8KUGLRn8JKuSPJbk6STPJflS135akm1JdiS5O8lRw5crqQ/j7PH3AJdU1dnAOuDyJOcDNwO3VNXpwJvAhuHKlNSnJYNfI//Vja7o/gq4BLi3a98MXD1IhZJ6N9Y5fpIjuifl7gYeBn4EvFVVe7tZdgInD1OipL6NFfyqeqeq1gFrgXOBM8ddQZKNSeaTzL/NnmWWKalPB3VVv6reAh4BLgCOS7LvabtrgV0HWGZTVc1V1dwKVk5UrKR+jHNV/6Qkx3XD7wMuA7Yz+gC4ppttPbBlqCIl9evIpWdhDbA5yRGMPijuqaoHkjwP3JXkL4AngdsGrFNSj5YMflU9A5yzSPtLjM73JR1mvHNPapDBlxpk8KUGGXypQQZfapDBlxpk8KUGGXypQQZfapDBlxpk8KUGGXypQQZfapDBlxo0zu/xm/WfD54+8b9x7BU7eqhE6pd7fKlBBl9qUJOH+n0cwi93XR7661DgHl9qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQaNHfzuUdlPJnmgGz8tybYkO5LcneSo4cqU1KeD2ePfwOhhmfvcDNxSVacDbwIb+ixM0nDGunMvyVrgU8BfAn+UJMAlwO90s2wG/hz42gA19u697p7zhzlqwbh7/K8CXwB+2Y1/CHirqvZ24zuBk3uuTdJAlgx+kiuB3VX1xHJWkGRjkvkk82+zZzn/hKSejXOofyHw6SRXAKuADwK3AsclObLb668Fdi22cFVtAjYBfDAnVC9VS5rIksGvqpuAmwCSXAz8cVV9Nsk3gGuAu4D1wJYB65waz8/Vgkm+x/8iowt9Oxid89/WT0mShnZQv8evqkeBR7vhl4Bz+y9J0tC8c09qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQYZfKlBBl9q0FhP0knyMvAz4B1gb1XNJTkBuBs4FXgZuLaq3hymTEl9Opg9/ieqal1VzXXjNwJbq+oMYGs3LukwMMmh/lXA5m54M3D15OVImoZxg1/Ad5I8kWRj17a6ql7rhl8HVvdenaRBjPu03IuqaleSXwUeTvKDhROrqpLUYgt2HxQbAVZx9ETFSurHWHv8qtrVve4G7mP0eOw3kqwB6F53H2DZTVU1V1VzK1jZT9WSJrJk8JMck+QD+4aBTwLPAvcD67vZ1gNbhipSUr/GOdRfDdyXZN/8f1dV307yOHBPkg3AK8C1w5UpqU9LBr+qXgLOXqT934FLhyhK0rC8c09qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQYZfKlBBl9qkMGXGmTwpQYZfKlBBl9q0FjBT3JcknuT/CDJ9iQXJDkhycNJXuxejx+6WEn9GHePfyvw7ao6k9HjtLYDNwJbq+oMYGs3LukwMM7Tco8FPg7cBlBVv6iqt4CrgM3dbJuBq4cqUlK/xtnjnwb8GPjbJE8m+Zvucdmrq+q1bp7XGT1VV9JhYJzgHwl8DPhaVZ0D/Jz9DuurqoBabOEkG5PMJ5l/mz2T1iupB+MEfyews6q2deP3MvogeCPJGoDudfdiC1fVpqqaq6q5Fazso2ZJE1oy+FX1OvBqko90TZcCzwP3A+u7tvXAlkEqlNS7I8ec7w+AO5McBbwE/B6jD417kmwAXgGuHaZESX0bK/hV9RQwt8ikS/stR9I0eOee1CCDLzXI4EsNMvhSgwy+1CCDLzXI4EsNyug2+ymtLPkxo5t9TgR+MrUVL+5QqAGsY3/W8W4HW8eHq+qkpWaaavD/d6XJfFUtdkNQUzVYh3XMqg4P9aUGGXypQbMK/qYZrXehQ6EGsI79Wce7DVLHTM7xJc2Wh/pSg6Ya/CSXJ3khyY4kU+uVN8ntSXYneXZB29S7B09ySpJHkjyf5LkkN8yiliSrkjyW5Omuji917acl2dZtn7u7/hcGl+SIrj/HB2ZVR5KXk3w/yVNJ5ru2WbxHptKV/dSCn+QI4K+B3wbOAq5LctaUVn8HcPl+bbPoHnwv8PmqOgs4H7i++z+Ydi17gEuq6mxgHXB5kvOBm4Fbqup04E1gw8B17HMDoy7b95lVHZ+oqnULvj6bxXtkOl3ZV9VU/oALgIcWjN8E3DTF9Z8KPLtg/AVgTTe8BnhhWrUsqGELcNksawGOBv4FOI/RjSJHLra9Blz/2u7NfAnwAJAZ1fEycOJ+bVPdLsCxwL/SXXsbso5pHuqfDLy6YHxn1zYrM+0ePMmpwDnAtlnU0h1eP8Wok9SHgR8Bb1XV3m6WaW2frwJfAH7ZjX9oRnUU8J0kTyTZ2LVNe7tMrSt7L+7x3t2DDyHJ+4FvAp+rqp/Oopaqeqeq1jHa454LnDn0OveX5Epgd1U9Me11L+KiqvoYo1PR65N8fOHEKW2XibqyPxjTDP4u4JQF42u7tlkZq3vwviVZwSj0d1bVt2ZZC0CNnor0CKND6uOS7OuHcRrb50Lg00leBu5idLh/6wzqoKp2da+7gfsYfRhOe7tM1JX9wZhm8B8Hzuiu2B4FfIZRF92zMvXuwZOE0aPItlfVV2ZVS5KTkhzXDb+P0XWG7Yw+AK6ZVh1VdVNVra2qUxm9H/6xqj477TqSHJPkA/uGgU8CzzLl7VLT7Mp+6Ism+12kuAL4IaPzyT+d4nq/DrwGvM3oU3UDo3PJrcCLwD8AJ0yhjosYHaY9AzzV/V0x7VqAjwJPdnU8C/xZ1/7rwGPADuAbwMopbqOLgQdmUUe3vqe7v+f2vTdn9B5ZB8x32+bvgeOHqMM796QGeXFPapDBlxpk8KUGGXypQQZfapDBlxpk8KUGGXypQf8Dn2Z1rLzKo0gAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f862de2ee10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "components = np.argmax(dict_learn.components_, axis=0) + 1\n",
    "atlas = components.reshape(atlas_gt.shape)\n",
    "plt.imshow(atlas)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.13"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}