Example of calculating median and 95th percentiles for NHANES PFAS data

pfas_example.R

devtools::install_github("silentspringinstitute/RNHANES")
library(tidyverse)
library(RNHANES)

# Data frame of variables for the analysis
variables <- tribble(
  ~column,   ~comment_column, ~file_name, ~cycle,
  "LBXPFNA", "LBDPFNAL",      "PFAS_H",   "2013-2014",
  "LBXPFHS", "LBDPFHSL",      "PFAS_H",   "2013-2014",
  "LBXPFDE", "LBDPFDEL",      "PFAS_H",   "2013-2014",
  "LBXMPAH", "LBDMPAHL",      "PFAS_H",   "2013-2014",
  "LBXPFBS", "LBDPFBSL",      "PFAS_H",   "2013-2014",
  "LBXPFHP", "LBDPFHPL",      "PFAS_H",   "2013-2014",
  "LBXPFUA", "LBDPFUAL",      "PFAS_H",   "2013-2014",
  "LBXPFDO", "LBDPFDOL",      "PFAS_H",   "2013-2014",
  "SSNPFOA", "SDNPFOAL",      "SSPFAS_H", "2013-2014",
  "SSBPFOA", "SDBPFOAL",      "SSPFAS_H", "2013-2014",
  "SSNPFOS", "SDNPFOSL",      "SSPFAS_H", "2013-2014",
  "SSMPFOS", "SDMPFOSL",      "SSPFAS_H", "2013-2014"
)

# Download NHANES data
nhanes_data <- nhanes_load_data(variables$file_name, variables$cycle, demographics = TRUE)

#
# Combine the variants of PFOA and PFOS into one variable
# 

# Add together the branched and linear variants of PFOA
nhanes_data$SSPFAS_H$PFOA <- nhanes_data$SSPFAS_H$SSBPFOA + nhanes_data$SSPFAS_H$SSNPFOA
# If both branched AND linear were non-detects, marked the combined value as a non-detect
nhanes_data$SSPFAS_H$PFOAL <- ifelse(nhanes_data$SSPFAS_H$SDBPFOAL == 1 & nhanes_data$SSPFAS_H$SDNPFOAL == 1, 1, 0)

nhanes_data$SSPFAS_H$PFOS <- nhanes_data$SSPFAS_H$SSMPFOS + nhanes_data$SSPFAS_H$SSNPFOS
nhanes_data$SSPFAS_H$PFOSL <- ifelse(nhanes_data$SSPFAS_H$SDNPFOSL == 1 | nhanes_data$SSPFAS_H$SDMPFOSL == 1, 1, 0)

# Add on those new variables to the list of variables for analysis
variables <- bind_rows(variables, tribble(
  ~column, ~comment_column, ~file_name,   ~cycle,
  "PFOA", "PFOAL",          "SSPFAS_H",   "2013-2014",
  "PFOS", "PFOSL",          "SSPFAS_H",   "2013-2014"
))

#
# Calculate 50th and 95th percentiles for women age greater than or equal to 18 in NHANES
#
quantiles <- nhanes_quantile(nhanes_data, variables,
                             quantiles = c(0.05, 0.25, 0.5, 0.75, 0.95),
                             filter = RIAGENDR == 2 & RIDAGEYR >= 18)

#
# Detection frequencies for same group
# 
detection_frequencies <- nhanes_detection_frequency(nhanes_data,
                                                    variables,
                                                    filter = RIAGENDR == 2 & RIDAGEYR >= 18)

#
# Means
#
survey_mean <- function(x, design) {
  survey::svymean(x, design, na.rm = TRUE)
}
means <- nhanes_survey(survey_mean, nhanes_data, variables, filter = RIAGENDR == 2 & RIDAGEYR >= 18)

quantiles
detection_frequencies
means