## 1. Introduction

The connection between public health and geography can be traced back to Hippocrates (c. 400 BC) who deduced that spatially varying factors such as climate, elevation, environmental toxins, ethnicity and race contributed to the spatial patterns of illness (Parchman *et al.*, 2002). The observations of Hippocrates still hold true today and these relationships between geography and disease have allowed geospatial methods to become valuable within the field of public health. Maps have long been a useful tool for visualizing patterns in health care. One of the earliest and most commonly cited examples is from the mid 1800s when Dr. John Snow deduced the source of a cholera outbreak in London based on a simple visualization of the incidents of cholera in relation to water pumps (Johnson, 2007). Although visualizing data geographically is still very valuable for uncovering patterns and associations over space, geospatial analysis has become more sophisticated over time.

One tool that can be used to apply advanced geospatial methods to health care problems is a geographic information system (GIS). The power of GIS lies in its ability to analyze, store and display large amounts of spatially referenced data. In a field where manual data analysis can become overwhelming, GIS is a valuable tool. The medical research applications of GIS are numerous and include finding disease clusters and their possible causes (Murray *et al.*, 2009; Srivastava *et al.*, 2009), improving deployment for emergency services (Ong *et al.*, 2009; Peleg & Pliskin, 2004) and determining if an area is being served adequately by health services (Cinnamon et al., 2009; Schuurman *et al.*, 2008). There have been several reviews and textbooks published in the past decade that focus on the application of GIS to different areas of health research (Cromley & McLafferty, 2002; Kurland & Gorr, 2009; McLafferty, 2003; Parchman *et al.*, 2002; Rushton, 2003).

While the use of GIS is gaining favor with health researchers, barriers remain in the uptake of more advanced geospatial methods by health care decision-makers. A recent qualitative study conducted in the UK found that although health care decision-makers see the value of using GIS in the decision making process, many (especially those in the community setting) still view GIS primarily as a visualization tool (Joyce, 2009). This chapter aims to go beyond a simple discussion of GIS and its standard capabilities (e.g. mapping, buffering, clipping, intersecting, etc.) in order to focus more broadly on geospatial methods (including spatial statistics) that can be used to describe and study the distribution of disease and health care delivery. It adds to current literature by taking an approach that compares and contrasts different geospatial methods with the goal of informing health care decision-makers about the strengths and weaknesses of each.

This chapter is written for an audience that has a basic knowledge of statistics and GIS. It does not give a comprehensive overview of all geospatial methods that are available for health research but instead highlights ten different geospatial methods that can be used to address problems posed when studying the distribution of disease and health service delivery. Five topic areas will be covered: representing populations, identifying disease clusters, identifying associations in a spatial context, identifying areas with access to health care and identifying the best location for a new health service. Although each method is discussed briefly, references will be provided where health researchers and decision-makers can find details. Each of these five sections will be further subdivided into two sections. ‘Representing populations’ will discuss the topics of *Geographically linked census data* and *Dasymetric modeling*. ‘Identifying disease clusters’ will discuss the topics of *Measures of spatial autocorrelation* and *Kernel density estimation*. ‘Identifying associations in a spatial context’ will cover *Spatial regression* and *Geographically weighted regression*. ‘Identifying areas with access to a health care service’ will discuss the topics of *Road network analysis* and *Interpolation*. Finally, ‘Identifying the best location for a new service’ will discuss the topics of *Coverage problems* and *Distance problems*.

## 2. Representing populations

Studying relationships between disease or health service delivery and environmental or population characteristics, requires researchers to represent populations in terms of their location in space. Two geospatial methods that can be used to represent population groups will be discussed in this section: geographically linked census data and dasymetric modeling.

### 2.1. Geographically linked census data

Studies in health research may sometimes be focused on the general population that is at risk for a disease or that has access to a specific health service in an area. This group can be represented by population counts in the form of census data. Census data can be represented at many different levels. For example, in Canada census data may be represented at the census block level, the dissemination area level, and at the census tract level (Apparicio et al., 2008). This method of representing population links census data tables to geographic boundary files using unique area level identifiers. The aggregation of populations into larger areal units can affect the results obtained from a spatial analysis. For example, when studying geographic access to urban amenities, Canadian studies have found that when accessibility was measured for smaller spatial units, the measures were less subject to aggregation error than when larger spatial units were used to represent populations (Apparicio et al., 2008; Hewko et al., 2002).

