Rural–Urban Disparities in Realized Spatial Access to General Practitioners, Orthopedic Surgeons, and Physiotherapists among People with Osteoarthritis in Alberta, Canada

Rural Canadians have high health care needs due to high prevalence of osteoarthritis (OA) but lack access to care. Examining realized access to three types of providers (general practitioners (GPs), orthopedic surgeons (Ortho), and physiotherapists (PTs)) simultaneously helps identify gaps in access to needed OA care, inform accessibility assessment, and support health care resource allocation. Travel time from a patient’s postal code to the physician’s postal code was calculated using origin–destination network analysis. We applied descriptive statistics to summarize differences in travel time, hotspot analysis to explore geospatial patterns, and distance decay function to examine the travel pattern of health care utilization by urbanicity. The median travel time in Alberta was 11.6 min (IQR = 4.3–25.7) to GPs, 28.9 (IQR = 14.8–65.0) to Ortho, and 33.7 (IQR = 23.1–47.3) to PTs. We observed significant rural–urban disparities in realized access to GPs (2.9 and IQR = 0.0–92.1 in rural remote areas vs. 12.6 and IQR = 6.4–21.0 in metropolitan areas), Ortho (233.3 and IQR = 171.3–363.7 in rural remote areas vs. 21.3 and IQR = 14.0–29.3 in metropolitan areas), and PTs (62.4 and IQR = 0.0–232.1 in rural remote areas vs. 32.1 and IQR = 25.2–39.9 in metropolitan areas). We identified hotspots of realized access to all three types of providers in rural remote areas, where patients with OA tend to travel longer for health care. This study may provide insight on the choice of catchment size and the distance decay pattern of health care utilization for further studies on spatial accessibility.


Introduction
Osteoarthritis (OA) is a significant cause of pain and disability affecting 1 in 8 (13%) Canadians [1]. The prevalence of OA is expected to increase, reaching 1 in 4 Canadians by 2040 due to ageing population and increasing prevalence of obesity [2,3]. Although there is no known cure for OA, early treatment and management are critical to slow the disease progression, reduce symptoms, and improve quality of life [4]. The vast majority of OA care takes place in community-based clinics and by primary care physicians. Pain is the well as the distance decay pattern of health care utilization. Health care utilization patterns vary across population groups, health care providers, and the spectrum of diseases [23,32]. Although numerous studies on spatial access have been published for different populations (general population [27], immigrants [33], older population [31]), different health services and specialties (hospitalizations [32], medical clinics [34], primary care [27], specialty [35]) and different diseases (hypertension and diabetes [36]), few studies have examined the issue comprehensively among people with OA across health care providers. Realized access refers to the actual utilization of services by patients overcoming all barriers (e.g., language, physical disabilities, finance, and shortage of physicians) [10,37]. One measurement of a barrier or facilitator to realized access is the driving time or straight-line distance between a patient's home address and the location of health services they utilize [22,23]. The measure of travel time and distance as a barrier or facilitator to access using empirical access behavior data is important to examine how far patients travel to seek health care and how travel distance and/or travel time interacts with health care utilization.
Examining realized access to three providers simultaneously provides a relatively comprehensive picture of gaps in access to needed OA care. It provides insight on the choice of catchment size and the distance decay pattern of health care utilization for accessibility assessment and health care resource allocation. It is important to provide information for improved identification of gaps in access to care that is provided by interprofessional teams and support evidence-informed decisions on coordinated and integrated care that is patient-centred [19]. In this study, we measured the realized access to three different types of health care providers for people with OA using utilized physician-seeking data. The aim of this study is three-fold: (1) to measure the rural-urban disparities in realized access to three types of health care providers (GPs, Ortho, and PTs); (2) to examine the spatial pattern of realized access at the local geographic area level; (3) to examine how the distance decay pattern of health care utilization varies along the rural-urban continuum.

