Delineations for Police Patrolling on Street Network Segments with p-Median Location Models

Abstract

Police patrolling intends to enhance traffic safety by mitigating the risks associated with vehicle crashes and accidents. From a view of operations, patrolling requires an effective distribution of resources and often involves area delineations for this distribution purpose. Given constraints such as budget and human resources for traffic safety, delineating geographic areas optimally for police patrol areas is an important agenda item. This paper considers two p-median location models using segments on a street network as observational units on which traffic issues such as vehicle crashes occur. It also uses two weight sets to construct an enhanced delineation of police patrol areas in the City of Plano, Texas. The first model for the standard p-median formulation gives attention to the cumulative number of motor vehicle crashes from 2011 to 2021 on the major transportation networks in Plano. The second model, an extension of this first p-median one, uses balancing constraints to achieve balanced spatial coverage across patrol areas. These two models are also solved with network kernel density count estimates (NKDCE) instead of crash counts. These smoothed densities on a network enable consideration of uncertainty affiliated with this aggregation. The analysis results of this paper suggest that the p-median models provide effective specifications, including their capability to define patrol areas that encompass the entire study region while minimizing distance costs. The inclusion of balancing constraints ensures a more equitable distribution of workloads among patrol areas, improving overall efficiency. Additionally, the model with NKDCE results in an improved workload balance among delineated areas for police patrolling activities, thus supporting more informed spatial decision-making processes for public safety.

1. Introduction

Emergency response vehicles, including ambulances, fire engines, and police cars, play a crucial role in swiftly addressing incidents to minimize casualties and property damage. Unlike other emergency services, police departments have unique responsibilities for prevention in addition to emergency response. Specifically, police must conduct routine patrols to maintain a visible presence and enforce laws, thereby preventing incidents before they occur [1]. However, studies in the literature on police patrolling primarily focus on optimizing dispatch with consideration of factors such as station locations, workforce sizes, and coverage areas based on call-for-service data [2,3], aiming to address the overall responsibilities of police departments. Despite the importance of routine patrolling, there has been limited exploration of location models specifically tailored for this function. Expanding research in this area is crucial for optimizing delineating for police patrolling and enhancing proactive law enforcement efforts.

This study investigates a location model approach for delineating police patrol areas based on vehicle crashes, as they are strongly linked to public safety and can be prevented through effective patrolling. Preventive measures include increasing drivers’ awareness of potential hazards at specific locations and times, while police enforcement plays a crucial role in deterring driving offenses and reducing accidents. Together, these efforts can significantly improve road safety [4,5]. Traffic accidents, in particular, pose one of the most serious threats to urban safety. According to the Centers for Disease Control and Prevention, accidents resulting in unintentional injuries are the fourth leading cause of death in the United States [6]. Since the top three causes are disease-related, accidents are the most prevalent non-health-related cause of death. Specifically, motor vehicle accidents rank as the second leading cause of death from unintentional injuries in the U.S. [7]. Given this, improving traffic safety management is critical for reducing fatalities, and enhanced patrolling strategies play a key role in achieving this goal.

The location model leverages the strengths of the p-median model for delineating police patrol areas, as it helps maintain homogeneity, ensures compactness, and incorporates factors such as road length, crime rates, and population density [8]. Bucarey et al. [9] improved the classical p-median model by setting bounds on demand for each sector, limiting deviations from the average demand to ensure homogeneity. Additionally, Salazar-Aguilar et al. [10] demonstrated that the model minimizes distances within districts and incorporates boundary measures, with compactness serving as a key objective.

To generate an adaptive police patrolling location model, p-median location models employed in this paper exhibit three distinct characteristics. First, with road segments on a street network as observational units, the location models aim to delineate police patrol areas by comprising p sets of coterminous street segments. Second, vehicle crashes are used to indicate a risk level at each segment. While the counts of vehicle crashes at each segment can serve as a surrogate for traffic volume, such counts might not precisely depict the spatial distribution of crashes within a segment and also between adjacent segments. Hence, the models utilize density-based weights, particularly network kernel density count estimates (NKDCE) [11]. Third, because this delineation is intended for police patrol reasons, maintaining the contiguity of road segments in a delineated area is a crucial modeling concern because non-continuous or fragmented delineation results may raise efficiency issues in practical patrolling situations.

In addition, to examine how the aforementioned characteristics affect the delineation of police patrol areas, the results of two p-median location specifications, each coupled with two different sets of weights, are explored and compared. Specifically, the four analyses are set with two model specifications and two weight sets. The standard p-median location model is solved with two different weight sets attached to road segments: vehicle crash counts and network kernel density estimated vehicle crash counts. The standard model is also extended to consider a balancing constraint across the sizes of patrol distances. The extended p-median location model with balancing factors is also solved with the same two different weight sets. These balancing constraints allow the investigation of any differences in the coverages of patrol areas and the computation of delineation results that minimize the differences among patrol areas to avoid workload imbalances.

2. Literature Review

The relevant existing literature references the employment of location models for police patrol designing, their considerations for patrol area delineation, and the use of vehicle crash counts for their weights.

2.1. Location Models for Police Patrolling

