Quentin Giboulot

Member, IEEE Rémi Cogranne

Multivariate Side-Informed Gaussian Embedding Minimizing Statistical Detectability

Steganography schemes based on a deflection criterion for embedding posses a clear advantage against schemes based on heuristics as they provide a direct link between theoretical detectability and empirical performance. However, this advantage depends on the accuracy of the cover and stego model underlying the embedding scheme. In this work we propose an original steganography scheme based on a realistic model of sensor noise, taking into account the camera model, the ISO setting and the processing pipeline. Exploiting this statistical model allows us to take correlations between DCT coefficients into account. Several types of dependency models are presented, including a very general lattice model which accurately models dependencies introduced by a large class of processing pipelines of interest. We show in particular that the stego signal which minimizes the KL divergence under this model has a covariance proportional to the cover noise covariance. The resulting embedding scheme achieves state-of-the-art performances which go well beyond the current standards in side-informed JPEG steganography.

I. INTRODUCTION

One of the stated goals of steganography is to design algorithms able to hide information in an innocuous medium, henceforth named a cover object. As of today, steganography as a discipline has concentrated its effort on using digital media as covers, particularly digital images which will be the focus of this work. Modern steganography is based on two main ingredients: a cost function and a coding scheme. On the one hand, coding schemes improve steganography performance by hiding more data with fewer changes. As of today, it is considered as a mostly solved problem since the Syndrome-Treillis Coding [START_REF] Filler | Minimizing Additive Distortion in Steganography Using Syndrome-Trellis Codes[END_REF] scheme is very close to Shannon efficiency. On the other hand, cost functions associate to each cover element a cost of modification. The goal of the steganographer is then to minimize the overall cost under the constraint of hiding messages of a given payload size. This allows for steganography to be adaptive; the underlying heuristic being that hiding in a smooth region should be more costly than hiding in a textured region because the former should be intuitively more detectable than the latter. It should be noted that cost functions currently used in the literature are almost always additive. This means that the cost of modifying one element does not depend on modifications on other elements. This key assumption allows for tremendous

The work presented in this paper was founded by the European Union's Horizon 2020 project UNCOVER (under grant 101021687), the ANR/DGA Project ALASKA (under grant ANR-18-ASTR-0009) and the ANR Project PACeS (grant ANR-21-CE39-0002).

Quentin Giboulot and Rémi Cogranne are with the LIST3N research unit, Troyes University of Technology. Patrick Bas is with the CNRS, École Centrale de Lille. simplifications when designing steganographic algorithms [START_REF] Filler | Minimizing Additive Distortion in Steganography Using Syndrome-Trellis Codes[END_REF]. However, we need to forego this assumption in this work to be able to leverage correlations between DCT coefficients.

The other side of the game -detecting objects in which information has been hidden -is named steganalysis. As of today, state-of-the-art methods are mainly based on supervised machine learning techniques. Until recently, the favored approach was the use of handcrafted, high-dimensional features specifically designed to capture artifacts introduced by steganographic schemes. These features were then fed to classifiers designed to handle such high-dimensional feature sets, the most popular being an ensemble of FLD classifiers [START_REF] Kodovsky | Ensemble Classifiers for Steganalysis of Digital Media[END_REF] and a fast ridge linear classifier [START_REF] Cogranne | Is ensemble classifier needed for steganalysis in high-dimensional feature spaces?[END_REF]. However, these techniques are now being largely superseded by deep neural networks. Those were at first specifically designed for steganography [START_REF] Xu | Structural Design of Convolutional Neural Networks for Steganalysis[END_REF], [START_REF] Boroumand | Deep Residual Network for Steganalysis of Digital Images[END_REF], but recent advances (see the results of the ALASKA2 competition [START_REF] Cogranne | ALASKA#2: Challenging Academic Research on Steganalysis with Realistic Images[END_REF]) showed that using neural networks pre-trained on ImageNet such as Efficient-Net [START_REF] Tan | EfficientNet: Rethinking Model Scaling for Convolutional Neural Networks[END_REF] can lead to similar or better performance than these specialized networks.

A. State of the art

As was alluded to earlier, steganography is now mainly concerned with the design of effective cost functions. There exists currently two main paths for their design: heuristic and statistical. The heuristic path has been by far the most popular. The approach is based on designing cost functions able to defeat the most effective steganalyzer available on a given standard dataset. The performance of these schemes is consequently only empirically validated. This approach gave rise to the most successful schemes both in the spatial domain [START_REF] Holub | Universal distortion function for steganography in an arbitrary domain[END_REF], [START_REF] Li | A new cost function for spatial image steganography[END_REF] and in the JPEG domain [START_REF] Holub | Universal distortion function for steganography in an arbitrary domain[END_REF], [START_REF] Guo | Using Statistical Image Model for JPEG Steganography: Uniform Embedding Revisited[END_REF]. Despite this success, this approach has several limitations. First, because these techniques are only validated empirically, one can observe significant differences in performance depending on how the steganalysis is performed [START_REF] Sedighi | Toss that BOSSbase, Alice![END_REF], [START_REF] Giboulot | Effects and solutions of Cover-Source Mismatch in image steganalysis[END_REF]. Secondly, and more generally, this approach uses cost functions that do not apriori have a clear link with theoretical or even empirical detectability. Consequently, this approach gives very little insight on why a strategy works and why another does not. Worse, it makes the approach incapable of giving any theoretical guarantees except under the precise setting for which it has been designed. The second path is based on minimizing a quantity which should be apriori linked to theoretical and/or empirical detectability. The most recent strategy following this path, adversarial embedding, directly tries to minimize empirical detectability on a given dataset by iteratively modifying the cost function to defeat a classifier which also updates at each step to defeat the steganographer. Though this approach has been proven to be quite successful [START_REF] Tang | CNN-Based Adversarial Embedding for Image Steganography[END_REF], [START_REF] Bernard | Explicit Optimization of min max Steganographic Game[END_REF] it has limitations of its own, namely its high computational cost and the need for large datasets, necessary to obtain a high-quality cost function. Furthermore, it does not solve the problem of theoretical guarantees for unseen datasets and there has not been, to the best of our knowledge, a systematic study of the impact cover-source mismatch for this strategy. For these reasons, we work with another strategy in this paper which works by bounding the theoretical power of an optimal detector. This approach first appeared in [START_REF] Fridrich | Multivariate gaussian model for designing additive distortion for steganography[END_REF] and culminated in the design of the MiPOD algorithm [START_REF] Sedighi | Content-Adaptive Steganography by Minimizing Statistical Detectability[END_REF]. The underlying framework is based on hypothesis testing theory. The idea is to cast the steganalysis problem as a simple test between two hypotheses: H 0 , the image under scrutiny is a cover, or H 1 , it is a stego. Using this theory, one can show under some conditions that an optimal detector exists and analytically compute its statistical performance. The goal of the steganographer is then to minimize the power of this detector under the constraint of embedding a payload of a given size. However, to be able to cast the problem in this setting, one must have a model of both the cover and the stego images. The limits of MiPOD and the MG algorithm which preceded it mostly came from their choice of the noise model. Indeed both of these algorithms modeled natural images as pixels corrupted by an independent, though not identical, Gaussian noise. Furthermore, it relied on the so-called fine quantization limit assumption, which states that the variances of the pixels are greater than one. These two assumptions are both erroneous in practice. First, it is known that the neighboring pixels are, more often than not, correlated due to the impact of the processing pipeline [START_REF] Taburet | Natural Steganography in JPEG Domain With a Linear Development Pipeline[END_REF], [START_REF] Giboulot | Detectability-Based JPEG Steganography Modeling the Processing Pipeline: The Noise-Content Trade-off[END_REF], [START_REF]Synchronization Minimizing Statistical Detectability for Side-Informed JPEG Steganography[END_REF]. Secondly, the fine quantization limit is often violated, especially in dark areas of an image due to the heteroscedastic nature of the noise [START_REF] Foi | Clipped noisy images: Heteroskedastic modeling and practical denoising[END_REF], [21, Chapter 5, Section 5] which implies that the variance in these zones is quite small. This assumption is even more problematic if one wants to extend the approach in the JPEG domain as the quantization of DCT coefficients can greatly lower the variance before rounding.

Despite these limitations, MiPOD still enjoys close to stateof-the-art performances in the spatial domain, hence demonstrating the merit of the approach.

B. Contributions and comparisons to current approaches

The work presented in this paper is the logical continuation of our previous works on Gaussian Embedding [START_REF] Giboulot | Detectability-Based JPEG Steganography Modeling the Processing Pipeline: The Noise-Content Trade-off[END_REF], [START_REF]Synchronization Minimizing Statistical Detectability for Side-Informed JPEG Steganography[END_REF]. This series of work has its roots in [START_REF] Giboulot | Effects and solutions of Cover-Source Mismatch in image steganalysis[END_REF] where it was observed that empirical steganalysis performance is mainly determined by only three factors: sensor, ISO, and processing pipeline. This motivated the construction of a statistical model of the noise that was only based on these three parameters. The importance of such a model was demonstrated by the success of the work on Natural Steganography culminating in [START_REF] Taburet | Natural Steganography in JPEG Domain With a Linear Development Pipeline[END_REF]. In particular, this work highlighted the importance of the covariance matrix of the cover noise for the security of a steganographic scheme. However, the idea behind Natural Steganography is to imitate the noise of an image if this image was taken at a higher ISO. Its methodology thus relies on cover generation whereas our approach is based on minimizing the statistical detectability of a given stego signal to provide security guarantees to the steganographer. To leverage the covariance of the noise, [START_REF] Giboulot | JPEG Steganography with Side Information from the Processing Pipeline[END_REF] builds a model of the sensor noise in the JPEG domain using a very general linear model of the processing pipeline. This allows improvements on MiPOD by allowing a better estimation of the variance maps. This work also innovated by modeling the stego signal in the continuous domain (DCT domain before rounding) as a Gaussian random variable while specifying the payload constraint in the discrete domain (DCT domain after rounding). This allows foregoing the fine quantization limit assumption entirely. Despite these improvements, this algorithm only used the covariance matrix to improve its cost estimation. It does not use a multivariate Gaussian signal which, as has been shown by Natural Steganography, has a tremendous impact on security. Furthermore, recent heuristic approaches [START_REF] Denemark | Improving Steganographic Security by Synchronizing the Selection Channel[END_REF], [START_REF] Li | A new cost function for spatial image steganography[END_REF], [START_REF] Hu | Model-based image steganography using asymmetric embedding scheme[END_REF] have shown that taking into account embedding modifications performed on neighboring cover elements can lead to a significant increase in security. However, until our work, there was no statistically founded approach to explain how to take neighboring modifications into account for imperfect steganography. The subsequent work [START_REF]Synchronization Minimizing Statistical Detectability for Side-Informed JPEG Steganography[END_REF] improved on [START_REF] Giboulot | Detectability-Based JPEG Steganography Modeling the Processing Pipeline: The Noise-Content Trade-off[END_REF] by allowing the use of a multivariate Gaussian stego signal. This is achieved by using the Cholesky decomposition of the covariance matrix which allows recasting the difficult multivariate problem as a sequence of simple univariate ones. This paper is an extension of our paper [START_REF]Synchronization Minimizing Statistical Detectability for Side-Informed JPEG Steganography[END_REF] and solves many of the limitations of this works, namely:

• Our previous works assumed that the steganographer has access to the RAW image to estimate the covariance matrix of the noise. We relax this assumption greatly by requiring only the knowledge of the camera and ISO used -information usually available in the EXIF metadataas well as black-box access to the processing pipeline. • Our previous work assumed macro-blocks of DCT coefficients to be independent. This means that dependencies between neighboring blocks not belonging to the same macro-block were not taken into account. We solve this limitation in this work using a lattice embedding strategy akin to what was used in Natural Steganography [START_REF] Taburet | Natural Steganography in JPEG Domain With a Linear Development Pipeline[END_REF].

• We extended our model to take saturation of pixels into account. • Last but not least, numerous new experimental results, including results using recent steganalysis tools such as Efficient-Net, are now given. Full proofs of all mathematical results are also now available in the appendices.

C. Notation and terminology

When referring to blocks of elements in an image, we always refer to 8 × 8 blocks of elements. When referring to "blocks of blocks" we use the term macro-block. Except when explicitly stated, we manipulate blocks in their vectorized representation arranged lexicographically. Vectors and matrices are always written in boldface. Vectors use small case letters while matrices use uppercase, except for collections of vectors which are denoted by bold small case letters without indices. When indexing individual elements of an image we use the letter i. When indexing vectors of elements of an image we use the letter k. Similarly, when referring to the number of individual elements in an image we use the letter n, whereas we use the letter m when referring to the number of blocks or macro-blocks of a specified size. For the sake of clarity, we use the same symbol for a random variable and its realization. The Gaussian distribution is designated with N . We do not differentiate the symbol between the univariate and multivariate cases as it can be inferred from the typesetting of the parameters. When referring to the distribution which is obtained by quantizing a Gaussian random variable with a uniform quantizer of step 1, we use the symbol N. We refer to the diagonal matrix constructed with individual elements from a vector σ as diag(σ). Similarly, we refer to the vector constructed with the diagonal of a matrix Σ as diag(Σ).

II. STATISTICAL MODEL OF THE NOISE IN THE DEVELOPED

DOMAIN

In this section, we derive a model of the noise of a developed image starting from the RAW image. We summarize each step of the model derivation in Figure 1. Note that we model the noise only up to the DCT transform of the JPEG compression but not including the rounding operation as we will mostly work with a model in the continuous domain throughout this paper.

A. Model of the RAW image

We begin by giving a model of the noise of the image in the RAW domain, that is, before any processing is applied. A widely adopted model of the sensor noise is the heteroscedastic Gaussian noise model studied by Foi et al. [START_REF] Foi | Clipped noisy images: Heteroskedastic modeling and practical denoising[END_REF]. Under this model, a RAW image is composed of n photo-sites where each photo-site x i follows a Gaussian distribution :

x i ∼ N µ i , σ 2 i , σ 2 i = c 1 µ i + c 2 , (1)
where µ i is the value the sensor would have registered in the absence of noise. The variance of the noise, σ 2 , depends linearly on µ i through two parameters c 1 and c 2 . These two parameters depend on the camera sensor and the ISO setting -see [START_REF] Foi | Practical Poissonian-Gaussian Noise Modeling and Fitting for Single-Image Raw-Data[END_REF], [START_REF] Thai | Statistical Model of Quantized DCT Coefficients: Application in the Steganalysis of Jsteg Algorithm[END_REF]. At his point, note that the noise is considered independent between photo-sites.

B. Model of the processing pipeline

We now go on to model the processing pipeline which takes a RAW image x as input and outputs a developed image y (in the continuous domain). Many operations in image processing can be modeled as convolutions (e.g demosaicking, sub-sampling, DCT transform) or sample-wise linear functions (e.g white balance, RGB to greyscale conversion). Even though non-linearities are also present, such as in the gamma correction, in order to keep our model computationally tractable, we model the processing pipeline as a stationary linear map acting on vectors of photo-sites. By stationary, we mean that the linear map is the same for every vector of elements of the input image. Formally, we represent the linear map as a full-rank matrix H ∈ R N ×M where M and N are freely chosen by the modeler depending on trade-off between the model accuracy and computational complexity: larger M and N allows for a model of more far-reaching dependencies between DCT coefficients but will be more computationally expensive due to the increased size of the processing pipeline matrix. Note that some values are optimal for certain kinds of pipelines, for example, the work in [START_REF] Taburet | Natural Steganography in JPEG Domain With a Linear Development Pipeline[END_REF] shows that choosing M = 26 2 and N = 24 2 is sufficient to model all the dependencies introduced by a bilinear demosaicking followed by a DCT transform. Finally, we require N to be a perfect square integer which is also a multiple of 8 since the DCT transform acts on 8 × 8 blocks of pixels.

C. Model of the developed image

To apply the processing pipeline matrix H as modeled in the previous subsection, we first write each macro-block of the RAW image as x k :

x k ∼ N (µ k , diag(σ k )) . (2)
By multiplying the matrix H with a macro-block x k containing M elements, we obtain a new developed macro-block y k which follows a multivariate Gaussian (MVG) distribution:

y k = Hx k , (3)
∼ N (H µ k , Σ k ) , (4)
with the covariance of y k simply obtained as:

Σ k = Hdiag(σ k )H T . (5)
Depending on the processing pipeline, the dependency structure between the y k might differ. In this paper we will only treat two models of dependency: 1) Independent model : Macro blocks of size

