Design of Convolutional Extreme Learning Machines for Vision-Based Navigation Around Small Bodies

Deep learning architectures such as convolutional neural networks are the standard in computer vision for image processing tasks. Their accuracy however often comes at the cost of long and computationally expensive training, the need for large annotated datasets, and extensive hyper-parameter searches. On the other hand, a different method known as convolutional extreme learning machine has shown the potential to perform equally with a dramatic decrease in training time. Space imagery, especially about small bodies, could be well suited for this method. In this work, convolutional extreme learning machine architectures are designed and tested against their deep-learning counterparts. Because of the relatively fast training time of the former, convolutional extreme learning machine architectures enable efficient exploration of the architecture design space, which would have been impractical with the latter, introducing a methodology for an efficient design of a neural network architecture for computer vision tasks. Also, the coupling between the image processing method and labeling strategy is investigated and demonstrated to play a major role when considering vision-based navigation around small bodies.

surface variations are only due to morphological characteristics (i.e. craters, boulders, etc.) which only vary under illumination conditions. This domain is fundamentally simpler than that of typical computer vision applications in urban environments, which need to account for a large variety of commonly used objects, where deep architectures are the state of the art [2]. Moreover, previous findings in [9] hint to the fact that the filtering capabilities CNNs on images of a small body are critical in pushing the performance compared to traditional methods.
Pivoting on these prior works and the existing gaps in the literature, in this paper the authors attempt to answer three fundamental questions, namely: 1) Although CNN are superior over complex scenery, can simpler methods perform better when it comes to analyzing images of small bodies?; 2) As suggested in [6], can a methodology be developed to bootstrap the training of CNN exploiting the capability of CELM to identify the most promising architectures?; 3) Which is the best labeling strategy and reference frame to be used when considering a visual-based navigation application around a small body?
In this work, CELM is thus investigated as a possible alternative to CNN for autonomous vision-based navigation systems around small bodies. This is done with an extensive analysis considering 4 different small body shapes and 5 different labeling strategies, for a total of 20 scenarios, each of which is examined with 4 different IP methods, resulting in tens of thousands of different architectures explored. Such efficient exploration of the architecture design space is possible thanks to the extremely fast training time of CELM, which is orders of magnitude faster than that of CNN [7]. It is also demonstrated that the labeling strategy and reference frame play a crucial role during training and significantly affect the performance of the methods considered. We believe that CELM can represent, with the proper labeling strategy, a promising alternative in space imagery, especially the one related to small bodies.
The rest of the paper is organized as follows. In Section II the three pillars sustaining the methodology are discussed in detail. These are the dataset generation, the preprocessing pipeline of the image-label pairs, and the description of the IP methods developed and used in this paper. The performance of these methods is then compared extensively in Section III while conclusions, and future works are discussed in Section IV.

A. Dataset generation
In this section, the dataset generation procedure is illustrated. Four small bodies are considered: 65803 Didymos (D), 103P/Hartley (H), 21 Lutetia (L) and 67P/Churyumov-Gerasimenko (P). Using the same methodology illustrated in [10], artificial morphological features such as boulders and craters are inserted into the shape models of these bodies.
For each body, a total of 17500 grayscale images are rendered in Blender * . These are split into training, validation, and test sets respectively made of 7500, 5000, and 5000 images. As illustrated in Fig. 1, these sets are randomly distributed across a cloud of points around the body in what is defined as the reference frame. This reference frame is an inertially fixed reference frame, centered on the Center of Mass (CoM) of the body, with the -axis oriented towards the projection of the Sun in the body equatorial plane and the -axis as the north pole of the body. Another reference frame that is used in this paper is the one, whose axes are fixed to the body surface and are obtained from a rigid rotation around the -axis of the reference frame. The points in Fig. 1  Since each dataset is composed of the same set of points across different bodies, uniform scaling is applied to each one to make sure that different body sizes are all filling the Field Of View (FOV) around 5 km. To generate the 8-bit grayscale images, a camera with a 10 ×10 deg FOV and a sensor of 1024 ×1024 pixel is considered.
All images are rendered assuming ideal pointing towards the CoM. From an optical-navigation point of view, images can be represented by different sets of labels. Five of these are considered in this work to investigate the influence of this choice on the performance of the IP methods. The first one is represented by optical observables linked with geometrical quantities directly detectable from images. In this work, these are: where is the difference in pixels between the Center of Brightness (CoB) and Center of Figure (CoF) of the body projected in the image plane and is the range from the CoM. These quantities can be used to generate a position estimate in the camera frame, which can be transformed in or by simulating on-board attitude determination from a star-tracker alongside the assumption of knowledge of the rigid rotation between the inertially fixed reference frames used by star-tracker and the or frames. The second and third sets are represented by the spacecraft position respectively in spherical and cartesian coordinates in the frame. Similarly, the fourth and fifth sets are represented by the same coordinates in the frame.
By combining the 4 different small bodies with these 5 different labeling strategies, a total of 20 cases are considered in this paper. The notation used to distinguish between these cases is summarized in Tab. 1. Table 1 Notation used for the datasets used in this work.

