A Method for Locational Risk Estimation of Vehicle–Children Accidents Considering Children’s Travel Purposes

The reduction in locational traffic accident risks through appropriate traffic safety management is important to support, maintain, and improve children’s safe and independent mobility. This study proposes and verifies a method to evaluate the risk of elementary school students-vehicle accidents (ESSVAs) at individual intersections on residential roads in Toyohashi city, Japan, considering the difference in travel purposes (i.e., school commuting purpose; SCP or non-school commuting purpose: NSCP), based on a statistical regression model and Empirical Bayes (EB) estimation. The results showed that the ESSVA risk of children’s travel in SCP is lower than that in NSCP, and not only ESSVAs in SCP but also most ESSVAs in NSCP occurred on or near the designated school routes. Therefore, it would make sense to implement traffic safety management and measures focusing on school routes. It was also found that the locational ESSVA risk structure is different depending on whether the purpose of the children’s travels is SCP or NSCP in the statistical model. Finally, it was suggested that evaluation of locational ESSVA risks based on the EB estimation is useful for efficiently extracting locations where traffic safety measures should be implemented compared to that only based on the number of accidents in the past.


Introduction
Road traffic accidents are one of the major causes of death and disability among children over the world [1]. The threat of being involved in a traffic accident contains not only the risk of direct physical damage but also the long-term impact on children's lifestyles. Namely, the greater the concern parents/guardians feel about traffic accidents involving their children, the less likely they are to allow their children to play freely outside or to go out independently without adults [2][3][4][5][6]. In fact, it has been reported that children's opportunities to walk and cycle, and play outside, as well as children's independent mobility (CIM), have decreased in the 21st century in some countries [7][8][9][10][11][12]. The decrease in daily walking and cycling opportunities leads to children's lack of physical activity [7]. Children's opportunities for physical activity are important for their physical maturation and the form of normal bone and behavior, and physical inactivity is widely recognized as one of the major risk factors in human health [8,13]. It has been argued that CIM is also important in fostering children's autonomy and in curbing car-independent transportation style which reduces time lost by parents taking their children as well as traffic congestion [10]. Therefore, the reduction in traffic accident risks and threats through appropriate traffic safety management is very important in urban planning to support, maintain, and improve children's safe and independent mobility (CSIM).
In traffic safety management, it is important to properly assess the risk of traffic accidents at each area or location, especially for identifying areas/locations where countermeasures should be implemented and verifying the effectiveness of those countermeasures. In previous studies on the evaluation of accident risk by areas/locations, statistical model analyses (so-called safety performance functions: SPF) have often been used to clarify the effects of various explanatory variables on the number of accidents [14]. Table 1 lists the main previous studies using a statistical modeling approach for the assessment of accident risk involving children by areas/locations [15][16][17][18][19][20]. The aggregation units for accident risk assessment are diverse, ranging from the microscopic level such as intersections and road segments to the macroscopic level such as prefectures. While these previous studies have made important contributions to the methodology for assessments of accident risks involving children by areas/locations, there have been issues to be dealt with. The first is the consideration of the children's travel purposes. The risk of children walking and cycling being involved in accidents may be different between traveling for school commuting purposes (SCP) and traveling for non-school commuting purposes (NSCP). In Japan, children are expected to be safer on the road in SCP than in NSCP, because many elementary schools have adopted a group walking system, called "Shudantogeko", to and from schools, and parents/guardians and district communities provide generous safety support [21,22]. In recent years, several countries, including the United Kingdom, also have begun to implement group walking initiatives, called "Walking bus", to school [21]. Thus, accident risk assessments that take into account the differences between children's travel in SCP and in NSCP are internationally important. Inada et al. [20] analyzed traffic accidents involving elementary school students and junior high school students by prefecture in Japan, dividing the accidents into those involving the students in SCP and NSCP, and found that the number of killed or seriously injured (KSI) cases during travels in SCP was approximately 30-40% lower than that during travels in NSCP. However, the study was conducted only at the macro level (i.e., prefecture-level), and, to the best of the author's knowledge, there have been no studies in which whether the children's purpose was SCP or NCSP was taken into account in a micro level (such as intersections or road segments level) analysis.
As the second issue, the accident risk of each location may be affected by a great variety of factors specific to the location. It is difficult to incorporate all of the factors as explanatory variables in a statistical model, and the effects of factors other than the explanatory variables are represented as the error term with a probability distribution. Therefore, in assessing accident risks at different locations, it is necessary to properly consider the balance between the average effects of explanatory variables and the effect of location-specific characteristics. While the Empirical Bayesian (EB) estimation method has been proposed and used in the context of a before-after study on the evaluation of traffic safety measures [23][24][25][26][27][28][29][30], it is also applicable to the purpose of extracting locations where countermeasures should be implemented by properly evaluating the accident risks at individual locations [31][32][33][34][35][36][37]. Although there have been some previous studies applying the EB method to evaluate the locational risks of accidents involving general pedestrian [30,33], no studies applying the EB method to assess the locational risks of accidents involving children have been seen.
Therefore, this study proposes and verifies a method to evaluate the risk of vehiclechildren accidents at individual intersections on residential roads in Toyohashi city, Japan, considering the difference in travel purposes (i.e., SCP or NSCP), based on a statistical regression model and EB estimation.

