Skip to content

Instantly share code, notes, and snippets.

@jiaweih
Created June 16, 2016 23:53
Show Gist options
  • Save jiaweih/6317f546858129e01d2e4f8bc67d6741 to your computer and use it in GitHub Desktop.
Save jiaweih/6317f546858129e01d2e4f8bc67d6741 to your computer and use it in GitHub Desktop.
{
"cells": [
{
"cell_type": "code",
"execution_count": 19,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"data <- read.csv('/ihme/forecasting/data/fertility/inputs/asfr_ldi_edu.csv')\n",
"data <- data[!(data$asfr==0),]\n",
"temp <- data[data$age_group_id==8,]\n",
"ldi.edu.lm <- lm('log(asfr) ~ ldi + edu',data=temp)"
]
},
{
"cell_type": "code",
"execution_count": 4,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/html": [
"<table>\n",
"<thead><tr><th></th><th scope=col>location_id</th><th scope=col>age_group_id</th><th scope=col>year_id</th><th scope=col>asfr</th><th scope=col>ldi</th><th scope=col>edu</th></tr></thead>\n",
"<tbody>\n",
"\t<tr><th scope=row>1</th><td>105</td><td>8</td><td>1951</td><td>0.1296</td><td>7.992734</td><td>5.097545</td></tr>\n",
"\t<tr><th scope=row>2</th><td>105</td><td>9</td><td>1951</td><td>0.2232</td><td>7.992734</td><td>5.097545</td></tr>\n",
"\t<tr><th scope=row>3</th><td>105</td><td>10</td><td>1951</td><td>0.2088</td><td>7.992734</td><td>4.822803</td></tr>\n",
"\t<tr><th scope=row>4</th><td>105</td><td>11</td><td>1951</td><td>0.1449</td><td>7.992734</td><td>4.822803</td></tr>\n",
"\t<tr><th scope=row>5</th><td>105</td><td>12</td><td>1951</td><td>0.1053</td><td>7.992734</td><td>4.474306</td></tr>\n",
"\t<tr><th scope=row>6</th><td>105</td><td>13</td><td>1951</td><td>0.0558</td><td>7.992734</td><td>4.474306</td></tr>\n",
"</tbody>\n",
"</table>\n"
],
"text/latex": [
"\\begin{tabular}{r|llllll}\n",
" & location_id & age_group_id & year_id & asfr & ldi & edu\\\\\n",
"\\hline\n",
"\t1 & 105 & 8 & 1951 & 0.1296 & 7.992734 & 5.097545\\\\\n",
"\t2 & 105 & 9 & 1951 & 0.2232 & 7.992734 & 5.097545\\\\\n",
"\t3 & 105 & 10 & 1951 & 0.2088 & 7.992734 & 4.822803\\\\\n",
"\t4 & 105 & 11 & 1951 & 0.1449 & 7.992734 & 4.822803\\\\\n",
"\t5 & 105 & 12 & 1951 & 0.1053 & 7.992734 & 4.474306\\\\\n",
"\t6 & 105 & 13 & 1951 & 0.0558 & 7.992734 & 4.474306\\\\\n",
"\\end{tabular}\n"
],
"text/plain": [
" location_id age_group_id year_id asfr ldi edu\n",
"1 105 8 1951 0.1296 7.992734 5.097545\n",
"2 105 9 1951 0.2232 7.992734 5.097545\n",
"3 105 10 1951 0.2088 7.992734 4.822803\n",
"4 105 11 1951 0.1449 7.992734 4.822803\n",
"5 105 12 1951 0.1053 7.992734 4.474306\n",
"6 105 13 1951 0.0558 7.992734 4.474306"
]
},
"execution_count": 4,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"head(data)"
]
},
{
"cell_type": "code",
"execution_count": 2,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/plain": [
"\n",
"Call:\n",
"lm(formula = \"log(asfr) ~ ldi + edu\", data = temp)\n",
"\n",
"Residuals:\n",
" Min 1Q Median 3Q Max \n",
"-3.9580 -0.3842 0.0986 0.4846 1.3845 \n",
"\n",
"Coefficients:\n",
" Estimate Std. Error t value Pr(>|t|) \n",
"(Intercept) -0.782794 0.051334 -15.25 <2e-16 ***\n",
"ldi -0.145784 0.007165 -20.34 <2e-16 ***\n",
"edu -0.139633 0.002493 -56.01 <2e-16 ***\n",
"---\n",
"Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1\n",
"\n",
"Residual standard error: 0.6883 on 12349 degrees of freedom\n",
"Multiple R-squared: 0.4563,\tAdjusted R-squared: 0.4562 \n",
"F-statistic: 5182 on 2 and 12349 DF, p-value: < 2.2e-16\n"
]
},
"execution_count": 2,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"summary(ldi.edu.lm)"
]
},
{
"cell_type": "code",
"execution_count": 3,
"metadata": {
"collapsed": false
},
"outputs": [
{
"name": "stderr",
"output_type": "stream",
"text": [
"Loading required package: Matrix\n"
]
}
],
"source": [
"library(lme4)"
]
},
{
"cell_type": "code",
"execution_count": 22,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"ldi.edu.lmer <- lmer('log(asfr) ~ 1 + ldi + edu + (1 | location_id)',data=temp)"
]
},
{
"cell_type": "code",
"execution_count": 23,
"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: temp\n",
"\n",
"REML criterion at convergence: 6320.3\n",
"\n",
"Scaled residuals: \n",
" Min 1Q Median 3Q Max \n",
"-5.6721 -0.4847 0.0381 0.5134 4.4876 \n",
"\n",
"Random effects:\n",
" Groups Name Variance Std.Dev.\n",
" location_id (Intercept) 0.40377 0.6354 \n",
" Residual 0.08923 0.2987 \n",
"Number of obs: 12352, groups: location_id, 193\n",
"\n",
"Fixed effects:\n",
" Estimate Std. Error t value\n",
"(Intercept) -1.662981 0.078984 -21.05\n",
"ldi -0.042999 0.008501 -5.06\n",
"edu -0.137635 0.001750 -78.66\n",
"\n",
"Correlation of Fixed Effects:\n",
" (Intr) ldi \n",
"ldi -0.810 \n",
"edu 0.541 -0.742"
]
},
"execution_count": 23,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"summary(ldi.edu.lmer)"
]
},
{
"cell_type": "code",
"execution_count": 24,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"mean.age <- read.csv('/ihme/forecasting/data/fertility/inputs/mean_fertile_age_by_asfr_ratio.csv')"
]
},
{
"cell_type": "code",
"execution_count": 28,
"metadata": {
"collapsed": false
},
"outputs": [
{
"data": {
"text/html": [
"<table>\n",
"<thead><tr><th></th><th scope=col>location_id</th><th scope=col>iso3</th><th scope=col>sex</th><th scope=col>year_id</th><th scope=col>mean_fertile_age</th></tr></thead>\n",
"<tbody>\n",
"\t<tr><th scope=row>1</th><td>6</td><td>CHN</td><td>2</td><td>1951</td><td>29.18958</td></tr>\n",
"\t<tr><th scope=row>2</th><td>6</td><td>CHN</td><td>2</td><td>1952</td><td>29.20997</td></tr>\n",
"\t<tr><th scope=row>3</th><td>6</td><td>CHN</td><td>2</td><td>1953</td><td>29.23125</td></tr>\n",
"\t<tr><th scope=row>4</th><td>6</td><td>CHN</td><td>2</td><td>1954</td><td>29.25348</td></tr>\n",
"\t<tr><th scope=row>5</th><td>6</td><td>CHN</td><td>2</td><td>1955</td><td>29.27673</td></tr>\n",
"\t<tr><th scope=row>6</th><td>6</td><td>CHN</td><td>2</td><td>1956</td><td>29.29278</td></tr>\n",
"</tbody>\n",
"</table>\n"
],
"text/latex": [
"\\begin{tabular}{r|lllll}\n",
" & location_id & iso3 & sex & year_id & mean_fertile_age\\\\\n",
"\\hline\n",
"\t1 & 6 & CHN & 2 & 1951 & 29.18958\\\\\n",
"\t2 & 6 & CHN & 2 & 1952 & 29.20997\\\\\n",
"\t3 & 6 & CHN & 2 & 1953 & 29.23125\\\\\n",
"\t4 & 6 & CHN & 2 & 1954 & 29.25348\\\\\n",
"\t5 & 6 & CHN & 2 & 1955 & 29.27673\\\\\n",
"\t6 & 6 & CHN & 2 & 1956 & 29.29278\\\\\n",
"\\end{tabular}\n"
],
"text/plain": [
" location_id iso3 sex year_id mean_fertile_age\n",
"1 6 CHN 2 1951 29.18958\n",
"2 6 CHN 2 1952 29.20997\n",
"3 6 CHN 2 1953 29.23125\n",
"4 6 CHN 2 1954 29.25348\n",
"5 6 CHN 2 1955 29.27673\n",
"6 6 CHN 2 1956 29.29278"
]
},
"execution_count": 28,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"head(mean.age)"
]
},
{
"cell_type": "code",
"execution_count": 26,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"library(plyr)"
]
},
{
"cell_type": "code",
"execution_count": 27,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"mean.age <- rename(mean.age,c('year' = 'year_id'))"
]
},
{
"cell_type": "code",
"execution_count": 29,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"data = merge.data.frame(x = data,y = mean.age,by = c('location_id','year_id'))"
]
},
{
"cell_type": "code",
"execution_count": 30,
"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>iso3</th><th scope=col>sex</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>CAN</td><td>2</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>CAN</td><td>2</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>CAN</td><td>2</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>CAN</td><td>2</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>CAN</td><td>2</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>CAN</td><td>2</td><td>27.82598</td></tr>\n",
"</tbody>\n",
"</table>\n"
],
"text/latex": [
"\\begin{tabular}{r|lllllllll}\n",
" & location_id & year_id & age_group_id & asfr & ldi & edu & iso3 & sex & mean_fertile_age\\\\\n",
"\\hline\n",
"\t1 & 101 & 1951 & 9 & 0.204566 & 9.36871 & 6.826699 & CAN & 2 & 27.82598\\\\\n",
"\t2 & 101 & 1951 & 13 & 0.03092 & 9.36871 & 6.414093 & CAN & 2 & 27.82598\\\\\n",
"\t3 & 101 & 1951 & 11 & 0.149179 & 9.36871 & 6.93765 & CAN & 2 & 27.82598\\\\\n",
"\t4 & 101 & 1951 & 10 & 0.2077373 & 9.36871 & 6.93765 & CAN & 2 & 27.82598\\\\\n",
"\t5 & 101 & 1951 & 14 & 0.0028385 & 9.36871 & 5.782532 & CAN & 2 & 27.82598\\\\\n",
"\t6 & 101 & 1951 & 12 & 0.0873706 & 9.36871 & 6.414093 & CAN & 2 & 27.82598\\\\\n",
"\\end{tabular}\n"
],
"text/plain": [
" location_id year_id age_group_id asfr ldi edu iso3 sex\n",
"1 101 1951 9 0.2045660 9.36871 6.826699 CAN 2\n",
"2 101 1951 13 0.0309200 9.36871 6.414093 CAN 2\n",
"3 101 1951 11 0.1491790 9.36871 6.937650 CAN 2\n",
"4 101 1951 10 0.2077373 9.36871 6.937650 CAN 2\n",
"5 101 1951 14 0.0028385 9.36871 5.782532 CAN 2\n",
"6 101 1951 12 0.0873706 9.36871 6.414093 CAN 2\n",
" mean_fertile_age\n",
"1 27.82598\n",
"2 27.82598\n",
"3 27.82598\n",
"4 27.82598\n",
"5 27.82598\n",
"6 27.82598"
]
},
"execution_count": 30,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"head(data)"
]
},
{
"cell_type": "code",
"execution_count": 32,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": [
"ldi.edu.age.lmer <- lmer('log(asfr) ~ 1 + ldi + edu + mean_fertile_age + (1 | location_id)',data=data)"
]
},
{
"cell_type": "code",
"execution_count": 33,
"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: 283414\n",
"\n",
"Scaled residuals: \n",
" Min 1Q Median 3Q Max \n",
"-5.7211 -0.4380 0.2407 0.6452 2.3737 \n",
"\n",
"Random effects:\n",
" Groups Name Variance Std.Dev.\n",
" location_id (Intercept) 0.4085 0.6392 \n",
" Residual 1.7105 1.3079 \n",
"Number of obs: 83715, groups: location_id, 188\n",
"\n",
"Fixed effects:\n",
" Estimate Std. Error t value\n",
"(Intercept) -0.825362 0.191263 -4.32\n",
"ldi -0.804326 0.013140 -61.21\n",
"edu 0.073875 0.002580 28.63\n",
"mean_fertile_age 0.158194 0.005268 30.03\n",
"\n",
"Correlation of Fixed Effects:\n",
" (Intr) ldi edu \n",
"ldi -0.571 \n",
"edu 0.293 -0.673 \n",
"mean_frtl_g -0.812 0.054 0.028"
]
},
"execution_count": 33,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"summary(ldi.edu.age.lmer)"
]
},
{
"cell_type": "code",
"execution_count": null,
"metadata": {
"collapsed": true
},
"outputs": [],
"source": []
},
{
"cell_type": "code",
"execution_count": 17,
"metadata": {
"collapsed": false
},
"outputs": [],
"source": [
"ldi.edu.lmer <- lmer('log(asfr) ~ (1 + ldi + edu | age_group_id)',data=data) "
]
},
{
"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 | age_group_id)\n",
" Data: data\n",
"\n",
"REML criterion at convergence: 157532\n",
"\n",
"Scaled residuals: \n",
" Min 1Q Median 3Q Max \n",
"-9.7681 -0.4219 0.0696 0.5225 5.8318 \n",
"\n",
"Random effects:\n",
" Groups Name Variance Std.Dev. Corr \n",
" age_group_id (Intercept) 0.28917 0.5377 \n",
" ldi 0.01901 0.1379 0.83 \n",
" edu 0.03392 0.1842 0.92 0.98\n",
" Residual 0.36259 0.6022 \n",
"Number of obs: 86319, groups: age_group_id, 7\n",
"\n",
"Fixed effects:\n",
" Estimate Std. Error t value\n",
"(Intercept) -0.46447 0.04826 -9.624"
]
},
"execution_count": 18,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"summary(ldi.edu.lmer)"
]
}
],
"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
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment