Vaccination site locations as of May 23rd 2021 were obtained through the NJ Department of Health's COVID-19 information hub (24) and geocoded with the ArcGIS World Geocoding Service (25). Daily available appointments at each site were used to represent the vaccination capacity of the individual sites. Socioeconomic data at the census tract level was compiled as Geographic Information System (GIS) maps and associated attribute tables containing variables from the US Census 2019 American Community Survey (ACS), including age and sex, race and ethnicity, income and poverty, housing characteristics, technology and internet availability, among others.

We used the two-step floating catchment area (2SFCA) method (26) to calculate accessibility to vaccine resources for each of the 318 census tracts in the two counties. The method computes a supply-to-population ratio to measure accessibility to facilities or resources such as healthcare facilities (27), food resources (28), and most recently COVID-19 testing sites (19) and hospital beds (21). To measure accessibility to vaccines, the first step is to assess the availability of supply at each vaccination site as the ratio of supply to the demand population located within a catchment area of each site (j). Delineation of catchment areas is based on a threshold travel distance (d₀) along road networks. A supply-to-demand ratio R_j is then computed for each site using Equation 1.

S_j is the available vaccination appointments of site j, P_k is the population of census tract k, whose centroid falls within the catchment of site j. The travel distance d_kj from k to j is no greater than a preset threshold driving distance d₀. In our study, we tested different thresholds in the delineation of catchments from 1–5 miles. It was observed that with 5-mile catchments, almost all areas in the two counties are covered. As vaccines can be readily administered at many pharmacies, accessibility to vaccination sites should be evaluated differently from accessibility to primary care doctors, for which the commonly used threshold is a 30 min travel time (19). With nearly 9 in 10 Americans living within 5 miles of a community pharmacy (29), we adopted the 5-mile threshold for calculating accessibility scores. One other note in the calculation of accessibility in our study is that the county mass vaccination sites are only open to the county's own residents. Thus we adjusted catchment delineation for county-operated mass sites using county boundaries.

One census tract may fall within the catchments of multiple vaccination sites. The second step of 2SFCA is to calculate the accessibility score A_i for each census tract i by summing up the supply of all nearby vaccination sites whose catchment areas contain the census tract i using Equation 2, where A_i is the accessibility score calculated for census tract i.

The accessibility scores of the census tracts were mapped with GIS. Cluster and outlier analysis was conducted using the Anselin Local Moran's I (30) to measure the concentration of high and low accessibilities in the study area with ArcGIS Pro (31).

In order to examine how accessibility varies among different populations, we selected a set of variables from the census attribute table. These include the percentage of population under poverty, percentage of 65 years and older, and percentage of different racial/ethnic groups. We calculated the average accessibility scores of each population group, following Lu et al. (19), using Equation 3. A^g is the average accessibility score of population group g, Pig is the percentage of population of group g in census tract i, and n is the total number of census tracts.

In order to test the significance of the differences in accessibility across population groups, we conducted spatial lag regressions with accessibility scores as the dependent variable and the socioeconomic variables as independent variables. A spatial lag of accessibility was added to the linear regression model as an independent variable to account for spatial autocorrelation between neighboring tracts because spatial accessibility measures are usually strongly spatially auto-correlated (19, 26). Specifically, let y_i be the accessibility score of a census tract i and the vector of covariates for tract i is x_i. The model is expressed as:

ε_i is a random error term. w_iy is the spatial lag, a weighted average of the spatial neighbors of census tract i, defined by a spatial contiguity matrix W. ρ represents the relationship between accessibility at a location with accessibilities of its neighbors and βis the vector of local regression coefficients associated with x_i. We adopted the first order Queen contiguity to define the spatial contiguity matrix (19, 32) and conducted the analyses using GeoDa 1.18 (33).

In order to examine how a certain socioeconomic factor is correlated with accessibility and how the correlation varies geographically in the two counties, we conducted local bivariate relationship analysis using local entropy maps (34). It allows for the quantification of spatial heterogeneity of the correlation between two variables. It uses a local entropy statistic to measure the amount of shared information between the two involved variables. Entropy can capture complex relationships including exponential, quadratic, and not just linear relationships like other statistics. The results will help us answer specific questions such as: if poverty has a significant correlation with accessibility to vaccination resources, how does the strength of the relationship vary across the different neighborhoods in the two NJ counties? The spatial variation of relationship types and strengths could give insights to resource disparity and guide future public health planning and responses.