Target Area
The target area was the entire Toyohashi City, Aichi Prefecture, Japan. There are 52 elementary schools in the city (see Figure 1), and the number of elementary school students is approximately 20,000 as of 2018. Almost all the students walk to school, except for about 10 students who cycle to school.

Target Area
The target area was the entire Toyohashi City, Aichi Prefecture, Japan. There are 52 elementary schools in the city (see Figure 1), and the number of elementary school students is approximately 20,000 as of 2018. Almost all the students walk to school, except for about 10 students who cycle to school.

Traffic Accident Data
The accident data used in this study were police-reported elementary school studentvehicle accidents (ESSVAs: injury accidents in which the first or second party was an elementary school student while walking or cycling) that occurred in Toyohashi City over a 12-year period from 2009 to 2020 (totally 638 ESSVAs). The individual accident data include the children's travel purpose (i.e., SCP or NSCP). Figure 1 shows the distribution of the locations of the ESSVAs by the children's travel purpose.
The ESSVAs which occurred within a radius of 50 m from the center of each community street intersection were extracted as "intersection ESSVAs", which were 360 ESSVAs (56.4% of all the ESSVAs). The community street intersections were defined as intersections with three or more legs of community streets other than "national expressways", "general national highways", "major regional roads", "general prefectural roads", and "major general roads" in the road network data of Zmap-Area II (ZENRIN corp.). Of the 7719 target intersections, one ESSVA occurred at 296 (3.8%) intersections, two ESSVAs at 24 (0.31%) intersections, three ESSVAs at four (0.052%) intersections, and four ESSVAs at one (0.013%) intersection, during the 12-year period.
We used the 299 ESSVAs in the first nine years as the "risk estimation period" to estimate the ESSVA risk at each intersection, and the 61 ESSVAs in the last three years as the "evaluation period" to evaluate the results of the estimation.

Traffic Accident Data
The accident data used in this study were police-reported elementary school studentvehicle accidents (ESSVAs: injury accidents in which the first or second party was an elementary school student while walking or cycling) that occurred in Toyohashi City over a 12-year period from 2009 to 2020 (totally 638 ESSVAs). The individual accident data include the children's travel purpose (i.e., SCP or NSCP). Figure 1 shows the distribution of the locations of the ESSVAs by the children's travel purpose.
The ESSVAs which occurred within a radius of 50 m from the center of each community street intersection were extracted as "intersection ESSVAs", which were 360 ESSVAs (56.4% of all the ESSVAs). The community street intersections were defined as intersections with three or more legs of community streets other than "national expressways", "general national highways", "major regional roads", "general prefectural roads", and "major general roads" in the road network data of Zmap-Area II (ZENRIN corp.). Of the 7719 target intersections, one ESSVA occurred at 296 (3.8%) intersections, two ESSVAs at 24 (0.31%) intersections, three ESSVAs at four (0.052%) intersections, and four ESSVAs at one (0.013%) intersection, during the 12-year period.
We used the 299 ESSVAs in the first nine years as the "risk estimation period" to estimate the ESSVA risk at each intersection, and the 61 ESSVAs in the last three years as the "evaluation period" to evaluate the results of the estimation.