There are limitations to geospatial analysis when populations are represented in this way. As discussed above, the aggregation errors differ depending on the spatial unit used. One study showed that for some amenities the aggregation to larger areas for abundant resources (e.g. playgrounds) should use smaller areal units more representative of where the population is distributed (Hewko et al., 2002). Another limitation occurs because census data is standardized based on population numbers. This causes rural areas to be represented by larger areal units that do not precisely represent locations of where the population is distributed. Finally, the modifiable areal unit problem (MAUP) should be considered when using different census boundaries to represent populations. The MAUP occurs because the results of any geographic analysis can be sensitive to the areal units that are being used for the study (Fotheringham & Wong, 1991; Oppenshaw, 1996). For example, results can differ depending on the level of census data used in analysis. Notable strengths of using census data to represent populations are: 1) in most developed countries this data is largely publically available and 2) this data can be grouped by age categories and linked with other census variables of interest to study associations (with, for example, income, employment, education, ethnicity and other socio-economic indicators).

### 2.2. Dasymetric modeling

Census data is often used by researchers to represent population numbers evenly within a specified area. For example, for analytical purposes the population within a census tract is often represented by a single point at the geographic centre of the polygon, known as the centroid of the area. A method used to represent population distribution more accurately is dasymetric modeling. Dasymetric modeling uses additional data sources to estimate the distribution of aggregated population data within a specified unit of analysis (Briggs et al., 2007a; Mennis, 2009; Poulsen & Kennedy, 2004). This technique essentially uses a combination of variables to represent where population is truly located within an area. Dasymetric models have used light emission and remotely sensed land cover data to represent populations at the small area level in Europe (Briggs et al., 2007b). One study has shown how population estimates in residential areas can be more accurately represented using grid based land cover data in combination with a mailing database to represent where residential populations are located within a census area (Langford *et al.*, 2008). Dasymetric modeling produces population numbers that are identical to those within a census area; it is the spatial distribution of the population that is modified within these spatial units. A major limitation of using this method for representing population distribution is that it is more time intensive and requires the use of additional data sources when compared to simply using the raw census data. The strength of this method is that it more accurately represents population distribution in large areas (e.g. rural census areas). One study found lower access to primary health care services using dasymetric modeling when compared to using point census data to represent populations in rural areas (Langford & Higgs, 2006). This shows that census methods of representing population data could be overestimating access in rural areas.

Both methods for representing population distribution have been shown to have strengths and limitations. In the absence of patient data or when the general population is of interest, it is possible to use census data at varying levels of aggregation. Using this method, populations can be represented as a centroid within the areal unit or as being evenly distributed throughout the area. Dasymetric mapping allows for a more accurate representation of the population distribution within areal units but at the cost of increased time and data collection effort. However, this method can be valuable when population data is represented for rural areas. Both methods suffer from the problem of the MAUP and this should be considered when applying any level of population data to analysis. Figure 1 highlights the differences in representing populations using two different representations of census data and dasymetric modeling.

## 3. Identifying disease clusters

The identification of clusters of disease can help to target health care delivery and identify where resource allocation efforts should be focused. Point pattern analysis is a broad category of geospatial methods that can be used to identify disease clusters. In this section we will discuss two specific methods: measures of spatial autocorrelation and kernel density estimates.

### 3.1. Measures of spatial autocorrelation

Spatial autocorrelation can generally be thought of as the association of a variable with itself over space. Statistical measures of spatial autocorrelation can be used to determine if the clustering of variables with similar attributes occurs over space. Moran’s I is a global measure for evaluating spatial autocorrelation. It is measured using the location of a given set of points and their associated attributes over an entire study area. The values of Moran’s I range from -1 to +1 (dispersed to clustered, respectively, with values close to 0 representing no spatial autocorrelation as might be observed in a random pattern) (Rogerson, 2010). One study conducted to evaluate the survival from Cardiac arrest in a US county used Moran’s I to determine that there was no clustering of survival rates for cardiac arrest in the study area (Warden *et al.*, 2007). The Moran’s I statistic can also be used to define the extent of spatial bands used for other spatial statistics as exemplified in a study that determined hot spots using the Getis-Ord Gi^{*} statistic (McEntee & Ogneva-Himmelberger, 2008). The value of the Moran’s I measure is that it allows for an overall assessment of whether spatial autocorrelation exists over a study area. The limitation of the Moran’s I measure is that it does not take into account variation at the local level and reports only one measure for the entire study area.

