# Pearson Coefficient Degradation in a Wasserstein/Gelbrich Ambiguity Set

Matthieu Borelle, Teodoro Alamo, Cristina Stoica, Sylvain Bertrand, Eduardo Camacho

► To cite this version:

Matthieu Borelle, Teodoro Alamo, Cristina Stoica, Sylvain Bertrand, Eduardo Camacho. Pearson Coefficient Degradation in a Wasserstein/Gelbrich Ambiguity Set. IEEE Control Systems Letters, 2024, 8, pp.1331-1336. 10.1109/LCSYS.2024.3410639. hal-04607918

HAL Id: hal-04607918

https://hal.science/hal-04607918v1

Submitted on 12 Nov 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Pearson coefficient degradation in a Wasserstein/Gelbrich ambiguity set

Matthieu Borelle, Teodoro Alamo, Cristina Stoica, Sylvain Bertrand, Eduardo F. Camacho

**Abstract**— This paper presents new results on the Gelbrich distance and the corresponding ambiguity sets, to analyze the correlation between two scalar random variables. A closed expression of the minimum disturbance in the Gelbrich metric necessary to achieve a specified correlation between two random variables is proposed. This expression allows us to analytically compute the minimal absolute correlation in a Gelbrich ball. This analysis facilitates the assessment of the robustness of the Pearson coefficient within an ambiguity set. Two numerical examples illustrate the validity of the proposed results.

## I. INTRODUCTION

The Pearson correlation coefficient [1] is a statistical measure aiming to assess the magnitude and direction of linear relationships between two variables. Its application extends across multiple statistical domains, encompassing areas like data analysis, variable selection, classification, clustering, decision-making, financial analysis, biological research, etc. However, the covariance and correlation coefficient can be very sensitive to outlier observations/measurement errors (e.g., bivariate case application in life sciences [2]). Moreover, relying solely on a correlation measure, typically computed empirically from a limited sample, could neglect crucial aspects of the data, such as the variability of a random sample data draw or its potential time-varying nature. Similarly, in state estimation, the lack of knowledge of the state-space model and of the noise distribution can cause distributional uncertainties. Classic estimators like the Kalman Filter cannot efficiently handle distributional uncertainties and can lead to poor estimation performance or even divergence in certain pathological cases [3]. This drawback motivates the development of *Distributionally Robust Optimization* (DRO) methods [4], [5], which seek solutions that are robust with respect to uncertainties in the underlying probability distribution. DRO approaches involve optimizing performance through a cost function evaluated under the worst-case scenario within a specified ambiguity

M. Borelle and C. Stoica are with Université Paris-Saclay, CNRS, CentraleSupélec, Laboratoire des signaux et systèmes, 91190, Gif-sur-Yvette, France (e-mail: {matthieu.borelle; cristina.stoica}@l2s.centralesupelec.fr).

M. Borelle and S. Bertrand are with Université Paris-Saclay, ONERA, Traitement de l'information et systèmes, 91123, Palaiseau, France (e-mail: {matthieu.borelle; sylvain.bertrand}@onera.fr).

T. Alamo and E.F. Camacho are with Department of Ingegniería de Sistemas y Automática, Universidad de Sevilla, Camino de los Descubrimientos, 41092 Sevilla, Spain (e-mail: talamo@us.es, eduardo@esi.us.es).

The authors acknowledge support from Agence de l'innovation de défense (grant 2022-65-0011-AID_ONERA), Agence Nationale de Recherche (grant ANR-21-CE48-0003), and the European Research Council under the OCONT SOLAR grant agreement ID: 789051. T. Alamo acknowledges support from grant PID2022-142946NA-100 funded by MCIN/AEI/10.13039/501100011033 and by ERDF, A way of making Europe.

set [6] of distributions, rather than assuming a fixed or known probability distribution.

The ambiguity set is often defined as a ball in the space of probability distributions that contains all distributions close to a nominal or a priori most likely distribution with respect to a given probability metric [6]. The ambiguity set radius therefore reflects the level of confidence in the nominal distribution.

Different notions of distance between distributions are available [7], e.g., Wasserstein distances [6] used to assess the difference between two probability measures. While in general, the Wasserstein metric is difficult to compute, the Gelbrich distance [8], [9] offers a simple expression to lower bound the Wasserstein distance, especially in the Gaussian case [10]. The Wasserstein and Gelbrich distances are versatile tools with broad machine learning and control applications, e.g., generative adversarial networks [11], outlier detection problem [12], Kalman filtering [13], [14], optimal control [15], predictive control [16]. Distributionally robust optimization problems are computationally challenging due to several factors, e.g., high dimensionality, complexity of the constraints related to the ambiguity set shape. Convex relaxations, dual representations or simplifying assumptions are generally used to address these challenges [9]. Finding tractable solutions to reduce the complexity and computational burden is of paramount importance for implementation and real-time applications.

The motivation of this work derives from the possible degradation of the Pearson coefficient in an ambiguity set. Therefore, in this paper, a distributionally robust optimization problem involving finding the minimal correlation in an ambiguity set is analytically solved for the case of bivariate distributions. In addition, this paper extends the results of minimum Gelbrich distance required to reach uncorrelation between two scalar random variables, proposed by the authors in [17], to the computation of the maximal degradation of the Pearson coefficient in a Gelbrich ball, which is a considerably more challenging problem. The proposed approach enables the identification of variables that exhibit both strong correlation and robustness to various data conditions.

This paper is organized as follows. Section II focuses on the Wasserstein and Gelbrich distances, and the Pearson correlation coefficient. A closed expression of the minimal Gelbrich distance to reach a given correlation index is provided and illustrated in Section III. The minimal absolute value of the Pearson correlation coefficient in a given Gelbrich ambiguity set is derived and illustrated in Section IV.