Materials and Methods
Cross-sectional study design was applied to examine the rural-urban disparities in realized access to multidisciplinary health care providers for adult people with OA. An adult OA-prevalent cohort (≥18 years of age at diagnosis) in fiscal year 2013 (1 April 2012-31 March 2013) was identified using the most current and validated case definition for administrative data [3,14].

Data Sources
The health care services provided by GPs, Orthos, and PTs in fiscal year 2013 were of interest in this study. As community-based physiotherapy is not an insured service, this study only captured those publicly funded outpatient PT visits which were supported by the Alberta Health Care Insurance Plan (AHCIP). A prevalence cohort of 359,638 adult OA patients (≥18 years of age at diagnosis) in fiscal year 2013 was identified using the most current and validated case definition for administrative data [3,14]. We obtained administrative data from five Alberta Health (AH) administrative databases: AHCIP population registry, the Physician Claims Database (Claims), the Discharge Abstract Database, and the Ambulatory Care Classification System/National Ambulatory Care Reporting System. All the OA-prevalent patients had to meet the following two criteria: (1) became an OA patient in Alberta between 1994 and 2013; (2) did not migrate out of the province or die between 1994 and 2013 [3]. A detailed description of the case definition was reported previously [3,14]. Unique deidentified patient IDs were used to link patients with health records and health care providers.
We obtained the provider type and the six-digit postal code of their practice location from the Claims Database and the patients' six-digit postal codes from the AHCIP population registry dataset. Patients were excluded if no health care was utilized for any of the three types of providers. Both patients and providers were geocoded using the 2013 Postal Code Translator File [38] created by AH.

Standard Geographic Areas
A set of standard geographic areas was jointly created by Alberta Health Services (AHS) and AH for planning, surveillance, monitoring, and reporting of population health, health outcomes, and health services across Alberta [39] (Figure 1), including five zones (North, Edmonton, Central, Calgary, and South) and 132 local geographic areas (LGA)-the lowest geographic level in Alberta. Rural-urban continuum areas were further developed to analyze, report, plan, monitor, and compare population health by rural and urban status which grouped 132 LGAs into seven distinct categories (Metro, Moderate Metro influence, Urban, Moderate Urban influence, Rural Centre, Rural, and Rural Remote) [39].
We obtained the provider type and the six-digit postal code of their practice location from the Claims Database and the patients' six-digit postal codes from the AHCIP population registry dataset. Patients were excluded if no health care was utilized for any of the three types of providers. Both patients and providers were geocoded using the 2013 Postal Code Translator File [38] created by AH.

Standard Geographic Areas
A set of standard geographic areas was jointly created by Alberta Health Services (AHS) and AH for planning, surveillance, monitoring, and reporting of population health, health outcomes, and health services across Alberta [39] (Figure 1), including five zones (North, Edmonton, Central, Calgary, and South) and 132 local geographic areas (LGA)the lowest geographic level in Alberta. Rural-urban continuum areas were further developed to analyze, report, plan, monitor, and compare population health by rural and urban status which grouped 132 LGAs into seven distinct categories (Metro, Moderate Metro influence, Urban, Moderate Urban influence, Rural Centre, Rural, and Rural Remote) [39].

Travel Time Calculations
Origin-destination network analysis was applied to calculate travel time in minutes from the origin of a patient's postal code to the destination of the physician's postal code. The travel time was calculated with the assumption that patients travelled by car under optimal conditions at posted speed limits, without delays due to rush hour, storms, road closures, detours, etc. AHS modeled a road network dataset from the DMTI Route Logistics Road Network [40] (2014), which classified roads into primary roads, secondary roads, and local roads. The Alberta Municipality Data Sharing Partnership (AMDSP) [41] road data (2013) were applied further to modify the DMTI road data and speed limits as AMDSP provided accurate speeds for the majority of municipalities and town in central and southern Alberta [42]. We adjusted the speed limits of road segments according to the AMDSP road data where they were available. Given the long winter time in Alberta, we applied a 20% reduction in the speed limit for travel distances of less than 50 km to obtain

