Risk-Aware Wasserstein Distributionally Robust Control of Vessels in Natural Waterways

In the realm of maritime transportation, autonomous vessel navigation in natural inland waterways faces persistent challenges due to unpredictable natural factors. Existing scheduling algorithms fall short in handling these uncertainties, compromising both safety and efficiency. Moreover, these algorithms are primarily designed for nonautonomous vessels, leading to labor-intensive operations vulnerable to human error. To address these issues, this study proposes a risk-aware motion control approach for vessels that accounts for the dynamic and uncertain nature of tide islands in a distributionally robust manner. Specifically, a model predictive control (MPC) method is employed to follow the reference trajectory in the time–space map while incorporating a risk constraint to prevent grounding accidents. To address uncertainties in tide islands, a novel modeling technique represents them as stochastic polytopes. Additionally, potential inaccuracies in waterway depth are addressed through a risk constraint that considers the worst case uncertainty distribution within a Wasserstein ambiguity set around the empirical distribution. Using sensor data collected in the Guadalquivir River, we empirically demonstrate the performance of the proposed method through simulations on a vessel. As a result, the vessel successfully navigates the waterway while avoiding grounding accidents, even with a limited dataset of observations. This stands in contrast to existing nonrobust controllers, highlighting the robustness and practical applicability of the proposed approach.


Introduction 1.Challenges in Navigating Tidal Inland Waterways
In today's interconnected and globalized society, maritime transportation and international logistics are vital for global economic development [1].In the case of inland ports connected to the sea by artificial channels or natural waterways, traffic management is the crucial factor for port efficiency.This issue has gained significant attention due to the substantial impact of transit channels on the global economy [2], concerns about the potential adverse impact of maritime transportation on climate change [3], and strong interest of maritime authorities in risk management processes to prevent shipping accidents [4].As a striking example, in this article, we consider the Guadalquivir River in the South of Spain, which is characterized by a changing depth.In this river, the spatial and temporal variation of the depth is mainly due to two natural phenomena: the irregular bathymetric profile of the waterway [5] and the effect of the tide [6].Due to its location in a highly protected area, the Doñana National Park, dredging operations [7] are not allowed, making the problem even more challenging.Additionally, due to the risk of grounding accidents caused by the varying depth of the river, expert river pilots are required to control the vessels in the waterway, and thus, navigation in this channel is prone to human errors.
In order to prevent grounding accidents, it is essential to create a map of time-space points, known as tide islands, which takes into account the depth at each location of the waterway and the draft of each vessel at every time instant.This map is produced by overlaying the depth map with a fixed cutting plane based on the vessel's draft.An example of such a map for a 7-meter draft vessel is displayed in Fig. 1, which was generated using interpolation techniques on data collected by a series of depth sensors located along the river.It is worth noting that the size of the tide island varies over time for a given location in the waterway.Additionally, different sources of uncertainty, such as measurement noise or the effect of the moon, make it difficult to determine the exact shape of the island and the depth of the waterway at a specific location, thereby increasing the risk of grounding accidents.Some works have paid particular attention to solving the scheduling problem in waterways affected by the tide, which complicates the problem due to the time-dependent depth of the waterway [8].However, despite all the efforts made to schedule and plan vessel trips in inland waterways influenced by tides, autonomous motion control of vessels in these channels remains a challenge, and expert river pilots are usually required to execute the optimal trip plans.Moreover, obtaining accurate dynamic models of time-varying depth for a specific waterway is a difficult task, if not impossible, due to their complexity.Meanwhile, interpolated models such as the one in Fig. 1 obtained from sensor readings are often used, which may not be accurate enough, jeopardizing security and increasing the risk of grounding accidents [9].

