Skip to content

Instantly share code, notes, and snippets.

@ikashnitsky
Last active October 12, 2016 19:03
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save ikashnitsky/a578eaef6b122aa2aa2e3469fd2dcbe7 to your computer and use it in GitHub Desktop.
Save ikashnitsky/a578eaef6b122aa2aa2e3469fd2dcbe7 to your computer and use it in GitHub Desktop.
Create a plot showing sex ratios at all ages for all countries from Human Mortality Database
################################################################################
#
# Sex ratios 11-10-2016
# Sex ratios VS age in all HMD countries
# Ilya Kashnitsky, ilya.kashnitsky@gmail.com
#
################################################################################
# load required packages
require(dplyr) # version 0.5.0
require(readr) # version 1.0.0
require(tidyr) # version 0.6.0
require(ggplot2) # version 2.1.0
require(ggthemes) # version 3.2.0
library(extrafont) # version 0.17
library(HMDHFDplus) # version 1.1.8
# Download Exposures to risk
# read HMD data directly into R
# please note! the arguments "ik_user_hmd" and "ik_pass_hmd" are my login credidantials
# at the website of Human Mortality Database, which are stored locally at my computer.
# In order to access the data, you need to create an account at
# http://www.mortality.org/
# and provide your own credidantials to the `readHMDweb()` function
country <- getHMDcountries()
exposures <- list()
for (i in 1: length(country)) {
cnt <- country[i]
exposures[[paste0(cnt,'_exp')]] <- readHMDweb(cnt,"Exposures_1x1",ik_user_hmd,ik_pass_hmd)
# to see the progress
print(paste(i,'out of',length(country)))
}
# plot Sex Ratio VS age as lines for all countries
sr_age <- list()
for (i in 1:length(exposures)) {
di <- exposures[[i]]
sr_agei <- di %>% select(Year,Age,Female,Male) %>%
filter(Year %in% 2012) %>%
select(-Year) %>%
transmute(country = gsub('_exp','',names(exposures)[i]),
age = Age, sr_age = Male / Female *100)
sr_age[[i]] <- sr_agei
}
sr_age <- bind_rows(sr_age)
# remove optional populations
sr_age <- sr_age %>% filter(!country %in% c("FRACNP","DEUTE","DEUTW","GBRCENW","GBR_NP"))
# summarize all ages older than 90 (to jerky)
sr_age_90 <- sr_age %>% filter(age %in% 90:110) %>%
group_by(country) %>% summarise(sr_age = mean(sr_age, na.rm = T)) %>%
ungroup() %>% transmute(country, age=90, sr_age)
sr_age_0090 <- bind_rows(sr_age %>% filter(!age %in% 90:110), sr_age_90)
# finaly - plot
gg_sr_age <- ggplot(sr_age_0090, aes(age, sr_age, color = country, group = country))+
geom_hline(yintercept = 1, color = 'grey50', size = 1)+
geom_line(size = 1)+
scale_y_continuous(limits = c(0,120), expand = c(0,0), breaks = seq(0,120,20))+
scale_x_continuous(limits = c(0,90), expand = c(0,0), breaks = seq(0,80,20))+
xlab('Age')+
ylab('Sex ratio, males per 100 females')+
facet_wrap(~country,ncol=6)+
theme_few(base_family = "DejaVu Sans Condensed", base_size = 15)+
theme(legend.position='none')
ggsave('sr_lines.png', gg_sr_age, width = 8, height = 12, type="cairo-png")
@ikashnitsky
Copy link
Author

@timriffe Thanks a lot for HMDHFDplus

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment