Self-Paced Deep Regression Forests with Consideration of Ranking Fairness

Deep discriminative models (DDMs), e.g. deep regression forests and deep decision forests, have been extensively studied recently to solve problems such as facial age estimation, head pose estimation, etc.. Due to a shortage of well-labeled data that does not have noise and imbalanced distribution problems, learning DDMs is always challenging. Existing methods usually tackle these challenges through learning more discriminative features or re-weighting samples. We argue that learning DDMs gradually, from easy to hard, is more reasonable, for two reasons. First, this is more consistent with the cognitive process of human beings. Second, noisy as well as underrepresented examples can be distinguished by virtue of previously learned knowledge. Thus, we resort to a gradual learning strategy -- self-paced learning (SPL). Then, a natural question arises: can SPL lead DDMs to achieve more robust and less biased solutions? To answer this question, this paper proposes a new SPL method: easy and underrepresented examples first, for learning DDMs. This tackles the fundamental ranking and selection problem in SPL from a new perspective: fairness. Our idea is fundamental and can be easily combined with a variety of DDMs. Extensive experimental results on three computer vision tasks, i.e., facial age estimation, head pose estimation, and gaze estimation, show our new method gains considerable performance improvement in both accuracy and fairness. Source code is available at https://github.com/learninginvision/SPU.

I. INTRODUCTION
D EEP discriminative models (DDMs) have been extensively studied recently to solve a variety of computer vision problems, with applications in image classification, action recognition, and object detection, to name a few. DDMs use deep neural networks (DNNs) [1]- [4] to formulate the input-output mapping in discriminative models. Due to a shortage of well-labeled training data, i.e. without noise and imbalanced distribution problems, learning DDMs is particularly challenging. Considerable literature has grown up around the theme of how to learn DDMs effectively [5]- [7].
One typical approach is to learn more discriminative features through rather deep neural networks and cost-sensitive discriminative functions [5], [8]- [11]. The other typical approaches are to reweight training samples according to esti-mation errors [12], [13] or gradient directions [14] (i.e. meta learning). Although these approaches can partially alleviate the noise and imbalanced distribution problem of training data, they are not consistent with the cognitive process of human beings; the difference is we learn new things based on previously learned knowledge. In addition, these approaches do not provide a way to distinguish noisy and underrepresented examples, and thus cannot fundamentally avoid noisy and biased solutions.
Intuitively, learning DDMs from easy to hard may be more reasonable because already-learned knowledge, in such a process, can be leveraged to learn DDMs. In addition to this, the noisy and underrepresented examples could be recognized by virtue of already-learned knowledge. Hence, we resort to a gradual learning strategy i.e., self-paced learning [15]- [18]. So far, however, there are few studies on learning DDMs in a self-paced manner. Then, a natural question arises: can SPL guide DDMs to achieve more robust and less biased solutions?
An underlying problem in SPL, which is firstly shown by this study, is that it assumes the distribution of training data is balanced, and thus the bias of solutions may be further exacerbated when such an assumption fails. The reason is that the underrepresented examples, due to large prediction error, would not be selected for training at early paces, thereby incurring unfair treatment. This means existing SPL methods have fairness issue.
To address the fairness issue in existing SPL methods as well as answer the above questions, this paper proposes a new self-paced learning method for learning DDMs. It tackles the fundamental ranking and selection problems in SPL from a new perspective: fairness. On one hand, our new method keeps SPL's advantage in robustness. On the other hand, it prevents seriously biased solutions induced by SPL. The major contributions of this work include: •   concerning some sensitive attributes. • For verification, we validate SPUDRFs on three computer vision tasks: (i) facial age estimation, (ii) head pose estimation, and (iii) gaze estimation. Extensive experimental results on the Morph II [19], FG-NET [20], BIWI [21] , BU3DFE [22] and MPIIGaze [23] datasets demonstrate the efficacy of our proposed new self-paced method. On all the above tasks, SPUDRFs almost achieve state-ofthe-art performance.
A preliminary version of this work was published in [24]. This paper significantly improves [24] in the following aspects. (i) We extend our model to incorporate capped likelihood, which could further promote robustness. The robust SPUDRFs can recognize and exclude examples with labeling noise, whereas the original is unable to do so, and only places more emphasis on reliable examples. (ii) We extend our model to combine with various weighting schemes, whereas [24] merely adopts a mixture weighting scheme. (iii) We define the fairness metric for regression in this work and show obvious fairness improvement of SPUDRFs against the original SP-DRFs. (iv) We also evaluate SPUDRFs on a new computer vision task, i.e., gaze estimation, and demonstrate its validity on the MPIIGaze dataset. (v) Both robustness and fairness evaluation results are added, along with more comprehensive discussions.
This paper is organized as follows. Sec. II introduces related work on SPL, and DDMs for facial age estimation, head pose estimation and gaze estimation. Sec. III details our proposed SPUDRFs. Sec. IV introduces robust SPUDRFs. Sec. V defines a fairness metric for regression problem. Sec. VI exhibits and analyzes the experimental results on three computer vision tasks and five related datasets. Sec. VII discusses the strengths and potential weaknesses of this work, followed by Sec. VIII concluding this paper with perspectives.
II. RELATED WORK This section reviews SPL methods and DDMs-based facial age estimation, head pose estimation and gaze estimation.

A. Self-Paced Learning
Self-paced learning, as a gradual learning method, builds on the intuition that, rather than considering all training samples simultaneously, the learning should be presented with the training data from easy to difficult [15], [25]. To date, a great deal of study on self-paced learning (SPL) has been undertaken, mostly about learning conventional discriminative models, e.g. SVM, logistic regressor. There also exists some literature in which SPL is employed for clustering problems [26]- [28].
Learning conventional discriminative models in a self-paced manner exhibits superiority over traditional methods. In [15], Kumar et al. proposed the original SPL framework in a regime that suggests processing the samples in order of an easy to hard order. In [16], Zhao et al. generalized the hard weighting scheme in SPL to a soft weighting scheme, which promoted discrimination accuracy. In [17], Ma et al. proposed self-paced co-training, where SPL is applied for multi-view or multimodality problems. In [29], Ren et al. introduced capped-norm into the objective function of SPL, so as to further exclude the interference of noisy examples. In fact, the above work can be cast as incorporating SPL with conventional classifiers, e.g., SVM, logistic regressors, etc.. To incorporate SPL with regression models, in [30], Han et al. learned a mixture of regressors in a self-paced manner, and added an 1,2 norm regularizer to avoid poorly conditioned linear sub-regressors.
In computer vision community, recently some researchers have realized that learning DDMs in a self-paced manner may lead to more robust solutions. One of the few studies is [18], where Li et al. sought to enhance the learning robustness of CNNs by SPL. However, they omitted one important problem in learning discriminative models: the imbalanced distribution of training data.
Our work is inspired by the existing studies [31], [32] which take the class diversity in sample selection of SPL into consideration. Jiang et al. [32] proposed to ensure the class diversity of samples at the early paces in self-paced training. Yang et al. [31] defined a metric called "complexity of image category" to measure the number of samples in each category, as well as jointly classifying difficulty, and used it to select samples. The two methods mentioned above find that a lack of class diversity in sample selection may result in seriously biased solutions. However, what causes this problem has not been studied. This work shows that the cause is fundamentally the unfairness of sample ranking in SPL, where underrepresented examples may often have large losses and not be selected at early paces. On the other hand, [31], [32] are only suitable for classification, but not for regression in which the output targets are continuous and high dimensional. This paper will go further in this direction, aiming to tackle the fundamental drawback in SPL-ranking unfairness-and to integrate SPL with DDMs.

B. Facial Age Estimation
Facial age estimation has been extensively studied for a decade in the computer vision community. In recent years, a large and growing body of literature [3], [8], [10], [33], [34] has been proposed for age estimation with varying degrees of success. Ordinal-based approaches [3], [10] are the bestknown and have demonstrated improved results. For example, Niu et al. proposed to estimating age through a set of sequential binary queries-each query refers to a comparison with a predefined age-to exploit the inter-relationship (ordinal information) among age labels. Furthermore, Gao et al. [8] proposed DLDL-v2 to explore the underlying age distribution patterns, so as to effectively accommodate age ambiguity. Shen et al. [33] used VGG-16 to extract facial features, and mapped the extracted features to age by deep regression forests (DRFs). Li et al. [34] proposed BridgeNet to predict age by mixing local regression results, where local regressors are learned on partitioned data space. Deng et al. [35] proposed an age classification method to deal with long-tailed distribution, where a variational margin is used to minimize the influence of head classes that misleads the prediction of tailed samples. Overall, these DDMs-based approaches have improved age estimation performance significantly. However, they plausibly ignore one problem: the interference with facial images duo to PIE (i.e. pose, illumination and expression) variation, occlusion, misalignment and so forth. More importantly, these approaches mostly ignore the fact that the training data is not distributed uniformally.

C. Head Pose Estimation
Head pose estimation has attracted vast attention from the computer vision community for a long time. Recently, an increasing amount of DDMs-based methods have been proposed for head pose estimation and have dramatically boosted estimation accuracy. Of these, Riegler et al. [36] utilized convolutional neural networks (CNNs) to extract patch features within facial images, and more precise results were achieved.
In [4], Huang et al. proposed multi-modal deep regression networks, with multi-layer perceptron (MLP) architecture, to fuse the RGB and depth images for head pose estimation. In [37], Wang et al. proposed a deep coarse-to-fine network to further boost pose estimation accuracy. In [7], Ruiz et al. used a large, synthetically expanded head pose dataset to train a rather deep multi-loss convolutional neural networks (CNNs), gaining obvious accuracy improvement. In [38], Kuhnke et al. proposed a domain adaptation method for head pose estimation where shared label spaces across different domains are assumed. Although accuracy has improved, handling noisy examples in the above models is not straightforward. Moreover, the methods mentioned above seldom deal with imbalanced training data, which may widely exist in visual problems.

D. Gaze Estimation
Appearance-based gaze estimation learns a mapping from facial or eye image to gaze. Recently, a considerable amount of DDMs based approaches have been proposed for personindependent gaze estimation. For example, in [39], Zhang et al. used spatial weights VGG-16 to encode face image to gaze direction. The spatial weights were used to flexibly suppress or enhance features from different facial regions. In the literature [40], Krafka et al. implemented a CNNs-based gaze tracker by mapping captured images from the left eye, right eye and face to gaze. The location and size of the head in the original image are provided for mapping networks. In [41], Fischer et al. adopted two backbone networks to extract features from two eye regions for gaze regression, and achieved improved robustness. Instead of predicting gaze directly, in [42], Park et al. proposed another method which maps the single eye image to an intermediate pictorial representation as additional supervision to simplify 3D gaze direction estimation. Wang et al. [43] proposed a Bayesian adversarial learning approach to address the appearance, head pose variations in gaze estimation, and overfitting problem. Xiong et al. [44] proposed mixed effect neural networks to combine global and person specific gaze estimation results together. Meanwhile, other studies show that facial images can help gaze estimation. In [45]- [47], multi-channel estimation architectures were adopted to utilize facial image and eye images jointly. These methods outperform the ones that only use eye images as input for gaze estimation. In addition to person-independent methods, there is some literature about person-specific gaze adaption methods [48]- [50], which have boosted performance dramatically when a few annotated examples are used for calibration. These methods indeed improve gaze estimation results; however, whether we should learn DDMs for gaze estimation in a gradual learning manner has not been discussed, nor has there been discussion of how to deal with imbalanced distribution of training data.

III. SELF-PACED DEEP REGRESSION FORESTS WITH CONSIDERATION ON UNDERREPRESENTED SAMPLES
This section first reviews the basic concepts in DRFs; and introduces the objective formulation for SPUDRFs, as well as variant weighting strategies and an underrepresented example augmentation method; and finally details the optimization algorithm.

A. Preliminaries
Deep regression forests (DRFs), as a deep regression model, connect deep neural netwoks (DNNs) to regression forests. We start by reviewing the basic concepts in DRFs [33]. Deep Regression Tree. DRFs usually consist of a number of deep regression trees, each of which, given input-output pairs , map extracted features through DNNs to target output by passing a regression tree. Further, a regression tree T consists of split (or decision) nodes N and leaf (or prediction) nodes L [33] (see Fig. 1). Specifically, each split node n ∈ N determines whether a sample x i goes to left or right subtree; each leaf node ∈ L represents a Gaussian distribution p (y i ) with mean µ l and variance σ 2 l -parameters of output distribution defined for each tree T . Split Node. Split node represents a split function s n (x i ; Θ) : where Θ denotes the parameters of DNNs, as in Fig. 1 , where σ(·) denotes a sigmoid function while ϕ(·) denotes an index function specifying the ϕ(n)-th element of f (x i ; Θ) in correspondence with the split node n, and f (x i ; Θ) denotes the extracted features through DNNs. An example to illustrate the sketch chart of the DRFs is shown in Fig. 1(c), where ϕ 1 and ϕ 2 are two index functions for two trees. Finally, the probability of the sample x i falling into the leaf node is given by: where [·] denotes an indicator function conditioned on the argument. In addition, L nleft and L nright correspond to the sets of leaf nodes owned by the sub-trees T nleft and T nright rooted at the left and right children n l and n r of node n. Leaf Node. For a tree T , each leaf node ∈ L defines a Gaussian distribution over output y i conditioned on µ l , σ 2 l . Since each example x i has the probability ω (x i |Θ) to reach leaf node , considering all leaf nodes, the predictive distribution over y i shall sum all leaf distribution weighted by ω (x i |Θ): where π represents the distribution parameters of all leaf nodes associated with tree T . Note that π varies along with tree T k , and thus we rewrite them in terms of π k . Forests of Regression Trees. Since a forest F consists of a number of deep regression trees {T 1 , ..., T k }, the predictive output distribution shall consider all trees: where K is the number of trees and Π = {π 1 , ...,

B. Objective
Underrepresented Examples.
Considering underrepresented examples in SPL is one of the main contributions of this work. In this work, underrepresented examples mean "minority". The underrepresented level could be measured by predictive uncertainty (i.e. entropy). In fact, underrepresented examples may incur unfair treatment in early paces due to imbalanced data distribution. Our proposed method tackles the ranking and selection problems in SPL from a new perspective: fairness.
Entropy. Given a sample x i , its predictive uncertainty is formulated by calculating the entropy of the conditional output distribution p F (y i |x i , Θ, Π): where H [·] denotes entropy, and the entropy which corresponds to the k th tree is: the more uncertain the prediction should be, i.e., the more underrepresented the sample x i is.
Because the integral of mixture of Gaussians in Eq. (5) is non-trivial, we resort to Monte Carlo sampling. However, because this has a large computational cost [51], we turn to use the lower bound of this integral to approximate the integration: Let v i denote a latent variable indicating whether the i th sample is selected (v i = 1) or not (v i = 0) for training. Our objective is to jointly maximize the log likelihood function with respect to DRFs' parameters Θ and Π, and learn the latent variables v = (v 1 , ..., v N ) T . Fig. 2 shows that the original self-paced method may miss the underrepresented examples at an early pace, resulting in ignorance.
Considering fairness, we prefer to select the underrepresented examples, which probably have higher predictive uncertainty (i.e. entropy), particularly at an early pace. Therefore, we maximize the likelihood function of DRFs coupled with the sample selection term meanwhile considering ranking fairness, where λ controls learning pace (λ ≥ 0), γ is a tradeoff between loss and uncertainty (λ ≥ 0), H i denotes the entropy of output distribution p F (y i |x i , Θ, Π). When γ decays to 0, the objective function is equivalent to the loglikelihood function with respect to DRFs' parameters Θ and Π. In Eq. (7), each sample is weighted by v i , and whether That is, the sample with high likelihood value or entropy may be selected. The optimal v * i is: where p F (y i |x i , Θ, Π) is denoted by p F i for simplicity, and we use this representation in all following parts. Iteratively increasing λ and decreasing γ, samples are dynamically involved for training DRFs, starting with easy and underrepresented examples and ending up with all samples. Thus, SPUDRFs are prone to achieve more robust and less biased solutions. One might argue that the noisy and hard examples tend to have high predictive uncertainty, with the result that they are selected in the early pace. In fact, from Eq. (8), we observe that whether a sample is selected is determined by both its predictive uncertainty and the log likelihood of being predicted correctly. The noisy and hard examples probably have relatively large loss i.e. low log likelihood log p F i , and would not be selected in the early pace.

C. Weighting Schemes
Weighting each sample according to its importance would be more reasonable in SPUDRFs. In the previous section, we adopted a hard weighting scheme in SPUDRFs, as defined in Eq. (7), where one sample selected or not is determined by v i . Such a weighting scheme appears to omit the importance of samples. SPUDRFs can easily incorporate other weighting schemes, including mixture weighting and soft weighting [52]. Mixture Weighting. Mixture weighting scheme [52] weights selected sample by its importance, i.e., 0 ≤ v i ≤ 1. The objective function under mixture weighting scheme is: where ζ is a parameter controlling the learning pace. We set ζ = 1 λ − 1 λ −1 , and λ > λ > 0 to construct a reasonable soft weighting formulation. The self-paced regularizer in Eq. (9) is convex with respect to v ∈ [0, 1]. Then, setting the partial gradient of Eq. (9) with respect to v i to be zero will lead the following equation: The optimal solution of v i is given by: If either the log likelihood or the entropy is too large, v * i equals 1. In addition, if the likelihood and entropy are both too small, v * i equals 0. Except for the above two cases, the soft weighting, i.e., the last line of Eq. (11), is adopted. Soft Weighting. A soft weighting scheme [52] weights a selected sample with respect to its output likelihood and entropy. Such a weighting scheme includes: linear weighting scheme and logarithmic weighting scheme. Linear weighting. This scheme linearly assigns weights to samples with respect to output likelihood and entropy. The objective of SPUDRFs under linear weighting scheme is formulated as: where We set the partial gradient of Eq. (12) with respect to v i to be zero, then the optimal solution for v i is: The larger the log likelihood and entropy are, the higher the weight associated with the i th sample should be. Logarithmic weighting. This scheme penalizes the output likelihood and entropy logarithmically. The objective of SPU-DRFs under a logarithmic weighting scheme can be formulated as: where

D. Underrepresented Example Augmentation
As previously explained, the SPUDRFs method places more emphasis on underrepresented examples and may achieve less biased solutions. Since the intrinsic reason for SPL's biased solutions is the ignorance of underrepresented examples, we further rebalance training data via distribution reconstruction. Specifically, we distinguish the underrepresented examples whose H i are larger than β, from regular examples at each pace and augment them. As the number of underrepresented examples increases through augmentation, the label distribution imbalance problem is alleviated.

Algorithm 1
The training process of SPUDRFs.

E. Optimization
To optimize SPUDRFs, as defined in Eq. (7), we propose a two-step alternative search strategy (ASS) algorithm: (i) For sample selection, update v with fixed Θ and Π (ii) update Θ and Π with current fixed sample weights v. In addition, with fixed v, our DRFs are learned by alternatively updating Θ and Π. In [33], the parameters Θ for split nodes (i.e. parameters for VGG) are updated through gradient descent since the loss is differentiable with respect to Θ. In comparison, the parameters Π for leaf nodes are updated by virtue of variational bounding [33] when fixing Θ.

EXAMPLES
As has already been discussed in Sec. III, SPL tends to place more emphasis on reliable examples to achieve more robust solutions. However, SPL's intrinsic selection scheme can not exclude noisy examples, especially the examples with labeling noise. To alleviate this drawback, we introduce capped-likelihood function in SPUDRFs, which can render an output likelihood with an especially small value as zero: where denotes the threshold and > 0.
The robust SPUDRFs can be optimized using the optimization method proposed in Sec. III-E. Note that, if p Fi ≤ , the capped likelihood p c Fi = 0, that is, log p c Fi = −∞. Since λ is a positive factor, maximizing Eq. (17) must yield v i = 0, meaning that the i th sample would not be selected. By adjusting , we can exclude a certain ratio of noisy examples to obtain more robust solutions.

V. FAIRNESS METRIC
In addition to accuracy, fairness should also be an essential metric to measure the performance of a regression model, particularly when the regression targets are related to sensitive attributes. For example, in age estimation, we expect our model to treat the younger and the elder fairly, i.e., not to result in large MAEs for the former but small MAEs for the latter. The notion of fairness was originally defined concerning a protected attribute such as gender, race or age. However, the term fairness has a range of potential definitions. Here, we adopt a notion of fairness that is for sensitive features [53]. Defining a fairness metric for regression is one of the contributions of this work.
The present studies [54]- [57] mostly refer to fairness constraints in regressions, for example, statistical parity or bounded group loss [54], [56], but seldom refers to fairness metrics. In this work, we define a new fairness metric for regression models. To be specific, the test dataset is divided into K subsets, each of which is denoted by D k = {x i , y i |y i ∈ G k }, and G k denotes the k th group. A fair model is expected to have the same performance on all subsets, which can be described mathematically as: denotes the expectation with respect to the loss L (·) over group D k ,ŷ is the predicted value of the model and y is the real target value. Motivated by p%−rule [57], which measure classification fairness, we evaluate fairness between two subsets D k and D l as follows: The loss expectation ratio characterizes the model's fairness on every two subsets. A small ratio means the performance on such two subsets is significantly distinct, reflecting model bias against a particular subset. On the contrary, a large ratio indicates the losses on D k and D l are similar, indicating the model is relatively fair for D k and D l . In particularly, when the ratio is equal to 1, the model satisfies Eq. (18), and D k and D l are treated fairly. Finally, the regression fairness is defined as the expectation of the fairness between any two subsets, In SPUDRFs, when we set γ>0, all samples are sorted by both likelihood and entropy. As a result, easy and underrepresented samples are selected first, which means samples from all subsets would be selected at early paces. Therefore, our model can alleviate the model's prejudice against different subsets and improve regression fairness. More extensive evaluation results can be found in Sec. VI-F.

A. Datasets and Experimental Setup
Dataset. There are five datasets used in our experiments, namely Morph II, FG-NET, BIWI, BU3DFE and MPIIGaze.   [24]. The hyper-parameters of DRFs were: tree number (5), tree depth (6), output unit number of feature learning (128). The iterations to update leaf node predictions was 20, and the number of images to update leaf node predictions was the whole training set for MPIIGaze and 50 times of batch size for other datasets. In the first pace, 50% samples were selected for training. Here, λ was set to guarantee the first 50% samples with large log p F i + γH i values were involved. λ was set to ensure that 10% ∼ 20% selected samples adopt soft weighting under a mixture weighing scheme. To promote efficiency, the samples selected in the previous pace were preserved and the rest were ranked for current sample selection. In addition, γ was initialized to be 15 for Morph II and BIWI, and 5 for FG-NET, BU-3DFE, and MPIIGaze. β was used to recognize underrepresented examples that need to be augmented twice in each pace (1000 for BIWI, 400 for BU3DFE, 2000 for MPIIGaze). The number of paces was empirically set to be 10 for Morph II, 6 for BIWI and MPIIGaze, 3 for FG-NET and BU-3DFE. Except for the first pace, an equal proportion of the rest data was gradually involved in each pace. Evaluation Details. To evaluate regression accuracy, we used the mean absolute error (MAE). MAE is defined as N i=1 |ŷ i − y i | /N , whereŷ i represents the predicted output for the i th sample, and N is the total number of images for testing. In addition, for age estimation, CS calculates the percentage of images sorted in the range of [y i − L, y i + L]:   measure regression fairness, we adopt our proposed fairness metric FAIR given in Sec.V. To calculate the above metrics, for Morph II, BIWI and BU3DFE, each dataset was divided into a training set (80%) and a testing dataset (20%). The random division was repeated 5 times, and the reported MAEs were averaged over 5 times. The leave-one-person-out protocol was used for FG-NET [33] and MPIIGaze [23], where one subject was used for testing and the left subjects for training.

B. Validity of SP-DRFs and SPUDRFs
This section validates the original SPL method as well as the new proposed SPL method for learning DDMs. In the following, SP-DRFs denotes self-paced deep regression forests, which learn DRFs using the original SPL method. Self-Paced Learning. Recall that learning DDMs in a gradual learning manner is more consistent with the cognitive process of human beings, and the noisy examples can be distinguished by virtue of learned knowledge. That means the learning can place more emphasis on "reliable" examples. To show this, Fig. 4(a) illustrates the representative face images at each learning pace of SP-DRFs on the Morph II dataset, sorted by increasing λ and decreasing γ. The two numbers below each image are the actual age (left) and predicted age (right). We observe that the training images involved in each pace become more confusing and noisy step by step, compared to images in early paces. Since the current model is initialized by the results obtained at the last pace, the updated model is adaptively calibrated by "reliable" examples rather than by confusing and noisy ones. Thus, SP-DRFs have improved MAE and CS compared to DRFs. We observe that SP-DRFs gain improvement by 0.33 in MAE (2.17 → 1.84), by 1.55% in CS (91.3% → 92.85%), as shown in Fig. 6(a).
Some representative eye images in the gradual learning sequence for MPIIGaze dataset are shown in Fig. 4(b). The images in each panel are sorted by increasing λ and decreasing γ. The two numbers below each pair of images are the actual gaze direction (left) and the predicted gaze direction (right). The same phenomena-the easy examples are prone to be selected in the early paces, while the confusing and noisy ones are prone to be selected in the later paces-can be observed. Since the updated model at each pace is adaptively calibrated by "reliable" examples rather than by confusing and noisy ones, thus, the MAE associated with each pace decreases step by step. Finally, the MAE of SP-DRFs decreases to be 4.57, whereas the MAE of DRFs is 4.62 (see Table. II). Considering Ranking Fairness. As was mentioned in Sec. III, the existing SPL methods may exacerbate the bias of solutioins. Fig. 5 visualizes the leaf node distributions of SP-DRFs in the gradual learning process. The Gaussian means µ l associated with the 160 leaf nodes, where each 32 leaf nodes are defined for 1 tree, are plotted in each sub-figures. Three paces, i.e. the 1 st , 3 rd , and 6 th pace, are randomly chosen for visualization. For clarity, only pitch and yaw angles are shown. Meanwhile, the leaf node distributions of SPUDRFs are also visualized in Fig. 5.
In Fig. 5, the comparison between SP-DRFs and SPUDRFs validates our proposed new SPL method. In SP-DRFs, because the ranking fairness is not considered, the leaf nodes (red points in the 2 rd row) are only concentrated in a small range that most samples are distributed over, thus leading to seriously biased solutions. The poor MAEs of SP-DRFs can serve as evidence for this. In contrast, because the ranking fairness is considered in SPUDRFs, the leaf nodes are distributed over a wide range that could cover underrepresented examples, thus improving performance. SPUDRFs, in the pitch and yaw directions, achieve the best performance with MAEs of 0.71 and 0.69, compared to SP-DRFs with MAEs of 1.05 and 1.23 (47.9% and 78.3% improvements).

C. Comparison with State-of-the-art Methods
We compare SPUDRFs with other state-of-the-art methods on three vision tasks: age estimation, head pose estimation and gaze estimation.
Results on Age Estimation. The comparison between our method and other baselines on the Morph II and FG-NET datasets are shown in Fig. 6(a) and Fig. 7(b). The baseline methods include: LSVR [60], RCCA [61], OHRank [62], OR-CNN [10], Ranking-CNN [3], DRFs [33], DLDL-v2 [8], and PML [35]. The results show some consistent trends. First, SPUDRFs have superior performance compared to con-ventional discriminative models, such as LSVR [60] and OHRank [62]. Second, SP-DRFs consistently outperform other DDMs. Compared to DRFs, our gains in MAE are 0.33 on Morph II and 0.12 on FG-NET, and in CS are 1.55% and 2.91%, respectively. The promotions demonstrate that learning DRFs in a self-paced manner is more reasonable. Third, SPUDRFs outperform SP-DRFs on both MAE and CS. Fig. 6(b) shows the CS curves of SPUDRFs, SP-DRFs and other baseline methods on the Morph II dataset, and Fig. 7(b) shows the CS curves of different methods on the FG-NET dataset. On both datasets above, SPUDRFs consistently outperform other DDMs. We observe, the CS of SPUDRFs reaches 93.34% at error level L = 5, significantly outperforming DRFs by 2.04% on the Morph II dataset. We also observe that SPUDRFs outperform DRFs by 4.09% in CS, on the FG-NET dataset. The CS increase clearly validates our proposed self-paced learning method.
Results on Head Pose Estimation. Recall that considering underrepresented examples in SPL is particularly crucial when the training data has an imbalanced distribution problem. Sec. VI-B has shown the effectiveness of SP-DRFs and SPUDRFs for head pose estimation tasks on the BIWI dataset. Fig. 5   ranking fairness in SPL can obtain more reasonable results.

Method
MAE↓ CS↑ IIS-LDL [63] 5.77 N/A LARR [64] 5.07 68.9% MTWGP [65] 4.83 72.3% DIF [66] 4.80 74.3% OHRank [62] 4.48 74.4% CAM [67] 4.12 73.5% DRFs [33] 2  Results on Gaze Estimation. For gaze estimation, the accuracy comparisons between our proposed SP-DRFs, SPUDRFs and other baselines on the MPIIGaze dataset are shown in Tab. II. Both the MAE and standard deviation across all persons are reported. Note that the standard deviation of RT-GENE [41], Pict-Gaze [42], and Ordinal Loss [72] are not reported because the original studies do not provide this information. As shown in Tab. II, we observe that SPUDRFs outperform all baseline methods in MAE, confirming the effectiveness of our proposed method. Because RT-GENE [41] directly maps the features extracted by VGG-16 to gaze through multiple FC layers, it has inferior performance relative to DRFs (4.62 → 4.80). For a fair comparison, we chose RT-GENE with 1 model. Further, because SPUDRFs method learn DRFs in a self-paced manner and take into account ranking fairness, it has a further gain in MAE over DRFs (4.62 → 4.45). Pict-Gaze [42] regresses an input image to an intermediate pictorial representation and then regresses the representation to the gaze direction. Ordinal Loss [72] utilizes ordinal loss with order regularization to solve the regression problem. The two methods mentioned above do not take into consideration the differences amongst samples; they train all examples simultaneously and thus have inferior MAEs. We observe that the MAE of Pict-Gaze [42] and Ordinal Loss [72] are 4.56 and 4.49 respectively, while ours is 4.45. It's noteworthy that SP-DRFs, when compared with DRFs, only promotes MAE slightly. This is probably due to the obvious distribution difference between training data and test data in the leave-one-out setting. Method MAE↓ MPIIGaze [23] 6.30±1.80 iTracker [40] 6.20±0.85 GazeNet+ [73] 5.40±0.67 MeNets [44] 4.90±0.59 RT-GENE [41] 4.80± − − Pict-Gaze [42] 4.56± − − Ordinal Loss [72]

D. Different Weighting Schemes
A potential concern for SP-DRFs and SPUDRFs is that different weighting schemes could affect the estimation performance on a variety of visual tasks. To evaluate this, we compare SP-DRFs/SPUDRFs under different weighting schemes, including hard, mixture, and soft weighting. Under all weighting schemes, λ was set as in Sec. VI-A. Under the mixture weighting scheme, λ was set to ensure that 10% selected samples adopt soft weighting for Morph II, and 20% samples for other datasets. Table. III shows the comparison results. The performances of SP-DRFs/SPUDRFs with different weighting schemes are only slightly different. We observe that the weights for a large proportion of examples are close to 1. For example, under the log weighting scheme, only 2% of examples have weights below 0.5. Under the mixture weighting scheme, the proportion of samples whose weights are 1 can be set manually. We chose to set this proportion to be 0.8 ∼ 0.9 on different tasks. The MAEs are not guaranteed to be better than other weighting schemes.

E. Robust SPUDRFs
The intuition for SPUDRFs to work better on different visual tasks is its improved robustness, i.e., emphasizing more on "reliable" examples. To further promote the robustness, we propose robust SPUDRFs in Sec. IV, which enable SPUDRFs to handle labeling noise. We added noise to labels to test the validity of robust SPUDRFs on the above datasets. Specifically, we chose 10% samples in Morph II, BIWI or MPIIGaze datasets, and added Gaussian noise N (0, 10) to their labels. To control the proportion (0% ∼ 20%) of samples whose likelihoods are capped to be 0, i.e., the portion of samples to be excluded, we set at variant values. Fig. 8 shows the MAE curves of SPUDRFs across variant capped proportions. We observe that, when no example's likelihood is capped at 0, due to the presence of noise, the corresponding MAE ais large for each dataset. As the capped proportion grows, the MAE gradually decreases. When the capped proportion changes to become 10%, the MAE in Fig. 8  may be discarded when the capped proportion becomes over 10%.

F. Fairness Improvement
This section discusses how SPUDRFs improve regression fairness on different visual tasks.
Taking into account of underrepresented examples, SPU-DRFs show improved accuracy for age, pose and gaze estimation, compared to SP-DRFs. In this section, we show how SPUDRFs further improve regression fairness. In Sec. V, we define FAIR as a regression fairness metric. To evaluate the regression fairness of DRFs, SP-DRFs, and SPUDRFs, we divided all the datasets mentioned above into different subsets and calculated FAIR on them. For Morph II, the entire data was divided into 7 groups, i.e., [10,20], [20,30], . . . , [70,80]. A similar division was conducted on the FG-NET dataset. For BIWI and BU-3DFE, the pitch, yaw, and roll directions were regarded as independent directions, and the interval for each regression group is 10 • . For MPIIGaze, the pitch and yaw directions were regarded as relevant angles. We set the group interval to be 5 • in each direction, resulting in a total of 60 groups.   Some pose estimation results of DRFs, SP-DRFs, and SPU-DRFs for BIWI are shown in Fig. 9. The estimation accuracy for different pose groups reveals some consistent trends. First, for most groups, SP-DRFs and SPUDRFs have superior MAEs over DRFs. In Sec. III-B, we discussed that these gains in fairness are due to SPL, which can guide DDMs to achieve more reasonable solutions. Second, SP-DRFs tend to have even worse MAEs for the underrepresented groups than DRFs, for example, [ ] in yaw direction. This also serves as an evidence that our proposed self-paced method alleviates the fairness problem in existing SPL methods.

VII. DISCUSSION
In the experiments above, we have evaluated SPUDRFs against other baseline methods on different visual tasks, such as facial age estimation, head pose estimation and gaze estimation (Sec. VI-C). We also evaluated the SPUDRFs method under different weighting schemes, its extension with capped likelihood formulation, and its performance improvement on fairness. On a number of tasks and datasets, SPUDRFs and SP-DRFs outperform other baseline methods. The advantage of considering ranking fairness in SPUDRFs is most obvious for BIWI. For Morph II and BU-3DFE, the performance improvement when considering ranking fairness in SPUDRFs is also observable.
Learning DDMs in a self-paced manner has some limitations. The most noticeable one is that, in a leave-one-out setting, SP-DRFs perform only slightly better than DRFs. We speculate that this is due to the data distribution difference between training set and test set. Therefore, whether SPL can improve the performance of DDMs may largely depend on the distribution divergence between training data and test data.

VIII. CONCLUSION
This paper explored how self-paced learning guided deep discriminative models (DDMs) to obtain robust and less biased solutions on different visual tasks, e.g. facial age estimation, head pose estimation, and gaze estimation. A novel self-paced method, which considers ranking fairness, was proposed. Specifically, a new ranking scheme that jointly considers loss and fairness was proposed in SPL. Such a method was combined with a typical DDM-deep regression forests (DRFs)and led to a new model, namely deep regression forests with consideration on underrepresented examples (SPUDRFs). In addition, SPUDRFs under different weighting schemes, their extension with capped likelihood formulation, and their performance improvement on fairness were discussed. Extensive experiments on three well-known computer vision tasks demonstrated the efficacy of our proposed new selfpaced method. The future work will include exploring how to incorporate such a method with other DDMs.


extfieldsofstudy:
Computer Science
Medicine

text:
A Tactical Conflict Resolution Proposal for U-Space Zu Airspace Volumes

Conflict management between UAVs is one of the key aspects in developing future urban aerial mobility (UAM) spaces, such as the one proposed in U-Space. In the framework of tactical conflict management, i.e., with the UAVs in flight, this paper presents PCAN (Prediction-based Conflict-free Adaptive Navigation). This relatively simple navigation technique predicts the occurrence of the conflict and avoids it by modifying the velocity vector of the UAVs involved. The performance evaluation carried out demonstrates its effectiveness compared to similar techniques, even in high-density scenarios, while proving a low overhead in flight time or in the distance traveled by the UAVs to reach their destinations.

Introduction
In recent years, interest in Unmanned Aerial Vehicles (UAVs) [1] has grown considerably, given the wide range of possibilities they offer. Within this revolution, it is more than likely that, in the short or medium term, the airspace of the metropolitan environment will be shared by "traditional" manned aerial vehicles and UAVs, which will be mainly electric and will have vertical takeoff and landing capabilities (eVTOL) [2]. These UAVs will fly autonomously [3] at low or very low altitude levels, covering a wide range of services, including the transport of goods or people, and contributing to reducing the surface and sub-surface congestion and the carbon footprint produced by our daily activity. These scenarios are called Urban Air Mobility/Advanced Air Mobility (UAM/AAM) [4].
Several public and private organizations have begun to restructure the airspace to integrate UAVs into it. In this sense, the concept of UTM (Unmanned Aircraft System Traffic Management) has been developed in the United States [5,6], while in Europe, this initiative has adopted the name U-Space [7,8]. Similar initiatives have appeared in other geographical areas, such as Korea [9] and Australia [10].
Conflict management is one of the many technical challenges to be solved in this scenario. Here, a "conflict" is understood as a situation in which two or more UAVs are at a distance less than a minimum separation predetermined by regulation. This is a classic problem in traditional Air Traffic Management (ATM) [11], where it is called Conflict Detection and Resolution (CD&R), and there are already countless resolution proposals [12,13]. In short, the aim is to avoid the occurrence of a collision between two aircraft at all costs.
One possible way of classifying conflict management techniques is to differentiate between those that obtain a priori a set of conflict-free routes or flight plans for the UAVs and those that detect conflicts during flight and resolve them by modifying the flight plan of the UAVs involved. Along the same lines, both UTM and U-space distinguish between strategic (pre-flight) and tactical (in-flight) conflict management.
Within the second group of techniques, there are strategies based on the calculation of the space of valid velocities that prevent any potential conflict between UAVs, such as

Tactical Conflict Management in U-Space
As stated, U-Space is focused on the safe integration of many UAVs into the airspace. To this end, it plans to introduce a set of services that will be deployed progressively, based on 4 phases (U1, U2, U3, and U4) with an incremental complexity and level of integration [8,20].
U-Space services related to tactical conflict management are in phase U2. Very briefly, a tracking service provides the estimated positions and trajectories of the UAVs. Then, the monitoring service uses the above information and, if applicable, alerts on potential conflicts. The tactical conflict resolution service is responsible for modifying the UAV flight plan. Finally, the emergency management service provides the UAVs involved with the action to be taken to mitigate the risk.
For UAV/UAS flight, U-Space defines a new Very Low Level (VLL) airspace, which coexists with the ICAO (International Civil Airspace Organization) A-G airspaces [21]. The VLL airspace is decomposed into three volumes, X, Y, and Z, which offer different services and have different access requirements. The main difference between these volumes lies precisely in the provision of conflict resolution services.
There are no such services in X-type volumes, and the pilot, in whose visual range the UAV is located (this is known as Visual Line of Sight, or VLOS), is responsible for ensuring a safe operation. In Y-type volumes, Beyond Visual Line of Sight (BVLOS) flights are allowed using remote pilots connected to U-Space, with strategic (but not tactical) conflict resolution services. Operation plans must be approved prior to takeoff. The drawback is that this can result in large separations between aircraft, both in space and time.
The Z-type volume allows for a higher density of operations, adding tactical conflict resolution services, which reduce the residual risk of strategic for in-and out-of-visualrange and autonomous flights. Specifically, and as discussed above, these services should check for any potential conflict in real-time, based on the current position of all aircraft (and, if possible, their flight plans), and provide the necessary avoidance instructions to the aircraft involved in the form of changes in velocity, altitude or course. If the airspace is controlled by the Air Traffic Service (ATS), it is called a Za volume. As we assume in this paper, the airspace is controlled by U-Space. It is referred to as a Zu volume.

Conflict Resolution Techniques for UAV
As discussed, UAV conflict management techniques have traditionally been classified into proactive (or strategic) and reactive (or tactical). Other possible classifications consider whether conflict control is centralized or distributed in each UAV or whether there is communication between UAVs sharing the airspace [22].
Within strategic management, the literature is abundant in collision-free path planning mechanisms, which have been widely explored for years in mobile robotics. Some of the most popular ones are visibility graphs, Probabilistic Road Maps (PRM), Rapidly-exploring Random Trees (RRT), A-Star (A*) and its multiple variants, Ant Colony Optimization (ACO), Tabu Search (TS), and Voronoi diagrams [23]. In general, all these techniques share a high computational complexity. In the following, we describe with more detail some of the techniques for tactical UAV conflict management.
ORCA [14] is a velocity-based conflict avoidance method. These methods are based on choosing the optimal velocity for a mobile agent (in our case, a UAV) from among all those velocities that avoid conflict with other mobile agents (called "valid velocities"). Each agent predicts its future position and nearby agents based on current velocities and, under this assumption, determines its new velocity according to some optimization criterion. In ORCA, each agent computes a Velocity Obstacle (VO) with respect to each agent flying in its environment. The VO includes all the prohibited velocities for the UAV. ORCA guarantees collision-free navigation in sparse scenarios. In dense scenarios it is possible not to find a valid velocity. Linear programming techniques are used to solve this problem, with a medium to high computational cost [24].
BBCA [15] is also a velocity-based mechanism with lower complexity than ORCA and, therefore, is more appropriate for real-time conflict management. Basically, in BBCA the set of valid velocities for each UAV is represented by a rectangular area. This area is modified in each execution of the algorithm (by very simple operations) to avoid conflict with each of the nearby UAVs. BBCA can avoid any conflict in scenarios with two UAVs, and most of them in denser scenarios, with a reasonably acceptable overhead on the distance traveled and the final flight time of the UAVs.
On the other hand, APF [16] is an example of a force-based conflict avoidance method. These methods simulate particle systems, where each particle exerts a certain force on the nearest ones. In the case of APF, the motion of each agent is determined by an attractive force towards its final destination, while the rest of the agents behave as obstacles exerting repulsive forces on it [25]. The main drawback of these methods is the existence of local minima that prevents the agents from moving towards their goal and a high computational cost.
There are also methods based on swarm intelligence, such as PSO [17]. They are inspired by the behavior of certain animal species, organizing agents into groups that work together to obtain an optimal solution to the problem of reaching their destinations (while avoiding conflicts between them). In PSO, obtaining this optimal solution is based on a continuous optimization problem. The distance of each agent to its destination is iteratively calculated and shared with the rest of the agents. These algorithms also involve a high computational cost, not being good candidates for real-time applications with multiple agents.
Another recent proposal is [26], where conflict management is performed in three separate stages. The first stage consists of strategic path planning using a PSO-type algorithm. Later, in a "pre-tactical" stage, it is proposed to delay the start of flights to avoid unresolved conflicts in the previous stage. Finally, residual conflicts that may appear in flight are solved using the "hovering" technique, which stops the UAV in the air for the time necessary to avoid the conflict.

Conflict-Free Navigation
This section presents our proposal for conflict-free navigation in UAM scenarios, which we have named PCAN (Prediction-based Conflict-free Adaptive Navigation). PCAN aims to avoid any conflict between a set of UAVs in flight, regardless of the airspace configuration. As the strategy will be based on modifying the velocity vector of the UAVs, the algorithm should try to introduce as little change as possible in the optimal velocity of each UAV. We will understand as optimal a solution in which the UAVs fly in a straight line towards their destination, traveling at the maximum possible velocity (which we call v max ). We call this strategy the "direct" method or algorithm.
We will base the choice of the new velocity vector for each UAV on the future state of the airspace. Considering the future position of a set of UAVs, we will decide which ones should modify their velocity vector and which ones should not to reach a conflict-free solution. PCAN employs two strategies to modify the velocities: (1) modify the direction of the velocity vector with a certain angle and (2) modify the modulus of the velocity vector, reducing the flight velocity below the maximum, but maintaining the direction and sense. Both modifications will be made on the direct velocity, trying to alter it as little as possible.
From now on, "airspace" refers to the region of the space where UAVs fly. We assume that all UAVs fly at the same altitude. Consequently, we tackle conflict detection in the 2D plane, working on an area where UAVs move from their initial position to their final (or destination) position. Nevertheless, it is possible to run the algorithm in parallel by using multiple layers at different heights. The criteria for the layer choice are left to the U-Space service provider and can be very varied (UAV heading, priority, maneuverability, airspace congestion...). The only requirement is that a UAV must participate in the conflict management of the layers below the deployment one during vertical takeoff and landing.
Before going into detail, we will introduce some starting definitions and assumptions.

UAV Dynamic Model
A comprehensive model of the dynamic behavior of a quadcopter-type UAV can be found in previous work [27]. In the present paper, it is sufficient to apply the simplified 2D model described below. Let x = p . p be the dynamic state of a UAV, where vectors p = p x p y ∈ P and . p = v x v y ∈ V, represent its current position and velocity, respectively. Let v = v x v y ∈ V be the commanded (or desired) velocity for that UAV. Then, its dynamic behavior is modeled by the following first order system: p /τ, being τ the response time of the system (time employed to achieve the 63% of the desired value). Note that, the first-order system applies to the velocity of the UAV. However, due to it has been described with respect to its position, a second-order derivative appears. Figure 1 shows an example.
Under these conditions, the position error experienced by the UAV due to the delay in following the commanded velocity is e = τ v − . p 0 . The maximum error would occur in a scenario with both vectors at maximum speed and opposite direction: e max = 2τv max . As we will see later, the conflict detection mechanism will consider the error due to the dynamic behavior of the UAVs, consequently increasing their safety radius.
Let A = {a 1,t , a 2,t . . . a n,t } be the set of aircraft in the airspace at time t. Each element a = [x R] ∈ A represents the status vector of a UAV, including its dynamic behavior and a route composed of several waypoints to its destination. Let R = {w 1 , w 2 . . . w m } be an ordered sequence of waypoints w i ∈ P.
Without loss of generality, we assume that all UAVs behave homogeneously, initially flying at a predetermined maximum velocity v max and they are considering a safety radius r (which delimits their "protected zone"). Figure 2 shows an example. As we will see later, the conflict detection mechanism will consider the error due to the dynamic behavior of the UAVs, consequently increasing their safety radius.
(a) (b) } be an ordered sequence of waypoints ∈ ℙ. Without loss of generality, we assume that all UAVs behave homogeneously, initially flying at a predetermined maximum velocity and they are considering a safety radius (which delimits their "protected zone"). Figure 2 shows an example. With this definition, the optimal behavior for each UAV in ideal airspace (without restrictions due, for example, to the presence of buildings or no-fly zones) would be to fly in a straight line towards the next waypoint in its route, at a velocity (what we have called the direct method). Table 1 shows the implementation of this navigation system, which is executed periodically. If the UAV would reach the current waypoint before the next execution (02), it switches to the next waypoint in route (03-04). If multiple UAVs behave in this way in the same airspace, they may collide in flight.  Without loss of generality, we assume that all UAVs behave homogeneously, initially flying at a predetermined maximum velocity and they are considering a safety radius (which delimits their "protected zone"). Figure 2 shows an example. With this definition, the optimal behavior for each UAV in ideal airspace (without restrictions due, for example, to the presence of buildings or no-fly zones) would be to fly in a straight line towards the next waypoint in its route, at a velocity (what we have called the direct method). Table 1 shows the implementation of this navigation system, which is executed periodically. If the UAV would reach the current waypoint before the next execution (02), it switches to the next waypoint in route (03-04). If multiple UAVs behave in this way in the same airspace, they may collide in flight. With this definition, the optimal behavior for each UAV in ideal airspace (without restrictions due, for example, to the presence of buildings or no-fly zones) would be to fly in a straight line towards the next waypoint in its route, at a velocity v max (what we have called the direct method). Table 1 shows the implementation of this navigation system, which is executed periodically. If the UAV would reach the current waypoint before the next execution (02), it switches to the next waypoint in route (03-04). If multiple UAVs behave in this way in the same airspace, they may collide in flight.  From now on, the expression "navigation computation" refers to the operation consisting in assigning to each UAV in the airspace a velocity with which it must move in order not to cause any conflict. A "valid velocity" for a UAV is a velocity that does not produce any conflict between it and the rest of the UAVs in the airspace. "Final velocity" is the valid velocity that each UAV will use to move after the execution of the navigation computation.
UAV movement through airspace occurs in given units of time (e.g., seconds). We assume that all displacements start simultaneously, but not all of them end at the same time. This will depend on the initial and final positions of the UAVs, their velocities during the travel, and the route they take. Navigation computation also occurs in units of time, with the same unit as UAV movement. We assume that navigation is executed every t units of travel. We call this parameter t nav . In other words, t nav defines how many time units a UAV can travel with the same navigation computation.

Conflict Prediction Mechanism
In this subsection, we will first detail a mechanism to check for the occurrence of a future conflict between two UAVs. From this, we will describe how to predict whether a given UAV will be involved in a conflict with any other UAV present in the airspace.
If we assume that a UAV a i maintains its velocity constant, then its position at a future time t (relative to the current time) is provided by Equation (1): Two UAVs, a 1 and a 2 , present a conflict if |p 1 p 2 | < 2r. The system of equations shown in Equation (2) allows us to determine the existence of a conflict between them, obtaining the instants and positions of its beginning and end.
(2) By solving this system, we obtain the expression shown in Equation (3): which results in a second-degree equation as a function of t. After being solved, several situations may occur: • No real root is obtained. This means that the protected zones of both UAVs do not contact at any time. No conflict situation arises and, therefore, there are not potential collisions. • A real root is obtained. This means that the protected zones of both UAVs contact each other without overlapping, which does not generate a conflict situation either. • Two real roots, t 1 and t 2 are obtained. These values indicate the instants of the beginning and end of the conflict. During this time interval, the protected zones of both UAVs are partially overlapped. If both roots are positive, the conflict is predicted in the future. If only one root is negative, a 1 and a 2 are currently in a conflict situation. Finally, if both roots are negative, the conflict was resolved in the past, or it never occurred (the navigation algorithm prevented it).
As an example, Figure 3 shows two UAVs whose trajectories lead to a conflict in the future. Finally, if both roots are negative, the conflict was resolved in the past, or it never occurred (the navigation algorithm prevented it).
As an example, Figure 3 shows two UAVs whose trajectories lead to a conflict in the future.  Table 2 presents the implementation of the conflict prediction operation between two UAVs. If the ConflictPrediction function is given a velocity, the prediction is performed using this value (lines 04-06). Otherwise, the velocity previously assigned to the UAV will be used (03). We will use this function later to check whether a velocity is suitable for resolving a conflict before assigning it to a UAV. Then, the system of equations is implemented in the variable (07-09) and its roots are computed (10). If any real root results (11), will be the smallest of them (12). If no real root exists, would be assigned to ∞ (14).
To predict whether a particular UAV will be involved in a conflict, we can use the conflict prediction operation just described (see Table 2), using as arguments the UAV  Table 2 presents the implementation of the conflict prediction operation between two UAVs. If the ConflictPrediction function is given a velocity, the prediction is performed using this value (lines 04-06). Otherwise, the velocity previously assigned to the UAV will be used (03). We will use this function later to check whether a velocity is suitable for resolving a conflict before assigning it to a UAV. Then, the system of equations is implemented in the variable t c (07-09) and its roots are computed (10). If any real root results (11), t c will be the smallest of them (12). If no real root exists, t c would be assigned to ∞ (14).
To predict whether a particular UAV will be involved in a conflict, we can use the conflict prediction operation just described (see Table 2), using as arguments the UAV under study and the rest of the UAVs present in the airspace. This check can be implemented through a loop, which will be interrupted as soon as the first conflict is detected, to replace the velocity of the UAV analyzed by a new valid velocity. If, on the other hand, this loop concludes without detecting any conflict, the current UAV velocity will be replaced by the direct velocity to the destination. We will discuss all this in more depth in later subsections. Table 3 implements this behavior using a Boolean function. The main loop (01-08) checks if there is a conflict between each pair of UAVs (05). If so, a true value is immediately returned (06). Obviously, we must omit this check for the current UAV (02-04). If, after considering all the UAVs in the airspace, no conflict has been found, a false value is provided (09), indicating that there is no conflict with the UAV under study. Table 3. Conflict prediction between a UAV and the rest.
return true 07 end if 08 end do 09 return false

