By Christian McDonald, Assistant Professor of Practice
School of Journalism and Media, Moody College of Communication
University of Texas at Austin

This is an analysis of the Primary Cesarean Delivery Rate, Uncomplicated, as defined by IQI 33 Primary Cesarean Delivery Rate, Uncomplicated in 2020:

“First-time Cesarean deliveries without a hysterotomy procedure per 1,000 deliveries. Excludes deliveries with complications (abnormal presentation, preterm delivery, fetal death, multiple gestation diagnoses, or breech procedure).”

This analysis has enough results to avoid defining as “per 1,000 deliveries”. I use a simple rate: Primary Cesarean / All uncomplicated deliveries.

AHRQ definition

This was updated for 2020.

Numerator

Discharges, among cases meeting the inclusion and exclusion rules for the denominator, with any- listed ICD-10-PCS procedure codes for Cesarean delivery (PRCSECP) and without any-listed ICD- 10-PCS procedure codes for hysterotomy (PRCSE2P).

Denominator

All deliveries, identified by any-listed ICD-10-CM diagnosis code for outcome of delivery (DELOCMD).

Exclude cases:

  • with any-listed ICD-10-CM diagnosis codes for abnormal presentation, preterm, fetal death, or multiple gestation (Appendix A: PRCSECD)
  • with any-listed ICD-10-CM diagnosis codes for previous Cesarean delivery (PRVBACD* )
  • with an ungroupable DRG (DRG=999)
  • with missing gender (SEX=missing), age (AGE=missing), quarter (DQTR=missing), year (YEAR=missing) or principal diagnosis (DX1=missing)
library(tidyverse)
library(janitor)
library(DT)
library(tigris)

# suppresses grouping warning
options(dplyr.summarise.inform = FALSE)

Import of deliveries

I start here with “All deliveries, identified by any-listed ICD-10-CM diagnosis code for outcome of delivery (DELOCMD)”. This was processed in 01-process-loop.

# set test flag to FALSE to run production data
test_flag <- F

### test data
path_test <- "data-test/ahrq_del_all_loop_test.rds"

### production data
path_prod <- "data-processed/ahrq_del_cleaned.rds"

### import based on flag
if (test_flag == T) del <- read_rds(path_test) else del <- read_rds(path_prod)

del %>% nrow()
## [1] 1457413

Set up various processing lists

These are lists from the IQI 33 reference used for filtering and such. See 01-process-lists for definitions.

# diagnostic columns
diag_cols <- read_rds("procedures-lists/cols_diag.rds") %>% .$diag

# surgical procedure columns
surg_cols <- read_rds("procedures-lists/cols_surg.rds") %>% .$surg

# appendix a complications list
prcsecd_list <- read_rds("procedures-lists/ahrq_prcsecd.rds") %>% .$prcsecd

# ICD-10-PCS procedure codes for Cesarean delivery
prcsecp_list <- read_rds("procedures-lists/ahrq_prcsecp.rds") %>% .$prcsecp

# Hysterotomy procedure codes: (PRCSE2P)
prcse2p_list <- read_rds("procedures-lists/ahrq_prcse2p.rds") %>% .$prcse2p

# Previous Cesarean delivery diagnosis codes: (PRVBACD)
prvbacd_list <- read_rds("procedures-lists/ahrq_prvbacd.rds") %>% .$prvbacd

Filter deliveries to set denominator

Filter out complicated deliveries

Start by filtering out rows “with any-listed ICD-10-CM diagnosis codes for abnormal presentation, preterm, fetal death, or multiple gestation (Appendix A: PRCSECD)”.

del_ucmp <- del %>% 
  filter_at(
    vars(all_of(diag_cols)),
    all_vars(
      !(. %in% prcsecd_list)
    )
  )

del_ucmp %>% nrow()
## [1] 1315670

Filter out previous Cesarean deliveries

For this analysis, I need to additionally exclude outcomes “with any-listed ICD-10-CM diagnosis codes for previous Cesarean delivery (PRVBACD)”. This prvbacd_list comes from the AHRQ specification.

del_ucmp <- del_ucmp %>% 
  filter_at(
    vars(all_of(diag_cols)),
    all_vars(
      !(. %in% prvbacd_list)
    )
  )

del_ucmp %>% nrow()
## [1] 1052162

Filter out ungroupable DRG (DRG=999)

