INTRODUCTION

Benzodiazepines (BZD) impair cognition1 and are associated with motor vehicle accidents,2 misuse and abuse,3 as well as falls4 and fractures5 among older adults. After opioids, BZDs are the second-most common medication class linked with pharmaceutical overdose deaths,6 the rate of which grew 13.6% per year from 1996 to 2013.7 The extent of BZD prescribing is high given the availability of safer, guideline-recommended treatments like selective serotonin reuptake inhibitors8,9 and cognitive behavioral therapy,10,11 though BZDs may have a role in treatment-resistant cases of anxiety or insomnia.8,9,12,13 Surprisingly, a significant amount of prescribing—perhaps even the majority—occurs to patients without insomnia or a clearly defined mental health disorder.14,15,16

General medical settings are the country’s primary mental health treatment setting.17 This is particularly true for older adults (i.e., ≥ 65 years), as psychiatrists are less likely than other medical specialties to accept Medicare.18 This may help account for the fact that among older adults—who have the highest prevalence of BZD use but the lowest use of specialty mental health care19—just 5.6% of those prescribed a BZD receive it from a psychiatrist.20

Within primary care, the links between geographic or prescriber characteristics and BZD prescribing have been subject to little research. Variation in BZD prescribing at the physician level may indicate a lack of knowledge about or consensus regarding their appropriate use. Case and Deaton recently suggested that growth in mortality among white Americans—partially driven by drug and alcohol poisonings and suicide—reflects an epidemic of distress.21 The growing disparities in life expectancy among US counties, largely driven by county-level race and socioeconomic status,22 suggest that such distress is concentrated within particular communities. If mental health disorders are not fully driving BZD prescribing, more extensive BZD prescribing in socioeconomically distressed counties might suggest that this distress—rather than a clinical indication such as panic disorder—is driving potentially inappropriate prescribing.

Understanding how BZD prescribing varies across primary care settings may help design and target interventions with the greatest potential to reduce prescribing, especially given the possible reluctance of patients23 and physicians24 to broadly reduce use. To our knowledge, no prior national analysis has examined BZD prescribing at the county and physician levels. The foundation of this analysis is the 2015 Medicare Part D Public Use Files, which contain virtually 100% of all Part D claims by Medicare beneficiaries in 2015. We pair this prescribing data with physician characteristics and county information from the Area Health Resources Files and County Health Rankings to determine the county and physician characteristics associated with high-intensity BZD prescribing.

METHODS

Primary Data Source

We used the 2015 Summary and Detailed Medicare Part D Public Use Files (PUF),25 which are derived from all Medicare Part D final-action claims from stand-alone and Medicare Advantage prescription drug plans.

The Summary PUF is organized by National Provider Identifier (NPI) and includes physician zip code, gender, and credentials; select patient panel characteristics, such as age, gender, race, low-income subsidy, and hierarchical condition category (HCC; a risk-adjustment score calculated by the Centers for Medicare & Medicaid Services [CMS]); and summary prescribing information, including total days prescribed overall and of specific medication classes, including opioids, antibiotics, and high-risk medications in the elderly.26 The Summary PUF includes information on over 57.7 billion days of prescription claims from 1.1 million prescribers; we limited this to 122,054 individual primary care physicians (PCPs), who as a group do the majority of BZD prescribing to older adults.20

The Detailed PUF is organized by NPI and drug name: each row includes a distinct combination of NPI + drug name and the total days’ supply for that medication, which we used to determine the total days of BZD prescribed by each PCP. To protect beneficiary privacy, the Detailed PUF excludes any record derived from ≤ 10 claims (e.g., if a PCP prescribed eight claims for lorazepam, this combination of NPI + lorazepam is not reported). After redaction, the 2015 Detailed PUF includes 1.2 billion prescription claims or 86.9% of total Part D claims.25

County- and Physician-Level Outcomes

Our county-level outcome is the intensity of BZD prescribing by a county’s PCPs defined as (days of BZDs prescribed by PCPs) × 100/(days of all Part D medications prescribed by PCPs). We used the Summary PUF to determine the total medication days prescribed by each physician; their BZD days were determined from the Detailed PUF by adding up all days prescribed by each NPI for each BZD.25 We derived the BZD intensity numerator and denominator by summing data for all PCPs within each zip code level and then using the U.S. Department of Housing and Urban Development zip-to-county crosswalk to generate county-level estimates.27 The final sample for the county-level analysis included 3038 counties with claims from the 122,054 individual PCPs.

As with the county-level outcome, we defined PCP-level BZD prescribing intensity as (days of BZD prescribed) × 100/(total days of all Part D medications prescribed).

Covariates for County-Level Analysis

We determined county-level demographic, health care system, socioeconomic status (SES), and health status characteristics primarily from the County Health Rankings28 (CHR; a joint initiative of the University of Wisconsin Population Health Institute and the Robert Wood Johnson Foundation) and the Area Health Resources Files (AHRF).29 CHR and AHRF incorporate data elements from a variety of sources, including the American Community Survey, the Behavioral Risk Factor Surveillance System, and CMS (eTable 1).

Demographic and SES variables used included percentage of female residents, percentage of residents ≥ 65 years, percentage of residents by race/ethnicity, percentage rural, percentage that completed any college, unemployment rate, and median income.

Health care system characteristics included primary care providers per 100,000 population, presence of a psychiatrist, federally qualified health center or community mental health center, and whether part of a county was designated a primary care or mental health provider shortage area. Finally, we used the Detailed PUF to determine the proportion of BZDs in each county prescribed by non-PCP physicians.

We used the following indicators of health status: proportion of disability-eligible Medicare beneficiaries, mean HCC, percentage reporting poor health in the past 30 days, days of poor health in the past 30 for physical or mental health reasons, binge alcohol use, and suicide rate.

Covariates for Physician-Level Analysis

We used the Summary PUF to categorize PCPs by specialty (internal medicine, family medicine, or geriatrics), MD or DO, gender, and years in practice (from the Medicare Physician Compare National Downloadable File).30

We used the opioid, antibiotic, and “high-risk medication” (i.e., inappropriate for adults 65+) claim totals provided for every provider in the Summary PUF to determine high-intensity prescribing for each of these three medication groups (i.e., high-intensity opioid PCP = top quartile of [opioid claims/total claims prescribed by that PCP]). We also identified overall high-volume and high-intensity prescribers (see eTable 2 for further details).

Finally, we included the following PCP panel characteristics from the Summary PUF: patient age, HCC score, percentage of white patients, and percentage of patients receiving the Part D low-income subsidy.

Analysis

We determined total BZD claims in 2015 and then limited the data to those prescribed by PCPs, who prescribed 61.9% of all Part D BZD days. We ranked BZDs by individual medication, grouped as short- and long-acting.20 We sorted counties into quartiles by intensity of BZD prescribing and compared quartiles on the county characteristics of interest using chi-square or ANOVA tests. To determine the characteristics that predicted high-intensity BZD prescribing, we grouped the bottom three quartiles together and tested the bivariate association between county characteristic and BZD prescribing intensity using logistic regressions (dependent variable: 1 = top quartile; 0 = bottom three quartiles). Next, we used multivariable logistic regression to predict top quartile BZD prescribing as a combined function of county demographic, SES, and health care system characteristics. The second and final county-level model incorporated indicators of health status to test whether other markers of poor health and distress were associated with BZD prescribing intensity.

For the PCP-level analysis, first we excluded the bottom decile of overall low-volume PCPs (n = 12,368; prescribed to a mean of 77 beneficiaries) because low-practice volume contributes to less reliable performance measurement,31 yielding n = 109,686 for the physician-level analysis. We split PCPs into quartiles by prescribing intensity and compared characteristics of interest across quartiles. For the model predicting high-intensity BZD prescribing, we again grouped the bottom three quartiles together. We used a multilevel logistic regression to predict the odds of a PCP being a high-intensity BZD prescriber. The model included PCP characteristics, patient panel characteristics, and county-level mental health treatment resources to account for local differences in availability of mental health specialty care as fixed effects, with a random effect for county. Analysis was completed in Stata 13.1 (College Station, TX). All statistical tests were two-tailed with α = 0.01 given the large sample size and multiple comparisons. Because this analysis only uses de-identified, aggregated patient data, it is not subject to IRB approval by the University of Michigan Medical School IRB.

RESULTS

Medicare Part D beneficiaries were dispensed over US$374 million worth and 1.2 billion total days of BZDs, which comprised 2.3% of all Part D medication days. PCPs prescribed over 728 million of all BZD days, accounting for 61.9% of BZD prescribing overall. Short-acting BZDs were prescribed most commonly (Table 1). On average, BZDs comprised 2.1% (standard deviation [SD] 0.9) of medications prescribed in a county (i.e., the intensity), with an interquartile range [IQR] of 1.5–2.6%. States with the highest county-level prescribing were all in the South. The top five were Louisiana, Florida, West Virginia, Tennessee, and Alabama (Fig. 1).

Table 1 Total Benzodiazepine Prescriptions to Part D Medicare Beneficiaries by Primary Care Physicians in 2015 (n = 122,054)
Figure 1
figure 1

County-level benzodiazepine prescribing intensity by primary care physicians across 50 states and the District of Columbia (n = 3038). The red horizontal line is the national county-level mean (2.1%); the dashed blue lines reflect the 25th and 75th percentiles (1.5 and 2.6%, respectively).

County-Level Analysis

The highest-prescribing county quartile had a mean intensity of 3.4%, 3.1 times higher than 1.1% in the low-intensity quartile (df = 3, F = 3293.8, p < 0.001). High-intensity counties were largely found in Appalachia and south-central states; New York counties were almost entirely low intensity (Fig. 2). The prescribing intensity quartiles varied on most county-level demographics (Table 2), with the largest difference being that high-intensity counties were more rural. High-intensity counties had fewer primary care providers per 100,000 population and fewer had a psychiatrist. High-intensity counties had lower SES and fared more poorly on every health status indicator with the exception of binge alcohol use.

Figure 2
figure 2

County-level Intensity of benzodiazepine prescribing by primary care physicians (n = 3038).

Table 2 Demographic, Healthcare System, Socioeconomic, and Health Status Characteristics Associated with County-Level Benzodiazepine Prescribing Intensity by Primary Care Physicians in 2015 (n = 3038)

In the first regression model that included county-level demographic, SES, and health care system characteristics, higher education and higher income were associated with lower odds of high-intensity prescribing. In the final model, which also included indicators of county-level health status, education and income were no longer significant. However, most characteristics suggesting worse health status were associated with increased odds of high-intensity prescribing, including the proportion of disability-eligible Medicare beneficiaries, days of poor mental health in the past 30, and suicide rate. Higher binge drinking, however, was still associated with lower odds of high-intensity BZD prescribing in the final model.

Physician-Level Analysis

BZD prescribing intensity ranged from 0.6% (SD 0.2) among the bottom PCP quartile to 3.9% (SD 1.8) in the top (df = 3, F = 63,910.2, p < 0.001), a 6.5-fold difference (Table 3). In a multivariable model that adjusted for the local availability of mental health services, females had significantly lower odds of being high-intensity prescribers, while additional years in practice were associated with higher odds. Being a high-intensity opioid prescriber was the single factor most strongly associated with high-intensity BZD prescribing, though every type of high-intensity prescribing was associated with high BZD prescribing intensity, along with simply being a high-volume prescriber. An older patient panel was associated with lower odds of high-intensity BZD prescribing, while higher proportions of patients who were white or receiving a low-income subsidy were associated with increased odds of high-intensity BZD prescribing.

Table 3 Primary Care Physician, Prescribing, and Panel Characteristics Associated with Physician-Level Benzodiazepine Prescribing Intensity, 2015 (n = 109,686)

Sensitivity Analyses

For the county-level model, we conducted sensitivity analyses using total BZD days as the numerator with two alternative denominators (eTable 3): (1) number of Medicare beneficiaries and (2) county population. For the provider-level model, we generated an alternative outcome using a denominator of Medicare beneficiaries seen (eTable 4). The association of characteristics with BZD prescribing intensity was generally consistent (i.e., the 99% CI for effect size overlapped), regardless of the denominator.

DISCUSSION

BZDs are widely prescribed to Medicare beneficiaries, and the intensity of prescribing varies markedly across US counties and individual primary care physicians. Among the top quartile of PCPs, nearly 4% of all Part D medications prescribed in 2015 were BZDs, a rate 6-fold higher than physicians in the lowest quartile.

To our knowledge, this is the first analysis to examine county-level variation in BZD prescribing. In our initial adjusted county-level model, counties with lower educational attainment and lower median income had high-intensity prescribing. In the final model, several other concerning county-level health status characteristics were associated with high-intensity BZD prescribing, including more days of poor mental health, a higher proportion of disability-eligible Medicare beneficiaries, and a higher suicide rate. Because a large amount of BZD prescribing occurs without a mental health or insomnia diagnosis,14,16,32,33 our findings fit Case and Deaton’s hypothesis, with prescribing potentially in response to distress in socioeconomically disadvantaged counties. It was unexpected that binge alcohol use was higher in low-intensity BZD counties, though the prevalence of binge alcohol use is higher in high-income households.34

