modeling_tfr

{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stderr",
     "output_type": "stream",
     "text": [
      "Loading required package: Matrix\n"
     ]
    }
   ],
   "source": [
    "library(lme4)\n",
    "library(plyr)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 2,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "data <- read.csv('/ihme/forecasting/data/fertility/inputs/asfr_ldi_edu.csv')\n",
    "data <- data[!(data$asfr==0),]\n",
    "mean.age <- read.csv('/ihme/forecasting/data/fertility/inputs/mean_fertile_age_by_asfr_ratio.csv')\n",
    "mean.age <- rename(mean.age,c('year' = 'year_id'))\n",
    "data <- merge.data.frame(x = data,y = mean.age,by = c('location_id','year_id'))\n",
    "data$sex <- NULL\n",
    "data$iso3 <- NULL"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 20,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/html": [
       "<table>\n",
       "<thead><tr><th></th><th scope=col>location_id</th><th scope=col>year_id</th><th scope=col>age_group_id</th><th scope=col>asfr</th><th scope=col>ldi</th><th scope=col>edu</th><th scope=col>mean_fertile_age</th></tr></thead>\n",
       "<tbody>\n",
       "\t<tr><th scope=row>1</th><td>101</td><td>1951</td><td>9</td><td>0.204566</td><td>9.36871</td><td>6.826699</td><td>27.82598</td></tr>\n",
       "\t<tr><th scope=row>2</th><td>101</td><td>1951</td><td>13</td><td>0.03092</td><td>9.36871</td><td>6.414093</td><td>27.82598</td></tr>\n",
       "\t<tr><th scope=row>3</th><td>101</td><td>1951</td><td>11</td><td>0.149179</td><td>9.36871</td><td>6.93765</td><td>27.82598</td></tr>\n",
       "\t<tr><th scope=row>4</th><td>101</td><td>1951</td><td>10</td><td>0.2077373</td><td>9.36871</td><td>6.93765</td><td>27.82598</td></tr>\n",
       "\t<tr><th scope=row>5</th><td>101</td><td>1951</td><td>14</td><td>0.0028385</td><td>9.36871</td><td>5.782532</td><td>27.82598</td></tr>\n",
       "\t<tr><th scope=row>6</th><td>101</td><td>1951</td><td>12</td><td>0.0873706</td><td>9.36871</td><td>6.414093</td><td>27.82598</td></tr>\n",
       "</tbody>\n",
       "</table>\n"
      ],
      "text/latex": [
       "\\begin{tabular}{r|lllllll}\n",
       "  & location_id & year_id & age_group_id & asfr & ldi & edu & mean_fertile_age\\\\\n",
       "\\hline\n",
       "\t1 & 101 & 1951 & 9 & 0.204566 & 9.36871 & 6.826699 & 27.82598\\\\\n",
       "\t2 & 101 & 1951 & 13 & 0.03092 & 9.36871 & 6.414093 & 27.82598\\\\\n",
       "\t3 & 101 & 1951 & 11 & 0.149179 & 9.36871 & 6.93765 & 27.82598\\\\\n",
       "\t4 & 101 & 1951 & 10 & 0.2077373 & 9.36871 & 6.93765 & 27.82598\\\\\n",
       "\t5 & 101 & 1951 & 14 & 0.0028385 & 9.36871 & 5.782532 & 27.82598\\\\\n",
       "\t6 & 101 & 1951 & 12 & 0.0873706 & 9.36871 & 6.414093 & 27.82598\\\\\n",
       "\\end{tabular}\n"
      ],
      "text/plain": [
       "  location_id year_id age_group_id      asfr     ldi      edu mean_fertile_age\n",
       "1         101    1951            9 0.2045660 9.36871 6.826699         27.82598\n",
       "2         101    1951           13 0.0309200 9.36871 6.414093         27.82598\n",
       "3         101    1951           11 0.1491790 9.36871 6.937650         27.82598\n",
       "4         101    1951           10 0.2077373 9.36871 6.937650         27.82598\n",
       "5         101    1951           14 0.0028385 9.36871 5.782532         27.82598\n",
       "6         101    1951           12 0.0873706 9.36871 6.414093         27.82598"
      ]
     },
     "execution_count": 20,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "head(data)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 6,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "regression <- function(formula) {\n",
    "    lm_lst <- list()\n",
    "    \n",
    "    for (age_group_id in 8:14) {\n",
    "        \n",
    "        temp <- data[data$age_group_id == age_group_id,]\n",
    "        lm_lst[[length(lm_lst) + 1]] <- formula(temp)\n",
    "        \n",
    "    }\n",
    "    return(lm_lst)\n",
    "}    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$log(ASFR) = \\alpha + \\beta_1 LDI + \\beta_2 \\text{matern_edu}$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 7,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "ldi.edu.formula <- function(data) {\n",
    "    formula <- lm('log(asfr) ~ ldi + edu',data=data)\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 8,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "ldi.edu.lm <- regression(ldi.edu.formula)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 9,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "Call:\n",
       "lm(formula = \"log(asfr) ~ ldi + edu\", data = data)\n",
       "\n",
       "Residuals:\n",
       "    Min      1Q  Median      3Q     Max \n",
       "-3.9380 -0.3764  0.0954  0.4697  1.3581 \n",
       "\n",
       "Coefficients:\n",
       "             Estimate Std. Error t value Pr(>|t|)    \n",
       "(Intercept) -0.800169   0.050520  -15.84   <2e-16 ***\n",
       "ldi         -0.141972   0.007044  -20.16   <2e-16 ***\n",
       "edu         -0.142461   0.002448  -58.19   <2e-16 ***\n",
       "---\n",
       "Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
       "\n",
       "Residual standard error: 0.6698 on 11972 degrees of freedom\n",
       "Multiple R-squared:  0.4705,\tAdjusted R-squared:  0.4704 \n",
       "F-statistic:  5320 on 2 and 11972 DF,  p-value: < 2.2e-16\n"
      ]
     },
     "execution_count": 9,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "summary(ldi.edu.lm[[1]])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$log(ASFR) = \\alpha_0 + \\alpha_{country} + \\beta_1 LDI + \\beta_2 \\text{matern_edu}$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 10,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "ldi.edu.country.formula <- function(data) {\n",
    "    formula <- lmer('log(asfr) ~ 1 + ldi + edu + (1 | location_id)',data=data)\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 11,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "ldi.edu.country.lmer <- regression(ldi.edu.country.formula)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 12,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Linear mixed model fit by REML ['lmerMod']\n",
       "Formula: log(asfr) ~ 1 + ldi + edu + (1 | location_id)\n",
       "   Data: data\n",
       "\n",
       "REML criterion at convergence: 6218.3\n",
       "\n",
       "Scaled residuals: \n",
       "    Min      1Q  Median      3Q     Max \n",
       "-5.6670 -0.4838  0.0367  0.5111  4.4597 \n",
       "\n",
       "Random effects:\n",
       " Groups      Name        Variance Std.Dev.\n",
       " location_id (Intercept) 0.37924  0.6158  \n",
       " Residual                0.08998  0.3000  \n",
       "Number of obs: 11975, groups:  location_id, 188\n",
       "\n",
       "Fixed effects:\n",
       "             Estimate Std. Error t value\n",
       "(Intercept) -1.695783   0.080025  -21.19\n",
       "ldi         -0.038773   0.008734   -4.44\n",
       "edu         -0.137771   0.001769  -77.87\n",
       "\n",
       "Correlation of Fixed Effects:\n",
       "    (Intr) ldi   \n",
       "ldi -0.822       \n",
       "edu  0.539 -0.731"
      ]
     },
     "execution_count": 12,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "summary(ldi.edu.country.lmer[[1]])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$log(ASFR) = \\alpha + \\beta_1 LDI + \\beta_2 \\text{matern_edu} + \\beta_3 \\text{mean_fertile_age}$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 13,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "ldi.edu.mean.age.formula <- function(data) {\n",
    "    formula <- lm('log(asfr) ~ ldi + edu + mean_fertile_age',data=data)\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 14,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "ldi.edu.mean.age.lm <- regression(ldi.edu.mean.age.formula)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 15,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "\n",
       "Call:\n",
       "lm(formula = \"log(asfr) ~ ldi + edu + mean_fertile_age\", data = data)\n",
       "\n",
       "Residuals:\n",
       "    Min      1Q  Median      3Q     Max \n",
       "-3.8569 -0.3503  0.0907  0.4391  1.3635 \n",
       "\n",
       "Coefficients:\n",
       "                  Estimate Std. Error t value Pr(>|t|)    \n",
       "(Intercept)       4.339577   0.128923   33.66   <2e-16 ***\n",
       "ldi              -0.158961   0.006572  -24.19   <2e-16 ***\n",
       "edu              -0.162916   0.002329  -69.94   <2e-16 ***\n",
       "mean_fertile_age -0.172521   0.004029  -42.82   <2e-16 ***\n",
       "---\n",
       "Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
       "\n",
       "Residual standard error: 0.6238 on 11971 degrees of freedom\n",
       "Multiple R-squared:  0.5409,\tAdjusted R-squared:  0.5407 \n",
       "F-statistic:  4700 on 3 and 11971 DF,  p-value: < 2.2e-16\n"
      ]
     },
     "execution_count": 15,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "summary(ldi.edu.mean.age.lm[[1]])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "$$log(ASFR) = \\alpha_0 + \\alpha_{country} + \\beta_1 LDI + \\beta_2 \\text{matern_edu} + \\beta_3 \\text{mean_fertile_age}$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 16,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "ldi.edu.country.mean.age.formula <- function(data) {\n",
    "    formula <- lmer('log(asfr) ~ 1 + ldi + edu + mean_fertile_age + (1 | location_id)',data=data)\n",
    "}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 17,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "ldi.edu.country.mean.age.lmer <- regression(ldi.edu.country.mean.age.formula)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 18,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "data": {
      "text/plain": [
       "Linear mixed model fit by REML ['lmerMod']\n",
       "Formula: log(asfr) ~ 1 + ldi + edu + mean_fertile_age + (1 | location_id)\n",
       "   Data: data\n",
       "\n",
       "REML criterion at convergence: 2578.6\n",
       "\n",
       "Scaled residuals: \n",
       "    Min      1Q  Median      3Q     Max \n",
       "-6.4446 -0.4733  0.0355  0.4847  4.9926 \n",
       "\n",
       "Random effects:\n",
       " Groups      Name        Variance Std.Dev.\n",
       " location_id (Intercept) 0.36115  0.601   \n",
       " Residual                0.06607  0.257   \n",
       "Number of obs: 11975, groups:  location_id, 188\n",
       "\n",
       "Fixed effects:\n",
       "                  Estimate Std. Error t value\n",
       "(Intercept)       3.509604   0.107345   32.69\n",
       "ldi              -0.047175   0.007511   -6.28\n",
       "edu              -0.144033   0.001521  -94.69\n",
       "mean_fertile_age -0.180393   0.002759  -65.38\n",
       "\n",
       "Correlation of Fixed Effects:\n",
       "            (Intr) ldi    edu   \n",
       "ldi         -0.541              \n",
       "edu          0.299 -0.729       \n",
       "mean_frtl_g -0.743  0.019  0.062"
      ]
     },
     "execution_count": 18,
     "metadata": {},
     "output_type": "execute_result"
    }
   ],
   "source": [
    "summary(ldi.edu.country.mean.age.lmer[[1]])"
   ]
  }
 ],
 "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": "3.2.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}