**Notation.** Consider $\mathbb{S}_+^d = \{ S \in \mathbb{S}^d : S \succeq 0 \}$, the set of symmetric positive semidefinite matrices in $\mathbb{R}^{d \times d}$. Denote
by $S^{\frac{1}{2}}$ the square root of $S \in \mathbb{S}_+^d$ satisfying $(S^{\frac{1}{2}})^2 = S$. The expectation of a random variable $x \in \mathbb{R}^n$ under a $n$-dimension probability distribution $\mathbb{Q}_S$ is denoted by $E_{\mathbb{Q}_S}\{x\}$.

## II. THEORETICAL BACKGROUND

This section provides the necessary elements on the Wasserstein and Gelbrich distances, and the Pearson coefficient.

### A. Wasserstein distance

The Wasserstein $p$-distance [6], [10] is a metric which provides a measure between two probability distributions, thus quantifying the effort necessary to transform one distribution into another.

**Definition 1 (Wasserstein $p$-distance [9], [10]):** Using, a metric space $(M, c)$, e.g., $M = \mathbb{R}^d$ and $c(x, y) = \|x - y\|_2$, $\forall x, y \in M$, define the Wasserstein $p$-distance $W_p(\mu, \nu)$ between two probability measures $\mu$ and $\nu$, with finite $p$-moments on $M$, as follows

$$W_p(\mu, \nu) = \left( \inf_{\pi \in \Pi(\mu, \nu)} \left\{ \int_{M \times M} c(x, y)^p \pi(dx, dy) \right\} \right)^{\frac{1}{p}}, \quad (1)$$

with $\Pi(\mu, \nu)$ the set of all probability distributions $\pi$ defined on $M \times M$ with marginals $\mu$ and $\nu$.

In general, the Wasserstein $p$-distance is difficult to compute. However, it is possible to derive simpler expressions to compute the Wasserstein 2-distance between Gaussian distributions, as detailed below.

### B. Gelbrich distance

The Gelbrich distance [8] provides a lower bound for the Wasserstein distance.

**Definition 2 (Gelbrich distance [8], [9]):** Consider the tuples $(\mu_\Sigma, \Sigma) \in \mathbb{R}^n \times \mathbb{S}_+^n$ and $(\mu_S, S) \in \mathbb{R}^n \times \mathbb{S}_+^n$ denoting the mean and covariance matrices of random distributions $\mathbb{Q}_\Sigma$ and $\mathbb{Q}_S$. The Gelbrich distance $G_D(\mathbb{Q}_\Sigma, \mathbb{Q}_S)$ between these distributions is given by

$$G_D(\mathbb{Q}_\Sigma, \mathbb{Q}_S) = \left( \|\mu_\Sigma - \mu_S\|_2^2 + \text{Tr} \left( S + \Sigma - 2(\Sigma^{\frac{1}{2}} S \Sigma^{\frac{1}{2}})^{\frac{1}{2}} \right) \right)^{\frac{1}{2}}, \quad (2)$$

with Tr denoting the trace of the considered matrix.

The Gelbrich distance provides a lower bound of the Wasserstein 2-distance [8], [9], as specified in Property 1.

**Property 1:** Consider $(\mu_\Sigma, \Sigma)$ and $(\mu_S, S)$ the mean and covariance matrices of distributions $\mathbb{Q}_\Sigma$ and $\mathbb{Q}_S$. Then, the Wasserstein 2-distance $W_2(\mathbb{Q}_\Sigma, \mathbb{Q}_S)$ between distributions $\mathbb{Q}_\Sigma$ and $\mathbb{Q}_S$ verifies

$$W_2(\mathbb{Q}_\Sigma, \mathbb{Q}_S) \geq G_D(\mathbb{Q}_\Sigma, \mathbb{Q}_S). \quad (3)$$

If $\mathbb{Q}_\Sigma$ and $\mathbb{Q}_S$ are Gaussian distributions, then the equality $W_2(\mathbb{Q}_\Sigma, \mathbb{Q}_S) = G_D(\mathbb{Q}_\Sigma, \mathbb{Q}_S)$ holds [9].

The next property (see [17] for a proof) offers a closed expression of the Gelbrich distance for the particular case of two-dimensional distributions.

**Property 2:** Defining $(\mu_\Sigma, \Sigma)$ and $(\mu_S, S)$ the mean and covariance matrices of the two-dimensional distributions $\mathbb{Q}_\Sigma$ and $\mathbb{Q}_S$, the Gelbrich

distance is computed as follows $G_D(\mathbb{Q}_\Sigma, \mathbb{Q}_S) = \left( \|\mu_\Sigma - \mu_S\|_2^2 + \text{Tr}(S + \Sigma) - 2\sqrt{\text{Tr}(S\Sigma) + 2\sqrt{\det(S\Sigma)}} \right)^{\frac{1}{2}}$.

### C. Ambiguity sets

The Wasserstein and Gelbrich distances can be used to define an ambiguity set.

**Definition 3 (Wasserstein and Gelbrich ambiguity sets [6]):** Consider the set $\mathcal{P}_d$ denoting all probability distributions defined on a $d$-dimensional space. Define the Gelbrich ambiguity set $\mathcal{A}_\alpha^G$ and the Wasserstein ambiguity set $\mathcal{A}_\alpha^W$, centred at $\mathbb{Q}_\Sigma \in \mathcal{P}_d$ with the radius $\alpha \ge 0$ as follows

$$\mathcal{A}_{\alpha}^{G}(\mathbb{Q}_{\Sigma}) = \{\mathbb{Q}_{S} \in \mathcal{P}_{d} \mid G_{D}(\mathbb{Q}_{\Sigma}, \mathbb{Q}_{S}) \le \alpha\},$$

$$\mathcal{A}_{\alpha}^{W}(\mathbb{Q}_{\Sigma}) = \{\mathbb{Q}_{S} \in \mathcal{P}_{d} \mid W_{2}(\mathbb{Q}_{\Sigma}, \mathbb{Q}_{S}) \le \alpha\}.$$

Since the Gelbrich distance is upper bounded by the Wasserstein 2-distance (3), we have $\mathcal{A}_\alpha^W(\mathbb{Q}_\Sigma) \subseteq \mathcal{A}_\alpha^G(\mathbb{Q}_\Sigma)$, $\forall \alpha \ge 0$.

### D. Pearson coefficient

The Pearson correlation coefficient is a statistical measure characterizing the degree of linear relationship between two variables. Consider the random scalar variables $x$ and $y$, with their corresponding bivariate distribution $\mathbb{Q}_S \in \mathcal{P}_2$. Denote by $S \in \mathbb{S}_+^2$ the covariance matrix

$$S = \begin{bmatrix} S_{xx} & S_{xy} \\ S_{xy} & S_{yy} \end{bmatrix} = E_{\mathbb{Q}_S} \left\{ \begin{bmatrix} x - \mu_x \\ y - \mu_y \end{bmatrix} \begin{bmatrix} x - \mu_x \\ y - \mu_y \end{bmatrix}^{\top} \right\}, \quad (4)$$

where $E_{\mathbb{Q}_S}\{\cdot\}$ is the expected value operator for $\mathbb{Q}_S$, with $\mu_x = E_{\mathbb{Q}_S}\{x\}$ and $\mu_y = E_{\mathbb{Q}_S}\{y\}$. The Pearson coefficient [1], also known as the correlation coefficient, for the distribution $\mathbb{Q}_S$ is defined as

$$\rho(\mathbb{Q}_S) = \frac{E_{\mathbb{Q}_S}\{(x - \mu_x)(y - \mu_y)\}}{\sqrt{E_{\mathbb{Q}_S}\{(x - \mu_x)^2\}E_{\mathbb{Q}_S}\{(y - \mu_y)^2\}}} = \frac{S_{xy}}{\sqrt{S_{xx}S_{yy}}}.$$

The Pearson coefficient takes values within $[-1, 1]$. A value close to the interval bounds indicates a strong linear relationship, while a value close to 0 is related to a weak or no linear relationship (no linear correlation). The sign indicates the direction of the relationship, while the magnitude represents the strength [1].

Given $x$, consider the following estimator for $y$: $\hat{y} = \mu_y + \lambda(x - \mu_x)$. It is well-known in the context of bivariate analysis, see e.g., chapter 2 in [18], that the variance of the estimation error satisfies

$$\sigma_e^2(\mathbb{Q}_S) = E_{\mathbb{Q}_S}\{(y - \hat{y})^2\} \geq S_{yy}(1 - \rho^2(\mathbb{Q}_S)), \forall \lambda \in \mathbb{R}, \quad (5)$$

where the equality is attained with the optimal choice $\lambda = \frac{S_{xy}}{S_{xx}}$, which assumes precise knowledge of the joint probability distribution. However, if the distribution belongs to an ambiguity set, the estimation problem becomes more challenging. It requires considering not only the nominal Pearson coefficient (center of the ambiguity set) but also its potential variations within the ambiguity set.
### III. MINIMAL GELBRICH PERTURBATION TO ATTAIN A GIVEN CORRELATION INDEX $\rho$

This section analyses the robustness of a given nominal Pearson correlation coefficient with respect to potential perturbations. We examine the minimum perturbation to a given bivariate distribution (measured by means of the Gelbrich distance) required to achieve a specified correlation value $\rho$.

#### A. Problem statement

Let us introduce the concept of the $S(\rho)$ set, which is instrumental in obtaining the main contribution of Section III.

**Definition 4 (S(ρ) set):** Given $\rho \in [-1, 1]$, define the set $S(\rho)$ as the set of bivariate distributions $Q_S$ with a covariance matrix $S \in \mathbb{S}_+^2$ that can be written as

$$S = \begin{bmatrix} \sigma_1^2 & \rho\sigma_1\sigma_2 \\ \rho\sigma_1\sigma_2 & \sigma_2^2 \end{bmatrix}. \quad (6)$$

Notice that the previous definition (4) does not constrain $\sigma_1$ or $\sigma_2$ to be positive values. Given a distribution $Q_\Sigma$ and $\rho \in [-1, 1]$, we focus now on the problem of computing

$$D_{\rho} = \min_{Q_S \in S(\rho)} G_D(Q_{\Sigma}, Q_S),$$

which is the minimal Gelbrich distance of $Q_\Sigma$ to the set $S(\rho)$. In [17], the authors solve this problem for the particular case $\rho = 0$ (i.e., minimal Gelbrich distance to uncorrelation $D_0$).

#### B. Minimal Gelbrich distance to S(ρ)

The next property derives a closed form for the minimal Gelbrich distance $D_\rho$ of $Q_\Sigma$ to the set $S(\rho)$, provided that the specified correlation $\rho$ is of the same sign as the original one (i.e., $\rho\Sigma_{xy} \ge 0$). We notice that if $\rho$ has a different sign, then $D_\rho$ would be larger than $D_0$.

**Property 3:** Suppose that $\mu_\Sigma \in \mathbb{R}^2$ and $\Sigma = \begin{bmatrix} \Sigma_{xx} & \Sigma_{xy} \\ \Sigma_{xy} & \Sigma_{yy} \end{bmatrix} \in \mathbb{S}_+^2$ represent the mean and covariance matrix of the bivariate distribution $Q_\Sigma$, respectively. Given $\rho \in [-1, 1]$, such that $\rho\Sigma_{xy} \ge 0$, then the minimal Gelbrich distance $D_\rho$ of $Q_\Sigma$ to the set $S(\rho)$ satisfies

$$D_{\rho}^{2} = \min_{Q_S \in S(\rho)} G_{D}^{2}(Q_{\Sigma}, Q_{S}) = \frac{\mathrm{Tr}(\Sigma) - \sqrt{\psi_{\rho}}}{2}, \quad (7)$$

$$\text{with } \psi_{\rho} = (\Sigma_{xx} - \Sigma_{yy})^2 + 4 \left( \rho \Sigma_{xy} + \sqrt{(1 - \rho^2) \det(\Sigma)} \right)^2.$$

**Proof:** Considering the closed-form expression for the Gelbrich distance from Property 2, the following result holds

