domain: info math

text:
Michele De Vita
email: michele.devita@stud.unifi.it

Carlotta Giannelli

Sofia Imperatore
email: sofia.imperatore@unifi.it

Sofia Imperatore 1⋆⋆

Angelos Mantzaflaris
email: angelos.mantzaflaris@inria.fr

Parameterization learning with convolutional neural networks for gridded data fitting

Keywords: deep learning, convolutional neural networks, parameterization, data fitting, THB-spline, tensor data

come

Parameterization learning with convolutional neural networks for gridded data fitting 1 Introduction

Gridded data consists of (measurement) values obtained at regularly spaced points that form a rectilinear grid. This kind of geometric data with a regular structure arises in different scientific application settings, such as geographic information systems, meteorology, and medical imaging, when, for example, measurements at regular space or time intervals are acquired through different types of technological instruments. The result is a set of ordered points that ⋆⋆ Corresponding author can represent curvilinear, surface, volumetric, or even higher dimension data. Since the possibility of suitably extending the measurement information to regions where no data is available is often required, efficient parameterization and fitting schemes are usually needed.

A data fitting scheme addresses the construction of curves and surfaces which properly approximate an arbitrary collection of geometric data, such as points or derivative vectors. This type of approximation problem is very common in different applications, where the generated points typically contain measurement errors or computational noise. The re-constructed geometric model should captures the shape of the data, while simultaneously avoiding artifacts and unwanted oscillations, which may easily arise when an (almost) exact match of all the data points is constructed. In the polynomial spline setting, a variety of different approximation schemes can be considered. For example, by focusing on spline approximations of a given structured point cloud, a non-linear optimization problem to minimize the fitting error could be performed with respect to the spline coefficients (e.g., the B-spline control points), the parameters, and the knot values. To avoid non-linear problems, the parameters and the knots are usually computed in advance. These two pre-processing steps are usually indicated as parameterization and knot placement problems, respectively.

In this paper we employ Convolutional Neural Networks (CNNs) for the parameterization learning problem to assign parameter values at an arbitrary number of points on a rectilinear grid for their subsequent use in multivariate polynomial spline approximation schemes. Convolutional neural networks are widely exploited in image processing, see e.g., [START_REF] Krizhevsky | ImageNet classification with deep convolutional neural networks[END_REF], but they also thrive when any kind of data involving spatial-correlated patterns, as for example audio files [START_REF] Stoller | Wave-u-net: A multi-scale neural network for end-to-end audio source separation[END_REF] and geometric data processing [START_REF] Wang | Dynamic graph CNN for learning on point clouds[END_REF], have to be processed. We train our model using synthetic data generated by sampling random tensor-product polynomial curves and surfaces by considering univariate and bivariate polynomial least squares approximations. In our experiments we also add noise to the test data to simulate measurement errors in the input of the proposed networks. Moreover, we use the resulting parameterization to fit point clouds by non-uniform spline constructions and show that the accuracy of the model properly scales under adaptive mesh refinement. The choice of CNNs allows us to obtain a method that is agnostic to the size of the input point cloud, and performs well with an input of a different size from the one considered in the training phase. The experiments highlight that our CNN model can also be effective when multivariate adaptive spline fitting with truncated hierarchical B-splines (THB-splines) is considered, see e.g., [START_REF] Giannelli | THB-splines: The truncated basis for hierarchical splines[END_REF][START_REF] Giannelli | Špeh THB-splines: An effective mathematical technology for adaptive refinement in geometric design and isogeometric analysis[END_REF][START_REF] Bracco | THBspline approximations for turbine blade design with local B-spline approximations[END_REF].

The structure of the paper is as follows. Section 2 presents a brief overview of related works on the parameterization of gridded data. The polynomial spline approximation problem for data belonging to R N is briefly summarized in Section 3 together with the parameterization problem in R D . The proposed convolutional parameterization learning model is presented in Section 4 in its general formulation. Details about the hyperparameters selection and the training procedure are provided in Section 5, with special focus on the curve and surface cases. A selection of numerical experiments is reported in Section 6. Finally, Section 7 contains few concluding remarks.

Related works

Over the last decades several approaches have been developed to handle the parameterization of gridded data. Standard methods rely either on a uniform distribution of parameter values over the parametric domain or on scaled configurations of the physical points on the parametric space, as for example the chord-length or centripetal parameterizations [START_REF] Lee | Choosing nodes in parametric curve interpolation[END_REF]. Alternatives to these standard methods were presented in the literature, with special focus on the univariate setting. In particular, variants of the scaled parameterizations were proposed in [START_REF] Shamsuddin | A hybrid parameterization method for NURBS[END_REF][START_REF] Fang | An improved parameterization method for B-spline curve and surface interpolation[END_REF], whereas an example of an iterative procedure of different kind can be found in [START_REF] Balta | Dynamic centripetal parameterization method for B-spline curve interpolation[END_REF]. In addition, the so-called universal parameterization method, closely related to the uniform parameterization, was presented in [START_REF] Lim | A universal parametrization in B-spline curve and surface interpolation[END_REF] in the context of univariate and bivariate spline interpolation.

The potential of exploiting machine learning techniques for B-spline curve approximation was originally outlined in [START_REF] Laube | Learnt knot placement in B-spline curve approximation using support vector machines[END_REF] and [START_REF] Laube | Deep learning parametrization for B-spline curve approximation[END_REF] by considering support vector machines and multilayer perceptrons, respectively. While the focus of [START_REF] Laube | Learnt knot placement in B-spline curve approximation using support vector machines[END_REF] was on identifying a suitable knot placement strategy, a scheme based on two interdependent deep neural networks to compute both the knot and parameter values was proposed in [START_REF] Laube | Deep learning parametrization for B-spline curve approximation[END_REF]. In particular, the neural network used to predict the parameters consists of a standard multilayer perceptron which takes in input a planar point sequence of fixed length equal to 100 and returns as output the corresponding parametric values. Therefore, if the given number of data points does not match 100, pre-and post-processing steps for super/sub-sampling the input data need to be performed. To overcome the computational limitations of this method, the authors in [START_REF] Scholz | Parameterization for polynomial curve approximation via residual deep neural networks[END_REF] proposed a method for univariate polynomial curve approximation based on a residual neural network. In this case the network takes in input a planar point sequence of length 2q + 1, where q ∈ N is the polynomial degree. The residual building blocks of the proposed networks are characterized by fully connected layers, which constrain the architecture to deal with an input of fixed size. Consequently, also this approach requires multiple network evaluations and a post-processing step to handle point sequences of arbitrary size. It should also be noted that, to the best of our knowledge, datadriven methods to address the parameterization problem of gridded data sets in a multivariate setting were not previously proposed. To properly process input datasets with varying dimension without requiring additional computations, we here propose a deep convolutional approach.

