Victoriya Kashtanova
email: victoriya.kashtanova@inria.fr

Mihaela Pop

Ibrahim Ayed

Patrick Gallinari

Maxime Sermesant
email: maxime.sermesant@inria.fr

Simultaneous Data Assimilation and Cardiac Electrophysiology Model Correction Using Differentiable Physics and Deep Learning

Keywords: Physics-based learning, Deep Learning, Cardiac electrophysiology, Simulations

Modelling complex systems, like the human heart, has made great progress over the last decades. Patient-specific models, called "digital twins", can aid in diagnosing arrhythmias and personalising treatments. However, building highly accurate predictive heart model requires a delicate balance between mathematical complexity, parameterisation from measurements, and validation of predictions. Cardiac electrophysiology models range from complex biophysical models to simplified phenomenological models. Complex models are accurate but computationally intensive and challenging to parameterise, while simplified models are computationally efficient but less realistic. In this paper, we propose a hybrid approach by leveraging deep learning to complete a simplified cardiac model from data. Our novel framework has two components, decomposing the dynamics into a physics-based and a datadriven term. This construction allows our framework to learn from data of different complexity, while simultaneously estimating model parameters. First, using in silico data, we demonstrate that this framework can reproduce the complex dynamics of cardiac transmembrane potential even in presence of noise in the data. Second, using ex vivo optical data of action potentials, we demonstrate that our framework can identify key physical parameters for anatomical zones with different electrical properties, as well as to reproduce the action potential wave characteristics obtained from various pacing locations. Our physicsbased data-driven approach may improve cardiac electrophysiology modeling by providing a robust biophysical tool for predictions.

Introduction

Computational cardiology is a multi-disciplinary field that has seen extensive progress in the past decade [START_REF] Lopez-Perez | Three-dimensional cardiac computational modelling: methods, features and applications[END_REF][START_REF] Richard H Clayton | An audit of uncertainty in multi-scale cardiac electrophysiology models[END_REF][START_REF] Bai | Editorial: Diagnosis, monitoring, and treatment of heart rhythm: new insights and novel computational methods[END_REF]. In particular, recent advances in numerical analysis and the development of virtual patient-specific models (known as 'digital twin') have allowed researchers to address critical challenges related to limitations of clinical methods routinely employed to diagnose arrhythmia, as well as to help planning the best therapy on an individual basis [START_REF] Anouk | From evidence-based medicine to digital twin technology for predicting ventricular tachycardia in ischaemic cardiomyopathy[END_REF]. However, in order to build such accurate predictive heart models, one needs to select the most suitable theoretical framework, balancing the degree of mathematical complexity needed for the specific problem studied, the correct parameterisation of model from measurements, and the validation of predictions.

Physiological and multi-physics phenomena characterising the heart function in normal and pathological conditions can be mechanistically described by mathematical models of different complexity [START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF]. Notably, several electrophysiology (EP) models are able to accurately reproduce the electrical behaviour of the heart at different scales (i.e. cell, tissue, organ). For instance, detailed biophysical models have been proposed in order to describe the dynamics of transmembrane potential, the flowing currents as well as the different ionic concentrations inside and outside the cardiac cell, such as the complex Ten Tusscher-Panfilov model [START_REF] Ten Tusscher | A model for human ventricular tissue[END_REF][START_REF] Ten Tusscher | Alternans and spiral breakup in a human ventricular tissue model[END_REF]. On the other hand, ionic models are not only intricate but also computationally expensive, requiring costly resources. In addition, these sophisticated models also employ numerous hidden variables that are impossible to be all measured, making difficult to accurately identify all model parameters.

An attractive alternative to complex computational models, is to use phenomenological models involving descriptions derived from simple biophysical models adapted to cardiac EP. Examples of such models are the FitzHugh-Nagumo [START_REF] Fitzhugh | Impulses and physiological states in theoretical models of nerve membrane[END_REF][START_REF] Nagumo | An active pulse transmission line simulating nerve axon[END_REF], Aliev-Panfilov [START_REF] Aliev | A simple two-variable model of cardiac excitation[END_REF][START_REF] Nash | Electromechanical model of excitable tissue to study reentrant cardiac arrhythmias[END_REF], and Mitchell-Schaeffer [START_REF] Colleen | A two-current model for the dynamics of cardiac membrane[END_REF] models, which only make use of a few variables and parameters. These models are used for rapid simulations of the action potential wave propagation at tissue level (e.g. on a 2D tissue layer or a 3D tissue slab), and organ level (through the entire heart). Nevertheless, being less realistic, these simplified models are also less accurate. Given the increasing availability of synthetically-generated data and observation data, an important research issue is to consider how machine learning (ML), and more specifically, deep learning (DL), could help complementing simple physical models in order to improve their accuracy.

Physics-aware deep learning is a recent field of research aiming at promoting the use of data-intensive methods for the modeling of complex physical phenomena [START_REF] George | Physics-informed machine learning[END_REF][START_REF] Willard | Integrating Scientific Knowledge with Machine Learning for Engineering and Environmental Systems[END_REF][START_REF] Wang | Physics-Guided Deep Learning for Dynamical Systems: A Survey[END_REF]. This research topic motivates works from different disciplines, ranging from climate to aeronautics and biology, by encompassing diverse objectives including: accelerating numerical simulations; improvement over physical models; building of emulators; solving differential equations in large variable spaces; discovering physical laws from data; etc. The methods developed for reaching these objectives are also multiple-folded: incorporating prior physical background in the loss functions [START_REF] Raissi | Physics-informed neural networks: A deep learning framework for solving forward and inverse problems involving nonlinear partial differential equations[END_REF][START_REF] Raissi | Hidden physics models: Machine learning of nonlinear partial differential equations[END_REF][START_REF] Yang | B-pinns: Bayesian physics-informed neural networks for forward and inverse pde problems with noisy data[END_REF] or as strict constraints in the NN architecture [START_REF] Chiyu | Enforcing physical constraints in cnns through differentiable pde layer[END_REF]; enforcing invariance or equivariance properties characteristic of physical laws in DL models [START_REF] Wang | Incorporating symmetry into deep dynamics models for improved generalization[END_REF][START_REF] Wang | Approximately equivariant networks for imperfectly symmetric dynamics[END_REF], combining DL and physical components into hybrid systems [START_REF] Yuan Yin | Augmenting physical models with deep networks for complex dynamics forecasting[END_REF].

Pure DL methods have been also increasingly used in order to learn dynamical models from data and make intelligent decision without human intervention. Notably, neural networks are capable to learn and to model various complex (e.g. nonlinear) relationships between input and output data [START_REF] Willard | Integrating physics-based modeling with machine learning: A survey[END_REF]. For example, a framework useful for automatic learning PDEs from data [START_REF] Long | PDE-net: Learning PDEs from data[END_REF][START_REF] Long | PDE-net 2.0: Learning PDEs from data with a numericsymbolic hybrid deep network[END_REF] has been proposed. Another group used an adjoint method to learn differential equations parameterised with neural networks [START_REF] Ricky | Neural ordinary differential equations[END_REF], while Ayed et al. [START_REF] Ayed | Learning dynamical systems from partial observations[END_REF] proposed a framework for learning models using a purely data-driven approach in partially observable settings.

However, despite achieving good progress and producing promising results in cardiac EP simulations [START_REF] Ayed | EP-Net: Learning cardiac electrophysiology models for physiology-based constraints in data-driven predictions[END_REF][START_REF] Kashtanova | EP-Net 2.0: Out-of-domain generalisation for deep learning models of cardiac electrophysiology[END_REF], these works suggested that data-driven models alone were not able to reproduce complex unseen dynamics (such as the repolarisation phase of an action potential) [START_REF] Kashtanova | EP-Net 2.0: Out-of-domain generalisation for deep learning models of cardiac electrophysiology[END_REF]; thus, the maximum forecasting horizon remains limited.

Subsequently, several researchers have started to use coupled physico-statistical approaches for cardiac EP simulations, in order to achieve a high precision at low cost. For example, one group designed a neural network that approximates the FitzHugh-Nagumo model [START_REF] Court | Design of the monodomain model by artificial neural networks[END_REF], while others used a physics-informed neural networks to construct cardiac activation maps by accounting for the underlying wave propagation dynamics [START_REF] Sahli Costabal | Physicsinformed neural networks for cardiac activation mapping[END_REF] and to estimate the cardiac fiber architecture of the human atria from catheter recordings of the electrical activation [START_REF] Herrera | Physics-informed neural networks to learn cardiac fiber orientation from multiple electroanatomical maps[END_REF]. Furthermore, another group proposed an approach to create a nonlinear reduced order model by employing deep learning algorithms (DL-ROM) designed for cardiac EP simulations [START_REF] Fresca | POD-enhanced deep learningbased reduced order models for the real-time simulation of cardiac electrophysiology in the left atrium[END_REF]. Finally, other researchers presented a physics-informed neural network for an accurate simulation of action potential wave and a correct estimation of the model parameters [START_REF] Herrero Martin | EP-PINNs: Cardiac electrophysiology characterisation using physics-informed neural networks[END_REF]. Unfortunately, the majority of these coupled approaches use a high-fidelity physical model as a core component of its structure. As a result, fitting those models to real data may not only be computationally expensive, but also difficult especially in order to properly deal with frequently observed large discrepancies between simulated and real data.

To address this critical limitation, here we propose a framework to Augment incomplete PHYsical models with a deep learning component for ideNtifying complex cardiac ElectroPhysiology dynamics (APHYN-EP) from data, based on a fast low-fidelity (or incomplete) physical model. This framework has two main components which decompose the model dynamics into a physical term and a data-driven term, respectively. The data-driven DL component is designed such that it captures only the information that cannot be modeled by the incomplete physical model. This construction allows our framework to learn from data of different complexity and to handle the different temporal and spatial resolution. The proposed framework closely follows the approach of Yin et al. [START_REF] Yuan Yin | Augmenting physical models with deep networks for complex dynamics forecasting[END_REF]. However, in contrast to this previous work (which considers fully-observable dynamics and simple test use cases), cardiac EP dynamics have a high complexity and represent simultaneously multiple underlying processes. Furthermore, most cardiac EP models lack measurements for some of the model variables, which requires inferring the dynamics from incomplete observations only, making a model partially-observable.