B. Preprocessing
After the dataset generation step, a total of 70000, 1024 × 1024 grayscale images with ideal pointing towards the CoM of each body are obtained. These images, however, cannot be used directly as input of the IP methods considered in this work for three main reasons.
First, the original image resolution is too high. Due to hardware limitations, image size needs to be reduced. This is typical of the data-driven IP methods used, which could encounter memory or processing saturation issues if working with images at native resolutions. Second, the ideal pointing assumed during rendering simplifies image generation but causes poor variability of the input-label relationship, which can cause poor generalization capability of the IP methods. Third, rendered images are ideal, thus far from realistic camera acquisitions.
All three issues are addressed together in a unique preprocessing pipeline which is a novel contribution of this work.
The pipeline transforms an image and its associated labels from a geometrical and rendering space referred to as S 0 into a new space S 2 that can be efficiently used by data-driven IP algorithms. A sketch of the pipeline is illustrated in Fig. 2 for clarity. The starting point of the pipeline is the image rendered in Blender, which by definition belongs to S 0 space. The image is binarized using the Otsu method [11] and a simple blob analysis is performed to compute the CoB (represented by the yellow point in Fig. 2) and the bounding box Γ (defined by two corner coordinates Γ 1 , Γ 2 , its width Γ 3 , and height Γ 4 ). All these quantities are computed in the frame [12], which is considered in the image plane as centered on the top-left corner of the image.
Cropping is performed around Γ to transform the image to a Γ 3 × Γ 4 snippet. Because the IP methods considered in this work assume fixed-size inputs, the variable-size snippet needs to be transformed into a fixed-size one.
The upper closest possible value between the largest size of Γ 3 × Γ 4 and 128, 256, 512, or 1024 is computed. Once the target size is determined, padding is performed using two parameters, , and , which are randomly generated.
These represent the necessary integer lengths that need to be added to Γ and are chosen respectively as 0 ≤ ≤ − Γ 3 and 0 ≤ ≤ − Γ 4 . After random padding, the image is said to be transformed in S 1 space. Note that the padding is performed by re-using the same pixels from the image defined in S 0 or by introducing zero-pixels values (using the nearest criteria) whenever this is not possible. When passing from S 0 to S 1 also the labels linked to the image need to be transformed, resulting in a change in the and coordinates: while the other labels remain unchanged: Noise is then added to the image, which is resized to a 128 × 128 resolution, transforming it in S 2 space. The noise step is optional for images in the datasets, while in an operational scenario it would be avoided. When passing from S 1 to S 2 the labels are transformed as follows: while 1 and 2 remain unchanged:

. Combining Equation 2 with Equation 3, Equation 4, and
Equation 5 it is possible to transform the image-label pairs for all the datasets in Tab. 1. All the IP methods used in this work use the image-label pairs in S 2 space, which need transformation to S 0 for on-board usage. The methods working on the datasets labeled with cartesian and polar coordinates either in the or reference frames generate directly a position estimate as: where reflects the reference frame used, and Ω is the transformation function from spherical to cartesian coordinates.
On the other hand, the methods working with the ( , ) labels need intermediate steps to generate a position estimate.
In inference, these methods generate the following optical observables: where S 2 and S 2 are output of the IP methods while S 0 and are parameters computed during the image preprocessing algorithm. In this work, the latter quantities are computed offline, however, it is noted that the algorithm can work also online during inference to prepare any incoming image with the proper format for the application of the IP method. The observable vector ,S 0 is transformed from the reference frame which express pixel coordinates on the image, to the Image Plane (ImP) reference frame using the inverse of the camera calibration matrix −1 [12]: The vector is then transformed into a Line of Sight (LoS) vector in the reference frame. Using the attitude quaternion of the spacecraft (which is assumed to be known from attitude determination from a Star-tracker) and assuming to know the rigid rotation between the inertial reference frame used by the star-tracker and a known asteroid frame, this LoS is transformed with the use of S 0 into a position estimate in reference frame: In this work, the frame is used for simplicity and no error is simulated on the attitude quaternion.

