Phenotype Discovery and Geographic Disparities of Late-Stage Breast Cancer Diagnosis across U.S. Counties: A Machine Learning Approach

Because of increased screening and awareness of symptoms in the past few decades, the rate of advanced stage breast cancer (for which prognosis is markedly poorer) has decreased dramatically (1). However, breast cancer remains the second leading cause of cancer death in U.S. women (1). In addition, segments of our population continue to bear a disproportionate burden and experience high rates of late-stage diagnosis.

Previous studies have linked place-based disparities in late-stage breast cancer (LSBC) to numerous contextual characteristics, including socioeconomic status (2, 3), neighborhood deprivation (4, 5), access to screening services (5–9), and availability of primary care physicians (10–12). Studies of LSBC traditionally have considered all predictor variables independently using a parametric regression framework (2–8, 10). However, such a framework does not lend itself to identifying homogeneous subgroups, nor is it efficient in detecting effect measure modification in predicting LSBC. For example, one variable might be important in predicting LSBC in a certain subgroup of a population but may not be as important in other subgroups of the population. Several studies have recognized this issue by stratifying their populations into subgroups by factors such as race/ethnicity (4, 11), income or poverty level (11, 12), and urban-rural status (4, 13). However, the decision on which attributes to stratify the data may be subjective, and the selection of threshold values often are arbitrary. In addition, the complexity of a model increases substantially when considering stratifications along multiple dimensions. Therefore, more advanced solutions are needed to explore the correlates of LSBC among different subgroups of the population.

In this study, we identified “phenotypes,” or combinations of county characteristics that predict the percentage of women with breast cancer presenting with late-stage disease. We applied the machine learning technique known as classification and regression tree (CART) with a broad range of county-level characteristics harvested from various sources. This resulted in the classification of counties into phenotypes based on the most important predictors of LSBC. We then examined the geographic distribution of counties with high-risk phenotypes. These findings were further examined using random forest.

Identifying specific clusters of characteristics associated with late-stage diagnosis acknowledges the complex relationships among selected drivers of cancer disparities. It also offers researchers and practitioners a better framework for addressing disparities across heterogeneous and more highly specified groups.

This study used cancer incidence data from the Surveillance, Epidemiology, and End Results (SEER) program, a resource from the NCI. Data in the SEER program cover 34.6% of the U.S. population with 97% completeness within SEER regions (14). These data cover geographically diverse regions of the country and are broadly representative of the U.S. population along the dimensions of poverty and education (15). The SEER*Stat software was used to query and extract the data. Given the deidentified nature of the data, the Case Western Reserve University Institutional Review Board determined that this work did not involve human subjects research and was thus exempted from review.

The study included 20 of the 21 registries in the SEER program. We excluded the Alaska Native Tumor Registry because it does not cover cancer cases of all demographic groups (16). The included registries cover the U.S. states of California (Greater California, Los Angeles, San Francisco-Oakland, and San Jose-Monterey registries), Connecticut, Georgia (Atlanta, Greater Georgia, and Rural Georgia registries), Hawaii, Idaho, Iowa, Kentucky, Louisiana, Massachusetts, New Jersey, New Mexico, New York, Utah, and metropolitan areas of Detroit and Seattle-Puget Sound. A total of 732 counties or equivalents (i.e., parishes in Louisiana; for convenience, “counties” is used in the rest of the text) from the 20 SEER registries were included in the study. The outcome of interest was the county-level percentage of LSBC among women diagnosed with invasive breast cancer during a 10-year period between 2009 and 2018. For Massachusetts, only cases from 2009 through 2017 were included because the stage variable in 2018 was not available at the time of the study. Further explanation of the Massachusetts data can be found in Supplementary Fig. S1 and Fig. S2. We used the “Combined Summary Stage” variable from the SEER*Stat software, which classifies tumors into five stages: in situ, local, regional, distant, and unknown stages (excluding in situ and unknown stage cases from our analysis). We collapsed the stage variable to “early stage” (which includes local stage only), and “late stage” (which includes regional and distant stages). To mitigate stochastic variations in percentages of LSBC, we excluded counties with fewer than 16 late-stage cases over the study period, a strategy also adopted by SEER in displaying cancer statistics, resulting in a total of 680 counties in the study. For individuals with multiple tumor records, we selected only the first record.

