Open Case Studies: Vaping Behaviors in American Youth

Disclaimer: The purpose of the Open Case Studies project is to demonstrate the use of various data science methods, tools, and software in the context of messy, real-world data. A given case study does not cover all aspects of the research process, is not claiming to be the most appropriate way to analyze a given data set, and should not be used in the context of making policy decisions without external consultation from scientific experts.

This work is licensed under the Creative Commons Attribution-NonCommercial 3.0 (CC BY-NC 3.0) United States License.

To cite this case study please use:

Wright, Carrie and Ontiveros, Michael and Meng, Qier and Jager, Leah and Taub, Margaret and Hicks, Stephanie. (2020). https://github.com/opencasestudies/ocs-bp-vaping-case-study. Vaping Behaviors in American Youth (Version v1.0.0).

To access the GitHub Repository with the data for this case study see here: https://github.com/opencasestudies/ocs-bp-vaping-case-study

You may also access and download the data using our OCSdata package. To learn more about this package including examples, see this link. Here is how you would install this package:

install.packages("OCSdata")

This case study is part of a series of public health case studies for the Bloomberg American Health Initiative.

The total reading time for this case study is calculated via koRpus and shown below:

Reading Time	Method
72 minutes	koRpus

Readability Score:

A readability index estimates the reading difficulty level of a particular text. Flesch-Kincaid, FORCAST, and SMOG are three common readability indices that were calculated for this case study via koRpus. These indices provide an estimation of the minimum reading level required to comprehend this case study by grade and age.

Text language: en

index	grade	age
Flesch-Kincaid	10	15
FORCAST	10	15
SMOG	12	17

Please help us by filling out our survey.

Motivation

In the United States, there have been significant and historical declines in cigarette smoking. In the 1970s, 75% of high school seniors were smoking, that number is below 10% now. This progress is largely due to the tobacco control movement and their focus on initiatives like ending advertising to children (like Joe Camel), passing indoor smoking laws, health communication, etc.

According to a recent report, overall tobacco/nicotine use increased in youths (middle school and high school students) in the United States in 2017 and 2018, despite previous years of declining use.

This major increase is attributed to an increase in the use of electronic cigarette (e-cigarette) products.

Forms of tobacco/nicotine include these categories:

Cigarette and other combustible tobacco (pipes, cigars, cigarillo, etc.)
E-cigarettes and vaporized tobacco/nicotine (hookah, Juul, etc.)
Other non-combustible, non-vapor tobacco/nicotine products (snus, chewing tobacco, etc.)

Youths are more likely to use e-cigarettes vs. combustible cigarettes these days, which is concerning because e-cigarettes are really efficient nicotine delivery devices that are reinforcing and easy to initiate. By contrast, it takes quite a while to become accustomed to cigarettes (eg. because of coughing) and become dependent. It is also harder for parents to detect e-cigarette use and intervene (eg. the smell is not as strong). This means that youths may be becoming physically dependent on nicotine more quickly than in past years, and that cessation services designed for youths will be needed.

Whereas in previous decades the focus was on advertising, the current era requires attention to the marketing broadly. Juul caught on through Instagram influencers. New policies that regulate these innovative marketing strategies will be critical.

E-cigarettes are referred to by many different names, including but not limited to:

Electronic nicotine delivery systems (ENDS)
Vapes
e-hookahs
vape pens
tanks
mods

The devices vary greatly:

[source]

See this CDC guide and the American Lung Association website for more information.

The report found that:

During 2017–2018, current use of any tobacco[/nicotine]\(^*\) product increased 38.3% (from 19.6% to 27.1%) among high school students and 28.6% (from 5.6% to 7.2%) among middle school students; e-cigarette use increased 77.8% (from 11.7% to 20.8%) among high school students and 48.5% (from 3.3% to 4.9%) among middle school students.

\(^*\) Note: we added “[/nicotine]” to this quote from the report because e-cigarettes deliver nicotine, but are not actually tobacco.

In 2018, the Federal Drug Administration (FDA) in the United States stated that e-cigarette usage use among youth reached:

“nothing short of an epidemic proportion of growth”

Additionally, you may learn more about e-cigarette or vaping use-associated lung injury (EVALI) from CDC’s Morbidity and Mortality Weekly Report (MMWR).

In this case study, we will be investigating the same data used in the report that generated the above findings. This data comes from the The National Youth Tobacco Survey (NYTS).

Gentzke, Andrea S., Melisa Creamer, Karen A. Cullen, Bridget K. Ambrose, Gordon Willis, Ahmed Jamal, and Brian A. King. “Vital Signs: Tobacco Product Use Among Middle and High School Students - United States, 2011-2018.” MMWR. Morbidity and Mortality Weekly Report 68 (6): 157–64 (2019).

Main Questions

Our main question:

How has tobacco and e-cigarette/vaping use by American youths changed since 2015?
How does e-cigarette use compare between males and females?
What vaping brands and flavors appear to be used the most frequently?
We will base this on the following survey questions:

“During the past 30 days, what brand of e-cigarettes did you usually use?”
“What flavors of tobacco products have you used in the past 30 days?”

Is there a relationship between e-cigarette/vaping use and other tobacco use?

Learning Objectives

In this case study, we will cover how to import data from multiple files efficiently, how to import data from excel files, and how to make a variety of visualizations to compare multiple groups across time. We will also demonstrate how to work with codebooks. We will cover the concept of survey weighting and introduce the srvyr package. We will discuss the difference between pooled cross-sectional data and panel data. We will especially focus on using packages and functions from the Tidyverse for wrangling data, such as tidyr and dplyr and for visualization, such as as ggplot2. The Tidyverse is a library of packages created by RStudio. While some students may be familiar with previous R programming packages, these packages make data science in R especially efficient.

The skills, methods, and concepts that students will be familiar with by the end of this case study are:

Data Science Learning Objectives:

Import data from Excel files
Merge data from multiple similar but not identical data structures
Create effective longitudinal data visualizations
Write functions in R
Apply functions across data subsets using purrr and dplyr functionality.

Statistical Learning Objectives:

Understanding of different types of longitudinal data
Usage of codebooks
Conceptual understanding of survey weighting
Implementing logistic regression with survey weighting

We will begin by loading the packages that we will need:

library(here)
library(readxl)
library(magrittr)
library(stringr)
library(purrr)
library(dplyr)
library(readr)
library(tidyr)
library(ggplot2)
library(scales)
library(viridis)
library(forcats)
library(naniar)
library(srvyr)
library(cowplot)
library(broom)
library(survey)
library(OCSdata)

Packages used in this case study:

Package	Use in this case study
here	to easily load and save data
readxl	to import the data in the excel files
magrittr	to use the compound assignment pipe operator `%<>%`
stringr	to manipulate the character strings within the data
purrr	to import the data in all the different excel and csv files efficiently
dplyr	to arrange/filter/select/compare specific subsets of the data
readr	to import the CSV file data
tidyr	to rearrange data in wide and long formats
ggplot2	to make visualizations with multiple layers
scales	to allow us to look at the colors within the viridis package
viridis	to make plots with a color palette that is compatible with color blindness
forcats	to allow for reordering of factors in plots
naniar	to make a visualization of missing data
syrvr	to use survey weights
cowplot	to allow plots to be combined
broom	to create nicely formatted model output
survey	to fit survey-weighted logistic regression
OCSdata	to access and download OCS data files

The first time we use a function, we will use the :: to indicate which package we are using. Unless we have overlapping function names, this is not necessary, but we will include it here to be informative about where the functions we will use come from.

Context

According to the cited Morbidity and Mortality Weekly Report this was what was already known about this topic and the implications of this study:

source

Importantly, the vapors used in e-cigarettes contain harmful chemicals:

source

E-cigarette usage has also been associated with lung injury:

source

See here for additional information about the potential health effects of e-cigarettes in teens and young adults.

Limitations

There are some important considerations regarding this data analysis to keep in mind:

The National Youth Tobacco Survey (NYTS) does not follow the same individual student respondents over time. A longitudinal study that does follow the same individuals over time collects data called panel data. The data in this study is called pooled cross-sectional data, and is obtained from random collection of observations across time.

According to Wikipedia:

Panel data differs from pooled cross-sectional data across time, because it deals with the observations on the same subjects in different times whereas the latter observes different subjects in different time periods

The data include percentages of student respondents reporting use of each particular tobacco product, but the survey questions did not ask the relative amount of use of one product compared to another. For example, the survey included questions like:“What flavors of tobacco products have you used in the past 30 days?” but did not ask how often one flavor was used by the same individual over another.

While gender and sex are not actually binary, the data used in this analysis only contain information for groups of individuals who answered the survey questions as male or female.

What are the data?

The data in this case study comes from the National Youth Tobacco Survey (NYTS) which is an annual survey that asks students in high school and middle school (grades 6-12) about tobacco usage in the United States of America.

The data for this survey is freely available online at this website with data from 1999, 2000, 2002, 2004, 2006, 2009, and 2011-2019. We will be using data from 2015-2019 due to the fact that these years are the most recent that asked questions regarding e-cigarette usage.

Each year includes documentation, such as a codebook and an excel file containing the data:

Therefore, since we are using data from 2015-2019, the data we are interested in are located in 5 different excel files, one for each year, each with its own codebook (this is the one for 2019).

The codebook contains information describing the data within the excel file.

