Chapter5GLM.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Chapter 5 Generalized Linear Models\n",
    "\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5.1 Logistic regression\n",
    "\n",
    "To model a binary outcome, $(X_i, y_i)$, $y_i \\in \\{0,1\\}$, $X_i \\in \\mathbb{R}^p$, we use \n",
    "$$\n",
    "{\\rm logit}(\\pi_i)=X_i^{T} \\beta,\n",
    "$$\n",
    "where $\\pi_i = p(y_i=1|X_i)$ and ${\\rm logit}(a) = \\log(a/(1-a))$. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "**Estimation** can be done with maximum likelihood estimation \n",
    "$$\n",
    "\\mathcal{L} = \\prod_{i=1}^n \\pi^{Y_i} (1-\\pi_i)^{1-Y_i}.\n",
    "$$\n",
    "Taking the logarithm yields \n",
    "$$\n",
    "\\log \\mathcal{L} = \\sum_{i=1}^n \\big[ Y_i \\log \\pi + (1-Y_i) \\log (1-\\pi_i)\\big].\n",
    "$$\n",
    "If ${\\rm logit}(\\pi)i)=\\beta_0$, we can derive that $\\pi_i = \\exp(\\beta_0)/[\\exp(\\beta_0)+1]$. \n",
    "Hence, $\\log \\pi_i = \\beta_0 - \\log[\\exp(\\beta_0)+1]$ and $\\log (1-\\pi_i) = 0 - \\log[\\exp(\\beta_0)+1]$. \n",
    "\n",
    "Plugging these into the log-likelihood yields \n",
    "$$\n",
    "\\log \\mathcal{L} = n_1 \\beta_0 - n\\log[\\exp(\\beta_0)+1].\n",
    "$$\n",
    "The first order condition yields \n",
    "$$\n",
    "\\frac{n_1}{n} = \\frac{\\exp(\\beta_0)  }{\\exp(\\beta_0)+1}.\n",
    "$$\n",
    "\n",
    "With covariates $X$, a closed form solution is not available. However, one can find the solution using optimization. \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Interpretation** \n",
    "\n",
    "A key concept for the logistic regression is the log odds ration, which is the interpretation of $\\beta$. Consider $\\beta_1$, we can see that \n",
    "$$\n",
    "\\beta_1= \\log \\left\\{ \\frac{ \\pi(X_{1}+1,\\cdots )/[1-\\pi(X_{1}+1,\\cdots )]}{\\pi(X_{1},\\cdots )/[1-\\pi(X_{1},\\cdots )]}    \\right\\}.\n",
    "$$\n",
    "\n",
    "The score function takes the form $U(\\beta)= \\partial \\log \\mathcal{L}/\\partial \\beta$. The MLE $\\hat{\\beta}$ satisfies that $U(\\hat{\\beta})=0$. The Fisher information is $I(\\beta)=\\mathbb{E} \\big[-\\partial^2 \\log \\mathcal{L}/\\partial \\beta^2 \\big]$. The confidence interval for each $\\beta_j$ is $\\hat{\\beta}_j \\mp z_{1-\\alpha/2} \\hat{\\rm se}(\\hat{\\beta}_j)$, where $z \\sim N(0,1)$. \n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "**Hypothesis testing**\n",
    "\n",
    "Consider the following null and alternative hypotheses $H_0: \\beta_1=0$ v.s. $H_1: \\beta_1=0$. There are three tests for this task. Note that all three tests are the same in linear regression, whereas they differ in the presence of nonlinearity. \n",
    "\n",
    "1. Likelihood ratio test. ${\\rm LR}=-2 \\big[ \\log \\mathcal{L}({\\rm reduced}) -\\log \\mathcal{L}({\\rm full}) \\big]$. We have that ${\\rm LR} \\sim \\chi^2_{K-k}$, under $H_0$ and large $n$. \n",
    "\n",
    "2. Score test. $S=U(\\beta)/I(\\beta) \\sim \\chi^2$ with $d.f.=1$.\n",
    "\n",
    "3. Wald test. Robust wald test. For $K-k=1$, we have $W=\\hat{\\beta}_1/{\\rm se}(\\hat{\\beta}_1) \\sim N(0,1)$ under $H_0$ and large $n$. The robustness comes from a robust standard error. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "\n",
    "\n",
    "**Visualization**\n",
    "\n",
    "Pearson residual: $r_{p_i}=(Y_i-\\hat{\\pi}_i)/[\\hat{\\pi}_i(1-\\hat{\\pi}_i) ]^{1/2}$\n",
    "\n",
    "Deviance residual: ${\\rm dev}_i = {\\rm sign}(Y_i-\\hat{\\pi}_i) \\sqrt{ -2[Y_i \\log \\hat{\\pi}_i +(1-Y_i)\\log  (1-\\hat{\\pi}_i) ]}$.\n",
    "\n",
    "One can define $H^2 = \\sum r_{p_i}^2$ and $G^2 = \\sum {\\rm dev}_i^2$, which are similar to ${\\rm SSE}$. They do not follow a $\\chi^2$ distribution, but we can use them anyway. \n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## 5.2 Generalized linear model \n",
    "\n",
    "The logistic regression is a special case of the generalized linear model. \n",
    "\n",
    "Suppose that the outcome $y$ follows a distribution that belongs to the exponential family \n",
    "$$\n",
    "p(y | \\theta, \\tau)= h(y,\\tau) \\exp\\left\\{ \\frac{b(\\theta) T(y)-A(\\theta)  }{d(\\tau)} \\right\\}.\n",
    "$$\n",
    "We usually use the canonical form \n",
    "$$\n",
    "p(y | \\theta, \\tau)= h(y,\\tau) \\exp\\left\\{ \\frac{\\theta T(y)-A(\\theta)  }{d \\tau} \\right\\}.\n",
    "$$\n",
    "\n",
    "For regular use, we specify a GLM in the following manner. \n",
    "\n",
    "1. Pick an exponential form of distribution.\n",
    "2. Set the predictor $X\\beta$.\n",
    "3. Pick a link function $g$ \n",
    "$$\n",
    "\\mathbb{E}[Y|X]=\\mu = g^{-1}(y).\n",
    "$$\n",
    "\n",
    "Note that steps (1) and (3) are often tied. \n"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "R",
   "language": "R",
   "name": "ir"
  },
  "language_info": {
   "codemirror_mode": "r",
   "file_extension": ".r",
   "mimetype": "text/x-r-source",
   "name": "R",
   "pygments_lexer": "r",
   "version": "4.0.3"
  },
  "latex_envs": {
   "LaTeX_envs_menu_present": true,
   "autoclose": false,
   "autocomplete": true,
   "bibliofile": "biblio.bib",
   "cite_by": "apalike",
   "current_citInitial": 1,
   "eqLabelWithNumbers": true,
   "eqNumInitial": 1,
   "hotkeys": {
    "equation": "Ctrl-E",
    "itemize": "Ctrl-I"
   },
   "labels_anchors": false,
   "latex_user_defs": false,
   "report_style_numbering": false,
   "user_envs_cfg": false
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}