Figure 1 presents the general scheme of our approach. The observed data, which represent only part of cardiac EP dynamics, is decomposed by the framework into two "physics-based" and "data-driven" parts. The "data-driven" part acts as a complement (or a correction) to the "physics-based" part, and handles quantities and dynamics which were neglected. The parameters of the physical model (inverse problem) and of the neural network (direct problem) were learned simultaneously during the training phase (see Fig. 2). After the complete training, the APHYN-EP framework can be used to forecast the learned dynamics at multiple horizons. We showed in previous papers [START_REF] Kashtanova | Deep learning for model correction in cardiac electrophysiological imaging[END_REF][START_REF] Kashtanova | Aphyn-ep: Physics-based deep learning framework to learn and forecast cardiac electrophysiology dynamics[END_REF] the preliminary results of the this framework application. In this current paper, we present in detail our novel framework, demonstrating results of its application on more complex settings. Furthermore, we also explore the generalisation capabilities of the framework.

In Section 2, we give a thorough overview of the method, with training details and a comparison with a traditional numerical method. In Section 3, we introduce the data used to test different abilities of the framework. In Section 4, we present numerical results obtained using our framework in comparison to the baseline methods. Finally, in Section 5, we discuss the possible framework improvements and future applications.

Learning Framework

The cardiac EP dynamics can be driven by an equation of the form:

dX t dt = F (X t ) (1)
defined over a finite time interval [0, T ], where the state X t is a spatio-temporal vector field over the domain Ω ∈ R k , with k ∈ {0, 2} in our experiments. We suppose that we have access to a set of observed trajectories:

D = {X : [0, T ] → A|∀t ∈ [0, T ], dX t dt = F (X t )},
where A is the set of vector field values X. In our case, the unknown F has A as domain and we only assume that F ∈ F, with (F, || • ||) a normed vector space. Since we consider that only part of the mesuared dynamics can be modeled via a family of PDEs characterised by their temporal evolution F p ∈ F p ⊂ F, APHYN-EP framework introduce a data-driven augmentation term F d ∈ F complementing F p . F being a vector space, thus we can write:

F = F p + F d , (2)
where F p represents a physical model (which is an incomplete description of the underlying phenomenon) and F d is a neural network complementing the physical model by capturing the information that cannot be modeled by the physics-described component. The physical and data-driven parameters defining F p and F d are unknown and need to be estimated from data by fitting the trajectories from D.

The decomposition (2) is not unique, as shown in [START_REF] Ayed | EP-Net: Learning cardiac electrophysiology models for physiology-based constraints in data-driven predictions[END_REF][START_REF] Kashtanova | EP-Net 2.0: Out-of-domain generalisation for deep learning models of cardiac electrophysiology[END_REF], all the measured dynamics could be captured by the F d component alone. This decomposition is therefore ill-defined, which complicates the interpretability and extrapolation possibilities of the framework. However, we aim to reduce the action of F d on the dynamics measured through its norm, and to find its optimal minimum. Therefore, in order to learn the cardiac EP dynamics X t , in this work, we solve following optimisation problem via our physics-based data-driven APHYN-EP framework:

min Fp∈Fp,F d ∈F ∥F d ∥ subject to ∀X ∈ D, ∀t ∈ [0, T ], dX t dt = F (X t ) = F p (X t ) + F d (X t ). (3)
Assuming that F p is a Chebyshev set, Propositions 1 and 2 from Yin et al. [START_REF] Yuan Yin | Augmenting physical models with deep networks for complex dynamics forecasting[END_REF] guarantee the existence and uniqueness of a minimising pair for [START_REF] Ayed | Learning dynamical systems from partial observations[END_REF].

Specifically, in our experiments the incomplete physical model is derived (F p (X t ) = v t ) from the two-variable (v, h) model [START_REF] Colleen | A two-current model for the dynamics of cardiac membrane[END_REF] for cardiac EP simulations, as described by equations (4). The variable v represents a normalised (v ∈ [0, 1]) dimensionless transmembrane potential, while the "gating" variable h controls the repolarisation phase (i.e., the gradual return to the initial resting state):

∂ t v = div (σI∇v) + hv 2 (1 -v) τ in - v τ out + J stim ∂ t h = 1-h τopen if v < v gate -h τ close if v > v gate (4)
where J stim is a transmembrane potential activation function, which is equal to 1 during the time the stimulus is applied (t stim ) in a certain stimulated area.

This physical model has been successfully used in patient-specific modelling [START_REF] Relan | Coupled personalization of cardiac electrophysiology models for prediction of ischaemic ventricular tachycardia[END_REF], covering general EP dynamics. Furthermore, in contrast to the very detailed ionic/cellular models, this model is flexible in terms of spatial and temporal steps set in the numerical analysis. Thus, assuming the initial conditions for this system (4) v(t = 0) = 0 and h(t = 0) = 1 we can compute an approximation of h for any time point t by employing a simple integration scheme.

In the experiments presented later (see Section 3), F p is the set of models spanned by the R.H.S. of the equations above for varying variables σ, τ in , τ out , τ close . This is a finite dimensional vector subspace which is indeed Chebyshev, thus falling under the assumption guaranteeing theoretical existence and uniqueness of a minimising pair.

The data-driven component (F d ) of the framework was implemented via a neural network. The choice of a neural network depends on the application problem and the dimension of the data. In this paper, we used a ResNet network [START_REF] He | Deep residual learning for image recognition[END_REF], because it could accurately reproduce complex cardiac EP dynamics [START_REF] Ayed | EP-Net: Learning cardiac electrophysiology models for physiology-based constraints in data-driven predictions[END_REF][START_REF] Kashtanova | EP-Net 2.0: Out-of-domain generalisation for deep learning models of cardiac electrophysiology[END_REF]. However, a simpler neural network could also be used for more rapid computations, as discussed in Section 5.

In APHYN-EP framework the physical (F p ) and the data-driven (F d ) components are trained simultaneously. This insures the finding of the best minimising pair for (3) determined by the set of parameters θ = (θ p , θ d ) based on observed measurements

X (i) h∆t , i ∈ [1, N ] and h ∈ [1, T /∆t].
The 'Loss function' (L) in training phase consisted of 2 parts: trajectory-based loss (L traj ) and loss on norm of F d , being represented as following:

L(θ) = λ * L traj (θ) + F θ d d = λ * N i=1 T /∆t h=1 ||X (i) h∆t - X(i) h∆t (θ)|| + F θ d d ( 5
)
where each state

X(i) h∆t (θ) = X (i) 0 +h∆t X (i) 0 (F θp p + F θ d d )(X s ) dX s was calculated from the initial state X (i)
0 via a differentiable ODE solver [START_REF] Ricky | Neural ordinary differential equations[END_REF][START_REF] Ricky | Learning neural event functions for ordinary differential equations[END_REF].

The key role of the λ coefficient is to balance the two parts of the loss. During training, we used λ in a dynamic state, as λ j+1 = λ j + γL traj (θ j+1 ) (j is an epoch number), in order to artificially increase the importance of L traj at the beginning of training and then to gradually decrease it, changing the focus of optimisation on the norm of F d .

This training algorithm, adapted from Yin et al. [START_REF] Yuan Yin | Augmenting physical models with deep networks for complex dynamics forecasting[END_REF], is presented below. Figure 2(a) is a schematic representation of the APHYN-EP training algorithm 1. This algorithm is similar to the traditional numerical method of Gradient Descent for solving Inverse and Ill-posed

Algorithm 1 APHYN-EP training

Initialisation: θ 0 , λ 0 ≥ 0, γ > 0; Problems, which is based on the solution of an adjoint problem to calculate the gradient of the Loss function [START_REF] Sergey | Inverse and Ill-posed Problems: Theory and Applications[END_REF] (Fig. 2(b, yellow part)). However, thanks to the automatic differentiation tools provided by the Pytorch library [START_REF] Paszke | Pytorch: An imperative style, highperformance deep learning library[END_REF], the gradient of ( 5) is calculated automatically inside of PyTorch framework (Fig. 2(b, blue part)).

for epoch = 1 : N epochs do for batch in 1 : B do L traj (θ j ) = N i=1 T /∆t h=1 ||X (i) h∆t - X(i) h∆t (θ j )|| θ j+1 = θ j -∇ λ j L traj (θ j ) + F
Additionally, in order to train simultaneously the physical and the data-driven components of APHYN-EP, we implemented the Laplace operator in (4) with a simple finite-difference scheme. Lastly, in order to avoid potential difficulties associated with the high time resolution required in this numerical scheme, we used two different time steps in the integration schemes for the computing the physical component and for computing the final forecast given by the framework, respectively.

Experiments

In order to test the performance of our APHYN-EP framework and to further show its capability to reproduce transmembrane potential dynamics of different complexities, we chose two types of experiments. First, using synthetic in silico data, we tested the ability of the framework to reproduce the complex dynamics of transmembrane potential including a case where noise is present in the data. Second, using optical fluorescence imaging data of action potentials recorded ex vivo on explanted porcine hearts, we aimed to show that our framework: a) can identify key physical parameters for different anatomical zones having abnormal electrical function (using 0D data from a heart with an ischemic region); b) is capable to reproduce the action potential wave characteristics obtained following heart stimulation from different locations (2D mapping data from a healthy heart).

The details of data collection and training settings used for the experiments are presented in detail below. 3.1 In silico data

Data collection

To evaluate our method, we used a dataset of transmembrane potential activations simulated by employing a monodomain reaction-diffusion equation and the Ten Tusscher -Panfilov ionic model [START_REF] Ten Tusscher | Alternans and spiral breakup in a human ventricular tissue model[END_REF], which represents twelve different transmembrane ionic currents. The simulations were performed using a spatial step of 0.2 mm and a time step of 1 ms (similarly to those used in the original model [START_REF] Ten Tusscher | Alternans and spiral breakup in a human ventricular tissue model[END_REF]), with the open-source finetwave software1 . For this, our computational domain was chosen to represent a 2D slab of cardiac tissue (isotropic), with 24 × 24 elements in size. For one data sample, in order to activate the transmembrane potential, an excitation pulse delivered via a stimulus was applied for 1 ms on a selected area. Each simulation represented 400 ms of a heart beat, and was intended to achieve a full depolarisation-repolarisation cycle (note that video examples showing simulated dynamics across time are available online2 ). This required 20 seconds of computation time on a 8-core Intel i7-7820HQ CPU per simulation.