$$D_{\rho}^{2} = \min_{\mu_{S} \in \mathbb{R}^{2}, S \in S(\rho)} \left( \|\mu_{\Sigma} - \mu_{S}\|_{2}^{2} + \mathrm{Tr}(S + \Sigma) - 2\sqrt{\mathrm{Tr}(S\Sigma) + 2\sqrt{\det(S\Sigma)}} \right),$$

leading to

$$D_{\rho}^{2} = \min_{S \in S(\rho)} \left( \mathrm{Tr}(S + \Sigma) - 2\sqrt{\mathrm{Tr}(S\Sigma) + 2\sqrt{\det(S\Sigma)}} \right).$$

Using the matrix $S$ from (6), the terms $\mathrm{Tr}(S + \Sigma)$, $\mathrm{Tr}(S\Sigma)$ and $\det(S\Sigma)$ can be developed as follows

$$\mathrm{Tr}(S + \Sigma) = \mathrm{Tr}(\Sigma) + \sigma_1^2 + \sigma_2^2 \quad (8)$$

$$\mathrm{Tr}(S\Sigma) = \sigma_1^2 \Sigma_{xx} + 2\sigma_1\sigma_2\rho\Sigma_{xy} + \sigma_2^2 \Sigma_{yy}$$

$$\det(S\Sigma) = \sigma_1^2 \sigma_2^2 (1 - \rho^2) \det(\Sigma). \quad (9)$$

Defining now $h(\sigma_1, \sigma_2) = \mathrm{Tr}(S\Sigma) + 2\sqrt{\det(S\Sigma)} = \sigma_1^2 \Sigma_{xx} + 2\sigma_1\sigma_2\rho\Sigma_{xy} + \sigma_2^2 \Sigma_{yy} + 2\sqrt{\sigma_1^2 \sigma_2^2 (1 - \rho^2) \det(\Sigma)}$, we rewrite the original optimization problem as

$$D_{\rho}^{2} = \min_{\sigma_{1} \in \mathbb{R}, \sigma_{2} \in \mathbb{R}} \left( \mathrm{Tr}(\Sigma) + \sigma_{1}^{2} + \sigma_{2}^{2} - 2\sqrt{h(\sigma_{1}, \sigma_{2})} \right). \quad (10)$$

By simple inspection, we notice that the signs of $\sigma_1$ and $\sigma_2$ have no effect on the functional to be minimized except for the term $\sigma_1\sigma_2\rho\Sigma_{xy}$ appearing in $h(\sigma_1, \sigma_2)$. Since $\rho\Sigma_{xy}$ is assumed to be non-negative, negative values for $\sigma_1\sigma_2$ would decrease $h(\sigma_1, \sigma_2)$ and thus increase the functional to be minimized. Thus, we conclude that $\sigma_1\sigma_2 \ge 0$ at the optimum that defines $D_\rho$. Thus, we can restrict the search of the minimum defined in (10) to the values of $\sigma_1$ and $\sigma_2$ satisfying $\sigma_1\sigma_2 \ge 0$, i.e.,

$$D_{\rho}^{2} = \min_{\sigma_{1} \in \mathbb{R}, \sigma_{2} \in \mathbb{R}, \sigma_{1}\sigma_{2} \ge 0} \left( \mathrm{Tr}(\Sigma) + \sigma_{1}^{2} + \sigma_{2}^{2} - 2\sqrt{h(\sigma_{1}, \sigma_{2})} \right).$$

Under the assumption $\sigma_1\sigma_2 \ge 0$, $h(\sigma_1, \sigma_2)$ has the next form

$$h(\sigma_1, \sigma_2) = \sigma_1^2 \Sigma_{xx} + 2\sigma_1\sigma_2\rho\Sigma_{xy} + \sigma_2^2 \Sigma_{yy} + 2\sigma_1\sigma_2\sqrt{(1 - \rho^2) \det(\Sigma)}.$$

Thus, under the assumption $\sigma_1\sigma_2 \ge 0$, $h(\sigma_1, \sigma_2)$ is a quadratic function that can be rewritten as

$$h(\sigma_1, \sigma_2) = \begin{bmatrix} \sigma_1 \\ \sigma_2 \end{bmatrix}^\top H(\rho) \begin{bmatrix} \sigma_1 \\ \sigma_2 \end{bmatrix},$$

with

$$H(\rho) = \begin{bmatrix} \Sigma_{xx} & \rho\Sigma_{xy} + \sqrt{(1 - \rho^2) \det(\Sigma)} \\ * & \Sigma_{yy} \end{bmatrix},$$

where “*” denotes the necessary term for the symmetry of $H(\rho)$. Thus, $D_\rho$ verifies

$$D_{\rho}^{2} = \min_{\sigma_{1} \in \mathbb{R}, \sigma_{2} \in \mathbb{R}, \sigma_{1}\sigma_{2} \ge 0} \left( \mathrm{Tr}(\Sigma) + \sigma_{1}^{2} + \sigma_{2}^{2} - 2\sqrt{\begin{bmatrix} \sigma_{1} \\ \sigma_{2} \end{bmatrix}^{\top} H(\rho) \begin{bmatrix} \sigma_{1} \\ \sigma_{2} \end{bmatrix}} \right).$$

Since $\rho\Sigma_{xy}$ is non-negative, the off-diagonal element of $H(\rho)$ is non-negative, and thus, negative values for $\sigma_1\sigma_2$ would increase the functional. Therefore, removing the constraint $\sigma_1\sigma_2 \ge 0$ does not change the optimal value (i.e., the constraint $\sigma_1\sigma_2 \ge 0$ is not active at the optimum). This leads to $D_\rho^2 = \min_{\sigma_1 \in \mathbb{R}, \sigma_2 \in \mathbb{R}} \left( \mathrm{Tr}(\Sigma) + \sigma_1^2 + \sigma_2^2 - 2\sqrt{\begin{bmatrix} \sigma_1 \\ \sigma_2 \end{bmatrix}^\top H(\rho) \begin{bmatrix} \sigma_1 \\ \sigma_2 \end{bmatrix}} \right)$. For a fixed value of $\tau^2 = \sigma_1^2 + \sigma_2^2$, the values of $\sigma_1, \sigma_2$ that minimize the functional
are the ones corresponding to the eigenvector of the largest eigenvalue of $H(\rho)$, that we denote by $\lambda_{max}$. Due to the Perron-Frobenious theorem [19], that states that the largest eigenvalue of a matrix with all its components non-negative is non-negative, we have that $\lambda_{max} \ge 0$. From all this, it can be inferred that

$$ 
\begin{aligned}
D_{\rho}^{2} &= \min_{\sigma_{1} \in \mathbb{R}, \sigma_{2} \in \mathbb{R}} \left( \operatorname{Tr}(\Sigma) + \sigma_{1}^{2} + \sigma_{2}^{2} \right. \\
&\qquad \left. - 2\sqrt{\lambda_{max}(\sigma_{1}^{2} + \sigma_{2}^{2})} \right) \\
&= \min_{\tau \in \mathbb{R}} \left( \operatorname{Tr}(\Sigma) + \tau^{2} - 2\sqrt{\tau^{2}\lambda_{max}} \right) \\
&= \min_{\tau \ge 0} \left( \operatorname{Tr}(\Sigma) + \tau^{2} - 2\tau\sqrt{\lambda_{max}} \right).
\end{aligned}
$$

The optimal value $\tau^*$ of $\tau$ satisfies $2\tau^* - 2\sqrt{\lambda_{max}} = 0$. Thus, $\tau^* = \sqrt{\lambda_{max}}$ and the minimal Gelbrich distance of $\mathbb{Q}_{\Sigma}$ to the set $S(\rho)$ verifies

$$ D_{\rho}^{2} = \min_{\tau \ge 0} \left( \operatorname{Tr}(\Sigma) + \tau^{2} - 2\tau\sqrt{\lambda_{max}} \right) = \operatorname{Tr}(\Sigma) - \lambda_{max}. \quad (11) $$

To determine $\lambda_{max}$, the characteristic polynomial of the matrix $H(\rho)$ can be written as

$$ 
\begin{aligned}
0 &= \det(\lambda I - H(\rho)) = \lambda^2 - \operatorname{Tr}(H(\rho))\lambda + \det(H(\rho)) \\
  &= \lambda^2 - a_1\lambda + a_2,
\end{aligned}
\quad (12) $$

where

$$ 
\begin{aligned}
a_1 &= \operatorname{Tr}(H(\rho)) = \operatorname{Tr}(\Sigma) \\
a_2 &= \Sigma_{xx}\Sigma_{yy} - (\rho\Sigma_{xy} + \sqrt{(1-\rho^2)\det(\Sigma)})^2.
\end{aligned}
$$

The largest eigenvalue of $H(\rho)$ is given by the largest root of (12)

$$ \lambda_{\max} = \frac{a_1 + \sqrt{a_1^2 - 4a_2}}{2}. \quad (13) $$

Denoting by $\psi_{\rho}$ the discriminant $a_1^2 - 4a_2$ leads to

$$ 
\begin{aligned}
\psi_{\rho} &= a_1^2 - 4a_2 = (\Sigma_{xx} + \Sigma_{yy})^2 \\
&\quad - 4 \left( \Sigma_{xx}\Sigma_{yy} - (\rho\Sigma_{xy} + \sqrt{(1-\rho^2)\det(\Sigma)})^2 \right) \\
&= (\Sigma_{xx} - \Sigma_{yy})^2 + 4 \left( \rho\Sigma_{xy} + \sqrt{(1-\rho^2)\det(\Sigma)} \right)^2.
\end{aligned}
$$

Thus, the discriminant is non-negative and the largest eigenvalue of matrix $H(\rho)$ is $\lambda_{max} = (\operatorname{Tr}(\Sigma) + \sqrt{\psi_{\rho}})/2$. From (11) we get $D_{\rho}^2 = \operatorname{Tr}(\Sigma) - \lambda_{max} = (\operatorname{Tr}(\Sigma) - \sqrt{\psi_{\rho}})/2$. ■

Notice that the set $S(\rho)$ consists of two disjoint parts corresponding to correlation $\rho$ and $-\rho$. Therefore, the minimum in (7) is always achieved at the distribution with the same sign of correlation. The assumption $\rho\Sigma_{xy} \ge 0$ is considered to simplify the proof of Property 3, and it does not interfere with the main result presented in the next section.

We observe that Property 3 enables the computation of the minimal Gelbrich perturbation required to attain a specified correlation $\rho$, provided it has the same sign as the original one (i.e., $\rho\Sigma_{xy} \ge 0$). In particular, setting $\rho$ equal to zero recovers the main result of [17], which was an

explicit expression for the minimal Gelbrich distance $D_{\mathcal{U}}$ to uncorrelation

$$ D_0^2 = D_{\mathcal{U}}^2 = \frac{\operatorname{Tr}(\Sigma) - \sqrt{(\Sigma_{xx} - \Sigma_{yy})^2 + 4 \det(\Sigma)}}{2}. \quad (14) $$

### C. Numerical example

In this numerical example a given bivariate covariance matrix $\Sigma = \begin{bmatrix} 2 & 1 \\ 1 & 1 \end{bmatrix}$ has been chosen with a Pearson correlation coefficient $\rho_0 = 1/\sqrt{2}$. The minimal distance to uncorrelation $D_{\mathcal{U}} = 0.618$ is computed with (14) and we randomly sample 3000 bivariate covariance matrices $\mathbb{Q}_S$ inside a Gelbrich ambiguity set centered on $\Sigma$, with the radius $\alpha = 1.5D_{\mathcal{U}}$. This choice of the ambiguity set radius permits to generate a bivariate covariance matrix with possible uncorrelation. For each sampled covariance matrix $\mathbb{Q}_S$, its Gelbrich distance to $\Sigma$ is plotted (blue dots) in Fig. 1 with respect to its Pearson correlation coefficient $\rho(\mathbb{Q}_S)$. The simulation results validate the Gelbrich distance lower bound (red points) given by (7) for the considered assumption $\rho\Sigma_{xy} \ge 0$, as all the sample covariance matrices (in blue) lie above the red dots. Expression (7) is a sharp lower bound and it is reached by some of the sampled distributions in Fig. 1.

