OpenSource/Haiti-WASH2016
1
0 Fork 0
mirror of https://github.com/flaneuse/Haiti-WASH2016.git synced 2026-10-26 17:35:24 +00:00
Code Issues Packages Projects Releases Wiki Activity
Files
master
Haiti-WASH2016/R/HT_05_improvedWater.R
T
flaneuse 7c6c7e7634 reset x limits
2016-10-26 13:46:43 -04:00

285 lines
13 KiB
R
Raw Permalink Blame History

# Haiti DHS Activity Design -----------------------------------------
#
# HT_04_improvedToilets: calculate improved latrine percentages
#
# Script to pull data from the Demographic and Health Surveys on variables
# related to water, sanitation, and hygiene (WASH)
#
# Data are from the 2012 DHS, available at http://dhsprogram.com/what-we-do/survey/survey-display-368.cfm
#
# Laura Hughes, lhughes@usaid.gov, 15 August 2016
# With Patrick Gault (pgault@usaid.gov) and Tim Essam (tessam@usaid.gov)
#
# Copyright 2016 by Laura Hughes via MIT License
#
# -------------------------------------------------------------------------
# Previous dependencies ---------------------------------------------------
# `HT_00_setupFncns.R`, `HT_01_importDHS_geo.R`, `HT_02_importDHS_hh.R` are meant to be run first.
# The following are dependencies in thoses files:
# * hh: dataframe with household-level indicators
# * admin1-2: shapefiles containing geographic polygons of Haiti
# Set colors and other basics ---------------------------------------------
colour_water = 'YlGnBu'
# classify improved/ not improved -----------------------------------------
# DHS claims to use WHO definitions for sanitation; similar to those provided by Dr. Elizabeth Jordan,
# WHO/UNICEF Joint Monitoring Programme for Water Supply and Sanitation: http://www.wssinfo.org/definitions-methods/watsan-categories/
# NOTE! DHS classifies bottled water (eau en bouteille, petit vendeur d'eau) as being improved.
# HOWEVER-- both WHO (JMP) and USAID WASH experts default to assuming the bottled water is unimproved unless proven otherwise.
# If the bottled water is consumed out of preference and the household has access to an improved water source for cooking and
# personal hygiene, the bottled water can be deemed "improved." In the absence of that issue, it's assumed to be unimproved.
# Unfortunately, the DHS didn't collect info on water sources for non-drinking water.
# Therefore, we are going AGAINST the DHS report and classifying bottled water as unimproved.
# -- WATER --
# Export type of water sources in survey.
water_types = attr(hh_raw$hv201, 'labels')
write.csv(water_types, '~/GitHub/Haiti-WASH2016/dataout/DHS_water_types.csv')
# read in Liz Jordan's classification of whether or not a toilet or water source is improved (see below in comments)
water_types = read.csv('~/GitHub/Haiti-WASH2016/dataout/DHS_water_classification.csv')
impr_water_codes = unlist(water_types %>% filter(improved == 1) %>% select(code))
unimpr_water_codes = unlist(water_types %>% filter(improved == 0) %>% select(code))
# To isolate the effect of the bottled water, create classificaiton if improved, unimproved, or bottled.
bottled_codes = c(71, 72)
unimpr_water_codes_bottled = setdiff(unimpr_water_codes, bottled_codes)
# To double check published stats, using the DHS-classification for bottled water (as improved).
impr_water_codes_dhs = c(impr_water_codes, bottled_codes)
unimpr_water_codes_dhs = setdiff(unimpr_water_codes, bottled_codes)
# Water sources are defined as being 'improved' if they are one of the following types and can be reached in <= 30 minutes
# watersource code improved
# piped water 10 1
# piped into dwelling 11 1
# piped to yard/plot 12 1
# piped from the neighbor 13 1
# public tap/standpipe 14 1
# tube well water 20 1
# tube well or borehole 21 1
# protected well to yard 31 1
# others protected well 32 1
# protected spring 41 1
# rainwater 51 1
# dug well (open/protected) 30 0
# unprotected well to yard 33 0
# public and others unprotected well 34 0
# surface water 40 0
# unprotected spring 42 0
# river/dam/lake/ponds/stream/canal/irrigation channel 43 0
# tanker truck 61 0
# cart with small tank 62 0
# bottled water 71 0
# sales company of water 72 0
# other 96 0
hh = hh %>%
mutate(
# -- straight classification of whether the source is improved --
improved_water = ifelse(water_source %in% impr_water_codes, 1,
ifelse(water_source %in% unimpr_water_codes, 0, NA)),
# -- improved source + <= 30 min. to acquire --
# time2water = 996 == "on premises"; assumed to be < 30 min.
impr_water_under30min = ifelse(is.na(time2water) | is.na(improved_water), NA,
ifelse(time2water <= 30 & improved_water == 1, 1, 0)),
# -- DHS classification --
impr_water_dhs = ifelse(water_source %in% impr_water_codes_dhs, 1,
ifelse(water_source %in% unimpr_water_codes_dhs, 0, NA)),
# -- bottled water --
impr_water_type = ifelse(water_source %in% impr_water_codes, 'improved',
ifelse(water_source %in% unimpr_water_codes_bottled, 'unimproved',
ifelse(water_source %in% bottled_codes, 'bottled',
NA)))
)
# Quick check bottled is correct:
# hh %>% group_by(impr_water_type) %>% summarise(n())
# Set up sampling weights --------------------------------------------------
DHSdesign = svydesign(id = ~prim_sampling_unit, strata = ~sample_strata, weights = ~sample_wt, data = hh)
summary(DHSdesign)
# Apply sampling weights --------------------------------------------------
# -- National-level stats (sans camps) --
water_natl = svymean(~impr_water_under30min, DHSdesign, na.rm = TRUE)
# Double-checking can replicate DHS numbers
# Off by ~ 0.1% because I removed NAs. Otherwise, seems to check out.
water_dhs = calcPtEst('impr_water_dhs', by_var = 'dhs_region', design = DHSdesign, df = hh)
# -- By Admin1 + Port-au-Prince --
water_admin1_PaP = calcPtEst('impr_water_under30min', by_var = 'region_name', design = DHSdesign, df = hh)
water_bottled_admin1_PaP = calcPtEst('impr_water_dhs', by_var = 'region_name', design = DHSdesign, df = hh)
# -- By Admin1 --
water_admin1 = calcPtEst('impr_water_under30min', by_var = 'admin1', design = DHSdesign, df = hh)
# -- By Admin2 --
water_admin2 = calcPtEst('impr_water_under30min', by_var = 'admin2', design = DHSdesign, df = hh)
# Merge in Admin1 names
water_admin2 = left_join(water_admin2, admin2_names, by = 'admin2')
# Merge in Admin1-level stats
water_admin1_sum = water_admin1 %>%
select(admin1, admin1_avg = avg, admin1_lb = lb, admin1_ub = ub)
water_admin2 = left_join(water_admin2, water_admin1_sum, by = 'admin1')
# Resort order of admin1 based on the admin1 averages. Best (most improved toilet) == top.
water_admin2$admin1 = factor(water_admin2$admin1, levels = water_admin1$admin1)
# -- By Admin3 --
water_admin3 = calcPtEst('impr_water_under30min', by_var = 'admin3', design = DHSdesign, df = hh)
# Admin2 dot plot ---------------------------------------------------------
pairGrid(water_admin2, savePlots = F, y_var = 'admin2')
# Admin 3 dot plot --------------------------------------------------------
pairGrid_admin3(water_admin3, savePlots = FALSE, y_var = 'admin3',
colorDot = colour_water, fill_limits = colour_limits) +
geom_vline(xintercept = water_natl[1], size = 0.25, colour = grey90K) +
annotate(geom = 'text', x = 0.25, y = 15, label = 'country average',
colour = grey90K,
family = 'Lato Light')
ggsave(filename = '~/Creative Cloud Files/MAV/Projects/Haiti_WASH-PAD_2016-09/exported_R/HTI_water_dotAdm3.pdf',
width = 10.75, height = 10, units = "in",
bg = "transparent",
scale = (8.9555/8.1219),
paper = "special", useDingbats = FALSE, compress = FALSE, dpi = 300)
# Admin1 map --------------------------------------------------------------
haiti_polygons = right_join(water_admin1_PaP, dhs_geo$df, by = c('region_name' = 'DHSREGFR'))
plotMap(haiti_polygons,
admin0 = hispaniola,
clipping_mask = admin0,
centroids = dhs_geo$centroids,
fill_var = 'impr_water_under30min',
fill_scale = colour_water,
fill_limits = colour_limits,
plot_base = FALSE,
exportPlot = TRUE,
fileName = '~/Creative Cloud Files/MAV/Projects/Haiti_WASH-PAD_2016-09/exported_R/HTI_imprwater30_adm1.pdf')
# Export data to krig surface ---------------------------------------------
# For krigging, unnecessary to apply weights to each EA (since same)
write_csv(hh, '~/Documents/USAID/Haiti/dataout/HT_DHS2012_imprsanitation_2016-08-29.csv')
# Admin2 map --------------------------------------------------------------
haiti_polygons = left_join(admin2$df, water_admin2, by = c('NAME_2' = 'admin2'))
# -- AVERAGE --
plotMap(haiti_polygons,
admin0 = hispaniola,
clipping_mask = admin0,
centroids = admin2$centroids,
fill_var = 'impr_water_under30min',
centroids_var = 'NAME_2',
fill_scale = colour_toilet,
fill_limits = colour_limits,
plot_base = TRUE,
exportPlot = FALSE)
# -- dot plot --
pairGrid(water_admin2,
y_var = 'admin2',
colorDot = colour_water,
savePlots = TRUE,
fill_limits = colour_limits,
file_name = '~/Creative Cloud Files/MAV/Projects/Haiti_WASH-PAD_2016-09/exported_R/HTI_imprwater_adm2dot.pdf')
# urban/rural splot -------------------------------------------------------
urb_rural_imprwater = calcPtEst('impr_water_under30min', by_var = 'region_urban', design = DHSdesign, df = hh)
urb_rural_imprwater = urb_rural_imprwater %>%
separate(region_urban, c('region', 'urban'), sep = '_and_') %>%
ungroup() %>%
arrange(region,urban) %>%
mutate(lagged = lag(avg))
water_avg = water_admin1_PaP %>% select(region = region_name, region_avg = avg)
urb_rural_imprwater = left_join(urb_rural_imprwater, water_avg, by = 'region')
ur_order_h2o = urb_rural_imprwater %>%
filter(urban == 'urban') %>%
arrange(region_avg)
urb_rural_imprwater$region = factor(urb_rural_imprwater$region, c(ur_order_h2o$region))
p = ur_pairGrid(urb_rural_imprwater, fill_scale = 'YlGnBu', fill_limits = colour_limits,
file_name = '~/Creative Cloud Files/MAV/Projects/HTI_WASH-PAD_2016-09/exported_R/HTI_UR_water30.pdf')
# export color bars -------------------------------------------------------
p + theme(legend.position = 'right', legend.direction = 'horizontal',
legend.text = element_text(colour = grey60K, family = 'Lato Light', size = 8)
) +
scale_fill_gradientn(colours = brewer.pal(9, colour_toilet),
limits = colour_limits,
labels = scales::percent)
ggsave(filename = '~/Creative Cloud Files/MAV/Projects/Haiti_WASH-PAD_2016-09/exported_R/HTI_RdPuscale.pdf',
width = 10.75, height = 9, units = "in",
bg = "transparent",
scale = (8.9555/8.1219),
paper = "special", useDingbats = FALSE, compress = FALSE, dpi = 300)
p + theme(legend.position = 'right', legend.direction = 'horizontal',
legend.text = element_text(colour = grey60K, family = 'Lato Light', size = 8)
) +
scale_fill_gradientn(colours = brewer.pal(9, colour_water),
limits = colour_limits,
labels = scales::percent)
ggsave(filename = '~/Creative Cloud Files/MAV/Projects/Haiti_WASH-PAD_2016-09/exported_R/HTI_YlGnBuscale.pdf',
width = 10.75, height = 9, units = "in",
bg = "transparent",
scale = (8.9555/8.1219),
paper = "special", useDingbats = FALSE, compress = FALSE, dpi = 300)
p + theme(legend.position = 'right', legend.direction = 'horizontal',
legend.text = element_text(colour = grey60K, family = 'Lato Light', size = 8)
) +
scale_fill_gradientn(colours = brewer.pal(9, colour_od),
limits = colour_limits,
labels = scales::percent)
ggsave(filename = '~/Creative Cloud Files/MAV/Projects/Haiti_WASH-PAD_2016-09/exported_R/HTI_YlOrBrscale.pdf',
width = 10.75, height = 9, units = "in",
bg = "transparent",
scale = (8.9555/8.1219),
paper = "special", useDingbats = FALSE, compress = FALSE, dpi = 300)
# Admin3 counts -----------------------------------------------------------
water_admin3_cts = hh %>%
filter(!is.na(impr_water_under30min)) %>%
group_by(admin3) %>% summarise(`households with improved sanitation` = sum(impr_water_under30min, na.rm = TRUE),
`total number of households` = n())
write.csv(water_admin3_cts, '~/Documents/USAID/Haiti/dataout/imprwater_admin3.csv')
Reference in new issue View Git Blame Copy Permalink