Risk Exposure Data
Accident risk exposure is an important factor for the statistical model analysis of locational accident risks, and some studies have made efforts to measure exposure for evaluating the locational risk of accidents involving pedestrians [38][39][40][41]. In this study, the exposure should be the number of passings or activities of elementary school students at each intersection. However, it is difficult to measure these data over the wide city area. As shown in Table 1, the previous studies have tried to estimate exposure by children population [18,20], by questionnaire surveys asking about travel routes [16,17], by an observation survey of travel mode to schools [19], and by a journey allocation model based on shortest routes [16].
On the other hand, the author and the Toyohashi City Board of Education have collaboratively georeferenced the designated school route data (see Figure 1) and the gathering point data, with the number of elementary school students for each school commuting group in the city. In this study, based on the designated school route data (2018 version), we calculated the following indices as surrogate variables for ESSVA risk exposure at the intersections [42]: Note that there are no major differences in the designated school route data in recent years.

Other Data
The ESSVA risk at an intersection is also may be greatly affected by the vehicle traffic volume [16,43,44]. Therefore, in order to take into account the effect of vehicle traffic volume at each target intersection, this study used probe vehicle data corrected by Pioneer Corp. for the years 2016-2019. The validity of the probe vehicle data has been confirmed in our previous studies [30]. The probe vehicle traffic volumes of all legs in the direction of inflow to the intersection were summed up to obtain the probe vehicle traffic volume for each target intersection, which was used as a variable to represent the vehicle traffic volume.
As other explanatory variables for road traffic environment conditions, we used the number of legs for each intersection, the presence or absence of traffic signals, the urban area land use type by 1/10 mesh (100 m mesh) of the National Land Information, whether each intersection is within the densely inhabited district (DID), and the distance to the nearest park.

Statistical Model
In this study, the locational accident risk is defined as the expected value of the number of ESSVAs by intersections in a specific period. Therefore, a count data model is used as the statistical model. The Poisson regression model is the simplest count data model. However, because the error structure (i.e., Poisson distribution) has a strong limitation that the variance equals the expected value, the negative binomial (NB) regression model, which is one of the extensions of the Poisson regression model, is commonly applied in recent accident count data modeling [14]. Additionally, in this study, the NB regression model is used.
To the best of the author's knowledge, the EB estimation has been found in the field of traffic accident risk assessment since Gipps [23], Abbess et al. [24], and Jarrett et al. [31] in the early 1980s, and later it was well-organized by Hauer [25] in the context of before/after study of traffic safety measures. In order to summarize the method for assessing accident risk considering location-specific factors using the EB estimation, it is first explained that the Poisson regression model can be extended to the NB regression model by assuming the property that the expected value in the Poisson regression model itself varies stochastically according to the gamma distribution. Then, the EB method is presented for the evaluation of the risk of accidents at each location, taking into account the location-specific effects that cannot be explained by the average effects (fixed effects) of the explanatory variables considered in the model.

Poisson Regression Model
Suppose that the number of accidents Y i at location i during a certain period of time follows a Poisson distribution with an expected value µ i that is a function of the explanatory variable vector x i and the corresponding parameter vector b: Thus, the probability that Y i accidents occur at location i during the period is and the expected value and variance of Y i are where