Valid Velocity Computation
After predicting a conflict, PCAN replaces the direct velocity of the UAV concerned by a new one that guarantees that no conflicts will occur. Each UAV must have a final velocity before moving during the following t nav time units, which guarantees that it will not collide with any other UAV in the airspace before the next execution of the navigation computation.
A set of candidate velocities "close" to the direct velocity are generated to obtain this final velocity, and these new velocities are checked for conflicts. Suppose any of them manages to resolve all the conflicts between the UAV under study and the rest of the UAVs in the airspace. In that case, this velocity will be a valid velocity, and it will be taken as the final velocity for the UAV. If none of these velocities leads to a conflict-free scenario, the algorithm proceeds to generate a new set of candidate velocity vectors. This process continues until a final velocity is found or a preset limit of iterations is exceeded. If, after exhausting all iterations and discarding all candidate velocity vectors, any valid velocity has been obtained, a decision is made to assign the UAV a zero velocity (assuming it is a rotary-wing UAV). This is done because, in this case, the airspace is very crowded, and no velocity in any direction would allow the UAV to proceed while guaranteeing a conflict-free scenario without increasing the time or distance traveled too much. The UAV should stop moving forward to leave its path clear for other UAVs to move in this situation.
To generate all these candidate velocities, a loop can be used that generates new velocities from the previous ones. In its first iteration, velocities generated are based on the direct one (v d ). In each iteration, three candidate velocities are generated. One of these velocities will have the same direction as the previous one but a lower modulus. The other two velocities will have a different direction but the same modulus. Through the coe f parameter, we can set how much the velocity is reduced when we are varying its modulus. In the same way, through the α parameter (angular displacement), we can control how much the direction will vary when we are obtaining two new velocities with different directions. These two velocities will have the same amount of direction variation, but one will be varied clockwise and the other counterclockwise.
We can see an example in Figure 4. As stated, the iterative process starts with v d . The three new candidate velocities generated in iteration i are called v j i , where j distinguishes between them. In the figure, we can see that v 1 1 and v 2 1 (orange color) present the same angular displacement with respect to v d , but in different directions. On the other hand, v 3 1 , which is the modified velocity in modulus, follows the direction of v d . If none of these three candidate velocities were valid, we would proceed to generate new ones (v j 2 ), following the same methodology, but this time starting from the velocities generated in the previous loop iteration.
three new candidate velocities generated in iteration are called , where distinguishes between them. In the figure, we can see that and (orange color) present the same angular displacement with respect to , but in different directions. On the other hand, , which is the modified velocity in modulus, follows the direction of . If none of these three candidate velocities were valid, we would proceed to generate new ones ( ), following the same methodology, but this time starting from the velocities generated in the previous loop iteration. Figure 4. Searching for a valid velocity from the direct one. Above, a UAV with direct velocity ( ) flying from the initial position ( ) to the destination ( ). In the center, three candidate velocity vectors have been generated in the first iteration. Below, new candidate velocity vectors have been generated, based on those of the previous iteration (if they did not lead to a conflict-free scenario). Table 4 provides the function that generates valid velocities for a UAV. As described above, three different velocities will be evaluated, which are initialized using the direct velocity (02). The computation is enclosed in a loop (03-17), which iterates until a valid velocity is found or a predetermined maximum number of iterations ( ) is reached. At each iteration, three modified velocities are obtained from the previous velocity vector. We use the parameter (a value less than 1) to modify the modulus of the velocity vector (04). The change in direction is performed by the Veer function (05-06), described below. These three candidate velocities are then checked for conflicts (07, 10, 13). If any of them did not cause any conflict, the corresponding velocity vector is returned (08, 11,14). If the loop concludes without finding any valid velocity, a null velocity is provided (18). Table 5 details the behavior of the Veer function. The argument is the velocity vector to vary. The argument defines the angular displacement between iterations (expressed in radians). In our study, this value has been assigned to the maximum allowed angular displacement divided by the maximum number of iterations ( ). The direction argument admits the values and , and represents the direction of the angular displacement (01-03). In (04), a new velocity vector with an angular displacement in a counterclockwise or clockwise direction, respectively, is obtained. Figure 4. Searching for a valid velocity from the direct one. Above, a UAV with direct velocity (v d ) flying from the initial position (p) to the destination (w). In the center, three candidate velocity vectors have been generated in the first iteration. Below, new candidate velocity vectors have been generated, based on those of the previous iteration (if they did not lead to a conflict-free scenario). Table 4 provides the function that generates valid velocities for a UAV. As described above, three different velocities will be evaluated, which are initialized using the direct velocity (02). The computation is enclosed in a loop (03-17), which iterates until a valid velocity is found or a predetermined maximum number of iterations (m) is reached. At each iteration, three modified velocities are obtained from the previous velocity vector. We use the coe f parameter (a value less than 1) to modify the modulus of the velocity vector (04). The change in direction is performed by the Veer function (05-06), described below. These three candidate velocities are then checked for conflicts (07, 10, 13). If any of them did not cause any conflict, the corresponding velocity vector is returned (08, 11,14). If the loop concludes without finding any valid velocity, a null velocity is provided (18). Table 5 details the behavior of the Veer function. The v argument is the velocity vector to vary. The α argument defines the angular displacement between iterations (expressed in radians). In our study, this value has been assigned to the maximum allowed angular displacement divided by the maximum number of iterations (m). The direction argument admits the values clockwise and counterclockwise, and represents the direction of the angular displacement (01-03). In (04), a new velocity vector with an angular displacement in a counterclockwise or clockwise direction, respectively, is obtained.