√ N × √ N are considered independent.
2) Lattice model : Each block is considered to be dependent on its neighboring block -including the diagonal ones. We illustrate each type of dependency model in Figure 2. . Finally, note that contrary to the independent model, which is heuristic in nature, the lattice model has been extensively studied and justified in the work of Natural Steganography [START_REF] Taburet | Natural Steganography in JPEG Domain With a Linear Development Pipeline[END_REF] as the optimal dependency model for certain types of linear pipeline.

III. ESTIMATION AND APPROXIMATION OF THE COVARIANCE MATRIX WITHOUT THE RAW FILE

Our previous work [START_REF] Giboulot | Detectability-Based JPEG Steganography Modeling the Processing Pipeline: The Noise-Content Trade-off[END_REF] gave a methodology to estimate the covariance matrix when the RAW file and processing pipeline are available. This section provides a novel method to estimate the covariance matrix without having access to the RAW file. Note that this section is an option However, we still assume access to the processing pipeline, at least as a black box, since the estimation method mostly relies on having a linear approximation of this pipeline. We will also assume that the steganographer knows the c 1 and c 2 parameters of the heteroscedastic model of the cover. This knowledge does not

              h 11 h 1M h N 1 h N M              

Lossy compression

Unprocessed image RAW image : x

x i ∼ N (µ i , c 1 µ i + c 2 ) x k ∼ N (µ k , diag (c 1 µ k + c 2 ))
Processing pipeline : H Ex : Demosaicking, White balance, DCT transform, . . .

Developed image

Precover : y

y k ∼ N (Hµ k , Σ k ) Compressed image Cover : z z k ∼ N (Hµ k , Σ k )
Fig. 1: Summary of the model of natural images that is described in Section II.

(a) Independent macroblock model: every nonoverlapping macro-blocks of a given size is considered to be independent with the others. In this figure,the chosen size is 24 × 24, hence every block is only dependent with blocks with which it shares its color. (b) Lattice model: every block is considered dependent only with its direct neighboring blocks. This induces four sets of blocks, called lattices, Λ1, Λ2, Λ3, Λ4 where each block in a given set is independent of every other block in the same set. Fig. 2: The two dependency models studied in this paper necessarily require access to the RAW image, it is sufficient to know the camera and ISO which were used to capture the cover; this information is usually available in the EXIF metadata of an image.

The estimation method presented in this section is based on estimating a generic correlation matrix which depends only on the processing pipeline. The idea is then to compute the variance of every DCT coefficient using an approximate heteroscedastic model of the DC coefficients before computing an estimation of the covariance of each block by scaling the correlation matrix using these estimated variances.

A. Heteroscedastic model of the DC coefficients

Our goal here is to show that, under some additional assumptions on the processing pipeline, the variance of the DC coefficients of each block is linear with respect to the expectation of this DC coefficient. Note that, for this section only, we distinguish between H DCT , the 64 × 64 matrix representing the DCT transform and H s the 64 × M matrix representing all operations performed in the RAW and spatial domain. We therefore have:

H = H DCT H s .
We assume, as usual, that the processing pipeline is both linear and stationary. Furthermore, we assume the RAW image is almost constant by block, that is, for all k we have

µ k = μ + e k with |e k,i | ≪ μ, 1 ≤ i ≤ M .
First of all, let us rewrite the models of the block of photosites and of DCT coefficients:

x k ∼ N (µ k , c 1 diag (µ k ) + c 2 ) , (6)
y k ∼ N H DCT (H s µ k -128), Hdiag (c 1 µ k +c 2 ) H T . (7)
Note here that, contrary to Eq 3, we take into account the fact that we subtract 128 to each pixel before the DCT transform as this has an impact on the estimation method of this section (whereas it does not when using the method which uses RAW file). We can express the first two moments of the DC coefficient y k,1 of each block:

E [y k,1 ] = M l=1 H 1,l (μ + e k,l ) -1024, (8)
Var [y k,1 ] = c 1 M l=1 H 2 1,l μ + H 2 1,l c 2 + c 1 M l=1 H 2 1,l e k,l . (9)
Let us note the following quantities:

Hi,j ≜ M l=1 H 2 i,l M l=1 H j,l , H(2) i,j ≜ M l=1 H 2 i,l M l=1 H 2 j,l . ( 10
)
The variance of y k,1 can be expressed using its expectation:

Var [y k,1 ] = c 1 H1,1 E [y k,1 ] + M l=1 H 2 1,l c 2 + 1024 • c 1 H1,1 + c 1 M l=1 H 2 1,l e k,l -H1,1 M l=1 H 1,l e k,l ≜ c DC 1 E [y k,1 ] + c DC 2 + error, (11)
with c DC 1 and c DC 2 defined as:

c DC 1 = c 1 H1,1 , (12)
c DC 2 = M l=1 H 2 1,l c 2 + 1024 • c 1 H1,1 . (13)
This shows that the variance of the DC coefficients is linear with their expectation up to an error which is small as long as the ratio H1,1 is small. As a particular case, note that if all the H 1,l are constant, for example if the only operation of the processing pipeline performed is the DCT transform, then the error is simply 0. In this paper, we assume that the error is negligible in practice.

Note however, that the H i,l usually alternate sign for i > 1 because of the structure of the DCT transform. This leads to the ratio Hi,l possibly exploding and the error can not be considered small anymore for AC coefficients in general. However we can still express the variance of the AC coefficients as a function of the variance of the DC coefficient:

Var [y k,i ] = c 1 μ M l=1 H 2 i,l + c 2 M l=1 H 2 i,l + c 1 M l=1 H 2 i,l e k,l ≃ H(2) i,1 Var [y k,1 ] , (14)
using the fact that |e k,l | is small compared to μ.

B. Processing pipeline and correlation matrix estimations

The first step is to estimate the processing pipeline matrix H since we only assumed a black-box access to the processing pipeline.

To do so, we first generate an image x so that:

xi ∼ N 0, σ 2 , ( 15
)
where the value of σ 2 can be freely chosen and does not impact the estimation.

Using the black-box access to the pipelines, we develop this image to obtain the developed image ȳ.

Using Eq (3), we know that for a given macro-block size M , the processing pipeline outputs a macro-block of size N and that the k-th macro-block of the image follows:

ȳk = Hx k , (16)
which can be solved using any type of linear regression method. For example, we can solve for H using a least-square estimation:

H = ȳx T xx T -1 . ( 17
)
where x and ȳ here correspond respectively to the M ×m and N × m matrices where each line corresponds to a vectorized macro-block. Note that, here, H is the closest solution (in the least-square sense) to the set of equations in Eq [START_REF] Sedighi | Content-Adaptive Steganography by Minimizing Statistical Detectability[END_REF] and not an exact solution.