The high-intensity counties share a striking overlap with county opioid use.35 The Part D PUF cannot determine whether a PCP is prescribing opioids and BZDs to the same patient. However, the fact that the same PCPs appear to be high-intensity prescribers of both medications is potential cause for concern. Our findings are consistent with a recent analysis of opioid and BZD prescribing to Medicaid beneficiaries with opioid use disorders, where prescribing rates varied significantly by provider and were higher in high-poverty counties.36 Some of the socioeconomic characteristics associated with BZD prescribing intensity have also been proposed as determinants in the opioid crisis37; it is unsurprising the same PCPs might be high-intensity prescribers of both.

A clear outlier is New York, which initiated a restrictive BZD triplicate prescription policy 30 years ago that does not permit refills: all but three counties were in the lowest intensity quartile. This cannot be because patients in New York are uniquely free from anxiety or insomnia, or that they have distress-free lives. It is most likely due to the restrictive policy in place, further suggesting that patient diagnosis is not the only factor to determine BZD prescribing. However, initiation of the New York policy also demonstrated that, in the face of pressure to reduce BZD prescribing, physicians may simply shift to even less desirable sedating medications.38

Over 90% of BZD days prescribed were for short-acting medications, such as alprazolam and lorazepam, consistent with prior work showing use of short-acting BZDs has increased.32 This may reflect providers’ response to earlier findings39 and recommendations40,41 that long-acting BZDs elevated risk of hip fracture relative to short-acting agents. However, subsequent analyses have suggested that the risk of short-acting agents is similarly high,5,42 and recent guidelines caution against all BZDs, regardless of half-life.43

The decision to prescribe a medication should depend on the presence of an underlying medical disorder. However, each physician may have their own BZD prescribing threshold and prescribing practice. One PCP characteristic associated with lower odds of high-intensity BZD prescribing was female gender. Other analyses have demonstrated that female clinicians may provider better care,44,44,46 though the reasons are unclear. Roter et al. have found that female physicians provide more psychosocial counseling,47 while “female primary care physicians engage in more communication that can be considered patient centered and have longer visits than their male colleagues”.48 For a patient in distress, it could be that a conversation with a female physician is perceived by both patient and physician as therapeutic in itself.

The single PCP factor that most strongly predicted high-intensity BZD prescribing was being a high-intensity opioid prescriber. However, being a high-intensity prescriber of antibiotics and high-risk medications for the elderly were also strongly associated. It is unclear why the fact that a PCP prescribes a lot of these medications—for unrelated conditions—would be associated with their BZD prescribing practice. However, as with the gender difference, some providers may be more comfortable using education and communication as part of treatment, which is useful whether addressing psychosocial distress or a viral upper respiratory infection.

Our analysis has several limitations. We do not have patient-level clinical data, which may account for some variation in physician- and county-level BZD prescribing, nor individual data on socioeconomic or demographic characteristics, which we analyzed based on county characteristics of the PCP practice location. While our analysis is a near-complete picture of BZD prescribing, because of data redaction the Medicare Public Use Files include the vast majority of—but not all—prescription claims. This means that our findings may be biased by the exclusion of PCPs that see few Medicare Part D patients or who have too few claims for a specific medication to be included in the Detailed PUF. We do not know about the dosing or chronicity of individual BZD prescriptions. Finally, we do not address prescribing by other types of prescribers—both other physician specialties as well as non-physician providers—who may vary widely in their BZD prescribing.

County and PCP characteristics significantly predict the intensity of BZD prescribing to Medicare beneficiaries, even when accounting for availability of mental health treatment resources. Given disparities in healthy aging reported by older adults49 and the non-specific nature of BZD prescribing, it is possible that some BZD prescribing to Medicare beneficiaries is in response to distress. Our county-level findings also demonstrate significant overlap with areas of high-intensity opioid prescribing,35 while high-intensity opioid prescribing is the physician characteristic most predictive of high-intensity BZD prescribing. Given the potential adverse outcomes associated with their use—particularly in combination with opioids—it is important to consider whether this variation in BZD prescribing may be contributing to the uneven distribution of gains in or even worsening of life expectancy in counties across the USA.22