0D data

In order to generate 0D databases, we selected the 2D data samples with the pacing stimulus applied as follows: on the left top corner for training datasets, and near the center for validation datasets, as shown in Figure 3. Next, we saved in separate files a time sequence for each pixel of the 2D cardiac slab, creating two databases: one for training and one for validation, respectively. To simplify the workflow for our framework, we removed the section of time sequences where the transmembrane potential was equal to zero (for non-activated at time 0 pixels), thus obtaining a time sequence of 350 ms per data sample, as shown in Figure 3(c). To test the ability of the framework to operate with noisy data, we added to each time sequence a 5 percent random noise characterised by a normal Gaussian distribution.

2D data

In order to generate 2D databases, we applied an excitation pulse at each grid point of 2D cardiac tissue slab. Next, we added a 5 percent random noise with normal Gaussian distribution directly on simulated data. Finally, we obtained a 2D database of around 500 training samples and 100 validation samples.

The data simulated via the Ten Tusscher -Panfilov model with added noise were considered here as the ground truth. The objective was then to learn the complex dynamics generated via this model using the APHYN-EP framework, by combining a simplified physics description with a deep learning component. We hypothesised that this approach will result in a low computational cost surrogate model of the computationally intensive, biophysically detailed Ten Tusscher -Panfilov model.

Training settings

The physical model (F p ) described by Eq. 4 was implemented with a standard finite-difference scheme for the Laplace operator, using a spatial resolution of 1 mm 2 pixels (absent for 0D data experiments) and an inner time resolution of 0.1 ms. We estimated only σ, τ in , τ out and τ close as unknown parameters in (4), since these control the major part of the model's dynamics (i.e., the main difference between the Mitchell-Schaeffer and the Ten Tusscher -Panfilov models in our simulations). The other Mitchell-Schaeffer model parameters were taken from the original paper [START_REF] Colleen | A two-current model for the dynamics of cardiac membrane[END_REF], as follows: τ open = 120, v gate = 0.13 and t stim = 1.

For the deep learning component (F d ) of the framework, we used: a ResNet [START_REF] He | Deep residual learning for image recognition[END_REF].with 4 input/output channels (assimilating the first milliseconds of dynamics). In order to demonstrate and to compare APHYN-EP framework's performance with simpler neural network, we also tested an MLP network (for 0D experiments) and a ConvNet (for 2D experiments) as a deep learning component (F d ) of the framework.

We used a time resolution of 1 ms to compute the forecast given by APHYN-EP framework. The training was performed using a horizon of 350 ms. More details on the training configurations can be found in the Supplementary material.

It is important to emphasise that while the framework took a few hours to train (e.g. approximately 3h on Nvidia Quadro M2200 GPU in case of 2D data), once this was done, the inference step was rapidly computed (i.e., less than 10 sec to compute 350 ms of 2D forecasting) and did not require any re-calibration. Furthermore, the major part of the inference time is taken by computing the integral of X(i)

T from X (i) 0 and could be reduced by using a larger time step.

Ex vivo data

Data collection

As mentioned above, we tested the APHYN-EP framework performance using ex vivo experimental datasets from optical fluorescence imaging of action potential. Briefly, the optical signals were recorded ex vivo from the epicardial surface of hearts explanted from juvenile swine ( 25kg in weight). Each heart was attached to a Langendorff perfusion system and a voltage sensitive dye (i.e., di-4-ANEPPS) was injected into the perfusate solution. In order to avoid cardiac motion artifacts during the optical recordings, the heart contraction was suppressed by a bolus of saline and Cytochalasin D, an electromechanical uncoupler. All optical images were acquired using a high-speed CCD camera (MICAM02, BrainVision Inc. Japan), with high-temporal resolution (3.7ms) as well as a high-spatial resolution (i.e., pixel size 0.7mm x 0.7mm). The action potential was then derived at each pixel from the relative change in the intensity of fluorescence signal. The ex vivo studies were described in more detail in Pop et al. [START_REF] Pop | Fusion of optical imaging and MRI for the evaluation and adjustment of macroscopic models of cardiac electrophysiology: A feasibility study[END_REF].

The recorded optical signals were saved and exported from the BV-Ana acquisition software. Notably, within each pixel, the relative change in optical signal intensity is directly proportional with the action potential.

0D data Our specific objective was to learn the complex dynamics of the measured action potential (AP), and then to identify its relevant physical parameters for different zones of the heart. Note that the optical data were considered as the ground truth. For this 0D experiment, we chose to use epicardial optical data recorded from a heart having an ischaemic region (artificially induced by occluding a blood vessel). For the analysis, we manually selected two rectangular regions of interest (ROIs) with different AP dynamics across time, as seen in Figure 4. For post-processing of the data, we normalised the optical signal intensity in order to obtain a [0, 1] min/max interval for the transmembrane potential, while keeping the noise in the data. Next, we took a first full cardiac cycle and removed the parts with zero potential, keeping only time-sequences of 300 ms per experiment. Then, we saved in separate files a time sequence for each pixel in the two selected ROI, creating two databases (ROI A and ROI B), each containing approximately 10 and 5 time-sequences for training and validation, respectively.

2D data For this experiment, we chose a healthy heart stimulated via a pacing electrode placed onto two different zones within the heart (see Fig. 5), having two recordings of depolarisation wave propagation. For one recording, we normalised the optical signal intensity in order to obtain a [0, 1] min/max interval for the transmembrane potential, slightly denoising 2D data. We manually selected a rectangular region of interest (see Fig. 5), where there were no changes in the tissue physical properties. Next, we divided this region of the heart into smaller 2D squares of 10 pixels. For each square, we took , where the excitation pulse was delivered), while the late depolarisation times are depicted in green-blue, color scale is represented in ms.

Training settings

The training settings for this experiment were similar to the ones described in section 3.1.2, except for training horizon (300 ms).

Results

4.1 In silico data

0D data

Results demonstrated that our proposed APHYN-EP framework was able to accurately reproduce several key features, the wave morphology and the electrical conduction properties of the transmembrane potential solution generated by the Ten Tusscher -Panfilov model (see Fig. 6), even in presence of noise in data. Figure 6(a) shows that in absence of learning of parameter τ close (controlling the repolarisation), the data-driven component (ResNet model) completed the dynamics generated by physical component. The value of τ close was fixed at 150 [START_REF] Colleen | A two-current model for the dynamics of cardiac membrane[END_REF]. In the presence of learning of τ close the error of dynamics corrected by data-driven component is minimal (Fig. 6(b)). Note that the predicted dynamic was generated via an Euler integration scheme, by assimilating only one first measurement of the transmembrane potential dynamics. Overall, the framework demonstrated robustness and was not sensitive to the noise in data, and, as a result, rapid changes in transmembrane potential activation were neglected, as observed during the first 40 ms of the action potential duration presented in Figure 6.

2D data

We include here our qualitative results obtained for the forecast over 8 ms, after assimilating only one first frame of dynamics (see Fig. 7). These first 8 ms (i.e., the AP upstroke) represent an important part of the cardiac dynamics, ranging from the earliest depolarisation phase to the full depolarisation phase. Importantly, one can observe a very good agreement between the ground truth and the forecast transmembrane potentials generated by APHYN-EP, as illustrated in Figure 8. The effect of the correction term introduced by F d is clearly visible. Moreover, in Figure 8(b), one can also observe that for the 2D data (which includes now a diffusion part) APHYN-EP framework achieves a good precision in transmembrane potential forecasting even when noise is present in the data. GT v [START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF][START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF] Prediction FW v [START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF][START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF] Prediction Fp v [START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF][START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF] Prediction Fd v [START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF][START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF] GT data v [START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF][START_REF] Camara | Inter-model consistency and complementarity: Learning from ex-vivo imaging and electrophysiological data towards an integrated understanding of cardiac physiology[END_REF] was created by the physical component of the framework (as seen in the second row of Fig. 9). The data-driven component was used only to correct the difference between the ground-truth dynamics and the physical part (Fig. 9 (third and fourth rows)). Table 1 presents the mean squared error (MSE) results for our framework on the training and validation data samples. Note that to calculate this error, for each data sample, we fed the model with only one initial test measurement, then allowed the model to predict 300 ms forward without any additional input information. Furthermore, for comparison, we also added two baseline models: the "incomplete" physical model (F p from APHYN-EP framework, trained alone) and a fully data-driven model (EP-Net 2.0 [START_REF] Kashtanova | EP-Net 2.0: Out-of-domain generalisation for deep learning models of cardiac electrophysiology[END_REF]) trained on the same dataset as APHYN-EP described in 3.1(2D). We clearly noticed that APHYN-EP captured the observed dynamics with good precision for a large time horizon (400 ms) and also outperformed the physical model for every dataset. At the same time, the pure data-driven model encountered difficulties to learn the proper dynamics.

Generalisation ability of APHYN-EP: Planar wave Since our objective was to train a model able to generalise to new conditions, outside of the training environment, we performed a test on out-of-domain data represented by planar wave dynamics (see Fig. 10). One can observe that APHYN-EP (trained on data from 3.1.1) has successfully generated the forecast of the new transmembrane potential wave dynamics. 4.2 Ex vivo data

0D data

Using optical imaging mapping data, our APHYN-EP framework was able to reproduce the observed action potential dynamics for different ROIs within the heart, identifying the 3 major physical dynamics parameters (τ in , τ out and τ close ). Figure 11 demonstrates that the framework correctly estimated the difference in value for the parameter τ close , which either increased APD or shortened it, respectively. Table 1 summarises the quantitative results for our framework forecasting on train and validation data samples, in comparison to baseline methods trained on the same data. The obtained MSE is relatively small for both ROI, and, despite the use of a limited dataset for training, the APHYN-EP framework achieved forecasting the dynamics with good accuracy for new data samples from the validation dataset. Furthermore, our framework clearly outperformed the physical model for every dataset, while the contribution of F d component was still minimal. Despite having a good results on ROI B, the pure data-driven model encountered difficulties to learn the dynamics from ROI A data.

