{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If this tutorial we are going to demonstrate the use and effect of different graph filtering (thresholding schemes)\n",
    "\n",
    "- - -\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Load and prepare data "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import sys"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "np.set_printoptions(precision=2, threshold=sys.maxsize)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "rng = np.random.RandomState(0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "all_conn_mtx = np.load(\"data/fmri_autism_conn_mtx.npy\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "conn_mtx = np.abs(all_conn_mtx[0, ])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Load all related modules and methods"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy.sparse.csgraph import minimum_spanning_tree"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from dyconnmap.graphs import threshold_omst_global_cost_efficiency"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Filter the matrix"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Minumum Spanning Tree "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "mtx_mst = minimum_spanning_tree(conn_mtx).todense()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Orthogonal Minimum Spanning Tree "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nCIJtree, CIJtree, _, _, global_cost_eff_max, cost_max, cost, global_cost_eff = threshold_omst_global_cost_efficiency(conn_mtx)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Plot"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as plt"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f, axes = plt.subplots(ncols=3, nrows=1, figsize=(10, 8), dpi=100, sharey=True, sharex=False)\n",
    "\n",
    "im = axes[0].imshow(conn_mtx, vmin=0.0, vmax=1.0, cmap=plt.cm.Spectral)\n",
    "axes[0].set_title('Input', fontsize=14)\n",
    "axes[0].set_ylabel('ROI', fontsize=14)\n",
    "axes[0].set_xlabel('ROI', fontsize=14)\n",
    "\n",
    "axes[1].imshow(mtx_mst, vmin=0.0, vmax=1.0, cmap=plt.cm.Spectral)\n",
    "axes[1].set_title('MST', fontsize=14)\n",
    "axes[1].set_xlabel('ROI', fontsize=14)\n",
    "\n",
    "axes[2].imshow(CIJtree, vmin=0.0, vmax=1.0, cmap=plt.cm.Spectral)\n",
    "axes[2].set_title('OMST', fontsize=14)\n",
    "axes[2].set_xlabel('ROI', fontsize=14)\n",
    "\n",
    "f.subplots_adjust(right=0.8)\n",
    "cbar_ax = f.add_axes([1.0, 0.35, 0.0125, 0.300])\n",
    "cb = f.colorbar(im, cax=cbar_ax)\n",
    "cb.set_label('Connectivity', fontsize=14)\n",
    "\n",
    "plt.tight_layout()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Sparsity "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nnz = 32 * 32\n",
    "nnz_mst = np.count_nonzero(mtx_mst == 0)\n",
    "nnz_omst = np.count_nonzero(CIJtree == 0)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.figure(figsize=(8, 6), dpi=100)\n",
    "plt.plot([nnz, nnz_mst, nnz_omst], 'o-')\n",
    "plt.axis('tight')\n",
    "plt.ylabel('Non-zero elements', fontsize=12)\n",
    "plt.xticks([0, 1, 2], ['Max', 'MST', 'OMST'], rotation=45, fontsize=12)\n",
    "plt.grid(alpha=0.25)\n",
    "plt.title('Sparsity', fontsize=14)\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "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.8.7"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}