C. Image Processing
In this section, the IP methods are described in detail with a standardized notation, illustrated in the architecture outlined in Fig. 3. The input is represented by a 128 × 128 × tensor while the output is a vector or matrix whose elements represent the specific labels associated with each image or tensor. To simplify the discussion, the architecture is conceptually divided into two main portions: the encoder and the head.
In the first one, a hierarchical sequence of convolutions, activation functions, and pooling operations are performed to extract spatial information which is synthesized in the last layer in a latent tensor represented in orange in Fig. 3.
Such 3D tensor is flattened into a 2D vector referred to as the fully connected layer which constitutes the first layer of the head of the architecture. In this portion, several hidden layers may or may not be introduced before the output layer. A network is therefore established to map the connection between the neurons from the fully connected layer, which embeds spatial information extracted from the image in the latent space, and the output layer, which expresses the desired output label of the architecture.
In this work, convolutional layers are driven by weights and biases that only influence the kernels used in the convolution operations. On the other hand, the weights and biases of the head are representing the influence of the neuron connections between the fully connected and output layer. In this work, it is important to consider these two sets of weight and biases separately, since depending on the training strategy, they are handled differently. The former will be referred to as and , the latter to and 0 . The set of weights and biases of the entire architecture defines the set of parameters referred to as = ( , , , 0 ). Another parameter used to describe the architecture is represented by which reflects the set of architectural and training choices often referred to as hyper-parameters.
Within this context, to find an IP method which performs well means to find the best desirable sets of and which optimize predefined metrics. Training an architecture can be reduced to a parameter estimation problem with a double nested optimization loop. During training, both architectural and training choices (global parameters) and inner weights and biases (local parameters) need to be found out through data. When these are found, they are frozen and at inference the architecture can be considered as a function (parameterized by and ) on the input tensor that generate an output vector : Now that a shared nomenclature is established, the characteristics of the four IP methods are described in details.

Convolutional Neural Network
In a CNN the architecture in Fig. 3 is treated as a neural network whose local parameters are found using GD All the GD-based methods are structured in a two-phase process. First, a forward pass is executed by a specific network defined by to generate an output. The comparison between such output and its ground-truth value is quantified by a loss metric, which in turn is fundamental in evaluating the performance of the network as well as to determine the necessary adjustments to for improvement. These are determined from the computation of the gradient, that is estimated in a backward pass in the network followed by an update of to be used in the next iteration. These two passages, especially the last, are computationally expensive, making the training a laborious procedure. A possible way to reduce the training time is to avoid the backward pass, that also require a different optimization scheme.

Convolutional Extreme Learning Machine
A CELM removes the need for a backward pass by solving the optimization problem on the connections between the last hidden layer and the output layer of the architecture.
In this work, the architecture and initialization of CELM are the same as the one of a CNN. However, greater importance is given in the CELM on the random distribution of the eightscnn and of the kernels and of 0 . Once they are randomly initialized when the architecture is generated, they are frozen and are not changed during training. This is not the same for eightscelm, that are the only ones changed during training. Given a set of true input-output samples ( , ), the forward pass of the input into the network generates a hidden layer output matrix right before the output layer: . . .

( )
where is the dimension of the hidden layer before the output one, represented by the number of neurons in the layer.
The training data target matrix is then defined as: To find the best set of weights that matches the matrix , the following optimization problem shall be solved [7]: Minimize : that is a regularized least square that depends on the coefficient. The inclusion of first term regarding the minimization of the weights vector increases the stability and improves generalization capabilities [4]. The minimization problem can be solved in an efficient way as: Since the remaining weights and biases ( , , and 0 ), are randomly fixed at initialization and are never changed, there is no need for a backward pass. Only the forward pass and the LS problem need to be processed, which makes the training of a CELM architecture much faster than the one of a CNN. Note that in this work 0 is not considered in the CELM architecture.

Hybrid Convolutional Extreme Learning Machine
By combining the design of the CNN and the one of the CELM, the capability of a hybrid IP method is also investigated.
In the Hybrid Convolutional Extreme Learning Machine (HCELM), transfer learning is used to get the weights and biases of the kernels from the encoder of a previously trained CNN architecture, while the head is trained using the LS method, as in the CELM paradigm. Two IP methods have been designed with this strategy: one that uses the encoder from the best CNN architecture for that dataset (referred to HCELM), and one that uses the encoder from the best CNN architecture for that reference frame (referred to Hybrid Convolutional Extreme Learning Machine 3 (HCELM 3 )). In

