{
 "cells": [
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    "# Tutorial 1: getting started with iDDN\n",
    "iDDN allows to integrate multi-omics data to study the common and differential networks under user-specified, biologically plausible constraints. iDDN is a natural extension of DDN3.0, which focus mainly on single omics data. \n",
    "\n",
    "In this tutorial, we will introduce basic usage of iDDN. First, we will discuss the preparation of data and constraints. Here we will use a synthetic data with known ground truth. Then we will run iDDN and check the results. The visualization part will be discussed in the next tutorial.\n",
    "\n",
    "We start by importing some basic packages. "
   ],
   "id": "63ded71417fe47dc"
  },
  {
   "cell_type": "code",
   "id": "initial_id",
   "metadata": {
    "collapsed": true,
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:39.414119Z",
     "start_time": "2024-08-26T14:56:38.229321Z"
    }
   },
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from IPython.display import Image\n",
    "\n",
    "from iddn import iddn, tools  # Run iDDN algorithm and process the output\n",
    "from iddn_data import load_data  # Load example data and images\n",
    "\n",
    "%load_ext autoreload\n",
    "%autoreload 2"
   ],
   "outputs": [],
   "execution_count": 1
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    "## Prepare data and specify constraints\n",
    "\n",
    "First, we need to prepare the data iDDN need. For each condition, iDDN needs a NumPy array whose rows are samples and columns are features (like genes and mRNAs).\n",
    "Here we provide sample data `dat1` and `dat2`, which are simulated from two conditions.\n",
    "If your own data is saved as `csv` files, you may use `pandas` to load them.\n",
    "We assume the data is already pre-processed (normalization, missing value imputation). It is also often suggested to perform log transform on the data. But it is not necessary to standardize each column for each condition, as iDDN will do this.\n",
    "The data from two conditions can have different sample sizes. However, the number and order of the columns (features) must be the same.\n"
   ],
   "id": "1482cd5110ad3271"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:39.525822Z",
     "start_time": "2024-08-26T14:56:39.442935Z"
    }
   },
   "cell_type": "code",
   "source": [
    "example = load_data.load_example(\"example.npz\")\n",
    "dat1 = example[\"dat1\"]\n",
    "dat2 = example[\"dat2\"]"
   ],
   "id": "62d1204d1dc7b5d6",
   "outputs": [],
   "execution_count": 2
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "The data just loaded are simulated from two networks with three layers. The figure below gives a conceptual illustration of the network. Here we assume there is a network among nodes in the middle layers (mRNAs), and nodes in the top and bottom layers (TFs and miRNA, respectively) will regulate the middle layer. In the simulation, The mRNAs form a scale free network, and there are 50 mRNAs, 50 TFs, and 50 miRNAs. Each condition contains 500 samples.",
   "id": "4ece02d936f35071"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:40.520169Z",
     "start_time": "2024-08-26T14:56:39.782124Z"
    }
   },
   "cell_type": "code",
   "source": "Image(filename=load_data.get_image_path(\"three_layers.png\"))",
   "id": "b17284cb4a562905",
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAkUAAAG3CAIAAAA1rrYZAABMAElEQVR42u3dCVwU9f8/cFSUS+4bQuXyRCEB0TzyxCstr0TT0vKuX4d5VFZiZaWYHf9v4o2FKabilSJoah7JpQKhorCAIMcCu7CciyL+Pzi10XII7Mx8Zndfz0dfHisuOzMf/O5rPzPveX/aPX78WAcAAEDNtUOeAQCABkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCeAQCAJkCe0SQpqfwjPv22KL9QWl4gLS8srigsrrI2N7A2N7Kx6Gxt0bmXq93zPi6WZoa09xRYJpaUX7qWcTtdXCgpLSqpKJBWFpVUWZkZ2FgYWpkZWVua9HKxHertbGvZmfaeAqgN5BkFeUVl56JT/4hL++uu+Pn+Nv3dDa3N9azNOlmb6Vma6EpKawpLqgtLHhQWV19Pq/rjmrhvd9vnfd1GDnS3tzKmve+gksyc4nOxqRfjRaJs6TAvGy9Xg39+9R0tjDtKyx4WyR4UFD8kv/qEdPnFG2JXJ4thPq4jB7h3czSnve8AQoc841XNo9odB6N3HIp76XmHYZ5mz3uZt2vXrvkfIb+gPxKKLyaWHP0jd+F034UzBup2aE/7OKDVZGXynYdjDpxOmvq8/TAvs+c8WpRPfyYXX0woCf8jd+Y4zwXT/EyN9WkfB4BwIc/4czAyYcfB2EH9LBdMtHWyMWjtj2cXVO08Kb6aJFk4Y8CMsV60jwZaYefh6F2H4qeMcHx9gp2VaafW/jiZtO0+lX/kfM4b030WTBtI+2gABAp5xoeHNbUffvtbVWXZwkkOXu4mqrxUQmrpjhO5BobGX733QkddTNSErkBa/vEPEVbGj96Y6OjqqNJ1UFFO5a6TOUVlHb54e7yNBa6rAShDnnEuK6/ko+9+6+diuGp2V7Zec+O+e0nplV+++0IXezPaxwdNunYr5+PvI6YMs1002YGt19xxIi/8j/wv3hnv3duR9vEBCAvyjFvRSVlrvjs1b4LT3LF27L5yaGT+nlPZ69+dMLBfF9pHCY04fuHWuh/PfrGo1/iBFuy+ckR00cfb76x9c/Tk4b1pHyWAgCDPOHRfLHv1w7AP53Qb42vFxeufiSv6am/mz18FPGNrSvtY4T/IzGzR2sNr5rlMHWbPxeuHX8xbvyd9+7ppmKUBKCDPODR/zf7nPY3njefkHY3x0+m8CwmlIetn0z5W+FeBtPy1j8LenuY0fiAnn2MYEdGSHw5n/fRlAK6lATCQZ1xZ92OEzqOKtfOdOd9QSIaObue1y8bRPmL426K1B3x7GC6cxNo1s6ZsP54bf7dy+7qZtI8YQBCQZ5w4ceH2kcjY3R958LO5179MnjJ2wKThvWgfN+jsPBybnpn+5SJ3fjb30fZUl24uC6YNoH3cAPQhzzgx5e2QlbOcWnjPrOr+TC4O2p995If5tI9b28nK5OMW79z7qaeKpfktJ8qpnPNZ4ultC3CrNQDyjH1hETeir6d893Z3Pjf67g93B/bvGTD+WdpHr9W+CTn/6GHJqlm8Vpxu3J/VoaPZ+/NH0D56AMqQZyyrffx43MIdQcvcPVW7b7q1ElNLV25JPb1jYfunNdACjmTmFL/8/t5Tm3za0AFEFUWyBxNWxP/6zRz0eAQthzxj2bFzNy/HJwUt5XVyxlgZfHeIT78XR/ahPQZaaveR2ALx/Q9eYe2u+Zb7+pd7NrbPvD4FV9FAqyHPWLYi6Ngwj46Th9jyv+njl8UXkx9uWvki7THQUvM+2r9okjVvF03r+zO5ePuJwj1fzqI9Bppm48aNq1evbv45M2fODAsLIw8CAgIOHDjQ8Aljx449ffo07UPRCsgzNpHBHDT7f6c2+ViY8HrGiSEtrTvvdHXfW0/t2Q+sE0vKp7/706Ut1JoFD10Wfei717BeGrtIUO3Zs0fxx7S0NJFIZG5uPmDAv1PhkSNHrlq1ijwYN25cZGSkq6urm5tb/Rfx9vZev3497UPRCsgzNl26lrHnyOVdq6l1IXpjw615U4YM9eb8pjdQcigy6catO+sXuNDagTU705/t3WP62H60R0KTMdO1puZbTJ5t2LCBiTfgH/KMTUEhF6wNZfMmPENrB/acul9Yabpy/nDaI6F1Pg8+09OhesZIDnvBNO/gubyUXL1Plo6hPRKaDHkmcMgzNq3efGJ4v47j/axp7UBETOGFpIcblk+iPRJa5+314VOHGg9/luXWwy134YY0/FLZD2um0h4JTYY8EzjkGZte/zhs2WQbn17UugPHp8i2HCvY/UUA7ZHQOrNXhn78qlPvbsa0duBmRvn60Kx9QXNpj4QmQ54JHPKMTZPf3PXDOz262bd67Wm2ZOZXvf393eP/e532SGid0Qu27fu0r425Hq0dKCiunv3ZX2d3LqY9EpoMeSZwyDM2PffK/858O8DIoAOtHaioejTmvdg/f3mL9khoHa9p313fPaR9e2qVpTU1tT4LriQcfpf2SGiyluRZw/rGt99+e8KECbT3XSsgz9hE8uz8D356ndrT2gHkGS1kfnYg0NPStCOtHcD8jActybOG39+6devixfi98AF5xqbJb+768b0eTrb0zjfmVb39/Z3jP75BeyS0zuyVP3/yWtdeXand/nUrs+yLn7P2Bb1KeyQ0Gc43ChzyjE2vf7z/zZdsvXvQqwe5LdtyXLz7C/SJ4Nvb6w9PG2b8vJclrR14Ut9Y+sOaabRHQpMhzwQOecam1d+cGNGv47iBVOv1Ex9seH8y7ZHQOp8HR/VyfDh9hB2tHXhy/1mnT5b60x4JTYY8EzjkGZuCdp+3NirF/dRa6FBU0o3klPULXWntwJodomc9ek73R38QDiHPBA55xqYn/a6u7FpNbZ3oNzbcnjdlMPpd8a+uf+N7P136kV7/xjejD32L/o3cQp4JHPKMTehHrM3mfbRv0SQbev31C/Z8OZv2GGg45JnAIc9Y9mS9mE6Th9jwv+njlwsuJj/AejG07D4SU5Cf/cEcCpPjr/dm2Ng5vT7Fj/YYaLht27YtWbJEsUCMEibP9u/fHxCABj10IM9YhvU8tdaT9alDT23ypbE+ddyv38zF+tSg5ZBnLKt9/Hjcwh1By9w93U343G5iaunKLamndyxsj5ON9HwTcu7Rw5JVs3hdonrj/nsdOpq9P38k7aMHoAx5xr6wiBvR11O+e5vXKdq7P9wd2L9nwPhnaR+9VpOVycct3rn3U09XR0N+tijKqZzzWeLpbQtMjfVpHz0AZcgzTkx5O2TlLCfeSgP+TC4O2p995If5tI8bdHYejknPzPhykTs/m/toe6pLN+cF03DlDAB5xo0TF24fiYzd/ZEHP5t7/cvkKWMHTBpO7T4BqG/R2gM+3Q0WTXbkekM7TuTG3ancvm4m7SMGEATkGVfW/Xha51H52vmcV7utC8nQ6dB57ZvjaB8x/K1AWv7ah/vfnu40nstOMRHRBT8cuv/TV7NsLHDPGUAd5BmH5q/ZN9zL5LVx9txtYk9E3h+JpSHrceORsFy7lbNo7aE181ynDuPktx9+MW/9HtH2ddO9e3M+CwRQF8gzDt0Xy179YP+Hc53H+Fpx8fpn4oo2HK4MemfIs92pdYyERpVX1azbnXA59uanc2zGD2T5ZsSI6MKPt99Z++boycNxbwbAv5Bn3IpOylrz/al5E5zm+rPcqTY0Kn/PqezRw4dU1+gGzu9J+0DhX4lpsuDjGZ31dcd6d968+8zU520XTmJtFrX9+P0jFwu+eGcCZmYASpBnnMvKK/nou9/6uRiums3abUkb991LSq/88t0XutibzVkfv3Kmu6cbtUVqQIFMy0Kjsq8kS+b6O431tdV5ci3t4+9PWpnUvjHxGRWL+EU5lbtO3i8qbf/FOxNxzQygIeQZHx7W1H743cmqipKFk57xUu0+64TU0h0ncgyMTL96d2JH3bqFsMm7Z/CxjL1rfGgfpbYj07KgA6merqZLX3TubKBb/692Hvpz1+FrU4Y7vD7RoQ3dQ4pkD3afzD1yIfeNad4Lpj9H+0ABBAp5xp+DpxN2HIod1NdiwQt2TjatXsM6u6Bq52/5V/+SLpw+YMY4r/p/tWJL8hhfa2ZCAPwj07JNYWn5xXIyLRvs0fiSnrIy+c5DVw9E/jX1ebthXhYtvDfxz+TiiwnS8D/yZ47tu2D6INw0DdAM5Bmvah7V7vj16o7D8S8Nsx/mZf68l/lTe+GTX9AfCeRNrfjoxbyF03wXvjxQt0N7pefkS+VLNycGL/e0s8D7Ha9IkpH5cWhU9tShDv6+NkrTsoYyc4rPxd65GJcuyi4e9qy1l6uhtbmejXlHMmmzMO4oLXtYWEL+e1BYXJ0gqrx4o9DVyXyYr8vIAT3QmxHgqZBnFOQVlZ27eveP+LS/Uguff9a6v3vdm5q1WSdrMz1LE11JaU1hSTXzpnY9tfKPG4V93a2f93EbOai7vZVxU68ZFJZKvq4M4KktBejUXdCqCDqQ2llfd+mLzq6ORq36WbGk/NK1jNuivEJpWVFJRYG0sqhEbmWmb2NhaGVmZG1h3MvVfqi3M9YzA2g55BlNkpLKP+LTb6flFhZXFBRXFMoeySofm+pXW5sb2ZiTN7XO5E3teR8XS7On1xGQicKSzQnr5vVq7RsrtIGi7qOF0zIA4AHyTEDCL+aKi6vJh/22/XhknPjIpbyty73a9uPQQopy/BUBbjjBCyAcyDMBCQpL7edqokpZB5miNVOSACpSTMvIZw4MMoDQIM8EZMWWZJJGqtxJJsqpWBGcfOQLdFtnH3NfBIkx8jvCCUYAAUKeCcic9fGblnqoeAorMCTF1dGIvOfSPhrNQT4lkGkZmZy1oe4DAHiDPBMK8nY5d/011adW5HWmfBwTusYbl3ZUh7oPADWCPBMKpviblWoO8hZMXg1NHVWEug8A9YI8E4rIOHGSqJStG8jmrI9HzUKb1Z+WTR3mQHt3AKBFkGdCoWKxvhI0dWyzZtowAoCQIc+EQvVifSWBISnkBTG9aDkyLSMfAkS5FbjnAUAdIc+EgsTPGF9rFt9G0dSxVSLjxC1vwwgAAoQ8EwpWivWVkNkGmXOgqWPzUI4PoBmQZ4LAVrF+w5dFU8dmoBwfQJMgzwSBxWJ9JcxqJmjq2BDK8QE0DPJMENgt1ldCpmhThtpjtU8FZloWFVew5MVuGBYAjYE8EwR2i/WVMIUhaOrIYKZlrg5GKMcH0DDIM0FgvVi/4euT926O8lJdMHUf+cXylTPdcUERQPMgzwSB9WJ9JVre1JEcflRcQfilXNR9AGgw5JkgcFGsryT8Yu6fydJNyzxoHyvfFHUfKMcH0GzIM/o4KtZvSNuaOqINI4BWQZ7Rx12xvhIyUwnck6IlhSFMG0Y7c32U4wNoCeQZfZwW6yvRhqaOzLQsUSRDG0YArYI8o4/TYn0lzLlNDW7qiDaMAFoLeUYf18X6SjR1tU+0YQTQcsgz+rgu1leieU0d0YYRAHSQZ0LAQ7G+Ek1a7RNtGAGAgTyjjLdifSUa0NRRMS2b6++k1gcCAKxAnlHGW7G+EnVv6siU43u6mqINIwAwkGeUkenFmbhCKtUZQWGp5KvarfZJpmXBxzJEuRUoxweA+pBnlPFZrK9E7Wr3yQ4zy7mh7gMAGkKeUcZzsb6SyDgxmR2qRVNH5sQs2jACQFOQZ5TxXKzfkPCbOqIcHwBaAnlGGf/F+koE3tQR5fgA0ELIM5qYZcnOfDOY7m6QOaKro9Fcfyfa4/EfimmZwKePACAQyDOaaBXrKxFgYQhT9+HpakpSFicYAaAlkGc0USzWVyKcpo5owwgAbYM8o4lisb4Spqnjypnunm6mFPcBdR8A0GbIM5roFusrodvUEXUfAKAi5BlN1Iv1lazYkkz2h+d8rT8t0+yFRgGAU8gzmqgX6ythmjqGrvHm7Vwf2jACAFuQZzSNef8K9WJ9Jbw1dUQbRgBgF/KMGoEU6yvhZ7XPyDgx2jACALuQZ9QIp1hfCadNHVGODwAcQZ5RI5xi/YbIFI3104AoxwcATiHPqBFUsb4SMotaEZzMYlNHlOMDANeQZ9QIrVi/4e6x0tRRMS0jLyXM8AYAzYA8o0ZoxfpKmF7JoWu8VdlDZlrm6mCEcnwA4BryjBoBFusrUaWpY75UvvVYZn6xfOVMd9R9AAAPkGd0CLNYv6E2rPZJJnZRcQXhl3JR9wEAfEKe0XElWXLkYh5HNfHs7mermjoq6j5Qjg8APEOe0SHkYn0lgSEpgzzMn1rKgTaMAEAX8owOIRfrK2GaOja/2ifThtHOXB/l+ABAC/KMDoEX6ysJPpZBpl+NNnVkpmWJIhnaMAIAXcgzOgRerK+kqaaOaMMIAMKBPKND+MX6Skh0HbmUpyjIRBtGABAa5BkF6lKsr4RM0aYMtR/sYYk2jAAgQMgzChLTZCQShF+sr4TE8Hv/+8vBWh9tGAFAgJBnFKhRsb4CU/dxKia/TzeTrxf1ob07AADKkGcUqFGxPoMpx/d0NX11rNOcL66p2NQRAIALyDMK1KhYn0zLgo9liHIrFOX4ZHKZJCoV4DKkAKDlkGcUqEuxflPl+G1o6ggAwDXkGQXCL9ZnKjCbasPInH5seVNHAAAeIM/4JvBi/fptGJspxw8MSennaoI+jQAgHMgzvgm5WF/RHf+p5fgk9uauv9Z8U0cAAD4hz/gmzGJ9xbSs5RfGgo9liKXVKAwBAIFAnvFNgMX6zLTM09V0rr9Ty/t9NNXUEQCACuQZ3wRVrK9iG8bWrvYJAMAd5BnfBFKs38K6j6dimjoKaroJANoJecY3IRTrt7zu46mY1T6PfOFH94gAAJBnvKJerF9/WsZWtX1QWCr52uhqnwAAvEGe8Ypusb6iDePSF51ZXOeFZOSUj2PQ1BEA6EKe8YpWsX7DNozsiowTn4krFOZNdQCgJZBnvKJSrN9UG0Z2oakjANCFPOMcmZP9mSwlQeLqaJSYJvN0Mx3jY83PqTkVy/FbhRxa4J4UFIYAAC3IM84xl5cafp/MlriroWCrHL9VAkNSSGrO9XfiYVsAAEqQZ3yYsz5eLK1W+mbg/J4cnZ1jsRy/VdDUEQAoQp7xIfxibvCxjPrfmTrMgYuqEMW0jEySqNzjTLYuyqlAU0cA4B/yjA9KpxxtLfS2Lvdi/RwgMy1zdTBitxy/tUeKpo4AQAXyjCeKU44kaQLn9fR0M2XxxUmKbApLyy+Wr5zpTj1I0NQRAKhAnvFEccqR3TONJMmi4grCL+XyWffxVCu2JI/xtUZTRwDgE/KMJ8wpRzJ52rTUo22pExSWqnQiUVH3wUM5fqswTR1D13gLJF8BQBsgz7giKan8Iz498U5egaSisLhCIqvQMXJuV11gYaRjbW5kZ9W5b3e7531cLM0MW/JqV5IlgSEpgz0smVILLtowsgtNHQGAZ8gzluUVlUVdvns2Jv1Ohti1q72ZqZmBob6hvp6BgYGurm5NTU1VVVWlvJr8r0RWIrqX18PZdrSfi/+Q7vZWxk29JkmvFcHJopwK8phM78jXoAOpdub6PJfjtwoKQwCAZ8gz1tQ8qt12IHpXeFxP9652tjZODjbt2rVr/kfI4GfnFuSLC1JS770x1XfxzIG6Hdo3fBqZiv0cmaX4I3PPsvA7S6GpIwDwCXnGjrBTidsOxZIMc3N1Nunc6hlJaXlFmiiDZNvi6QMCJnjW/yvmWhSZ7ii+8+rYLurSg4NM0dQiegFAAyDPVPWwpnbFppPZ4oqe3V1trCxUeamCImnKXZGTrdGmFRM76v49UVuxJTlRJFN65tblXmpxHk+UU7EiOBlNHQGAB8gzlWTllSwPOmVqYty3T2+2XvOvm7dkpWWbV07oYm/GlIE0fA4JM4qLgrYKmjoCAD+QZ20XnZS16ptTnr3dXV26sfvKovTMxFup6/5v3K7IIsVd2EYGHchXO3N9Wws9W3M94dxt1jw0dQQAfiDP2ui+WPbK6rBn+/Xp5mTPxetnZufdSL77ydJxzg6mJMDUIrqagqaOAMAD5FkbzVoVZmFhSaZm3G0iLT1DIinavzGA9rGyYM76+JUz3dnt8gUAUB/yrC1Wb44UF8u9+nFeiZ6QlGxnYfD1e/60j1hVaOoIAFxDnrXaiQu3tx+6Nur5Qfxs7vc/ri6a7j1peC/ax62qwJCUQR7maOoIABxBnrXaxGV7evfq6Whnzc/mcvILb95OObVlHu3jVhVzIx0KQwCAI8iz1vnlZMKJP1IH+vbnc6PRcdcnPe/+ykT1KNBvBpo6AgB3kGetUPv48eg3dg708bKxVum+6dYqKJRGxyec3bWg/dMaaAkcmjoCAHeQZ61w7NzNfadvDeJ3csa4Gnd99rjeL47sQ3sMVBUZJz5yKU9dbgYHADWCPGuFt7483r6jiUs3Cq0u0jOzax+W/u+jybTHgAVo6ggAXECetRQZKL9Z/5s6cZS+vh7/W5fLq8NP/h6z/62n9uwXPjR1BAAuIM9a6tK1jG9Co4cPHkBrBy5ciX1/7sCh3s60R4IFQWGpdhb6aOoIACxCnrXU17suiHIf9OnpSmsHbqaIXB06ffDGcNojwYLyqpopH8eErvFG7T4AsAV51lLvbvitVtfEpYsjrR1Iz8ppX1P63eoXaI8EO8Iv5iaJStHUEQDYgjxrqdmrD3RxcrazpVbFkF8gycrK2LdhJu2RYM2c9fFLX3RGYQgAsAJ51lLjl4QM9O1vatKZ1g7IysqjY69HbJ1PeyRYk5gmCzqQiqaOAMAK5FlLDZz947QXRnfsSG3dlocPaw7/djZ635u0R4JNgSEp/VxNpg5zoL0jAKD2kGct5Tfrfy+/6N+hQwdaO0Dy7NCJMzH736I9EmzCap8AwBbkWUuNW7x7sJ+PsTG1Rk2y0vKrsddOb3ud9kiwLPhYBkk1NHUEABUhz1pq9uqwLl1d7fjt3FhfvliSlZ2+b4MmLO9ZH5o6AgArkGct9e7Xvz3WNXHuSrVe/2Hpdx9oSL1+fVeSJaFR2WjqCACqQJ611Fc7z4tyH3j0cqO1A8m309wc9TTjfuqGyBRtylB7rPYJAG2GPGupS9cyNv18dcQQal0Hz1+OWfHqIM3od9UQs9onmjoCQJshz1oK/Yi5FhSW2tlAd+mLmhnYAMA15FkrvLn+eLuOJm7OFLropmVkP35Y+uMaTVgvpilo6ggAqkCetcKxczf3RtwaPIDCep5XYq/PGa8J63k2LzJOfCaucNMyD9o7AgDqB3nWCrWPH498fcdzvs/a8Fu1X1Ao/TPuxrndC9tr7slGBTR1BIC2QZ61TuiJGycupj43wJvPjf4Ze23SMPe5k56lffR8SEyTBe5JQWEIALQW8qzVJizb06dXT0c7a342l5NfePN2yqkt82gfN3/Q1BEA2gB51monLtzedvDa6OGD+Nnc2QtXF8/wnjS8F+3j5o9SU0cyYyuX1+AMJAA0D3nWFqs3R+ZJq/p79uV6Q9cT/7K3MNiwfCztI+ZbaFS2KKdirr8TeXAlWfLq2C7kMe2dAgBBQ561UcDK/ZaW1u6uHN4sdSdVVFwsCQuaRftYKWBur9Z5MlcjX6cOc8B9aQDQPORZG90Xy2avCuvv2aebkz0Xr5+ZnZeUVjjCr8f04U69uhrTPlxekQlZ8LEMsbRa8Z3BHpaB83vS3i8AEDTkWdtFJ2Wt3hzh0cvN3aUbu6+cmp6ZfDttw/LxmUUdIqLzbcz1xg2wHeXNUwUKXUs2J4hyKpS+iTwDgKdCnqkkK6/kvY0nTU1N+/VhrV4j6eZtmUz27aqJXezNmO/8kVAUESvOLqga72c7foCtpWkn2sfNIRJmZHKWKJLV/6anqylusgaA5iHPVPWwpvb9TSez8sp793SzsVLpPuuCIumtFFEXe6NvVkzsqNte6W9T75dHxBRExOSP7G9Ngs3D2YT2oXMo/GIuSTXFH20t9Pau8aG9UwAgaMgzdvzyW8LOw3GO9lbd3V1NOrd6XcrS8oq7qaKcvKIF03xfeaG5ZcDKq2oiYsTkP1OjjuMH2vr72NA+dK6QI90UlnYlWUIedzbQxR3WANA85Blrah7VBodFhxyJ6+7WxcHO1snB5qm98MngZ+cW5OaL76ZlzZ/iuzRgoG6H9i3bms6VvyQRsQVk0jbBz5YEm40Zha7/PGAmasgzAHgq5BnL8orKIi7d/T0mLTWj0LmrnYWZuYGhvqG+ngGh36lK/qCqqqpSXk3+Jy0pzriX7+5sPcrPbfzQ7vZWbSliTM+riIgWk2Ab7GEx3s/Wy82U9gCwj0zUSKStDHCnvSMAIGjIM65ISiovxKcnpuSJJeUF0ori0sqK6g7Gxoam+o9sLIxsLTt79rQf7uNiaWao+rbkD2ojYvIjYgr0OrafMNCWBBvto1cJGbFL1zJuiQryCksLiyuKSipKSqvMTAyszIyszY3srU16u9oM9XYmY0h7TwFAQJBn/AkMSXF1NOK0z0X0LWlEjPhmRul4P7vxA20dLNVpIbHMnOLfY1J/j07PzJG6dLU3NjYhc1pDAz1DA319vU7yajK1lVdUVZMJbnl5mSgzt5ujxaiBLqP83Ls5mtPedwCgD3nGnznr4zct9eBhscoscRWZrp2KEfv0MJvgZ+fdw4z2oT+FrEy+/WDMr5FJvbt3tbGxaWGv55z8woKCglt37708tt+iGX6mxuoU3gDAOuQZT64kS0Kjsrcu91L9pVqo5tHjU9FiEmzt2rUb72dDZmy6HYS4fNrWX2NCwuP69e7m4uxMZmKt/fEquTw9IyPpVub8qb5LXkbNCID2Qp7xhCnSo9JUN/5OSUR0/vVU2bgBNhMG2jnZGNAejL8VSMtXf3u66kG7nu6uZqYq9fQqkZWlpIoMOj3e8N44GwtcV2sR8v/9y9czryTcExeVk9+FRFZZUlpFvm9mbGBpZkiG0daq82CvrkP6d3tqpS6AECDPeMLbycam5BbJT8WISbD1dTUdP8DGrzevS2w3dO1WzurNEd1du3R3c2XrNVPT0u+I7m1YPt67tyPdoxOy2sePT5y/dSZaFJt0z97W0t7WupOenqG+voFB3VfyhMoqeVVVdWW1/IFcnl9YlJsvGdCv65iBrpNG9NaGFdJBfSHP+JCYJgs+nsHnycamkPey00+ajDyoecx0z9Lr1NI73lh07Pytz7acHTG4/zOOLHdzzs7JvXDlxqfLRr84ojf/xyV8YRGJ2w/GWluZW1tZPWNvo6//9NsW5fLq+3kFhUVFhUXFi2YMCBjvSfsgABqHPONDaFQ2+SqoFbwS0mSnosVXb0rrUs3Pxtm+1T1N2ozMzBatPezn49HdpSsXr383/V5MfPL2ddMwS6vv2Llb2w/FGhsZujg721i3ZXZeUChNz8goq6hcPMNv8ggtWmAW1AXyjA9z1sevm9fL1ZG/zGihfKn8dGxBRLS4u1NnEmzPeXB+ErJAWj57VZinRw8nRwfutkJmaYnJd/ZtDMC1NMYn/+9Mwh1xzx5uLawdbUZOfmHKnTSvHraf/98Y2ocF8B/IM86JciqCDqQK4WRjMyJj65qMlFU+rLtxzc/GSF+Xow3N/fBXYxPz7m4uXB/R3VRRWVlJ6Fcvc70hgbsvlq38JkJPz8CrH5sLFNxISq6WV25aMeEZWw1sSQNqCnnGOQGebGxKUrrsdEzB+RuFT5qM2LmxPaHccTj299jsgT48RXt0fMKoAU4Lpw3gZ3MCxCzR18PNuXcP1opuFG6miO6kpW98f8LAfl1oHyhAHeQZ5+asj185091TfTorFskeMN2zutoajB9oN6yfJSsvKyuT+y/cOX70EBVL81uuRFYWcfZy1I4F2nmrdVZeyZwPwnyf9XBiu+hGITM770bSzdCvZirW6gOgCHnGLbU42diUs9cKSbCReHvSwt/OxFClk5Bf7riQnlvu0ZvXOoLkW7ddHDp/tHA4nxsVgoc1tbNWhdnb2rqyvXi6kruiDLFYHBY0q+GKfQA8Q55xKzQqu7yqZumLzrR3pO1uZZaRVIuMLRj/5CRkzy5tqbDIzCmesXzv1BdGtqEDiCqq5PLw384d3DxH23o8vrn+uPxhu759+LhpIenmbX3d2i0fT6Z90KDtkGfcUruTjU0pKX/IdM+ys9CfMNB2xLOtK5PbcTj2apK4nweFe8KSkm8N6merVVfRDkYm7T99a9hzvrxt8cLluNkTer88th/tQwethjzjkCinYu2e23vX+NDeETZdSCiKiBbnSuTjBtiMH2hnYdyxJT81a1XYM07dVC8Wb4Oc/ML72Zn7Nwbwv2kqah7Vjlm4a5Cvl40Vfy1gCoqkV2NvnNm5oOUL0gKwDnnGIQ042diUu9nlETHiU9H5o31sJvjZ9nE2aebJYkn5S2//NPOlcbT29sDR00d/eE1L1kv79ufLN+5K+nuyWZ3fEtcSkvv3sHzv1SG0BwC0F/KMQyu2JE8ZZj/Yg536QAEiaf3kJKTYwrjTOD+bMT42jT7tYGTSsYvpPl59ae1nfMJfLw5zmaEFZ8PyisomLtn90sQRJp35vnm/tLzi6MnzJ7e+3raV1gFUhzzjiiinYkVw8pEvtGIFk8t/SSKixel5Fczt2NZm/+kK+OmPZwtL2/Vw60Zr9+6kZVqbPP7szdG0x4lzPx29djYu25vSR4fricmjfJ557SVv2sMAWgp5xpXwi7ni4mqNPNnYFFFOecST7lnDPK3GDbBRVMEsDjxiZGbbxdGW1o5l5YgrSsTbAqfQHiHOzf3woLWtI62hJuNcKM4J/WoG7WEQlm3bti1ZsqT+d8zNzf39/d96660hQ5RPz1pYWBQXF5MnnD9/3tNTufVzQEDAgQMHLl261PAHicTERC+vuluDxo4de/r0adrHTQHyjCsaf7KxKZXVj07H1J2ENNTvQKZrJNimvbe3d89elhbUbrktkpbcTrl9+Ns5tMeGW5KSyvFLds+aOo7WcmXkzWR/+OmIra9bmhnSHgwB2bhx4+rVq0lEDRjwd5FtZGQk82D//v0kouo/WfG78/HxiYuLU3qpcePGkZ9tKs+WLVsWHBzMPBaJRC4unHeVExrkGSeYk42ha7w7G3DVCFH4rt6UklS7fa9Mkp8+enBfQwNqTToqq+RR569cCFlEe0i4dejMX4fP3vXjq51Yo2LiE6aN7j59DLVrpQLE5Fn9OZNEIvnkk0+Y7FEKHibPSPiRWVrDtGs+z8jcTufJHI688oYNG1atWkX70PmGPONE+MVcUW7FygB32jtCX2Z+5Utvbn9t5gsU1ziura09HJVgZedkpK9rqNeBTBwN6x60f/JY10if+U4HQ72/Hxvodah78ORvdTuozQqWn/7vTGFp+x7u3Sjuw53UTGuT2s/eQuv9fzXMM50nkWZlZaXTYIrG/N+EpBEzpUtNTbW0/PccTzN5FhYWNmvWrKVLl77wwgsTJ050dXVNS0ujfeh8Q55xQmtPNjZq+PztY0cOMWjB0pEcYeZnJ4Nfr5Q/qqyurZTXkAcV1Y+YB3WPyYO6Pz5SfK345690ddsZ6ekaNsg5o/9GIPmOgV57o3/CUq8jhduwlnx21NDEuoujHa1x1qm7hJZfWVq49dOXKO4D65gUIW+Vly9ffu+99+Lj48k3Z86c+eGHH3p6ejKTLRInZEZFUuSHH36YMGFC/R9vNM90/omurVu3Ll68WOmbZFtubm5k6qY0zWomz5hLaydPniRbZy7CNTWN02DIM/aVV9XMXX9Ny0821jf13dA+vftYmlNrkiKRlty8fSv8u7lt+Fn5AyYCn4RfXc79nYJ1D6rrHlRV1/4dfvVSsPbx47q54N+Z10ExKVSOQ/KcTu3/fvDkr1Qb5729e/a0sqTZ2atIUnwrJSX8O426VMmkyEcfffTll18yl8FiY2OZqo3jx49PnjyZPCZxReZDJIHI85lQUfx4o3mWnp5Owo88SEhIqF/3ocizU6dOkWmWzn9PSDaVZ8yrkf2RSqU6/1xII3O1LVu20B48XiHP2Bd+MTdJVBo4vyftHRGKRYHhxuZ2Tg406xvLS/K3B07lbYsPa2r/nghW/zP/YyaC/0bgI+W/rX4kf1D77wnP+hH43xQ06NTh3zjUr3tm+/Z1b4LD528bM3yIkSHNxQQqKuVnLly+ELJY9ZcSDiZFyAOSEJ9//rmlpSXJDx8fHxJjOk+udSlqEZkgUTrX1zDPSFatXbuWzPMaRo4izxTbJRNBMvmrvycN84wpoVS8GpOFinjTHsgz9gWGpIzxtcbJRoWP/98ZSXmHHm5dae3AnbRMy861Xwh+PWXy/0XFtE8xEVRM+/6eCDJnROU1T55Q92TyQL8TicD2eXniKSN7UbxOqfPkUuXPv55MOPwu7bFkE5MiStmzZs0aZrpWv7BeMeuq/77K5JnSa5If/OCDDxqWbNTPM0X9vSLAmsozX19fko7154XMKceGFSWaDXnGsvKqmiWbE7Yu98LJRoUDp/869kfagP7U2nPEXk968Xm3meM0tuiuqrou9qa9E/LkOiXmZyxrNEWavypW/8lK9frMVI/84C+//FK/1qP+jyvelpkJn6J2v9E9YWJPaTbGxG39uZ02QJ6xDCcbG3rSv/HnmS+NpbUDB45GHv3hVY3v31h3/axXbysLmos5aPD1MxXzTPFMiURCXpBMpxq9w0wpz8iT3d3dFTOtRveEiS6l6aNibqdVN6Ihz1iGk42NCli536mLM63++tlZGWFBs2iPAeee1DfaUOzDoqPR9Y1s5ZlOg5Rq+OMNT1e6urrGxMS88sorDfeEObVI0lFptsdMBJXqJzUb8oxNONnYlG0HY6/+le/Vtw//m0746+agvnaLZ2j++mdrfzwjlrXvSfv+Mzuzx2uXaVSrTNbzTOefSVXDKVrDPCMUtfvnzp1T2hNFGWRTtOpGNOQZm64kS87EFeJkY0OZOcXTl4dOe2EU/+tTH/7t90Ob52rD+tQC6Q8ydVT3Gf4adamSizx7ar1+/ddU1CuSH4mPj6//4sxtZ42W5itOOSptQoMhz9iEk43N+GLb+fTccp6XqE5KvuXi0PnjxSNoHz0f0L+RI1zkmU4TZZON5plOvXsG6r+4oslIU4nF1D1qz41oyDPWMCcbNy31sLOgWWAmWLIyuf+iXeNHDTYz5Wl9rBJZWcTvV6K2v2FqrC2/EfTX5wJHecZ0qCIPioqKFJe+msozxWSr/oszt501c0aReYL23IiGPGMNTjY+1fZDsWdjs5/z5emE2J9xCaMHOC2arvlXzhR+Onrt9/j7/C9OzdDU9c8azTMmKhoWxDPVGfUnTEyeNTpJYp5c/74xkmdNxY9i3RnFizM71kzFh6LwREt6XyHPWBMYktLP1WTqMAfaOyJocz781cjItFcPN643dPuuqKK8ZO9XL9M+Yl5hfWrQZsgzdjA9G4OXe+JkY/MKpOWzVoZ59Orh3JXD4M+4l5N8++7+oAAbCw2/56yhb3++fOOuhP8p2rWE5Gd7WC5/VfMnASBYyDN24GRjy127lbPw00MDfft2d+GkA9bd9HvRcX/t+Gy6d29H2sdKQc2j2jELdw3y9bKxsuBtowVF0quxN87sXKDbgcLCAgAM5Bk7AkNSXB2N5vo70d4R9XD0/M3Pfvx96EBP567PsPvKZGZ2KTrh0zdHvTSCwr1uAnEwMmn/6VvDnvPlbYsXLsfNntD75bHUWpoB6CDPWIHKxjaIv3l/9ebTbs5devVwZes1b91JFWXc3/j+eO2cmdX35vrj1TXtPXr34mFbScm39Ds+3vLxZNoHDdoOecaCK8mS0Kjsrctp3seqjgqk5Ss2nap62KFPD1cVi/hLZGU376QZdKzdtGKCFl4za+hhTe2slfvt7OzcXLpxuqG7aRniAnFY0KyOujjTCJQhz1gQFJZKZmY42dg2/2/f1dDj1/r06NrdzaUN3UOq5PK7aek379ybO9n7/2YPon00ApKVVzL3o4NkiubezZ6jTWRm511PTN77dUAXezPahwuAPGPDnPXxONmoClmZPPhA9KEzST1cuzrY27awbXFOfmFunviO6N70Mf2WzhyoPTdNt1BimmzTgTs5ufmudvp9erJ2UlfhZoooJTU9aMWEgf260D5WgDrIM1XhZCNbMnOKI67cvRCbfi9X6tzV3szE1MDAwMhAz4DM2vQ6yasfVFbJK6uqq6qqSmSyjKy8rg4Wwwe4jB/cXRt6M7YKSbLg4xminAry+I0JjkdPR3fqpN/fi82eitcSkh48kH+zYsIztjSXpwGoD3mmKhJm5CtONrJILCm/eC0j8U5+gaSsqKRCUlIpK6syNTawNDO0MjOysTT27GE3zNtZ49czawOSZEcu5ZHPWIrvbFrq4elm+vEPUQl3C3v1cFV9yR4yLb6dkubZ3Wb9O/60DxfgP5BnqsLJRhCCfKmcJFlUXEF5VU3974eu8Wb+cYb/fmvHoVhjI0M3F2cb67bcmlZQKE1LzyirqFw4fcDUUbz2lQZoCeSZSkQ5FUEHUnGyEagLjcr+OTJL6Zu2Fnp71/jU/87+iMTtv8ZYWJjZ29o8Y2+jr6/31FeWy6vv5xXkiQuk0pJFL/vNGq8VK4+AOkKeqQQnG0E4EtNkgXtS6s/PPF1NNy1TbnxV+/jx0d+Tz0anxydn2dlY2FpbGxjqG+rpGxjoGRrUzeQq5fKqqmrytbKiqkAiyRdLfDy6jh7o8tKoPu0prUQD0BLIs5YiA3X5emZ0UlZeYVmBtEIqq5SWVBhauus9llgad7CxMLK3Nh7Yr8uQ/t1orT4FkC+VL92cqIg0f1+blQHuTT2Z+Sd96XpmvqSisLiC/HsuKa0k3zczMbQwM7I2N7KzNBravxv+SYO6QJ49Bfkwe+L8rd9jM6ITMx1tLRzsrDt01DPUf/Jh9snNUpVVTz7MVssfPXyQKy7MyZcM7Nd1lJ/LpBG98WEWeBYalU0izdXBKPhYBvnj1GEOS190pr1TADxBnjUnLCJxx6FYaytzK0urVl1sKJIUFRYVL5w+IAAXG4AvV5IlJMa2LvfqbKDLXNmd6++E1dJBeyDPGnfiwu2tv8aYdDZy7tatzcVgGZmZpeUVS172mzScjzZ6oM2YJqJkNlY/wMg3SbbR3jUAniDPGvHp/84k3ino3p2dm3Xupor6udt8/n9jaB8WaLLAkBRbCz2cXQRthjz7j/ti2arNEXp6Bp592VwOMTEpWV5dFfT+eDRTAC6EX8yNii/YtNQDszHQZsizf0UnZa365pRXH3cX526sv3iaKCPptmjj++PR7A7YlS+VrwhOXjnT3dMNn5ZAqyHP/lbXjPzDA77P9nnGgatm5Nk5efE3kkPRjBzYU15VExiS8pyHxdRhDrT3BYAy5FmdusWiVu13dLB37tqV0w2JMjLz8sX7NwZgsShgRfjF3CRRaeD8nrR3BIA+5Fmdt76sW8y3Ty8+qhD/unVbr0MtFvMF1YlyKtbuuY32oQAM5JnOwciksMjbQwf5qP5SLXTxz7iAcb1fHtuP9qGDGiuvqlkRnIw7zAAUtD3Pah7V+i/c5efjaWvN35tCQZH0alzCmR1v6HbAWUdoI6YDCAr0ARS0Pc9+2Hvl+p0ir35sVue3xPXE5P49rN6dO5j2AIBaqt8KhPa+AAiFVudZXlHZhCW7p0wcadLZkOdNl5ZXHD15/uTW1+2tjGkPA6gZphUICvQBlGh1nv3y242zsdmefftQ2TqZoo0e4PTq5P60hwHUDFqBADRKq/PstTUHLa0duzjaUtl6Vo64qCDn5y9n0B4GUCdoBQLQFO3NM0lJ5YSlIQFTxtJa24mM/P7w0xFbX7c04/tsJ6gppkB/3bxero5GtPcFQHC0N8/CzyaHn7vr8yzN9VxirydOG+k+bUxf2oMBagCtQACap7159sXW33OkOj3cu1HchzupmXZmj9cuG017MEANhEZlk/kZWoEANEV782zZF0f1O1t3cbSjuA9ZOflVZYXBn7xEezBA6BLTZEEHUlGgD9AM7c2zGcv39eje3crSnOI+FEmKb9+5c/jbV2gPBggaWoEAtIT25tmI17ePfn6wkSHNxncVlfIzFy5fCFlMezBA0NAKBKAltDfPvKZ999rMF2gVNzJqa2t//vVkwuF3aQ8GCFdknPjIpTwU6AM8lfbmGZmf+Y8YbKCP+RkIF9bqBGg57c2z6e/90rNnLysLmm8TddfPUu4c/g7Xz6BxgSEpro5Gc/2daO8IgBrQ3jxb9vlRfWMbWs1BGFk5+ZWlhVs/RX0jNCL8Yu6fydLA+T1xphGgJbQ3z9ZtOZtX0q4n7fvPbExq1701hvZggOCgFQhAa2lvnoWfTT70+90B/Wn2B4mJT5g6qscMf75XqwGBYwr0/X1s0AoEoOW0N88kJZXjl+yeNXUc+jeC0AQfyyCRtjLAnfaOAKgT7c0z4tWPDlrZ0OyvXyi+H/rVy7SHAYQFa3UCtI1W51no8etn47Kf5X1xasa1hORRvs/Me8mb9jCAgJBp2dz111YEuKEVCEBraXWe5RWVTVyy+6WJI0w6833JHetTQ6OwVidAm2l1nhHf/nz5+p0iby++V2yJv5H0bHfr9+cNoT0AICBYqxNAFdqeZzWPascs3DXI18vGyoK3jRYUSa/GJpzZ+YZuh/a0BwCEgmkFggJ9gDbT9jwjfo1MCjt9e9hzPrxt8cLl2IBxvQPG07xVAAQFa3UCqA55VmfpF8eqH7bv59Gbh20lJt/S7/g4+OPJtA8aBARrdQKoDnlW52FNbcDK/Xa2tu6u3F6HT0lNLygoOLBpVkddnGmEv2GtTgBWIM/+lpVXMvfDA/09+3R9xp6jTWRm511PvLn365ld7M1oHy4IBdbqBGAL8uxf0UlZq76J6OHm3KenK+svfjNFdDcjb9wIbzOTuqv95F2souoR+Vr3WF73deVMdxQCaCGs1QnAFuTZf9wXy1ZsiujUSb8/qxX81xKSHjyQf7Niws5TeVeSJQ2f4OlqumkZujhqHbQCAWAR8qwRH/8QlXC3sFcPV0c7axVfKie/8HZKmmd3m/Xv+Ov8c3JJlFOh9LTA+T1xuknbYK1OAHYhzxoX/vutHYdijY0M3VycbazbcmtaQaE0LT2jrKJy4fQBU0f9WznJNDRizjQqnPlmMO0jBr5hrU4AdiHPmrM/InH7rzEWFmb2tjbP2Nvo6+s99Ufk8ur7eQV54gKptGTRy36zGrvJLDFNRj6YK/5I3tTIjG3qMIcpQ+3tLPRpHzTwAWt1ArAOefYUtY8fH/09+Wx0enxylp2Nha21tYGhvqGevoGBnqFBXfZUyuVVVdXka2VFVYFEki+W+Hh0HT3Q5aVRfdo3vRJNZJx4U1gaeUDezkLXeJPp2pn4wp8jswZ7WC55sRtSTbNhrU4ALiDPWooM1OXrmZeuZ+ZLKgqLK6QlFSWlle06dDQ16mhhZmRtbmRnaTS0f7ch/bu1cEG1oLDUqLiCV8d2qX/GiXxsD7+Ua2euP8bXeqwvnYVsgFNoBQLAEeRZ2zHFHYHzerZ5OkXe1xqdjV1JloRGZZMHU4baI9U0DPnN5kvlWKsTgHXIs7Yjb0w/R2YpTbBYlJgmi4ovIHM4XFrTGOR3Gnw8Ax30AbiAPGsjMjlbsjlBLK0mb0zByz25CxtRTkX4pVykmgZAKxAATiHP2ij8Yi7T2YHw97Xh4fQRc2nNzaHzIA9znIRUR0FhqeTTD1qBAHAEedYWiskZ80eup2j1MalGtohLa+oFa3UCcA151haKansFfqZoCleSJX8mS5mTkPi8L3xoBQLAA+RZqzXVsyp0jTfPF7dwaU1drNiSTJIMrUAAOIU8a7WGkzMGz1O0+sIv5oZGZXu6mo7xtUatgdAwrUDQbxqAa8izVktMk+UXy5nHW49lkg/dRgYdmIsidLMEl9YECK1AAHiDPFPJnPXxm5Z6COpE35VkyZm4wrTcchKuuLRGFwr0AfiEPFOJAPOMgUtrQoC1OgH4hDxTiWDzjEHmByTSyLsqmR+QVENxHZ+wVicAz5BnKhF4nimgzTHPmFXuAuf1xGcIAN4gz1SiLnnGIDOGIxfz8ovlU4c6oLk7pwJDUmwt9HCmEYBPyDOVLNmcsHKmu3qVruHSGtfQCgSACuSZStQxzxj5UvmRS3nkndff12auvxNSjS1MKxAU6APwD3mmEvXNMwVcWmMR1uoEoAh5phINyDNGZJw4NCob92KriIyhKKcicH5P2jsCoI2QZyrRmDxjYAVRVZDRCzqQigJ9AFqQZyrRsDxjoGCkDdAKBIA65JlKNDLPFLCCaMuhFQgAdcgzlWh2njHQ5vipmKuPONMIQBfyTCUrttSdYtKGHhBYQbQpWKsTQCCQZyrRnjxj4NJaQ4EhKWSCjrU6AahDnqlE2/JMIfxiLtoc6/yzVmfg/J440whAHfJMJVqbZwwtv7SGtToBBAV5phItzzOGos0xiTTtOe2GViAAQoM8UwnyTEHbLq2FRmXnS+UrA9xp7wgA/A15phLkmRLyFn8mvvDnyKzBHpZkZDT1RBzW6gQQIOSZSpBnTdHgNsflVTVLNicsfdEZrUAABAV5ppLAkBTyfo33taaQeUxoVDZ5oEkFI1irE0CYkGcqQZ61hCa1OcZanQCChTxTCfKs5RQriKpvqqEVCICQIc9UgjxrAzVtc4wCfQCBQ56pBHnWZmq3giiJ4SRRKdbqBBAs5JlKkGcqqn9pba6/k2AvSjGtQDYt9VDH06QAWgJ5phLkGSsEfi821uoEUAvIM5UEhaX2czVRi9NlaiH8Ym5oVLanq6mgPiVgrU4AtYA8UwnyjAuCanOMViAA6gJ5phLkGXdIkJyJK0wUyfx9bWjNjZhWICjQB1ALyDOVIM+4RvfSGtbqBFAjyDOVIM94w/8KomgFAqBekGcqQZ7xjLdLa1irE0DtIM9UgjyjQrGC6NShDlx060ArEAB1hDxTCfKMIu4urYVGZZMXRysQAPWCPFNJ8LEMW3M9fIqnSLGCqL+vzVx/J9VTLTFNFnw8A5fNANQO8kwlyDPhYGUFUbQCAVBfyDOVIM+ERsUVRNEKBEB9Ic9UgjwTpratIBoZJz5yKQ9nGgHUFPJMJcgzIWtVwQjW6gRQd8gzlSDP1EJLVhBFKxAAdYc8UwnyTI00cy82+as/k6WB83viTCOA+kKeqQR5pnauJEtIdDEnIZm6D7QCAdAMyDOVMKV0OEmldhSX1iY9Z5eUXjrBzxYfSgDUHfJMJcgzdff+luTb98oG9DTnrc0xAHAEeaYS5JlaU6zVSSZqwllBFADaBnmmEuSZ+iqvqpm7/tqKADdFKxBmBdG03HISafidAqgd5JlKkGfqKzAkxdZCr2ErELoriAJAmyHPVII8U1NPXauTzN5IpDEriJLfL0ofAYQPeaYS5Jk6YlqBtLBAn5U2xwDAA+SZSpBnaqdta3UqVhAlv2ukGoAwIc9UQj68i4ur0Y5djaiyVicurQEIGfJMJcgz9ZKYJgs6kLp1uZcqfa3ypfIjl/LIrx6pBiAoyDOVIM/UCOtrdeLSGoCgIM9UgjxTIxyt1RkZJw6Nysa92ADUIc9UgjxTF0zqqHimsRn1VxAlU0D06QfgH/JMJcgztcDbWp0oGAGgCHmmEuSZWuB/rU7FCqJjfK3ZulwHAM1DnqkkMk6cJCpdGeBOe0egSRTX6mxmBVEAYB3yTCXIM4ETwlqdzAqi5D9/XxtM5QG4gzxTCfJMyNrWCoQjuLQGwDXkmUqQZ0IWGpWdL5UL7bcTfjGXaXOMFUQB2IU8UwnyTLBYaQXCHVxaA2Ad8kwlyDNhYr0VCEcUbY6xgiiA6pBnKkGeCVNTa3UKEy6tAbACeaYS8vn6TFxh25q1A0eeulanMOVL5WfiC3+OzCJzyiUvdkOqAbQW8kwlyDOh4a0VCHfQ5higbZBnrSCWlF+6lpGcllsoLS8qqSgsriytbG9iamTd+ZGVmZGNpXEfV/uh3s62lp1p76k6aTiq0hK5hZm+tblhG0Z1xZZkkmQacC2KfFRiVotlq2CE3XEGECDk2dNl5hSfjblzIU6Ucb/Y08O0a5cOpia6df+ZduhspFteUVNaWlMieyQrrbmX9SgxWeb8jPlwX9fRfj26OZrT3nfh4mJUycwmSVSqSdPl+m2O23ZpDf96QXsgz5ojK5NvPxR9MDJpsJ9pnz6GvXsYtuSnbt2pvHmz8kqMbMbYfoumDzQ1xoWQ/+BoVJlWIJuWemjelae2rSCKf72gbZBnTdp28ErIkRtDB5mPHNHZ1LjVlQWysppz58svXS2eP+XZxTMG0z4aoeBoVNWlQF9FijbHgzzMmz8JiX+9oIWQZ40okJZ/+P1JXf0q/5Gm9nadVHmpvPwHUedkNXKDr96ZaGOh1VcmOB1VjtbqFKbmVxDFv17QWsgzZddu5Xz43cnn/Iz9R5uw9ZpRZ8v+jCn96t2J3r0daR8fHZyOqrxWn+SZYFuBcKT+pTVFkONfL2gz5Nl/HDt/c92PZ197xda3vzG7rxx3vfSnXwoC3xozeXhv2kfJN65HdeSIwXPHqnGBvirq34vdWbfs25Dz+NcLWgt59i/y2Xbhp4dmvWw12M+Mi9e/ElOy/9eiHZ9N16rPuRhVfnz/6+1Dv9+dPKoc4wxaC3n2twJp+ZwP9k0ab+rdn8PrBNeuV5yIKNn79WwtuRqBUeUHxhlAB3mm8Pqn+7t0a+c/iuUTNQ1FnS3Nuqez+7NZtI+YDxhVfmCcAXSQZ4wdh2Pi795+bbY1P5v7aV+hT/deC6f50T5ubmFU+YFxBmAgz+puOx27eOfKdxxVLG5uubz8B0Hf50RuW6DBN6tiVPmBcQZQQJ7pfL37rKRCPGUyJ1fRm3LkeImlke0Hr4+mffRcwajyA+MMoKDteZaZU/zy+6HrPu7ahh4KqpCV1az94t6v38zVyC55GFV+YJwB6tP2PNsZHp2SLZr6Eq8fbxnhR0t6OrkumDqQ9hiwD6PKD4wzQH3anmdzPvxl+IhOLWzVyq5bdyovnH+w96tXaI8B+zCq/MA4A9Sn1XkmlpRPfXdP0BfUmv6t/Dgj/Lt5GrbiFEaVHxhnACVanWcHI5POJyTOCbCgtQN7w6QjvDxnjO1HeyTYhFHlB8YZQIlW51lgcKSBmWToc9T6/l36U1ZVYhm4dCztkWATRpUfGGcAJVqdZ8s+P9Svf20/DyNaO5CUXJF0vf2WT6bTHgk2YVT5gXEGUKLVeTZzxc9TXjLq4kTtttB7WdVHj5Uf2PQq7ZFgE0aVHxhnNSWRSHbt2uXh4TFhwgTa+6JptDrPRr6xdcU79mam1BbNKpHVbPo+79yuJbRHgk0YVX5gnNXUxo0bV69eTR7Uf++1sLAoLi6u/zQfH58ZM2a88cYblpb/WXI9ICDgwIED5MH+/fvJY6UX37Zt25IlSzZs2LBq1apGt6748UuXLg0ZMoT2YLBMq/PMa9p3/9vk3r49tR2oqdF5e1VqwuF3aY8EmzCq/MA4q6mwsLBZs2aRuIqLi1N8s127djpPMoxJr7S0NJFIxHzn9OnT9SNt3LhxkZGR5IG5uXl8fLyLi0v9F2fCsqk8S09Pd3V1ZR4vXbp0y5YttAeDZVqdZ+QT7qr3HEyMO9DaAY38hItR5QfGWZMweVZ/znT58uXJkyeTSZtS8DB5RsKs4V/pPC3PmNnbRx999OWXX5JXSE1NVZr8qTutzrO6KxBTOnd5Ro/WDmRly48cLT+w6TXaI8EmjCo/MM6apGGeEWvWrCHBQ2ZUZLqm+CaTZySxmJOWSj/SfJ65ubmJnpg5cyaZ2zV6xlKtaXWeLfv8kKd3bd8+NCvEEq+3D9asCjGMKj8wzrQwmUGCxMHBgcx1mMtRPj4+69atY0o8yBO2b99OYoPMgT744AOlaCETr6FDh5JECQsLU3yz0TxjNkReRCqVKr7J5Bl55r59+4KDg8eOHXv69GmlH2k0zxITE728vJjznMxETelnNYBW51nglkgDC8nQQTTv4JGXWK1d6k97JNiEUeUHxpkWJjOWLl1KAqm4uJikguJy18mTJ9euXUumPiQ2yB/JA/KVZN769esVP87kmVKWNJpny5YtI4mllHyKPCNpSrZCdqD+NKuZPGNebevWrYsXL1ZcSCO7rXQFTq1pdZ6hwwIXMKr8wDjToihQJHFCJmckDyQSySuvvMKUaejUqzxkSj90/hsbLckzkjeHDh1iJmfnz5/39PRUPFORZ+SZzJ6QZIqJiWGuhDWTZ0wJpWJPfH19Sdwy8UZ7RFmj1XmGDnhcwKjyA+NMC5MZSpWHp06dmjhxok6DMnrmkpVSoUdTeaaEPIckU/0w0/lvnpEc9fPzI6+vCLCm8ozZvfpFlcwpR6WLc+pOq/NMBx3KuYFR5QfGmYpGM4NJKZ3/3lWmU698Q/HkZvKMqddnzl6SmdmWLVsa1mvUzzOdf6aAitr9pvKMue2s/myMZKGVlRV5kJCQoBSZ6kvb8wwrSHEBo8oPjDMV3OWZUkrpPLkgp9RGRCnPFN9havcb3TdFdCldLWNCTpNuRNP2PMvMKZ6xPPSzTyis8Pvp5/cObtbMFX4xqvzAOFPBQ54pttLwfGDDPFNsmsy0yF813Dfm1CKZww0YMKD+SykmgvXrJ9WatucZ8fXus5IK8ZTJvH7IPXK8xNLQ9oM3RtM+eq5gVPmBceYfP3mmmFQpTdEa5pnOP7WL5DVHjhzZcN+Y0o9mjkhjbkRDnpEPm/Kxi3eufMfR3q4TP1vMy38Q9H1O5LYFpsbUmslyDaPKD4wz//jJM52n1evXf2Z6ejpTu09eVmlzzG1nOk2U5jOnHJU2ob6QZ3V2HI6Jv3v7tdnW/Gzup32FPt17LZzmR/u4uYVR5QfGmWe85ZniNetHUaN5plPvLgKi/uaYJiNNJZaiLFMzbkRDnv3t9U/3d+mi4z/GhOsNRZ0ty7r3ePdns2gfMR8wqvzAOPOJtzzT+afcv/4d2U3lmaJ2X+e/ecbcdtbMGUXmCZpxIxry7G8F0vI5H+ybNN7Muz+HDYSuXS8/ESHb+/VsGwutuGsHo8oPjDOfGs0z5rRew9qKhoXyzfS7ajhJYko56t83xuRZo0X2ihbGis0x22q+4oM5HM3ofYU8+9e1WzkLPz0062WrwX6cXF2/ElOy/9eiHZ9N9+7tSPtY+YNR5QfGGQB59h9Hz//12Y/n5r9i592f5U+g165XhPyS9+mbI18a0Zf2UfINo8oPjDNoOeSZMvI598PvTj7nZ+I/2pit14w6U/pnbNlX707U2s+2GFV+YJxBmyHPGlEgLf/w+5O6+lX+I01VLIPOy38QdU5WIzf46p2JWn7VAaPKD4wzaC3kWZO2HbwScuTG0EHmI0d0bkP/BVlZzbnz5ZeuFs+f8uziGYNpH41QYFT5gXEGLYQ8a46sTL79UPTByKTBfqZ9+hi2sPHrrTuVN29WXomRzRjbb9H0gVp722lTMKr8wDiDtkGePV1mTvHZmDsX4kQZ94s9PUy7dulgaqJrZtrBxES3s5FueUWNTPZIVlpD/ruX9SgxWeb8jPlwX9fRfj20s7tdC2FU+YFxBu2BPGsFsaT80rWM5LTcQml5UUlFYXGltERuYaZvbW5oZWZkY2ncx9V+qLezFq4IpQqMKj8wzqDxkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJkGcAAKAJ/j8OyVpj420VtwAAAABJRU5ErkJggg==",
      "text/plain": [
       "<IPython.core.display.Image object>"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "execution_count": 3
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    "In many problems, there are some biological knowledge about the (non)existence of edges between nodes and iDDN is able to incorporate this information during the optimization process. Here we use a pre-specified constraint matrix `dep_mat` (here `dep` means dependency). Since we have 150 nodes, the shape of the constraints matrix is 150 by 150. If element `[i,j]` in the constraint matrix is 1, we allow edges from node `i` to node `j`; if the value of that element is 0, we do not allow such edge. \n",
    "More details about the constraints will be given in the next tutorial."
   ],
   "id": "2ef2b387b2540a94"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:41.428139Z",
     "start_time": "2024-08-26T14:56:41.081075Z"
    }
   },
   "cell_type": "code",
   "source": [
    "dep_mat = example[\"dep_mat\"]\n",
    "plt.imshow(dep_mat)"
   ],
   "id": "f440e06194702720",
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<matplotlib.image.AxesImage at 0x15624d507c0>"
      ]
     },
     "execution_count": 4,
     "metadata": {},
     "output_type": "execute_result"
    },
    {
     "data": {
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ],
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAakAAAGhCAYAAADbf0s2AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAi8UlEQVR4nO3df3DU9Z3H8deGhSUFdkNis5vVRGKHHqgUkWgaYXq27DQig3DSWpiUcpYxZ5uoEAchrcFrq0Y4f1CQkurc+GMKtXVGUJkRJw02qWMIkEhbEQOOKaTQTdrG7JJgQsh+7o8OO7cSJNEN+9n4fMx8Z8z3+9kv788Iec6SL4nDGGMEAICFUhI9AAAA50OkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWSlikNm/erEmTJmns2LHKz8/X3r17EzUKAMBSCYnUb37zG5WVlemBBx5QU1OTpk+frsLCQrW3tydiHACApRyJ+Aaz+fn5uu666/Tkk09KkiKRiLKzs3XXXXdpzZo1F3x9JBLRiRMnNGHCBDkcjuEeFwAQZ8YYnTx5Un6/Xykp53+/5LyIM0mSTp8+rcbGRpWXl0fPpaSkKBAIqL6+fsDX9Pb2qre3N/rx8ePHdeWVVw77rACA4dXa2qrLLrvsvNcveqT+8Y9/qL+/X16vN+a81+vVe++9N+BrKisr9ZOf/OSc80ebJsk9nmc/ACDZhLsiuvzav2jChAmfuO6iR+rTKC8vV1lZWfTjcDis7OxsucenyD2BSAFAsrrQl2wueqQuueQSjRo1Sm1tbTHn29ra5PP5BnyNy+WSy+W6GOMBACxy0d+GjBkzRjNnzlRNTU30XCQSUU1NjQoKCi72OAAAiyXkr/vKysq0bNky5eXl6frrr9eGDRvU3d2t22+/PRHjAAAslZBIfec739Hf//53rV27VsFgUNdcc4127dp1zsMUAIDPt4T8O6nPKhwOy+Px6MPDV/DgBAAkofDJiCZ++QOFQiG53e7zruMzPADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1op7pCorK3XddddpwoQJyszM1MKFC9Xc3ByzpqenRyUlJcrIyND48eO1aNEitbW1xXsUAECSi3ukamtrVVJSoj179qi6ulp9fX365je/qe7u7uialStX6tVXX9WLL76o2tpanThxQrfeemu8RwEAJDmHMcYM5y/w97//XZmZmaqtrdXXvvY1hUIhffGLX9S2bdv0rW99S5L03nvvaerUqaqvr9dXv/rVC94zHA7L4/How8NXyD2Bv7EEgGQTPhnRxC9/oFAoJLfbfd51w/4ZPhQKSZLS09MlSY2Njerr61MgEIiumTJlinJyclRfXz/gPXp7exUOh2MOAMDIN6yRikQiWrFihWbNmqWrr75akhQMBjVmzBilpaXFrPV6vQoGgwPep7KyUh6PJ3pkZ2cP59gAAEsMa6RKSkr0zjvv6IUXXvhM9ykvL1coFIoera2tcZoQAGAz53DduLS0VDt37lRdXZ0uu+yy6Hmfz6fTp0+rs7Mz5t1UW1ubfD7fgPdyuVxyuVzDNSoAwFJxfydljFFpaam2b9+u3bt3Kzc3N+b6zJkzNXr0aNXU1ETPNTc369ixYyooKIj3OACAJBb3d1IlJSXatm2bXn75ZU2YMCH6dSaPx6PU1FR5PB4tX75cZWVlSk9Pl9vt1l133aWCgoJBPdkHAPj8iHuktmzZIkm68cYbY84/88wz+s///E9J0hNPPKGUlBQtWrRIvb29Kiws1C9+8Yt4jwIASHLD/u+khgP/TgoAkps1/04KAIBPi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGCtYY/UI488IofDoRUrVkTP9fT0qKSkRBkZGRo/frwWLVqktra24R4FAJBkhjVS+/bt0y9/+Ut95StfiTm/cuVKvfrqq3rxxRdVW1urEydO6NZbbx3OUQAASWjYItXV1aWioiI9/fTTmjhxYvR8KBTS//7v/+rxxx/XN77xDc2cOVPPPPOM3nrrLe3Zs2fAe/X29iocDsccAICRb9giVVJSonnz5ikQCMScb2xsVF9fX8z5KVOmKCcnR/X19QPeq7KyUh6PJ3pkZ2cP19gAAIsMS6ReeOEFNTU1qbKy8pxrwWBQY8aMUVpaWsx5r9erYDA44P3Ky8sVCoWiR2tr63CMDQCwjDPeN2xtbdU999yj6upqjR07Ni73dLlccrlccbkXACB5xP2dVGNjo9rb23XttdfK6XTK6XSqtrZWGzdulNPplNfr1enTp9XZ2Rnzura2Nvl8vniPAwBIYnF/JzVnzhz9+c9/jjl3++23a8qUKVq9erWys7M1evRo1dTUaNGiRZKk5uZmHTt2TAUFBfEeBwCQxOIeqQkTJujqq6+OOTdu3DhlZGREzy9fvlxlZWVKT0+X2+3WXXfdpYKCAn31q1+N9zgAgCQW90gNxhNPPKGUlBQtWrRIvb29Kiws1C9+8YtEjAIAsJjDGGMSPcRQhcNheTwefXj4Crkn8J2dACDZhE9GNPHLHygUCsntdp93HZ/hAQDWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwVkJ+nhQwWP0movb+U+qIjEr0KADiqOt0ZFDriBSs1t5/SnP23qlRDW4p6X7yGYDz6e/tkfSjC64jUrBaR2SUnPVuZT32VqJHARBHZ0yfmgexjq9JAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFrDEqnjx4/ru9/9rjIyMpSamqpp06Zp//790evGGK1du1ZZWVlKTU1VIBDQkSNHhmMUAEASi3ukPvzwQ82aNUujR4/Wa6+9pnfffVePPfaYJk6cGF2zfv16bdy4UVVVVWpoaNC4ceNUWFionp6eeI8DAEhiznjfcN26dcrOztYzzzwTPZebmxv9b2OMNmzYoPvvv18LFiyQJD3//PPyer3asWOHFi9eHO+RAABJKu7vpF555RXl5eXp29/+tjIzMzVjxgw9/fTT0estLS0KBoMKBALRcx6PR/n5+aqvrx/wnr29vQqHwzEHAGDki3ukPvjgA23ZskWTJ0/W66+/rh/84Ae6++679dxzz0mSgsGgJMnr9ca8zuv1Rq99XGVlpTweT/TIzs6O99gAAAvFPVKRSETXXnutHn74Yc2YMUPFxcW64447VFVV9anvWV5erlAoFD1aW1vjODEAwFZxj1RWVpauvPLKmHNTp07VsWPHJEk+n0+S1NbWFrOmra0teu3jXC6X3G53zAEAGPniHqlZs2apubk55tzhw4d1+eWXS/rXQxQ+n081NTXR6+FwWA0NDSooKIj3OACAJBb3p/tWrlypG264QQ8//LBuu+027d27V0899ZSeeuopSZLD4dCKFSv04IMPavLkycrNzVVFRYX8fr8WLlwY73EAAEks7pG67rrrtH37dpWXl+unP/2pcnNztWHDBhUVFUXX3Hffferu7lZxcbE6Ozs1e/Zs7dq1S2PHjo33OACAJOYwxphEDzFU4XBYHo9HHx6+Qu4JfGenkezg6Y+0ZNO9ynrsrUSPAiCOzpg+/V4vKxQKfeJzBnyGBwBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWnGPVH9/vyoqKpSbm6vU1FR96Utf0s9+9jMZY6JrjDFau3atsrKylJqaqkAgoCNHjsR7FABAkot7pNatW6ctW7boySef1KFDh7Ru3TqtX79emzZtiq5Zv369Nm7cqKqqKjU0NGjcuHEqLCxUT09PvMcBACQxZ7xv+NZbb2nBggWaN2+eJGnSpEn69a9/rb1790r617uoDRs26P7779eCBQskSc8//7y8Xq927NihxYsXx3skAECSivs7qRtuuEE1NTU6fPiwJOmPf/yj3nzzTc2dO1eS1NLSomAwqEAgEH2Nx+NRfn6+6uvrB7xnb2+vwuFwzAEAGPni/k5qzZo1CofDmjJlikaNGqX+/n499NBDKioqkiQFg0FJktfrjXmd1+uNXvu4yspK/eQnP4n3qAAAy8X9ndRvf/tbbd26Vdu2bVNTU5Oee+45Pfroo3ruuec+9T3Ly8sVCoWiR2traxwnBgDYKu7vpFatWqU1a9ZEv7Y0bdo0HT16VJWVlVq2bJl8Pp8kqa2tTVlZWdHXtbW16Zprrhnwni6XSy6XK96jAgAsF/d3UqdOnVJKSuxtR40apUgkIknKzc2Vz+dTTU1N9Ho4HFZDQ4MKCgriPQ4AIInF/Z3U/Pnz9dBDDyknJ0dXXXWV3n77bT3++OP6/ve/L0lyOBxasWKFHnzwQU2ePFm5ubmqqKiQ3+/XwoUL4z0OACCJxT1SmzZtUkVFhX74wx+qvb1dfr9f//Vf/6W1a9dG19x3333q7u5WcXGxOjs7NXv2bO3atUtjx46N9zgAgCTmMP//W0EkiXA4LI/How8PXyH3BL6z00h28PRHWrLpXmU99laiRwEQR2dMn36vlxUKheR2u8+7js/wAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYa8iRqqur0/z58+X3++VwOLRjx46Y68YYrV27VllZWUpNTVUgENCRI0di1nR0dKioqEhut1tpaWlavny5urq6PtNGAAAjz5Aj1d3drenTp2vz5s0DXl+/fr02btyoqqoqNTQ0aNy4cSosLFRPT090TVFRkQ4ePKjq6mrt3LlTdXV1Ki4u/vS7AACMSM6hvmDu3LmaO3fugNeMMdqwYYPuv/9+LViwQJL0/PPPy+v1aseOHVq8eLEOHTqkXbt2ad++fcrLy5Mkbdq0STfffLMeffRR+f3+z7AdAMBIEtevSbW0tCgYDCoQCETPeTwe5efnq76+XpJUX1+vtLS0aKAkKRAIKCUlRQ0NDQPet7e3V+FwOOYAAIx8cY1UMBiUJHm93pjzXq83ei0YDCozMzPmutPpVHp6enTNx1VWVsrj8USP7OzseI4NALBUUjzdV15erlAoFD1aW1sTPRIA4CKIa6R8Pp8kqa2tLeZ8W1tb9JrP51N7e3vM9TNnzqijoyO65uNcLpfcbnfMAQAY+eIaqdzcXPl8PtXU1ETPhcNhNTQ0qKCgQJJUUFCgzs5ONTY2Rtfs3r1bkUhE+fn58RwHAJDkhvx0X1dXl95///3oxy0tLTpw4IDS09OVk5OjFStW6MEHH9TkyZOVm5uriooK+f1+LVy4UJI0depU3XTTTbrjjjtUVVWlvr4+lZaWavHixTzZBwCIMeRI7d+/X1//+tejH5eVlUmSli1bpmeffVb33Xefuru7VVxcrM7OTs2ePVu7du3S2LFjo6/ZunWrSktLNWfOHKWkpGjRokXauHFjHLYDABhJHMYYk+ghhiocDsvj8ejDw1fIPSEpnv3Ap3Tw9EdasuleZT32VqJHARBHZ0yffq+XFQqFPvE5Az7DAwCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALDWkCNVV1en+fPny+/3y+FwaMeOHdFrfX19Wr16taZNm6Zx48bJ7/fre9/7nk6cOBFzj46ODhUVFcntdistLU3Lly9XV1fXZ94MAGBkGXKkuru7NX36dG3evPmca6dOnVJTU5MqKirU1NSkl156Sc3Nzbrlllti1hUVFengwYOqrq7Wzp07VVdXp+Li4k+/CwDAiOQc6gvmzp2ruXPnDnjN4/Gouro65tyTTz6p66+/XseOHVNOTo4OHTqkXbt2ad++fcrLy5Mkbdq0STfffLMeffRR+f3+T7ENAMBINOxfkwqFQnI4HEpLS5Mk1dfXKy0tLRooSQoEAkpJSVFDQ8OA9+jt7VU4HI45AAAj37BGqqenR6tXr9aSJUvkdrslScFgUJmZmTHrnE6n0tPTFQwGB7xPZWWlPB5P9MjOzh7OsQEAlhi2SPX19em2226TMUZbtmz5TPcqLy9XKBSKHq2trXGaEgBgsyF/TWowzgbq6NGj2r17d/RdlCT5fD61t7fHrD9z5ow6Ojrk8/kGvJ/L5ZLL5RqOUQEAFov7O6mzgTpy5Ih+97vfKSMjI+Z6QUGBOjs71djYGD23e/duRSIR5efnx3scAEASG/I7qa6uLr3//vvRj1taWnTgwAGlp6crKytL3/rWt9TU1KSdO3eqv78/+nWm9PR0jRkzRlOnTtVNN92kO+64Q1VVVerr61NpaakWL17Mk30AgBhDjtT+/fv19a9/PfpxWVmZJGnZsmX67//+b73yyiuSpGuuuSbmdW+88YZuvPFGSdLWrVtVWlqqOXPmKCUlRYsWLdLGjRs/5RYAACPVkCN14403yhhz3uufdO2s9PR0bdu2bai/NADgc4bv3QcAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWGHKm6ujrNnz9ffr9fDodDO3bsOO/aO++8Uw6HQxs2bIg539HRoaKiIrndbqWlpWn58uXq6uoa6igAgBFuyJHq7u7W9OnTtXnz5k9ct337du3Zs0d+v/+ca0VFRTp48KCqq6u1c+dO1dXVqbi4eKijAABGOOdQXzB37lzNnTv3E9ccP35cd911l15//XXNmzcv5tqhQ4e0a9cu7du3T3l5eZKkTZs26eabb9ajjz46YNQAAJ9Pcf+aVCQS0dKlS7Vq1SpdddVV51yvr69XWlpaNFCSFAgElJKSooaGhgHv2dvbq3A4HHMAAEa+uEdq3bp1cjqduvvuuwe8HgwGlZmZGXPO6XQqPT1dwWBwwNdUVlbK4/FEj+zs7HiPDQCwUFwj1djYqJ///Od69tln5XA44nbf8vJyhUKh6NHa2hq3ewMA7BXXSP3hD39Qe3u7cnJy5HQ65XQ6dfToUd17772aNGmSJMnn86m9vT3mdWfOnFFHR4d8Pt+A93W5XHK73TEHAGDkG/KDE59k6dKlCgQCMecKCwu1dOlS3X777ZKkgoICdXZ2qrGxUTNnzpQk7d69W5FIRPn5+fEcBwCQ5IYcqa6uLr3//vvRj1taWnTgwAGlp6crJydHGRkZMetHjx4tn8+nf/u3f5MkTZ06VTfddJPuuOMOVVVVqa+vT6WlpVq8eDFP9gEAYgz5r/v279+vGTNmaMaMGZKksrIyzZgxQ2vXrh30PbZu3aopU6Zozpw5uvnmmzV79mw99dRTQx0FADDCDfmd1I033ihjzKDX/+UvfznnXHp6urZt2zbUXxoA8DnD9+4DAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYK24/vj4i+Xsz7MKd0USPAmGW9fpiPp7e3TG9CV6FABxdEb/+jN9oZ9P6DBD+QmGlvjrX/+q7OzsRI8BAPiMWltbddlll533elJGKhKJqLm5WVdeeaVaW1vldrsTPdKnFg6HlZ2dndT7YA92GAl7kEbGPtjDhRljdPLkSfn9fqWknP8rT0n5130pKSm69NJLJUlutztpfxP8fyNhH+zBDiNhD9LI2Ad7+GQej+eCa3hwAgBgLSIFALBW0kbK5XLpgQcekMvlSvQon8lI2Ad7sMNI2IM0MvbBHuInKR+cAAB8PiTtOykAwMhHpAAA1iJSAABrESkAgLWIFADAWkkbqc2bN2vSpEkaO3as8vPztXfv3kSPdF6VlZW67rrrNGHCBGVmZmrhwoVqbm6OWdPT06OSkhJlZGRo/PjxWrRokdra2hI08YU98sgjcjgcWrFiRfRcMuzh+PHj+u53v6uMjAylpqZq2rRp2r9/f/S6MUZr165VVlaWUlNTFQgEdOTIkQROfK7+/n5VVFQoNzdXqamp+tKXvqSf/exnMd+o07Z91NXVaf78+fL7/XI4HNqxY0fM9cHM29HRoaKiIrndbqWlpWn58uXq6uqyYg99fX1avXq1pk2bpnHjxsnv9+t73/ueTpw4YdUeLrSPj7vzzjvlcDi0YcOGmPMXcx9JGanf/OY3Kisr0wMPPKCmpiZNnz5dhYWFam9vT/RoA6qtrVVJSYn27Nmj6upq9fX16Zvf/Ka6u7uja1auXKlXX31VL774ompra3XixAndeuutCZz6/Pbt26df/vKX+spXvhJz3vY9fPjhh5o1a5ZGjx6t1157Te+++64ee+wxTZw4Mbpm/fr12rhxo6qqqtTQ0KBx48apsLBQPT09CZw81rp167RlyxY9+eSTOnTokNatW6f169dr06ZN0TW27aO7u1vTp0/X5s2bB7w+mHmLiop08OBBVVdXa+fOnaqrq1NxcfHF2sIn7uHUqVNqampSRUWFmpqa9NJLL6m5uVm33HJLzLpE70G68P+Ls7Zv3649e/bI7/efc+2i7sMkoeuvv96UlJREP+7v7zd+v99UVlYmcKrBa29vN5JMbW2tMcaYzs5OM3r0aPPiiy9G1xw6dMhIMvX19Ykac0AnT540kydPNtXV1ebf//3fzT333GOMSY49rF692syePfu81yORiPH5fOZ//ud/ouc6OzuNy+Uyv/71ry/GiIMyb9488/3vfz/m3K233mqKioqMMfbvQ5LZvn179OPBzPvuu+8aSWbfvn3RNa+99ppxOBzm+PHjF232sz6+h4Hs3bvXSDJHjx41xti3B2POv4+//vWv5tJLLzXvvPOOufzyy80TTzwRvXax95F076ROnz6txsZGBQKB6LmUlBQFAgHV19cncLLBC4VCkqT09HRJUmNjo/r6+mL2NGXKFOXk5Fi3p5KSEs2bNy9mVik59vDKK68oLy9P3/72t5WZmakZM2bo6aefjl5vaWlRMBiM2YPH41F+fr41e5CkG264QTU1NTp8+LAk6Y9//KPefPNNzZ07V1Ly7OOswcxbX1+vtLQ05eXlRdcEAgGlpKSooaHhos88GKFQSA6HQ2lpaZKSZw+RSERLly7VqlWrdNVVV51z/WLvI+m+C/o//vEP9ff3y+v1xpz3er167733EjTV4EUiEa1YsUKzZs3S1VdfLUkKBoMaM2ZM9DfzWV6vV8FgMAFTDuyFF15QU1OT9u3bd861ZNjDBx98oC1btqisrEw/+tGPtG/fPt19990aM2aMli1bFp1zoN9btuxBktasWaNwOKwpU6Zo1KhR6u/v10MPPaSioiJJSpp9nDWYeYPBoDIzM2OuO51OpaenW7mnnp4erV69WkuWLIl+B/Fk2cO6devkdDp19913D3j9Yu8j6SKV7EpKSvTOO+/ozTffTPQoQ9La2qp77rlH1dXVGjt2bKLH+VQikYjy8vL08MMPS5JmzJihd955R1VVVVq2bFmCpxu83/72t9q6dau2bdumq666SgcOHNCKFSvk9/uTah8jVV9fn2677TYZY7Rly5ZEjzMkjY2N+vnPf66mpiY5HI5EjyMpCR+cuOSSSzRq1Khznhpra2uTz+dL0FSDU1paqp07d+qNN96I+UmUPp9Pp0+fVmdnZ8x6m/bU2Nio9vZ2XXvttXI6nXI6naqtrdXGjRvldDrl9Xqt30NWVpauvPLKmHNTp07VsWPHJCk6p+2/t1atWqU1a9Zo8eLFmjZtmpYuXaqVK1eqsrJSUvLs46zBzOvz+c55MOrMmTPq6Oiwak9nA3X06FFVV1fH/BymZNjDH/7wB7W3tysnJyf65/zo0aO69957NWnSJEkXfx9JF6kxY8Zo5syZqqmpiZ6LRCKqqalRQUFBAic7P2OMSktLtX37du3evVu5ubkx12fOnKnRo0fH7Km5uVnHjh2zZk9z5szRn//8Zx04cCB65OXlqaioKPrftu9h1qxZ5zz6f/jwYV1++eWSpNzcXPl8vpg9hMNhNTQ0WLMH6V9Pkn38J5mOGjVKkUhEUvLs46zBzFtQUKDOzk41NjZG1+zevVuRSET5+fkXfeaBnA3UkSNH9Lvf/U4ZGRkx15NhD0uXLtWf/vSnmD/nfr9fq1at0uuvvy4pAfuI+6MYF8ELL7xgXC6XefbZZ827775riouLTVpamgkGg4kebUA/+MEPjMfjMb///e/N3/72t+hx6tSp6Jo777zT5OTkmN27d5v9+/ebgoICU1BQkMCpL+z/P91njP172Lt3r3E6neahhx4yR44cMVu3bjVf+MIXzK9+9avomkceecSkpaWZl19+2fzpT38yCxYsMLm5ueajjz5K4OSxli1bZi699FKzc+dO09LSYl566SVzySWXmPvuuy+6xrZ9nDx50rz99tvm7bffNpLM448/bt5+++3ok2+Dmfemm24yM2bMMA0NDebNN980kydPNkuWLLFiD6dPnza33HKLueyyy8yBAwdi/pz39vZas4cL7WMgH3+6z5iLu4+kjJQxxmzatMnk5OSYMWPGmOuvv97s2bMn0SOdl6QBj2eeeSa65qOPPjI//OEPzcSJE80XvvAF8x//8R/mb3/7W+KGHoSPRyoZ9vDqq6+aq6++2rhcLjNlyhTz1FNPxVyPRCKmoqLCeL1e43K5zJw5c0xzc3OCph1YOBw299xzj8nJyTFjx441V1xxhfnxj38c88nQtn288cYbA/4ZWLZs2aDn/ec//2mWLFlixo8fb9xut7n99tvNyZMnrdhDS0vLef+cv/HGG9bs4UL7GMhAkbqY++DnSQEArJV0X5MCAHx+ECkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWv8HlqYq5xWvU40AAAAASUVORK5CYII="
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "execution_count": 4
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    "## Run iDDN\n",
    "Now we will apply iDDN to estimate the common and differential networks. The common network consists of edges that occurs in both conditions, and the differential network consists of edges that exist only in one of the conditions. iDDN estimate the networks under two conditions with sparsity penalties and difference penalties, and then find the common and differential networks.\n",
    "\n",
    "We first need to specify two parameters: $\\lambda_1$ and $\\lambda_2$. The first one determines the overall sparsity of the network (i.e., the number of edges), and the second one controls the similarity of the estimated networks between the two conditions. The choice of these two parameters depends on the biological understanding of the problems: how many edges are expected, and how similar the two conditions should be. Alternatively, we can use a data driven approach to choose the two parameters using cross validation, which will be covered in the third tutorial. However, it is still better to utilize the understanding of each problem to experiment and select the parameters. \n",
    "\n",
    "Here we provide the values of $\\lambda_1=0.15$ and $\\lambda_2=0.02$ that works well in this synthetic data (as we already know hte ground truth). These values may or may not work well for real applications, though they can serve as a reasonable starting point. "
   ],
   "id": "a12a8933a7b7e044"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:41.519895Z",
     "start_time": "2024-08-26T14:56:41.460054Z"
    }
   },
   "cell_type": "code",
   "source": [
    "lambda1 = example[\"lambda1\"]\n",
    "lambda2 = example[\"lambda2\"]"
   ],
   "id": "e151c4b383d75303",
   "outputs": [],
   "execution_count": 5
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "Actually, iDDN allows the specification of $\\lambda_1$ and $\\lambda_2$ for each edge, and we need to provide two matrices encoding the penalties for each edge. For simplicity, we use the same value for all edges. If we do have some prior information, we may wish to set different penalties for different types of edges. For example, we can have stronger penalties for edges among mRNAs than edges between regulators and mRNAs.",
   "id": "7f7902f35abab29e"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:41.611649Z",
     "start_time": "2024-08-26T14:56:41.553804Z"
    }
   },
   "cell_type": "code",
   "source": [
    "lambda1_mat = np.copy(dep_mat) * lambda1\n",
    "lambda2_mat = np.copy(dep_mat) * lambda2"
   ],
   "id": "d14f6e37fbb53c48",
   "outputs": [],
   "execution_count": 6
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    "Now we can run iDDN. It requires the following information:\n",
    "1. Input data under the first condition. NumPy array: number_sample_condition_1 x feature.\n",
    "2. Input data under the second condition. NumPy array: number_sample_condition_2 x feature.\n",
    "3. $\\lambda_1$ penalty for each pair of edge.\n",
    "4. $\\lambda_2$ penalty for each pair of edge.\n",
    "5. Constraints matrix.\n",
    "\n",
    "As we already have all of these, we can run iDDN. Note that the first time we run iDDN, it will take one or two seconds to compile part of the code *just in time*, as we use `Numba` to accelerate part of the code. If you call iDDN for the second time, you can find it to be much faster (especially when the data is small)."
   ],
   "id": "acfc91c7edd55a13"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:44.629432Z",
     "start_time": "2024-08-26T14:56:41.645558Z"
    }
   },
   "cell_type": "code",
   "source": [
    "out_iddn = iddn.iddn(\n",
    "    dat1,\n",
    "    dat2,\n",
    "    dep_mat=dep_mat,\n",
    "    lambda1=lambda1_mat,\n",
    "    lambda2=lambda2_mat,\n",
    ")"
   ],
   "id": "3a786c40b9dacc5f",
   "outputs": [],
   "execution_count": 7
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    "iDDN will provide a NumPy array of shape $2\\times 150 \\times 150$. Here `out_ddn[0]` is the estimated dependency matrix for condition 1, and `out_ddn[1]` for condition 2. Due to numerically precisions and convergence threshold, there could be some very small values that need to be cleaned. In addition, these two matrices should be converted to the common and differential networks we need.\n",
    "\n",
    "For convenience, we provide a pipeline that can run iDDN and provide these information. Besides, we do not need to create `lambda1_mat` and `lambda2_mat` ourselves.\n",
    "The output of this function is a dictionary with five items:\n",
    "1. \"comm\": estimated common network matrix\n",
    "2. \"diff\": estimated differential network matrix\n",
    "3. \"g1\": estimated network for the first condition\n",
    "4. \"g2\": estimated network for the first condition\n",
    "5. \"out_ddn\": the raw output of iDDN"
   ],
   "id": "3ef8e06844e53151"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:44.960964Z",
     "start_time": "2024-08-26T14:56:44.693453Z"
    }
   },
   "cell_type": "code",
   "source": [
    "res = tools.iddn_basic_pipeline(dat1, dat2, dep_mat, lambda1, lambda2)\n",
    "out_iddn_1 = res[\"out_iddn\"]"
   ],
   "id": "544e59e746427edc",
   "outputs": [],
   "execution_count": 8
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "We can compare the raw output of iDDN with reference outputs and their difference should be very small.",
   "id": "8b08788f65121aa"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:45.054867Z",
     "start_time": "2024-08-26T14:56:44.992996Z"
    }
   },
   "cell_type": "code",
   "source": [
    "out_iddn_expected = example[\"out_ddn\"]\n",
    "assert np.sum(np.abs(out_iddn_expected - out_iddn_1)) < 1e-8\n",
    "assert np.sum(np.abs(out_iddn_expected - out_iddn)) < 1e-8"
   ],
   "id": "6ca3805349d03e73",
   "outputs": [],
   "execution_count": 9
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    "## Analysis of the results\n",
    "As a simple example, we calculate the node degrees and find hubs nodes in common and differential networks. This can be achieved by simply taking summations of each row in the common and differential network matrices."
   ],
   "id": "d7013a8dca9300f0"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:45.149840Z",
     "start_time": "2024-08-26T14:56:45.086994Z"
    }
   },
   "cell_type": "code",
   "source": [
    "comm_est = res[\"comm\"]\n",
    "diff_est = res[\"diff\"]\n",
    "deg_comm = np.sum(comm_est, axis=1)\n",
    "deg_diff = np.sum(diff_est, axis=1)"
   ],
   "id": "b66b21c59e286056",
   "outputs": [],
   "execution_count": 10
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "If the node names are available, we can display the names of hub nodes.",
   "id": "e7fab07e02d56f74"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:45.245554Z",
     "start_time": "2024-08-26T14:56:45.182749Z"
    }
   },
   "cell_type": "code",
   "source": [
    "node_names = example[\"node_names\"]\n",
    "print(\"Common network hub nodes: \", node_names[deg_comm>5])\n",
    "print(\"Differential network hub nodes: \", node_names[deg_diff>5])"
   ],
   "id": "c3219f05c5a3d3f7",
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Common network hub nodes:  ['mrna1' 'mrna2' 'mrna4' 'mrna32' 'mrna35' 'mrna45']\n",
      "Differential network hub nodes:  ['mrna1' 'mrna2' 'mrna24' 'mrna39']\n"
     ]
    }
   ],
   "execution_count": 11
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    "## Save the networks\n",
    "We can convert the detected networks to Pandas data frames."
   ],
   "id": "d2723a4354235dde"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:45.353545Z",
     "start_time": "2024-08-26T14:56:45.278466Z"
    }
   },
   "cell_type": "code",
   "source": "df_edge_comm, df_edge_diff, _, _ = tools.iddn_output_to_csv(out_iddn, node_names)",
   "id": "4f15f33b2a7a2cc5",
   "outputs": [],
   "execution_count": 12
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "We can save them to CSV files if needed.",
   "id": "2e6337166a56add"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:45.480665Z",
     "start_time": "2024-08-26T14:56:45.409843Z"
    }
   },
   "cell_type": "code",
   "source": [
    "# df_edge_comm.to_csv(\"./common_network.csv\")\n",
    "# df_edge_diff.to_csv(\"./differential_network.csv\")"
   ],
   "id": "f54aab1d2941ea09",
   "outputs": [],
   "execution_count": 13
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    "## Optinoal: compare with ground truth\n",
    "As our example is from a synthetic data with known ground truth, we can evaluate the accuracy of our results.\n",
    "First we load the ground truth network matrices."
   ],
   "id": "be38cc860119520b"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:45.638568Z",
     "start_time": "2024-08-26T14:56:45.536658Z"
    }
   },
   "cell_type": "code",
   "source": [
    "comm_gt = example[\"comm_gt\"]  # common network ground truth\n",
    "diff_gt = example[\"diff_gt\"]  # differential network ground truth\n",
    "g1_gt = example[\"g1_gt\"]  # first condition network ground truth\n",
    "g2_gt = example[\"g2_gt\"]  # second condition network ground truth"
   ],
   "id": "1767cfc85df80894",
   "outputs": [],
   "execution_count": 14
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "We provide a function to calculate the precision, recall, and F1 score for each network. Below shows the results for common and differential networks.",
   "id": "189d27f441ec0ac1"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:45.749288Z",
     "start_time": "2024-08-26T14:56:45.671497Z"
    }
   },
   "cell_type": "code",
   "source": [
    "tools.evaluate_metrics(comm_gt, comm_est)\n",
    "tools.evaluate_metrics(diff_gt, diff_est)"
   ],
   "id": "42d28b846b503034",
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "recall=0.9929577464788732, precision=1.0, F1=0.9964664310954063\n",
      "recall=0.8990825688073395, precision=0.98989898989899, F1=0.9423076923076924\n"
     ]
    }
   ],
   "execution_count": 15
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    " The sample size is this simulaiton is pretty large, so we expect that the performance should be quite good. In real applications, the sample sizes could be much smaller, and the accuracy of iDDN, especially for differential networks, may decrease. \n",
    " \n",
    "We can also draw the difference between ground truth network and estimated network matrices. First we show common network estimation errors."
   ],
   "id": "95dd6eedde1367ab"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:45.971693Z",
     "start_time": "2024-08-26T14:56:45.783198Z"
    }
   },
   "cell_type": "code",
   "source": [
    "plt.figure()\n",
    "plt.imshow(comm_gt!=comm_est)\n",
    "plt.show()"
   ],
   "id": "953655487a867eef",
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ],
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAakAAAGhCAYAAADbf0s2AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAA9hAAAPYQGoP6dpAAAjPUlEQVR4nO3df3BU9b3/8deGDZsIZENis5vVRNIOcwNKEQmmEaaXltwGYRCutF6clHItY6pNVIgXIfcavL1VI9QqBSmpTkftXKytMwKFueKkQZIyhhAS8AdixGkGInST9tLsQjAhZD/3j37Zb1eCEtywn43Px8yZMed89vD+jMBzlpwkDmOMEQAAFkqI9QAAAFwMkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWCtmkdq4caPGjRunpKQk5efna9++fbEaBQBgqZhE6je/+Y3Ky8v1yCOPqKWlRZMnT1ZRUZE6OztjMQ4AwFKOWHyD2fz8fE2bNk3PPPOMJCkUCikrK0v33XefVq1a9ZmvD4VCOnHihMaMGSOHwzHU4wIAoswYo1OnTsnn8ykh4eLvl5xXcCZJ0tmzZ9Xc3KyKiorwuYSEBBUWFqqhoWHA1/T29qq3tzf88fHjxzVx4sQhnxUAMLTa29t17bXXXvT6FY/UX/7yF/X398vj8USc93g8ev/99wd8TVVVlX70ox9dcH6G5sipxCGZEwAwdM6pT3v0PxozZsynrrvikbocFRUVKi8vD38cDAaVlZUlpxLldBApAIg7/+8TTZ/1KZsrHqmrr75aI0aMUEdHR8T5jo4Oeb3eAV/jcrnkcrmuxHgAYsj55XHqzv2SEs6GlPROu/o7eJjqi+6KP903cuRITZ06VbW1teFzoVBItbW1KigouNLjALCI/58yNe3R/brq4RM6c9N1sR4HFojJP/eVl5dryZIlysvL080336x169apu7tbd911VyzGAWCJvtEOzXa/rbHOM3otOTvW48ACMYnUv/zLv+jPf/6zVq9eLb/frxtvvFE7d+684GEKAF8s3oZu3T/yB3Keka451Kn+WA+EmIvJ10l9XsFgUG63WzM1nwcngOHm/CfS4++vJgzCOdOn3dqmQCCglJSUi66Li6f7AHyBECf8Hb7BLADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKwV9UhVVVVp2rRpGjNmjDIyMrRgwQK1trZGrOnp6VFpaanS09M1evRoLVy4UB0dHdEeBQAQ56Ieqbq6OpWWlmrv3r2qqalRX1+fvvWtb6m7uzu8Zvny5dq+fbteeeUV1dXV6cSJE7r99tujPQoAIM45jDFmKH+BP//5z8rIyFBdXZ2+/vWvKxAI6Etf+pJeeuklffvb35Ykvf/++5owYYIaGhr0ta997TPvGQwG5Xa7NVPz5XQkDuX4AIAhcM70abe2KRAIKCUl5aLrhvxzUoFAQJKUlpYmSWpublZfX58KCwvDa3Jzc5Wdna2GhoYB79Hb26tgMBhxAACGvyGNVCgU0rJlyzR9+nTdcMMNkiS/36+RI0cqNTU1Yq3H45Hf7x/wPlVVVXK73eEjKytrKMcGAFhiSCNVWlqqd999Vy+//PLnuk9FRYUCgUD4aG9vj9KEAACbOYfqxmVlZdqxY4fq6+t17bXXhs97vV6dPXtWXV1dEe+mOjo65PV6B7yXy+WSy+UaqlEBAJaK+jspY4zKysq0ZcsW7dq1Szk5ORHXp06dqsTERNXW1obPtba26tixYyooKIj2OACAOBb1d1KlpaV66aWXtG3bNo0ZMyb8eSa3263k5GS53W4tXbpU5eXlSktLU0pKiu677z4VFBRc0pN9AIAvjqhHatOmTZKkmTNnRpx//vnn9a//+q+SpKeffloJCQlauHChent7VVRUpJ///OfRHgUAEOeG/OukhgJfJwUA8c2ar5MCAOByESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKw15JF64okn5HA4tGzZsvC5np4elZaWKj09XaNHj9bChQvV0dEx1KMAAOLMkEaqqalJv/jFL/TVr3414vzy5cu1fft2vfLKK6qrq9OJEyd0++23D+UoAIA4NGSROn36tIqLi/Xcc89p7Nix4fOBQEC//OUv9dRTT+mb3/ympk6dqueff15vvvmm9u7dO+C9ent7FQwGIw4AwPA3ZJEqLS3V3LlzVVhYGHG+ublZfX19Eedzc3OVnZ2thoaGAe9VVVUlt9sdPrKysoZqbACARYYkUi+//LJaWlpUVVV1wTW/36+RI0cqNTU14rzH45Hf7x/wfhUVFQoEAuGjvb19KMYGAFjGGe0btre364EHHlBNTY2SkpKick+XyyWXyxWVewEA4kfU30k1Nzers7NTN910k5xOp5xOp+rq6rR+/Xo5nU55PB6dPXtWXV1dEa/r6OiQ1+uN9jgAgDgW9XdSs2bN0jvvvBNx7q677lJubq5WrlyprKwsJSYmqra2VgsXLpQktba26tixYyooKIj2OACAOBb1SI0ZM0Y33HBDxLlRo0YpPT09fH7p0qUqLy9XWlqaUlJSdN9996mgoEBf+9rXoj0OACCORT1Sl+Lpp59WQkKCFi5cqN7eXhUVFennP/95LEYBAFjMYYwxsR5isILBoNxut2ZqvpyOxFiPAwAYpHOmT7u1TYFAQCkpKRddx/fuAwBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkMSqePHj+u73/2u0tPTlZycrEmTJmn//v3h68YYrV69WpmZmUpOTlZhYaGOHDkyFKMAAOJY1CP117/+VdOnT1diYqJee+01vffee/rpT3+qsWPHhtesXbtW69evV3V1tRobGzVq1CgVFRWpp6cn2uMAAOKYM9o3XLNmjbKysvT888+Hz+Xk5IT/2xijdevW6eGHH9b8+fMlSb/61a/k8Xi0detWLVq0KNojAQDiVNTfSf3ud79TXl6evvOd7ygjI0NTpkzRc889F77e1tYmv9+vwsLC8Dm32638/Hw1NDQMeM/e3l4Fg8GIAwAw/EU9Un/84x+1adMmjR8/Xq+//rruvfde3X///XrxxRclSX6/X5Lk8XgiXufxeMLXPqmqqkputzt8ZGVlRXtsAICFoh6pUCikm266SY8//rimTJmikpIS3X333aqurr7se1ZUVCgQCISP9vb2KE4MALBV1COVmZmpiRMnRpybMGGCjh07Jknyer2SpI6Ojog1HR0d4Wuf5HK5lJKSEnEAAIa/qEdq+vTpam1tjTj3wQcf6LrrrpP0t4covF6vamtrw9eDwaAaGxtVUFAQ7XEAAHEs6k/3LV++XLfccosef/xx3XHHHdq3b5+effZZPfvss5Ikh8OhZcuW6dFHH9X48eOVk5OjyspK+Xw+LViwINrjAADiWNQjNW3aNG3ZskUVFRX6r//6L+Xk5GjdunUqLi4Or3nooYfU3d2tkpISdXV1acaMGdq5c6eSkpKiPQ4AII45jDEm1kMMVjAYlNvt1kzNl9ORGOtxAACDdM70abe2KRAIfOpzBnzvPgCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArBX1SPX396uyslI5OTlKTk7WV77yFf34xz+WMSa8xhij1atXKzMzU8nJySosLNSRI0eiPQoAIM5FPVJr1qzRpk2b9Mwzz+jw4cNas2aN1q5dqw0bNoTXrF27VuvXr1d1dbUaGxs1atQoFRUVqaenJ9rjAADimDPaN3zzzTc1f/58zZ07V5I0btw4/frXv9a+ffsk/e1d1Lp16/Twww9r/vz5kqRf/epX8ng82rp1qxYtWhTtkQAAcSrq76RuueUW1dbW6oMPPpAkvfXWW9qzZ49uvfVWSVJbW5v8fr8KCwvDr3G73crPz1dDQ8OA9+zt7VUwGIw4AADDX9TfSa1atUrBYFC5ubkaMWKE+vv79dhjj6m4uFiS5Pf7JUkejyfidR6PJ3ztk6qqqvSjH/0o2qMCACwX9XdSv/3tb7V582a99NJLamlp0Ysvvqgnn3xSL7744mXfs6KiQoFAIHy0t7dHcWIAgK2i/k5qxYoVWrVqVfhzS5MmTdLRo0dVVVWlJUuWyOv1SpI6OjqUmZkZfl1HR4duvPHGAe/pcrnkcrmiPSoAwHJRfyd15swZJSRE3nbEiBEKhUKSpJycHHm9XtXW1oavB4NBNTY2qqCgINrjAADiWNTfSc2bN0+PPfaYsrOzdf311+vAgQN66qmn9P3vf1+S5HA4tGzZMj366KMaP368cnJyVFlZKZ/PpwULFkR7HABAHIt6pDZs2KDKykr98Ic/VGdnp3w+n37wgx9o9erV4TUPPfSQuru7VVJSoq6uLs2YMUM7d+5UUlJStMcBAMQxh/n7bwURJ4LBoNxut2ZqvpyOxFiPAwAYpHOmT7u1TYFAQCkpKRddx/fuAwBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWoOOVH19vebNmyefzyeHw6GtW7dGXDfGaPXq1crMzFRycrIKCwt15MiRiDUnT55UcXGxUlJSlJqaqqVLl+r06dOfayMAgOFn0JHq7u7W5MmTtXHjxgGvr127VuvXr1d1dbUaGxs1atQoFRUVqaenJ7ymuLhYhw4dUk1NjXbs2KH6+nqVlJRc/i4AAMOSwxhjLvvFDoe2bNmiBQsWSPrbuyifz6cHH3xQ//Zv/yZJCgQC8ng8euGFF7Ro0SIdPnxYEydOVFNTk/Ly8iRJO3fu1Jw5c/TRRx/J5/N95q8bDAbldrs1U/PldCRKkpxfHif/P2Wqb7RD3je75Wh463K3BQAYYudMn3ZrmwKBgFJSUi66Lqqfk2pra5Pf71dhYWH4nNvtVn5+vhoaGiRJDQ0NSk1NDQdKkgoLC5WQkKDGxsYB79vb26tgMBhxfFJ37pc0s6RRT977nE58fZTkcERzawCAGIhqpPx+vyTJ4/FEnPd4POFrfr9fGRkZEdedTqfS0tLCaz6pqqpKbrc7fGRlZV2wJuFsSEdOZ+jAx+Pk/DgauwEAxFpcPN1XUVGhQCAQPtrb2y9Yk/ROu4JPZGnHI9+U7/VO6fL/FRMAYAlnNG/m9XolSR0dHcrMzAyf7+jo0I033hhe09nZGfG6c+fO6eTJk+HXf5LL5ZLL5frUX7u/o1Ou1/523/7L3QAAwCpRfSeVk5Mjr9er2tra8LlgMKjGxkYVFBRIkgoKCtTV1aXm5ubwml27dikUCik/Pz+a4wAA4tyg30mdPn1aH374YfjjtrY2HTx4UGlpacrOztayZcv06KOPavz48crJyVFlZaV8Pl/4CcAJEyZo9uzZuvvuu1VdXa2+vj6VlZVp0aJFl/RkHwDgi2PQkdq/f7++8Y1vhD8uLy+XJC1ZskQvvPCCHnroIXV3d6ukpERdXV2aMWOGdu7cqaSkpPBrNm/erLKyMs2aNUsJCQlauHCh1q9fH4XtAACGk8/1dVKxMtDXSQEA4kdMvk4KAIBoIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYC0iBQCwFpECAFhr0JGqr6/XvHnz5PP55HA4tHXr1vC1vr4+rVy5UpMmTdKoUaPk8/n0ve99TydOnIi4x8mTJ1VcXKyUlBSlpqZq6dKlOn369OfeDABgeBl0pLq7uzV58mRt3LjxgmtnzpxRS0uLKisr1dLSoldffVWtra267bbbItYVFxfr0KFDqqmp0Y4dO1RfX6+SkpLL3wUAYFhyGGPMZb/Y4dCWLVu0YMGCi65pamrSzTffrKNHjyo7O1uHDx/WxIkT1dTUpLy8PEnSzp07NWfOHH300Ufy+Xyf+esGg0G53W7N1Hw5HYmXOz4AIEbOmT7t1jYFAgGlpKRcdN2Qf04qEAjI4XAoNTVVktTQ0KDU1NRwoCSpsLBQCQkJamxsHPAevb29CgaDEQcAYPgb0kj19PRo5cqVuvPOO8Ol9Pv9ysjIiFjndDqVlpYmv98/4H2qqqrkdrvDR1ZW1lCODQCwxJBFqq+vT3fccYeMMdq0adPnuldFRYUCgUD4aG9vj9KUAACbOYfipucDdfToUe3atSvi3xu9Xq86Ozsj1p87d04nT56U1+sd8H4ul0sul2soRgUAWCzq76TOB+rIkSP6/e9/r/T09IjrBQUF6urqUnNzc/jcrl27FAqFlJ+fH+1xAABxbNDvpE6fPq0PP/ww/HFbW5sOHjyotLQ0ZWZm6tvf/rZaWlq0Y8cO9ff3hz/PlJaWppEjR2rChAmaPXu27r77blVXV6uvr09lZWVatGjRJT3ZBwD44hj0I+i7d+/WN77xjQvOL1myRP/5n/+pnJycAV/3xhtvaObMmZL+9sW8ZWVl2r59uxISErRw4UKtX79eo0ePvqQZeAQdAOLbpT6C/rm+TipWiBQAxDdrvk4KAIDLRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgLSIFALAWkQIAWItIAQCsRaQAANYiUgAAaxEpAIC1iBQAwFpECgBgrUFHqr6+XvPmzZPP55PD4dDWrVsvuvaee+6Rw+HQunXrIs6fPHlSxcXFSklJUWpqqpYuXarTp08PdhQAwDA36Eh1d3dr8uTJ2rhx46eu27Jli/bu3Sufz3fBteLiYh06dEg1NTXasWOH6uvrVVJSMthRAADDnHOwL7j11lt16623fuqa48eP67777tPrr7+uuXPnRlw7fPiwdu7cqaamJuXl5UmSNmzYoDlz5ujJJ58cMGoAgC+mqH9OKhQKafHixVqxYoWuv/76C643NDQoNTU1HChJKiwsVEJCghobGwe8Z29vr4LBYMQBABj+oh6pNWvWyOl06v777x/wut/vV0ZGRsQ5p9OptLQ0+f3+AV9TVVUlt9sdPrKysqI9NgDAQlGNVHNzs372s5/phRdekMPhiNp9KyoqFAgEwkd7e3vU7g0AsFdUI/WHP/xBnZ2dys7OltPplNPp1NGjR/Xggw9q3LhxkiSv16vOzs6I1507d04nT56U1+sd8L4ul0spKSkRBwBg+Bv0gxOfZvHixSosLIw4V1RUpMWLF+uuu+6SJBUUFKirq0vNzc2aOnWqJGnXrl0KhULKz8+P5jgAgDg36EidPn1aH374YfjjtrY2HTx4UGlpacrOzlZ6enrE+sTERHm9Xv3DP/yDJGnChAmaPXu27r77blVXV6uvr09lZWVatGgRT/YBACIM+p/79u/frylTpmjKlCmSpPLyck2ZMkWrV6++5Hts3rxZubm5mjVrlubMmaMZM2bo2WefHewoAIBhzmGMMbEeYrCCwaDcbrdmar6cjsRYjwMAGKRzpk+7tU2BQOBTnzPge/cBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsBaRAgBYi0gBAKxFpAAA1iJSAABrESkAgLWIFADAWkQKAGAtIgUAsFZUf3z8lXL+R2CdU58Udz8NCwBwTn2S/v/f5xcTl5E6deqUJGmP/ifGkwAAPo9Tp07J7XZf9Hpc/mTeUCik1tZWTZw4Ue3t7Z/6Ux1tFwwGlZWVFdf7YA92GA57kIbHPtjDZzPG6NSpU/L5fEpIuPhnnuLynVRCQoKuueYaSVJKSkrc/ib4e8NhH+zBDsNhD9Lw2Ad7+HSf9g7qPB6cAABYi0gBAKwVt5FyuVx65JFH5HK5Yj3K5zIc9sEe7DAc9iANj32wh+iJywcnAABfDHH7TgoAMPwRKQCAtYgUAMBaRAoAYC0iBQCwVtxGauPGjRo3bpySkpKUn5+vffv2xXqki6qqqtK0adM0ZswYZWRkaMGCBWptbY1Y09PTo9LSUqWnp2v06NFauHChOjo6YjTxZ3viiSfkcDi0bNmy8Ll42MPx48f13e9+V+np6UpOTtakSZO0f//+8HVjjFavXq3MzEwlJyersLBQR44cieHEF+rv71dlZaVycnKUnJysr3zlK/rxj38c8Y06bdtHfX295s2bJ5/PJ4fDoa1bt0Zcv5R5T548qeLiYqWkpCg1NVVLly7V6dOnrdhDX1+fVq5cqUmTJmnUqFHy+Xz63ve+pxMnTli1h8/axyfdc889cjgcWrduXcT5K7mPuIzUb37zG5WXl+uRRx5RS0uLJk+erKKiInV2dsZ6tAHV1dWptLRUe/fuVU1Njfr6+vStb31L3d3d4TXLly/X9u3b9corr6iurk4nTpzQ7bffHsOpL66pqUm/+MUv9NWvfjXivO17+Otf/6rp06crMTFRr732mt577z399Kc/1dixY8Nr1q5dq/Xr16u6ulqNjY0aNWqUioqK1NPTE8PJI61Zs0abNm3SM888o8OHD2vNmjVau3atNmzYEF5j2z66u7s1efJkbdy4ccDrlzJvcXGxDh06pJqaGu3YsUP19fUqKSm5Ulv41D2cOXNGLS0tqqysVEtLi1599VW1trbqtttui1gX6z1In/3/4rwtW7Zo79698vl8F1y7ovswcejmm282paWl4Y/7+/uNz+czVVVVMZzq0nV2dhpJpq6uzhhjTFdXl0lMTDSvvPJKeM3hw4eNJNPQ0BCrMQd06tQpM378eFNTU2P+8R//0TzwwAPGmPjYw8qVK82MGTMuej0UChmv12t+8pOfhM91dXUZl8tlfv3rX1+JES/J3Llzzfe///2Ic7fffrspLi42xti/D0lmy5Yt4Y8vZd733nvPSDJNTU3hNa+99ppxOBzm+PHjV2z28z65h4Hs27fPSDJHjx41xti3B2Muvo+PPvrIXHPNNebdd9811113nXn66afD1670PuLundTZs2fV3NyswsLC8LmEhAQVFhaqoaEhhpNdukAgIElKS0uTJDU3N6uvry9iT7m5ucrOzrZuT6WlpZo7d27ErFJ87OF3v/ud8vLy9J3vfEcZGRmaMmWKnnvuufD1trY2+f3+iD243W7l5+dbswdJuuWWW1RbW6sPPvhAkvTWW29pz549uvXWWyXFzz7Ou5R5GxoalJqaqry8vPCawsJCJSQkqLGx8YrPfCkCgYAcDodSU1Mlxc8eQqGQFi9erBUrVuj666+/4PqV3kfcfRf0v/zlL+rv75fH44k47/F49P7778doqksXCoW0bNkyTZ8+XTfccIMkye/3a+TIkeHfzOd5PB75/f4YTDmwl19+WS0tLWpqarrgWjzs4Y9//KM2bdqk8vJy/fu//7uampp0//33a+TIkVqyZEl4zoF+b9myB0latWqVgsGgcnNzNWLECPX39+uxxx5TcXGxJMXNPs67lHn9fr8yMjIirjudTqWlpVm5p56eHq1cuVJ33nln+DuIx8se1qxZI6fTqfvvv3/A61d6H3EXqXhXWlqqd999V3v27In1KIPS3t6uBx54QDU1NUpKSor1OJclFAopLy9Pjz/+uCRpypQpevfdd1VdXa0lS5bEeLpL99vf/labN2/WSy+9pOuvv14HDx7UsmXL5PP54mofw1VfX5/uuOMOGWO0adOmWI8zKM3NzfrZz36mlpYWORyOWI8jKQ4fnLj66qs1YsSIC54a6+jokNfrjdFUl6asrEw7duzQG2+8oWuvvTZ83uv16uzZs+rq6opYb9Oempub1dnZqZtuuklOp1NOp1N1dXVav369nE6nPB6P9XvIzMzUxIkTI85NmDBBx44dk6TwnLb/3lqxYoVWrVqlRYsWadKkSVq8eLGWL1+uqqoqSfGzj/MuZV6v13vBg1Hnzp3TyZMnrdrT+UAdPXpUNTU1ET+HKR728Ic//EGdnZ3Kzs4O/zk/evSoHnzwQY0bN07Sld9H3EVq5MiRmjp1qmpra8PnQqGQamtrVVBQEMPJLs4Yo7KyMm3ZskW7du1STk5OxPWpU6cqMTExYk+tra06duyYNXuaNWuW3nnnHR08eDB85OXlqbi4OPzftu9h+vTpFzz6/8EHH+i6666TJOXk5Mjr9UbsIRgMqrGx0Zo9SH97kuyTP8l0xIgRCoVCkuJnH+ddyrwFBQXq6upSc3NzeM2uXbsUCoWUn59/xWceyPlAHTlyRL///e+Vnp4ecT0e9rB48WK9/fbbEX/OfT6fVqxYoddff11SDPYR9UcxroCXX37ZuFwu88ILL5j33nvPlJSUmNTUVOP3+2M92oDuvfde43a7ze7du82f/vSn8HHmzJnwmnvuucdkZ2ebXbt2mf3795uCggJTUFAQw6k/298/3WeM/XvYt2+fcTqd5rHHHjNHjhwxmzdvNldddZX57//+7/CaJ554wqSmpppt27aZt99+28yfP9/k5OSYjz/+OIaTR1qyZIm55pprzI4dO0xbW5t59dVXzdVXX20eeuih8Brb9nHq1Clz4MABc+DAASPJPPXUU+bAgQPhJ98uZd7Zs2ebKVOmmMbGRrNnzx4zfvx4c+edd1qxh7Nnz5rbbrvNXHvttebgwYMRf857e3ut2cNn7WMgn3y6z5gru4+4jJQxxmzYsMFkZ2ebkSNHmptvvtns3bs31iNdlKQBj+effz685uOPPzY//OEPzdixY81VV11l/vmf/9n86U9/it3Ql+CTkYqHPWzfvt3ccMMNxuVymdzcXPPss89GXA+FQqaystJ4PB7jcrnMrFmzTGtra4ymHVgwGDQPPPCAyc7ONklJSebLX/6y+Y//+I+Ivwxt28cbb7wx4J+BJUuWXPK8//u//2vuvPNOM3r0aJOSkmLuuusuc+rUKSv20NbWdtE/52+88YY1e/isfQxkoEhdyX3w86QAANaKu89JAQC+OIgUAMBaRAoAYC0iBQCwFpECAFiLSAEArEWkAADWIlIAAGsRKQCAtYgUAMBaRAoAYK3/A7gzXXMDqf16AAAAAElFTkSuQmCC"
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "execution_count": 16
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": "Then we show the errors in differential network. There are more errors compared with the common network case, which is expected. ",
   "id": "edbbcbc79b6ce124"
  },
  {
   "metadata": {
    "ExecuteTime": {
     "end_time": "2024-08-26T14:56:46.193100Z",
     "start_time": "2024-08-26T14:56:46.004605Z"
    }
   },
   "cell_type": "code",
   "source": [
    "plt.figure()\n",
    "plt.imshow(diff_gt!=diff_est)\n",
    "plt.show()"
   ],
   "id": "ed77df8d58412459",
   "outputs": [
    {
     "data": {
      "text/plain": [
       "<Figure size 640x480 with 1 Axes>"
      ],
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAakAAAGhCAYAAADbf0s2AAAAOXRFWHRTb2Z0d2FyZQBNYXRwbG90bGliIHZlcnNpb24zLjguMiwgaHR0cHM6Ly9tYXRwbG90bGliLm9yZy8g+/7EAAAACXBIWXMAAA9hAAAPYQGoP6dpAAArrElEQVR4nO3de3hU9YH/8c+ZTDIJl0wImAyRBCLLCihFBIwRfi2uqaAsQsVafFJKLStVg3JxvbBbsN2qEeqFohSq28fL1ntXsPCsuDQolJ8hQCIqFwP+RAjgJGDITBJIMpk5vz9cZzsS5DZhvpO8X89znic553sOn29N8uk5c2aOZdu2LQAADOSIdQAAAE6GkgIAGIuSAgAYi5ICABiLkgIAGIuSAgAYi5ICABiLkgIAGIuSAgAYi5ICABgrZiW1dOlS9evXT8nJycrLy9PmzZtjFQUAYKiYlNRrr72muXPn6sEHH1RFRYWGDh2qsWPHqqamJhZxAACGsmLxAbN5eXkaOXKknn76aUlSKBRSdna27rrrLj3wwAOn3D8UCunQoUPq3r27LMtq77gAgCizbVv19fXKysqSw3Hy8yXnecwkSWppaVF5ebnmzZsXXudwOFRQUKDS0tI292lublZzc3P4+4MHD2rw4MHtnhUA0L6qqqrUp0+fk24/7yV15MgRBYNBZWZmRqzPzMzUJ5980uY+xcXF+tWvfnXC+tG6Xk4ltktOAPHNkezS0UlDdWRESN32JujCP+1Va3XnfknBSkySf+JQ1eTZ6nIwQdlvHlDr/oNnfbyEXr30xU39Vf93QaVvc6jnmx8rdOz4ae3bqoA26r/UvXv3bx133kvqbMybN09z584Nf+/3+5WdnS2nEuW0KCkAJ0ro4taXYxO0/dplunnPjbI3ZMpZczTWsWLKkZSiI1e79MHEJ3VX1fU6XJotZ9XZF3dCrwvUc3Ktyga/riG9fy7PO24Fj7ee3s7/80LTqV6yOe8l1atXLyUkJKi6ujpifXV1tTweT5v7uFwuuVyu8xHPaM6L+qlx4AVytISU/HGVgp38/xUC38ZubVXyviTddeD72rX7Qg0+dlihWIeKtVBIXaoSNOvAOL2/u78GNtQreC7Hc1hKTAiqiyNJlsOWrOjfi3fe7+5LSkrS8OHDVVJSEl4XCoVUUlKi/Pz88x0nrni/31sjH9qqLr84pGOX9411HMBoocZjuujFgzp0dz8NWlyn4MEvYh0p5kLNzer7apWqZ/XVwEX1sj/bH+tIpxSTy31z587VtGnTNGLECF1xxRVavHixGhsbdeutt8YiTtwIdLM0zv2RejiP6e2UnFjHAcwWCqp17z5pr87tbKEjsW217quS9lVF53+TYEi+pmTtb21QsDlBsqN/rhqTkvrRj36kw4cPa8GCBfJ6vbrsssu0Zs2aE26mQCRPaaPuTvq5nMekC3fU8IsHILa+PCrrtYs1tuw+XbgjKLuhMer/REzeJ3Wu/H6/3G63xmhi57tx4usXGePvPxuAjugs/ya12gG9p7fk8/mUmpp60nFxcXcf/gblBMAk7fw3iQ+YBQAYizOpOOLo0kVW925SMKiQr152oCXWkaLO0b27HF27yG5uUdDnl0K88gZ0ZpRUHGn8/qU68IOg1ODUgD8elzZ9FOtIUeXo0kXeqUPU+N0GJXzSTRf9+z61Hjj7d8MDiH9c7osjtQOdWnP1Ei0oWKGG7JRYx4k6y+WS74pm7fo/z2vA1Z8plP7tH5cCoOPjTCqOpH4e0k0f/JMa/Cnq/0XzqXeIM3ZLi7pud+l7npvk3Z6hi/2H+IQAdHoJF/+dakf2khWU0jd98dV7vzoRSiqOpL29Uz22pEuhRoVqjnS4P+ChY8eU/cIeaWU3uY/v52OfAEk1oy/Qj+b8t44Eumtj8ZXqRknBVEG/X/L7Yx2j/di2gocPS4cPxzoJYA6H1MXRoi6OFtmd8PF5lBQAGCzj/Vo9F/pHOVqlXhXVne6TZigpADBYcEeleu74n69jGyUmuLsPAGAsSgoAYCxKCgBgLEoKAGAsSgoAYCxKCgBgLEoKAGAsSgoAYCxKCgBgLEoKAGAsSgoAYCxKCgBgLEoKAGAsSgoAYCwe1QEABrKcTslyyA4GpVBnfEjHVygpADBMQppbNZMHq26glL5D6vXmjq+ezN0JcbkPAAxjuVPVMqFO5bc8oYZ/rJeV2j3WkWKGkgIA07QE1HAgVY8czlfzgW6yA4FYJ4oZLvcBgGGCR2p18fI0lf/n5bq45qhCtXWxjhQzlBQAGMYOtMje/okSJHXeWya+wuU+AICxKCkAgLEoKQCAsSgpAICxKCkAgLEoKQCAsSgpAICxKCkAgLGiXlLFxcUaOXKkunfvroyMDE2aNEmVlZURY5qamlRUVKSePXuqW7dumjx5sqqrq6MdBQAQ56JeUuvXr1dRUZE2bdqktWvXKhAI6Nprr1VjY2N4zJw5c7Rq1Sq98cYbWr9+vQ4dOqQbb7wx2lEAAHHOsm3bbs9/4PDhw8rIyND69ev13e9+Vz6fTxdccIFefvll3XTTTZKkTz75RIMGDVJpaamuvPLKUx7T7/fL7XZrjCbKaSW2Z3wAQDtotQN6T2/J5/MpNTX1pOPa/TUpn88nSUpPT5cklZeXKxAIqKCgIDxm4MCBysnJUWlpaZvHaG5ult/vj1gAAB1fu5ZUKBTS7NmzNWrUKF166aWSJK/Xq6SkJKWlpUWMzczMlNfrbfM4xcXFcrvd4SU7O7s9YwMADNGuJVVUVKTt27fr1VdfPafjzJs3Tz6fL7xUVVVFKSEAwGTt9qiOmTNnavXq1dqwYYP69OkTXu/xeNTS0qK6urqIs6nq6mp5PJ42j+VyueRyudorKgDAUFE/k7JtWzNnztSKFSu0bt065ebmRmwfPny4EhMTVVJSEl5XWVmp/fv3Kz8/P9pxAABxLOpnUkVFRXr55Zf11ltvqXv37uHXmdxut1JSUuR2uzV9+nTNnTtX6enpSk1N1V133aX8/PzTurMPANB5RL2kli1bJkkaM2ZMxPrnnntOP/3pTyVJTz75pBwOhyZPnqzm5maNHTtWv/vd76IdBQAQ59r9fVLtgfdJAUB8M+Z9UgAAnC1KCgBgLEoKAGCsdnuf1PnmvKifvN/vrUA3S573G2WVfhje5ujSRY3XXqrai51K/TyktLd3KshHKwGA8TrMmVTjwAs0ZkaZHrvjWR36blfJssLbrO7ddOAHrXrjjsfUcItPVq/0GCYFAJyuDnMm5WgJaU9DhjKS+sl5/Bsbg0HZjU5tacpRY32yFGps8xgAALN0mJJK/rhK/kf7anVKX2XtqFHwb+6sD/nq9fcvNun3707WRV80K1RzJIZJAQCnq8OUVLC6Rq63a776+hvb7ECLtOkjddv01feh8xsNAHCWOsxrUgCAjoeSAgAYi5ICABiLkgIAGIuSAgAYi5ICABiLkgIAGIuSAgAYi5ICABiLkgIAGIuSAgAYi5ICABiLkgIAGIuSAgAYi5ICABiLkgIAGIuSAgAYq8M8mRdA9CT06KHjV/RXc5pT7kq/Qh/ukmw71rHQCXEmBeAEof4X6visOk3/5Urt/UGarKSkWEdCJ0VJATiB7XQoo2uDhiRXKZhsy7KsWEdCJ8XlPgAncH5erep/v0j/lDFLOR80K9QSiHUkdFKUFIATtHqrlfYf1UqLdRB0elzuAwAYi5ICABiLkgIAGIuSAgAYi5ICABiLkgIAGIuSAgAYi5ICABiLkgIAGKvdS+rRRx+VZVmaPXt2eF1TU5OKiorUs2dPdevWTZMnT1Z1dXV7RwEAxJl2LaktW7bo97//vb7zne9ErJ8zZ45WrVqlN954Q+vXr9ehQ4d04403tmcUAEAcareSamhoUGFhoZ599ln16NEjvN7n8+kPf/iDnnjiCf3DP/yDhg8frueee07vv/++Nm3a1Oaxmpub5ff7IxYAQMfXbiVVVFSk8ePHq6CgIGJ9eXm5AoFAxPqBAwcqJydHpaWlbR6ruLhYbrc7vGRnZ7dXbACAQdqlpF599VVVVFSouLj4hG1er1dJSUlKS0uLWJ+ZmSmv19vm8ebNmyefzxdeqqqq2iM2AMAwUX9UR1VVlWbNmqW1a9cqOTk5Ksd0uVxyuVxRORaAM2RZkuWQ7BCPkMd5F/UzqfLyctXU1Ojyyy+X0+mU0+nU+vXrtWTJEjmdTmVmZqqlpUV1dXUR+1VXV8vj8UQ7DoBz4UhQy7XDtX9+nmruzJezL5facX5FvaSuueYaffzxx9q2bVt4GTFihAoLC8NfJyYmqqSkJLxPZWWl9u/fr/z8/GjHAXAOrIQEHRyTqD/d+rhG/ORDtfTrFetI6GSifrmve/fuuvTSSyPWde3aVT179gyvnz59uubOnav09HSlpqbqrrvuUn5+vq688spoxwFwLuyQknyW/uQbroqaC9XreGusE6GTicnj45988kk5HA5NnjxZzc3NGjt2rH73u9/FIgqAb2EHg8pZUaN3Pxql9PpWJez+TMFYh0KnYtl2/L0S6vf75Xa7NUYT5bQSYx0HAHCGWu2A3tNb8vl8Sk1NPek4PrsPAGAsSgoAYKyYvCYFAG1xZvdRc/8MOQIhOSurFDzyZawjIcY4kwJgjOprs9X30d1qebBOx4fnxjoODEBJATBGcw9Lt2Ws1w1ZHynQLSHWcWAALvcBMIZny3H99I8z5TxmKeeTL7ndHZQUAHM4NmxTv//71RlUMEhFgZICYBLblt3Kp1rgf/GaFADAWJ3mTMrRpYssl0t2S4tCx47xyAEAiAOdoqQc3bvLO3WI6kY2q/t2l/q8uEfBw4djHQsAcAqd4nKfo2sXHftegz699ll1/75XcneLdSQAwGnoFCVlN7fI+qSbfvDp9areniHreHOsIwEATkOnuNwX9PnV/9/3K/CfabrYf0jB6ppYRwIAnIZOUVIKBdVadUCqkkKxzgIAOG2d4nIfACA+UVIAAGNRUgAAY1FSAABjUVIAAGNRUgAAY1FSAABjUVIAAGNRUgAAY1FSAABjUVIAAGNRUgAAY1FSAABjUVIAAGNRUgAAY1FSAABjUVIAAGNRUgAAY1FSAABjUVIAAGNRUgAAY1FSAABjUVIAAGO1S0kdPHhQP/7xj9WzZ0+lpKRoyJAh2rp1a3i7bdtasGCBevfurZSUFBUUFGjPnj3tEQUAEMeiXlJHjx7VqFGjlJiYqLfffls7d+7U448/rh49eoTHLFq0SEuWLNHy5ctVVlamrl27auzYsWpqaop2HABAHHNG+4ALFy5Udna2nnvuufC63Nzc8Ne2bWvx4sX6xS9+oYkTJ0qSXnzxRWVmZmrlypWaMmVKtCMBAOJU1M+k/vznP2vEiBH64Q9/qIyMDA0bNkzPPvtsePvevXvl9XpVUFAQXud2u5WXl6fS0tI2j9nc3Cy/3x+xAAA6vqiX1GeffaZly5ZpwIABeuedd3THHXfo7rvv1gsvvCBJ8nq9kqTMzMyI/TIzM8Pbvqm4uFhutzu8ZGdnRzs2AMBAUS+pUCikyy+/XI888oiGDRumGTNm6LbbbtPy5cvP+pjz5s2Tz+cLL1VVVVFMDAAwVdRLqnfv3ho8eHDEukGDBmn//v2SJI/HI0mqrq6OGFNdXR3e9k0ul0upqakRCwCg44t6SY0aNUqVlZUR63bv3q2+fftK+uomCo/Ho5KSkvB2v9+vsrIy5efnRzsOACCORf3uvjlz5uiqq67SI488optvvlmbN2/WM888o2eeeUaSZFmWZs+erYceekgDBgxQbm6u5s+fr6ysLE2aNCnacQAAcSzqJTVy5EitWLFC8+bN07/9278pNzdXixcvVmFhYXjMfffdp8bGRs2YMUN1dXUaPXq01qxZo+Tk5GjHAQDEMcu2bTvWIc6U3++X2+3WGE2U00qMdRwAwBlqtQN6T2/J5/N9630GfHYfAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFhRL6lgMKj58+crNzdXKSkp6t+/v37961/Ltu3wGNu2tWDBAvXu3VspKSkqKCjQnj17oh0FABDnol5SCxcu1LJly/T0009r165dWrhwoRYtWqSnnnoqPGbRokVasmSJli9frrKyMnXt2lVjx45VU1NTtOMAAOKYM9oHfP/99zVx4kSNHz9ektSvXz+98sor2rx5s6SvzqIWL16sX/ziF5o4caIk6cUXX1RmZqZWrlypKVOmRDsSACBORf1M6qqrrlJJSYl2794tSfrwww+1ceNGXXfddZKkvXv3yuv1qqCgILyP2+1WXl6eSktL2zxmc3Oz/H5/xAIA6Piifib1wAMPyO/3a+DAgUpISFAwGNTDDz+swsJCSZLX65UkZWZmRuyXmZkZ3vZNxcXF+tWvfhXtqAAAw0X9TOr111/XSy+9pJdfflkVFRV64YUX9Nhjj+mFF14462POmzdPPp8vvFRVVUUxMQDAVFE/k7r33nv1wAMPhF9bGjJkiPbt26fi4mJNmzZNHo9HklRdXa3evXuH96uurtZll13W5jFdLpdcLle0owIADBf1M6ljx47J4Yg8bEJCgkKhkCQpNzdXHo9HJSUl4e1+v19lZWXKz8+PdhwAQByL+pnUhAkT9PDDDysnJ0eXXHKJPvjgAz3xxBP62c9+JkmyLEuzZ8/WQw89pAEDBig3N1fz589XVlaWJk2aFO04AIA4FvWSeuqppzR//nzdeeedqqmpUVZWln7+859rwYIF4TH33XefGhsbNWPGDNXV1Wn06NFas2aNkpOTox0HABDHLPtvPwoiTvj9frndbo3RRDmtxFjHAQCcoVY7oPf0lnw+n1JTU086js/uAwAYi5ICABiLkgIAGIuSAgAYi5ICABiLkgIAGIuSAgAYi5ICABiLkgIAGIuSAgAYK+qf3Yf4k9Cjh0IXZclOcMi5v0at3upYRwIASZxJQdLxK/qr7tfNcv3msL4syI11HAAIo6Sg5jSnbsvdqF/krNLxC6LwI+FIkJWYJMvplCzr3I8HoNPich/krvTrN6/fqGCyrb7bms7tYJallu8P08ExiUqqs5SzskbByk+jExRAp0NJQaEPdyn3E9dXX7cEzu1glkPevCStuOVxvVJ3hTZ+dKVclVEICaBToqQg2bZCTed4BhU+VkjJX0qLqwtUXt1HPRtao3NcAJ0SJYXosm1lra7S59v/Xhccb1XC7s8UjHUmAHGLkkLUte6rkmNflSRRUADOCXf3AQCMRUkBAIxFSQEAjEVJAQCMRUkBAIxFSQEAjEVJAQCMRUkBAIxFSQEAjEVJAQCMRUkBAIxFSQEAjEVJAQCMRUkBAIxFSQEAjEVJAQCMRUkBAIxFSQEAjEVJAQCMRUkBAIx1xiW1YcMGTZgwQVlZWbIsSytXrozYbtu2FixYoN69eyslJUUFBQXas2dPxJja2loVFhYqNTVVaWlpmj59uhoaGs5pIgCAjueMS6qxsVFDhw7V0qVL29y+aNEiLVmyRMuXL1dZWZm6du2qsWPHqqmpKTymsLBQO3bs0Nq1a7V69Wpt2LBBM2bMOPtZAAA6JMu2bfusd7YsrVixQpMmTZL01VlUVlaW7rnnHv3zP/+zJMnn8ykzM1PPP/+8pkyZol27dmnw4MHasmWLRowYIUlas2aNrr/+eh04cEBZWVmn/Hf9fr/cbrfGaKKcVuLZxkcMOLP7qPrabDX3sOTZclyODduks/8RBBCnWu2A3tNb8vl8Sk1NPem4qL4mtXfvXnm9XhUUFITXud1u5eXlqbS0VJJUWlqqtLS0cEFJUkFBgRwOh8rKyto8bnNzs/x+f8SC+NTcP0NDbtuu5Xc+rQNXp8hKSIh1JAAGi2pJeb1eSVJmZmbE+szMzPA2r9erjIyMiO1Op1Pp6enhMd9UXFwst9sdXrKzs6MZG+eRIxDS5/Xp2nL8IiU0SXaIsygAJxcXd/fNmzdPPp8vvFRVVcU6Es6Ss7JKeuwC/WnBWOWs+lIKBWMdCYDBnNE8mMfjkSRVV1erd+/e4fXV1dW67LLLwmNqamoi9mttbVVtbW14/29yuVxyuVzRjIoYCR75UknvfKkkSdQTgFOJ6plUbm6uPB6PSkpKwuv8fr/KysqUn58vScrPz1ddXZ3Ky8vDY9atW6dQKKS8vLxoxgEAxLkzPpNqaGjQp59+Gv5+79692rZtm9LT05WTk6PZs2froYce0oABA5Sbm6v58+crKysrfAfgoEGDNG7cON12221avny5AoGAZs6cqSlTppzWnX0AgM7jjEtq69atuvrqq8Pfz507V5I0bdo0Pf/887rvvvvU2NioGTNmqK6uTqNHj9aaNWuUnJwc3uell17SzJkzdc0118jhcGjy5MlasmRJFKYDAOhIzul9UrHC+6QAIL7F5H1SAABEEyUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwVlQf1RELjuRkWSkpsltbFWo8xvOJAKADieuSshKTVHvzMNV8L6CUz5OU+x8H1Pr5/ljHAgBESVxf7rMSE3Q4L6SPxz6tvPEfqzXDHetIAIAoiuuSsoMhdft/Cbp5z41678OBSqhvinUkAEAUxfXlPrulRdmvfKbQf6dp8LHDCh78ItaRAABRFNclJdtW6xde6QuvQrHOAgCIuri+3AcA6NgoKQCAsSgpAICxKCkAgLEoKQCAsSgpAICxKCkAgLEoKQCAsSgpAICxKCkAgLEoKQCAsSgpAICxKCkAgLEoKQCAsSgpAICxKCkAgLEoKQCAsSgpAICxKCkAgLEoKQCAsSgpAICxKCkAgLEoKQCAsc64pDZs2KAJEyYoKytLlmVp5cqV4W2BQED333+/hgwZoq5duyorK0s/+clPdOjQoYhj1NbWqrCwUKmpqUpLS9P06dPV0NBwzpMBAHQsZ1xSjY2NGjp0qJYuXXrCtmPHjqmiokLz589XRUWF3nzzTVVWVuqGG26IGFdYWKgdO3Zo7dq1Wr16tTZs2KAZM2ac/SwAAB2SZdu2fdY7W5ZWrFihSZMmnXTMli1bdMUVV2jfvn3KycnRrl27NHjwYG3ZskUjRoyQJK1Zs0bXX3+9Dhw4oKysrFP+u36/X263W2M0UU4r8WzjAwBipNUO6D29JZ/Pp9TU1JOOa/fXpHw+nyzLUlpamiSptLRUaWlp4YKSpIKCAjkcDpWVlbV5jObmZvn9/ogFANDxtWtJNTU16f7779ctt9wSbkqv16uMjIyIcU6nU+np6fJ6vW0ep7i4WG63O7xkZ2e3Z2wAgCHaraQCgYBuvvlm2batZcuWndOx5s2bJ5/PF16qqqqilBIAYDJnexz064Lat2+f1q1bF3G90ePxqKamJmJ8a2uramtr5fF42jyey+WSy+Vqj6gAAINF/Uzq64Las2eP/vKXv6hnz54R2/Pz81VXV6fy8vLwunXr1ikUCikvLy/acQAAceyMz6QaGhr06aefhr/fu3evtm3bpvT0dPXu3Vs33XSTKioqtHr1agWDwfDrTOnp6UpKStKgQYM0btw43XbbbVq+fLkCgYBmzpypKVOmnNadfQCAzuOMb0F/7733dPXVV5+wftq0afrlL3+p3NzcNvd79913NWbMGElfvZl35syZWrVqlRwOhyZPnqwlS5aoW7dup5WBW9ABIL6d7i3o5/Q+qVihpAAgvhnzPikAAM4WJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADBWu3x2X2eSkJoqy50qBQIKfnlUdqAl1pEAoMOgpM6B5XTq8E2XqGmCT40Humvg8jQFd1TGOhYAdBhc7jsXlkNHB0kbR/5BPxxdppaMrrFOBAAdSqc7k0q4+O9UM/oCyZIy3v9SwZ27z/pYdjCo9B1S/qbb1Hygmy6uOapgFLMCQGfX6UqqdmQv/WjOfyvZatWLweuVvvMcDhYKqtebO3TBX7rLDtQoVFsXrZgAAHXCkrKC0pFAd3VxtMgRhdOeoN8v+f3nfiAAwAk6XUmll3n110evlCypZ3k1l+cAwGCdrqRaP/tc3T/7XJIoKAAwHHf3AQCMRUkBAIxFSQEAjEVJAQCMRUkBAIxFSQEAjEVJAQCMRUkBAIxFSQEAjEVJAQCMRUkBAIxFSQEAjEVJAQCMRUkBAIwV94/qsBKTZCU6pVBIoeZmybZjHQkAECVxXVJWYpL8ky+Xd5StrlUJynmtSq37qmIdCwAQJXF9uc9KTJB3tK3ySU9q6KSdas1Mi3UkAEAUxXVJKRRSl6oEzTowTu/v7q+EhuZYJwIARFFcX+4LNbeo76tVqn63rwY21Mv+bH+sIwEAoiiuS0q2/dVrUPuqFIx1FgBA1MX35T4AQIdGSQEAjEVJAQCMdcYltWHDBk2YMEFZWVmyLEsrV6486djbb79dlmVp8eLFEetra2tVWFio1NRUpaWlafr06WpoaDjTKACADu6MS6qxsVFDhw7V0qVLv3XcihUrtGnTJmVlZZ2wrbCwUDt27NDatWu1evVqbdiwQTNmzDjTKACADu6M7+677rrrdN11133rmIMHD+quu+7SO++8o/Hjx0ds27Vrl9asWaMtW7ZoxIgRkqSnnnpK119/vR577LE2Sw0A0DlF/TWpUCikqVOn6t5779Ull1xywvbS0lKlpaWFC0qSCgoK5HA4VFZW1uYxm5ub5ff7IxYAQMcX9ZJauHChnE6n7r777ja3e71eZWRkRKxzOp1KT0+X1+ttc5/i4mK53e7wkp2dHe3YAAADRbWkysvL9dvf/lbPP/+8LMuK2nHnzZsnn88XXqqq+BBZAOgMolpSf/3rX1VTU6OcnBw5nU45nU7t27dP99xzj/r16ydJ8ng8qqmpidivtbVVtbW18ng8bR7X5XIpNTU1YgEAdHxR/VikqVOnqqCgIGLd2LFjNXXqVN16662SpPz8fNXV1am8vFzDhw+XJK1bt06hUEh5eXnRjAMAiHNnXFINDQ369NNPw9/v3btX27ZtU3p6unJyctSzZ8+I8YmJifJ4PLr44oslSYMGDdK4ceN02223afny5QoEApo5c6amTJnCnX0AgAhnfLlv69atGjZsmIYNGyZJmjt3roYNG6YFCxac9jFeeuklDRw4UNdcc42uv/56jR49Ws8888yZRgEAdHCWbcff89b9fr/cbrfGaKKcVmKs4wAAzlCrHdB7eks+n+9b7zOI70d1SEq44AKpR6r09d2EoZB05KiCR4/GNhgA4JzFdUlZLpcO3TJAqeO/UIIjJElqaHYp4dWL5f7jphinAwCcq/guqQSHfANb9f4lr6mLI0mStL+1QeM23Se3ZUnxdyUTAPA34rqk7EBQvbYmaEjqz2U5viqkYHOC+uwMUlAA0AHEeUm1qNdrH+mCVV0i19fXKxSjTACA6InrkpKkUGOj1NgY6xgAgHbAk3kBAMaipAAAxqKkAADGoqQAAMaipAAAxqKkAADGoqQAAMaipAAAxqKkAADGistPnPj6EVitCkh8RB8AxJ1WBST979/zk4nLkqqvr5ckbdR/xTgJAOBc1NfXy+12n3R7XD6ZNxQKqbKyUoMHD1ZVVdW3PtXRdH6/X9nZ2XE9D+Zgho4wB6ljzIM5nJpt26qvr1dWVpYcjpO/8hSXZ1IOh0MXXnihJCk1NTVufwj+VkeYB3MwQ0eYg9Qx5sEcvt23nUF9jRsnAADGoqQAAMaK25JyuVx68MEH5XK5Yh3lnHSEeTAHM3SEOUgdYx7MIXri8sYJAEDnELdnUgCAjo+SAgAYi5ICABiLkgIAGIuSAgAYK25LaunSperXr5+Sk5OVl5enzZs3xzrSSRUXF2vkyJHq3r27MjIyNGnSJFVWVkaMaWpqUlFRkXr27Klu3bpp8uTJqq6ujlHiU3v00UdlWZZmz54dXhcPczh48KB+/OMfq2fPnkpJSdGQIUO0devW8HbbtrVgwQL17t1bKSkpKigo0J49e2KY+ETBYFDz589Xbm6uUlJS1L9/f/3617+O+KBO0+axYcMGTZgwQVlZWbIsSytXrozYfjp5a2trVVhYqNTUVKWlpWn69OlqaGgwYg6BQED333+/hgwZoq5duyorK0s/+clPdOjQIaPmcKp5fNPtt98uy7K0ePHiiPXncx5xWVKvvfaa5s6dqwcffFAVFRUaOnSoxo4dq5qamlhHa9P69etVVFSkTZs2ae3atQoEArr22mvV2NgYHjNnzhytWrVKb7zxhtavX69Dhw7pxhtvjGHqk9uyZYt+//vf6zvf+U7EetPncPToUY0aNUqJiYl6++23tXPnTj3++OPq0aNHeMyiRYu0ZMkSLV++XGVlZeratavGjh2rpqamGCaPtHDhQi1btkxPP/20du3apYULF2rRokV66qmnwmNMm0djY6OGDh2qpUuXtrn9dPIWFhZqx44dWrt2rVavXq0NGzZoxowZ52sK3zqHY8eOqaKiQvPnz1dFRYXefPNNVVZW6oYbbogYF+s5SKf+b/G1FStWaNOmTcrKyjph23mdhx2HrrjiCruoqCj8fTAYtLOysuzi4uIYpjp9NTU1tiR7/fr1tm3bdl1dnZ2YmGi/8cYb4TG7du2yJdmlpaWxitmm+vp6e8CAAfbatWvt733ve/asWbNs246POdx///326NGjT7o9FArZHo/H/s1vfhNeV1dXZ7tcLvuVV145HxFPy/jx4+2f/exnEetuvPFGu7Cw0LZt8+chyV6xYkX4+9PJu3PnTluSvWXLlvCYt99+27Ysyz548OB5y/61b86hLZs3b7Yl2fv27bNt27w52PbJ53HgwAH7wgsvtLdv32737dvXfvLJJ8Pbzvc84u5MqqWlReXl5SooKAivczgcKigoUGlpaQyTnT6fzydJSk9PlySVl5crEAhEzGngwIHKyckxbk5FRUUaP358RFYpPubw5z//WSNGjNAPf/hDZWRkaNiwYXr22WfD2/fu3Suv1xsxB7fbrby8PGPmIElXXXWVSkpKtHv3bknShx9+qI0bN+q6666TFD/z+Nrp5C0tLVVaWppGjBgRHlNQUCCHw6GysrLznvl0+Hw+WZaltLQ0SfEzh1AopKlTp+ree+/VJZdccsL28z2PuPsU9CNHjigYDCozMzNifWZmpj755JMYpTp9oVBIs2fP1qhRo3TppZdKkrxer5KSksI/zF/LzMyU1+uNQcq2vfrqq6qoqNCWLVtO2BYPc/jss8+0bNkyzZ07V//yL/+iLVu26O6771ZSUpKmTZsWztnWz5Ypc5CkBx54QH6/XwMHDlRCQoKCwaAefvhhFRYWSlLczONrp5PX6/UqIyMjYrvT6VR6erqRc2pqatL999+vW265JfwJ4vEyh4ULF8rpdOruu+9uc/v5nkfclVS8Kyoq0vbt27Vx48ZYRzkjVVVVmjVrltauXavk5ORYxzkroVBII0aM0COPPCJJGjZsmLZv367ly5dr2rRpMU53+l5//XW99NJLevnll3XJJZdo27Ztmj17trKysuJqHh1VIBDQzTffLNu2tWzZsljHOSPl5eX67W9/q4qKClmWFes4kuLwxolevXopISHhhLvGqqur5fF4YpTq9MycOVOrV6/Wu+++qz59+oTXezwetbS0qK6uLmK8SXMqLy9XTU2NLr/8cjmdTjmdTq1fv15LliyR0+lUZmam8XPo3bu3Bg8eHLFu0KBB2r9/vySFc5r+s3XvvffqgQce0JQpUzRkyBBNnTpVc+bMUXFxsaT4mcfXTievx+M54cao1tZW1dbWGjWnrwtq3759Wrt2bcRzmOJhDn/9619VU1OjnJyc8O/5vn37dM8996hfv36Szv884q6kkpKSNHz4cJWUlITXhUIhlZSUKD8/P4bJTs62bc2cOVMrVqzQunXrlJubG7F9+PDhSkxMjJhTZWWl9u/fb8ycrrnmGn388cfatm1beBkxYoQKCwvDX5s+h1GjRp1w6//u3bvVt29fSVJubq48Hk/EHPx+v8rKyoyZg/TVnWTffJJpQkKCQqGQpPiZx9dOJ29+fr7q6upUXl4eHrNu3TqFQiHl5eWd98xt+bqg9uzZo7/85S/q2bNnxPZ4mMPUqVP10UcfRfyeZ2Vl6d5779U777wjKQbziPqtGOfBq6++artcLvv555+3d+7cac+YMcNOS0uzvV5vrKO16Y477rDdbrf93nvv2V988UV4OXbsWHjM7bffbufk5Njr1q2zt27daufn59v5+fkxTH1qf3t3n22bP4fNmzfbTqfTfvjhh+09e/bYL730kt2lSxf7j3/8Y3jMo48+aqelpdlvvfWW/dFHH9kTJ060c3Nz7ePHj8cweaRp06bZF154ob169Wp779699ptvvmn36tXLvu+++8JjTJtHfX29/cEHH9gffPCBLcl+4okn7A8++CB859vp5B03bpw9bNgwu6yszN64caM9YMAA+5ZbbjFiDi0tLfYNN9xg9+nTx962bVvE73lzc7MxczjVPNryzbv7bPv8ziMuS8q2bfupp56yc3Jy7KSkJPuKK66wN23aFOtIJyWpzeW5554Ljzl+/Lh955132j169LC7dOli/+AHP7C/+OKL2IU+Dd8sqXiYw6pVq+xLL73Udrlc9sCBA+1nnnkmYnsoFLLnz59vZ2Zm2i6Xy77mmmvsysrKGKVtm9/vt2fNmmXn5OTYycnJ9kUXXWT/67/+a8QfQ9Pm8e6777b5OzBt2rTTzvvll1/at9xyi92tWzc7NTXVvvXWW+36+noj5rB3796T/p6/++67xszhVPNoS1sldT7nwfOkAADGirvXpAAAnQclBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAwFiUFADAWJQUAMBYlBQAw1v8HQz+oO3IS+ikAAAAASUVORK5CYII="
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "execution_count": 17
  },
  {
   "metadata": {},
   "cell_type": "markdown",
   "source": [
    "## Next steps\n",
    "In the following tutorials, we will learn more details about the constraints."
   ],
   "id": "af714ca497053de3"
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 2
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython2",
   "version": "2.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 5
}