<figure><img src="image_2.png" /><figcaption>Fig. 1. Illustrating the minimal Gelbrich distance to reach a correlation.</figcaption></figure>

## IV. MINIMAL CORRELATION IN A GELBRICH BALL

This section proposes a closed-form expression of the minimal correlation in a Gelbrich ball.

### A. Problem statement

Denote by $\mathcal{P}_2$ the set of all the probability distributions of dimension 2. Given the central distribution $\mathbb{Q}_{\Sigma} \in \mathcal{P}_2$ and the radius $\alpha \ge 0$, we concentrate now on determining the smallest absolute value that attains the Pearson correlation coefficient in the ambiguity set $\mathcal{A}_{\alpha}^{\text{G}}(\mathbb{Q}_{\Sigma}) = \{\mathbb{Q}_S \in \mathcal{P}_2 | G_D(\mathbb{Q}_{\Sigma}, \mathbb{Q}_S) \le \alpha\}$. Therefore, we address the problem

$$ \min_{\mathbb{Q}_S} |\rho(\mathbb{Q}_S)| \quad (15) $$

$$ \text{s.t.} \quad G_D(\mathbb{Q}_{\Sigma}, \mathbb{Q}_S) \le \alpha. \quad (16) $$
## B. Closed-form expression for the minimal correlation

The next property provides a closed-form expression for the solution of optimization problem (15). This allows us to determine the minimum Pearson correlation coefficient (absolute value) in an ambiguity set defined by the Gelbrich distance.

**Property 4:** Given the distribution $Q_Σ \in P_2$, with covariance matrix $\Sigma = \begin{bmatrix} \Sigma_{xx} & \Sigma_{xy} \\ \Sigma_{xy} & \Sigma_{yy} \end{bmatrix}$, where $\Sigma_{xx} > 0$ and $\Sigma_{yy} > 0$, consider the Gelbrich ambiguity set $A_α^G(Q_Σ)$, where $α > 0$. Then, using the notations

$$ 
\begin{aligned}
\rho_0 &= \frac{\Sigma_{xy}}{\sqrt{\Sigma_{xx}\Sigma_{yy}}}, \\
D_U^2 &= \frac{\text{Tr}(\Sigma) - \sqrt{(\Sigma_{xx} - \Sigma_{yy})^2 + 4 \det(\Sigma)}}{2}, \\
\beta &= \frac{(\text{Tr}(\Sigma) - 2\alpha^2)^2 - (\Sigma_{xx} - \Sigma_{yy})^2}{4\Sigma_{xx}\Sigma_{yy}}, \\
\eta &= 1 + \frac{\beta - 1}{\rho_0^2},
\end{aligned}
$$

the Pearson coefficient $\rho_α$ with minimum absolute value in $A_α^G(Q_Σ)$ is

$$ \rho_{\alpha} = \begin{cases} 0 & \text{if } \alpha \ge D_{\mathcal{U}}, \\ \rho_0 (\sqrt{\beta} - \sqrt{\beta - \eta}) & \text{otherwise.} \end{cases} \quad (17) $$

**Proof:** From (14), $D_U$ represents the minimum Gelbrich distance to uncorrelation. Thus, it leads to $\rho_α = 0$ if $\alpha \ge D_U$, and $\rho_α^2 > 0$ otherwise. Focus now on the case $\alpha < D_U$ (or equivalently, $\rho_α^2 > 0$). It is not difficult to see that, in this case, constraint (16) must be active. If the constraint is satisfied with strict inequality, then it would be possible to slightly reduce the absolute value of $S_{xy}$ (see (2)), the off-diagonal element of the covariance matrix $S$ of $Q_S$ in (4), to obtain a strictly smaller value of $\rho(Q_S)$ without violating the constraint (16). This means that, provided $\alpha < D_U$, the minimum absolute value of the correlation coefficient is attained at a point that is at a Gelbrich distance equal to $\alpha$ from the central distribution $Q_Σ$. From this, and using Definition 4, it is inferred that the value of $\rho_α$ corresponds to the value for which the distance of the central distribution to set $S(\rho_α)$ is exactly equal to $\alpha$. Therefore, the optimal value $\rho_α$ satisfies

$$ \alpha = \min_{Q_S \in S(\rho_\alpha)} G_D(Q_\Sigma, Q_S) = D_{\rho_\alpha}, \quad (18) $$

where $D_{\rho_α}$ represents the minimal the Gelbrich distance required to reach the set of matrices with correlation coefficient $\rho_α$. Since $\alpha$ is smaller than the minimum Gelbrich distance to uncorrelation, $\rho_α$ and $\rho_0$ have the same sign, which implies that $\rho_α\Sigma_{xy} \ge 0$. Thus, we can use the closed expression of $D_{\rho_α}$ presented in Property 3. We obtain that $\rho_α$ is the solution with the smallest absolute value to the equation

$$ \alpha^2 = \frac{\text{Tr}(\Sigma) - \sqrt{\psi(\rho_\alpha)}}{2}, \quad (19) $$

with

$$ \psi(\rho_{\alpha}) = (\Sigma_{xx} - \Sigma_{yy})^2 + 4 \left( \rho_{\alpha} \Sigma_{xy} + \sqrt{(1 - \rho_{\alpha}^2) \det(\Sigma)} \right)^2. $$

Rewrite (19) as follows $\sqrt{\psi(\rho_α)} = \text{Tr}(\Sigma) - 2\alpha^2$. Since we are considering the case $\alpha^2 < D_U^2$, we have from the expression for $D_U$ that $\text{Tr}(\Sigma) - 2\alpha^2 \ge 0$. This means that (19) can be further rewritten as $\psi(\rho_α) = (\text{Tr}(\Sigma) - 2\alpha^2)^2$. From the expression of $\psi(\rho_α)$, we obtain