Our Contributions
In the proposed framework, the vessel control problem is viewed as a trajectory-tracking problem, where the optimal reference trajectory is provided by a planner (e.g., [10]).Due to the dynamic and uncertain nature of tide islands, it is crucial to design a motion controller that efficiently follows the reference trajectory while avoiding grounding accidents, especially when operating with a limited dataset of sensor readings.
Since the main source of uncertainty is the depth of the waterway at each point and time instance, we first propose a novel formulation to model tide islands as stochastic uni-dimensional polytopes in the time-space map where the constraint imposed by the depth of the waterway is violated.The only information available about these uncertainties is the limited set of data obtained from noisy measurements taken by depth sensors placed along the waterway.Using this dataset, an interpolated map of the time-space regions of depth constraint violation is obtained.Based on this information, a time-varying empirical distribution of the width of each tide island is derived.
Then, we design a motion controller that assesses the risk of accidents via the conditional valueat-risk (CVaR), which measures the average worst-case loss.CVaR is a coherent risk measure in the sense of Artzner [11] that can distinguish tail events and is considered an adequate risk measure for robotics applications [12].However, in order to properly assess the risk of accidents, it is required to have accurate information about the probability distribution of the underlying uncertainties.Given that we have only a limited dataset of sensor readings, we may end up with undesired behaviors and a lack of safety.To address this issue, we utilize the distributionally robust risk-aware model predictive control (DR-MPC) method introduced in [13], which is robust against distributional inaccuracies within a specific set of distributions called the ambiguity set and can result in high-performance and safe control actions even with inaccurate information.To the best of our knowledge, this marks the pioneering application of distributionally robust optimization (DRO) and distributionally robust control (DRC) in the field of maritime transportation.
In this work, the ambiguity set is designed as a statistical ball around the given nominal distribution, where the distance between distributions is measured via the Wasserstein distance.Unfortunately, the proposed DR-MPC problem is intractable as it constitutes an infinite-dimensional problem.To overcome this issue, a Kantorovich duality-based strategy [14] is employed to reformulate the original problem into a finite-dimensional nonconvex optimization problem that can be solved using a standard commercial solver.
Finally, using depth sensor data collected in the Guadalquivir River, we empirically demonstrate the performance of our risk-aware Wasserstein DR-MPC approach that solves the motion control problem of vessels in the presence of inaccuracies in the empirical distribution.Our findings indicate that the approach successfully controls the vessel along the waterway without any collisions, even with a limited dataset of only six observations, while the standard non-robust, risk-aware model predictive controllers (MPC) result in collision rates of at least 63%. 1he remainder of this paper is organized as follows.In Section 2, we detail the preliminaries needed for solving the problem.Section 3 presents the adapted version of the Wasserstein DR-MPC approach, explaining how tide islands are modeled as stochastic uni-dimensional polytopes and how the DR-MPC problem is formulated.In Section 4, we apply and test the proposed motion control strategy in a real-world scenario in the Guadalquivir River, comparing the results to those obtained when using the empirical distribution directly.Finally, in Section 5, we summarize the main conclusions of our work.

System and Obstacles Modeling
In this work, we treat each vessel as a water surface robot whose dynamics can be modeled by a discrete-time system where x(t) ∈ R nx is the state of the system (e.g., location, velocity), and u(t) ∈ R nu is the control input (e.g., thrust, acceleration).The system output y(t) ∈ R ny represents the position of the vessel in a n y -dimensional space.The functions f : R nx × R nu → R nx and h : R nx × R nu → R ny are general nonlinear functions representing the system dynamics and the output mapping, respectively.Perfect knowledge of the state and noiseless output measuring is assumed.Moreover, the system state and the control input are constrained as where X ⊆ R nx and U ⊆ R nu are convex compact sets, designed based on the physical constraints.Vessels must safely sail through the waterway while avoiding obstacles, which are the regions in the time-space map where the depth constraint is violated.Due to the nonconvex nature of these regions, each of them is approximated by the ball of minimum radius containing all points in the region, as it is depicted in Fig. 2, where three fictitious tide islands have been represented2 .Here, the dark blue areas represent the tide islands, while ω l,t,k represents the distance from the extreme of the ball to the center of each obstacle l at time t + k.
Once the convex region containing each tide island is found, each of the obstacles to be avoided can be seen as a one-dimensional ball whose center is always located in the same position of the waterway but whose radii periodically change in time due to the effect of the tide.In a general formulation of the collision avoidance problem, we consider a total of L size-changing obstacles (tide islands) to be avoided by vessels while sailing in the waterway.The region occupied by each where O c l is the center of the ball, ω l,t,k its radius at time In this work, we consider ω l,t,k to be a random variable whose true probability distribution is unknown.The obstacles' stochastic nature arises from various factors, including measurement noise from depth sensors, information loss during depth map interpolation, and the impact of unquantifiable natural phenomena on waterway depth (e.g., the effect of suspended sediments on the density of the river water and the tidal amplification effect [15]).
The procedure of obtaining an empirical distribution for ω l,t,k is detailed in Section 2.2.Note that these obstacles are given by the intersection of the ball containing all points where the depth constraint is violated with a plane of constant value given by the time instance considered.This is shown in Fig. 2, where the red lines represent the obstacles to be avoided at some particular time t * .
Having the definition of each l-th obstacle, we define its safe region, Y l (t+k), as the complement of the region occupied by the obstacle, i.e., where For collision avoidance, vessels should always sail through the intersection of the safe regions associated with all obstacles, i.e., in order to guarantee safe navigation.