The study of police patrolling location modeling encompasses two primary themes: the deployment networks for police cars and officers, as well as the delineation of police patrol areas. The early development of police deployment network models occurred in the 1970s. Larson [12] utilized the hypercube queuing model, which considers metrics such as the number of emergency service resources, arrival rates, service rates, and queue capacities in urban emergency (911 calls) dispatch services. This model enables the analysis of intricate queuing scenarios where multiple factors interact to influence overall system performance. Chaiken and Dormont [13] created a model for police car allocation, optimizing the assignment of police patrol resources. These models represented a significant advancement beyond hand calculations or simple computer algorithm output. Subsequently, these initial models have been extended to incorporate multiple dispatch queuing [14] and various queuing model variations [15,16].

Previous research has provided insights into deployment network models (

Table 1. Studies on police patrol area designs and optimization models.

Study	Model	Contribution
Mitchell [18]	p-median model	Designed patrol beats to reduce response time and balance workloads.
Vlćek et al. [20]	Contiguous p-median problem	Introduced preceding adjacent basic areas, ensuring compactness and contiguity in police service districting.
Larson [17]	p-median based approach, maximal covering location problem	Minimized travel distance for police service calls in New York City.
Curtin et al. [2,21]	Maximal covering location problem	Optimized patrol centers in Dallas, focusing on crime and service call priorities.
D’Amico et al. [22]	Heuristic approach (simulated annealing)	Partitioned police jurisdictions using a patrol car allocation model.

), yet there has been relatively little focus on partitioning and delineating police patrol areas in the literature [2]. On the one hand, Larson [17] developed an algorithm aimed at minimizing the total weighted travel distance for responding to service calls received by the New York City Police Department. Two optimization approaches commonly utilized in this context are the p-median location problem and the maximal covering location problem. Mitchell [18] introduced a p-median-based approach to designing police patrol beats, with the goal of reducing the average overall response time while highlighting disparities in beat workloads and response times. This approach demonstrated that the p-median location model yields a solution where the marginal improvement in weighted total (or average) cost or distance decreases monotonically with the addition of facilities [19]. Furthermore, the application of p-median objectives to districting problems has not been widely demonstrated despite their potential for broad applicability. There remains room for improvement in delineating police patrol areas, as shown by [20], who recently introduced the concept of preceding adjacent basic areas (BAs) to the p-median problem, ensuring compactness and contiguity in police service districting.

On the other hand, Ref. [2,21] introduced a maximal covering location problem (MCLP) to determine the optimal number of patrol centers required to cover their research area in Dallas, Texas. Their focus was on prioritizing crime incidents and service calls by leveraging the MCLP to solve the challenge of placing police cars in each area. Their objective was to maximize the number of incidents that could be responded to within a specified service distance or response time rather than minimizing weighted total (or average) cost as in the p-median location model. Additionally, they implemented a maximal backup coverage model with a multi-objective function to identify alternative patrol areas. D’Amico et al. [22] proposed a heuristic approach using simulated annealing to achieve a suitable partitioning of a police jurisdiction by adapting a variant of the patrol car allocation model and evaluating the effectiveness of district design.

As this preceding literature reveals, researchers have made significant progress in developing location models for police patrolling by integrating locationally optimized methodologies, with a particular emphasis on deployment networks and delineating patrol areas. However, these models often encompass a broad spectrum of events related to all police tasks rather than focusing on specific topics (i.e., vehicle crashes). There is still an opportunity to delve deeper into detailed topics within this field. From a practical perspective, by employing models tailored to specific topic goals, authorities in police departments can implement flexible plans by incorporating detailed strategies that target specific aspects of their tasks. In this context, the p-median approach is advantageous because it can generate optimized solutions that designate each patrol area to specialize in a specific task, potentially leading to cost minimization.

2.2. Model Considerations for Patrol Areas Delineation

The districting problem involves grouping preferably coterminous small geographic areal units into larger, somewhat homogeneous clusters known as districts or regions [23,24]. When addressing this problem, it is essential to consider three key criteria: contiguity, balance (especially workload), and compactness (Table 2), as outlined by [23]. A significant application of the districting problem is in delineating police patrol areas, which involves aggregating small geographical areas (e.g., road segments) within a police department’s jurisdiction into coherent patrol districts. For example, Zhang and Brown [25] emphasized the importance of contiguity and compactness in shaping police patrol areas, as well as reducing discrepancies in patrol area sizes to achieve workload balance. Furthermore, the choice of geographic units used as the building blocks for small geographic areas greatly influences the fulfillment of these criteria.

First, contiguity is a key geometric characteristic desired in the delineation of patrol areas, referring to the ability to travel within a district without exiting its boundaries and eliminating interstices under other police officers’ authority. However, incorporating contiguity constraints into location models presents challenges due to the sizable number of constraints it involves when solving its optimization problem, the sheer number possibly making solutions intractable [26]. To address this issue, contiguity can be defined using proximity grids, such as a Gabriel graph, a minimum spanning tree, or a Voronoi diagram, which provides a practical means of conceptualizing contiguity in location models [23].

Second, ensuring balanced workloads among districts is essential to prevent delays or failures in providing efficient and effective police patrol services [12]. Various criteria can be employed to quantify balance. One criterion is to consider the relative deviation of district sizes, such as the sum of distances within each district, from the mean district size [27,28,29]. This criterion allows for a measurement of the deviation from an ideal workload balance. Additionally, the range of district sizes can be used as a measure of balance [30]. Another criterion to consider is the volume of service calls. Because their response time depends on location, calls may require different amounts of service time; this criterion is closely linked to achieving balanced workloads [31]. Effectively dispatching service calls can play a significant role in enhancing workload balance. Curtin et al. [2,21] utilized this metric to define police patrol areas. By considering these criteria, location models can be prescribed for balanced workloads across delineated districts, thereby promoting efficient and equitable police patrol service areas.