As you can see the excel file contains very short variable names and numeric values, and it is not clear what they mean without the codebook:

The codebook explains what the variables (the columns) are:

And the codebook explains what the values for each variable are:

We will explain more later about what the values on the right indicate.

The reason that there are codebooks for each year is because the questions asked each year varied slightly.

The data in this survey are pooled cross-sectional data. In this study design, different subsets of student respondents are surveyed each year and it is not clear which, if any, individuals participate more than once from one year to the next.

Data Import

Reading in the excel files

Since these excel files are so large (each has roughly 20,000 rows), it takes a bit of time for the data to load. To make the process faster, we previously imported these files, selected only the columns relevant to our questions of interest, and saved these data subsets as comma-separated (.csv) files.

You may download these files using the OCSdata package:

# install.packages("OCSdata")

library(OCSdata)
simpler_import_data("ocs-bp-vaping-case-study", outpath = getwd())

You may also obtain these files from the GitHub Repository. If you have trouble accessing the GitHub Repository, these CSV files can be downloaded from here.

Click here for details on how the data were originally imported

If you have trouble accessing the GitHub Repository, you can follow the links to download the excel files for 2015, 2016, 2017, 2018, and 2019.

First we created a vector of file names of all the different excel files. Using the here() function of the here package, we looked in all the directories of the project. The list.files() function looked for all files with .xlsx within these sub-directories.

excel_files <- list.files(here::here(), recursive = TRUE,
                  pattern = "*.xlsx")
excel_files

[1] "data/raw/2015-nyts-dataset-and-codebook-microsoft-excel/nyts2015.xlsx"
[2] "data/raw/2016-nyts-dataset-and-codebook-microsoft-excel/nyts2016.xlsx"
[3] "data/raw/2017-nyts-dataset-and-codebook-microsoft-excel/nyts2017.xlsx"
[4] "data/raw/2018-nyts-dataset-and-codebook-microsoft-excel/nyts2018.xlsx"
[5] "data/raw/2019-nyts-dataset-and-codebook-microsoft-excel/nyts2019.xlsx"

All the files were read using read_excel() of the readxl package. Using the map() function of the purrr package this was done efficiently for all of the excel files in the vector using one command. The . is used to indicate that we want to apply the read_excel() function to each element of the data that we just piped into the map() function.

Here we also used the %>% pipe which can be used to pass the input from one function to another for later sequential steps. You can use pipes if you load the dplyr package or the magrittr package.

This created a single list of tibbles (one for each file).

tbl_files <- excel_files %>%
       map(~ readxl::read_excel(.))

The elements of this list are in the same order as the elements of the excel_files vector, so we can extract the data for a particular file (year) by selecting a certain element of the list. However, it is safer to be able to select the data from this list for a specific year using a name based on the original vector of file names. We extracted a name from each Excel file name using the str_extract() function of the stringr package. Here we are keeping occurrences of the character string “nyts201” followed by a “5”, “6”, “7”, “8”, or “9”.

tbl_names <- excel_files %>%
  str_extract("nyts201[5-9]")

tbl_names

[1] "nyts2015" "nyts2016" "nyts2017" "nyts2018" "nyts2019"

These names became the names of the tibbles in the list of tibbles.

names(tbl_files) <- tbl_names

Specific columns were selected using the select() function of dplyr from each of the tibbles using the variable name, as identified in the codebook as being of interest for our analysis. In some cases functions like starts_with() of the dplyr package were used to select several variables at once. Most of the survey questions about tobacco use start with an "E" or a "C" according to the codebooks.

tbl_files[["nyts2015"]] <- tbl_files[["nyts2015"]] %>%
    dplyr::select(psu, # Primary Sampling Unit
                  finwgt, # Analysis Weight
                  stratum, # Sampling stratum
                  Qn1, # Age
                  Qn2, # Sex
                  Qn3, # Grade
                  starts_with("E",
                              ignore.case = FALSE),
                  starts_with("C",
                              ignore.case = FALSE),
                  )


tbl_files[["nyts2016"]] <- tbl_files[["nyts2016"]] %>%
    dplyr::select(psu,
                  finwgt,
                  stratum,
                  Q1, # Age
                  Q2, # Sex
                  Q3, # Grade
                  starts_with("E",
                              ignore.case = FALSE),
                  starts_with("C",
                              ignore.case = FALSE),
                  Q50A, # Menthol # What flavors of tobacco products have you used in the past 30 days? (Select one or more)
                  Q50B, # Clove or spice
                  Q50C, # Fruit
                  Q50D, # Chocolate
                  Q50E, # Alcoholic Drink
                  Q50F, # Candy/Desserts/Other Sweets
                  Q50G  # Some Other Flavor Not Listed Here
                  )

tbl_files[["nyts2017"]] <- tbl_files[["nyts2017"]] %>%
    dplyr::select(psu,
                  finwgt,
                  stratum,
                  Q1, # Age
                  Q2, # Sex
                  Q3, # Grade
                  starts_with("E",
                              ignore.case = FALSE),
                  starts_with("C",
                              ignore.case = FALSE),
                  Q50A, # Menthol # What flavors of tobacco products have you used in the past 30 days? (Select one or more)
                  Q50B, # Clove or spice
                  Q50C, # Fruit
                  Q50D, # Chocolate
                  Q50E, # Alcoholic Drink
                  Q50F, # Candy/Desserts/Other Sweets
                  Q50G  # Some Other Flavor Not Listed Here
                  )

tbl_files[["nyts2018"]] <- tbl_files[["nyts2018"]] %>%
    dplyr::select(psu,
                  finwgt,
                  stratum,
                  Q1, # Age
                  Q2, # Sex
                  Q3, # Grade
                  starts_with("E",
                              ignore.case = FALSE),
                  starts_with("C",
                              ignore.case = FALSE),
                  Q50A, # Menthol #What flavors of tobacco products have you used in the past 30 days? (Select one or more)
                  Q50B, # Clove or spice
                  Q50C, # Fruit
                  Q50D, # Chocolate
                  Q50E, # Alcoholic Drink
                  Q50F, # Candy/Desserts/Other Sweets
                  Q50G  # Some Other Flavor Not Listed Here
                  )

tbl_files[["nyts2019"]] <- tbl_files[["nyts2019"]] %>%
    dplyr::select(psu,
                  finwgt,
                  stratum,
                  Q1, # Age
                  Q2, # Sex
                  Q3, # Grade
                  starts_with("E",
                              ignore.case = FALSE),
                  starts_with("C",
                              ignore.case = FALSE),
                  Q40, # Brand, e-cigarettes
                  Q62A, # Menthol # What flavors of tobacco products have you used in the past 30 days? (Select one or more)
                  Q62B, # Clove or spice
                  Q62C, # Fruit
                  Q62D, # Chocolate
                  Q62E, # Alcoholic Drink
                  Q62F, # Candy/Desserts/Other Sweets
                  Q62G, # Some Other Flavor Not Listed Here
                  )

A directory called reduced was created for the new .csv files using the base dir.create() function. A csv file was created for each of the tibbles in the list using the write_csv() function of the readr package. This was done all at once using the map() function.

dir.create("data/reduced")

names(tbl_files) %>% map(~ write_csv(tbl_files[[.]], 
                                     path = paste0("data/reduced/", ., ".csv")))

Now we will show how to read in the data from the five csv files that were created from the five different Excel files.

Reading in the CSV files

Using the here() function of the here package, we looked in all the directories of the project. The here package automatically starts looking for files based on where you have a .Rproj file which is created when you start a new RStudio project.

Click here to see more about creating new projects in RStudio.

You can create a project by going to the File menu of RStudio like so:

You can also do so by clicking the project button:

See here to learn more about using RStudio projects.

The list.files() function looked for all files with .csv within the data/reduced sub-directories.

All the files were read using read_csv() of the readr package. Using the map() function of the purrr package this was done efficiently for all of the CSV files in the vector using one command. The . is used to indicate that we want to apply the read_csv() function to each element of the data that we just piped into the map() function. For more information about the map() function, see here.

Here we also used the %>% pipe which can be used to pass the input from one function to another for later sequential steps.

nyts_data <- list.files(here::here(), recursive = TRUE,
                        pattern = "*.csv") %>%
  stringr::str_subset(., "wrangled", negate = TRUE) %>%
  map(~ read_csv(.))


nyts_data_names <- list.files(recursive = TRUE,
                              pattern = "*.csv") %>%
  stringr::str_subset(., "wrangled", negate = TRUE) %>%
  str_extract("nyts201[5-9]")

names(nyts_data) <- nyts_data_names

We can save our imported data as an rda file (stands for R data file) using the save() function.

save(nyts_data, file = here::here("data", "imported", "imported_data.rda"))

Data Exploration and Wrangling

If you have been following along but stopped, we could load our imported data from the “data” directory like so:

load(here::here("data", "imported", "imported_data.rda"))

If you skipped the data import section click here.

First you need to install and load the OCSdata package:

install.packages("OCSdata")
library(OCSdata)

Then, you may load the imported data using the following code:

imported_data("ocs-bp-vaping-case-study", outpath = getwd())
load(here::here("OCSdata", "data", "imported", "imported_data.rda"))