Generalisation ability of APHYN-EP: Fast personalisation Having noted that the obtained correction of the DL component is similar for both ROIs (Fig. 11), we assumed that it is possible to obtain a single DL correction term (F d ) suitable for each selected ROI in the whole heart.

To prove our hypothesis, we first performed a test where we simply replaced the F d component of APHYN-EP trained on data in one ROI by the the F d component of APHYN-EP trained on data from another ROI (and vice versa), and obtained confirming results.

Second, we preformed a series of personalisation experiments. Since we can obtain an estimation of physical parameters specific for each ROI (e.g. by training only the F p component of the framework), we fixed the estimated θ p parameters (for ROI A and ROI B, respectively) in the APHYN-EP framework and trained only the F d component (separately for each ROI). Next, we replaced in the framework the F d component trained on ROI A data (F d A (A)) by the F d component trained on ROI B data (F d B (B)) (and vice versa), and achieved similar results in forecasting (see Fig. 12131415).

Figures 13 and15 visually confirm our initial hypothesis, suggesting that in order to obtain an appropriate forecast of 0D ex vivo dynamics for any region selected on the epicardium, it is sufficient to train the whole APHYN-EP framework using only one selected ROI.

2D data

Finally, we observed that APHYN-EP framework was able to reproduce the features of action potential wave from 2D optical mapping data (see Fig. 16 and Fig. 17). Table 1 shows the quantitative results for our framework forecasting on train, validation and out-of-domain test data samples in comparison to baseline methods trained on the same data. One can notice that the Physical model and a simple Data-driven model outplay the APHYN-EP framework on training and validation datasets, but they have worse performance on a test dataset, which indicates either their possible overfitting or disability to generalise to new conditions. It is also important to note that APHYN-EP has the best results of forecasting for first 150 ms (see Table 1) and retains persistent performance on out-of-domain test data.

Generalisation ability of APHYN-EP: LV pacing As shown in Figure 17, the framework keeps the capability to generalise to unseen conditions (e.g. LV pacing). The absolute error was slightly larger than that on images obtained with RV pacing (used for training); however, this error was still acceptable. The quantitative results are provided in Table 1.

Discussion

In this work, we developed and tested a robust learning framework that is able to assimilate cardiac electrophysiology dynamics. Overall, our results suggest that automated learning of cardiac EP dynamics is feasible and has great potential in predicting the features of action potential wave. In particular, we demonstrated that the APHYN-EP framework can master the EP data of different origin (i.e, data simulated via the Ten Tusscher -Panfilov ionic model, as well as real optical mapping data), even in the presence of noise in the data. This framework may be useful for applications concerning fast parameterisation of computational models that predict electrical wave propagation through the heart.

In the section 4.2.1, we demonstrated that APHYN-EP framework can identify the major physical dynamics parameters for different heart regions. Alongside, we also demonstrated an example of fast personalisation for the framework.

Moreover, in sections 4.1.2 and 4.2.2 we illustrated the generalisation abilities of our framework, by applying the framework on new data without the need of further re-training or a change in its parameters.

However, we acknowledge that the performance of the framework depends on many hyper-parameters, such as a choice of DL component of the framework, or of the training strategy, etc. For instance, due to the data-driven architecture, our framework training could lead to a local minimum for the physical component parameters and additional involvement of the DL component. Several solutions for this problem may be given by the rigid boundaries on physical model parameters or a more strict training protocol such as separate sequential training of physical and DL components. Furthermore, we observed that the usage of simple neural networks for the DL component of our APHYN-EP framework, is able to improve the results (see Table 1); however, the choice of such models depends entirely on the learning data configurations such as data dimension, presence/absence of noise, etc. In this work we have focused mostly on results generated by the framework with DL component represented by a ResNet network, that is because of its robustness and capability of resting stable during different simulations. Additionally, thanks to its residual connections, ResNet can accurately reproduce complex cardiac EP dynamics [START_REF] Ayed | EP-Net: Learning cardiac electrophysiology models for physiology-based constraints in data-driven predictions[END_REF][START_REF] Kashtanova | EP-Net 2.0: Out-of-domain generalisation for deep learning models of cardiac electrophysiology[END_REF].

Conclusion

In this article we successfully demonstrated the ability of a novel APHYN-EP framework to learn cardiac EP dynamics from data of different complexity. The main advantage of our proposed framework is its coupled architecture, which allowed us to use a simplified low-fidelity EP model as a physical component of the framework, along with a neural network as a data-driven correction mechanism for the EP model. Our original framework opens up several possibilities in order to introduce prior knowledge in deep learning approaches through explicit equations, as well as to correct the physical model errors from assimilated data.

Figure 1 :

 1 Figure 1: General APHYN-EP framework. Correction of a PDE-based model with an additional term learned from the data.

Figure 2 :

 2 Figure 2: (a) General APHYN-EP framework training scheme. (b) Calculation of the Loss function gradient via Pytorch framework (blue part) vs. Adjoint model (yellow part).

Figure 3 :

 3 Figure 3: (a,b) Example of a selected 2D myocardial tissue slab with the transmembrane potential activation (yellow) and resting phase (dark blue), for train and validation dataset, respectively. (c) Typical temporal sequence for the simulation experiment (without noise), with the normalised amplitude of the transmembrane potential being represented on the Y-axis and the time (ms) on the X-axis.

Figure 4 :

 4 Figure 4: Example of mapping data from the optical experiment: (left) explanted porcine heart, and (right) tracings of denoised action potential waves recorded from the epicardial surface of the left ventricle, LV. Note that ROI B represents an ischaemic region characterised by a shorter action potential duration (APD) compared to the normal APD recorded from ROI A.

Figure 5 :

 5 Figure 5: Example of depolarisation maps selected from the optical data recorded ex vivo in a healthy porcine heart, with the stimulating electrode (used for pacing the heart) placed onto the right ventricle (a), and left ventricle (b), respectively. Red areas correspond to early activation times (i.e., where the excitation pulse was delivered), while the late depolarisation times are depicted in green-blue, color scale is represented in ms.

Figure 6 :

 6 Figure 6: Validation results of the trained framework with learning of: (a) 2 (τ in and τ out ) and (b) 3 (τ in , τ out and τ close ) physical parameters. Legend: ground truth (GT), prediction of the framework (Prediction FW), decomposition of prediction on physical (F p ) and DL (F d ) components.

Figure 7 :

 7 Figure 7: APHYN-EP predicted dynamics for the transmembrane potential diffusion. Figure shows a period of 8 ms of the forecast. Red point is the reference point for Figure 8.

Figures 8

 8 Figures 8 and 9 present the performance of different components of APHYN-EP and associated contribution to the final result. We can observe which part of the generated transmembrane potential

Figure 8 :

 8 Figure 8: Transmembrane potential at point (5,5) in the cardiac slab (red point, see Fig. 7): (a) Original, (b) Zoom-in of first 40 ms. Legend: ground truth (GT), APHYN-EP, physical (F p ) and data-driven (F d ) component of APHYN-EP.

Figure 9 :

 9 Figure 9: Exemplary results illustrating: the dynamics of the transmembrane potential diffusion predicted by the APHYN-EP physical component (second row); the error with ground-truth diffusion for this physical component of APHYN-EP (third row); and, the trained APHYN-EP data-driven component contribution (bottom row).

Figure 10 :

 10 Figure 10: APHYN-EP predicted dynamics for the transmembrane potential diffusion of plannar wave. The frames show a period of 8 ms of forecast obtained without re-training the APHYN-EP framework.

Figure 11 :

 11 Figure 11: Validation results of the framework trained on: (a) ROI A data and (b) ROI B data, respectively. Ground truth (GT) data, prediction of the framework (Prediction FW), decomposition of prediction on physical (F p ) and DL (F d ) components.

Figure 12 :Figure 13 :

 1213 Figure 12: (a) Validation results of the F p trained alone on ROI A data to obtain θ p (A) parameters. (b) Validation results for APHYN-EP framework (with fixed θ p (A) parameters) trained on ROI A data. (c) Validation results of the APHYN-EP framework (with fixed θ p (A) parameters) with F d component trained on ROI B data. Legend: ground truth (GT) data, prediction of the framework (Prediction FW), decomposition of prediction on physical (F p ) and DL (F d ) components.

Figure 14 :Figure 15 :

 1415 Figure 14: (a) Validation results of the F p trained alone on ROI B data to obtain θ p (B) parameters. (b) Validation results of the APHYN-EP framework (with fixed θ p (B) parameters) trained on ROI B data. (c) Validation results of the APHYN-EP framework (with fixed θ p (B) parameters) with F d component trained on ROI A data. Legend: ground truth (GT) data, prediction of the framework (Prediction FW), decomposition of prediction on physical (F p ) and DL (F d ) components.

Figure 16 :

 16 Figure 16: APHYN-EP predicted dynamics for the transmembrane potential diffusion, RV pacing. The frames illustrate time sequences from a 30 ms period of forecast.

Table 1 :

 1 Mean-squared error, MSE (x 10 -3 ) of the normalised transmembrane potential (adimensional) forecasting (forecasting horizon of 400 ms for In silico experiments and of 300 ms for Ex vivo experiments). Baseline models: the Physical model (4) and a fully data-driven model (EP-Net 2.0[START_REF] Kashtanova | EP-Net 2.0: Out-of-domain generalisation for deep learning models of cardiac electrophysiology[END_REF]) trained on the same dataset as APHYN-EP. Out-of-domain tests: Plannar wave for In silico 2D dataset, LV pacing for Ex vivo 2D dataset.

	Dataset	Method	Training	Validation	Out-of-
				data (||F d || 2 )	data (||F d || 2 )	domain
						test (||F d || 2 )
	In silico	2D	APHYN-EP ResNet	2.54 (0.47)	2.54 (0.472)	4.2 (0.4)
			APHYN-EP ConvNet	2.5 (0.81)	2.5 (0.8)	4.4 (1)
			Physical model	5.7	5.6	4.3
			Data-driven model	10	10	100
	Ex vivo	0D ROI A APHYN-EP ResNet	9.12 (0.16)	5.72 (0.08)
			APHYN-EP MLP	9 (0.0785)	5.37 (0.077)
			Physical model	14	10
			Data-driven model	20	9.78
		ROI B APHYN-EP ResNet	10 (0.08)	8 (0.07)
			APHYN-EP MLP	8.79 (0.18)	7 (0.21)
			Physical model	14.5	9.3
			Data-driven model	7.78	6.79
		2D RV	APHYN-EP ResNet	54.5 (92.6)	53.26 (94.2)	49.9 (134.6)
			APHYN-EP ConvNet	60.3 (105.6)	58.9 (108.1)	42 (128.6)
			Physical model	39.2	37.5	106.4
			Data-driven model	33.89	31.7	93
		RV	APHYN-EP ResNet	8.9 (110.85)	7.66 (109.9)	8.93 (135.8)
		150	APHYN-EP ConvNet	8.26 (101.16)	6.93 (107.5)	7.47 (124.4)
			Physical model	20	18.42	32.5
			Data-driven model	16.84	15.23	14.9

https://github.com/TiNezlobinsky/Finitewave

https://doi.org/10.6084/m9.figshare.21648752.v3

Acknowledgments

The authors are grateful to the OPAL infrastructure from Université Côte d'Azur for providing resources and support and to Dr. A. Panfilov and his PhD student T. Nezlobinsky for showing us how to use their software FiniteWave to easily simulate the EP data.


domain: info phys stat

text:
Bou-Laouz Moujahid

Vadaine Rodolphe

Hajduch Guillaume

Fablet Ronan

Automated Ship Detection and Characterization in Sentinel-2 Images: A Comprehensive Approach

Keywords: Deep Neural Network, Sentinel-2 Images, Ship Detection, Ship Characterization

The automatic detection and characterization of ships in optical remote sensing images is a key challenge for maritime surveillance applications. This paper presents an automated system specically designed for ship detection in medium-resolution Sentinel-2 images. The proposed approach relies on a deep learning model trained on a dataset comprising over 6000 annotated Sentinel-2 images. It achieves a detection rate of 93%, with an average of 2.1 to 3.9 false alarms per Sentinel-2 image. Besides the detection task, it also addresses the estimation of ship lengths as well as ship headings. It yields a mean error of 15.36m ± 19.57m for ship lengths, and estimates ship headings with an accuracy of 93%. This contribution signicantly enhances the performance of ship detection and characterization systems in optical remote sensing imagery.

I. INTRODUCTION

Spaceborne remote sensing imagery conveys invaluable information for the monitoring of maritime activities, especially the maritime trafc. This is of critical importance for both surveillance and defense issues [START_REF] Emsa | Emsa outlook 2023[END_REF]. We may cite among others the monitoring of maritime borders, the identication of illegal maritime behaviours, the ght against illegal shing and smuggling activities, search and rescue operations... We can distinguish two main categories of satellite imagery for maritime surveillance topics: Synthetic Aperture Radar (SAR) imagery [START_REF] Iodice | Maritime Surveillance with Synthetic Aperture Radar[END_REF] and optical imagery [START_REF]Satellite Imaging for Maritime Surveillance of the European Seas[END_REF]. Currently, SAR imagery is more widespread due to its applicability both day and night and in all weather conditions (e.g., cloudy conditions). Ship detection in SAR images [START_REF] Crisp | The state-of-the-art in ship detection in synthetic aperture radar imagery[END_REF] relies on relatively simple low-level image processing schemes for ships which possess a metallic structure highly responsive to radar signals. This may lead to detection ambiguities along the seashore, where other metallic structures (buoys, pontoons, etc.) cause high-amplitude patterns in SAR images [START_REF] Zhai | Segmentation-based ship detection in harbor for sar images[END_REF], as well as in the case of sea clutter [START_REF] Stasolla | A comparative study of operational vessel detectors for maritime surveillance using satellite-borne synthetic aperture radar[END_REF]. Additionally, inatable and wooden boats like Zodiacs and plastic boats can hardly be identied in SAR images. By contrast, optical imagery offers a practical alternative for ship detection, increasing not only detection capabilities in complex environments but also providing the ability to detect all types of ships.

In recent years, ship detection from optical imagery has seen increased research effort, with a focus on deep learning methods [START_REF] Li | Ship detection and classication from optical remote sensing images: A survey[END_REF]. Most studies focus on ship detection from very-high-resolution (VHR) (0.3 m resolution) [START_REF] Voinov | Towards automated vessel detection and type recognition from vhr optical satellite images[END_REF] and high-resolution (HR) (0.4 -2 m resolution) [START_REF] Chen | Mssdet: Multi-scale shipdetection framework in optical remote-sensing images and new benchmark[END_REF] imagery. However, these satellites images are typically accessible exclusively through tasking modes, resulting in high acquisition costs and a limited interest for continuous monitoring and surveillance tasks. By contrast, medium-resolution (MR) optical imagery (typically, a 10 m resolution), as deployed on Sentinel-2, delivers freely available remote sensing images for numerous locations on Earth with a revisit time ranging between 5 and 10 days. This seems particularly adapted to maritime surveillance tasks. Yet, only few studies [START_REF] Li | Ship detection and classication from optical remote sensing images: A survey[END_REF] have addressed the automated detection and characterization of ships in MR optical satellite imagery.

Here, we address these challenging issues and present a deep learning approach. We rst collect a representative groundtruthed dataset comprising more than 12000 ship examplars from 6000 Sentinel-2 images. Our multi-task deep learning scheme relies on a Faster R-CNN. Our numerical experiments explore data splitting strategies during the training phase to account for class imbalance. We also asses the impact of the neural backbone of the Faster R-CNN architecture. Overall, we report state-of-the-art performance with a detection rate above 93% including for small ships with a ship length below 20 m. These results support the relevance MR optical imagery for maritime monitoring besides HR optical and SAR imagery.

This paper is organized as follows. In Section II, we provide an overview of related work and analyze their drawbacks. The details of the proposed approach are presented in Section III. Section IV substantiates the relevance of our proposed methods through experiments. Discussion is covered in Section V, and the conclusion is presented in Section VI.

II. RELATED WORK AND MOTIVATION

Recent review papers [START_REF] Li | Ship detection and classication from optical remote sensing images: A survey[END_REF], [START_REF] Kanjir | Vessel detection and classi-cation from spaceborne optical images: A literature survey[END_REF] provide surveys on ship detection and characterization in optical satellite images. They compile over a hundred research articles dating back to 1978, providing a comprehensive overview of the subject. The majority of these studies have investigated ship detection using High-Resolution (HR) and Very High-Resolution (VHR) images, employing deep learning approaches such as Faster R-CNN [START_REF] Zhang | R-cnn-based ship detection from high resolution remote sensing imagery[END_REF], YOLO [START_REF] Tang | Hyolo: A single-shot ship detection approach based on region of interest preselected network[END_REF] and U-net [START_REF] Karki | Ship detection and segmentation using unet[END_REF]. However, our specic challenge is different from working with high or very high-resolution (<5 m) images. In these cases, identifying ships among other objects is relatively straightforward due to their visually distinct features. As a result, the challenge of ship detection in high-resolution images primarily revolves around maximizing the detection rate and closely aligns with the established problem of object detection in computer vision, involving issues such as adjacent and small object detection [START_REF] Zhang | R-cnn-based ship detection from high resolution remote sensing imagery[END_REF].

In medium-resolution imagery (10 m in our case), as illustrated in Figure 1, ships, especially smaller ones of 50-meterlong or below, the visual detection and characterization may be complex: for instance, a large ship can resemble parts of a platform; a small island may have a topology similar to a ship; a very small ship might appear almost identical to a bright spot on rough waters, or even be mistaken for a small cloud. As a result, a key challenge is to maximize the detection rate while minimizing false alarms. This may question how above cited studies for HR images apply to MR optical imagery. Only a limited number of studies have specically explored ship detection in MR imagery. Most of these works tend to concentrate on specic aspects. For example, [START_REF] Heiselberg | A direct and fast methodology for ship recognition in sentinel-2 multispectral imagery[END_REF] delves into ship detection and characterization in favorable conditions. [START_REF] Ciocarlan | Ship detection in sentinel 2 multi-spectral images with self-supervised learning[END_REF] addresses the issue of scarce annotated ship data and introduces a 'self-supervised learning' approach. Additionally, [START_REF] Kanjir | Detecting migrant vessels in the mediterranean sea: Using sentinel-2 images to aid humanitarian actions[END_REF] introduces a method for identifying particular ship shapes associated with migrant activity. Few articles adopt a comprehensive approach [START_REF] Štepec | Automated system for ship detection from medium resolution satellite optical imagery[END_REF], which is our primary area of interest.

This scarcity of studies is, in part, attributable to the absence of freely available ship datasets. In contrast to ship detection at high resolution, which benets from an established reference dataset [START_REF] Kızılkaya | Vhrships: An extensive benchmark dataset for scalable deep learning-based ship detection applications[END_REF], the only publicly available Sentinel-2 ship datasets are [START_REF] Vieilleville | Sentinel-2 dataset for ship detection[END_REF] and [START_REF] Ciocarlan | Ship detection in sentinel 2 multi-spectral images with self-supervised learning[END_REF]. These datasets comprise 31 (2000 ship exemplars) and 16 (1053 ship exemplars) Sentinel-2 images respectively, which may be limited to deploy stateof-the-art learning-based frameworks.

Only [START_REF] Štepec | Automated system for ship detection from medium resolution satellite optical imagery[END_REF] addresses the detection of ships in mediumresolution satellite optical images within a relatively general framework. In this article, the analyzed images are sourced from both the Sentinel-2 mission satellites and the Planet Labs Dove satellite constellation. The images are divided into 800 × 800 pixels patches, which are then used as input for the Faster R-CNN [START_REF] Ren | Faster r-cnn: Towards real-time object detection with region proposal networks[END_REF] object detection model. The annotation process relies on colocating the satellite images with AIS (Automatic Identication System) data [START_REF] Milios | Automatic fusion of satellite imagery and ais data for vessel detection[END_REF]. AIS data comprise ship identiers and locations to groundtruth bounding boxes in the considered dataset. As most small ships are not equipped with AIS systems, this study only addresses large ships (> 100 m). Overall, this study reports a 85% detection rate but does not document the false alarm rate. These results do not appear fully conclusive compared with the performance reported for HR and SAR imagery.

