Community-Detection.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Table of Contents\n",
    "* [Graph Community Detection](#Graph-Community-Detection)\n",
    "\t* [Setting](#Setting)\n",
    "\t* [Motivation](#Motivation)\n",
    "\t* [Proto-Algorithm](#Proto-Algorithm)\n",
    "\t* [Unions of $V_p$'s](#Unions-of-$V_p$'s)\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Graph Community Detection"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Setting"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let m be the number of labeled nodes\n",
    "\n",
    "$ w_j = \\sum_{i=1}^m \\unicode{x1D7D9}$(i reachable within r+1 hops) $\\forall j \\in [1,2,3...,n]$ \n",
    "\n",
    "Let the set of communities be $V_i$. WLOG, choose the target community as $V_1$. To honor the side information constraint, the node weights $W_i$ must follow the below condition:\n",
    "\n",
    "\n",
    "Intuitively, we can see the the quantity (p-q) will get lower bounded.\n",
    "\n",
    "$ \\frac{{(p-q)}^2}{p\\sqrt{p}} \\geq \\frac{\\alpha_{max}}{{\\alpha_{min}}^2\\sqrt n} $\n",
    "\n",
    "where, $\\alpha_i$ denotes the fraction of nodes in community i, i.e. $\\alpha_i=\\frac{n_i}{n}$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Motivation"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Since the graph is random, the adjacency matrix of graph behaves randomly as well\n",
    "* Suppose, we knew the expectation of the adjacency matrix a priori. The expected value of a column j, $\\mathbb{E}[X_j]$ would be then 'p' for the points belonging to the community j belongs to, and 'q' for all other points. \n",
    "> i.e $\\mathbb{E}[X_j]_i=p\\unicode{x1D7D9}($i in same community as j)$+q\\unicode{x1D7D9}($i in diff. community)\n",
    "* Define *community membership vectors*, $\\mu_1,\\mu_2,....\\mu_k$, with $\\mu_i \\in \\mathbb{R}^n$ in the same spirit as above, with:\n",
    "> $\\mathbb{E}[X_j]=\\mu_{c_j}$ where $\\mu_{c_j}$, where $c_j$ is the community j belongs to\n",
    "* Now, the whole problem boils down to estimating $\\mu_1$ from the graph. Let $\\bar{w_j} :=\\mathbb{E}[w_j]$. Define A and B matrices as below."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Proto-Algorithm"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* The main idea now is that we can still recover $\\mu_1$ by “whitening” B using A (All this given that we know, or can estimate $\\mathbb{E}[X]$ closely).\n",
    "\n",
    "* Do rank k SVD of A, to get $A=UDU^T$, let $W:=UD^{-\\frac{1}{2}}$.\n",
    "> Observe, $W^TAW = I_k = \\sum_{i=1}^k \\tilde{\\mu_i}\\tilde{\\mu_i}^T$, where $\\tilde{\\mu_i} := \\sqrt{\\alpha_i}W^T\\mu_i$. Also, note that addition of k terms of type $\\tilde{\\mu_i}\\tilde{\\mu_i}^T$, results in $I_k$, which happens only if corresponding $\\mu_i$ are *orthonormal* vectors in $\\mathbb{R}^k$\n",
    "\n",
    "* Hence we have obtained $\\tilde{\\mu_i}$ which are whitened versions of $\\mu_i$. Now compute $R := W^TBW = \\sum_{i=1}^k w_i\\tilde{\\mu_i}\\tilde{\\mu_i}^T$\n",
    "> Now, since  $\\tilde{\\mu_i}$ are orthonormal, the above equation represents an eigenvalue decomposition of the k × k size matrix R, with eigenvectors $\\tilde{\\mu_i}$ and corresponding eigenvalues $w_i$ . Thus, $\\tilde{\\mu_1}$ – the whitened vector corresponding to the target community – is now the leading eigenvector of R, because $w_1>w_i, \\forall i \\neq 1$ \n",
    "\n",
    "* Find $\\tilde{\\mu_1}$ by setting it to be the leading eigenvector of R. Finally recover $\\mu_1$ from $\\tilde{\\mu_1}$ in the foll. two steps.\n",
    "    * First compute $z := UD^{\\frac{1}{2}}\\tilde{\\mu_1}=\\sqrt{\\alpha_1}\\mu_1$\n",
    "    * Next compute $m_1 := \\frac{1}{n}\\sum_{j=1}^n\\mathbb{E}[X_j] = \\sum_{i=1}^n \\alpha_iu_i.$ Using this, recover $\\sqrt{\\alpha_1}=\\tilde{\\mu_1}^TW^Tm_1$\n",
    "    * Divide the obtained z, by recovered $\\sqrt{\\alpha_1}$ to get $\\mu_1$ and hence the nodes in $V_1$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Unions of $V_p$'s"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Recall that we partitioned graph into 4 parts, $P_1,P_2P_3,P_4$. Under conditions (1)-(3) discussed we can guarantee that the estimated community $\\hat{V}_{P_1}$ has at most $o(\\alpha_1n)$ erroneous nodes. We now want to make a similar statement for $\\hat{V}_{P_i}, i \\in [2,3,4]$. Then, we can take unions of $\\hat{V}_{P_i}$'s and recover all nodes belonging to $V_1$\n",
    "\n",
    "Define $\\forall j \\in P_2, d_j(\\hat{V}_{P_1})$ be the number of edges node j shares with nodes in $\\hat{V}_{P_1}$. \n",
    "\n",
    "Let $v_1 = | \\hat{V}_{P_1}\\cap V_1 |$ be the number of correctly detected nodes in $\\hat{V}_{P_1}$, and $e_1 = | \\hat{V}_{P_1}\\cap V_1^c |$ be the number of erroneous nodes in $\\hat{V}_{P_1}$. \n",
    "\n",
    "From earlier discussions we have with high probability $v_1 = \\Theta(\\alpha_1n)$, $e_1 = o(\\alpha_1n)$. Let $v_0=|V_1\\cap P_1| = \\Theta(\\alpha_1 n)$. Note that $v_1 \\geq v_0 − e_1$ . \n",
    "\n",
    "Now $ \\forall j \\in V_1 \\cap P_2$ using Chernoff bound we get with high probability:\n",
    "\n",
    "$$ d_j(\\hat{V}_{P_1}) \\geq v_1p−\\sqrt{v_1plog(n)}+e_1q−\\sqrt{e_1qlog(n)} \\implies d_j(\\hat{V}_{P_1}) \\geq v_0\\frac{p+q}{2} $$ \n",
    "\n",
    "Similarly $\\forall j \\in V_1^c \\cap P_2$ we get:\n",
    "\n",
    "$$ d_j(\\hat{V}_{P_1}) \\leq v_1q−\\sqrt{v_1qlog(n)}+e_1p−\\sqrt{e_1plog(n)} \\implies d_j(\\hat{V}_{P_1}) < v_0\\frac{p+q}{2} $$ \n",
    "\n",
    "Therefore using a threshold $v_0\\frac{(p+q)}{2}$ we can correctly recover all nodes in $V_1 \\cap P_2$. Observe that this threshold is independent of $P_2$, and hence the same analysis will hold for $P_3$ and $P_4$ as well. Thus, we can detect all the nodes in $V_1$ using a same threshold, i.e. $\\tau = min(\\sqrt{\\alpha_1}\\frac{(p + q)}{2},v_0\\frac{(p+q)}{2})$.\n",
    "\n",
    "$$\\tau = \\sqrt{\\alpha_1}\\frac{(p + q)}{2}$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Attempt: 100\n",
      "Attempt: 242\n",
      "Attempt: 447\n"
     ]
    }
   ],
   "source": [
    "import igraph as ig\n",
    "from IPython.display import SVG, display\n",
    "import numpy as np\n",
    "import random\n",
    "from community_search import community_search \n",
    "n=1000\n",
    "k=5\n",
    "alpha_min=0.1\n",
    "q=0.01\n",
    "p=0.5\n",
    "\n",
    "#get Pref Matrix, i*p \n",
    "Pref_mat=np.matrix(p*np.identity(k)+q*(np.ones(shape=(k,k))-np.identity(k)))\n",
    "#Get block vector by num=amin*n+random(n-amin)\n",
    "min_commsize=int(n*alpha_min)\n",
    "block_vec=[min_commsize]\n",
    "sum_cnt=min_commsize\n",
    "left_elem=k-1\n",
    "for i in range(k-2):\n",
    "\t#print (i)\n",
    "\twhile(1):\n",
    "\t\tprint (\"Attempt: \"+str(sum_cnt))\n",
    "\t\trandom_gen=random.randint(min_commsize,int((n-sum_cnt)/left_elem))\n",
    "\t\tsum_cnt=sum_cnt+random_gen\n",
    "\t\tif((n-sum_cnt)<(left_elem*min_commsize)):\n",
    "\t\t\tsum_cnt=sum_cnt-random_gen\n",
    "\t\t\tcontinue\n",
    "\t\telse:\n",
    "\t\t\tblock_vec.append(random_gen)\n",
    "\t\t\tleft_elem=left_elem-1\n",
    "\t\t\tbreak\n",
    "block_vec.append(n-sum_cnt)\n",
    "#print (block_vec)\n",
    "#print (Pref_mat.tolist())\n",
    "g=ig.Graph.SBM(n,Pref_mat.tolist(),block_vec)\n",
    "\n",
    "#Label m vertices as blue, other as pink\n",
    "m=10\n",
    "r=3 #5 Hops\n",
    "'''Wt assignment\n",
    "#print (g.get_shortest_paths(g.vs[100],g.vs[0:10]))\n",
    "#print ([1 if len(x)<4 else 0 for x in g.get_shortest_paths(g.vs[100],g.vs[0:10])])\n",
    "#print (sum([1 if len(x)<4 else 0 for x in g.get_shortest_paths(g.vs[100],g.vs[0:10])]))\n",
    "less than r hops\n",
    "'''\n",
    "w=[sum([1 if len(x)<r else 0 for x in g.get_shortest_paths(g.vs[i],g.vs[0:m])]) for i in range(n)]\n",
    "#print(w[::4])\n",
    "\n",
    "g.vs[\"side_info_wt\"]=w\n",
    "X_list=g.get_adjacency()\n",
    "X_arr=np.array(X_list[0])\n",
    "for i in range(1,n):\n",
    "\tX_arr=np.vstack([X_arr,X_list[i]])\n",
    "X=np.matrix(X_arr)\n",
    "#g.vs[\"color\"]=[\"red\" if i>(m-1) else \"pink\" for i in range(n)]\n",
    "#display(SVG(ig.plot(g,bbox = (300,300))._repr_svg_()))\n",
    "\n",
    "\n",
    "colors=[\"blue\",\"green\",\"yellow\",\"black\",\"grey\",\"blue\",\"green\",\"yellow\",\"black\",\"grey\"]\n",
    "color_vec=[]\n",
    "for count in range(1):\n",
    "\tnodes_in_V1=community_search(X,k,w,0.0012,n)\n",
    "\tif(count==0):\n",
    "\t\tcolor_vec=[colors[count] if nodes_in_V1[i]==1 else \"red\" for i in range(nodes_in_V1.size)]\n",
    "\t\tg.vs[\"color\"]=color_vec\n",
    "\t\tig.plot(g)\n",
    "\telse:\n",
    "\t\tfor i in range(nodes_in_V1.size):\n",
    "\t\t\tif((nodes_in_V1[i]==1) and (color_vec[i]==\"red\")):\n",
    "\t\t\t\tcolor_vec[i]=colors[count]\n",
    "\n",
    "g.vs[\"color\"]=color_vec\n",
    "#display(SVG(ig.plot(g,bbox = (300,300))._repr_svg_()))\n",
    "#ig.summary(g)\n",
    "\n",
    "\n",
    "#g=ig.Graph.GRG(100,0.2)\n",
    "#display(SVG(ig.plot(g)._repr_svg_()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "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.4.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}