Uncovering spatial variation in maternal healthcare service use at subnational level in Jimma Zone, Ethiopia

Background: Analysis of disaggregated national data suggest uneven access to essential maternal healthcare services within countries. This is of concern as it hinders equitable progress in health outcomes. Mounting an effective response requires identification of subnational areas that may be lagging behind. This paper aims to explore spatial variation in maternal healthcare service use at health centre catchment, village and household levels. Spatial correlations of service use with household wealth and women’s education levels were also assessed. Methods: Using survey data from 3758 households enrolled in a cluster randomized trial geographical variation in the use of maternity waiting homes (MWH), antenatal care (ANC), delivery care and postnatal care (PNC) was investigated in three districts in Jimma Zone. Correlations of service use with education and wealth levels were also explored among 24 health centre catchment areas using choropleth maps. Global spatial autocorrelation was assessed using Moran’s I. Cluster analyses were performed at village and household levels using Getis Ord Gi*and Kulldorf spatial scan statistics to identify cluster locations. Results: Significant global spatial autocorrelation was present in ANC use (Moran’s I = 0.15, p value = 0.025), delivery care (Moran’s I = 0.17, p value = 0.01) and PNC use (Moran’s I = 0.31, p value < 0.01), but not MWH use (Moran's I = -0.005, p value = 0.94) suggesting clustering of villages with similarly high (hot spots) and/or low (cold spots) service use. Hot spots were detected in health centre catchments in Gomma district while Kersa district had cold spots. High poverty or low education catchments generally had low levels of service use, but there were exceptions. At village level, hot and cold spots were detected for ANC, delivery care and PNC use. Household-level analyses revealed a primary cluster of elevated MWH-use not detected previously. Further investigation of spatial heterogeneity is warranted. Conclusions: Sub-national variation in maternal healthcare services exists in Jimma Zone. There was relatively higher poverty and lower education in areas where service use cold spots were identified. Re-directing resources to vulnerable sub-groups and locations lagging behind will be necessary to ensure equitable progress in maternal health. Ethiopia is composed of nine regional states and two city administrations; these are sub-divided into woredas (districts) which comprise several kebeles (villages). This analysis focuses on three rural woredas (Gomma, Seka Chekorsa and Kersa) in Jimma Zone, located in south-western Ethiopia in Oromiya region as shown in Fig. 1. Populations in the woredas range from 180,000 to 270,000 in 2016 [16]. Coffee production is an important source of revenue for residents of Gomma while Seka Chekorsa and Kersa residents engage mainly in small-scale cereal production [16]. As part of the tiered health care system in Ethiopia, Jimma Zone has two general hospitals, six district hospitals, 122 health centres and 566 health posts [17]. The lowest health system level functions at the woreda level and consists of primary health care units (PHCUs). Each PHCU includes a health centre serving around 25,000 people and community-based health posts responsible for a population of 3000-5000. Regional level data indicate that in 2016, 48.6% of women surveyed in Oromiya Region in the DHS received no antenatal care and just 18.8% of women delivered at a health facility, placing the region among the poorer performing areas in the country [3]. Study area map showing PHCU boundaries and locations of health centres within PHCUs. The figure was generated in ArcMap 10.6.1 using study GPS data shapefiles obtained from the Jimma Zone Health Office Data from a cross-sectional household survey, conducted prior to intervention roll-out (baseline survey) between October 2016 and January 2017 as part of a cluster-randomized controlled trial (ClinicalTrials.gov; {"type":"clinical-trial","attrs":{"text":"NCT03299491","term_id":"NCT03299491"}}NCT03299491), was used for this analysis [18]. Women living in the catchment area of 24 PHCUs within the study’s three districts and who reported a pregnancy outcome (livebirth, stillbirth, miscarriage or spontaneous abortion) within 12 months of the baseline survey were eligible to participate. Sample size (24 PHCU clusters with 160 women each) was calculated based on 80% power to detect an absolute difference of 0.17 in the proportion of institutional births (primary trial outcome) assuming a control arm proportion of 0.4 and using a two-sided alpha of 0.025 to account for two pairwise comparisons. Verbal informed consent was sought from all participants owing to low levels of literacy. Approximately 98.5% (n = 3784) of the targeted women were successfully enrolled and interviewed in either Afaan Oromo or Amharic by trained research assistants. Data were collected on sociodemographic characteristics and maternal healthcare service utilization using tablet computers programmed with surveys using Open Data Kit software. GPS locations of households and health centres were also collected on tablet computers; GPS locations were available for 3758 (97.9%) of the enrolled households. The analysis focused on four services across the continuum of maternal healthcare: antenatal care, maternity waiting home use, delivery care at health facilities, and postnatal care. As women’s education levels and household wealth have previously been linked to access inequities [19], we adopted them for our analyses. Operational definitions of all the variables of interest are described in Table 1. Operational definitions of analysis variables used to describe maternal healthcare service use among study women Principal components analysis was used to create the household wealth variable using ownership of assets and animals, utilization of health insurance, presence of electricity supply, type of drinking water source, type of toilet present, and type of materials used for floor construction of the home. Using methods described by Vyas & Kumaranayake [20], socio-economic scores were generated for each household and then categorized into quintiles. The first quintile corresponds to poorest households while the fifth quintile corresponds to the least poor households. Frequencies and percentages of education and wealth levels in survey clusters were generated in STATA version 15 as part of descriptive analyses. Differences between survey clusters in the percentages of wealthy and educated households were then compared using a chi square test. P values less than 0.05 considered to be statistically significant. Euclidean distances between households and the health centre within the PHCU catchment were calculated in kilometres for each PHCU using the Generate Near Table tool in ArcMap version 10.6.1. Distances were summarized using means, standard deviations and ranges and classified into three categories (< 2 km, 2-5 km and > 5 km) to facilitate comparison between PHCUs using a chi square test. Geographical variation in maternal healthcare service use was explored at the PHCU-level, the level at which maternal and other primary care services are coordinated and at the kebele-level, which is the smallest administrative unit and where health posts, staffed by two health extension workers (HEWs) are located. To this end, the proportion of women who reported using a given service during their last pregnancy within each geographical unit (PHCU or kebele) in the study area was calculated. Additionally, household-level binary data was used to explore variation in service use without constraining analyses within administrative boundaries. Instead, a maximum radius of 5 km, at which we hypothesize interpersonal and social factors affecting service use could operate, was used. Other studies in similar settings have estimated a one-hour walking distance to be between 3 km to 5 km depending on the season and terrain [21, 22]. District-level analyses were not feasible given the limited number of districts. At the PHCU-level, choropleth maps were generated; choropleth maps present interval data superimposed over geographic units using colours and symbols to distinguish between intervals [23]. To explore utilization levels among marginalized population sub-groups (women from poor households and women with no education), education and household wealth were symbolized as separate layers and superimposed to visualize correlations. Intervals for all variables were manually created based on what best suited the distribution of each variable, meaning that only qualitative comparisons can be made between services as intervals differ between services. Each maternal health care service was assigned a different colour with darker shades indicating higher extents of service. Wealth and education indicators are overlayed using proportional symbols where larger circles correspond to higher levels. Administrative boundaries for PHCUs were created by dissolving boundaries of villages known to fall within the catchment area of the PHCU using the Dissolve tool in ArcMap. PHCU size varied from 53 km2 (sq km) to 186 sq. km; the mean area was 108.5 sq. km (standard deviation 39.6 sq. km). To examine variation in service use at the kebele level, the presence of spatial autocorrelation for each service was first examined using the global Moran’s I statistic. The results of this test indicate whether or not the spatial patterns observed in the data are random by looking at how each kebele deviates from the mean among neighbouring kebeles. Statistically significant and positive Moran’s I indices indicate the presence of clustering (i.e a high degree of similarity in levels of service use between neighbouring kebeles) while negative values specify dispersion (i.e neighbouring kebeles have dissimilar values) [24]. An inverse-distance-based spatial relationship conceptualization was selected for the process; this means that kebeles outside the threshold distance are not included in the computations. The threshold distance at which every kebele had at least one neighbouring kebele was determined to be 7.6 km. Of the 100 kebeles present in the study area, data were available for 96. An average of 48 women were enrolled per kebele. One kebele with less than five women enrolled was excluded as this number was considered to be too low for the analysis. The next step was to pinpoint the location of clusters. The Optimized Hot Spot Analysis tool in ArcMap was used to identify where clustering was occurring. This tool relies on the Getis Ord Gi* spatial statistic to uncover statistically significant hot spots (clusters of high service use) and cold spots (clusters of low service use). Clusters are considered statistically significant only if they are surrounded by similarly high or low values [25]. The optimal fixed distance band is determined by the tool for each service outcome by determining at what distance clustering is maximized. This was found to be 10 km for ANC, 7 km for MWH use and 19 km for delivery care. For PNC, a threshold distance of 7 km was used as no optimal distance was found using the clustering intensity method. The resulting output is a map of statistically significant hot (red) and cold (blue) spots presented at 99% (+/− 3 Gi* values), 95% (+/− 2 Gi* values) and 90% significance (+/− 1 Gi* values) and corrected for multiple testing. Finally, household-level data on women’s reported service use were analysed using the Kulldorf spatial scan statistic in SaTScan 9.6. This spatial tool functions by running circular “scanning windows” of multiple sizes across the study space. The events observed within the window are compared to those expected under the null hypothesis of no difference inside and outside the window; a relative risk (RR) is generated which represents the ratio of events observed and expected within the window to those observed and expected in the study area. The formula used in SaTScan [26] is: where c is the observed number of events within a potential cluster (window). E[c] is the expected number of events within a potential cluster (window). C is the total number of events in the dataset. The overall proportion of events in the study area were 84.3% for ANC, 6.6% for MWH use, 48.5% for delivery care and 39.0% for PNC and are used to calculate the expected number of events within a cluster (E[c]). A relative risk greater than one would suggest proportion observed within the window is higher than expected while relative risks less than one suggest observed proportions are lower than expected. A Bernoulli model was used due to the binary nature of the outcomes and the window radius set to a maximum of 5 km. The most likely cluster is identified as the window with the maximum likelihood, meaning it was least likely to have occurred by chance. Monte Carlo hypothesis testing is used to generate p-values where ranks of maximum likelihood of the observed dataset is compared to random datasets. Secondary clusters are also identified and ranked according to the Likelihood Ratio Test statistic [26].