There are two Diagnosis Related Groups in the data, MS_DRG: Centers for Medicare and Medicaid Services (CMS) Diagnosis Related Group (DRG), as assigned for hospital payment for Medicare beneficiaries; and APR_DRG: All Patient Refined (APR) Diagnosis Related Group (DRG) as assigned by 3M APR-DRG Grouper). I filter out records with “999” in either of them.

del_ucmp <- del_ucmp %>% 
  filter(
    (APR_DRG != "999"),
    (MS_DRG != "999")
  )

del_ucmp %>% nrow()
## [1] 1052162

This completes our “denominator” set of data. This is saved in del_ucmp.

Code records for numerator

To review: Discharges, among cases meeting the inclusion and exclusion rules for the denominator, with any- listed ICD-10-PCS procedure codes for Cesarean delivery (PRCSECP) and without any-listed ICD- 10-PCS procedure codes for hysterotomy (PRCSE2P).

NOTE: In 2019, MS-DRG codes for Cesarean delivery changed and that is not reflected in the AHRQ June 2019 specification. See 01-process-lists for more information.

How I handle this:

Create column PRCSECP for cesarean deliveries

del_ucmp <- del_ucmp %>% 
  mutate(
    PRCSECP = case_when(
      PRINC_SURG_PROC_CODE %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_1 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_2 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_3 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_4 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_5 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_6 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_7 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_8 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_9 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_10 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_11 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_12 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_13 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_14 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_15 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_16 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_17 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_18 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_19 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_20 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_21 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_22 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_23 %in% prcsecp_list ~ TRUE,
      OTH_SURG_PROC_CODE_24 %in% prcsecp_list ~ TRUE,
      TRUE                     ~  FALSE
    )
  )

del_ucmp %>% 
  count(PRCSECP)

Create column based on hysterotomy

… “without any-listed ICD-10-PCS procedure codes for hysterotomy (PRCSE2P).”

It is not surprising that there are zero results here as these are “Abortion of Products of Conception” procedures . Abortions are filtered out in the method used to define deliveries.

I did test the code against unfiltered data to ensure that it does catch the PRCSE2P codes.

This is a code block that I would like to refactor as a loop through the surg_cols list.

#case_when version
del_ucmp <- del_ucmp %>% 
  mutate(
    PRCSE2P = case_when(
      PRINC_SURG_PROC_CODE %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_1 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_2 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_3 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_4 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_5 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_6 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_7 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_8 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_9 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_10 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_11 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_12 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_13 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_14 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_15 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_16 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_17 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_18 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_19 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_20 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_21 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_22 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_23 %in% prcse2p_list ~ TRUE,
      OTH_SURG_PROC_CODE_24 %in% prcse2p_list ~ TRUE,
      TRUE                     ~  FALSE
    )
  )

# peek at results
del_ucmp %>% 
  count(PRCSE2P)

Create Primary Cesarean indicator column

Created the PCSEC column and set as TRUE only if both PRCSE2G is TRUE and PRCSE2P is FALSE.

del_csec <- del_ucmp %>% 
  mutate(
    PCSEC = if_else((PRCSECP == T & PRCSE2P == F), T, F)
  )

# peek at data
del_csec %>% count(PCSEC)

This completes our “numerator” set of data. This is saved in del_csec.

Create Medicaid indicator column

Creates an MC column to indicate TRUE if FIRST_PAYMENT_SRC value is “MC”, which is the indication for Medicaid Also changes NA values to FASLE since those would be unknown.

ahrq_pcsec <- del_csec %>% 
  mutate(
    MC = if_else(FIRST_PAYMENT_SRC == "MC", T, F),
    MC = if_else(is.na(MC), F, MC)
  )

This completes the processing of the data. ahrq_pcsec is our working dataframe.


Table: Primary Cesarean rate by hospital, combined years

Get Primary Cesarean rate by hospital from all years combined. Hospitals with fewer than 120 deliveries over the data set excluded.

ahrq_pcsec_rate_hosp <- ahrq_pcsec %>% 
  group_by(THCIC_ID, PROVIDER_NAME, PCSEC) %>% 
  summarize(CNT = n()) %>% 
  pivot_wider(names_from = PCSEC, values_from = CNT) %>% 
  rename(
    NPCSEC_CNT = "FALSE",
    PCSEC_CNT = "TRUE"
  ) %>% 
  mutate(
    TOTAL = NPCSEC_CNT + PCSEC_CNT
  ) %>% 
  filter(TOTAL >= 120) %>% 
  mutate(
    PCRATE = round_half_up((PCSEC_CNT / TOTAL) * 100,1)
  ) %>%
    arrange(PCRATE %>% desc())

ahrq_pcsec_rate_hosp %>% datatable()

Table: Primary Cesarean rates, by hospital and year

Hospitals with fewer than 30 total deliveries within a given year are excluded.

ahrq_pcsec_rate_hosp_yr <- ahrq_pcsec %>% 
  group_by(YR, THCIC_ID, PROVIDER_NAME, PCSEC) %>% 
  summarize(CNT = n()) %>% 
  pivot_wider(names_from = PCSEC, values_from = CNT) %>% 
  rename(
    NPCSEC_CNT = "FALSE",
    PCSEC_CNT = "TRUE"
  ) %>% 
  mutate(
    TOTAL = NPCSEC_CNT + PCSEC_CNT
  ) %>% 
  filter(TOTAL >= 30) %>% 
  mutate(
    PCRATE = round_half_up((PCSEC_CNT / TOTAL) * 100,1)
  ) %>%
    arrange(PCRATE %>% desc())

ahrq_pcsec_rate_hosp_yr_table <- ahrq_pcsec_rate_hosp_yr %>% 
  arrange(YR) %>% 
  select(YR, THCIC_ID, PROVIDER_NAME, PCRATE) %>% 
  pivot_wider(names_from = YR, values_from = PCRATE) %>% 
  arrange(`2019` %>% desc())

ahrq_pcsec_rate_hosp_yr_table %>% datatable()

Laredo hospitals

Primary Cesarean rate combined years, Laredo hospitals

The TRUE value below is the percentage of deliveries involving all low-risk mothers (without prior Cesareans) that result in a Cesarean delivery.

The numbers in parenthesis are counts.

ahrq_pcsec %>% 
  filter(
     str_detect(PROVIDER_NAME, "Laredo")
  ) %>% 
  tabyl(PROVIDER_NAME, PCSEC) %>% 
  adorn_percentages() %>% 
  adorn_pct_formatting() %>% 
  adorn_ns()

Primary Cesarean rate by year, Laredo hospitals

ahrq_pcsec_rate_hosp_yr_table %>% 
  filter(str_detect(PROVIDER_NAME, "Laredo"))

Chart: Primary Cesarean Rates by year, Laredo

ahrq_pcsec_rate_hosp_yr %>% 
  filter(str_detect(PROVIDER_NAME, "Laredo")) %>% 
  ggplot(aes(YR, PCRATE)) +
  geom_line(aes(group = PROVIDER_NAME, color = PROVIDER_NAME)) +
  expand_limits(y = c(0,40)) +
  theme(legend.position="bottom", legend.box = "vertical") +
  labs(title = "Primary Cesarean rate by year, Laredo", x = "YEAR", y = "Primary Cesarean Rate")

Chart: Primary Cesarean rates by quarter, Laredo

ahrq_pcsec_rate_hosp_qr <- ahrq_pcsec %>% 
  group_by(DISCHARGE, PROVIDER_NAME, PCSEC) %>% 
  summarize(CNT = n()) %>% 
  pivot_wider(names_from = PCSEC, values_from = CNT) %>% 
  rename(
    NPCSEC_CNT = "FALSE",
    PCSEC_CNT = "TRUE"
  ) %>% 
  mutate(
    TOTAL = NPCSEC_CNT + PCSEC_CNT
  ) %>% 
  mutate(
    PCRATE = round_half_up((PCSEC_CNT / TOTAL) * 100,1)
  ) %>%
    arrange(PCRATE %>% desc())

ahrq_pcsec_rate_hosp_qr %>% 
  filter(str_detect(PROVIDER_NAME, "Laredo")) %>% 
  ggplot(aes(DISCHARGE, PCRATE)) +
  geom_line(aes(group = PROVIDER_NAME, color = PROVIDER_NAME)) +
  expand_limits(y = c(0,40)) +
  theme(legend.position="bottom", legend.box = "vertical", axis.text.x=element_text(angle = -45, hjust = 0)) +
  labs(title = "Primary Cesarean rate by quarter, Laredo", x = "QUARTER", y = "Primary Cesarean Rate")

