logistic_regression.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "### Statistical Classification\n",
    "We want to classify survining the Titanic. This inference is based the set of features shown below. \n",
    "\n",
    "* $x_1$: Age\n",
    "* $x_2$: Class\n",
    "* $x_3$: Sex\n",
    "\n",
    "The dataset will be split (80/20) into training and testing data. Two classification models will be used: logistic regression (parametric model) and k-nearest-neighbors (non-parametric model)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 356,
   "metadata": {},
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 357,
   "metadata": {},
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "/Library/Frameworks/Python.framework/Versions/3.9/lib/python3.9/site-packages/pandas/core/frame.py:4906: SettingWithCopyWarning: \n",
      "A value is trying to be set on a copy of a slice from a DataFrame\n",
      "\n",
      "See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy\n",
      "  return super().drop(\n"
     ]
    }
   ],
   "source": [
    "TRAINING_TEST_SPLIT = 0.8\n",
    "CLASSIFER_THRESHOLD = 0.5\n",
    "\n",
    "data = pd.read_csv('data/titanic.csv')\n",
    "data.dropna(inplace=True)\n",
    "split_index = int(np.ceil(TRAINING_TEST_SPLIT * len(data)))\n",
    "y = data['Survived']\n",
    "X = data[['Sex', 'Age', 'Pclass']]\n",
    "\n",
    "# Transforming text based binary data to intergers.\n",
    "sex = pd.get_dummies(X['Sex'])\n",
    "sex = sex['male']\n",
    "X.drop(columns=['Sex'], inplace=True)\n",
    "X.insert(1, 'Sex', sex)\n",
    "\n",
    "# Adding bias feature.\n",
    "constants = np.ones(len(X))\n",
    "X.insert(0, 'consant', constants)\n",
    "\n",
    "# numpy conversion.\n",
    "y = np.array(y)\n",
    "X = np.array(X)\n",
    "\n",
    "# Split data into training and test (80/20)\n",
    "training_X = X[0:split_index]\n",
    "training_y = y[0:split_index]\n",
    "test_X = X[split_index:]\n",
    "test_y = y[split_index:]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Logistic Regression\n",
    "\n",
    "We begin by applying logistic regression. Logistic regression is a transformed version of linear regression, which limits the output to $[0,1]$, which represents a probability of being in a class or not. Note, we are only making inference on one class.\n",
    "\n",
    "For logistic regression, we transform the $y_i$ using the Sigmoid function defined as: $f(x) \\triangleq \\frac{1}{1 + e^x}$. Thus we obtain the following model for our classification problem:\n",
    "$$\n",
    "\\mathbb{P}({y_i}=1| \\mathbf{x}, \\mathbf{\\theta}) = \\frac{1}{1 + \\exp \\{\\theta_0 + \\theta_1 x_1 + \\ldots + \\theta_6 x_6\\}} = \\frac{1}{1 + \\exp \\{\\mathbf{x}^{\\text{T}}\\mathbf{\\theta}\\}}\n",
    "$$\n",
    "\n",
    "The sample of observations $\\{y_i\\}_{i=1}^n$ can now be intepreted as a sequence of Bernoulli random variables. Thus we can obtain the following log-likelihood function rather easily:\n",
    "\n",
    "$$\n",
    "-\\ell(\\bold{X}, \\bold{y} | \\mathbf{\\theta}) = \\sum_{i=1}^n \\big[(y_i-1)\\mathbf{\\theta}^{\\text{T}}\\mathbf{x}_i - \\ln(1+\\exp\\{\\mathbf{\\theta}^{\\text{T}}\\mathbf{x}_i\\})\\big]\n",
    "$$\n",
    "\n",
    "Next, we fit the model to the training data using the log-likelihood function above. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 358,
   "metadata": {},
   "outputs": [],
   "source": [
    "def sigmoid(x: np.array) -> float:\n",
    "    return 1 / (1 + np.exp(-x))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 359,
   "metadata": {},
   "outputs": [],
   "source": [
    "def log_loss(params: np.array, y: np.array, X: np.array) -> float:\n",
    "    assert(len(y) == len(X))\n",
    "    n_observations = len(y)\n",
    "    epsilon = 1e-5   \n",
    "    log_loss = 0\n",
    "    for i in range(0, n_observations):\n",
    "        log_loss = log_loss + y[i] * np.log(sigmoid(params.T @ X[i])+epsilon) + (1 - y[i]) * np.log(1-sigmoid(params.T @ X[i])+epsilon)\n",
    "    return log_loss"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 360,
   "metadata": {},
   "outputs": [],
   "source": [
    "from scipy.optimize import minimize\n",
    "\n",
    "# Solving minimzing loss function.\n",
    "inital_solution = np.zeros(4)\n",
    "solution = minimize(log_loss, x0=inital_solution, args=(training_y, training_X))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 361,
   "metadata": {},
   "outputs": [],
   "source": [
    "def classifier(data: np.array, params: np.array, threshold: float = CLASSIFER_THRESHOLD) -> int:\n",
    "    prob = 1 / (1 + np.exp(params.T @ data))\n",
    "    if prob > threshold:\n",
    "        return 1\n",
    "    else:\n",
    "        return 0"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 362,
   "metadata": {},
   "outputs": [],
   "source": [
    "def classify_data(data: np.array, params: np.array) -> np.array:\n",
    "    n_observations = len(data)\n",
    "    prediction = np.zeros(n_observations)\n",
    "    for i in range(0, n_observations):\n",
    "        prediction[i] = classifier(data[i], params)\n",
    "    return prediction"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 363,
   "metadata": {},
   "outputs": [],
   "source": [
    "def calculate_accuracy(prediction: np.array, validation: np.array) -> float:\n",
    "    assert(len(prediction)==len(validation))\n",
    "    return np.sum([1 for i in prediction-validation if i==0]) / len(validation)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 364,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Self-implemented logistic regression accuracy: 72.2 %\n"
     ]
    }
   ],
   "source": [
    "prediction = classify_data(test_X, solution.x)\n",
    "validation = test_y\n",
    "accuracy = calculate_accuracy(prediction, validation)\n",
    "print('Self-implemented logistic regression accuracy:', str(accuracy*100)[:4], '%')   "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Comparison\n",
    "\n",
    "To test the model accuracy two comparisons are made with models from sklearn. Firstly with sklearn's logistic regression package, and secondly with sklearn's kNN (k-nearest-neighbors) algoritm. In doing so we test both a parametric and a non-parametric model against the model implemented above."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 365,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Sklearn-implemented logistic regression accuracy: 75.0 %\n"
     ]
    }
   ],
   "source": [
    "from sklearn.linear_model import LogisticRegression\n",
    "\n",
    "# sklearn logistic regression.\n",
    "lr = LogisticRegression().fit(training_X, training_y)\n",
    "prediction = lr.predict(test_X)\n",
    "validation = test_y\n",
    "accuracy = calculate_accuracy(prediction, validation)\n",
    "print('Sklearn-implemented logistic regression accuracy:', str(accuracy*100)[:4], '%')   "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 366,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Sklearn-implemented kNN accuracy: 66.6 %\n"
     ]
    }
   ],
   "source": [
    "from sklearn.neighbors import KNeighborsClassifier\n",
    "\n",
    "# sklearn kNN algorithm.\n",
    "knn = KNeighborsClassifier().fit(training_X, training_y)\n",
    "prediction = knn.predict(test_X)\n",
    "validation = test_y\n",
    "accuracy = calculate_accuracy(prediction, validation)\n",
    "print('Sklearn-implemented kNN accuracy:', str(accuracy*100)[:4], '%')   "
   ]
  }
 ],
 "metadata": {
  "interpreter": {
   "hash": "aee8b7b246df8f9039afb4144a1f6fd8d2ca17a180786b69acc140d282b71a49"
  },
  "kernelspec": {
   "display_name": "Python 3.9.7 64-bit",
   "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.9.7"
  },
  "orig_nbformat": 4
 },
 "nbformat": 4,
 "nbformat_minor": 2
}