If the package does not work for you, alternatively, an RDA file (stands for R data) of the data can be found in our GitHub repository or slightly more directly here. Download this file and then place it in your current working directory within a subdirectory called “imported” within a subdirectory called “data” to copy and paste our code. We used an RStudio project and the here package to navigate to the file more easily.

load(here::here("data", "imported", "imported_data.rda"))

Click here to see more about creating new projects in RStudio.

You can create a project by going to the File menu of RStudio like so:

You can also do so by clicking the project button:

See here to learn more about using RStudio projects and here to learn more about the here package.

Our goal in this section is to adjust or “wrangle” the data from each year into a common format so that we can combine the datasets across years for our analysis, and so that we have values in our variables that are correct and easy to interpret. We will need to understand what is the same and what is different across the data from different years, rename and recode the variables (e.g., by replacing the numbers 1 and 2 with the values “Male” and “Female” for the Sex variable), and combine the data. We will walk through these steps below.

First, let’s take a look at our data. We can get a good sense of it using the glimpse() function of the dplyr package.

dplyr::glimpse(nyts_data[["nyts2015"]])

Rows: 17,711
Columns: 29
$ psu        <chr> "015438", "015438", "015438", "015438", "015438", "015438",…
$ finwgt     <dbl> 216.7268, 324.9620, 324.9620, 397.1552, 264.8745, 264.8745,…
$ stratum    <chr> "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", "BR…
$ Qn1        <dbl> 10, 9, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, 10, …
$ Qn2        <dbl> 2, 1, 1, 1, 2, 2, 1, 2, 1, 2, 2, 2, 1, 2, 2, 1, 2, 2, 2, 2,…
$ Qn3        <dbl> 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 7, 5, 5,…
$ ECIGT      <dbl> 2, 1, 2, 1, 2, 1, 1, 1, 2, 2, 2, 2, 1, 2, 2, 1, 2, 1, 2, 1,…
$ ECIGAR     <dbl> 1, 1, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 1, 1, 2,…
$ ESLT       <dbl> 2, 2, 2, 2, 2, 2, 1, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ EELCIGT    <dbl> 2, 1, 2, 1, 2, 1, 1, 1, 2, 2, 2, 1, 2, 2, 2, 1, 2, 2, 2, 1,…
$ EROLLCIGTS <dbl> 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2,…
$ EFLAVCIGTS <dbl> 2, 2, 2, 1, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ EBIDIS     <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ EFLAVCIGAR <dbl> 2, 1, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 1, 2, 1, 1, 2,…
$ EHOOKAH    <dbl> 2, 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ EPIPE      <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ ESNUS      <dbl> 2, 2, 2, 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ EDISSOLV   <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CCIGT      <dbl> 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CCIGAR     <dbl> 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CSLT       <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CELCIGT    <dbl> 2, 2, 2, 1, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CROLLCIGTS <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CFLAVCIGTS <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CBIDIS     <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CHOOKAH    <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CPIPE      <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CSNUS      <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…
$ CDISSOLV   <dbl> 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2,…

Updating the set of variables and their names

The easiest way of making it so that the data from the different years can be combined is by making sure that the same type of data in different datasets share the same names. In addition, giving the columns informative names will help make our code more readable. Currently, it isn’t very clear what most of the variables indicate since the variable names are uninformative on their own, without the codebook.

We want to rename variables like Qn1 to something more meaningful like Age.

To do this we will use the rename() function of the dplyr package. The new name is always listed first before the =. This function will replace the old variable names with the new ones, i.e., after running the code below, there will no longer be a Qn1 variable in the dataset, but there will be an Age variable instead. We will start working with the 2015 data, and then move on to the other years down below.

nyts_data[["nyts2015"]] <- nyts_data[["nyts2015"]] %>%
  dplyr::rename(Age = Qn1,
                Sex = Qn2,
                Grade = Qn3)

Ultimately we will be combining the data from each year using the bind_rows() function of the dplyr package, which will fill in NA values for any columns that do not exist for one of the years.

The data for 2016-2018 have many common attributes, so we will want to write code that can be applied to all three datasets. To do this, we will use a function in R, which is basically a piece of code that can be applied to similar but different objects in R (e.g., the data tibbles from each of these three years). Recall that you are using functions from packages all the time, like the rename() function of the dplyr package. Now you are going to write your own function! For more information on functions, see here.

These next 3 years have the same structure for many of the questions we are interested in. For example, they all have flavor questions, but not a brand question. Moreover, their variable names are consistent across the years; for each year, we want to replace the vague question name Q50A with the value menthol in all three datasets, and the same is true for the other flavor variables.