County-level characteristics were harvested from the Census-American Community Survey (ACS; ref. 17), County Health Rankings & Roadmaps (CHR; ref. 18), Area Health Resources Files (AHRF; ref. 19), Behavioral Risk Factor Surveillance System (BRFSS; ref. 20), and FDA (21), as well as from SEER*Stat. In total, 53 variables were included in our models (Table 1). These variables were selected on the basis of several domains of health care resources, behavioral risk factors, population health status, demographic compositions, and other measures of social determinants of health including income, education, occupation, housing, transportation, and neighborhood safety. Their relationships with diagnosis and treatment of breast cancer were explained under the Conceptual Reason for Inclusion column in Table 1. We temporally harmonized the data sources by selecting the years of these variables that overlapped with or were closest to the mid-years of the breast cancer data, with the assumption that no significant secular trends would substantively change any of the factors described in Table 1.

Given the nature of the study, all measures were aggregated at the county level; therefore, we did not account for variables at the individual level.

As described in detail below, machine learning methods (including CART and random forest) and geographic information systems were used to accomplish the objectives of this study. The outcome for all models was county-level percentage of LSBC and all variables in Table 1 were included as candidate predictors.

CART uses conditional inference that recursively partitions data into smaller, more homogeneous groups characterized by combinations of predictors (22, 23). At each split, the data are divided into two homogeneous groups according to a threshold value of one of the predictors, a predictor that results in the two groups with greatest difference in the outcome (23). The splitting procedure is repeatedly applied for each of the split groups by selecting one of the predictors that holds the lowest P-value based on Pearson correlation test if the predictor variable is numeric, or Kruskal–Wallis test if the predictor variable is categorical (22). This procedure continues until all possible splits are exhausted or until some stopping criteria is met. In this study, we set the following stopping criteria: a maximum tree depth of 6 splits, a minimum number of 80 counties in a terminal node, and lack of statistical significance for variable splits (α > 0.05). We also conducted a sensitivity analysis of the CART model with a minimum number of 20 counties in a terminal node.

CART was used to identify phenotypes associated with differing levels of LSBC. We defined a phenotype of LSBC as the combination of characteristics along a top-down path of a tree to a terminal node (a node without any further split) which includes a group of homogenous counties with similar percentages of LSBC. Conceptually, this results in the identification of the combinations of county characteristics that predict the percentage of women with breast cancer presenting with late-stage disease.

Next, we visualized the identified phenotypes using geographic information systems and examined their distribution among regions. To optimize interpretability, the minimum number of counties in a terminal node of the CART model was set to 80 to limit the number of phenotypes presented on the map.

Random forest analysis was used to determine whether our CART model captured the most important variables in predicting the percentage of LSBC. While CART has advantages in variable identification and group classification, a major disadvantage of this single-tree model is that it is sensitive to changes in the data. Hence, the entire tree could be altered if, for example, additional counties are included in the model. In contrast, random forest analysis is more “stable” (24). It uses the same algorithm as CART, but instead of relying on only one tree, the algorithm creates and aggregates an ensemble of trees using random variable selection and bootstrap sampling (24). It then takes an overall average of these tree models' outputs as a prediction. Next, the mean decrease in accuracy was used to calculate variables' relative importance in predicting the outcome. We created 200,000 trees with all predictor variables included in the analysis. The number of variables randomly sampled as candidates at each split was set equal to the number of splits in the results of the CART model.

Because of the nature of CART and random forest, in that both algorithms select only one of many variables at each split of their trees, they can handle highly correlated variables. However, the associations among candidate predictors remain unknown. For highly correlated predictors, while CART selects only the one that most significantly splits the group, it does not suggest that the rest of the predictors are not predictive of the outcome. Random forest partially addresses this issue with the rankings of variable importance. To further explore the associations among candidate predictors, we conducted a correlation analysis using Pearson correlation coefficients among all splitting variables in the CART model and top 10 variables in the variable importance plot of the random forest, as well as the variable representing the proportion of women in race/ethnic minorities.

SAS v9.4 and R v3.6.1 were used for the analyses, and ArcGIS Pro v2.7.0 and Tableau v2021.1 were used for mapping and visualization. R packages “rpart,” “partykit,” and “randomForest” were used for conducting machine learning analyses.

Within the 680 included counties between 2009 and 2018, there were 812,048 women diagnosed with invasive breast cancer, among whom 276,305 (34.0%) were diagnosed at a late stage. The median percentage of LSBC among the counties was 35.4%. The geographic distribution of counties by percentage of LSBC is presented in Fig. 1.

We observed that counties in the Northeast states (New York, Massachusetts, Connecticut, and New Jersey) had lower percentages of LSBC compared with the majority of those in the south and west. The box plots showing the county distribution of LSBC by region are included in Supplementary Fig. S2.

The results of the CART analysis are shown in Fig. 2. Each path down to a terminal node of the tree represents a phenotype of LSBC with corresponding characteristics. P values on the splitting nodes (<0.05) suggest that the groups of counties split by the thresholds of the variables are statistically significantly different from each other in terms of percentage of LSBC. Counties within the same terminal node have similar percentages of LSBC and belong to the same phenotype. Phenotypes are classified into low-risk (LR), medium-risk (MR), and high-risk (HR) by their median percentage of LSBC among counties, with LR 1 having the lowest percentage and HR 3 having the highest percentage.

The results show that among all candidate predictors, CART selected five variables as splitting nodes, with percentage of uninsured women ages 18–64 on the top, followed sequentially by percentage of mammography use among women aged 67–69 enrolled in Medicare, Area Deprivation Index (ADI), Urban Influence Code (UIC), percentage of people under poverty, and per capita income. The ADI is an index of social deprivation calculated from census variables, which incorporates 17 separate factors covering domains of education, employment, income, housing (costs and crowding), and transportation access (25). The ADI ranges from 30.5 to 154.5 among all counties in the study area (mean: 98.8, median: 98.5). A higher ADI indicates that the county is more deprived. Counties with an ADI greater than 99.7, as shown in the CART output, means that they were more deprived than 52.2% of all counties. Counties with a higher UIC were more rural. Specifically, counties with UIC less than or equal to 3 were metropolitan counties or micropolitan counties adjacent to a large metropolitan county.

HR 3 is the phenotype associated with the highest median percentage of LSBC (40.1%). It includes 89 counties and is characterized as having a higher percentage of uninsured middle-aged women (>11.6%), greater area deprivation (ADI > 99.7), and more people under poverty (>26.1%). HR 2 is the phenotype associated with the second highest median percentage of LSBC (38.4%). HR 2 has the same levels of percentage of uninsured middle-aged women and area deprivation compared with HR 1 but has a lower poverty rate (≤26.1%) and higher per capita income (>32,946 U.S. dollars). HR 1 has a slightly lower median percentage of LSBC than HR 2 (37.0% vs. 38.4%), and their only difference is in per capita income.

MR 1 has the largest group of counties (n = 139) among all phenotypes and has a median percentage of LSBC that is close to that of the overall study area (35.5% vs. 35.4%). Counties of MR 1 has greater rates of uninsured middle-aged women (>11.6%) but lower area deprivation (≤99.7).

Phenotypes LRs 1, 2, and 3, have better outcomes of LSBC compared with HRs 1, 2, 3, and MR 1. The key difference between LR phenotypes and MR or HR phenotypes is the top splitting variable (uninsured middle-aged women), with a threshold of 11.6% that separates the tree into two large branches. The variable that separates LR 1 from LRs 2 and 3 is mammography use among Medicare beneficiaries ages 67–69 years (>68.1%). When mammography use is at a lower level (≤68.1%), UIC comes into play and differentiates LR 2 (urban) with LR 3 (rural).

The sensitivity analysis of CART (with a minimum of 20 counties in a terminal node) presents additional splits compared with the main model, which include availability of obstetrics and gynecologists, access to exercise opportunities, breast cancer incidence rate, availability of primary care physicians, women in professional occupations, and Medicare eligibility (Supplementary Fig. S3). Note that UIC no longer appears in this model, suggesting potential correlations between this and other variables.

Figure 3 shows the geographic distribution of the phenotypes identified in the main CART model. We observed strong regional variations in the composition of phenotypes. For example, Massachusetts is the only region that has only one phenotype, which is LR 1, the most favorable phenotype of LSBC, while California, Kentucky, and Georgia tend to have the greatest variability of phenotypes. LR 1 also covers a large number of counties in Upstate New York and Iowa. LR 2 is common in metropolitan areas of San Francisco Bay Area (CA), New York (NY, NJ, and CT), Seattle, Detroit, Louisville (KY), Lexington (KY), metropolitan areas in Upstate New York, and urban areas in Hawaii and Iowa. LR 3 is mostly found in rural areas in Iowa, Kentucky, and Upstate New York. No counties of any LR phenotypes appear in Georgia, Louisiana, and Idaho, and only one county in New Mexico (Los Alamos County). MR 1 covers most counties in Utah, Idaho, metropolitan Atlanta in Georgia, and Coastal California. The three HR phenotypes (with higher rates of uninsured women, greater area deprivation, and worst outcome of LSBC) are found mostly in rural counties of the southern and western states of Georgia, Louisiana, New Mexico, Idaho, and California. We also summarized these descriptions by phenotype in Table 2 and visualized the distribution of phenotypes by region in Supplementary Fig. S4.

The ranking of variables by their importance for predicting percentages of LSBC based on the random forest analysis is shown in Fig. 4. This dot plot ranks the variables in descending order relative to the most important predictor. Among all variables, uninsured rate among middle-aged women and mammography use among Medicare beneficiaries ages 67 to 69 are the two most important predictors of percentage of LSBC. Other important predictors are adult obesity, percentage of children in poverty, ADI, percentage of people under 200% of poverty, percentage of female-headed households, and percentage of people with poor or fair health. Detailed information regarding the calculation of variable importance with intermediate results, as well as an alternative measure of mean decrease in node impurity of the random forest analysis are included in Supplementary Table S1.

The correlation matrix in Supplementary Table S2 suggests that some of the candidate predictors were highly correlated with absolute values of correlation coefficient (cc) greater than 0.7. Notably, ADI was highly correlated with multiple variables, including women with high school degree or above, per capita income, teen birth, percent people with poor of fair health, poverty, child poverty, and population below 200% of the poverty level (cc: −0.77, −0.74, 0.81, 0.85, 0.88, 0.89, and 0.90, respectively). The correlation between the percentage of women in race/ethnic minority groups and other variables were moderate, with the largest cc at 0.65 with female-headed households.

Using CART analysis, our study identified low, medium, and high risk phenotypes of LSBC consisting of county-level characteristics that were predictive of LSBC. These phenotypes were defined by combinations of indicators of uninsured rate, mammography use, area deprivation, urban-rural status, poverty rate, and per capita income. Among those, the importance of uninsured rate and mammography use was further evidenced by their top rankings in the random forest analysis.