Navigation Computation. The PCAN Algorithm
We have seen the above two mechanisms to predict a conflict between a UAV and the rest of the UAVs in the airspace and provide a new valid velocity for a UAV in that situation. This subsection will describe the general mechanism to compute the final velocity for all the UAVs in the airspace (what we have called navigation computation). Table 6 details the general behavior of the PCAN algorithm. First, all UAVs are assigned direct velocity to the destination (02-06). Then, for each UAV we check if it produces any conflict with the airspace (08). If there are no conflicts, we assign the direct velocity as the final velocity for the considered UAV (15). However, if a conflict is predicted, a valid velocity is computed for that UAV (09).
if not AirspaceConflicts(v 1 , a 1 , A) 08 return if not AirspaceConflicts(v 2 , a 1 , A) 11 return if not AirspaceConflicts(v 3 , a 1 , A) 14 return end while 18 return v = 0 If the ValidVelocity function returns a null velocity (10), we must recalculate the final velocities for the previously processed UAVs (11). If the velocity is not null, we assign the obtained velocity as the final velocity for the UAV (15).
It is important to clarify that the process must be carried out in sequence, processing one UAV after another until all of them have a final velocity. When a UAV is about to calculate its final velocity, conflict prediction will be made using the final velocities of the UAVs that already have calculated their final velocity and the direct velocities of the UAVs that have not been processed yet. In the case of the first UAV processed, conflict prediction will be made from the direct velocities of the rest of the UAVs, while the last UAV processed will perform the prediction from the final velocities of the rest of the UAVs in the airspace. Since a final velocity for a UAV will only be accepted if it does not lead to a conflict with any other UAV, the order in which the final velocities are obtained will affect the solution for each of the UAVs. The UAV whose final velocity is calculated first will have to fight more conflicts than the last UAV processed, which will not have to fight any conflict, since the previous UAVs treated have already prevented the conflict with it.
For this work, the order in which the UAVs are processed is given by their number within the airspace (ordinal within the set). This implies a priority which, although not chosen, is necessary since there must be an order. However, as future work, we plan to implement a priority system based, for example, on UAV categories (as established by U-space), on the urgency of the service provided, or on any other criteria. In this way, those UAVs with higher priority will minimally modify the optimal trajectory offered by the direct method, while those UAVs with lower priority will have to avoid a greater number of conflicts, obtaining for them a solution farther away from the optimal one. Finally, note that assigning a null velocity to a UAV implies a higher computational cost because it invalidates the final velocities of UAVs that have been processed before it. By assigning a null velocity, the prediction made by another UAV for which the final velocity was already calculated is invalidated and must be performed again. This process only occurs in really crowded airspaces.