The objective of our work is to create an automated system for detecting and characterizing ships in Sentinel-2 images.

III. PROPOSED APPROACH

We present a multi-stage approach for ship detection and characterization in Sentinel-2 images. In the initial stage, we employ a sliding window technique to systematically cover the image. Each window is of size 100 × 100 pixels, representing 1 km × 1 km, and overlaps with neighboring windows by 25%. This overlap ensures that if a ship extends across two windows, a signicant portion of it remains detectable in at least one of them. A window is considered valid if it contains at least 5% of sea pixels, determined using the land/sea mask. Subsequently, these valid windows are categorized as either "Ship" or "No ship" using a Resnet-type [START_REF] He | Deep residual learning for image recognition[END_REF] classier. For the windows classied as "Ship", we utilize a Faster R-CNN detector to obtain a bounding box around the detected ships (only the coordinates of the ships are necessary). This detector incorporates a dedicated branch to estimate the ship's heading. Figure 2 provides an overview of our ship detection and heading estimation system. When different adjacent patches detect parts of the same ship, we apply the non-maximum suppression (NMS) algorithm. However, in our approach, the suppression is not based on network condence but rather on the size of the bounding box. Our primary objective in this context is to retain the largest bounding box to obtain the most precise ship coordinates. Finally, once the coordinates are identied, we create a 50 × 50 image centered on these coordinates to estimate the ship length using a Resnet-type network (see Figure 3).

The motivation behind proposing a two-stage strategy is as follows. In current state-of-the-art schemes, candidate proposals are typically classied based on their internal features. Consequently, a region containing a very small vessel may exhibit the same features as a region containing a portion of rough sea. The context becomes crucial in this scenario, and the initial classication step is necessary to capture the context within sufciently large patches that encompass the entire context while remaining small enough to avoid missing very small vessels.

In our approach, each individual component undergoes distinct training phases, following a traditional methodology that includes training, validation, and testing. It is strongly recommended to include a substantial number of complete Sentinel-2 images during the model development phase. This approach not only allows for the identication and resolution of potential methodological limitations but also facilitates a comprehensive examination of the test dataset's distribution.

A. Classication phase

Developing a "Ship" and "No ship" category classier involves two main phases. Firstly, we focus on generating datasets and verifying their distribution. Secondly, we work on the model aspect, which includes selecting models and tuning parameters.

1) Considered datasets: We generate our dataset by extracting patches with dimensions of 100 × 100 pixels from Sentinel-2 images. Specically, we utilize only the Red, Green, and Blue (RGB) bands (10 m resolution). This selection aims to leverage the advantages of a pre-trained neural network which is a powerful tool for increasing the accuracy and the robustness of our classier [START_REF] Hendrycks | Using pre-training can improve model robustness and uncertainty[END_REF]. These images are generated using the vessel detection reports provided by CLS analysts. These vessel detection reports are routinely generated by CLS in the context of its commercial activities and contains for each Sentinel-2 satellite image analysed the position of the detected vessels along its characteristics such as length and heading when these are measurable. Using the vessel detection reports available, we process over 6000 Sentinel 2 images, extracting approximately 12000 ship images and 20000 non-ship images. To our knowledge, this represents the rst instance where a signicantly large dataset has been employed for ship detection in medium-resolution optical imagery. The "No ship" category receives particular attention, as it must encompass all the scenarios encountered in optical images, including sea, land, and clouds. Images consisting solely of land are excluded, as their presence does not impact performance due to the availability of a land mask that eliminates these false alarms. Furthermore, it is essential to ensure that our dataset adequately represented objects resembling ships, such as small clouds, small islands, rough seas, and platforms. Manual verication is performed extensively to conrm the presence of all these "sub-categories" within our dataset.

2) Verication of Distribution: The partitioning of training, validation, and test sets may appear straightforward in the majority of machine learning problems. However, when dealing with automation applied to a wide range of diverse images under various scenarios, a careful partitioning becomes crucial. Let us assume that our network is trained on images featuring rough seas but is not tested on samples from this specic scenario to assess its "understanding". This can lead to vulnerabilities during the deployment phase: results on the test set may appear very promising, but performance on real-world images may be signicantly poorer. Indeed, our issue involves binary classication. However, given the diversity of scenarios, we consider working with subcategories. In our work, We divide the "Ship" category into 11 sub-categories: Ship 10, Ship 20, Ship 30, Ship 40, Ship 50, Ship 100, Ship 150, Ship 200, Ship 300, Ship 400, and Ship 458. In this terminology, "Ship x" signies a ship with a length less than "x". Figure 4 illustrates the distribution of "Ship" category in our dataset based on their size. We also subdivide the "No ship" category into 6 subcategories: Others compries patches devoid of vessels, randomly sampled from 6000 Sentinel-2 images, Artifacts contains patches with unwanted visual anomalies or distortions that can occur during image acquisition or processing, Cloud contains various cloud formations, Island includes patches featuring islands and underwater reliefs, Land predominantly consists of coastal areas, ports, and a few images with only land terrain, lastly, Swell contains various forms of swell. Figure 5 illustrates the distribution of "No ship" category based on this subdivision. Figure 6 shows some exemples of each sub-category.

To reduce the random variability caused by the train/test split, we further divide each of these subclasses into 5 clusters using an automatic classication method, specically Kmeans [START_REF] Macqueen | Some methods for classication and analysis of multivariate observations[END_REF]. This clustering is based on the histogram values of the image's 3 color channels. Each of these clusters is subsequently divided into training, evaluation, and test sets. This approach ensures that the neural network is trained and tested on a wider range of data, resulting in a model that is not only more accurate but also more robust. Our approach is summarized in Diagram 7.

We also extensively employ data augmentation techniques to ensure a balanced representation of all our sub-cateogies. For the 'Ship' category, we apply vertical and horizontal ips as well as rotations. Regarding the 'No ship' category, we employ the same augmentation techniques in addition to elastic deformation to diversify the shapes of islands and clouds. It's worth noting that conducting a comprehensive enumeration of all the natural sub-categories is of paramount importance. The more comprehensive the enumeration, the more effective the network will be in practical applications. Carrying out this enumeration during the data generation phase can save considerable time and effort, as it reduces the need for frequent ne-tuning of the network.

3) Classication Network: We conduct extensive experiments using various pretrained models on the Imagenet dataset [START_REF] Deng | Imagenet: A large-scale hierarchical image database[END_REF], including InceptionV3 [START_REF] Szegedy | Rethinking the inception architecture for computer vision[END_REF], Resnet18, Resnet34, and Resnet152 [START_REF] He | Deep residual learning for image recognition[END_REF], while ne-tuning different parameters. Our ndings consistently favore the Resnet architecture, which demonstrates superior performance on our dataset.

In our training process, binary cross-entropy serves as our loss function. We utilize the Adam Optimizer with a batch size of 128 and a learning rate set at 10 -4 .

To further enhance our model's ability to reduce false alarms, we implemented an ensemble approach using two classiers: 'Resnet152' and 'Resnet34.' These two models exhibit comparable results and show signicant improvements in false alarm reduction, albeit with a slight decrease in ship detection performance.

B. Detection phase

The core of our system is the classication component, which plays a pivotal role. The main goal of the detection component is to assign bounding boxes to ships that Fig. 7: Our K-means method for forming the training, validation, and test sets. our classier identies within patches labeled as 'Ship'. Additionally, it may eliminate some false alarms if no ship is detected in the patch by the detector. To build our detection network, we employe the Faster R-CNN framework [START_REF] Ren | Faster r-cnn: Towards real-time object detection with region proposal networks[END_REF]. Despite its introduction in 2015, Faster R-CNN, along with adapted versions, continues to demonstrate its state-of-the-art performance on established reference datasets [START_REF] He | Rethinking imagenet pre-training[END_REF], [START_REF] Mahendrakar | Performance study of yolov5 and faster r-cnn for autonomous navigation around non-cooperative targets[END_REF]. However, it's important to note that we do not explore other detector types in this work, including YOLO and its variants [START_REF] Redmon | You only look once: Unied, real-time object detection[END_REF] and detection transformer (DETR) and its variants [START_REF] Carion | End-to-end object detection with transformers[END_REF].

We utilize a modied version of Faster R-CNN for detection, based on the code referenced in [START_REF] Trzynadlowski | Faster r-cnn[END_REF], implemented on the TensorFlow [START_REF] Abadi | TensorFlow: Large-scale machine learning on heterogeneous systems[END_REF] framework. The primary motivation for these modications is to enhance the detection of small vessels. The Faster R-CNN architecture can be divided into three key components: the Backbone network, the Region Proposal Network (RPN), and the Detector network. Since we are dealing with a single object type, we can use only the RPN along with the backbone network to detect ships. However, by employing the entire structure, we achieve signicantly improved results. This improvement is attributed to the Detector network, which incorporates Fully Connected layers, rening the outcomes of the RPN, which relies solely on convolutional operations.

Our patches are sized at 100 × 100 pixels and contain sometimes very small vessels, representing just 1 pixel of the image, and sometimes large ships, occupying more than 40 pixels. To accommodate this size range, we congure the Region Proposal Network (RPN) to propose different regions of interest, even for the smallest objects. To achieve this, we utilize anchor boxes of various sizes based on the statistics of the ship length in our dataset, we choose 4 × 4, 7 × 7, 10 × 10, 14 × 14, 18 × 18, and 28 × 28 for square-shaped anchor boxes. Additionally, we maintain the original aspect ratios from the Faster R-CNN paper [START_REF] Ren | Faster r-cnn: Towards real-time object detection with region proposal networks[END_REF], which are 1, 0.5, and 2. Furthermore, the input feature map for the RPN is set to a size of 25 × 25 (with a stride length of 4) to preserve crucial information from very small vessels, especially in challenging scenarios. This conguration allows the RPN to effectively handle objects of different sizes within our patches.