Figure 1 shows the distribution of accessibility scores with different types of vaccination sites: medical/pharmacy sites (Figure 1A), county mass vaccination sites (Figure 1B), and FEMA site (Figure 1C). There are notable variations of accessibility depending on the type of facilities. The county sites map visualizes a large cluster in Union county with high accessibility scores compared to the lower accessibilities with only medical/pharmacy sites. In Essex county, high accessibility clusters appear on the northwest side for both county mass sites and medical/pharmacy sites but not the southeast. It is the FEMA community site established in Newark that improved accessibility to the southeast portion of the county. By May 2021, all three types of vaccination sites were in operation. Spatial heterogeneity and local clustering are still notable on the accessibility map combining all vaccine resources (Figure 1D).

The significance of spatial clusters is also confirmed by Anselin Local Moran's I measures. Figure 1E shows significant local clusters include high value clusters (HH) and low value clusters (LL). Significant outliers are neighborhoods with high accessibility surrounded by those with low values (HL), and vice versa (LH). Two small high accessibility score clusters are identifiable in part of Essex (Newark) and Union (Union Township). Two large clusters of low accessibility are located near Montclair in Essex and Springfield in Union. Contrasting accessibility in neighboring communities can be found in many areas such as part of Belleville with high accessibility tracts surrounding low accessibility neighborhoods. These results indicate significant spatial heterogeneity of accessibility to vaccine resources.

Table 1 lists the computed average accessibilities of different socioeconomic population groups. In Essex County, Black and Native Americans' accessibilities are lower than average with county mass sites only. So is the accessibility of the population under poverty. Adding medical/pharmacy sites increased overall accessibilities just slightly and did not change the pattern. It was the FEMA site that increased accessibilities greatly and especially for Black, Native, Hispanic populations and the population under poverty. Union county has overall higher accessibilities than Essex. Black and Native Americans' accessibilities are always higher than average starting with county mass sites. So is the accessibility of the population under poverty. Asian and Hispanic populations' accessibilities are slightly lower than average. The FEMA site did not change the situation much as most Union residents live more than 5 miles away from Newark. With all vaccine resources combined, most minority population groups and populations under poverty have higher than average accessibilities in the two counties. Elderly population has slightly lower than average accessibility in both counties. These results indicate that the low accessibility clusters in Union County are not necessarily associated with unfavorable social vulnerability in the aspect of minority population and only corresponds slightly to elderly population.

We fitted spatial lag regression models with two sets of variables to further explore the relationships between the socioeconomic variables and vaccine accessibility. The first set includes the percentage of the total minority population, percentage of population under poverty, percentage of people aged 65 years or older, percentage of the population with no computer and smartphone access, and population density. We conducted spatial lag regression with a second set of variables to break down the percentage of minority population into individual racial/ethnic groups. We performed analyses for accessibility in multiple scenarios: with all vaccination sites combined, with medical sites only, with county mass sites only, and with county mass sites plus medical sites combined. Complete results of this analysis can be found in the Supplementary Material.

The pseudo R² values for most spatial lag models range between 0.7 and 0.8 in the various scenarios for the two counties, indicating a moderately strong predictive power. The most significant explanatory variable, however, is the spatial lag variable. This is expected as spatial accessibility measures are often strongly spatially autocorrelated, since neighboring tracts are generally located in close proximity to the same facilities (19). The results show that for Essex County with all vaccination sites combined, most variability of accessibility are explained by such spatial auto-correlation and not most of the other explanatory variables except for the percentage of minority. A positive coefficient, though, suggests a positive correlation, meaning that the higher the minority population, the higher accessibility. This positive relationship is not present with only medical sites, county mass sites, and when the two are combined. Results with the minority population broken down to individual groups show that the percentages of Black and Hispanic populations have significant correlations with accessibility. The correlations, once again, are both positive, indicating the higher these minority population percentages, the higher accessibility to vaccine resources.