Empirical Distribution Generation
In practice, the true distribution µ l,t,k of the uncertain parameter ω l,t,k , representing the width of each obstacle l at time t + k, is unknown, and the only information available about its distribution is a limited set of historical data.In this work, we assume that the tidal wave is a perfect sine of constant period.Therefore, we assume that the depth at time instances separated by intervals given by the tide period are samples drawn from the same probability distribution.Assumption 1.Let the depth at some location of the waterway be a periodic wave [16].Let the period of this wave be given by a constant value, which is reasonable for the considered time period, given the small variation caused by the gravitational effect of the moon.Also, let T be the period of the tidal wave.Then the following expression holds: Under Assumption 1, a series of sets of all possible values of the width of each obstacle l at each time instance t has been obtained using the map of depth depicted in Fig. 1.Each of these datasets is a collection of different observations of ω l,t,iT for all i = 1, ..., N , which is represented in Fig. 3.Then, these datasets are employed to obtain an empirical distribution of ω l,t,k for each obstacle and for each time instance.
To obtain an estimation of the real probability distribution, the only information available is a limited set of N possibly noisy data samples, D l,t,k = {ω 1 l,t,k , ..., ωN l,t,k }.A possible option would be to obtain an empirical distribution defined as , where δ ω is the Dirac delta function concentrated at ω.This way of obtaining an empirical distribution is, however, not capable of providing a reliable estimation of the safety risk, especially when the size of the sample data set N is small.

Wassertein Distributionally Robust Motion Control
Taking into account the changing depth and all operational and regulatory constraints for preventing collisions at sea [17], a reference trajectory can be obtained, for example, following the novel approach recently presented in [10].Once the optimal reference trajectory for a vessel is determined, we adopt the Wasserstein DR-MPC control approach [13] and apply it to track this trajectory while guiding the vessel to its final destination safely.This approach limits the risk of grounding accidents caused by the waterway's shifting depths, which are treated as uncertainties subject to an unknown probability distribution. 3