Table: Considering Medicaid

The questions we wanted to answer here is: What is the percentage of patients on Medicaid that get a primary cesarean? How does that compare to patients not on medicare?

MEDI_W_PCSEC is the percentage of Medicaid (as first payment option) patients that get primary cesareans. NOT_MEDI_W_PCSEC is the percentage of patients NOT using Medicaid that get primary cesareans. Printed table filters out hospitals with less than 20 Medicaid births.

mc_compare <- ahrq_pcsec %>% 
  group_by(THCIC_ID, PROVIDER_NAME) %>% 
  summarize(
    MCT = sum(MC == T), # medicaid count
    MCPCT = sum(MC == T & PCSEC == T), # medicaid with primary csection
    MEDI_W_PCSEC = ((MCPCT/MCT) * 100) %>% round_half_up(1),
    NMCT = sum(MC != T), # not medicare count
    NMCPCT = sum(MC != T & PCSEC == T), #not medicaid with primary csection
    NOT_MEDI_W_PCSEC = ((NMCPCT / NMCT) * 100) %>% round_half_up(1),
    PDIFF =  MEDI_W_PCSEC - NOT_MEDI_W_PCSEC
  )

# select fewer cols for easy comparison
mc_compare %>% 
  filter(MCT >= 20) %>% # at lest 20 medicaid patients
  select(THCIC_ID, PROVIDER_NAME, MEDI_W_PCSEC, NOT_MEDI_W_PCSEC, PDIFF) %>%
  arrange(PDIFF %>% desc()) %>% 
  datatable()

### Considering Medicaid, Laredo hospitals

mc_compare %>% 
  filter(str_detect(PROVIDER_NAME, "Laredo"))

Patient County

A look at Primary cesarean by the county of the patient. Note the hospital could be elsewhere.

tx_fips <- fips_codes %>% 
  filter(state == "TX") %>% 
  select(county_code, county)
# peek
tx_fips %>% head(2)

Primary cesarean by county, full data

Currently not filtered for any minimum caseload since it is a 3+ year period.

ahrq_pcsec_rate_county <- ahrq_pcsec %>% 
  group_by(PAT_COUNTY, PCSEC) %>% 
  summarize(CNT = n()) %>% 
  ### join for county names
  left_join(tx_fips, by=c( "PAT_COUNTY" = "county_code")) %>% 
  ungroup() %>% 
  rename(COUNTY = county) %>% 
  pivot_wider(names_from = PCSEC, values_from = CNT) %>% 
  rename(
    NPCSEC_CNT = "FALSE",
    PCSEC_CNT = "TRUE"
  ) %>% 
  mutate(
    TOTAL = NPCSEC_CNT + PCSEC_CNT
  ) %>% 
  ### filter by cases
  # filter(TOTAL >= 10) %>%
  mutate(
    PCRATE = round_half_up((PCSEC_CNT / TOTAL) * 100,1)
  ) %>%
    arrange(PCRATE %>% desc())

ahrq_pcsec_rate_county %>% 
  select(-PAT_COUNTY) %>% 
    datatable()

Primary cesarean by county, by year

This table excludes records with fewer than 10 cases in a county in a given year, which reduces the number of available counties.

ahrq_pcsec_rate_county_yr <- ahrq_pcsec %>% 
  group_by(YR, PAT_COUNTY, PCSEC) %>% 
  summarize(CNT = n()) %>% 
  ### join for county names
  left_join(tx_fips, by=c( "PAT_COUNTY" = "county_code")) %>% 
  ungroup() %>% 
  rename(COUNTY = county) %>% 
  pivot_wider(names_from = PCSEC, values_from = CNT) %>% 
  rename(
    NPCSEC_CNT = "FALSE",
    PCSEC_CNT = "TRUE"
  ) %>% 
  mutate(
    TOTAL = NPCSEC_CNT + PCSEC_CNT
  ) %>% 
  ### filter by caseload
  filter(TOTAL >= 10) %>%
  mutate(
    PCRATE = round_half_up((PCSEC_CNT / TOTAL) * 100,1)
  ) %>%
    arrange(PCRATE %>% desc())