Travel Time Calculations
Origin-destination network analysis was applied to calculate travel time in minutes from the origin of a patient's postal code to the destination of the physician's postal code. The travel time was calculated with the assumption that patients travelled by car under optimal conditions at posted speed limits, without delays due to rush hour, storms, road closures, detours, etc. AHS modeled a road network dataset from the DMTI Route Logistics Road Network [40] (2014), which classified roads into primary roads, secondary roads, and local roads. The Alberta Municipality Data Sharing Partnership (AMDSP) [41] road data (2013) were applied further to modify the DMTI road data and speed limits as AMDSP provided accurate speeds for the majority of municipalities and town in central and southern Alberta [42]. We adjusted the speed limits of road segments according to the AMDSP road data where they were available. Given the long winter time in Alberta, we applied a 20% reduction in the speed limit for travel distances of less than 50 km to obtain a conservative estimate of winter travel time [22]. Ferry routes, one-way roads, and trails were used as restrictions within the analysis. A snap distance of 5 km was applied to capture postal codes within 5 km of the roads. Time-based turn restrictions (such as prohibited left turn from 6:00 to 18:00) and delays at signals/stop signs were not included in the analysis as point data locations and delay times were not available. The road network was validated against the nonurgent emergency medical services (EMS) trips covering Alberta in 2013 [42], showing a difference of 0.26 min between the calculated travel time and the EMS travel time. Detailed information on analysis setting is available from the AHS document [42].

Descriptive Statistics
We grouped patients with OA by sex, age, rural-urban continuum, and type of health care providers (GPs, Orthos, and PTs). The numbers of patients/visits were presented for each subgroup. The percentage of patients/visits was calculated as the number of patients/visits within a subgroup divided by the number of patients/visits in Alberta. Boxplots of travel time in minutes by type of providers were provided to capture the distribution of travel times. Due to the positive skewness of travel time, as shown in Appendix A, the data on the boxplots were truncated at the 95th percentile [43] of travel time.

Hot Spot Analysis at the LGA Level
We wanted to determine if there was a spatial pattern in realized access at the LGA level to GPs, Orthos, and PTs among people with OA. Median travel time was used to represent the level of realized access within each LGA due to the skewed distribution of travel time values. Hot spot analysis based on the Getis-Ord Gi* statistic identified hot spots of realized access (LGAs with a significantly long median travel time surrounded by other LGAs with similar values) and cold spots of realized access (LGAs with a significantly short median travel time surrounded by other LGAs with similar values) [44][45][46]. The underlying null hypothesis of hot spot analysis was that the median travel time of each LGA followed a random distribution across Alberta. The alternative hypothesis was that it was nonrandom. Using the Getis-Ord Gi* routine, we calculated the Z-score and the P-value for each LGA to determine statistical significance of hot spots and cold spots [45]. The results of the hot spot analysis are presented in Section 3.4.
We captured the spatial interaction among neighboring LGAs using a spatial weights matrix that is based on the number of nearest neighbors or the distance between two LGAs [46]. Spatial weights are equally assigned to the LGAs within a defined neighborhood, indicating an equal influence on the LGA of interest [46]. We applied hot spot analysis using three different spatial weight matrices (eight nearest neighbors, fixed distance of 40 km with at least eight nearest neighbors, and queen's continuity with at least eight nearest neighbors) to examine the sensitivity of the identified hot spots to different models of spatial interaction [17]. Detailed results are reported in Appendix B. Based on our previous study of OA prevalence patterns at the LGA level [17] as well as for maintaining consistency in the results of the hot spot analysis, we only reported results using the eight nearest neighbors spatial weights matrix.

Distance Decay Pattern of Health Care Utilization
Our null hypothesis was that there was no relationship between distance and health care utilization. Our alternative hypothesis was that there was a relationship and we tested a variety of alternatives represented by distance decay functions. The distance decay pattern depicted the principle that health care utilization decreases along with travel time. In this analysis, the cumulative probability approach [32,47] was applied to describe the distance decay behavior of health care visits to three types of providers along the rural-urban continuum. Health care visits were grouped by rural-urban status and then cumulated inversely by travel time in minutes. The distance decay function could be estimated by fitting a continuous relationship of the cumulative proportion of health care utilization with corresponding travel times, which was formulated as follows: where d is the travel time in minutes from a patient's residential postal code to the provider's practice postal code; Y d is the cumulative probability of health care utilization with the travel time greater than d; µ is a scalar; f (d) is a set of three popular distance decay functions found in the literature, including exponential function (e −βd ), log-logistic function 1/[1 + (d/θ)ˆβ], and power function (d −β ) [26,32]. Parameters β and θ are to be estimated for distance decay functions, indicating how fast the utilization of health care decreases with travel time.
The nonlinear least square estimator was used to estimate the parameters, and pseudo-R 2 and the Akaike information criterion (AIC) were calculated to identify the best-fitting models [32]. To aggregate the number of health care visits to providers by travel time, we defined 68 travel time datapoints including 31 points within the 0-30 min range (every minute), 18 points from 31 to 120 min (every 5 min), 12 points from 121 to 240 min (every 10 min), six points from 241 min to 600 min (every 60 min), and the maximum (1 point). The model for Alberta and each rural-urban continuum had 68 observations to fit the distance decay functions.
Travel time calculation and hot spot analysis were conducted using ArcMap 10.8. Distance decay functions were estimated using R 4.0.2, the "aomisc" package. Ethics approval for this project was provided by the Conjoint Health Research Ethics Board at the University of Calgary (REB13-0100).  Table 1). The women utilized approximately 60% of the visits with each type of providers. The patients aged between 55-74 years accounted for approximately 50% of the visits with each type of health care providers. The patients residing in metropolitan areas accounted for 50% of the health care visits, compared to 20% in rural/remote areas. The number and percentage of patients and visits by provider type, age, sex, and comorbidities across the rural-urban continuums are shown in Appendix A.