To approximate the covariance matrices, we will also use a correlation matrix ρ based on the processing pipeline. First we compute the sample covariance matrix of ȳ:

Σ ȳ = 1 m -1 ȳk ȳT k , (18)
and obtain the correlation matrix:

ρ = diag ϱ 2 ȳ -1 2 Σ ȳdiag ϱ 2 ȳ -1 2 . ( 19
)
where

ϱ 2 ȳ ≜ diag (Σ ȳ).

C. Approximation of the covariance matrix

With the use of the heteroscedastic parameters c DC 1 and c DC 2 , it is now possible to compute an approximation of the covariance matrix using the correlation matrix ρ computed in Section III-B. In this subsection, we denote Var [y k,i ] as ϱ 2 k,i . We first compute an approximation of the true variance map of the DCT coefficients and use it to scale the correlation matrix:

1) Compute the variance ϱ 2 k,1 of the DC coefficient of the k-th block y k as:

ϱ 2 k,1 = c DC 1 y k,1 + c DC 2 . ( 20
)
2) Compute the variances of the l-th coefficients of the k-th block simply as:

ϱ 2 k,i = H(2) i,1 ϱ 2 k,1 . (21)
. 3) Compute the covariance matrix of the k-th block by scaling the correlation matrix ρ with the variances of the DCT coefficients of the block using the standard formula:

Σk = diag(ϱ k ) ρ diag(ϱ k ), (22)
where ϱ k is the vector of standard deviation of the k-th block of DCT coefficients and Σk the resulting estimation of the covariance matrix.

IV. OPTIMAL DETECTOR The goal of the steganographer is to find the optimal prestego signal to use in order to evade the steganalyst. By prestego, we refer to the stego signal in the continuous domain, which will be the only domain of interest in this section. Section V addresses the transition to the quantized domain in detail.

To find the form of an optimal prestego signal, three ingredients are needed:

1) A model of both the precover, in our case, the cover after the DCT transform but before rounding -see y k in Figure 1 and the prestego, 2) A criterion of optimality, 3) An optimal detector built from the two first items. Once these ingredients are available, the optimal prestego signal is the signal which minimizes the power of the optimal detector under a given payload size constraint.

In this section we provide these three elements, beginning with a choice of prestego model. The optimal detector is then derived using the precover and the prestego model under the Neyman-Pearson of optimality. Simple hypotheses

H 0 = {ξ k ∼ N (µ k , Σ k )} H 1 = {ξ k ∼ N (µ k , Σ k + E k )}
Likelihood Ratio Test Fixed false alarm rate : P F A Maximizes power :

P D Λ(ξ, Σ, E) = qΣ,E(ξ) pΣ(ξ) ≶ H0 H1 τ
Fig. 3: Summary of the construction of the Likelihood Ratio Test for two simple hypotheses H 0 , the image is cover, H 1 , the image is stego.

A. Prestego model

Despite the fact that the precover model has been fully defined in Section II, the steganographer still needs a model of the prestego signal.

We will assume for the rest of the paper that the steganographer uses a centered multivariate Gaussian signal:

s k ∼ N (0, E k ) , (23)
where E k is the covariance of the prestego signal s k . The signal s k is then simply added to the k-th macro-block of the precover, creating the prestego macro-block γ k :

γ k = y k + s k , ∼ N (µ k , Σ k + E k ) . (24)
This choice of distribution, though we do not prove its optimality, is motivated by three nice properties of this distribution. First, the sum of two Gaussian random variables is also a Gaussian random variable which obviously facilitates the derivations in this paper. Secondly, for a given expectation and variance, the Gaussian distribution is the maximum entropy distribution in the continuous domain [START_REF] Cover | Elements of information theory[END_REF], we can thus expect the distribution to maximize the embedded payload for a given KL-divergence. Finally, using a multivariate Gaussian distribution allows us to construct extremely efficient algorithms to compute the stego signal is the quantized domain -see Section V.

B. Optimal detector

Now that both the precover and prestego model are defined, it is possible to construct the optimal detector. To do so, we use the Neyman-Pearson criterion of optimality [28, Chapter 3, Section 2]. In this setting the steganalyst constructs a test δ : R n → {H 0 , H 1 } which maximizes the power of all possible tests P D ≜ P (δ (x) = H 1 |H 1 ) whose false-alarm probability

P F A ≜ P (δ (x) = H 1 |H 0 ) is upper bounded by a chosen α 0 .
We also assume a worst-case adversary for the steganographer and as such consider that the steganalyst has access to all the model parameters : E = (E 1 , E 2 , . . . , E m ), µ = (µ 1 , µ 2 , . . . , µ m ) and Σ = (Σ 1 , Σ 2 , . . . , Σ m ). Finally, the image under scrutiny by the steganalyst, which can be either a precover or a prestego, is written as ξ.

The problem of the steganalyst reduces to a choice between two hypotheses: H 0 , the image under consideration is cover, H 1 , the image under consideration is stego:

H 0 = {ξ k ∼ N (µ k , Σ k )} , H 1 = {ξ k ∼ N (µ k , Σ k + E k )} . ( 25
)
Because under our model the two hypotheses are simple, the Neyman-Pearson Lemma states that the most-powerful test is the likelihood ratio test (LRT) which maximizes the power P D for a given P F A ≤ α 0 , defined, in our case as follows:

Λ(ξ, Σ, E) = q Σ,E (ξ) p Σ (ξ) H0 ≶ H1 τ, (26)
where p Σ and q Σ,E are the joint pdf of all macro-blocks of the precover and of the prestego respectively. The threshold τ is fixed in advance by the steganalyst depending on her choice of

P F A = P (δ(x) > τ | H 0 ) = α 0 .
The study of the performance of the LRT is detailed in Appendix A and B.

V. DESIGN OF A SIDE-INFORMED MULTIVARIATE STEGANOGRAPHIC SCHEME

A. Problem specification and optimal stego signal At this point of the paper, the steganographer has access to a model of both the precover and of the prestego as well as to the detector of the worst-case adversary. His goal is then to design a steganographic scheme which minimizes the power of this detector while still hiding a secret message of a given size in the cover.

Until now, all of our models have been specified in the continuous domain. However, the payload has to be embedded in the cover in the quantized domain. We make the assumption that minimizing the power of the LRT in the continuous domain also minimizes the power of the optimal detector in the quantized domain. As such we specify the problem of the steganographer as minimizing the power P D of the optimal detector in the continuous domain under a payload constraint R in the quantized domain:

       min E P D (E) R = n i=0 j∈Z β (j) i log β (j) i (27)
where

β (j) i
is the probability of modifying the i-th DCT coefficient by +j. Note that the payload constraint is written under the assumption that we embed at the optimal rate. Pre-cover P Covariances Σ k

Sample white Gaussian noise w

Compute Cholesky decomposition

Σ k = L k L T k -Eq (31)
Compute conditionnal Gaussians parameters ηi ,ε i -Eq (33 ,34)

Compute prestego signal s k = L k w k -Eq [START_REF] Van Trees | Detection estimation and modulation theory[END_REF] Scale signal and parameters by α Compute embedding probabilities β

(j) i -Eq 36 Compute entropy H(β). Compare to payload constraint R Stego image If H(β) = R If H(β) ̸ = R
Fig. 4: Principle of all variants of Multivariate Gaussian Embedding (MGE) presented in this paper. The main idea is to generate a multivariate Gaussian stego signal with covariance proportional to the covariance of the cover noise. The Cholesky decomposition of the cover noise covariance matrix is used to generate the stego signal as well as to compute the entropy of the signal.