ahrq_pcsec_rate_county_yr_table <- ahrq_pcsec_rate_county_yr %>% 
  arrange(YR) %>% 
  select(YR, COUNTY, PCRATE) %>% 
  pivot_wider(names_from = YR, values_from = PCRATE) %>% 
  arrange(`2019` %>% desc())

ahrq_pcsec_rate_county_yr_table %>%  datatable()

Summaries

Primary Cesarean rate statewide

This compares all the cases and gives us a percentage of those that are Primary Cesarean births. The TRUE value is the percentage of deliveries involving all low-risk mothers (without prior Cesareans) that result in a Cesarean delivery.

ahrq_pcsec %>% 
  tabyl(PCSEC) %>% 
  rename(count = n) %>% 
  adorn_pct_formatting()
ahrq_pcsec_rate_tx_summary <- ahrq_pcsec %>% 
  group_by(PCSEC) %>% 
  summarize(CNT = n()) %>% 
  pivot_wider(names_from = PCSEC, values_from = CNT) %>% 
  rename(
    NPCSEC_CNT = "FALSE",
    PCSEC_CNT = "TRUE"
  ) %>% 
  mutate(
    TOTAL = NPCSEC_CNT + PCSEC_CNT
  ) %>% 
  mutate(
   SUMMARY = "TX",
   CATEGORY = "PRIMARY_CESAREAN",
   MEASUREMENT = "RATE",
   VALUE = round_half_up((PCSEC_CNT / TOTAL) * 100,1) # RATE
  ) %>%
  select(SUMMARY, CATEGORY, MEASUREMENT, VALUE)

ahrq_pcsec_rate_tx_summary

Primary Cesarean rate statewide by year

Excludes hospitals with fewer than 30 deliveries a year.

ahrq_pcsec_rate_tx_yr_summary <- ahrq_pcsec %>% 
  group_by(PCSEC, YR) %>% 
  summarize(CNT = n()) %>% 
  pivot_wider(names_from = PCSEC, values_from = CNT) %>% 
  rename(
    NPCSEC_CNT = "FALSE",
    PCSEC_CNT = "TRUE"
  ) %>% 
  mutate(
    TOTAL = NPCSEC_CNT + PCSEC_CNT
  ) %>% 
  filter(TOTAL >= 30) %>% 
  mutate(
   SUMMARY = "TX",
   CATEGORY = "PRIMARY_CESAREAN",
   MEASUREMENT = "RATE",
   VALUE = round_half_up((PCSEC_CNT / TOTAL) * 100,1) # RATE
  ) %>%
  select(YR, SUMMARY, CATEGORY, MEASUREMENT, VALUE)

ahrq_pcsec_rate_tx_yr_summary

Primary Cesarean rate by hospital: Averaged, by year

ahrq_pcsec_rate_hosp_yr_summary <- ahrq_pcsec_rate_hosp_yr %>% 
  ungroup() %>% 
  group_by(YR) %>% 
  summarize(
    SUMMARY = "HOSPITAL",
    CATEGORY = "PRIMARY_CESAREAN",
    MEASUREMENT = "MEAN_OF_RATE",
    # MEDIAN = median(PCRATE),
    VALUE = round_half_up(mean(PCRATE),1) #MEAN
  )

ahrq_pcsec_rate_hosp_yr_summary

Write files

Writing out aggregate files. Here is a list of CSVs exports:

And some summary data, including:

# if (test_flag == F) ahrq_pcsec %>% write_rds("data-processed/ahrq_pcsec.rds")

if (test_flag == F) ahrq_pcsec_rate_tx_yr_summary %>%
  write_rds("data-processed/ahrq_pcsec_rate_tx_yr_summary.rds")

if (test_flag == F) ahrq_pcsec_rate_hosp_yr_summary %>%
  write_rds("data-processed/ahrq_pcsec_rate_hosp_yr_summary.rds")

if (test_flag == F) ahrq_pcsec_rate_hosp_yr %>%
  write_rds("data-processed/ahrq_pcsec_rate_hosp_yr.rds")

if (test_flag == F) ahrq_pcsec_rate_hosp_yr %>% 
  arrange(PROVIDER_NAME,YR) %>%
  write_csv("exports/ahrq_pcsec_rate_hosp_yr.csv")

if (test_flag == F) ahrq_pcsec_rate_county_yr %>%
  write_csv("exports/ahrq_pcsec_rate_county_yr.csv")

# A klaxon to indicate the processing is complete
beepr::beep(4)