Third, compactness is another crucial geometric consideration when delineating police patrol areas because it relates to the shapes of the formed districts. Ideally, districts should have round, undistorted shapes without any holes, resembling a circle, not an ellipse or a doughnut. Evaluating compactness commonly involves two local measures. One measure is the ratio of a district’s area to the area of its smallest circumscribing circle [32].

Table 2. Studies on districting problem considerations.

Criterion	Study	Methods for Addressing
Contiguity	Kalcsics [23]	Use of proximity grids with Gabriel graphs, minimum spanning trees, or Voronoi diagrams to define contiguity in models
Balanced workloads	Camacho-Collados et al. [31], Curtin et al. [21]	Use metrics with deviations of district sizes, ranges of district sizes, and service call volumes to define balanced workloads
Compactness	Garfinkel and Nemhauser [33], Young [32]	Measure compactness using the area-to-circle ratio, Schwarzberg compactness score, or elongation index to assess district shape

Table 2. Studies on districting problem considerations.

Criterion	Study	Methods for Addressing
Contiguity	Kalcsics [23]	Use of proximity grids with Gabriel graphs, minimum spanning trees, or Voronoi diagrams to define contiguity in models
Balanced workloads	Camacho-Collados et al. [31], Curtin et al. [21]	Use metrics with deviations of district sizes, ranges of district sizes, and service call volumes to define balanced workloads
Compactness	Garfinkel and Nemhauser [33], Young [32]	Measure compactness using the area-to-circle ratio, Schwarzberg compactness score, or elongation index to assess district shape

This provides an assessment of how closely a district’s shape approximates a circle, with a higher ratio indicating a more compact and circular shape. Another measure is the Schwarzberg compactness score, which calculates the ratio of a district’s perimeter length to the circumference of a circle with an equivalent area [33]. This score helps determine how well a district’s perimeter aligns with a circular shape. Additionally, Griffith et al. [34] analyzed the shape of geographical areas using the elongation index, calculated from the ratio of the orthogonal axes, to evaluate how circular a polygon is. By utilizing these three measures of compactness, the delineation of police patrol areas can create ideally shaped districts (e.g., regular hexagons or squares), avoiding distorted or irregular boundaries. This design promotes efficient patrol coverage and facilitates effective law enforcement efforts. However, for this study, assessing compactness presents challenges due to using road segments shaped by the road network rather than continuous space polygons as geographic units. Consequently, this paper primarily emphasizes contiguity and balance.

Additionally, when it comes to location modeling for police patrol areas, geographical units for the purpose of this paper can be defined generally in two ways: road segments and areal-based surfaces (e.g., beats, sectors, and grids). Previous studies have predominantly utilized road segments in their models because these units align with the practical understanding of police officers who primarily navigate and patrol on roads. For example, Rosser et al. [35] suggested that street-based districting models offer a better depiction and prediction of crime risks. Furthermore, criminology research conducted by [36,37,38] demonstrated that analyzing crime patterns using a street network base map is more beneficial due to its intuitive representation and understanding of patrol areas by sworn officers. Therefore, incorporating road network segments into models ensures more realistic results than outcomes using areal-based surfaces when addressing location model applications, specifically the p-median location models that help formulate patrolling areas.

Moreover, the literature on network-based spatial optimization for police patrol districting has recently advanced with several new approaches in mathematical modeling. Integrating road network segments into models yields more realistic and pragmatic results. For example, Chen et al. [39] incorporated street segments and factors such as crime risk rate and area size in a mixed-integer programming approach, and Kong et al. [40] employed center-based clustering for mixed-integer linear programming to consider balance workload, demonstrating the effectiveness of road-focused modeling. These approaches are often combined with how a shift toward data-driven, network-aligned patrol optimization, enhancing police responsiveness and adaptability across urban environments.

2.3. Vehicle Crash Counts as Weights for the Location Models

Vehicle crash counts can be suitable weights for a location model due to the strong correlation between prevention efforts and the occurrence and volume of accidents. Police services must provide reasonable service levels to ensure public safety and security as stipulated in the laws and practices of their jurisdictions. Traffic police fulfill two major functions: enforcing traffic laws and traffic accident/crime prevention [41]. On the one hand, prevention of these traffic troubles is carried out by increasing drivers’ awareness of the potential dangers in driving at different locations and times. On the other hand, police enforcement is an important factor in deterring driving offenses, reducing accidents, and thus improving road safety [4,5].

Weights are represented in two ways: aggregated vehicle crash counts or network kernel density vehicle crash count estimates (NKDCEs). These two measurements may generate different modeling solution outcomes. Vehicle crash counts are aggregated numbers of vehicle crashes based on road segment units. At the same time, NKDCEs are computed as a smoothed surface that takes into account the uncertainties of the crash locations. Each approach has different strengths and weaknesses for modeling purposes. First, vehicle crash counts have the advantage of reflecting the actual locations of vehicle crashes, which enables researchers to fine-tune and add to the sophistication that accompanies the utilization of location models. However, when formulating the models using road segments for a p-median model as a form of discrete optimization, it becomes necessary to aggregate crash locations based on these segments. While aggregating points into road segments, inherent uncertainties and biases arise from the use of spatially aggregated interacting masses [42,43]. For instance, the segmentation of road sections over a certain distance can result in a low variance in the number of crashes across these sections [44]. However, the methods to determine an appropriate segment width (i.e., size, starting/ending point for the segmentation process) remain unclear. According to [45], although an arbitrary spatial unit of analysis can be defined to create homogeneity across an entire region to facilitate comparisons and taxonomic exercises, its careful planning does not necessarily compensate for its arbitrary nature.