When a smaller number of counties was allowed in a terminal node, additional characteristics came into play in the CART model. Surprisingly, the percentage of racial and ethnic minorities, a factor frequently emphasized in previous studies of LSBC (3, 4, 7, 10–13, 26–28), did not appear in the results of the CART models and was only ranked 13th in the random forest plot. Our study suggests that other contextual factors might have played more critical roles than the constructs of race and ethnicity themselves.

We also recognized the correlations among variables. Some predictors, although did not appear in CART, could still have important implications due to their strong correlations with the outcome and splitting variables. For example, several measures of population socioeconomic status were highly correlated with ADI, most of which also had relatively high rankings in the random forest plot. In addition, despite the absence of adult obesity as a splitting variable in CART, it was identified as a top-ranking variable in predicting LSBC by the random forest analysis, which is notable because it was the only physiologic factor associated with LSBC.

The spatial distribution of phenotypes shows that LR phenotypes, or those comprised of counties with a relatively lower percentage of LSBC, were prevalent in northeastern states, Iowa, and select metropolitan areas, whereas HR phenotypes were mostly observed in the southern states of Georgia, Louisiana, New Mexico, and some rural areas in other states. The unbalanced distribution of phenotypes suggests that there were geospatial disparities in LSBC, and these disparities were strongly associated with population characteristics along multiple dimensions.

The geographic clustering of phenotypes suggests that the association between LSBC and various socioeconomic characteristics may be mediated by geographic-related factors. We observed strong differences in the phenotype composition of states, such as northeastern states versus southern and western states. These state-level disparities could be related to the differences in state-specific culture and policies. For example, states with stringent eligibility criteria for Medicaid enrollment may observe higher rates of uninsured among people with low socioeconomic status compared with states that do not. This was especially the case after 2014 when some states expanded their Medicaid programs under the Patient Protection and Affordable Care Act, while others had not by the end of the study period (29).

We also noted within-state variations, especially in California, Kentucky, and Georgia, where counties in large metro areas generally had more favorable phenotypes than their rural neighbors. This may be due to the large differences in demographic and socioeconomic characteristics and availability of resources between rural and urban areas as indicated in earlier studies (9, 13). However, not all urban areas outperformed rural areas in LSBC. For example, several rural counties in northern California had LR phenotypes, while Los Angeles County and San Diego County, the state's two most populous counties, were of the MR phenotype. This suggests that there was a complex relationship between LSBC and urbanicity, which could involve associations with other variables, such as uninsured rates and mammography use as shown in the CART models.

To our knowledge, this is the first study that combines machine learning and geographic methods to explore the association between LSBC and population characteristics using cancer registry data. Prior studies have investigated associations between LSBC and various demographic, socioeconomic, and behavioral factors using geospatial analyses and parametric regression models (2–8, 10–12, 30–32). However, these studies did not evaluate the roles of these factors in specific subgroups of the population. For example, one indicator may be important to a certain group of a population but less important to other subgroups of the population. As indicated in Fig. 2, when mammography use was less than or equal to 68.1%, UIC was able to significantly distinguish LR 2 and LR 3 based on the percentage of LSBC. However, it was no longer used as a splitting node when mammography use was higher. Similarly, when uninsured rate was greater than 11.6%, neither mammography use nor rural-urban status was as important as the other variables that appeared in the right branch of the tree. By using the CART-defined phenotypes, our study identifies correlates that are specific to various subgroups of the population that share common characteristics.

There are several methodologic strengths that lend confidence to the study results. First, compared with parametric regression models (such as logistic regression), CART and random forest methods can deal with a greater number of predictor variables simultaneously without concerns about outliers, multicollinearity, heteroscedasticity, or distributional error structures that affect parametric procedures. Second, both CART and random forest methods are able to handle highly correlated data due to their variable selection and bootstrap sampling strategies. Finally, the identified phenotypes capture both the outcome and top predictor variables, allowing the examination of geographic patterns of LSBC from multiple aspects.

