Code/assumption_checking_tyler.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Checking Assumptions"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Library/Python/2.7/site-packages/rpy2/robjects/packages.py:63: UserWarning: Note: the specification for S3 class “family” in package ‘MatrixModels’ seems equivalent to one from package ‘lme4’: not turning on duplicate class definitions for this class.\n",
      "\n",
      "  return _reval(expr)\n"
     ]
    }
   ],
   "source": [
    "import itertools\n",
    "import csv\n",
    "import numpy as np\n",
    "from scipy import linalg\n",
    "from scipy.stats import cumfreq\n",
    "import matplotlib.pyplot as plt\n",
    "import matplotlib as mpl\n",
    "from sklearn import mixture\n",
    "from sklearn import cluster\n",
    "from rpy2 import robjects as ro\n",
    "from rpy2.robjects.packages import importr\n",
    "mvn = importr('MVN')\n",
    "%matplotlib inline\n",
    "np.random.seed(1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Load data"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "data loaded\n"
     ]
    }
   ],
   "source": [
    "X = np.load('SynapseFeatures.npy')\n",
    "print 'data loaded'"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Cluster data into 17 clusters using k-means"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "k-means complete\n"
     ]
    }
   ],
   "source": [
    "kmeans = cluster.MiniBatchKMeans(n_clusters=17)\n",
    "ClusterIdx = kmeans.fit_predict(X)\n",
    "print 'k-means complete'"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "from rpy2.robjects.packages import importr\n",
    "mvn = importr('MVN')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "ro.r.assign('ClusterIdx',ClusterIdx.tolist())\n",
    "ro.r('ClusterIdx <- unlist(ClusterIdx)')\n",
    "ro.r.assign('Xr',ro.FloatVector(np.ravel(X)))\n",
    "nr = X.shape[0]\n",
    "ro.r.assign('nr',nr)\n",
    "ro.r('Xr <- matrix(Xr, nrow = nr, byrow = TRUE)')\n",
    "ro.r('''\n",
    "     for (i in seq(1,max(ClusterIdx)))\n",
    "     {pdf(paste(\"qq_plot_\",i,\".pdf\"))\n",
    "     qqnorm(Xr[ClusterIdx==i,1:ncol(Xr)], main = \"Normal Q-Q Plot\",\n",
    "         xlab = \"Theoretical Quantiles\", ylab = \"Sample Quantiles\",\n",
    "         plot.it = TRUE)\n",
    "     qqline(Xr[ClusterIdx==i,1:ncol(Xr)])\n",
    "     dev.off()\n",
    "     hzTestResult <- hzTest(Xr, qqplot = FALSE)\n",
    "     print(hzTestResult)\n",
    "     }\n",
    "     ''')\n",
    "#ro.r('pdf(\"qq_plot.pdf\")')\n",
    "#ro.r('''qqnorm(Xr, main = \"Normal Q-Q Plot\",\n",
    "#   xlab = \"Theoretical Quantiles\", ylab = \"Sample Quantiles\",\n",
    "#   plot.it = TRUE)''')\n",
    "#ro.r('''qqline(Xr)''')\n",
    "#ro.r('dev.off()')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "ro.r('hzTestResult <- hzTest(Xr, qqplot = FALSE)')\n",
    "hzTestStat = ro.r('hzTestResult')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "  Henze-Zirkler's Multivariate Normality Test \n",
      "--------------------------------------------- \n",
      "  data : Xr \n",
      "\n",
      "  HZ      : 1 \n",
      "  p-value : 0 \n",
      "\n",
      "  Result  : Data are not multivariate normal. \n",
      "--------------------------------------------- \n",
      "\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<RS4 - Python:0x10d022638 / R:0x7fc67ef02db8>"
      ]
     },
     "execution_count": 11,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ro.r('print(hzTestResult)')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 19,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "[1] 100000\n",
      "[1] 8897\n"
     ]
    },
    {
     "data": {
      "text/plain": [
       "<IntVector - Python:0x10d020ab8 / R:0x7fc67da26aa8>\n",
       "[    8897]"
      ]
     },
     "execution_count": 19,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "ro.r('print(length(ClusterIdx))')\n",
    "ro.r('print(nrow(Xr))')"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 2",
   "language": "python",
   "name": "python2"
  },
  "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.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}