Distribution of Patients and Providers
We identified 2312 GP practices (54% in metropolitan areas vs. 16% in rural/remote areas) accessed by the 170,342 patients residing in 45,614 postal code areas. We identified 243 Ortho practices which provided services to the 47,370 patients from 23,845 postal code areas. We identified 87 PT clinics (14 in metropolitan areas vs. 55 in rural/remote areas) that served the 39,923 patients from 20,039 postal code areas.

Travel Time in the Rural-Urban Continuum
The median travel time in Alberta is 11.6 min (IQR = 4.3-25.7) to visit a GP, 28.

Hot Spots of Median Travel Time at the LGA Level
The median travel time to a GP at the LGA level within the rural and remote areas is not homogeneously distributed (Figure 3a). The top 20% longest median travel times were observed mainly in northern remote, moderate metropolitan, and southwest rural areas. We identified hot spots of median travel time in northern rural remote areas and cold spots in western rural areas (Figure 3a).
The median travel time to an Ortho was shortest in urban and metropolitan areas, which increased along with the increase in rurality. We observed a 114.7-fold difference in travel time to an Ortho at the LGA level (4.2 min in urban Red Deer vs. 481.6 min in rural remote High Level). We identified nine hot spots in northern remote areas and five hot spots in eastern central rural areas. All the cold spots were identified in Metro-Edmonton and Metro-Calgary (Figure 3b).
Similarly, the median travel time to PT services was not evenly distributed within the rural and remote areas. Hot spot analysis identified three hot spots in northern rural remote areas, six hot spots in southern rural areas and the urban Lethbridge and the surrounding areas. All the cold spots were identified in the urban Red Deer and the surrounding areas (Figure 3c).   The median travel time to a GP at the LGA level within the rural and remote areas is not homogeneously distributed (Figure 3a). The top 20% longest median travel times were observed mainly in northern remote, moderate metropolitan, and southwest rural areas. We identified hot spots of median travel time in northern rural remote areas and cold spots in western rural areas (Figure 3a). The median travel time to an Ortho was shortest in urban and metropolitan areas, which increased along with the increase in rurality. We observed a 114.7-fold difference in travel time to an Ortho at the LGA level (4.2 min in urban Red Deer vs. 481.6 min in rural remote High Level). We identified nine hot spots in northern remote areas and five hot spots in eastern central rural areas. All the cold spots were identified in Metro-Edmonton and Metro-Calgary (Figure 3b).
Similarly, the median travel time to PT services was not evenly distributed within the rural and remote areas. Hot spot analysis identified three hot spots in northern rural remote areas, six hot spots in southern rural areas and the urban Lethbridge and the surrounding areas. All the cold spots were identified in the urban Red Deer and the surrounding areas (Figure 3c).