The main limitation of the current study is that individual regions in the study area were disconnected in geography, preventing the use of geospatial analyses designed for contiguous regions, such as the spatial scan statistics (33), the Local Indicators of Spatial Association (LISA; ref. 34), and the recently developed geographically weighted random forest (35, 36). A geospatial analysis on a contiguous study area, such as the contiguous United States, could help us understand a broad scope of the disparities in LSBC and provide insights into why neighboring areas present similar or different patterns. However, this limitation was tempered by the wide distribution of the study area across the country which covered diverse populations that were comparable with the overall United States. Another limitation is that we were not able to incorporate additional individual-level characteristics in the analysis. Nevertheless, the findings of the study are still valid as the nature of our study was to discover community-level drivers of LSBC disparities.

In summary, our study shows that the use of machine learning and geographic methods is a promising avenue for future disparities research. The findings of our study suggest that the disparities of LSBC are associated with multiple characteristics of the population, and these associations vary greatly across geographies. Local interventions to reduce late-stage diagnosis of breast cancer, such as community education and outreach programs, should consider the characteristics of their communities; thus translational and implementation researchers should consider phenotype-tailored interventions.

W. Dong reports grants from the NCI during the conduct of the study; and grants from Cleveland Clinic Foundation outside the submitted work. U. Kim reports grants from National Center for Advancing Translational Sciences, National Institute of General Medical Sciences, and PhRMA Foundation during the conduct of the study. S.M. Koroukian reports grants from Centers for Disease Control and Prevention, Ohio Medicaid Technical Assistance and Policy Program, NIH, American Cancer Society, and Contracts from Cleveland Clinic Foundation, including a subcontract from Celgene Corporation outside the submitted work; and her spouse has stock ownership in American Renal Associates. No disclosures were reported by the other authors.

W. Dong: Conceptualization, resources, data curation, software, formal analysis, validation, investigation, visualization, methodology, writing–original draft, writing–review and editing. W.P. Bensken: Conceptualization, resources, investigation, methodology, writing–original draft, writing–review and editing. U. Kim: Conceptualization, resources, investigation, methodology, writing–original draft, writing–review and editing. J. Rose: Conceptualization, resources, funding acquisition, investigation, methodology, writing–original draft, project administration, writing–review and editing. N.A. Berger: Conceptualization, resources, writing–original draft, writing–review and editing. S.M. Koroukian: Conceptualization, resources, supervision, funding acquisition, investigation, methodology, writing–original draft, project administration, writing–review and editing.

The authors thank Will Penman of Composition Coaching for his advice and assistance in the earlier planning of this article. This study was funded by a grant from the NCI, Case Comprehensive Cancer Center (P30 CA043703) to S.M. Koroukian, J. Rose, and W. Dong. W. Dong and S.M. Koroukian are supported by a grant from the American Cancer Society (132678-RSGI-19-213-01-CPHPS) and by contracts from Cleveland Clinic Foundation, including a subcontract from Celgene Corporation. S.M. Koroukian is also supported by grants from the NIH (R15 NR017792, and UH3-DE025487) and the American Cancer Society (RWIA-20-111-02 RWIA). W.P. Bensken is supported by a grant from the National Institute on Minority Health and Health Disparities (F31MD015681). U. Kim is supported by grants from the National Institute of General Medical Sciences (5T32GM007250), National Center for Advancing Translational Sciences (5TL1TR002549), and the PhRMA Foundation (PDHO18). J. Rose is supported by grants from the National Institute of Dental and Craniofacial Research (1UH2DE025487-01), the National Heart Lung and Blood Institute (R01 HL153175), and the American Cancer Society (RWIA-20-111-02 RWIA). N.A. Berger is supported by grants from the NCI (2P50 CA150964, 2U54CA163060, P20CA233216, R25CA221718, and R25CA225461).

This study was funded by a grant from the NCI, Case Comprehensive Cancer Center (P30 CA043703, to S.M. Koroukian, J. Rose, and W. Dong).

The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.