The feature map is constructed from the classication network mentioned earlier. In fact, the classication network demonstrates excellent results in distinguishing between "Ship" and "No ship" images. As a result, we believe that the feature maps generated by this network encompass the essential information required to serve as the backbone of our detector. While the use of the initial 25 × 25 feature map from the classication network aids in proposing smaller regions that may contain small ships, it remains relatively shallow for capturing information about larger vessels and other objects within the image. To address this, we integrate intermediate feature maps of sizes 13 × 13 and 7 × 7. By employing upsampling techniques, we combine these feature maps to create the feature map for region proposals, ensuring a comprehensive representation of potential objects of interest [START_REF] Cao | An improved faster r-cnn for small object detection[END_REF].

The region proposal network (RPN) architecture, after creating the feature maps, remains unchanged from the original Faster R-CNN. It consists solely of a convolutional layer (512 channels) and two additional convolutional layers, one for classication and the other for regression.

During the training of the RPN, an anchor is considered positive only if its Intersection over Union (IOU) with a ground truth box exceeds 0.7. Conversely, it is labeled as negative if its IOU falls below 0.3. Anchors falling within this range are not used for training purposes. To address the class imbalance issue, we employ a weighted loss function, which helps balance the small number of positive anchors per image.

After proposing the regions of interest, we apply ROI pooling to transform these regions into a xed size of 6 × 6, which is determined based on the lengths of the ships in our datasets. Subsequently, these resized regions are forwarded to the detector network. The structure of our modied Faster R-CNN is illustrated within the green rectangle in Figure 2.

Only the RPN part and the Detector part are trained with a learning rate of 10 -4 . The weights of the backbone layers are frozen. We believe that since the classier performs exceptionally well on a large dataset containing both "Ship" and "No ship" images, the network only needs to learn how to express the precise position of the ship, which is what the RPN and Detector are designed to accomplish.

C. Characterization of ships

In this section, we introduce the deep learning component for the estimation of the heading of the ship by adding an additional branch to the Faster R-CNN architecture. Additionally, we use a ResNet type model to estimate the size of the ship.

1) Headings: While it was feasible to train a bounding box with the ship's orientation, as demonstrated in [START_REF] Xie | Oriented r-cnn for object detection[END_REF], we have access to a dataset containing 3000 ships annotated with their heading. This information proves to be more valuable than a simple rotated bounding box for our specic task.

Concerning the training process, we keep all the weights of our previously discussed Faster R-CNN model frozen. Then, we introduce an additional branch after the ROI pooling layer, as illustrated within the sky blue rectangle in Figure 2, to handle ship orientation estimation. In this branch, we utilize both the cosine and sine of the ship's angle as model outputs. This approach is adopted to convey circular information effectively to the network.

2) Length Estimation: Since regions of interest (ROIs) are resized during the process, some size information is inevitably lost. To address this issue, one approach is to recover size information using ship orientation and the scaling factor. However, given the current resolution, attempts to introduce an additional branch for length estimation, similar to what we do for ship heading, do not yield satisfactory results. Therefore, we decided to employ a separate network for this purpose.

We utilize the Resnet50 network, with inputs sized at 50×50 pixels and centered on the ship. While smaller patches can be considered (46 since the biggest ship we can meet is 458), during the deployment of the entire ship detection system, the inputs for the length estimation model are constructed around the center of the ship's bounding box predicted by our detector. This bounding box may be slightly affected by ship wakes or irregular ship shapes. Hence, we opt for this patch size, taking into account potential small errors in the detection phase.

It is essential to note that for ensuring the robustness of our approach, we systematically exclude all images featuring partial cropping and inadequately labeled lengths. However, we do not exclude images of ships partially obscured by clouds if their lengths are still measurable; otherwise, they are excluded. This renement leads to a dataset consisting of 7460 images, distributed randomly across training (80%), validation (10%), and test (10%) sets.

IV. EXPERIMENTS AND RESULTS

In this section, we present the outcomes of our system, along with experiments demonstrating the effectiveness of the methods used in both the classication and detection phases. We commence by introducing the experiments and results of the classication phase, followed by those of the detection phase. Subsequently, we delve into the characterization phase, and nally, we present the results of our complete ship detection system. To assess the performance of our system, we utilize a distinct dataset comprising 60 Sentinel-2 images, herein referred to as the Evaluation set. It is essential to note that this dataset is distinct from the training, validation, and test datasets used for training each component of the system

A. Description of the Evaluation Set

The evaluation set consists of 60 Sentinel-2 images captured from approximately 20 different Earth locations, as indicated by red stars in Figure 8, at various acquisition times. These images encompass a wide range of scenarios, including coastlines, turbulent seas, artifacts, images with a high density of ships, low-density ship images, and cloud-covered scenes. In total, the 60 images collectively contain 878475 valid patches. A patch is considered valid if it contains at least 5% of sea pixels. In total, the 60 images collectively contain around 1147 ships. The image annotations are carried out by CLS analysts. It's important to note that very small ships may not have been annotated due to the inherent difculty and resolution limitations. Our evaluation set is available in [START_REF] Bou-Laouz | Sentinel-2 dataset for ship detection and characterization[END_REF].

B. Classication Experimentations 1) Evaluating Distribution Verication Effectiveness:

In this experiment, we analyze and compare the results achieved by employing various data splitting methods. We utilize the pre-trained Resnet152 network on the Imagenet dataset, as it exhibits the best performance. We examine three distinct data splitting congurations:

• In the rst approach, we implement a random split, as depicted in the rst phase of the diagram 7. This split is referred to as the "Random split."

• In the second approach, we employ a split based on subclasses, as shown in the second phase of the diagram 7. This split is referred to as the "Subclass split."

• In the third approach, we create a split based on clusters within each subclass, represented in the third phase of the diagram 7. This split is referred to as the "Subclass k-means split." We present the performance metrics for the test sets in Table I, each comprising approximately 6500 images, with 2500 ship images and 4000 no-ship images. The metrics used to assess our classiers include the False Positive Number (FP) and Recall. We opted for FP instead of precision due to the minimal occurrence of false positives. These metrics are compared across different condence thresholds: 0.50, 0.75, and 0.95.

Table I shows comparative results; however, we cannot draw conclusions without testing the methods on a common dataset: the evaluation set (constructed using 60 Sentinel-2 images, as previously mentioned).

Table II This outcome underscores the signicance of employing veri-cation techniques and methodological renements to reduce the random variability inherent in splitting datasets, particularly in complex scenarios with a wide range of image characteristics.

2) Ensemble network classication: Throughout the remainder of our results, we exclusively employ the dataset split derived from the k-means method. To minimize the False Positive Number, we conduct extensive experiments with various models pretrained on Imagenet. Table III Our aim is to create the optimal ensemble, prioritizing FP reduction while tolerating a slight drop in Recall. To achieve this, we employ ensemble models and manually selected aggregation rules, guided by a criterion: False Positives should not include any images clearly devoid of vessels. Consequently, we choose Resnet152 and Resnet54 with an aggregation rule stipulating that both models must output a condence score exceeding 0.95 to classify an image as containing a ship.

The Ensemble network yields a Recall of 93% with only 4 false positives on the test set. Figure 9 illustrates images of these four false positives. It is not evident whether these images contain a ship. Our MRS-Faster R-CNN demonstrates promising performance in comparison to a standard version of Faster R-CNN. Even though Resnet152 and Resnet34 yield superior classication results (Table III), the feature maps constructed from Resnet18 produce the best outcomes. We hypothesize that the comparatively shallower architecture of Resnet18 allows the intermediate feature maps which are of a relatively large size to encompass a richer set of information. It's worth noting that utilizing a backbone constructed from Resnet152 or Resnet34 for detection will help minimize image processing time in our system, albeit with a minor trade-off in performance.

Detection and False Alarm rates are inuenced by two key factors. Firstly, our dataset comprises numerous small vessels with wakes, leading to instances where wakes are considered as part of the ship due to resolution limitations. However, our detector consistently distinguishes between wakes and vessels on larger ships, as depicted in Figure 10. Secondly, our model encounters challenges in detecting closely positioned ships, particularly small vessels, as illustrated in Figure 11. Rening the Intersection over Union (IOU) condition to 0.25, Detection rate reaches up to 95%, and False alarm rate decreases to 5%.

2) Direct application of the Faster R-CNN: In this experiment, we directly apply our detector to the evaluation set, excluding the use of the classication component, following the approch of [START_REF] Štepec | Automated system for ship detection from medium resolution satellite optical imagery[END_REF]. This direct application yields a Detection Rate of 96% and an Average False Alarm number per image of 50. The expected outcome aligns with the inherent design of the Faster R-CNN, which lacks contextual information. To ensure a fair experiment with direct applications, additional enhancements are essential. Incorporating branches to encode contextual information, as demonstrated in [START_REF] Leng | Context learning network for object detection[END_REF], represents a potential avenue for improvement. Notably, such enhancements were not explored in our work.

D. System Performance evaluation

We assess the performance of our ship detection system, emphasizing the detection rate and false alarm count. Figure 12 presents the detection rates for various ship sizes, showcasing an overall detection rate of 93.2%. The system records an average of 2.1 to 3.9 false alarms per Sentinel-2 image. The uncertainty in false alarm numbers stems from the presence of 106 detections resembling small ships (< 20m) but not annotated by CLS analysts due to resolution limitations. Notably, our system excels in detecting large vessels (> 100m), achieving a detection rate exceeding 97%. Table V presents detailed results for each image, wherein false alarms are predominantly attributed to harbor facilities, offshore platforms, and swell (Figure 13). These results demonstrate signicant promise for achieving complete automation in ship detection on Sentinel-2 images and, more broadly, on optical mediumresolution images.

E. Ship Characterization Results