$$ (\text{Tr}(\Sigma) - 2\alpha^2)^2 = (\Sigma_{xx} - \Sigma_{yy})^2 + 4 \left( \rho_{\alpha} \Sigma_{xy} + \sqrt{(1 - \rho_{\alpha}^2) \det(\Sigma)} \right)^2 $$

or, equivalently,

$$ \left( \rho_{\alpha} \Sigma_{xy} + \sqrt{(1 - \rho_{\alpha}^2) \det(\Sigma)} \right)^2 = \frac{(\text{Tr}(\Sigma) - 2\alpha^2)^2}{4} - \frac{(\Sigma_{xx} - \Sigma_{yy})^2}{4}. $$

Dividing by $\Sigma_{xx}\Sigma_{yy}$ and using the notation

$$ \beta = \frac{(\text{Tr}(\Sigma) - 2\alpha^2)^2 - (\Sigma_{xx} - \Sigma_{yy})^2}{4\Sigma_{xx}\Sigma_{yy}}, \quad (20) $$

leads to

$$ \frac{1}{\Sigma_{xx}\Sigma_{yy}} \left( \rho_{\alpha}\Sigma_{xy} + \sqrt{(1 - \rho_{\alpha}^2)\det(\Sigma)} \right)^2 = \beta. $$

Simple calculations show that $\alpha < D_U$ implies $\beta \ge 0$. Thus,

$$ \frac{1}{\sqrt{\Sigma_{xx}\Sigma_{yy}}} \left( \rho_{\alpha}\Sigma_{xy} + \sqrt{(1 - \rho_{\alpha}^2)\det(\Sigma)} \right) = \pm \sqrt{\beta}. $$

Therefore, since $\rho_α\Sigma_{xy} \ge 0$, the unique possible solution is

$$ \frac{1}{\sqrt{\Sigma_{xx}\Sigma_{yy}}} \left( \rho_{\alpha}\Sigma_{xy} + \sqrt{(1 - \rho_{\alpha}^2)\det(\Sigma)} \right) = \sqrt{\beta}. $$

With the notation $\rho_0 = \frac{\Sigma_{xy}}{\sqrt{\Sigma_{xx}\Sigma_{yy}}}$, we obtain

$$ \rho_{\alpha}\rho_0 + \sqrt{(1 - \rho_{\alpha}^2)(1 - \rho_0^2)} = \sqrt{\beta}, $$

which can be rewritten as

$$ \rho_{\alpha}^2 - 2\rho_{\alpha}\rho_0\sqrt{\beta} + \rho_0^2\eta = 0, \quad (21) $$

with $\eta = 1 + (\beta - 1)/\rho_0^2$. One can prove that the discriminant<fnref n="1" /> of this quadratic equation is positive $4\rho_0^2\beta - 4\rho_0^2\eta \ge 0$, and thus, the solutions of (21) are $\rho_α = \rho_0\sqrt{\beta} \pm \sqrt{\rho_0^2\beta - \rho_0^2\eta}$. Since the solution of the original optimization problem (15) considers the minimal absolute value of $\rho_α$, we have

$$ \rho_{\alpha} = \begin{cases} \rho_0\sqrt{\beta} - \sqrt{\rho_0^2\beta - \rho_0^2\eta} & \text{if } \rho_0 > 0, \\ \rho_0\sqrt{\beta} + \sqrt{\rho_0^2\beta - \rho_0^2\eta} & \text{otherwise,} \end{cases} $$

which is equivalent to $\rho_α = \rho_0(\sqrt{\beta} - \sqrt{\beta - \eta})$. ■ We notice that the inclusion $A_α^W(Q_Σ) \subseteq A_α^G(Q_Σ)$, $\forall \alpha \ge 0$, see Subsection II-C, implies

$$ \min_{Q_S \in A_\alpha^W(Q_\Sigma)} |\rho(Q_S)| \ge \min_{Q_S \in A_\alpha^G(Q_\Sigma)} |\rho(Q_S)|. $$

Thus, Property 4 can be used to obtain a lower bound of the smallest absolute value of the Pearson coefficient in $A_α^W(Q_Σ)$.

<footnote n="1">¹For $\alpha < D_U$, it can be proven that $\beta \ge \eta$, or equivalently $(\beta - 1)(\rho_0^2 - 1) \ge 0$. Since $\rho_0 \in [-1, 1]$, it remains to prove that $\beta - 1 \le 0$. Using (20) we obtain $\beta - 1 = \alpha^2(\alpha^2 - \Sigma_{xx} - \Sigma_{yy})/(\Sigma_{xx}\Sigma_{yy})$. From (19), we also have $\alpha^2 \le \Sigma_{xx} + \Sigma_{yy}$. Thus, we get $\beta - 1 \le 0$.</footnote>
### C. Numerical example

We consider the same bivariate covariance matrix $\Sigma$ as in subsection III-C, with a Pearson correlation coefficient $\rho_0 = 1/\sqrt{2}$. The minimal distance to uncorrelation $D_U = 0.618$ is computed with equation (14) and 3000 bivariate covariance matrices are also randomly drawn inside a Gelbrich ambiguity set centered on $\Sigma$, with the radius $\alpha = 1.5D_U$. We compute and display their absolute Pearson correlation coefficient as a function of their Gelbrich distance from $\Sigma$ and we compare it to the lower bound of the absolute Pearson correlation coefficient given by (17). Sampled bivariate covariance matrix are represented in blue dots and the red dots stand for the correlation absolute lower bound given by the equation (17). Notice that this lower bound seems relevant as all the sampled covariance matrices that have a Gelbrich distance less or equal to a red dot lie above the red dots, i.e., they have a greater or equal absolute Pearson correlation coefficient. When the Gelbrich radius is greater than $D_U = 0.618$, the minimal absolute correlation is equal to 0. Denote by $\mu_S$ and $S = \begin{bmatrix} S_{xx} & S_{xy} \\ S_{yx} & S_{yy} \end{bmatrix}$ the mean and the variance of the distribution $Q_S$. From (5), we have that the variance of the error of any linear estimator of $y$ in terms of $x$ is larger or equal to $S_{yy}(1 - \rho^2(Q_S))$. When $Q_S$ equals the nominal distribution $Q_\Sigma$, and the estimator is $\hat{y} = \mu_y + \frac{S_{xy}}{S_{xx}}(x - \mu_x)$, then the variance of the error is $S_{yy}(1 - \rho^2(Q_\Sigma)) = S_{yy}(1 - \rho_0^2) = 0.5$. However, if $\alpha \ge D_U = 0.618$, then there is at least an element $Q_S \in A_\alpha^G(Q_\Sigma)$ such that the Pearson correlation is zero and thus the corresponding error variance is larger or equal than $S_{yy} = 1$, indicating that the use of $x$ as a covariate to estimate $y$ may not be advisable.