3 The parameterization problem for data fitting schemes Let P = {p I ∈ R N : I ∈ Γ } be a gridded point cloud, with a rectilinear grid structure along D ≤ N directions, defined by the set of multi-indices Γ , so that each item of the point cloud belongs to R N . Let d = 1, . . . , D be the index associated to the parametric direction and let m d be the point sequence length along each direction d. We define the multi-index set

Γ = {I = (i 1 , . . . , i D ) : i d = 1, . . . , m d , d = 1, . . . , D} , (1)
where each multi-index I uniquely identifies a point p I of the point cloud P . Consequently, the point cloud can be rewritten as

P = p I ∈ R N : I = (i 1 , . . . , i D ) , i d = 1, . . . , m d , d = 1, . . . , D . (2)
We denote the cardinality of the point cloud as m := D d=1 m d . Fig. 1 shows an example of a point cloud P , characterized by D = 2 parametric directions and elements belonging to R N , with N = 3. The point sequences in the two directions have lengths m 1 = m 2 = 4. Each item of the point cloud p I is uniquely determined by the multi-index I = (i 1 , i 2 ), with i 1 and i 2 varying between 1 and 4.

The least squares spline fitting method addresses the construction of a spline model on a certain parametric space Ω ⊆ R D , which minimizes the residual approximation error. More precisely, by choosing the polynomial multi-degree q ∈ N D and a domain tessellation T (Ω), the objective of the method consists in finding a spline s :

Ω → R N , s(x) = n j=1 c j ϕ j (x), x ∈ Ω, (3)
where ϕ j : Ω → R for i = 1, . . . , n are D-variate splines, which minimizes the residual min c1,...,cn

I∈Γ ∥p I -s(x I )∥ 2 2 = min c ∥P -Bc∥ 2 2 (4)
for certain choice of the parameters x I ∈ Ω, for I ∈ Γ , where P are the samples of the given point cloud, B := ϕ j (x I ), for j = 1, . . . , n and I ∈ Γ , is the system collocation matrix, c := c 1 , . . . c n T is the matrix of unknowns containing the coefficients of the reconstructed geometry. Note that in equation ( 4), p I ∈ R N , for I ∈ Γ , are the Cartesian coordinates of the physical points in R N , whereas x I ∈ Ω for I ∈ Γ are their associated parametric values.

To construct the gridded data approximation, a parameterization of the point cloud P is needed. This problem requires to identify a set of suitable parametric values

X = {x I ∈ Ω : I = (i 1 , . . . , i D ) , i d = 1, . . . , m d , d = 1, . . . , D} , (5)
so that any point p I ∈ P is associated to the parameter x I ∈ X. Without loss of generality, we set Ω = [0, 1] D . In the tensor-product polynomial spline setting, once the univariate parameterization along each parametric direction is computed, the multivariate parameterization can be simply obtained by averaging across all the other parametric directions [START_REF] Piegl | The NURBS Book[END_REF]. The identification of a suitable parametrization method for the definition of the set X of parameter values is a challenging open problem within the polynomial spline fitting framework. In particular, it can naturally influence the quality of the final approximation leading to sub-optimal approximation error results if not carefully handled, see e.g., [START_REF] Floater | On the deviation of a parametric cubic spline interpolant from its data polygon[END_REF][START_REF] Floater | Parameterization for curve interpolation[END_REF] which address the case of curve interpolation. Among other optimization approaches, in the next sections, we design a parameterization scheme based on convolutional neural network to address this point.

The convolutional parameterization model

The data connectivity information on proximity and distance plays a fundamental role in geometric data processing. Within the data learning framework, modeling the interactions between data relays on their weighted sums, computed via neural networks. Convolutional operators have become a building block for deep learning architectures to process data with a grid-like topology, for example when reconstructing a free-form model starting from a set of gridded observations p I ∈ R N , for I ∈ Γ . Architectures which involve convolutional operators within their layers are addressed as convolutional neural networks (CNNs). Thanks to the nature of convolutional operators, these networks are characterized by three unique properties: sparse connectivity, parameter sharing, and shift equivariant representations [START_REF] Goodfellow | Deep Learning[END_REF]. Sparse connectivity means that the networks have local receptive fields, which collect information jointly from spatially neighbouring inputs. Parameter sharing is achieved by requiring the receptive field to assume the same values on different instances of the input, and, consequently, the same parameters are involved to compute each output unit. Finally, shift equivariant representations follow by the definition of the convolutional operations and from the previous two properties. Fig. 2 illustrates the locality and parameter sharing properties for the one (D = 1) and the two (D = 2) dimensional convolutional operator. In both cases, the input is processed by a filter through a matrix convolution operator to generate the output feature map, usually given in input to the successive convolutional operator. In particular, each output item is given by the sum of the element wise product between different local values of the input (green numbers in Fig. 2) and the corresponding values of the filter (red numbers in Fig. 2).

The architecture we designed is a pure CNN, hence consisting only of convolutional blocks. This choice has been made to take advantage of the locality of convolutional operators in order to be able to support variable input sizes without any additional effort and/or pre-or post-processing of the data. In particular, since the convolutional neural network is characterized by local operators, the particular size of the considered point clouds does not play a fundamental role, whereas the filter dimension does. Consequently, even if a convolutional model is trained on point cloud of fixed size, it can be easily generalized to point cloud of different dimensions. In the following we detail the proposed architecture, while