1) Heading Results: Figure 14 displays the heading results for 600 ship images. We assess accuracy in two ways. In the left part of the gure, a predicted heading is considered accurate if the absolute difference between the predicted heading (h p ) and the ground truth heading (h gt ) is less than a specied heading error (h error ), with h error taking values of 15 • , 30 • , and 45 • . In the right part of the gure, accuracy is determined by verifying if the absolute difference between the predicted heading and the actual heading modulo 180 • is less than h error . This relaxation of the rst condition helps accurately assess the model's capability, especially in situations where the direction is unclear due to resolution limitations, as depicted in Figure 15. The accuracy for predicting headings within a range of ±30 degrees is 93%. 2) Length estimation Results: The results on length estimation are based on 746 ship images with dimensions of 50 × 50. The mean length error is reported as 15.36m ± 19.57m, accompanied by an R-squared value of 92 (Figure 16). The relatively high standard deviation error is attributed to instances where wake shapes or partial cloud coverage make the length measurement task challenging (Figure 17). By excluding these challenging instances, the dataset now consists of 695 images, yielding a reduced mean length error of 11.29m ± 10.17m. Despite these challenges, our results demonstrate commendable performance, particularly considering the resolution limitations of Sentinel-2 data.

V. DISCUSSION

The detection of ships in medium resolution optical images remains an underexplored eld, with only a few studies adopting a comprehensive approach to this task. The scarcity of research in this area is not indicative of its insignicance but rather stems from the limited availability of data and annotated resources. In an attempt to bridge this gap, we The rst number denotes the ground truth length, while the predicted length is enclosed in parentheses make available 60 annotated Sentinel-2 images (1147 ship exemplars) which include length, and heading information when theses are mesurable [START_REF] Bou-Laouz | Sentinel-2 dataset for ship detection and characterization[END_REF].

Our study demonstrates a signicant improvement in performance compared to existing results: our detection rate is 93%, compared to 75% in [START_REF] Ciocarlan | Ship detection in sentinel 2 multi-spectral images with self-supervised learning[END_REF]. Notably, our detection rate on large vessels (> 100 m) is 97% compared to 85% in [START_REF] Štepec | Automated system for ship detection from medium resolution satellite optical imagery[END_REF]. It is crucial to highlight that this improvement is primarily associated with the substantial difference in the dataset sizes rather than the methodology employed. Importantly, it should be noted that reproducing our results does not necessitate an equivalent volume of data. Our dataset contains a substantial amount of redundant information, suggesting that achieving comparable results may be feasible with a more carefully selected dataset.

It's crucial to highlight that our primary research objective centers around establishing an effective ship detection system, with ship heading and length estimation emerging as subsequent considerations. While there is room for improvement, particularly in length estimation (our mean error: 15.36m ± 19.57m on 746 ship images), in [START_REF] Heiselberg | A direct and fast methodology for ship recognition in sentinel-2 multispectral imagery[END_REF], the reported mean error is 15.48m ± 10.94m on 34 ships, only 3 of which contain ship wakes. Our higher standard deviation error arises from the presence of ships partially obscured under clouds and a signicant number of images with ship wakes. Upon removing these challenging images, the mean length error is reduced to 11.29m ± 10.17m.

The results for heading estimation have proven highly satisfactory. We report a 93% accuracy in heading estimation, where correct heading estimation is dened as ground truth heading ±30 degrees. The remaining 7% discrepancy is at-tributed to instances where the heading is visually unclear due to resolution limitations. It is noteworthy that, to the best of our knowledge, we did not nd any studies presenting results in heading estimation specically on medium-resolution (MR) optical imagery.

The demonstrated results hold promise for achieving complete task automation. However, it is essential to note that full conrmation of these ndings awaits deployment in production settings (Figure 18), considering the diversity and occasional complexity of scenes encountered. Our system exhibits a time limitation, requiring approximately 10 minutes to process a complete Sentinel-2 image on an 8-Core CPU. This processing time is attributed to the multiple stages involved, which, with the application of alternative ideas, could potentially be consolidated into a single stage without compromising performance.

One potential avenue for system enhancement involves introducing a context branch to the Faster R-CNN architecture [START_REF] Leng | Context learning network for object detection[END_REF]. This modication holds promise in addressing issues such as false alarms triggered by large vessels, which may generate noise in their surroundings, producing bright points that our detector may misinterpret as small ships. To overcome the time limitations inherent in our system, the utilization of a one-stage detector, such as YOLO [START_REF] Terven | A Comprehensive Review of YOLO: From YOLOv1 and Beyond[END_REF], could be considered. Additionally, we acknowledge the challenge of accurately detecting very close ships, reecting a broader issue in general object detection [START_REF] Kaur | A comprehensive review of object detection with deep learning[END_REF]. Future research efforts may explore innovative strategies to mitigate these challenges and further rene the performance of ship detection in optical images.

VI. CONCLUSION

In this research, we introduced an automated system designed to detect and characterize ships in Sentinel-2 images through a multi-stage process. Our system comprises distinct components focusing on classication, detection, and characterization. Within the classication component, we highlighted the importance of validating distribution in complex environments and diverse scenarios. State-of-the-art models, specically Resnet152 and Resnet34, were employed to minimize false positives. During the detection phase, modications were made to the original Faster R-CNN to effectively handle our specic dataset. In the characterization phase, a dedicated branch was integrated into the Faster R-CNN to estimate ship headings, and Resnet50 was separately used for length estimation. The results of our system are promising, indicating the potential for complete automation of ship detection and characterization in Sentinel-2 images. Additionally, the presented approach is adaptable to other satellite constellations. However, our system encounters challenges in differentiating small ships from their wakes and struggles with detecting closely positioned ships. Acknowledging these limitations, further research could explore these specic points, ensuring a continuous improvement in automated ship detection and characterization methodologies.

VII. ACKNOWLEDGEMENTS

The dataset employed in this paper was provided by Collecte Localisation Satellites (CLS). Our sincere appreciation goes to Vincent Kerbaol for his pivotal support, making this project possible. We express our gratitude to the CLS analysts for their meticulous annotation of the datasets.

Special recognition is extended to Etienne Gauchet, Mohammed Benayad, and Maxime Vandevoorde for their contributions to data cleaning and generation. We offer a heartfelt thanks to Chaïmae Sriti and Aurélien Colin for their comprehensive review of the paper and the insightful comments they provided.

Fig. 1 :

 1 Fig. 1: Illustration of complex ship images, presenting challenges in visual detection and characterization.

Fig. 2 :

 2 Fig. 2: An Overview of Our Ship Detection and Heading Estimation System on Sentinel-2 Images. The Faster R-CNN model comprises a backbone constructed from the Resnet18 classier trained on Sentinel-2 data. The layer names in black used to build detection feature map are derived from keras applications [23].

Fig. 3 :

 3 Fig. 3: Ship Length Estimation Component.

Fig. 4 :

 4 Fig. 4: Distribution of ship sub-categories.

Fig. 5 :

 5 Fig. 5: Distribution of No ship sub-categories.

Fig. 6 :

 6 Fig. 6: Examples of No ship sub-categories images.

Fig. 8 :

 8 Fig. 8: Location of the 60 Sentinel 2 Images Constituting the Evaluation Set.

Fig. 9 :

 9 Fig. 9: The only four false positives produced by the Resnet152 and Resnet34 ensemble.

Fig. 10 :

 10 Fig. 10: Model detections exemples : Green bounding boxes denote annotations, while yellow indicates model predictions.

Fig. 11 :

 11 Fig. 11: Model struggles with close ship differentiation, generating a single bounding box for closely positioned vessels.

Fig. 12 :

 12 Fig. 12: System Performance evaluation : Overall Detection rate : 93.2%, Potential False Alarms (Range): 128-234, Average False Alarms for One Image (Range): 2.1-3.9.

Fig. 13 :

 13 Fig. 13: Examples of False Alarms.

Fig. 14 :

 14 Fig. 14: Heading Prediction Results for 600 Ship Images. Left: A correct prediction is indicated when |h p -h gt | ≤ h error . Right: A correct prediction occurs when |h p -h gt | mod 180 ≤ h error .

Fig. 15 :

 15 Fig. 15: Examples of Ambiguous Heading.

Fig. 16 :

 16 Fig. 16: Predicted vs. True ship lengths with an R-squared of 0.92.

Fig. 17 :

 17 Fig. 17: Examples of Ship Length Prediction. The rst number denotes the ground truth length, while the predicted length is enclosed in parentheses

Fig. 18 :

 18 Fig. 18: Integration of the Ship Detection System in CLS Processing Pipeline -Initial System Trials.

  displays the results of the three distinct models, evaluated on the Evaluation set. The FP metric represents the aggregate count of false positives observed across all 60 evaluation images, encompassing approximately 877470

	Split / Threshold		0.5		0.75		0.95
		FP	Rec.(%)	FP	Rec.(%)	FP	Rec.(%)
	Random	45	98.1	34	97.7	15	96.0
	Subclass	23	98.3	12	97.9	5	96.0
	Subclass k-means	25	98.1	16	97.5	8	95.9

TABLE I :

 I Evaluation of Data Splitting Methods on Their Respective Test Sets Using Resnet152.

	Split / Threshold		0.5		0.75		0.95
		FP	Rec.(%)	FP	Rec.(%)	FP	Rec.(%)
	Random	2529	98.2	1818	97.6	1003	96.1
	Subclass	1050	98.4	840	97.9	605	96.2
	Subclass k-means	546	98.2	406	97.1	240	95.5

TABLE II :

 II Evaluation of Data Splitting Methods on The Evaluation Set Using Resnet152. valid patches devoid of ships. Remarkably, the 'k-means split' method outperforms both the 'random split' and 'subclass split' approaches, resulting in the lowest occurrences of false positives, averaging approximately 4 per Sentinel-2 image.

  summarizes the outcomes of these models in terms of False Positive Number (FP) and Recall metrics.

	Model / Threshold		0.5		0.75		0.95
		FP	Rec.(%)	FP	Rec.(%)	FP	Rec.(%)
	InceptionV3	36	97.6	25	97.0	11	95.0
	Resnet18	40	97.0	28	96.9	13	94.0
	Resnet34	33	97.5	24	97.2	9	94.8
	Resnet152	25	98.1	16	97.5	8	95.9
	Ens. Resnet34/Resnet152	14	96.5	10	95.3	4	93.0

TABLE III :

 III Evaluation of Various Pretrained Models' Performance on the Test Dataset.

TABLE IV :

 IV Performance Metrics for Different detection Models.

TABLE V :

 V Detection System's performance on the Evaluation Set[START_REF] Bou-Laouz | Sentinel-2 dataset for ship detection and characterization[END_REF] with Detection Rate and False Alarm Count metrics.