The National Survey on Drug Use and Health (NSDUH) is a national study run by the DHHS that collects information on several health-related issues including tobacco, alcohol, drug use, and mental health in the United States. The population of interest is the civilian, non-institutionalized population aged 12 or older at the time of the survey, i.e., excluding active military personnel, people living in institutional group quarters, and homeless individuals. Our analysis only uses data on adults ages 18 and older from the 2020 survey (27,170 observations).

Participants are selected for the survey using an independent, multistage area probability sample within each state and the District of Columbia. Because geographic identifiers such as state are not included in the 2020 public use file, it is only possible to make national estimates. The weights in this file are approximately the inverse probability of selection for each record.

The response variables of interest from this survey include level of psychological distress over the past 30 days, worst psychological distress over the past 30 days, predicted probability of serious mental health illness, and mental illness category as well as indicators for the past month of serious psychological distress, serious suicidal thoughts, suicidal plans, suicide attempt, lifetime major depressive episode, and major depressive episode in last year. Available demographic and social characteristics include gender, educational status, marital status, age, race, employment status, income, and poverty level.

Due to the COVID-19 pandemic, Quarter 4 of 2020 was the first time that the NSDUH utilized web-based interviewing. However, a concerningly high number of adults provided usable information on their substance use but did not complete the mental health or later questions (i.e., “break-offs”). The DHHS created additional analysis weights to analyze the unimputed outcomes starting from the mental health and subsequent sections of the questionnaire to account for respondents who broke off the interview before completing these sections. Unfortunately, the public use file does not contain these “break-off analysis weights,” so estimates from 2020 should not be compared longitudinally to other years of survey data.

The Behavioral Risk Factor Surveillance System (BRFSS) is a nationwide system run by the CDC of health-related telephone surveys that provide state-level information about health-related risk behaviors, chronic health conditions, and the use of preventive services among U.S. residents. The population of interest is the noninstitutionalized adult population (18 years or older) residing in private residences or college housing in the United States or participating areas who have a working cellular telephone. Participants are selected through an overlapping, dual-frame landline and cell phone sample. The BRFSS samples high and medium-density strata to obtain a probability sample of all households with telephones. State health departments may directly collect data from their state residents, or they may use a contractor. Person-level analysis weights involve two components: design-based weights and weight adjustment factors. Design weights reflect probabilities of selection at the sample stage. Weight adjustments are performed using the generalized exponential model (GEM) developed by Folsom and Singh [7], which calibrated the design-based weights to reduce non-response bias, poststratify to known population control totals, and control for extreme weights when necessary.

The response variables of interest from this survey include an indicator for whether individuals were ever told they had a depressive disorder and a self-assessment of general health. Relevant demographic and social variables include health insurance status, age, education, and income.

The main package we used to specify complex survey designs and produce unbiased summary statistics was survey: Analysis of Complex Survey Samples [9], written by Thomas Lumley. This package creates svydesign objects that ensure that the design information cannot be separated or used with the wrong data.

The fundamental concept that underlies design-based inference is that an individual sampled with sampling probability π i represents 1 π i individuals in the population. 1 π i is referred to as the sampling weight. All of the analyses run with this package use the Horvitz-Thompson estimator to estimate the population total [8]. For a sample of size n, the Horvitz-Thompson estimator T X ˆ for the population total of X is (2.1) T X ˆ = ∑ i = 1 n 1 π i X i = ∑ i = 1 n X i ˆ . The variance estimate is (2.2) v a r [ T X ˆ ] ˆ = ∑ i , j X i X j π i j − X i π i X j π j , which can also be written as (2.3) v a r [ T X ˆ ] ˆ = N − n N x N 2 x v a r [ X ] n , for a population of size N. The svymean() function used to produce mean estimates in the app is estimated by dividing the estimated total by the population size N. The variance estimate is estimated by dividing the variance estimate for the total by N 2 . The Horvitz-Thompson estimator of the population size is (2.4) ∑ i = 1 n 1 π i .

For estimates in a subpopulation, the survey package handles the computational details of domain estimation and sets the sampling weights to 0 for observations outside of the subpopulation.

The visualizations on the National tab are implemented with the package, ggsurvey [1] that simplifies ggplot2 functions for svydesign objects. The package’s functions call “ggplot2” to make bar charts, histograms, boxplots, and hexplots of survey objects to accurately represent the weighted sample distributions.

The interactive app is implemented using the shiny package [4] and it is deployed on shinyapps.io in its own protected environment with SSL-encrypted access.

We begin by visualizing national-level statistics from the NSDUH data set through boxplots and histograms for continuous response variables and bar charts for categorical responses. The visualizations pipe the current input values into the ggsurvey functions to dynamically render the histogram and box plots in response to the user’s selected input values. With the log transform Y radio button, the plots remain in the original scale for interpretability, but the breaks in the y-axis for the boxplot and the x-axis for the histogram are spaced according to the log scale. This allows the user to analyze and compare the centers and spreads of the skewed distributions with greater precision.

Figure 1 shows the distribution of psychological distress within the past month without the log transformation.

The table of national-level summaries appears on the National tab and includes data on the estimated population mean standard error, and confidence intervals for the population distributions represented in the visualizations. Logistically, we first find the estimated survey mean, standard error, and confidence interval for the response using the Horvitz-Thompson estimators (see Equations (2.1) and (2.2)) as calculated in the survey package. The total number of individuals represented in the estimated population for each subgroup is calculated using Equation (2.4), i.e., adding the sampling weights of the individuals in the survey that fit into each category. Formatted data table objects may be sorted by any of the columns, allowing users to order estimates.

Figure 2 shows the table view of the national distribution of the categorical mental illness variable by gender.

The state-level choropleths visually represent the geographic distribution of response variables of interest in the BRFSS across levels of social factors. In order to increase the speed at which the outputs are rendered, we preprocessed the BRFSS survey object to produce lists of tibbles that represent the distribution of the General Health and Depressive Diagnosis variables across health insurance status, education, age, and race (only for General Health). These tibbles are then merged with a simplified shapefile from the tigris [14] package that has data on the primary governmental divisions of the 50 states in the United States, as well as the District of Columbia, Puerto Rico, American Samoa, the Commonwealth of the Northern Mariana Islands, Guam, and the U.S. Virgin Islands.

Figure 3 shows the geographic distribution of depressive disorder prevalence by health insurance status.

After the variables are specified, a function is applied over the corresponding tibble to render a choropleth with leaflet [5], a JavaScript library used to build web mapping applications. Each state on the leaflet has hover text including the state name, estimated mean of the selected response variable, and estimated standard error.

The state table organizes the data visualized in the choropleths into rows that represent the estimated mean of the selected response variable for each state across the different levels of the selected demographic variable (in columns). We combined the reactive list of tibbles produced by the input calls by variable name using a binary function and then spread the rows by the demographic variable to show the state-level mean estimates for the response variables for all states across all levels of the demographic variable.

Figure 4 shows the table view of the state-level estimates of the prevalence of depressive diagnoses as well as the marginal distribution by health insurance status.