Accident Risk Evaluation That Considers Location-Specific Factors
Assuming that the expected value of the Poisson regression model described above is not completely determined only by the function of explanatory variables and parameters, but is also affected by factors specific to location i, and then assuming as an alternative to the expected value µ i the random variable λ i which follows a gamma distribution with the shape parameter α = φ and the rate parameter β = φµ −1 i : The probability that λ i be a certain value r i is and the expected value and variance of λ i are The structure that the expected value of the number of accidents λ i is determined following the gamma distribution, and the realization y i of the number of accidents Y i occurs following the Poisson distribution with λ i as the parameter can be interpreted as a hierarchical model. Therefore, under the information that y i accidents actually occurred at location i, the probability that λ i is a certain value r i can be obtained by Bayes' theorem: Here, µ i and φ are the prior parameters. Furthermore, the gamma distribution is the natural conjugate distribution of the Poisson distribution, thus it is known that the posterior distribution is also gamma, so that: and its expected value is The expected value of λ i can be used as an accident risk evaluation index that takes into account location-specific factors.

EB Estimation through an NB Regression Model
The accident risk evaluation index in Equation (11) depends on the prior parameters µ i and φ. Therefore, in this study, these prior parameters are estimated by the EB method. That is, we use µ i and φ which maximize the marginal likelihood: Expanding Equation (12) further and setting p i = φ φ+µ i , we obtain: indicating that Y i follows an NB distribution with parameters µ i and φ: The expected value and variance of Y i are Here, φ is called the dispersion parameter. The larger φ, the smaller the variance (residual deviation) of Y i . When φ = ∞, the variance is equal to the expected value, which means that the NB distribution is consistent with the Poisson distribution. From the above, the EB estimation of λ i is possible by estimating the NB regression model using the actual number of accidents y i by locations and the explanatory variable vector x i , and then by using the variance parameter φ NB and the expected number of accidents by location µ NB i obtained from the estimated NB regression model as prior parameters. Denoting the EB estimator of λ i as λ EB i , it can be shown that, Therefore, it can be found that this EB estimator is a weighted average of the two variables: (1) the expected number of accidents µ NB i , which is calculated based on the average effects of the explanatory variables in the estimated NB regression model; and (2) the number of accidents y i that actually occurred. It can be also found that the larger the dispersion parameter φ NB of the NB regression model is (i.e., the smaller the residual deviation by the NB regression model is), the greater the weight is given to µ NB i . According to Ghosh and Meeden [45], "one of the main features of EB analysis is to borrow strength from the ensemble; that is, use information from similar sources in constructing estimators and predictors in addition to the most directly available source of information." Therefore, this study first estimated the parameter b k and the dispersion parameter φ NB of the NB regression model shown in Equations (4) and (14), using the actual number of ESSVAs at each intersection i (explained in Section 2.2.1) as the objective variable y i and the risk exposure indices (explained in Section 2.2.2) and other road traffic environment conditions (explained in Section 2.2.3) as explanatory variables x ik . Next, using the explanatory variables x ik and the estimated model parameter b k , the NB estimate of the expected number of ESSVAs at each intersection, µ NB i , was calculated by Equation (4). Then, using the NB estimate of the expected number of ESSVAs, the actual number of ESSVAs y i , and the parameter φ NB , the EB estimate of the expected number of ESSVAs at each intersection, λ EB i , was calculated by Equation (17). This EB estimate of the expected number of ESSVAs was the risk evaluation value for each intersection.

Evaluation of Efficiency to Extract Risky Locations
To evaluate the usefulness of the EB estimation from the view of efficient extraction of high ESSVA risk locations, we calculated efficiency criteria by the following procedure. The concept of efficiency here is that "the higher the number of ESSVAs in the evaluation period per one extracted ESSVA risk intersection, the higher the efficiency".

1.
Each location was ranked in order of riskiness based on the three criteria, respectively: a. Actual number criterion, i.e., the higher the actual number of ESSVAs y i in the risk estimation period (the first nine years), the higher the risk is; b.
NB estimate criterion, i.e., the higher the expected number of ESSVAs µ NB i calculated by the NB regression model, the higher the risk is; c.
EB estimate criterion, i.e., the higher the expected number of ESSVAs λ EB i corrected by the EB estimation method, the higher the risk is.

2.
Based on the rankings of the three criteria, the cumulative sums of the actual number of ESSVAs in the evaluation period (the last three years) up to the rank were, respectively calculated.

