Although the traditional kernel method has mostly assumed the underlying data is continuous, in many cases, especially in transportation problems, the data might include categorical or binary attributes. Traditionally, in case of having both continuous and categorical data, the frequency approach would be implemented, where the data would be categorized into cells being assumed for categorical data, and then the density approach would be implemented on the continuous attributes in each cell.

However, the frequency approach is expected to perform unsatisfactorily as it would result in efficiency loss due to the employment of the sample splitting method, especially in the case of having a large sample size. Therefore, a kernel-based estimator could be employed for any type of dataset, categorical or a mixture of both continuous and categorical, to address the shortcoming of the frequency method as it does not rely on sample splitting [17].

We start with a univariate kernel density function for continuous and discrete variables, then expand it to our methods. The kernel function for the continuous attribute could be written as:

(1)

Where k is the kernel function, X^c_ij are univariate independent and identically distributed samples drawn from some unknown distributions at a point of x_j^c. h is called smoothing parameter, bandwidth or window width. It is intuitive from the above equation that when the smoothing parameter of h is too large, the important features would be obscured [18].

On the other hand, the kernel for discrete variables could be written as:

(2)

Where λ_q are smoothing parameters for discrete attributes. While the boundary of 0 < h_j < ∞, the boundary of λ_j is 0 ≤ λ_j ≤ (r_j-1)/r_j, where r_j is number of categories of a discrete variable. So, in the case of having a binary attribute, the maximum of λ_j would be 0.5. In other words, while in case of being an unrelated attribute for a continuous attribute, h_j → ∞ for categorical variable and in case of binary predictor λ_j → 0.5, making a binary attribute to be an unrelated attribute.

While using any nonparametric method it should be noted that the inclusion of the parameters that are independent of the response should be removed before conducting any statistical analysis as not excluding the irrelevant attributes would degrade the parameters’ estimate accuracy and prediction. The problem of irrelevant attributes in the nonparametric setting could be alleviated by choosing smoothing parameters. The cross-validation technique could be implemented by simultaneously choosing smoothing parameters for attributes and thus removing irrelevant explanatory components [9]. It would be implemented by extending the irrelevant attributes’ bandwidths to infinity and thus eliminating irrelevant attributes. The cross-validation could be based on the simplest form of minimization of least-square for selection of the smoothing parameters [19] or other methods.

This subsection is based on the work in the literature [9]. The objective function is to identify smoothing parameters that minimize the MISE as:

(3)

Where h_p and λ_q are bandwidths for the continuous components and smoothing parameters for the discrete components, respectively. Recall, 0 < h_p < ∞, and 0 ≤ λ_q ≤ (r_j - 1)/r_j, where x^d takes values of 0, 1,…, r_j-1. x^c is p-variate, while x^d is q-variate. m is a marginal density and m(x) = m ^c(X ^c| X ^d) P (X ^d = x ^d). It is intuitive from Equation 3 that the objective function is mainly based on minimization of the expected variation between the expected conditional density. g ^ (y | x), and g (y | x).

g (y | x) highlights the density of y conditional on x, which could be written as g (y | x) = f (x,y) / m(x) while g ^ (y | x) = f ^ (x,y) / m ^ (x). Now for estimation of g ^ (y | x), in Equation 3, the estimated of the two parts in the numerator and denominator, f ^ (x, y) and m ^ (x) , could be written as:

(4)

(5)

Where K and L are generalized kernels. Comparably, K for considering both continuous and categorical would be written as:

(6)

X, as the whole explanatory variables, which could be split into X= (x^c, x^d), where x ^c and x ^d accommodate continuous and discrete predictors, respectively. L, which was used in 4, would be changed accordingly based on y and Y in L (y, Y_i)

(7)

For this method, cross-validation would be used to overcome challenges of incorporating attributes that seem to be irrelevant and included in the model. It would be implemented by assigning large smoothing parameters and shrinking the parameters to the uniform distribution so they would not be considered in model performance.

Obtaining the optimal values of smoothing parameters for various attributes, the MISE needs to be minimized. The MISE in 3, and after considering the above equations, could be extended to [9]:

(8)

Where X is a function including dominant terms of MISE, where “a” and “b” are functions of h and λ, and n is the number of observations. For instance, h _j = a_j n^{(-1 / p + 5)}, or λ_j = a_j n^{(-1 / p + 5)}. In summary, the parameters, such as smoothing, would be identified by minimizing X or MISE, by having constrained the parameters to be non-negative. It should be noted, for instance, that X, in addition to being a function of all a_p and b_q, it is also a function of the second derivative of f (x^c, x^d, y) with respect to y.

