{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Binary Pattern Dictionary Learning"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We present an image processing pipeline which accepts a large number of images, containing spatial expression information for thou- sands of genes in Drosophila imaginal discs. We assume that the gene activations are binary and can be expressed as a union of a small set of non-overlapping spatial patterns, yielding a compact representation of the spatial activation of each gene. This lends itself well to further automatic analysis, with the hope of discovering new biological relation- ships. Traditionally, the images were labeled manually, which was very time consuming. The key part of our work is a binary pattern dictionary learning algorithm, that takes a set of binary images and determines a set of patterns, which can be used to represent the input images with a small error."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      ":0: FutureWarning: IPython widgets are experimental and may change in the future.\n"
     ]
    }
   ],
   "source": [
    "%matplotlib inline\n",
    "%load_ext autoreload\n",
    "%autoreload 2\n",
    "import os, sys\n",
    "import numpy as np\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 -> /mnt/30C0201EC01FE8BC/TEMP/atomicPatternDictionary_v0/datasetFuzzy_raw\nloading (True) <- /mnt/30C0201EC01FE8BC/TEMP/atomicPatternDictionary_v0/dictionary/atlas.png\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAEICAYAAAB/KknhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAEjFJREFUeJzt3X2QXXV9x/H3xxASnkII0BASxmBB\nMHUk2JWHkaISeSyV/MEwICNrJ87OtLaFaiugo1Mcx0E7I9rRod0CElsUEMVgxgdipKPTSmAj4TEi\nUYIkTVgQ0gCtIcFv/zi/rZftPpzdex528/u8Zu7sPeeee883ufu5v9/v3LO/o4jAzPLyurYLMLPm\nOfhmGXLwzTLk4JtlyME3y5CDb5YhB3+Kk7RZ0rvbrmM0kj4q6Ya267CJcfBbIOnfJL0gadaw9TdL\n+lRbdU1GRHw6Ij7Q9H4lfVfSSx23VyQ93HQd05WD3zBJi4E/AgJ4T6vFTGMRcW5EHDh0A/4D+Hrb\ndU0XDn7zLgPuBW4GeodWSuoDLgU+klqwbw9/oqSTJP1E0g5J2yR9UdK+6TFJuk7SoKSdkh6W9Ob0\n2HslPTRaQcOHE5L+TtK/pvuLJYWkXkm/kvScpI+NtG1afp+kpyT9WtLHOl97eI9G0jslbelYPlLS\nNyQ9K+lJSX9V5j+048P0K2W2Nwe/DZcBt6Tb2ZLmA0REf1r32dSK/ckIz30V+GvgMOBUYBnw5+mx\ns4DTgTcCBwMXAb9Or/3ViHhLl3WfBhyX9vkJSW8avoGkJcD1wPuAI4FDgUVlXlzS64BvAw8CC9N+\nrpB0domnXwb8OCI2l9mXOfiNknQa8Hrg9ohYD/wCeG/Z50fE+oi4NyL2pF/yfwLekR7eDRwEHA8o\nIjZGxLYKy78mIv4nIh6kCOcJI2xzIbA6In4UEbuAjwO/Lfn6bwMOj4hPRsQrEfFL4J+Bi0s89zKK\nHpSV5OA3qxe4OyKeS8tfpaO7Px5Jb5S0WtJ2STuBT1O0/kTED4EvAl8CBiX1S5pTYe3bO+7/N3Dg\nCNscCTw9tBARL5N6HSW8HjgyDWN2SNoBfBSYP9aT0ofpEcAdJfdjOPiNkbQfRff7HSm42ym67SdI\nGmo9x/tTyeuBnwHHRsQcimBo6MGI+IeI+ENgCUWX/29LlvcysH/H8hElnzfcNuCooQVJ+1N098vs\n52ngyYiY23E7KCLOG2efvcA3I+KlSdacJQe/OcspxuhLgKXp9ibgxxRdVYBngDeM8RoHATuBlyQd\nD/zZ0AOS3ibpZEkzKQL2G1I3W9L7JW0e43U3ABdLmimph6LLPhl3AOdLOi0ddPwkr/0d2wCcJ2me\npCOAKzoeuw94UdKVkvaTNEPSmyW9bbSddXyY3jzJerPl4DenF/hyRPwqIrYP3Si655dK2ge4EViS\nurrfGuE1/obimMCLFOPf2zoem5PWvQA8RdHF/vv02FHAv49R28eB30/PvYZiCDJhEfEo8MH0/G3p\n9bZ0bPIvFMcHNgN3d9YfEa8C51N8ID4JPAfcQHGgcjTLgR3APZOpN2fyRBx7P0l3A5dHxMYW9r0Z\n+EBE/KDpfdvo9mm7AKtfRJzVdg02tbirb5Yhd/XNMtRViy/pHEmPS9ok6aqqijKzek26xZc0A/g5\ncCbFkdv7gUsi4rHRnrOvZsVsDpjU/sxsfL/hZV6JXRpvu24O7p0EbEqnViLpVuACYNTgz+YATtay\nLnZpZmNZF2tLbddNV38hHadnUrT6C4dvJKlP0oCkgd3s6mJ3ZlaV2o/qR0R/RPRERM9MZo3/BDOr\nXTfB30rHedkUf365tbtyzKwJ3QT/fuBYSUen87IvBu6qpiwzq9OkD+5FxB5JfwF8H5gB3JTO1Taz\nKa6rU3Yj4jvAdyqqxcwa4lN2zTLk4JtlyME3y5CDb5YhB98sQw6+WYYcfLMMOfhmGXLwzTKU5WSb\n//WdY0ptd/B5m2quxKwdbvHNMuTgm2Uom65+2e79WM9x19/2Fm7xzTLk4JtlyME3y5CDb5YhB98s\nQw6+WYYcfLMMOfhmGXLwzTLk4JtlKJtTdi0P5z664zXL3/2DuS1VMrW5xTfL0LjBl3STpEFJj3Ss\nmydpjaQn0s9D6i3TzKqkiBh7A+l04CXgKxHx5rTus8DzEXGtpKuAQyLiyvF2Nkfz4mQtq6DsiZvM\nX+eNxX+pV7/h3fZu5dDtXxdr2RnPa7ztxm3xI+JHwPPDVl8ArEz3VwLLJ1yhmbVmsgf35kfEtnR/\nOzB/tA0l9QF9ALPZf5K7M7MqdX1UPyJC0qjjhYjoB/qh6Op3u782uXtfr6q79ja6yR7Vf0bSAoD0\nc7C6ksysbpMN/l1Ab7rfC6yqphwza0KZr/O+BvwEOE7SFkkrgGuBMyU9Abw7LZvZNDHuGD8iLhnl\noXa+l7O9SpPj+s595fDV3lh85p5Zhhx8swxl80c6nV/F+RJa7ZkqX9nl/sc8bvHNMuTgm2XIwTfL\nUDZj/E5Vj92f6zu10tebiMP6f9LavsuaKuN6+x23+GYZcvDNMjTuRBxVanMijiq02aWfjKk4DJgO\n3f7p/NVeZRNxmNnex8E3y1CWR/XLmm5d++HGqn8yw4Dv/+eG0tuefeTSCb++NcctvlmGHHyzDDn4\nZhnyGH8Mw8fB03nMP9mv9iYyrh/teR7vTz1u8c0y5OCbZchd/TFM5679cMP/LVPxrD5rjlt8sww5\n+GYZcvDNMuQxfiY8ph/bdP6LvMlwi2+WoTKX0DpK0j2SHpP0qKTL0/p5ktZIeiL9PKT+cs2sCmW6\n+nuAD0fETyUdBKyXtAZ4P7A2Iq6VdBVwFXBlfaXaRLXZve88W286TL6Rm3Fb/IjYFhE/TfdfBDYC\nC4ELgJVps5XA8rqKNLNqTejgnqTFwInAOmB+RGxLD20H5o/ynD6gD2A2+0+2TjOrUOmDe5IOBL4B\nXBEROzsfi2LivhEn74uI/ojoiYiemczqqlgzq0apFl/STIrQ3xIR30yrn5G0ICK2SVoADNZVZFv2\npr/Oa9pof53n8f7UUOaovoAbgY0R8bmOh+4CetP9XmBV9eWZWR3KtPhvB94HPCxp6GP8o8C1wO2S\nVgBPARfVU6KZVc3z6lcgl0to1THZZltd/731TD3Pq29mo3LwzTLkrn6Lqp73frqru9u/t3bvO7mr\nb2ajcvDNMuTgm2XIE3G0KMdx/FiGj8HLjvlzGLtXzS2+WYYcfLMMuatvU5a78PVxi2+WIQffLEMO\nvlmGHHyzDDn4Zhly8M0y5OCbZcjBN8uQg2+WIQffLEMOvlmGHHyzDDn4Zhly8M0y5OCbZajMtfNm\nS7pP0oOSHpV0TVp/tKR1kjZJuk3SvvWXa2ZVKNPi7wLOiIgTgKXAOZJOAT4DXBcRxwAvACvqK9PM\nqjRu8KPwUlqcmW4BnAHckdavBJbXUqGZVa7UGF/SjHSl3EFgDfALYEdE7EmbbAEWjvLcPkkDkgZ2\ns6uKms2sS6WCHxGvRsRSYBFwEnB82R1ERH9E9EREz0xmTbJMM6vShI7qR8QO4B7gVGCupKHJOhcB\nWyuuzcxqUuao/uGS5qb7+wFnAhspPgAuTJv1AqvqKtLMqlVmeu0FwEpJMyg+KG6PiNWSHgNulfQp\n4AHgxhrrNLMKjRv8iHgIOHGE9b+kGO+b2TTjM/fMMuTgm2XIwTfLkINvliEH3yxDDr5ZhnyZ7Dqc\n8pbf3b/3ofbqMBuFW3yzDDn4ZhnKs6vf2RUfrrNrPtZ2k92Xu/42BbjFN8uQg2+WIQffLEP5jPHL\njterGNeXfX2P960lbvHNMuTgm2Vo7+3q191lr4K/6rOWuMU3y5CDb5YhB98sQw6+WYYcfLMMOfhm\nGXLwzTLk4JtlqHTw06WyH5C0Oi0fLWmdpE2SbpO0b31lmlmVJnLm3uUUF8uck5Y/A1wXEbdK+kdg\nBXB9xfVN3vCz4KbimXw+U89aUqrFl7QI+GPghrQs4AzgjrTJSmB5HQWaWfXKdvU/D3wE+G1aPhTY\nERF70vIWYOFIT5TUJ2lA0sBudnVVrJlVY9zgSzofGIyI9ZPZQUT0R0RPRPTMZNZkXsLMKlZmjP92\n4D2SzgNmU4zxvwDMlbRPavUXAVvrK7MCVUyiWfVrNGjT504pve0xH7q3xkpsKhi3xY+IqyNiUUQs\nBi4GfhgRlwL3ABemzXqBVbVVaWaV6uZ7/CuBD0naRDHmv7GaksysboqIxnY2R/PiZC1rbH+5m0j3\nfjTu9k8v62ItO+N5jbedz9wzy5CDb5ahvXfOvUxV0b0f7fXc7d97uMU3y5CDb5YhB98sQw6+WYYc\nfLMMOfhmGXLwzTLk4JtlyME3y5CDb5YhB98sQw6+WYYcfLMMOfhmGXLwzTLk4JtlyME3y5CDb5Yh\nB98sQw6+WYYcfLMMOfhmGSo1vbakzcCLwKvAnojokTQPuA1YDGwGLoqIF+op08yqNJEW/10RsTQi\netLyVcDaiDgWWJuWzWwa6KarfwGwMt1fCSzvvhwza0LZ4Adwt6T1kvrSuvkRsS3d3w7MH+mJkvok\nDUga2M2uLss1syqUvYTWaRGxVdLvAWsk/azzwYgISSNedjci+oF+KK6W21W1ZlaJUi1+RGxNPweB\nO4GTgGckLQBIPwfrKtLMqjVu8CUdIOmgofvAWcAjwF1Ab9qsF1hVV5FmVq0yXf35wJ2Shrb/akR8\nT9L9wO2SVgBPARfVV6aZVUkRzQ2752henKxlje0vd1VcMtuXxp5e1sVadsbzGm87n7lnliEH3yxD\nDr5Zhsp+j2/TUOf4fCLjfY/r935u8c0y5OCbZchd/Uy4+26d3OKbZcjBN8uQg2+WIQffLEMOvlmG\nHHyzDDn4Zhly8M0y5OCbZcjBN8uQg2+WIQffLEMOvlmGHHyzDDn4Zhly8M0y5OCbZcjBN8tQqeBL\nmivpDkk/k7RR0qmS5klaI+mJ9POQuos1s2qUbfG/AHwvIo4HTgA2AlcBayPiWGBtWjazaaDM1XIP\nBk4HbgSIiFciYgdwAbAybbYSWF5XkWZWrTIt/tHAs8CXJT0g6YZ0uez5EbEtbbOd4qq6/4+kPkkD\nkgZ2s6uaqs2sK2WCvw/wVuD6iDgReJlh3fooLrk74mV3I6I/Inoiomcms7qt18wqUCb4W4AtEbEu\nLd9B8UHwjKQFAOnnYD0lmlnVxg1+RGwHnpZ0XFq1DHgMuAvoTet6gVW1VGhmlSt7JZ2/BG6RtC/w\nS+BPKT40bpe0AngKuKieEs2saqWCHxEbgJ4RHlpWbTlm1gSfuWeWIQffLEMOvlmGHHyzDDn4Zhly\n8M0ypOJs24Z2Jj1L8Z3/YcBzje14ZFOhBnAdw7mO15poHa+PiMPH26jR4P/fTqWBiBjpvICsanAd\nrqOtOtzVN8uQg2+WobaC39/SfjtNhRrAdQznOl6rljpaGeObWbvc1TfLkINvlqFGgy/pHEmPS9ok\nqbFZeSXdJGlQ0iMd6xqfHlzSUZLukfSYpEclXd5GLZJmS7pP0oOpjmvS+qMlrUvvz21p/oXaSZqR\n5nNc3VYdkjZLeljSBkkDaV0bvyONTGXfWPAlzQC+BJwLLAEukbSkod3fDJwzbF0b04PvAT4cEUuA\nU4APpv+DpmvZBZwREScAS4FzJJ0CfAa4LiKOAV4AVtRcx5DLKaZsH9JWHe+KiKUd35u38TvSzFT2\nEdHIDTgV+H7H8tXA1Q3ufzHwSMfy48CCdH8B8HhTtXTUsAo4s81agP2BnwInU5whts9I71eN+1+U\nfpnPAFYDaqmOzcBhw9Y1+r4ABwNPkg6611lHk139hcDTHctb0rq2lJoevC6SFgMnAuvaqCV1rzdQ\nTJK6BvgFsCMi9qRNmnp/Pg98BPhtWj60pToCuFvSekl9aV3T70tXU9lPhA/uMfb04HWQdCDwDeCK\niNjZRi0R8WpELKVocU8Cjq97n8NJOh8YjIj1Te97BKdFxFsphqIflHR654MNvS9dTWU/EU0Gfytw\nVMfyorSuLa1MDy5pJkXob4mIb7ZZC0AUV0W6h6JLPVfS0DyMTbw/bwfeI2kzcCtFd/8LLdRBRGxN\nPweBOyk+DJt+Xxqbyr7J4N8PHJuO2O4LXEwxRXdbGp8eXJIoLkW2MSI+11Ytkg6XNDfd34/iOMNG\nig+AC5uqIyKujohFEbGY4vfhhxFxadN1SDpA0kFD94GzgEdo+H2JJqeyr/ugybCDFOcBP6cYT36s\nwf1+DdgG7Kb4VF1BMZZcCzwB/ACY10Adp1F00x4CNqTbeU3XArwFeCDV8QjwibT+DcB9wCbg68Cs\nBt+jdwKr26gj7e/BdHt06Hezpd+RpcBAem++BRxSRx0+ZdcsQz64Z5YhB98sQw6+WYYcfLMMOfhm\nGXLwzTLk4Jtl6H8BEGZDb3ifE5cAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f384f1fc090>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "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)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlMAAACDCAYAAACp4J7uAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAG0tJREFUeJzt3V2MI9lVB/D/sd3uD0/PdM9Mumdm\nZ0M2EIHyQhJQJESCIhAobEd0eIkCLwtCmhciJRIRGcFLHnggSATyQtCQBK1QpCUkQb1iJGBZJRJP\nUTbf7C6b3SyJMjsfPTM93e222x9lXx7se+e4pqpddvmjfP3/Sda6y2XX3Tmuqutzbt0SYwyIiIiI\naDi5aTeAiIiIaJaxM0VERESUAjtTRERERCmwM0VERESUAjtTRERERCmwM0VERESUQqrOlIi8X0Re\nEZHXROTqqBpF08F4+oOx9Avj6Q/G0k8y7DxTIpIH8EMAvwngBoBvAvg9Y8xLo2seTQrj6Q/G0i+M\npz8YS3+lyUy9G8BrxpjXjTENAM8A2B5Ns2gKGE9/MJZ+YTz9wVh6Kk1n6jEAP1V/3+guo9nEePqD\nsfQL4+kPxtJThXFvQESuALjS/fOXYtZ55LmI9DyPY8uUxpie5+HXqcMYE/+P2QdjmS1pYgkkiydN\nzrj3TZoc7pt+SRLPNJ2pNwA8rv6+3F0WbsQ1ANcAQESMPbHmcp2kWC6XQz6f7zSmUECxWAQALC4u\nYmlpyT2362itVgv1eh0AUKvV3PNGo4EgCNw67XYbANBut3lyjtc3nozlzBhq35xM02gIA++bk2sa\nDYj7JqJ/VIeXRf0wjxM+F0zj3JCmM/VNAG8TkSfQ+TJ8GMDv93uTPQHbf5x8Po9CodOMYrGIxcVF\nAMCpU6ewtrbmntsTs9ZoNHB0dAQA2N/fj92m/cfUAeHJ9xEDx5OxzKyh9k3KLMbTH3Mfy3AnqV+n\nSa9jf7hb9sd1uANlGWPce8d9nhi6M2WMCUTkIwD+A0AewBeMMS+OrGU0UYynPxhLvzCe/mAs/TX0\n1AhDbUzE5HI52AcALCwsuAzG0tISTp06BQBYW1vDhQsXAADr6+tuHa1er+PBgwcAgNu3b7uMxtHR\nEWq1mlun2WwC6PRidU/WmqfMRtpavsVYTt+oYgn4WUqYNaPcN0fxOTQ87puPSjKeNioDlcvlEmWm\n7HM9HGRUJb9xj5kaii4JAZ1/HFsaWlpacuWgCxcuYH19HQBw+fJlt1yLKwcFQeDG2TSbTbetdrs9\nsZTfPGAsiYjoJFFjoXTnSJ8/4jpT9rySz+d7jvutVgtA73lCs50qu/44zxm8nQwRERFRChPPTAG9\nPdV8Pu96pouLi640tL6+jsuXLwPolIm2tx+d12xnZ6fnbzuAuVwu4/j42H2+7bGKSM8AZmY00mMs\nadr6Xemj8XtCNDlxpT2bedLnjEKhEJmBAuAuWlpYWHDrt1otN+yj0Wi4K8DDwtkp24ZRHwum0pkC\nov+R8/l8z+X0thy0vb2NUqkEAFhdXUW5XHbL7Un47t277r3hQAxysKXBMZY0KWnjH34/O1dEoxfX\niQJ6S3v5fB4LCwsAOueJqE4TADe1zsrKilun0WigWq0CgDuPRJnUmFqW+YiIiIhSmFpmalCrq6sA\n4LIaAFCpVPq+Tw9q0wPfdA+VZaLJYixpGCdlKAeZ1E8v40UMRKMVN49UVGlPXwG+vLzszg06AwXA\nlf/OnTvXk5m6f/9+ZBv0HTT0Ph73fBRmpjNlhS+N1Ow/8srKikv/BUHg6qr5fN7VT3nSnT7GkvqJ\n6zQlmexPi0v1cxJYotFIMj7K/leX9paXlwF0fmSfP38eQG+nSVtfX3cdrrjSXqvVcmNrgyDoKReO\nc/wUy3xEREREKcxcZirOm970Jjea/+joyPVaj4+PXc80CIKxpvloNBhLiss66ZJB1PM49nuh7+mo\nf6XqdZjpJBrMSaU9nZEC4kt758+fx7lz5wAAm5ubbnnY1tYWAOD69es9y20Fo1KpuLKgzYBZ4xyM\nPpOdKRuUUqnkLrPf2dnB3bt3AXRKRHYdXgmWbYwlnUTHWXea9LgLPQt/1PdCz47cbrfdRH/2Nbuc\nHXKi9ML3a9WdKODk0t7m5qZbbjtNYXas7dbWlutQPXjwACsrK+71qHLeuMdPscxHRERElMLMZab0\nL8tKpeLmJtrf33eloUaj4daZ1/u2zQLGkqKEs0s666Qn9bO/dPP5vFsel5myA1L1RH8i4pbncrme\nX7PMUhENJq4cb/fNQUp7W1tbLgN1+vTpyO3pK8DX19fRaDQeWSduMPo4ziUz05nSI/ejRvHfvXvX\n3Sj36OjInYxbrVZkWj/8nCaHsaQo/W6EWigU3MGwWCy6g7Oe7E9fuWO1Wi13oK3X6+7z9MFXj8ED\nOH6KKImoHy/hmxfbzpTdR5OU9kqlEs6ePQugc0W33q/tOaNUKg00fqpQKETu17yaj4iIiCgDppaZ\nisoqhH9B7u/vA+gMSLaDk3Vqb2dnx61Tr9d77udWq9UAdH5x2vR9XFaD0mEsaZROuhpI36LIDjhd\nW1tzJYGouWkajYb7rtnvmH4N6P0ORQ1eJaLkbGle77+2LL+ystK3tHf27Fm3fy8vL/fcj88u104a\njG7L+kEQ9NzbddSTPk+lM6Ubq0s3+iRqyzxW+Ea4QOfAeOPGDbe+fW+9Xnf/aM1m0x0cjTGxMyHT\ncBhLGgddJoi6pHplZcXd73FjY8OVBOzrWr1ex97eXuR29JU+nASWaPT0lX1A5wePLvnFlfbsFX/6\nBsgi4pYDcOeGuPFT1WrVvZbP53s6eeGpUdJimY+IiIgohYlnpuycDjaDkcvlXO+yVqs9koa37LxD\nWr1ed1mP27dvu/fWajX3mfqKMZ3N4C/P9BhLGpW4ecP0HFL5fN5lntbW1rCxsQGgUxK4dOkSAODM\nmTOPfMbBwUHkZ+sr+/R2wr9YeWUf0fD00A+gU1q3maNyuezKc1tbW+7KPT0/lb7wRFtYWHAZq/Bg\ndHsuiZuncBymVubTk+XpS5WjHB0dxY6FsOWg/f39ntKQ/Uw943H4YMiDY3qMJY1a+Gog28kpFAru\nu6NLApcuXXKdKDseTwuXlW3av1wuu4Nxs9nsexNtIkomapJcexyvVqvuBsXFYjF2pnMt6m4HhULB\njY0KgqCn1DcNLPMRERERpTCVMp/+r74qK8yWd8rlcuz8MXYOolqt1jPRo56kT2+TvzRHh7GkSYga\njK7nmTpz5ozLSNmrgVZXV918NNvb2y47df/+/Z45qfRnE9Fg9FxNepk+Ttvjty3t6aEbPpn61Ahx\nI+p1EI6Pj2NnNrYn8CAIei6BtMt1aShq+5QeY0lZYUsGp06dcsumnf4n8tUwk3YuLS25ZYPSHbUg\nCNzknIeHh0N93ij1LfOJyOMi8jUReUlEXhSRj3aXnxWR50Tk1e5/18ffXEqLsfQH902/MJb+4L45\nf5J0DwMAf2KM+baIrAL4log8B+APADxvjPlLEbkK4CqATyTZaDiTEJXR0IOZo8pCls5a6Akdo+Yj\nYgYDAGPpk5Hvm76o1+uRtyrKOMbSH3O9byadtHN9PVlfMuoq7na73TPPlL0qsFwuu5Ji3L1dx6Fv\nZ8oYcwvAre7zsoi8DOAxANsA3tdd7WkAX0eCL4W9nD7qf0yXcfT4iPB9s/RnWXHjacJja+b9JGyM\n+Xb3v4zljBv1vumTZrPpyns7OztueoR6ve4OtOHv2bSNct+k6eK++dAgk3ba8p+eADo8Hssujxuj\n++DBA3e1YLVadT/k9Q/zUU/YCQw4ZkpE3gLgnQC+AWCz+4UBgNsANmPecwXAleGbSOPAWPqF8fQH\nY+kXxnM+JO5MicgpAF8B8DFjzGGot2hEJPKnnTHmGoBr3c8w3WWR87jo5SLSc3uHOFGZirjsRRZ+\nfWYBY+mXUcZz1tnS3urqqruyb5au5mMs/cJ4DjZpZ5Tl5eXIiZur1WrswHM9ebTdVhAEfa8GT3Ne\nSdSZEpEFdL4QXzTGfLW7+I6IXDTG3BKRiwB2B9mwLgHF/U9FnaT7fV7U+jzxPsRY+mUc8cwafQDU\nB2Q7fcbBwYHrLMXdRHsWynzzEMt5Ms/xTDppp3b9+vXYDpW9ak9PqXB4eOj28+vXr7sfUrqzpvdx\nPRZ3HD/Mk1zNJwA+D+BlY8yn1UvPAniq+/wpAI/evZayiLH0BPdN7zCWnuC+OX+kX69MRN4D4L8B\n/ACAHbX1Z+jUf78E4M0AfgLgQ8aY6FuzP/ysEzcWl2pPWhpKspwAMJa+eC8mtG+mkSTmeh0J3ULG\nXgW0uLjo5o86d+7cQPfmu3nzJgBgb28Pu7udZMD9+/cjb1sUBMFYf8H2MfOxJAAzsm+OisTcAmph\nYcHtv0tLSwA6ZTs7H9z58+fdlX2bm5tu+dbWlivTnz592pXzwnPG6Sv47ty5A6CzX9+7dw9A59Zk\nNqtVq9V6bk02yD5ujOk7DqBvZ2qUBvlSDDOGgTtzf0m+FEkwltM3qlgCgx+wxznGSB+Mc7mcOxgX\ni0V3QF5ZWcHa2hoAYGNjw92nz86KrtXrdeztdc5Xu7u77iba1WoVtVoNQKc0EHegndR3cRr7Znf9\nUWw20rzux9PcN6dBf4fC988LT42wuLiI5eVlAJ2xjefPnwfQ+YG0ubnplkeV/GznybKlvTt37rjS\n4b1799zy4+Njt483m83Y+7z2+54miSfvzUdERESUwtRuJ9PPvP6i8RFjOfsGzV6E1+83t1j476hB\n581mM7YdNv1vr9TTGo2Ge12n/XVpL2uD0cdpkrGMuyiF/KXjfdL9WvsJZ6GATibqwYMH7m870FyX\n9srlMo6PjwF09nF9JeA4y/aZ7UwR0fSFxzKd9Pqg700y7YWeXE9E3MFTa7VaLq0fdxNt+756vd5z\nE209AeA0SnuTNK1YhiddJP/EXbGdZnJM3WmyGo2GK+cFQeD262q12lPas/t4s9ns2cfH+YOJZT4i\nIiKiFJiZIqIecVmLqAkuk2Qz9IBUzRjTkw2K+hWrs0R6jplGo+GWN5tNdxuKuO3Y97ZarZ7Snl2u\nt+PT5LBZi+Ug883RbIkr7Q6TnapWq7Ele311nh4CoLPPeh/X93wd5z7OzhQROXEnV32FXfh5HH1F\nnl4/ahI9PbYifDK2z3O5XGTnJ5fLuYNnvxO9nkQwPKFf+D2zLmux1OtzLJW/0pT87HenUqnEluzt\nvq5/UOkfRkEQJCrfj/r7xzIfERERUQrMTBERAESWfnQWwv5SzOfzPZmKqAHGupSUz+d7ynC6bBd1\npV4QBO7zwxkqvS2dsep3hZrOoEQ91+v4IKux1BkLlvz8dVIGMi5Dpcvx+r6Z4XXs90xPrhverydV\n2tPYmSIiJ+nke/Z5oVCIPQFbi4uLPTcXjrq/XvhKPXtQzeVyPQdffcC02wh3sqLETbsQdVD15eSe\n5ViyzDcfBin56c6U7YBHiRqbFzcJ57hLexrLfEREREQpMDNFRD10OciWdPL5vMtILC4uutu2FIvF\nnoGiNlOhl62srLj7bBWLRZe1qFQq7tYuWjiroctEcb8yB8lMRS3zNUuS5VhazFL5LZzhDGeoosrJ\nSTNT/eaqm+Q+zs4UEQF4eDCLKgcVi0V30tX3xSuVSj2XMNuTp162sLDQc+88Ww6y98uLElUmOmny\nx2EPlL6exGctlhw/5b9+V3Lq70GSK/+STPob9fe4sMxHRERElAIzU0QE4OGcQ/pKL10OWllZAQCs\nra1hY2MDAHD27FmX5QDgMhV6GdC5IzwAnDlzBgcHB4nbpOeSARA5GH0Q85L5mIdY0mwKZyLTzD3W\nb/1Jfq+YmSIiIiJKgZkpIgLwcJyNvmw+amzNxsaGGzdz6dIlnDlzxn2GzVToZdr29jZ2dnYSt0n/\ncg3PU6TXoV6zGkuaD1HTbqT9Hkz7e8TOFBE5+nYh+XzenYDD5aBLly4B6Jxot7e33fvtyVUv00ql\nknvtpBOxngBSTyoZNTCVV4NFYywp68Kx7ndV7knvnTaW+YiIiIhSYGaKiHoupbfZg0Kh4AYtl0ql\nyHLQ9va2m3fI/m3XX11dfWQ7CwsLPetGZTQqlQrK5bJrg76Jsb5NStZ+mWYFY0mzapa/B+xMEVEP\nfZKLmptIl4P0iRmAO+nqMTp6Ar52u92zvi4T3b9/323LblefdGlwjCXRZLDMR0RERJRC4s6UiORF\n5Dsi8m/dv58QkW+IyGsi8s8iUuz3GZQNjKVfphnP1dVVd0uSUqmE06dP4/Tp0yiVSm4AdC6XQ7FY\nRLFYxNLSkls/qnQ07xhLv/BYOz8GyUx9FMDL6u9PAfgbY8zPAXgA4I9G2TAaK8bSL1OLZ7lcRr1e\nR71eRxAEaLfbaLfbqFQqODw8xOHhISqVCmq1mnvY9e1YGurBWPqFx9o5kagzJSKXAWwB+Fz3bwHw\n6wC+3F3laQAfHEcDabQYS78wnv5gLP3CeM6XpAPQ/xbAnwKwudxzAPaNMUH37xsAHot6o4hcAXAl\nTSNppBhLv4w8nvqO63aOoEaj4W4vcnBwEDkHkS71lMtlVCoVAI9eDWazGJVKxX3OwcGB+/xGo+G2\nG74zvOcYS7/wWDtH+namROQDAHaNMd8SkfcNugFjzDUA17qfNVd7UkYxlp4Y9b5pT3Z2MsUgCNBo\nNAB0TpZ7e3uPfMbOzo47CdsTrn4N6Jykw6/Z1+0s2zdv3nSfX6lU3HZtuanb3p7OgYcYS0/wvDl/\nkmSmfhXA74jIkwCWAJwG8BkAayJS6PayLwN4Y3zNpBFiLP3BfdMvjKU/uG/OGRnkV0G3h/1xY8wH\nRORfAHzFGPOMiPw9gO8bY/6uz/vZw54yY4wAjKUPbCyB0cSzUCggl8u5yRjtFVtA8vu5af3u7XZw\ncICbN28CAPb29rC7uwsA2N/fR7VaBQDUajWX2Wg2my6z0W63vctsjHLfZCyna9T75nhbS/3oeMZJ\nM2nnJwA8IyJ/AeA7AD6f4rNouhhLvwwVT1t60fdSSzLJop2gMcyOmznpdVsO2t3dxf7+PgCgWq26\n9zabzXkecwMwlr7hsdZTA3WmjDFfB/D17vPXAbx79E2iSWAs/cJ4+oOx9AvjOR94OxkiAoCeActA\n5/YftiwTp1KpuHu+hdn3nvS6Hcysy0H1et29t9VquWyGj+WgcWEsiSaLnSkiAoBHTm5BELjXwidi\ne1Isl8vu3mthdp2TXrefayd/tNuy7w2CoKdd+sTLk3A8xpJosnhvPiIiIqIUmJkioh62PJTL5SIz\nGsYYNJtNAEChUIgd2GyzDSe9bj+/1Wq5z2y1Wm55uFwV/mw6GWNJNBnsTBGRY4xxJ8x2u41crpO8\nDpdo7PKTrhJLcgLWl8fbclCr1XLv1SdfnnQHw1gSTQ7LfEREREQpMDNFRAAeZjJ0FsJmE3RGQmc5\nksxd1G+b9jOjnut1ws8pHmNJNFnsTBGRo0tD4ZOdPhlHnZiH3V74edyEjjz5DoaxJJoclvmIiIiI\nUmBmioh66KxBVGYj6vVRbCtqGTMY6TCWRJPBzhQRxdJjbiZZruGJd/QYS6LxYZmPiIiIKAVmpoio\nr3B2IW1JqN/n0/gwlkSjx84UEQ2MJ0x/MJZE6bHMR0RERJQCO1NEREREKbAzRURERJQCO1NERERE\nKbAzRURERJQCO1NEREREKbAzRURERJRCos6UiKyJyJdF5H9F5GUR+RUROSsiz4nIq93/ro+7sZQe\nY+kXxtMfjKVfGM/5kjQz9RkA/26M+QUAvwjgZQBXATxvjHkbgOe7f1P2MZZ+YTz9wVj6hfGcJ8aY\nEx8AzgD4PwASWv4KgIvd5xcBvJLgswwfU38wlv48uG/69WAs/Xlw3/To0S9GxphEmaknANwF8I8i\n8h0R+ZyIlABsGmNudde5DWAzwWfR9DGW/uC+6RfG0h/cN+dMks5UAcC7AHzWGPNOABWEUpOm0302\nUW8WkSsi8oKIvJC2sTQSjKU/uG/6hbH0B/fNeZMgxXgBwI/V3+8FcB1MV87qg7H058F9068HY+nP\ng/umR4+RlPmMMbcB/FREfr676DcAvATgWQBPdZc9BWCn32dRJjCWnuC+6R3G0hPcN+ePdHu+J68k\n8g4AnwNQBPA6gD9Ep0T4JQBvBvATAB8yxuz1+Zz+G6NxeycYSy8YY4T7ple4b3qC+6ZfjDHSb51E\nnalR4Zdi+pJ8KZJgLKdvVLEEGM8s4L7pD+6bfkkSz8IkGqLcQ2cg3r0Jb7ef88hem4DRt+tnRvhZ\nWY0lkM14ZjmWAHCEzniOrMliLIFsx5P75mCyHEuA++agphLPiWamAEBEXjDG/PJEN9pHFtsEZLdd\nVlbbl8V2ZbFNWlbbx3YNJ6vty2K7stgmLavtY7t68d58RERERCmwM0VERESUwjQ6U9emsM1+stgm\nILvtsrLaviy2K4tt0rLaPrZrOFltXxbblcU2aVltH9ulTHzMFBEREZFPWOYjIiIiSmFinSkReb+I\nvCIir4nI1f7vGFs7HheRr4nISyLyooh8tLv8kyLyhoh8t/t4csLt+rGI/KC77Re6y86KyHMi8mr3\nv+uTbNNJshDPrMay24aZiSdj2bdtjOXg7WA8RyAL8WQsE0pyz5m0DwB5AD8C8FZ0ZoP9HoC3T2Lb\nEW25COBd3eerAH4I4O0APgng49NoU7ctPwZwPrTsrwBc7T6/CuBT02pfFuOZ1VjOUjwZS8aS8WQ8\nGcv0j0llpt4N4DVjzOvGmAaAZwBsT2jbPYwxt4wx3+4+LwN4GcBj02hLAtsAnu4+fxrAB6fYFi0T\n8ZyxWALZjCdjORzG8gSM50hkIp6MZTKT6kw9BuCn6u8byEAwROQt6NwP6xvdRR8Rke+LyBemkOY1\nAP5TRL4lIle6yzaNMbe6z28D2Jxwm+JkLp4ZiyUwO/FkLPtjLFNgPIeWuXgylvHmdgC6iJwC8BUA\nHzPGHAL4LICfBfAOALcA/PWEm/QeY8y7APw2gD8WkV/TL5pOzpKXXkbIYCwBxnMojKVfGE9/MJYn\nm1Rn6g0Aj6u/L3eXTYWILKDzpfiiMearAGCMuWOMaRlj2gD+AZ0U68QYY97o/ncXwL92t39HRC52\n23wRwO4k23SCzMQzi7HstmFW4slY9sFYDofxTC0z8WQs+5tUZ+qbAN4mIk+ISBHAhwE8O6Ft9xAR\nAfB5AC8bYz6tll9Uq/0ugP+ZYJtKIrJqnwP4re72nwXwVHe1pwDsTKpNfWQinlmMZXf7sxRPxvLk\ndjGWQ2A8RyIT8WQsE5rEKPdOpg1PonMVwI8A/PmkthvRjvegk/b7PoDvdh9PAvgnAD/oLn8WwMUJ\ntumt6Fyp8T0AL9p/HwDnADwP4FUA/wXg7LT+3bIYzyzGchbjyVgylown48lYpntwBnQiIiKiFOZ2\nADoRERHRKLAzRURERJQCO1NEREREKbAzRURERJQCO1NEREREKbAzRURERJQCO1NEREREKbAzRURE\nRJTC/wMN9uRYKVVqVQAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f384e8adf10>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(10, 2))\n",
    "for i in range(5):\n",
    "    plt.subplot(1, 5, i + 1), plt.imshow(list_imgs[i], cmap=plt.cm.Greys_r)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## BPDL"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "from bpdl import dictionary_learning as dl\n",
    "from bpdl import pattern_atlas as ptn_dict\n",
    "\n",
    "init_atlas_msc = ptn_dict.init_atlas_mosaic(list_imgs[0].shape, 6)\n",
    "# init_encode_rnd = ptn_weigth.initialise_weights_random(len(imgs), np.max(atlas))\n",
    "\n",
    "atlas, w_bins, deforms = dl.bpdl_pipeline(list_imgs, init_atlas=init_atlas_msc, max_iter=9)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**show the estimated components - dictionary**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAP4AAAEICAYAAAB/KknhAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAEiBJREFUeJzt3X+QXWV9x/H3xxASfocIjSFhDBYU\nU0eCXfkxUlQiiJRK/mAYlJG1E2dnWtuCtRXQ0SmO46CdEe3o0G4BiS0KCGgwY5UY6ei0Ekgk/IxI\nhCBJE1aENEBrSPDbP86zzmW7P87ee37s7vN5zdzZe8499z7f5O7nPs9z7tlzFBGYWV5e1XYBZtY8\nB98sQw6+WYYcfLMMOfhmGXLwzTLk4E9xkrZKelfbdYxF0sclXdt2HTY5Dn4LJP27pOckzRmx/gZJ\nn2mrrm5ExGcj4kNNtyvp3yS90HF7SdKDTdcxXTn4DZO0BPgjIID3tlrMNBYR74mIg4dvwH8C32y7\nrunCwW/excDdwA1A//BKSQPARcDHUg/2nZFPlHSSpJ9I2iVph6QvS9o/PSZJV0sakrRb0oOS3pQe\ne7+kB8YqaOR0QtLfSfrXdH+JpJDUL+mXkp6R9InRtk3LH5D0pKRfS/pE52uPHNFIeoekbR3LR0m6\nTdKvJD0h6a/K/Id2fJh+rcz25uC34WLgxnR7t6QFABExmNZ9PvVifzLKc18GPgIcAZwKLAf+PD12\nFnA68HrgMOAC4Nfptb8eEW/use7TgDekNj8l6Y0jN5C0FLgG+ABwFPBqYHGZF5f0KuA7wP3AotTO\npZLeXeLpFwM/joitZdoyB79Rkk4DXgvcEhEbgV8A7y/7/IjYGBF3R8S+9Ev+T8Db08N7gUOA4wFF\nxOaI2FFh+VdGxP9GxP0U4TxhlG3OB9ZExI8iYg/wSeC3JV//rcCREfHpiHgpIh4H/hm4sMRzL6YY\nQVlJDn6z+oE7I+KZtPx1Oob7E5H0eklrJO2UtBv4LEXvT0T8EPgy8BVgSNKgpEMrrH1nx/3/AQ4e\nZZujgKeGFyLiRdKoo4TXAkelacwuSbuAjwMLxntS+jB9DXBryXYMB78xkg6gGH6/PQV3J8Ww/QRJ\nw73nRH8qeQ3wM+C4iDiUIhgafjAi/iEi/hBYSjHk/9uS5b0IHNix/JqSzxtpB3D08IKkAymG+2Xa\neQp4IiLmddwOiYhzJmizH7g9Il7osuYsOfjNWUExR18KLEu3NwI/phiqAjwNvG6c1zgE2A28IOl4\n4M+GH5D0VkknS5pNEbDfkIbZkj4oaes4r7sJuFDSbEl9FEP2btwKnCvptLTT8dO88ndsE3COpPmS\nXgNc2vHYPcDzki6TdICkWZLeJOmtYzXW8WF6Q5f1ZsvBb04/8NWI+GVE7By+UQzPL5K0H3AdsDQN\ndb89ymv8DcU+gecp5r83dzx2aFr3HPAkxRD779NjRwP/MU5tnwR+Pz33SoopyKRFxMPAh9Pzd6TX\n29axyb9Q7B/YCtzZWX9EvAycS/GB+ATwDHAtxY7KsawAdgF3dVNvzuQTccx8ku4ELomIzS20vRX4\nUET8oOm2bWz7tV2A1S8izmq7BptaPNQ3y5CH+mYZ6qnHl3S2pEclbZF0eVVFmVm9uu7xJc0Cfg6c\nSbHn9l7gfRHxyFjP2V9zYi4HddWemU3sN7zIS7FHE23Xy869k4At6dBKJN0EnAeMGfy5HMTJWt5D\nk2Y2nvWxrtR2vQz1F9FxeCZFr79o5EaSBiRtkLRhL3t6aM7MqlL7Xv2IGIyIvojom82ciZ9gZrXr\nJfjb6Tgum+LPL7f3Vo6ZNaGX4N8LHCfpmHRc9oXAHdWUZWZ16nrnXkTsk/QXwPeBWcD16VhtM5vi\nejpkNyK+C3y3olrMrCE+ZNcsQw6+WYYcfLMMOfhmGXLwzTLk4JtlyME3y5CDb5YhB98sQ1mebPP7\n/7Wp1HbvPmpZzZWYtcM9vlmGHHyzDGUz1C87vB/vOR7620zhHt8sQw6+WYYcfLMMOfhmGXLwzTLk\n4JtlyME3y5CDb5YhB98sQw6+WYayOWTX8vDMwKmvWD5i8CctVTK1ucc3y9CEwZd0vaQhSQ91rJsv\naa2kx9LPw+st08yqpIgYfwPpdOAF4GsR8aa07vPAsxFxlaTLgcMj4rKJGjtU8+NkLa+g7Mnr5q/z\nxuO/1KvfyGF7r3IY9q+PdeyOZzXRdhP2+BHxI+DZEavPA1al+6uAFZOu0Mxa0+3OvQURsSPd3wks\nGGtDSQPAAMBcDuyyOTOrUs979SMiJI05X4iIQWAQiqF+r+21ycP7elU9tLexdbtX/2lJCwHSz6Hq\nSjKzunUb/DuA/nS/H1hdTTlm1oQyX+d9A/gJ8AZJ2yStBK4CzpT0GPCutGxm08SEc/yIeN8YD7Xz\nvZzNKE3O6zvbyuGrvfH4yD2zDDn4Zhma8Mi9KrV55F4nX0KrPVP1K7uZMvSv7Mg9M5t5HHyzDDn4\nZhnK8kQcVc/dt1x9SqWvNxnHfuTu1toua6rO63PmHt8sQw6+WYay/DqvW20O6bsxFacB02HYP52/\n2vPXeWY2JgffLENZ7tUva7oN7Ucar/6pOA2w5rjHN8uQg2+WIQffLEOe449j5Dx4Os/5a5nTn/Lm\nctvd/UD1bVtP3OObZcjBN8uQh/rjmM5D+5FG/lu6GvqXHdrblOce3yxDDr5Zhhx8swx5jp8JH6I7\nvun8F3ndcI9vlqEyl9A6WtJdkh6R9LCkS9L6+ZLWSnos/Ty8/nLNrAplhvr7gI9GxE8lHQJslLQW\n+CCwLiKuknQ5cDlwWX2l2mS1OrzvOFpvOpx8IzcT9vgRsSMifpruPw9sBhYB5wGr0margBV1FWlm\n1ZrUzj1JS4ATgfXAgojYkR7aCSwY4zkDwADAXA7stk4zq1DpnXuSDgZuAy6NiN2dj0Vx4r5RT94X\nEYMR0RcRfbOZ01OxZlaNUj2+pNkUob8xIm5Pq5+WtDAidkhaCAzVVWRbZtJf51Vi5F/ZjXcIb8dj\nnV+Veb4/NZTZqy/gOmBzRHyh46E7gP50vx9YXX15ZlaHMj3+24APAA9KGr7M7MeBq4BbJK0EngQu\nqKdEM6uaz6tfAV9Cq3ttDf1n6pF6Pq++mY3JwTfLkIf6LfJ571+p7mH/TB3ed/JQ38zG5OCbZcjB\nN8uQ5/g2ZZWd8+cwdy/Lc3wzG5ODb5Yhn3PPpiwP4evjHt8sQw6+WYYcfLMMOfhmGXLwzTLk4Jtl\nyME3y5CDb5YhB98sQw6+WYYcfLMMOfhmGXLwzTLk4JtlyME3y1CZa+fNlXSPpPslPSzpyrT+GEnr\nJW2RdLOk/esv18yqUKbH3wOcEREnAMuAsyWdAnwOuDoijgWeA1bWV6aZVWnC4EfhhbQ4O90COAO4\nNa1fBayopUIzq1ypOb6kWelKuUPAWuAXwK6I2Jc22QYsGuO5A5I2SNqwlz1V1GxmPSoV/Ih4OSKW\nAYuBk4DjyzYQEYMR0RcRfbOZ02WZZlalSe3Vj4hdwF3AqcA8ScMn61wMbK+4NjOrSZm9+kdKmpfu\nHwCcCWym+AA4P23WD6yuq0gzq1aZ02svBFZJmkXxQXFLRKyR9Ahwk6TPAPcB19VYp5lVaMLgR8QD\nwImjrH+cYr5vZtOMj9wzy5CDb5YhB98sQw6+WYYcfLMMOfhmGfJlsmtw7sPP/e7+mj84vMVKzEbn\nHt8sQw6+WYayHOp3DsVH6hyaj7ddt2156G9TgXt8sww5+GYZcvDNMpTNHL/sfL2KeX3Z1/d839ri\nHt8sQw6+WYZm7FC/7iF7FfxVn7XFPb5Zhhx8sww5+GYZcvDNMuTgm2XIwTfLkINvliEH3yxDpYOf\nLpV9n6Q1afkYSeslbZF0s6T96yvTzKo0mSP3LqG4WOahaflzwNURcZOkfwRWAtdUXF/XRh4FNxWP\n5PORetaWUj2+pMXAHwPXpmUBZwC3pk1WASvqKNDMqld2qP9F4GPAb9Pyq4FdEbEvLW8DFo32REkD\nkjZI2rCXPT0Va2bVmDD4ks4FhiJiYzcNRMRgRPRFRN9s5nTzEmZWsTJz/LcB75V0DjCXYo7/JWCe\npP1Sr78Y2F5fmb2r4iSaVb9Gk/77u8eW3vawc7bUWIlNBRP2+BFxRUQsjoglwIXADyPiIuAu4Py0\nWT+wurYqzaxSvXyPfxnw15K2UMz5r6umJDOrmyKiscYO1fw4Wcsbay93kxnej8XD/ullfaxjdzyr\nibbzkXtmGXLwzTI0Y8+5l6sqhvdjvZ6H/TOHe3yzDDn4Zhly8M0y5OCbZcjBN8uQg2+WIQffLEMO\nvlmGHHyzDDn4Zhly8M0y5OCbZcjBN8uQg2+WIQffLEMOvlmGHHyzDDn4Zhly8M0y5OCbZcjBN8uQ\ng2+WoVKn15a0FXgeeBnYFxF9kuYDNwNLgK3ABRHR3ZUkzaxRk+nx3xkRyyKiLy1fDqyLiOOAdWnZ\nzKaBXob65wGr0v1VwIreyzGzJpQNfgB3StooaSCtWxARO9L9ncCC0Z4oaUDSBkkb9rKnx3LNrApl\nL6F1WkRsl/R7wFpJP+t8MCJC0qiX3Y2IQWAQiqvl9lStmVWiVI8fEdvTzyHgW8BJwNOSFgKkn0N1\nFWlm1Zow+JIOknTI8H3gLOAh4A6gP23WD6yuq0gzq1aZof4C4FuShrf/ekR8T9K9wC2SVgJPAhfU\nV6aZVWnC4EfE48AJo6z/NbC8jqKse52Xsq7iktm+NPbM5CP3zDLk4JtlyME3y1DZ7/FtGup2vu95\n/cznHt8sQw6+WYY81M+Eh+/WyT2+WYYcfLMMOfhmGXLwzTLk4JtlyME3y5CDb5YhB98sQw6+WYYc\nfLMMOfhmGXLwzTLk4JtlyME3y5CDb5YhB98sQw6+WYYcfLMMlQq+pHmSbpX0M0mbJZ0qab6ktZIe\nSz8Pr7tYM6tG2R7/S8D3IuJ4istpbQYuB9ZFxHHAurRsZtNAmavlHgacDlwHEBEvRcQu4DxgVdps\nFbCiriLNrFplevxjgF8BX5V0n6Rr0+WyF0TEjrTNToqr6v4/kgYkbZC0YS97qqnazHpSJvj7AW8B\nromIE4EXGTGsj4gAYrQnR8RgRPRFRN9s5vRar5lVoEzwtwHbImJ9Wr6V4oPgaUkLAdLPoXpKNLOq\nTRj8iNgJPCXpDWnVcuAR4A6gP63rB1bXUqGZVa7slXT+ErhR0v7A48CfUnxo3CJpJfAkcEE9JZpZ\n1UoFPyI2AX2jPLS82nLMrAk+cs8sQw6+WYYcfLMMOfhmGXLwzTLk4JtlSMXRtg01Jv2K4jv/I4Bn\nGmt4dFOhBnAdI7mOV5psHa+NiCMn2qjR4P+uUWlDRIx2XEBWNbgO19FWHR7qm2XIwTfLUFvBH2yp\n3U5ToQZwHSO5jleqpY5W5vhm1i4P9c0y5OCbZajR4Es6W9KjkrZIauysvJKulzQk6aGOdY2fHlzS\n0ZLukvSIpIclXdJGLZLmSrpH0v2pjivT+mMkrU/vz83p/Au1kzQrnc9xTVt1SNoq6UFJmyRtSOva\n+B1p5FT2jQVf0izgK8B7gKXA+yQtbaj5G4CzR6xr4/Tg+4CPRsRS4BTgw+n/oOla9gBnRMQJwDLg\nbEmnAJ8Dro6IY4HngJU11zHsEopTtg9rq453RsSyju/N2/gdaeZU9hHRyA04Ffh+x/IVwBUNtr8E\neKhj+VFgYbq/EHi0qVo6algNnNlmLcCBwE+BkymOENtvtPerxvYXp1/mM4A1gFqqYytwxIh1jb4v\nwGHAE6Sd7nXW0eRQfxHwVMfytrSuLaVOD14XSUuAE4H1bdSShtebKE6Suhb4BbArIvalTZp6f74I\nfAz4bVp+dUt1BHCnpI2SBtK6pt+Xnk5lPxneucf4pwevg6SDgduASyNidxu1RMTLEbGMosc9CTi+\n7jZHknQuMBQRG5tuexSnRcRbKKaiH5Z0eueDDb0vPZ3KfjKaDP524OiO5cVpXVtaOT24pNkUob8x\nIm5vsxaAKK6KdBfFkHqepOHzMDbx/rwNeK+krcBNFMP9L7VQBxGxPf0cAr5F8WHY9PvS2Knsmwz+\nvcBxaY/t/sCFFKfobkvjpweXJIpLkW2OiC+0VYukIyXNS/cPoNjPsJniA+D8puqIiCsiYnFELKH4\nffhhRFzUdB2SDpJ0yPB94CzgIRp+X6LJU9nXvdNkxE6Kc4CfU8wnP9Fgu98AdgB7KT5VV1LMJdcB\njwE/AOY3UMdpFMO0B4BN6XZO07UAbwbuS3U8BHwqrX8dcA+wBfgmMKfB9+gdwJo26kjt3Z9uDw//\nbrb0O7IM2JDem28Dh9dRhw/ZNcuQd+6ZZcjBN8uQg2+WIQffLEMOvlmGHHyzDDn4Zhn6P0BQVOzB\nVN2HAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f384e79ebd0>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.imshow(atlas, interpolation='nearest')\n",
    "_ = plt.title('Atlas; unique %i' % len(np.unique(atlas)))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**particular coding of each sample**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAaQAAADUCAYAAADTGVdAAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAAG+VJREFUeJzt3Xm8X3V95/HXm7CkRBQliKxCNfgY\nRECaBi10CqIQIgOdGZektS5Dm+qIg7Zji1MHWjozdak6Kgw0SmQpgrjgpBoJKaLIKJgQ2cLWNEVJ\nYIgssjRFknvf/eOcqz8uv+0uv9859/d7Px+P87hn+Z7v+f5y783nfr/nu8g2ERERVduh6gJERERA\nAlJERNREAlJERNRCAlJERNRCAlJERNRCAlJERNRCAlJEREyYpOWStki6o8V1SfqMpA2SbpN0ZKc8\nE5AiImIyLgIWtrl+EjCv3JYC53fKMAEpIiImzPb1wKNtkpwKXOLCjcDukvZul2cCUkRE9MK+wP0N\nx5vKcy3t2NPitCBpIfBpYBbwedsfGXd9F+AS4NeAR4C32r6vU747axfPZs70FzgiZqSDD9tadRGm\nzX33b+PhR0c0XfmdeNwcP/LoSMvrN9/28/XA0w2nltleNl3Pb6bvAUnSLOA84A0UEXONpBW272xI\ndhrwmO2XS1oMfBR4a6e8ZzOHo3R8L4odETPQqlW3VF2EabPgxPs7J5qAhx/dzvevbl1hmb3PPz1t\ne/4UHrEZ2L/heL/yXEtVNNktADbY3mj7GeAKirbGRqcCF5f7XwGOlzRtfxlERAw7A6O45TYNVgBv\nL3vbvQZ43PaD7W6oosmuWbviUa3S2N4u6XFgD+DhvpQwImLAGbPNrZvsOpF0OXAsMFfSJuBsYCcA\n2xcAK4FFwAZgK/CuTnlW8g5pOklaStGlkNnsWnFpIiJmjqnUhGwv6XDdwHsnkmcVTXbdtCv+Io2k\nHYEXUHRueA7by2zPtz1/J3bpQXEjIgaPgW2MttyqUEVAWgPMk3SQpJ2BxRRtjY1WAO8o998EfNtZ\nSTAiYtoYGLFbblXoe5Nd+U7odGAVRbfv5bbXSzoHWGt7BXAhcKmkDRQDrxb3u5wREYPMmG3T03lh\n2lTyDsn2SooXXo3nzmrYfxp4c7/LFRExLGzYVq94NPM7NURExGSIEeo1miYBKSJiCBnY5gSkiIio\nmCE1pIiIqF5RQ6rX/NoJSBERQ8iIkZot+JCAFBExhFJD6rGDD9s6ULP7nrjPEVUXIYbMqgcG5/cH\nBut36F43naxmCsRIAlJERFStmDpoVtXFeJYEpIiIIWSLbU5AioiIihXdvtNkFxERFTNim+sVAupV\nmoiI6JuRms3U0Pf6mqT9JV0n6U5J6yWd0STNsZIel3RLuZ3VLK+IiJicsRpSq60KVTx1O/DHttdJ\n2g24WdJq23eOS/c92ydXUL6IiIGXd0iA7QeBB8v9JyXdBewLjA9IERHRI0UNKb3sfkHSgcCrgZua\nXH6tpFuBB4D/ant9izyWAksBZrPrQA2Ei4jolXT7biDpecBXgffbfmLc5XXAS20/JWkR8HVgXrN8\nbC8DlgE8Xy+q2XJTERH1VCxhXq8mu0pKI2knimB0me2vjb9u+wnbT5X7K4GdJM3tczEjIgbWWJNd\nq60Kfa8hSRJwIXCX7U+2SPMS4CHblrSAInBO90ROERFDbeg7NQBHA78H3C5pbCbH/wYcAGD7AuBN\nwHskbQf+BVhsO81xERHTJJ0aANs3QPtlCm2fC5zbnxJFRAwfA6NTeIckaSHwaWAW8HnbHxl3/QDg\nYmD3Ms2Z5SuYljJTQ0TEEJpKLztJs4DzgDcAm4A1klaMG0/6YeBK2+dLOgRYCRzYLt8EpIiIIVQs\n0DfpJrsFwAbbGwEkXQGcyrPHkxp4frn/AoohPG0lINVYFkuLfsv3aJh0XKBvrqS1DcfLymE2UExm\ncH/DtU3AUePu/3PgGknvA+YAr+9UogSkiIgh1EUN6WHb86fwiCXARbY/Iem1wKWSDrU92uqGBKSI\niCFkxOjkZ/veDOzfcLxfea7RacBCANs/kDQbmAtsaZVpvTqhR0REX9hMZWDsGmCepIMk7QwsBlaM\nS/MT4HgASf8GmA38tF2mqSFFRAypydaQbG+XdDqwiqJL93Lb6yWdA6y1vQL4Y+Bzkj5A0UL4zk7j\nSROQIiKG0FQHxpZjilaOO3dWw/6dFBMhdC0BKSJiCBmxfXTIZ2qIiIh6GG0/aU7fDVRAOviwraxa\nNThjdzImJCJ6xYZtqSEVJN0HPAmMANvH93cvZwX/NLAI2ErxQmxdv8sZETGIptjtuyeqriEdZ/vh\nFtdOoliUbx7FCODzee5I4IiImAQD22u2QF/VAamdU4FLym6CN0raXdLeth+sumAREYNgKrN990KV\npTHFPEc3S1ra5HqzuZL27UvJIiIGnC22e4eWWxWqrCEdY3uzpBcDqyXdbfv6iWZSBrOlAAfsW+cK\nX0REfRjYPpoaEgC2N5dftwBXUUxn3qibuZKwvcz2fNvz99yjXj1GIiLqbNRquVWhkoAkaY6k3cb2\ngROAO8YlWwG8XYXXAI/n/VFExPQwabIbsxdwVdGzmx2BL9q+WtK7AWxfQDElxSJgA0W373d1yvTe\n23bN2J0aG7T1nQZRfn+GiCc/l12vVBKQylUGD29y/oKGfQPv7We5IiKGRR3fIaUXQETEEMrA2IiI\nqI0OS5j3XQJSRMQQsmEkTXYREVG9NNlFREQNmNSQIiKiDlw029XJQAWkrIdUb4P2eTKuKmYyk04N\nERFRC3mHFBERNTE6moAUEREVS7fviIiojbp1auh7eJT0Ckm3NGxPSHr/uDTHSnq8Ic1Z/S5nRMQg\nM2J0dIeWWyeSFkq6R9IGSWe2SPMWSXdKWi/pi53y7HsNyfY9wBEAkmZRrHF0VZOk37N9cj/LFhEx\nTCZbQSr/7z4PeAPFat5rJK2wfWdDmnnAh4CjbT9WLsbaVlc1JEkfk/R8STtJulbSTyW9bXIf5VmO\nB/7R9o+nIa+IiOiWwaNquXWwANhge6PtZ4ArgFPHpfkD4Dzbj8EvFmNtq9smuxNsPwGcDNwHvBz4\nYJf3trMYuLzFtddKulXStyS9chqeFRERDWy13DrYF7i/4XhTea7RwcDBkv6fpBslLeyUabdNdmPp\n3gh82fbj5eJ6kyZpZ+AUiirdeOuAl9p+StIi4OvAvBb5LAWWAhyw72D10Ri0gZeDNjB20D5PDBfT\nsdv3XElrG46X2V42gUfsSPH/9rHAfsD1kl5l+2ftbujGNyTdDfwL8B5JewJPT6BgzZwErLP90PgL\nZW1sbH+lpP8jaa7th5ukXQYsA5h/+Oya9RmJiKipssmujYdtz29xbTOwf8PxfuW5RpuAm2xvA/5J\n0r0UAWpNqwd21WRn+0zgN4D5ZeZbeW574UQtoUVznaSXqKyCSVpQlvORKT4vIiIauc3W3hpgnqSD\nytauxcCKcWm+TlE7QtJciia8je0y7bZTw67AfwbOL0/tA7SKnN3kN4eid8bXGs69W9K7y8M3AXdI\nuhX4DLC4XNI8IiKmResODZ06NdjeDpwOrALuAq60vV7SOZJOKZOtAh6RdCdwHfBB220rFt022X0B\nuJmilgRF1ezLwDe6vP9ZbP8zsMe4cxc07J8LnDuZvCMiogumm84LrW+3VwIrx507q2HfwB+VW1e6\n7WX3MtsfA7aVD9oK1GsSpIiImBir9VaBbmtIz0j6FcqWRUkvA37es1JFRETv1exFSLcB6WzgamB/\nSZcBRwPv7FWhIiKixwzMxNm+ba+WtA54DUVT3RnNumBX7d7bds3YkIj4hUEay7fgxK3TnqdHpz3L\nKWkbkCQdOe7Ug+XXAyQdYHtdb4oVERE9N8MW6PtEm2sGXjeNZYmIiH4xaCbVkGwf16+CREREP1XX\nm66Vrt4hSZpNMTD2GIqa0feAC2xPdfqgiIioykyqITW4BHgS+Gx5/DvApcCbe1GoiIjogxna7ftQ\n24c0HF9XTgcREREzkUE16/bd7UwN6yS9ZuxA0lHA2jbpIyKi7iY/uWpPdFtD+jXg+5J+Uh4fANwj\n6XaKKYsO60npYqAM0piQQTVo4/gG6fPc235e0knRDG2y67jSXzOSllOsMrvF9qHluRcBXwIOpFh9\n9i1jS9yOu/cdwIfLw/9h++LJlCEiIpqo4UwN3a6H9GPgCeAFFLN07wHsYfvH5bVWLuK5wexM4Frb\n84Bry+NnKYPW2cBRFGu3ny3phd2UNSIiujQTm+wk/SXF3HX/yC+L2nFgrO3rJR047vSplIs2ARcD\n3wH+dFyaE4HVth8tn7+aIrA1XdAvIiImbkYNjG3wFoolKJ6ZhmfuZXtsCqL/D+zVJM2+wP0Nx5vK\nc88haSmwFGA2u05D8SIihkTN3iF128vuDmD36X54uYDTlP5JbC+zPd/2/J3YZZpKFhEx2FR2+261\nVaHbGtJfAT+SdAcN6yDZPqX1LS09JGlv2w9K2hvY0iTNZn7ZrAewH0XTXkRETJOZ2mR3MfBR4Ham\nPtnECuAdwEfKr/+3SZpVwP9q6MhwAvChKT43IiIa1azJrtuAtNX2ZyaauaTLKWo6cyVtoug59xHg\nSkmnAT+meD+FpPnAu23/vu1Hy44Ua8qszhnr4DBMMm6n3gZpjEsMoZk223eD70n6K4raTWOTXdv1\nkGwvaXHp+CZp1wK/33C8HFjeZfkiImKiZmgN6dXl19c0nMt6SBERM9iMnKkh6yJFRAygmRiQACS9\nEXglMHvsnO1zelGoiIjosRq+Q+pqHJKkC4C3Au8DRLEO0kt7WK6IiOghUQSkVlvH+6WFku6RtEHS\nc6aAa0j3HyW57LjWVrcDY3/D9tuBx2z/BfBa4OAu742IiDqa5Fx2kmYB5wEnAYcASyQd0iTdbsAZ\nwE3dFKfbgDS2VPlWSfsA24G9u7w3IiLqxlOqIS0ANtjeWE4pdwXFPKXj/SXFGNanm1x7jm7fIf2d\npN2BjwPrio/C57q8t28OPmwrq1YNztidQRvnknFV0W+D9DO34MSt05/p5Ds1NJtv9KjGBJKOBPa3\n/U1JH+wm024D0t3AiO2vltWyI4Gvd3lvRETUUIea0FxJjSuDL7O9rKt8pR2AT1KsEtG1bgPSf7f9\nZUnHUIw9+mvgfMZFxIiImCE6vyt62Harjgibgf0bjvcrz43ZDTgU+I4kgJcAKySdUk6C0FS375BG\nyq9vBD5n+5vAzl3eGxERNTSFd0hrgHmSDpK0M7CYYiYfAGw/bnuu7QNtHwjcCLQNRtB9QNos6W8o\nun6vlLTLBO6NiIgakltv7djeDpxOMRH2XcCVttdLOkfSZFaBACa2QN9C4K9t/6xcNqLtSypJy4GT\ngS22Dy3PfRz4d8AzFKvPvsv2z5rcex/wJEXNbHubamNEREyGmdLaDbZXAivHnTurRdpju8mzq1qO\n7a22v2b7H8rjB21f0+G2iyiCWKPVwKG2DwPupf2SEsfZPiLBKCJi+onJ15B6pWfNbravBx4dd+6a\nsqoHRZvifr16fkREtDc0AakL/wn4VotrBq6RdLOkpX0sU0TE8Bhts1Wg68lVp5OkP6OY7eGyFkmO\nsb1Z0ouB1ZLuLmtczfJaCiwFmM2uAzeYNCImb5D+P7jXj0xvhhXWhFrpew1J0jspOjv8ru2m/xy2\nN5dftwBXUUxT0ZTtZbbn256/E7v0oMQREYNpKpOr9kJfA5KkhcCfUPRHbzoPhqQ55YR8SJoDnADc\n0b9SRkQMiUlOrtorPQtIki4HfgC8QtImSacB51KM4F0t6ZZyWQsk7SNprPvgXsANkm4Ffgh80/bV\nvSpnRMRQmtrkqj3Rs3dItpc0OX1hi7QPAIvK/Y3A4b0qV0RE/HI9pDqppFNDRETUQM06NSQgRUQM\nI4NG6xWREpAiIoZU3bp9JyBF3wzSmBAYrMXfxgza9yjayzukiIioh9SQIiKick4NKSIiaqDo9l2v\nKlICUkTEkEqnhoiIqJ5BI1UX4tkSkCIihlVqSBERUbkMjO2tgw/byqpVgzM2JGNC6i3fn5jp6vYO\nqZezfS+XtEXSHQ3n/lzS5nKm71skLWpx70JJ90jaIOnMXpUxImJYjU2uWqfZvnu5HtJFwMIm5z9l\n+4hyWzn+oqRZwHnAScAhwBJJh/SwnBERw8dGo623KvQsIJVLjj86iVsXABtsb7T9DHAFcOq0Fi4i\nIoZngb42Tpd0W9mk98Im1/cF7m843lSea0rSUklrJa396SM168MYEVFXBo245VaFfgek84GXAUcA\nDwKfmGqGtpfZnm97/p57zJpqdhERw2MKNaRO7/ol/ZGkO8sKyLWSXtopz74GJNsP2R6xPQp8jqJ5\nbrzNwP4Nx/uV5yIiYhpN9h1Sl+/6fwTMt30Y8BXgY53K09eAJGnvhsN/D9zRJNkaYJ6kgyTtDCwG\nVvSjfBERw0RuvXXQ8V2/7etsby0Pb6SoXLTVs3FIki4HjgXmStoEnA0cK+kIigrhfcAflmn3AT5v\ne5Ht7ZJOB1YBs4Dlttd388x7b9t1oMaGDNp6O4P0vYmZYZB+hxacuLVzoglQ54GxcyWtbTheZntZ\nud/sXf9RbfI6DfhWpzL1LCDZXtLk9IUt0j4ALGo4Xgk8p0t4RERMnw6dFx62PX/Kz5DeBswHfqtT\n2oGaqSEiIro0te7dXb3rl/R64M+A37L9806ZVtHtOyIiKjelgbEd3/VLejXwN8Aptrd0U6LUkCIi\nhpUnV0Vq9a5f0jnAWtsrgI8DzwO+LAngJ7ZPaZdvAlJExDByx3dI7W9v8q7f9lkN+6+faJ4JSBER\nw6pms30nIEVEDCmNVjStdwsJSDWWcTv1NkhjXMYM2s/cIH2ee/3ItOYnVzdnXSsJSBERw2qSnRp6\nJQEpImIYGUgNKSIi6kCpIUVERPUMw9KpQdJy4GRgi+1Dy3NfAl5RJtkd+Jnt57x1lHQf8CQwAmyf\njvmUIiKigRmqd0gXAecCl4ydsP3WsX1JnwAeb3P/cbYf7lnpIiKG3ND0srN9vaQDm11TMY/EW4DX\n9er5ERHRhoGRejXZVTW56m8CD9n+hxbXDVwj6WZJS9tlJGmppLWS1m6j42SyEREBgIsmu1ZbBarq\n1LAEuLzN9WNsb5b0YmC1pLttX98sYblg1DKA5+tF9ap/xkAbpEGXYwZtsO8gfo+mVc06NfS9hiRp\nR+A/AF9qlcb25vLrFuAqiuVyIyJiuhgYdeutAlU02b0euNv2pmYXJc2RtNvYPnACcEcfyxcRMQQM\noyOttwr0LCBJuhz4AfAKSZsknVZeWsy45jpJ+0gam8Z8L+AGSbcCPwS+afvqXpUzImIo1bCG1Mte\ndktanH9nk3MPAIvK/Y3A4b0qV0RElGr2DikzNUREDCMbRqppmmslASkiYlgN0UwNERFRW67dwNgE\npBrLmJCI6BmDnYAUERF1kBpSRERUzkO0/ERERNSba9bLrqrJVSMiokouOzW02jqQtFDSPZI2SDqz\nyfVdJH2pvH5Tq9UfGiUgRUQMK4+23tqQNAs4DzgJOARYIumQcclOAx6z/XLgU8BHOxUnASkiYgjZ\nxiMjLbcOFgAbbG+0/QxwBXDquDSnAheX+18Bji/XwmspASkiYkh51C23DvYF7m843lSea5rG9naK\nFcL3aJfpQHVqeJLHHv57f+XHfXjUXKDny6vP2rvXT/iFvnwe2ND7RxT69Hn6pm+fp08/c338/vTl\nZ65fn+el05nZkzy26u9Hr5zbJslsSWsbjpeV68/1zEAFJNt79uM5ktbant+PZ/VDPk+95fPU20z9\nPLYXTuH2zcD+Dcf7leeapdlUroP3AuCRdpmmyS4iIiZqDTBP0kGSdqZYVmjFuDQrgHeU+28Cvm23\nnzxvoGpIERHRe7a3SzodWAXMApbbXi/pHGCt7RXAhcClkjYAj1IErbYSkCanp+2oFcjnqbd8nnob\ntM/TFdsrgZXjzp3VsP808OaJ5KkONaiIiIi+yDukiIiohQSkCeg0VcZMI2m5pC2S7qi6LFMlaX9J\n10m6U9J6SWdUXaapkDRb0g8l3Vp+nr+oukzTQdIsST+S9I2qyzIdJN0n6XZJt4zrIh2TkCa7LpVT\nZdwLvIFiENgaYIntOyst2BRI+rfAU8Altg+tujxTIWlvYG/b6yTtBtwM/PZM/f6UI9rn2H5K0k7A\nDcAZtm+suGhTIumPgPnA822fXHV5pkrSfcB824M07q0yqSF1r5upMmYU29dT9H6Z8Ww/aHtduf8k\ncBfPHTk+Y7jwVHm4U7nN6L8eJe0HvBH4fNVliXpKQOpeN1NlRA2Uswq/Grip2pJMTdm8dQuwBVht\ne0Z/HuB/A38C1GsRnqkxcI2kmyUtrbowM10CUgwUSc8Dvgq83/YTVZdnKmyP2D6CYhT8AkkztllV\n0snAFts3V12WaXaM7SMpZr1+b9kMHpOUgNS9bqbKiAqV71q+Clxm+2tVl2e62P4ZcB0wlaleqnY0\ncEr5zuUK4HWS/rbaIk2d7c3l1y3AVRRN+zFJCUjd62aqjKhI2QngQuAu25+sujxTJWlPSbuX+79C\n0Znm7mpLNXm2P2R7P9sHUvzufNv22you1pRImlN2oEHSHOAEYMb3WK1SAlKXyunTx6bKuAu40vb6\naks1NZIuB34AvELSJkmnVV2mKTga+D2Kv7xvKbdFVRdqCvYGrpN0G8UfQ6ttD0RX6QGyF3CDpFuB\nHwLftH11xWWa0dLtOyIiaiE1pIiIqIUEpIiIqIUEpIiIqIUEpIiIqIUEpIiIqIUEpIiakfR+SbtW\nXY6Ifku374iayQzSMaxSQ4qBIOntkm4r1w+6VNKBkr5dnrtW0gFluosknS/pRkkbJR1brgt1l6SL\nGvJ7StKnyrWIrpW0Z3n+iPLe2yRdJemF5fnvSPpouYbRvZJ+szw/S9LHJa0p7/nD8vyx5T1fkXS3\npMtU+C/APhSDYq8r779I0h3lujsf6PM/bUTfJCDFjCfplcCHgdfZPhw4A/gscLHtw4DLgM803PJC\n4LXAByimf/oU8ErgVZKOKNPMAdbafiXwXeDs8vwlwJ+W+d7ecB5gR9sLgPc3nD8NeNz2rwO/DvyB\npIPKa68u0x4C/CpwtO3PAA8Ax9k+DjgC2Nf2obZfBXxhCv9UEbWWgBSD4HXAl8eauGw/ShFwvlhe\nvxQ4piH937loq74deMj27bZHgfXAgWWaUeBL5f7fAsdIegGwu+3vlucvBhpndx6b0PXmhnxOAN5e\nLiNxE7AHMK+89kPbm8pn39JwT6ONwK9K+qykhcCMnsE8op0EpBhGPy+/jjbsjx3v2OKebl62juU1\n0pCPgPfZPqLcDrJ9zbj04+/55UPtx4DDge8A7yaL28UAS0CKQfBt4M2S9gCQ9CLg+xSzSgP8LvC9\nCea5A/Cmcv93gBtsPw48NvZ+iGIy1+82u7nBKuA95dIYSDq4nBm6nSeBsVmk5wI72P4qRbPkkRP8\nHBEzRqu/BiNmDNvrJf1P4LuSRoAfAe8DviDpg8BPgXdNMNt/plgU78MUK7a+tTz/DuCCslv2xi7y\n/TxFU9y6comMnwK/3eGeZcDVkh6geMf0BUljfzx+aIKfI2LGSLfviCYkPWX7eVWXI2KYpMkuIiJq\nITWkiIiohdSQIiKiFhKQIiKiFhKQIiKiFhKQIiKiFhKQIiKiFhKQIiKiFv4V1ol6DcHM7GkAAAAA\nSUVORK5CYII=\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f384eb82610>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(7, 3))\n",
    "plt.imshow(w_bins[:20, :], interpolation='nearest', aspect='auto'), plt.colorbar()\n",
    "_= plt.xlabel('components'), plt.ylabel('samples')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**backword reconstruction from encoding and dictionary**"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [],
   "source": [
    "img_rct = ptn_dict.reconstruct_samples(atlas, w_bins)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAXQAAAD7CAYAAAB68m/qAAAABHNCSVQICAgIfAhkiAAAAAlwSFlz\nAAALEgAACxIB0t1+/AAAADl0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uIDIuMS4yLCBo\ndHRwOi8vbWF0cGxvdGxpYi5vcmcvNQv5yAAABR5JREFUeJzt3dtu4jAUQNF61P//Zc/TSCMEJYAL\nyc5ar0gtkqWdE5PLmHN+AXB8fz79BQBYQ9ABIgQdIELQASIEHSBC0AEiBB0gQtABIgQdIOL7nf9s\njOG21J2Yc45Vf8u67od1bdq6riZ0gAhBB4gQdIAIQQeIEHSACEEHiHjrZYuvuPYijjGWXaHFh1hX\nWOcQE/qttyrNOW9+xv5ZV1jrEEEH4D5BB4gQdIAIQQeIEHSAiETQXRHRZF3hMYmgAxAJuhtRmqwr\nPCYRdAAEHSDjEEG/deo9xnBafmDWFdY6RNABuO8wT1s0sTVZV1jHhA4QIegAEYIOEHGYPXSAV10+\nTqL2G44JHSDChA7kneVBbyZ0gAhBB4gQdOCUaj+Ifn0JOnAC154PVNxXF3SACFe5AKdR3Gb5nwkd\nIELQASIEHSBC0AEiBB0gQtABIgQdIELQASIEHSBC0AEiBB0gQtABIgQdIELQASIEHSBC0AEiBB0g\nQtABIgQd+Kg5Z/KFzZ/gnaLAcpeBHmPcjfacM//Oz99mQgeIEHRgqWuT+NYtFdsvrxF0YImVMRb1\n5wg6QISgA0QIOkCEoANECDqwxBhj2XXkrkd/jqADS92L8U+frzwonJGgA0S49R9Y7pUpneeZ0AEi\nBB0gQtABIgQdIELQASIEHSBC0AEiBB0gQtABIgQdIELQASIEHSBC0AEiBB0gQtABIgQdIELQASIE\nHSBC0AEiBB0gQtABIgQdIELQASIEHSBC0AEiBB0gQtABIgQdIELQASIEHSBC0AEiBB0gQtABIsac\n89PfAYAFTOgAEYIOECHoABGCDhAh6AARgg4QIegAEYIOECHoABGCDhAh6AARgg4QIegAEYIOECHo\nABGCDhAh6AARgg4QIegAEd/v/GdjDC8w3Yk551j1t6zrfljXpq3rakIHiBB0gAhBB4gQdIAIQQeI\nEHSACEEHiBB0gAhBB4gQdIAIQQeIEHSACEEHiBB0+CVzelgh7/XWx+dCzb1o//t8jGVPtYWbTOgA\nEYIOT7Klwt4IOkCEoMOD5pwPT+emed5B0OEBr4T5mQMBPELQASJctggbmKw5AhM6bDDGWHYtuYMD\nv0XQASIEHSDCHjpwGpfbXbVHMpjQASJM6EDeWX6INqHDm9VO89kPQQeIsOUCG5zllP1MimdKJnTY\nYOWNRbzftfUrHqQFHSDClgs8YIzx9GRnwv+8+hqY0AEiBB0e9Mx+en0yZB9sucCTbm2/iDefYkIH\niDChwwtM4+yJCR0gQtABIgQdIELQASIEHSBC0AEiBB0gQtABIgQdIELQASIEHSBC0AEiBB0gIhn0\nOWfyBbAAPznM43MvA73l3Y5zTo83BU4jOaEDnNEhgn5tEt+6pWL7BTiLXQd9ZYxFHajbddAB2E7Q\nASIEHSBC0AEidh30Mcay68hdjw7U7Tro/9yL8U+frzwoAOzZIYIOwH2HufX/lSkd4AxM6AARgg4Q\nIegAEYIOECHoABGCDhAh6AARgg4QIegAEYIOECHoABGCDhAh6AARgg4QIegAEYIOECHoABGCDhAh\n6AARgg4QIegAEYIOECHoABGCDhAh6AARgg4QIegAEYIOECHoABGCDhAh6AARgg4QIegAEWPO+env\nAMACJnSACEEHiBB0gAhBB4gQdIAIQQeIEHSACEEHiBB0gAhBB4gQdIAIQQeIEHSACEEHiBB0gAhB\nB4gQdIAIQQeIEHSACEEHiBB0gAhBB4gQdICIv7KGyBoKdZesAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f38496a0a50>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.figure(figsize=(6, 4))\n",
    "for i in range(6):\n",
    "    plt.subplot(2, 3, i + 1), plt.imshow(img_rct[i], cmap=plt.cm.Greys_r), plt.axis('Off')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    ""
   ]
  }
 ],
 "metadata": {
  "signature": "sha256:18bf5511bac6f35a7a840ff6d9cb0cb31a12a102abf2ff3461ff081aad59342c"
 },
 "nbformat": 4,
 "nbformat_minor": 0
}