Abstract
Objectives
Evaluating population health initiatives at the community level necessitates valid counterfactual communities, which includes having similar population composition, health care access, and health determinants. Estimating appropriate county counterfactuals is challenging in states with large intercounty variation. We describe an application of K-means cluster analysis for determining county-level counterfactuals in an evaluation of an intervention, a county perinatal system of care for Medicaid-insured pregnant women.
Methods
We described counties by using indicators from the American Community Survey, Area Health Resources Files, University of Wisconsin Population Health Institute County Health Rankings, and vital records for Michigan Medicaid-insured births for 2009, the year the intervention began (or the closest available year). We ran analyses of 1000 iterations with random starting cluster values for each of a range of number of clusters from 3 to 10 with commonly used variability and reliability measures to identify the optimal number of clusters.
Results
The use of unstandardized features resulted in the grouping of 1 county with the intervention county in all solutions for all iterations and the frequent grouping of 2 additional counties with the intervention county. Standardized features led to no solution, and other distance measures gave mixed results. However, no county was ideal for all subpopulation analyses.
Practice Implications
Although the K-means method was successful at identifying comparison counties, differences between the intervention county and comparison counties remained. This limitation may be specific to the intervention county and the constraints of a within-state study. This method could be more useful when applied to other counties in and outside Michigan.
Complex, multifactorial public health problems necessitate population-wide solutions and, thus, community-level interventions are frequently used. The gold standard for evaluating community-level interventions, such as those at the county level, is the group-randomized trial. 1 However, group-randomized trials are often not feasible or appropriate in a small number of communities, particularly in community-based participatory research and other scenarios in which the intervention is initiated by community members or in system-of-care interventions that depend on coordinating local resources. It is not clear how best to evaluate non–randomly assigned community interventions under these circumstances. 2 -4 As in individual-level quasi-experimental analyses, the challenge is to identify an appropriate counterfactual group or communities that are as similar as possible to the intervention community to approximate group randomization, 3,5 but little guidance exists on selection criteria or methods.
We conducted a demonstration project during 2009-2015 to determine whether enhancements to a county perinatal system of care improved population outcomes for Medicaid-insured pregnant women and their infants in Kent County, Michigan, a mixed urban/rural county containing the second-largest city in the state, Grand Rapids. Because our intervention targeted the entire county, the ideal counterfactual community or communities would have the same sociodemographic characteristics and level of complexity in service delivery as the intervention county. Complexity is influenced by many factors, including the numbers of residents, health care agencies, health care providers, and enhanced prenatal care programs.
Because we wanted to account for the broader context and administrative features of an entire county, we explored a clustering method, K-means cluster analysis, to identify 1 or more comparison counties in Michigan. K-means cluster analysis has been used increasingly in the public health and health services research literature to create clusters of exposures or individuals, most frequently for sources of air pollution 6 and for segmentation of health care populations. 7 Where the method was used for county-level clustering, it grouped all counties nationwide into clusters for descriptive health rankings rather than for intervention evaluation. 8,9
The objective of this study was to describe and discuss an application of K-means cluster analysis for determining county-level counterfactuals for evaluation of a county perinatal system of care for Medicaid-insured pregnant women. To our knowledge, this is the first use of this method to find the most similar counterfactual county.
Materials and Methods
Michigan has 83 counties. We selected county-level features from 3 publicly available data sources: the American Community Survey, 10 the Area Health Resources Files, 11 and the University of Wisconsin Population Health Institute County Health Rankings. 12 The fourth data source was a limited dataset of Michigan resident live birth records for Medicaid-insured births only (hereinafter called Michigan Medicaid vital records [MMVR]), including county identifiers. We retrieved the dataset from the Michigan Department of Health and Human Services (MDHHS) Health Services Data Warehouse through an honest broker (a neutral third party, who is not part of the research team in any way).
We considered the following categories of features: population composition, social and economic factors, health care access, health outcomes, health behaviors, and physical environment. We selected 35 features across these broad categories, some chosen from Rettenmaier and Wang 13 and others based on the evaluation team’s experience with the Michigan Medicaid Maternal and Infant Health Program and Strong Beginnings, a federal Healthy Start program. 14 -16 These features represented a mixture of compositional (ie, aggregated from individual-level characteristics) and contextual (ie, pertaining to the county environment) measures. 17 Whenever possible, we used data from 2009 because the enhanced county perinatal system of care was implemented in that year. For data elements without information from 2009, we used the following years as proxies, in priority order: 2010, 2008, or 2011. We did not have data on 4 features in some counties (number of counties ranged from 1435); for these counties, we imputed mean values. We described the 35 features with the data source, year, mean, and SD (Table 1).
aUS Census Bureau. 10
bHealth Resources and Services Administration, US Department of Health and Human Services. 11
cMMVR is a limited dataset of live birth records among Michigan residents for Medicaid-insured births only, retrieved from the Michigan Department of Health and Human Services Health Services Data Warehouse.
dUniversity of Wisconsin Population Health Institute. 12
eEighteen counties missing data for this variable; values replaced with mean values.
fData for 2009 and 2010 not available.
gThirty-five counties missing data for this variable; values replaced with mean values.
hEighteen counties missing data for this variable; values replaced with mean values.
iFourteen counties missing data for this variable; values replaced with mean values.
Analysis
We used K-means cluster analyses to identify the comparison county or counties that would be classified most frequently in the same cluster as the intervention county. Cluster analysis, also called data segmentation, groups or segments a collection of objects (counties in our study) into subgroups or clusters such that objects (counties) in the same cluster are more similar to each other than objects (counties) assigned to other clusters. 18 In a cluster analysis of counties, counties are described by a set of characteristics (features, measurements) that are preselected according to the purpose of an evaluation; in our study, we selected characteristics relevant to the evaluation of a pregnancy-centered health program. Counties in the same cluster may have different degrees of similarity or dissimilarity with each other and to a specific county (Kent County in our study). Among the most widely used clustering methods is the K-means (or median) clustering, for which the within-cluster variation is as small as possible. 6,7
We tested solutions for 3-10 clusters for 1000 iterations each with randomly assigned starting cluster values. To identify the best solution, we used scree plots to visualize a kink in the curve generated from the within sum of squares (WSS) or its logarithm (log[WSS]) for all cluster solutions, that is, the point at which reduction in WSS or log(WSS) is not appreciable enough to warrant increasing the number of clusters.
19
The WSS can be thought of as the error sum of squares in regression analysis, where the total sum of squares is defined when all objects are in 1 cluster. We also used another criterion, the η2 coefficient:
Nineteen of the 35 features are proportions that have a natural range (0-100), and the other features are continuous measures that may have wide SDs. If we standardized these features, some natural variation (eg, the size of the population in each county) could be masked. In K-means analyses, if all variables are standardized, then clustering based on correlation (similarity) is equivalent to clustering based on squared distance (dissimilarity). 21 Therefore, we ran all analyses twice, with and without a z score standardizing all variables. We then conducted several robustness checks on both unstandardized and standardized analyses, including using the 5 standardization methods in Schaffer and Green 22 and the mean absolute deviation standardization to calculate the Euclidean distance, L1 distance, Canberra distance, and 1 – correlation distance of the 35 features between Kent County and other counties. We used Stata version 15 (StataCorp LLC) for all analyses. The evaluation project was reviewed by the Michigan State University, MDHHS, and Spectrum Health institutional review boards and determined not to be human subjects research.
Results
By visual inspection of the scree plot, when features were not standardized, the optimal number of clusters was 4, where WSS could be seen at a kink point and the reduction in PRE was the largest (Figure 1). On the other hand, when features were standardized, we found no clear optimal solution of the number of clusters (Figure 2). The gap statistic suggested a 5-cluster solution when features were not standardized and no solution when features were standardized (details available upon request).
Scree plot of within sum of squares (WSS), log(WSS), η 2 coefficient, and proportional reduction of error (PRE) to demonstrate the optimal number of clusters for 10 cluster solutions when county-level features for characterizing Michigan counties are not standardized, Michigan, 2008-2011. Visual inspection of the plot demonstrates that the optimal number of clusters is 4.
Scree plot of within sum of squares (WSS), log(WSS), η 2 coefficient, and proportional reduction of error (PRE) to demonstrate the optimal number of clusters for 10 cluster solutions when county-level features used for characterizing Michigan counties are standardized, Michigan, 2008-2011. Visual inspection of the plot demonstrates no clear optimal solution of the number of clusters.
Because our goal was to find another county or other counties most suitable as a comparison with Kent County, we viewed higher frequency of a county in the same cluster as Kent County to be an indication of similarity to Kent County. The use of unstandardized features grouped 3 counties with Kent County in 1000 iterations when the number of clusters varied from 3 to 6, the range around the optimal number of 4 or 5. Macomb County clustered with Kent County in all iterations of all solutions, Oakland County clustered with Kent County in 979-1000 iterations, and Genesee County clustered with Kent County in 127-1000 iterations, depending on the number of clusters. In contrast, when we standardized all features by z score, no county was grouped with Kent County in a majority of iterations (details on standardized features and robustness checks available upon request).
Given the lack of an optimal solution under z score standardization and that the counties grouped together were not similar, we turned to different distance and standardization methods (Table 2). In terms of Euclidean distance (the distance measure used in K-means clustering), unstandardized features led to the same 3 counties, whereas Ingham County had the smallest distance to Kent County when we used 4 of 6 standardized measures. For the most part, the use of L1 and Canberra distance resulted in the same counties that were obtained by using Euclidean distance, while 1 – correlation distance produced different results.
Top 3 counties as nearest neighbors of Kent County based on various distance and standardization methods, Michigan, 2008-2011
a x represents the value for each observation for a given feature.
bStandardized by dividing each value x by the maximum value x for that feature.
cStandardized by dividing each value x by the difference between the minimum value x and the maximum value x for that feature.
dStandardized by dividing each value x by the sum of the values of x for that feature.
eStandardized using the rank of each value x for that feature.
fStandardized by taking the absolute value of the difference between each value x and the mean value x for that feature.
The 3 most similar counties to Kent County, using unstandardized features, were Macomb, Oakland, and Genesee counties, and the 3 most similar counties, using standardized features, were Ingham, Kalamazoo, and Muskegon counties (Table 3) under most of the difference measures. In examining values individually, Macomb and Kalamazoo counties were within 2 percentage points of Kent County for 10 of 19 percentage-based features, with nearly identical median annual household income and numbers of Medicaid-eligible women and Medicaid-eligible people aged <21. Ingham and Muskegon counties were within 2 percentage points of Kent County for 7 of 19 percentage-based features, and Oakland and Genesee counties were within 2 percentage points of Kent County for 5 of 19 percentage-based features; all 7 counties had similar median annual household incomes. We found broad variation in access-to-care features across counties. However, we found 2 related features for which all comparison counties had smaller values than Kent County had, most of which were considerably smaller: the 2 Latinx population features, “percentage of Latina females aged 15-54” (Kent County, 6.7%; comparison counties, 1.5%-5.0%) and “percentage of Latinx people aged ≥25 with <high school diploma” (Kent County, 46.3%; comparison counties, 16.9%-35.7%).
aUS Census Bureau. 10
bHealth Resources and Services Administration, US Department of Health and Human Services. 11
cMMVR is a limited dataset of live birth records among Michigan residents for Medicaid-insured births only, retrieved from the Michigan Department of Health and Human Services Health Services Data Warehouse.
dUniversity of Wisconsin Population Health Institute. 12
Discussion
K-means clustering achieved identification of 1-3 comparison counties when preselected features were not standardized in our study, whereas other distance measures identified comparison counties when features were standardized across methods. Overall, the evidence suggested that the selection of counties based on a single distance measure is sensitive to the choice of standardization and the choice of distance measure.
However, the features that differed between Kent County and similar counties were concerning. Compared with the best-matched county by the K-means analysis, Kent County had a greater percentage of Latina women and a greater percentage of people with <high school education. We were also cautious about interpreting differences in county maternal and infant health programming and infrastructure that were not reflected in the data sources. For example, one of the selected counties includes the city of Flint, which received an influx of maternal and child health resources during the study period because of the city’s water crisis. Ultimately, we opted against using the method described here and instead selected individual women from across the state (outside Kent County), using propensity score-weighted analyses, matched exactly on race, and using census tract–level and block group–level data, which are more granular than county-level data. 23 This approach was an option for us given our intervention design and available data; however, it would not be appropriate for other evaluation designs involving primary data collection and/or more limited resources.
Our results are consistent with the descriptive rankings of Wallace et al, 8 who used K-means clustering and national data to create clusters of counties with similar sociodemographic features. They noted considerable within-state heterogeneity, such that most states (including Michigan) contained counties from many different clusters. Of note, the 6 comparison counties that we identified were in the same cluster as Kent County in Wallace’s schema as well. However, they also indicated that the counties that were sociodemographically heterogeneous would benefit from more granular subcounty data, which we found applicable to our diverse intervention county.
Our findings are also congruent with recent discussion of how machine learning algorithms can reproduce biases against racial/ethnic minority groups. 24 -26 Although previous literature has focused on racial/ethnic disparities in predictive algorithms, our analysis demonstrates how a racial/ethnic subgroup can be obscured by a clustering method. Thus, the same cautions that should be extended to ensure equity in predictive modeling should be applied to cluster analyses as well.
Limitations
General limitations of the K-means method to select counterfactual communities include the following. First, as in any selection method, the results are only as good as the chosen features. In our example, we selected features through discussion with the research team and community partners with substantive knowledge of factors influencing perinatal health in the community. We captured data on a wide range of complex indicators, including health care determinants and social determinants of maternal and infant health outcomes, triangulated across 4 data sources. That said, our selection was likely not the optimal set, and it is probable that no single optimal set exists. A more data-driven approach than the one we took could be applied to selecting features. 8 Relatedly, given our interest in the Latinx population, adding more measures specific to this population would weight their importance more heavily than they were weighted in our study. However, this approach would not address the problem that the counties whose other features are most similar to Kent County have much smaller Latinx populations and would have led to no better solution.
Second, no matter how practically important a certain measure might be (in our case, important to the evaluation team and community partners), a variable could be unrelated to the classification of counties, if a natural grouping exists. Inclusion of such variables, typically called “masking variables,” 27 can substantially affect the accuracy of the solution. Two broad approaches for mitigating the effect of masking variables are variable weighting and variable selection. 28 However, the literature is vast, and no consensus exists when the true cluster is unknown in observational data. 29,30 Collaborations between practitioners and methodologists become especially important in complex real-world intervention evaluation.
Third, our method involved the assumption that if a counterfactual community is comparable at baseline, it will remain so throughout the study. This assumption is shared by any method for identifying a comparison at baseline, 3 but it may not hold because of differential changes in the communities, unrelated to the intervention, that can confound the effect of the intervention. The assumption was particularly salient in our study given the Flint water crisis and the resulting resource allocation. Careful selection of an analytic method can help to mitigate biases introduced by extrinsic temporal trends. 31 Moreover, in our data sources, not all features were measured during the optimal baseline year, necessitating the use of proxy years. However, the intervention evaluated was complex and involved many changes during the study period, and even the postbaseline proxy years preceded the years in which changes were implemented.
Finally, the use of K-means clustering in general comes with its own set of assumptions and limitations. 21 One issue of the K-means clustering method is that the resulting assignments depend on the random starting point. The K-means algorithm gives local minima and does not guarantee to give the global minimum, so the starting points should be varied to examine the end partitioning. Another issue of the K-means algorithm is that a variable with high variability can dominate the cluster analysis. This problem was evident in our data as the features “number of Medicaid-eligible females” and “number of Medicaid-eligible people aged <21 years” drove the solutions for our unstandardized analyses. A common solution is to standardize variables, but standardizing could also hide the true groupings in the data. Steinley 29 pointed out several diverging viewpoints on standardization, including critiques of z score standardization, that justify our use of different standardization methods. On the other hand, in our study, because the true cluster structure was unknown, not standardizing features may have given undue weights to variables with high variability, with the potential consequence of selecting a comparison county similar in these statistically dominating features and missing other possibly more important confounding variables. For example, when ranking the features by the standardized difference between Kent and Macomb counties, we found that the feature “number of Medicaid-eligible people aged <21 years” had the third smallest standardized difference, whereas the same feature ranked fifth and ninth in the other counties. If this variable was in fact not an important factor (which cannot be ascertained unequivocally by researchers), then perhaps another county should have been selected. This decision is a case-by-case decision and depends on the type of data and the nature of the groups. The number of clusters K is the tuning parameter that can be chosen by cross-validation if the number of independent and identically distributed observations is large, which was not possible in our example.
Given these limitations, one might question whether K-means clustering is an improvement over more subjective heuristic selection methods. We argue that it is an improvement. Some limitations of K-means analysis are shared with other selection methods (eg, selection is only as good as the attributes selected for; external temporal trends may confound the effect). In addition, any county-level selection method could mask subpopulations of interest when counties are heterogeneous. However, compared with a subjective heuristic method, the K-means clustering method allows selection to be based on a greater number of features, more efficiently and with greater transparency. That the choice of a comparison has ramifications for the apparent effectiveness of an intervention is certain; K-means clustering is a data-driven option for selecting an appropriate comparison.
Practice Implications
K-means clustering provides a rigorous, data-driven approach to selecting 1 or more counterfactual groups for a community-level intervention, representing an improvement over subjective selection of comparison communities. Although it provided a consistent comparison county solution, it did not result in an optimal comparison for the intervention county in our study. The method could be more useful at other geographic levels or for other counties in or outside Michigan. Nonrandomly assigned community interventions are common in population health; our study can stimulate conversation about how best to select appropriate comparisons. In particular, our study demonstrates that selection of counterfactual communities should be objective, transparent, and examined critically through a health equity lens.
Footnotes
Acknowledgments
The authors thank the Michigan Department of Health and Human Services for providing access to its Health Services Data Warehouse and consultation from the Maternal and Infant Health Division and the Division for Vital Records and Health Statistics.
Declaration of Conflicting Interests
The authors declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.
Funding
The authors disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This project was supported by the Agency for Healthcare Research and Quality (grant R18-HS020208) and the Spectrum Health–Michigan State University Alliance Corporation. The content is solely the responsibility of the authors and does not necessarily represent the official views of either funding organization.