In order to deal with the limitation of the global measure of Moran’s I, a local measure of spatial association was introduced by Anselin (Anselin, 1995). This local indicator of spatial association (LISA) breaks down the global indicator to allow for an understanding of how each observation contributes to the overall measure. There are three main purposes that the LISA statistic can serve: 1) to indicate local areas of spatial non-stationarity, 2) to evaluate how different locations within an area are contributing to the overall global spatial autocorrelation statistic and 3) to evaluate if spatial outliers exist (Anselin, 1995). A positive LISA value for an area indicates that the feature is surrounded by features with similar values (a cluster). A negative LISA value indicates that the feature is surrounded by features with dissimilar values (an outlier) (Anselin, 1995). LISA has been used to uncover clusters of overweight and obesity in Canada (Pouliou & Elliott, 2009) and clusters of chronic ischemic heart disease in relation to aerosol air pollution in the Eastern United States (Hu & Rao, 2009). As with most geospatial methods, the above methods of cluster analysis are limited by the quality and scale of the data that is available (Shi, 2009).

### 3.2. Kernel density estimation

Quadrat analysis is a method that evaluates the points within polygon areas to assess whether a pattern is random, dispersed or clustered in space. Standard quadrat analysis is limited by edge effects, where the density of points between adjoining areas can change abruptly from one value to another. Kernel density estimates address the limitations of simple quadrat analysis though the use of a kernel function that smoothes the values of interest over space by forming a count of events per unit area within a moving quadrat or ‘window’. This produces a more spatially ‘smooth’ estimate of variation in spatial autocorrelation than can be obtained by using a fixed grid of quadrats (Gatrell et al., 1996). A detailed description of the kernel function is outlined in several resources (de Smith *et al.*, 2007a; Gatrell et al., 1996).

One important consideration when using kernel density estimation is that the degree of smoothing produced is dependent on the size of the bandwidth. It is often useful to examine several density plots of the same samples, all smoothed by different amounts, in order to gain greater insight into the data (Farwell, 1999). Guidelines for choosing bandwidths vary from simple visual inspection, to specific parameters such as specified resolution limits or nearest neighbor distances (Levine, 2007). One study by Lerner *et al.* used kernel density estimation in ArcGIS with default bandwidth settings to identify out of hospital cardiac arrests in Rochester, NY. By integrating census data with the data from a cluster analysis using kernel density estimation, they were able to discover that out of hospital cardiac arrest clusters were associated with communities with lower income, higher number of African-Americans residents and more residents without a high school diploma than the population of the city in general (Lerner *et al.*, 2005). The identification of clusters of events could lead to an increased understanding of where health care resources should be allocated. One of the major limitations of kernel density estimation is that the choice of size of the ‘moving window’ (bandwidth of the kernel function) can alter the results of the cluster analysis (de Smith *et al.*, 2007a). This method is superior to a simple quadrat analysis because it takes into account edge effects and creates a more realistic representation of spatial variation through a smoothing process (Gatrell et al., 1996).

The two methods discussed to identify disease clusters are similar in that they both allow for a statistical measure of spatial association within areas. Using the global and local measures of spatial association can give a quantitative measure of whether global or local spatial autocorrelation (reflecting clustering) exists in a dataset while the kernel density estimates are a way of using data locally within windows to create a grid surface of where clustering is occurring. The spatial measures of association that show whether spatial autocorrelation exists in a dataset can indicate what type of further analysis should be undertaken. The presence of spatial autocorrelation suggests that spatial regression techniques should be used in place of simple linear regression. While the LISA method can be used to measure spatial association within a defined area, the method of kernel density estimation allows for minimizing edge effects between areal units with sharp contrasts. One recent study has shown how LISA and kernel density estimation can be used together to uncover hot spots of diarrhea in Thailand (Chaikaew et al., 2009). Figure 2 shows how kernel density estimation has been used to identify areas of high spousal violence incidents in Calgary, Canada.

## 4. Identifying associations in a spatial context

The purpose of many studies on the distribution of health services is to determine if inequity exists. These studies strive to determine if those who are considered disadvantaged are served as equally by health services as those who are not (Meade & Emch, 2010). Studies may also seek to understand relationships between disease prevalence and various environmental variables. Spatial regression and geographically weighted regression are two spatial statistical methods that can be used to study these types of associations between dependent and independent variables.