Due to the difficulty of obtaining an exact expression of the power of the LRT in our setting -see Appendix Awe simplify this optimization problem by minimizing the KLdivergence between the precover and prestego instead:

       min E D KL (p Σ |q Σ,E ) R = n i=0 j∈Z β (j) i log β (j) i ( 28
)
This simplification is justified by a data-processing inequality stating that the KL-divergence between the precover and the prestego provides an upper bound on the power of the LRT [START_REF] Cachin | An Information-Theoretic Model for Steganography[END_REF]. Consequently, we sacrifice the guarantee of optimality while conserving a guarantee of security. We now use the key result of this paper -proven in Appendix Bwhich states that the prestego signal which minimizes the KL-divergence between the precover and the prestego under a given entropy in the continuous domain has the following form:

s k ∼ N (0, αΣ k ) , ∀k (29)
with α ∈ R + . In other words, the prestego signal has a covariance proportional to the covariance of the precover noise. In particular, α is the same for all macro-blocks of the precover. However, the reader should be aware that the form of the matrix in Eq (29) has only been shown to be optimal -in the KL-divergence sense -for an entropy constraint in the continuous domain. We have to assume this result to translate for an entropy constraint in the quantized domain. This is not always guaranteed. Indeed, in the case where we have αε i < 0.5 for some i (see Eq (33) for a definition of εi ), the expression of the entropy is a complicated function of the rounding errors r i , εi and ηi (again see Eq (34) for a definition of ηi ) and might not be strictly increasing.

In practice, these cases are rare enough in a single image, for payloads of interest, that we have not observed their impact. Therefore, we assume Eq [START_REF] Cachin | An Information-Theoretic Model for Steganography[END_REF] to hold for payload constraint in the quantized domain and, in particular, that

n i=0 j∈Z β (j) i log β (j) i
is strictly increasing in α. Notice that, in the particular case where all the εi are large enough, the entropy is well approximated by

n i=1 1 2 log 2πe εi + 1 12 2 .
Following this assumption, the system in Eq ( 28) is easily solved by a simple bisection search on α.

However, to be able to compute the value of the payload size for a given α it is necessary to compute all individual β (j) i which is not an obvious task since we want to compute an individual value in the quantized domain from a multivariate signal in the continuous domain. The next subsection addresses this difficulty and presents a solution that does not rely on expensive Monte-Carlo simulations.

Note however, that for the rest on this paper, we only discuss methods related to simulating the pre-stego signal as simply and effectively as possible. These methods can be adapted to practical uses with a multi-layered STC by using the rejection sampling method used in [17, Section V.C].

B. Computing embedding probabilities

First, let s k be a pre-stego signal macro-block. It follows a centered multivariate Gaussian (MVG) random variable (rv) of N elements with a full-rank covariance E k . Denote the i-th element of s k as s k,i . Then we have, for all i:

s k,i |s k,1 , s k,2 , ...s k,i-1 ∼ N ηk,i , ε2 k,i . (30)
In other words, every element of an MVG rv conditioned on all its previous elements follows a univariate Gaussian distribution with mean ηk,i and variance ε2 k,i . Secondly, since E k is a full-rank covariance matrix, it is symmetric positive definite. Consequently it can factorized uniquely by a lower triangular matrix L k with positive diagonal entries [30, Corollary 7.2.9]:

E k = L k L T k , (31)
which corresponds to the Cholesky decomposition of the covariance matrix.

Finally, let w k be a vector of M univariate standard Gaussian rvs. We can correlate white noise by multiplying it by the Cholesky decomposition of a chosen covariance matrix:

L k w k ∼ N (0, E k ) . ( 32
)
The parameters η and ε can be computed efficiently using the Cholesky decomposition and the realization of s using the following equations:

εk = diag(L k ), (33) ηk = (L k w k -diag(L k )w k ) . ( 34
)
Now observe that all of this methodology can be applied to every macro-block of the prestego signal s k as defined in Eq [START_REF] Cachin | An Information-Theoretic Model for Steganography[END_REF]. If we apply Eq (33 ,34) to every s k , we obtain a vector s of n elements such that:

si ∼ N (η i , εi ) . ( 35
)
Finally, using the chain rule of probability on each β

(j)
i , the embedding probabilities β (j) i are obtained by:

β (j) i = Φ j -r i -ηi + 0.5 εi -Φ j -r i -ηi -0.5 εi , (36)
where Φ(•) represents the cumulative distribution function of the standard normal distribution and r i = y i -[y i ] denotes the rounding error of i-th DCT coefficient. In practice, we perform a 2q +1-ary embedding. Consequently, the alphabet size of the embedding scheme is finite; j is thus constrained to a finite range q and the β (j) i must be normalized accordingly:

β (j) i,normalized = β (j) i q j=-q β (j) i . ( 37
)
To understand Equation (36), observe that after adding the prestego signal, the new value will either stay in the same integer bin after rounding or fall into a neighboring one. The probability of falling into one bin or another depends on the original rounding error r i , the conditioned mean ηi and variance ε2

i of the prestego signal -see Figure 5. Computing β (j) i is then simply a matter of computing the probability of falling into each bin which is simply the area under the curve of the prestego signal centered at the original rounding error inside each integer bin.

Note that for the STC implementation as described in [17, Section V.C], one would have to convert the embedding probabilities into costs using the following equation:

β (j) i = e -λρ (j) i 1 + j̸ =0 e -λρ (j) i , (38)
where ρ

(j)
i is the cost of modifying the i-th DCT coefficient by +j and λ is the Lagrange multiplier which allows solving Eq [START_REF] Lehmann | Testing statistical hypotheses[END_REF].

VI. PRACTICAL IMPLEMENTATION

We are now ready to give the algorithms for the simulators of the steganographic scheme outlined in Section V. One algorithm is given for each of the two dependency models defined in Section II-C, namely the independent macro-block model and the lattice model. are computed with a prestego signal centered at the original rounding error r i . The dashed boxes represent the integer bins after rounding the DCT coefficient whereas the colored parts represent the area under the curve of the pre-stego signal which will be taken into account for each β (j) i .

β 0 i β 1 i β 2 i β -1 i β -2 i r i + ηi -2 -1 0 1 2
All the algorithms specified in this section assume that the steganographer has access to a precover y, its corresponding RAW file x and full knowledge of the processing pipeline. It is assumed that the covariance matrices of the cover are estimated using the method given in [START_REF] Giboulot | Detectability-Based JPEG Steganography Modeling the Processing Pipeline: The Noise-Content Trade-off[END_REF] which will not be repeated here due to space considerations.

The effect of pixel saturation, which is not directly taken into account by our model, is discussed a the end of the section.

A. Independent macro-block model

The first dependency model assumes that all the elements in the same √ N × √ N macro-block are possibly dependent while elements in two different macro-blocks are considered independent. As shown in Figure 2a, only non-overlapping macro-blocks are considered.

To build the macro-blocks in this model, one just has to split the whole image y into non-overlapping macro-blocks y k of the same size so that every element y i of the precover belongs to one and only one macro-block.

Consequently we have each y k independent from every other y k and following a MVG distribution:

y k ∼ N (µ k , Σ k ) . ( 39
)
To simulate the stego signal, we follow Eq [START_REF] Lehmann | Testing statistical hypotheses[END_REF]. The steganographer first samples the prestego signal according to the precover distribution and secondly, performs a bisection search on α in order to scale the signal so that it matches the payload constraint. The signal is finally added to the precover and quantized.

B. Lattice model

The second dependency model we study in this paper is the lattice model. In this model, we assume dependencies between DCT coefficients within the same block as well as among DCT coefficients with neighboring blocks.

First, we introduce some notation that will be used throughout this subsection. We write Λ to denote a set of blocks, i is sampled conditionally on each of its drawn neighboring blocks. which we call a lattice, such that all blocks in the set are independent from one another.

For a given block y k we write y cardinal k with cardinal ∈ {N, S, E, W, N E, N W, SE, SW } to designate the block which is respectively above, below, right of, left of, etc. . . of y k . Now, from Figure 2b, observe that the lattice model assumptions lead to a natural decomposition of the image into four lattices