Since we want to perform the same modifications on the data from all three years, rather than repeating the same somewhat messy piece of code three times, we can do this more efficiently if we create a function to do all of these steps at once. Then we can use the map_at() function of the purrr package, which is an extension of the map() function that we used in the [Data Import] (https://www.opencasestudies.org/ocs-bp-vaping-case-study/#Data_Import){target=“_blank”} section to apply the function we just created (for renaming variables etc.) specifically to the data from 2016-2018 within the nyts_data. By using vars() inside of the map_at() function we can specify what tibbles within our nyts_data list we want to include or exclude.

Functions

So how do you create a function? You first need to specify that you are creating a function by using the function() base function. Yes, that’s right it is a function for creating functions called function!

First we specify our input within the parentheses of function(). Thus if our function will apply something to an input called x then we would use function(x). Theoretically our input can be named whatever we want, but we we need to refer to it consistently within our function to indicate what we want done with the input data. We can actually have more than one input as well, we would indicate two inputs like this: function(x, b). Here we would be using both x and b to do something in our function.

The next part of a function is defined within curly brackets {}. This is where we write what we want done to our inputs.

Click here to see a simple example

Our function will be called simple_function, which will add 2 to our input. It will take the input vector_data (note that we usually want to use an informative input argument like vector_data).

simple_function <- function(vector_data) {
  vector_data + 2
}

We can use an input that is a vector called x.

x = c(1, 2, 3, 4)
x

[1] 1 2 3 4

Now we specify using the argument vector_data = to indicate that x is our input that we want to perform the function on.

simple_function(vector_data = x)

[1] 3 4 5 6

This function also works for other vector input. For example, vector y is now our input that we want to perform the function on.

y = c(5, 10, 15, 20)

simple_function(vector_data = y)

[1]  7 12 17 22

In our case we will be applying our function to the variable names for the dataset for each year. The output of our function is the result of renaming these variables for each year.

update_survey <- function(dataset) { dataset %>%
           rename(Age = Q1,
                  Sex = Q2,
                Grade = Q3,
              menthol = Q50A,
          clove_spice = Q50B,
                fruit = Q50C,
            chocolate = Q50D,
      alcoholic_drink = Q50E,
 candy_dessert_sweets = Q50F,
                other = Q50G)
}

# options to apply the function to the data:
# nyts_data <- nyts_data %>% map_at(vars(nyts2016, nyts2017, nyts2018), Update_survey)
nyts_data <- nyts_data %>% map_at(vars(-nyts2015, -nyts2019), update_survey)

The final year, 2019, has a slightly different data structure compared to these earlier datasets. It actually has a brand_ecig variable already and different question numbers correspond to our flavor questions of interest. So we will rename the variables in this dataset individually. We could also write this as a function, but since we are only applying this one time, there is no need to. Functions are really helpful for repeating the same task repeatedly using different data inputs.

nyts_data[["nyts2019"]] <- nyts_data[["nyts2019"]] %>%
    rename(brand_ecig = Q40,
                  Age = Q1,
                  Sex = Q2,
                Grade = Q3,
              menthol = Q62A,
          clove_spice = Q62B,
                fruit = Q62C,
            chocolate = Q62D,
      alcoholic_drink = Q62E,
 candy_dessert_sweets = Q62F,
                other = Q62G)

Now let’s take a look at the variable names for each of the years using the map function from purrr.

map(nyts_data, names)

$nyts2015
 [1] "psu"        "finwgt"     "stratum"    "Age"        "Sex"       
 [6] "Grade"      "ECIGT"      "ECIGAR"     "ESLT"       "EELCIGT"   
[11] "EROLLCIGTS" "EFLAVCIGTS" "EBIDIS"     "EFLAVCIGAR" "EHOOKAH"   
[16] "EPIPE"      "ESNUS"      "EDISSOLV"   "CCIGT"      "CCIGAR"    
[21] "CSLT"       "CELCIGT"    "CROLLCIGTS" "CFLAVCIGTS" "CBIDIS"    
[26] "CHOOKAH"    "CPIPE"      "CSNUS"      "CDISSOLV"  

$nyts2016
 [1] "psu"                  "finwgt"               "stratum"             
 [4] "Age"                  "Sex"                  "Grade"               
 [7] "ECIGT"                "ECIGAR"               "ESLT"                
[10] "EELCIGT"              "EHOOKAH"              "EROLLCIGTS"          
[13] "EFLAVCIGAR"           "EPIPE"                "ESNUS"               
[16] "EDISSOLV"             "EBIDIS"               "CCIGT"               
[19] "CCIGAR"               "CSLT"                 "CELCIGT"             
[22] "CHOOKAH"              "CROLLCIGTS"           "CPIPE"               
[25] "CSNUS"                "CDISSOLV"             "CBIDIS"              
[28] "menthol"              "clove_spice"          "fruit"               
[31] "chocolate"            "alcoholic_drink"      "candy_dessert_sweets"
[34] "other"               

$nyts2017
 [1] "psu"                  "finwgt"               "stratum"             
 [4] "Age"                  "Sex"                  "Grade"               
 [7] "ECIGT"                "ECIGAR"               "ESLT"                
[10] "EELCIGT"              "EHOOKAH"              "EROLLCIGTS"          
[13] "EPIPE"                "ESNUS"                "EDISSOLV"            
[16] "EBIDIS"               "CCIGT"                "CCIGAR"              
[19] "CSLT"                 "CELCIGT"              "CHOOKAH"             
[22] "CROLLCIGTS"           "CPIPE"                "CSNUS"               
[25] "CDISSOLV"             "CBIDIS"               "menthol"             
[28] "clove_spice"          "fruit"                "chocolate"           
[31] "alcoholic_drink"      "candy_dessert_sweets" "other"               

$nyts2018
 [1] "psu"                  "finwgt"               "stratum"             
 [4] "Age"                  "Sex"                  "Grade"               
 [7] "ECIGT"                "ECIGAR"               "ESLT"                
[10] "EELCIGT"              "EHOOKAH"              "EROLLCIGTS"          
[13] "EPIPE"                "ESNUS"                "EDISSOLV"            
[16] "EBIDIS"               "CCIGT"                "CCIGAR"              
[19] "CSLT"                 "CELCIGT"              "CHOOKAH"             
[22] "CROLLCIGTS"           "CPIPE"                "CSNUS"               
[25] "CDISSOLV"             "CBIDIS"               "menthol"             
[28] "clove_spice"          "fruit"                "chocolate"           
[31] "alcoholic_drink"      "candy_dessert_sweets" "other"               

$nyts2019
 [1] "psu"                  "finwgt"               "stratum"             
 [4] "Age"                  "Sex"                  "Grade"               
 [7] "ECIGT"                "ECIGAR"               "ESLT"                
[10] "EELCIGT"              "EHOOKAH"              "EROLLCIGTS"          
[13] "EPIPE"                "ESNUS"                "EDISSOLV"            
[16] "EBIDIS"               "EHTP"                 "CCIGT"               
[19] "CCIGAR"               "CSLT"                 "CELCIGT"             
[22] "CHOOKAH"              "CROLLCIGTS"           "CPIPE"               
[25] "CSNUS"                "CDISSOLV"             "CBIDIS"              
[28] "CHTP"                 "brand_ecig"           "menthol"             
[31] "clove_spice"          "fruit"                "chocolate"           
[34] "alcoholic_drink"      "candy_dessert_sweets" "other"

It’s looking better! The data that overlap across years have the same variable names.

Updating Values

Now that we have made some progress on the selection and names of the variables themselves, we will work on the values contained in the different variables.

We can start with updating the values for Age and Grade, so that they are more understandable.

Recall from the codebook for this year’s dataset that Age isn’t listed in the way one might expect, i.e., it is not just a number of years, but a numerically valued categorical variable.

The same is true for Grade:

This is why it is so important to always check the codebook!!

We also want to replace the value of 19 for Age to be ">18" and the value of 13 for Grade to be replaced with "Ungraded/Other" Also according to the codebooks, numeric values of 1 indicate a survey answer of FALSE, while a value of 2 indicates TRUE. Sex also needs to be recoded. If we take a look at the codebooks carefully (make sure you look at the questions that we pulled, not the recoded values), we will see that males are indicated by 1 and females are indicated by 2. Finally some values are indicated with "*" or"**" when they are missing. We want to replace these with NA.

Let’s create a function to make all these updates. We will use the mutate function of the dplyr package to modify these variables. This function can also be used to create new variables. We will also use the recode() function of the dplyr package to replace specific values of certain variables.

update_values <- function(dataset) {
  dataset %>%
    dplyr::mutate(Age = as.numeric(Age) + 8,
                  Grade = as.numeric(Grade) + 5) %>%
    mutate(Age = as.factor(Age),
           Grade = as.factor(Grade),
           Sex = as.factor(Sex)) %>%
    mutate(Sex = recode(Sex,
                        `1` = "male",
                        `2` = "female")) %>%
    dplyr::mutate_all( ~ replace(., . %in% c("*", "**"), NA)) %>%
    mutate(Age = recode(Age,
                        `19` = ">18"),
           Grade = recode(Grade,
                          `13` = "Ungraded/Other")) %>%
    dplyr::mutate_at(vars(
      starts_with("E", ignore.case = FALSE),
      starts_with("C", ignore.case = FALSE)
    ),
    list( ~ recode(
      .,
      `1` = TRUE,
      `2` = FALSE,
      .default = NA,
      .missing = NA
    )))
}

nyts_data <- nyts_data %>% map(., update_values)

Now if we wanted to check that everything is expected we could do something like this to check the Sex variable using the count() function of the dplyr package. It is advisable to check your data frequently to make sure that it is as expected!

According to the codebook, we should have:

8,958 males in 2015
10,438 males in 2016
8,881 males in 2017
10,069 males in 2018
9,803 males in 2019

count_sex <- function(dataset) {dataset %>% count(Sex)}
nyts_data %>% map(., count_sex)

$nyts2015
# A tibble: 3 × 2
  Sex        n
  <fct>  <int>
1 male    8958
2 female  8622
3 <NA>     131

$nyts2016
# A tibble: 3 × 2
  Sex        n
  <fct>  <int>
1 male   10438
2 female 10082
3 <NA>     155

$nyts2017
# A tibble: 3 × 2
  Sex        n
  <fct>  <int>
1 male    8881
2 female  8815
3 <NA>     176

$nyts2018
# A tibble: 3 × 2
  Sex        n
  <fct>  <int>
1 male   10069
2 female  9920
3 <NA>     200

$nyts2019
# A tibble: 3 × 2
  Sex        n
  <fct>  <int>
1 .N       116
2 male    9803
3 female  9099

Looks good!

The years (2016-2019) that have flavors also need the flavor data to be logical (meaning TRUE or FALSE):

In this case we also are setting missing values to FALSE because then the TRUE values will represent those who reported using a specific flavor out of all users, rather than those that used a specific flavor compared to those who used a different flavor.

update_flavors <- function(dataset) {dataset %>%
   mutate_at(vars(menthol:other),
              list(~ recode(.,
                           `1` = TRUE,
                      .default = FALSE,
                      .missing = FALSE))) }

nyts_data  <- nyts_data  %>% map_at(vars(-nyts2015), update_flavors)

Now there are just a few changes needed that are specific to 2019. Specifically, some of the 2019 questions use the values “.N”, “.S”, and “.Z” to indicate different types of missing data (see for example Q2 of the 2019 codebook); we just want them to be replaced with NA values.

nyts_data[["nyts2019"]] <- nyts_data[["nyts2019"]] %>%
  mutate_all(~ replace(., . %in% c(".N", ".S", ".Z"), NA)) %>%
  mutate(psu = as.character(psu)) %>%
  mutate(brand_ecig = recode(brand_ecig,
                                             `1` = "Other", # levels 1,8 combined to `Other`
                                             `2` = "Blu",
                                             `3` = "JUUL",
                                             `4` = "Logic",
                                             `5` = "MarkTen",
                                             `6` = "NJOY",
                                             `7` = "Vuse",
                                             `8` = "Other"))

Great! Now our values don’t need to be handled any differently for any of the years, thus we can combine the data across years.

Even though we have different numbers of variables for each year, we can coerce the data to be combined into one tibble by using the bind_rows() function of dplyr. Importantly, this function does not require that the columns be the same. This will create NA values for any variable that is not present in given data frame but is present in one of the other data frames that is being combined. Note that the bind_cols() function does expect that the rows match. The .id argument will create a new variable with values to link each row to its original data frame. For more information see here.

nyts_data <- nyts_data %>%
  map_df(dplyr::bind_rows, .id = "year")

glimpse(nyts_data)

Rows: 95,465
Columns: 40
$ year                 <chr> "nyts2015", "nyts2015", "nyts2015", "nyts2015", "…
$ psu                  <chr> "015438", "015438", "015438", "015438", "015438",…
$ finwgt               <dbl> 216.7268, 324.9620, 324.9620, 397.1552, 264.8745,…
$ stratum              <chr> "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", …
$ Age                  <fct> 18, 17, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 1…
$ Sex                  <fct> female, male, male, male, female, female, male, f…
$ Grade                <fct> 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 1…
$ ECIGT                <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ ECIGAR               <lgl> TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, FAL…
$ ESLT                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, T…
$ EELCIGT              <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ EROLLCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, F…
$ EFLAVCIGTS           <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, TRUE, FA…
$ EBIDIS               <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ EFLAVCIGAR           <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, FA…
$ EHOOKAH              <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ EPIPE                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ ESNUS                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, F…
$ EDISSOLV             <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CCIGT                <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, F…
$ CCIGAR               <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, F…
$ CSLT                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CELCIGT              <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, F…
$ CROLLCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CFLAVCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CBIDIS               <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CHOOKAH              <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CPIPE                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CSNUS                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CDISSOLV             <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ menthol              <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ clove_spice          <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ fruit                <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ chocolate            <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ alcoholic_drink      <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ candy_dessert_sweets <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ other                <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ EHTP                 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ CHTP                 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ brand_ecig           <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…

We will want to do some of our analysis split by year, so we would like to be sure we have one variable that has the correct value for year. It looks like we just need to remove "nyts" from the year variable that we created from the names of the tibbles in our list and we should be all set. We will use another function from the stringr package to do this. The str_remove() function takes a string followed by a pattern and removes the pattern from the string.

nyts_data <- nyts_data %>%
  mutate(year = as.numeric(str_remove(year, "nyts")))

Here is our clean and wrangled data:

glimpse(nyts_data)

Rows: 95,465
Columns: 40
$ year                 <dbl> 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2…
$ psu                  <chr> "015438", "015438", "015438", "015438", "015438",…
$ finwgt               <dbl> 216.7268, 324.9620, 324.9620, 397.1552, 264.8745,…
$ stratum              <chr> "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", …
$ Age                  <fct> 18, 17, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 1…
$ Sex                  <fct> female, male, male, male, female, female, male, f…
$ Grade                <fct> 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 1…
$ ECIGT                <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ ECIGAR               <lgl> TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, FAL…
$ ESLT                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, T…
$ EELCIGT              <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ EROLLCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, F…
$ EFLAVCIGTS           <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, TRUE, FA…
$ EBIDIS               <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ EFLAVCIGAR           <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, FA…
$ EHOOKAH              <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ EPIPE                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ ESNUS                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, F…
$ EDISSOLV             <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CCIGT                <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, F…
$ CCIGAR               <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, F…
$ CSLT                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CELCIGT              <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, F…
$ CROLLCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CFLAVCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CBIDIS               <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CHOOKAH              <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CPIPE                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CSNUS                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CDISSOLV             <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ menthol              <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ clove_spice          <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ fruit                <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ chocolate            <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ alcoholic_drink      <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ candy_dessert_sweets <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ other                <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ EHTP                 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ CHTP                 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ brand_ecig           <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…

Note that there are several variables where there are similar names, but with a C compared to an E in the variable name. Those starting with C are related to questions about current usage (last 30 days), while those with an E are related to usage across the student respondent’s whole life (“ever” usage). We will discuss these groups further below.

Now we will save our wrangled data. We will save it as an rda file for ourselves and as csv files, as this is often a good option for collaborators. We will save this file in a directory called “wrangled” within our “data” directory of our project.

save(nyts_data, file = here::here("data", "wrangled", "wrangled_data.rda"))

write_csv(nyts_data, path = here::here("data", "wrangled", "nyts_data.csv"))

Variable Table

Click here to see a table about the final variables in our data set.

Variable	Details
year	the year that the survey results from a particular student respondent were acquired
psu	the primary sampling unit for the survey weighting
finwgt	the final analysis weight for the survey weighting
stratum	the stratum used for variance estimation for the survey weighting
Age	the age of the student when they took the survey
Sex	the sex of the student when they took the survey
Grade	the grade of the the student when the took the survey
ECIGT	student reported having ever tried cigarette smoking, even one or two puffs
ECIGAR	student reported having ever tried cigar, cigarillo, or little cigar smoking, even one or two puffs
ESLT	student reported having ever tried chewing tobacco, snuff, or dip
EELCIGT	student reported having ever tried e-cigarettes
EROLLCIGTS	student reported having ever tried roll-your-own cigarettes
EFLAVCIGTS	(2015 only) based on answer to “Which of the following tobacco products that you used in the past 30 days were flavored?”
EBIDIS	student reported having ever tried bidis (small brown cigarettes wrapped in a leaf)
EFLAVCIGAR	student reported having ever tried a flavored cigar (2015-2016)
EHOOKAH	student reported having ever smoked tobacco from a hookah or a waterpipe
EPIPE	student reported having ever smoked tobacco from a pipe (not hookah)
ESNUS	student reported having ever used snus, such as Camel or Malboro Snus
EDISSOLV	student reported having ever tried dissolvable tobacco products such as Ariva, Stonewall, Camel orbs, Camel sticks, Marlboro sticks, or Camel strips
CCIGT	student reported they smoked cigarettes on >= 1 of the past 30 days
CCIGAR	student reported they smoked cigars on >= 1 of the past 30 days
CSLT	student reported they used chewing tobacco, snuff, or dip on >= 1 of the past 30 days
CELCIGT	student reported they used electronic cigarettes or e-cigarettes one or more days in the past 30
CROLLCIGTS	student reported they smoked roll-your-own cigarettes during the past 30 days
CFLAVCIGTS	(2015 only) based on answer to “Which of the following tobacco products that you used in the past 30 days were flavored?”
CBIDIS	student reported they smoked bidis during the past 30 days
CHOOKAH	student reported they smoked tobacco in a hookah on >= 1 of the past 30 days
CPIPE	student reported they smoked tobacco in a pipe (not hookah) during the past 30 days
CSNUS	student reported they used snus during the past 30 days
CDISSOLV	student reported they used dissolvable tobacco products such as Ariva, Stonewall, Camel orbs, Camel sticks, Marlboro sticks, or Camel strips during the past 30 days
brand_ecig	student answer to “During the past 30 days, what brand of e-cigarettes did you usually use?”
menthol	student selected Menthol or mint as the answer to “What flavors of tobacco products have you used in the past 30 days? (select one or more)”
clove_spice	student selected clove or spice as the answer to “What flavors of tobacco products have you used in the past 30 days? (select one or more)”
fruit	student selected fruit as the answer to “What flavors of tobacco products have you used in the past 30 days? (select one or more)”
chocolate	student selected chocolate as the answer to “What flavors of tobacco products have you used in the past 30 days? (select one or more)”
alcoholic_drink	student selected alcoholic drink (such as wine, cognac, margarita, or other cocktails) as the answer to “What flavors of tobacco products have you used in the past 30 days? (choose one or more)”
candy_dessert_sweets	student selected candy, desserts or other sweets as the answer to “What flavors of tobacco products have you used in the past 30 days? (choose one or more)”
other	student selected some other flavor not listed as the answer to “What flavors of tobacco products have you used in the past 30 days? (choose one or more)”
EHTP	student reported having ever tried heated (also known as “heat-not-burn”) tobacco products
CHTP	student reported they used heated tobacco products during the past 30 days

Data Visualization

If you have been following along but stopped, we could load the wrangled data from the “data” directory like so:

load(here::here("data", "wrangled", "wrangled_data.rda"))

If you skipped the data exploration and wrangling section click here.

First you need to install and load the OCSdata package:

install.packages("OCSdata")
library(OCSdata)

Then, you may load the wrangled data using the following code:

wrangled_rda("ocs-bp-vaping-case-study", outpath = getwd())
load(here::here("OCSdata", "data", "wrangled", "wrangled_data.rda"))

If the package does not work for you, alternatively, an RDA file (stands for R data) of the data can be found in our GitHub repository or slightly more directly here. Download this file and then place it in your current working directory within a subdirectory called “wrangled” within a subdirectory called “data” to copy and paste our code. We used an RStudio project and the here package to navigate to the file more easily.

load(here::here("data", "wrangled", "wrangled_data.rda"))

Click here to see more about creating new projects in RStudio.

You can create a project by going to the File menu of RStudio like so:

You can also do so by clicking the project button:

See here to learn more about using RStudio projects and here to learn more about the here package.

Recall that our main questions were:

How has tobacco and e-cigarette/vaping use by American youths changed since 2015?
How does e-cigarette use compare between males and females?
What vaping brands and flavors appear to be used the most frequently?
We will base this on the following survey questions:
> “During the past 30 days, what brand of e-cigarettes did you usually use?”
> “What flavors of tobacco products have you used in the past 30 days?”
Is there a relationship between e-cigarette/vaping use and other tobacco use?

We are now going to create data visualizations to explore each of these questions.

For many of these questions we will be interested in both current and ever users, so we will want to create a variable for labeling individuals who are current users of any tobacco product (or not, i.e., who do not currently use a tobacco product) and a variable for labeling individuals who are “ever users” of any tobacco product (or not, i.e., who have never used a tobacco product).

We define these two groups as follows:

current = students who used a product for >=1 day in the past 30 days
ever = students who report having used or tried a product at any point in time

All current users are therefore ever users but not all ever users are current users. Thus, current users are a subset of ever users.

To add these grouping variables to our data we will do a bit more wrangling using the mutate() function again of the dplyr package. As discussed above, our dataset contains a set of questions that relate to whether the student has ever used the particular tobacco product (questions that start with the letter “E”), and questions that relate to whether the student currently uses the particular tobacco product (questions that start with the letter “C”).

Here are some examples for these data entries:

EPIPE: Students who reported they have smoked tobacco from a pipe (not hookah).
CPIPE: Students who reported they smoked tobacco in a pipe (not hookah) during the past 30 days.
EROLLCIGTS: Students who reported they have tried smoking roll-your-own cigarettes.
CROLLCIGTS: Students who reported they smoked roll-your-own cigarettes during the past 30 days.

Based on many questions like this:

In the past 30 days, which of the following products have you used on at least one day? (Select one or more)
A. Roll-your-own cigarettes
B. Pipes filled with tobacco (not hookah or waterpipe)
C. Snus, such as Camel, Marlboro, or General Snus
D. Dissolvable tobacco products such as Ariva, Stonewall, Camel orbs, Camel sticks, Marlboro sticks, or Camel strips
E. Bidis (small brown cigarettes wrapped in a leaf)
F. I have not used any of the products listed above in the past 30 days

Which of the following tobacco products have you ever tried, even just one time? (Select one or more)
A. Roll-your-own cigarettes
B. Pipes filled with tobacco (not hookah or waterpipe)
C. Snus, such as Camel, Marlboro, or General Snus
D. Dissolvable tobacco products such as Ariva, Stonewall, Camel orbs, Camel sticks, Marlboro sticks, or Camel strips
E. Bidis (small brown cigarettes wrapped in a leaf)
F. I have never tried any of the products listed above

We will sum across the variables that relate to ever or current tobacco usage questions to determine if the student answered yes to any of the ever or current questions. To do this we will use the base rowSums function.

We will then use the case_when() function of the dplyr package to convert the sum values to TRUE or FALSE based on the threshold of zero. If the sum is greater than zero, then we know the student answered yes to at least one question.

nyts_data %<>%
  mutate(tobacco_sum_ever = rowSums(select(., starts_with("E", 
                                    ignore.case = FALSE)), na.rm = TRUE),
      tobacco_sum_current = rowSums(select(., starts_with("C", 
                                    ignore.case = FALSE)), na.rm = TRUE)) %>%
      mutate(tobacco_ever = case_when(tobacco_sum_ever > 0 ~ TRUE,
                                      tobacco_sum_ever == 0 ~ FALSE),
          tobacco_current = case_when(tobacco_sum_current > 0 ~ TRUE,
                                      tobacco_sum_current == 0 ~ FALSE))

glimpse(nyts_data)

Rows: 95,465
Columns: 44
$ year                 <dbl> 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2…
$ psu                  <chr> "015438", "015438", "015438", "015438", "015438",…
$ finwgt               <dbl> 216.7268, 324.9620, 324.9620, 397.1552, 264.8745,…
$ stratum              <chr> "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", …
$ Age                  <fct> 18, 17, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 1…
$ Sex                  <fct> female, male, male, male, female, female, male, f…
$ Grade                <fct> 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 1…
$ ECIGT                <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ ECIGAR               <lgl> TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, FAL…
$ ESLT                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, T…
$ EELCIGT              <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ EROLLCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, F…
$ EFLAVCIGTS           <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, TRUE, FA…
$ EBIDIS               <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ EFLAVCIGAR           <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, FA…
$ EHOOKAH              <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ EPIPE                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ ESNUS                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, F…
$ EDISSOLV             <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CCIGT                <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, F…
$ CCIGAR               <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, F…
$ CSLT                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CELCIGT              <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, F…
$ CROLLCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CFLAVCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CBIDIS               <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CHOOKAH              <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CPIPE                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CSNUS                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CDISSOLV             <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ menthol              <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ clove_spice          <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ fruit                <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ chocolate            <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ alcoholic_drink      <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ candy_dessert_sweets <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ other                <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ EHTP                 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ CHTP                 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ brand_ecig           <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ tobacco_sum_ever     <dbl> 1, 4, 0, 3, 0, 2, 8, 4, 0, 0, 0, 1, 1, 0, 0, 4, 0…
$ tobacco_sum_current  <dbl> 0, 2, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
$ tobacco_ever         <lgl> TRUE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE,…
$ tobacco_current      <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, FALSE, FALSE, FA…

We are also interested in e-cigarette/vaping product usage compared to other tobacco products, so we will create some variables related to the sum of all e-cigarette usage question variables and the sum of all tobacco usage question variables excluding those that are about e-cigarettes. There is only one variable about e-cigarette usage ever (EELCIGT) and one about current usage (CELCIGT).

nyts_data <- nyts_data %>% 
  mutate(ecig_sum_ever = rowSums(select(., EELCIGT), na.rm = TRUE),
      ecig_sum_current = rowSums(select(., CELCIGT), na.rm = TRUE),
     non_ecig_sum_ever = rowSums(select(., starts_with("E", 
                                           ignore.case = FALSE), 
                                           -EELCIGT), na.rm = TRUE),
  non_ecig_sum_current = rowSums(select(., starts_with("C",
                                           ignore.case = FALSE), 
                                           -CELCIGT), na.rm = TRUE)) %>%
      mutate(ecig_ever = case_when(ecig_sum_ever > 0 ~ TRUE,
                                   ecig_sum_ever == 0 ~ FALSE),
          ecig_current = case_when(ecig_sum_current > 0 ~ TRUE,
                                   ecig_sum_current == 0 ~ FALSE),
         non_ecig_ever = case_when(non_ecig_sum_ever > 0 ~ TRUE,
                                   non_ecig_sum_ever == 0 ~ FALSE),
      non_ecig_current = case_when(non_ecig_sum_current > 0 ~ TRUE,
                                   non_ecig_sum_current == 0 ~ FALSE))

glimpse(nyts_data)

Rows: 95,465
Columns: 52
$ year                 <dbl> 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2…
$ psu                  <chr> "015438", "015438", "015438", "015438", "015438",…
$ finwgt               <dbl> 216.7268, 324.9620, 324.9620, 397.1552, 264.8745,…
$ stratum              <chr> "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", …
$ Age                  <fct> 18, 17, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, 1…
$ Sex                  <fct> female, male, male, male, female, female, male, f…
$ Grade                <fct> 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 1…
$ ECIGT                <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ ECIGAR               <lgl> TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, FAL…
$ ESLT                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, T…
$ EELCIGT              <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ EROLLCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, F…
$ EFLAVCIGTS           <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, TRUE, FA…
$ EBIDIS               <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ EFLAVCIGAR           <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, FA…
$ EHOOKAH              <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ EPIPE                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ ESNUS                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, F…
$ EDISSOLV             <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CCIGT                <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, F…
$ CCIGAR               <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, F…
$ CSLT                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CELCIGT              <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, F…
$ CROLLCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CFLAVCIGTS           <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CBIDIS               <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CHOOKAH              <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CPIPE                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CSNUS                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CDISSOLV             <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ menthol              <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ clove_spice          <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ fruit                <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ chocolate            <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ alcoholic_drink      <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ candy_dessert_sweets <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ other                <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ EHTP                 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ CHTP                 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ brand_ecig           <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, N…
$ tobacco_sum_ever     <dbl> 1, 4, 0, 3, 0, 2, 8, 4, 0, 0, 0, 1, 1, 0, 0, 4, 0…
$ tobacco_sum_current  <dbl> 0, 2, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
$ tobacco_ever         <lgl> TRUE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE,…
$ tobacco_current      <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, FALSE, FALSE, FA…
$ ecig_sum_ever        <dbl> 0, 1, 0, 1, 0, 1, 1, 1, 0, 0, 0, 1, 0, 0, 0, 1, 0…
$ ecig_sum_current     <dbl> 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
$ non_ecig_sum_ever    <dbl> 1, 3, 0, 2, 0, 1, 7, 3, 0, 0, 0, 0, 1, 0, 0, 3, 0…
$ non_ecig_sum_current <dbl> 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0…
$ ecig_ever            <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ ecig_current         <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, F…
$ non_ecig_ever        <lgl> TRUE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE,…
$ non_ecig_current     <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, F…

Finally, we are also interested in grouping students that only use e-cigarettes and those that only use other forms of tobacco.

Recall that current users are a subset of ever users, thus students would typically answer yes to having tried vaping products if they had used them one or more days in the past 30 days.

First we will make a small toy dataset called test to show what we will do with the larger dataset:

test <- tibble(ecig_ever = c("TRUE", "TRUE", "TRUE", "TRUE", "FALSE",
                             "FALSE", "TRUE", "FALSE", "FALSE"),
               non_ecig_ever = c("TRUE", "FALSE", "FALSE", "FALSE", "FALSE",
                                 "FALSE", "TRUE", "TRUE", "TRUE"),
               ecig_current = c("TRUE", "FALSE", "FALSE", "TRUE", "TRUE", 
                                "FALSE", "FALSE", "FALSE", "FALSE"),
               non_ecig_current = c("TRUE", "FALSE", "TRUE", "FALSE", "TRUE",
                                    "FALSE", "FALSE", "FALSE", "TRUE"))

test

# A tibble: 9 × 4
  ecig_ever non_ecig_ever ecig_current non_ecig_current
  <chr>     <chr>         <chr>        <chr>           
1 TRUE      TRUE          TRUE         TRUE            
2 TRUE      FALSE         FALSE        FALSE           
3 TRUE      FALSE         FALSE        TRUE            
4 TRUE      FALSE         TRUE         FALSE           
5 FALSE     FALSE         TRUE         TRUE            
6 FALSE     FALSE         FALSE        FALSE           
7 TRUE      TRUE          FALSE        FALSE           
8 FALSE     TRUE          FALSE        FALSE           
9 FALSE     TRUE          FALSE        TRUE

Now, let’s look at identifying students who have tried e-cigarettes, but are not current users, and who have never tried other tobacco products (and are therefore not current users). We will again use the case_when() and the mutate function to create new variables with specific values when certain conditions are met. In this case, we will specify that several conditions must be met by using the & symbol. For a value of TRUE for the new ecig_only_ever variable, all of the conditions combined with & must be met. If any of the conditions are not met then the ecig_only_ever value will be FALSE based on the last line TRUE ~ FALSE.

test <- test %>% mutate(ecig_only_ever = case_when(ecig_ever == TRUE &
                                               non_ecig_ever == FALSE &
                                                ecig_current == FALSE &
                                            non_ecig_current == FALSE ~ TRUE,
                                                         TRUE ~ FALSE))
test

# A tibble: 9 × 5
  ecig_ever non_ecig_ever ecig_current non_ecig_current ecig_only_ever
  <chr>     <chr>         <chr>        <chr>            <lgl>         
1 TRUE      TRUE          TRUE         TRUE             FALSE         
2 TRUE      FALSE         FALSE        FALSE            TRUE          
3 TRUE      FALSE         FALSE        TRUE             FALSE         
4 TRUE      FALSE         TRUE         FALSE            FALSE         
5 FALSE     FALSE         TRUE         TRUE             FALSE         
6 FALSE     FALSE         FALSE        FALSE            FALSE         
7 TRUE      TRUE          FALSE        FALSE            FALSE         
8 FALSE     TRUE          FALSE        FALSE            FALSE         
9 FALSE     TRUE          FALSE        TRUE             FALSE

We can see from the second row, that the ecig_only_ever is TRUE when we would expect it to be. We can also see from the fourth row, that even though the student reported yes to ever trying e-cigarettes, because they also reported yes to currently using e-cigarettes the value for only ever trying e-cigarettes is FALSE. Additionally we can see from the seventh row that similarly even though the student reported yes to ever trying e-cigarettes, they also reported yes to ever trying other products, and the value for only ever trying e-cigarettes is FALSE. Importantly, we can see from the 6th row, that if all responses are negative than the value is FALSE.

Now we will expand this to the other possible categories. In this case we note that since current users are a subset of ever users, it doesn’t matter if a user reports yes to ever trying e-cigarettes, they can still be a current user.

test <- test %>%
         mutate(ecig_only_ever = case_when(ecig_ever == TRUE &
                                       non_ecig_ever == FALSE &
                                        ecig_current == FALSE &
                                    non_ecig_current == FALSE ~ TRUE,
                                                        TRUE ~ FALSE),
          ecig_only_current = case_when(ecig_current == TRUE &
                                       non_ecig_ever == FALSE &
                                    non_ecig_current == FALSE ~ TRUE,
                                                        TRUE ~ FALSE),
        non_ecig_only_ever = case_when(non_ecig_ever == TRUE &
                                           ecig_ever == FALSE &
                                        ecig_current == FALSE &
                                    non_ecig_current == FALSE ~ TRUE,
                                                        TRUE ~ FALSE),
  non_ecig_only_current = case_when(non_ecig_current == TRUE &
                                           ecig_ever == FALSE &
                                        ecig_current == FALSE ~ TRUE,
                                                        TRUE ~ FALSE),
                    no_use = case_when(non_ecig_ever == FALSE &
                                           ecig_ever == FALSE &
                                        ecig_current == FALSE &
                                    non_ecig_current == FALSE ~ TRUE,
                                                        TRUE ~ FALSE))
glimpse(test)

Rows: 9
Columns: 9
$ ecig_ever             <chr> "TRUE", "TRUE", "TRUE", "TRUE", "FALSE", "FALSE"…
$ non_ecig_ever         <chr> "TRUE", "FALSE", "FALSE", "FALSE", "FALSE", "FAL…
$ ecig_current          <chr> "TRUE", "FALSE", "FALSE", "TRUE", "TRUE", "FALSE…
$ non_ecig_current      <chr> "TRUE", "FALSE", "TRUE", "FALSE", "TRUE", "FALSE…
$ ecig_only_ever        <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ ecig_only_current     <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, …
$ non_ecig_only_ever    <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ non_ecig_only_current <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ no_use                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, FALSE, …

Take a minute to check that the values are what we would expect.

OK, now we are going to make a Group variable based on the new variables we just made to classify students into one of four mutually exclusive and exhaustive categories. In this case we will have a particular value based on one condition or another. This or conditional is specified by the | symbol. Only one of the conditions needs to exist for that particular value, whereas when we used the & symbol, all of the conditions had to be met.

If a student has ever tried or currently uses e-cigarettes, but has never tried other tobacco products, the value will be Only e-cigarettes. If a student has ever tried or is a current user of other tobacco products, but has never tried e-cigarettes, the value will be Only other products. If the value of the no_use variable is simply TRUE, then the Group variable value will be Neither. Finally, if a student has tried or currently uses both e-cigarettes and other tobacco products the Group variable value will be Combination of products. Thus in this case the values for the usage of the variables based on only using e-cigarettes or only other products will all be FALSE.

test <- test %>%
  mutate(Group = case_when(ecig_only_ever == TRUE |
                        ecig_only_current == TRUE ~ "Only e-cigarettes",
                       non_ecig_only_ever == TRUE |
                    non_ecig_only_current == TRUE ~ "Only other products",
                                   no_use == TRUE ~ "Neither",
                           ecig_only_ever == FALSE &
                        ecig_only_current == FALSE &
                       non_ecig_only_ever == FALSE &
                    non_ecig_only_current == FALSE &
                                   no_use == FALSE ~ "Combination of products"))


test %>% count(Group)

# A tibble: 4 × 2
  Group                       n
  <chr>                   <int>
1 Combination of products     4
2 Neither                     1
3 Only e-cigarettes           2
4 Only other products         2

glimpse(test)

Rows: 9
Columns: 10
$ ecig_ever             <chr> "TRUE", "TRUE", "TRUE", "TRUE", "FALSE", "FALSE"…
$ non_ecig_ever         <chr> "TRUE", "FALSE", "FALSE", "FALSE", "FALSE", "FAL…
$ ecig_current          <chr> "TRUE", "FALSE", "FALSE", "TRUE", "TRUE", "FALSE…
$ non_ecig_current      <chr> "TRUE", "FALSE", "TRUE", "FALSE", "TRUE", "FALSE…
$ ecig_only_ever        <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ ecig_only_current     <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, …
$ non_ecig_only_ever    <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ non_ecig_only_current <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ no_use                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, FALSE, …
$ Group                 <chr> "Combination of products", "Only e-cigarettes", …

OK, now that we have seen how this works with our toy dataset, we will apply our code to our nyts_data.

nyts_data %<>%
             mutate(ecig_only_ever = case_when(ecig_ever == TRUE &
                                           non_ecig_ever == FALSE &
                                            ecig_current == FALSE &
                                        non_ecig_current == FALSE ~ TRUE,
                                                             TRUE ~ FALSE),
              ecig_only_current = case_when(ecig_current == TRUE &
                                           non_ecig_ever == FALSE &
                                        non_ecig_current == FALSE ~ TRUE,
                                                            TRUE ~ FALSE),
            non_ecig_only_ever = case_when(non_ecig_ever == TRUE &
                                               ecig_ever == FALSE &
                                            ecig_current == FALSE &
                                        non_ecig_current == FALSE ~ TRUE,
                                                            TRUE ~ FALSE),
      non_ecig_only_current = case_when(non_ecig_current == TRUE &
                                               ecig_ever == FALSE &
                                            ecig_current == FALSE ~ TRUE,
                                                            TRUE ~ FALSE),
                        no_use = case_when(non_ecig_ever == FALSE &
                                               ecig_ever == FALSE &
                                            ecig_current == FALSE &
                                        non_ecig_current == FALSE ~ TRUE,
                                                            TRUE ~ FALSE)) %>%
                 mutate(Group = case_when(ecig_only_ever == TRUE |
                                       ecig_only_current == TRUE ~ "Only e-cigarettes",
                                      non_ecig_only_ever == TRUE |
                                   non_ecig_only_current == TRUE ~ "Only other products",
                                                  no_use == TRUE ~ "Neither",
                                          ecig_only_ever == FALSE &
                                       ecig_only_current == FALSE &
                                      non_ecig_only_ever == FALSE &
                                   non_ecig_only_current == FALSE &
                                                  no_use == FALSE ~ "Combination of products"))

Lastly, it can be very helpful to have the total number of students surveyed each year. We can easily add a variable for this by using the add_count() function of the dplyr package. This will create a variable called n which will show the total number of survey responses for that year.

nyts_data %<>% dplyr::add_count(year)

glimpse(nyts_data)

Rows: 95,465
Columns: 59
$ year                  <dbl> 2015, 2015, 2015, 2015, 2015, 2015, 2015, 2015, …
$ psu                   <chr> "015438", "015438", "015438", "015438", "015438"…
$ finwgt                <dbl> 216.7268, 324.9620, 324.9620, 397.1552, 264.8745…
$ stratum               <chr> "BR3", "BR3", "BR3", "BR3", "BR3", "BR3", "BR3",…
$ Age                   <fct> 18, 17, 18, 18, 18, 18, 18, 18, 18, 18, 18, 18, …
$ Sex                   <fct> female, male, male, male, female, female, male, …
$ Grade                 <fct> 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, 12, …
$ ECIGT                 <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRU…
$ ECIGAR                <lgl> TRUE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, FA…
$ ESLT                  <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, …
$ EELCIGT               <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRU…
$ EROLLCIGTS            <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, …
$ EFLAVCIGTS            <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, TRUE, F…
$ EBIDIS                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ EFLAVCIGAR            <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, TRUE, F…
$ EHOOKAH               <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ EPIPE                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ ESNUS                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, TRUE, …
$ EDISSOLV              <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ CCIGT                 <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CCIGAR                <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ CSLT                  <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ CELCIGT               <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, …
$ CROLLCIGTS            <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ CFLAVCIGTS            <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ CBIDIS                <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ CHOOKAH               <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ CPIPE                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ CSNUS                 <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ CDISSOLV              <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ menthol               <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ clove_spice           <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ fruit                 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ chocolate             <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ alcoholic_drink       <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ candy_dessert_sweets  <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ other                 <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ EHTP                  <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ CHTP                  <lgl> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ brand_ecig            <chr> NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, NA, …
$ tobacco_sum_ever      <dbl> 1, 4, 0, 3, 0, 2, 8, 4, 0, 0, 0, 1, 1, 0, 0, 4, …
$ tobacco_sum_current   <dbl> 0, 2, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
$ tobacco_ever          <lgl> TRUE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ tobacco_current       <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, FALSE, FALSE, F…
$ ecig_sum_ever         <dbl> 0, 1, 0, 1, 0, 1, 1, 1, 0, 0, 0, 1, 0, 0, 0, 1, …
$ ecig_sum_current      <dbl> 0, 0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
$ non_ecig_sum_ever     <dbl> 1, 3, 0, 2, 0, 1, 7, 3, 0, 0, 0, 0, 1, 0, 0, 3, …
$ non_ecig_sum_current  <dbl> 0, 2, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, …
$ ecig_ever             <lgl> FALSE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRU…
$ ecig_current          <lgl> FALSE, FALSE, FALSE, TRUE, FALSE, FALSE, FALSE, …
$ non_ecig_ever         <lgl> TRUE, TRUE, FALSE, TRUE, FALSE, TRUE, TRUE, TRUE…
$ non_ecig_current      <lgl> FALSE, TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ ecig_only_ever        <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ ecig_only_current     <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ non_ecig_only_ever    <lgl> TRUE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, …
$ non_ecig_only_current <lgl> FALSE, FALSE, FALSE, FALSE, FALSE, FALSE, FALSE,…
$ no_use                <lgl> FALSE, FALSE, TRUE, FALSE, TRUE, FALSE, FALSE, F…
$ Group                 <chr> "Only other products", "Combination of products"…
$ n                     <int> 17711, 17711, 17711, 17711, 17711, 17711, 17711,…

Now let’s save our wrangled data for future use. We will save it as an r compatible file (a .rda file), as well as csv files, as these are useful to give to collaborators. We will use the write_csv() function of the readr package to do this.

save(nyts_data, file = here::here("data", "wrangled", "wrangled_data_with_var_for_plots.rda"))
readr::write_csv(nyts_data, 
                 path = here::here("data", "wrangled", "nyts_data_for_plots.csv"))

Question 1

Click here if you skipped the previous sections and want to start here.

First you need to install and load the OCSdata package:

install.packages("OCSdata")
library(OCSdata)

Then, you may load the wrangled data for plots using the following code:

wrangled_rda("ocs-bp-vaping-case-study", outpath = getwd())
load(here::here("OCSdata", "data", "wrangled", "wrangled_data_with_var_for_plots.rda"))

If the package does not work for you, alternatively, an RDA file (stands for R data) of the data can be found in our GitHub repository or slightly more directly here. Download this file and then place it in your current working directory within a subdirectory called “wrangled” within a subdirectory called “data” to copy and paste our code. We used an RStudio project and the here package to navigate to the file more easily.

load(here::here("data", "wrangled", "wrangled_data_with_var_for_plots.rda"))

Recall that we are interested in investigating how vaping product use has compared with other tobacco use over time. To answer this, we first want to get a sense of how tobacco use has changed in general since 2015.

To create a visualization of how tobacco usage has changed over time, we will first convert the usage data to a percent value for each year, telling us what percent of student respondents fall into a particular usage category. To do this we will use the group_by() and summarize() functions of the dplyr package. This will create new variables which we will name Ever and Current based on the percentages of TRUE values for tobacco_ever and tobacco_current for each year. In this case the mean() function is used to calculate the percentages based on an automatic conversion that R does where for TRUE/FALSE variables, TRUE is given a value of one and FALSE is given a value of zero. The mean of a 0-1 binary variable is just the percent of the time the value is 1. NA values do not contribute to the total count when we include the argument na.rm = TRUE to our function call.

Click here to see a toy example:

# the test data has 3 TRUE values and 7 FALSE values
test <- tibble("var1" = c("TRUE", "TRUE", "TRUE", rep("FALSE", 7)))
test %<>% mutate(var1 = as.logical(var1))
test

# A tibble: 10 × 1
   var1 
   <lgl>
 1 TRUE 
 2 TRUE 
 3 TRUE 
 4 FALSE
 5 FALSE
 6 FALSE
 7 FALSE
 8 FALSE
 9 FALSE
10 FALSE

test %>% summarize(Percentage = mean(var1) * 100)

# A tibble: 1 × 1
  Percentage
       <dbl>
1         30

# the test data has 3 TRUE values, 3 FALSE values, and 4 NA value
test <- tibble("var1" = c("TRUE", "TRUE", "TRUE", rep("FALSE", 3), rep("NA", 4)))
test %<>% mutate(var1 = as.logical(var1))
test

# A tibble: 10 × 1
   var1 
   <lgl>
 1 TRUE 
 2 TRUE 
 3 TRUE 
 4 FALSE
 5 FALSE
 6 FALSE
 7 NA   
 8 NA   
 9 NA   
10 NA

test %>% summarize(Percentage = mean(var1, na.rm = TRUE) * 100)

# A tibble: 1 × 1
  Percentage
       <dbl>
1         50

And now back to our data:

nyts_data %>%
  dplyr::group_by(year) %>%
  dplyr::summarize(Ever = (mean(tobacco_ever, na.rm = TRUE) * 100),
                   Current = (mean(tobacco_current, na.rm = TRUE) * 100))

# A tibble: 5 × 3
   year  Ever Current
  <dbl> <dbl>   <dbl>
1  2015  36.8    18.0
2  2016  33.4    14.0
3  2017  31.8    14.4
4  2018  34.7    18.7
5  2019  39.7    22.4

We will use the pivot_longer function from the tidyr package to take all columns except year (in this case the Ever and Current columns), to create a column called User that will contain the current column name information and a column called Percentage of students which will contain the mean percentage values that we just calculated. This converts our data into a format called “long” format.

nyts_data %>%
  dplyr::group_by(year) %>%
  dplyr::summarize("Ever \n (any lifetime use)" = (mean(tobacco_ever, na.rm = TRUE) * 100),
                   "Current \n (any past-30-day use)" = (mean(tobacco_current, na.rm = TRUE) * 100)) %>%
  tidyr::pivot_longer(cols = -year,
                      names_to = "User",
                      values_to = "Percentage of students")

# A tibble: 10 × 3
    year User                               `Percentage of students`
   <dbl> <chr>                                                 <dbl>
 1  2015 "Ever \n (any lifetime use)"                           36.8
 2  2015 "Current \n (any past-30-day use)"                     18.0
 3  2016 "Ever \n (any lifetime use)"                           33.4
 4  2016 "Current \n (any past-30-day use)"                     14.0
 5  2017 "Ever \n (any lifetime use)"                           31.8
 6  2017 "Current \n (any past-30-day use)"                     14.4
 7  2018 "Ever \n (any lifetime use)"                           34.7
 8  2018 "Current \n (any past-30-day use)"                     18.7
 9  2019 "Ever \n (any lifetime use)"                           39.7
10  2019 "Current \n (any past-30-day use)"                     22.4

You may have noticed that our data is longer than it used to be! Hence the name of the function pivot_longer(). Data is often easier to plot when it is in this format.

Now we will use this data to create a plot using the ggplot2 package.

The first thing we do to create a plot is specify what data we are using for our x axis and y axis with theaes() argument of the ggplot() function. Then we add layers to our plot that specify what type of plot we would like to create. We can use the geom_line() function to create lines for each type of user. We specify that we want to use different line types for each user using aes(). We will also add points to our lines using the geom_point() function. We can add additional layers to specify colors and details about labels and legends etc.

plot1 <- nyts_data %>%
    dplyr::group_by(year) %>%
    dplyr::summarize("Ever \n (any lifetime use)" = (mean(tobacco_ever, na.rm = TRUE) * 100),
                     "Current \n (any past-30-day use)" = (mean(tobacco_current, na.rm = TRUE) * 100)) %>%
  pivot_longer(cols = -year, names_to = "User", values_to = "Percentage of students") %>%
  ggplot(aes(x = year, y = `Percentage of students`)) +
  geom_line(aes(linetype = User)) +
  geom_point(show.legend = FALSE, size = 2) +
  # this allows us to choose what type of line we want for each line
  scale_linetype_manual(values = c(1, 2), 
                        breaks = c("Ever \n (any lifetime use)", 
                                   "Current \n (any past-30-day use)")) +
  # this allows us to specify how the y-axis should appear
  scale_y_continuous(breaks = seq(0, 70, by = 10),
                     labels = seq(0, 70, by = 10),
                     limits = c(0, 70)) +
  # this adjusts the background style of the plot
    theme_linedraw() +
  # this moves the legend to the bottom of the plot and removes the x axis title
    theme(legend.position = "bottom",
          axis.title.x = element_blank()) +
    labs(title = "How has tobacco use varied over the years?",
         y = "% of students")

plot1 + theme(text = element_text(size = 15))