3.
The cumulative sums of the actual number of ESSVAs in the evaluation period for the three criteria were, respectively divided by the number of intersections up to the rank, which is the efficiency in this study.  Figure 3 shows the distribution of the occurrence time of ESSVAs by the children's travel purposes, indicating that many of the ESSVAs occurred in the time period of 4:00-5:00 p.m. in NSCP. The high number of ESSVAs occurring in NSCP at this time means that elementary school students are mostly involved in an accident when being out to play or go somewhere (e.g., friends' houses or shops) after returning home on weekdays or on holidays.  Figure 4 shows the distribution of the occurrence of ESSVAs regarding the distance from the nearest school route by the children's travel purposes. As expected, most (roughly 80%) of the ESSVAs in SCP occurred on the school routes since children basically walk to school and walk back home along the designated school routes. On the other hand, it is a surprise that most ESSVAs in NSCP also occurred on or close to (roughly within 100 m of) the school routes. This could be an important fact and will be discussed, along with the statistical model estimation results below, in detail in Section 4.   Figure 3 shows the distribution of the occurrence time of ESSVAs by the children's travel purposes, indicating that many of the ESSVAs occurred in the time period of 4:00-5:00 p.m. in NSCP. The high number of ESSVAs occurring in NSCP at this time means that elementary school students are mostly involved in an accident when being out to play or go somewhere (e.g., friends' houses or shops) after returning home on weekdays or on holidays.  Figure 4 shows the distribution of the occurrence of ESSVAs regarding the distance from the nearest school route by the children's travel purposes. As expected, most (roughly 80%) of the ESSVAs in SCP occurred on the school routes since children basically walk to school and walk back home along the designated school routes. On the other hand, it is a surprise that most ESSVAs in NSCP also occurred on or close to (roughly within 100 m of) the school routes. This could be an important fact and will be discussed, along with the statistical model estimation results below, in detail in Section 4.  Figure 4 shows the distribution of the occurrence of ESSVAs regarding the distance from the nearest school route by the children's travel purposes. As expected, most (roughly 80%) of the ESSVAs in SCP occurred on the school routes since children basically walk to school and walk back home along the designated school routes. On the other hand, it is a surprise that most ESSVAs in NSCP also occurred on or close to (roughly within 100 m of) the school routes. This could be an important fact and will be discussed, along with the statistical model estimation results below, in detail in Section 4.  Table 2 shows the estimation results of two types of final NB regression models using the numbers of ESSVAs in SCP and the numbers of ESSVAs in NSCP, at the target intersections in the risk estimation period as the objective variable. They indicated that the statistically significant explanatory variables were different between the SCP model and the NSCP model. The parameter of the constant term in the SCP model was strongly lesser than those in the other models, indicating that children's risk to be involved in an accident during school commuting is very lower than during other travels.  Table 2 shows the estimation results of two types of final NB regression models using the numbers of ESSVAs in SCP and the numbers of ESSVAs in NSCP, at the target intersections in the risk estimation period as the objective variable. They indicated that the statistically significant explanatory variables were different between the SCP model and the NSCP model. The parameter of the constant term in the SCP model was strongly lesser than those in the other models, indicating that children's risk to be involved in an accident during school commuting is very lower than during other travels.

Estimation Results of the NB Regression Model
Regarding the surrogate variables for the ESSVA risk exposure, the more the number of elementary school students on the designated school route passing within the intersections (NSSR15), the higher the ESSVA risk was in the SCP model. In the NSCP model, on the other hand, the nearer the distance from the intersections to the nearest designated school route, the higher the ESSVA risk was. These results suggest that the designated school route data could be useful as surrogate variables for ESSVA risk exposure indicators.
As for the other explanatory variables, the main findings were as follows: • ESSVA risk increased with vehicle traffic volume in both models, but the impact was less in the SCP model than in the NSCP model; • ESSVA risk was higher at four-leg intersections than at three-leg intersections in the NCSP model, but there was no such effect of road structure in the SCP model; • The SCP model showed a lower ESSVA risk at intersections near parks, but the NCSP model did not show such an effect; • ESSVA risk was higher in areas with building use in both models, but was particularly strong in the SCP model; • Intersections within the DID had lower ESSVA risk in the SCP model, but the effect is not shown in the NSCP model. Figure 5 shows the scatter charts plotting the ESSVA risks calculated by the NB regression model (µ NB i ) on the horizontal axes and those corrected by the EB estimation (λ EB i ) on the vertical axes, according to the number of ESSVAs actually occurred during the risk estimation period shown at the top of the charts (0-3). The 45-degree line is also plotted. It can be seen that in the chart for the intersections where the actual numbers of ESSVAs were zero (the leftmost chart in Figure 5), all points were below the 45-degree line. This means that the expected numbers of ESSVAs were adjusted downward through the EB estimation because the actual numbers of ESSVAs had been smaller than the expected numbers estimated by the NB model. On the other hand, in the chart for the intersections where the actual numbers of ESSVAs were two or three (the right two charts in Figure 5), all points were above the 45-degree line. This means that the expected numbers of ESSVAs were adjusted upward through the EB estimation because the actual numbers of ESSVAs had been greater than the expected numbers estimated by the NB model. Additionally, in the chart for the intersections where the actual numbers of ESSVAs were one (the second left chart in Figure 5), the downward and upward adjustments described above were mixed. It can also be seen that even for the intersections with the same actual numbers of ESSVAs, the further to the right, the greater the difference between each point and the 45-degree line tended to be. This means that the higher the ESSVA risks calculated by the NB regression model were, the larger the magnitude of the adjustment by the EB estimation was. Figure 6 shows the efficiencies calculated by the procedure explained in Section 2.4 based on the three criteria. It can be seen that the actual number criterion had the highest efficiency up to the 6th ESSVA risk rank, whereas the EB estimate criterion had the highest efficiency from the 7th up to about the 60th risk rank. This means that, in this study case, while using the actual number of ESSVAs was more efficient in extracting the top six risky locations, using the EB estimation was more efficient in extracting the risky locations beyond these six locations.  Figure 6 shows the efficiencies calculated by the procedure explained in Section 2.4 based on the three criteria. It can be seen that the actual number criterion had the highest efficiency up to the 6th ESSVA risk rank, whereas the EB estimate criterion had the highest efficiency from the 7th up to about the 60th risk rank. This means that, in this study case, while using the actual number of ESSVAs was more efficient in extracting the top six risky locations, using the EB estimation was more efficient in extracting the risky locations beyond these six locations.

Discussion and Conclusions
The results of the aggregation analyses and the statistical model analysis showed that ESSVAs in SCP are less frequent than those in NSCP. This tendency is similar throughout Japan [20], and this may be due in large part to the effect of the cultural safety systems such as group school commuting (called "shudan-togeko") and the safety supports such   Figure 6 shows the efficiencies calculated by the procedure explained in Section 2.4 based on the three criteria. It can be seen that the actual number criterion had the highest efficiency up to the 6th ESSVA risk rank, whereas the EB estimate criterion had the highest efficiency from the 7th up to about the 60th risk rank. This means that, in this study case, while using the actual number of ESSVAs was more efficient in extracting the top six risky locations, using the EB estimation was more efficient in extracting the risky locations beyond these six locations.

Discussion and Conclusions
The results of the aggregation analyses and the statistical model analysis showed that ESSVAs in SCP are less frequent than those in NSCP. This tendency is similar throughout Japan [20], and this may be due in large part to the effect of the cultural safety systems such as group school commuting (called "shudan-togeko") and the safety supports such

Discussion and Conclusions
The results of the aggregation analyses and the statistical model analysis showed that ESSVAs in SCP are less frequent than those in NSCP. This tendency is similar throughout Japan [20], and this may be due in large part to the effect of the cultural safety systems such as group school commuting (called "shudan-togeko") and the safety supports such as crossing guard by district communities when traveling in SCP, as suggested by Waygood et al. [21] and Matsuo et al. [22] in macroscopic levels. Rothman et al. [19] showed a different result that intersections with school crossing guards had higher collision risk. However, there is a possibility that the locations where school crossing guards are located reflect the areas where many children walk, namely the risk exposure is higher resulting in higher accident frequencies. The result in this study is considered reasonable because the amount of risk exposure in SCP used in this is valid as described below.
The results of the aggregation analysis revealed that not only ESSVAs in SCP but also most ESSVAs in NSCP occurred on or near the designated school routes. The statistical model analysis showed that the number of elementary students commuting to school on each designated school route passing within a radius of 15 m from the center of each target intersection (NSSR15) was significant in the SCP model. This can be interpreted as the number of children passing through the designated school routes functioning as a risk exposure variable for ESSVAs in SCP. This is precisely because the number of children commuting to school on the school route itself is appropriate as exposure for ESSVA risk in SCP. Lee et al. [17] showed similar results using the number of children on the route to school based on a questionnaire survey as an exposure variable. On the other hand, in the NSCP model, the distance from each intersection to the nearest school route (DSR) was significant, indicating that the structure regarding ESSVA risk exposure differs depending on the children's travel purposes. This result suggests that elementary school students tend to use the designated school routes even for NSCP because they are familiar with the school routes due to daily use in SCP, resulting that the number of ESSVAs on or near school routes would be higher. Therefore, even though the fact described above is that the majority of the ESSVAs occurred while walking and cycling for NSCP, it would make sense to implement traffic safety management and measures focusing on school routes as the Japanese national and regional government currently makes efforts [46].
As Bennet and Yiannakoulias [16] showed, ESSVA risk increased with vehicle traffic volume. However, the impact was less in the SCP model than in the NSCP model. This may be because there is a cultural system that contributes to the safety of children traveling to and from school in Japan, as mentioned above, and also because crossing guards are present especially in locations with heavy vehicle traffic. On the other hand, it suggests that traffic calming is more important for the safety of children traveling in NCSP. Many previous studies showed the safety effect of traffic calming devices on residential roads [47,48].
Although some previous studies used the number of legs of intersections but could not find a significant effect [17,18], this study found that ESSVA risk was higher at fourleg intersections than at three-leg intersections for NCSP travels. While the result that ESSVA risk was higher in areas with building use in both models is similar to the previous studies [15][16][17]19], this study additionally found that the effect was particularly strong for SCP travels. Inada et al. [20] showed that locations of DID have the effect of an increase in children's accident risk through the prefectural level macroscopic analyses; however, this study revealed that the intersections within the DID have lower ESSVA risk for SCP travels. This may be because the size of group commuting and safety support by district communities, as described above, are more extensive in densely populated areas. As mentioned so far, the locational ESSVA risk structure, including the exposure, varies depending on whether the purpose of the children's travels is SCP or NSCP, which is considered an important finding of this study.
Regarding the EB estimation results, it was shown that in this study case, using the EB estimate criterion was more efficient when extracting more than six risky locations. This suggested that the evaluation of locational accident risks based on the EB estimation is useful also for efficiently extracting locations where traffic safety measures for ESSVAs should be implemented compared to that only based on the number of accidents in the past, as argued in the previous studies applying EB method for other types of accidents [33][34][35][36][37]. In practice, it would be better to use both the number of accidents and the EB estimation method to identify risky locations.
Finally, this study has some limitations and further issues. The first is that the exposure data for locational ESSVA risks in NSCP were insufficient. While we used designated school route data as a surrogate variable for exposure in this study, the spatial distribution of children's travels and activities in NSCP should be measured and used in some way (e.g., questionnaire/diary survey [15,17], observation survey [19], GPS logger survey [49], etc.).
The second is to be able to assess in detail the safety effects of group commuting and safety support by district communities. For this purpose, it is necessary to make comparisons between districts with different approaches to safety support, and between elementary school districts with and without group commuting.
Thirdly, although the same statistical modeling and predicting approach (i.e., NB regression model and EB estimation method) can be used for other regions and countries, the road traffic environment conditions and data availability differ among regions or countries. Therefore, it is necessary to consider an objective variable and a variety of explanatory variables including risk exposure indicators by tailoring them to those conditions.