Distance Decay Pattern of Health Care Utilization
The pseudo-R 2 and AIC values for the three distance functions are reported in Table  2. The log-logistic function fitted best with higher pseudo-R 2 and lower AIC values compared to the exponential function and the power function for most rural-urban groups in this study. As shown in Table 2, all parameters β and θ of the log-logistic models were statistically significant (p < 2 * 10 −16 ). These models produced low residual standard errors (RSE: 0.02-0.09 for GP visits; 0.01-0.09 for Ortho visits; 0.02-0.09 for PT visits) with a good curve fit [32,47] (Figure 4).

Distance Decay Pattern of Health Care Utilization
The pseudo-R 2 and AIC values for the three distance functions are reported in Table 2. The log-logistic function fitted best with higher pseudo-R 2 and lower AIC values compared to the exponential function and the power function for most rural-urban groups in this study. As shown in Table 2, all parameters β and θ of the log-logistic models were statistically significant (p < 2 * 10 −16 ). These models produced low residual standard errors (RSE: 0.02-0.09 for GP visits; 0.01-0.09 for Ortho visits; 0.02-0.09 for PT visits) with a good curve fit [32,47] (Figure 4).   The distance decay patterns of health care utilization were captured by the fitted log-logistic functions, highlighting the association of increasing travel time with decreasing health care utilization. As shown in Figure 4a, most GP visits in urban/metropolitan areas occurred within 30 min from a patient's residence (90%), whereas in rural/remote areas, the percentage dropped to 67%. Although 60% of the rural and remote patients travelled to a GP within 15 min, 24% travelled over 60 min. The remote patients showed the strongest distance decay effect within 15 min travel time, but also the weakest distance decay effect on GP visits when the travel time was longer than 30 min. As shown in Figure 4b, 90% of the patients in metropolitan and urban areas travelled up to 45 min for Ortho visits, which increased to 240 min in rural/rural centre areas and to about 510 min in remote areas. The effects of increasing travel time on Ortho visits were strongest in metropolitan/urban areas, while rural remote areas had the slowest decline in the number of Ortho visits with travel time. As shown in Figure 4c, urban areas had the strongest distance effects on PT visits with the steepest decline in the number of PT visits within 15 min travel time. Most patients had access to a PT within 60 min (90% in metropolitan areas, 70% in rural/rural centre, 60% in remote areas).

Discussion
Using administrative data, we observed significant rural-urban disparities in spatial access to health care providers along the rural-urban continuum. The rural remote patients had the shortest median travel time to visit a GP, while urban patients had the easiest access to an Ortho and a PT. Hot spot analysis revealed that northeastern rural remote areas experienced greatest geographic barriers to all the three types of providers. The rural remote patients were the least sensitive to the increase in travel time when seeking health care.
We observed that the median travel time to GPs in metropolitan areas was 12.6 min in Alberta, which is consistent with the literature. As reported in the literature, various travel threshold times to a GP, ranging from 9 to 30 min, have been used in the two-step floating catchment area method (2SFCA) [24][25][26] and other accessibility models, including the 11 min threshold to a GP in Toronto [30] and the maximum of 30 min catchment size in metropolitan Australia [48].
The significant rural-urban disparities in realized access to health care providers suggested that varying the catchment size may be required to capture the interaction between patients and health care providers. Catchment size delimits how far geographically services delivering health care to patients and at the same time determines how far patients travel to access the services. Distance and geographical isolation are the foremost health care access barriers [49,50]. As population dispersion increases, greater travel times are required to access nearby health care providers in more dispersed areas. Varying catchment size captures diverse travel patterns by rural-urban status and provides accurate estimates of access to care [48]. The results shed light on the selection of appropriate catchment sizes along the rural-urban continuum, which will inform health care planners and policymakers on health resource allocation by defining health services areas based on travel time thresholds.
Our study utilized empirical data on actual physician-seeking behavior, which provided evidence on the choice of catchment sizes and distance decay association when modeling spatial accessibility. The measure of spatial accessibility can be broadly spilt into realized accessibility and potential accessibility [22]. Potential accessibility refers to potential availability of health services to potential users [19], which is usually measured with floating catchment area (FCA) metrics, such as the 2SFCA [24][25][26] and the enhanced 2SFCA [27]. To date, many applications of the FCA method have not been verified against empirical access behavior data [24][25][26]. A key limitation in modeling spatial accessibility is the lack of available empirical data on the actual physician-seeking behavior [28][29][30][31]. The FCA method assumes that patients visit the nearest health care facility within a specific catchment area, which ignores patient and system barriers to health care such as language, physical disabilities, and health literacy [10]. Furthermore, the FCA method requires assumptions of catchment size [33] and distance decay association between travel time/distance and health care utilization [32], which may vary by rural-urban status and type of health services [32,48]. However, the choice of reasonable catchment sizes and appropriate distance decay effect by rural-urban status and type of health services cannot be settled without empirical utilized health care data [30,31,48,51]. Our study sheds light on the key components of spatial accessibility modeling-choice of catchment sizes and distance decay effect on health care utilization, which may be applied to other regions that have a shortage of physicians in rural and remote areas due to misalignment with population health care needs.
We observed large variations in travel even within the same rural-urban continuum as shown in Figure 2 and Appendix B. Almost half of the health care visits to a GP in rural remote areas happened within the local postal code, while 95% of GP visits took up to 361.8 min. Urban areas had the shortest median travel time to an Ortho and a PT; however, 95% travelled up to 352.7 min to Ortho and up to 422.6 min to a PT, which was the second longest travel time to visit health care providers. As shown in the Figure 3, both the easiest (cold spot) and the most difficult (hot spot) access to a GP were identified in northern rural remote areas. The northern urban centre had significant geographic barriers to Ortho when compared to southern urban areas. Further studies at the local level as well as the individual level may help understand the most important factors driving this travel pattern.
The findings of the hot spot analysis shed light on health resource allocation to reduce rural-urban disparities in access to care. The government of Alberta initiated the Alberta Surgical Initiative in 2020 [52], emphasizing the importance of expanding telephone and electronic medical advice programs through which primary care providers can receive timely consultations from specialists for patients under care. Within northeastern remote areas where the patients had geographic barriers to all the three types of providers, health policies may be targeted to retain and recruit physicians and health care professionals (e.g., nurse practitioners and PTs) [53], as well as promote innovative methods of delivering health care (such as virtual care, TeleHealth, travel clinics, and integrated care) in remote areas [54,55]. For the northern remote areas that had a relatively easy access to GPs and PTs and, at the same time, a significantly difficult access to Ortho, integrated care may be targeted to increase Ortho capacity by having primary practitioners, PTs, and/or nurse practitioners to deliver timely advice to patients.
Physiotherapists play an integral role in optimizing OA patient outcomes before and after surgery [56,57], reducing health care costs [58], and increasing patient satisfaction [59]. However, it is challenging to evaluate rehabilitation services across Alberta due to variation in practice and funding resources [60]. The model of care was zone-specific. There was no provincial mechanism to collect, collate, and report data/outcomes related to physiotherapy services [60]. A previous study reported that reimbursement was linked to services rendered and the available human resources, showing that the ability or willingness to pay is a barrier to access to private community physiotherapy services in Alberta [61]. Health system efficiency and resource allocation are needed for rehabilitation services.
This study has numerous strengths. Firstly, we utilized empirical data to investigate the realized access to different health care providers within the rural-urban continuum for people with OA. This analysis is informative and will assist health planners and health service delivery by defining appropriate catchment sizes of health services, relocating health resources, and reducing the rural-urban disparities in access to health and health outcomes. Secondly, this study fills a gap in the modeling of spatial accessibility by providing evidence on the choice of catchment areas and distance decay effect on health care utilization along the rural-urban continuum. Thirdly, given the commonly reported outcomes of pain and reduced physical activity among people with OA, people with OA may present a different health care-seeking behavior as compared to other diseases [62]. Our study specifically focused on OA, providing evidence on the travel pattern of people with OA. Fourthly, we applied a validated algorithm based on administrative data to identify physiotherapy services that are coved by AHCIP. Lastly, we obtained the patients' addresses at the six-digit postal code level from the administrative dataset, the finest spatial scale available for spatial analysis, which allows for the greatest detail regarding travel times.
We acknowledge the following limitations. Firstly, the travel time was calculated between the population-weighted centroids of two postal codes instead of geocoded addresses. However, due to privacy and confidentiality, six-digit postal codes at the patient level constitute the most detailed spatial scale for spatial analysis. Among the 57,874 postal codes identified with at least one OA case in Alberta, the average nearest distance between the postal codes ranged from 0.06 km in metropolitan/urban areas to 6.62 km in rural remote areas [14]. Given the large variation in the sizes of postal code areas, the travel time of rural patients visiting providers in their home postal code area may be longer compared to that of their urban counterparts. Secondly, the location and number of practices do not equal the number of providers. The capacity of each clinic was not accounted in this analysis. Thirdly, the travel time was calculated assuming optimal travel conditions. However, to make the best estimates, we utilized the impedance travel time which accounted for real-world conditions by taking a 20% reduction in speed limit [63]. Fourthly, the findings of the rural-urban disparities in access to OA care provide insights to the current health care delivery. Between 2014 and 2018, Alberta had fast growth in the number of doctors in urban Alberta (14% increase) compared to rural Alberta (6% increase). In particular, the number of doctors in rural Alberta in 2018 was down 3.6% from 2017. Though we saw a growth in the number of physicians in Alberta, physician service shortages continue in rural and remote areas of the province as well as in some underserved urban areas. Further, PT services in this study included only publicly funded outpatient visits, while those visits not supported by AHCIP were excluded. This study focused on the pattern of spatial access to three types of providers for people with OA. Confounding factors such as socioeconomic status on the pattern of realized travel time was not accounted for in this analysis due to unavailability of individual-level data. As a referral is required to visit an Ortho, the Ortho utilization pattern may be conditionally associated with the distribution and referral pattern of GPs. Lastly, as this analysis focused on the key providers relevant to OA management including GPs, Ortho, and PTs, patients who did not seek care from any of the three types of providers were not captured in this study. Further studies on the geospatial pattern of those OA patients not included in this analysis may be used as a weight for estimating the overall burden of geographic residence in assessing care.

Conclusions
This study examined the rural-urban disparities in realized access to health care providers for patients with OA. It is important to provide health planners with evidence on the thresholds of travel time to different health care providers along the rural-urban continuum. These findings may inform policies on accessibility evaluation and health resource allocation in Alberta. Further studies on the driving factors of travel patterns will be of interest. In the context of the COVID-19 pandemic where virtual care has been applied for many health services, the interpretation and meaning of "accessibility" may change accordingly. Further questions may be asked to answer to what extent the provision of virtual care may improve access to care for rural patients and what factors may affect equitable access to virtual care.  Data Availability Statement: Restrictions apply to the availability of these data. Data was obtained from Alberta Health Services and are available from the corresponding author with the permission of Alberta Health Services.

Conflicts of Interest:
The authors declare no conflict of interest.