Input Features D-conv ReLU ×B D-conv Softmax CumSum Output
Fig. 3: The CNN architecture.

additional information on the data generation, the hyperparameters selection, and the training process are provided in the next section.

Learning model architecture

The proposed architecture, presented in Fig. 3, is a sequential CNN which consists of a feature extraction layer, B convolutional building blocks (characterized by convolutional layers coupled with ReLU activation functions), a final block of one convolutional layer, the softmax activation function, and a cumulative sum layer connected to the output. We now give a detailed description of each component of our architecture to illustrate how to obtain the output parameterization result.

Input The learning model takes in input a point cloud P of gridded data, organized in a tensor structure, as in [START_REF] Bracco | THBspline approximations for turbine blade design with local B-spline approximations[END_REF].

L d out =     L d in + 2 • pad d -dil d • ks d -1 s d + 1     . (7)
In ReLU Activation functions are non linear transformations which have been introduced in order to obtain nonlinear learning models [START_REF] Goodfellow | Deep Learning[END_REF]. The rectified linear unit (ReLU) is a function which is applied element-wise to the output of all except the last convolutional hidden layer. It is defined as ρ(t) := max {0, t} and it is used for practical applications among most feedforward neural networks due both to its simple piecewise linear structure and its high performance [START_REF] Yarotsky | Error bounds for approximations with deep ReLU networks[END_REF]. Softmax As a final activation function we apply the softmax function along each direction d = 1, . . . , D to generate parameter intervals that form a partition of unity along each direction d. In particular, the output of the softmax are the elements

∆ d I ∈ [0,

Loss function

Once the parameterization as in ( 5) is computed with the convolutional parameterization model, we estimate the control points of the approximating geometry by solving the polynomial least squares minimization problem in Bernstein-Bézier form. The reconstructed model s : [0, 1] D → R N is then used to defined the loss function of our unsupervised learning model.

In particular, let L : R m×N × [0, 1] m×D → R be the mean squared error (MSE) between the data P and the corresponding values of the polynomial approximation model P = p I ∈ R N : I ∈ Γ , so that s(x I ) = p I for each multi-index I. Therefore, the loss function is the residual

L (P , X) := 1 m I∈Γ ∥p I -p I ∥ 2 2
of the least squared fitting problem corresponding to the parameters X in output from the network. Note that our method falls into the unsupervised learning class, since the real parametric values are not considered to guide the optimization problem and the label role is played by the input points.

Curve and surfaces parameterization learning

In this section, we provide the information concerning the training of our models. The method we propose can be implemented for any spatial dimension N ≥ 1 and parametric dimension D ≤ N , thanks to the general definition of the convolutional operator. Note that we just set up and train a different model for any distinct pair of values N and D, independently of the considered point cloud.

In particular, we focus on two cases: the parameterization of point sequences, when N = 2 or N = 3 and D = 1, and the parameterization of spatial gridded data, when N = 3 and D = 2, to subsequently compute polynomial and spline approximations in terms of univariate curves and bivariate surfaces, respectively.

Data generation

Deep learning has its roots in data-driven artificial intelligence and the considered datasets naturally play a fundamental role for the performance of the model. In view of the lack of public datasets for the problem at hand, the points used in the training phase have been randomly generated.

We generated multiple synthetic data by sampling different kinds of parametric objects, either curves (D = 1, N = 2, 3) or surfaces (D = 2, N = 3). The proposed algorithm for data generation is the following.

1. Generate a random parametric object s :

[0, 1] D → R N . 2. Sample u d = {u i d } m d i d =1
parameters according to the random uniform distribution Uniform(0, 1), for d = 1, . . . , D. N .

Normalize the parametric values u

i d in [0, 1], so that u 1 = 0, u m d = 1, i.e.
Note that the geometry s : [0, 1] D → R N can belong to any dimension and since the input processed by the convolutional layers are the edges computed in (6), our parameterization method is affine invariant.

The type of the parametric object in step 1 naturally characterizes the generated point cloud. In particular, for N = 2, 3 and D = 1 we considered random polynomial curves of fixed degree q in Bernstein-Bézier form, where ϕ j for j = 1, . . . , q in (3) are defined in terms of Bernstein polynomials and c j ∈ R N for j = 1, . . . , q are the corresponding control points sampled from the standard normal distribution N (0, 1). For N = 3 and D = 2, we generate an initial polynomial tensor-product surface with control points c j = [c 1 j , c 2 j , c 3 j ] T , where [c 1 j , c 2 j ] T are chosen as the tensor product Greville abscissae for degree q and c 3 j are sampled randomly accordingly to the uniform distribution in [0, 1]. The resulting surface is then also randomly rotated.

Hyperparameters selection and training process

Training our neural network consists in using our synthetic data to minimize the loss function, so that the model weights tend to be optimal for both seen and unseen input point clouds. In order to optimize the performance of the network, we set up an extensive grid search on optimizers, learning rates, filter sizes, and number of convolutional building block when applying the learning parameterization model in the context of planar polynomial curves, namely for N = 2 and D = 1. In particular, we tested three optimizers (Adam, SGD, RMSprop) with learning rates γ = 1e -i, i = 1, . . . , 6. The filter size is set to be twice the polynomial degree plus one, as in [START_REF] Scholz | Parameterization for polynomial curve approximation via residual deep neural networks[END_REF], and the number of building blocks is B = 4. In addition, we randomly select the point sequence length of each training batch in range m ∈ [2d+1, 100]. We obtained the best results using RMSprop as optimizer, learning rate equals to 1e -5, and momentum 0.9. We then used these hyperparameter configurations also to train the model designed for higher dimensions, in particular for N = 3 and D = 2 to address the surface fitting problem. In this case we increase B to 5 and we fix the point cloud size to m = 100. Note that we train a different model for any value of N and D without a fixed number of steps. By following the so-called early stopping criteria, the training phase continues until the validation loss stagnates, a strategy commonly used to avoid overfitting.

Numerical results

In this section we present a selection of experiments to analyze the performance of the introduced learning parameterization model and its generalization capabilities with respect to a variety of structured data and different spline approximation schemes. The quality of the parameterization model is defined in terms of the final geometric approximant accuracy, reported in terms of mean squared error (MSE).

We show the effectiveness of our method starting from simple planar polynomial curve approximations (for N = 2 and D = 1) to arrive at complex adaptive spline surface approximation with THB-splines (for N = 3 and D = 2). In particular, Section 6.1 presents a comparison with state-of-the-art neural networks for polynomial curve approximation in the case of polynomial data. Trigonometric and noisy data approximation is addressed in Section 6.2, in the context of adaptive spline curve approximation. The performance of the proposed CNN model on benchmark datasets for spline curve approximation is then illustrated in Section 6.3. By moving to the surface case, Section 6.4 shows the results on polynomial data, together with the model robustness to noise. Subsequently, in Section 6.5, we provide numerically evidence of the generalization capabilities of the proposed parameterization model with respect to datasets significantly different (in nature and size) from the synthetic data used during the training phase. Finally, in Section 6.6 and Section 6.7 we investigate the performance of the learning parameterization model for tensor-product B-spline and THB-spline fitting schemes, respectively.

Polynomial curve approximation

In this experiment we exploit the convolutional parameterization model on univariate polynomial data for different degrees q and variable point sequence lengths m. For each degree q = 2, . . . , 5, and length m = 20, 50, 80, we generate 100 polynomial point sequences P = p i ∈ R 2 : i = 1, . . . , m as detailed in Section 5.1. Subsequently, we compute X = {x i ∈ [0, 1] : i = 1, . . . , m} with the proposed CNN approach and with standard parameterization schemes of the form

x i+1 = x i + ∆p i ∆p for i = 1, . . . , m -1,
with

x 1 = 0, ∆p i := ||p i+1 -p i || α , ∆p = m-1 j=1 ∆p i .
In particualr, the common uniform, chord-length, and centripetal parameterizations can be recovered by setting α = 0, 1, 1/2 respectively. We then compare the averaged MSE error obtained building the Bézier fitting model at the parametric values given by the best standard technique (STD), the results presented in [16, methods and the RES model in terms of mean squared error. In particular, the MSE obtained with CNN is reduced up to one and two orders of magnitude as regards RES and STD, respectively. The advantage of our convolutional parameterization learning model lies not only in the better performance shown in Table 1, but also in the novelty of its self-contained nature. For each input point sequence of any length, our convolutional learning method directly computes the corresponding point parameterization, without any additional computation. The approach proposed in [START_REF] Scholz | Parameterization for polynomial curve approximation via residual deep neural networks[END_REF] instead requires requires multiple network evaluations and a post processing step to handle arbitrary sequence lengths.

Generalization to spline curve approximation for datasets of different nature

In this example we investigate the generalization capabilities of our convolutional parameterization model to different datasets and adaptive B-spline curve fitting of various degrees (q = 2, 3, 4, 5), as well as its robustness to noise. By considering the algorithms presented in Section 5.1, we generate two distinct datasets of size 100 by sampling trigonometric and polynomial curves. The trigonometric curves are defined as c(t) = a 0 + r j=1 a j cos(jt) + i r j=1 b j sin(jt), for a fixed degree r, where a 0 ∈ C and a j , b j ∈ R are values sampled from the standard normal distribution N (0, 1), for j = 1, . . . , r. In the polynomial case, we added a factor of 1e-2 random Gaussian noise to the data points.

We test the CNN parameterization by considering adaptive B-spline least squares approximation obtained by performing dyadic refinement of the knot intervals which exhibit the highest error values. The resulting MSE is shown in in Table 2 and Table 3 for trigonometric data of degree r = 5 and noisy polynomial data, respectively, for different input sequence lengths (m = 20, 50, 80) and adaptive B-spline approximations with a final number of 5 interior knots. The CNN parameterization model always obtains better results than the best standard methods (STD). In particular, in the case of trigonometric data (Table 2), it gains up to two order accuracy. For the results on noisy polynomial data (Ta-Table 2: MSE error for adaptive B-spline approximation with 5 interior knots of degree q = 2, 3, 4, 5 on trigonometric data for different input lengths m = 20, 50, 80 in Section 6.2.

m = 20 m = 50 m = 80 STD CNN STD CNN STD CNN q = 2 3.01e-4
1.56e-5 2.17e-4 6.08e-6 1.71e-4 5.82e-6 q = 3 1.35e-4 1.63e-5 1.10e-4 1.05e-5 3.16e-5 7.84e-6 q = 4 1.74e-4 5.53e-6 2.48e-5 3.66e-6 3.14e-5 4.74e-6 q = 5 3.09e-5 2.54e-6 1.17e-5 1.88e-6 4.03e-5 1.92e-6

Table 3: MSE error for adaptive B-spline approximation with 5 interior knots of degree q = 2, 3, 4, 5 on polynomial data with random gaussian noise for different input lengths m = 20, 50, 80 in Section 6.2.

m = 20 m = 50 m = 80 STD CNN STD CNN STD CNN q = 2 1.27e-4 7.69e-5 1.32e-4 9.46e-5 1.54e-4 1.02e-4 q = 3
9.54e-5 3.36e-5 1.21e-4 6.33e-5 1.12e-4 7.24e-5 q = 4 7.37e-5 2.68e-5 1.66e-4 8.89e-5 8.41e-5 4.04e-5 q = 5 9.02e-5 2.72e-5 1.68e-4 7.64e-5 1.04e-4 6.11e-5

ble 3), the CNN model accuracy is improved up to one order of magnitude when compared to STD.

To better understand the generalization capabilities of the network when an increasing number of knots is considered, we present a final test with fixed input length equals to m = 20. We then perform a dyadic adaptive B-spline curve fitting with a varying number of interior knots, starting from the polynomial case of 0 internal knots up to 10 internal knots. The results in terms of MSE for the centripetal (CEN), chordlenght (CHL), uniform (UNI) and CNN parameterizations for q = 3, 4, 5 are shown in Fig. 4 (top) for trigonometric data of degree r = 5 and in Fig. 4 (bottom) for noisy polynomial data. We may observe that in Fig. 4 (top) for q = 3, 5 the MSE error gap between CNN and STD is maintained or even increased with respect to the number of interior knots inserted. In Fig. 4 (top) for q = 4 and in Fig. 4 (bottom) for all the degrees, the MSE error gap between CNN and STD tends to reduce with the number of interior knots inserted, but the CNN parameterization model always achieves the best results.

Spline curve approximation of benchmark datasets

In this experiment, we test our parameterization model on point sequences sampled by the benchmark presented in [START_REF] Ohrhallinger | Muthuganapathy. 2D points curve reconstruction survey and benchmarks[END_REF]. By considering B-spline approximation of degree 4 with at most 15 adaptive interior knots, Fig. 5 presents the results Fig. 4: MSE error for adaptive B-spline approximations with increasing number of knots tested on m = 20 trigonometric (top) and noisy (bottom) data for degree q = 3 (left), q = 4 (center) and q = 5 (right) in Section 6.2. obtained employing the standard (UNI, CHL, CEN) and the CNN parameterizations on an increasing number of input points, namely m = 20, 50, 80. An analogous comparison on a sequence of m = 65 points, approximated with degree q = 5 B-spline and adaptive knot vectors with a varying number (4, 7, and 10) of interior knots is shown in Fig. 6. Table 4: MSE error for polynomial least squares approximation of polynomial (left) and noisy (right) data of bidegree (q, q), for q = 2, 3, 4, 5 in Section 6. [START_REF] Floater | On the deviation of a parametric cubic spline interpolant from its data polygon[END_REF]

Polynomial surface approximation

In this experiment we illustrate the performance of the proposed parameterization model on both exact and noisy data sampled from polynomial patches in the bivariate case. In particular, we collect structured point clouds by sampling Bézier surfaces of bi-degree (q, q). More precisely, for each q = 2, 3, 4, 5, we generate 100 point clouds consisting of m = 144 gridded points following the algorithm described in Section 5.1. In addition we also create a noisy version of this datasets by perturbing each point cloud with Gaussian noise of factor ϵ = 5e-3. The results of the parametric polynomial fitting of bidegree (q, q) of the aforementioned generated point clouds are reported in Table 4 in terms of MSE for the best standard parameterization technique (STD) and the proposed convolutional model (CNN). The error values give numerically evidence of the generalization capabilities of the trained model with respect to datasets which are somehow similar but at the same time independent from the ones used in the training phase. In particular, the proposed model avoids potential over-fitting problems, it is robust to noise, and it always increases the accuracy of the polynomial approximation more than 50% over STD parameterization techniques in all the considered configurations.

Generalization to point clouds sampled from geometric models

In this experiment, we analyse the generalization capabilities of the parameterization learning model for data which are different both in dimension and nature from the one used in the training phase. More precisely, we sample point clouds of dimension either m = 1089 or m = 3025 from 3 different B-spline geometries of bidegree (q, q), representing a car, a ship hull, and a wind turbine. The corresponding point clouds are illustrated in Fig. 7 for the case of 1089 (top) and 3025 (bottom) samples. Moreover, for this experiment we have also stored the true parametric values using during the sampling phase and we will refer to them as the ground truth (GT) in the comparisons. Once the data have been parametrized with the convolutional method (CNN), we performed polynomial approximation of bi-degree (q, q) and compare the MSE and MAX obtained with GT and CNN parameterizations. Note that the GT parameters come from a specific B-spline geometry, therefore if used for different spline models they not necessarily lead to zero error at the evaluation sites. The results are reported in Table 5. They numerically show the generalization capabilities of the proposed model with respect to the size of the input point clouds and to the nature of the data. For the car geometry, in all the considered configurations, the MSE error of the final reconstructed polynomial geometry is smaller in the case of the CNN parameterization, with a gain in accuracy up to more than 50% of for q = 3. When the ship hull geometry is considered, the best performances are achieved if the CNN parameterization is applied in combination with lower polynomial bidegrees, namely q = 2, 3, where for q = 2 the final reconstructed geometry of the bigger point cloud gains more than 60% of accuracy. As concerns the final wind turbine geometry, the best advantage is registered on the point cloud of size m = 1089 among all the different bidegrees, when the MSE is improved up to more than 60% for bidegree [START_REF] Fang | An improved parameterization method for B-spline curve and surface interpolation[END_REF][START_REF] Fang | An improved parameterization method for B-spline curve and surface interpolation[END_REF]. The same bidegree [START_REF] Fang | An improved parameterization method for B-spline curve and surface interpolation[END_REF][START_REF] Fang | An improved parameterization method for B-spline curve and surface interpolation[END_REF] leads also to the best error results for the bigger wind turbine point cloud.

Generalization to tensor product B-splines approximation

In this example we show the generalization capabilities of the convolutional parameterization for tensor product B-spline models. In particular, we consider the point clouds obtained from the car shell geometry introduced in Section 6.5 for which we perform a tensor product B-spline approximation for bidegree (3, 3) and 3 levels of uniform refinement, both considering the given parameterization (GT) and the convolutional parameterization (CNN).

Table 5: MSE for polynomial least squares surface approximation with different bidegrees (q, q) of the car shell, the ship hull, and the wind turbine point clouds of size m = 1089 and m = 3025 in Section 6.5. Concerning the point cloud of dimension 1089, the MSE error registered using GT parameterization is 8.18e-6. If the CNN parameterization is considered, the MSE is more than halved and reduces to 3.54e-6. Analogous results are achieved when considering the point clouds of dimension 3025 for which the GT parameterization leads to MSE of 8.21e-6, whereas by using the CNN parameterization the accuracy improves more than 40% obtaining an MSE equals to 4.31e-6.

The scaled error distributions of the final approximation plotted on the parametric and physical domain, for both GT and CNN parameterizations, are illustrated in Fig. 8 for the point cloud of size 3025. More precisely, plots (a) and (b) in Fig. 8 are the error distributions for the GT parameterization, whereas the corresponding results for the CNN parameterization are shown in plots (c) and (d) of the same figure.

Hierarchical spline approximation of a ship-hull

In this experiment we exploit the performance of the convolutional parameterization model when hierarchical B-spline fitting with THB-splines is considered. In this case the final approximation is obtained with the adaptive THB-spline fitting scheme presented in [START_REF] Kiss | Adaptive CAD model (re-)construction with THB-splines[END_REF], enriched by a step of parameter correction at each adaptive iteration. In particular, given a point cloud and a suitable parameterization, the considered adaptive fitting scheme designs a hierarchical spline model of the form (3), where the basis function are bivariate truncated hierarchical B-splines (THB-splines). The final THB-spline geometric model is obtained by the iterative alternation of two steps: the computation of the control points and an adaptive refinement of the spline space. At each iteration of the adaptive

Closure

We presented a general parameterization learning method for polynomial and spline approximation of gridded data point by exploiting CNN architectures. In particular, the learning model takes in input a structured point cloud and gives in output suitable parameters to be used for the geometric modeling of the input data via spline constructions. Our network is characterized by a feature extraction layer, which yields to translation invariance of the proposed method. This layer is then followed by a sequence of CNN layers which encode and propagate proximity information to successive layers. Finally, the parametric values for the input point cloud can be recovered by applying the softmax activation function and a cumulative sum layer. For the numerical experiments, we applied the proposed method by considering its univariate and bivariate configuration. The implemented schemes compares favorably with respect to standard techniques and alternative deep parameterization methods. The results also confirm the flexibility of the convolutional approach, naturally able to approximate gridded point clouds of variable size. In particular, we highlighted the generalization capabilities of our model for spline approximation on unseen datasets of different kind, including more general (noisy) data configurations and available datasets. We numerically demonstrate that the accuracy of the spline approximation based on the CNN parameterization properly scales with respect to (non-)uniform spline configurations, while still comparing favorably with respect to standard choices. Finally, we also show that the proposed method can be suitably combined with THB-spline approximation for geometric modeling applications. The potential of exploiting the convolutional spline learning model both for optimal knot placement and mesh refinement, as well as for control points estimation, are interesting directions for future research.

Fig. 1 :

 1 Fig. 1: Example of a gridded point cloud with items in R N , with N = 3, and D = 2 parametric directions. The feature extraction step derives the edge vectors e ij of the input grid.

Fig. 2 :

 2 Fig. 2: Examples of convolutional operators for D = 1 (a) and D = 2 (b).

  Features The feature extraction layer produces a relational representation of the point cloud that is translation invariant. It consists in the computation of m d -1 edges in R N between each pair of consecutive points along each direction d = 1, . . . , D. Fig. 1 illustrates the feature extraction process for a structured point cloud with N = 3 and D = 2. In particular, the edge associated to the element p 22 along the direction d = 1 is e 1 22 := p 32 -p 22 (swown in red in Fig. 1), whereas the edge associated to the element p 32 along the direction d = 2 is e 2 32 := p 33 -p 32 (shown in blue in Fig. 1). More precisely, for any D ≥ 1, let I = (i 1 , . . . , i D ) for i d = 1, . . . , m d -1, d = 1, . . . , D, and let I + d be the multi-index I augmented by 1 along the d-th direction, for d = 1, . . . , D, namely I + d := (i 1 , . . . , i d + 1, . . . , i D ). For all I = (i 1 , . . . , i D ) with i d = 1, . . . , m d -1 and d = 1, . . . , D, we define the edges e d I as e d I := p I+d -p I . Note that e d I ∈ R N . By concatenating the edges e I = e 1 I , . . . , e D I , the extracted features are defined as the collection of edges E = e I ∈ R N ×D : I = (i 1 , . . . , i D ), i d = 1, . . . , m d -1, d = 1, . . . , D . (6) The extracted features are then passed as input to the convolutional building blocks. D-conv The convolutional layers perform a D-convolution over the input. Each block is characterized by discrete convolutional operators whose hyperparameters along each direction d = 1, . . . , D are the number of input and output channels c d in , c d out ∈ N \ {0}, the size of the convolving filter ks d ∈ 2N + 1, the amount of padding added to the boundaries of the input pad d ∈ N, the stride of the convolutional operator s d ∈ N \ {0}, and the spacing between the filter elements dil d ∈ N. By denoting with u d ∈ R c d in ×L d in the input instance of any convolutional layer along the direction d = 1, . . . , D and with y d ∈ R c d out ×L d out its output, the relationship between the input (L d in ) and the output (L d out ) sizes is

  order to fully convert the input along any direction d = 1, . . . , D through the filter of dimension ks d and obtain an output with the same sequence length, namely L d in = L d out , we specify s d = 1 and dil d = 0 and suitably pad the input of each convolutional layer by adding ks d /2 elements around the boundary of each item along any direction d = 1, . . . , D. In particular, we repeat the first and the last input element (reflect padding mode) ks d /2 times at the beginning and at the end of each input item, along each direction d = 1, . . . , D. Moreover, we set the hyperparameter values of the architecture along each direction d = 1, . . . , D as ks d = 2k d -1, where k d is the spline order along direction d, c d in = m (amount of point items) in the first convolutional layer and c d out = D in the last convolutional layer. The input/output channels of the hidden layers are fixed to 100. In particular, their values depends on the user choice, with the only constrain that c d out = c d in for two consecutive layers.

  1] associated to each edge e d I computed in the feature extraction layer, for I = (i 1 , . . . , i D ) with i d = 1, . . . , m d -1, d = 1, . . . , D. CumSum In conclusion, the parametric values are recovered by applying the cumulative sum operation along each direction d = 1, . . . , D. More precisely, we define the parametric values asx d (1,i2,...,i D ) = 0, x d I+d = x d I + ∆ d I ,where I + d is the multi-index I augmented by 1 along the d-th direction, for d = 1, . . . , D. Thanks to the application of the softmax activation function, for each d = 1, . . . , D it is also ensured that the remaining boundary points are mapped to the boundary of the parameter domain, namely x d (i1,i2,...,i D-1 ,m D ) = 1. Output The output of the learning model consists of the parametric values X as defined in[START_REF] Floater | Parameterization for curve interpolation[END_REF], associated to the input point cloud P .

u 4 .

 4 i d = u i d -min j=1,...,m d u j max j=1,...,m d u j -min j=1,...,m d u j . Define a tensor-product mesh of samples U := D d=1 u d , so that U = u I ∈ [0, 1] D : u I , I = (i 1 , . . . , i D ) , i d = 1, . . . , m d , d = 1, . . . , D . 5. Evaluate the curve on the parametric values U , so that p I = s(u I ) and collect them in P . 6. Normalize P in [0, 1]

  97e-3 5.3e-4 3.21e-5 1.79e-3 7.3e-4 2.23e-5 1.46e-3 7.9e-4 4.13e-5 q = 3 2.20e-3 6.7e-4 1.31e-4 1.43e-3 8.7e-4 8.80e-5 1.41e-3 8.5e-4 1.09e-4 q = 4 1.66e-3 3.3e-4 7.24e-5 1.06e-3 4.3e-4 5.30e-5 9.48e-4 4.1e-4 6.13e-5 q = 5 1.35e-3 4.8e-4 3.71e-5 7.65e-4 5.8e-4 2.50e-5 5.09e-4 5.6e-4 3.55e-5

Fig. 5 :

 5 Fig. 5: B-spline approximation of degree q = 4 with 15 interior knots for 20 (left), 50 (center), and 80 points (right) in Section 6.3.

Fig. 6 :

 6 Fig. 6: B-spline approximation of 65 points with q = 5 and 4 (left), 7 (center), 10 (right) interior knots in Section 6.3.

Fig. 7 :

 7 Fig. 7: Point clouds obtained from a car (left), a ship hull (center), and a wind turbine (right) model of 1089 (top) and 3025 (bottom) items.

Fig. 8 :

 8 Fig. 8: Point-wise error on the parametric (a-c) and geometric (b-d) domain of the tensor product B-spline approximation of bi-degree (3, 3) for the car shell point cloud of size 3025 using the given GT (a-b) and the CNN (c-d) parameterizations in Section 6.6.

Fig. 9 :

 9 Fig. 9: Hierarchical mesh (left) and geometric spline model (right) of the ship hull point cloud of dimension 3025 via the adaptive scheme with THB-splines in Section 6.7.

Table 2 ]

 2 (RES), and the ones obtained with our convolutional model (CNN). The outcome of this analysis is reported in Table1. The CNN outperforms STD

Table 1 :

 1 MSE error for Bézier approximation of degree q = 2, 3, 4, 5 on polynomial data for different input lengths m = 20, 50, 80 in Section 6.1.

  .

		exact		noisy
		STD	CNN	STD	CNN
	q = 2	1.20e-3	5.21e-4	1.26e-3	5.94e-4
	q = 3	6.07e-4	2.25e-4	6.45e-4	2.57e-4
	q = 4	3.91e-4	1.46e-4	4.62e-4	1.94e-4
	q = 5	2.19e-4	8.97e-5	2.75e-4	1.35e-4

Acknowledgements CG acknowledges the contribution of the National Recovery and Resilience Plan, Mission 4 Component 2 -Investment 1.4 -CN 00000013 "CENTRO NAZIONALE HPC, BIG DATA E QUANTUM COMPUTING", spoke 6. CG and SI are members of the INdAM group GNCS. The INdAM-GNCS support is gratefully acknowledged. CG, SI, AM also acknowledge the PHC GALILEE project 47786N and AM acknowledges H2020 Marie Sklodowska-Curie grant GRAPES No. 860843.


Marie Collet
email: marie.collet@inrae.fr

Nicolas Rodriguez

Selma Baati

Alain Husson

Eddy Renaud

Kevin Caillaud

Yves Le Gat

Kévin Caillaud

Understanding the leakage process for multi-scale water infrastructure asset management: necessity for a dialogue between sociological and data sciences

Keywords: Data, Human/non-human actors, Infrastructure Asset Management, Pipe leakage model, Sociology, Sociotechnical translation

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

Water supply systems are designed to provide drinkable water with a flow and pressure that should be sufficient for each point of consumption. Water losses occur mainly due to pipe leaks, which are caused by deteriorating infrastructures. They may affect the water distribution network performance and therefore, may lead to an increase in the abstraction of raw water from drinkable water sources. One of the most pressing issues in the management of water utilities is to reduce leaks in order to ensure a sustainable use of water resources [START_REF] Hafskjold | SIP -Future Rehabilitation Strategies for Physical Infrastructure: COORDINATION OF REHABILITATION PLANNING AND MEASURES -CO-INFRASTRUCTURE INTERACTIONS[END_REF][START_REF] Renaud | Gestion patrimoniale des réseaux d'eau potable Volume II :Optimiser ses pratiques pour un service durable et performant[END_REF]. Thus, water systems operators are moving towards developing strategies to control leakages in pipe networks. Apart from service pressure reduction or modulation, such strategies are focused on two major types of action: the implementation of research campaigns which aim at leaks detection and repair as soon as possible; and rehabilitation/renewal of parts of the distribution network which are characterized by the highest frequencies of leakage occurrence [START_REF] Renaud | Réduction des pertes d'eau des réseaux de distribution d'eau potable: guide pour l'élaboration du plan d'actions[END_REF](Renaud et al., , 2012;;[START_REF] Rokstad | An Approach to planning Data Collection for IAM in VMW[END_REF]. The cost of those shares, especially the expense of replacing pipes, may raise the operating cost of water utility, hence water price. Thus, the implementation of an efficient leakage control strategy could be based on tools that assist in the selection of actions that are better suited for water utility management and its goals for short and long terms. An example could be a leak prediction model for water distribution networks, for which an efficient and well-documented water leakage (or failure) detection and repair process generating accurate data is essential.

Data from many sources and in varying formats are collected, stored, and used by water utilities for most of their administrative and business processes (e.g., data on the hydraulic characteristics of the network, its structure, data on the network environment, failure history). Infrastructure Asset Management (IAM) relies on detailed and structured technical information on the water infrastructures and operations involved in supplying water to users, but even more on the understanding of the rationale behind the organization of water utilities and the work practice of its many employees. Yet, these data, when available, are not tailored for the construction of probabilistic models. This is often due to their production and exploitation methods (processing to which they are subject, archiving methods, recording in databases, etc.). Therefore, the pipe leakage state may not be obvious in the raw data, leading to uncertainties.

Based on the experience from an applied research project, called GePaME (Multi-scale IAM), led in collaboration between Société Wallonne des Eaux (SWDE) and INRAE-Bordeaux, this paper presents a set-up of research tasks designed to extract pipe maintenance (e.g., inspection, repair) data embedded in SAP (Enterprise Resource Planning software) records and Geographical Information Systems (GIS) data from SWDE in order to infer a pipe leakage probability model. The goal of this model is to help SWDE improve its long-term IAM process. The originality of this work lies particularly in the link between pipe leakage detection and repairs, as well as in the inter-disciplinary work involved to process the available raw data.

Methods

In fact, the raw data provided by SWDE cannot be used to supply a probabilistic model. In order to make better use of these data, we propose to study and analyse its production and exploitation process in it whole sociotechnical dimension, meaning: actors, tools and technical artefacts (software, databases, instruments, etc.), and current practices in leaks detection and repair, in order to complete and transform existing raw data into information suitable for IAM. This analysis is be partly based on planned interviews with actors involved in generating the raw data at different levels.

For this, analyses mobilizing four disciplines are combined: (i) Data science, to cleanse and connect network segment description and localisation data to leakage detection operations and pipe leak repairs; (ii) Drinking water engineering science, to exploit DMA monitoring data; (iii) Sociology (sociology of Sciences and Technology, sociology of work) to understand how DMA, and pipes within a DMA, are practically selected for leakage inspection, and how inspection data are produced and reported; (iv) Probability and statistics, for building a modelling tool that could account for this cumulated knowledge and the associated uncertainties.

Mainly three problems can be a hindrance to exploiting the vast amount of data that water utilities generate on a daily basis for IAM [START_REF] Okwori | Data-driven asset management in urban water pipe networks: a proposed conceptual framework[END_REF][START_REF] Rokstad | An Approach to planning Data Collection for IAM in VMW[END_REF]. First, water utilities may not always store the crucial metadata describing how, when, and by whom the data were recorded (or updated) and validated. Therefore, the data can contain redundancies and inconsistencies, making the accuracy and priority given to the information not always clear. In addition, various technical terms (e.g., District Metered Area (DMA), pipe section, leak, and failure) can have different meanings for different operators involved in the data production. Second, many technical details linked with work practice and experience are not reported due to the administrative procedures that historically framed the data collection (data overwritten or not archived). Again, extracting information useful for IAM from raw data is difficult, even with the most advanced data mining techniques. Finally, the data may not be structured as a unique and consistent database [START_REF] Rokstad | An Approach to planning Data Collection for IAM in VMW[END_REF]. Ideally, all technical objects (such as pipes) and technical tasks (such as pipe leak repairs) are represented as entities, that is, tables with primary keys (IDs which are not duplicated nor null), attributes (characteristics) which are categorized (not free text), and foreign keys (links to other entities). This structure allows a seamless use of the data by merging different tables to create information useful for probabilistic modelling.

We therefore develop advanced data processing algorithms that overcome these difficulties inherent in the data, with the R software language (R Core Team, 2020).

We eventually develop a statistical model of the pipe leakage state probability based on the Linear Extension of the Yule Model with selective survival (LEYP, Le [START_REF] Gat | Extending the Yule process to model recurrent pipe failures in water supply networks[END_REF] and feed it with the processed data combining the network characteristics, DMA, and the leakage detection and pipe leak repair (technical and sociological data).

Results

This multidisciplinary approach allows a better understanding of the technical and social context in which SWDE's leak detection strategy is employed. The production and exploitation of the data are indeed oriented by the finalities of leak inspection process (regulatory, environmental, economic, administrative, etc.). This process implies many tools/instruments used to generate data and individuals/organizations to implement these tools and interpret these data. Taken into consideration, the production and exploitation data can be presented as a "chain of translations" which links human and non-human entities [START_REF] Akrich | Sociologie de la traduction: textes fondateurs[END_REF]. By moving from one entity to another during production and operation, data undergoes transformations and takes on meaning. Actually, these issues of finalities and the attribution meaning are not known very well.

From the clarification of the chain of translations related to leakage and break detection and repair, we expect a better understanding of how to filter, cleanse and interconnect available raw data in order to identify which DMA, and which pipes inside these DMA's, were inspected and when, and which pipes among these were found to be subject to leakage or breakage. This will then make it possible to build a leakage probability model, as well as a break intensity model, both usable within a long term simulation framework devoted to assess the relevance of pipe inspection and renewal strategies. This mobilization of inter-disciplinary skills requires a common vocabulary and thus sharing different epistemological cultures, from the field agents to the data scientists, via sociology of work. The study of the translation chains along the data course adds essential information that was not recorded in the current database. Assessing the adequacy of this database design and identifying its areas of improvement on a short and long-term perspective for IAM heavily depends on a clear understanding of the organization of SWDE on both on a governance level and on an operational level, which pertains again to sociological practice.

Finally, the shift of temporalities (from daily basis for operational work on the pipe to longterm strategies of network renewal) is taken into account to make sure that the SWDE will be able to make use of and adapt the tools produced and integrate them in their daily practice of data acquisition.

Conclusion

This study thus highlights the need to better understand the processes that underlie the leakage detection operations both technically and sociologically. On the one hand, the leakage detection operations could answer multiple objectives (administrative, financial, technical, regulatory…) that may influence the technics used and/or the data stored. On the other hand, these detection operations are carried out by human field operators, and the question of the impact of their perception and understanding of where and how to implement an operation remains open. The dialogue between sociology, data science, water supply engineering and probability is also at the core of the comprehension of the mutual influence between the technical tools and the realization of the leakage detection operations, as well as their evolutions in order to take into account possible biases in the statistical model.