Session 4 - Centrality, Spectral.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# EigenCentrality\n",
    "\n",
    "We have a network with N nodes in it. N by N Matrix $X$ describes strengths of the connections between nodes.\n",
    "Our goal is to define some notion of \"centrality\" or \"importance\" of a node in this network.\n",
    "\n",
    "Idea: \"central\" nodes are strongly connected with other \"central\" nodes. Denote centrality of $i$th node as $v_i$:\n",
    "\n",
    "$$ v_i = \\frac{1}{\\lambda} \\sum_{j = 1}^N{X_{ij} v_j} $$\n",
    "\n",
    "Pack all $v_i$ in a single vector $\\mathbf{v}$ and write the above line with matrix notation:\n",
    "\n",
    "$$ \\mathbf{v} = \\frac{1}{\\lambda} X \\mathbf{v} $$\n",
    "\n",
    "$$ X \\mathbf{v} = \\lambda \\mathbf{v} $$\n",
    "\n",
    "Eigen equation!\n",
    "However, we know it has many solutions. But only non-negative centralities make sense. \n",
    "\n",
    "*Perron-Frobenius Theorem*: If all elements of matrix $X$ are non-negative and there is a path from any node to any node, then there is only one eigenvector with all real non-negative components and it corresponds to the highest by absolute value real positive eigenvector.\n",
    "\n",
    "\n",
    "# Example: Slice of Six Degrees of Francis Bacon Dataset"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from mpl_toolkits import mplot3d\n",
    "import csv\n",
    "%matplotlib notebook\n",
    "plt.rcParams['figure.dpi'] = 150"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open('quakers_nodelist.csv', 'r') as nodecsv: # Open the file\n",
    "    nodereader = csv.reader(nodecsv) # Read the csv\n",
    "    # Retrieve the data \n",
    "    nodes = dict([(n[0], n[1:]) for n in nodereader][1:])\n",
    "    \n",
    "ids = dict([(name, i) for (i, name) in enumerate(nodes.keys())])\n",
    "ids_names = dict(enumerate(nodes.keys()))\n",
    "N = len(ids)\n",
    "\n",
    "with open('quakers_edgelist.csv', 'r') as edgecsv: # Open the file\n",
    "    edgereader = csv.reader(edgecsv) # Read the csv\n",
    "    edges = [tuple(e) for e in edgereader][1:] # Retrieve the data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(list(nodes.keys())[:10])\n",
    "print(edges[:10])\n",
    "print(N)\n",
    "print(len(edges))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X = np.zeros((N, N))\n",
    "\n",
    "for (p1, p2) in edges:\n",
    "    id1 = ids[p1]\n",
    "    id2 = ids[p2]\n",
    "    \n",
    "    X[id1, id2] = 1.\n",
    "    X[id2, id1] = 1.\n",
    "    \n",
    "X.max(axis=1) # Each person has at least one connection"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "eigvals, eigvecs = np.linalg.eig(X)\n",
    "sorted_idx = np.argsort(eigvals)[::-1]\n",
    "eigvals, eigvecs = np.real(eigvals), np.real(eigvecs)\n",
    "\n",
    "    \n",
    "eigvals, eigvecs = eigvals[sorted_idx], eigvecs[:,sorted_idx]\n",
    "eigvals[:10]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "eigvecs[:, 0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "epsilon = 1e-5\n",
    "off_cluster_ids = np.array(range(N))[np.abs(eigvecs[:, 0]) < epsilon]\n",
    "off_cluster_ids"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(ids_names[1], nodes[ids_names[1]])\n",
    "print(X[1])\n",
    "print(ids_names[4], nodes[ids_names[4]])\n",
    "print(X[4])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "most_popular_ids = np.argsort(eigvecs[:, 0])[::-1]\n",
    "print(ids_names[most_popular_ids[0]], nodes[ids_names[most_popular_ids[0]]])\n",
    "print(ids_names[most_popular_ids[1]], nodes[ids_names[most_popular_ids[1]]])\n",
    "print(ids_names[most_popular_ids[2]], nodes[ids_names[most_popular_ids[2]]])\n",
    "print(ids_names[most_popular_ids[3]], nodes[ids_names[most_popular_ids[3]]])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Spectral Clustering\n",
    "\n",
    "If a network is divided into disconnected (or weakly connected) clusters, we can identify them.\n",
    "\n",
    "Several definitions: degree matrix $D$ of $X$ is a diagonal matrix with $D_{ii}$ being a degree of $i$th node:\n",
    "\n",
    "$$ D_{ii} = \\sum_j{X_{ij}} $$\n",
    "\n",
    "Laplacian of $X$ is the matrix \n",
    "\n",
    "$$ L = D - X $$\n",
    "\n",
    "It has the following property, for any vector $\\mathbf{v}$:\n",
    "\n",
    "$$ \\mathbf{v}^T L \\mathbf{v} = \\frac{1}{2} \\sum_{i, j = 1}^N {X_{ij} (v_i - v_j)^2} $$\n",
    "\n",
    "Obviously, it has $\\mathbf{1}$, the vector consisting as all ones, as an eigenvector, which corresponds to zero eigenvalue.\n",
    "\n",
    "Moreover, if the network has several disconnected components, indicator vectors for each component will also be an eigenvectors corresponding to zero eigenvalue.\n",
    "\n",
    "Number of connected components is exactly the multiplicity of zeroth eigenvalue."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "with open('spiral.txt', 'r') as csv_file: # Open the file\n",
    "    data = csv.reader(csv_file, delimiter='\\t') # Read the csv\n",
    "    # Retrieve the data \n",
    "    data = np.array([d[:2] for d in data]).astype(np.float)\n",
    "\n",
    "# number of points\n",
    "N = len(data)\n",
    "    \n",
    "fig = plt.figure()\n",
    "plt.scatter(data[:, 0], data[:, 1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# X is pairwise similarity between points\n",
    "X = np.zeros((N, N))\n",
    "\n",
    "for i in range(N):\n",
    "    for j in range(N):\n",
    "        X[i, j] = np.exp(-((data[i] - data[j])**2).sum() / 2.)\n",
    "        \n",
    "D = np.diag(np.einsum(\"ii -> i\", X))\n",
    "X = D - X"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "eigvals, eigvecs = np.linalg.eig(X)\n",
    "sorted_idx = np.argsort(np.abs(eigvals))\n",
    "eigvals, eigvecs = np.real(eigvals), np.real(eigvecs)\n",
    "\n",
    "eigvals, eigvecs = eigvals[sorted_idx], eigvecs[:,sorted_idx]\n",
    "print(eigvals[:10])\n",
    "\n",
    "relevant_eigencomponents = eigvecs[:, :3]\n",
    "relevant_eigencomponents"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#from sklearn.cluster import KMeans\n",
    "\n",
    "#kmeans = KMeans(n_clusters=3).fit(relevant_eigencomponents)\n",
    "\n",
    "labels1 = np.abs(relevant_eigencomponents[:, 0]) > 0.5e-2\n",
    "labels2 = np.abs(relevant_eigencomponents[:, 1]) > 0.5e-2\n",
    "labels3 = np.abs(relevant_eigencomponents[:, 2]) > 0.5e-2\n",
    "\n",
    "\n",
    "fig = plt.figure()\n",
    "plt.scatter(data[labels1, 0], data[labels1 == 1, 1])\n",
    "plt.scatter(data[labels2, 0], data[labels2 == 1, 1])\n",
    "plt.scatter(data[labels3, 0], data[labels3 == 1, 1])"
   ]
  }
 ],
 "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.7.6"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}