### 4.1. Spatial regression

When exploring variables over space it is often found that spatial autocorrelation exists. Spatial dependence violates the assumption of independence needed for standard linear regression (Legendre, 1993). Spatial regression is a valuable method in situations where spatial autocorrelation exists because it considers the spatial dependence between variables by giving weights to them based on distance. The weighting matrix used in the spatial regression process accounts for the spatial dependency in the dependent variable (Fotheringham et al., 2005b). In order to determine accurately the weights used within this matrix it is necessary to first determine the range (distance) of the spatial dependence on the dependent variable using a semivariogram, a graph that visually provides a description of how the data are related (correlated) with distance. Beyond the range defined by the semivariogram there is no longer spatial dependence (Fotheringham et al., 2005b). This incorporation of spatial dependence is important because spatial autocorrelation can reduce the efficiency of the Ordinary Least Squares estimator used in linear regression. Since most features over space have an element of spatial dependence associated with them, it suggests that spatial regression may be a better choice than simple linear regression when modeling features over space, although spatial dependence can be removed using a variety of methods including detrending procedures where spatial trends are evident.

The benefit of using the spatial regression technique is its ability to create a model that can be used to make generalized inferences in an area while accounting for the effects of spatial dependence. A study by Zenk *et al*. examined the association between neighborhood racial composition, neighborhood poverty and spatial accessibility to supermarkets using spatial regression and found that racial segregation places African-Americans in more disadvantaged neighborhoods and at further distances from supermarkets (Zenk *et al.*, 2005). Although spatial regression can account for spatial dependence in a dataset, a recent review states that there is always the possibility of ecological bias when area level data are used to represent individual characteristics within studies (Wakefield, 2007).

### 4.2. Geographically weighted regression

Geographically weighted regression (GWR) is a local regression model that creates unique local coefficients in order to describe the variation in the dependent variable. This creates local models that are well suited to explore the variation in the data graphically (Jetz *et al.*, 2005). Instead of incorporating the spatial dependence into the model as is done through the process of spatial regression, GWR “measures the relationships inherent in the model at each point” (Fotheringham et al., 2005a). The weighting matrix used in GWR is dependent on the location of individual points and thus must be recalculated at each point. This is in contrast to the spatial regression method where a single contiguity weighting matrix is used for the entire study area (based on the extent of spatial dependence on the dependent variable). In GWR, kernels can be used which associate aspects of density of data into the process (Fotheringham et al., 2005a). GWR also has the benefit of being able to evaluate changes in strength of relationships over different scales, which is difficult to do using spatial regression (Jetz et al., 2005). Previous publications have summarized GWR and local analysis in further detail (Brunsdon et al., 1996; Fotheringham et al., 1998).

GWR has been used in a recent study to evaluate the relationship between cold surges and cardiovascular disease (CVD) mortality. This study found that the CVD mortality rates increased significantly after cold surges and that the tolerance of different populations to these weather changes differed over space (Yang *et al.*, 2009). Another study compared the findings from linear regression and GWR when studying access to parks and physical activity sites and found that while OLS regression found weak relationships between park density and physical activity site density, GWR indicated disparities in accessibility that varied over space (Maroko *et al.*, 2009). This study shows the importance of using the appropriate method of analysis for the data that is under study. Yet another study used GWR as a tool to indicate the need for spatial coefficients when linking health and environmental data in geographical analysis (Young *et al.*, 2009).

The main difference between spatial regression and GWR is that one is a global model and one is a local model. Overall, GWR should be used when there is spatial non-stationarity or when relationships between variables differ over space (Fotheringham et al., 2005a). A simple map of the residuals of a model indicate whether a spatial regression method should be employed over a traditional linear regression method; when the residuals are clustered into high and low values over space it is likely necessary to use a spatial regression method. The residuals of a spatial regression or GWR can also be mapped to assess the areas where the models are under or over estimating (Fotheringham et al., 2005b).

## 5. Identifying areas with access to a health care service

The concept of access to health care is a central focus of many health research studies. Spatial accessibility is one aspect of access to health care that can be measured using geospatial methods. Geographic access to a service can be represented as straight line distance, network distance and travel time based on travel along a road network (Apparicio et al., 2008); travel effort can also be represented as a continuous travel time surface over a specified area. The two methods that we will focus our discussion on in this section are road network analysis and interpolation.

