class04_clustcoeff_python3.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# CSX46 - Class Session 4 - Transitivity (Clustering Coefficients)\n",
    "\n",
    "In this class session we are going to compute the local clustering coefficient of all vertices in the undirected human\n",
    "protein-protein interaction network (PPI), in two ways -- first without using `igraph`, and the using `igraph`. We'll obtain the interaction data from the Pathway Commons SIF file (in the `shared/` folder), we'll make an \"adjacency forest\" representation of the network, and we'll manually compute the local clustering coefficient of each vertex (protein) in the network using the \"enumerating neighbor pairs\" method described by Newman. Then we'll run the same algorithm using the `transitivity_local_undirected` function in `igraph`, and we'll compare the results in order to check our work.  Grad students:  you should also group vertices by their \"binned\" vertex degree k (bin size 50, total number of bins = 25) and plot the average local clustering coefficient for the vertices within a bin, against the center k value for the bin, on log-log scale (compare to Newman Fig. 8.12)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {},
   "outputs": [],
   "source": [
    "from igraph import Graph\n",
    "from igraph import summary\n",
    "import pandas\n",
    "import numpy\n",
    "import timeit\n",
    "from pympler import asizeof\n",
    "import bintrees"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 1:  load in the SIF file (refer to Class 3 exercise) into a data frame `sif_data`, using the `pandas.read_csv` function, and name the columns `species1`, `interaction_type`, and `species2`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [],
   "source": [
    "sif_data = pandas.read_csv(\"shared/pathway_commons.sif\",\n",
    "                           sep=\"\\t\", names=[\"species1\",\"interaction_type\",\"species2\"])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 2:  restrict the interactions to protein-protein undirected (\"in-complex-with\", \"interacts-with\"), by using the `isin` function and then using `[` to index rows into the data frame.  Call the returned ata frame `interac_ppi`."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [],
   "source": [
    "interaction_types_ppi = set([\"interacts-with\",\n",
    "                             \"in-complex-with\"])\n",
    "interac_ppi = sif_data[sif_data.interaction_type.isin(interaction_types_ppi)]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 3: restrict the data frame to only the unique interaction pairs of proteins (ignoring the interaction type), and call that data frame `interac_ppi_unique`. Make an igraph `Graph` object from `interac_ppi_unique` using `Graph.TupleList`, `values`, and `tolist`. Call `summary` on the `Graph` object. Refer to the notebooks for the in-class exercises in Class sessions 3 and 6."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "IGRAPH UN-- 17020 491784 -- \n",
      "+ attr: name (v)\n"
     ]
    }
   ],
   "source": [
    "for i in range(0, interac_ppi.shape[0]):\n",
    "    if interac_ppi.iat[i,0] > interac_ppi.iat[i,2]:\n",
    "        temp_name = interac_ppi.iat[i,0]\n",
    "        interac_ppi.set_value(i, 'species1', interac_ppi.iat[i,2])\n",
    "        interac_ppi.set_value(i, 'species2', temp_name)\n",
    "        \n",
    "interac_ppi_unique = interac_ppi[[\"species1\",\"species2\"]].drop_duplicates()        \n",
    "\n",
    "\n",
    "ppi_igraph = Graph.TupleList(interac_ppi_unique.values.tolist(), directed=False)\n",
    "summary(ppi_igraph)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 4: Obtain an adjacency list representation of the graph (refer to Class 2 notebooks), using the `igraph.Graph.get_adjlist` method."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {},
   "outputs": [],
   "source": [
    "ppi_adj_list = ppi_igraph.get_adjlist()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 5: You can get a count *n* of the vertices in the graph by calling `len` on `ppi_adj_list`. Call this adjacency forest, `ppi_adj_forest`. Make an \"adjacency forest\" data structure as a list of *n* `AVLTree` objects (refer to Class 2c exercise). "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {},
   "outputs": [],
   "source": [
    "def get_bst_forest(theadjlist):\n",
    "    g_adj_list = theadjlist\n",
    "    n = len(g_adj_list)\n",
    "    theforest = []\n",
    "    for i in range(0,n):        \n",
    "        itree = bintrees.AVLTree()\n",
    "        for j in g_adj_list[i]:\n",
    "            itree.insert(j,1)\n",
    "        theforest.append(itree)\n",
    "    return theforest\n",
    "\n",
    "def find_bst_forest(bst_forest, i, j):\n",
    "    return j in bst_forest[i]\n",
    "\n",
    "ppi_adj_forest = get_bst_forest(ppi_adj_list)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Step 6: Compute the local clustering coefficient (C<sub>i</sub>) values of the first 100 vertices (do timing on this operation using the `timeit` package) as a `numpy.array`; for any vertex with degree=1, it's C<sub>i</sub> value can be numpy `NaN`. You'll probably want to have an outer `for` loop for vertex ID `n` going from 0 to 99, and then an inner `for` loop iterating over neighbor vertices of vertex `n`. Store the clustering coefficients in a list `civals`. Print out how many seconds it takes to perform this calculation."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "78.19 s\n"
     ]
    }
   ],
   "source": [
    "N = len(ppi_adj_list)\n",
    "civals = numpy.zeros(100)\n",
    "civals[:] = numpy.NaN\n",
    "start_time = timeit.default_timer()\n",
    "for n in range(0, 100):\n",
    "    neighbors = ppi_adj_list[n]\n",
    "    nneighbors = len(neighbors)\n",
    "    if nneighbors > 1:\n",
    "        nctr = 0\n",
    "        for i in range(0, nneighbors):\n",
    "            for j in range(i + 1, nneighbors):\n",
    "                if neighbors[j] in ppi_adj_forest[neighbors[i]]:\n",
    "                    nctr += 1\n",
    "        civals[n] = nctr/(nneighbors*(nneighbors-1)/2)\n",
    "        \n",
    "ci_elapsed = timeit.default_timer() - start_time\n",
    "print(\"%0.2f s\" % ci_elapsed)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 7: Calculate the local clustering coefficients for the first 100 vertices using \n",
    "the method `igraph.Graph.transitivity_local_undirected` and save the results as a list `civals_igraph`. Do timing on the call to `transitivity_local_undirected`, using `vertices=` to specify the vertices for which you want to compute the local clustering coefficient."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "0.95 s\n"
     ]
    }
   ],
   "source": [
    "start_time = timeit.default_timer()\n",
    "civals_igraph = ppi_igraph.transitivity_local_undirected(vertices=list(range(0,100)))\n",
    "ci_elapsed = timeit.default_timer() - start_time\n",
    "print(\"%0.2f s\" % ci_elapsed)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Why do you think the `igraph` function `transitivity_local_undirected` is so much faster than our Python code above (which uses the adjacency forest representation of the graph)?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 8: Compare your C<sub>i</sub> values to those that you got from `igraph`, using a scatter plot where `civals` is on the horizontal axis and `civals_igraph` is on the vertical axis."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYIAAAEGCAYAAABo25JHAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3deZgU1bnH8e8LKIqCG6gIjOwiokQZwB1XRFzQq0Y0eqNRuRrRGE2McYNgjLjEJZGrQYNrIhoTFSOIouDCFZ1BQGWihjUMBEQwiiLIMO/9o4qhu5kZepbq6uX3eZ556HPqdPU7DcxvTlX1KXN3RESkcDWJuwAREYmXgkBEpMApCERECpyCQESkwCkIREQKXLO4C6ir1q1be8eOHeMuQ0Qkp8ycOfNzd29T3bacC4KOHTtSWloadxkiIjnFzBbXtE2HhkREClykQWBmg8zsEzObZ2bXVbO9yMymmtksM/vAzAZHWY+IiGwpsiAws6bAGOBEoCdwjpn1TBl2I/CMux8IDAX+N6p6RESkelHOCPoB89x9gbt/B4wHhqSMcaBV+HgnYFmE9YiISDWiDIJ2wJKEdnnYl2gkcJ6ZlQMTgSuq25GZDTOzUjMrXblyZRS1iogUrLhPFp8DPOru7YHBwBNmtkVN7j7W3YvdvbhNm2qvfhIRkXqKMgiWAh0S2u3DvkQXAc8AuPs7wHZA6whrEhGRFFEGQQnQzcw6mdm2BCeDJ6SM+RdwLICZ7UsQBDr2IyKSYNXX67n7lU/4at2GSPYf2QfK3L3CzIYDk4GmwDh3n2tmo4BSd58AXAM8ZGY/JThxfIHrBgkiIlXumvwJ90+dB8AR3dvQt+Oujf4akX6y2N0nEpwETuy7OeFxGXBYlDWIiOSi8i/WcvjtU6vaVx/fPZIQgBxcYkJEJN/94tkPeLp080WXs28+np1bbBvZ6ykIRESyxCfL13DCvW9WtW89vRc/6L935K+rIBARiZm7c8EjJbzxaXCtTPNmTZh980C237ZpRl5fQSAiEqOZi1dzxgPvVLUf+MFBnLh/24zWoCAQEYnBxkrnpN+9xcfL1wBQtGsLXrtmANs0zfznfBUEIiIZNvXjz7jw0ZKq9p8v7s+hXeP7LK2CQEQkQ9Zt2Mght73GF2uDD4b167gr44cdTJMmFmtdCgIRkQz42/vlXP3MnKr23684nF7tdoqxos0UBCIiEVqzbgP7j3ylqn3yAW35/TkHYhbvLCCRgkBEJCIPv7WAX7/0j6r21J8dRafWO8RYUfUUBCIijWzlmvX0vXVKVfvCwzoy4pT9YqyodgoCEZFGdNukf/CHNxZUtd+7/lh2b7VdjBVtnYJARKQRLFm9liPu2LxI3LWD9uHHR3WNsaL0KQhERBro6qdn87dZm++7NWfEQHbafpsYK6obBYGISD2VLfuKwb97q6p9+xn7c3bfohgrqh8FgYhIHbk75/3xXabPWwXAjs2bUXrjcWy3TWYWiWtskQaBmQ0C7iO4Q9nD7j46Zfs9wNFhswWwu7vvHGVNIiIN8d7C1Xz/D5sXifvD+X04Yb89Y6yo4SILAjNrCowBjgfKgRIzmxDelQwAd/9pwvgrgAOjqkdEpCEqNlYy6L63mPfZ1wB0brMDr1x1JM1iWCSusUU5I+gHzHP3BQBmNh4YApTVMP4cYESE9YiI1MurZSu45PHSqvbTww6mf+fdYqyocUUZBO2AJQntcqB/dQPNbG+gE/B6DduHAcMAiopy70SMiOSmdRs20vfXU1izvgKAQ7vsxp8u7p9Vy0M0hmw5WTwUeNbdN1a30d3HAmMBiouLPZOFiUhheqZ0Cdc++0FVe+KVR9Bzr1YxVhSdKINgKdAhod0+7KvOUODyCGsREUnLl99uoPevNi8Sd9r39uLeofl9+jLKICgBuplZJ4IAGAqcmzrIzHoAuwDvpG4TEcmkB9+Yz+hJH1e13/z50RTt1iLGijIjsiBw9wozGw5MJrh8dJy7zzWzUUCpu08Ihw4Fxru7DvmISCw++2od/X7zWlV72JGduX7wvjFWlFmRniNw94nAxJS+m1PaI6OsQUSkNrf8vYw/vr2wql1yw3G0adk8xooyL1tOFouIZNSiz7/hqLumVbVvGLwvlxzZOb6CYqQgEJGCc8VTs3hxzrKq9gcjB9Jqu9xZJK6xKQhEpGB8tPRLTv7921Xtu87qzZl92sdYUXZQEIhI3qusdIY+NIP3Fq4GYOcW2zDjl8fm7CJxjU1BICJ57f/mf865D71b1f7jD4s5dt89Yqwo+ygIRCQvbdhYyXF3v8HiVWsB6LFnS1668giaNsmv5SEag4JARPLOyx/9m0uffL+q/eylh1DccdcYK8puCgIRyRvffreRA295hXUbKgE4snsbHruwb94tEtfYFAQikhf+/O6/uP65D6vak686kn32bBljRblDQSAiOW3xqm8YcOe0qvaZfdpz11m94ysoBykIRCRn9bhpUtVhIIC3rj2aDrvm/yJxjU1BICI5Z/aS/3DamOlJfYtGnxRTNblPQSAiOaXjdS8ltV/96ZF020PnAhpCQSAiOWFK2QouTrhvcNGuLXjz2qNjrCh/KAhEJKu5O51+mbSaPe9efyx7tNouporyj4JARLLWkzMWc+PzH1W1j96nDY9c2C/GivJTpEFgZoOA+wjuUPawu4+uZsz3gZGAA3PcfYvbWYpIYanYWEnXGyYl9X04ciAtC3ip6ChFFgRm1hQYAxwPlAMlZjbB3csSxnQDfgkc5u5fmNnuUdUjIrnh9pc/5oFp86vaFx7WkRGn7BdjRfkvyhlBP2Ceuy8AMLPxwBCgLGHMJcAYd/8CwN0/i7AeEclia7+roOfNk5P6/nnriWzTtElMFRWOKIOgHbAkoV0O9E8Z0x3AzKYTHD4a6e4vp+7IzIYBwwCKiooiKVZE4vPjP81k4ofLq9o3n9yTHx3eKcaKCkvcJ4ubAd2Ao4D2wJtmtr+7/ydxkLuPBcYCFBcXe6aLFJFofP71eop/PSWpb+Ftg7VIXIZFGQRLgQ4J7fZhX6Jy4F133wAsNLNPCYKhJMK6RCQLnHDPm3yyYk1V+4EfHMSJ+7eNsaLCFWUQlADdzKwTQQAMBVKvCHoeOAd4xMxaExwqWhBhTSISswUrv+aY376R1KflIeIVWRC4e4WZDQcmExz/H+fuc81sFFDq7hPCbQPNrAzYCPzc3VdFVZOIxCt1eYi/XnYIffbWDWPiZu65dci9uLjYS0tLtz5QRLLGzMWrOeOBd5L6NAvILDOb6e7F1W2L+2SxiOS51FnAa9cMoEubHWOqRqqjIBCRSKTeN7jr7jsy5eoBMVYkNVEQiEijqm6RuJIbjqNNy+YxVSRboyAQkUbzyPSF/OrFzYsHDNpvTx48v0+MFUk6FAQi0mDrKzayz43JiwKUjTqBFtvqR0wu0N+SiDTIsb+dxvyV31S1Lx3QhetO7BFjRVJXCgIRqZcvvvmOA295Nalv3q0n0kyLxOUcBYGI1FnqJaFn9WnPnWf1jqkaaSgFgYikrbrlIbRIXO5TEIhIWlJnAded2INLB3SJqRppTAoCEanVjAWrGDp2RlKflofILwoCEalR6izgwfP6MKjXnjFVI1FREIjIFv46s5xr/jInqU+zgPylIBCRJKmzgBcuP4zeHXaOqRrJBAWBiABw1+RPuH/qvKQ+zQIKg4JApMBVVjqdr09eJG76dcfQbuftY6pIMi3SjwCa2SAz+8TM5pnZddVsv8DMVprZ7PDr4ijrEZFkwx4vTQqB5s2asGj0SQqBAhPZjMDMmgJjgOMJblJfYmYT3L0sZejT7j48qjpEZEvrNmykx03Ji8R9OHIgLbfbJqaKJE5RHhrqB8xz9wUAZjYeGAKkBoGIZNCht73Gsi/XVbX7dtyFv1x6aIwVSdyiDIJ2wJKEdjnQv5pxZ5jZkcCnwE/dfUk1Y0SkgVauWU/fW6ck9WmROIH4Txa/CDzl7uvN7H+Ax4BjUgeZ2TBgGEBRUVFmKxTJA6mXhJ5/8N7cclqvmKqRbBNlECwFOiS024d9Vdx9VULzYeCO6nbk7mOBsQDFxcXeuGWK5K9PV6xh4D1vJvVpkThJFWUQlADdzKwTQQAMBc5NHGBmbd3932HzVOAfEdYjUlBSZwG/OnU/fnhox3iKkaxWaxCY2SHAecARQFvgW+Aj4CXgSXf/sqbnunuFmQ0HJgNNgXHuPtfMRgGl7j4BuNLMTgUqgNXABQ3/lkQK21v/XMn5f3wvqU8fDJPamHv1R1rMbBKwDHgBKAU+A7YDugNHA6cAd4c/0DOmuLjYS0tLM/mSIjkjdRYw7oJijumxR0zVSDYxs5nuXlzdttpmBOe7++cpfV8D74dfvzWz1o1Uo4g0wO0vf8wD0+Yn9WkWIOmqMQiqCYF6jRGRaKXOAp665GAO6bJbTNVILtrqyWIz+y/gdmB3wMIvd/dWEdcmIrU47+F3eXte8u9imgVIfaRz1dAdwCnurit6RLJAxcZKut4wKanvjZ8fxd677RBTRZLr0gmCFQoBkezQ46ZJrNtQmdSnWYA0VI1BEB4SAig1s6eB54H1m7a7+98irk1EQl+t28ABI19J6vtg5EBaaZE4aQS1zQhOSXi8FhiY0HZAQSCSAakng7dt1oRPf31iTNVIPqrtqqELM1mIiCRbvOobBtw5LalPi8RJFNK5aqgzcB9wMMFM4B3gKndfGHFtIgUrdRZweNfWPHlxdYv3ijRcOieL/0xwg5nTw/ZQYDzVLyktIg0wY8Eqho6dkdSnk8EStXSCoIW7P5HQftLMfh5VQSKFKnUWcNlRXfjFoB4xVSOFJJ0gmBTeb3g8waGhs4GJZrYrgLuvjrA+kbz3dMm/+MVfP0zq0yxAMimdIPh++Of/pPQPJQiGzo1akUgBSZ0F3HN2b04/sH1M1Uih2moQuHunTBQiUkgufOQ9pn6yMqlPswCJS1o3pjGzXkBPgmWoAXD3x6MqSiRfuTudfjkxqe+xH/VjQPc2MVUkkt7loyOAowiCYCJwIvA2oCAQqYPuN07iuwotDyHZJ50ZwZlAb2CWu19oZnsAT0Zblkj++K6iku43Ji8S9/o1A+jcZseYKhJJlk4QfOvulWZWYWatCO5U1mFrTwIws0EEH0ZrCjzs7qNrGHcG8CzQ1911+zHJG6kng0GzAMk+6QRBqZntDDwEzCS4S9k7W3uSmTUl+CDa8UA5UGJmE9y9LGVcS+AnwLt1rF0kay3/ch0H3/ZaUp8WiZNstbWb1xtwm7v/B3jQzF4GWrn7B2nsux8wz90XhPsaDwwBylLG3UJw4xt9SE3ygmYBkmtqDQJ3dzObCOwfthfVYd/tgCUJ7XJSlqUws4OADu7+Um2fVjazYcAwgKKiojqUIJI5JYtWc9aDyZPl+b8ZTNMmFlNFIulJ59DQ+2bW191LGvOFzawJcDdwwdbGuvtYYCxAcXGxN2YdIo1BswDJZekEQX/gB2a2GPiGzfcsPmArz1tK8knl9mHfJi2BXsC04AgUewITzOxUnTCWXPHEjMXc9PxHSX0KAMk16QTBCfXcdwnQzcw6EQTAUODcTRvd/Uug9aa2mU0DfqYQkFyhpaIlX6QTBGvS7Evi7hVmNhyYTHD56Dh3n2tmo4BSd59Qt1JFssPVT8/mb7OWJvVpFiC5LK1zBASHeL4gOCy0M7DczFYAl7j7zJqe6O4TCT6NnNh3cw1jj0qzZpHYpM4Crj6+O1ce2y2makQaRzpB8CrwrLtPBjCzgcAZwCPA/6Ib1EgB6HvrFFauWZ/Up1mA5It0bn568KYQAHD3V4BD3H0G0DyyykSygLvT8bqXkkJg3AXFCgHJK+nMCP5tZr8guDENBDemWRF+criy5qeJ5DZdEiqFIp0gOBcYATwftqeHfU3ZfNMakbzxzfoK9hsxOalPi8RJPkvnxjSfA1fUsHle45YjEi/NAqQQ1RgEZnavu19lZi8S3JIyibufGmllIhm0YOXXHPPbN5L65owYyE7ba5E4yX+1zQieCP+8KxOFiMRFswApdDUGwabPB7j7GzWNEcllr8xdzrAnkj8Gs+A3g2miReKkwNR2aOhFgoXeXnb3DSnbOhMsFrfI3cdFWqFIBDQLENmstkNDlwBXA/ea2WpgJcHN6zsC84H73f2FyCsUaUSjXixj3PSFSX0KACl0tR0aWg5cC1xrZh2BtsC3wKfuvjYj1Yk0otRZQO8OO/PC5YfFVI1I9kjncwSbbkizKNJKRCJy9F3TWPj5N0l9mgWIbJZWEIjkqtRZwEWHd+Kmk3vGVI1IdlIQSF7SyWCR9CkIJK9UVjqdr09a+ZwHzzuIQb3axlSRSPar7fLRZ9z9+2b2IcmfLE73VpUiGaVZgEj91DYj+En458n13bmZDQLuI1ig7mF3H52y/VLgcmAj8DUwzN3L6vt6Upg+/3o9xb+ektT32jUD6KJF4kTSUlsQ7GBmh7n79MROMzsMWL61HYfLVI8BjgfKgRIzm5Dyg/7P7v5gOP5U4G5gUB2/BylgmgWINFxtN6a5F/iqmv6vwm1b0w+Y5+4L3P07gvsZDEkc4O6J+9+Baha3E6nO2//8fIsQmDNioEJApB5qmxHs4e4fpna6+4fhB8y2ph2wJKFdTjW3tTSzywk+wbwtcEx1OzKzYcAwgKKiojReWvKZZgEijau2INi5lm3bN1YB7j4GGGNm5wI3Aj+sZsxYgnWPKC4u1qyhQI144SMee2dxUt/C2wZjpkXiRBqitiAoNbNL3P2hxE4zuxiYWcNzEi0FOiS024d9NRkPPJDGfqUAaRYgEp3aguAq4Dkz+wGbf/AXExzCOT2NfZcA3cysE0EADCW4xWUVM+vm7v8MmycB/0QkQa8Rk/l6fUVSnwJApHHVtujcCuBQMzsa6BV2v+Tur6ezY3evMLPhwGSCy0fHuftcMxsFlLr7BGC4mR0HbAC+oJrDQlK4NAsQyQxzz61D7sXFxV5aWhp3GRIhBYBI4zOzme5eXN222i4fFcm41BAo3nsXhYBIxLTWkGQFzQJE4qMgkFh9V1FJ9xsnJfXdfHJPfnR4p5gqEik8CgKJjWYBItlBQSAZt2T1Wo64Y2pS39+vOJxe7XaKqSKRwqYgkIzSLEAk+ygIJCMmzFnGlU/NSuqb+6sT2KG5/gmKxE3/CyVymgWIZDcFgUTmJ+Nn8cLsZUl9CgCR7KMgkEhoFiCSOxQE0qj2vellvt2wMalPASCS3RQE0mg0CxDJTQoCaTAFgEhu06Jz0iCpIXBw510VAiI5RjMCqRfNAkTyh4JA6uTr9RX0GjE5qe/nJ+zD5Ud3jakiEWmoSIPAzAYB9xHcoexhdx+dsv1q4GKgAlgJ/MjdF2+xI8kKmgWI5KfIgsDMmgJjgOOBcqDEzCa4e1nCsFlAsbuvNbPLgDuAs6OqSepn9pL/cNqY6Ul9Lw4/nP3ba5E4kXwQ5YygHzDP3RcAmNl4YAhQFQTunrgE5QzgvAjrkXrQLEAk/0UZBO2AJQntcqB/LeMvAiZVt8HMhgHDAIqKihqrPqnF6Ekf8+Ab85P6ykadQIttdVpJJN9kxf9qMzsPKAYGVLfd3ccCYyG4eX0GSytImgWIFJYog2Ap0CGh3T7sS2JmxwE3AAPcfX2E9chW9Lt1Cp+tSf4rUACI5L8og6AE6GZmnQgCYChwbuIAMzsQ+AMwyN0/i7AW2QrNAkQKV2RB4O4VZjYcmExw+eg4d59rZqOAUnefANwJ7Aj8xcwA/uXup0ZVk2xJASAikZ4jcPeJwMSUvpsTHh8X5etL7RQCIgJZcrJYMksBICKJtOhcAXH3LUKg3c7bKwRECpxmBAVCswARqYlmBHlu5Zr1W4TAT47tphAQkSqaEeQxzQJEJB0Kgjw0Yc4yrnxqVlLf3684nF7ttEiciGxJQZBnNAsQkbpSEOSJoWPfYcaC1Ul9H98yiO22aRpTRSKSKxQEeUCzABFpCAVBDlMAiEhj0OWjOUohICKNRTOCHKMAEJHGphlBDlEIiEgUNCPIAQoAEYmSZgRZrLJyy0XimphCQEQal2YEWUqzABHJlEhnBGY2yMw+MbN5ZnZdNduPNLP3zazCzM6MspZcMX/l11uEwE+P664QEJHIRDYjMLOmwBjgeKAcKDGzCe5eljDsX8AFwM+iqiOXaBYgInGI8tBQP2Ceuy8AMLPxwBCgKgjcfVG4rTLCOrLeo9MXMvLFsqS+l686gh57toqpIhEpJFEGQTtgSUK7HOgf4evlJM0CRCRuOXGy2MyGAcMAioqKYq6mcZz0u7eYu+yrpL55t55Is6a6kEtEMivKIFgKdEhotw/76szdxwJjAYqLi73hpcVLswARySZRBkEJ0M3MOhEEwFDg3AhfL+spAEQkG0V2HMLdK4DhwGTgH8Az7j7XzEaZ2akAZtbXzMqBs4A/mNncqOqJm0JARLJVpOcI3H0iMDGl7+aExyUEh4zylgJARLKdzkxGSCEgIrkgJ64ayjUKABHJJZoRNKINGyu3CIFTe++lEBCRrKYZQSPRLEBEcpWCoIHKv1jL4bdPTer7/TkHckrvvWKqSESkbhQEDaBZgIjkAwVBPUz9+DMufLQkqe/9m45n1x22jakiEZH6UxDUkWYBIpJvFARpenLGYm58/qOkvoW3DcbMYqpIRKRxKAjSoFmAiOQzBUEtfvynmUz8cHlSnwJARPKNgqAGqbOAHnu25OWrjoypGhGR6CgIUugwkIgUGi0xkSA1BM47uEghICJ5TzMCYNC9b/Lx8jVJfQoAESkUBR0EGyudLtcn3S6Bv152KH323iWmikREMq9gg+CSx0t5tWxFUp9mASJSiCINAjMbBNwHNAUedvfRKdubA48DfYBVwNnuvijKmio2VtL1hklJfXNuHshOLbaJ8mVFRLJWZCeLzawpMAY4EegJnGNmPVOGXQR84e5dgXuA26OqB2BK2YqkEDjpgLYsGn2SQkBEClqUM4J+wDx3XwBgZuOBIUBZwpghwMjw8bPA/WZm7u6NXUzJotVc/HgpAP077cpTlxxMkyZaHkJEJMogaAcsSWiXA/1rGuPuFWb2JbAb8HniIDMbBgwDKCoqqlcxbXZsTt+OuzDilP3o1W6neu1DRCQf5cTnCNx9rLsXu3txmzZt6rWPjq134C+XHqoQEBFJEWUQLAU6JLTbh33VjjGzZsBOBCeNRUQkQ6IMghKgm5l1MrNtgaHAhJQxE4Afho/PBF6P4vyAiIjULLJzBOEx/+HAZILLR8e5+1wzGwWUuvsE4I/AE2Y2D1hNEBYiIpJBkX6OwN0nAhNT+m5OeLwOOCvKGkREpHY5cbJYRESioyAQESlwCgIRkQKnIBARKXCWa1drmtlKYHE9n96alE8tZxHVVj+qrX5UW/3kcm17u3u1n8jNuSBoCDMrdffiuOuojmqrH9VWP6qtfvK1Nh0aEhEpcAoCEZECV2hBMDbuAmqh2upHtdWPaqufvKytoM4RiIjIlgptRiAiIikUBCIiBS4vg8DMBpnZJ2Y2z8yuq2Z7czN7Otz+rpl1zKLajjSz982swszOzFRdadZ2tZmVmdkHZvaame2dRbVdamYfmtlsM3u7mvtjx1ZbwrgzzMzNLCOXH6bxnl1gZivD92y2mV2cibrSqS0c8/3w39tcM/tzttRmZvckvGefmtl/sqi2IjObamazwv+ng9Pasbvn1RfBktfzgc7AtsAcoGfKmB8DD4aPhwJPZ1FtHYEDgMeBM7PsfTsaaBE+vizL3rdWCY9PBV7OltrCcS2BN4EZQHE21AVcANyfqX9jdaytGzAL2CVs754ttaWMv4Jgif2sqI3ghPFl4eOewKJ09p2PM4J+wDx3X+Du3wHjgSEpY4YAj4WPnwWONbNM3Ml+q7W5+yJ3/wCozEA9da1tqruvDZszCO46ly21fZXQ3AHI1FUQ6fx7A7gFuB1Yl2V1xSGd2i4Bxrj7FwDu/lkW1ZboHOCpjFSWXm0OtAof7wQsS2fH+RgE7YAlCe3ysK/aMe5eAXwJ7JYltcWlrrVdBEyKtKLN0qrNzC43s/nAHcCV2VKbmR0EdHD3lzJUU1p1hc4IDyE8a2YdqtkehXRq6w50N7PpZjbDzAZlUW0AhIdGOwGvZ6AuSK+2kcB5ZlZOcC+YK9LZcT4GgUTMzM4DioE7464lkbuPcfcuwC+AG+OuB8DMmgB3A9fEXUs1XgQ6uvsBwKtsniVng2YEh4eOIvit+yEz2znWirY0FHjW3TfGXUiCc4BH3b09MJjgDpBb/Tmfj0GwFEj8zaZ92FftGDNrRjCFWpUltcUlrdrM7DjgBuBUd1+fTbUlGA+cFmlFm22ttpZAL2CamS0CDgYmZOCE8VbfM3dflfB3+DDQJ+Ka0q6N4LfdCe6+wd0XAp8SBEM21LbJUDJ3WAjSq+0i4BkAd38H2I5gMbraZeIkRya/CH6TWEAwZdt0QmW/lDGXk3yy+JlsqS1h7KNk9mRxOu/bgQQnq7pl4d9pt4THpxDcFzsraksZP43MnCxO5z1rm/D4dGBGtrxnwCDgsfBxa4JDIrtlQ23huB7AIsIP5WbR+zYJuCB8vC/BOYKt1piRbyDTXwRTok/DH1o3hH2jCH6LhSAl/wLMA94DOmdRbX0Jfhv6hmCWMjeLapsCrABmh18Tsqi2+4C5YV1Ta/thnOnaUsZmJAjSfM9uC9+zOeF71iNb3jPACA6plQEfAkOzpbawPRIYnama6vC+9QSmh3+ns4GB6exXS0yIiBS4fDxHICIidaAgEBEpcAoCEZECpyAQESlwCgIRkQKnIJC8Y2Z7mtl4M5tvZjPNbKKZdTezvczs2Rqe09bM/p7pWmtjZo9ubQVaM7vLzI7JVE2SnxQEklfCxQOfA6a5exd37wP8EtjD3Ze5e00/WK8GHspUnY3o90CNS1+LpENBIPnmaGCDuz+4qcPd57j7W2bW0cw+quF5ZwAvQ9U6/c+b2atmtsjMhof3YpgVLoC2q5l1MbP3Nz3ZzLolthP6u5rZFDObE95noosF7jSzj8J7KJwdjjUzuz9cb34KsHvCfvqY2RvhDGeymbUNv7fFwG5mtmcjvHdSoBQEkm96ATPr8gQz6xBrmacAAAIbSURBVAR84clrJ/UC/ovgk963Amvd/UDgHeC/3X0+8KWZfS8cfyHwSDW7/xPBcsq9gUOBf4f7/R7QGzgOuDP8wX46sA/Bp0P/OxyPmW1D8Jv/meEMZ1xY0ybvA4fV5XsWSdQs7gJEskBbYGVK31R3XwOsMbMvCVbqhGC5gwPCxw8DF5rZ1cDZBOvFVzGzlkA7d38OwN3Xhf2HA095sGrlCjN7gyBwjkzoX2Zmm5Y33ocgmF4Nb5vRlCBQNvkM2KsB378UOAWB5Ju5QF1v8fktwfpTiRJnB5UJ7Uo2/7/5KzCCYD36me4e1Qq2RrDm1CE1bN+O4HsQqRcdGpJ88zrQ3MyGbeowswPM7IhanvMpwS1C6yT8DX8y8ADVHBYKZxTlZnZaWEdzM2sBvAWcbWZNzawNwUzgPYJbWW7qb0twvgPgE6CNmR0S7mcbM9sv4aW6AzWd+xDZKgWB5BUPVlE8HTguvHx0LsEqm8trec43wHwz61qPl/wTwSzhlRq2nw9caWYfAP8H7ElwVdMHBCtEvg5c6+7Lw/5/Eqy4+TjB+Qg8uC3hmcDtZrZpVcnE8wddgdJ61C4CoNVHRQDM7HSgj7vX6c5mZvYzYCd3vymayrb6+qcDB8X1+pIfdI5ABHD358ysTvetNrPngC5AnB/oagb8NsbXlzygGYGISIHTOQIRkQKnIBARKXAKAhGRAqcgEBEpcAoCEZEC9/8h9nPdS7urlgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "import matplotlib.pyplot\n",
    "matplotlib.pyplot.plot(civals, civals_igraph)\n",
    "matplotlib.pyplot.xlabel(\"Ci (my code)\")\n",
    "matplotlib.pyplot.ylabel(\"Ci (igraph)\")\n",
    "matplotlib.pyplot.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 9:  scatter plot the average log(C<sub>i</sub>) vs. log(k) (i.e., local clustering coefficient vs. vertex degree) for 25 bins of vertex degree, with each bin size being 50 (so we are binning by k, and the bin centers are 50, 100, 150, 200, ...., 1250)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "scrolled": false
   },
   "outputs": [
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYoAAAEKCAYAAAAMzhLIAAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAADh0RVh0U29mdHdhcmUAbWF0cGxvdGxpYiB2ZXJzaW9uMy4xLjIsIGh0dHA6Ly9tYXRwbG90bGliLm9yZy8li6FKAAAgAElEQVR4nO3dd3ib5b038O8tyZK8JO+peI8MJ84iicNIQigklIQRyii0h1PKaFmnpT2l1+nb6z2jV3vO+3LeLnYodABlhUOgaUkhJCQhIXuQacd7xVuesiX5fv+Q7HjIsuxYejS+n+vy5UR6Hj0/Q6Sv7/Hct5BSgoiIaCIqpQsgIiL/xqAgIiK3GBREROQWg4KIiNxiUBARkVsMCiIickujdAHekJCQILOyspQug4goYBw+fLhFSpno6rmgDIqsrCwcOnRI6TKIiAKGEKJqoufY9URERG4xKIiIyC0GBRERucWgICIitxgURETkFoOCiIjcYlBQQKlp68X5i11Kl0EUUoLyPgoKLuY+K7adbMCWI7U4WNkOALg6PwGPrsnD8px4hasjCn4MCvJLNvsgdpe24N0jtdh++iIGbIPISYzED28ohFolsHl3Oe58cT+WZcXhsbV5uCovAUIIpcsmCkoMCvIrp+s7seVILf7nWD1auvsRExGGu66YhU2LTVhgMg6HwX0rs/DnA9V4flc5vvHyARTPisGja/Jw3ZwkBgbRDBPBuBXq0qVLJZfwCBxNXRZsPVaPd4/U4UxDJ8LUAmsKk3DbYhOunZ0ErWbiobR+mx1bjtTh2Z1lqGnrw+yUaDx6bR7WF6VCrWJgEHlKCHFYSrnU5XMMClKCxWrH309fxJYjtfistAX2QYlikxG3LTZhQ3Ea4iK1U3o9m30QW4/X45lPy3ChuQc5iZF4ZHUebl6YBo2aczaIJsOgIL8gpcShqnZsOVKLD080oMtiQ4pBj1sXp2PT4nTkJUVf9jXsgxJ/+7IRv9lRirONXZgVF47vrMrDpiXp0GnUM/BTEAUnBgUpqrq1F1uO1mLLkTpUt/UiPEyN9UUp2LTEhBU58V7pIpJS4pMzTfjNp2U4XtOBFIMeD63KwV1XZCBcy8AgGitkgkIIsQHAhry8vAdKS0uVLiekdVqs2HaiAVuO1OFAZRuEAEpy4nHbYhPWF6UgUuebeRRSSuwpa8FvdpThQEUbEqK0+PbVObh3RSaifFQDUSAImaAYwhaFMmz2Qewua8GWI3XYfqoR/c4prZsWm3DLonSkx4QrWt8X5a347adl2F3aAmN4GL51ZTbuW5kFY0SYonUR+QMGBXnVmYZLU1qbuxxTWjcWp+G2xSYUj5jS6i+O1XTgtzvK8PGZi4jSafDNkkzcf1U24qN0SpdGpBgGBc245q5+vH+sDluO1OF0Qyc0KoE1s5OwabEJa2YnBsTA8en6TjyzswzbTjZAr1Hj68sz8OA1OUg26JUujcjnGBQ0IyxWOz4+cxFbjtRh1/lm2AclFpiM2DTNKa3+oqypG8/uLMP7x+qhVgncsdSEh1flwhQboXRpRD7DoKBpk1LicFU73j1Shw9P1I+a0nrbonTkJ1/+lFZ/Ud3ai+d2XcA7h2sgJXDronR8d00eshMilS6NyOsYFDRlNW292HKkDluO1qKq9dKU1tsWm1CS650prf6iwdyHF3aV440D1bDaB7GhOA2PXZs3I/d5EPkrBgUNk1Ki02LDxU4LGswWXDRb0Dj05xHf23oGADimtG5aYsK6opSQm07a3NWPzbvL8cf9Veiz2vHV+al47Np8FKYwMCj4MChChH1QoqW7H43mSx/4jZ0WNJodX0NB0Ge1jzs3LlKLFIMeKUY9kg16ZCdE4Mb5qeynB9DWM4DNu8vx+88r0TNgx43zU/DYtfmYk2pQujSiGcOgCAIWq93xgT/0wT8iABo7HSHQ1NUP++Do/58alUCyMwBSDI4QSDXqkez8e6pRjySDLiBmKSmtvWcAv9tbgVf3VqKr34Yb5iXjsWvzUZRuVLo0osvGoPBjUkp09FodH/wufvsfahV09FrHnRul0yDZoEOqMdwZBjqkGMMdLQNnOMRHaqEK4vEEJZh7rfjd3gr8bm8Fuiw2XDcnCY+vzccCU4zSpRFNG4NCITb7IJq7+0eNBTSOGRNoNFvQbxscdZ4QQHykzvHBbwh3ftdfCgGjDskGPaL1vKNYSZ0WK36/txKb91TA3GfFmsJEPL42H4syYpUujWjKGBRe0DtgG9X146o10NLdjzE9QdCqVUge9cGvc3YHhQ8HQFK03u0eDORfuixW/GFfFTbvLkd7rxXXFCTiibV5WJIZp3RpRB5jUEyBlBJtPQOjfvu/6BwcHhoLaDBb0GWxjTvXoNcMDwanDo0JGJ1jAs7uoLhIrd8taUEzo7vfhj/tr8JLn5WjtWcAV+bF44m1BViWzcAg/8eg8NDtz32OE7VmDNhHdwWpBJAYrXM5GJwy4nuENrSmj5JrvQM2vLa/Gi98dgEt3QNYkROHJ9YWYEVOHH9JIL/FoPDQf//9PPpt9lGDwSlGPRKjdNwljaasb8CO1w9U4/ldF9Dc1Y9lWXF44rp8rMyNZ2CQ32FQECnIYrXjzYM1eG7nBTR2WrAkMxYPXZOD1YXu9wMn8iUGBZEfsFjtePtQDZ7deQENZgsMeg3WF6Vi48I0r+30R+QpBgWRH7HaB7GntAVbj9dj+6lG9AzYkRClw00LUrGhOA2LM2LYNUU+x6Ag8lMWqx07zjZh67F67DjXhAHbINJjwrGhOA0bi9MwJzWaoUE+waAgCgBdFiu2n7qIrcfrsaesBfZBibykKGwsTsOG4jQud05exaAgCjBtPQPYdrIBW4/X42BlG6QE5qcbsbE4DTcVpyLVqOz+4xR8GBREAazB3Ie/nHCExolaMwBgWVYcbl6UhjuWzkIYp27TDGBQEAWJypYefHC8HluP16O0qRvXzUnCb7++GPowrv5Ll8ddUPBXEaIAkpUQicfW5mP7967Bv99ShE/ONuG+Vw6gu3/8kjJEM4VBQRSAhBD4xopM/PLOhThY2Y6vv7R/eFdCopnGoCAKYDcvTMeL31iCc41duPOFfWg0W5QuiYIQg4IowK2dk4zff2sZGswW3P7856hq7VG6JAoyDAqiILAiJx6vP7AcPf023P78Ppxt7FS6JAoiDAqiILHAFIO3HiqBWgjc+cJ+HKluV7okChIMCqIgkp8cjbcfLkFMRBju3fwF9pS2KF0SBQEGBVGQmRUXgbcfLkFGXAS+9epB/O3LRqVLogDHoCAKQknRerz5YAnmpRvw3dcO453DtUqXRAGMQUEUpIwRYfjT/cuxMjcBP3j7OF7ZW6F0SRSgGBREQSxSp8HL9y3Funkp+NcPTuOXH59HMC7bQ94VVEEhhNgghHjRbDYrXQqR39Bp1Pjt1xfh9iUm/PLjUvzbh6cxOMiwIM8FVVBIKT+QUj5oNBqVLoXIr2jUKvzXpgX41pXZeGVvJZ7acoItC/KYRukCiMg3VCqB/3XTHOjCVHhu5wWsL0rFmtlJSpdFASCoWhRE5J4QAt//SgGy4iPwn387Czu7oMgDDAqiEBOmVuHJ6wtxtrEL7x+rU7ocCgAMCqIQ9NX5qZifbsTT28/DYrUrXQ75OQYFUQhSqQSeWj8bdR19+NP+KqXLIT/HoCAKUVfmJeDq/AQ882kZOi1WpcshP8agIAphP1o3G+29Vry4q1zpUsiPMSiIQlhRuhEbi9OweU85mjq5Ox65xqAgCnFPXl8Am13iV5+UKl0K+SkGBVGIy4yPxD3LM/DngzUob+5WuhzyQwwKIsJja/Oh16jwf7efU7oU8kMMCiJCQpQOD1yTg20nG3GspkPpcsjPMCiICADw7atzEB+pxS/+eoYLBtIoDAoiAgBE6TR4fG0+9pe3Ydf5ZqXLIT/CoCCiYXcvy0BGXAR+8dez3LOChjEoiGiYVqPCk9cXOBYMPM4FA8mBQUFEo2xYkIZ5aQY8vf08+m1cMJAYFEQ0xtCCgbXtfXhtf7XS5ZAfYFAQ0ThX5yfiyrx4/GZHKRcMJAYFEbk2tGDgS59xwcBQx6AgIpcWmGJw04JUbN5dgaYuLhgYyhgURDShH1xfCKt9EL/mgoEhjUFBRBPKSojE3csy8MaBGlS09ChdDimEQUFEbj22Ng86jQo/evcEatt7lS6HFMCgICK3kqL1+N8b5+FEbQfWPr0LT28/h55+m9JlkQ8xKIhoUncsnYVPnlyNdUUp+M2OMlz79E68e7iWy3yECAYFEXkkPSYcv7prEd79zkqkGMPx5NvHccuze3Gosk3p0sjLGBRENCVLMmPx3ndW4v/dWYymzn7c/vw+PPbGUY5fBDEGBRFNmUolcOsiE3b8YBUeX5uP7acaOX4RxBgURDRtEVoNvv+VAuz4wWrcMO/S+MX7x7jybDBhUBDRZUuPCcev716Ed79TghSDHk/8+RiOVLcrXRbNEAYFEc2YJZlxeP2BFTCGh+HZTy8oXQ7NEAYFEc2oSJ0G/3hlFj4+cxFnGzuVLodmAIOCiGbcfSuzEKlV47mdl9+qqO/og5S8X0NJDAoimnExEVrcuyITHxyvR+VlrBG1u7QZK3+xAzvPNc9gdTRVDAoi8or7r8qGRq3CC59Nr1UxOCjx821nAQA7zzXNZGk0RQwKIvKKJIMedyw14Z3DtWg0T30/iw9O1ON0Qyei9RrsK2/1QoXkKQYFEXnNQ9fkYlACL+2e2i55A7ZBPL39POakGvDwqlycv9iNlu5+L1VJk2FQEJHXzIqLwM0L0/D6F9Vo6xnw+Lw3DlSjuq0XP1pXiJW58QCA/WxVKIZBQURe9d3VubDY7Hhlb4VHx3f32/DrT0pRkhOPVQWJmJ9uRKRWjX0XGBRKYVAQkVflJUXjhrkpePXzSnRZrJMev3l3OVp7BvCj9bMhhIBGrcIV2XFsUSiIQUFEXvfImjx0WWz40/5qt8e1dPfjpc/KceP8FCycFTP8eElOPC4096Cpc+qD4nT5GBRE5HXzTUZcU5CIl/eUw2K1T3jcb3eUwWIbxA+uLxz1eIlznIKzn5TBoCAin3hkdS5augfw5sEal89Xt/bitS+qcOcVs5CTGDXquXlpRkTrNex+UgiDgoh8Yll2HJZmxuKFXRcwYBsc9/zTfz8HtUrgibX5455TqwSWZ8dxQFshUwoKIcRSIYTWW8UQUfASQuCRa/NQb7aM26/iyzoz3j9Wj/uvykayQe/y/BU58ahs7UWDuc8X5dIIHgeFECIVwOcAvua9cogomK0uSMTcVAOe23kB9sFLC/3910fnEBMRhodW5U547vA4BVsVPjeVFsU/APg9gG97qRYiCnJCCDyyJg/lLT3425eNAIDPy1rw2flmPLomDwZ92ITnzkkxwBgexnEKBUwlKL4B4McAdEKIiWOfiMiNdUUpyEmMxDOflmFwUOIXfzuLNKMe967IdHueamicgkHhcx4FhRBiDYCzUsoWAK8CuN+bRRFR8FKrBL6zKhenGzrx4y0ncaLWjO9fXwh9mHrSc0ty41HT1ofa9l4fVEpDPG1RfAvAy84//xnA14QQnDFFRNNyy6J0pMeE481DNShMjsati9I9Oo/jFMqY9MNeCBEDoATAXwFAStkJYD+AG71bGhEFqzC1Cg+vygEA/PO6QqhVwqPzCpKiERepZfeTj2kmO0BK2QEgb8xj3/BaRUQUEu5ZnonFmbGYl2b0+ByVSmBFThz2X2iFlBJCeBYwSjlQ0YbClGgYwycepA8EblsUQgiD83ucqy/flEhEwUilElMKiSElOfGoN1tQ0+bf91N0Wqy468V9eP0L9+tbBYLJWhSvA7gJwGEAEsDI+JYAcrxU1zAhRA6AfwFglFLe7u3rEZF/u7TuUwsy4jMUrmZiVS29GJRAU1fgL2TotkUhpbzJ+T1bSpnj/D70NWlICCF+J4RoEkJ8OebxdUKIc0KIMiHEU5PUUC6l5CwrIgIA5CZGISFK5/cD2pWtPQCA9ils2OSvJut6ukEIMe63eCHEJiHEVzx4/VcBrBtzrhrAMwDWA5gL4G4hxFwhxHwhxIdjvpI8/kmIKCQI4Rin2FfuGKfwV5UtjqBoDfagAPBTALtcPL4LwL9N9uJSys8AtI15eBmAMmdLYQCO6bY3SylPSilvGvPV5MHPAAAQQjwohDgkhDjU3Nzs6WlEFIBKcuNxsbMfFc4P44k0dVnwL++dVGQfi8pWx70e7b3BHxQ6KeW4T13njXeR07xmOoCR6wzXOh9zSQgRL4R4HsAiIcSPJzpOSvmilHKplHJpYmLiNEsjokBQkjP5/hT9Njse/uNhvPZFNf6wr8pXpQ2rGu56mnxXP383WVAYhBDjBryFEGEAwr1T0mhSylYp5cNSylwp5c99cU0i8m/ZCZFINkw8TiGlxE/e+xJHqjuQZtRj6/F6n3dTDbUo2kKg62kLgJeEEMOtByFEFIDnnc9NRx2AWSP+bnI+RkTkESEESnLisb+8zWUAvLK3Em8frsXja/PxT9cVoLqtF8drzT6rr7vfhpbufhj0GvRZ7egbmHhXv0AwWVD8BMBFAFVCiMNCiMMAKgA0O5+bjoMA8oUQ2c69Le4CsHWar0VEIaokNx4t3f240Nw96vE9pS342bYzuH5uMv5pbT5uKEqBVq3C1mP1PqttqNtpUUYsgMAfp5hseqxNSvkUHC2A+5xfGVLKp6SUk3a8CSHeALAPQKEQolYIcb+U0gbgUQAfATgD4C0p5anL+zGIKNSU5CQAGL3uU1VrDx55/QhyEyPx33cuhEolYAwPw6rCRHx4on7UHhjTUdbUjVue2YvW7n63x1U5u50WZcQACPzuJ0/WeooAUOCclXRSStknhMgQQky6ipeU8m4pZaqUMkxKaZJSvux8fJuUssA57vCzGfg5iCjEzIoLR5pRPzyg3d1vwwN/OAQhgM3fvAJRukvDqxuL09DU1Y8vKi7v3ottJxtwrKYD+8vHTuYcbWg21lCLIuiDAoAVwJaR4xQANgNI9U5JRESTE0JgRa5jnMI+KPG9N4/hQnMPnvn6YmTER4w69ro5yYjQqvHB8cvrfjpY6QiIk3XuxzuqWnuQGK2DKdYx5yeou54AwNnF9B6AOwBACJEBIFFKecjLtRERuVWSE4+2ngE8/sZR/P30Rfzkq3NwZV7CuOPCtWp8ZW4ytp1sxIBtcFrXstkHcaSqHYBjj293Klt7kRUfgbgILYDQaFEAjhbEPzr//E0Ar3innMsjhNgghHjRbPbd7AYiUs7Quk9/OdmAO5aacN/KrAmP3VicBnOfFXvKpndD7umGTvQM2BEfqcXJOrPb6bZVrT3IjI+EMTwMKhH4y3h4FBRSyrMAhBCiAI5ZSn/0alXTJKX8QEr5oNE49RUpiSjwmGIjkJ8UhaWZsfj3W4rcLjt+dX4ijOFh0579dKDC0e10z/IMmPusqG13vXpt74ANFzv7kRUfAZVKIDZCG/DLeExll7qX4WhZnJRStnupHiKiKXn3uyvx2gPLodO430pVq1Hhxvkp2H764rTuazhQ0YZZceFYOycZwMTdT9VtjhlPmfGOYd3YSG3wj1GM8BaAYlzaEpWISHEGfdikITFkQ3Eaegfs+OTsxSldQ0qJQ1XtWJYVj8KUaGhUYsIB7coWR1BkOYMiLkIbMmMUkFL2SimNUsqPvVkQEZG3LM+OR1K0bsrdTxeau9HWM4Bl2bHQh6lRkBw9cVA4b7bLTHDMvIqNDAv49Z6m0qIgIgpoapXATQvSsPNcM8x9nn94H6hw9LZfkeXY2HN+uhFfTjCgXdXag/hILQx6x/ancZFatIVQ1xMRUcDbuDANA/ZBfHSq0eNzDla2ISFKi+wER3dSkcmI9l4r6jrGD2hXtvQic8R9HLERWrT3DPj13hmTYVAQUUgpNhmRGR8xpZvvDlS04YqsuOFZVUVpBgCuB7SrWnuGxycAR4vCNijRabFdZuXKYVAQUUgRQmDDgjTsLWtBc5f7NZsAoL6jD3UdfcPdTgAwJ9UAtUrgy7rOUcdarHbUmy3DM54AR1AAgX0vBYOCiELOxoVpGJSOtZsmM7Rsx7LsS0GhD1MjPylq3IB2jXNqbFbCiK4nZ1AE8jgFg4KIQk5BcjQKk6M96n76oqINUToN5qQaRj3uakB7aLOiUS2KCLYo/AqX8CAiT21cmIZDVe0uB6RHOljRhiWZsVCrRt/1Pd9kRGvPABrMl/bjrnSuGpvtouspkO+lCKqg4BIeROSpjcVpEAJ4YdeFCY9p7xlAaVP3qG6nIfPSHJ8zI7ufKlt7EBMRBmNE2PBjQ11PgXx3dlAFBRGRp2bFReAfSrLwx/1VOFzlen+JofGJkQPZQ+amGqASo2c+VbX2jup2AoBIrRpajSqg13tiUBBRyPrhDYVIM4bjR++eRL9t/PpPByvboFWrsMA0vpciXKtGflL0qKCobO1B1pi9MIQQiHPeSxGoGBREFLIidRr8x61FKGvqxrOfju+COlDZjuJZRujDXK8lVZRuxMm6Tkgp0W+zo76jb1yLAnB0P7UF8DIeDAoiCmlrCpNwy8I0PLuzDOcvdg0/3jtgw6k6s8tupyFF6Qa0dPfjYmc/atv7MCgxrkUBAHGRYRyjICIKZD/dMA/R+jD88zsnYB90THc9Wt0B26DEFS4GsofMT780oF01tBigqxYFu56IiAJbXKQWP71pLo7VdOAP+yoBOJbtEAJYkhk74Xlz0xwD2ifrzMPLiw+tBzX29TmYTUQU4G5emIbVhYn4Px+dQ217Lw5UtGFuqmF4FVhXIrQa5CZG4cs6MypbexCt1yA2YvzxcZFamPussNmnt1+30hgURERwzE76j1uKAAA/3nISR2va3Y5PDJmfbnS0KFp7kRUf6XI71qGb7jqmsLS5P2FQEBE5mWIj8MMbCrG7tAUW66DLG+3GmpduRHNXP45Wt49aXnyk2ABfxiOogoJLeBDR5fpmSRYWZcQAAJZmTTw+MWRoQLvLYhu1vPhIgb6MR1AFBZfwIKLLpVYJPHvPYjx7z2IkResnPX5emgFDvU2TtigCdIpsUAUFEdFMSDWG48b5qR4dG6nTIMc50ynLxYwn4FKLIlBnPjEoiIgu01D300RdT7GRjplQgTpGoVG6ACKiQHfLonRY7RIJUVqXz+s0akTpNAG7jAeDgojoMq0uTMLqwiS3x8QG8DIe7HoiIvKBuAgtZz0REdHEHCvIMiiIiGgCbFEQEZFbcZFajlEQEdHEYiO16B2ww2Idv5Oev2NQEBH5wNBNd4HYqmBQEBH5wNAyHq7GKb735jG8f6zO1yV5LKiCgosCEpG/mmhhwOrWXrx3tA6/+qQUUkolSptUUAUFFwUkIn8V51zGY2xQ7C5rBgCUN/fgaE2Hz+vyRFAFBRGRv5poT4o9pS1IitYhPEyNtw/VKlHapBgUREQ+EBOhhRBAW++l9Z7sgxKfX2jF6sJErC9KwYfH6/1yVhSDgojIB9QqgZjwsFEtipN1Zpj7rLgqPxG3LzGhq9+Gj041KlilawwKIiIfGbuMx+7zjvGJK3PjsSInHukx4XjnsP91PzEoiIh8ZOwyHrvLWlCUbkB8lA4qlcCmJSbsKWtBXUefglWOx6AgIvKR2BHLePT023C0uh1X5SUOP3/7YhOkBN474l+tCgYFEZGPjGxRfFHRCqtd4ur8hOHnM+IjsDw7Du8crvWreyoYFEREPhIX5WhRSCmxu7QFOo0KSzJjRx3ztaWzUNnai0NV7QpVOR6DgojIR+IitLDaJbr7bdhd2oJl2XHQh6lHHbO+KAURWjXePlSjUJXjMSiIiHwk1rmMx+n6TpQ1deOa/MRxx0TqNPjq/FT85UQDegdsvi7RJQYFEZGPDC3jsfV4PQDgqhHjEyPdttiEngE7dp5r9llt7jAoiIh8ZGgZj20nG5AQpcPslGiXx12RFQuDXoNPzzb5srwJMSiIiHzk0p4UVlyVFw8hhMvjNGoVrilIxM7zzRgcVH72E4OCiMhHhoICAK52MT4x0prCJDR39eNUfae3y5pUUAUF96MgIn8WpdMgTO1oRUw0PjFkVWEihAA+Pad891NQBQX3oyAifyaEQGyEFgXJUUg26N0emxClwwJTDIOCiCjU3L0sAw9dk+vRsWsKE3GspsPl9qm+xKAgIvKh732lAJuWmDw6dk1hEqQEdp1XtlXBoCAi8lPz041IiNLi07PK3k/BoCAi8lMqlcCqgiTsOt8Mu4LTZBkURER+bM3sRJj7rDhardwigQwKIiI/dnV+ItQqoejsJwYFEZEfM4aHYUlmrKLjFAwKIiI/t6YwCacbOtFotihyfQYFEZGfu3Z2EgBgp0LdTxpFrkpERB4rSI6CKTYcP9t2BqfqO3H7EhMWmIwTLio409iiICLyc0IIvPTNpVg7OwlvHarBzc/sxQ2//Gx4XwtvY1AQEQWAOakG/PKuRTj4k+vw89vmQ0rgybeOYcA26PVrMyiIiAKIQR+Gu5dl4NFr82C1S1S09Hj9mgwKIqIAVJDs2B3v3MUur1+LQUFEFIByEiOhVgmcb2RQEBGRCzqNGtkJkWxREBHRxAqSo1DKoJgaboVKRKGkIDkaVW296Buwe/U6QRUU3AqViEJJYXI0pATKmrq9ep2gCgoiolCS75z5dN7L3U8MCiKiAJUVHwGtWsWgICIi1zRqFXKTorw+84lBQUQUwBwznzhGQUREEyhIjkZdRx+6LFavXYNBQUQUwAqHB7S916pgUBARBbACH8x8YlAQEQUwU2w4wsPUDAoiInJNpRIoSI5iUBAR0cQKkqNxrpFjFERENIGC5Gi0dPejrWfAK6/PoCAiCnAFKd4d0GZQEBEFuEIvz3xiUBARBbhkgw7Reg2DgoiIXBNCoDA5Gue9NKDNoCAiCgIFKdE4d7ELUsoZf20GBRFREJibakBClBadFtuMv7Zmxl+RiIh87t4Vmbh3RaZXXpstCiIicotBQUREbjEoiIjILQYFERG5xaAgIiK3gioohBAbhBAvms1mpUshIgoaQRUUUsoPpJQPGo1GpUshIgoaQRUUREQ084Q3bvdWmhCiGUDVNE83Arjcvgg4FHQAAAMKSURBVKvLeY3pnDuVczw9NgFAyxTrCFYz8W/CW3xZmzeuNVOv6cv3nLeOV/o9lymlTHT5jJSSXyO+ALyo5GtM59ypnOPpsQAOKf3/wl++ZuLfRDDU5o1rzdRr+vI9563j/fk9x66n8T5Q+DWmc+5UzpmJny/U+PN/M1/W5o1rzdRr+vI95+3j/U5Qdj3R5RNCHJJSLlW6DqJQ4c/vObYoaCIvKl0AUYjx2/ccWxREROQWWxREROQWg4KIiNxiUBARkVvc4Y48IoS4BcBXARgAvCyl3K5wSURBSwgxB8ATcNyE94mU8jkl62GLIoQJIX4nhGgSQnw55vF1QohzQogyIcRTACCl/B8p5QMAHgZwpxL1EgWyKb7fzkgpHwZwB4Arlah3JAZFaHsVwLqRDwgh1ACeAbAewFwAdwsh5o445CfO54loal7FFN5vQoiNAP4CYJtvyxyPQRHCpJSfAWgb8/AyAGVSynIp5QCAPwO4WTj8J4C/SimP+LpWokA3lfeb8/itUsr1AO7xbaXjcYyCxkoHUDPi77UAlgN4DMB1AIxCiDwp5fNKFEcUZFy+34QQqwHcBkAHP2hRMCjII1LKXwP4tdJ1EIUCKeVOADsVLmMYu55orDoAs0b83eR8jIhmXkC83xgUNNZBAPlCiGwhhBbAXQC2KlwTUbAKiPcbgyKECSHeALAPQKEQolYIcb+U0gbgUQAfATgD4C0p5Skl6yQKBoH8fuOigERE5BZbFERE5BaDgoiI3GJQEBGRWwwKIiJyi0FBRERuMSiIiMgtBgWRDwghssYuL00UKBgURETkFoOCyMeEEDlCiKNCiCuUroXIE1w9lsiHhBCFcOw5cJ+U8rjS9RB5gkFB5DuJAN4HcJuU8rTSxRB5il1PRL5jBlAN4CqlCyGaCrYoiHxnAMCtAD4SQnRLKV9XuiAiTzAoiHxIStkjhLgJwN+dYeF3ew8QjcVlxomIyC2OURARkVsMCiIicotBQUREbjEoiIjILQYFERG5xaAgIiK3GBREROQWg4KIiNz6/++BxKda2fqWAAAAAElFTkSuQmCC\n",
      "text/plain": [
       "<Figure size 432x288 with 1 Axes>"
      ]
     },
     "metadata": {
      "needs_background": "light"
     },
     "output_type": "display_data"
    }
   ],
   "source": [
    "civals_igraph = numpy.array(ppi_igraph.transitivity_local_undirected())\n",
    "deg_igraph = ppi_igraph.degree()\n",
    "deg_npa = numpy.array(deg_igraph)\n",
    "deg_binids = numpy.rint(deg_npa/50)\n",
    "binkvals = 50*numpy.array(range(0,25))\n",
    "civals_avg = numpy.zeros(25)\n",
    "for i in range(0,25):\n",
    "    civals_avg[i] = numpy.mean(civals_igraph[deg_binids == i])\n",
    " \n",
    "matplotlib.pyplot.loglog(\n",
    "     binkvals,\n",
    "     civals_avg)\n",
    "matplotlib.pyplot.ylabel(\"<Ci>\")\n",
    "matplotlib.pyplot.xlabel(\"k\")\n",
    "matplotlib.pyplot.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Compare the above plot to the same type of plot (from Figure 2 of the Potapov et al. article, left panel) for the human gene regulatory network \"TN_RN\". Do they have similar asymptotic (large-*k*) behavior?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Step 10: Now try computing the local clustering coefficient for *all* vertices, this time using a \"list of hashtables\" approach. Estimate how long the same calculation would take using an `AVLTree` adjacency forest, and compare. Which is faster?  Why do you think that is so?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "49.74\n"
     ]
    }
   ],
   "source": [
    "civals = numpy.zeros(len(ppi_adj_list))\n",
    "civals[:] = numpy.NaN\n",
    "\n",
    "ppi_adj_hash = []\n",
    "for i in range(0, len(ppi_adj_list)):\n",
    "    newhash = {}\n",
    "    for j in ppi_adj_list[i]:\n",
    "        newhash[j] = True\n",
    "    ppi_adj_hash.append(newhash)\n",
    "\n",
    "    \n",
    "start_time = timeit.default_timer()\n",
    "for n in range(0, len(ppi_adj_list)):\n",
    "    neighbors = ppi_adj_hash[n]\n",
    "    nneighbors = len(neighbors)\n",
    "    if nneighbors > 1:\n",
    "        nctr = 0\n",
    "        for i in neighbors:\n",
    "            for j in neighbors:\n",
    "                if (j > i) and (j in ppi_adj_hash[i]):\n",
    "                    nctr += 1\n",
    "        civals[n] = nctr/(nneighbors*(nneighbors-1)/2)\n",
    "        \n",
    "ci_elapsed = timeit.default_timer() - start_time\n",
    "print(\"%0.2f\" % ci_elapsed)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So the built-in python dictionary type gave us fantastic performance.  But is this coming at the cost of huge memory footprint? Let's check the size of our adjacency \"list of hashtables\", in MB:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {},
   "outputs": [
    {
     "data": {
      "text/plain": [
       "68.834888"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "asizeof.asizeof(ppi_adj_hash)/1000000"
   ]
  }
 ],
 "metadata": {
  "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.6.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}