Λ 1 , Λ 2 , Λ 3 , Λ 4 .
The steganographer's goal is, first, to sample a MVG signal s with the same covariance as the noise of the precover. To do so, we use the fact that the pdf of the whole image can be decomposed using the chain rule of probability so that the pdf of s is:

p (s) =p s Λ1 p s Λ2 |s Λ1 p s Λ3 |s Λ1 , s Λ2 p s Λ4 |s Λ1 , s Λ2 , s Λ3 . (40)
Consequently, we can first sample the prestego signal in Λ 1 by sampling from p s Λ1 , then sample the prestego signal in Λ 2 according to p s Λ2 |s Λ1 and so on for the two other lattices. Now let Λ * 1 , Λ * 2 , Λ * 3 , Λ * 4 be such that:

Λ * 1 = Λ 1 , (41)
Λ * 2 = {[y N E , y N W , y SE , y SW , y]|y ∈ Λ 2 ]}, (42)
Λ * 3 = {[y N , y S , y E , y W , y]|y ∈ Λ 3 ]}, (43)
Λ * 4 = {[y N , y S , y E , y W , y N E , y N W , y SE , y SW , y]|y ∈ Λ 4 ]}. (44)
In other words, each Λ * contains vectors built from each block in Λ along with the blocks from previous lattices on which it is dependent. See Figure 6 for a graphical representation of the dependency structure of each lattice used to construct each Λ * . In what follows, we refer to the covariance matrix of the k-th vector of the i-th lattice Λ * i as Σ

Λ * i k
By using exactly the same reasoning as in Section V-B, adapted to a multivariate setting, we can easily sample the blocks in each lattice by using the Cholesky decomposition of the covariance of each vector in each Λ * i : Once every s Λi k is computed, the algorithm to sample the pre-stego signal is exactly the same as for the independent macro-block model.

s Λi k = L Λ * i k w Λi k , ∀i ∈ {1, 2, 3, 4} (45) with L

C. Impact of pixel saturation

An important feature of the noise model in the RAW domain, which was not taken into account in Section II, is the saturation of photo-sites. It is known from the work of Foi et al. [START_REF] Foi | Practical Poissonian-Gaussian Noise Modeling and Fitting for Single-Image Raw-Data[END_REF] that the heteroscedastic model breaks near saturating values of camera sensor. The variance decreases sharply in this regime due to the clipping of the photo-sites values -see Figure 7 for an example.

In order to solve this problem while keeping the simplicity of our model, we have experimentally found that forbidding some embedding locations is counter-productive as it destroys the covariance structure of the block where the changes are forbidden.

A better solution is to modify directly the covariance matrix in the RAW domain during the estimation phase. In particular, we empirically found that a relevant heuristic is to set the variance to 1 when the mean value of the photo-site µ i is greater than 0.95S where S is the saturating value of the camera sensor.

Similarly, it is important to set a threshold for variances near zero so that the Cholesky decomposition can still be numerically computed. As long as the value is small we have not found the exact value of the threshold to matter; as a consequence, we fixed it to 10 -5 , that is we fix the variance σ 2 i to 10 -5 if σ 2 i < 10 -5 .

VII. NUMERICAL EVALUATIONS AND COMPARISONS

In this section, we study the performance of our different extensions of Gaussian Embedding in the JPEG domain. To have access to a precise estimation of the covariance matrix, we use the estimation method described in our previous work [START_REF] Giboulot | Detectability-Based JPEG Steganography Modeling the Processing Pipeline: The Noise-Content Trade-off[END_REF]Section II]. Consequently, we use the knowledge of the RAW file and the processing pipeline in those cases and prefix the name of the embedding scheme with Σ. When

GE

Minimizes the power of of the MP detector in the continuous domain supposing every DCT coefficients to be independent as described in [START_REF] Giboulot | Detectability-Based JPEG Steganography Modeling the Processing Pipeline: The Noise-Content Trade-off[END_REF].

MGE Intra Only

Minimizes the KL divergence detector in the continuous domain supposing 8 × 8 DCT blocks to be independent. MGE Intra+Inter Minimizes the KL divergence in the continuous domain supposing 24 × 24 DCT macro-blocks to be independent. Lat MGE Minimizes the power of the MP detector in the continuous domain using lattice embedding as described in Section VI-B.

SI-UNIWARD

Side informed distortion based schemes as described in [START_REF] Holub | Universal distortion function for steganography in an arbitrary domain[END_REF].

using the estimation method described in Section III which does not require access to the RAW dataset, we add the suffix noRAW at the end of the embedding scheme name. We use the BOSS RAW dataset excluding the M9 camera because of the peculiar distribution of its photonic noise (see [START_REF] Denemark | Natural Steganography in JPEG Compressed Images[END_REF], Fig. 2) which would lead to an imprecise estimation of the covariance matrix. From this dataset comprising 7240 RAW images taken with 6 different cameras, we produce two new datasets using two different processing pipelines: a linear processing pipeline and a processing pipeline close to the original BOSSBase. Both these pipelines output 264 × 264 greyscale JPEG images. The details are exposed in Table I.

Note that for cropping we used an algorithm -Edge crop 1 -1 Available at https://alaska.utt.fr/#material which selects crops containing the greatest number of edges -that is the zone that should contain the most textured areas. This choice of cropping was made because it is known that SI-UNIWARD does not perform well on smooth images [18, Section V.A]; using such a cropping algorithm allows us to compare our embedding schemes in a situation where SI-UNIWARD is not disfavored due to a suboptimal content choice.

The parameters of the photonic noise c 1 and c 2 were estimated as described in [18, Section II] using the algorithm described in [START_REF] Thai | Statistical Model of Quantized DCT Coefficients: Application in the Steganalysis of Jsteg Algorithm[END_REF]. The H matrix is estimated once for each processing pipeline using the method described in Section III-B. The different embedding schemes used as well as their parameters are described in Table II. Since we performed optimal simulations of the embedding, which finds the optimal stego-signal in the continuous domain, we always chose the smallest alphabet necessary to encode all amplitudes of the stego-signal once discretized.

Steganalysis was performed with Efficient-Net-b3 [START_REF] Tang | CNN-Based Adversarial Embedding for Image Steganography[END_REF] modified so that the stride of the stem is equal to 1. The network was trained by training from the highest payload for 30 epochs, then 10 epochs for other payloads. The model for the highest payload was initialized with ImageNet weights. The base learning rate was fixed at 0.0005 and divided by 2 on loss plateau. The batch size was fixed to 24. This configuration of Efficient-Net-b3 was the one used during the ALASKA2 competition. The rest of the parameters are initialized in the same way as in the original paper.

A. Performance evaluation

We now proceed to discuss the performance of the embedding schemes developed in this paper with SI-UNIWARD on the two datasets presented at the beginning of the section. The results are presented in Figure 8a-9c.

Beginning with QF100, regardless of the processing pipeline, there is a clear gap between the performance of SI-UNIWARD and the different variants of MGE. For the linear pipeline, the gap is, on average, of 6% and 9% in terms of absolute P E for Σ-MGE Intra only and Σ-MGE intra + inter respectively. However, it is of 22% on average when using Σ-Lat MGE, showing the importance of using the most precise model for this pipeline. The difference between the different models is less pronounced for the Bosslike pipeline with an average gain of 13.5% in terms of absolute P E wrt SI-UNIWARD for the MGE schemes whatever the chosen correlation model. This is due to the fact that the downsampling operation removes most correlations between blocks due to the removal of neighboring pixels before the DCT transform.