### 5.1. Road network analysis

Network analysis has been considered by many researchers to be a more accurate method for evaluating accessibility in terms of travel time and distance when the focus is on travel via roads (Brabyn & Skelly, 2002; Christie & Fone, 2003). These researchers argue that although straight line Euclidean distance and the evaluation of catchments/service areas using Thiessen polygons are commonly used in accessibility analysis, these methods only allow for a crude measure of accessibility (Brabyn & Skelly, 2002). Network analysis allows for a more realistic representation of travel times along a road network by accounting for the different travel speeds along the various road links. Network analysis has been used in health care research to study access to video lottery terminals (Robitaille & Herjean, 2008), access to parks (Comber et al., 2008) and access to palliative care (Cinnamon et al., 2008). The measures of access obtained through network analysis can then be used in combination with census data to assess the population numbers with access to a facility within a specified distance or travel time (Christie & Fone, 2003; Nallamothu *et al.*, 2006). The travel times and distances can also be used to evaluate associations between measures of access and various health indicator (socio-economic) variables (Hare & Barcus, 2007). One UK study found that travel time estimates based on network analysis may be preferred to reported travel times for modeling purposes because the reported times can reflect unusual circumstances and reporter error (Haynes et al., 2006). The strength of road network analysis is that it allows for the estimation of travel times between different point locations along a road network. It does not allow for continuous travel time surfaces to be mapped. In order to do this we can use a method referred to as interpolation.

### 5.2. Interpolation

The method of interpolation is based on Tobler’s first law of geography: things that are near are more similar than things that are further apart (Tobler, 2004). The process of interpolation uses mathematical modeling to ‘fill in the gaps’ between point data values in order to create a continuous surface of values. There are a variety of resources available to researchers that discuss in detail the theoretical basis of the different interpolation methods (Lam, 1983; Mitas & Mitasova, 1999; Waters, 1997a, 1997b). The advantage of interpolation over simple network analysis is that it can allow for a comparison between different modes of transportation to health care facilities. To do this it would be necessary to have point data representing travel times via different modes for this type of analysis. For example, using interpolation, a study conducted by Lerner *et al.* considered how GIS could be used as a tool to help make transport decisions for trauma patients (Lerner *et al.*, 1999). The goal of that study was to create a map that would show where air or ground ambulance should be preferred to transport patients to a trauma center from a given patient starting location. One Canadian study used this method to compare where ground, helicopter or fixed-wing transport was faster when transporting patients to a cardiac catheterization facility (Patel *et al.*, 2007). Figure 3 highlights how interpolation can be used to model travel times via different modes of emergency transport over an area.

There are strengths and limitations to using interpolation as part of a geospatial analysis. One limitation of using this method is that an accurate interpolation requires many points that are evenly distributed over the study area. Another limitation is that depending on the type of interpolation used (e.g. global or local, exact or inexact) the results could vary (Lam, 1983; Mitas & Mitasova, 1999). In spite of these limitations, interpolation is a powerful method to compare access via different modes of travel. This is achieved by evaluating those areas that are served faster by different modes of travel using grid cell comparisons. Both network analysis and interpolation can be used to estimate travel time values that can be used as a measures of geographic access. These measures can then be used to evaluate the associations between spatial access to services and other variables (e.g. socio-economic variables). These access measures can also be linked with census or patient data to study the populations with access to certain services (Nallamothu *et al.*, 2006; Schuurman *et al.*, 2008). It is important to note that because the method of interpolation focuses on creating continuous surfaces from point data, it can also be used to model the concentration of environmental variables from a point source or the spread of an infectious agent over boundaries (Kistemann *et al.*, 2002).

## 6. Identifying the best location for a new service

The spatial analysis of health services is based on the principle that populations and their need for healthcare vary across space. People are not located randomly about the earth and it is often observed that different areas are populated by groups with differing characteristics (e.g. socio-economic status, race, age, income, education, etc.); these varying characteristics of a group often influence their need for health services (McLafferty, 2003). Recent studies have used GIS and the principles of location analysis to evaluate optimal locations for pre-hospital helicopter emergency medical services (Schuurman *et al.*, 2009) and to evaluate how patient access can change depending on where cardiac services are located (Pereira *et al.*, 2007). The operations research literature is a strong source of methods papers and includes literature on how models can be used to maximize service area (Mahmud & Indriasari, 2009), how models can simulate different starting point locations for ambulance (Ingolfsson et al., 2003) and how stochastic approaches can be superior to deterministic approaches when planning health service delivery (Harper et al., 2005). It is beyond the scope of this chapter to discuss the mathematical models that are used for location analysis; this information is available elsewhere (de Smith *et al.*, 2007b; ReVelle & Eiselt, 2005). This section will discuss two different approaches for identifying the best location for a new service: coverage and distance problems.

