frankdavenport — Feb 10, 2014, 2:36 PM
#-------------------Base Setup--------------------------------------------------
rm(list=ls())
svnDir<-'/Users/frankdavenport/Education/R_Work/SVN/'
rcodeDir<-paste0(svnDir,'fd_FEWS/R_code/PAA_bweights/')
latexDir<-paste0(svnDir,'LaTex/')
rdatDir<-'~/Dropbox/FEWS_clim_low_bweight/'
dhsDir<-'~/Google Drive/dhs/'
chirpsDir<-'~/Google Drive/chirps/'
datDir<-paste0(dhsDir,'all/')
setwd(datDir)
library(plyr)
library(reshape2)
library(ggplot2)
library(scales)
library(car)
library(stringr)
library(lme4)
Loading required package: lattice
Loading required package: Matrix
Attaching package: 'lme4'
The following object is masked from 'package:ggplot2':
fortify
#library(dplyr)
#===============================================================================
#-----------------------------------Load the 'final' data.frame
load('06_childWithprecipandtemp.Rdata')
d<-subset(d,in_res_long=='yes')
#--------Subset by Women with >1 child, where at least 1 child is lbw
d$lbw<-ifelse(d$bw<2500,'low','normal')
#----Get cases where the woman has more than 1 child under 5 (a record in the child recode)
#PASTE ID TO CASE TO MAKE SURE THEY ARE UNIQUE
nchil<-as.data.frame(table(d$case))
nchil<-subset(nchil,Freq>1)
dm<-subset(d,case %in% nchil$Var1)
#----Now get cases where at least 1 child is lbw------
nlb<-droplevels(d[,c('.id','case','lbw')])
nlb<-melt(nlb,id.vars=c('.id','case'))
nlb<-dcast(nlb,case~value)
Aggregation function missing: defaulting to length
nlb<-subset(nlb,low>0 & normal>0) #at least 1 child in each category
nlb2<-subset(nlb,low==1 & normal==1) #exaclty 1 child in each category
#--
dm<-subset(dm,case %in% nlb2$case)
dm<-droplevels(dm)
#===============================================================
#----Examine Difference in Means Across Trimester----
tvars<-paste0('tmax.t',0:3)
pvars<-paste0('precip.t',0:3)
pd<-dm[,c('country','lbw',pvars,tvars)]
#------Overall----NOT MUCH DIFFERENCE
pdm<-melt(pd,id.vars=c('country','lbw'))
pdmm<-dcast(pdm,lbw~variable,mean)
pdms<-dcast(pdm,lbw~variable,sd)
pdmm #Overall Means
lbw precip.t0 precip.t1 precip.t2 precip.t3 tmax.t0 tmax.t1 tmax.t2
1 low 156.5 220.8 219.2 293.9 29.46 43.57 42.57
2 normal 148.2 219.9 223.7 305.0 28.72 43.13 42.87
tmax.t3
1 57.41
2 56.77
pdms #Overall Standard Devations
lbw precip.t0 precip.t1 precip.t2 precip.t3 tmax.t0 tmax.t1 tmax.t2
1 low 187.7 261.0 260.8 309.7 11.12 14.60 15.20
2 normal 183.5 250.9 270.0 331.0 11.17 15.34 14.83
tmax.t3
1 17.66
2 18.01
#--By Country-------
pdmm<-dcast(pdm,country+lbw~variable,mean)
pdms<-dcast(pdm,country+lbw~variable,sd)
#---Merge the 2
pdmm<-melt(pdmm,id.vars=c('country','lbw'),value.name='mean')
pdms<-melt(pdms,id.vars=c('country','lbw'),value.name='sd')
pd<-merge(pdmm,pdms)
pd$low<-pd$mean-(1.96*pd$mean)
pd$high<-pd$mean+(1.96*pd$mean)
pd$treatment<-str_sub(pd$variable,end=-4)
#-----All Countries
p0<-ggplot(pd,aes(x=mean,y=variable,color=lbw))+geom_point()+theme_bw()+facet_wrap(~treatment,nrow=1,scales='free')
p0
#-------Plot Temperature
pdt<-subset(pd,treatment=='tmax')
p1<-ggplot(data=pdt,aes(x=mean,y=variable,color=lbw))+geom_point()+theme_bw()+facet_wrap(~country)
p1<-p1+labs(x='Mean Count of Days Above Tmax',y='Terms',title='Tmax By Country')+scale_color_discrete(name='Birthweight')
p1
p1a<-ggplot(data=pdt,aes(x=mean,y=country,color=lbw))+geom_point()+theme_bw()+facet_wrap(~variable)
p1a<-p1a+labs(x='Mean Count of Days Above Tmax',y='Country',title='Tmax By terms')+scale_color_discrete(name='Birthweight')
p1a
#---Plot Precip-------
pdp<-subset(pd,treatment=='precip')
p2<-ggplot(data=pdp,aes(x=mean,y=variable,color=lbw))+geom_point()+theme_bw()+facet_wrap(~country)
p2<-p2+labs(x='Mean Total Precip',y='Country',title='Precip By Country')+scale_color_discrete(name='Birthweight')
p2
pdp<-subset(pd,treatment=='precip')
p2a<-ggplot(data=pdp,aes(x=mean,y=country,color=lbw))+geom_point()+theme_bw()+facet_wrap(~variable)
p2a<-p2a+labs(x='Mean Total Precip',y='Country',title='Precip By terms')+scale_color_discrete(name='Birthweight')
p2a
#------print them out----
setwd('~/Dropbox/Birthweight_temp_precip/')
pdf(file='MatchingonMothers_figs.pdf',onefile=T,width=7.5,height=8.5)
print(p0)
print(p1)
print(p1a)
print(p2)
print(p2a)
dev.off()
pdf
2