Equation 8 is also based on “u^d” . “u^d”, compared with x ^d, being only different from attributes in jth components, in case of the component having no information regarding the response. In other words, in case of irrelevant information, the jth component would be removed from x^c, or x^d. For estimation of MISE, and as expected, due to the multivariate nature of the algorithm, in the process of X every pairwise consideration of x^c, u^d given y, f (x^c, u^d, y), m (x^c, u^d) and f (x, y) would be considered. The interested readers are referred to the literature review for more description [9].

The following paragraphs will outline a single index model or semiparametric estimation. Besides nonparametric method, the semiparametric method of the single-index model was considered in this study. The description of the semiparametric single index model in this subsection is based on the work in the literature [10]. The method, the functional form of the choice probability function, is characterized by an index, and thus the model makes no assumption regarding the distribution or generating the disturbance. The method has been discussed for the binary choice model only as:

(9)

where ϑ in ϑ (x; θ_o) is a known function with x as a vector of exogenous variables, and θ_o as an unknown parameter vector, and u_o is the random disturbance. The model is called the index restriction model, as E (y | x) is restricted to E (y | ϑ (x; θ_o)), where E is the conditional expectation, and ϑ (x; θ_o) is an index. u_o from above could be written as u_o ≡ s(x; β_o). ε, where s is the positive scaling function known, and other parameters were defined earlier. u_o highlights the known heteroscedasticity of a known form of s(x; β_o) and heteroscedasticity of unknown ε being independent of x.

Now, maximum likelihood could be employed to estimate θ_o in 9 as:

(10)

It should be noted that Equation 10 looks similar to the binary logit model with the difference that ln (P_i^*) and ln (1 - P_i^*) would be estimated. Where P_i^* could be written as:

(11)

Where F_u|x is a known cumulative density function (CDF). To simplify further, while assuming u and x are independent, F_u|x could be written as F_u. Now, F_u would be maximized by θ jointly. However, still, the distribution function might be replaced by its maximum likelihood estimate. Thus P_i^*(θ) would be replaced by P(θ), which is tractable and could locally approximate P_i^*(θ).

After doing some algebra, the objective function process would be transformed to estimate θ by minimizing the Kullback-Leibler information criteria (KLIC), or discrepancy between θ^(p^) and θ^(p) which would be equivalent to maximizing the quasi-likelihood function as follows:

(12)

The above is based on the fact that the known probability function of P performs asymptotically similar as Q(P^(θ)) or conditional probability function, where Q stands for quasi-likelihood.

In summary, P_i^*(θ) in 11 is replaced with tractable P_i(θ) in 12, which could be locally approximated. That starts from constructing P_i(θ) by using C as the event u < v (x; θ), so P_i^*(θ) in 11 would be written as P^*[v(x; 0); 0] = Pr[C]_gv|c(v; θ) / g(v; θ), which could be replaced by C_o, which is an event C at θ_o, or u_o < ϑ(x; θ_o) equivalent for y=1 in Equation 9.

As a result, we could write p(v_o; θ_o) based on probability of y=1, conditioned on various explanatory variables, Pr[y=1|x]. So, the analysis would be employed on the proportion of pedestrian crashes that result in severe crashes.

The models described above were applied to a case study of the pedestrian crash dataset in Wyoming. The data consists of a sample of pedestrian crashes that occurred in Wyoming during 2010-2019. The dataset was directly obtained from the Wyoming Department of Transportation (WYDOT). The considered attributes in Table 1 are based on their significance in parametric analysis. Table 1 highlights that across all considered variables, only posted speed limit is continuous. The characteristics are all related to pedestrian crashes and the mean highlights the distribution of those pedestrians involved in crashes. For instance, the pedestrian’s gender could be divided into almost equal proportions, mean=0.43.

This section is conducted to compare the performances of the implemented methods. As semi- and nonparametric methods are not based on log-likelihood, standard measures such as the Akaike information criterion (AIC) are not suitable for comparison. As a result, the performances are evaluated by means of confusion matrices. To address the concern that the methods have seen the used data, we split our data into two subsets of the test and training dataset and the confusion matrices of both datasets are presented in Table 2. It should be noted that only the correct classification rate (CCR) of the test data is considered and presented in Table 2.