At QF95, the difference between the different schemes becomes less pronounced for the Linear pipeline. However, there is still a gap in performance between these schemes and SI-UNIWARD. Compared with Σ-Lat MGE, there is an average gain of 10.5% wrt to SI-UNIWARD for the Linear pipeline and 4% for the BOSS pipeline. At QF75, every scheme performs approximately the same irrespective of the pipeline. This is most likely due to the fact that, at such a low QF, most of the covariances are now close to 0 -see [START_REF] Giboulot | Detectability-Based JPEG Steganography Modeling the Processing Pipeline: The Noise-Content Trade-off[END_REF][Section VI.B]. As such most of the performance of the steganography is likely due to the side-information related to the rounding errors.

We also note that the implementation of Σ-Lat MGE which does not use the RAW file has a performance virtually identical with the original implementation showing that the assumptions on the processing pipeline we used in Section III are quite practical for a steganography context.

Finally, note that when using Σ-GE, which does not use any correlation between DCT coefficient, there is always a small gain with respect to SI-UNIWARD for both QF100 and QF95. However, its performance are always subpar even compared to Σ-MGE intra only, showing the importance of taking these correlations into account.

One interesting thing to note here is the impact of the processing pipeline on the correlation structure necessary to obtain good performances. In the Bosslike case, the downsampling operation has made most inter-block dependencies very small compared to intra-block dependencies. Consequently, using a more sophisticated model of dependencies does not bring any gain in performance. On the other hand, in the Linear pipeline case where all dependencies are preserved until the end of the pipeline since no downsampling is performed, there is quite a substantial gain when using the most sophisticated model.

B. Impact of the estimation of the covariance matrix

In this subsection we study the impact of errors on the estimation of the covariance matrix. Three main types of error can occur on the estimation: errors on the heteroscedastic parameters c 1 and c 2 , errors on the pipeline matrix H and errors due to the saturation of the pixels (see Section VI-C).

Errors on the heteroscedastic parameters can be studied analytically. Let c 1 and c 2 be the true parameters; ĉ1 and ĉ2 are the estimated parameters. The variances σ i and estimated variances σi in the RAW domain are given by σ i = c 1 µ i + c 2 and σi = ĉ1 µ i + ĉ2 respectively. Let α ′ = ĉ1 c1 be relative estimation error on the c 1 parameter. We have that σi = α ′ σ i + C where C = ĉ2α ′ c 2 is a constant which does not depend on the photo-site. Consequently, the resulting estimated covariance Σk in the developed domain is given by Σk = α ′ Σ k + CHH T . As a particular case, note that if we have ĉ1 ĉ2 = c1 c2 , then the estimation error has no impact on the MGE steganography schemes, because a multiplicative error on the covariance does not change the optimal solution. If we take the Bosslike pipeline as an example, we computed HH T for each camera and found that the average absolute value of the non-diagonal entries is of the order of 10 -7 and 10 -4 for the diagonal elements. If we set ĉ2 = 0, the relative estimation error would then have to be very high (i.e, in the order of at least 100) to begin to have an impact. To validate this analysis on the Bosslike pipeline, we repeated the experiments of Section VII-A except that the covariance matrices were estimated by fixing c 1 = 0.5 and c 2 = 0 for every image. The results are given in Figure 10a. As expected they are extremely close to the original Σ-Lat-MGE results on this pipeline.

The two other types of error are an error on the estimation of the processing pipeline and the error due to saturation of the pixels. In this paper, we always assumed the steganographer had access to the correct processing pipeline for the estimation of the covariance matrix, even in the "no RAW" case presented in Section III. It would then be interesting to see the impact on security when using mismatched pipelines for the covariance estimation. To do so, we repeated the experiments of Section VII-A, but this time using the covariance matrices estimated with the Boss pipeline on the Linear dataset and vice versa. Similarly, we also repeated the experiments when using the heuristic presented in Section VI-C and when not using it. The results are presented in Figure 10a-10b.

The results are clear: using strongly mismatched pipelines or not taking the saturation of pixels into account leads to near useless schemes, showing the importance of having a good model of the pipeline in the first place. Note, however, that such a mismatch is quite significant since one pipeline uses cropping and the other downsampling, leading to very different correlation structures. A complete study of this phenomenon should find a good measure of mismatch using the covariance matrices and link this measure to empirical detectability but it is out of the scope of this paper.

VIII. CONCLUSION

This paper is the conclusion to a series of work on Gaussian Embedding. Through the development of a multivariate Gaussian model of the sensor noise and the processing pipeline, we derived the optimal detector when the sensor noise follows this very general model. We showed that in this setting, the signal which minimizes the KL-divergence under a given payload constraint has a covariance matrix proportional to the covariance matrix of the cover noise. With these results, we designed the general form of Gaussian Embedding using a lattice embedding strategy. This allows us to use an embedding scheme taking into account correlations between every neighboring block. Furthermore, this algorithm does not need access to a RAW file. This yields an algorithm that beats the previous state-of-the-art by an important margin, even when our model does not match exactly the real distribution of the sensor noise.

To do so, we begin by showing the equality for the trace term:

trace (Σ + E) -1 Σ = trace L T (Σ + E) -1 L = trace L T (Σ + E) -1 L -1 -1 = trace L -1 (Σ + E) L -1 T -1 = trace I + L -1 E L -1 T -1
= trace (I + E w ) -1 .

(

) 59
and secondly for the determinant term:

|Σ + E | • |Σ -1 | = |I + E Σ -1 | = |I + E L -1 (L -1 ) T | = |I + E w |, (60)
with the last line obtained using Sylvester's Law of determinant. This validates the invariance of the KL-divergence to our change of basis. Using Eq (58), we rewrite the system in Eq (54) as: With a change of variable such that k ′ i = log (k i ) we obtain the following optimization problem:

   min k ′ i 1 2 n i=1 1 e k ′ i +1 + log(e k ′ i + 1) -n R ′ = n i=1 k ′ i . ( 63
)
where R ′ = -2R + log (2πe|E w |)n log (2πe).

Since the objective function is a sum of identical convex functions and the constraint is the sum of k ′ i , the problem is convex and consequently admits a single, global minimum which is attained when all k ′ i are equal. From this observation, it follows that E w is a a matrix where all eigenvalues are equal. Since, E w is positive definite (hence normal), this implies that it is proportional to the identity matrix:

E w = αI, (64)
for some α > 0.

Going back to the original system in Eq (54), the solution is obtained by:

E = L E w L T = LαIL T = αΣ.
(65)

  Observations ξ = {ξ 0 , . . . , ξ m } Precover parameters Denoised image: µ = (µ 1 , . . . , µ m ) Precover noise covariances: Σ = (Σ 1 , . . . , Σ m ) Prestego parameters Prestego signal covariances: E = (E 1 , . . . , E m )

Fig. 5 :

 5 Fig. 5: Figure explaining how β (j) i

4 Fig. 6 :

 46 Fig. 6: Graphical representation of the blocks used to construct each Λ * i . Note that the central block of each figure, corresponding to the last block of each Λ *i is sampled conditionally on each of its drawn neighboring blocks.

kFig. 7 :k

 7 Fig. 7: Variance σ 2 10000 samples of N (µ, c 1 µ + c 2 ) as function of µ with c 1 = 1 and c 2 = 1000. The curve in red is obtained when the samples are clipped at 0 and 4096 while the curve in black is obtained when the samples are not clipped..

Fig. 8 :Fig. 9 :

 89 Fig.8: P E as function of payload size for BossBase developed with the BOSS pipeline using Efficient-Net-b3.

Fig. 10 :

 10 Fig.10: P E as function of payload size for different the different heuristics and design decisions made in this paper.

minD

  Ew KL (p I || q I,Ew ) R = 1 2 log (2πe|E w |)log 2πe|Σ -1 | .(61)Now, let k i be the i-th eigenvalues of E w , we have that:D KL (p I || q I,Ew )