Second, NKDCEs provide an alternative approach when uncertainties in spatial network data are of concern, such as inaccuracies of events, particularly in representing the spatial distribution of vehicle crashes. NKDCEs have an advantage over the raw vehicle crash counts approach because they offer a smooth surface that reduces uncertainties in spatial patterns when visualizing and exploring hotspots [46,47] or modeling traffic accidents [48,49]. Unlike aggregated vehicle crash counts, using NKDCE mitigates the problem with the use of an arbitrary spatial unit and provides a nuanced understanding of crash patterns. The rationale for using NKDCEs instead of conventional Kernel Density Estimates (KDEs) lies in the limitations of applying KDEs to a road network. KDEs rely on Euclidean distance in a continuous, homogeneous, and isotropic two-dimensional (2-D) space, which is not suitable for networks [46,50]. Therefore, the use of NKDCE, which constrains the functional space of interest to a road network, is recommended to better reflect the distribution of vehicle crashes.

To sum up, despite the potential benefits of using spatial optimization models for police area delineation problems, especially when analyzing the geographic pattern of crashes on road segments in a transportation network, the p-median location model approach has not been widely employed in previous studies for designing them specifically for traffic safety. This paper aims to fill this gap by presenting p-median location models to construct optimal patrolling areas considering patrolling distances and vehicle crashes as weights.

3. Study Area and Methods

This paper seeks to compare models with different constraints and refinement methods for model weights related to vehicle crashes. To facilitate practical comparisons, the models are constructed using vehicle crash data and the road network in the City of Plano, Texas. These models are structured around the road network where vehicle crashes occur. The p-median location model is particularly well-suited for this task, as it helps maintain homogeneity, ensures compactness, and incorporates vehicle crash data. Weights are assigned based on the frequency of vehicle crashes, with the assumption that road segments with higher crash counts pose a greater risk and, therefore, require increased police patrolling. Figure 1 presents the process of model construction. Additionally, NKDE methods are adaptable to other geographic areas.

3.1. Study Area

The study area encompasses the City of Plano, located in the northeastern part of the Dallas-Fort Worth metropolitan area. The Plano Police Department has established eight safety goals for the city, including (i) achieving a low crime rate, (ii) ensuring safe travel, and (iii) providing timely responses [51]. This paper focuses on contributing to the goal of ensuring safe travel by utilizing models to delineate patrol areas based on vehicle crashes that have occurred on the city’s primary road network. It includes tollways, highways, and arterial roads, which account for approximately 92% of the total recorded crashes (see

Table 3. Major causes of vehicle crashes on the main roads in Plano, TX, from 1/2011 to 6/2021.

Cause	Total	Ratio
Intersection Related Crashes	16,955	51.33%
Distracted Driving Crashes	9815	29.71%
Work Zone Crashes	2157	6.53%
Alcohol Involved Crashes	1914	5.79%
Speed Related Crashes	1271	3.85%
Driving Under the Influence of Drugs	337	1.02%
Cell Phone Use Involved Crashes	320	0.97%
Motorcycle Related Crashes	262	0.79%
Total	33,031	100.00%

), data obtained from the Crash Report Information System (CRIS) of the Texas Department of Transportation, January 2011 to June 2021 [52]. During this period, a total of 33,031 crashes were reported. Among these incidents, intersection-related crashes accounted for 16,955 (51.3%), while crashes attributable to distracted driving constituted 9815 (29.7%) of the total accidents. Together, these two categories comprise approximately 80% of the overall crash occurrences during the studied time period.

Figure 2 displays the geographic distribution of aggregated vehicle crash counts across a total of 686 study-created road sections. These road sections were constructed as roughly 300 m intervals between road network intersection boundaries to optimize traffic flow and maintain simulation accuracy. According to [53], setting dual-lane road segments to be at least 300 m from an intersection and 200 m long results in significantly less deviation in traffic simulations. This choice ensures a balance between computational efficiency and the fidelity of the traffic models. Data manipulation for the p-median location modeling involved road segments at an intersection being merged into a single unit (‘+’ shape at a four-way intersection) under an assumption that all such segments are observable from the center of their intersection. This segmentation process supports solving the p-median problem in a discrete space using constructed weight values that are the number of crashes on each segment. Distances separating segments were calculated as shortest paths through the road network between their center locations (i.e., the midpoint of each segment). Note that the counts of vehicle crashes have a potential limitation as weights because unreported vehicle crashes may not be included and recording inaccuracy issues regarding GPS errors and incorrect address information plague the data.

3.2. The p-Median Location Model

The p-median location model uses crash counts that have occurred in each road segment and the distances between the midpoints of these road segments. The formulation of this p-median location model is adapted from [54,55] as an Integer Programming specification. Its formulation is as follows:

$$minimize\sum _{i=1}^n\sum _{j=1}^n{a}_i{d}_{ij}{x}_{ij}$$

(1)

Subject to:

$$\sum _{j=1}^n{x}_{ij}=1,\forall i$$

(2)

$$\sum _{j=1}^n{x}_{jj}=p$$

(3)

$$x_{ij}-x_{jj}\le 0\forall i,jandi\ne j$$

(4)

$$x_{ij}=\left\{\mathrm{0,1}\right\}\forall i,j$$

(5)

where

• $n$ denotes the number of road segments (n = 686),
• $i$ denotes an index of road segments as demand points,
• $j$ denotes an index of patrol areas, which are the centers of segments (j = 1, 2, …, n),
• $a_i$ denotes the number of vehicle crashes on road segment $i$,
• $d_{ij}$ denotes the shortest path distance between road segments $i$ and $j$,
• $x_{ij}$ is 1 if the $i$th road segment is assigned to the $j$th road segment as a patrol area center and 0 otherwise (i.e., the areas partition the road segments into mutually exclusive and collectively exhaustive groupings), and
• $p$ denotes the number of patrol areas to be delineated.

The goal of the model is to determine the allocation of road segments i to the center of a patrol area j considering the weights (a_i) and the distances separating road segments. Note that the center of a patrol area does not indicate where police cars will be stationed but rather the place from which all allocated segments are reachable with the lowest cost. The objective function (1) minimizes the sum of weighted distances, where the weights are based upon the number of crashes on road segment i. Constraints (2) are to ensure that each segment is allocated to one and only one patrol area (i.e., collectively exhaustive). Constraints (3) specify the number of patrol areas. Constraints (4) are the Balinski constraints to enforce that a demand point $i$ must have only one $j$ as the center of a patrol area to avoid any multiple allocations of a road segment to patrol areas (i.e., mutually exclusive). Finally, constraints (5) indicate a binary variable identifying which road segment is assigned to each patrol area. The p-median location model is solved using a commercial optimization solver, the IBM CPLEX 20.1, for a sequence of p values to ensure the identification of the most desirable optimal solution.

3.3. The Extended p-Median Location Model with Workload Balancing Constraints

The p-median model is extended here to reduce workload differences among the delineated patrol areas. This extension incorporates the constraint of capacity restrictions, which is specified by following Equation (6):

$${\sum }_{i=1}^n{d}_{ij}{x}_{ij}\le T,\forall j$$

(6)

Note that ${\sum }_{i=1}^n{d}_{ij}{x}_{ij}$ relates to the patrolling distance in the form of the sum of the network distance between a patrol area ($j$) center and assigned road segments ($i)$. Here, T represents a predetermined threshold that serves to cap the allowable patrolling distance. The inclusion of balancing constraints guarantees that the total patrolling distance within each delineated district remains below the set threshold T. This threshold stands as a maximal requirement, maintaining the connectivity of the delineated area for each value of p. Table S2 reports T values for each p instance studied in this paper.

3.4. NKDE: Estimation of Crash Counts for the p-Median Location Model

According to [56] (p. 183), the general formula of network kernel density estimation (NKDE) can be expressed as

$$\mathbb{K}\left(\gamma \right)={\displaystyle \frac{1}{n}}\sum {\mathbb{K}}_q\left(\gamma \right)$$

(7)

