Tuto-GUDHI-persistence-diagrams.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "# TDA  with Python using the Gudhi Library  \n",
    "\n",
    "# Persistent homology and persistence diagrams"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "import gudhi as gd  \n",
    "from sklearn import manifold\n",
    "from pylab import *\n",
    "%matplotlib inline\n",
    "# import sklearn_tda"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Persistent homology and persistence diagram"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "One strong motivation for defining homology groups is the fact that two shapes can be distinguished by examining their holes.\n",
    "\n",
    "Persistent homology is a powerful tool to compute, study and encode efficiently multiscale topological features  of nested families of simplicial complexes and topological spaces. It encodes the evolution of the homology groups of the nested complexes (holes) across the scales.\n",
    "\n",
    "In the example below we consider the filtration given by a union of growing balls centered on the finite set of points $C$. \n",
    "<img src=\"Images/persistence.png\"  style=\"width: 600px;\"/>"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "a) For the radius r = 0, the union of balls is reduced to the initial finite set of point, each of them corresponding to a 0-dimensional feature, i.e. a connected component; an interval is created for the birth for each of these features at r = 0.  \n",
    "\n",
    "b) Some of the balls started to overlap resulting in the death of some connected components that get merged together; the persistence diagram keeps track of these deaths, putting an end point to the corresponding intervals as they disappear.    \n",
    "\n",
    "c) New components have merged giving rise to a single connected component and, so, all the intervals associated to a 0-dimensional feature have been ended, except the one corresponding to the remaining components; two new 1-dimensional features, have appeared resulting in two new intervals (in blue) starting at their birth scale.  \n",
    "\n",
    "d) One of the two 1-dimensional cycles has been filled, resulting in its death in the filtration and the end of the corresponding blue interval.     \n",
    "\n",
    "e) all the 1-dimensional features have died, it only remains the long (and never dying) red interval. \n",
    "\n",
    "The final barcode can also be equivalently represented as **a persistence diagram** where every interval (a,b) is represented by the the point of coordinate (a,b) in $\\mathbb R^2$. \n",
    "\n",
    "Intuitively the longer is an interval in the barcode or, equivalently the farther from the diagonal is the corresponding point in the diagram, the more persistent, and thus relevant, is the corresponding homological feature across the filtration. \n",
    " \n",
    "          "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Protein binding dataset\n",
    "\n",
    "The data we study in this notebook pertains to configurations of protein binding. This example is borrowed from the paper of Kovacev-Nikolic et al. [[1]](https://arxiv.org/pdf/1412.1394.pdf).\n",
    "\n",
    "The paper compares closed and open forms of the maltose-binding protein (MBP), a large biomolecule containing $370$ amino-acid residues. The analysis is not based on geometric distances in $\\mathbb R^3$ but on a metric of *dynamical distances* defined by\n",
    "\n",
    "$$ D_{ij} = 1 - |C_{ij}|, $$\n",
    "\n",
    "where $C$ is the correlation matrix between residues.\n",
    "\n",
    "Correlation matrices between residues can be found at this [link](https://www.researchgate.net/publication/301543862_corr). We are greatful to the authors for sharing data.\n",
    "\n",
    "We start from the Vietoris-Rips filtrations of the Protein binding distance matrices. See this [previous tutorial](https://github.com/GUDHI/TDA-tutorial/blob/master/Tuto-GUDHI-simplicial-complexes-from-distance-matrix.ipynb) for more details on how to build the Rips complex."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "path_file = \"./datasets/Corr_ProteinBinding/\"\n",
    "files_list = [\n",
    "    '1anf.corr_1.txt', \n",
    "    '1ez9.corr_1.txt', \n",
    "    '1fqa.corr_2.txt', \n",
    "    '1fqb.corr_3.txt', \n",
    "    '1fqc.corr_2.txt', \n",
    "    '1fqd.corr_3.txt', \n",
    "    '1jw4.corr_4.txt', \n",
    "    '1jw5.corr_5.txt', \n",
    "    '1lls.corr_6.txt', \n",
    "    '1mpd.corr_4.txt', \n",
    "    '1omp.corr_7.txt', \n",
    "    '3hpi.corr_5.txt', \n",
    "    '3mbp.corr_6.txt', \n",
    "    '4mbp.corr_7.txt'\n",
    "]\n",
    "\n",
    "corr_list = [pd.read_csv(path_file + u,\n",
    "                         header = None,\n",
    "                         delim_whitespace = True) for u in files_list]\n",
    "\n",
    "dist_list = [ 1- np.abs(c) for c in corr_list]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We first consider the first distance matrix:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "D0 = dist_list[0]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Persistence barcodes and persistence diagrams"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We build the Vietoris-Rips complex from the distance matrix `D0`. See this [notebook](https://github.com/GUDHI/TDA-tutorial/blob/master/Tuto-GUDHI-simplicial-complexes-from-distance-matrix.ipynb) for more details."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [],
   "source": [
    "skeleton_protein0 = gd.RipsComplex(\n",
    "    distance_matrix = D0.values, \n",
    "    max_edge_length = 0.8\n",
    ") \n",
    "\n",
    "Rips_simplex_tree_protein0 = skeleton_protein0.create_simplex_tree(max_dimension = 2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we can compute persistence on the simplex tree structure using the `persistence()` method of the simplex tree class:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "BarCodes_Rips0 = Rips_simplex_tree_protein0.persistence()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The object `BarCodes_Rips0` is the list of barcodes: each element in the list is a tuple (dim,(b,d)) where dim is a dimension, b is birth parameter  and d is death parameter.\n",
    "\n",
    "Let's print the 20 first elements in the list:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "(1, (0.07963602000000003, 0.35798637))\n",
      "(1, (0.12677510000000003, 0.39508646999999997))\n",
      "(1, (0.26003449999999995, 0.5273952))\n",
      "(1, (0.07943339999999999, 0.31429881000000004))\n",
      "(1, (0.08248586999999996, 0.30429980999999995))\n",
      "(1, (0.11378021999999999, 0.31171713999999995))\n",
      "(1, (0.07726765000000002, 0.26078758))\n",
      "(1, (0.09107215000000002, 0.25065161999999996))\n",
      "(1, (0.0709843, 0.22765623000000001))\n",
      "(1, (0.09347402000000005, 0.24999733000000002))\n",
      "(1, (0.07013614000000001, 0.22504734000000004))\n",
      "(1, (0.08752541000000003, 0.20355559))\n",
      "(1, (0.21541215000000002, 0.32814707))\n",
      "(1, (0.06835270000000004, 0.17527247))\n",
      "(1, (0.08857625000000002, 0.19539684000000002))\n",
      "(1, (0.08241111999999995, 0.18353136999999997))\n",
      "(1, (0.10362273, 0.20264340000000003))\n",
      "(1, (0.09289979999999998, 0.19181444000000003))\n",
      "(1, (0.09581541999999998, 0.19310879999999997))\n",
      "(1, (0.09541275000000005, 0.18175165000000004))\n"
     ]
    }
   ],
   "source": [
    "for i in range(20):\n",
    "    print(BarCodes_Rips0[i])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "These 20 topolological features have dimension 1, they corresponds to holes of dimension 1."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We have access to persistence_intervals per dimension using the `persistence_intervals_in_dimension()` method, for instance for dimension 0:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([[0.        , 0.01498817],\n",
       "       [0.        , 0.01614179],\n",
       "       [0.        , 0.01618629],\n",
       "       [0.        , 0.01964641],\n",
       "       [0.        , 0.02002598],\n",
       "       [0.        , 0.02019621],\n",
       "       [0.        , 0.02128167],\n",
       "       [0.        , 0.02139837],\n",
       "       [0.        , 0.02189414],\n",
       "       [0.        , 0.02211046],\n",
       "       [0.        , 0.02279135],\n",
       "       [0.        , 0.02305466],\n",
       "       [0.        , 0.02376384],\n",
       "       [0.        , 0.02398958],\n",
       "       [0.        , 0.02427977],\n",
       "       [0.        , 0.02431161],\n",
       "       [0.        , 0.02514183],\n",
       "       [0.        , 0.0252621 ],\n",
       "       [0.        , 0.025566  ],\n",
       "       [0.        , 0.02572414],\n",
       "       [0.        , 0.02668906],\n",
       "       [0.        , 0.02703624],\n",
       "       [0.        , 0.02723942],\n",
       "       [0.        , 0.02730715],\n",
       "       [0.        , 0.02737215],\n",
       "       [0.        , 0.02743006],\n",
       "       [0.        , 0.02758413],\n",
       "       [0.        , 0.02769299],\n",
       "       [0.        , 0.0277112 ],\n",
       "       [0.        , 0.02815389],\n",
       "       [0.        , 0.02889452],\n",
       "       [0.        , 0.0291826 ],\n",
       "       [0.        , 0.02937602],\n",
       "       [0.        , 0.02954253],\n",
       "       [0.        , 0.02960991],\n",
       "       [0.        , 0.02965513],\n",
       "       [0.        , 0.02989476],\n",
       "       [0.        , 0.0299387 ],\n",
       "       [0.        , 0.03005481],\n",
       "       [0.        , 0.03019355],\n",
       "       [0.        , 0.0302212 ],\n",
       "       [0.        , 0.03044412],\n",
       "       [0.        , 0.03061732],\n",
       "       [0.        , 0.03066776],\n",
       "       [0.        , 0.03093647],\n",
       "       [0.        , 0.03093904],\n",
       "       [0.        , 0.03096135],\n",
       "       [0.        , 0.03119385],\n",
       "       [0.        , 0.03127303],\n",
       "       [0.        , 0.03137071],\n",
       "       [0.        , 0.03139056],\n",
       "       [0.        , 0.03141663],\n",
       "       [0.        , 0.03147105],\n",
       "       [0.        , 0.03165437],\n",
       "       [0.        , 0.03197555],\n",
       "       [0.        , 0.0320141 ],\n",
       "       [0.        , 0.03203755],\n",
       "       [0.        , 0.0320416 ],\n",
       "       [0.        , 0.03204383],\n",
       "       [0.        , 0.03205409],\n",
       "       [0.        , 0.03209123],\n",
       "       [0.        , 0.03216506],\n",
       "       [0.        , 0.03251925],\n",
       "       [0.        , 0.03252582],\n",
       "       [0.        , 0.0326297 ],\n",
       "       [0.        , 0.03320018],\n",
       "       [0.        , 0.03320555],\n",
       "       [0.        , 0.03322223],\n",
       "       [0.        , 0.03329767],\n",
       "       [0.        , 0.03348698],\n",
       "       [0.        , 0.03385555],\n",
       "       [0.        , 0.03414255],\n",
       "       [0.        , 0.03438535],\n",
       "       [0.        , 0.03452845],\n",
       "       [0.        , 0.03464101],\n",
       "       [0.        , 0.0351205 ],\n",
       "       [0.        , 0.03512744],\n",
       "       [0.        , 0.03513682],\n",
       "       [0.        , 0.03527348],\n",
       "       [0.        , 0.03541263],\n",
       "       [0.        , 0.03559076],\n",
       "       [0.        , 0.03564263],\n",
       "       [0.        , 0.03568637],\n",
       "       [0.        , 0.03572321],\n",
       "       [0.        , 0.03593942],\n",
       "       [0.        , 0.03595195],\n",
       "       [0.        , 0.03602799],\n",
       "       [0.        , 0.03609681],\n",
       "       [0.        , 0.03610733],\n",
       "       [0.        , 0.03616306],\n",
       "       [0.        , 0.03644183],\n",
       "       [0.        , 0.03654824],\n",
       "       [0.        , 0.03661895],\n",
       "       [0.        , 0.03724071],\n",
       "       [0.        , 0.03733498],\n",
       "       [0.        , 0.03768468],\n",
       "       [0.        , 0.03774493],\n",
       "       [0.        , 0.03798143],\n",
       "       [0.        , 0.03815909],\n",
       "       [0.        , 0.03847974],\n",
       "       [0.        , 0.03851494],\n",
       "       [0.        , 0.03857502],\n",
       "       [0.        , 0.03892696],\n",
       "       [0.        , 0.03907595],\n",
       "       [0.        , 0.03913963],\n",
       "       [0.        , 0.0392337 ],\n",
       "       [0.        , 0.03926192],\n",
       "       [0.        , 0.03947112],\n",
       "       [0.        , 0.03964138],\n",
       "       [0.        , 0.0397164 ],\n",
       "       [0.        , 0.03975846],\n",
       "       [0.        , 0.03984852],\n",
       "       [0.        , 0.03997199],\n",
       "       [0.        , 0.04015835],\n",
       "       [0.        , 0.04017049],\n",
       "       [0.        , 0.04023692],\n",
       "       [0.        , 0.04046704],\n",
       "       [0.        , 0.04072276],\n",
       "       [0.        , 0.04074775],\n",
       "       [0.        , 0.04088823],\n",
       "       [0.        , 0.04095619],\n",
       "       [0.        , 0.04096195],\n",
       "       [0.        , 0.04110269],\n",
       "       [0.        , 0.04123369],\n",
       "       [0.        , 0.04127687],\n",
       "       [0.        , 0.04131072],\n",
       "       [0.        , 0.04140789],\n",
       "       [0.        , 0.04165062],\n",
       "       [0.        , 0.04176154],\n",
       "       [0.        , 0.04181859],\n",
       "       [0.        , 0.04191496],\n",
       "       [0.        , 0.04196476],\n",
       "       [0.        , 0.0420939 ],\n",
       "       [0.        , 0.04227029],\n",
       "       [0.        , 0.04236735],\n",
       "       [0.        , 0.04254986],\n",
       "       [0.        , 0.04260688],\n",
       "       [0.        , 0.04261734],\n",
       "       [0.        , 0.04266441],\n",
       "       [0.        , 0.04296581],\n",
       "       [0.        , 0.04303701],\n",
       "       [0.        , 0.04306739],\n",
       "       [0.        , 0.04317269],\n",
       "       [0.        , 0.04346509],\n",
       "       [0.        , 0.0437037 ],\n",
       "       [0.        , 0.04372465],\n",
       "       [0.        , 0.04383413],\n",
       "       [0.        , 0.04391404],\n",
       "       [0.        , 0.04398658],\n",
       "       [0.        , 0.04400293],\n",
       "       [0.        , 0.04407297],\n",
       "       [0.        , 0.04409098],\n",
       "       [0.        , 0.04435811],\n",
       "       [0.        , 0.04439728],\n",
       "       [0.        , 0.04443404],\n",
       "       [0.        , 0.0444774 ],\n",
       "       [0.        , 0.04451368],\n",
       "       [0.        , 0.04482899],\n",
       "       [0.        , 0.04495304],\n",
       "       [0.        , 0.04547506],\n",
       "       [0.        , 0.04547828],\n",
       "       [0.        , 0.0455478 ],\n",
       "       [0.        , 0.04615431],\n",
       "       [0.        , 0.04619047],\n",
       "       [0.        , 0.04625533],\n",
       "       [0.        , 0.04628039],\n",
       "       [0.        , 0.04678061],\n",
       "       [0.        , 0.04679402],\n",
       "       [0.        , 0.04680115],\n",
       "       [0.        , 0.04704526],\n",
       "       [0.        , 0.04719844],\n",
       "       [0.        , 0.04747816],\n",
       "       [0.        , 0.04767326],\n",
       "       [0.        , 0.04768366],\n",
       "       [0.        , 0.04792482],\n",
       "       [0.        , 0.04797942],\n",
       "       [0.        , 0.04802758],\n",
       "       [0.        , 0.04811274],\n",
       "       [0.        , 0.04824104],\n",
       "       [0.        , 0.04824237],\n",
       "       [0.        , 0.04833753],\n",
       "       [0.        , 0.04848755],\n",
       "       [0.        , 0.0485553 ],\n",
       "       [0.        , 0.04857819],\n",
       "       [0.        , 0.04860536],\n",
       "       [0.        , 0.04885771],\n",
       "       [0.        , 0.04897734],\n",
       "       [0.        , 0.04901851],\n",
       "       [0.        , 0.04944188],\n",
       "       [0.        , 0.04948698],\n",
       "       [0.        , 0.04950527],\n",
       "       [0.        , 0.049558  ],\n",
       "       [0.        , 0.04972988],\n",
       "       [0.        , 0.04981983],\n",
       "       [0.        , 0.04982386],\n",
       "       [0.        , 0.04984923],\n",
       "       [0.        , 0.04986276],\n",
       "       [0.        , 0.05027105],\n",
       "       [0.        , 0.0504428 ],\n",
       "       [0.        , 0.05055695],\n",
       "       [0.        , 0.05066297],\n",
       "       [0.        , 0.05070345],\n",
       "       [0.        , 0.0509538 ],\n",
       "       [0.        , 0.05103037],\n",
       "       [0.        , 0.0510853 ],\n",
       "       [0.        , 0.0511247 ],\n",
       "       [0.        , 0.05119141],\n",
       "       [0.        , 0.05125808],\n",
       "       [0.        , 0.0514161 ],\n",
       "       [0.        , 0.05146225],\n",
       "       [0.        , 0.05157315],\n",
       "       [0.        , 0.05181604],\n",
       "       [0.        , 0.0518984 ],\n",
       "       [0.        , 0.05228542],\n",
       "       [0.        , 0.05245802],\n",
       "       [0.        , 0.05249016],\n",
       "       [0.        , 0.05256426],\n",
       "       [0.        , 0.05261345],\n",
       "       [0.        , 0.052687  ],\n",
       "       [0.        , 0.05269135],\n",
       "       [0.        , 0.05286959],\n",
       "       [0.        , 0.05308782],\n",
       "       [0.        , 0.05312981],\n",
       "       [0.        , 0.05328802],\n",
       "       [0.        , 0.05388745],\n",
       "       [0.        , 0.05406629],\n",
       "       [0.        , 0.05424359],\n",
       "       [0.        , 0.05432115],\n",
       "       [0.        , 0.05434062],\n",
       "       [0.        , 0.05447437],\n",
       "       [0.        , 0.05475445],\n",
       "       [0.        , 0.0547687 ],\n",
       "       [0.        , 0.05480814],\n",
       "       [0.        , 0.05501963],\n",
       "       [0.        , 0.05515901],\n",
       "       [0.        , 0.0551603 ],\n",
       "       [0.        , 0.05521654],\n",
       "       [0.        , 0.05539072],\n",
       "       [0.        , 0.05547035],\n",
       "       [0.        , 0.05559308],\n",
       "       [0.        , 0.05608519],\n",
       "       [0.        , 0.05611108],\n",
       "       [0.        , 0.05611354],\n",
       "       [0.        , 0.05622496],\n",
       "       [0.        , 0.05627249],\n",
       "       [0.        , 0.05637189],\n",
       "       [0.        , 0.057297  ],\n",
       "       [0.        , 0.05729971],\n",
       "       [0.        , 0.05746997],\n",
       "       [0.        , 0.05753706],\n",
       "       [0.        , 0.05846568],\n",
       "       [0.        , 0.05849421],\n",
       "       [0.        , 0.05855491],\n",
       "       [0.        , 0.05858433],\n",
       "       [0.        , 0.05881015],\n",
       "       [0.        , 0.05891398],\n",
       "       [0.        , 0.05910734],\n",
       "       [0.        , 0.05927985],\n",
       "       [0.        , 0.05930276],\n",
       "       [0.        , 0.05937792],\n",
       "       [0.        , 0.0595791 ],\n",
       "       [0.        , 0.05959093],\n",
       "       [0.        , 0.05960105],\n",
       "       [0.        , 0.06001639],\n",
       "       [0.        , 0.06091708],\n",
       "       [0.        , 0.06102111],\n",
       "       [0.        , 0.0611472 ],\n",
       "       [0.        , 0.061409  ],\n",
       "       [0.        , 0.06150065],\n",
       "       [0.        , 0.06153714],\n",
       "       [0.        , 0.06154331],\n",
       "       [0.        , 0.06164122],\n",
       "       [0.        , 0.0619127 ],\n",
       "       [0.        , 0.06200894],\n",
       "       [0.        , 0.06212457],\n",
       "       [0.        , 0.0624782 ],\n",
       "       [0.        , 0.06275173],\n",
       "       [0.        , 0.06289645],\n",
       "       [0.        , 0.06309263],\n",
       "       [0.        , 0.06329755],\n",
       "       [0.        , 0.06339085],\n",
       "       [0.        , 0.06340604],\n",
       "       [0.        , 0.0634087 ],\n",
       "       [0.        , 0.06365607],\n",
       "       [0.        , 0.06370852],\n",
       "       [0.        , 0.06376314],\n",
       "       [0.        , 0.06452474],\n",
       "       [0.        , 0.06464462],\n",
       "       [0.        , 0.06539584],\n",
       "       [0.        , 0.0655788 ],\n",
       "       [0.        , 0.06560913],\n",
       "       [0.        , 0.06565494],\n",
       "       [0.        , 0.06628835],\n",
       "       [0.        , 0.0664005 ],\n",
       "       [0.        , 0.06651979],\n",
       "       [0.        , 0.0665555 ],\n",
       "       [0.        , 0.06667502],\n",
       "       [0.        , 0.0667084 ],\n",
       "       [0.        , 0.06679813],\n",
       "       [0.        , 0.0668217 ],\n",
       "       [0.        , 0.06688167],\n",
       "       [0.        , 0.06728802],\n",
       "       [0.        , 0.06736142],\n",
       "       [0.        , 0.06764264],\n",
       "       [0.        , 0.06825241],\n",
       "       [0.        , 0.06869392],\n",
       "       [0.        , 0.06877336],\n",
       "       [0.        , 0.06892966],\n",
       "       [0.        , 0.06897429],\n",
       "       [0.        , 0.06918575],\n",
       "       [0.        , 0.06954789],\n",
       "       [0.        , 0.06996146],\n",
       "       [0.        , 0.07021497],\n",
       "       [0.        , 0.07103847],\n",
       "       [0.        , 0.07105814],\n",
       "       [0.        , 0.07122811],\n",
       "       [0.        , 0.07162656],\n",
       "       [0.        , 0.07171408],\n",
       "       [0.        , 0.07273788],\n",
       "       [0.        , 0.07326412],\n",
       "       [0.        , 0.07355539],\n",
       "       [0.        , 0.07385251],\n",
       "       [0.        , 0.0738902 ],\n",
       "       [0.        , 0.07394715],\n",
       "       [0.        , 0.07403273],\n",
       "       [0.        , 0.0746566 ],\n",
       "       [0.        , 0.07496114],\n",
       "       [0.        , 0.07510989],\n",
       "       [0.        , 0.07526754],\n",
       "       [0.        , 0.0761947 ],\n",
       "       [0.        , 0.07620041],\n",
       "       [0.        , 0.07628743],\n",
       "       [0.        , 0.07870332],\n",
       "       [0.        , 0.07881973],\n",
       "       [0.        , 0.07893213],\n",
       "       [0.        , 0.07893759],\n",
       "       [0.        , 0.07933292],\n",
       "       [0.        , 0.08004957],\n",
       "       [0.        , 0.08036911],\n",
       "       [0.        , 0.08080199],\n",
       "       [0.        , 0.08109694],\n",
       "       [0.        , 0.081755  ],\n",
       "       [0.        , 0.08176781],\n",
       "       [0.        , 0.08248695],\n",
       "       [0.        , 0.08258634],\n",
       "       [0.        , 0.08299805],\n",
       "       [0.        , 0.08464198],\n",
       "       [0.        , 0.08729047],\n",
       "       [0.        , 0.08744106],\n",
       "       [0.        , 0.08878104],\n",
       "       [0.        , 0.08969459],\n",
       "       [0.        , 0.09321352],\n",
       "       [0.        , 0.09341871],\n",
       "       [0.        , 0.09479173],\n",
       "       [0.        , 0.09541855],\n",
       "       [0.        , 0.10237913],\n",
       "       [0.        , 0.10285434],\n",
       "       [0.        , 0.10715871],\n",
       "       [0.        , 0.10830089],\n",
       "       [0.        , 0.10954903],\n",
       "       [0.        , 0.11366507],\n",
       "       [0.        , 0.11958183],\n",
       "       [0.        , 0.14547953],\n",
       "       [0.        , 0.18964299],\n",
       "       [0.        , 0.39963213],\n",
       "       [0.        , 0.40540978],\n",
       "       [0.        , 0.4251596 ],\n",
       "       [0.        , 0.50723615],\n",
       "       [0.        , 0.66869636],\n",
       "       [0.        ,        inf]])"
      ]
     },
     "execution_count": 7,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "Rips_simplex_tree_protein0.persistence_intervals_in_dimension(0)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The last bars (0.0, inf) die at infinity."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally we can plot the points (birth, death) in the so-called **persistence diagram**:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAZIAAAEaCAYAAAA7YdFPAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAA+gklEQVR4nO3deXjU5bnw8e8za0JCGIIgIiAEa2WTEIK7XTTYxVNbLVhtj6d6qqBlURRBbWvfesSyqFWrIkFaFSxCwiKuFcSqdFFWFUQUwk5YkyH7rM/7x/wGJmGSmSSTWZL7c125kvktz9wZwtzz7EprjRBCCNFSpkQHIIQQIrVJIhFCCNEqkkiEEEK0iiQSIYQQrSKJRAghRKtIIhFCCNEqkkhEXCilcpRSM5VSWim1QSk11fiaq5QqUkrltdHzzlVKzW2LshNNKeUwXjvd4Hi7/Z1FclIyj0TEk1JqJzBXaz0r5FgOsBMYoLUuifHzFQBorVdHce1YrXVhLJ8/HpRSWmutQh5H/TsLEQuWRAcgOiRn6AOtdYlSaiMwDRgXyydq5pvpiFg+d6JIAhHxJk1bQgBKqalAdqLjSBZKKUeiYxCpQ2okIuGMpq084HbjsQN4AFgHjARWaa1XG002M4H1wCrgZ8BirXWxkQg2Ag5glHGdA5gHlGmtRxllh7sux3iePOO8M9jEFUUsJUCwP2IUsDO0eazB/U7AobUubqrsKF6vmQSaAstoULsz+poa/s55xu9IuOcxXv9xIXHsBMqUUtnG79fYa95ouY28PnlGvKuNn7OBEVrrmNZCRQJoreVLvuL2ReBNqggYDYwFphqP8xpc4wh5vCH42LhnJ4FEkGd8jQVGh1w/OlgeUEDgDY4I140GihqJt7FYRhvnc4zHjsB/qXr3bwjzu+VFKruJ129Vg/Jywjznyd85pNzQ37s8zO+YE3LvhgbnT3vNoyw37OsDFDT4fUY3/D3lK7W+pEYiEmGdNj6VN6SUGk2gRuAMObweuAEoJPApvMQ4v9G4JxuYa3xf3VjZGJ+Oo7gu2licEOjjMb47lVKh9xcQeBPdGHL/GK31xijKDhdPXsPydKB/qalfI/ic9QYxKKUcIc+dTeB1xfjuaHD/aa95lOU6gzEa34OvT+g9TqRJMeVJIhHJJgdOjTwyFHH6m89JOtDUNI5A88xco+P+qgZv0lFf18xYThtlFvJmmtPwfEgSiKbshvJp8LtHqSzYZMephJEdUlYhgVrHLAJNV9PClBHueSOVC+F/n7JGfhYpShKJSDYlBJp3oh55pJQqMK4Pts/P5NQbY0uuC7b7NzuWML9LThPnmlt2U+U1ZQNGTQhAKTXP+B5MeDuBYL/PH5tIrM0tV3QQMmpLJBWjuckZ8mYenMzY1ITFvAaf7Be34LrQN+kcrXVJlLE02ixjJImS0OcM3t+S3zOkvLyQe5qcyGmczw55s3eEnA7GNYJA5/zqaBNAlOVC5GYradZqB6RGIuIiZGRQDvAzpRQ6ZFJiA1cBDyil1hmPnfrUSKlxQL7RpFJovPE5gRyj3wECiWCW8WY3zbh+bGPXQaDJSSlVYlwXTSzBsnOCsRAYgQUwUyk10+gbuMp4nIMxyiqkFhK27AgvZfCeHEKakJRSRQRGveWE/s5a60Kl1JKQ0WoY143jVCLdAJQrpYLllRCYNFrY2GtuvF6NlhvN60Mg4eQZ15S0ouYnEkxmtgvRgRlv+AW6/koDDuA94PYGAwWECEuatoTo2PI5fUCAk8AIMiGiIjUSITq4kOa8MgJ9Fg5gozQ1iWhJIhFCCNEq0rQlhBCiVdr1qK0zzjhD9+vXL9FhCCFEStmwYcMxrXX3aK9v14mkX79+rF8vfYZCCNEcSqk9zblemraEEEK0iiQSIYQQrSKJRAghRKtIIhFCCNEq7bqzPZ5cmzZTW1SEd+8+LH37kD5mDPbhuYkOSwgh2lxS1EhCFtFLSa5Nm6mc/Ri+snJMPXviKyuncvZjuDZtTnRoQgjR5pIikTRGKeVQShWELlEdPJbAsE5TW1SEysrClJWFMpkC37OyqC0qSnRoQgjR5pK2actIHjcQWDxunlJqLoG1gHK01sUhmxQ1vG8sgc2KyMzMZNy4cW0ea92a91Gd0oHQLU81uqaWtOPH2vz5hRAikRK21paxn0LoxjobgueC+yCEJgql1AZgcYPlrsMmk6D8/HwdjwmJzgd/E2jWyso6ecxfUYE5uyuOR6e3+fMLIUQsKaU2aK3zo70+YU1bxg50hVrrQmBV8GfjcTiLCWyI5ICTNZak2O85fcwYdEUF/ooKtN8f+F5RQfqYMYkOTQgh2lzSNm0Zu9AFO+GzjR3vigk0c60jsJtcY0knruzDc+G+KfVHbY0bK6O2hBAdQtImEji5f3fo4xIgKT/m24fnSuIQQnRISTFqq2HCEEIIkTqSIpEIIYRIXXFLJEqp0cackKlhzuUppXYqpTYYXzMbu6epcoQQQsRfXBJJsNPcGKrrDDOhMFtrPUBrPQK4HZgb7p4oyhFCCNFCpaWl3Hzzzc2+L141kpFAifFzCZAXerLBXJAco1M93D1NlgOBCYlKqfVKqfVHjx6N3W8ghBDtlMfj4YknnmDo0KF07dq12ffHK5E4GjzuFu4ipdTYkI73cPdELMeYi5Kvtc7v3j3qnSKFEKJD+sc//sHw4cN58803WbBgQYtWA4nX8F8nkB3FdaOA4NyQcPdEW44QQogmHDhwgClTpvDPf/6TBx54gCuuuAKAlqx2Eq8ayTpO1SZygFUNLwhdmLGJeyKWI4QQonFut5vZs2czbNgwzjzzTJYtW3YyibRUXBKJ0VyVE+wcD/aJKKVCE0E2IUuehLunsXKEEEJE9t577zFs2DBWrVrFggULuP3220lLS2t1uQlbtDEe4rVoI8jGVkKI5LVv3z7uvfdePvnkEx544AEuu+yyJq8fOnRoaiza2J7IxlZCiGTkcrmYMWMGw4cPp0+fPixdujRiEmmJpF5rK1WEbmwFoLKy8BvHpVYihEiEd999l4kTJ9K/f39eeeUVevfu3aKO9GhIIokB7959mHr2rHdMZWbi3bsvQREJITqqPXv2cM8997Bp0yYeeOABLr30UrTWbZZEQJq2YsLStw+6qqreMV1VhaVvnwRFJIToaOrq6njkkUfIy8sjJyeH4uJiLrnkkjZNIEGSSGJANrYSQiTS22+/zdChQ/nPf/7Dq6++yi233ILNZovb80vTVgzIxlZCiETYtWsXkydPZsuWLTz44INcdNFFbd6MFY4kkhiRja2EEPFSW1vLrFmz+POf/8xtt93G73//e6xWa9wTSJAkEiGESCGvv/46d999N4MHD+bVV1/lrLPOSlgCCZJEIoQQKWDnzp3cddddfPXVV/z2t7/lwgsvTEgzVjjS2S6EEEmspqaGhx56iIsuuogLLriAJUuWMHLkyKRIIEGSSIQQIglprVmxYgWDBw/m888/Z8mSJfziF7/AYkm+hqTki0gIITq4r7/+mkmTJrF7927+8Ic/MGLEiKRpxgpHaiRCCJEkqqur+c1vfsMll1zCiBEjWLRoEXl5eUmbQIIkkQghRIJprVm6dCmDBg1i+/btFBUVcdNNNyVlM1Y4qRGlEEK0U9u3b2fixIkcOHCA6dOnM3z48KRuxgpHaiRCCJEAVVVV3H///Vx++eVccsklvPLKK+Tm5qZUAgmSRCKEEHGktWbJkiUMHDiQkpISiouLueGGG1KmGSuc1I1cCCFSzLZt25gwYQJHjhxh5syZDBs2LOWascKJW41EKTVaKVWglJrayPk845rRTd0TqRwhhEg2lZWV3HfffXzrW9/i29/+NgsXLuSCCy5I+QQSFJdEEkwOWuvVgFMpVRDmsnFa62IgRymVE+6eKMsRQoikoLVm0aJFDBw4kH379rF06VKuv/56zGZzokOLqXg1bY0EFhs/lwB5wOrgSaXUWGCDUipHaz3LODYuzD3dmionpKyxAH379m2L30UIISLasmULEyZMoLy8nMcee4yhQ4e2mxpIQ/Fq2nI0eNytweMBxleZUmquUsrRyD2RykFrXai1ztda53fv3r2l8QohRIucOHGCyZMn893vfpeCggIWLFjAkCFD2m0SgfjVSJxAdoRrdmqtnUqpDQRqFOHuiaYcIYRoU1v3O3lz80FKy2s5q2s61+T2YtDZXVi4cCHTpk3jyiuvZNmyZXTt2jXRocZFvBLJOk7VJnKAVWHOBxOEg0DCKAlzjyNCOUII0aa27ndSuGYHmWkWemTZOVHjZvqLb7F95Z/RXhdPPvkkgwcPbtc1kIbi0rQV0oleYDxeDaCUWhVy3hFyvjDcPY2VI4QQ8fLm5oNkplnonGalrrqS91+ezcrHJnL5laN48cUXGTRoUIdKIhDHeSTBTvQGx0aFOb86wj2nHRNCiFgK13Q1uLcDgNLyWs7ItLLhvRX8fcGTXHj5d3hj5Wt0dXTBZFKJDTxBZEKiEEKECNd0VbhmB2OvPJfBvR34y3Yz55GHsSrN00//mdwLhmDp4GuESCIRQnRoDWsfZVUuMtMseH2aT/c5qXH5sJoV8975FM+mJby6pIhf/O+d3HLT9dgs7Ws+SEtJIhFCdFjhah+b9pTT/4xOHHTWYTGbSLPA1/96gzfeeoHrrr2G119bQVZWVqJDTyqSSIQQHVZoxzlA5zQrmXYLXx+uoksnG1UHvmJD0Z+wWxTPzXmOEdKMFZYkEiFEh1VaXkuPLHu9YwN6ZPLBpyWUvPEyh7f+k/+9YyK3/Ow60mzSjNUYSSRCiA7rrK7pnKhxn6yR+H0+tvxjOZ8uepYLv1XAX99YSdcuWZg76GisaEkiEUK0a00N5b0mtxeFa3YAUL5nGyuefwS7zcpLf5nHeeedh9UsCSQakkiEEO1WaGe61az4eMcx3tt6iOHndOXmy/szuLeDMblduee+aWxb9wETJk3m5z+9FiX5o1kkkQghUlpTNY5gZ/q+49XsPFKF1qCAzXvKqa5z4zi4lrlPzuS6667j6T+8TmZmZkJ/l1QliUQIkbIiTR4sLa+l1u1hx+EqAExK4dea/ds/Zf1zcznD0ZkXXniBAQMGdLhlTWJJEokQhqY+2YrkFG74bvD44N4Ozuqazt8/LQcCScRd7WTPqvlUlGzgjomTuXn0j7BZTJJEWkkSiRBE/mQr2l5LEnnD4bt7jlXx9aFK6jx+Nu0u49vnd8fl9YHfx8F1b7D/w79x8ZU/4JE/r6RLVmcZjRUjkkiEIPInW9G2GibyfceqeWDxZnp0TuO8XlmNJpXQ4bt7jlWxdX8FKE2a1US1y8vy9QdwHdjGttf/TJcuXfjLX+Yz8JvfwC6zCmNKXk0hCHyyzbTX/1yVabdQWl6boIg6ltBE7qzxsOd4NQBVLu/J2uHW/c7T7rsmtxdVdV4q6zx8fagSjcbn07i8fsqPH2VL0Qy+KHqUO8fezmuL/kru4G9KEmkDUiMRgtMnpkHgTeysrukJjKrjCG2i2nO8GovZhNWkOFHrYeeRKipqPDyyYgu//ckQBvd2sHW/kwVrd7HtQAVur490m4Ualw8UWLSfQ5+8zu5//I3LR/2Ih595nS6ZGTKktw1JIokR16bN1BYV4d27D0vfPqSPGYN9eG6iwxJRCp2Ylmm3UOXyUlXn5eeX9ktsYB1EaCKvcflIt5mocXmpc3tJt5nJTDNTUeuhcM0OCoacyYoNBygtr8FuNZNuM1Pn8aGByj2fs/vtZ+nW7QxefuklzhvQH5vUQNqcvMIx4Nq0mcrZj+ErK8fUsye+snIqZz+Ga9PmRIcmojS4t4OxV55Ll042jlS46NLJJh3tcRTaRJVuM1Hr9lFZ56WT3YLNbMLng6x0K5lpFv72rz04q92k2SzYLebAUu415exYPpOSFbMYf+edLF34AkO+OUCSSJxIjSQGaouKUFlZmIylpVVWFn7juNRKUsfg3g5JHAkSTORvbj5I5zQrtW4fVpMiw27G4/Xj8fkZcGYmmXYL5dVurGYTnexm/D4vX31QzJerFnLF93/M/yv8I5mdOslorDiLWyJRSo0GnEBeuO1ylVLlQAmwWms9rbF7IpWTCN69+zD17FnvmMrMxLt3X4IiEiL1hCbyrfudPLJiCxW1HrLSrQw4M5PsTDuVdR66Zthwe/0c/HIDny17ih5n9mThgpfJyemHzSw1kESISyIx3vzRWq9WSuUopQq01qsbXDYm9Fi4ewBHFOXEnaVvH3xl5aiQzW50VRWWvn0SGJUQqWtwbwe//cmQk0OCM+0WKus8VNV5uXqAjf976DccKfmcyffex09+MAq7VRJIIsXr1R9JoLaB8T0vzDUOpVROhHsilqOUGquUWq+UWn/06NFYxB5R+pgx6IoK/BUVaL8/8L2igvQxY+Ly/EK0Rw37rTKsYNn+Jvfd/ENGjTyfVW+tZMx/XS1JJAnEq2nL0eBxtzDXZANlSqm5WutxjdwTsRytdSFQCJCfnx+XdQ/sw3Phvin1R22NGyv9I0K0UrC5a/Xq1UycOJE+ffqwYMEC+vbtK8uaJJF4JRIngUTRKCMBoJRyhvSDNLwnYjmJYh+eK4lDiBjbt28f99xzD+vXr+f+++/nsssuA5AkkmTiVSdcx6naRA6wKvSk0Rw12nh4vIl7mixHCNE+uFwu/vjHPzJ8+HD69OlDcXHxySQikk9cEonWuhgIdpgT7CBXSgUTwRLAGXK+ONw9jZUjhGg/3nnnHYYOHcqHH37IK6+8wq9+9SvsdnvkG0XCqPZcRczPz9fr169PdBhCiCjs3r2byZMn89lnn3H//fdz6aWXShNWggwdOnSD1jo/2utlQqIQIiZaup9LXV0ds2fP5qmnnuLWW2/ld7/7HTabTZJICpFxc0KIVgsuA3+ixl1vP5dwK/aGevPNNxkyZAgff/wxixYt4pZbbsFms8UnaBEzUiMRQrRaU/u5BL+H1lTS3WXcfffdbNu2jQcffJCLLroIrbXUQlKU1EiEEK3W2H4uX5VW1KupHHNW8L8TpzIifySDBg1iyZIlXHjhhZJAUpzUSIQQrdbYfi7VLi9ndkkj025h2yf/4I35M8g5byB/W/QqfXv3kgTSTkgiEUK0WmP7uWTYLLjKDrL0hRk4D+/j9w89xLcvvxSLSSYVtifStCWEaLVw+7ncfEkvdq15kTlTf8FlF+XzxoplXHlFIImI9kVqJEKImAiui6W1ZsWKFYweNZmBQy5gafESzurZU/YIacckkQghYuarr75i0qRJ7N27lz/84Q/k5+dLE1YHIIlEiHaupRMFm6O6uprp06dTWFjIHXfcwYwZM7BYLJJEOghprRSiHWvpRMFoaa0pLi5m4MCBfPXVVxQVFXHjjTdischn1I5E/rVjxLVpc/39SMaMkWXlRcI1NVEwtFbSVK1l634nC9buYtuBCgAGnp3FzZf3x1x1iIkTJ3Lw4EGmT59OXl6e1EA6KEkkMeDatJnK2Y+hsrIw9eyJr6ycytmPwX1TJJmIhNm638k/tx/F5fXh82ssJhNdMqz0ze5EaXltvesK1+zA79ccqaxje2kFa7cf4YfDzmL3sRrWl5Th8fnJtJuxW81s/LqUN//6OAfWv8PECeO57rrrpBmrg4tJIlFK9dNa745FWamotqgokESMPdtVVhZ+47gkEpEIJ5OD1tS4vZiUwuPzY66DzXudZNotTF6wgbO6plNW5cLv1+w4UoXL48OvNT6/5qWPdmM2gdZgMimqXF6OfP4Ptr7xPINyR7KkqJiePc5I9K8qkkCLE4lSKivk4TTgztaHk5q8e/dh6tmz3jGVmYl3774ERSTaq2g7zoNNWmYTKBQmpdBaU1XnxevzYzMrrGbFxzuOcfhELSYFWiusFhPKSDoAfqOSUXV4F7venoPFU81js2Zx0cg87DIhRBia9ZeglLpdKeVXSvkIbHtbbnwfG/vQUoelbx90VVW9Y7qqCkvfPgmKSLRHzek4D6595dcKR4YVs0nh0xqXx4/WUFHn5bN9TjRgs5hx+8CnNRp9MokowF1XzZ53C/nipWn8+AejeHP5Yr51Sb4kEVFPc2skA4CuWusToQeVUjNiF1LqSR8zhhMPPYTn8zL8Ljcmuw1TdjaZ4zp0fhUxFm3HOYDdamL9rjIqaz0oBek2My6vD6XAbjHh9mlq3T7SrGY6p1mocfsA8Hj9+HVgNNaxLe+zb/V8huRdzIsrl9P9jDOwyKRCEUZzE8m6hknE8MdYBJPK/FXV+I8fR7vc+O02sMnWoCK2Sstr6ZFV/+8q026p13EOgZrLkRN1VNZ58OPH5dYnEwVAndePUoEmr+NVLpRSRtMW+DTUHd1NyVvPYvHV8ac/Pc5FI4ZLDUQ0qbl/HeVKqX5hjt8e6Ual1GilVIFSamqE62Y2dU+05cRT1XNz8JeXY+qajblPH0xds/GXl1P13JxEhybakbO6plPl8tY7VuXyclbX9JOPt+538siKLewvr6HO7QWt6tUi7JZA0vD7AwnF49P4fH7SbSb8rmr2rZrLtgX3M/rHP+TtFYv51kUjJImIiJqskSil/t7wEDBCKVUGlBjHugH9gceaKGc0gNZ6tVIqRylVoLVeHea6AiCnsXsARzTlxJv7082otDRUcGc3mw2lNe5PNyc0rliJx8xoEVljK+z+/NJ+QGgfigefX2MyBRKAxRzoaA/2hyg0FXW+kJI1+ze8y+5VfyHvost5eeUKunXLlrWxRNQiNW11I9Bs5WziGgVEqh2MBBYbP5cAeUC9BKCUyuFUcmrsnm5RlDMWo/O/b9++EcKKkcaGz7eDYfXBN6fMNEu9Dt6xV54rySTOgivshib1n1/a7+S/Q7APJauTldJyL1aLCb9fU+cJJA2/hmqXl9DpHrWHStj19jOkmXw8++enGJE7DKtZEohonkiJZJrW+r1IhSiljke4xNHgcbcw1+QYNY2m7olYjta6ECgEyM/Pj8tbuTV3GO5169FKgdUKHg+6rg7byPx4PH2bak4Hr2h7wRV2wwn2oZzTLYNDzlq8Pj+g8evApz2LGXy+wOcbb10VB95/mfJtH/HL2+/kV78YTSe7NWy5QkTSZONnpCQS7C/RWm+K8DxOILuJcsI1UYW7p8lyEqXzr3+NuXdvUApdUwNKYe7dm86//nWiQ2u1xrZQbdjBKxIvtA8lw27B7fVT59GYFJhMYFYmTPg5tvldPnvuds7ponjr9de485YbJYmIVmnWqC2l1G1a6xfqH1JXAf0bHG9oHadqEznAqgbny0L6QHKUUnmN3OOIUE5C2Ifn4vi/h9vlWluNbaEa2sErksM1ub144u0vKS2vwW41k51h43iVGwU4OtmoKf2az4ufpJNFMeeZZ8gdNhS7RZqxROs1d/ivI/SB1noXsEspdVtTN2mti5VSU41kQbD2oZRapbUepbXeaDwey6kO9cbuOe1YMrAPz20XiaOhSB28InkM7u2ge2c7zmo3Hp+fzulWLGZFXXUFX618nsNbPuL2Oydw4+ifkG6TZfZE7KhIC60ppW4HRgBdCXRub2xwSQ5QorX+WZtE2Ar5+fl6/fr1iQ4j5cmordQxecEGemTZUUrh9/v58K0i/vHqM4y87EqmPzgZR5csTDIaS0QwdOjQDVrrqDt5I34s0VrPA+YZ8zYUsKTBJSVR9JGIFNZUB69ILsGmyBP7tvPa3OnYrSZemPs8gwYNxCajsUQbibp+q7WepZS6KppRXEKIxLjsnDTumPRb9nz6Eb+eeBdjrvsRdos50WGJdq5ZDaUNk4hSqgtQAOzUWm+OYVyiFaQpquPx+XzMmzeP3//+94z6/g957v9eo0tWZ2nGEnHRoh43Ywn54DDcDcA4YHOMYhKtIBMIO56PP/6Y8ePHY7PZeP755znvvPNkkykRV80d/nsVUMSpOdvK+DniWlsiPmQCYXKIR63w6NGjPPDAA7z11lvce++9fO973wOQJCLirrmrsRVorbO11t2AsVrrbAKjtpwxj0y0iEwgTLzm7BvSEj6fj2effZbBgwdjNptZvnz5ySQiRCI0t2krdCxtVwCt9YmQZU1EgskEwsRry1rhv/71L8aPH09GRgYvvPACAwYMkBqISLhmrw+tlLre+LGrUmqY8XNe7EISrXFNbi+q6rxU1nnQWlNZ56Gqzss1ub0SHVqH0Ra1wsOHD3PLLbcwZswYfvnLX/L888+Tk5MjSUQkheYmkhLgQWONrUJgqbFg44BYByZaJrhCbJdONo5UuOjSySYd7XEWzb4h0fJ6vTz99NMMGTKE9PR0li9fTkFBQaxCFSImmjv8dxMQOtvxXKXUcJmQmFxkAmFixWpZmY8++ogJEybgcDiYP3++NGOJpNXs4b/G0N8bgNVa6920i103RKpKxjkzkfYNiaS0tJSpU6fy/vvvM3XqVK688kpARmOJ5NWS4b8zCazMWwLsBk4opa7UWq+JfXhCNC6Z58y0pFbo8Xh45plnePTRR7nhhhtYvnw56ekySEIkv+bWSEYFF/Iykgpa613Gsu9CxFV7mjPzwQcfMH78eLp3786LL75Iv379pAYiUkZzE8knjRyXv3gRd8EdAUOl2pyZgwcPMmXKFNauXcvUqVP5zne+A0gzlkgtzU0kFxp7iFRiJA9jBNeFwLIYxyZaqDX9BsnY59CYVJ4z4/F4eOqpp5gxYwY33XQTS5culWYskbKaO/z3j8AmpdQ6YKZS6msCS6Y8GvPIRIu0ZlZ1W8/IjrVUnTOzZs0ahg0bxjvvvMOCBQu44447JImIlNbc4b8nCAz5HQ30J7AXydI2iUy0SLDfwOvTfLrPSY3Lh9WsWLB2FzNuHB7VvanS59Da0VHxtn//fu69914+/vhjpk2bxhVXXAFIM5ZIfS1a/VdrXRzrQERslJbXYjUrviytwGI2kW4z4fb62bS7nK37nU2+yaZin0MqzJlxu9386U9/Yvbs2dx8880UFxeTlpaW6LCEiJmomraUUv2UUnOUUl8rpXzG98UhS6REU8ZopVSBsdNiuPMFxtfMpu6JVE5Hd1bXdHYeqcJiNmEzm1AoTCgy0yy8uflgxHtjNSNbBLz77rsMHTqUNWvWsHDhQm677TZJIqLdiVgjMfZsnwusJrAsStAA4H2l1HSt9eMRyhgNoLVerZTKUUoVaK1Xh5zPIzC0eJpSappSKgdj/a7QewBHU+WIQL/Be1sPkWEzo7XG69N4fH6+eVbniDWLWM3Ibi9aM/Bg7969TJ48mY0bN3L//fdz+eWXo7WWZizRLjVZIzHmiowDumqtr9Zazw75usNYRv4ipdSVEZ5nJIEJjBjf68070VpvNJKIg0C/S0kj9zRZjhHzWKXUeqXU+qNHj0YIq/0Z3NvB8HO6opSi1uPDZjVxfq8sbFZzxJqFrNN1SksHHrhcLqZPn05eXh79+vWjuLiYyy67TBKIaNci1UhGBycgNkZrfYNSag7Q1Mx2R4PH3Rq5Lp9Te5uEuydiOVrrQoyaU35+fof833vz5f1Pzvhubs0iFfoc4qElAw/efvttJk2axDe+8Q3+9re/0atXco8eEyJWIiWSaDcaORHhvJNTW/M2ymiyGmM0hYW7J6pyOrpUG82UjJoz8GD37t3cfffdbNmyhWnTpnHppZdKDUR0KJESyc4oyzke4fw6TtUmcoBVoSeNDvadRm3CSSBZhLvH0VQ54hSpWbRONJMd6+rqmDVrFk8//TS33norDz30EDabTZKI6HAijdqK9n9Ek9cZw4WDHeYEO8iVUsFEMBcoCXaoa60Lw93TWDlCxFqkyY5vvPEGgwcP5pNPPmHRokX88pe/xGazJThqIRJDNfXpyZjBvpim92R3AD/TWo+MaWQxkJ+fr9evXx/5QiHCCDdqK91dxl133cX27duZNm0aF198sdRARLszdOjQDZH6x0NFatoagbE3ewT9o31CIVJFaPNgbW0tM2bM4Nlnn+VXv/oVDz/8MFarVZKIEEROJLO01vdHKkQpNSNG8QiRVLTWrFy5krvvvpsLLriAV199lZ49eyY6LCGSSqREMjfKcqK9TrSxVFq9N9nt2LGDSZMmsXPnTh566CEuvPBCqYEIEUaTne1a613RFBLtdaJtpdrqvcmqpqaG3/72t1x88cXk5uayePFiRo4cKUlEiEa0aNFGkZyamkQHsGDtLrYdqABg4NlZ3Hx5/5SvrcSyBqa1Zvny5UyePJm8vDyWLFlCjx49YhuwEO2QJJJ2pLFJdF+VVrC9tILS8hrsVjMK+Gyvk8ff3Ma91wxM2WQSyz3bt2/fzqRJk9i7dy8PP/ww+fn5UgMRIkrN3dhKJLHGVu+tdnlxVrtJs1mwW8zYLGbSrWZO1HoirgiczEJrYEopOqdZo1rlOFRVVRX3338/l112GSNHjmTRokWMGDFCkogQzSCJpB1pbBJdhs2C2+vHaj614o3FrPD4/Em910gkpeW1ZNrrV6qj3T9Fa01RURGDBg1i586dFBcXc+ONN2KxSCVdiOaS/zXtSGNrbL25+SDlNW48Po3NSCZen8ZqNqX0XiORljFprP9k27ZtTJw4kUOHDvHoo48yfPhwqYEI0QqSSNqZxtbYCvaRaKOPpM7jp2eXtKTf37wpTe2fEq7/5Jk3P6VuUzFvLH2VX//611x//fWYzWZJIkK0kiSSdiDSyKXBvR3c84Pz643auqCvI+VHbTW1yvGsN7442X+itWbnJ6t466+PkX/hxRQXF9OtW2M7GQghmksSSYrbut/JE29/ibPajdvrZ+/xaraXVnDPD84/LZnMuHF44gKNkXBJc+p/DTrtuuAItsN7d7By7nQ8tZX86fHHyBueW6+vSAjRepJIUtyCtbsoLa8hzWahk92Mx6cpLa9hwdpdCU0cbTHDvjnDfbvafayYN4Mv1r7N/469g5tvHEOaTf7chWgLMmorxW07UIHdasZmNqFQ2Mwm7FbzySasRGirGfbRDPfVWrNw4UKenfQT0rxVvLZ8Kb/675skiQjRhuR/VzvQsKEm0Q03LdmmNhqRdi387LPPmDBhAhUVFTz15BMMGTKU6LfUEUK0lNRIUtzAs7Oo8/jxeP1orfF4/dR5/Aw8OythMbVmfkdTGptw6bB6uOuuuygoKODqq6/m5ZdfZsiQIUgSESI+pEaS4m6+vD+Pv7mNE7Ueaj0+rGYTPbukcfPlgS1iErEacDTb1DYULk6g3rGBvTqzesthIJCYKmvdrF+zkk0rnud7VxewbNkyHI62/d2EEKdrcofEVNdRdkhsLFmEdk6HzrMIdk63VZKJ9LzRXH/oRB3arzmra3q9MgqGnMm2g5V8unkz//zbbDpb4Xe/eYDBgwfLfBAhYiTWOyTGjFJqNIEte/O01rManHMABcbDkVrraY3d01Q5HVVjkxAjrQYcqwUPw8XT2PyOoNAkduhELd0ybfXi/PJgYLDAeWdl1Yt9w1f7Obb2FVYXFTFp0iSuvfZalFKSRIRIoLgkEuPNH631aqVUjlKqQGu9OuSSG4zzhUqpkUqpsUBZw3sI7A/fVDkiRLBzurzazZ7j1Zyo8VDn9uLxad7efJDMNDODznacHAEFre8QDwpNbsGk8cL7O+s1UQWT2PbSCqrrvHSyWcjODHSmu71+VMioAb/fz5cfvc47C55i9I+vYfny5WRlJa4fSAhxSrxqJCOBxcbPJUAecDIBaK0LQ67NIbDj4rgw93RrqhwAIwmNBejbt28sf4cmuTZtpraoCO/efVj69iF9zBjsw3Pj9vzhnNU1nX3HqtlzvBqtodrlxWO8QXt8Pk7U+Pl0bznD+nala4YtYod4NP0YDZvHws39mPf+Tvp263QyeWV1slLj8rKvrOZkIrFZTo0DObBjK6/NnY5FaQrnPEvuBUOkBiJEEolXInE0eBx2fQqlVA5QprUuMZq7Gt4TsRwjKRVCoI+kBbE2m2vTZipnP4bKysLUsye+snIqZz8G901JaDK5JrcXDyzeDECN24vPF0giVpPCD5hMCrfXz57j1XTNsDXZIf76xn3Me38nXp+mc7oFr8/PE29/ebIfo7HmsXDNa16f5lilm97ZGQCc0y2DT/eWc6Kmlqo6LzaLCbvVhHJVs+SZ/8eOde9z58RJ3HDdj7FbZW0sIZJNvBKJE8iO4rrRWutxTdwTbTlxVVtUFEgiRlOLysrCbxxPVCIJ1h6qXV6UUrg8flBgUQqv1vj94PNrzCaFs9rNupLjVLm8DD+nK1v3O0+rVcx7fycAndMteHzaqOVobBbzyX4Mr0+zv6yGBxd/ymXf7M41ub3Czv3onG6hss5T75jPr7GYFUqB9vvZ8++3+OLt+Vz2nQKefP01HF26YEr0BBkhRFjxSiTrOFWbyAFWNbxAKTU6pEO9oJF7HJHKSQTv3n2Yevasd0xlZuLduy8h8YQ2J3XvnEaN24vX58fvB69fowjURvx+jcenOVHjwWo2cUHvLphN6mStAgI1in9uP8qJGg9dM6zG7PnAO3pZpYusToGfy6vdbDt4ArNJodEnayfpNjNVLm+9ocA9OqdR666mss5Dpt3CjsOVWM0mhvTuQvXBr1g5dzp2q5k5z81h6JCBUgMRIsnFJZForYuVUlONBEGwg1wptUprPco4PlMp9YBxy7Qm7jntWKJZ+vbBV1aOCun81VVVWPr2afPnDtdvEdqc1LdbJ748WEG6zYyzOlAL0ASWEgl+wLeaA7PhQzvdF6zdRa3bR2aaBU2g5lJW7aZbpsJuMWM1K6Mc2Ly3nMPOOpSCDLuZzDTryXJ8fk1VXWASYXAYr8mkuP27A9h2sJLS8lq8Pk3/LD//ePFRtq//gIl3TeanP/lR4DkkiQiR9GQeSQyE9pGozEx0VRW6ooLObdxH0th8jeo6LwPOzEQZw57KqlzsPV7DvuPV2K1mat0+lAK7xYTH58diNtEtw47NamJY365orfnPjuMM6dMFr0+zcXcZNS4vPh1IOJ3sZmwWE16fxu/XZKRZOFHjRimFz68ZcnYXzumeidaaIxUubvvugEY75H0+H2MmPczfF/6Zgu9dw5RJd9KlS5Y0YwmRQEk7j6Q9sw/Phfum1B+1NW5sm/ePNDZP5EhFXb3mpOxMO1aLCaUgp0cm63Yep8btw+X149fg136OVblQCvp2y6DG5aWi1sOGkjLqvD7MJoXf+LyhAZfHj8er6XdGOlmdbByvcnOixoPXpwHNjiNVZHWyYTErzuqaXq/jvbS89uQ8lsr9XzJ+/HgwW5g/bx7nn38eFskgQqQcSSQxYh+eG/eO9cYWMcywW05rTgrsHHgOi/69h8q6+utV+XVg3obNYmLTnnJ8fj82i6LS5cHvP33FKqVgQI8Myqo9DDy7Exl2K+XVLmrdfsxK4fL4+HxfOWd17RR2t8KDpYe4Ydb9HPtqHVPvm8L3vve9tnyZhBBtTBJJCmtsTavzzsrimtxeLFi7i4+2H6XW7cNigr3HqjhSEah5aF0/QWjA6/OjIDD0VoHP3/hzH6ty4fMHnm/P8Wo62a2kWf1U1noxmyDdZqF7Z3u93Qo7WRX/fmsRaxY/z1Xfu4bCGcvpkpXZVi+PECJOJJGksKb2LAc4WunC69f4tZ8ql6bK5UMDNjN4fKfKUQQSSc8u6RyvduHy+snOtFNR6z0t4QQ3F1RKkWY1UVXnpaLGQ2aaGZ82kWG3cH6vLLpm2DhS4QICNafaA1/w18LpOLpk8eJf5vONc8/FLGtPC9EuSCJJYeHWtLr43G6nhuzWuvH6/IAKDMvVGp8Gd0gSMZsAHUgMzlo3Lk9g0qLfmGPi9wfuCfJpAm1hGrIzbIy98lweWbGFiloPdosZkwm+LK3EalZ8o2dnDh8+zIcvPszOT//DPfdO4Uc//L5sdStEOyOJJMU1XNMq2Bfh1xq3MQlR60BSCOdk85XW+FyBDOPXcKzShdWscDeskgAo6JZp47xeWQzu7eDnl57Dc6u/5mhFHRazIsNuodrt4/1lL1E4fiE/ue56nn7kNbIyM06OJBNCtB+SSGIkGdbaCh3FZTEHlkEJJgGfX0fc5il43qIC10OgLyVIAVazopPdzPFqN9fk9mLrfiertxzGrBRWiwmfX3Ng+0b2vjOHM7ufwQvz/8I3vzFA5oMI0Y5JIomBZFlr66vSCirrPFTUeqmsrb8ESejbuOL0SkYok0nh9Wm8/vpXKSPBpNssdMuo35FuMpnI1E4+WzmH4yWfcu+9U/jxNd+XSYVCdADS3RkDoWttKZMp8D0ri9qiorjFsHW/kyMVddS6fXh8fpqqf0R6W/f4Tr8i2DLWOd2KSSnO6xWYxV9aXku6WbPvn0W8O+NWhp17Nu+8sZLrfyR9IUJ0FFIjiYFkWGvrzc0H6dutE3uP1eD2BtbVCjIpsJgVHm/k5i2on2hMgJ9Av4lJQWWtBw0nl5Cv2f85C5/7P87scSavLHyZAf37Y5EEIkSHIokkBhK51lZQsFnL7fPj8/vrJQOb0XfREiqYSQgkEwUMP6crDlXNTTfdyQcffsTku+7h+98rwCbjeYXokOR/fgykjxmDrqjAX1GB9vsD3ysqSB8zJi7PH9qsZTap0yYS1nn8eHzR1UYaCpaljC/t83D4P8UMGzaMM844g5WvreDaH14tSUSIDkxqJDFgH56L59prqS6ci+/oMczdzyBj7Lg2X7AxdM/z7AwrB8rrTlv+JFY0ULlrE/venYOrX19efvnluO5AKYRIXpJIYsC1aTN1K1diOe+bWPNGoKuqqFu5Eus3z2uTZNJw7artpRVofWq59lhznTjKwffm4Tmyg989cD9Xfefb9bbCFUJ0bJJIYqC2qAjt8+Hb/iX+qmpMmRmYepzZZjskNlz1N6uTlYPlNS1qumqK3+vm0H+Wceg/y/jZTb9g/NjHyeyUFuNnEUKkOkkkMeDeshXP119DTQ34/fhNJjh0GF1b2ybPV1pei9Ws2Ly3nBqXD7MJ3N7YphHnjvXs/fvz9Ol7DkWLF9H/nL6yR4gQIixJJDHg278fKioCM/YAfD6oqAgcbwN2q4nP9paTZrOQbjNR4/JFvilKLudh9rw7F8+xXfzuwQe4+spvY5dmLCFEEySRxIC/rCzwQ4MZ3CePtwmFCuwjRbXLG3LUCKWZpfm9bkr/VczhT1bwo5/exL3jnyQrIz1WwQoh2rG4JRKl1GjACeRprWeFOV9AYK/2UU3dE6mchPA20snd2PFWcnn8DD47i/3ltVS7vKfNEWluEnF+/Ql73plDvwHn8lzRq/Tp3Vt2KhRCRC0uicR480drvVoplaOUKtBarw69xjg3ral7AEekchKisbWk2miNKbvVxNeHKgPb5Xp8uL3+JmsiIXMK66krP8Sevz+Pv3w/f/j9b7jy21eQJs1YQohmite7xkigxPi5BMhr4T0Ry1FKjVVKrVdKrT969Girgk42W/c7uf/VTazbeTywL3utmzpPYLOq4Fc4DZOI3+Ni/wcL+WL+JK6+LI933ljOD6/6tiQRIUSLxKtpy9HgcbcW3hOxHK11IVAIkJ+f326WnQ3OHdlfVkNWJytllS48OrBsSXOUb/8Pe9+dS843vsnc4iWcfXYvacYSQrRKvBKJE8iOwT0tKaddCM4d8fg06TYTymQizQx1xp65kfZYrys7yJ6/P4+uKGX6Hx7i21dcKqOxhBAxEa9Eso5TtYkcYFUL73G0oJx2obS8lh5ZdjrZzbi9fixmhc/nR+tAEmmsO8bnqaN07RKObHiD0Tf9D5PG3UJmuj2+wQsh2rW4fCTVWhcDwQ5zgh3kSqmTicDoXM8P6WQ/7Z7GyukIzuqaTpXLyzndMvD6/KRbzSc3nvKHaeLSWlP25b/4fM44sj2HWLGsmKmTxkoSEULEXNyG/4Ybqhs61NdIEsVR3JMcQ37j7JrcXifX1+rlSOOLgxVhN6ACqDt+gN3vzEFVH2XGIw9zxWUXSzOWEKLNyITEFDG4t4OCIWfy1w92cehEbdhOdp+7joNrX+Xoxrf42X/fwvjb/kdqIEKINieJJAW8vnEff/1gF0cq6oDT54VorSn/8p/sfbeQ84cM4y/Ll9Kz55mYZTSWECIOJJEkuedXb+eVf+1ptBmr9tg+dr8zB3NdGY/NnM4lF10ozVhCiLiSRJLEtu53svg/e9FhhmT53LUc+PBvHNv8d/77ltsYe8vPpRlLCJEQkkiS1Nb9Th5ZsYVaT/3JIVpryr74iL2rChmSm8+Lry3nzB7dpRlLCJEwkkiS0Nb9Tp54+0sOltffz6Tm6B72vD0Hq7eSp56YzcgRedKMJYRIOEkkSWjB2l2UltdgNissPnDV1bD/g1co+3w1//O/t/Or/7mJzDRbosMUQghAEklS2nagArvVjNWs+GrTe5S8+wLDRlzIy68tp0f3M6QZSwiRVCSRJKnK0hI+W/40Jk81zzz5OCOG52KTZiwhRBKSRJJkTpw4wcH3Ctny0Vv8auyd/PLnN5Bul38mIUTyko+4SUJrzYIFCxg4cCB9M7y8uXIFY395kyQRIUTSk3epJPDpp58yYcIEqqqqeOKJJxg6dGjYuSNCCJGMpEaSQE6nk0mTJnH11Vdz9dVX89JLLzFkyBBJIkKIlCKJJAH8fj8vvvgiAwcOpKysjGXLlnHddddhMsk/hxAi9UjTVpxt2rSJ8ePH43K5eOqppxg8eLDUQIQQKU0+AsdJeXk548eP5/vf/z7/9V//xYsvvsigQYMkiQghUp4kkjbm9/uZP38+AwcOpKqqiuXLl3PttdeilEwqFEK0D9K01YY+dbu5/pJL0FrzzDPPMHDgQKmBCCHanbjVSJRSo5VSBUqpqdGej/ZYomVMmlTvcbnfx/3Ocm5113H99dczf/58zj//fEkiQoh2KS6JRCk1GkBrvRpwKqUKIp2P9lg84o/EMe0+MiZNwt8li4W1tXz32FFMBVex4p13uOaaa6QZSwjRrsWrRjISKDF+LgHyojgf7bF6lFJjlVLrlVLrjx49GrNfIJLt3/kW12d0YsU3v8HcxYuZ+uijdO7cOW7PL4QQiRKvROJo8LhbFOejPVaP1rpQa52vtc7v3r17s4JsiaNHj3Lbbbfxk5/8hBtuuIF58+Zx3nnnSTOWEKLDiFcicQLZzTwf7bGE8Pl8PPfccwwePBiA5cuX84Mf/CDBUQkhRPzFa9TWOk7VJnKAVVGcd0R5LO7+/e9/M378eNLS0igsLOTcc89NRBhCCJEU4lIj0VoXAznBznGjsxyl1KrGzkd7LB7xBx05coRbb72Vn/70p/ziF7+QJCKEEMRxHonWelaYY6MinI/qWFvzer3MmTOHhx9+mOuuu47ly5eTkZER7zCEECIpyYTECNauXcv48ePJyspi/vz55OTkJDokIYRIKpJIGnHo0CGmTp3KmjVruPfeexk1alTkm4QQogOStbYa8Hq9PPnkkwwdOpROnTqxfPlySSJCCNEEqZGE+OCDD5gwYQLZ2dn85S9/oX///okOSQghkp4kEuDgwYPcd999fPjhh0yZMoWrrroq0SEJIUTK6NBNWx6Ph8cff5wLLrgAh8PBsmXLJIkIIUQzddgayZo1a5gwYQI9e/bkpZde4pxzzkl0SEIIkZI6XCLZv38/U6ZM4d///jf33Xcf3/nOdxIdkhBCpLQO07TldruZNWsWubm59OjRg6VLl0oSEUKIGOgQNZJVq1YxceJE+vTpw4IFC+jTp0+iQxJCiHajXScSt9vN6NGj2bBhA1OnTuVb3/qWLO8uhBAx1q6btr744gt69+5NcXExV1xxhSQRIYRoA6o9v7kqpY4Ce+L8tGcAx+L8nK2RSvGmUqyQWvGmUqyQWvGmUqwQiDdDax31zoDtOpEkglJqvdY6P9FxRCuV4k2lWCG14k2lWCG14k2lWKFl8bbrpi0hhBBtTxKJEEKIVpFEEnuFiQ6gmVIp3lSKFVIr3lSKFVIr3lSKFVoQr/SRCCGEaBWpkQghhGgVSSRCCCFaRRJJKyilRiulCpRSU1tyPp6iiLVAKbUq3nE1pql4lVIO4/xopdTMRMTXIJ5oXtuCZIgVov+7TIZ4o3hty5VSG5IhVogq3rzg3268Ywsnwv+zPKXUTuP1bfI1lkTSQsE/BK31asCplCpozvl4iiYW41xSiCLeG4BsrXWxcf3YOId4UhR/B3nAKON8nlIqJwFhhsYT1d+lcTwVYh2jtR6htZ4W3+hOF2W844y/25wU+FvI1loP0FqPAG4H5jZWliSSlhsJlBg/lwB5zTwfT8kUSzSajFdrXai1Do4syQESmQQjxbpRaz1NKeUASrTWJSRWxL8F4w0u0XFCdH+3jkS/IYdoMl7jA88GpVSO1npWsv8tNPhwmdNUvJJIWs7R4HG3Zp6PJ0eDx4mMJRqOBo/Dxmu8gZQl+D+ko8Hjxl7bfMDZppFEx9Hgcbh4m3zTiCNHg8fhYs0GypRSjX5ajiNHg8cN4x1gfJUppeYaHy4SqeHzN/b/bGyw9t8YSSQt5yTwR9zS8/HkJHliiYaT6OIdrbUe18axROIkiliNT3eOJGgbd9JEvEqpgiRq5nQS4bU1aqdOAk0zSf3aGnYa8W4AEtYka3AS3f+zUZEukETScus4ldFzgIYd1ZHOx1MyxRKNiPEqpUZrrWcZPyes/4kIsSqlZob04ThJfEKP9NqWGZ2vowm04yeyGTTSazs2JHkcj2NcjYnmPSHIQeJrqNH8P3M0PBaOJJIWCukwKzAerwYIjnxq7Hwyxmr8PBrIT4JPdRHjNY7PDI4mSVykUb22c4ES47wjpG8nIaL4u91oHMvm9KaPuIritV1CSCdxpOaXthble4Ij5HxS/y0YsoGySGXJzHYhhBCtIjUSIYQQrSKJRAghRKtIIhFCCNEqkkiEEEK0iiQSIYQQrWJJdABCpBJjNv04YCqwEVhsnOpGYCz+3NCh3sEZ1y2ZOBkc5kxg9n7ESWFCJIoM/xWiBZRSO4GZoXMBjMlb5ZxapJHmzCMylqIobHBsNIGF/iSRiKQlTVtCtIyz4QFj6YvVBGoswWOrmzEZdUQ0zyNEspFEIkRs5VN/KYyoGPtBJHr5FCFaRPpIhIgBo+9kJlAYsgZYHjAPo48jpM9jPYF1jX5GoI/FSWBJ7zwjoTjDNHEVEFiyZCSBhf8SuryGEKEkkQjRciNC1ibLJlATcQZPaq03KqWmAdOMx6uNzvfgsZKQ6xzGz7PCPE8+xl4mSqnVBPphJJGIpCGJRIiW29BwoUCl1Cql1IgmRmmVEUgKTgKjvqJxcs8VrbVTKdXigIVoC9JHIkRsFRF5nwlnM8s8bZOpJNgUSYiTJJEIkWSUUgnfz1uI5pBEIkTLOBo5Po6W9V+UEJjQCKdvdSujuURSkz4SIZohZGZ7DjAmpImpG5AHrNZaTzOuzSPQqZ5v7JJYYtybb4zOCm4TG+xwLwnZTTH0/pzg9cADxumZSqmZSbK3uujgZGa7EEKIVpGmLSGEEK0iiUQIIUSrSCIRQgjRKpJIhBBCtIokEiGEEK0iiUQIIUSrSCIRQgjRKpJIhBBCtMr/B+JnO+tF/xBlAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "gd.plot_persistence_diagram(BarCodes_Rips0);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "In this representation, $0$-dimensional features are the red points (holes of dimension $0$, namely connected components). The last connected component dies at infinity in the filtration (red point at the top). The $1$-dimensional features are represented in blue. \n",
    "\n",
    "The most persistent topological features are those points that are far from the diagonal. Further in the tutorial we give statistical methods to identify significant topological features.\n",
    "\n",
    "Note that this representation does not say which points are \"at the origin\" of a given feature. Moreover, a given topological feature (namely a homology class that corresponds to a hole of dimension $d$) is by definition a class of cycles defined on the point cloud and thus it can be represented by several cycles."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Bottleneck distance"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To exploit the topological information and topological features inferred from persistent homology, one needs to be able to compare persistence diagrams.\n",
    "\n",
    "We see a persistence diagram as the union of its points and of the diagonal, where the point of the diagonal are counted with infinite multiplicity. \n",
    "\n",
    "A matching  between two diagrams $\\mathrm{dgm}_1$ and $\\mathrm{dgm}_2$ is a subset $m \\subseteq \\mathrm{dgm}_1 \\times \\mathrm{dgm}_2$ such that every points in $\\mathrm{dgm}_1 \\setminus \\Delta$ and $\\mathrm{dgm}_2 \\setminus \\Delta$ appears exactly once in $m$. \n",
    "\n",
    "<img src=\"Images/MatchingDiag.png\"  style=\"width: 400px;\"/>\n",
    "\n",
    "The **Bottleneck distance** between $\\mathrm{dgm}_1$ and $\\mathrm{dgm}_2$ is then defined by\n",
    "\n",
    "$$\n",
    "\\mathrm{d}_b (\\mathrm{dgm}1, \\mathrm{dgm}_2) =\\inf_{\\scriptsize{\\mbox{matching}} \\ m} \\max_{(p,q) \\in m} \\| p-q \\|_\\infty.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let us compute the Rips complex filtration for another configuration of protein"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [],
   "source": [
    "D1 = dist_list[1]\n",
    "\n",
    "skeleton_protein1 = gd.RipsComplex(\n",
    "    distance_matrix = D1.values, \n",
    "    max_edge_length = 0.8\n",
    ") \n",
    "\n",
    "Rips_simplex_tree_protein1 = skeleton_protein1.create_simplex_tree(max_dimension = 2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "and the barcode for this filtration:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {},
   "outputs": [],
   "source": [
    "BarCodes_Rips1 = Rips_simplex_tree_protein1.persistence()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The bottleneck distance between the two persistence diagrams can be computed using the <code>bottleneck_distance()</code> function. The bottleneck distance is computed per dimension (for dimension 1 in the example below). We can give in argument of the function the persistence_intervals for a given dimension, which can be computed using the `persistence_intervals_in_dimension()` function. \n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.05052142999999998"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "I0 = Rips_simplex_tree_protein0.persistence_intervals_in_dimension(1)\n",
    "I1 = Rips_simplex_tree_protein1.persistence_intervals_in_dimension(1)\n",
    "\n",
    "gd.bottleneck_distance(I0, I1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By default, the function uses an expensive algorithm to compute the exact bottleneck distance. It is also possible to compute an approximated bottleneck distance (additional approximation error given in argument), which is usually a lot faster to compute. See the [documention](https://gudhi.inria.fr/python/latest/bottleneck_distance_user.html#gudhi.bottleneck_distance) for more details."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "0.05061900451758505"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "gd.bottleneck_distance(I0, I1, 0.01)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### MDS on the Bottleneck distances"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "A simple approach for comparing the MPB configurations is to use the Bottleneck distances between configurations. We first compute the Rips complex filtration for each MPB:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0\n",
      "1\n",
      "2\n",
      "3\n",
      "4\n",
      "5\n",
      "6\n",
      "7\n",
      "8\n",
      "9\n",
      "10\n",
      "11\n",
      "12\n",
      "13\n"
     ]
    }
   ],
   "source": [
    "persistence_list0 = []\n",
    "persistence_list1 = []\n",
    "\n",
    "i=0\n",
    "for d in dist_list:\n",
    "    print(i)\n",
    "    rips_complex = gd.RipsComplex(\n",
    "        distance_matrix = d.values, \n",
    "        max_edge_length = 0.8\n",
    "    )\n",
    "    simplex_tree = rips_complex.create_simplex_tree(max_dimension = 2)\n",
    "    diag = simplex_tree.persistence()\n",
    "    persistence_list0.append(simplex_tree.persistence_intervals_in_dimension(0))\n",
    "    persistence_list1.append(simplex_tree.persistence_intervals_in_dimension(1))\n",
    "    i = i + 1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Next, we compute the matrix of Bottleneck distances for dimensions $0$ and $1$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {},
   "outputs": [],
   "source": [
    "l = len(files_list)\n",
    "B0 = np.zeros((l, l))\n",
    "B1 = np.zeros((l, l))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in range(l):\n",
    "    for j in range(i):\n",
    "        B0[i,j] = gd.bottleneck_distance(persistence_list0[i], persistence_list0[j])\n",
    "        B1[i,j] = gd.bottleneck_distance(persistence_list1[i], persistence_list1[j])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {},
   "outputs": [],
   "source": [
    "B0 = B0 + B0.transpose()\n",
    "B1 = B1 + B1.transpose()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We apply a Multidimensional Scaling method (from the [<code>scikit-learn</code>](http://scikit-learn.org/stable/) library)  to visualize the Bottleneck distances:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {},
   "outputs": [],
   "source": [
    "mds = manifold.MDS(\n",
    "    n_components = 2, \n",
    "    max_iter = 3000, \n",
    "    eps = 1e-9, \n",
    "    dissimilarity = \"precomputed\", \n",
    "    n_jobs = 1\n",
    ")\n",
    "pos = mds.fit(B0).embedding_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAD3CAYAAADxJYRbAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAARmUlEQVR4nO3dz24a1//G8edjNVJkKZJ/dr2LDJmsunXGd0DvgKgL9uQObFVZW1V8B2EdL9r6DuAOgtl2VZpE3vlrB6mSFalSPr/FMDZgMOBh+HfeL8nCzGGGw2Ty+Mw5zBlzdwEAwrCx6AoAAOaH0AeAgBD6ABAQQh8AAkLoA0BACH0ACMgPi67AOD/++KMXi8VFVwMAVsb5+fn/3H13WNnSh36xWFSz2Vx0NQBgZZjZ51FldO8AQEAIfQAICKEPAAEh9AEgIJlD38zKZlYys8MR5SUzqw8s+2pm52b2Luv7A8AinZ5KxaK0sZE8np4uukYPyxT6ZlaWJHdvSOqYWWnwNd2yQa/d/ZW7H2V5fwBYpNNTqVqVPn+W3JPHanW5gz9rS/9AUrv7e1vS/oTrbZlZlPG9AWCh3r6Vbm76l93cJMuXVdbQ3xp4vjPhetuSrs3s/bBCM6uaWdPMmpeXl1nqBwC5+fJluuXLIGvod5QE+FTcvebuHSVdQuUR5bG7x7u7Qy8qA4CF29ubbvkyyBr6H3XX2o8k1Ue/NNFtxadBf5Xx/YHHW7UROCyd42Npc7N/2eZmsnxZZQp9dz+TFKUDuOmgbe+3dboBH/cE/R/qGfTtbgOYr1UcgcPSqVSkWk0qFCSz5LFWS5YvK1v2e+TGcezMvYOZKxaToB9UKEifPs27NsBMmdm5u8fDyrg4C2FaxRE4YAYIfYRpFUfggBkg9BGmVRyBA2aA0EeYVnEEDpiBpb+JCpCbSoWQR3Bo6QNAQAh9AAgIoQ8AASH0ASAghD4ABITQB4CAEPoAEBBCHwACQugDQEAIfQAICKEPAAEh9AEgIIQ+AASE0AeAgBD6ABAQQh8AAkLoA0BAMoe+mZXNrGRmhyPKS2ZWn2YdAEA+MoW+mZUlyd0bkjpmVhp8TbdsqnWAZXR6KhWL0sZG8nh6uugaAdPL2tI/kNTu/t6WtJ/TOsBCnZ5K1ar0+bPknjxWqwQ/Vk/W0N8aeL4zi3XMrGpmTTNrXl5ePrJqwOy8fSvd3PQvu7lJluMOZ0PLL2vodyRtz3odd6+5e+zu8e7u7iOrBszOly/TLQ8RZ0OrIWvof9Rdyz2SVB/90kzr4LFoes3E3t50y0PE2dBqyBT67n4mKUoHY9NB295v63QHbuOeAdyh6yAHNL1m5vhY2tzsX7a5mSxHgrOh7ObSRnP3pf559eqVT+3DB/dCwd0sefzwYfptrINCwT2J+/6fQmHRNVtJHFYP43DL5sMH983N/n23ufm440xS00dkqiXlyyuOY282m5OvkLZue88zNzelWk2qVGZfwWW2sZEcO4PMpO/f518frDX+62VTLCYn44MKBenTp+m2ZWbn7h4PK1u/K3LpWLxDRzTmqFJJAr5QSNoVhQKBP415dY+tX+jTsXiHjmjMWaWStEq/f08eCfzJzauNtn6hT+v2Dk0vYGXMq422fqFP67YfTS9gJcyrjfbDbDe3BNI99PZt0qWzt5cEPmEHYMlVKvlH1fqFvjSfPQcAK2j9uncAACMR+gAQEEIfAAJC6ANAQAh9AAgIoQ8AASH0ASAghD4ABITQB4CAEPoAEBBCHwACQugDQEAIfQAICKEPAAEh9AEgIIQ+AAQkc+ibWdnMSmZ2OGm5mX01s3Mze5f1/QEAk8sU+mZWliR3b0jqmFlpwvLX7v7K3Y+yvD8AYDpZW/oHktrd39uS9ics3zKzKON7AwCmlDX0twae70xYvi3p2szeD9uomVXNrGlmzcvLy4xVBACksoZ+R0mAT1Xu7jV37yjp8imPKI/dPd7d3c1YRQBAKmvof9Rdaz6SVB9X3m3Fp0F/lfH9AQBTyBT67n4mKUoHaLsDtjKz+gPlf6hnULf7GgDAHJi7L7oOD4rj2JvN5qKrAQArw8zO3T0eVsbFWQAQEEIfAAJC6ANAQAh9AAgIoQ8AASH0ASAghD4ABITQB4CAEPoAEBBCHwACQugDQEAIfQAICKEPAAEh9AEgIIQ+AASE0AeAgBD6ABAQQh8AAkLoA0BACH0ACAihDwABIfQBICA/ZN2AmZUldSTtu/vJJOXj1gEA5CNTS78b3nL3hqSOmZXGlY9bBwCQn6zdOweS2t3f25L2Jygftw4AICdZQ39r4PnOBOXj1pGZVc2saWbNy8vLLPUDAPTIGvodSdtTlo9bR+5ec/fY3ePd3d0s9QMA9Mg6kPtRdy33SFJ9gvKtMesAAHKSqaXv7meSonQwtjs4KzOrjyoftQ4AIH/m7ouuw4PiOPZms7noagDAyjCzc3ePh5VxcRYABITQB4CAEPoAEBBCHwACQugDQEAIfQAICKEPAAEh9AEgIIQ+AASE0AeAgBD6ABAQQh8AAkLoA0BACH0ACAihDwABIfQBICCEPgAEhNAHgIAQ+gAQEEIfAAJC6ANAQAh9AAhI5tA3s7KZlczscNJyM/tqZudm9i7r+wMAJpcp9M2sLEnu3pDUMbPShOWv3f2Vux9leX8AyMPpqVQsShsbyePp6aJrNDtZW/oHktrd39uS9ics3zKzKON7A8DMnZ5K1ar0+bPknjxWq+sT/FlDf2vg+c6E5duSrs3s/bCNmlnVzJpm1ry8vMxYRQCY3Nu30s1N/7Kbm2T5POR9lvHDuBeYWXXI4nbaZaMkwEcZWu7ute62O2ZWdvezIeU1SYrj2MfVEQBm5cuX6ZbPUnqWkf7RSc8yJKlSmc17jA39NKBH+Ki71nwkqT6uvPtH5Lob9FfTVBYA8ra3l4TtsOV5e+gsY1ahn6l7pxvcUTpA2239y8zqD5T/oZ5B3cFWPgAs0vGxtLnZv2xzM1met3mcZYxt6Y/j7idDlv08qtzdO5Ia3acNAcASSVvUb98mYbu3lwT+rFraD5nHWQYXZwHAgEpF+vRJ+v49eZxH4EvzOcsg9AFgSVQqUq0mFQqSWfJYq832j07m7h0AwOxUKvmeWdDSB4CAEPoAEBBCH8B8rfPENiuAPn0A8zOPS07xIFr6AOZn0RPbgNAHMEeLnNgGkgh9APM06tLSeUxsA0mEPoB5WuTENpBE6AOYp3lccooH8e0dAPOV9yWneBAtfQAICKEPAAEh9AEgIIQ+AASE0AeAgBD6ABAQQh8AAkLoA0BACH0ACAihDwAByRz6ZlY2s5KZHY4oL5lZfZp1AMwAd6jCEJlC38zKkuTuDUkdMysNvqZbNtU6ADJK71D1+bPkfneHKoI/eFlb+geS2t3f25L2c1oHwDS4QxVGyBr6WwPPd2axjplVzaxpZs3Ly8tHVg0IGHeowghjp1Y2s+qQxe20e0bS9pTvOXYdd69JqklSHMc+5fYB7O0lXTrDliNoY0O/G8CjfNRdyz2SVB/90kzrAJjG8XHSh9/bxcMdqqCM3TvufiYpSgdj00Hb3m/rdAdu454B3KHrAJgh7lCFEcx9uXtP4jj2ZrO56GoAwMows3N3j4eVcXEWAASE0AeAgBD6ABAQQh8AAkLoA0BACH0ACAihDwABIfQBICCEPgAEhNAHgICMnXBtVfz333+6uLjQt2/fFl2VuXn69KmeP3+uJ0+eLLoqAFbE2oT+xcWFnj17pmKxKDNbdHVy5+66urrSxcWFXrx4sejqhOH0NLkJyZcvyRTFx8dMYIaVszbdO9++fdPOzk4QgS9JZqadnZ2gzmwWitsPYk2sTehLCibwU6F93oXi9oNYE2sV+kBuuP0g1kS4oX96KhWL0sZG8jjj0/RGo6Gff/55Jts6OjrSycnJTLaFRxp1m0FuP4gVE2boz6F/tlQqaWtraybb+uWXX2ayHWRwfJzcbrAXtx/ECgoz9HPqnz05OVGr1dLZ2dm9slqtplarpVotueVwq9VSo9G4/el0Ojo5OVGj0bh9Tfq80eCOkgvH7QexJsIM/Rz6Z8/OzhRFkfb39/X777/3lZ2cnCiOY+3v7yuKItVqtdvXlEolRVGk3377TaVSSaVSSefn56rVatrf379dhiVQqUifPknfvyePBD5WUJihn0P/bL1eVxRFkqQ///zzXlna1RNFker1un799VfV63W9fPlSnU5HrVZL19fXarVaevPmTd/2AGBWwgz9HPpnX758qevra0lSp9PpK9vf31e73ZYktdttHRwcqNFo6N27dzo/P+8b9E3PBg4ODtRqtSTpdrsAkNXaXJE7lfS0fIZXVx4eHuro6EhSEvpRFKnVaqndbuvdu3e3375ptVp9r5WkcrmsKIr6vqFzeHh4O0bQarVUr9dVrVZnNjgMIEzm7tk2YFaW1JG07+73vldoZiVJR+7+c8+yr5LakhrufjS4Tq84jr3ZbI6tx19//aWffvppytqvvlA/N4DRzOzc3eNhZZm6d7qBL3dvSOp0A75Pt2zQa3d/NS7wAQCzlbVP/0BJi13dx/0J19syM0YpAWDOsob+1sDznQnX25Z0bWbvhxWaWdXMmmbWvLy8zFI/AECPsQO5ZlYdsriddukoCfCpuHutu+2OmZXd/WxIeU1K+vSn3T4AYLixoZ8G9AgfddfajyTVx22v+0fkuhv0VxPUEQAwI5m6d7rBHaUDuOmgrZndhn93sDdOB30l/aGeQd/BVj4AID+ZL85y9xN3b/SeEfR+PdPdz9z9/9Jwd/dO9/WNYV/xnJecJ9kEgKUU5MVZ6SSb6Zxr6SSbEtOpAFhvQU7DkNdNkAZn0mw0Gnr16tXtzJnpbJmDM2qm0zA0Go3bq3ABIA9Bhn4eN0EaNpNmqVTS9va2SqWSqtWq3rx5I0n3ZtQslUq6vr5WqVRSuVy+N0snAMxKkKGfx02Qhs2kKalvrpwoitRut+/NqJmWAfPAeFbYggz9PG6CNGwmTal/xs12u60oiu7NqAnMyxxuGoclF+RAbg6TbA6dSVPSbYu+2Wzq/fvkAuR0Bs3U9vb27YycjUZDrVZLnU6HGTUxcw+NZ/ElhjBknmUzb6s+y+br16/v3VRllpb1c2M5bWwkLfxBZskNwbAecptlEw9LW+1ptw+waHmMZ2G1EPo5KpVK+vvvv+m3x+PkMOKax3gWVstahf6yd1XNWmifNyg5jbhWKlKtJhUKSZdOoZA8pz8/HGvTp//PP//o2bNn2tnZkZnNoWaL5e66urrSv//+qxcvXiy6Opi1YjEJ+kGFgvTp07xrgxXzUJ/+2nx75/nz57q4uFBI8+8/ffpUz58/X3Q1kIc8riAEtEah/+TJE1q8WB97e8Nb+oy4IqO16tMH1gYjrsgJoQ8sI0ZckZO16d4B1k6lQshj5mjpA0BACH0ACAihDwABIfQBICCEPgAEZOmnYTCzS0npVSo/SvrfAquzbNgf97FP+rE/+oWyPwruvjusYOlDv5eZNUfNJxEi9sd97JN+7I9+7A+6dwAgKIQ+AARk1UK/tugKLBn2x33sk37sj37B74+V6tMHAGSzai19AEAGhD4ABGTpQ9/MymZWMrPDEeUlM6sPLPtqZudm9m4+tZyfR+6PB9dZZRPsj3vl63R8PPLzczys6fEwiaUOfTMrS5K7NyR1zKw0+Jpu2aDX7v7K3Y/yruM8PWZ/TLLOqhr32R4oX4vj4zGfn+NhfY+HSS116Es6kNTu/t6WtD/heltmFuVTpYV6zP547D5cBeM+26jydTk+HvP5OR7ul6/L8TCRZQ/9rYHnOxOuty3p2szez7Y6C7c18HyS/fGYdVbF1sDzwc82qnxdjo+tgeeTfP5x66yyrYHnoR0PE1n4nbPMrDpkcTs9BVPyDzIVd691t90xs7K7n2Wr5fzksD8es87SyLg/hpav8vExoKPpP/+4dVZZR2EfDxNZeOinO3yEj7r76xxJqo9+aaIbEtfdf7irzBWcs1nvj0euszQy7o975at+fAyY+vN3nz+0zioL/XiYyFJ373T/IaJ0wCUdpOz9dkp3cCZOB2kk/aGeQZp1+qv9mP0xap11MG5/jChfm+PjMZ+f42F9j4dJcUUuAARkqVv6AIDZIvQBICCEPgAEhNAHgIAQ+gAQEEIfAAJC6ANAQP4fNGnzUkQ67pEAAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(pos[0:7, 0], pos[0:7, 1], color = 'red' , label = \"closed\")\n",
    "plt.scatter(pos[7:l, 0], pos[7:l, 1], color = 'blue', label = \"open\")\n",
    "plt.legend(loc = 3, borderaxespad = 1);"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {},
   "outputs": [],
   "source": [
    "pos = mds.fit(B1).embedding_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAX0AAAD3CAYAAADxJYRbAAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjMuMSwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy/d3fzzAAAACXBIWXMAAAsTAAALEwEAmpwYAAAR6klEQVR4nO3dv24bV/rG8ecVfl5oBQQgpAhbxJDkcZWWGd0BDewFMEnBeulmaxmBayOQ6i3W7FkE8R2Qd2CKxTapzESGt1hoZRAIoDUcwO+vmKFkUaRIavhnOOf7AQyR5xWp147y6OjMzBlzdwEAwrCx6gYAAMtD6ANAQAh9AAgIoQ8AASH0ASAghD4ABOT/Vt3AJF9++aUfHBysug0AWBunp6f/dffdUbXch/7BwYE6nc6q2wCAtWFmZ+NqLO8AQEAIfQAICKEPAAHJvKZvZlVJfUlldz+Zpm5mlbT8xN2fZe0BADCdTDP9NNDl7m1J/c/CfGzdzMpKwr4tqWxmUZYeACSaTengQNrYSD42m6vuCHmUdXnnUFIvfdyTVJ5Ud/euuz8zs5Kknrv3BCCTZlOq16WzM8k9+VivE/y4LWvol4ae78xQj5Us+9xiZnUz65hZ5/z8PEt/QBCeP5cuL2+OXV4m48DnsoZ+X9L2ferp8k5psAQ0VGu4e+zu8e7uyOsLMAcsBxTH27ezjSNcWUP/ta5n85Gk1qS6mR2bWT0d6+vuHxpYEJYDimVvb7ZxhCtT6Lv7K0nR4ABuOnuXmbXuqL+U1EvHSu7eyNID7oflgGJ58ULa2ro5trWVjAOfs7zfLjGOY2cbhvnb2Ehm+MPMpE+flt8Psms2kx/ab98mM/wXL6RabdVdYRXM7NTd41G13O+9g8XY20uWdEaNYz3VaoQ8JuOK3ECxHACEidAPVK0mNRrS/n6ypLO/nzxnpggUG8s7AWM5AAgPM30ACAihDwABIfQBICCEPgAEhNAHgIAQ+gAQEEIfAAJC6ANAQAh9AAgIoQ8AASH0ASAghD4ABITQB4CAEPoAEBBCHwACQugDQEAIfQAICKEPAAEh9AEgIIQ+gJVrNqWDA2ljI/nYbK66o+LixugAVqrZlOp16fIyeX52ljyXpFptdX0VVeaZvplVzaxiZkfT1M2slI5Vzew469cHsN6eP78O/IHLy2Qc85cp9M2sKknu3pbUN7PKFPXvJG27+6v0c+pZegCw3t6+nW0c2WSd6R9K6qWPe5LKk+ru3nD3RjoWSWpn7AHAGtvbm20c2WQN/dLQ851p62YWSXrv7r2hz5GZ1c2sY2ad8/PzjC0CyLMXL6StrZtjW1vJOOYva+j3JW3fs15196ejCulvA7G7x7u7u9k6BJBrtZrUaEj7+5JZ8rHR4CDuomQ9e+e1rmfzkaTWNHUzq7r7Sfq4kq75AwhUrUbIL0ummX56MDYaHMAdhLeZtcbV08fHZnZqZqeZugcAzCTzefqDGfvQ2JNx9fQHw+OsXxcAMDuuyAWAgBD6ABAQQh8AAkLoA0BACH0ACAihDwBDirzVM1srA8Bnir7VMzN9APhM0bd6JvSBeyryEkDIir7VM6EP3MNgCeDsTHK/XgIg+Ndf0bd6JvSBeyj6EkDIir7VM6EP3EPRlwBCVvStnjl7B7iHvb1kSWfUONZfkbd6ZqYP3EPRlwBQXIQ+cA9FXwJAcbG8A9xTkZcAUFzM9AEgIIQ+AASE0AeAgBD6ABAQQh8AAkLoLwq7cQHIIU7ZXISib8gNYG0x018EduMCkFOE/iKwGxeAnMoc+mZWNbOKmR1NW0+ft7J+7dwq+obcANZWptA3s6okuXtbUt/MKtPU0+fFxW5cAHIq60z/UFIvfdyTVJ6xXkzsxgUgp7KevVMaer4zY30kM6tLqkvS3rouibAbF4AcyjrT70vazlAfyd0b7h67e7y7u3vP1gAAw7KG/mtdz+YjScMHZyfVAQBLlCn03f2VpGj4AO3gzJw76lVJ8eBALwBgOczdV93DneI49k6ns+o2AGBtmNmpu8ejalycBQABIfQBICCEPgAEhNAHgIAQ+gAQEEIfAAJC6ANAQAh9AMiTBd9qldslAkBeLOFWq8z0ASAvlnCrVUIfAPJiCbdaJfQBIC+WcKtVQh8A8mIJt1ol9AEgL5Zwq1XO3gGAPFnwrVaZ6QNAQAh9AAgIoQ8AASH0ASAghD4ABITQB4CAEPoAEBBCHwACQugDQEAIfQAISOZtGMysKqkvqezuJ9PUJ70GALAYmWb6aXjL3duS+mZWmVSf9BoAwOJkXd45lNRLH/cklaeoT3oNAGBBsoZ+aej5zhT1Sa+RmdXNrGNmnfPz8yz9AQA+kzX0+5K2Z6xPeo3cveHusbvHu7u7WfoDAHwm64Hc17qeuUeSWlPUSxNeAwBYkEwzfXd/JSkaHIxND87KzFrj6uNeAwBYPHP3VfdwpziOvdPprLoNAFgbZnbq7vGoGhdnAUBACH0ACAihDwABIfQBICCEPgAEhNAHgIAQ+gAQEEIfAAJC6ANAQAh9AAgIoQ8AASH0ASAgxQz9ZlM6OJA2NpKPzeaqOwKAXMh8Y/TcaTalel26vEyen50lzyWpVltdXwCQA8Wb6T9/fh34A5eXyTgABK54of/27WzjABCQ4oX+3t5s4wAQkOKF/osX0tbWzbGtrWQcAAJXvNCv1aRGQ9rfl8ySj40GB3EBQEU8e0dKAp6QB4BbijfTBwCMRegDQEAIfQD5xdX1c1fMNX0A64+r6xeCmT6AfOLq+oXIHPpmVjWzipkdTVtPn7eyfm0ABcbV9QuRKfTNrCpJ7t6W1DezyjT19DkAjMfV9QuRdaZ/KKmXPu5JKs9YB4DRuLp+IbKGfmno+c6M9ZHMrG5mHTPrnJ+f37M1AGuNq+sXYmLopwE8/GewjNOXtH3HyyfVR3L3hrvH7h7v7u7O+vKl4WwyYMFqNem336RPn5KPBH5mE0/ZdPfGHeXXup7NR5KGD85Oqq8tziYDsI4yLe+4+ytJ0fAB2sGZOXfUq5LiwYHedcTZZADWkbn7qnu4UxzH3ul0Vt3GLRsb0qh/OrPkN1EAWBUzO3X3eFSNi7PuibPJAKwjQv+eOJsMwDoi9O+Js8kArCM2XMuAe7UAWDfM9AEgIIQ+AASE0AeAgBD6ABAQQh8AAkLoA0BACH0ACAihDwABIfQBICCEPgAEhNAHgIAQ+gAQEEIfAAJC6ANAQAh9AAgIoQ8AASH0ASAghD6KqdmUDg6kjY3kY7O56o6AXOB2iSieZlOq16XLy+T52VnyXOL+lggeM30Uz/Pn14E/cHmZjAOBI/RRPG/fzjYOBCTz8o6ZVSX1JZXd/WRS3cxKkipp+dDdn2XtAbhhby9Z0hk1DgQu00w/DXS5e1tS38wqU9S/k7Tt7q/Sz6ln6QG45cULaWvr5tjWVjIOBC7r8s6hpF76uCepPKnu7g13b6RjkaT28JuaWd3MOmbWOT8/z9giglOrSY2GtL8vmSUfGw0O4gLKvrxTGnq+M23dzCJJ7929N/Q5Sn8oNCQpjmPP2CNCVKsR8sAIE0N/zPJLb7BkI2n7jpffVa+6+9NJXx8AMD8Tl3cGyzFDfwZLMq91PZuPJLWGXj6ybmbVwUHf4eMAAJAbBbzIL9Pyjru/MrOjQXAPfhiYWcvdn4yqp4+PzeyH9G3mcvbOH3/8oXfv3unDhw/zeLu1sLm5qYcPH+rBgwerbgUonoJe5Gfu+V4yj+PYO53OxM/79ddf9cUXX2hnZ0dmtoTOVsvddXFxod9//12PHj1adTtA8RwcjD71d39f+u23ZXczEzM7dfd4VK0wF2d9+PAhmMCXJDPTzs5OUL/ZAEtV0Iv8ChP6koIJ/IHQ/r7AUo27mG/NL/IrVOgDwNwU9CK/cEN/wUfl2+22njx5Mpf3evbsmU5Obu1wAWCRCnqRX5ihPzgqf3YmuV8flZ9j8FcqFZVKpbm81/fffz+X9wEwo1otOWj76VPycc0DXwp1P/27tt7N8B/15ORElUpFvV5P1Wr1Rq3RaCiOY3U6HdXrdXW7Xb1///6qHsexGo2GyuWyer2e6vW6Tk5OVC6X1e12790TAHwuzJn+Ao7Kv3r1SlEUqVwu66effrpROzk5URzHKpfLiqJIjUbj6nMqlYqiKNKPP/6oSqWiSqWi09PTqx8AgzEAmIcwQ38BR+VbrZaiKJIk/fzzz7dqg6WeKIrUarX0ww8/qNVq6fHjx+r3+1cz/263q6dPn954PwCYlzBDfwFH5R8/fny1XNPv92/UBks2ktTr9XR4eKh2u63j42Odnp7eOOg7+G3g8PDwalnn82UgAMgizDX9wbr98+fJks7eXhL4Gdbzj46O9OxZsqNEv99XFEXqdrvq9Xo6Pj6+Ovum2+3e+FxJqlariqLoxhk6R0dHOjk5UbfbVbfbVavVUr1en9vBYQBhKsw2DL/88ou+/vrrJXSUL6H+vQGMF8Q2DACAyQh9AAgIoQ8AASH0ASAghD4ABITQB4CABBv6ozbZvLiQ/vUvqdNJPl5crLpLAJivIEN/1Cabf/ub9M9/Sh8/Jp/z8WMyTvADKJIgQ3/UJpv/+5/0j3/cHPv0Sfr3v6d/30ajoW63q0ajISnZU/+bb75Ru91Wo9FQu92WlFyxe3JycmP8yZMnarfbV1fhAsAiBBn64zbT/M9/bo8NZv6TjNpJs1KpaHt7W5VKRfV6XU+fPpWkWztqVioVvX//XpVKRdVq9dYunQAwL0GG/rjNNP/yl9tjf/rTdO85aidNSTf2yomiSL1e79aOmoMaACxakKE/apPNP/9Z+vvfb45tbEhffTXde47aSVO6ueNmr9dTFEW3dtQEgGUJcpfNcZts/vWvyRr+x4/JDP+rr6Sdnenec9ROmpKuZvSdTkcvX76UdL2D5sD29vbVjpztdlvdblf9fp8dNQHMHbtsLti3335766Yq85TXvzeA1VnoLptmVjWzipkdTVtPn1fM7Djr18+zwax9sOwDAKuWKfTNrCpJ7t6W1DezyqS6mZUlPUnHymZW2EXtSqWiN2/esG6fV6Ou0AMKLutM/1DSYBrbk1SeVHf3rrs/M7OSpJ6735oGm1ndzDpm1jk/P5+6mbwvVc1baH/fuRp1hV69TvCj8LKGfmno+fBhz7vqsaT+qDd194a7x+4e7+7uTtXI5uamLi4ugglCd9fFxYU2NzdX3cp6GnWF3uVlMg4U2MSzd8ysPmK4N1iykbR9x8vH1t29bWbfmlnV3V9N0eudHj58qHfv3mmW3wzW3ebmph4+fLjqNtbTuCv0xo0DBTEx9N29cUf5ta5n85Gk1qR6evD2Tfq+fd39Q2NqDx480KNHj+bxVgjB3l6ypDNqHCiwTMs76Qw9GhzATWf/MrPWHfWXknrpWGnCDxVgMUZdobe1lYwDBVaY8/SBmTWbt6/QG1y5B6yxu87TD/KKXEBSEvCEPAIT5N47ABAqQh8AAkLoA0BACH0ACAihDwAByf0pm2Z2LmnEVTRL9aWk/664h1msU7/r1KtEv4tGv/Ox7+4j97DJfejngZl1xp3zmkfr1O869SrR76LR7+KxvAMAASH0ASAghP501m1/oHXqd516leh30eh3wVjTB4CAMNMHgIAQ+gAQEEJ/BDOrpjdxP5q2nj6vpDeJWapZ+zWzUjpWXYd+07HK4D4NOeztztcsUp7/Laft5676Gn+vriQLpkHoDzGzqnR1w5f+4AYwd9XNrCzpSTpWNrMoz/1K+k7S9uA2lWNuiZmnfq9u0JO33ia9Jm/9fvZ86UL4Xl1lFkyL0L/tUFIvfdyTVJ5Ud/euuz8zs5KS+wf3tDz36bfx2R3LIknLDIGZ+11SX9N87VH1det3lQr/vbriLJgKoX9baej5zgz1WMl9f5epNPR86n7TWcj7JX9jju1nyvoiTfrao+qTXrNIk772pPqylYaeF/l7dRVZMJUg75w15lfE3uDXNN19s/axdXdvm9m3ZlYd/Do6D4vqV1LV3Z9m6+62Bfa7aJO+9qj6pNcs0qSvPam+bH3l7Ht1gr5ylgXzEGToT7gZ+2td/wSPJA0f8LpVTw/YvEnft685/482736lZD3S3U/Sx5V5rvMuot8luU9vpQmvWaQ8/1uOkrvv1QlylwXzwPLOkPSncjR80GtwtsOY+ktJvXSsNCH0Vt5v+vjYzE7N7HRZvd6337RelRQPDp7lpbdxr1mGPP9bzqvfNfxeXVkWTIsrcgEgIMz0ASAghD4ABITQB4CAEPoAEBBCHwACQugDQEAIfQAIyP8D9csiuGzpfC4AAAAASUVORK5CYII=\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "plt.scatter(pos[0:7, 0], pos[0:7, 1], color = 'red' , label = \"closed\")\n",
    "plt.scatter(pos[7:l, 0], pos[7:l, 1], color = 'blue', label = \"open\")\n",
    "plt.legend(loc = 3, borderaxespad = 1);"
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "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.7.4"
  },
  "toc": {
   "base_numbering": 1,
   "nav_menu": {},
   "number_sections": true,
   "sideBar": true,
   "skip_h1_title": false,
   "title_cell": "Table of Contents",
   "title_sidebar": "Contents",
   "toc_cell": false,
   "toc_position": {},
   "toc_section_display": true,
   "toc_window_display": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}