Airspace Bounding
As mentioned before, valid velocities are generated from the airspace state. However, in airspaces with a high density of UAVs, which is expected to happen in UAM scenarios, the computational cost of calculating these valid velocities increases. Moreover, considering UAVs far enough away from a given UAV not to cause a conflict with it in the short term does not make sense and leads to worse solutions. For this reason, an immediate improvement of PCAN would consist of analyzing only the portion of the airspace containing those UAVs that could collide with the UAV understudy before the next execution of the navigation computation.
In short, we define a boundary radius (BR), establishing a circular region around the UAV so that only the UAVs inside that area will be considered for the calculation of the valid velocity. This radius should be greater than the threshold (BR min ) defined by the worst possible situation, consisting of two UAVs flying in opposite directions, according to: (4) Figure 5 shows an example. Given the airspace A = {a 1 . . . a 6 }, we proceed to calculate a valid velocity for a 1 . Consequently, only the UAVs within the circular region defined by BR, plotted in green, are considered for the valid velocity calculation. Thus, a 5 and a 6 do not influence the final velocity of a 1 .
to: Figure 5 shows an example. Given the airspace = { … }, we proceed to calculate a valid velocity for . Consequently, only the UAVs within the circular region defined by , plotted in green, are considered for the valid velocity calculation. Thus, and do not influence the final velocity of .

Performance Evaluation
In this section, we present the performance evaluation results of the proposed algorithm (PCAN), with the airspace bounding improvement described in Subsection 3.5. This evaluation has been carried out using simulation techniques.
The simulation tool employed has been developed in Matlab R2020b [28] by making use of the object-oriented programming language offered by this platform. The simulator allows the configuration of the parameters of the proposed algorithm and the generation of random scenarios of any size. It also provides result reports, generating plots such as the ones shown below.
In total, more than 150,000 scenarios have been simulated, with different configurations. As we will see next, the PCAN algorithm has prevented the occurrence of conflicts in all cases.

Two-UAV Study
First, we will analyze the behavior of PCAN when = { , }. Figure 6 shows the different scenarios evaluated. Given a circumference of 4 km in diameter, in all the scenarios starts its flight from the west (point ), crosses the circumference passing

Performance Evaluation
In this section, we present the performance evaluation results of the proposed algorithm (PCAN), with the airspace bounding improvement described in Section 3.5. This evaluation has been carried out using simulation techniques.
The simulation tool employed has been developed in Matlab R2020b [28] by making use of the object-oriented programming language offered by this platform. The simulator allows the configuration of the parameters of the proposed algorithm and the generation of random scenarios of any size. It also provides result reports, generating plots such as the ones shown below.
In total, more than 150,000 scenarios have been simulated, with different configurations. As we will see next, the PCAN algorithm has prevented the occurrence of conflicts in all cases.

Two-UAV Study
First, we will analyze the behavior of PCAN when A = {a 1 , a 2 }. Figure 6 shows the different scenarios evaluated. Given a circumference of 4 km in diameter, in all the scenarios a 1 starts its flight from the west (point p), crosses the circumference passing through its center, and ends at the opposite point of the circumference (w). On the other hand, a 2 starts its flight from a different position in each scenario. In particular, we have considered 18 different initial positions for a 2 (p 1 . . . p 18 ) resulting in 18 relative angles between both UAVs (between 0 • and 170 • , in 10 • intervals). Each initial position (p 1 . . . p 18 ) has an associated destination (w 1 . . . w 18 ), which is also reached by passing through the center of the circumference. In each configuration, a 1 and a 2 start the flight simultaneously, with v max = 13.9 m/s (50 km/h), finally colliding at the center of the circumference if the navigation algorithm does not prevent it. We have set r = 45 m, t nav = 3 s and τ = 0.3 s, and PCAN parameters as BR = ∞, m = 50, coe f = 0.95, and α = π m rad. Moreover, the safety radius r has been increased 8.34 m, as Figure 1 shown explains. between both UAVs (between 0° and 170°, in 10° intervals). Each initial position ( … ) has an associated destination ( … ), which is also reached by passing through the center of the circumference. In each configuration, and start the flight simultaneously, with = 13.9 m/s (50 km/h), finally colliding at the center of the circumference if the navigation algorithm does not prevent it. We have set = 45 m, = 3 s and = 0.3 s, and PCAN parameters as = ∞, = 50, = 0.95, and = . Moreover, the safety radius has been increased 8.34 m, as Figure 1 shown explains.  Figure 7 shows a scenario number corresponding to each of the 18 configurations described above on its abscissa axis. The ordinate axis represents the number of conflicts produced in each scenario. The two series represent the results obtained by the direct and PCAN algorithms. As expected, the direct algorithm produces conflicts (leading to collisions) in all scenarios. In contrast, the PCAN algorithm manages to avoid them in all the scenarios analyzed. As detailed, PCAN avoids conflicts by causing UAVs to deviate from the optimal trajectory to the destination or by reducing their velocity (but maintaining the optimal trajectory). This can lead to increases in the distance traveled by the UAVs and the time required to reach their destinations. Figure 8a shows the average distance traveled by the  Figure 7 shows a scenario number corresponding to each of the 18 configurations described above on its abscissa axis. The ordinate axis represents the number of conflicts produced in each scenario. The two series represent the results obtained by the direct and PCAN algorithms. As expected, the direct algorithm produces conflicts (leading to collisions) in all scenarios. In contrast, the PCAN algorithm manages to avoid them in all the scenarios analyzed.
between both UAVs (between 0° and 170°, in 10° intervals). Each initial position ( 1 … 18 ) has an associated destination ( 1 … 18 ), which is also reached by passing through the center of the circumference. In each configuration, 1 and 2 start the flight simultaneously, with = 13.9 m/s (50 km/h), finally colliding at the center of the circumference if the navigation algorithm does not prevent it. We have set = 45 m, = 3 s and = 0.3 s, and PCAN parameters as = ∞, = 50, = 0.95, and = . Moreover, the safety radius has been increased 8.34 m, as Figure 1 shown explains.  Figure 7 shows a scenario number corresponding to each of the 18 configurations described above on its abscissa axis. The ordinate axis represents the number of conflicts produced in each scenario. The two series represent the results obtained by the direct and PCAN algorithms. As expected, the direct algorithm produces conflicts (leading to collisions) in all scenarios. In contrast, the PCAN algorithm manages to avoid them in all the scenarios analyzed. As detailed, PCAN avoids conflicts by causing UAVs to deviate from the optimal trajectory to the destination or by reducing their velocity (but maintaining the optimal trajectory). This can lead to increases in the distance traveled by the UAVs and the time required to reach their destinations. Figure 8a shows the average distance traveled by the As detailed, PCAN avoids conflicts by causing UAVs to deviate from the optimal trajectory to the destination or by reducing their velocity (but maintaining the optimal trajectory). This can lead to increases in the distance traveled by the UAVs and the time required to reach their destinations. Figure 8a shows the average distance traveled by the UAVs, expressed in meters. The direct algorithm indicates the minimum distance between the initial and the final position. The increase due to the deviations made by PCAN to avoid conflict is negligible in this plot. In the case of the BBCA algorithm, an increase in the distance traveled is clearly visible.
To better analyze these results, Figure 8b shows the increase in the distance traveled by the UAVs when using PCAN and BBCA relative to the minimum distance provided by the direct algorithm. As we can see, the penalty in the distance for BBCA is appreciable in some scenarios. In PCAN, this increase is negligible again (about 0.11% in the worst case, but about 0.06% in most of the scenarios). the distance traveled is clearly visible.
To better analyze these results, Figure 8b shows the increase in the distance traveled by the UAVs when using PCAN and BBCA relative to the minimum distance provided by the direct algorithm. As we can see, the penalty in the distance for BBCA is appreciable in some scenarios. In PCAN, this increase is negligible again (about 0.11% in the worst case, but about 0.06% in most of the scenarios). We have also analyzed the impact on the flight time due to the avoidance maneuvers performed by the UAVs. Figure 9a shows the average flight time of UAVs, expressed in seconds. In the case of PCAN, it can be observed that the increase in time is greater than in the distance. This is because PCAN tends to modify the velocity modulus (instead of the UAV trajectory) to avoid the conflict, thereby slightly increasing the time the UAV remains in flight.
Analogous to Figure 8b, Figure 9b shows the relative increase (with respect to the direct algorithm) in the UAV flight time when using PCAN and BBCA. For PCAN, we can observe an increase of less than 2% in all the scenarios. One more time, the penalty introduced by BBCA in some scenarios is notable.  We have also analyzed the impact on the flight time due to the avoidance maneuvers performed by the UAVs. Figure 9a shows the average flight time of UAVs, expressed in seconds. In the case of PCAN, it can be observed that the increase in time is greater than in the distance. This is because PCAN tends to modify the velocity modulus (instead of the UAV trajectory) to avoid the conflict, thereby slightly increasing the time the UAV remains in flight.
case, but about 0.06% in most of the scenarios). We have also analyzed the impact on the flight time due to the avoidance maneuvers performed by the UAVs. Figure 9a shows the average flight time of UAVs, expressed in seconds. In the case of PCAN, it can be observed that the increase in time is greater than in the distance. This is because PCAN tends to modify the velocity modulus (instead of the UAV trajectory) to avoid the conflict, thereby slightly increasing the time the UAV remains in flight.
Analogous to Figure 8b, Figure 9b shows the relative increase (with respect to the direct algorithm) in the UAV flight time when using PCAN and BBCA. For PCAN, we can observe an increase of less than 2% in all the scenarios. One more time, the penalty introduced by BBCA in some scenarios is notable.  Analogous to Figure 8b, Figure 9b shows the relative increase (with respect to the direct algorithm) in the UAV flight time when using PCAN and BBCA. For PCAN, we can observe an increase of less than 2% in all the scenarios. One more time, the penalty introduced by BBCA in some scenarios is notable.