Measuring the Risk of Unsafety Using CVaR
The risk-aware control approach proposed in this work employs the notion of safety risk introduced in [18].In particular, the loss of safety regarding obstacle l be the distance between the vessel's location at time t + l and its associated safe region, i.e., the vessel's penetration into the obstacle measured as dist(y where ∥•∥ 2 is the standard Euclidean norm.The objective is to drive the vessel so that the loss of safety is zero.Otherwise, if this quantity takes a positive value, the vessel would penetrate into the obstacle, and a grounding accident would occur.However, the resulting actions will be overly conservative if a deterministic constraint is imposed on the penetration distance.As an alternative, we propose quantifying the unsafety associated with the obstacle l using a risk measure. In particular, we measure the risk of unsafety by computing the CVaR of the loss of safety defined in (2).The CVaR of a random loss X following a distribution γ is defined as the conditional expectation of the loss within the (1 − α) worst-case quantile.As defined in [19], it can be expressed as where (x) + = max{x, 0}.Thus, the safety risk defined by CVaR measures the conditional expectation of the distance between the vessel's position y(t + k) and the safe region Y l (t + k) within the (1 − α) worst-case quantile of the safety loss distribution.Then, the risk constraint imposed on each vessel can be written as: where δ l ≥ 0 is a user-defined parameter that can be employed to adjust the risk tolerance of the vessel.It can be chosen to adjust the risk tolerance of the vessel, usually chosen to obtain the desired trade-off between robustness and optimal performance in terms of sailing times.

Motion Control via DR-MPC
The limitation imposed by the lack of information may lead to unsafe decision-making without respecting the original risk constraint, which could lead to a grounding accident as the safety risk is poorly assessed.To address this issue, the approach proposed in [13] is adopted, where the safety risk is evaluated under the worst-case distribution of ω l,t,k belonging to a given set of distributions D l,t,k called the ambiguity set, leading to the following distributionally robust risk constraint: By upper-bounding the worst-case risk value that the vessel can bear, we ensure that the resulting control action is robust against distribution errors.The construction of the ambiguity set is a crucial factor that impacts the trade-off between robustness and performance.Similar to [13], we model the ambiguity set as a statistical ball, which includes all probability measures whose distance from the empirical one is no greater than the given radius θ > 0, i.e., D l,t,k := {µ ∈ P(W) | W p (µ, ν l,t,k ) ≤ θ}, where P(W) denotes the set of Borel probability measures on the support W ⊆ R. The distance between two measures is computed via the Wasserstein metric of order p [20], which measures the minimum amount of work for transporting a certain mass from one distribution to another through non-uniform perturbations.
Definition 1.The p-order Wasserstein metric between two probability measures µ ∈ P(W) and ν ∈ P(W) is defined as where κ is the transportation plan, with Π i κ being its ith marginal, with ∥•∥ representing an arbitrary norm on R n .
Though the Wasserstein metric has often been viewed in the context of the optimal transport problem [21], recent works in DRO and DRC have demonstrated the application and utility of Wasserstein ambiguity sets for promoting robustness in the distributional sense while providing substantial performance guarantees.Therefore, in this work, we aim to develop a control framework that is robust against distributional errors characterized by the Wasserstein ambiguity set.
Remark 1.Other types of ambiguity sets can be employed in the DR-MPC formulation.A popular approach in the literature of DRO is the use of moment-based ambiguity sets [22], which are often overly conservative, requiring a large amount of data to reliably estimate moment information.Another approach involves the use of statistical-distance-based ambiguity sets employing ϕ-divergence [23] and total variation distance [24].In contrast to these methods, Wasserstein ambiguity sets [25] contain a richer set of relevant probability distributions.The associated DRO problem has been proven to offer superior finite-sample performance guarantees and address the closeness between two points in the support [14], especially useful in sequential decision-making problems [26].Moreover, as demonstrated in [27], Wasserstein DRO can alleviate overfitting due to its strong relation to regularization techniques used in machine learning algorithms.All these features make ambiguity sets based on Wasserstein distance an attractive approach for the proposed motion control tool.
Taking this into account, the discrete-time Wasserstein DR-MPC problem for trajectory tracking can be formulated as inf u,x,y x k ∈ X (4e) where all constraints hold for k = 1, ..., K and l = 1, ..., L, except for (4b) and (4f) which hold for k = 0, ..., K − 1.Here, K is the prediction horizon and u := (u 0 , ..., u K−1 ), x := (x 0 , ..., x K ), and y := (y 0 , ..., y K ), are the control input, state, and output sequences over the time horizon considered, respectively.The equality constraints (4b) and (4c) represent the system dynamics and the measurement of the predicted output, respectively.Polytopic constraints on state and control input are imposed via (4e) and (4f), while the initial state is set via (4d).Finally, (4g) bounds the distributionally robust risk at a given tolerance level.The cost function (4a) is chosen to penalize the deviation from the reference trajectory τ k with weighting matrices Q ⪰ 0 and P ⪰ 0, as well as penalize large deviations of the control action ∆u k := u k − u k−1 , which is included with a weighting matrix R ≻ 0. The Wasserstein DR-MPC is solved in a receding horizon fashion, where only the first component u * 0 of the optimal solution u * is applied to the system at time t, i.e., u(t) = u * 0 .However, solving this DR-MPC control problem is a challenging task since it involves an infinitedimensional optimization problem over the space of distributions due to the risk constraint (4g).To alleviate the intractability, we perform reformulation, detailed in the following subsection.
Remark 2. In the proposed DR-MPC problem formulation, we consider uncertainty in the obstacles' dimensions.This differs from [13] in that the latter considers uncertainty in the obstacle's motion.The proposed formulation is also different from the one studied in [28], where authors consider uncertainties in the system.Several DRC problems solved via dynamic programming approaches can be found in the literature, as is the case of the problems studied in [29][30][31], among others.

Finite-dimensional Reformulation via the Kantorovich Duality
In this subsection, we reformulate the DR-MPC problem to tackle the issue of infinite dimensionality and develop a tractable algorithm for solving it.For ease of explanation, the obstacle index and time subscripts will be omitted.First, we derive a simple expression for the loss of safety (2).Given that the region occupied by an obstacle O is a convex polytope, the corresponding safe region can be represented by the union of half-spaces as Then, we use the following lemma to reformulate the loss of safety using the half-space definition of the safe region.( The proof of this lemma can be found in Appendix A.1. Next, we find a conservative approximation of the distributionally robust risk constraint given by (3).Specifically, as shown in [13], the distributionally robust safety risk can be upper-bounded by This new risk formulation is still difficult to evaluate since it involves a maximization problem over a continuous set of distributions.Fortunately, when the distance between two distributions is measured via the Wasserstein distance of order one, we can reformulate the risk constraint using the Kantorovich duality [14], as proposed in [13], to transform it into a finite-dimensional optimization problem as follows.
Proposition 1.Consider the Wasserstein distance of order one.Let the uncertainty set be the entire space W = R and ∥ • ∥ represent the Euclidean norm.Then, the following equality holds: where all constraints hold for i = 1, ..., N .Here, e m ∈ R m is a vector whose elements are all equal to one.
The proof of this proposition can be found in Appendix A.2.It is important to note that although Proposition 1 is rooted in [13, Proposition 1], the key difference lies in the decoupling of disturbances and the obstacle's position due to the way disturbances influence the system.This decoupling allows for constraint reduction when reformulating the DR-MPC problem.Specifically, the original DR-MPC problem (4) can be reformulated as inf u,x,y, z,λ,s,ρ where all constraints hold for k = 1, ..., K, l = 1, ..., L and i = 1, ..., N , except for constraints (7m) and (7l), which hold for k = 0, ..., K − 1 and k = 0, ..., K, respectively.The reformulated above problem is a finite-dimensional problem.However, it is still nonconvex due to the bilinearity of constraint (7f).A locally optimal solution could be obtained by means of nonlinear programming algorithms such as interior point methods or sequential quadratic programming [32].If the model of the system dynamics and output equations are affine, relaxation techniques can be employed to find a globally optimal solution [33].
A notable advantage of the Wasserstein DR-MPC approach ( 7) is that for a carefully designed radius θ, it assures a probabilistic finite-sample performance guarantee, meaning that the safety risk constraint is satisfied with a probability no less than a certain threshold, even when evaluated under a set of new samples chosen independently of the training data [13].
Remark 3. The deviation of the vessel from its reference trajectory could be such that maneuvers among vessels may not occur in the expected section of the waterway according to the optimal trip plan.To solve this, a multi-agent control approach could be adopted.The proposed method can be used for centralized control of multiple vessels by considering a single (possibly large) MPC problem considering all vessels and their associated operational constraints.This scenario is out of the scope of this work, as we aim to propose a general formulation of the DR-MPC problem for a single vessel in natural waterways.

Case Study
In this section, we carry out a series of simulations in which the proposed control approach is applied to a real case scenario in the Guadalquivir River.We use the DR-MPC method as a tracking controller for the vessel of interest with its corresponding reference trajectory.For training purposes, it is assumed that only a set of N = 11 possibly noisy values of past records of ω l,t,k is available.A series of experiments are carried out to analyze the performance and robustness of the proposed approach, measuring the cost of the DR-MPC problem and the minimum distance between the location of the vessel y(t) and the contour of the safe region Y.In the first one, for a given value of δ l , the risk tolerance level, we analyze the performance and robustness of the system as a function of θ, the radius of the Wasserstein ball defining the ambiguity set centered at the empirical distribution.Then, for a given value of θ, we analyze the performance and safety of the vessel as a function of the number of samples N employed to generate the empirical distribution.
All the results are compared to two baselines: a standard chance-constrained MPC (CC-MPC) [34] and a risk-aware MPC (SAA-MPC) [35], where in both cases the obstacle avoidance constraint solely relies on the empirical distribution.We maintained the same probability level α as used in DR-MPC for the baselines.All the experiments have been solved in AMD RYZEN 4000 CPU @ 3GHz using Gurobi solver, which incorporates a spatial branch and bound algorithm for solving bilinear programs [36].

Setup
The Guadalquivir River is divided into a total of 24 sections, representing the number of boundaries.Each of these sections p is characterized by its length, minimum width, and the maximum speed at which vessels can sail through it.These parameter values can be found in [10].The maximum time we allow each vessel to spend waiting and sailing through the waterway is set to 12.In the scheduling problem, a total of four vessels are considered.The maximum speed at which these vessels can sail is set to v i = 25 km/h, and the beam of each vessel is randomly selected.The draft of all vessels is set to 7 m.Thus, the depth map depicted in Fig. 1 is employed for all vessels, which shows that L = 5 tide islands must be avoided by each vessel.The planning objective is then to minimize the time vessels spend sailing through the waterway and waiting to access it.We assume that the times when vessels arrive at the waterway are sampled from a Poisson distribution.

Vessel Motion Control
Given the optimal reference, we now demonstrate the utility of the proposed DRC method.The objective is to drive a single vessel from one extreme of the waterway to the other (i.e., track the reference trajectory) while avoiding the obstacles representing the tide islands.Because the analysis of vessel maneuverability is not the objective of this work, the model of the vessel is given by the linear system where x k ∈ R and u k ∈ R are the position and speed of the vessel at time instance k, respectively, and T s is the sampling time, which in this case is set to 1 h/10 (6 minutes).Note that in this case, the measured state is the output of the system.We set the MPC prediction horizon to K = 10, while the weighting factors of the cost function are set to Q = P = 1 and R = 0.1.

Effect of Wasserstein Radius
In the first experiment, we analyze the effect of the Wasserstein radius θ on the control performance and safety of the DR-MPC method.A total of 100 simulations are performed for different values of θ.We generate the empirical distribution assuming that, at each time step, only the N last values of ω l,t,k for each time step t are available.These values are randomly selected from a validation data set of 1200 values of ω l,t,k .
Examples of the different trajectories of the vessel controlled by the baselines and our DR-MPC with parameters δ l = 0.02, α = 0.95, and N = 11 for a single uncertainty realization are shown in Fig. 4. Note that, in this case, only the first and last obstacles affect the navigation due to their distances to the reference trajectory.For all the scenarios generated, the average value of the total costs, the minimum distance between the vessel and any of the obstacles (safety margin), as well as the collision rates are reported in Table 1.Notably, larger values of θ correspond to higher costs for DR-MPC and larger safety margins.In contrast, both SAA-MPC and CC-MPC prioritize cost over safety, leading to negative safety margins and collisions occurring in 63% and 89% of cases, respectively.Remarkably, all trajectories generated by DR-MPC with θ ≥ 0.0015 are collision-free.Even for DR-MPC with θ = 0.001, the collision rate is lower than that for the baselines, underscoring the necessity of adopting the distributionally robust approach to ensure collision-free navigation.Moreover, our DR-MPC method achieves a positive average safety margin for θ ≥ 0.00125, indicating no collisions with obstacles.In contrast, both SAA-MPC and CC-MPC consistently show negative safety margins, signifying vessel penetration of obstacles.

Effect of Sample Size
In the second experiment, we examine how, for some fixed value of the Wasserstein radius, the control performance and safety are affected by the amount of data samples.It is assumed that only the last N samples of ω l,t,k for each time step t are available.These values are randomly selected from a validation data set of 1200 values of ω l,t,k for each time step t.The results of adopting the DR-MPC approach are compared to those obtained when the baselines are applied for 100 randomly generated scenarios using the same configuration as in the first experiment.Assuming the best trade-off between cost and safety is achieved at θ = 0.0015, we present the average values of total cost, safety margin, and collision rate for DR-MPC, SAA-MPC, and CC-MPC in Table 2.These results were obtained for different values of N , with δ l = 0.02 and α = 0.95.As anticipated, trajectories generated by SAA-MPC and CC-MPC consistently result in collisions, regardless of the number of data samples considered.Moreover, the collision rate for vessels controlled by both SAA-MPC and CC-MPC rises as the number of samples decreases while the safety margin remains negative.In contrast, our DR-MPC method ensures collision-free navigation with a larger safety margin, even with a small dataset of N = 6 samples.

Conclusions
In this work, a DRC approach for vessels in natural inland waterways, paying particular attention to the case of the Guadalquivir River, where navigation is critically affected by the effect of the tide.A risk-aware DR-MPC method was proposed for sailing through the waterway while avoiding stochastic grounding accidents.For that purpose, we modeled the regions of depth constraint violation as balls whose radii belong to an unknown probability distribution.We employ the dataset of past measurements of the tide to build the empirical distribution of the newly modeled obstacles' radii.To hedge against inaccuracies in the resulting distributional information, the DR-MPC method limits the safety risk by allowing perturbation within a Wasserstein ambiguity set.By doing so, the proposed approach resolves the issue related to the inexact empirical distribution obtained from the limited amount of available data.After performing a series of experiments, the results obtained show that even with a very small sample size, the vessel driven by the Wasserstein DR-MPC successfully avoids the obstacles that randomly change in size over time.We also show how, for different values of the Wasserstein radius and the number of data samples, the proposed DR-MPC approach outperforms SAA-MPC and CC-MPC approaches, where only the empirical distribution is considered.The proof is concluded using the fact that by definition |A j | = 1 and taking ω outside the minimum.

A.2 Proof of Proposition 1
Proof.Since the loss function ( 5) is affine in ω, we first apply the results in [14][Corollary 5.1], so that the left-hand side of (6) equals to Finally, it is required to show that the original problem ( 6) is equivalent to the reformulated problem (8).Let J * and J ′ be the optimal values of the original and reformulated problems, respectively.Let (λ * , s * ) and (λ ′ , s ′ , ρ ′ ) denote the corresponding optimal solutions.By the feasibility of the point s * i and by using ( 9), we have that min j=1,2 where ρ * is the optimal solution to the right-hand side in (9).Since all the other constraints are the same in both problems, we conclude that the optimal solution (λ * , s * ) to the original problem is feasible for the reformulation one, i.e., J ′ ≤ J * .To show that J ′ ≥ J * , we first see that ω(i) + min j=1,2 {A j (O c − y)} ≤ ω(i) + (ρ ′ ) ⊤ A(O c − y) where the inequalities follow from the feasibility of ρ ′ .By verifying all the remaining constraints, we conclude that the optimal solution (λ ′ , s ′ , ρ ′ ) to the reformulated problem is feasible for the original one, i.e., J ′ = J * .Thus, the two problems are equivalent.

Figure 1 :
Figure 1: Depth map example for a 7-meter draft vessel.

Figure 2 :
Figure 2: Fictitious time-space map with random size-changing obstacles.

Figure 3 :
Figure 3: Obtaining the uncertainty data set for ω l,t,k at t for different observations of the same obstacle l.

Lemma 1 .
Consider the obstacle O defined in (1) and the associated safe region Y.Then, the loss of safety can be expressed as dist(y, Y) = ω + min j=1,2 {A j (O c − y)} + .

Figure 4 :
Figure 4: Trajectories of the vessels controlled by the SAA-MPC, CC-MPC, and DR-MPC for different values of θ.The dotted black lines around each tide island represent the maximum value of the uncertainty in the observation dataset.

A. 1 2 ω
Proof of Lemma 1 Proof.Applying the steps in the proof of [13, Lemma 1], it can be shown that dist(y, Y) = min j=1,+ A j (O c − y) |A j | + .

Table 1 :
Average total cost, safety margin, and collision rate for the DR-MPC, SAA-MPC, and CC-MPC approaches for different values of θ with δ l = 0.02, α = 0.95, and N = 11 computed over 100 simulations.