Training
In this section, the strategy to train the IP methods is described. The architecture of each method is generated following a rigid procedural methodology to define the parameters of Θ. A schematic of the procedure is illustrated in Fig. 5.
To begin with, a great number of CELM architectures are initialized following the same rules for each dataset. The encoder of each architecture is designed with a hierarchical structure [6]. While going deeper from input to the fully connected layer, the starting 128 × 128 × 1 tensor is squeezed; its size is halved as function of the depth level as = 2 7− while its depth is doubled starting from an arbitrary value of 16. Each depth level is made up of the consecutive application of a convolution, activation function, and pooling operation.
The convolutions are performed with 3 × 3 kernels while the number of kernels used at each depth is set to increase exponentially as: The number of weights and biases at each depth level of the encoder can thus be determined: Store β(C) on the depth level as = 2 , the number of and 0 of the head is computed as:

CNN -Training and Validation
where 0 is the number of neurons composing the output. For improved handling of the labels, whenever the ( 1 , 2 , ) labels are considered the azimuth angle 1 is transformed in the adimensional ( 1 , 1 ) pair. 0 is thus equal to 3 in all datasets but those for which spherical coordinates are used as labels, for which 0 = 4. An example of a 5 layers architecture is illustrated in Tab. 2 with a parameter count divided per depth. Having defined a procedural set of rules to generate each architecture, these are generated with the hyper-parameters Θ summarized in Tab. 3.
By considering the possible combination between , , , and a total of 120 different CELM architectures are considered. For each kernel distribution, a random initialization is executed 3 times, meaning that a total of 360 CELM networks are generated for each dataset in Tab. 1.
Once the forward pass of the CELM is executed, for each architecture the LS optimization problem is run 9 different times during training of the CELM with different regularization terms . The training set is used to determine all possible values of depending on , while the validation set is used to determine the best value of . The combination of determined from the training set and determined from the validation set is used in inference on the test set to produce the estimated labels. Because the training time of CELM is orders of magnitude faster than the one of the CNN (∼ 1 compared to ∼ 600 , on average), with the use of CELM networks it is possible to explore the architecture space to find those that seems inherently more suitable for the task at hand. As suggested in [6], a great portion of the performance of the CNN seems to be given by the architectural choices, which are often neglected nor sufficiently  For each dataset, the best CELM network is defined as the one achieving the minimum positioning error (defined in Sec. III) on the test set while its parameters are saved.
For each dataset, a CNN architecture is then initialized with the hyper-parameters found from the CELM architectures.
This CNN is then trained by varying and and performing 3 runs for each combination, for a total of 45 cases. For each of them, the CNN is trained for 300 epochs while the best value on the validation loss is used to instantiate the weights and biases of the best possible realization of the CNN.
In this way, we provide the CNN with an architecture that has proven to work properly even with random kernels and investigate whether the re-arrangement of the weights and biases of these kernels together with the ones in the head portion can further improve the baseline performances by extracting additional information or not. The set of for the best CELM and CNN architectures is summarized in Tab frame. This means that in the first case 20 encoders are used, one for each architecture, then 3 encoders are used, each of them being the representative of the CNN architecture working best in that reference frame.
All IP methods are considered to work with normalized input and output, which has been observed to improve the overall performance. The labels are normalized over the maximum and minimum values found in each training set.
While CELM are trained with the entire dataset at once, CNN architectures are trained with batches of images. In the