Multi-UAV Study
Next, we study the behavior of PCAN when it is used in airspaces with multiple UAVs. We have considered a 5 × 5 km flight region and with different UAV densities. In particular, |A| = {10, 20,30,40,50,60,70,80,90, 100}. Figure 10 shows an example. For each density, 24 random route configurations have been generated for the UAVs, considering that the initial and the final position cannot be less than 100 m from the edges of the region and less than 1 km away. Each route configuration has been simulated with the direct, BBCA, APF, and PCAN algorithms, resulting in 960 runs. In all cases, we have set v max = 50 km/h, r = 40 m, t nav = 3 s and τ = 0.3 s, and PCAN parameters as BR = 3v max t nav + 2r, m = 50, coe f = 0.95, and α = π m rad. Moreover, the safety radius r has been increased 8.34 m, as Figure 1 explains.
For each density, 24 random route configurations have been generated for the UAVs, considering that the initial and the final position cannot be less than 100 from the edges of the region and less than 1 km away. Each route configuration has been simulated with the direct, BBCA, APF, and PCAN algorithms, resulting in 960 runs. In all cases, we have set = 50 km/h , = 40 m , = 3 s and = 0.3 s , and PCAN parameters as = 3 + 2 , = 50, = 0.95, and = . Moreover, the safety radius has been increased 8.34 m, as Figure 1 explains. Figure 10. Random scenario with 80 UAVs and PCAN. Each number represents the initial position of a UAV; each circle represents the current position of the UAV and its safety radius; each line represents the trajectory followed by the UAV. Units on the X and Y axes are meters.
The abscissa axis in the plot of Figure 11a indicates the number of UAVs in the flight region. In contrast, the ordinate axis shows the mean number of conflicts produced and the corresponding standard deviation. For the direct method, conflicts increase exponentially with the number of UAVs in the region. The APF algorithm manages to resolve isolated conflicts, but its performance decreases as the UAV density increases. BBCA significantly reduces the number of conflicts but does not eliminate them. Finally, as expected, PCAN successfully resolves all the conflicts. Figure 11b shows the same results, but now from the point of view of the percentage of solved conflicts. We can clearly see that PCAN and BBCA outperform APF. Obviously, the direct method is not plotted since it does not avoid any conflict. The abscissa axis in the plot of Figure 11a indicates the number of UAVs in the flight region. In contrast, the ordinate axis shows the mean number of conflicts produced and the corresponding standard deviation. For the direct method, conflicts increase exponentially with the number of UAVs in the region. The APF algorithm manages to resolve isolated conflicts, but its performance decreases as the UAV density increases. BBCA significantly reduces the number of conflicts but does not eliminate them. Finally, as expected, PCAN successfully resolves all the conflicts. Figure 11b shows the same results, but now from the point of view of the percentage of solved conflicts. We can clearly see that PCAN and BBCA outperform APF. Obviously, the direct method is not plotted since it does not avoid any conflict.  Figure 12a shows the distance traveled by the UAVs, in absolute terms, as a function of the UAV density. As expected, when applying the direct method, the distance traveled by each UAV does not vary with its number, while the standard deviation decreases. On the contrary, the detours introduced by BBCA increase the distance linearly. Finally, PCAN also increases the distance traveled by the UAVs, but to a lesser extent since it  Figure 12a shows the distance traveled by the UAVs, in absolute terms, as a function of the UAV density. As expected, when applying the direct method, the distance traveled by each UAV does not vary with its number, while the standard deviation decreases. On the contrary, the detours introduced by BBCA increase the distance linearly. Finally, PCAN also increases the distance traveled by the UAVs, but to a lesser extent since it resolves the conflicts more efficiently. Figure 12a shows the distance traveled by the UAVs, in absolute terms, as a function of the UAV density. As expected, when applying the direct method, the distance traveled by each UAV does not vary with its number, while the standard deviation decreases. On the contrary, the detours introduced by BBCA increase the distance linearly. Finally, PCAN also increases the distance traveled by the UAVs, but to a lesser extent since it resolves the conflicts more efficiently.
In Figure12b, we can see the relative increase in the distance traveled by the UAVs regarding the direct algorithm. We can observe the linearity in this increase, which in the case of PCAN progresses from 0% to about 4% for the 100-UAV configuration. If we compare PCAN to BBCA, we can see that the increase in distance is marginal, even avoiding all conflicts. APF performs slightly better than PCAN, but at the cost of not being conflictfree.
(a) (b)  Finally, Figure 13b compares the flight time for BBCA, APF, and PCAN with respect to the direct algorithm. In the case of PCAN, in the highest density scenarios, the flight time is increased by about 6%, outperforming, in any case, the other two techniques. In Figure 12b, we can see the relative increase in the distance traveled by the UAVs regarding the direct algorithm. We can observe the linearity in this increase, which in the case of PCAN progresses from 0% to about 4% for the 100-UAV configuration. If we compare PCAN to BBCA, we can see that the increase in distance is marginal, even avoiding all conflicts. APF performs slightly better than PCAN, but at the cost of not being conflict-free. Figure 13a shows the flight time of UAVs, in absolute terms, as a function of the UAV density. For the direct method, flight time does depend on the number of UAVs. As shown, PCAN offers better results than the rest of conflict management algorithms.

Computation Time
To conclude the analysis, we studied the time used by the direct, PCAN, BBCA, and APF algorithms to offer their solutions. Figure 14 shows the time consumed to produce an output for every UAV each time the navigation computation is executed, expressed in milliseconds. Results were generated with an Intel i9-10900KF@3.7GHz processor, using one core.
As shown, PCAN works very well in low-density scenarios exhibiting better behavior than BBCA from 10 to 50 UAVs. As the density of UAVs increases, PCAN requires more time to solve all the potential conflicts produced, which were shown in Figure 11a, series direct. In any case, the computation time does not represent a dramatic bottleneck. Finally, Figure 13b compares the flight time for BBCA, APF, and PCAN with respect to the direct algorithm. In the case of PCAN, in the highest density scenarios, the flight time is increased by about 6%, outperforming, in any case, the other two techniques.

Computation Time
To conclude the analysis, we studied the time used by the direct, PCAN, BBCA, and APF algorithms to offer their solutions. Figure 14 shows the time consumed to produce an output for every UAV each time the navigation computation is executed, expressed in milliseconds. Results were generated with an Intel i9-10900KF@3.7GHz processor, using one core.
respect to the direct algorithm.

Computation Time
To conclude the analysis, we studied the time used by the direct, PCAN, BBCA, and APF algorithms to offer their solutions. Figure 14 shows the time consumed to produce an output for every UAV each time the navigation computation is executed, expressed in milliseconds. Results were generated with an Intel i9-10900KF@3.7GHz processor, using one core.
As shown, PCAN works very well in low-density scenarios exhibiting better behavior than BBCA from 10 to 50 UAVs. As the density of UAVs increases, PCAN requires more time to solve all the potential conflicts produced, which were shown in Figure 11a, series direct. In any case, the computation time does not represent a dramatic bottleneck. For the 100-UAV configuration, the time employed is less than 0.9 ms. Note that these are extremely dense scenarios in which avoiding all conflict is complex.

Conclusions and Future Works
This paper proposes the PCAN algorithm for conflict-free navigation among a set of UAVs flying over urban airspace according to a set of predetermined flight plans. In each run of the algorithm, assuming that all UAVs are heading in a straight line and at maximum velocity towards their destinations, the algorithm predicts future conflicts between them and proposes modifications in their velocities to prevent their occurrence. The analysis carried out shows that PCAN results, in the worst case, in an increase in the distance traveled by the UAVs of about 4% and an increase in flight time of approximately 6%. This makes our conflict avoidance proposal suitable for the deployment of tactical conflict resolution services in the framework of the future U-Space UAM space. As shown, PCAN works very well in low-density scenarios exhibiting better behavior than BBCA from 10 to 50 UAVs. As the density of UAVs increases, PCAN requires more time to solve all the potential conflicts produced, which were shown in Figure 11a, series direct. In any case, the computation time does not represent a dramatic bottleneck. For the 100-UAV configuration, the time employed is less than 0.9 ms. Note that these are extremely dense scenarios in which avoiding all conflict is complex.

Conclusions and Future Works
This paper proposes the PCAN algorithm for conflict-free navigation among a set of UAVs flying over urban airspace according to a set of predetermined flight plans. In each run of the algorithm, assuming that all UAVs are heading in a straight line and at maximum velocity towards their destinations, the algorithm predicts future conflicts between them and proposes modifications in their velocities to prevent their occurrence. The analysis carried out shows that PCAN results, in the worst case, in an increase in the distance traveled by the UAVs of about 4% and an increase in flight time of approximately 6%. This makes our conflict avoidance proposal suitable for the deployment of tactical conflict resolution services in the framework of the future U-Space UAM space.
As future work, we plan to revise the decision-making of PCAN to reduce the overhead involved in conflict avoidance. Instead of performing an intensive search of modified velocities, the algorithm can start by exploring the most promising velocities according to the airspace state. Moreover, as discussed in Section 3.4, we can implement a priority system based on categories. In this way, those UAVs with higher priority will see their trajectories modified to a lesser extent. Finally, we must consider the existence of geofences in the airspace that define no-fly zones.