In Union County, there is a significant correlation between minority and accessibility and it is also a positive relationship. This relationship starts to be present when we add county mass sites to medical sites. Breaking down to individual minority groups, the significant correlations are present with both Black and Hispanic populations. And similar to Essex county, the correlations are positive, indicating that these minority groups have higher accessibilities.

Based on our analyses, a few conclusions can be drawn. First, the COVID-19 pandemic exposed health and socioeconomic inequities in our communities that need to be addressed at the various stages from testing to vaccination. Public health decision making at different levels (county, state, and federal) are all critical in mitigating the disparities. In the case of the two diversely populated New Jersey counties, there was clear spatial heterogeneity in the accessibility to vaccine resources with the existing infrastructure (medical/pharmacy sites). Accessibility patterns changed with the introduction of county mass sites and the FEMA community site. In particular, county mass sites were set up in Union County targeting densely populated areas and minority populations, so was the FEMA community site in Essex County. These resulted in notable changes in the accessibility landscape. It shows how both the local health department and the federal government play an important role in mitigating pre-existing inequalities.

Second, social determinants of health need to be addressed (13) and taken into explicit consideration when distributing resources and evaluating resource accessibility. In our study, the county mass sites in Union County greatly increased accessibilities for populations of minority and poverty. The county mass sites in Essex County, on the other hand, were set up evenly across space, and did not show the same effect. This suggests that the setup of mass vaccination sites should not just cover the geographic space, but need to address the attribute space defined by variables including population density and socioeconomic factors, especially when some socioeconomic groups have been found to be disproportionately affected by the pandemic. Spatial distributions of such socioeconomic variables should direct the vaccine site selections.

Third, utilizing appropriate tools and technology can help improve public health response and decision making. The quick accessibility measure and geospatial analysis used in this study could be applied prospectively during an ongoing mass vaccination campaign to provide essential and timely evaluation of resource accessibility to guide decision making. County health departments can repeatedly run such analyses with various scenarios to guide site selection. Continued analyses could be done to monitor the change of accessibility landscape as vaccine shipments and daily operations of vaccination sites change from time to time.

The accessibility measure used in this study can be quickly calculated with readily available data. Improvements can be made to increase the accuracy of accessibility evaluation, such as considering distance decay (21) or using dynamic catchments (27). Catchment areas could also be delineated based on travel time and transportation modes such as driving, public transportation, and walking could be separately considered. Ecological fallacy is another factor that may affect the accuracy of our accessibility analysis. As we use aggregated data at the census tract level, the inferences about the group may differ from the real experience of individuals (35). Additional analysis is needed to provide more in-depth investigations to individual's experience and perception of accessibility. Additionally, if detailed vaccination record is available, gravity models could be employed to study the flow of population to specific vaccine sites and examine other influencing factors of vaccination such as characteristics of vaccine site locations (17). Lastly, regardless of the setup of the county mass sites and FEMA site prioritizing accessibility for populations disproportionately affected by the pandemic, there is still disparity in vaccination rate among populations and communities. Our study only focused on accessibility to vaccine resources and did not consider one major factor that influences vaccination rates, which is vaccine acceptance/hesitancy (36–39). More in-depth studies on this front could provide valuable insights to innovative solutions in mass vaccination campaigns. The findings could also inform the constraints built in optimization-based methods (18) regarding individual choices.

This pandemic has shown how marginalized and minority groups have and will suffer disproportionately due to the inequities in society perpetuated by systematic practices. This study is timely to investigate one aspect of health equity, raise consciousness, and use tools and resources to confront inequities. Vaccine distribution system design can greatly influence equity and accessibility at the community level (18). As public health decision makers at different levels set up vaccination sites, evaluating accessibility helps to inform policy and improve coverage and accessibility for disproportionately affected populations. During a vaccination campaign, health equity should be addressed through incorporating social factors in the distribution of vaccine resources to mitigate disparities and prioritize accessibility of the disproportionately affected groups. The results from this study indicate that improving accessibility of these groups can be achieved when site selection considers explicitly the socioeconomic landscape of the population. The methodology employed in this study provides a tool for quick and timely assessment of resource accessibility to make necessary public health decision adjustments during the pandemic.