TABLE I :

 I Names and operations of the processing pipelines used in the experiments. The operations are performed in the order they are presented in the table

	Pipeline name	Demosaicking	White Balance	RGB to grey Downsampling method
	Linear Pipeline Bilinear	No	Yes	Edge crop, 264 × 264
	BOSS Pipeline	PPG	Yes, Camera	Yes	Resize from 792 × 792 (Edge crop) to 264 × 264, Lanczos kernel

TABLE II :

 II Nomenclature of the embedding schemes

	Name	Meaning

APPENDIX

In this appendix, we provide proofs and derivations of the main results of the paper. First, we derive the optimal detector of the hypothesis testing problem given in Section IV. Secondly we simplify the power minimization problem as a minimization of the KL-divergence of the likelihood-ratio and show that the optimal covariance in this setting is a scaling of the covariance of the cover noise.

A. Derivation of the optimal detector

To be as general as possible, we will not consider any particular dependency model in this appendix. As such we will consider that the sample ξ is the whole image under scrutiny. Consequently we try to derive the optimal test which discriminate between the two following hypotheses:

Notice that we consider here the covariance of the whole image Σ and not the covariance of the blocks. Now recall that the optimal test under our setting is the likelihood-ratio test given by Eq [START_REF] Thai | Statistical Model of Quantized DCT Coefficients: Application in the Steganalysis of Jsteg Algorithm[END_REF]. To compute the power of this test, we need to compute the distribution of the likelihood ratio under both hypotheses. To do so, we use the fact that the statistic of the LRT can be written as a quadratic form -see [START_REF] Van Trees | Detection estimation and modulation theory[END_REF][Chapter 3, Section 3] for a derivation:

This test can be simplified as:

by putting the contribution of the constant values in the threshold.

We will show that this test can be rewritten as a sum of weighted independent standard chi-squared rv.s. We will develop only the case under H 0 ; the case for H 1 follows exactly from the same method.

Let:

and decompose Σ using the symmetric square-root matrix :

Using the spectral theorem we can write:

where U is an orthogonal matrix and K a diagonal matrix. Finally, let ξ u = U T Σ -1/2 ξ and note that it follows a centered multivariate Gaussian distribution with covariance the identity matrix (by the orthogonality of U).

Now write:

The statistic in Eq (48) is thus a realization of a sum of independent standard chi-squared random variable weighted by the eigenvalues of AΣ and A (Σ + E) under H 0 and H 1 respectively.

B. Form of the covariance matrix which minimizes the KLdivergence

Computing the power of the LRT when the distribution under both hypotheses are weighted sum of chi-square distribution is difficult in the general case. Furthermore, we cannot appeal to the central limit theorem since we have no information on the distribution of the weights K ii .

If we cannot give the form of the optimal pre-stego signal due to this fact, we can however give security guarantees for a given pre-stego signal. To do so, we use the result of Cachin's work on steganographic security [START_REF] Cachin | An Information-Theoretic Model for Steganography[END_REF], which shows, using a simple data-processing inequality, that the performance of an optimal detector is upper bounded by the KL divergence between the distributions of the two hypotheses.

Instead of minimizing the power of LRT under an entropy constraint we will thus minimize the KL-divergence D KL (p||q) under an entropy constraint:

where R is the payload constraint in the continuous domain and |.| the matrix determinant. Recall that we assume here Σ and E to be both positive definite matrices. Note that the KL divergence between p and q is given by:

The main idea of this proof is that we can simplify the problem by working in a basis where the covariance of the cover noise is the identity matrix. This is done by showing that the KL-divergence is invariant when changing to this basis. We then express the Lagrangian of the system as a function of the eigenvalues of E and show that there exists a unique minimum for our system.

First, let L be the Cholseky decomposition of Σ such that:

Let us also define E w as:

We now want to show the following equality:

where I is the identity matrix (with relevant dimensions).


Aurélien De Jode
email: aurelien.dejode@imbe.fr

Sophie Dubois

Zinovia Erga

Dorian Guillemain

Romain David

Frédéric Zuberer

Jean-Pierre Feral

Anne Chenuil

Understanding ecological functionning of ecological habitats, and building New Indicators based on Genetic Tools to assess their GES (Good Environmental Status)

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

More precisely, it is defined by the bioconcretions of corallines red algae that grows at dim light conditions. Corallines are the main builders, or enginners species of this ecosystems but their diversity and their ability to colonize new sites is poorly understood. The genus of calcareous red algae Lithophyllum is one of the most important builders in the area of Marseilles, but some marine invertebrates, such as Myriapora truncata, directly contribute to the framework of the habitats because they build their own calcareous skeleton (Hong 1980). However coralligenous habitats are built by various organisms and provide very different ecosystems services from recreational diving to shelter or nurseries for baby fishes and other species.

The aim of this project is to improve our understanding of its functioning and resilience capacity, and establish biodiversity and connectivity patterns among localities and ecological profiles, providing information for rationalizing the design of Marine Protected Area networks and monitoring methods. To reach this goal, a multidisciplinary approach will be used, combining population genetics, community ecology, microbial ecology and physical oceanography. This requests four tasks, partly using the same data, but different analysis methods. 1-Establishment of the taxonomic composition of the engineer corallinale algae. 2-Establishment of the species composition using metabarcoding (together with traditional taxonomy and photo quadrate validation for some samples) 3-Population structure and phylogeography of two selected taxa: a red alga Lithophyllum spp. and a bryozoan Myriapora truncata, on the CIGESMED sampling network. 4-Synthesis and Comparisons: Propose new GES (good environmental status) indicators for the coralligenous and biodiversity management rules.

Preliminary studies of the two engineers species revealed genetic differentiation in the Bay of Marseille.

Our first results suggests that the engineers species are actually composed of cryptic species that may have live in different ecological niches. It emphasizes our lack of knowledge about this endemic Mediterranean ecosystem and stresses our requirement for a better understanding of its functioning in order to rationalize its monitoring and protection.

In this project different methods DNA extraction, PCR, metabarcoding and disciplines will be used : population genetics, phylogeography, community ecology and hydrological modeling.

Genetic data obtained after

Principle of metabarcoding

Wall rugosity

Orientation Cover Slope

The coralligenous is a typical Mediterranean marine habitat, Paradoxically, unlike Posidonia meadows, it is poorly studied. Like the famous tropical coral reefs, it is a complex biogenic habitat, based on encrusting red algae (Ballesteros 2006).

Heterogeneous coralligenous habitats

Paramuricea clavata Eunicella singularis Eunicella cavolinii

Peyssonnelia squamaria

Acknowledgements

This work is a contribution to the Labex OT-Med (n° ANR-11-LABX-0061) funded by the French Government «Investissements d'Avenir» program of the French National Research Agency (ANR) through the A*MIDEX project (n° ANR-11-IDEX-0001-02). *Other Contributors to the CIGESMED program : Açik Çinar S., Andral B., Arvanitidis C., Aurelle D., Aysel V.,Bakir K., Bellan G., Bellan-Santini D., Bouchoucha M., Cahill A., Celik C., Chatzigeorgiou G.,Chatzinikolaou E., Chenesseau S., Çinar M.E., Dağli E., Dailianis T., DimitriadisC., D'Iribarne C., Doğan A., Dounas C., Egea E., Emery E., Evcen A., Faulwetter S.,Gatti G., Gerovasileiou V., Güçver S.M., Issaris Y., Katağan T., Keklikoglou K., Kirkim F.,Koçak F., Koutsoubas D., Marschal C., Önen M., Önen S., Öztürk B., Panayiotidis P., Panteri E.,Pavloudi C., Pergent G., Pergent-Martini C., Poursanidis D., Ravel C., Reizopoulou S., Rocher C.,Ruitton S., Sakher S., Salomidi M., Sarropoulou E., Sartoretto S., Selva M., Sini M., Sourbes L., Simboura N.,Taşkin E., Thierry de Ville d'Avray L., Tschanz L., Vacelet J., Valavanis V., Vasileiadou A., Verlaque M.