### 6.1. Coverage problems

There are two types of coverage problems that are commonly referred to in the operations research literature: the location set covering problem (LSCP) and the maximal covering location problem (MCLP). The goal of LSCP is to “find the minimum number of facilities and their locations such that all neighborhoods are covered within the maximal distance or time standard” (Church & Murray, 2009a). In this model, each demand point is covered at least once. The goal of MCLP is to “find a prespecified number of facilities such that coverage within a maximal service distance or time is maximized” (Church & Murray, 2009a). Given these constraints this model seeks to provide the best possible coverage with the limited available resources and does not guarantee service. Geospatial tools available within GIS software can be used to assemble the needed data to address both types of coverage problems. The three essential data components that need to be represented in a GIS in order to solve this problem are the demand points that are to be served by the facilities, the set of possible facility locations and the service area capabilities of these facilities (Church & Murray, 2009a). Ambulance and fire department response are two examples of services that require complete coverage of an area and thus a LSCP approach to coverage. Other services may strive for as much coverage as possible based on the resources available. One example is the placement of public health care facilities. By creating catchment areas within GIS, the best possible locations for new facilities can be evaluated (based on the most area and population covered) (Murad, 2008). The strength of using these types of coverage models is that they provide a spatial-standard based framework for siting new health service facilities that considers population needs and available resources. While these frameworks can be used to study the best locations based on the coverage of a facility over a given area, in other instances a distance based approach may be required.

### 6.2. Distance problems

There are two common approaches to finding an optimal location based on distance for a new facility: p-median and p-centre approaches. The p-median approach focuses on minimizing the average distance/cost to or from patients, while the p-centre approach focuses on minimizing the maximum distance/cost to or from patients. The limitation of solving these problems discretely is that when the number of facilities to be placed increases, the computational intensity of the calculations quickly increases (Church & Murray, 2009b). In order to deal with this limitation, heuristic algorithms are employed that use systematic procedures to trade off the quality of the solution with processing speed (de Smith *et al.*, 2007b). These algorithms strive to find an optimal solution based on specified search strategies (Church & Murray, 2009b). While these methods are founded theoretically in the mathematical modeling and operations research literature, there are a number of components to these strategies that can be evaluated using GIS (Church, 2002). For example, previous studies have used geospatial methods to evaluate the suitability of selected sites as possible locations for a new facility (Cinnamon et al., 2009; Schuurman et al., 2008) and to measure the distance between populations and proposed services (Pereira et al., 2007).

Both the coverage approach and the distance based approach to the placement of new facilities are important. They provide a framework for optimizing the solution to the best location for a new facility. The coverage model focuses on ensuring the population is best served in terms of the placement of a new facility. The location-allocation approach is more focused on the tradeoffs in terms of minimizing the average or maximum distance that patients have to travel in order to access a service. The incorporation of patient service utilization data for both coverage and distance models will add value to the optimal solutions; these data appear to be one component that is lacking in current research that uses GIS and location-allocation approaches to recommend the placement of new services.

## 7. Conclusions

This chapter has shown that the geospatial methods that can be used with GIS extend beyond the simple visualization of patterns. There are powerful geospatial methods available to address problems that are commonly posed by epidemiologists and health services researchers. This chapter goes beyond a simple discussion of the basic methods of data creation and mapping with GIS, to introduce researchers having a basic understanding of GIS and statistical methods to more robust methods for spatial data analysis. It is hoped that a broader understanding of the tools available for spatial analysis will ensure that the most appropriate method to address a spatial problem is used while a consideration is made of the strengths and limitations of the geospatial methods employed in health research. While sophisticated tools have been introduced to deal with spatial data issues such as spatial autocorrelation, the mapping of spatial data remains a useful way of conveying the results of any spatial analysis. In future, health researchers should aim to work in collaborative teams with geographers and operations research specialists to ensure that the most appropriate geospatial methods are being used to study the distribution of disease and health service delivery.