III. Results
The four IP methods trained over all datasets of Tab where ,S 0 , ,S 0 are respectively the estimated and true position in the reference frame while S 0 is the true range from the body CoM, all evaluated in S 0 . Also, for the architectures that work with the ( , ) labels, the following additional metrics are defined: From Fig. 6, Fig. 7, and Fig. 8 global performance in terms of are summarized for all cases considered. From these plots, it is possible to draw important considerations on the importance of labeling strategy, reference frame, and training method.
In Fig. 6 the box plot are organized from top-down in groups of 5 by macro-categories based on the labeling strategy and then in clusters of 4 based on the IP method, ordered from top-down as CELM, CNN, HCELM, and HCELM 3 . It is possible to immediately appreciate that all the methods trained with the labeling strategy based on ( , ) outperform all others considered. This is especially remarkable for the CELM method, which is performing similarly to the other IP methods when considering this labeling strategy. On the other hand, in all other cases, CELM-based methods perform poorly. This is especially true when compared with the CNN ones, which outperform all methods considered in all possible combinations. Focusing on the importance of the reference frame, it is observed that overall better performances are achieved when the labels are expressed in the reference frame than in the one, while the choice of the coordinate system (cartesian or polar) seems to cause only minor differences.
An even more concise representation is visible in the semi-log plot of the mean error clustered by body in Fig. 7 Fig. 7 it is also possible to observe a trend depending on the shape considered: simpler, regular shapes such as the one of Didymos (D) are better exploited for navigation than highly irregular ones such as Lutetia (L), 67P (P), and Hartley (H).
Finally, Fig. 8 reports in a stacked histogram plot the share for each dataset in which an IP method performs better than all others considered. While the CNN is always the best candidate, HCELM consistently scores as the second best method across all datasets considered. Moreover, the CELM is considered the best third option only in the case of the ( , ) labeling strategy. This graphic representation is possible thanks to the fact that the geometric points considered across the test sets are the same.
Since it is observed that the labeling strategy based on ( , ) works better than the others, the performance is further investigated in this case. In Fig. 9 the histograms of the and metrics in S 2 space are illustrated. When considering D, it is commented that the CELM shows a much larger mean error in than the other methods while only a smaller variability both in terms of means and variance is observed in . The fundamental failure mechanism for which the CELM method performs worse than the CNN one with D is thus represented by a larger mean error in the estimated CoF coordinates. When considering L and P, the difference between the histograms of the various methods is more subtle since only small variations in the mean and variance are observed across the different methods. On the other hand, when considering H, it is possible to see that CNN is much more accurate both in the CoF and range estimates than all other methods. The hybrid CELM performs better in the range estimate (with similar variance than the CNN) but does not perform at the same level as the CNN in the CoF estimate.
In Fig. 10 the error ellipses of the CoF coordinate in S 2 space are illustrated together with the error ellipse (dashed dark ellipse) that would have been obtained by not correcting the CoB with a data-driven scattering law implemented by the IP. It is also observed that different than all other ellipses, the CELM ellipse with D is not centered in the zero error point thus introducing a bias that is causing the larger mean error already observed in Fig. 9.
In Fig. 11 a representative case of the position error of the CELM in reference frame for the L body is illustrated.
It is possible to see that the position error is one order of magnitude higher in the boresight direction of the camera than on the other axis. This result is expected from optical-based navigation systems. This means that the error in the range estimate from the body is the major contributor to the positioning performances when using this labeling strategy.
Eventually, a visualization of a sample of 150 images in S 2 and the CoF prediction by the different IP methods is illustrated in Fig. 12.

IV. Conclusions
In this work, CELMs are investigated as a possible alternative to CNNs for autonomous vision-based navigation systems around small bodies. This is done in an extensive analysis considering 4 different small body shapes and 5 different labeling strategies, for a total of 20 scenarios, each of which is examined with 4 different IP methods, resulting in tens of thousands of different architectures explored.
It is demonstrated that the coupling between IP method and labeling strategy plays a fundamental role to achieve the desired performance. In particular, the ( , ) labeling strategy is found to be the best one. In this case, the navigation labels are given by two geometrical quantities which are estimated directly from the image, and a third, more difficult, label associated with the range from the body. The error on such a label is the one driving the error in the position estimate.
It is also noted that amongst the other labeling strategies, the ones adopting the reference frame perform better than those using the one. This hints at the fact that in the former the estimated position relies on the correlation with global geometric properties of the body which are not dependent on the specific rotational state but are rather dependent on the illumination conditions. On the contrary, in the latter, the methods need to learn illumination invariant features that correlate to the specific rotational state of the body. This seems to be much harder to do, at least when 7500 images are considered for training. Future work will be focused on changes in the method to allow larger datasets to be used. Irrespective of the reference frame, it is commented that the mapping from images to polar or cartesian coordinates seems to degrade the performance. In particular, while the ( , ) and ( 1 , 2 , ) labeling strategies both share as a label, the latter performs considerably worse than the former. This may hint at a more difficult mapping from images into two angles in spherical coordinates than the mapping to projective quantities in the image plane, the latter correlating much easily and directly from image.
Although the CNN outperform all other methods considered, the CELM seems a credible alternative, since it performs in the same order of magnitude as the CNN method when considering the best labeling strategy.
It is also commented that the usage of CELM for exploration of the architecture space generated CNN architectures that experienced robust training. This bootstrap training strategy could be deployed as a standard and efficient way to explore the global architecture space as suggested in [6]. Moreover, consistently with the findings in [4][5][6][7], a preference is observed for orthogonalization of the weights and biases of the kernels, as it is possible to see from the best hyper-parameters in Tab. 4.
Future works will address the addition of dropout into the networks to generate uncertainty estimates, algorithmic design choices to adapt the CELM for small body images, larger datasets, and the usage of segmentation maps as input in place of grayscale images, the former representing easier input to work with.