<figure><img src="image_3.png" /><figcaption>Fig. 2. Illustrating the minimal (absolute) correlation in a Gelbrich ball.</figcaption></figure>

## V. CONCLUSION AND PERSPECTIVES

This paper provides new developments on the Wasserstein/Gelbrich distance and the characterization of their bivariate distribution ambiguity sets, easing the task to analyze correlations between scalar random variables and facilitating robust decision making. Through closed expressions, the paper establishes the minimum Gelbrich perturbation needed

to reach a specified correlation between two random scalar variables. Additionally, it allows to compute the minimum absolute Pearson correlation coefficient in a Gelbrich ambiguity set. The obtained results have potential application in the field of variable selection. Numerical examples validate the two proposed lower bounds closed expressions. Ongoing efforts are concentrated on expanding the obtained results to encompass new applications in estimation and control.

## REFERENCES

[1] J. Lee Rodgers and W. A. Nicewander, "Thirteen ways to look at the correlation coefficient," *The American Statistician*, vol. 42, no. 1, pp. 59-66, 1988.

[2] E. Saccenti and M. H.W.B. Hendriks and A. K. Smilde, "Corruption of the Pearson correlation coefficient by measurement error and its estimation, bias, and correction under different error models," *Scientific Reports*, vol. 10, pp. 1-19, 2020.

[3] D. Simon, "The discrete-time Kalman filter," in *Optimal State Estimation: Kalman, H∞, and Nonlinear Approaches*, ch. 5, pp. 121-148, John Wiley & Sons, Ltd, 2006.

[4] J. M. Nadales, A. Hakobyan, D. Muñoz de la Peña, D. Limon, and I. Yang, "Risk-aware Wasserstein distributionally robust control of vessels in natural waterways," *IEEE Transactions on Control Systems Technology*, pp. 1-0, 2024.

[5] P. Esfahani and D. Kuhn, "Data-driven distributionally robust optimization using the Wasserstein metric: Performance guarantees and tractable reformulations," *Mathematical Programming*, vol. 171, 05 2015.

[6] D. Kuhn, P. M. Esfahani, V. A. Nguyen, and S. Shafieezadeh-Abadeh, "Wasserstein distributionally robust optimization: Theory and applications in machine learning," in *Operations Research & Management Science in the Age of Analytics*, pp. 130-166, INFORMS, 2019.

[7] Y. Cai and L.-H. Lim, "Distances between probability distributions of different dimensions," *IEEE Transactions On Information Theory*, vol. 68, no. 6, pp. 4020-4031, 2022.

[8] M. Gelbrich, "On a formula for the L2 Wasserstein metric between measures on Euclidean and Hilbert spaces," *Mathematische Nachrichten*, vol. 147, no. 1, pp. 185-203, 1990.

[9] V. A. Nguyen, S. Shafieezadeh-Abadeh, D. Kuhn, and P. Mohajerin Esfahani, "Bridging Bayesian and minimax mean square error estimation via Wasserstein distributionally robust optimization," *Mathematics of Operations Research*, vol. 48, no. 1, pp. 1-37, 2023.

[10] C. R. Givens and R. M. Shortt, "A class of Wasserstein metrics for probability distributions.," *Michigan Mathematical Journal*, vol. 31, no. 2, pp. 231-240, 1984.

[11] M. Arjovsky, S. Chintala, and L. Bottou, "Wasserstein generative adversarial networks," in *Proceedings of the 34th International Conference on Machine Learning* (D. Precup and Y. W. Teh, eds.), vol. 70, pp. 214-223, PMLR, 2017.

[12] R. Chen and I. C. Paschalidis, "A distributionally robust optimization approach for outlier detection," in *57th IEEE Conference on Decision and Control*, pp. 352-357, 2018.

[13] S. Shafieezadeh Abadeh, V. A. Nguyen, D. Kuhn, and P. M. Mohajerin Esfahani, "Wasserstein distributionally robust Kalman filtering," *Advances in Neural Information Processing Systems*, vol. 31, 2018.

[14] D.-J. Xin, L.-F. Shi, and X. Yu, "Distributed Kalman filter with faulty/reliable sensors based on Wasserstein average consensus," *IEEE Transactions on Circuits and Systems II: Express Briefs*, vol. 69, no. 4, pp. 2371-2375, 2022.

[15] F. A. Taha, S. Yan, and E. Bitar, "A distributionally robust approach to regret optimal control using the Wasserstein distance," in *62nd IEEE Conference on Decision and Control*, pp. 2768-2775, 2023.

[16] L. Aolaritei, M. Fochesato, J. Lygeros, and F. Dörfler, "Wasserstein Tube MPC with exact uncertainty propagation," 2023.

[17] M. Borelle, T. Alamo, C. Stoica, S. Bertrand, and E. F. Camacho, "Minimal Gelbrich distance to uncorrelation," *IEEE Control Systems Letters*, vol. 8, pp. 61-66, 2023.

[18] J. Cohen, P. Cohen, S. G. West, and L. S. Aiken, *Applied Multiple Regression/Correlation Analysis for the Behavioral Sciences*. Lawrence Erlbaum Associates, Inc., Third ed., 2003.

[19] G. Birkhoff and S. Mac Lane, *A survey of modern algebra*. A.K. Peters, Ltd., 1997.