About 14% (No=100) of balanced observations randomly were assigned for the test dataset, while the remaining (No=711) were used to train the model. The numbers are chosen due to limited observations and a balanced vision about the model performance on the test dataset. As can be seen from Table 2, the accuracy for the training dataset is 75%, 68%, and 67% for nonparametric, semi-parametric, and parametric methods, respectively.

However, as it is possible that the methods were overfitted, it is more reliable to look at the models’ performance for the test dataset. The model performance on the test dataset highlights that the nonparametric method outperforms the other method by 79% accuracy compared with 71% and 77% for semi-parametric and parametric methods, respectively. The severe crash prediction is of crucial importance for policymakers ofthe state as reduction of injury or severe crashes is aparamount priority. As can be seen from Table 2, CCR (1) has the superior performance for semi-and nonparametric methods, 76% and 68%, respectively. So, semi- and nonparametric methods especially work better than the parametric method on the category of pedestrian severe crashes.

The performance of the included models is in line with the work in the literature review that the nonparametric method would outperform the semi- and parametric counterparts [20]. Finally, it should be noted that another study was conducted to investigate the pedesterian crash severity.

Although nonparametric was able to adjust itself for both test and training datasets, both semi- and parametric methods perform largely better than the parametric method over the training dataset. More studies are needed, especially in traffic studies, to confirm the findings.

As the nonparametric method is not dealing with parameters’ estimates, distributions, and related bandwidths, and parameters’ estimates could not be obtained. As a result, the point estimates of semi- parametric and parametric are only included in Table 3. Although most parameters are significant at the 0.05 significance level based on the binary logistic model, it should be noted that the two variables were significant at 0.1 level only, and the p-value was .19 for an attribute of pedestrian gender for that model.

For both models, the estimations are in line and have expected signs. For semi-parametric, based on design, the first attribute is normalized to 1. Comparison across the parameters’ estimates magnitudes is not recommended as the parametric method incorporates an intercept, and the first attribute of the semi-parametric is normalized to 1. However, across both models, alcohol and drug involvement, driving on non-level grade, and bad lighting conditions are some of the factors that increase the likelihood of pedestrian crash severity.

To have a vision about the bandwidth estimates of the model attributes and to see if any irrelevant variables have been removed from the analysis by assigning upper bandwidth close to the maximum value, Table 4 is provided. Table 4 presents the maximum and assigned bandwidth values for continuous and categorical variables. For the nonparametric method, cross-validation automatically identifies irrelevant and relevant components and diverges the irrelevant to infinity so the distribution would be virtually uniform on the real line [9].

In simple words, the process would be summarized as follows: first uniform priors, or initial parameters values, would be given to the parameters, where λ is close to 0.5. Then the probability estimates would be updated by seeing observations, wherein case of λ being close to 0, the prior would be updated based on the samples, while in case of λ → 0.5, the uniform prior, or initial value, would be the estimated value of the parameter. Recall that for categorical attributes λ_j ≤ (r_j - 1) / r_j.

It is clear that a larger variance is associated with those λ being closer to 0, while a smaller variance, or larger bias, is related to those λ being closer to 0.5. When bandwidth is close to a maximum value of λ, there is proof that the parameter distribution is flat with a variance of almost zero, so the parameter is irrelevant in the model.

Two models were included in Table 4. First, the model, which we discussed in the results of previous subsection. We also added some noisy variables, such as driver age, in a continuous format to the second model. From the first model in Table 4, the nonparametric method even took advantage of those variables that were based on the parametric method found to be associated with uncertainty and not being different than 0: for instance, the bandwidth of pedestrian gender is 0.27 compared with a maximum value of 0.5. So, the attribute was kept in the model.

Moving to the second model, we incorporated noisy attributes. As can be seen from the second model in Table 4, a wider bandwidth closer to the maximum values is given to driver age. A wider bandwidth has been given also to the drug involved attribute. Previously based on Table 3, there was an uncertainty for this attribute based on the logit model. Although the attribute was kept in the first model, when adding some noisy attribute, this predictor was removed from the dataset. The other variable is a discrete attribute of pedestrian age. The nonparametric method has done a good job in giving a wide bandwidth close to the maximum to this value, so the distribution for this attribute would be flat, and thus this attribute would be considered irrelevant, see the second model in Table 4.