classify_with_res.ipynb

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "\n",
    "from sklearn import cross_validation\n",
    "from sklearn.cross_validation import StratifiedKFold as KFold\n",
    "from sklearn.metrics import classification_report\n",
    "from sklearn.ensemble import RandomForestClassifier as RF\n",
    "import xgboost as xgb\n",
    "import pandas as pd\n",
    "import scipy as sp\n",
    "from glob import glob\n",
    "from skimage import measure"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 45,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "def getRegionFromMap(slice_npy):\n",
    "    thr = np.where(slice_npy > np.mean(slice_npy),0.,1.0)\n",
    "    label_image = measure.label(thr)\n",
    "    labels = label_image.astype(int)\n",
    "    regions = measure.regionprops(labels)\n",
    "    return regions\n",
    "\n",
    "def getRegionMetricRow(fname = \"nodules.npy\"):\n",
    "    # fname, numpy array of dimension [#slices, 1, 512, 512] containing the images\n",
    "    seg = np.load(fname)\n",
    "    nslices = seg.shape[0]\n",
    "    \n",
    "    #metrics\n",
    "    totalArea = 0.\n",
    "    avgArea = 0.\n",
    "    maxArea = 0.\n",
    "    avgEcc = 0.\n",
    "    avgEquivlentDiameter = 0.\n",
    "    stdEquivlentDiameter = 0.\n",
    "    weightedX = 0.\n",
    "    weightedY = 0.\n",
    "    numNodes = 0.\n",
    "    numNodesperSlice = 0.\n",
    "    # crude hueristic to filter some bad segmentaitons\n",
    "    # do not allow any nodes to be larger than 10% of the pixels to eliminate background regions\n",
    "    maxAllowedArea = 0.10 * 512 * 512 \n",
    "    \n",
    "    areas = []\n",
    "    eqDiameters = []\n",
    "    #for slicen in range(nslices):\n",
    "    regions = getRegionFromMap(seg[:,:])\n",
    "    if len(regions) > 0:\n",
    "        for region in regions:\n",
    "            #if region.area > maxAllowedArea:\n",
    "                #continue\n",
    "            totalArea += region.area\n",
    "            areas.append(region.area)\n",
    "            #avgEcc += region.eccentricity\n",
    "            #avgEquivlentDiameter += region.equivalent_diameter\n",
    "            #eqDiameters.append(region.equivalent_diameter)\n",
    "            #weightedX += region.centroid[0]*region.area\n",
    "            #weightedY += region.centroid[1]*region.area\n",
    "            numNodes += 1\n",
    "            \n",
    "        #weightedX = weightedX / totalArea \n",
    "        #weightedY = weightedY / totalArea\n",
    "        #avgArea = totalArea / numNodes\n",
    "        avgArea = totalArea\n",
    "        #avgEcc = avgEcc / numNodes\n",
    "        #avgEquivlentDiameter = avgEquivlentDiameter / numNodes\n",
    "        #stdEquivlentDiameter = np.std(eqDiameters)\n",
    "\n",
    "        #maxArea = max(areas)\n",
    "\n",
    "\n",
    "        #numNodesperSlice = numNodes*1. / nslices\n",
    "    \n",
    "    \n",
    "        #return np.array([avgArea,maxArea,avgEcc,avgEquivlentDiameter,\\\n",
    "        #                 stdEquivlentDiameter, weightedX, weightedY, numNodes, numNodesperSlice])\n",
    "        print (\"numNodes: \", numNodes)\n",
    "        print (\"aveArea: \", avgArea)\n",
    "        return np.array([numNodes, avgArea])\n",
    "    \n",
    "    return np.zeros(2)\n",
    "\n",
    "\n",
    "def createFeatureDataset(nodfiles=None):\n",
    "    if nodfiles == None:\n",
    "        # directory of numpy arrays containing masks for nodules\n",
    "        # found via unet segmentation\n",
    "        noddir = \"training_set/\" \n",
    "        nodfiles = glob(noddir +\"*npy\")\n",
    "    # dict with mapping between training examples and true labels\n",
    "    # the training set is the output masks from the unet segmentation\n",
    "    truthdata = pd.read_csv('data/stage1_labels.csv')\n",
    "    \n",
    "    numfeatures = 2\n",
    "    feature_array = np.zeros((len(nodfiles),numfeatures))\n",
    "    truth_metric = np.zeros((len(nodfiles)))\n",
    "    \n",
    "    for i,nodfile in enumerate(nodfiles):\n",
    "        patID = nodfile.split(\"/\")[1].split(\".\")[0]\n",
    "        truth_metric[i] = truthdata[truthdata.id == patID].cancer\n",
    "        feature_array[i] = getRegionMetricRow(nodfile)\n",
    "    \n",
    "    np.save(\"dataY.npy\", truth_metric)\n",
    "    np.save(\"dataX.npy\", feature_array)\n",
    "\n",
    "\n",
    "def logloss(act, pred):\n",
    "    epsilon = 1e-15\n",
    "    pred = sp.maximum(epsilon, pred)\n",
    "    pred = sp.minimum(1-epsilon, pred)\n",
    "    ll = sum(act*sp.log(pred) + sp.subtract(1,act)*sp.log(sp.subtract(1,pred)))\n",
    "    ll = ll * -1.0/len(act)\n",
    "    return ll\n",
    "\n",
    "\n",
    "def classifyData():\n",
    "    X = np.load(\"dataX.npy\")\n",
    "    Y = np.load(\"dataY.npy\")\n",
    "\n",
    "    print \"Random Forest\"\n",
    "    kf = KFold(Y, n_folds=3)\n",
    "    y_pred = Y * 0\n",
    "    for train, test in kf:\n",
    "        X_train, X_test, y_train, y_test = X[train,:], X[test,:], Y[train], Y[test]\n",
    "        clf = RF(n_estimators=100, n_jobs=3)\n",
    "        clf.fit(X_train, y_train)\n",
    "        y_pred[test] = clf.predict(X_test)\n",
    "    print classification_report(Y, y_pred, target_names=[\"No Cancer\", \"Cancer\"])\n",
    "    print(\"logloss\",logloss(Y, y_pred))\n",
    "\n",
    "    # All Cancer\n",
    "    print \"Predicting all positive\"\n",
    "    y_pred = np.ones(Y.shape)\n",
    "    print classification_report(Y, y_pred, target_names=[\"No Cancer\", \"Cancer\"])\n",
    "    print(\"logloss\",logloss(Y, y_pred))\n",
    "\n",
    "    # No Cancer\n",
    "    print \"Predicting all negative\"\n",
    "    y_pred = Y*0\n",
    "    print classification_report(Y, y_pred, target_names=[\"No Cancer\", \"Cancer\"])\n",
    "    print(\"logloss\",logloss(Y, y_pred))\n",
    "\n",
    "    # try XGBoost\n",
    "    print (\"XGBoost\")\n",
    "    kf = KFold(Y, n_folds=3)\n",
    "    y_pred = Y * 0\n",
    "    for train, test in kf:\n",
    "        X_train, X_test, y_train, y_test = X[train,:], X[test,:], Y[train], Y[test]\n",
    "        clf = xgb.XGBClassifier(objective=\"binary:logistic\")\n",
    "        clf.fit(X_train, y_train)\n",
    "        y_pred[test] = clf.predict(X_test)\n",
    "    print classification_report(Y, y_pred, target_names=[\"No Cancer\", \"Cancer\"])\n",
    "    print(\"logloss\",logloss(Y, y_pred))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 46,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 259621.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 259621.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 259621.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 259621.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260994.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 259621.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 259677.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 259792.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 259792.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260543.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 259621.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260723.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 260462.0)\n",
      "('numNodes: ', 1.0)\n",
      "('aveArea: ', 261303.0)\n",
      "Random Forest\n",
      "             precision    recall  f1-score   support\n",
      "\n",
      "  No Cancer       0.86      0.95      0.90        64\n",
      "     Cancer       0.67      0.38      0.48        16\n",
      "\n",
      "avg / total       0.82      0.84      0.82        80\n",
      "\n",
      "('logloss', 5.6125811490766164)\n",
      "Predicting all positive\n",
      "             precision    recall  f1-score   support\n",
      "\n",
      "  No Cancer       0.00      0.00      0.00        64\n",
      "     Cancer       0.20      1.00      0.33        16\n",
      "\n",
      "avg / total       0.04      0.20      0.07        80\n",
      "\n",
      "('logloss', 27.631660793872687)\n",
      "Predicting all negative\n",
      "             precision    recall  f1-score   support\n",
      "\n",
      "  No Cancer       0.80      1.00      0.89        64\n",
      "     Cancer       0.00      0.00      0.00        16\n",
      "\n",
      "avg / total       0.64      0.80      0.71        80\n",
      "\n",
      "('logloss', 6.9077552789821342)\n",
      "XGBoost\n",
      "             precision    recall  f1-score   support\n",
      "\n",
      "  No Cancer       0.82      0.98      0.89        64\n",
      "     Cancer       0.67      0.12      0.21        16\n",
      "\n",
      "avg / total       0.79      0.81      0.76        80\n",
      "\n",
      "('logloss', 6.4760305690136279)\n"
     ]
    }
   ],
   "source": [
    "createFeatureDataset()\n",
    "classifyData()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "kernelspec": {
   "display_name": "Python [default]",
   "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.12"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}