$${\mathbb{K}}_q\left(\gamma \right)=\left\{\begin{array}{cc}k\left(d_S\left(q,\gamma \right)\right)& \hspace{1em}\hspace{1em}\mathrm{f}\mathrm{o}\mathrm{r}0\le d_S\left(q,\gamma \right)<2d_s\left(q,v_i\right)-h\\ k\left(d_S\left(q,\gamma \right)\right)-{\displaystyle \frac{n_i-2}{n_i}}k(2{\mathrm{d}}_{\mathrm{S}}\left(q,v_i\right)-d_S\left(q,\gamma \right))& \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\mathrm{f}\mathrm{o}\mathrm{r}{2d}_S\left(q,v_i\right)-h\le d_S\left(q,\gamma \right)<d_S\left(q,v_i\right)\\ {\displaystyle \frac{2}{n_i}}k\left(d_S\left(q,\gamma \right)\right)& \mathrm{f}\mathrm{o}\mathrm{r}{d}_S\left(q,v_i\right)\le d_S\left(q,\gamma \right)<h\end{array}\right.$$

(8)

where

• ${\mathbb{K}}_q\left(\gamma \right)$ denotes the network kernel density estimator function at kernel center q on nondirected connected network $\tilde{L}$ in the spatial coverage of the NKDE,
• $\gamma$ denotes an arbitrary point on a subnetwork of $\tilde{L}$ that is defined as a network between two nodes,
• $k\left(d_S\left(q,\gamma \right)\right)$ denotes a base kernel density function using the shortest-path distance, $d_S\left(q,\gamma \right)$,
• $n_i$ denotes the total numbers of kernel centers in the buffer network of the ith node, $v_i$, in the set of nodes in $\tilde{L}$, and
• $h$ denotes the bandwidth that determines the range of the kernel function used to estimate the density.

NKDE is an extension of continuous surface KDE, a construct that is designed to estimate data within a network applicable to vehicle crashes on road systems [57]. Like KDE, NKDE approximates the intensity of a point pattern using a smoothed continuous surface that visualizes the variation of point event density but is computed only along a network segment [46]. In KDE, the density value is calculated based on observational points and a kernel weighting function that considers distance decay effects. However, traditional KDE in a 2D space is not suitable for vehicle crashes that only occur on a network of road segments. To address this, NKDE utilizes the non-Euclidean shortest-path distance $d_S\left(q,\gamma \right)$ along a road network [58].

Equation (7) represents the main NKDE formula, where $\mathbb{K}\left(\gamma \right)$ is the overall density estimate at a point, $\gamma$ on the network, calculated by averaging vehicle crashes across $n$ kernel centers. Each kernel center, $q$, contributes to calculating this density estimation through its own localized density function, ${\mathbb{K}}_q\left(\gamma \right)$, as defined in Equation (8). The formula for ${\mathbb{K}}_q\left(\gamma \right)$ is divided into three cases, based on distance from $q$: (1) For $0\le d_S\left(q,\gamma \right)<2d_s\left(q,v_i\right)-h$, where the distance is within the bandwidth range $d_S$. (2) For ${2d}_S\left(q,v_i\right)-h\le d_S\left(q,\gamma \right)<d_S\left(q,v_i\right)$, where the density function adapts, accounting for proximity to the network boundaries, and (3) For $d_S\left(q,v_i\right)\le d_S\left(q,\gamma \right)<h$, where the point $\gamma$ lies near or at the bandwidth edge, limiting further influence from the kernel center.

Figure 3 illustrates the spatial distribution of vehicle crash counts estimated using the NKDE with a bandwidth of 1000 m and a cell size of 200 m, processed with SANET 4.1 [57]. The parameter values are adapted from [53], considering a 300 m length for road segments. The 1000 m bandwidth covers three adjacent segments, while the 200 m cell size allows for an approximate overlapping of two adjacent segments. When compared to the crash counts map depicted in Figure 2, the NKDE map shows the smoothed weights surface. Specifically, the concentration of crashes at the intersections is spread out (i.e., smoothed) to their neighboring segments. This modification seeks to mitigate the segmentation issue, which can be related to the modifiable area unit problem and the aggregation of points. Specifically, it may reflect the distribution of crashes within a segment better than simple counts do. NKDEs provide more refined estimates of vehicle crash counts, resulting in continuous point patterns from the kernel function across a road network.

This capability helps reduce uncertainties associated with the aggregation process, unlike the aggregated observed crash counts by road segments Furthermore, NKDEs can help alleviate uncertainties and biases that arise during the somewhat arbitrary discretization process, particularly at intersections where roads intersect in multiple directions.

3.5. Model Results Efficiency Measurement

The p-median location model aims to minimize the objective function values. In general, the objective function values of optimal solutions decrease as p increases, because the distance from each demand point to its corresponding facility location gets shorter with an increasing number of facilities. However, it is important to note that an increase in p can also entail higher management costs, requiring more patrol vehicles to serve the increased number of delineated patrol areas. This characteristic often necessitates the measurement of efficiency for the system in operation because of p, serving as the basis for decision-making to determine an appropriate p number. Thus, an efficiency measure is used, which is defined as the decreasing rate of total objective function values between the current and previous p cases (i.e., between $p$ and $p-1$). Specifically, the efficiency of the two models is evaluated using the percentage change, ${\delta }_p$(%), calculated with the objective function values as follows:

$${\delta }_p(\%)={\displaystyle \frac{\left(ob{j}_p-ob{j}_{p-1}\right)}{ob{j}_{p-1}}}\times 100$$

(9)

where

• ${\delta }_p$ denotes the percentage change between two consecutive objective function values at $p$ and $p-1$,
• ${obj}_p$ denotes the value of the objective function at $p$, and
• $ob{j}_{p-1}$ denotes the value of the objective function at $p-1.$

4. Results

This section summarizes the assortment of p-median location problem solutions computed for this paper’s analysis.

4.1. Standard p-Median Location Model Findings

The p-median location model suitability was tested with consecutive integer values from p = 2 to 30. This range was carefully selected to thoroughly evaluate the impact of each model setting on the solutions. It covers a spectrum from simple scenarios to extreme cases, examining how varying the number of police patrol areas influences the global partitioning outcomes. Figure 4 presents the efficiency measure for the standard p-median location model defined in Section 3.4. When the ${\delta }_p$ value for a specific p shows a notable decrease compared with other ${\delta }_p$ values, then its p can be interpreted as a target for the number of patrol areas in terms of system-wide management. The p values between 3 and 5 exhibit the greatest efficiency gain (from −21.6% to −11.9%). However, despite these favorable ${\delta }_p$ values, a small number of patrol areas (i.e., small values of p) might not practically cover the entire study area in an adequate manner. Therefore, ideal numbers of patrol areas can be chosen from ${\delta }_p$ larger than −10%, excluding the smaller p = 3 to 5 cases. Among these, potential cutoff points to determine the ideal number of patrol areas can be considered at p = 6, 11, 14, and 25 (denoted as filled black circles in Figure 4) because their ${\delta }_p$ values are relatively lower than those for their subsequent $p$. Detailed model results appear in Table S1.

Nevertheless, the findings recounted here also disclose the potential for imbalanced workloads due to noticeable variations in the sizes of the delineated patrol areas. Excessively large patrol areas could lead to increased workloads, resulting in longer patrol hours and reduced patrol frequency within a given timeframe in practical scenarios. For instance, in Figure 5, particularly for instances with p values of 6, 11, 14, and 25, it is evident that the southeastern area of the city tends to feature a patrol district with the longest patrolling distance among all districts, potentially leading to disproportionately large patrol areas in this region. To address this issue and achieve a more equitable distribution of workloads, it is imperative to adjust the size of patrol areas in relation to their patrolling distances within this southeastern region. This adjustment would contribute to a more balanced and effective implementation of police patrolling strategies, although it still needs to be conditional upon the number of traffic accidents.

4.2. The Results of the Extended p-Median Location Model with Workload Balancing Constraints

Figure 6 displays the differences in patrol area sizes between two p-median location models: the original standard model and the extended model with balancing constraints. Detailed descriptive statistics of patrol area sizes are reported in Table S2. The ratios between the maximum and minimum patrolling distances (blue and orange lines) and their max-min differences (black bars) are used to assess the improvement in workload balance. For example, the difference is 1.43 (black) at p = 22, calculated by subtracting 5.94 (orange) from 7.37 (blue). A difference greater than zero indicates that the balancing constraints improve workload balance by reducing maximum patrol distances and/or increasing minimum patrol distances to minimize disparities. Noticeable reductions were observed for p = 22 (1.43), p = 20 (1.25), and p = 19 (0.73).

The ratios from both models increase as the number of patrol areas increases, indicating that disparities in patrol area sizes become larger with more patrol areas. However, the ratio from the extended model relieves the disparities because the workload balancing constraints reassign some road segments from the largest patrol areas to neighboring areas, achieving a more equitable but efficient distribution of patrol resources.

Additionally, the inclusion of balancing constraints does not result in a completely different arrangement of delineated patrol areas; instead, it leads to minor adjustments, primarily affecting the area with the maximal patrolling distance. Consequently, the rearranged patrol areas remain largely consistent with the optimal solutions of the standard model. This approach minimizes discrepancies in the sizes of each patrol area while maintaining geographical connectivity. A detailed description of this effect is provided in Section 4.4.

4.3. The Results of the Standard p-Median Location Model with NKDCE

The standard p-median location model using network kernel density counts estimates (NKDCE) was executed for p = 2 to 30. Detailed results are reported in Table S3. The model results using NKDCE yield larger objective function values compared to those with aggregated vehicle crash counts because the NKDCE method accounts for the densities of crash events after smoothing counts across the road network, thereby increasing certain counts. Changes in efficiency, ${\delta }_p$ are also observed, as shown in Figure 7. Noticeable increases in their objective function values are found for p = 6, 9, 11, 13, and 21, while the ${\delta }_p$ of the standard model with aggregated crash counts monotonically decreases (as shown in Figure 4). Notably, the NKDCE implementations suggest fewer patrol areas (13 instead of 14, and 21 instead of 25, respectively) and require lower costs compared to the preceding standard model.

Comparing Figure 5a,b with Figure 8a,b reveals absolute changes in spatial patterns. Different delineations affect the positions of the largest and smallest patrol areas at p = 6 and 11. This indicates that smoothed NKDCEs, which address uncertainties from arbitrarily segmenting roads, influence the spatial patterns in the optimal solutions.

When comparing the ratios between maximum and minimum values using the blue line in Figure 9 and two lines in Figure 6, the ratios of the standard model using NKDCE are consistently lower than models using aggregated vehicle crash counts, except in six instances (i.e., p = 4, 5, 10, 16, 17, and 20). Specifically, for p = 6, the ratio between the maximum and minimum distances is smaller with vehicle crash count NKDCE (2.02) than with aggregated vehicle crash counts (2.86). Similarly, for p = 11, the ratio with NKDCE (1.89) is lower than its counterpart (3.25). Additionally, for p = 13 and 21, the ratios with NKDCE are smaller than their counterparts, with values of 3.34 to 1.93 and 5.73 to 4.96, respectively.

4.4. The Results of the p-Median Model with NKDCE and Balancing Constraints

As the final extended model, we can consider the one incorporating NKDCE along with workload balancing constraints. This investigation allows us to assess how effectively the constraints induce any improvements in workload distribution. These improvements are observed for all instances except seven cases (p = 7, 16, 17, 18, 19, 20, and 26) as shown in Figure 9. Detailed results of delineated patrol areas appear in Table S4. Notably, for $p\ge 23$, differences in NKDCEs (black bars in Figure 9) are larger than those of the models using aggregated vehicle crash counts (black bars in Figure 6). Reflecting this outcome, the average max-min difference value for all p instances is 0.34 for models with NKDCE, and 0.18 for models with the aggregated vehicle crash counts. The difference between these two averages indicates that the extended model with NKDCE can largely improve workload distribution compared to the models using aggregated vehicle crash counts.

In summary, when comparing the results of the two implementations of the two model specifications, a pair of key features emerges based on whether aggregated vehicle crash counts or vehicle crash count NKDCE are used, and whether balancing constraints are applied. First, models formulated with NKDCE significantly decrease the mean values of patrolling distances for patrol areas. This result is due to the unbiased and smoothed nature of vehicle crash count NKDCE, which notably reduces maximum patrolling distances as an entailed effect. Second, the introduction of balancing constraints also results in smaller maximum patrolling distances for policing districts. There are no complete changes in spatial patterns, as reassignments occur only among neighboring road segments within each patrol area. However, this reduction is not uniform across all instances, with some cases showing no significant reduction. Consequently, the utilization of the NKDCE method appears to be a pivotal factor in effectively delineating police patrol areas while achieving improved workload balance.

For the p = 22 instance, both the standard model and the extended model using NKDCE as weights (Figure 10c,d) demonstrate a substantial decrease in the difference in patrol area sizes. Specifically, the max-to-min ratio is reduced from 4.83 to 3.69 compared to models using vehicle crash counts (Figure 10a,b), which show a reduction from 7.37 to 5.94. Additionally, the two extended specifications with balancing constraints (Figure 10b,d) show a slight decrease in the difference in patrol area sizes compared to the standard models (Figure 10a,c). This result stems from a decrease in the size of the largest area (78.3 km to 59.8 km of patrolling distances for models with NKDCEs and 137.6 km to 107.6 km for models with aggregated vehicle crash counts) after workload balancing constraints reassigned some road segments from the largest patrol areas to their neighboring areas. This observation underscores the equitable contribution of incorporating NKDCE and balancing constraints to achieve more equitable workloads, which is evident in the reduction in mean values when considering both NKDCE and balancing constraints.

In terms of computational complexity and solution times, however, the models incorporating NKDCE generally required more time to achieve optimal solutions across p instances, except for p = 21. The average solution times for p ranging from 2 to 30 were 137.68 s for the standard model with aggregated vehicle crash counts and 532.77 s for the standard model with NKDCE. The maximum solution time was 1484.28 s at p = 16 for the NKDCE under the machine (e.g., Intel^® Core™ i7-9700 CPU @ 3.0GHz, 16GB DDR4 RAM) and IBM CPLEX 20.1. These details are presented in Tables S2 and S4. The variation in solution times can be attributed to the differences in spatial patterns of vehicle crashes, as NKDCEs exhibit similar values in adjacent road segments, allowing the CPLEX to explore more candidate optimal solutions, generating larger numbers of iterations and increased computational complexity in solving MIP problems.

5. Discussion and Conclusions

This paper proposes a pair of police patrol area delineation models for Plano, TX, using the p-median location modeling approach. The utilized p-median model in our study is prescribed based upon road segments as observational areal units, with vehicle crash counts serving as weights, an aspect that has received little attention in the literature. This model is further extended by introducing workload balancing constraints. These two p-median models are also examined with network kernel density count estimates (NKDCE) for spatial smoothing to evaluate their effectiveness as an alternative weights option. These specifications and implementations lead to different results in cost efficiency and geographic coverage for each patrol area. The selection of these options allows for the creation of adaptable models that align with different patrolling scenarios, taking into account the specific features highlighted in this paper when executing a spatial location model.

The findings demonstrate the effectiveness of the p-median location problem approach specifically for delineating police patrol areas in the context of vehicle crashes. The common advantage of the p-median among location models is its ability to encompass patrol areas that cover an entire study region, a capability not shared by other location models such as those for the maximal covering or location set covering problem. These models leverage vehicle crash counts as weights and shortest-path distances on road networks to optimize patrol area allocations. Moreover, the inclusion of balancing constraints allows for a more equitable distribution of workloads among patrol areas, enhancing overall efficiency.

Additionally, refining vehicle crash weights by introducing network kernel density estimation (NKDE) results in significant improvements in delineating patrol areas. This refinement considers the uncertainties associated with vehicle crash locations and contributes to a more balanced workload for police patrolling operations. The smoothed surface estimates lead to a more balanced distribution of workload for police patrolling activities. The combination of p-median modeling with NKDCE furnishes valuable options for optimizing police patrol area delineations and enhancing workload balance for law enforcement efforts. However, future studies exploring models with different parameters for the NKDE would be illuminating.

The principal contribution of this study is twofold. First, it outlines a procedure for delineating Plano police patrol areas using the $p$-median location problem approach. This strategy can provide a flexible methodology to improve the delineation problems in a discretized network space. Second, the p-median problem model is examined with NKDCE so that a spatial pattern of point events on a network space is approximated and reflected in smoothed weights. These smoothed values may mitigate a potential misrepresentation of a spatial pattern caused by the somewhat arbitrary discretization of observation units. That is, it can better represent crashes as weights and avoid a potential bias resulting from irregular distributions of vehicle crashes across road network segments.

However, there is a limitation in model comparison due to the different weight sets used, despite the improvement from arbitrary discretization. The objective function values in this experiment depend on these specific weight sets. Consequently, while both of the two considered model specifications provide optimal results with their respective weight sets, a direct comparison based on numerical measures may not be feasible. This limitation, common in many optimization problems, arises because optimal solutions are conditional on their input values. Another limitation is that the model cannot consider the practical methods employed by the Plano Police Department because they declined to share open data regarding police patrolling. The Plano Police Department cited the strong connection between this information and their operational strategies aimed at ensuring citizen safety.

This study intends to underscore the contributions of the p-median model, particularly when it is combined with supportive elements, such as balancing constraints and network kernel density estimates, which have not been adequately addressed in the current literature. This study emphasizes the contribution of the p-median model in a practical context rather than from a purely theoretical modeling perspective.

However, this research is limited in that it lacks diverse evaluation metrics. Future research can address this issue with cross-comparisons with alternative models to assess the strengths and weaknesses of each.

Also, exploring a model that incorporates service calls and/or response times, along with resources such as the number of patrol cars and police officers in Plano, could enhance its practical applicability. Additionally, integrating temporal dynamic components such as time of day and days of the week could improve alignment with real operational practices.