{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "import os\n",
    "os.chdir('/Users/toshan/dev/pysal/pysal/contrib/spint')\n",
    "from count_model import CountModel\n",
    "import numpy as np\n",
    "import pandas as pd\n",
    "from gravity import Gravity, Production, Attraction, Doubly, BaseGravity\n",
    "import statsmodels.formula.api as smf\n",
    "from statsmodels.api import families\n",
    "os.chdir('/Users/toshan/dev/pysal/pysal/contrib/glm')\n",
    "from glm import GLM\n",
    "from family import Poisson, QuasiPoisson\n",
    "\n",
    "import pysal\n",
    "import os\n",
    "os.chdir('/Users/toshan/dev/pysal/pysal/contrib/spint')\n",
    "from dispersion import alpha_disp, phi_disp\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "rec = pd.read_csv('/Users/toshan/Documents/RecreationDemand.csv')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<div style=\"max-height:1000px;max-width:1500px;overflow:auto;\">\n",
       "<table border=\"1\" class=\"dataframe\">\n",
       "  <thead>\n",
       "    <tr style=\"text-align: right;\">\n",
       "      <th></th>\n",
       "      <th>Unnamed: 0</th>\n",
       "      <th>trips</th>\n",
       "      <th>quality</th>\n",
       "      <th>ski</th>\n",
       "      <th>income</th>\n",
       "      <th>userfee</th>\n",
       "      <th>costC</th>\n",
       "      <th>costS</th>\n",
       "      <th>costH</th>\n",
       "    </tr>\n",
       "  </thead>\n",
       "  <tbody>\n",
       "    <tr>\n",
       "      <th>0</th>\n",
       "      <td> 1</td>\n",
       "      <td> 0</td>\n",
       "      <td> 0</td>\n",
       "      <td> yes</td>\n",
       "      <td> 4</td>\n",
       "      <td> no</td>\n",
       "      <td> 67.59</td>\n",
       "      <td> 68.620</td>\n",
       "      <td> 76.800</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>1</th>\n",
       "      <td> 2</td>\n",
       "      <td> 0</td>\n",
       "      <td> 0</td>\n",
       "      <td>  no</td>\n",
       "      <td> 9</td>\n",
       "      <td> no</td>\n",
       "      <td> 68.86</td>\n",
       "      <td> 70.936</td>\n",
       "      <td> 84.780</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>2</th>\n",
       "      <td> 3</td>\n",
       "      <td> 0</td>\n",
       "      <td> 0</td>\n",
       "      <td> yes</td>\n",
       "      <td> 5</td>\n",
       "      <td> no</td>\n",
       "      <td> 58.12</td>\n",
       "      <td> 59.465</td>\n",
       "      <td> 72.110</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>3</th>\n",
       "      <td> 4</td>\n",
       "      <td> 0</td>\n",
       "      <td> 0</td>\n",
       "      <td>  no</td>\n",
       "      <td> 2</td>\n",
       "      <td> no</td>\n",
       "      <td> 15.79</td>\n",
       "      <td> 13.750</td>\n",
       "      <td> 23.680</td>\n",
       "    </tr>\n",
       "    <tr>\n",
       "      <th>4</th>\n",
       "      <td> 5</td>\n",
       "      <td> 0</td>\n",
       "      <td> 0</td>\n",
       "      <td> yes</td>\n",
       "      <td> 3</td>\n",
       "      <td> no</td>\n",
       "      <td> 24.02</td>\n",
       "      <td> 34.033</td>\n",
       "      <td> 34.547</td>\n",
       "    </tr>\n",
       "  </tbody>\n",
       "</table>\n",
       "</div>"
      ],
      "text/plain": [
       "   Unnamed: 0  trips  quality  ski  income userfee  costC   costS   costH\n",
       "0           1      0        0  yes       4      no  67.59  68.620  76.800\n",
       "1           2      0        0   no       9      no  68.86  70.936  84.780\n",
       "2           3      0        0  yes       5      no  58.12  59.465  72.110\n",
       "3           4      0        0   no       2      no  15.79  13.750  23.680\n",
       "4           5      0        0  yes       3      no  24.02  34.033  34.547"
      ]
     },
     "execution_count": 3,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "rec.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 4,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "y = rec['trips'].values.reshape((-1,1))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 5,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "X = rec[['quality', 'income', 'costC', 'costS', 'costH']].values"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 22,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "test = CountModel(y, X, constant=False)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 23,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "glm_results = test.fit(framework='glm')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 24,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([  7.30811593e+00,   2.71035909e+00,   3.36051997e-03])"
      ]
     },
     "execution_count": 24,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "phi_disp(glm_results)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 26,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([  6.30811593e+00,   2.71035909e+00,   3.36051997e-03])"
      ]
     },
     "execution_count": 26,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "alpha_disp(glm_results)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 27,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([  1.55402055e+00,   3.38253708e+00,   3.59097912e-04])"
      ]
     },
     "execution_count": 27,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "alpha_disp(glm_results, lambda x: x**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 28,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#Prepare some test data - columbus example\n",
    "db = pysal.open(pysal.examples.get_path('columbus.dbf'),'r')\n",
    "y = np.array(db.by_col(\"HOVAL\"))\n",
    "y = np.reshape(y, (49,1))\n",
    "X = []\n",
    "#X.append(np.ones(len(y)))\n",
    "X.append(db.by_col(\"INC\"))\n",
    "X.append(db.by_col(\"CRIME\"))\n",
    "X = np.array(X).T\n",
    "\n",
    "poisson_y = np.round(y).astype(int)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 29,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "test = CountModel(poisson_y, X)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 30,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "glm_results = test.fit(framework='glm')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 31,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 5.39968689,  2.3230411 ,  0.01008847])"
      ]
     },
     "execution_count": 31,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "phi_disp(glm_results)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 32,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 4.39968689,  2.3230411 ,  0.01008847])"
      ]
     },
     "execution_count": 32,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "alpha_disp(glm_results)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 33,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "array([ 0.10690133,  2.24709978,  0.01231683])"
      ]
     },
     "execution_count": 33,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "alpha_disp(glm_results, lambda x:x**2)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "<class 'family.Poisson'>\n",
      "<class 'family.Poisson'>\n",
      "<class 'family.Poisson'>\n"
     ]
    }
   ],
   "source": [
    "model1 = GLM(y, X, constant=False, family=Poisson()).fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "<class 'family.QuasiPoisson'>\n",
      "<class 'family.QuasiPoisson'>\n",
      "<class 'family.QuasiPoisson'>\n"
     ]
    }
   ],
   "source": [
    "model2 = GLM(y, X, constant=False, family=QuasiPoisson()).fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "1.0\n",
      "<class 'family.QuasiPoisson'>\n",
      "7.02573401193\n"
     ]
    }
   ],
   "source": [
    "print model1.scale\n",
    "print model2.scale"
   ]
  }
 ],
 "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.9"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}