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

HAL
open science

# On the Shift Invariance of Max Pooling Feature Maps in Convolutional Neural Networks

Hubert Leterme, Kévin Polisano, Valérie Perrier, Karteek Alahari

► To cite this version:

Hubert Leterme, Kévin Polisano, Valérie Perrier, Karteek Alahari. On the Shift Invariance of Max Pooling Feature Maps in Convolutional Neural Networks. 2023. hal-03779434v2

HAL Id: hal-03779434

https://hal.science/hal-03779434v2

Preprint submitted on 24 Oct 2023

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_2.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# On the Shift Invariance of Max Pooling Feature Maps in Convolutional Neural Networks

Hubert Leterme, Kévin Polisano, Valérie Perrier, Karteek Alahari

► To cite this version:

Hubert Leterme, Kévin Polisano, Valérie Perrier, Karteek Alahari. On the Shift Invariance of Max Pooling Feature Maps in Convolutional Neural Networks. 2025. hal-03779434v3

HAL Id: hal-03779434

https://hal.science/hal-03779434v3

Preprint submitted on 18 Apr 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
# On the Shift Invariance of Max Pooling Feature Maps in Convolutional Neural Networks*

Hubert Leterme†, Kévin Polisano‡, Valérie Perrier‡, and Karteek Alahari†

**Abstract.** This paper focuses on improving the mathematical interpretability of convolutional neural networks (CNNs) in the context of image classification. Specifically, we tackle the instability issue arising in their first layer, which tends to learn parameters that closely resemble oriented band-pass filters when trained on datasets like ImageNet. Subsampled convolutions with such Gabor-like filters are prone to aliasing, causing sensitivity to small input shifts. In this context, we establish conditions under which the max pooling operator approximates a complex modulus, which is nearly shift invariant. We then derive a measure of shift invariance for subsampled convolutions followed by max pooling. In particular, we highlight the crucial role played by the filter's frequency and orientation in achieving stability. We experimentally validate our theory by considering a deterministic feature extractor based on the dual-tree complex wavelet packet transform, a particular case of discrete Gabor-like decomposition.

**Key words.** deep learning, image processing, shift invariance, max pooling, dual-tree complex wavelet packet transform, aliasing

**MSC codes.** 42C40, 68U10, 94A08

**1. Introduction.** Understanding the mathematical properties of deep convolutional neural networks (CNNs) [22] remains a challenging issue today. On the other hand, wavelet and multi-resolution analysis are built upon a well-established mathematical framework. They have proven to be efficient for tasks such as signal compression and denoising [48], and have been widely used as feature extractors for signal, image and texture classification [17, 21, 37, 52]. There is a broad literature revealing strong connections between these two paradigms, as discussed in subsections 1.1 and 1.2. Inspired by this line of research, the present paper extends existing knowledge about CNN properties. Specifically, we assess the shift invariance of max pooling feature maps through both theoretical and empirical approaches in the context of image classification, by leveraging the properties of oriented band-pass filters.

**1.1. Motivations and Main Contributions.** CNNs rely on convolutions and nonlinear pooling operations to transform input images into high-level feature vectors, which are in turn processed for the task at hand. In the context of image classification, the feature vectors are fed into a linear classifier. In order to achieve high classification accuracy, a convolutional network is expected to retain discriminative image components while reducing intra-class variability [9, 23]. A key property that is often desired in CNNs is their ability to remain invariant to small

<footnote>*This work has been partially supported by the LabEx PERSYVAL-Lab (ANR-11-LABX-0025-01) funded by the French program Investissement d'avenir, as well as the ANR grant MIAI (ANR-19-P3IA-0003). Most of the computations presented in this paper were performed using the GRICAD infrastructure (https://gricad.univ-grenoble-alpes.fr), which is supported by Grenoble research communities.</footnote>

<footnote>†Univ. Grenoble Alpes, CNRS, Inria, Grenoble INP, LJK, 38000 Grenoble, France (hubert.leterme@univ-grenoble-alpes.fr, karteek.alahari@inria.fr).</footnote>

<footnote>‡Univ. Grenoble Alpes, CNRS, Grenoble INP, LJK, 38000 Grenoble, France (kevin.polisano@univ-grenoble-alpes.fr, valerie.perrier@univ-grenoble-alpes.fr).</footnote>
# On the Shift Invariance of Max Pooling Feature Maps in Convolutional Neural Networks

**Hubert Leterme**

Université Caen Normandie, ENSICAEN, CNRS, Normandie Univ, GREYC UMR 6072, F-14000 Caen, France

Univ. Grenoble Alpes, CNRS, Inria, Grenoble INP, LJK, 38000 Grenoble, France

hubert.leterme@ensicaen.fr

**Kévin Polisano**

Univ. Grenoble Alpes, CNRS, Grenoble INP, LJK, 38000 Grenoble, France

kevin.polisano@univ-grenoble-alpes.fr

**Valérie Perrier**

Univ. Grenoble Alpes, CNRS, Grenoble INP, LJK, 38000 Grenoble, France

valerie.perrier@univ-grenoble-alpes.fr

**Karteek Alahari**

Univ. Grenoble Alpes, CNRS, Inria, Grenoble INP, LJK, 38000 Grenoble, France

karteek.alahari@inria.fr

## Abstract

This paper focuses on improving the mathematical interpretability of convolutional neural networks (CNNs) in the context of image classification. Specifically, we tackle the instability issue arising in their first layer, which tends to learn parameters that closely resemble oriented band-pass filters when trained on datasets like ImageNet. Subsampled convolutions with such Gabor-like filters are prone to aliasing, causing sensitivity to small input shifts. In this context, we establish conditions under which the max pooling operator approximates a complex modulus, which is nearly shift invariant. We then derive a measure of shift invariance for subsampled convolutions followed by max pooling. In particular, we highlight the crucial role played by the filter's frequency and orientation in achieving stability. We experimentally validate our theory by considering a deterministic feature extractor based on the dual-tree complex wavelet packet transform, a particular case of discrete Gabor-like decomposition.

## 1 Introduction

Understanding the mathematical properties of deep convolutional neural networks (CNNs) (LeCun et al., 2015) remains a challenging issue today. On the other hand, wavelet and multi-resolution analysis are built upon a well-established mathematical framework. They have proven to be efficient for tasks such as signal compression and denoising (Vetterli, 2001), and have been widely used as feature extractors for signal, image and texture classification (Laine & Fan, 1993; Pittner & Kamarthi, 1999; Yen, 2000; Huang & Aviyente, 2008). There is a broad literature revealing strong connections between these two paradigms, as discussed in Sections 1.1 and 1.2. Inspired by this line of research, the present paper extends existing knowledge about CNN properties. Specifically, we assess the shift invariance of max pooling feature maps through both theoretical and empirical approaches in the context of image classification, by leveraging the properties of oriented band-pass filters.

### 1.1 Motivations and Main Contributions

CNNs process input images through convolutions and nonlinear pooling operations, transforming them into high-level feature vectors that are subsequently used for the task at hand. In image classification, the feature vectors are fed into a linear classifier. To achieve high classification accuracy, a convolutional network must
input transformations, such as translations, rotations, distortions, or scaling [6, 9, 26, 43, 51]. Since perfect invariance is seldom achieved, we shall also use the term *stability* to refer to this behavior. This paper targets translations, also called shifts.

Furthermore, we focus on a configuration that is commonly observed in CNNs when trained on image datasets: many convolution kernels in the first layer resemble band-pass oriented waveforms [38, 53], referred to as *Gabor-like filters*. Whether extracted features are stable to translations is partly addressed by [2, 56]. These papers point out that strided convolution and pooling operators may greatly diverge from shift invariance, due to aliasing when subsampling high-frequency signals. In response, recent works [56, 59] introduced an antialiasing method based on low-pass filtering. They managed to increase both stability and predictive power of CNNs, despite the resulting loss of information.

In the current paper, we show that, under specific conditions that we establish, the max pooling operator can actually partially restore shift invariance. We unveil a connection between the output of the first max pooling layer and the modulus of complex Gabor-like coefficients, which is known to be nearly shift invariant. This work led us to develop a method for improving shift invariance in CNNs which, unlike the previously-mentioned papers, preserves high-frequency information [25].

**1.2. Related Work.** Analyzing the invariance properties of CNNs is critical as it enables to identify their shortcomings and provides an opportunity to enhance their performance. In recent years, several works focused on this topic.

**1.2.1. Wavelet Scattering Networks.** Most notably, Bruna and Mallat [9] developed a family CNN-like architectures, named *wavelet scattering networks* (ScatterNets), based on a succession of complex convolutions with wavelet filters followed by nonlinear modulus pooling. They produce translation-invariant image representations which are stable to deformation and preserve high-frequency information [28, 29]. A variation has been proposed by Sifre and Mallat [43] to include rotational invariance. ScatterNets achieve strong performance on handwritten digits and texture datasets, but do not scale well to more complex ones. To overcome this, Oyallon et al. [32, 33] introduced hybrid ScatterNets, where the scattering coefficients are fed into a standard CNN architecture, showing that the network complexity can be reduced while keeping competitive performance. Derived models include ScatterNets built upon the dual-tree complex wavelet transform [44], learnable and parametric ScatterNets [10, 14], geometric ScatterNets operating on Riemanian manifolds [36], and graph ScatterNets [13, 58]. Also worth mentioning, Czaja and Li [11, 12] studied ScatterNets based on uniform covering frames, i.e., frames splitting the frequency domain into windows of roughly equal size, much like DT-CWPT frames (as used in the present paper). Other works by Zarka et al. [54, 55] proposed to sparsify wavelet scattering coefficients by learning a dictionary matrix, to learn 1 × 1 convolutions between feature maps of scattering coefficients and to apply soft thresholding to reduce within-class variability.

ScatterNets are specifically designed to meet some desired properties. As deep learning architectures with well-established mathematical properties, they are sometimes used as explanatory models for standard, freely-trained networks. However, whether their properties are transferable to a broader class of models is unclear, because the former rely on complex-valued convolutions whereas more conventional architectures exclusively employ real-valued kernels.
preserve discriminative image features while reducing intra-class variability (LeCun et al., 1998; Bruna & Mallat, 2013). An important and often desired property of CNNs is their ability to remain invariant to small input transformations, such as translations, rotations, distortions, or scaling (Liao & Peng, 2010; Bruna & Mallat, 2013; Sifre & Mallat, 2013; Bietti & Mairal, 2017; Wiatowski & Bölcskei, 2018; Cahill et al., 2024).

In particular, an image in which the main subject is slightly shifted from its original position should retain its initial classification. The absence of this property, known as *translation invariance* or *shift invariance*, can negatively impact the model's predictive performance. Consequently, a reliable model must learn this property through extensive training on large datasets. The main focus of this paper is therefore to assess whether shift invariance is inherently guaranteed by the model's architecture. Since perfect invariance is rarely achieved, we also the term *stability* to refer to this behavior.

More specifically, we focus on commonly observed phenomenon in CNNs when trained on image datasets: many convolution kernels in the first layer resemble band-pass oriented waveforms (Yosinski et al., 2014; Rai & Rivas, 2020), referred to as *Gabor-like filters*. Whether the features extracted by these filters remain stable under translations has been partly addressed by Azulay & Weiss (2019); Zhang (2019), who highlight that strided convolution and pooling operations can significantly diverge from shift invariance, due to aliasing effects when subsampling high-frequency signals. In response, Zhang (2019); Zou et al. (2023) introduced antialiasing methods based on low-pass filtering, improving both stability and predictive performance—albeit at the cost of some loss of information.

In the current paper, we show that, under specific conditions that we establish, the max pooling operator can partially restore shift invariance. We unveil a connection between the output of the first max pooling layer and the modulus of complex Gabor-like coefficients, which is known to be nearly shift invariant. This work offers a promising direction for improving shift invariance in CNNs while preserving high-frequency information—unlike the previously-mentioned approaches.

Before proceeding further, we emphasize that our study is not limited to purely convolutional architectures. In recent years, self-attention mechanisms have gained significant interest in computer vision due to their ability to model complex, long-range dependencies in image representations. Notably, the *vision transformer* (ViT) was introduced by Dosovitskiy et al. (2021) by adapting the transformer architecture, initially developed for natural language processing (NLP) (Vaswani et al., 2017), to computer vision tasks. Unlike CNNs, the ViT operates without any convolutional layers. Instead, input images are partitioned into fixed-size patches, which serve as inputs to the first self-attention module. However, more recent research has explored hybrid architectures that integrate self-attention with convolutional components (Wu et al., 2021; Yuan et al., 2021; Hassani et al., 2022; Li et al., 2023; Yin et al., 2024). This approach allows for reducing the amount of labeled data required, while achieving faster training and improving generalizability. In particular, the first layers of a CNN can be used as a “convolutional token embedding,” replacing the naive patch extraction used in the original ViT. The theoretical framework presented in this paper also applies to such hybrid architectures, providing a better understanding of the invariance properties of the inputs to self-attention modules.

## 1.2 Related Work

Analyzing the invariance properties of CNNs is critical as it enables to identify their shortcomings and provides an opportunity to enhance their performance. In recent years, several works focused on this topic.

### 1.2.1 Wavelet Scattering Networks

Most notably, Bruna & Mallat (2013) developed a family CNN-like architectures, named *wavelet scattering networks* (ScatterNets), based on a succession of complex convolutions with wavelet filters followed by nonlinear modulus pooling. They produce translation-invariant image representations which are stable to deformation and preserve high-frequency information (Mallat, 2012; 2016; Czaja et al., 2024). A variation has been proposed by Sifre & Mallat (2013) to include rotational invariance. ScatterNets achieve strong performance on handwritten digits and texture datasets, but do not scale well to more complex ones. To overcome this, Oyallon (2017); Oyallon et al. (2018) introduced hybrid ScatterNets, where the scattering
Moreover, the modulus operator is used as an activation and pooling layer in ScatterNets, whereas standard CNNs implement pointwise nonlinear operators such as ReLU and spatial pooling layers such as max pooling. This limitation has been pointed out by Tygert et al. [47] as an argument in favor of complex-valued CNNs. In this context, our work seeks evidence that properties established for complex-valued networks are—to some extent—embedded in standard architectures.

### 1.2.2. Invariance Studies in CNNs.
Wiatowski and Bölcskei [51] considered a wide variety of feature extractors involving convolutions, Lipschitz-continuous non-linearities and pooling operators. The paper shows that outputs become more translation invariant with increasing network depth. However, these results do not fully extend to the discrete framework, because subsampled convolutions with band-pass real-valued filters can introduce aliasing artifacts, resulting in instability to translations [2, 56]. The current paper specifically addresses this issue.

Another line of work is focused on modeling and studying CNNs from the point of view of convolutional kernel networks [5, 7, 8, 40]. These authors showed that certain classes of CNNs are contained into the reproducing kernel Hilbert space (RKHS) of a multilayer convolutional kernel representation. As such, stability metrics are estimated, based on the RKHS norm which is difficult to control in practice. Kernel representations do not seem to suffer from aliasing effects; this can be explained by the Gaussian pooling layers that have been employed instead of max pooling: by discarding high-frequency information, shift invariance is preserved. Finally, some papers studied stability of CNNs in a broader sense, measured in terms of Lipschitz continuity [3, 35, 45, 49, 57]. However, the Lipschitz bounds, which have been obtained theoretically, are generally several orders of magnitude higher than empirical results. This discrepancy may be due to the fact that these bounds were obtained for generic situations and represent overly conservative worst-case scenarios, rather than typical real-world situations. Furthermore, the specific case of convolutions with band-pass Gabor-like filters have been overlooked, except for Pérez et al. [35].

In summary, we have identified the following blind spots in the literature, regarding the topic of studying shift invariance in CNNs.

* The effect of the max pooling operator on network stability under small input shifts has not been investigated, particularly when used in combination with Gabor-like convolutions.

* While the shift invariance of CNNs tends to increase with network depth in the continuous framework, in the discrete case, the presence of subsampled convolutions with oriented band-pass filters can lead to aliasing artifacts. To our knowledge, the literature lacks theoretical studies that take these aliasing effects into account.

* Although extensive studies have been conducted on complex-valued convolutions followed by modulus, a link is missing to extend these results to standard CNNs, which implement real-valued convolutions and spatial pooling operators.

All these points have been tackled in the present paper, from both theoretical and empirical perspectives.

## 1.3. Paper Outline.
In what follows, $l_{\mathbb{R}}^2(Z^2)$ and $l_{\mathbb{C}}^2(Z^2)$ represent the discrete spaces of square-summable two-dimensional sequences with values in $\mathbb{R}$ and $\mathbb{C}$, respectively. Let
coefficients are fed into a standard CNN architecture, showing that the network complexity can be reduced while keeping competitive performance. Derived models include ScatterNets built upon the dual-tree complex wavelet transform (Singh & Kingsbury, 2017), learnable and parametric ScatterNets (Cotter & Kingsbury, 2019; Gauthier et al., 2022), geometric ScatterNets operating on Riemanian manifolds (Perlmutter et al., 2020), and graph ScatterNets (Gama et al., 2019; Zou & Lerman, 2020). Also worth mentioning, Czaja & Li (2019; 2020) studied ScatterNets based on uniform covering frames, i.e., frames splitting the frequency domain into windows of roughly equal size, much like DT-CWPT frames (as used in the present paper). Other works by Zarka et al. (2020; 2021) proposed to sparsify wavelet scattering coefficients by learning a dictionary matrix, to learn 1 × 1 convolutions between feature maps of scattering coefficients and to apply soft thresholding to reduce within-class variability.

ScatterNets are specifically designed to meet some desired properties. As deep learning architectures with well-established mathematical properties, they are sometimes used as explanatory models for standard, freely-trained networks. However, whether their properties are transferable to a broader class of models is unclear, because the former rely on complex-valued convolutions whereas more conventional architectures exclusively employ real-valued kernels. Moreover, the modulus operator is used as an activation and pooling layer in ScatterNets, whereas standard CNNs implement pointwise nonlinear operators such as ReLU and spatial pooling layers such as max pooling. This limitation has been pointed out by Tygert et al. (2016) as an argument in favor of complex-valued CNNs. In this context, our work seeks evidence that properties established for complex-valued networks are—to some extent—embedded in standard architectures.

1.2.2 Invariance Studies in CNNs

Wiatowski & Bölcskei (2018) considered a wide variety of feature extractors involving convolutions, Lipschitz-
continuous non-linearities and pooling operators. The paper shows that outputs become more translation
invariant with increasing network depth. Additionally, Cahill et al. (2024) designed a family of operators
called *max filters*, which encompass a wide variety of operators including, in specific cases, the max pooling
operator. Stability with respect to diffeomorphisms were established, following the ideas developed for scat-
tering networks. However, these results do not fully extend to the discrete framework, because subsampled
convolutions with band-pass real-valued filters can introduce aliasing artifacts, resulting in instability to
translations (Azulay & Weiss, 2019; Zhang, 2019). The current paper specifically addresses this issue.

Another line of work is focused on modeling and studying CNNs from the point of view of convolutional
kernel networks (Bietti & Mairal, 2019a;b; Scetbon & Harchaoui, 2020; Bietti, 2022; Chen, 2023). These
authors showed that certain classes of CNNs are contained in the reproducing kernel Hilbert space (RKHS)
of a multilayer convolutional kernel representation. As such, stability metrics are estimated, based on the
RKHS norm which is difficult to control in practice. Kernel representations do not seem to suffer from
aliasing effects; this can be explained by the Gaussian pooling layers that have been employed instead of
max pooling: by discarding high-frequency information, shift invariance is preserved.

Finally, some papers studied stability of CNNs in a broader sense, measured in terms of Lipschitz continuity (Szegedy et al., 2014; Balan et al., 2018; Virmaux & Scaman, 2018; Pérez et al., 2020; Zou et al., 2020; Gupta et al., 2022; Zühlke & Kudenko, 2024). However, the Lipschitz bounds, which have been obtained theoretically, are generally several orders of magnitude higher than empirical results. This discrepancy may be due to the fact that these bounds were obtained for generic situations and represent overly conservative worst-case scenarios, rather than typical real-world situations. Furthermore, the specific case of convolutions with band-pass Gabor-like filters have been overlooked, except for Pérez et al. (2020).

In summary, we have identified the following blind spots in the literature, regarding the topic of studying
shift invariance in CNNs.

* The effect of the max pooling operator on network stability under small input shifts has not been investigated, particularly when used in combination with Gabor-like convolutions.

* While the shift invariance of CNNs tends to increase with network depth in the continuous framework, in the discrete case, the presence of subsampled convolutions with oriented band-pass filters
$W \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ denote a two-dimensional band-pass, oriented and analytic Gabor-like filter, for which a formal definition will be provided in (2.5). We first consider an operator, referred to as real-max-pooling (RMax), which computes the subsampled cross-correlation between an input image $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and the real part of W; then calculates the maximum value over a sliding discrete grid:

$$ (1.1) \qquad U_{m,q}^{\max}[W] : X \mapsto \mathrm{MaxPool}_q((X * \overline{\mathrm{Re} W}) \downarrow m), $$

where $m \in \mathbb{N} \setminus \{0\}$ denotes a subsampling factor, $\overline{V}$ denotes the “flipped” sequence for any given $V \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ or $l_{\mathbb{C}}^2(\mathbb{Z}^2)$, satisfying, for any $n \in \mathbb{Z}^2$,

$$ (1.2) \qquad \overline{V}[n] := V[-n], $$

and $*$, $\downarrow$ respectively refer to the convolution and subsampling operations, defined by

$$ (1.3) \qquad (X * \overline{V})[n] := \sum_{\mathbf{p} \in \mathbb{Z}^2} X[\mathbf{p}] \overline{V}[n - \mathbf{p}] \quad \text{and} \quad (Y \downarrow m)[n] := Y[mn]. $$

In the above expression, $\mathrm{MaxPool}_q$ selects the maximum value over a sliding grid of size $(2q + 1) \times (2q + 1)$, with a subsampling factor of 2. More formally, for any $Y \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and any $n \in \mathbb{Z}^2$,

$$ (1.4) \qquad \mathrm{MaxPool}_q(Y)[n] := \max_{\|\mathbf{p}\|_{\infty} \le q} Y[2n + \mathbf{p}]. $$

On the other hand, we consider an operator, referred to as complex-modulus (CMod), computing the modulus of subsampled cross-correlation between X and W:

$$ (1.5) \qquad U_m^{\mathrm{mod}}[W] : X \mapsto |(X * \overline{W}) \downarrow (2m)|. $$

First, we show that, under the Gabor hypothesis, CMod is stable with respect to small input shifts. We then establish conditions on the filter's frequency and orientation under which CMod and RMax produce comparable outputs:

$$ (1.6) \qquad U_m^{\mathrm{mod}}[W](X) \approx U_{m,q}^{\max}[W](X). $$

We deduce a measure of shift invariance for RMax operators, which benefits from the stability of CMod. Next, we extend our results to multichannel operators (i.e., applied on RGB input images), such as implemented in conventional CNN architectures. Our framework therefore provides a theoretical grounding to study these networks.

**Remark 1.1.** In the above definitions, cross-correlations are computed with a subsampling factor which is twice larger for CMod, compared to RMax. However, since max pooling is also computed with subsampling, both operators have the same subsampling factor of $2m$.

We assess our theoretical findings on a deterministic setting based on the dual-tree complex wavelet packet transform (DT-CWPT), a particular case of discrete Gabor-like decomposition with perfect reconstruction properties [4]. DT-CWPT spawns a set of convolution kernels
can lead to aliasing artifacts. To our knowledge, the literature lacks theoretical studies that take these aliasing effects into account.

* Although extensive studies have been conducted on complex-valued convolutions followed by modulus, a link is missing to extend these results to standard CNNs, which implement real-valued convolutions and spatial pooling operators.

All these points have been tackled in the present paper, from both theoretical and empirical perspectives.

## 1.3 Paper Outline

In what follows, $l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and $l_{\mathbb{C}}^2(\mathbb{Z}^2)$ represent the discrete spaces of square-summable two-dimensional sequences with values in $\mathbb{R}$ and $\mathbb{C}$, respectively. Let $W \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ denote a two-dimensional band-pass, oriented and analytic Gabor-like filter, for which a formal definition will be provided in (11). We first consider an operator, referred to as real-max-pooling ($\mathbb{R}$Max), which computes the subsampled cross-correlation between an input image $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and the real part of $W$; then calculates the maximum value over a sliding discrete grid:

$$ U_{m,q}^{\max}[W] : X \mapsto \text{MaxPool}_q \left( (X * \overline{\text{Re} W}) \downarrow m \right), \quad (1) $$

where $m \in \mathbb{N} \setminus \{0\}$ denotes a subsampling factor, $\overline{V}$ denotes the “flipped” sequence for any given $V \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ or $l_{\mathbb{C}}^2(\mathbb{Z}^2)$, satisfying, for any $\mathbf{n} \in \mathbb{Z}^2$,

$$ \overline{V}[\mathbf{n}] := V[-\mathbf{n}], \quad (2) $$

and $*$, $\downarrow$ respectively refer to the convolution and subsampling operations, defined by

$$ (X * \overline{V})[\mathbf{n}] := \sum_{\mathbf{p} \in \mathbb{Z}^2} X[\mathbf{p}] \overline{V}[\mathbf{n} - \mathbf{p}] \quad \text{and} \quad (Y \downarrow m)[\mathbf{n}] := Y[m\mathbf{n}]. \quad (3) $$

In the above expression, $\text{MaxPool}_q$ selects the maximum value over a sliding grid of size $(2q+1) \times (2q+1)$, with a subsampling factor of 2. More formally, for any $Y \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and any $\mathbf{n} \in \mathbb{Z}^2$,

$$ \text{MaxPool}_q(Y)[\mathbf{n}] := \max_{\|\mathbf{p}\|_{\infty} \le q} Y[2\mathbf{n} + \mathbf{p}]. \quad (4) $$

On the other hand, we consider an operator, referred to as complex-modulus (CMod), computing the modulus of subsampled cross-correlation between X and W:

$$ U_m^{\text{mod}}[W] : X \mapsto |(X * \overline{W}) \downarrow (2m)|. \quad (5) $$

First, we show that, under the Gabor hypothesis, CMod is stable with respect to small input shifts (Section 2). We then establish conditions on the filter's frequency and orientation under which CMod and $\mathbb{R}$Max produce comparable outputs (Section 3):

$$ U_m^{\text{mod}}[W](X) \approx U_{m,q}^{\max}[W](X). \quad (6) $$

We deduce a measure of shift invariance for $\mathbb{R}$Max operators, which benefits from the stability of CMod (Section 4). Next, we extend our results to multichannel operators (i.e., applied on RGB input images), such as implemented in conventional CNN architectures (Section 5). Our framework therefore provides a theoretical grounding to study these networks.

**Remark 1.** In the above definitions, cross-correlations are computed with a subsampling factor which is twice larger for CMod, compared to $\mathbb{R}$Max. However, since max pooling is also computed with subsampling, both operators have the same subsampling factor of $2m$.

Finally, in Section 6, we assess our theoretical findings on a deterministic setting based on the dual-tree complex wavelet packet transform (DT-CWPT), a particular case of discrete Gabor-like decomposition with perfect reconstruction properties (Bayram & Selesnick, 2008). DT-CWPT spawns a set of convolution kernels which tile the Fourier domain into square regions of identical size. Such kernels possess characteristics that are comparable to those found in the first convolution layer of CNNs after training with image datasets such
which tile the Fourier domain into square regions of identical size. Such kernels possess
characteristics that are comparable to those found in the first convolution layer of CNNs after
training with image datasets such as ImageNet [39]. More specifically, given an input image,
we compute the mean square error between the outputs of CMod and RMax, for each wavelet
packet filter. We then observe that shift invariance, when measured on RMax feature maps,
is nearly achieved when they remain close to CMod outputs. We therefore establish a domain
of validity for shift invariance of the RMax operator.

Prior to this work, we presented a preliminary study [24], where we experimentally showed
that an operator based on the real part of DT-CWPT can mimic the behavior of the first
convolution layer with fewer parameters, while keeping the network's predictive power. Our
model was solely based on real-valued filters, which are prone to aliasing. Yet, it produced
relatively stable outputs when compared with other models based on the standard, poorly-
oriented wavelet packet transform. At the same time, we became aware of a preliminary
work by Waldspurger [50, pp. 190-191], suggesting a potential connection between the com-
binations "real wavelet transform → max pooling" on the one hand and "complex wavelet
transform → modulus" on the other hand. Following this idea, we decided to study whether
invariance properties of complex moduli could somehow be captured by the max pooling oper-
ator. As shown in the present paper, Waldspurger's work does not fully extend to the discrete
framework. We address this issue by adopting a probabilistic point of view.

**2. Shift Invariance of CMod Outputs.** The primary goal of this paper is to theoretically establish conditions for near-shift invariance at the output of the first max pooling layer. In this section, we start by proving shift invariance of CMod operators. Then, in section 3, we establish conditions under which RMax and CMod produce closely related outputs. Finally, in section 4, we derive a probabilistic measure of shift invariance for RMax.

**2.1. Notations.** The complex conjugate of any number $z \in \mathbb{C}$ is denoted by $z^*$. For any $p \in \mathbb{R}_{>0} \cup \{\infty\}$, $x \in \mathbb{R}^2$ and $r \in \mathbb{R}_+$, we denote by $B_p(x, r) \subset \mathbb{R}^2$ the closed $l^p$-ball with center $x$ and radius $r$. When $x = 0$, we write $B_p(r)$.

*Continuous Framework.* Considering a measurable subset $E$ of $\mathbb{R}^2$, we denote by $L_{\mathbb{C}}^2(E)$ the Hilbert space of square-integrable functions $F: E \to \mathbb{C}$. Whenever we talk about equality in $L_{\mathbb{C}}^p(E)$ or inclusion in $E$, it shall be understood as “almost everywhere with respect to the Lebesgue measure.” Besides, we denote by $L_{\mathbb{R}}^2(E) \subset L_{\mathbb{C}}^2(E)$ the subspace of real-valued functions. For any $F \in L_{\mathbb{C}}^2(\mathbb{R}^2)$, $\bar{F}$ denotes its flipped version: $\bar{F}(x) := F(-x)$.

The 2D Fourier transform of any $F \in L_{\mathbb{C}}^2(\mathbb{R}^2)$ is denoted by $\hat{F} \in L_{\mathbb{C}}^2(\mathbb{R}^2)$, such that

$$
(2.1) \qquad \forall \nu \in \mathbb{R}^2, \hat{F}(\nu) := \iint_{\mathbb{R}^2} F(x)e^{-i\langle \nu, x \rangle} \, d^2x.
$$

For any $\varepsilon > 0$ and $\nu \in \mathbb{R}^2$, we denote by $\mathcal{V}(\nu, \varepsilon) \subset L_{\mathbb{C}}^2(\mathbb{R}^2)$ the set of functions whose Fourier transform is supported in a square region of size $\varepsilon \times \varepsilon$ centered in $\nu$:

$$
(2.2) \qquad \mathcal{V}(\nu, \varepsilon) := \left\{ \Psi \in L_{\mathbb{C}}^2(\mathbb{R}^2) \mid \operatorname{supp} \hat{\Psi} \subset B_{\infty}(\nu, \varepsilon/2) \right\}.
$$

$\nu$ and $\varepsilon$ are respectively referred to as characteristic frequency and bandwidth. Finally, for any $h \in \mathbb{R}^2$, we consider the translation operator, denoted by $\mathcal{T}_h$, defined by

$$
(2.3) \qquad \mathcal{T}_h F : x \mapsto F(x - h).
$$
as ImageNet (Russakovsky et al., 2015). More specifically, given an input image, we compute the mean square error between the outputs of CMod and RMax, for each wavelet packet filter. We then observe that shift invariance, when measured on RMax feature maps, is nearly achieved when they remain close to CMod outputs. We therefore establish a domain of validity for shift invariance of the RMax operator.

This work builds upon an idea sketched by Waldspurger (2015, pp. 190–191), which suggests a potential connection between the combinations “real wavelet transform → max pooling” on the one hand and “complex wavelet transform → modulus” on the other hand. Building on this idea, we investigated whether the invariance properties of complex moduli could be captured by the max pooling operator. However, as shown in this paper, these principles do not fully extend to the discrete framework. To address this limitation, we adopted a Bayesian point of view.

## 2 Shift Invariance of CMod Outputs

The primary goal of this paper is to theoretically establish conditions for near-shift invariance at the output of the first max pooling layer. In this section, we start by proving shift invariance of CMod operators. Then, in Section 3, we establish conditions under which RMax and CMod produce closely related outputs. Finally, in Section 4, we derive a probabilistic measure of shift invariance for RMax.

### 2.1 Notations

The complex conjugate of any number $z \in \mathbb{C}$ is denoted by $z^*$. For any $p \in \mathbb{R}_{>0} \cup \{\infty\}$, $x \in \mathbb{R}^2$ and $r \in \mathbb{R}_+$, we denote by $B_p(x, r) \subset \mathbb{R}^2$ the closed $l^p$-ball with center $x$ and radius $r$. When $x = 0$, we write $B_p(r)$.

**Continuous Framework** Considering a measurable subset $E$ of $\mathbb{R}^2$, we denote by $L_{\mathbb{C}}^2(E)$ the Hilbert space of square-integrable functions $F : E \to \mathbb{C}$. Whenever we talk about equality in $L_{\mathbb{C}}^p(E)$ or inclusion in $E$, it shall be understood as “almost everywhere with respect to the Lebesgue measure.” Additionally, we denote by $L_{\mathbb{R}}^2(E) \subset L_{\mathbb{C}}^2(E)$ the subspace of real-valued functions. For any $F \in L_{\mathbb{C}}^2(\mathbb{R}^2)$, $\bar{F}$ denotes its flipped version: $\bar{F}(x) := F(-x)$.

The 2D Fourier transform of any $F \in L_{\mathbb{C}}^2(\mathbb{R}^2)$ is denoted by $\hat{F} \in L_{\mathbb{C}}^2(\mathbb{R}^2)$, such that

$$ \forall \nu \in \mathbb{R}^2, \hat{F}(\nu) := \iint_{\mathbb{R}^2} F(x)e^{-i\langle \nu, x \rangle} d^2x. \qquad (7) $$

For any $\varepsilon > 0$ and $\nu \in \mathbb{R}^2$, we denote by $\mathcal{V}(\nu, \varepsilon) \subset L_{\mathbb{C}}^2(\mathbb{R}^2)$ the set of functions whose Fourier transform is supported in a square region of size $\varepsilon \times \varepsilon$ centered in $\nu$:

$$ \mathcal{V}(\nu, \varepsilon) := \{ \Psi \in L_{\mathbb{C}}^2(\mathbb{R}^2) \mid \text{supp} \hat{\Psi} \subset B_{\infty}(\nu, \varepsilon/2) \}. \qquad (8) $$

$\nu$ and $\varepsilon$ are respectively referred to as characteristic frequency and bandwidth. Finally, for any $h \in \mathbb{R}^2$, we consider the translation operator, denoted by $\mathcal{T}_h$, defined by

$$ \mathcal{T}_h F : x \mapsto F(x - h). \qquad (9) $$

**Discrete Framework** We denote by $l_{\mathbb{C}}^2(\mathbb{Z}^2)$ the space of 2D complex-valued square-summable sequences, represented by straight capital letters. Indexing is made between square brackets: $\forall X \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$, $\forall n \in \mathbb{Z}^2$, $X[n] \in \mathbb{C}$, and we denote by $l_{\mathbb{R}}^2(\mathbb{Z}^2) \subset l_{\mathbb{C}}^2(\mathbb{Z}^2)$ the subset of real-valued sequences. For any $V \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$, $\bar{V}$ denotes its “flipped” version as defined in (2). The convolution and subsampling operators, respectively denoted by $*$ and $\downarrow$, are defined in (3). 2D images, feature maps and convolution kernels are considered as elements of $l_{\mathbb{C}}^2(\mathbb{Z}^2)$. Additionally, multichannel arrays of 2D sequences are denoted by bold straight capital letters, for instance: $\mathbf{X} := (X_k)_{k \in \{0, ..., K-1\}}$. Note that indexing starts at 0 to comply with practical implementations.

The 2D discrete-time Fourier transform of any $X \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$, denoted by $\hat{X} \in L_{\mathbb{C}}^2([-\pi, \pi]^2)$, is defined by

$$ \forall \theta \in [-\pi, \pi]^2, \hat{X}(\theta) := \sum_{n \in \mathbb{Z}^2} X[n]e^{-i\langle \theta, n \rangle}. \qquad (10) $$
**Discrete Framework.** We denote by $l_{\mathbb{C}}^2(\mathbb{Z}^2)$ the space of 2D complex-valued square-summable sequences, represented by straight capital letters. Indexing is made between square brackets: $\forall X \in l_{\mathbb{C}}^2(\mathbb{Z}^2), \forall n \in \mathbb{Z}^2, X[n] \in \mathbb{C}$, and we denote by $l_{\mathbb{R}}^2(\mathbb{Z}^2) \subset l_{\mathbb{C}}^2(\mathbb{Z}^2)$ the subset of real-valued sequences. For any $V \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$, $\overline{V}$ denotes its “flipped” version as defined in (1.2). The convolution and subsampling operators, respectively denoted by $*$ and $\downarrow$, are defined in (1.3). 2D images, feature maps and convolution kernels are considered as elements of $l_{\mathbb{C}}^2(\mathbb{Z}^2)$. Besides, multichannel arrays of 2D sequences are denoted by bold straight capital letters, for instance: $\mathbf{X} := (X_k)_{k \in \{0..K-1\}}$. Note that indexing starts at 0 to comply with practical implementations.

The 2D discrete-time Fourier transform of any $X \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$, denoted by $\hat{\mathbf{X}} \in L_{\mathbb{C}}^2([-\pi, \pi]^2)$, is defined by

$$ (2.4) \qquad \forall \theta \in [-\pi, \pi]^2, \hat{\mathbf{X}}(\theta) := \sum_{n \in \mathbb{Z}^2} X[n]e^{-i\langle\theta, n\rangle}. $$

For any $\kappa \in ]0, 2\pi]$ and $\theta \in B_{\infty}(\pi)$, we denote by $\mathcal{J}(\theta, \kappa) \subset l_{\mathbb{C}}^2(\mathbb{Z}^2)$ the set of 2D sequences whose Fourier transform is supported in a square region of size $\kappa \times \kappa$ centered in $\theta$:

$$ (2.5) \qquad \mathcal{J}(\theta, \kappa) := \{W \in l_{\mathbb{C}}^2(\mathbb{Z}^2) \mid \text{supp } \hat{W} \subset B_{\infty}(\theta, \kappa/2)\}. $$

As in the discrete framework, $\theta$ and $\kappa$ are respectively referred to as characteristic frequency and bandwidth. The elements of $\mathcal{J}(\theta, \kappa)$ are designated as Gabor-like filters.

**Remark 2.1.** The support $B_{\infty}(\theta, \kappa/2)$ actually lives in the quotient space $[-\pi, \pi]^2/(2\pi \mathbb{Z}^2)$. Consequently, when $\theta$ is close to an edge, a fraction of this region is located at the far end of the frequency domain. From now on, the choice of $\theta$ and $\kappa$ is implicitly assumed to avoid such a situation.

**2.2. Intuition.** In many CNNs for computer vision, input images are first transformed through subsampled (or strided) convolutions. For instance, in AlexNet, convolution kernels are of size $11 \times 11$ and the subsampling factor is equal to 4. Figure 1 displays the corresponding kernels after training with ImageNet. This linear transform is generally followed by rectified linear unit (ReLU) and max pooling.

We can observe that many kernels display oscillating patterns with well-defined orientations (Gabor-like filters). We denote by $V \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ one of these “well-behaved” filters. Its Fourier spectrum roughly consists in two bright spots which are symmetric with respect to the origin.<fnref n="1" /> Now, we consider a complex-valued companion $W \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ such that

$$ (2.6) \qquad \hat{W}(\omega) := (1 + \operatorname{sgn}\langle\omega, u\rangle) \cdot \hat{V}(\omega) \quad \forall \omega \in [-\pi, \pi]^2, $$

where $u$ denotes a unit vector orthogonal to the filter's orientation.

We can show that $V$ is the real part of $W$, and that $W = V + i\mathcal{H}(V)$, where $\mathcal{H}$ denotes the two-dimensional Hilbert transform as introduced by Havlicek et al. [16]. It satisfies

$$ (2.7) \qquad \hat{\mathcal{H}}(V)(\omega) := -i \operatorname{sgn}\langle\omega, u\rangle \hat{V}(\omega). $$

<footnote n="1">¹Actually, the Fourier transform of any real-valued sequence is centrally symmetric: $\hat{V}(-\omega) = \hat{V}(\omega)^*$. The specificity of well-oriented filters lies in the concentration of their power spectrum around two precise locations.</footnote>
<figure><img src="image_2.png" /><figcaption>Figure 1. Spatial (left) and Fourier (right) representations of convolution kernels in the first layer of AlexNet, after training with ImageNet ILSVRC 2012-2017 (Russakovsky et al., 2015). Each kernel connects the 3 RGB input channels to one of the 64 output channels.</figcaption></figure>

For any $\kappa \in ]0, 2\pi]$ and $\theta \in B_{\infty}(\pi)$, we denote by $\mathcal{J}(\theta, \kappa) \subset l_{\mathbb{C}}^2(\mathbb{Z}^2)$ the set of 2D sequences whose Fourier transform is supported in a square region of size $\kappa \times \kappa$ centered in $\theta$:

$$ \mathcal{J}(\theta, \kappa) := \{ W \in l_{\mathbb{C}}^2(\mathbb{Z}^2) \mid \text{supp} \hat{W} \subset B_{\infty}(\theta, \kappa/2) \}. \quad (11) $$

As in the continuous framework, $\theta$ and $\kappa$ are respectively referred to as characteristic frequency and bandwidth. The elements of $\mathcal{J}(\theta, \kappa)$ are designated as Gabor-like filters.

**Remark 2.** The support $B_{\infty}(\theta, \kappa/2)$ actually lives in the quotient space $[-\pi, \pi]^2/(2\pi\mathbb{Z}^2)$. Consequently, when $\theta$ is close to an edge, a fraction of this region is located at the far end of the frequency domain. From now on, the choice of $\theta$ and $\kappa$ is implicitly assumed to avoid such a situation.

## 2.2 Intuition

In many CNNs for computer vision, input images are first transformed through subsampled (or strided) convolutions. For instance, in AlexNet, convolution kernels are of size 11 × 11 and the subsampling factor is equal to 4. Figure 1 displays the corresponding kernels after training with ImageNet. This linear transform is generally followed by rectified linear unit (ReLU) and max pooling.

We can observe that many kernels display oscillating patterns with well-defined orientations (Gabor-like filters). We denote by $V \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ one of these “well-behaved” filters. Its Fourier spectrum roughly consists in two bright spots which are symmetric with respect to the origin.<fnref n="1" /> Next, we consider a complex-valued companion $W \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ such that

$$ \hat{W}(\omega) := (1 + \operatorname{sgn}\langle\omega, u\rangle) \cdot \hat{V}(\omega) \quad \forall \omega \in [-\pi, \pi]^2, \quad (12) $$

where $u$ denotes a unit vector orthogonal to the filter’s orientation.

We can show that $V$ is the real part of $W$, and that $W = V + i\mathcal{H}(V)$, where $\mathcal{H}$ denotes the two-dimensional Hilbert transform as introduced by Havlicek et al. (1997). It satisfies

$$ \hat{\mathcal{H}}(\hat{V})(\omega) := -i \operatorname{sgn}\langle\omega, u\rangle \hat{V}(\omega). \quad (13) $$

As a consequence, $\hat{W}$ is equal to $2\hat{V}$ on one half of the Fourier domain, and 0 on the other half. Therefore, only one bright spot remains in the spectrum. We refer the reader to fig. 2 for visual example of complex-valued Gabor-like filter. It turns out that such complex filters with high frequency resolution produce stable signal representations, as we will see in Section 2. In the subsequent sections, we then wonder whether this property is kept when considering the max pooling of real-valued convolutions.

<footnote n="1">¹Actually, the Fourier transform of any real-valued sequence is centrally symmetric: $\hat{V}(-\omega) = \hat{V}(\omega)^*$. The specificity of well-oriented filters lies in the concentration of their power spectrum around two precise locations.</footnote>
<figure><img src="image_3.png" /><figcaption>Figure 1. Spatial (left) and Fourier (right) representations of convolution kernels in the first layer of AlexNet, after training with ImageNet ILSVRC 2012-2017 [39]. Each kernel connects the 3 RGB input channels to one of the 64 output channels.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 2. (a), (b): Real and imaginary parts of a Gabor-like filter W as defined in (2.6). (c), (d): Magnitude spectra (modulus of the Fourier transform) of V and W, respectively.</figcaption></figure>

As a consequence, $\hat{W}$ is equal to $2\hat{V}$ on one half of the Fourier domain, and 0 on the other half. Therefore, only one bright spot remains in the spectrum. We refer the reader to Figure 2 for visual example of complex-valued Gabor-like filter. It turns out that such complex filters with high frequency resolution produce stable signal representations, as we will see in section 2. In the subsequent sections, we then wonder whether this property is kept when considering the max pooling of real-valued convolutions.

In what follows, W will be referred to as a discrete Gabor-like filter, and the coefficients resulting from the convolution with W will be referred to as discrete Gabor-like coefficients. The aim of this section is to show that, under the Gabor hypothesis on the convolution kernels $W \in l_2^2(\mathbb{Z}^2)$, CMod is nearly shift-invariant. To clarify, we establish that

$$ (2.8) \qquad U_m^{\text{mod}}[W](X) \approx U_m^{\text{mod}}[W](\mathcal{T}_u X), $$

for “small” translation vectors $u \in \mathbb{R}^2$, where a formal definition of the translation operator will be defined in (2.34). This result is hinted by Kingsbury and Magarey [20] but not formally proven.
<figure><img src="image_3.png" /><figcaption>Figure 2. (a), (b): Real and imaginary parts of a Gabor-like filter W as defined in (12). (c), (d): Magnitude spectra (modulus of the Fourier transform) of V and W, respectively.</figcaption></figure>

In what follows, W will be referred to as a discrete Gabor-like filter, and the coefficients resulting from the convolution with W will be referred to as discrete Gabor-like coefficients. The aim of this section is to show that, under the Gabor hypothesis on the convolution kernels $W \in L_{\mathbb{C}}^2(\mathbb{Z}^2)$, $\mathbb{C}$Mod is nearly shift-invariant. To clarify, we establish that

$$U_m^{\text{mod}}[W](X) \approx U_m^{\text{mod}}[W](\mathcal{T}_u X), \qquad (14)$$

for “small” translation vectors $u \in \mathbb{R}^2$, where a formal definition of the translation operator will be defined in (28). This result is hinted by Kingsbury & Magarey (1998) but not formally proven.

## 2.3 Continuous Framework

We introduce several results regarding functions defined on the continuous space $\mathbb{R}^2$. Near-shift invariance on discrete 2D sequences will then be derived from these results by taking advantage of sampling theorems. Lemma 1 below is adapted from Waldspurger (2015, pp. 190–191).

**Lemma 1.** Given $\varepsilon > 0$ and $\nu \in \mathbb{R}^2$, let $\Psi \in \mathcal{V}(\nu, \varepsilon)$ denote a complex-valued filter such as defined in (8). Next, for any real-valued function $F \in L_{\mathbb{R}}^2(\mathbb{R}^2)$, we consider the complex-valued function $F_0 \in L_{\mathbb{C}}^2(\mathbb{R}^2)$ defined by

$$F_0: \boldsymbol{x} \mapsto (F * \bar{\Psi})(\boldsymbol{x}) e^{i\langle\boldsymbol{\nu}, \boldsymbol{x}\rangle}. \qquad (15)$$

Then $F_0$ is low-frequency. Specifically,

$$\text{supp } \widehat{F}_0 \subset B_{\infty}(\varepsilon/2). \qquad (16)$$

*Proof.* See Appendix A.1. □

On the other hand, the following proposition provides a shift invariance bound for low-frequency functions such as introduced above.

**Proposition 1.** For any $F_0 \in L_{\mathbb{R}}^2(\mathbb{R}^2)$ such that $\text{supp } \widehat{F}_0 \subset B_{\infty}(\varepsilon/2)$, and any $h \in \mathbb{R}^2$,

$$\| \mathcal{T}_h F_0 - F_0 \|_{L^2} \le \alpha(\varepsilon h) \| F_0 \|_{L^2}, \qquad (17)$$

where we have defined

$$\alpha: \boldsymbol{\tau} \mapsto \frac{\|\boldsymbol{\tau}\|_1}{2}. \qquad (18)$$

*Proof.* See Appendix A.2. □

## 2.4 Adaptation to Discrete 2D Sequences

Given $\kappa \in ]0, 2\pi]$ and $\theta \in B_{\infty}(\pi)$, let $W \in \mathcal{J}(\theta, \kappa)$ denote a discrete Gabor-like filter such as defined in (11). For any image $X \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ with finite support and any subsampling factor $m \in \mathbb{N} \setminus \{0\}$, we express $(X * \bar{W}) \downarrow m$ using the continuous framework introduced above, and derive an invariance formula.

For any sampling interval $s \in \mathbb{R}_{>0}$, let $\Phi^{(s)} \in L_{\mathbb{R}}^2(\mathbb{R}^2)$ denote the Shannon scaling function parameterized by $s$, such that

$$\widehat{\Phi}^{(s)} := s \mathbb{1}_{B_{\infty}(\pi/s)}. \qquad (19)$$
**2.3. Continuous Framework.** We introduce several results regarding functions defined on the continuous space $\mathbb{R}^2$. Near-shift invariance on discrete 2D sequences will then be derived from these results by taking advantage of sampling theorems. Lemma 2.2 below is adapted from Waldspurger [50, pp. 190–191].

**Lemma 2.2.** Given $\varepsilon > 0$ and $\nu \in \mathbb{R}^2$, let $\Psi \in \mathcal{V}(\nu, \varepsilon)$ denote a complex-valued filter such as defined in (2.2). Now, for any real-valued function $F \in L_{\mathbb{R}}^2(\mathbb{R}^2)$, we consider the complex-valued function $F_0 \in L_{\mathbb{C}}^2(\mathbb{R}^2)$ defined by

$$ (2.9) \qquad F_0 : \boldsymbol{x} \mapsto (F * \bar{\Psi})(\boldsymbol{x}) e^{i\langle\boldsymbol{\nu}, \boldsymbol{x}\rangle}. $$

Then $F_0$ is low-frequency. Specifically,

$$ (2.10) \qquad \operatorname{supp} \widehat{F_0} \subset B_{\infty}(\varepsilon/2). $$

*Proof.* Applying the Fourier transform on (2.9) yields, for any $\boldsymbol{\xi} \in \mathbb{R}^2$,

$$ (2.11) \qquad \widehat{F_0}(\boldsymbol{\xi}) = (\widehat{F * \Psi})(\boldsymbol{\xi} - \boldsymbol{\nu}) = \mathcal{T}_{\nu}(\widehat{F * \Psi})(\boldsymbol{\xi}). $$

By hypothesis on $\Psi$, we have

$$ (2.12) \qquad \operatorname{supp}(\widehat{F * \Psi}) \subset \operatorname{supp} \widehat{\Psi} \subset B_{\infty}(-\boldsymbol{\nu}, \varepsilon/2). $$

The translation operator $\mathcal{T}_{\nu}$ shifts the support with respect to $\nu$, which yields (2.10). ■

On the other hand, the following proposition provides a shift invariance bound for low-frequency functions such as introduced above.

**Proposition 2.3.** For any $F_0 \in L_{\mathbb{R}}^2(\mathbb{R}^2)$ such that $\operatorname{supp} \widehat{F_0} \subset B_{\infty}(\varepsilon/2)$, and any $\boldsymbol{h} \in \mathbb{R}^2$,

$$ (2.13) \qquad \|\mathcal{T}_{\boldsymbol{h}} F_0 - F_0\|_{L^2} \le \alpha(\varepsilon \boldsymbol{h}) \|F_0\|_{L^2}, $$

where we have defined

$$ (2.14) \qquad \alpha : \boldsymbol{\tau} \mapsto \frac{\|\boldsymbol{\tau}\|_1}{2}. $$

*Proof.* Using the 2D Plancherel formula, we compute

$$ \begin{align*} \|\mathcal{T}_{\boldsymbol{h}} F_0 - F_0\|_{L^2}^2 &= \frac{1}{4\pi^2} \left\| \widehat{\mathcal{T}_{\boldsymbol{h}} F_0} - \widehat{F_0} \right\|_{L^2}^2 \\ &= \frac{1}{4\pi^2} \iint_{B_{\infty}(\varepsilon/2)} |\widehat{F_0}(\boldsymbol{\xi})|^2 |e^{-i\langle\boldsymbol{h}, \boldsymbol{\xi}\rangle} - 1|^2 d^2\boldsymbol{\xi} \\ &= \frac{1}{4\pi^2} \iint_{B_{\infty}(\varepsilon/2)} |\widehat{F_0}(\boldsymbol{\xi})|^2 (2 - 2 \cos \langle\boldsymbol{h}, \boldsymbol{\xi}\rangle) d^2\boldsymbol{\xi} \\ &\le \frac{1}{4\pi^2} \iint_{B_{\infty}(\varepsilon/2)} |\widehat{F_0}(\boldsymbol{\xi})|^2 |\langle\boldsymbol{h}, \boldsymbol{\xi}\rangle|^2 d^2\boldsymbol{\xi}, \end{align*} $$
This 2D function is a tensor product of scaled and normalized sinc functions. For any $n \in \mathbb{Z}^2$, we denote by $\Phi_n^{(s)}$ a shifted version of $\Phi^{(s)}$, satisfying

$$ \Phi_n^{(s)}(\boldsymbol{x}) := \Phi^{(s)}(\boldsymbol{x} - s\boldsymbol{n}). \qquad (20) $$

Then, $\{\Phi_n^{(s)}\}_{n \in \mathbb{Z}^2}$ is an orthonormal basis of

$$ \mathcal{V}^{(s)} := \{F \in L_{\mathbb{C}}^2(\mathbb{R}^2) \mid \text{supp } \hat{F} \subset B_{\infty}(\pi/s)\}. \qquad (21) $$

Then, using the notation introduced in (8), we have $\mathcal{V}^{(s)} = \mathcal{V}(0, 2\pi/s)$.

We now consider the following lemma.

**Lemma 2.** Let $s > 0$. For any $F \in \mathcal{V}^{(s)}$ and any $\boldsymbol{\xi} \in B_{\infty}(\pi/s)$, we have

$$ \hat{F}(\boldsymbol{\xi}) = s \hat{X}(s\boldsymbol{\xi}), \qquad (22) $$

where $X \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ is a uniform sampling of $F$, defined such that $X[n] := sF(sn)$, for any $n \in \mathbb{Z}^2$. Moreover, we have the following norm equality:

$$ \|F\|_{L^2} = \|X\|_2. \qquad (23) $$

*Proof.* See Appendix A.3. □

We then get the following proposition, which draws a bond between the discrete and continuous frameworks.

**Proposition 2.** Let $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ denote an input image with finite support, and $W \in \mathcal{J}(\boldsymbol{\theta}, \kappa)$. Considering a sampling interval $s \in \mathbb{R}_{>0}$, we define $F_X$ and $\Psi_W \in \mathcal{V}^{(s)}$ such that

$$ F_X := \sum_{n \in \mathbb{Z}^2} X[n] \Phi_n^{(s)} \quad \text{and} \quad \Psi_W := \sum_{n \in \mathbb{Z}^2} W[n] \Phi_n^{(s)}. \qquad (24) $$

Then,

$$ \Psi_W \in \mathcal{V}(\boldsymbol{\theta}/s, \kappa/s). \qquad (25) $$

Moreover, for all $n \in \mathbb{Z}$,

$$ X[n] = s F_X(sn); \quad W[n] = s \Psi_W(sn), \qquad (26) $$

and, for a given subsampling factor $m \in \mathbb{N} \setminus \{0\}$,

$$ ((X * \bar{W}) \downarrow m)[n] = (F_X * \bar{\Psi}_W)(msn). \qquad (27) $$

*Proof.* See Appendix A.4. □

Proposition 2 introduces a latent subspace of $L_{\mathbb{R}}^2(\mathbb{R}^2)$ from which input images are uniformly sampled. This allows us to define, for any $u \in \mathbb{R}^2$, a translation operator $T_u$ on discrete sequences, even if $u$ has non-integer values:

$$ T_u X[n] := s T_{su} F_X(sn), \qquad (28) $$

where $F_X$ is defined in (24). We can indeed show that this definition is independent from the choice of sampling interval $s > 0$. Moreover, given $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$, we have

$$ \forall \boldsymbol{p} \in \mathbb{Z}^2, \quad \mathcal{T}_p X[n] = X[n - \boldsymbol{p}]; \qquad (29) $$

$$ \forall \boldsymbol{u}, \boldsymbol{v} \in \mathbb{R}^2, \quad \mathcal{T}_u(\mathcal{T}_v X) = \mathcal{T}_{u+v} X, \qquad (30) $$

which shows that $T_u$ corresponds to the intuitive idea of a translation operator. Expressions (29) and (30) are direct consequence of the following lemma, which bonds the shift operator in the discrete and continuous frameworks.
because $\cos t \ge 1 - \frac{l^2}{2}$. Note that the integral is computed on a compact domain because, according to Lemma 2.2, $\text{supp} \hat{F}_0 \subset B_\infty(\varepsilon/2)$. Now, we use the Cauchy-Schwarz inequality to compute:

$$ \begin{aligned} \forall \xi \in B_\infty(\varepsilon/2), \quad & |\langle h, \xi \rangle| \le \|h\|_1 \cdot \|\xi\|_\infty \\ & \le \frac{\varepsilon}{2} \|h\|_1. \end{aligned} $$

Therefore,

$$ (2.15) \qquad \| \mathcal{T}_h F_0 - F_0 \|_{L^2}^2 \le \frac{\varepsilon}{4} \|h\|_1^2 \|F_0\|_{L^2}^2, $$

which yields the result.

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

**2.4. Adaptation to Discrete 2D Sequences.** Given $\kappa \in ]0, 2\pi]$ and $\theta \in B_\infty(\pi)$, let $W \in \mathcal{J}(\theta, \kappa)$ denote a discrete Gabor-like filter such as defined in (2.5). For any image $X \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ with finite support and any subsampling factor $m \in \mathbb{N}\setminus\{0\}$, we express $(X*\overline{W}) \downarrow m$ using the continuous framework introduced above, and derive an invariance formula.

For any sampling interval $s \in \mathbb{R}_{>0}$, let $\Phi^{(s)} \in L_{\mathbb{R}}^2(\mathbb{R}^2)$ denote the Shannon scaling function parameterized by $s$, such that

$$ (2.16) \qquad \widehat{\Phi^{(s)}} := s \mathbb{1}_{B_\infty(\pi/s)}. $$

This 2D function is a tensor product of scaled and normalized sinc functions. For any $n \in \mathbb{Z}^2$, we denote by $\Phi_n^{(s)}$ a shifted version of $\Phi^{(s)}$, satisfying

$$ (2.17) \qquad \Phi_n^{(s)}(x) := \Phi^{(s)}(x - sn). $$

Then, $\{\Phi_n^{(s)}\}_{n \in \mathbb{Z}^2}$ is an orthonormal basis of

$$ (2.18) \qquad \mathcal{V}^{(s)} := \{F \in L_{\mathbb{C}}^2(\mathbb{R}^2) \mid \text{supp } \hat{F} \subset B_\infty(\pi/s)\}. $$

Then, using the notation introduced in (2.2), we have $\mathcal{V}^{(s)} = \mathcal{V}(\mathbf{0}, 2\pi/s)$.

We now consider the following lemma.

**Lemma 2.4.** Let $s > 0$. For any $F \in \mathcal{V}^{(s)}$ and any $\xi \in B_\infty(\pi/s)$, we have

$$ (2.19) \qquad \hat{F}(\xi) = s \hat{X}(s\xi), $$

where $X \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ is a uniform sampling of $F$, defined such that $X[n] := sF(sn)$, for any $n \in \mathbb{Z}^2$. Besides, we have the following norm equality:

$$ (2.20) \qquad \|F\|_{L^2} = \|X\|_2. $$
**Lemma 3.** For any $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and any $\mathbf{u} \in \mathbb{R}^2$,

$$F_{\mathcal{T}_{\mathbf{u}}X} = \mathcal{T}_{\mathbf{u}}F_{X}. \qquad (31)$$

*Proof.* See Appendix A.5. □

We now consider the following corollary to Proposition 2.

**Corollary 1.** For any shift vector $\mathbf{u} \in \mathbb{R}^2$, we have

$$((\mathcal{T}_{\mathbf{u}}X * \overline{W}) \downarrow m)[\mathbf{n}] = (\mathcal{T}_{\mathbf{u}}F_{X} * \overline{\Psi}_{W})(ms\mathbf{n}). \qquad (32)$$

*Proof.* Applying (27) in Proposition 2 with $X \leftarrow \mathcal{T}_{\mathbf{u}}X$, we get

$$((\mathcal{T}_{\mathbf{u}}X * \overline{W}) \downarrow m)[\mathbf{n}] = (F_{\mathcal{T}_{\mathbf{u}}X} * \overline{\Psi}_{W})(ms\mathbf{n}), \qquad (33)$$

and Lemma 3 concludes the proof. □

## 2.5 Shift Invariance in the Discrete Framework

We consider the $\mathbb{C}$Mod operator defined in (5). For the sake of conciseness, in what follows we will write $U_m^{\text{mod}}$ instead of $U_m^{\text{mod}}[W]$, when no ambiguity is possible. First, we state the following lemma.

**Lemma 4.** For any input image $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ with finite support, and any Gabor-like filter $W \in \mathcal{J}(\theta, \kappa)$, we consider the low-frequency function

$$F_0 : \mathbf{x} \mapsto (F_X * \overline{\Psi}_W)(\mathbf{x}) e^{i\langle\theta/s, \mathbf{x}\rangle}, \qquad (34)$$

with $F_X$ and $\overline{\Psi}_W$ satisfying (24). If $\kappa \le \pi/m$, then

$$F_0 \in \mathcal{V}^{(s')}. \qquad (35)$$

Moreover, for any $\mathbf{h} \in \mathbb{R}^2$,

$$\sum_{\mathbf{n} \in \mathbb{Z}^2} |\mathcal{T}_{\mathbf{h}} F_0(s' \mathbf{n}) - F_0(s' \mathbf{n})|^2 = \frac{1}{s'^2} \|\mathcal{T}_{\mathbf{h}} F_0 - F_0\|_{L^2}^2, \qquad (36)$$

where we have denoted $s' := 2ms$. Finally,

$$\|U_m^{\text{mod}}X\|_2 = \frac{1}{s'} \|F_0\|_{L^2}. \qquad (37)$$

*Proof.* See Appendix A.6. □

We are now ready to state the main result about shift invariance of $\mathbb{C}$Mod outputs.

**Theorem 1** (Shift invariance of $\mathbb{C}$Mod). Let $W \in \mathcal{J}(\theta, \kappa)$ denote a discrete Gabor-like filter and $m \in \mathbb{N} \setminus \{0\}$ denote a subsampling factor. Then, under the following condition:

$$\kappa \le \pi/m, \qquad (38)$$

we have, for any input image $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ with finite support and any translation vector $\mathbf{u} \in \mathbb{R}^2$,

$$\|U_m^{\text{mod}}(\mathcal{T}_{\mathbf{u}}X) - U_m^{\text{mod}}X\|_2 \le \alpha(\kappa\mathbf{u}) \|U_m^{\text{mod}}X\|_2, \qquad (39)$$

where $\alpha$ has been defined in (18).

*Proof.* See Appendix A.7. □
*Proof.* Since $F \in \mathcal{V}^{(s)}$, the two-dimensional version of Shannon's sampling theorem [27, Theorem 3.11, p. 81] yields

$$ (2.21) \qquad F = \sum_{n \in \mathbb{Z}^2} X[n] \Phi_n^{(s)}, \quad \text{and} \quad \hat{F} = \sum_{n \in \mathbb{Z}^2} X[n] \widehat{\Phi_n^{(s)}}. $$

Besides, using (2.16), we can show that, for any $\xi \in B_\infty(\pi/s)$,

$$ (2.22) \qquad \widehat{\Phi_n^{(s)}}(\xi) = \widehat{\Phi^{(s)}}(\xi) e^{-i\langle s\xi, n \rangle} = s e^{-i\langle s\xi, n \rangle}. $$

Therefore, plugging (2.22) into (2.21) proves (2.19).

Then, by combining (2.19) with the Plancherel formula, we get

$$ \begin{aligned} \|F\|_{L^2}^2 &= \frac{1}{4\pi^2} \|\hat{F}\|_{L^2}^2 \\ &= \frac{1}{4\pi^2} \iint_{B_\infty(\pi/s)} |\hat{F}(\xi)|^2 d^2\xi \\ &= \frac{1}{4\pi^2} \iint_{B_\infty(\pi/s)} |s\hat{X}(s\xi)|^2 d^2\xi. \end{aligned} $$

The integral is performed on $B_\infty(\pi/s)$ because $F \in \mathcal{V}^{(s)}$. Then, by applying the change of variable $\xi' \leftarrow s\xi$, we get

$$ \begin{aligned} \|F\|_{L^2}^2 &= \frac{1}{4\pi^2} \iint_{B_\infty(\pi)} |\hat{X}(\xi')|^2 d^2\xi' \\ &= \frac{1}{4\pi^2} \|\hat{X}\|_{L^2}^2 = \|X\|_2^2, \end{aligned} $$

hence (2.20), which concludes the proof.

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

We then get the following proposition, which draws a bond between the discrete and continuous frameworks.

**Proposition 2.5.** Let $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ denote an input image with finite support, and $W \in \mathcal{J}(\theta, \kappa)$. Considering a sampling interval $s \in \mathbb{R}_{>0}$, we define $F_X$ and $\Psi_W \in \mathcal{V}^{(s)}$ such that

$$ (2.23) \qquad F_X := \sum_{n \in \mathbb{Z}^2} X[n] \Phi_n^{(s)} \quad \text{and} \quad \Psi_W := \sum_{n \in \mathbb{Z}^2} W[n] \Phi_n^{(s)}. $$

Then,

$$ (2.24) \qquad \Psi_W \in \mathcal{V}(\theta/s, \kappa/s). $$

Moreover, for all $n \in \mathbb{Z}$,

$$ (2.25) \qquad X[n] = s F_X(sn); \quad W[n] = s \Psi_W(sn), $$

and, for a given subsampling factor $m \in \mathbb{N} \setminus \{0\}$,

$$ (2.26) \qquad ((X * \overline{W}) \downarrow m)[n] = (F_X * \overline{\Psi}_W)(msn). $$
Interestingly, the reference value used in Theorem 1, i.e., $||U_m^{\text{mod}}X||_2$, is fully shift-invariant, as stated in the following proposition.

**Proposition 3.** Let $W \in \mathcal{J}(\theta, \kappa)$ and $m \in \mathbb{N} \setminus \{0\}$. Under condition (38), we have, for any $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and any $u \in \mathbb{R}^2$,

$$ ||U_m^{\text{mod}}(\mathcal{T}_u X)||_2 = ||U_m^{\text{mod}}X||_2. \qquad (40) $$

*Proof.* See Appendix A.8. □

# 3 From CMod to RMax

CMod operators are found in ScatterNets and complex-valued convolutional networks (Tygert et al., 2016). However, they are absent from conventional, freely-trained CNN architectures. Therefore, Theorem 1 cannot be applied as is. Instead, the first convolution layer contains real-valued kernels, and is generally followed by ReLU and max pooling. As shown in Section 5, this process can be described with RMax operators, such as defined in (1).

As explained in Section 1.1, an important number of trained convolution kernels exhibit oscillating patterns with well-defined frequencies and orientations. To elaborate, let $V \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ denote such a trained kernel, and consider $W \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ as the complex-valued companion of $V$ satisfying (12). Then, $W$ has its energy concentrated in a small region of the Fourier domain. We thus state the hypotheses that $W \in \mathcal{J}(\theta, \kappa)$ (11) for a certain value of $\theta \in [-\pi, \pi]^2$ and $\kappa \in ]0, 2\pi]$. For the sake of conciseness, from now on we write $U_{m,q}^{\text{max}}$ instead of $U_{m,q}^{\text{max}}[W]$, when no ambiguity is possible. In what follows, we establish conditions on $W$ under which CMod (5) and RMax (1) operators produce comparable outputs. The final goal, achieved in Section 4, is to provide a shift invariance bound for RMax.

To give an intuition about why RMax may act as a proxy for CMod, we place ourselves in the continuous framework. Consider the real-valued wavelet transform output $\text{Re} F_1 := F * \text{Re} \bar{\Psi}$, employed in RMax, as the real part of the complex-valued wavelet transform output $F_1 := F * \bar{\Psi}$, used in CMod. At a given location $x \in \mathbb{R}^2$, the corresponding imaginary part may carry a large amount of information, which somehow needs to be retrieved. The key idea is that, if $\Psi$ is sufficiently localized in the Fourier domain, then only the phase of $F_1$ significantly varies in the vicinity of $x$, whereas its magnitude remains nearly constant. Therefore, finding the maximum value of $\text{Re} F_1$ within a local neighborhood around $x$ is nearly equivalent to shifting the phase of $F_1(x)$ towards 0. The resulting value then approximates $|F_1(x)|$. To put it differently, max pooling pushes energy towards lower frequencies, in a similar way as the modulus does for complex-valued transforms (Bruna & Mallat, 2013). This result is hinted in Section 3.1.

Regretfully, things do not work so smoothly in the discrete case. At first glance, this is surprising because Shannon's sampling theorem allows to cast discrete problems into the continuous framework, as done in Section 2.4. However, as explained in Section 3.2, max pooling operates over a discrete grid instead of a continuous window. Consequently, in some situations, the maximum value may fall far away from any zero-phase coefficient. Taking into account this behavior, we adopt a probabilistic point of view, as detailed in Section 3.4. Then, we provide in Section 3.5 an upper bound for the expected gap between CMod and RMax outputs.

## 3.1 Continuous Framework

This section, inspired from Waldspurger (2015, pp. 190–191), provides an intuition about resemblance between RMax and CMod in the continuous framework. As will be highlighted in Section 3.2, adaptation to discrete 2D sequences is not straightforward and will require a probabilistic approach.

We consider an input function $F \in L_{\mathbb{R}}^2(\mathbb{R}^2)$ and a band-pass filter $\Psi \in \mathcal{V}(\nu, \varepsilon)$. Let us also consider

$$ G : (x, h) \mapsto \cos(\langle \nu, h \rangle - H(x)), \qquad (41) $$
Proof. First, $F_X$ and $\Psi_W$ are well defined because $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and $W \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$. By construction, $F_X$ and $\Psi_W \in \mathcal{V}^{(s)}$. Therefore, according to Shannon's sampling theorem [27, Theorem 3.11, p. 81],

$$ (2.27) \qquad F_X := s \sum_{n \in \mathbb{Z}^2} F_X(sn) \Phi_n^{(s)} \quad \text{and} \quad \Psi_W := s \sum_{n \in \mathbb{Z}^2} \Psi_W(sn) \Phi_n^{(s)}. $$

By uniqueness of decompositions in an orthonormal basis, we get (2.25). Moreover, using (2.19) in Lemma 2.4, we get, for any $\xi \in B_\infty(\pi/s)$,

$$ (2.28) \qquad \hat{\Psi}_W(\xi) = s \hat{W}(s\xi). $$

Since $\hat{\Psi}_W(\xi) = 0$ outside $B_\infty(\pi/s)$, (2.28) is true for any $\xi \in \mathbb{R}^2$. Therefore, by hypothesis on $W$,

$$ (2.29) \qquad \operatorname{supp} \hat{\Psi}_W \subset B_\infty(\theta/s, \kappa/(2s)), $$

which yields (2.24).

We now prove (2.26). For $n \in \mathbb{Z}^2$, we compute:

$$ \begin{aligned} (F_X * \bar{\Psi}_W)(msn) &= \iint_{\mathbb{R}^2} F_X(msn - x) \bar{\Psi}_W(x) \, d^2x \\ &= \iint_{\mathbb{R}^2} \sum_{p \in \mathbb{Z}^2} X[p] \Phi_p^{(s)}(msn - x) \bar{\Psi}_W(x) \, d^2x \\ &= \sum_{p \in \mathbb{Z}^2} X[p] \iint_{\mathbb{R}^2} \Phi_p^{(s)}(msn - x) \bar{\Psi}_W(x) \, d^2x. \end{aligned} $$

The sum-integral interchange is possible because $X$ has a finite support. Then:

$$ (2.30) \qquad (F_X * \bar{\Psi}_W)(msn) = \sum_{p \in \mathbb{Z}^2} X[p] \iint_{\mathbb{R}^2} \bar{\Psi}_W(x) \Phi^{(s)}(s(mn - p) - x) \, d^2x $$

$$ (2.31) \qquad = \sum_{p \in \mathbb{Z}^2} X[p] (\bar{\Psi}_W * \Phi^{(s)})(s(mn - p)) $$

Since $\{\Phi_n^{(s)}\}_{n \in \mathbb{Z}^2}$ is an orthonormal basis of $\mathcal{V}^{(s)}$, the definition of $\Psi_W$ in (2.23) implies, for any $p' \in \mathbb{Z}^2$,

$$ (2.32) \qquad \overline{W}[p'] = \langle \Psi_W, \Phi_{-p'}^{(s)} \rangle = (\overline{\Psi}_W * \Phi^{(s)})(sp'), $$

because $\Phi^{(s)}$ is even. Therefore, plugging (2.32) with $p' \leftarrow (mn - p)$ into (2.31) yields

$$ (2.33) \qquad (F_X * \bar{\Psi}_W)(msn) = \sum_{p \in \mathbb{Z}^2} X[p] \overline{W}[mn - p] = (X * \overline{W})[mn], $$

hence the result.

<figure><img src="image_7.png" /></figure>
where $H: \mathbb{R}^2 \to [0, 2\pi]$ denotes the phase of $F * \bar{\Psi}$. Lemma 1 introduced low-frequency functions $F_0$, with slow variations. In a nutshell, since $\text{supp } F_0 \subset B_\infty(\varepsilon/2)$, we can write

$$ \|h\|_2 \ll \lambda_{F_0} \implies F_0(x + h) \approx F_0(x), \quad (42) $$

where we have defined $\lambda_{F_0} := 2\pi/\varepsilon$. Therefore, according to Proposition 4 below, we get the following approximation of $F * \text{Re}\bar{\Psi}$ in a neighborhood around any point $x \in \mathbb{R}^2$:

$$ \|h\|_2 \ll \lambda_{F_0} \implies (F * \text{Re}\bar{\Psi})(x + h) \approx |(F * \bar{\Psi})(x)| G(x, h). \quad (43) $$

**Proposition 4.** For any $h \in \mathbb{R}^2$,

$$ |(F * \text{Re}\bar{\Psi})(x + h) - |(F * \bar{\Psi})(x)| G(x, h)| \le |F_0(x + h) - F_0(x)|. \quad (44) $$

*Proof.* See Appendix A.9. □

On the one hand, we consider a continuous equivalent of the $\mathbb{C}\text{Mod}$ operator $U_m^{\text{mod}}[W]$ as introduced in (5). Such an operator, denoted by $U^{\text{mod}}[\Psi]$, is defined, for any $F \in L_{\mathbb{R}}^2(\mathbb{R}^2)$, by

$$ U^{\text{mod}}[\Psi](F) : x \mapsto |(F * \bar{\Psi})(x)|. \quad (45) $$

On the other hand, we consider the continuous counterpart of $\mathbb{R}\text{Max}$ as introduced in (1). It is defined as the maximum value of $F * \text{Re}\bar{\Psi}$ over a sliding spatial window of size $r > 0$. This is possible because $F$ and $\text{Re}\bar{\Psi}$ both belong to $L_{\mathbb{R}}^2(\mathbb{R}^2)$, and therefore $F * \text{Re}\bar{\Psi}$ is continuous. Such an operator, denoted by $U_r^{\text{max}}[\Psi]$, is defined, for any $F \in L_{\mathbb{R}}^2(\mathbb{R}^2)$, by

$$ U_r^{\text{max}}[\Psi](F) : x \mapsto \max_{\|h\|_\infty \le r} (F * \text{Re}\bar{\Psi})(x + h). \quad (46) $$

For the sake of conciseness, the parameter between square brackets is ignored from now on. If $r \ll \lambda_{F_0}$, then (43) is valid for any $h \in B_\infty(r)$. Then, using (45) and (46), we get

$$ r \ll \lambda_{F_0} \implies U_r^{\text{max}} F(x) \approx U^{\text{mod}} F(x) \max_{\|h\|_\infty \le r} G(x, h). \quad (47) $$

Using the periodicity of $G$, we can show that, if $r \ge \frac{\pi}{\|\nu\|_2}$, then $h \mapsto G(x, h)$ necessarily reaches its maximum value (i.e., 1) on $B_\infty(r)$. We therefore get

$$ \frac{\pi}{\|\nu\|_2} \le r \ll \frac{2\pi}{\varepsilon} \implies U_r^{\text{max}} F(x) \approx U^{\text{mod}} F(x). \quad (48) $$

## 3.2 Adaptation to Discrete 2D Sequences

As in Section 2.4, we consider an input image $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$, a complex, analytic convolution kernel $W \in \mathcal{J}(\theta, \kappa)$, a subsampling factor $m \in \mathbb{N} \setminus \{0\}$ and an integer $q \in \mathbb{N} \setminus \{0\}$, referred to as a *half-size*, such that max pooling operates on a grid of size $(2q+1) \times (2q+1)$. We seek a relationship between

$$ \text{Y}^{\text{max}} := U_{m,q}^{\text{max}}[W](X) \quad \text{and} \quad \text{Y}^{\text{mod}} := U_m^{\text{mod}}[W](X), \quad (49) $$

where $U_{m,q}^{\text{max}}[W](\mathbb{R}\text{Max})$ and $U_m^{\text{mod}}[W](\mathbb{C}\text{Mod})$ have been respectively defined in (1) and (5). As before, in what follows we omit the parameter between square brackets.

We now use the sampling results from Proposition 2. Let $F_X$ and $\Psi_W \in \mathcal{V}^{(s)}$ denote the functions satisfying (24). Recall that the continuous versions of $\mathbb{C}\text{Mod}$ and $\mathbb{R}\text{Max}$ operators have been defined in (45) and (46), respectively. On the one hand, we apply (27) with $m \leftarrow 2m$ to $\text{Y}^{\text{mod}}$. For any $n \in \mathbb{Z}^2$,

$$ U_m^{\text{mod}}X[n] = (F_X * \bar{\Psi}_W)(x_n) \quad (50) $$

$$ = U^{\text{mod}} F_X(x_n), \quad (51) $$
**Proposition 2.5** introduces a latent subspace of $L_{\mathbb{R}}^2(\mathbb{R}^2)$ from which input images are uniformly sampled. This allows us to define, for any $\mathbf{u} \in \mathbb{R}^2$, a translation operator $\mathcal{T}_{\mathbf{u}}$ on discrete sequences, even if $\mathbf{u}$ has non-integer values:

$$ (2.34) \qquad \mathcal{T}_{\mathbf{u}}X[n] := s \mathcal{T}_{su}F_X(sn), $$

where $F_X$ is defined in (2.23). We can indeed show that this definition is independent from the choice of sampling interval $s > 0$. Besides, given $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$, we have

$$ (2.35) \qquad \forall \mathbf{p} \in \mathbb{Z}^2, \quad \mathcal{T}_{\mathbf{p}}X[n] = X[n - \mathbf{p}]; $$

$$ (2.36) \qquad \forall \mathbf{u}, \mathbf{v} \in \mathbb{R}^2, \quad \mathcal{T}_{\mathbf{u}}(\mathcal{T}_{\mathbf{v}}X) = \mathcal{T}_{\mathbf{u}+\mathbf{v}}X, $$

which shows that $\mathcal{T}_{\mathbf{u}}$ corresponds to the intuitive idea of a translation operator. Expressions (2.35) and (2.36) are direct consequence of the following lemma, which bonds the shift operator in the discrete and continuous frameworks.

**Lemma 2.6.** For any $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and any $\mathbf{u} \in \mathbb{R}^2$,

$$ (2.37) \qquad F_{\mathcal{T}_{\mathbf{u}}X} = \mathcal{T}_{su}F_X. $$

*Proof.* Let $\mathbf{u} \in \mathbb{R}^2$. By definition of $F_{\mathcal{T}_{\mathbf{u}}X}$ and $\mathcal{T}_{\mathbf{u}}X$,

$$ (2.38) \qquad F_{\mathcal{T}_{\mathbf{u}}X} = s \sum_{n \in \mathbb{Z}^2} \mathcal{T}_{su}F_X(sn) \Phi_n^{(s)}. $$

On the other hand, $F_X \in \mathcal{V}^{(s)}$ by construction. Therefore, $\mathcal{T}_{su}F_X \in \mathcal{V}^{(s)}$. Then, according to Shannon's sampling theorem [27, Theorem 3.11, p. 81], we get

$$ (2.39) \qquad \mathcal{T}_{su}F_X = s \sum_{n \in \mathbb{Z}^2} \mathcal{T}_{su}F_X(sn) \Phi_n^{(s)}, $$

which concludes the proof.

We now consider the following corollary to Proposition 2.5.

**Corollary 2.7.** For any shift vector $\mathbf{u} \in \mathbb{R}^2$, we have

$$ (2.40) \qquad ((\mathcal{T}_{\mathbf{u}}X * \overline{W}) \downarrow m)[n] = (\mathcal{T}_{su}F_X * \overline{\Psi}_W)(msn). $$

*Proof.* Applying (2.26) in Proposition 2.5 with $X \leftarrow \mathcal{T}_{\mathbf{u}}X$, we get

$$ (2.41) \qquad ((\mathcal{T}_{\mathbf{u}}X * \overline{W}) \downarrow m)[n] = (F_{\mathcal{T}_{\mathbf{u}}X} * \overline{\Psi}_W)(msn), $$

and Lemma 2.6 concludes the proof.
with $\boldsymbol{x}_n := 2msn$. On the other hand, we postulate that

$$U_{m,q}^{\max} X[n] = U_r^{\max} F_X(\boldsymbol{x}_n) \qquad (52)$$

for a certain value of $r \in \mathbb{R}_{>0}$. Then, (48) implies $Y^{\text{mod}} \approx Y^{\text{max}}$. However, as explained hereafter, (52) is not satisfied, due to the discrete nature of the max pooling grid. According to (1) and (4), we have

$$U_{m,q}^{\max} X[n] = \max_{\|\boldsymbol{p}\|_{\infty} \le q} \mathrm{Re}((X * \bar{W}) \downarrow m) [2\boldsymbol{n} + \boldsymbol{p}]. \qquad (53)$$

Therefore, according to (27) in Proposition 2, we get

$$U_{m,q}^{\max} X[n] = \max_{\|\boldsymbol{p}\|_{\infty} \le q} (F_X * \mathrm{Re} \bar{\Psi}_W)(\boldsymbol{x}_n + \boldsymbol{h}_p), \qquad (54)$$

with

$$\boldsymbol{x}_n := 2msn \quad \text{and} \quad \boldsymbol{h}_p := ms\boldsymbol{p}. \qquad (55)$$

By considering $r_q := ms(q + \frac{1}{2})$, we get a variant of (52) in which the maximum is evaluated on a discrete grid of $(2q + 1)^2$ elements, instead of the continuous region $B_{\infty}(r_q)$, as defined in (46) with $r \leftarrow r_q$. As a consequence, (47) is replaced in the discrete framework by

$$q \ll 2\pi/(m\kappa) \quad \Rightarrow \quad U_{m,q}^{\max} X[n] \approx U_m^{\text{mod}} X[n] \max_{\|\boldsymbol{p}\|_{\infty} \le q} G_X(\boldsymbol{x}_n, \boldsymbol{h}_p), \qquad (56)$$

where we have introduced, similarly to (41),

$$G_X : (\boldsymbol{x}, \boldsymbol{h}) \mapsto \cos(\langle \boldsymbol{\nu}, \boldsymbol{h} \rangle - H_X(\boldsymbol{x})), \qquad (57)$$

with

$$\boldsymbol{\nu} := \boldsymbol{\theta}/s \quad \text{and} \quad H_X := \angle(F_X * \bar{\Psi}_W), \qquad (58)$$

where $\angle : \mathbb{C} \to [0, 2\pi]$ denotes the phase operator. Unlike the continuous case, even if the window size $r_q$ is large enough, the existence of $\boldsymbol{p} \in \{-q, \dots, q\}^2$ such that $G_X(\boldsymbol{x}_n, \boldsymbol{h}_p) = 1$ is not guaranteed, as illustrated in fig. 3 with $q = 1$. Instead, we can only seek a probabilistic estimation of the normalized mean squared error between $Y^{\text{max}}$ and $Y^{\text{mod}}$.

Approximation (56) implies

$$q \ll 2\pi/(m\kappa) \quad \Rightarrow \quad \|U_m^{\text{mod}} X - U_{m,q}^{\max} X\|_2 \approx \|\delta_{m,q} X\|_2, \qquad (59)$$

where $\delta_{m,q} X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ is defined such that, for any $\boldsymbol{n} \in \mathbb{Z}^2$,

$$\delta_{m,q} X[n] := U_m^{\text{mod}} X[n] \left( 1 - \max_{\|\boldsymbol{p}\|_{\infty} \le q} G_X(\boldsymbol{x}_n, \boldsymbol{h}_p) \right). \qquad (60)$$

Expression (59) suggests that the difference between the left and right terms can be bounded by a quantity which only depends on the product $m\kappa$ (subsampling factor $\times$ frequency localization) and the grid half-size $q$. In what follows, we establish a bound characterizing this approximation, which will be provided in Proposition 5.

For the sake of conciseness, we introduce the following notations:

$$A_X : (\boldsymbol{x}, \boldsymbol{h}) \mapsto (F_X * \mathrm{Re} \bar{\Psi}_W)(\boldsymbol{x} + \boldsymbol{h}); \qquad (61)$$

$$\tilde{A}_X : (\boldsymbol{x}, \boldsymbol{h}) \mapsto |(F_X * \bar{\Psi}_W)(\boldsymbol{x})| G_X(\boldsymbol{x}, \boldsymbol{h}). \qquad (62)$$

We now consider, for any $\boldsymbol{n} \in \mathbb{Z}^2$, the vectors $\boldsymbol{h}_n^{\max}$ and $\boldsymbol{h}'_n^{\max} \in ms\{-q, \dots, q\}^2$ achieving the maximum value of $A_X(\boldsymbol{x}_n, \boldsymbol{h}_p)$ and $\tilde{A}_X(\boldsymbol{x}_n, \boldsymbol{h}_p)$ over the max pooling grid, respectively. They satisfy

$$A_X^{\max}(\boldsymbol{x}_n) := A_X(\boldsymbol{x}_n, \boldsymbol{h}_n^{\max}) = \max_{\|\boldsymbol{p}\|_{\infty} \le q} A_X(\boldsymbol{x}_n, \boldsymbol{h}_p); \qquad (63)$$

$$\tilde{A}_X^{\max}(\boldsymbol{x}_n) := \tilde{A}_X(\boldsymbol{x}_n, \boldsymbol{h}'_n^{\max}) = \max_{\|\boldsymbol{p}\|_{\infty} \le q} \tilde{A}_X(\boldsymbol{x}_n, \boldsymbol{h}_p). \qquad (64)$$
**2.5. Shift Invariance in the Discrete Framework.** We consider the $\mathbb{C}\text{Mod}$ operator defined in (1.5). For the sake of conciseness, in what follows we will write $U_m^{\text{mod}}$ instead of $U_m^{\text{mod}}[W]$, when no ambiguity is possible. First, we state the following lemma.

**Lemma 2.8.** For any input image $X \in l_1^2(\mathbb{Z}^2)$ with finite support, and any Gabor-like filter $W \in \mathcal{J}(\theta, \kappa)$, we consider the low-frequency function

$$ (2.42) \qquad F_0: \boldsymbol{x} \mapsto (F_X * \bar{\Psi}_W)(\boldsymbol{x}) e^{i\langle\theta/s, \boldsymbol{x}\rangle}, $$

with $F_X$ and $\Psi_W$ satisfying (2.23). If $\kappa \le \pi/m$, then

$$ (2.43) \qquad F_0 \in \mathcal{V}^{(s')}. $$

Moreover, for any $\boldsymbol{h} \in \mathbb{R}^2$,

$$ (2.44) \qquad \sum_{\boldsymbol{n} \in \mathbb{Z}^2} \left| \mathcal{T}_{\boldsymbol{h}} F_0(s' \boldsymbol{n}) - F_0(s' \boldsymbol{n}) \right|^2 = \frac{1}{s'^2} \| \mathcal{T}_{\boldsymbol{h}} F_0 - F_0 \|_{L^2}^2, $$

where we have denoted $s' := 2ms$. Finally,

$$ (2.45) \qquad \| U_m^{\text{mod}} X \|_2 = \frac{1}{s'} \| F_0 \|_{L^2}. $$

*Proof.* Let us write:

$$ (2.46) \qquad \sum_{\boldsymbol{n} \in \mathbb{Z}^2} \left| \mathcal{T}_{\boldsymbol{h}} F_0(s' \boldsymbol{n}) - F_0(s' \boldsymbol{n}) \right|^2 = \sum_{\boldsymbol{n} \in \mathbb{Z}^2} \left| F_1(s' \boldsymbol{n}) \right|^2 = \frac{1}{s'^2} \| X_1 \|_2^2, $$

where we have denoted, for any $\boldsymbol{n} \in \mathbb{Z}^2$,

$$ (2.47) \qquad F_1 := \mathcal{T}_{\boldsymbol{h}} F_0 - F_0 \quad \text{and} \quad X_1[\boldsymbol{n}] := s' F_1(s' \boldsymbol{n}). $$

According to Proposition 2.5 (2.24), $\Psi_W \in \mathcal{V}(\theta/s, \kappa/s)$. Therefore, according to Lemma 2.2,

$$ (2.48) \qquad \operatorname{supp} \widehat{F}_0 \subset B_\infty \left( \frac{\kappa}{2s} \right). $$

Moreover, by hypothesis, $\kappa \le \pi/m$; thus,

$$ (2.49) \qquad B_\infty \left( \frac{\kappa}{2s} \right) \subset B_\infty \left( \frac{\pi}{s'} \right), $$

which yields (2.43), and $F_1 \in \mathcal{V}^{(s')}$. Then, according to Lemma 2.4 (2.20) with $X \leftarrow X_1$, $F \leftarrow F_1$ and $s \leftarrow s'$,

$$ (2.50) \qquad \| X_1 \|_2 = \| F_1 \|_{L^2} = \| \mathcal{T}_{\boldsymbol{h}} F_0 - F_0 \|_{L^2}. $$

Therefore, plugging (2.50) into (2.46) yields (2.44).
Then, according to (50) and (54), we get, for any $\mathbf{n} \in \mathbb{Z}^2$,

$$A_X^{\max}(\mathbf{x}_n) = U_{m,q}^{\max} X[n]; \qquad (65)$$

$$\tilde{A}_X^{\max}(\mathbf{x}_n) = U_m^{\text{mod}} X[n] \max_{\|\mathbf{p}\|_{\infty} \le q} G_X(\mathbf{x}_n, \mathbf{h}_p), \qquad (66)$$

and (56) becomes

$$q \ll 2\pi/(m\kappa) \quad \implies \quad A_X^{\max}(\mathbf{x}_n) \approx \tilde{A}_X^{\max}(\mathbf{x}_n). \qquad (67)$$

**Remark 3.** Expression (43) implies that, if $q \ll 2\pi/(m\kappa)$, then $A_X(\mathbf{x}_n, \mathbf{h}_p) \approx \tilde{A}_X(\mathbf{x}_n, \mathbf{h}_p)$ for all $\mathbf{p} \in \{-q, \dots, q\}^2$. However, this property does not guarantee that $A_X$ and $\tilde{A}_X$ reach their maximum in the same exact location; i.e., that $\mathbf{h}_n^{\max} = \mathbf{h}_n'^{\max}$.

The following lemma provides a bound for approximation (67).

**Lemma 5.** For any $\mathbf{x} \in \mathbb{R}^2$,

$$|A_X^{\max}(\mathbf{x}_n) - \tilde{A}_X^{\max}(\mathbf{x}_n)| \le \max_{\mathbf{h} \in \{\mathbf{h}_n^{\max}, \mathbf{h}_n'^{\max}\}} |F_0(\mathbf{x}_n + \mathbf{h}) - F_0(\mathbf{x}_n)|. \qquad (68)$$

*Proof.* See Appendix A.10. $\square$

Before stating Proposition 5, we consider the following hypothesis:

**Hypothesis 1.** There exists $\mathbf{h}_0 \in \mathbb{R}^2$ with $\|\mathbf{h}_0\|_2 = \sqrt{2qms}$, such that

$$\sum_{\mathbf{n} \in \mathbb{Z}^2} \max_{\mathbf{h} \in \{\mathbf{h}_n^{\max}, \mathbf{h}_n'^{\max}\}} |F_0(\mathbf{x}_n + \mathbf{h}) - F_0(\mathbf{x}_n)|^2 \le \sum_{\mathbf{n} \in \mathbb{Z}^2} |F_0(\mathbf{x}_n + \mathbf{h}_0) - F_0(\mathbf{x}_n)|^2. \qquad (69)$$

The underlying idea is explained as follows. The absolute difference between $F_0(\mathbf{x}_n + \mathbf{h})$ and $F_0(\mathbf{x}_n)$ is more likely to increase with the norm of $\mathbf{h}$. For any given $\mathbf{n} \in \mathbb{Z}^2$, we have, by construction, $\|\mathbf{h}_n^{\max}\|_2 \le \sqrt{2qms}$ and $\|\mathbf{h}_n'^{\max}\|_2 \le \sqrt{2qms}$. Therefore, we can expect to observe

$$\max_{\mathbf{h} \in \{\mathbf{h}_n^{\max}, \mathbf{h}_n'^{\max}\}} |F_0(\mathbf{x}_n + \mathbf{h}) - F_0(\mathbf{x}_n)|^2 \le |F_0(\mathbf{x}_n + \mathbf{h}_0) - F_0(\mathbf{x}_n)|^2. \qquad (70)$$

While this might occasionally not be true, Hypothesis 1 postulates that, when summing over all the data-points, the inequality holds.

We now formally state the result characterizing approximation (59).

**Proposition 5.** We assume that condition (38) is satisfied: $\kappa \le \pi/m$. Then, under Hypothesis 1,

$$\|U_m^{\text{mod}}X - U_{m,q}^{\max}X\|_2 \le \|\delta_{m,q}X\|_2 + \beta_q(m\kappa) \|U_m^{\text{mod}}X\|_2, \qquad (71)$$

where $\beta_q : \mathbb{R}_+ \to \mathbb{R}_+$ is defined by

$$\beta_q : \kappa' \mapsto q\kappa'. \qquad (72)$$

*Proof.* See Appendix A.11. $\square$

We now seek a probabilistic estimation of $\|\delta_{m,q}X\|_2$. For this purpose, we first reformulate the problem using the unit circle $\mathbb{S}^1 \subset \mathbb{C}$, before introducing a probabilistic framework in Section 3.4.
Besides, according again to Lemma 2.4,

$$ (2.51) \qquad \|F_0\|_{L^2}^2 = \|X_0\|_2^2, $$

where we have defined, for any $\mathbf{n} \in \mathbb{Z}^2$,

$$ (2.52) \qquad X_0[\mathbf{n}] := s'F_0(s'\mathbf{n}). $$

Then,

$$ \|X_0\|_2^2 = s'^2 \sum_{\mathbf{n} \in \mathbb{Z}^2} |(F_X * \bar{\Psi}_W)(s'\mathbf{n})|^2 \qquad \text{(acc. to (2.42))} $$

$$ = s'^2 \sum_{\mathbf{n} \in \mathbb{Z}^2} |(X * \bar{W}) \downarrow (2m)[\mathbf{n}]|^2 \quad \text{(acc. to Proposition 2.5 with } m \leftarrow 2m \text{)} $$

$$ = s'^2 \|U_m^{\text{mod}} X\|_2^2. \qquad \text{(acc. to (1.5))} $$

Finally, plugging this result into (2.51) yields (2.45) and concludes the proof. ■

We are now ready to state the main result about shift invariance of CMod outputs.

**Theorem 2.9 (Shift invariance of CMod).** Let $W \in \mathcal{J}(\theta, \kappa)$ denote a discrete Gabor-like filter and $m \in \mathbb{N} \setminus \{0\}$ denote a subsampling factor. Then, under the following condition:

$$ (2.53) \qquad \kappa \le \pi/m, $$

we have, for any input image $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ with finite support and any translation vector $\mathbf{u} \in \mathbb{R}^2$,

$$ (2.54) \qquad \|U_m^{\text{mod}}(\mathcal{T}_u X) - U_m^{\text{mod}} X\|_2 \le \alpha(\kappa \mathbf{u}) \|U_m^{\text{mod}} X\|_2, $$

where $\alpha$ has been defined in (2.14).

**Proof.** As in Lemma 2.8, we consider the low-frequency function $F_0$ satisfying (2.42), and denote $s' := 2ms$. We can write

$$ (2.55) \qquad |F_X * \bar{\Psi}_W| = |F_0| \quad \text{and} \quad |\mathcal{T}_{su} F_X * \bar{\Psi}_W| = |\mathcal{T}_{su} F_0|. $$

Recall that $U_m^{\text{mod}} X = |(X * \bar{W}) \downarrow (2m)|$, such as defined in (1.5). According to Proposition 2.5 (2.26) and Corollary 2.7 (2.40) with $m \leftarrow 2m$, we therefore get

$$ (2.56) \qquad U_m^{\text{mod}} X[n] = |F_0(s'n)|; $$

$$ (2.57) \qquad U_m^{\text{mod}}(\mathcal{T}_u X)[n] = |(\mathcal{T}_{su} F_0)(s'n)|. $$

Then, using (2.56), (2.57) and the reverse triangle inequality,

$$ \begin{aligned} \|U_m^{\text{mod}}(\mathcal{T}_u X) - U_m^{\text{mod}} X\|_2^2 &= \sum_{\mathbf{n} \in \mathbb{Z}^2} \left| |(\mathcal{T}_{su} F_0)(s'n)| - |F_0(s'n)| \right|^2 \\ &\le \sum_{\mathbf{n} \in \mathbb{Z}^2} \left| (\mathcal{T}_{su} F_0)(s'n) - F_0(s'n) \right|^2. \end{aligned} $$
<figure><img src="image_4.png" /><figcaption>Figure 3. Search for the maximum value of $h \mapsto G_X(x, h)$ over a discrete grid of size $3 \times 3$, i.e., $q = 1$. This figure displays 3 examples with different frequencies $\nu := \theta/s$ and phases $H_X(x)$. Hopefully the result will be close to the true maximum (left), but there are some pathological cases in which all points in the grid fall into pits (middle and right).</figcaption></figure>

### 3.3 Notations on the Unit Circle

In what follows, for any $z \in \mathbb{C} \setminus \{0\}$, we denote by $\angle z \in [0, 2\pi[$ the argument of $z$. For any $z, z' \in \mathbb{S}^1$, the angle between $z$ and $z'$ is given by $\angle(z^*z')$. We then denote by $[z, z']_{\mathbb{S}^1} \subset \mathbb{S}^1$ the arc on the unit circle going from $z$ to $z'$ counterclockwise:

$$[z, z']_{\mathbb{S}^1} := \{z'' \in \mathbb{S}^1 \mid \angle(z^*z'') \le \angle(z^*z')\}. \qquad (73)$$

We remind readers that $x_n$ and $h_p \in \mathbb{R}^2$ have been defined in (55). By using the relation $\cos \alpha = \text{Re}(e^{i\alpha})$, (57) becomes, for any $n \in \mathbb{Z}^2$ and any $p \in \{-q, \dots, q\}^2$,

$$G_X(x_n, h_p) = \text{Re}(Z_X^*(x_n) Z_p(m\theta)), \qquad (74)$$

where we have defined the following functions with outputs on the unit circle:

$$Z_X : x \mapsto e^{i H_X(x)} \quad \text{and} \quad Z_p : \omega \mapsto e^{i \langle \omega, p \rangle}, \qquad (75)$$

where $H_X$ denotes the phase of $F_X * \bar{\Psi}_W$ as introduced in (58). On the one hand, $Z_X(x_n)$ is the phase (represented on the unit circle $\mathbb{S}^1$) of the complex wavelet transform $F_X * \bar{\Psi}_W$ at location $x_n$. On the other hand, $Z_p(m\theta)$ approximates the phase shift between any two evaluations of $F_X * \bar{\Psi}_W$ at locations $x, x'$ such that $x' - x = h_p$. This however is only true if we assume that $\Psi_W$ exhibits slow amplitude variations. Then, $G_X(x_n, h_p)$ approximates the cosine of the phase of $F_X * \bar{\Psi}_W$ at location $x_n + h_p$.

According to (56), $\max_{\|p\|_{\infty} \le q} G_X(x_n, h_p)$ approximates the ratio between $\mathbb{R}$Max and $\mathbb{C}$Mod outputs at discrete location $n \in \mathbb{Z}^2$. The intuition behind this is that max pooling seeks a point in a discrete grid around $x_n$ where the phase of $F_X * \bar{\Psi}_W$ is the closest to 1, thereby maximizing the amount of energy on the real part of the signal. Assuming slow amplitude variations of $\Psi_W$, the result therefore approximates the modulus of the complex coefficients.

To get an estimation of $\delta_{m,q}X[n]$ (60), we will exploit the following property. If the phases $Z_p(m\theta)$ for $p \in \{-q, \dots, q\}^2$ are well distributed on the unit circle, then the values of $G_X(x_n, h_p)$ are evenly spread out on $[-1, 1]$. Therefore, its maximum value is more likely to be close to 1, and (60) becomes

$$\delta_{m,q} X[n] \ll U_m^{\text{mod}} X[n] \quad \forall n \in \mathbb{Z}^2. \qquad (76)$$

Let $n_q := (2q + 1)^2$ denote the number of evaluation points for the max pooling operator. For any $\omega \in \mathbb{R}^2$, we consider a sequence of values on $\mathbb{S}^1$, denoted by $(Z_i^{(q)}(\omega))_{i \in \{0, \dots, n_q - 1\}}$, obtained by sorting $\{Z_p(\omega)\}_{p \in \{-q, \dots, q\}^2}$ (75) in ascending order of their arguments:

$$0 = H_0^{(q)}(\omega) \le \dots \le H_{n_q-1}^{(q)}(\omega) < 2\pi, \qquad (77)$$
Since condition (2.53) is satisfied, we can use Lemma 2.8 (2.44) with $h \leftarrow su$:

$$ (2.58) \qquad \|U_m^{\text{mod}}(\mathcal{T}_u X) - U_m^{\text{mod}}X\|_2^2 \le \frac{1}{s'^2} \|\mathcal{T}_{su}F_0 - F_0\|_{L^2}^2 $$

Now, according to Proposition 2.3 with $\varepsilon \leftarrow \kappa/s$ and $h \leftarrow su$, we then get the following bound:

$$ (2.59) \qquad \|U_m^{\text{mod}}(\mathcal{T}_u X) - U_m^{\text{mod}}X\|_2^2 \le \frac{\alpha(\kappa u)^2}{s'^2} \|F_0\|_{L^2}^2. $$

Finally, using Lemma 2.8 (2.45) yields (2.54), which completes the proof.

Interestingly, the reference value used in Theorem 2.9, i.e., $\|U_m^{\text{mod}}X\|_2$, is fully shift-invariant, as stated in the following proposition.

**Proposition 2.10.** Let $W \in \mathcal{J}(\theta, \kappa)$ and $m \in \mathbb{N} \setminus \{0\}$. Under condition (2.53), we have, for any $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and any $u \in \mathbb{R}^2$,

$$ (2.60) \qquad \|U_m^{\text{mod}}(\mathcal{T}_u X)\|_2 = \|U_m^{\text{mod}}X\|_2. $$

*Proof.* Let $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ and $s > 0$. We consider $F_0 \in L_{\mathbb{C}}^2(\mathbb{R}^2)$ as the “low-frequency” function satisfying (2.42). Again, we introduce $s' := 2ms$ and $X_0 \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ satisfying (2.52). Moreover, for any $Y \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$, we denote by $F_Y^{(s')}$ the Shannon interpolation of Y parameterized by $s'$, analogously to (2.23):

$$ (2.61) \qquad F_Y^{(s')} := \sum_{n \in \mathbb{Z}^2} Y[n] \Phi_n^{(s')}. $$

On the one hand, Lemma 2.8 provides (2.45). On the other hand, we seek a similar result with $X \leftarrow \mathcal{T}_u X$. For this purpose, (2.57) can be rewritten

$$ (2.62) \qquad U_m^{\text{mod}}(\mathcal{T}_u X)[n] = |\mathcal{T}_{s'u'} F_0(s'n)|, $$

with $u' := u/(2m)$. Besides, according to Lemma 2.8 (2.43), $F_0 \in \mathcal{V}^{(s')}$. Therefore, Shannon's sampling theorem [27, Theorem 3.11, p. 81] with $s \leftarrow s'$ yields

$$ \begin{aligned} F_0 &= s' \sum_{n \in \mathbb{Z}^2} F_0(s'n) \Phi_n^{(s')} \\ &= \sum_{n \in \mathbb{Z}^2} X_0[n] \Phi_n^{(s')} = F_{X_0}^{(s')}, \end{aligned} $$

where we have used the notations introduced in (2.52) and (2.61). Then, using Lemma 2.6 with $X \leftarrow X_0$, $u \leftarrow u'$ and $s \leftarrow s'$, we get

$$ (2.63) \qquad F_{\mathcal{T}_{u'}X_0}^{(s')} = \mathcal{T}_{s'u'}F_{X_0}^{(s')} = \mathcal{T}_{s'u'}F_0. $$

Besides, (2.25) (from Proposition 2.5) with $X \leftarrow \mathcal{T}_{u'}X_0$ and $s \leftarrow s'$ becomes

$$ (2.64) \qquad \mathcal{T}_{u'}X_0[n] = s' F_{\mathcal{T}_{u'}X_0}^{(s')}(s'n), $$
where $H_i^{(q)}(\omega)$ denotes the phase of $Z_i^{(q)}(\omega)$. In addition, we close the loop with $H_{n_q}^{(q)}(\omega) := 2\pi$ and $Z_{n_q}^{(q)}(\omega) := 1$. Then, we split $\mathbb{S}^1$ into $n_q$ arcs delimited by $Z_i^{(q)}(\omega)$:

$$ \mathfrak{A}_i^{(q)}(\omega) := \begin{cases} [Z_i^{(q)}(\omega), Z_{i+1}^{(q)}(\omega)]_{\mathbb{S}^1} & \text{if } H_{i+1}^{(q)}(\omega) - H_i^{(q)}(\omega) < 2\pi; \\ \mathbb{S}^1 & \text{otherwise.} \end{cases} \qquad (78) $$

Finally, for any $i \in \{0, \dots, n_q - 1\}$, we denote by

$$ \delta H_i^{(q)} : \omega \mapsto H_{i+1}^{(q)}(\omega) - H_i^{(q)}(\omega) \qquad (79) $$

the function computing the angular measure of arc $\mathfrak{A}_i^{(q)}(\omega)$, for any $\omega \in \mathbb{R}^2$.

## 3.4 Probabilistic Framework

From now on, input X is considered as a discrete 2D stochastic process. In order to “randomize” $F_X$ introduced in (24), we define a continuous stochastic process from X, denoted by $F_X$, such that

$$ \forall \mathbf{x} \in \mathbb{R}^2, F_X(\mathbf{x}) := \sum_{n \in \mathbb{Z}^2} X[n] \Phi_n^{(s)}(\mathbf{x}). \qquad (80) $$

Next, we consider the following stochastic processes, which are parameterized by X:

$$ M_X := |F_X * \bar{\Psi}_W|; \quad H_X := \angle(F_X * \bar{\Psi}_W); \quad Z_X := e^{iH_X}, \qquad (81) $$

and, for any $\mathbf{p} \in \{-q, \dots, q\}^2$,

$$ G_{X, \mathbf{p}} := \mathrm{Re}(Z_X^* Z_{\mathbf{p}}(m\theta)); \quad G_X^{\mathrm{max}} := \max_{\|\mathbf{p}\|_{\infty} \le q} G_{X, \mathbf{p}}, \qquad (82) $$

where the deterministic function $Z_{\mathbf{p}}$ has been defined in (75).

**Remark 4.** By continuous extension, $H_X(\mathbf{x})$ and $Z_X(\mathbf{x})$ are uniquely defined at $\mathbf{x}$ such that $M_X(\mathbf{x}) = 0$.

For any $\mathbf{x} \in \mathbb{R}^2$, $F_X(\mathbf{x})$ (24) and $H_X(\mathbf{x})$ (58) are respectively drawn from $F_X(\mathbf{x})$ and $H_X(\mathbf{x})$. Then, $Z_X(\mathbf{x})$ (75) is a realization of $Z_X(\mathbf{x})$. Consequently, according to (74), $G_X(\mathbf{x}, \mathbf{h}_p)$ is a realization of $G_{X, \mathbf{p}}(\mathbf{x})$. Furthermore, according to the definition of $\mathrm{CMod}$ in (5) and $\mathbf{x}_n$ in (55), Proposition 2 with $m \leftarrow 2m$ implies that

$$ M_X(\mathbf{x}_n) = U_m^{\mathrm{mod}} X[n]. \qquad (83) $$

We remind that $\theta \in [-\pi, \pi]^2$ and $\kappa \in ]0, 2\pi]$ respectively denote the center and size of the Fourier support of the complex kernel $W \in \mathcal{J}(\theta, \kappa)$. To compute the expected discrepancy between $Y^{\mathrm{max}}$ and $Y^{\mathrm{mod}}$ such as introduced in (49), we assume that

$$ \|\theta\|_2 \gg 2\pi/N; \qquad (84) $$

$$ \|\theta\|_2 \gg \kappa, \qquad (85) $$

where $N \in \mathbb{N} \setminus \{0\}$ denotes the support size of input images. These assumptions exclude low-frequency filters from the scope of our study. We then state the following hypotheses, for which a justification is provided in Appendix B.

**Hypothesis 2.** For any $\mathbf{x} \in \mathbb{R}^2$, $Z_X(\mathbf{x})$ is uniformly distributed on $\mathbb{S}^1$.

**Hypothesis 3.** For any $n \in \mathbb{N} \setminus \{0\}$ and $\mathbf{x}, \mathbf{y}_0, \dots, \mathbf{y}_{n-1} \in \mathbb{R}^2$, the random variables $M_X(\mathbf{y}_i)$ for $i \in \{0, \dots, n-1\}$ are jointly independent of $Z_X(\mathbf{x})$.
and inserting (2.63) into (2.64) yields

$$ (2.65) \qquad \mathcal{T}_{u'} X_0[n] = s' \mathcal{T}_{s'u'} F_0(s'n). $$

Therefore, (2.62) and (2.65) imply

$$ (2.66) \qquad \|U_m^{\text{mod}}(\mathcal{T}_u X)\|_2 = \frac{1}{s'} \|\mathcal{T}_{u'} X_0\|_2. $$

Moreover, since $F_0 \in \mathcal{V}^{(s')}$, and according to (2.65), we can use Lemma 2.4 with $s \leftarrow s'$, $F \leftarrow \mathcal{T}_{s'u'} F_0$ and $X \leftarrow \mathcal{T}_{u'} X_0$. We get

$$ (2.67) \qquad \|\mathcal{T}_{u'} X_0\|_2 = \|\mathcal{T}_{s'u'} F_0\|_{L^2} = \|F_0\|_{L^2}, $$

and plugging (2.67) into (2.66) yields

$$ (2.68) \qquad \|U_m^{\text{mod}}(\mathcal{T}_u X)\|_2 = \frac{1}{s'} \|F_0\|_{L^2}. $$

Finally, considering Lemma 2.8 (2.45) concludes the proof.

**3. From CMod to RMax.** CMod operators are found in ScatterNets and complex-valued convolutional networks [47]. However, they are absent from conventional, freely-trained CNN architectures. Therefore, Theorem 2.9 cannot be applied as is. Instead, the first convolution layer contains real-valued kernels, and is generally followed by ReLU and max pooling. As shown in section 5, this process can be described with RMax operators, such as defined in (1.1).

As explained in subsection 1.1, an important number of trained convolution kernels exhibit oscillating patterns with well-defined frequencies and orientations. To elaborate, let $V \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ denote such a trained kernel, and consider $W \in l_{\mathbb{C}}^2(\mathbb{Z}^2)$ as the complex-valued companion of $V$ satisfying (2.6). Then, $W$ has its energy concentrated in a small region of the Fourier domain. We thus emit the hypotheses that $W \in \mathcal{J}(\theta, \kappa)$ (2.5) for a certain value of $\theta \in [-\pi, \pi]^2$ and $\kappa \in ]0, 2\pi]$. For the sake of conciseness, from now on we write $U_{m,q}^{\text{max}}$ instead of $U_{m,q}^{\text{max}}[W]$, when no ambiguity is possible. In what follows, we establish conditions on $W$ under which CMod (1.5) and RMax (1.1) operators produce comparable outputs. The final goal, achieved in section 4, is to provide a shift invariance bound for RMax.

To give an intuition about why RMax may act as a proxy for CMod, we place ourselves in the continuous framework. Consider the real-valued wavelet transform output $\text{Re} F_1 := F * \text{Re} \bar{\Psi}$, employed in RMax, as the real part of the complex-valued wavelet transform output $F_1 := F * \bar{\Psi}$, used in CMod. At a given location $x \in \mathbb{R}^2$, the corresponding imaginary part may carry a large amount of information, which somehow needs to be retrieved. The key idea is that, if $\Psi$ is sufficiently localized in the Fourier domain, then only the phase of $F_1$ significantly varies in the vicinity of $x$, whereas its magnitude remains nearly constant. Therefore, finding the maximum value of $\text{Re} F_1$ within a local neighborhood around $x$ is nearly equivalent to shifting the phase of $F_1(x)$ towards 0. The resulting value then approximates $|F_1(x)|$. To put it differently, max pooling pushes energy towards lower frequencies, in a similar way as the modulus does for complex-valued transforms [9]. This result is hinted in subsection 3.1.
### 3.5 Expected Quadratic Error between $\mathbb{R}Max$ and $\mathbb{C}Mod$

In this section, we propose to estimate the expected value of the stochastic quadratic error $\tilde{\mathbb{P}}_X^2$, defined such that

$$ \tilde{\mathbb{P}}_X := \|U_m^{\text{mod}}X - U_{m,q}^{\text{max}}X\|_2 / \|U_m^{\text{mod}}X\|_2. \quad (86) $$

According to (49), this is an estimation of the relative error between $Y^{\text{mod}}$ and $Y^{\text{max}}$.

First, let us reformulate $\delta_{m,q}X$, introduced in (60), using the probabilistic framework. According to (74) and (82), we have, for any $n \in \mathbb{Z}^2$,

$$ \delta_{m,q}X[n] := U_m^{\text{mod}}X[n] (1 - G_X^{\text{max}}(x_n)). \quad (87) $$

We now consider the stochastic process

$$ Q_X := 1 - G_X^{\text{max}}, \quad (88) $$

and the random variable

$$ \tilde{Q}_X := \|\delta_{m,q}X\|_2 / \|U_m^{\text{mod}}X\|_2. \quad (89) $$

The next steps are as follows: (1) at the pixel level, show that $\mathbb{E}[Q_X(x)^2]$ depends on the subsampling factor $m$ and the filter frequency $\theta$, and remains close to zero with some exceptions; (2) at the image level, show that the expected value of $\tilde{Q}_X^2$ is equal to the latter quantity; (3) use Proposition 5, which implies that $\tilde{\mathbb{P}}_X \approx \tilde{Q}_X$, to deduce an upper bound on the expected value of $\tilde{\mathbb{P}}_X^2$.

The first point, which is established in Proposition 6 below, is a key result of this paper. It will be used to prove Theorem 2, which corresponds to the two remaining points.

**Proposition 6.** *Assuming Hypothesis 2, the expected value of $Q_X(x)^2$ is independent from the choice of $x \in \mathbb{R}^2$, and*

$$ \mathbb{E}[Q_X(x)^2] = \gamma_q(m\theta)^2, \quad (90) $$

where we have defined

$$ \gamma_q : \omega \mapsto \sqrt{\frac{3}{2} + \frac{1}{4\pi} \sum_{i=0}^{n_q-1} \left( \sin \delta H_i^{(q)}(\omega) - 8 \sin \frac{\delta H_i^{(q)}(\omega)}{2} \right)}, \quad (91) $$

with $\delta H_i^{(q)}(\omega) \in [0, 2\pi]$ (79) being the length of arc $\mathfrak{A}_i^{(q)}(\omega)$.

*Proof.* For the sake of readability, in this proof we omit the argument of functions $Z_p$ (75), $Z_i^{(q)}$, $H_i^{(q)}$ (77), $\mathfrak{A}_i^{(q)}$ (78), and $\delta H_i^{(q)}$ (79); we assume they are evaluated at $\omega \leftarrow m\theta$. We consider the “Lebesgue” Borel $\sigma$-algebra on $\mathbb{S}^1$ generated by $\{[z, z']_{\mathbb{S}^1} \mid z, z' \in \mathbb{S}^1\} \cup \{\mathbb{S}^1\}$, on which we have defined the angular measure $\vartheta$ such that $\vartheta(\mathbb{S}^1) := 2\pi$, and

$$ \forall z, z' \in \mathbb{S}^1, \vartheta([z, z']_{\mathbb{S}^1}) := \angle(z^*z'). \quad (92) $$

For any $p \in \mathbb{N} \setminus \{0\}$, we compute the $p$-th moment of $G_X^{\text{max}}(x)$ defined in (82). By considering

$$ g_{\text{max}} : \mathbb{S}^1 \to [-1, 1] \\ z \mapsto \max_{\|p\|_{\infty} \le q} \text{Re}(z^* Z_p), \quad (93) $$

we get $G_X^{\text{max}}(x) = g_{\text{max}}(Z_X(x))$. A visual representation of $g_{\text{max}}$ is provided in fig. 4, for two different values of $\theta$. According to Hypothesis 2, $Z_X(x)$ follows a uniform distribution on $\mathbb{S}^1$. Therefore,

$$ \mathbb{E}[G_X^{\text{max}}(x)^p] = \frac{1}{2\pi} \int_{\mathbb{S}^1} g_{\text{max}}(z)^p \, d\vartheta(z), \quad (94) $$
Regretfully, things do not work so smoothly in the discrete case. At first glance, this is surprising because Shannon's sampling theorem allows to cast discrete problems into the continuous framework, as done in subsection 2.4. However, as explained in subsection 3.2, max pooling operates over a discrete grid instead of a continuous window. Consequently, in some situations, the maximum value may fall far away from any zero-phase coefficient. Taking into account this behavior, we adopt a probabilistic point of view, as detailed in subsection 3.4. Then, we provide in subsection 3.5 an upper bound for the expected gap between CMod and RMax outputs.

**3.1. Continuous Framework.** This section, inspired from Waldspurger [50, pp. 190–191], provides an intuition about resemblance between RMax and CMod in the continuous framework. As will be highlighted in subsection 3.2, adaptation to discrete 2D sequences is not straightforward and will require a probabilistic approach.

We consider an input function $F \in L_{\mathbb{R}}^2(\mathbb{R}^2)$ and a band-pass filter $\Psi \in V(\nu, \varepsilon)$. Let us also consider

$$ (3.1) \qquad G: (x, h) \mapsto \cos(\langle \nu, h \rangle - H(x)), $$

where $H: \mathbb{R}^2 \to [0, 2\pi]$ denotes the phase of $F * \bar{\Psi}$. Lemma 2.2 introduced low-frequency functions $F_0$, with slow variations. In a nutshell, since supp $F_0 \subset B_\infty(\varepsilon/2)$, we can write

$$ (3.2) \qquad \|h\|_2 \ll \lambda_{F_0} \implies F_0(x+h) \approx F_0(x), $$

where we have defined $\lambda_{F_0} := 2\pi/\varepsilon$. Therefore, according to Proposition 3.1 below, we get the following approximation of $F * \mathrm{Re}\bar{\Psi}$ in a neighborhood around any point $x \in \mathbb{R}^2$:

$$ (3.3) \qquad \|h\|_2 \ll \lambda_{F_0} \implies (F * \mathrm{Re}\bar{\Psi})(x+h) \approx |(F * \bar{\Psi})(x)| G(x, h). $$

**Proposition 3.1.** For any $h \in \mathbb{R}^2$,

$$ (3.4) \qquad |(F * \mathrm{Re}\bar{\Psi})(x+h) - |(F * \bar{\Psi})(x)| G(x, h)| \le |F_0(x+h) - F_0(x)|. $$

*Proof.* Let us write:

$$ \begin{align*} (F * \mathrm{Re}\bar{\Psi})(x+h) - |(F * \bar{\Psi})(x)| G(x, h) & \\ &= \mathrm{Re}((F * \bar{\Psi})(x+h)) - |(F * \bar{\Psi})(x)| \mathrm{Re}(e^{-i\langle\nu, h\rangle} e^{H(x)}) \\ &= \mathrm{Re}((F * \bar{\Psi})(x+h)) - \mathrm{Re}(|(F * \bar{\Psi})(x)| e^{H(x)} e^{-i\langle\nu, h\rangle}) \\ &= \mathrm{Re}((F * \bar{\Psi})(x+h)) - \mathrm{Re}((F * \bar{\Psi})(x) e^{-i\langle\nu, h\rangle}) \\ &= \mathrm{Re}((F * \bar{\Psi})(x+h) - (F * \bar{\Psi})(x) e^{-i\langle\nu, h\rangle}). \end{align*} $$

Therefore,

$$ \begin{align*} |(F * \mathrm{Re}\bar{\Psi})(x+h) - |(F * \bar{\Psi})(x)| G(x, h)| & \le |(F * \bar{\Psi})(x+h) - (F * \bar{\Psi})(x) e^{-i\langle\nu, h\rangle}| \\ &= |F_0(x+h) e^{-i\langle\nu, x+h\rangle} - F_0(x) e^{-i\langle\nu, x+h\rangle}|, \end{align*} $$

which yields (3.4) and concludes the proof.

<figure><img src="image_8.png" /></figure>
<figure><img src="image_5.png" /><figcaption>Figure 4. Top: 2D representation of $h \mapsto G_X(x_n, h)$ (57), for two different values of $\theta \in \mathbb{R}^2$, $q=1$ and arbitrary values of $m \in \mathbb{N} \setminus \{0\}$ and $s \in \mathbb{R} \setminus \{0\}$. Assuming the plots are centered around $h=0$, each point materializes a location $h_p$ in the max pooling grid, for $p \in \{-q, \dots, q\}^2$. The desirable situation occurs when one of these locations falls near a ridge (bright areas), in which case the outputs produced by $\mathbb{R}Max$ and $\mathbb{C}Mod$ are similar—see (56). Each number $i \in \{0, \dots, 8\}$ represents the rank of $Z_p \in \mathbb{S}^1$ (75), when these values are sorted by ascending order of their arguments (77). If location $h_p$ gets ranked $i$, then we have $Z_p = Z_i^{(q)}$. Bottom: polar representations of $g_{\max}: \mathbb{S}^1 \to [-1, 1]$ (93), corresponding to the same settings. The closer the curve is from the outer ring, the more likely some points $h_p$ will fall near a ridge of $G_X$. (a) Case where the values $Z_p$ are roughly evenly distributed on $\mathbb{S}^1$. (b) Case where these values are concentrated in a small portion of the unit circle. The most extreme cases occurs when $Z_p = 1$ for any $p$. Figure 3 (middle and right) depicts two such situations.</figcaption></figure>

which proves that $\mathbb{E}[G_X^{\max}(x)^p]$ does not depend on $x$. Let us split the unit circle $\mathbb{S}^1$ into the arcs $\mathfrak{A}_0^{(q)}, \dots, \mathfrak{A}_{n_q-1}^{(q)}$ such as introduced in (78):

$$ \mathbb{E}[G_X^{\max}(x)^p] = \frac{1}{2\pi} \sum_{i=0}^{n_q-1} \int_{\mathfrak{A}_i^{(q)}} g_{\max}(z)^p \, d\vartheta(z). \qquad (95) $$

Let $i \in \{0, \dots, n_q - 1\}$. We show that

$$ \forall z \in \mathfrak{A}_i^{(q)}, g_{\max}(z) = \max \left( \mathrm{Re}(z^* Z_i^{(q)}), \mathrm{Re}(z^* Z_{i+1}^{(q)}) \right). \qquad (96) $$

Let $z \in \mathfrak{A}_i^{(q)}$ and $i' \notin \{i, i+1\}$. We prove that

$$ \mathrm{Re}(z^* Z_{i'}^{(q)}) \le \mathrm{Re}(z^* Z_i^{(q)}) \quad \text{or} \quad \mathrm{Re}(z^* Z_{i'}^{(q)}) \le \mathrm{Re}(z^* Z_{i+1}^{(q)}). \qquad (97) $$

On the one hand, we assume that $\angle(z^* Z_{i'}^{(q)}) \le \pi$. By design of $(Z_i^{(q)})_{i \in \{0, \dots, n_q-1\}}$, we have

$$ Z_{i+1}^{(q)} \in [z, Z_{i'}^{(q)}]_{\mathbb{S}^1}. \qquad (98) $$
On the one hand, we consider a continuous equivalent of the CMod operator $U_m^{\text{mod}}[W]$ as introduced in (1.5). Such an operator, denoted by $U^{\text{mod}}[\Psi]$, is defined, for any $F \in L_{\mathbb{R}}^2(\mathbb{R}^2)$, by

$$ (3.5) \qquad U^{\text{mod}}[\Psi](F) : \boldsymbol{x} \mapsto |(F * \bar{\Psi})(\boldsymbol{x})|. $$

On the other hand, we consider the continuous counterpart of $\mathbb{R}Max$ as introduced in (1.1). It is defined as the maximum value of $F * \text{Re}\bar{\Psi}$ over a sliding spatial window of size $r > 0$. This is possible because $F$ and $\text{Re}\bar{\Psi}$ both belong to $L_{\mathbb{R}}^2(\mathbb{R}^2)$, and therefore $F * \text{Re}\bar{\Psi}$ is continuous. Such an operator, denoted by $U_r^{\text{max}}[\Psi]$, is defined, for any $F \in L_{\mathbb{R}}^2(\mathbb{R}^2)$, by

$$ (3.6) \qquad U_r^{\text{max}}[\Psi](F) : \boldsymbol{x} \mapsto \max_{\|\boldsymbol{h}\|_{\infty} \le r} (F * \text{Re}\bar{\Psi})(\boldsymbol{x} + \boldsymbol{h}). $$

For the sake of conciseness, the parameter between square brackets is ignored from now on. If $r \ll \lambda_{F_0}$, then (3.3) is valid for any $\boldsymbol{h} \in B_{\infty}(r)$. Then, using (3.5) and (3.6), we get

$$ (3.7) \qquad r \ll \lambda_{F_0} \implies U_r^{\text{max}} F(\boldsymbol{x}) \approx U^{\text{mod}} F(\boldsymbol{x}) \max_{\|\boldsymbol{h}\|_{\infty} \le r} G(\boldsymbol{x}, \boldsymbol{h}). $$

Using the periodicity of $G$, we can show that, if $r \ge \frac{\pi}{\|\boldsymbol{v}\|_2}$, then $\boldsymbol{h} \mapsto G(\boldsymbol{x}, \boldsymbol{h})$ necessarily reaches its maximum value (i.e., 1) on $B_{\infty}(r)$. We therefore get

$$ (3.8) \qquad \frac{\pi}{\|\boldsymbol{v}\|_2} \le r \ll \frac{2\pi}{\varepsilon} \implies U_r^{\text{max}} F(\boldsymbol{x}) \approx U^{\text{mod}} F(\boldsymbol{x}). $$

### 3.2. Adaptation to Discrete 2D Sequences.
As in subsection 2.4, we consider an input image $X \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$, a complex, analytic convolution kernel $W \in \mathcal{J}(\boldsymbol{\theta}, \kappa)$, a subsampling factor $m \in \mathbb{N} \setminus \{0\}$ and an integer $q \in \mathbb{N} \setminus \{0\}$, referred to as a *half-size*, such that max pooling operates on a grid of size $(2q + 1) \times (2q + 1)$. We seek a relationship between

$$ (3.9) \qquad \text{Y}^{\text{max}} := U_{m,q}^{\text{max}}[W](X) \quad \text{and} \quad \text{Y}^{\text{mod}} := U_m^{\text{mod}}[W](X), $$

where $U_{m,q}^{\text{max}}[W]$ ($\mathbb{R}Max$) and $U_m^{\text{mod}}[W]$ (CMod) have been respectively defined in (1.1) and (1.5). As before, in what follows we omit the parameter between square brackets.

We now use the sampling results from Proposition 2.5. Let $F_X$ and $\Psi_W \in \mathcal{V}^{(s)}$ denote the functions satisfying (2.23). Recall that the continuous versions of CMod and $\mathbb{R}Max$ operators have been defined in (3.5) and (3.6), respectively. On the one hand, we apply (2.26) with $m \leftarrow 2m$ to $\text{Y}^{\text{mod}}$. For any $\boldsymbol{n} \in \mathbb{Z}^2$,

$$ (3.10) \qquad U_m^{\text{mod}}X[n] = (F_X * \bar{\Psi}_W)(x_n) $$

$$ (3.11) \qquad = U^{\text{mod}} F_X(x_n), $$

with $x_n := 2msn$. On the other hand, we postulate that

$$ (3.12) \qquad U_{m,q}^{\text{max}}X[n] = U_r^{\text{max}}F_X(x_n) $$
Therefore, by definition of arcs on the unit circle (73), we get

$$ \angle(z^*Z_{i+1}^{(q)}) \le \angle(z^*Z_{i'}^{(q)}). \qquad (99) $$

Then, since cos is non-increasing on $[0, \pi]$, we get

$$ \cos \angle(z^* Z_{i+1}^{(q)}) \ge \cos \angle(z^* Z_{i'}^{(q)}), \qquad (100) $$

which yields the right part of (97). On the other hand, if $\angle(z^*Z_{i'}^{(q)}) \ge \pi$, a similar reasoning yields the left part of (97). Then, (96) holds.

Next, we show that, as observed in fig. 4, $g_{\max}$ is piecewise-symmetric with respect to the center value of each arc $\mathfrak{A}_i^{(q)}$, denoted by

$$ \bar{Z}_i^{(q)} := \sqrt{Z_i^{(q)} Z_{i+1}^{(q)}}. \qquad (101) $$

Let $z_1, z_2 \in \mathfrak{A}_i^{(q)}$ which are symmetric with respect to $\bar{Z}_i^{(q)}$. Therefore, there exists $z' \in \mathbb{S}^1$ such that $z_1 = \bar{Z}_i^{(q)} z'$ and $z_2 = \bar{Z}_i^{(q)} z'^*$. We now prove that

$$ g_{\max}(z_1) = g_{\max}(z_2). \qquad (102) $$

A simple calculation yields

$$ z_1^* Z_{i+1}^{(q)} = z'^* \tilde{Z}_i^{(q)} \quad \text{and} \quad z_2^* Z_i^{(q)} = (z'^* \tilde{Z}_i^{(q)})^*, \qquad (103) $$

with

$$ \tilde{Z}_i^{(q)} := (Z_i^{(q)*} \bar{Z}_i^{(q)}) = (\bar{Z}_i^{(q)*} Z_{i+1}^{(q)}). \qquad (104) $$

Therefore,

$$ \mathrm{Re}(z_1^* Z_{i+1}^{(q)}) = \mathrm{Re}(z_2^* Z_i^{(q)}). \qquad (105) $$

Since $z_1, z_2$ both belong to $\mathfrak{A}_i^{(q)}$, $g_{\max}(z_1)$ and $g_{\max}(z_2)$ satisfy (96). Then, by symmetry, (105) implies (102). One can observe from fig. 4 that $g_{\max}$ reaches its local minimum at the center of arc $\mathfrak{A}_i^{(q)}$, i.e., $\bar{Z}_i^{(q)}$. This corresponds to a point where $g_{\max}$ is non-differentiable.

We denote by $\bar{\mathfrak{A}}_i^{(q)} := [Z_i^{(q)}, \bar{Z}_i^{(q)}]_{\mathbb{S}^1}$ the first half of arc $\mathfrak{A}_i^{(q)}$. Then,

$$ \forall z \in \bar{\mathfrak{A}}_i^{(q)}, g_{\max}(z) = \mathrm{Re}(z^* Z_i^{(q)}). \qquad (106) $$

As a consequence, using symmetry, we get

$$ \int_{\mathfrak{A}_i^{(q)}} g_{\max}(z)^p \mathrm{d}\vartheta(z) = 2 \int_{\bar{\mathfrak{A}}_i^{(q)}} g_{\max}(z)^p \mathrm{d}\vartheta(z) \\ = 2 \int_{\bar{\mathfrak{A}}_i^{(q)}} \mathrm{Re}(z^* Z_i^{(q)})^p \mathrm{d}\vartheta(z). $$

By using the change of variable formula (Athreya & Lahiri, 2006, p. 81) with $z \leftarrow e^{i\eta}$, we get

$$ \int_{\mathfrak{A}_i^{(q)}} g_{\max}(z)^p \mathrm{d}\vartheta(z) = 2 \int_{H_i^{(q)}}^{\bar{H}_i^{(q)}} \cos^p(\eta - H_i^{(q)}) \mathrm{d}\eta, \qquad (107) $$

where $\bar{H}_i^{(q)} := (H_i^{(q)} + H_{i+1}^{(q)})/2$ denotes the argument of $\bar{Z}_i^{(q)}$. Then, the change of variable $\eta' \leftarrow \eta - H_i^{(q)}$ yields

$$ \int_{\mathfrak{A}_i^{(q)}} g_{\max}(z)^p \mathrm{d}\vartheta(z) = 2 \int_{0}^{\delta H_i^{(q)}/2} \cos^p \eta' \mathrm{d}\eta'. \qquad (108) $$
for a certain value of $r \in \mathbb{R}_{>0}$. Then, (3.8) implies $Y^{\text{mod}} \approx Y^{\text{max}}$. However, as explained hereafter, (3.12) is not satisfied, due to the discrete nature of the max pooling grid. According to (1.1) and (1.4), we have

$$ (3.13) \qquad U_{m,q}^{\max} X[n] = \max_{\|\boldsymbol{p}\|_{\infty} \le q} \mathrm{Re}((X * \bar{W}) \downarrow m) [2\boldsymbol{n} + \boldsymbol{p}]. $$

Therefore, according to (2.26) in Proposition 2.5, we get

$$ (3.14) \qquad U_{m,q}^{\max} X[n] = \max_{\|\boldsymbol{p}\|_{\infty} \le q} (F_X * \mathrm{Re}\,\bar{\Psi}_W)(\boldsymbol{x}_n + \boldsymbol{h}_p), $$

with

$$ (3.15) \qquad \boldsymbol{x}_n := 2ms\boldsymbol{n} \quad \text{and} \quad \boldsymbol{h}_p := m\boldsymbol{sp}. $$

By considering $r_q := ms(q + \frac{1}{2})$, we get a variant of (3.12) in which the maximum is evaluated on a discrete grid of $(2q + 1)^2$ elements, instead of the continuous region $B_{\infty}(r_q)$, as defined in (3.6) with $r \leftarrow r_q$. As a consequence, (3.7) is replaced in the discrete framework by

$$ (3.16) \qquad q \ll 2\pi/(m\kappa) \implies U_{m,q}^{\max} X[n] \approx U_m^{\text{mod}} X[n] \max_{\|\boldsymbol{p}\|_{\infty} \le q} G_X(\boldsymbol{x}_n, \boldsymbol{h}_p), $$

where we have introduced, similarly to (3.1),

$$ (3.17) \qquad G_X : (\boldsymbol{x}, \boldsymbol{h}) \mapsto \cos(\langle \boldsymbol{\nu}, \boldsymbol{h} \rangle - H_X(\boldsymbol{x})), $$

with

$$ (3.18) \qquad \boldsymbol{\nu} := \boldsymbol{\theta}/s \quad \text{and} \quad H_X := \angle (F_X * \bar{\Psi}_W), $$

where $\angle : \mathbb{C} \to [0, 2\pi[$ denotes the phase operator. Unlike the continuous case, even if the window size $r_q$ is large enough, the existence of $\boldsymbol{p} \in \{-q..q\}^2$ such that $G_X(\boldsymbol{x}_n, \boldsymbol{h}_p) = 1$ is not guaranteed, as illustrated in Figure 3 with $q = 1$. Instead, we can only seek a probabilistic estimation of the normalized mean squared error between $Y^{\text{max}}$ and $Y^{\text{mod}}$.

Approximation (3.16) implies

$$ (3.19) \qquad q \ll 2\pi/(m\kappa) \implies \|U_m^{\text{mod}}X - U_{m,q}^{\max}X\|_2 \approx \|\delta_{m,q}X\|_2, $$

where $\delta_{m,q}X \in l_2^2(\mathbb{Z}^2)$ is defined such that, for any $\boldsymbol{n} \in \mathbb{Z}^2$,

$$ (3.20) \qquad \delta_{m,q} X[n] := U_m^{\text{mod}} X[n] \left( 1 - \max_{\|\boldsymbol{p}\|_{\infty} \le q} G_X(\boldsymbol{x}_n, \boldsymbol{h}_p) \right). $$

Expression (3.19) suggests that the difference between the left and right terms can be bounded by a quantity which only depends on the product $m\kappa$ (subsampling factor $\times$ frequency localization) and the grid half-size $q$. In what follows, we establish a bound characterizing this approximation, which will be provided in Proposition 3.5.
<figure><img src="image_6.png" /><figcaption>Figure 5. $\gamma(m\theta)^2$ as a function of the kernel characteristic frequency $\theta \in [-\pi, \pi]^2$. According to Theorem 2, this quantity provides an approximate bound for the expected quadratic error between $\mathbb{R}Max$ and $\mathbb{C}Mod$ outputs. The subsampling factor $m$ has been set to 2 as in ResNet (left), and 4 as in AlexNet (right). The bright regions correspond to frequencies for which the two outputs are expected to be similar. However, in the dark regions, pathological cases such as illustrated in fig. 3 are more likely to occur.</figcaption></figure>

Next, we insert (108) into (95), and compute $\mathbb{E}[G_X^{\max}(\boldsymbol{x})^p]$ for $p \leftarrow 1$ and $p \leftarrow 2$:

$$ 
\begin{aligned}
\mathbb{E}[G_X^{\max}(\boldsymbol{x})] &= \frac{1}{\pi} \sum_{i=0}^{n_q-1} \sin \frac{\delta H_i^{(q)}}{2}; \\
\mathbb{E}[G_X^{\max}(\boldsymbol{x})^2] &= \frac{1}{2} + \frac{1}{4\pi} \sum_{i=0}^{n_q-1} \sin \delta H_i^{(q)}.
\end{aligned}
$$

We recall that $Q_X := 1 - G_X^{\max}$. By linearity of $\mathbb{E}$, we get

$$ \mathbb{E}[Q_X(\boldsymbol{x})^2] := \frac{3}{2} + \frac{1}{4\pi} \sum_{i=0}^{n_q-1} \left( \sin \delta H_i^{(q)} - 8 \sin \frac{\delta H_i^{(q)}}{2} \right), \quad (109) $$

which concludes the proof. $\square$

We consider an ideal scenario where $(Z_i^{(q)}(m\theta))_{i \in \{0, \dots, n_q-1\}}$ are evenly spaced on $\mathbb{S}^1$. Then, an order 2 Taylor expansion yields

$$ \gamma_q(m\theta) = o(1/q^2), \quad (110) $$

providing an order-two-polynomial decay rate for $Q_X(\boldsymbol{x})$, when the grid half-size $q$ increases. Figure 5 displays $\theta \mapsto \gamma_q(m\theta)^2$ for $\theta \in [-\pi, \pi]^2$, with $m = 4$ and $q = 1$ as in AlexNet. We notice that, for the major part of the Fourier domain, $\gamma_q$ remains close to 0. However, we observe a regular pattern of dark regions, which correspond to pathological frequencies where the repartition of $(Z_i^{(q)}(m\theta))_{i \in \{0, \dots, n_q-1\}}$ is unbalanced.

So far, we established a result at the pixel level. Before stating Theorem 2, which extends the result to the image level, we need the following intermediate statement.

**Proposition 7.** We consider the random variable

$$ \tilde{S}_X := \|U_m^{\text{mod}}X\|_2. \quad (111) $$

Under Hypothesis 3, for any $\boldsymbol{x} \in \mathbb{R}^2$,

* $Z_X(\boldsymbol{x})$ is independent of $\tilde{S}_X$;
For the sake of conciseness, we introduce the following notations:

$$ (3.21) \qquad A_X : (\mathbf{x}, \mathbf{h}) \mapsto (F_X * \mathrm{Re} \, \bar{\Psi}_W)(\mathbf{x} + \mathbf{h}); $$

$$ (3.22) \qquad \tilde{A}_X : (\mathbf{x}, \mathbf{h}) \mapsto |(F_X * \bar{\Psi}_W)(\mathbf{x})| G_X(\mathbf{x}, \mathbf{h}). $$

We now consider, for any $\mathbf{n} \in \mathbb{Z}^2$, the vectors $\mathbf{h}_n^{\max}$ and $\mathbf{h}_n'^{\max} \in \mathrm{ms}\{-q..q\}^2$ achieving the maximum value of $A_X(\mathbf{x}_n, \mathbf{h}_p)$ and $\tilde{A}_X(\mathbf{x}_n, \mathbf{h}_p)$ over the max pooling grid, respectively. They satisfy

$$ (3.23) \qquad A_X^{\max}(\mathbf{x}_n) := A_X(\mathbf{x}_n, \mathbf{h}_n^{\max}) = \max_{\|\mathbf{p}\|_{\infty} \le q} A_X(\mathbf{x}_n, \mathbf{h}_p); $$

$$ (3.24) \qquad \tilde{A}_X^{\max}(\mathbf{x}_n) := \tilde{A}_X(\mathbf{x}_n, \mathbf{h}_n'^{\max}) = \max_{\|\mathbf{p}\|_{\infty} \le q} \tilde{A}_X(\mathbf{x}_n, \mathbf{h}_p). $$

Then, according to (3.10) and (3.14), we get, for any $\mathbf{n} \in \mathbb{Z}^2$,

$$ (3.25) \qquad A_X^{\max}(\mathbf{x}_n) = U_{m,q}^{\max} X[n]; $$

$$ (3.26) \qquad \tilde{A}_X^{\max}(\mathbf{x}_n) = U_m^{\mathrm{mod}} X[n] \max_{\|\mathbf{p}\|_{\infty} \le q} G_X(\mathbf{x}_n, \mathbf{h}_p), $$

and (3.16) becomes

$$ (3.27) \qquad q \ll 2\pi/(m\kappa) \quad \implies \quad A_X^{\max}(\mathbf{x}_n) \approx \tilde{A}_X^{\max}(\mathbf{x}_n). $$

*Remark 3.2.* Expression (3.3) implies that, if $q \ll 2\pi/(m\kappa)$, then $A_X(\mathbf{x}_n, \mathbf{h}_p) \approx \tilde{A}_X(\mathbf{x}_n, \mathbf{h}_p)$ for all $\mathbf{p} \in \{-q..q\}^2$. However, this property does not guarantee that $A_X$ and $\tilde{A}_X$ reach their maximum in the same exact location; i.e., that $\mathbf{h}_n^{\max} = \mathbf{h}_n'^{\max}$.

The following lemma provides a bound for approximation (3.27).

**Lemma 3.3.** For any $\mathbf{x} \in \mathbb{R}^2$,

$$ (3.28) \qquad \left| A_X^{\max}(\mathbf{x}_n) - \tilde{A}_X^{\max}(\mathbf{x}_n) \right| \le \max_{\mathbf{h} \in \{\mathbf{h}_n^{\max}, \mathbf{h}_n'^{\max}\}} \left| F_0(\mathbf{x}_n + \mathbf{h}) - F_0(\mathbf{x}_n) \right|. $$

*Proof.* We apply Proposition 3.1 with $\mathbf{h} \leftarrow \mathbf{h}_n^{\max}$ and $\mathbf{h} \leftarrow \mathbf{h}_n'^{\max}$, respectively:

$$ (3.29) \qquad A_X^{\max}(\mathbf{x}_n) \le \tilde{A}_X(\mathbf{x}_n, \mathbf{h}_n^{\max}) + |F_0(\mathbf{x}_n + \mathbf{h}_n^{\max}) - F_0(\mathbf{x}_n)|; $$

$$ (3.30) \qquad \tilde{A}_X^{\max}(\mathbf{x}_n) \le A_X(\mathbf{x}_n, \mathbf{h}_n'^{\max}) + |F_0(\mathbf{x}_n + \mathbf{h}_n'^{\max}) - F_0(\mathbf{x}_n)|. $$

By construction, we have, for any $\mathbf{p} \in \{-q..q\}^2$,

$$ (3.31) \qquad \tilde{A}_X(\mathbf{x}_n, \mathbf{h}_p) \le \tilde{A}_X^{\max}(\mathbf{x}_n) \quad \text{and} \quad A_X(\mathbf{x}_n, \mathbf{h}_p) \le A_X^{\max}(\mathbf{x}_n). $$

Moreover, by definition, there exists $\mathbf{p}$ and $\mathbf{p}' \in \{-q..q\}^2$ such that $\mathbf{h}_n^{\max} = \mathbf{h}_p$ and $\mathbf{h}_n'^{\max} = \mathbf{h}_p'$. Therefore, (3.29) and (3.30) yield, respectively,

$$ (3.32) \qquad A_X^{\max}(\mathbf{x}_n) \le \tilde{A}_X^{\max}(\mathbf{x}_n) + |F_0(\mathbf{x}_n + \mathbf{h}_n^{\max}) - F_0(\mathbf{x}_n)|; $$

$$ (3.33) \qquad \tilde{A}_X^{\max}(\mathbf{x}_n) \le A_X^{\max}(\mathbf{x}_n) + |F_0(\mathbf{x}_n + \mathbf{h}_n'^{\max}) - F_0(\mathbf{x}_n)|, $$

which yields (3.28) and concludes the proof.

<figure><img src="image_9.png" /></figure>
• $Z_X(x)$, $M_X(x)$ are conditionally independent given $\tilde{S}_X$.

*Proof.* See Appendix A.12. □

Finally, Propositions 6 and 7 yield the following theorem. It provides an upper bound on the expected value of the normalized mean squared error $\tilde{P}_X^2$, such as defined in (86).

**Theorem 2** (MSE between $\mathbb{C}Mod$ and $\mathbb{R}Max$). Let $W \in \mathcal{J}(\theta, \kappa)$ denote a discrete Gabor-like filter, $m \in \mathbb{N} \setminus \{0\}$ a subsampling factor and $q \in \mathbb{N} \setminus \{0\}$ a grid half-size. We consider a stochastic process $X$ whose realizations are elements of $l_{\mathbb{R}}^2(\mathbb{Z}^2)$. We assume that condition (38) is satisfied: $\kappa \le \pi/m$. Then, under Hypotheses 1 to 3,<fnref n="2" />

$$ \mathbb{E}[\tilde{P}_X^2] \le (\beta_q(m\kappa) + \gamma_q(m\theta))^2, \quad (112) $$

where $\tilde{P}_X^2$ (86) denotes the stochastic quadratic error between $\mathbb{C}Mod$ and $\mathbb{R}Max$ outputs. We remind that $\beta_q$ and $\gamma_q$ have been introduced in (72) and (91), respectively.

*Proof.* See Appendix A.13. □

Let us analyze the bound obtained in (112). The first term, $\beta_q(m\kappa)$, accounts for the localized property of the convolution filter $W$. This term decreases linearly with the product $m\kappa$. In the limit case where $\kappa = 0$ (infinite, nonlocal filter), we get $\beta_q(m\kappa) = 0$. Note that a smaller subsampling factor $m$ allows for a larger bandwidth $\kappa$. Moreover, $\beta_q(m\kappa)$ increases linearly with the size of the max pooling grid, which is characterized by $q$. The second term, $\gamma_q(m\theta)$, accounts for the discrete nature of the max pooling grid. It strongly depends on the characteristic frequency $\theta$, as illustrated in fig. 5. According to (110), this term has a polynomial decay when $q$ increases. However, increasing the size of the max pooling grid also results in increasing the term $\beta_q(m\kappa)$, as explained above. Therefore, a tradeoff must be found to get an optimal bound.

# 4 Shift Invariance of $\mathbb{R}Max$ Outputs

In this section, we present the main theoretical claim of this paper. Based on the previous results, we provide a probabilistic measure of shift invariance for $\mathbb{R}Max$ operators. First, we state the following lemma.

**Lemma 6.** If Hypotheses 2 and 3 are satisfied, then they are also true with $X \leftarrow \mathcal{T}_u X$, for any $u \in \mathbb{R}^2$.

*Proof.* See Appendix A.14. □

We are now ready to state the main result about shift invariance of $\mathbb{R}Max$ outputs.

**Theorem 3** (Shift invariance of $\mathbb{R}Max$). We assume that the requirements stated in Theorem 2 are satisfied. Additionally, given a translation vector $u \in \mathbb{R}^2$, we consider the following random variable:

$$ \tilde{\mathbb{R}}_{X, u} := \|U_{m, q}^{\max}(\mathcal{T}_u X) - U_{m, q}^{\max} X\|_2 / \|U_m^{\text{mod}} X\|_2. \quad (113) $$

Then, under condition (38), we have

$$ \mathbb{E}[\tilde{\mathbb{R}}_{X, u}] \le 2 (\beta_q(m\kappa) + \gamma_q(m\theta)) + \alpha(\kappa u), \quad (114) $$

where $\alpha$, $\beta_q$ and $\gamma_q$ are defined in (18), (72) and (91), respectively.

*Proof.* See Appendix A.15. □

<footnote n="2">²We can easily prove that these properties are independent from the choice of sampling interval $s > 0$.</footnote>
Before stating Proposition 3.5, we consider the following hypothesis:

**Hypothesis 3.4.** There exists $h_0 \in \mathbb{R}^2$ with $\|h_0\|_2 = \sqrt{2qms}$, such that

$$ (3.34) \quad \sum_{n \in \mathbb{Z}^2} \max_{h \in \{h_n^{\text{max}}, h_n^{\prime \text{max}}\}} |F_0(\mathbf{x}_n + \mathbf{h}) - F_0(\mathbf{x}_n)|^2 \le \sum_{n \in \mathbb{Z}^2} |F_0(\mathbf{x}_n + \mathbf{h}_0) - F_0(\mathbf{x}_n)|^2. $$

The underlying idea is explained as follows. The absolute difference between $F_0(\mathbf{x}_n + \mathbf{h})$ and $F_0(\mathbf{x}_n)$ is more likely to increase with the norm of $\mathbf{h}$. For any given $\mathbf{n} \in \mathbb{Z}^2$, we have, by construction, $\|h_n^{\text{max}}\|_2 \le \sqrt{2qms}$ and $\|h_n^{\prime \text{max}}\|_2 \le \sqrt{2qms}$. Therefore, we can expect to observe

$$ (3.35) \quad \max_{h \in \{h_n^{\text{max}}, h_n^{\prime \text{max}}\}} |F_0(\mathbf{x}_n + \mathbf{h}) - F_0(\mathbf{x}_n)|^2 \le |F_0(\mathbf{x}_n + \mathbf{h}_0) - F_0(\mathbf{x}_n)|^2. $$

While this might occasionally not be true, **Hypothesis 3.4** postulates that, when summing over all the datapoints, the inequality holds.

We now formally state the result characterizing approximation (3.19).

**Proposition 3.5.** We assume that condition (2.53) is satisfied: $\kappa \le \pi/m$. Then, under **Hypothesis 3.4**,

$$ (3.36) \quad \|U_m^{\text{mod}}X - U_{m,q}^{\text{max}}X\|_2 \le \|\delta_{m,q}X\|_2 + \beta_q(m\kappa) \|U_m^{\text{mod}}X\|_2, $$

where $\beta_q: \mathbb{R}_+ \to \mathbb{R}_+$ is defined by

$$ (3.37) \qquad \beta_q : \kappa' \mapsto q\kappa'. $$

*Proof.* Let us write:

$$ \begin{align*} \|U_m^{\text{mod}}X - U_{m,q}^{\text{max}}X\|_2^2 &= \sum_{n \in \mathbb{Z}^2} \left( U_m^{\text{mod}}X[n] - U_{m,q}^{\text{max}}X[n] \right)^2 \\ &= \sum_{n \in \mathbb{Z}^2} \left( U_m^{\text{mod}}X[n] - U_m^{\text{mod}}X[n] \max_{\|\mathbf{p}\|_{\infty} \le q} G_X(\mathbf{x}_n, \mathbf{h}_p) \right. \\ & \qquad \left. + U_m^{\text{mod}}X[n] \max_{\|\mathbf{p}\|_{\infty} \le q} G_X(\mathbf{x}_n, \mathbf{h}_p) - U_{m,q}^{\text{max}}X[n] \right)^2 \\ &= \sum_{n \in \mathbb{Z}^2} \left( \delta_{m,q}X[n] + \tilde{A}_X^{\text{max}}(\mathbf{x}_n) - A_X^{\text{max}}(\mathbf{x}_n) \right)^2, \end{align*} $$

according to (3.20), (3.25) and (3.26). Then, using the triangle inequality, we get

$$ (3.38) \quad \|U_m^{\text{mod}}X - U_{m,q}^{\text{max}}X\|_2 \le \|\delta_{m,q}X\|_2 + \left( \sum_{n \in \mathbb{Z}^2} \left( \tilde{A}_X^{\text{max}}(\mathbf{x}_n) - A_X^{\text{max}}(\mathbf{x}_n) \right)^2 \right)^{1/2}. $$
In the bound established in (114), the sum $\beta_q(m\kappa) + \gamma_q(m\theta)$ accounts for the discrepancy between $\mathbb{R}Max$ and $\mathbb{C}Mod$ outputs, as stated in Theorem 2, whereas the term $\alpha(\kappa\mathbf{u})$ characterizes the stability of $\mathbb{C}Mod$ outputs, as stated in Theorem 1. If $\kappa$ is sufficiently small, then $\alpha(\kappa\mathbf{u})$ and $\beta_q(m\kappa)$ become negligible with respect to $\gamma_q(m\theta)$, and the bound can be approximated by $2\gamma_q(m\theta)$. Theorem 3 therefore provides a validity domain for shift invariance of $\mathbb{R}Max$ operators, as illustrated in fig. 5 with $q=1$.

**Remark 5.** The stochastic discrepancy introduced in (113) is estimated relatively to the $\mathbb{C}Mod$ output. This choice is motivated by the perfect shift invariance of its norm, as shown in Proposition 3.

**Remark 6.** In practice, most of the time max pooling is performed on a grid of size $3 \times 3$; therefore $q=1$. For the sake of conciseness, we shall sometimes drop $q$ in the notations, which implicitly means $q=1$.

# 5 Adaptation to Multichannel Convolution Operators

In this section, we adapt Theorems 1 to 3 to multichannel inputs (e.g., RGB images), employed in conventional CNNs such as AlexNet or ResNet.

First, we define multichannel $\mathbb{R}Max$ and $\mathbb{C}Mod$ operators relatively to (1) and (5). We denote by $K$ and $L \in \mathbb{N} \setminus \{0\}$ the number of input and output channels, respectively. Additionally, we consider a multichannel convolution tensor

$$ \mathbf{W} := (W_{lk})_{l \in \{0, \dots, L-1\}, k \in \{0, \dots, K-1\}} \in (l_{\mathbb{C}}^2(Z^2))^{L \times K}. \quad (115) $$

Multichannel $\mathbb{R}Max$ and $\mathbb{C}Mod$ operators take as input images, denoted by

$$ \mathbf{X} := (X_k)_{k \in \{0, \dots, K-1\}} \in (l_{\mathbb{R}}^2(Z^2))^K. \quad (116) $$

They are defined, for any given output channel $l \in \{0, \dots, L-1\}$, by

$$ U_{m, q, l}^{\max}[\mathbf{W}] : \mathbf{X} \mapsto \mathrm{MaxPool}_q \left( \sum_{k=0}^{K-1} (X_k * \mathrm{Re} \, \overline{W}_{lk}) \downarrow m \right); \quad (117) $$

$$ U_{m, l}^{\mathrm{mod}}[\mathbf{W}] : \mathbf{X} \mapsto \left| \sum_{k=0}^{K-1} (X_k * \overline{W}_{lk}) \downarrow (2m) \right|, \quad (118) $$

where $m, q \in \mathbb{N} \setminus \{0\}$ respectively denote a subsampling factor and the max pooling grid half-size. Analogously to (49) for single-channel inputs, we now consider

$$ Y_l^{\max} := U_{m, q, l}^{\max}[\mathbf{W}](\mathbf{X}) \quad \text{and} \quad Y_l^{\mathrm{mod}} := U_{m, l}^{\mathrm{mod}}[\mathbf{W}](\mathbf{X}). \quad (119) $$

Again, in what follows we omit the parameter between square brackets. To apply Theorems 1 to 3 to the current setting on the $l$-th output channel, we need the following hypotheses.

**Hypothesis 4** (Monochorome filters). Let

$$ \widetilde{W}_l := \frac{1}{K} \sum_{k=0}^{K-1} W_{lk} \quad (120) $$

denote the mean kernel of the $l$-th output channel. Then, there exists $\mu_l \in \mathbb{R}^K$ such that

$$ \forall k \in \{0, \dots, K-1\}, W_{lk} = \mu_{lk} \widetilde{W}_l. \quad (121) $$

**Hypothesis 5** (Gabor-like filters). There exists a bandwidth $\kappa > 0$ satisfying $\kappa \le \pi/m$ and a frequency vector $\theta_l \in [-\pi, \pi]^2$ such that

$$ \widetilde{W}_l \in \mathcal{J}(\theta_l, \kappa). \quad (122) $$

Note that the bandwidth $\kappa$ is not indexed by $l$, because it shall later be assumed to be shared across the output channels. Then, under Hypothesis 4, $Y_l^{\max}$ and $Y_l^{\mathrm{mod}}$ are the outputs of single-channel $\mathbb{R}Max$ and $\mathbb{C}Mod$ operators, as introduced in (1) and (5):

$$ Y_l^{\max} = U_{m, q}^{\max}[\widetilde{W}_l](X_l^{\mathrm{lum}}) \quad \text{and} \quad Y_l^{\mathrm{mod}} = U_m^{\mathrm{mod}}[\widetilde{W}_l](X_l^{\mathrm{lum}}), \quad (123) $$
Furthermore, Lemma 3.3 yields

$$ (3.39) \quad \sum_{n \in \mathbb{Z}^2} \left( \tilde{A}_X^{\max}(\mathbf{x}_n) - A_X^{\max}(\mathbf{x}_n) \right)^2 \le \sum_{n \in \mathbb{Z}^2} \max_{h \in \{h_n^{\max}, h_n^{\max}\}} |F_0(\mathbf{x}_n + h) - F_0(\mathbf{x}_n)|^2 $$

$$ (3.40) \quad \le \sum_{n \in \mathbb{Z}^2} |F_0(\mathbf{x}_n + h_0) - F_0(\mathbf{x}_n)|^2, $$

according to Hypothesis 3.4. Now, since (2.53) is satisfied, we can use Lemma 2.8 (2.44) with $h \leftarrow h_0$. We get

$$ \begin{aligned} \sum_{n \in \mathbb{Z}^2} \left( \tilde{A}_X^{\max}(\mathbf{x}_n) - A_X^{\max}(\mathbf{x}_n) \right)^2 &\le \frac{1}{4m^2s^2} \|T_{h_0}F_0 - F_0\|_{L^2}^2 \\ &\le \alpha(\kappa h_0/s)^2 \frac{1}{4m^2s^2} \|F_0\|_{L^2}^2 \quad (\text{acc. to Proposition 2.3}) \\ &= \alpha(\kappa h_0/s)^2 \|U_m^{\text{mod}}X\|_2^2. \quad (\text{acc. to Lemma 2.8 (2.45)}) \end{aligned} $$

Since, according to Hypothesis 3.4, $\|h_0\|_2 = \sqrt{2}qms$, it comes that $\|h_0\|_1 = 2qms$. Therefore,

$$ (3.41) \quad \alpha(\kappa h_0/s)^2 = \frac{\kappa^2 \|h_0\|_1^2}{4s^2} = (qm\kappa)^2, $$

which yields

$$ (3.42) \quad \sum_{n \in \mathbb{Z}^2} \left( \tilde{A}_X^{\max}(\mathbf{x}_n) - A_X^{\max}(\mathbf{x}_n) \right)^2 \le \beta_q(m\kappa)^2 \|U_m^{\text{mod}}X\|_2^2. $$

Finally, plugging (3.42) into (3.38) concludes the proof.

We now seek a probabilistic estimation of $\|\delta_{m,q}X\|_2$. For this purpose, we first reformulate the problem using the unit circle $S^1 \subset \mathbb{C}$, before introducing a probabilistic framework in subsection 3.4.

**3.3. Notations on the Unit Circle.** In what follows, for any $z \in \mathbb{C} \setminus \{0\}$, we denote by $\angle z \in [0, 2\pi[$ the argument of $z$. For any $z, z' \in S^1$, the angle between $z$ and $z'$ is given by $\angle(z^*z')$. We then denote by $[z, z']_{S^1} \subset S^1$ the arc on the unit circle going from $z$ to $z'$ counterclockwise:

$$ (3.43) \quad [z, z']_{S^1} := \{z'' \in S^1 \mid \angle(z^*z'') \le \angle(z^*z')\}. $$

We remind readers that $\mathbf{x}_n$ and $h_p \in \mathbb{R}^2$ have been defined in (3.15). By using the relation $\cos \alpha = \text{Re}(e^{i\alpha})$, (3.17) becomes, for any $n \in \mathbb{Z}^2$ and any $p \in \{-q..q\}^2$,

$$ (3.44) \quad G_X(\mathbf{x}_n, h_p) = \text{Re}(Z_X^*(\mathbf{x}_n) Z_p(m\theta)), $$

where we have defined the following functions with outputs on the unit circle:

$$ (3.45) \quad Z_X : \mathbf{x} \mapsto e^{i H_X(\mathbf{x})} \quad \text{and} \quad Z_p : \omega \mapsto e^{i \langle \omega, p \rangle}, $$
where $X_l^{\text{lum}} \in l_\mathbb{R}^2(Z^2)$ ("luminance" image) is defined as the following linear combination:

$$X_l^{\text{lum}} := \sum_{k=0}^{K-1} \mu_{lk} X_k. \qquad (124)$$

The results established for single-channel inputs can therefore be extended to multichannel operators. Specifically, we get the following corollaries to Theorems 1 to 3.

**Corollary 2 (Shift invariance of CMod).** For a given output channel $l \in \{0, \dots, L-1\}$, we postulate Hypotheses 4 and 5. Then, for any input image $\mathbf{X} \in (l_\mathbb{R}^2(Z^2))^K$ with finite support and any translation vector $\mathbf{u} \in \mathbb{R}^2$,

$$\|U_{m,l}^{\text{mod}}(\mathcal{T}_{\mathbf{u}}\mathbf{X}) - U_{m,l}^{\text{mod}}\mathbf{X}\|_2 \le \alpha(\kappa\mathbf{u}) \|U_{m,l}^{\text{mod}}\mathbf{X}\|_2, \qquad (125)$$

where $\alpha$ has been defined in (18).

**Corollary 3 (MSE between CMod and RMax).** As in Corollary 2, we postulate Hypotheses 4 and 5. Again, we assume that condition (38) is satisfied: $\kappa \le \pi/m$. Additionally, we consider $\mathbf{X}$ as a stack of $K$ discrete stochastic processes, and assume Hypotheses 1 to 3 with $X \leftarrow X_l^{\text{lum}}$ and $W \leftarrow \tilde{W}_l$. Then,

$$\mathbb{E}[\tilde{\mathbf{P}}_{\mathbf{X},l}^2] \le (\beta_q(m\kappa) + \gamma_q(m\theta_l))^2, \qquad (126)$$

where we have defined the following random variable:

$$\tilde{\mathbf{P}}_{\mathbf{X},l} := \|U_{m,l}^{\text{mod}}\mathbf{X} - U_{m,l}^{\text{max}}\mathbf{X}\|_2 / \|U_{m,l}^{\text{mod}}\mathbf{X}\|_2. \qquad (127)$$

**Corollary 4 (Shift invariance of RMax).** We assume that the requirements stated in Corollary 3 are satisfied. Then, for any translation vector $\mathbf{u} \in \mathbb{R}^2$,

$$\mathbb{E}[\tilde{\mathbf{R}}_{\mathbf{X}, \mathbf{u}, l}] \le 2(\beta_q(m\kappa) + \gamma_q(m\theta_l)) + \alpha(\kappa\mathbf{u}), \qquad (128)$$

where we have defined the following random variable:

$$\tilde{\mathbf{R}}_{\mathbf{X}, \mathbf{u}, l} := \|U_{m,l}^{\text{max}}(\mathcal{T}_{\mathbf{u}}\mathbf{X}) - U_{m,l}^{\text{max}}\mathbf{X}\|_2 / \|U_{m,l}^{\text{mod}}\mathbf{X}\|_2. \qquad (129)$$

**Remark 7.** In the above results, we used a translation operator on multichannel tensors, obtained by applying $\mathcal{T}_{\mathbf{u}}$, as defined in (28), to each channel $X_k$.

# 6 A Case Study Implementing the Dual-Tree Complex Wavelet Packet Transform

In this section, we experimentally validate the results stated in Theorems 1 to 3. To this end, we consider a fully-deterministic scenario implementing the dual-tree complex wavelet packet transform (DT-CWPT), which exhibit characteristics akin to those observed in the initial convolution layer of freely-trained CNNs such as AlexNet or ResNet. In particular, as stated in Section 6.1, DT-CWPT achieves subsampled convolutions with oriented band-pass filters tiling the Fourier domain into overlapping square windows. As such, it provides a convenient framework to experimentally validate our theoretical findings in a controlled environment. Then, in Section 6.2, we build CMod and RMax operators based on DT-CWPT convolution kernels.

## 6.1 Main Properties

In what follows, we outline the principal characteristics of DT-CWPT. A detailed description of the transform itself is provided in Appendix C.1, whereas the results presented hereafter are formally established in Appendices C.2 and C.3.

For a given decomposition depth $J \in \mathbb{N} \setminus \{0\}$, DT-CWPT achieves subsampled convolutions with $4 \times 4^J$ oriented band-pass filters that tile the Fourier domain into overlapping square windows of size

$$\kappa_J := \pi/m_J, \quad \text{with} \quad m_J := 2^{J-1}. \qquad (130)$$
<figure><img src="image_10.png" /><figcaption>Figure 3. Search for the maximum value of $h \mapsto G_X(x, h)$ over a discrete grid of size $3 \times 3$, i.e., $q = 1$. This figure displays 3 examples with different frequencies $\nu := \theta/s$ and phases $H_X(x)$. Hopefully the result will be close to the true maximum (left), but there are some pathological cases in which all points in the grid fall into pits (middle and right).</figcaption></figure>

where $H_X$ denotes the phase of $F_X * \bar{\Psi}_W$ as introduced in (3.18). On the one hand, $Z_X(x_n)$ is the phase (represented on the unit circle $S^1$) of the complex wavelet transform $F_X * \bar{\Psi}_W$ at location $x_n$. On the other hand, $Z_p(m\theta)$ approximates the phase shift between any two evaluations of $F_X * \bar{\Psi}_W$ at locations $x, x'$ such that $x' - x = h_p$. This however is only true if we assume that $\Psi_W$ exhibits slow amplitude variations. Then, $G_X(x_n, h_p)$ approximates the cosine of the phase of $F_X * \bar{\Psi}_W$ at location $x_n + h_p$.

According to (3.16), $\max_{\|p\|_{\infty} \le q} G_X(x_n, h_p)$ approximates the ratio between $\mathbb{R}Max$ and $\mathbb{C}Mod$ outputs at discrete location $n \in \mathbb{Z}^2$. The intuition behind this is that max pooling seeks a point in a discrete grid around $x_n$ where the phase of $F_X * \bar{\Psi}_W$ is the closest to 1, thereby maximizing the amount of energy on the real part of the signal. Assuming slow amplitude variations of $\Psi_W$, the result therefore approximates the modulus of the complex coefficients.

To get an estimation of $\delta_{m,q}X[n]$ (3.20), we will exploit the following property. If the phases $Z_p(m\theta)$ for $p \in \{-q..q\}^2$ are well distributed on the unit circle, then the values of $G_X(x_n, h_p)$ are evenly spread out on $[-1, 1]$. Therefore, its maximum value is more likely to be close to 1, and (3.20) becomes

$$ (3.46) \qquad \delta_{m,q}X[n] \ll U_m^{\text{mod}}X[n] \quad \forall n \in \mathbb{Z}^2. $$

Let $n_q := (2q + 1)^2$ denote the number of evaluation points for the max pooling operator. For any $\omega \in \mathbb{R}^2$, we consider a sequence of values on $S^1$, denoted by $(Z_i^{(q)}(\omega))_{i \in \{0..n_q-1\}}$, obtained by sorting $\{Z_p(\omega)\}_{p \in \{-q..q\}^2}$ (3.45) in ascending order of their arguments:

$$ (3.47) \qquad 0 = H_0^{(q)}(\omega) \le \dots \le H_{n_q-1}^{(q)}(\omega) < 2\pi, $$

where $H_i^{(q)}(\omega)$ denotes the phase of $Z_i^{(q)}(\omega)$. Besides, we close the loop with $H_{n_q}^{(q)}(\omega) := 2\pi$
More specifically, considering an input image $X \in l_2^2(\mathbb{Z}^2)$, it produces a set of $4 \times 4^J$ output feature maps

$$ \mathbf{D}^{(J)} := (D_l^{\lambda^{(J)}}, D_l^{\lambda^{(J)}}, D_l^{\kappa^{(J)}}, D_l^{\kappa^{(J)}})_{l \in \{0, \dots, 4^J - 1\}}, \qquad (131) $$

where each arrow points to the Fourier quadrant where the feature map's energy is concentrated. Moreover, as stated in Proposition 10, for any $l \in \{0, \dots, 4^J - 1\}$, there exists $W_l^{\lambda^{(J)}} \in l_2^2(\mathbb{Z}^2)$ such that

$$ D_l^{\lambda^{(J)}} = (X * \overline{W}_l^{\lambda^{(J)}})^* \downarrow 2^J. \qquad (132) $$

An interesting property is that each kernel $W_l^{\lambda^{(J)}}$ approximately satisfies

$$ W_l^{\lambda^{(J)}} \in \mathcal{J}(\theta_l^{(J)}, \kappa_J) \qquad (133) $$

for a certain characteristic frequency $\theta_l^{(J)} \in [0, \pi]^2$. In other words, it approximately behaves as a Gabor-like filter in the discrete framework (11). Moreover, each kernel corresponds to a different frequency, thereby covering the top-right quadrant of the Fourier domain. Similar results can be established for the other three Fourier quadrants. Graphical representations of $\mathbf{W}^{\lambda^{(J)}} := (W_l^{\lambda^{(J)}})_{l \in \{0, \dots, 4^J - 1\}}$ and $\mathbf{W}^{\lambda^{(J)}} := (W_l^{\lambda^{(J)}})_{l \in \{0, \dots, 4^J - 1\}}$ are provided in fig. 6 with $J = 2$ (fig. 6a, 32 filters) and $J = 3$ (fig. 6b, 128 filters).

<figure><img src="image_7.png" /><figcaption>Figure 6. Real part of the convolution kernels $\mathbf{W}^{\lambda^{(J)}}$, $\mathbf{W}^{\lambda^{(J)}}$, with $J = 2$ (32 filters, $m_J = 2$) and $J = 3$ (128 filters, $m_J = 4$), respectively. The kernels have been computed using Q-shift orthogonal QMFs of length 10 (Kingsbury, 2003). The kernels have been respectively cropped to size $11 \times 11$ and $19 \times 19$, for the sake of legibility. Note that the filters displayed in (a) and (b) share similarities with those found in, respectively, ResNet ($m = 2$) and AlexNet ($m = 4$), after training with ImageNet.</figcaption></figure>

The $\mathbb{R}$Max and $\mathbb{C}$Mod operators implemented in our experiments respectively satisfy (1) and (5) with $W \leftarrow W_l^{\lambda^{(J)}}$ or $W_l^{\lambda^{(J)}}$, and $m \leftarrow m_J$. Note that increasing the decomposition depth $J$, and therefore the subsampling factor $m_J$, results in a decreased Fourier support size $\kappa_J$, therefore matching the condition stated in (38) $\kappa \leftarrow \kappa_J$ and $m \leftarrow m_J$.
and $Z_{n_q}^{(q)}(\omega) := 1$. Then, we split $\mathbb{S}^1$ into $n_q$ arcs delimited by $Z_i^{(q)}(\omega)$:

$$ (3.48) \qquad \mathfrak{A}_i^{(q)}(\omega) := \begin{cases} [Z_i^{(q)}(\omega), Z_{i+1}^{(q)}(\omega)]_{\mathbb{S}^1} & \text{if } H_{i+1}^{(q)}(\omega) - H_i^{(q)}(\omega) < 2\pi; \\ \mathbb{S}^1 & \text{otherwise.} \end{cases} $$

Finally, for any $i \in \{0..n_q - 1\}$, we denote by

$$ (3.49) \qquad \delta H_i^{(q)} : \omega \mapsto H_{i+1}^{(q)}(\omega) - H_i^{(q)}(\omega) $$

the function computing the angular measure of arc $\mathfrak{A}_i^{(q)}(\omega)$, for any $\omega \in \mathbb{R}^2$.

**3.4. Probabilistic Framework.** From now on, input X is considered as discrete 2D stochastic processes. In order to “randomize” $F_X$ introduced in (2.23), we define a continuous stochastic process from X, denoted by $F_X$, such that

$$ (3.50) \qquad \forall \mathbf{x} \in \mathbb{R}^2, F_X(\mathbf{x}) := \sum_{n \in \mathbb{Z}^2} X[n] \Phi_n^{(s)}(\mathbf{x}). $$

Now, we consider the following stochastic processes, which are parameterized by X:

$$ (3.51) \qquad M_X := |F_X * \bar{\Psi}_W|; \quad H_X := \angle(F_X * \bar{\Psi}_W); \quad Z_X := e^{iH_X}, $$

and, for any $\mathbf{p} \in \{-q..q\}^2$,

$$ (3.52) \qquad \mathbf{G}_{X,p} := \mathrm{Re}(\mathbf{Z}_X^* \mathbf{Z}_p(m\theta)); \qquad \mathbf{G}_X^{\max} := \max_{\|\mathbf{p}\|_{\infty} \le q} \mathbf{G}_{X,p}, $$

where the deterministic function $Z_p$ has been defined in (3.45).

**Remark 3.6.** $H_X(\mathbf{x})$ is ill-defined if $M_X(\mathbf{x}) = 0$. To overcome this, it is designed to follow a uniform conditional probability distribution on $[0, 2\pi[, \text{ given } M_X(\mathbf{x}) = 0$. Moreover, we impose the following conditional independence, for any $n \in \mathbb{N} \setminus \{0\}$ and $\mathbf{x}, \mathbf{y}_0, \dots, \mathbf{y}_{n-1} \in \mathbb{R}^2$:

$$ (3.53) \qquad H_X(\mathbf{x}) \perp \mathbf{M} \mid M_X(\mathbf{x}) = 0, \quad \text{with} \quad \mathbf{M} := (M_X(\mathbf{y}_0), \dots, M_X(\mathbf{y}_{n-1}))^\top. $$

Finally, we impose the following relationship between $H_{T_u X}$ and $H_X$, for any $\mathbf{u} \in \mathbb{R}^2$:

$$ (3.54) \qquad M_{T_u X}(\mathbf{x}) = 0 \implies H_{T_u X}(\mathbf{x}) = T_{su} H_X(\mathbf{x}). $$

For any $\mathbf{x} \in \mathbb{R}^2$, $F_X(\mathbf{x})$ (2.23) and $H_X(\mathbf{x})$ (3.18) are respectively drawn from $F_X(\mathbf{x})$ and $H_X(\mathbf{x})$. Then, $Z_X(\mathbf{x})$ (3.45) is a realization of $Z_X(\mathbf{x})$. Consequently, according to (3.44), $G_X(\mathbf{x}, \mathbf{h}_p)$ is a realization of $G_{X,p}(\mathbf{x})$. Besides, according to the definition of $\mathbb{C}Mod$ in (1.5) and $\mathbf{x}_n$ in (3.15), Proposition 2.5 with $m \leftarrow 2m$ implies that

$$ (3.55) \qquad M_X(\mathbf{x}_n) = U_m^{\mathrm{mod}} X[n]. $$
<figure><img src="image_8.png" /><figcaption>Figure 7. Detailed illustration of the RMax (a) and CMod (b) operators based on DT-CWPT, with $J = 3$ decomposition stages. The numbers between modules correspond to the number of feature maps, height and width. The orange modules represent subsampled convolutions using one of the four 2D filter banks $G^{[0-3]}$, such as introduced in (243). The FB index is indicated between square brackets. The RMax model (a) only computes the real part of the dual-tree coefficients, and the last stage of decomposition is performed without subsampling (red modules). Additionally, the blue modules represent linear combinations of feature maps such as described in (248).</figcaption></figure>

**Remark 8.** Because $X$ is real-valued, the feature maps $D_l^{\perp(J)}$ and $D_l^{\perp(J)}$ are the respective complex conjugates of $D_l^{\perp(J)}$ and $D_l^{\perp(J)}$, and thus do not need to be explicitly computed. Then, we can easily show that $W_l^{\perp(J)}$ and $W_l^{\perp(J)}$ are also the complex conjugates of $W_l^{\perp(J)}$ and $W_l^{\perp(J)}$, respectively.

## 6.2 DT-CWPT-Based RMax and CMod Operators

According to (130), (132), and (133), we can apply Theorems 1 to 3 to the dual-tree framework. More precisely, for any output channel $l \in \{0, \dots, 4^J - 1\}$, we consider the following RMax and CMod operators:

$$ U_l^{\max\perp} : X \mapsto \mathrm{MaxPool} \left( \left( X * \mathrm{Re} \overline{W_l^{\perp(J)}} \right) \downarrow 2^{J-1} \right); \qquad (134) $$

$$ U_l^{\mathrm{mod}\perp} : X \mapsto \left| \left( X * \overline{W_l^{\perp(J)}} \right) \downarrow 2^J \right|. \qquad (135) $$

Using the notations introduced in (5) and (1), we have

$$ U_l^{\max\perp} = U_{m_J}^{\max} [W_l^{\perp(J)}] \quad \text{and} \quad U_l^{\mathrm{mod}\perp} := U_{m_J}^{\mathrm{mod}} [W_l^{\perp(J)}], \qquad (136) $$

where we have defined $m_J := 2^{J-1}$. Note that, following Remark 6, we have omitted the grid half-size $q$, which is equal to 1 (max pooling operates on a grid of size $3 \times 3$). Furthermore, for the sake of brevity, we have omitted the depth $J$ in the above notations.

**Remark 9.** Both $U_l^{\max\perp}$ and $U_l^{\mathrm{mod}\perp}$ are implemented using DT-CWPT with $J$ decomposition stages. However, in (134), the subsampling factor is equal to $2^{J-1}$, instead of $2^J$, as stated in (132). In order to accommodate this property of RMax operators, the last stage of DT-CWPT decomposition is carried out without subsampling, resulting in higher redundancy. This is similar to the concept of stationary wavelet transform as described by Nason & Silverman (1995). Furthermore, only the real component of the wavelet feature maps is preserved. On the other hand, $U_l^{\mathrm{mod}\perp}$ implements a fully-decimated wavelet packet transform, and keeps both real and imaginary parts. Figure 7 illustrates these technical details.
We remind that $\theta \in [-\pi, \pi]^2$ and $\kappa \in ]0, 2\pi]$ respectively denote the center and size of the Fourier support of the complex kernel $W \in \mathcal{J}(\theta, \kappa)$. To compute the expected discrepancy between $Y^{\text{max}}$ and $Y^{\text{mod}}$, we assume that

$$ (3.56) \qquad \|\theta\|_2 \gg 2\pi/N; $$

$$ (3.57) \qquad \|\theta\|_2 \gg \kappa, $$

where $N \in \mathbb{N} \setminus \{0\}$ denotes the support size of input images. These assumptions exclude low-frequency filters from the scope of our study. We then state the following hypotheses, for which a justification is provided in Appendix A.

**Hypothesis 3.7.** For any $\boldsymbol{x} \in \mathbb{R}^2$, $Z_X(\boldsymbol{x})$ is uniformly distributed on $\mathbb{S}^1$.

**Hypothesis 3.8.** For any $n \in \mathbb{N} \setminus \{0\}$ and $\boldsymbol{x}, \boldsymbol{y}_0, \dots, \boldsymbol{y}_{n-1} \in \mathbb{R}^2$, the random variables $M_X(\boldsymbol{y}_i)$ for $i \in \{0 \dots n-1\}$ are jointly independent of $Z_X(\boldsymbol{x})$.

**3.5. Expected Quadratic Error between $\mathbb{R}^{\text{Max}}$ and $\mathbb{C}^{\text{Mod}}$.** In this section, we propose to estimate the expected value of the stochastic quadratic error $\tilde{P}_X^2$, defined such that

$$ (3.58) \qquad \tilde{P}_X := \|U_m^{\text{mod}}X - U_{m,q}^{\text{max}}X\|_2 / \|U_m^{\text{mod}}X\|_2. $$

According to (3.9), this is an estimation of the relative error between $Y^{\text{mod}}$ and $Y^{\text{max}}$.

First, let us reformulate $\delta_{m,q}X$, introduced in (3.20), using the probabilistic framework. According to (3.44) and (3.52), we have, for any $\boldsymbol{n} \in \mathbb{Z}^2$,

$$ (3.59) \qquad \delta_{m,q}X[\boldsymbol{n}] := U_m^{\text{mod}}X[\boldsymbol{n}] (1 - G_X^{\text{max}}(\boldsymbol{x}_n)). $$

We now consider the stochastic process

$$ (3.60) \qquad Q_X := 1 - G_X^{\text{max}}, $$

and the random variable

$$ (3.61) \qquad \tilde{Q}_X := \|\delta_{m,q}X\|_2 / \|U_m^{\text{mod}}X\|_2. $$

The next steps are as follows: (1) at the pixel level, show that $\mathbb{E}[Q_X(\boldsymbol{x})^2]$ depends on the subsampling factor $m$ and the filter frequency $\theta$, and remains close to zero with some exceptions; (2) at the image level, show that the expected value of $\tilde{Q}_X^2$ is equal to the latter quantity; (3) use Proposition 3.5, which implies that $\tilde{P}_X \approx \tilde{Q}_X$, to deduce an upper bound on the expected value of $\tilde{P}_X^2$.

The first point is established in Proposition 3.9 below, and the two remaining ones are the purpose of Theorem 3.11.

**Proposition 3.9.** *Assuming Hypothesis 3.7, the expected value of $Q_X(\boldsymbol{x})^2$ is independent from the choice of $\boldsymbol{x} \in \mathbb{R}^2$, and*

$$ (3.62) \qquad \mathbb{E}[Q_X(\boldsymbol{x})^2] = \gamma_q(m\theta)^2, $$
## 6.3 Experiments and Results

We implemented the $\mathbb{R}Max$ and $\mathbb{C}Mod$ operators $U_l^{\max\uparrow}$ and $U_l^{\mod\uparrow}$, as introduced in (134) and (135), with both $J = 2$ and 3 stages of wavelet packet decomposition. To cover the whole frequency plane, we also implemented similar operators, denoted by $U_l^{\max\downarrow}$ and $U_l^{\mod\downarrow}$. They are associated with the convolution filters $W_l^{\downarrow(J)}$, introduced in Proposition 10, with energy being located in the bottom-right quadrant. However, as explained in Remark 8, we did not need to deal with the two other quadrants (negative $x$-values). Using the validation set of ImageNet-1K (Russakovsky et al., 2015), ($N := 50000$ images), we measured the mean discrepancy between $\mathbb{R}Max$ and $\mathbb{C}Mod$ outputs, and evaluated the shift invariance of both models. Dual-tree decompositions have been performed with Q-shift orthogonal filters of length 10 (Kingsbury, 2003).

### 6.3.1 MSE between $\mathbb{R}Max$ and $\mathbb{C}Mod$

Each image $n \in \{0, ..., N-1\}$ in the dataset was converted to grayscale, from which a center crop of size $224 \times 224$ was extracted. We denote by $X_n \in l_{\mathbb{R}}^2(\mathbb{Z}^2)$ the resulting input feature map. For any $l \in \{0, ..., 4^J - 1\}$, we denote by

$$ Y_{nl}^{\max\uparrow} := U_l^{\max\uparrow}(X_n) \quad \text{and} \quad Y_{nl}^{\mod\uparrow} := U_l^{\mod\uparrow}(X_n) \quad (137) $$

the outputs of the $l$-th $\mathbb{R}Max$ and $\mathbb{C}Mod$ operators as defined in (134) and (135), respectively. We adopt similar notations for the bottom-right Fourier quadrant. Then, the normalized mean squared error between $Y_{nl}^{\mod\uparrow}$ and $Y_{nl}^{\max\uparrow}$ was computed. It is defined by the square of

$$ \rho_{nl}^{\uparrow} := \| Y_{nl}^{\mod\uparrow} - Y_{nl}^{\max\uparrow} \|_2 / \| Y_{nl}^{\mod\uparrow} \|_2. \quad (138) $$

Finally, for each output channel $l$, an empirical estimate for $E[\tilde{\rho}_X^2]$, introduced in (86), was obtained by averaging $\tilde{\rho}_{nl}^{\uparrow 2}$ over the whole dataset. We denote by $\tilde{\rho}_l^{\uparrow 2}$ the corresponding quantity.

Since $U_l^{\max\uparrow}$ and $U_l^{\mod\uparrow}$ are parameterized by $W_l^{\uparrow(J)}$, it follows that $\tilde{\rho}_l^{\uparrow 2}$ depends on the filter's characteristic frequency $\theta_l^{(J)}$ (133). According to Proposition 11, these frequencies form a regular grid in the top-right quadrant of Fourier domain. This provides a visual representation of $\tilde{\rho}_l^{\uparrow 2}$, as shown in fig. 8. This figure also displays $\tilde{\rho}_l^{\uparrow 2}$, corresponding to the bottom-right quadrant. The half-plane of negative $x$-values has simply been symmetrized, following Remark 8. We can observe a regular pattern of dark spots. More precisely, high discrepancies between max pooling and modulus seem to occur when the energy of $W_l^{\uparrow(J)}$ or $W_l^{\downarrow(J)}$ overlaps a dark region of fig. 5. This result corroborates Theorem 3, which states that high discrepancies are expected for certain pathological frequencies, due to the search for a maximum value over a discrete grid.

### 6.3.2 Shift invariance

For each input image previously converted to grayscale, two crops of size $224 \times 224$ were extracted, such that the corresponding sequences $X_n$ and $X'_n$ are shifted by one pixel along the $x$-axis. From these inputs, the following quantity was then computed:

$$ \rho_{nl}^{\max\uparrow} := \| Y_{nl}^{\max\uparrow} - Y_{nl}^{\max\downarrow} \|_2 / \| Y_{nl}^{\mod\uparrow} \|_2, \quad (139) $$

where $Y_{nl}^{\max\uparrow}$ satisfies (136) with $X_n \leftarrow X'_n$. Finally, for each output channel $l \in \{0, ..., 4^J - 1\}$, an empirical estimate for $E[\tilde{\rho}_{X, u}]$, satisfying (113) with $u = (1, 0)^T$, was obtained by averaging $\tilde{\rho}_{nl}^{\max\uparrow}$ over the whole dataset. We denote by $\tilde{\rho}_l^{\max\uparrow}$ the corresponding quantity. We point out that shift invariance is measured relatively to the norm of the $\mathbb{C}Mod$ output, as explained in Remark 5.

On the other hand, the same procedure was applied to the $\mathbb{C}Mod$ operators:

$$ \rho_{nl}^{\mod\uparrow} := \| Y_{nl}^{\mod\uparrow} - Y_{nl}^{\mod\downarrow} \|_2 / \| Y_{nl}^{\mod\uparrow} \|_2, \quad (140) $$
where we have defined

$$ (3.63) \qquad \gamma_q : \boldsymbol{\omega} \mapsto \sqrt{\frac{3}{2} + \frac{1}{4\pi} \sum_{i=0}^{n_q-1} \left( \sin \delta H_i^{(q)}(\boldsymbol{\omega}) - 8 \sin \frac{\delta H_i^{(q)}(\boldsymbol{\omega})}{2} \right)}, $$

with $\delta H_i^{(q)}(\boldsymbol{\omega}) \in [0, 2\pi]$ (3.49) being the length of arc $\mathfrak{A}_i^{(q)}(\boldsymbol{\omega})$.

*Proof.* For the sake of readability, in this proof we omit the argument of functions $Z_p$ (3.45), $Z_i^{(q)}$, $H_i^{(q)}$ (3.47), $\mathfrak{A}_i^{(q)}$ (3.48), and $\delta H_i^{(q)}$ (3.49); we assume they are evaluated at $\boldsymbol{\omega} \leftarrow m\boldsymbol{\theta}$. We consider the “Lebesgue” Borel $\sigma$-algebra on $\mathbb{S}^1$ generated by $\{[z, z']_{\mathbb{S}^1} \mid z, z' \in \mathbb{S}^1\} \cup \{\mathbb{S}^1\}$, on which we have defined the angular measure $\vartheta$ such that $\vartheta(\mathbb{S}^1) := 2\pi$, and

$$ (3.64) \qquad \forall z, z' \in \mathbb{S}^1, \vartheta([z, z']_{\mathbb{S}^1}) := \angle(z^*z'). $$

For any $p \in \mathbb{N} \setminus \{0\}$, we compute the $p$-th moment of $G_X^{\max}(\boldsymbol{x})$ defined in (3.52). By considering

$$ (3.65) \qquad \begin{gathered} g_{\max} : \mathbb{S}^1 \to [-1, 1] \\ z \mapsto \max_{\|\boldsymbol{p}\|_{\infty} \le q} \mathrm{Re}(z^* Z_{\boldsymbol{p}}), \end{gathered} $$

we get $G_X^{\max}(\boldsymbol{x}) = g_{\max}(Z_X(\boldsymbol{x}))$. A visual representation of $g_{\max}$ is provided in Figure 4, for two different values of $\boldsymbol{\theta}$. According to Hypothesis 3.7, $Z_X(\boldsymbol{x})$ follows a uniform distribution on $\mathbb{S}^1$. Therefore,

$$ (3.66) \qquad \mathbb{E} [G_X^{\max}(\boldsymbol{x})^p] = \frac{1}{2\pi} \int_{\mathbb{S}^1} g_{\max}(z)^p \, d\vartheta(z), $$

which proves that $\mathbb{E}[G_X^{\max}(\boldsymbol{x})^p]$ does not depend on $\boldsymbol{x}$. Let us split the unit circle $\mathbb{S}^1$ into the arcs $\mathfrak{A}_0^{(q)}, \dots, \mathfrak{A}_{n_q-1}^{(q)}$ such as introduced in (3.48):

$$ (3.67) \qquad \mathbb{E} [G_X^{\max}(\boldsymbol{x})^p] = \frac{1}{2\pi} \sum_{i=0}^{n_q-1} \int_{\mathfrak{A}_i^{(q)}} g_{\max}(z)^p \, d\vartheta(z). $$

Let $i \in \{0..n_q - 1\}$. We show that

$$ (3.68) \qquad \forall z \in \mathfrak{A}_i^{(q)}, g_{\max}(z) = \max \left( \mathrm{Re}(z^* Z_i^{(q)}), \mathrm{Re}(z^* Z_{i+1}^{(q)}) \right). $$

Let $z \in \mathfrak{A}_i^{(q)}$ and $i' \notin \{i, i+1\}$. We prove that

$$ (3.69) \qquad \mathrm{Re}(z^* Z_{i'}^{(q)}) \le \mathrm{Re}(z^* Z_i^{(q)}) \quad \text{or} \quad \mathrm{Re}(z^* Z_{i'}^{(q)}) \le \mathrm{Re}(z^* Z_{i+1}^{(q)}). $$

On the one hand, we assume that $\angle(z^*Z_{i'}^{(q)}) \le \pi$. By design of $(Z_i^{(q)})_{i \in \{0..n_q-1\}}$, we have

$$ (3.70) \qquad Z_{i+1}^{(q)} \in [z, Z_{i'}^{(q)}]_{\mathbb{S}^1}. $$
<figure><img src="image_9.png" /><figcaption>Figure 8. Empirical estimates of the normalized mean squared error between $\mathbb{R}$Max and CMod outputs, computed on ImageNet-1K (validation set). For each channel $l \in \{0, \dots, 4^J - 1\}$, $\tilde{\rho}_l^{\wedge 2}$ is plotted as a grayscale pixel centered in $\theta_l^{(J)}$ such as introduced in eq. (133) (top-right quadrant). Similarly, $\tilde{\rho}_l^{\times 2}$ is plotted in the bottom-right quadrant. Finally, the bottom- and top-left quadrants ($\tilde{\rho}_l^{\wedge \times 2}$ and $\tilde{\rho}_l^{\times \wedge 2}$) are simply obtained by symmetrizing the figures. Since the subsampling factor $m_J$ is equal to $2^{J-1}$, these experimental results can be compared with the left and right parts of fig. 5. Note that the low-pass filters have been discarded because they are outside the scope of this study.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Figure 9. Shift invariance of $\mathbb{R}$Max and CMod outputs, computed on ImageNet 2012 (validation set). For each $l \in \{0, \dots, 4^J - 1\}$, $\tilde{\rho}_l^{\max \wedge}$ (fig. 9a) and $\tilde{\rho}_l^{\max \times}$ (fig. 9b) are plotted by applying the same procedure as in fig. 8.</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>Figure 4. Top: 2D representation of $h \mapsto G_X(x_n, h)$ (3.17), for two different values of $\theta \in \mathbb{R}^2$, $q = 1$ and arbitrary values of $m \in \mathbb{N} \setminus \{0\}$ and $s \in \mathbb{R} \setminus \{0\}$. Assuming the plots are centered around $h = 0$, each point materializes a location $h_p$ in the max pooling grid, for $p \in \{-q..q\}^2$. The desirable situation occurs when one of these locations falls near a ridge (bright areas), in which case the outputs produced by $\mathbb{R}Max$ and $\mathbb{C}Mod$ are similar—see (3.16). Each number $i \in \{0..8\}$ represents the rank of $Z_p \in \mathbb{S}^1$ (3.45), when these values are sorted by ascending order of their arguments (3.47). If rank $i$ is affected to location $h_p$, then we have $Z_p = Z_i^{(q)}$. Bottom: polar representations of $g_{\max}: \mathbb{S}^1 \to [-1, 1]$ (3.65), corresponding to the same settings. The closer the curve is from the outer ring, the more likely some points $h_p$ will fall near a ridge of $G_X$. (a) Case where the values $Z_p$ are roughly evenly distributed on $\mathbb{S}^1$. (b) Case where these values are concentrated in a small portion of the unit circle. The most extreme cases occurs when $Z_p = 1$ for any $p$. Figure 3 (middle and right) depicts two such situations.</figcaption></figure>

Therefore, by definition of arcs on the unit circle (3.43), we get

$$ (3.71) \qquad \angle(z^*Z_{i+1}^{(q)}) \le \angle(z^*Z_{i'}^{(q)}). $$

Then, since cos is non-increasing on $[0, \pi]$, we get

$$ (3.72) \qquad \cos \angle(z^*Z_{i+1}^{(q)}) \ge \cos \angle(z^*Z_{i'}^{(q)}), $$
and $\tilde{\rho}_l^{\text{mod}\perp}$ was obtained as before by averaging $\rho_{nl}^{\text{mod}\perp}$ over the whole dataset.

A visual representation of $\tilde{\rho}_l^{\text{max}\perp}$ and $\tilde{\rho}_l^{\text{mod}\perp}$ are provided in fig. 9 (as well as the other Fourier quadrants). Two observations can be drawn here. (1) When the filter is horizontally oriented, the corresponding output is highly stable with respect to horizontal shifts. This can be explained by noticing that such kernels perform low-pass filtering along the x-axis. The exact transposed phenomenon occurs for vertical shifts. (2) Elsewhere, we observe that high discrepancies between $\mathbb{R}Max$ and $\mathbb{C}Mod$ outputs (fig. 8) are correlated with shift instability of $\mathbb{R}Max$ (fig. 9, top). This is in line with (112) and (114) in Theorems 2 and 3. Note that $\mathbb{C}Mod$ outputs are nearly shift invariant regardless the characteristic frequency $\theta_l^{(J)}$ (fig. 9, bottom), as predicted by Theorem 1 (39).

# 7 Conclusion

In this paper, we investigated the shift invariance properties captured by the max pooling operator, when applied to a convolution layer with Gabor-like kernels. We established a validity domain for near-shift invariance and confirmed our predictions through an experimental setting based on the dual-tree complex wavelet packet transform. Our results indicate that the $\mathbb{C}Mod$ operator can serve as a stable proxy for $\mathbb{R}Max$, extracting comparable features—except at certain filter frequencies, where potential degeneracies may arise after max pooling. In this context, a follow-up study applying these principles to real-world architectures was published as a conference paper (Leterme et al., 2024).

A link was missing between real- and complex-valued convolutions in CNNs. By comparing the outputs of $\mathbb{C}Mod$ and $\mathbb{R}Max$ operators, we established a connection between these two worlds, creating opportunities for extensions of the results obtained for complex wavelet transforms. To paraphrase Tygert et al. (2016), the correspondence between standard real-valued CNNs (using max pooling) and complex wavelets is no longer “just a vague analogy.”

# References

* Krishna B. Athreya and Soumendra N. Lahiri. *Measure Theory and Probability Theory*, volume 19. Springer, 2006.

* Aharon Azulay and Yair Weiss. Why do deep convolutional networks generalize so poorly to small image transformations? *Journal of Machine Learning Research*, 20(184):1–25, 2019.

* Radu Balan, Maneesh Singh, and Dongmian Zou. Lipschitz Properties for Deep Convolutional Networks. *Contemporary Mathematics*, 706:129–151, 2018.

* Ilker Bayram and Ivan W. Selesnick. On the Dual-Tree Complex Wavelet Packet and M-Band Transforms. *IEEE Transactions on Signal Processing*, 56(6):2298–2310, 2008.

* Alberto Bietti. Approximation and Learning with Deep Convolutional Models: A Kernel Perspective. In *International Conference on Learning Representations*, 2022.

* Alberto Bietti and Julien Mairal. Invariance and Stability of Deep Convolutional Representations. In *Advances in Neural Information Processing Systems*, pp. 1622, 2017.

* Alberto Bietti and Julien Mairal. Group invariance, stability to deformations, and complexity of deep convolutional representations. *The Journal of Machine Learning Research*, 20(1):876–924, 2019a.

* Alberto Bietti and Julien Mairal. On the Inductive Bias of Neural Tangent Kernels. In *Advances in Neural Information Processing Systems*, volume 32. Curran Associates, Inc., 2019b.

* Joan Bruna and Stéphane Mallat. Invariant Scattering Convolution Networks. *IEEE Transactions on Pattern Analysis and Machine Intelligence*, 35(8):1872–1886, 2013.

* Jameson Cahill, Joseph W. Iverson, Dustin G. Mixon, and Daniel Packer. Group-Invariant Max Filtering. *Foundations of Computational Mathematics*, 2024.
which yields the right part of (3.69). On the other hand, if $\angle(z^*Z_{i'}^{(q)}) \ge \pi$, a similar reasoning yields the left part of (3.69). Then, (3.68) holds.

Now, we show that, as observed in Figure 4, $g_{\max}$ is piecewise-symmetric with respect to the center value of each arc $\mathfrak{A}_i^{(q)}$, denoted by

$$ (3.73) \qquad \bar{Z}_i^{(q)} := \sqrt{Z_i^{(q)} Z_{i+1}^{(q)}}. $$

Let $z_1, z_2 \in \mathfrak{A}_i^{(q)}$ which are symmetric with respect to $\bar{Z}_i^{(q)}$. Therefore, there exists $z' \in \mathbb{S}^1$ such that $z_1 = \bar{Z}_i^{(q)}z'$ and $z_2 = \bar{Z}_i^{(q)}z'^*$. We now prove that

$$ (3.74) \qquad g_{\max}(z_1) = g_{\max}(z_2). $$

A simple calculation yields

$$ (3.75) \qquad z_1^* Z_{i+1}^{(q)} = z'^* \tilde{Z}_i^{(q)} \quad \text{and} \quad z_2^* Z_i^{(q)} = (z'^* \tilde{Z}_i^{(q)})^*, $$

with

$$ (3.76) \qquad \tilde{Z}_i^{(q)} := (Z_i^{(q)*} \bar{Z}_i^{(q)}) = (\bar{Z}_i^{(q)*} Z_{i+1}^{(q)}). $$

Therefore,

$$ (3.77) \qquad \mathrm{Re}(z_1^* Z_{i+1}^{(q)}) = \mathrm{Re}(z_2^* Z_i^{(q)}). $$

Since $z_1, z_2$ both belong to $\mathfrak{A}_i^{(q)}$, $g_{\max}(z_1)$ and $g_{\max}(z_2)$ satisfy (3.68). Then, by symmetry, (3.77) implies (3.74). One can observe from Figure 4 that $g_{\max}$ reaches its local minimum at the center of arc $\mathfrak{A}_i^{(q)}$, i.e., $\bar{Z}_i^{(q)}$. This corresponds to a point where $g_{\max}$ is non-differentiable.

We denote by $\bar{\mathfrak{A}}_i^{(q)} := [Z_i^{(q)}, \bar{Z}_i^{(q)}]_{\mathbb{S}^1}$ the first half of arc $\mathfrak{A}_i^{(q)}$. Then,

$$ (3.78) \qquad \forall z \in \bar{\mathfrak{A}}_i^{(q)}, g_{\max}(z) = \mathrm{Re}(z^* Z_i^{(q)}). $$

As a consequence, using symmetry, we get

$$ \int_{\mathfrak{A}_i^{(q)}} g_{\max}(z)^p \mathrm{d}\vartheta(z) = 2 \int_{\mathfrak{A}_i^{(q)}} g_{\max}(z)^p \mathrm{d}\vartheta(z) \\ = 2 \int_{\bar{\mathfrak{A}}_i^{(q)}} \mathrm{Re}(z^* Z_i^{(q)})^p \mathrm{d}\vartheta(z). $$

By using the change of variable formula [1, p. 81] with $z \leftarrow e^{i\eta}$, we get

$$ (3.79) \qquad \int_{\mathfrak{A}_i^{(q)}} g_{\max}(z)^p \mathrm{d}\vartheta(z) = 2 \int_{H_i^{(q)}}^{\bar{H}_i^{(q)}} \cos^p(\eta - H_i^{(q)}) \mathrm{d}\eta, $$
Zhengdao Chen. Multi-Layer Neural Networks as Trainable Ladders of Hilbert Spaces. In *Proceedings of the 40th International Conference on Machine Learning*, pp. 4294–4329, 2023.

Fergal Cotter and Nick Kingsbury. A Learnable Scatternet: Locally Invariant Convolutional Layers. In *2019 IEEE International Conference on Image Processing (ICIP)*, pp. 350–354, 2019.

Wojciech Czaja and Weilin Li. Analysis of time-frequency scattering transforms. *Applied and Computational Harmonic Analysis*, 47(1):149–171, 2019.

Wojciech Czaja and Weilin Li. Rotationally Invariant Time-Frequency Scattering Transforms. *Journal of Fourier Analysis and Applications*, 26(1):4, 2020.

Wojciech Czaja, Brandon Kolstoe, and David Koralov. Translation-Invariant Behavior of General Scattering Transforms. arXiv:2410.21030, 2024.

Alexey Dosovitskiy, Lucas Beyer, Alexander Kolesnikov, Dirk Weissenborn, Xiaohua Zhai, Thomas Unterthiner, Mostafa Dehghani, Matthias Minderer, Georg Heigold, Sylvain Gelly, Jakob Uszkoreit, and Neil Houlsby. An Image is Worth 16x16 Words: Transformers for Image Recognition at Scale. In *International Conference on Learning Representations*. arXiv, 2021.

Fernando Gama, Alejandro Ribeiro, and Joan Bruna. Stability of Graph Scattering Transforms. In *Advances in Neural Information Processing Systems*, volume 32. Curran Associates, Inc., 2019.

Shanel Gauthier, Benjamin Thérien, Laurent Alsène-Racicot, Muawiz Chaudhary, Irina Rish, Eugene Belilovsky, Michael Eickenberg, and Guy Wolf. Parametric Scattering Networks. In *Proceedings of the IEEE/CVF Conference on Computer Vision and Pattern Recognition*, pp. 5749–5758, 2022.

Kavya Gupta, Fateh Kaakai, Beatrice Pesquet-Popescu, Jean-Christophe Pesquet, and Fragkiskos D. Malliaros. Multivariate Lipschitz Analysis of the Stability of Neural Networks. *Frontiers in Signal Processing*, 2, 2022.

Paul R. Halmos. *Measure Theory*. Springer, 2013. ISBN 978-1-4684-9440-2.

Ali Hassani, Steven Walton, Nikhil Shah, Abulikemu Abuduweili, Jiachen Li, and Humphrey Shi. Escaping the Big Data Paradigm with Compact Transformers. arXiv:2104.05704, 2022.

J.P. Havlicek, J.W. Havlicek, and A.C. Bovik. The analytic image. In *Proceedings of International Conference on Image Processing*, volume 2, pp. 446–449 vol.2, 1997.

Ke Huang and Selin Aviyente. Wavelet feature selection for image classification. *IEEE Transactions on Image Processing*, 17(9):1709–1720, 2008.

Nick Kingsbury. Complex wavelets for shift invariant analysis and filtering of signals. *Applied and computational harmonic analysis*, 10(3):234–253, 2001.

Nick Kingsbury. Design of Q-shift complex wavelets for image processing using frequency domain energy minimization. In *Proceedings International Conference on Image Processing*, volume 1, pp. I–1013, 2003.

Nick Kingsbury and Julian Magarey. Wavelet Transforms in Image Processing. In *Signal Analysis and Prediction*, Applied and Numerical Harmonic Analysis, pp. 27–46. Birkhäuser, Boston, MA, 1998. ISBN 978-1-4612-1768-8.

Andrew Laine and Jian Fan. Texture Classification by Wavelet Packet Signatures. *IEEE Transactions on Pattern Analysis and Machine Intelligence*, 15(11):1186–1191, 1993.

Yann LeCun, Léon Bottou, Yoshua Bengio, and Patrick Haffner. Gradient-based learning applied to document recognition. *Proceedings of the IEEE*, 86(11):2278–2323, 1998.

Yann LeCun, Yoshua Bengio, and Geoffrey Hinton. Deep learning. *Nature*, 521(7553):436–444, 2015.
where $\overline{H}_i^{(q)} := (H_i^{(q)} + H_{i+1}^{(q)})/2$ denotes the argument of $\overline{Z}_i^{(q)}$. Then, the change of variable $\eta' \leftarrow \eta - H_i^{(q)}$ yields

$$ (3.80) \qquad \int_{\mathbb{Q}_i^{(q)}} g_{\max}(z)^p d\vartheta(z) = 2 \int_0^{\delta H_i^{(q)}/2} \cos^p \eta' d\eta'. $$

Now, we insert (3.80) into (3.67), and compute $\mathbb{E}[G_X^{\max}(\boldsymbol{x})^p]$ for $p \leftarrow 1$ and $p \leftarrow 2$:

$$ \mathbb{E}[G_X^{\max}(\boldsymbol{x})] = \frac{1}{\pi} \sum_{i=0}^{n_q-1} \sin \frac{\delta H_i^{(q)}}{2}; $$

$$ \mathbb{E}[G_X^{\max}(\boldsymbol{x})^2] = \frac{1}{2} + \frac{1}{4\pi} \sum_{i=0}^{n_q-1} \sin \delta H_i^{(q)}. $$

We recall that $Q_X := 1 - G_X^{\max}$. By linearity of $\mathbb{E}$, we get

$$ (3.81) \qquad \mathbb{E}[Q_X(\boldsymbol{x})^2] := \frac{3}{2} + \frac{1}{4\pi} \sum_{i=0}^{n_q-1} \left( \sin \delta H_i^{(q)} - 8 \sin \frac{\delta H_i^{(q)}}{2} \right), $$

which concludes the proof.

We consider an ideal scenario where $(Z_i^{(q)}(m\boldsymbol{\theta}))_{i \in \{0..n_q-1\}}$ are evenly spaced on $\mathbb{S}^1$. Then, an order 2 Taylor expansion yields

$$ (3.82) \qquad \gamma_q(m\boldsymbol{\theta}) = o(1/q^2), $$

providing an order-two-polynomial decay rate for $Q_X(\boldsymbol{x})$, when the grid half-size $q$ increases. Figure 5 displays $\boldsymbol{\theta} \mapsto \gamma_q(m\boldsymbol{\theta})^2$ for $\boldsymbol{\theta} \in [-\pi, \pi]^2$, with $m = 4$ and $q = 1$ as in AlexNet. We notice that, for the major part of the Fourier domain, $\gamma_q$ remains close to 0. However, we observe a regular pattern of dark regions, which correspond to pathological frequencies where the repartition of $(Z_i^{(q)}(m\boldsymbol{\theta}))_{i \in \{0..n_q-1\}}$ is unbalanced.

So far, we established a result at the pixel level. Before stating Theorem 3.11, which extends the result to the image level, we need the following intermediate statement.

**Proposition 3.10.** We consider the random variable

$$ (3.83) \qquad \tilde{S}_X := \|U_m^{\text{mod}}X\|_2. $$

Under Hypothesis 3.8, for any $\boldsymbol{x} \in \mathbb{R}^2$,

* $Z_X(\boldsymbol{x})$ is independent of $\tilde{S}_X$;

* $Z_X(\boldsymbol{x})$, $M_X(\boldsymbol{x})$ are conditionally independent given $\tilde{S}_X$.

*Proof.* We suppose that Hypothesis 3.8 is satisfied and we consider $\boldsymbol{x} \in \mathbb{R}^2$. For a given $n \in \mathbb{N} \setminus \{0\}$, we introduce the random variable

$$ (3.84) \qquad \tilde{S}_{X, n} := \sqrt{\sum_{\|\boldsymbol{p}\|_{\infty} \le n} M_X(\boldsymbol{x}_p)^2}. $$
Hubert Leterme, Kévin Polisano, Valérie Perrier, and Karteek Alahari. From CNNs to Shift-Invariant Twin Models Based on Complex Wavelets. In *2024 32nd European Signal Processing Conference (EUSIPCO)*, 2024.

Kunchang Li, Yali Wang, Junhao Zhang, Peng Gao, Guanglu Song, Yu Liu, Hongsheng Li, and Yu Qiao. UniFormer: Unifying Convolution and Self-Attention for Visual Recognition. *IEEE Transactions on Pattern Analysis and Machine Intelligence*, 45(10):12581-12600, 2023.

Bin Liao and Fen Peng. Rotation-invariant texture features extraction using Dual-Tree Complex Wavelet Transform. In *2010 International Conference on Information, Networking and Automation (ICINA)*, volume 1, 2010.

Stéphane Mallat. *A Wavelet Tour of Signal Processing: The Sparse Way*. Academic Press, 2009.

Stéphane Mallat. Group invariant scattering. *Communications on Pure and Applied Mathematics*, 65(10):1331-1398, 2012.

Stéphane Mallat. Understanding deep convolutional networks. *Philosophical Transactions of the Royal Society A: Mathematical, Physical and Engineering Sciences*, 374(2065):20150203, 2016.

Yves Meyer. Principe d'incertitude, bases hilbertiennes et algèbres d'opérateurs. In *Séminaire Bourbaki*, volume 662. Société Mathematique de France, 1985.

G. P. Nason and B. W. Silverman. The Stationary Wavelet Transform and some Statistical Applications. In *Wavelets and Statistics*, Lecture Notes in Statistics, pp. 281-299. Springer, 1995.

Edouard Oyallon. *Analyzing and Introducing Structures in Deep Convolutional Neural Networks*. Doctoral thesis, Paris Sciences et Lettres, 2017.

Edouard Oyallon, Eugene Belilovsky, Sergey Zagoruyko, and Michal Valko. Compressing the Input for CNNs with the First-Order Scattering Transform. In *Proceedings of the European Conference on Computer Vision (ECCV)*, 2018.

Kun Il Park and M. Park. *Fundamentals of Probability and Stochastic Processes with Applications to Communications*. Springer, 2018.

Juan C. Pérez, Motasem Alfarra, Guillaume Jeanneret, Adel Bibi, Ali Thabet, Bernard Ghanem, and Pablo Arbeláez. Gabor Layers Enhance Network Robustness. In *Proceedings of the European Conference on Computer Vision (ECCV)*, 2020. ISBN 978-3-030-58545-7.

Michael Perlmutter, Feng Gao, Guy Wolf, and Matthew Hirn. Geometric Wavelet Scattering Networks on Compact Riemannian Manifolds. In *Proceedings of The First Mathematical and Scientific Machine Learning Conference*, pp. 570-604. PMLR, 2020.

Stefan Pittner and Sagar V. Kamarthi. Feature extraction from wavelet coefficients for pattern recognition tasks. *IEEE Transactions on Pattern Analysis and Machine Intelligence*, 21(1):83-88, 1999.

Mehang Rai and Pablo Rivas. A Review of Convolutional Neural Networks and Gabor Filters in Object Recognition. In *2020 International Conference on Computational Science and Computational Intelligence (CSCI)*, pp. 1560-1567, 2020.

Olga Russakovsky, Jia Deng, Hao Su, Jonathan Krause, Sanjeev Satheesh, Sean Ma, Zhiheng Huang, Andrej Karpathy, Aditya Khosla, Michael Bernstein, Alexander C. Berg, and Li Fei-Fei. ImageNet Large Scale Visual Recognition Challenge. *International Journal of Computer Vision*, 115(3):211-252, 2015.

Meyer Scetbon and Zaid Harchaoui. Harmonic Decompositions of Convolutional Networks. In *International Conference on Machine Learning*, pp. 8522-8532. PMLR, 2020.

Ivan W. Selesnick, Richard Baraniuk, and Nick Kingsbury. The dual-tree complex wavelet transform. *IEEE Signal Processing Magazine*, 22(6):123-151, 2005.
<figure><img src="image_12.png" /><figcaption>Figure 5. $\gamma(m\theta)^2$ as a function of the kernel characteristic frequency $\theta \in [-\pi, \pi]^2$. According to Theorem 3.11, this quantity provides an approximate bound for the expected quadratic error between $\mathbb{R}Max$ and $\mathbb{C}Mod$ outputs. The subsampling factor $m$ has been set to 2 as in ResNet (left), and 4 as in AlexNet (right). The bright regions correspond to frequencies for which the two outputs are expected to be similar. However, in the dark regions, pathological cases such as illustrated in Figure 3 are more likely to occur.</figcaption></figure>

According to Hypothesis 3.8, $Z_X(\boldsymbol{x})$ is jointly independent of $M_X(\boldsymbol{x}_p)$ for $p \in \{-n..n\}^2$. Therefore, by composition, $Z_X(\boldsymbol{x})$ is also independent of $\tilde{S}_{X,n}$. Moreover, according to (3.55) and (3.83), $\tilde{S}_{X,n}$ converges almost surely towards $\tilde{S}_X$, which proves independence between $Z_X(\boldsymbol{x})$ and $\tilde{S}_X$.

Now, we prove conditional independence between $Z_X(\boldsymbol{x})$ and $M_X(\boldsymbol{x})$ given $\tilde{S}_X$. According to Hypothesis 3.8,

$$ (3.85) \qquad (M_X(\boldsymbol{x}), \tilde{S}_{X,n}) \perp Z_X(\boldsymbol{x}), $$

where $\perp$ stands for independence. This is because $\tilde{S}_{X,n}$ only depends on a finite number of $M_X(\boldsymbol{x}_p)$. Therefore,

$$ (3.86) \qquad Z_X(\boldsymbol{x}) \perp M_X(\boldsymbol{x}) \mid \tilde{S}_{X,n}. $$

Finally, since $\tilde{S}_{X,n}$ converges almost surely towards $\tilde{S}_X$, it comes that $Z_X(\boldsymbol{x})$ and $M_X(\boldsymbol{x})$ are conditionally independent given $\tilde{S}_X$.

Finally, Propositions 3.9 and 3.10 yield the following theorem. It provides an upper bound on the expected value of the normalized mean squared error $\tilde{P}_X^2$, such as defined in (3.58).

**Theorem 3.11 (MSE between $\mathbb{C}Mod$ and $\mathbb{R}Max$).** Let $W \in \mathcal{J}(\theta, \kappa)$ denote a discrete Gabor-like filter, $m \in \mathbb{N} \setminus \{0\}$ a subsampling factor and $q \in \mathbb{N} \setminus \{0\}$ a grid half-size. We consider a stochastic process $X$ whose realizations are elements of $l_{\mathbb{R}}^2(\mathbb{Z}^2)$. We assume that
C.E. Shannon. Communication in the Presence of Noise. *Proceedings of the IRE*, 37(1):10-21, 1949.

Laurent Sifre and Stéphane Mallat. Rotation, Scaling and Deformation Invariant Scattering for Texture Discrimination. In *Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition*, pp. 1233-1240, 2013.

Amarjot Singh and Nick Kingsbury. Dual-Tree wavelet scattering network with parametric log transformation for object classification. In *ICASSP, IEEE International Conference on Acoustics, Speech and Signal Processing - Proceedings*, 2017.

Christian Szegedy, Wojciech Zaremba, Ilya Sutskever, Joan Bruna, Dumitru Erhan, Ian Goodfellow, and Rob Fergus. Intriguing properties of neural networks. *arXiv:1312.6199 [cs]*, 2014.

Antonio Torralba and Aude Oliva. Statistics of natural image categories. *Network: Computation in Neural Systems*, 14(3):391-412, 2003.

Mark Tygert, Joan Bruna, Soumith Chintala, Yann LeCun, Serkan Piantino, and Arthur Szlam. A Mathematical Motivation for Complex-Valued Convolutional Networks. *Neural Computation*, 28(5):815-825, 2016.

Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N. Gomez, Lukasz Kaiser, and Illia Polosukhin. Attention is all you need. In *Advances in Neural Information Processing Systems*, pp. 5998-6008, 2017.

Martin Vetterli. Wavelets, approximation, and compression. *IEEE Signal Processing Magazine*, 18(5):59-73, 2001.

Aladin Virmaux and Kevin Scaman. Lipschitz regularity of deep neural networks: Analysis and efficient estimation. In *Advances in Neural Information Processing Systems*, 2018.

Irène Waldspurger. *Wavelet Transform Modulus : Phase Retrieval and Scattering*. Doctoral thesis, Ecole normale supérieure, Paris, 2015.

Thomas Wiatowski and Helmut Bölcskei. A Mathematical Theory of Deep Convolutional Neural Networks for Feature Extraction. *IEEE Transactions on Information Theory*, 64(3):1845-1866, 2018.

Haiping Wu, Bin Xiao, Noel Codella, Mengchen Liu, Xiyang Dai, Lu Yuan, and Lei Zhang. CvT: Introducing Convolutions to Vision Transformers. In *Proceedings of the IEEE/CVF International Conference on Computer Vision*, pp. 22-31, 2021.

Gary G. Yen. Wavelet packet feature extraction for vibration monitoring. *IEEE Transactions on Industrial Electronics*, 47(3):650-667, 2000.

Lirong Yin, Lei Wang, Siyu Lu, Ruiyang Wang, Youshuai Yang, Bo Yang, Shan Liu, Ahmed AlSanad, Salman A. AlQahtani, Zhengtong Yin, Xiaolu Li, Xiaobing Chen, and Wenfeng Zheng. Convolution-Transformer for Image Feature Extraction. *Computer Modeling in Engineering & Sciences*, 141(1):87-106, 2024.

Jason Yosinski, Jeff Clune, Yoshua Bengio, and Hod Lipson. How transferable are features in deep neural networks? In *Advances in Neural Information Processing Systems*, pp. 3320-3328, 2014.

Kun Yuan, Shaopeng Guo, Ziwei Liu, Aojun Zhou, Fengwei Yu, and Wei Wu. Incorporating Convolution Designs Into Visual Transformers. In *Proceedings of the IEEE/CVF International Conference on Computer Vision*, pp. 579-588, 2021.

John Zarka, Louis Thiry, Tomás Angles, and Stéphane Mallat. Deep Network Classification by Scattering and Homotopy Dictionary Learning. In *International Conference on Learning Representations*, 2020.

John Zarka, Florentin Guth, and Stéphane Mallat. Separation and Concentration in Deep Networks. In *International Conference on Learning Representations*, 2021.
condition (2.53) is satisfied: $\kappa \le \pi/m$. Then, under Hypotheses 3.4, 3.7, and 3.8,<fnref n="2" />

$$ (3.87) \qquad \mathbb{E}[\tilde{\mathbb{P}}_X^2] \le (\beta_q(m\kappa) + \gamma_q(m\theta))^2, $$

where $\tilde{\mathbb{P}}_X^2$ (3.58) denotes the stochastic quadratic error between $\mathbb{C}Mod$ and $\mathbb{R}Max$ outputs. We remind that $\beta_q$ and $\gamma_q$ have been introduced in (3.37) and (3.63), respectively.

*Proof.* We consider $\mathbf{n} \in \mathbb{Z}^2$. By construction, $\mathbb{Q}_X(\mathbf{x}_n) := 1 - \mathbb{G}_X^{\max}(\mathbf{x}_n)$ only depends on $\mathbb{Z}_X(\mathbf{x}_n)$. Therefore, under Hypothesis 3.8, Proposition 3.10 implies

$$ (3.88) \qquad \mathbb{Q}_X(\mathbf{x}_n) \perp \mathbb{M}_X(\mathbf{x}_n) | \tilde{\mathbb{S}}_X^2 \quad \text{and} \quad \mathbb{Q}_X(\mathbf{x}_n) \perp \tilde{\mathbb{S}}_X^2. $$

Besides, we introduce

$$ (3.89) \qquad \tilde{\Delta}_X := \|\delta_{m,q}X\|_2, $$

where $\delta_{m,q}X$ is defined in (3.59). Then, using the linearity of $\mathbb{E}$, we get

$$ \begin{aligned} \mathbb{E}[\tilde{\Delta}_X^2 | \tilde{\mathbb{S}}_X^2 = \sigma] &= \sum_{\mathbf{n} \in \mathbb{Z}^2} \mathbb{E}[\delta_{m,q}[\mathbf{n}]^2 | \tilde{\mathbb{S}}_X^2 = \sigma] \\ &= \sum_{\mathbf{n} \in \mathbb{Z}^2} \mathbb{E}[U_{m,l}^{\text{mod}} X[\mathbf{n}]^2 (1 - \mathbb{G}_X^{\max}(\mathbf{x}_n))^2 | \tilde{\mathbb{S}}_X^2 = \sigma] \\ &= \sum_{\mathbf{n} \in \mathbb{Z}^2} \mathbb{E}[\mathbb{M}_X(\mathbf{x}_n)^2 \mathbb{Q}_X(\mathbf{x}_n)^2 | \tilde{\mathbb{S}}_X^2 = \sigma] \quad (\text{acc. to (3.55) and (3.60)}) \\ &= \sum_{\mathbf{n} \in \mathbb{Z}^2} \mathbb{E}[\mathbb{M}_X(\mathbf{x}_n)^2 | \tilde{\mathbb{S}}_X^2 = \sigma] \mathbb{E}[\mathbb{Q}_X(\mathbf{x}_n)^2] \quad (\text{acc. to (3.88)}). \end{aligned} $$

According to (3.55) and (3.83), we have

$$ (3.90) \qquad \sum_{\mathbf{n} \in \mathbb{Z}^2} \mathbb{M}_X(\mathbf{x}_n)^2 = \|U_m^{\text{mod}}X\|_2^2 = \tilde{\mathbb{S}}_X^2. $$

Therefore, using again the linearity of $\mathbb{E}$, we get

$$ \begin{aligned} \mathbb{E}[\tilde{\Delta}_X^2 | \tilde{\mathbb{S}}_X^2 = \sigma] &= \mathbb{E}[\tilde{\mathbb{S}}_X^2 | \tilde{\mathbb{S}}_X^2 = \sigma] \mathbb{E}[\mathbb{Q}_X(\mathbf{x}_n)^2] \\ &= \sigma \cdot \mathbb{E}[\mathbb{Q}_X(\mathbf{x}_n)^2]. \end{aligned} $$

Under Hypothesis 3.7, Proposition 3.9 yields

$$ (3.91) \qquad \mathbb{E}[\tilde{\Delta}_X^2 | \tilde{\mathbb{S}}_X^2 = \sigma] = \sigma \cdot \gamma_q(m\theta)^2. $$

Besides, we can reformulate $\tilde{\mathbb{Q}}_X$ such as defined in (3.61): $\tilde{\mathbb{Q}}_X = \tilde{\Delta}_X/\tilde{\mathbb{S}}_X$. Therefore,

$$ (3.92) \qquad \mathbb{E}[\tilde{\mathbb{Q}}_X^2 | \tilde{\mathbb{S}}_X^2 = \sigma] = \frac{1}{\sigma} \mathbb{E}[\tilde{\Delta}_X^2 | \tilde{\mathbb{S}}_X^2 = \sigma] = \gamma_q(m\theta)^2. $$

<footnote n="2">²We can easily prove that these properties are independent from the choice of sampling interval $s > 0$.</footnote>
Richard Zhang. Making Convolutional Networks Shift-Invariant Again. In *International Conference on Machine Learning*, 2019.

Dongmian Zou and Gilad Lerman. Graph convolutional neural networks via scattering. *Applied and Computational Harmonic Analysis*, 49(3):1046–1074, 2020.

Dongmian Zou, Radu Balan, and Maneesh Singh. On Lipschitz Bounds of General Convolutional Neural Networks. *IEEE Transactions on Information Theory*, 66(3):1738–1759, 2020.

Xueyan Zou, Fanyi Xiao, Zhiding Yu, Yuheng Li, and Yong Jae Lee. Delving Deeper into Anti-Aliasing in ConvNets. *International Journal of Computer Vision*, 131(1):67–81, 2023.

Monty-Maximilian Zühlke and Daniel Kudenko. Adversarial Robustness of Neural Networks From the Perspective of Lipschitz Calculus: A Survey. *ACM Comput. Surv.*, 2024.

# A Appendix – Proofs

## A.1 Proof of Lemma 1

*Proof.* Applying the Fourier transform on (15) yields, for any $\xi \in \mathbb{R}^2$,

$$ \widehat{F}_0(\xi) = (\widehat{F} * \widehat{\Psi})(\xi - \nu) = \mathcal{T}_{\nu}(\widehat{F} * \widehat{\Psi})(\xi). \qquad (141) $$

By hypothesis on $\Psi$, we have

$$ \mathrm{supp}(\widehat{F} * \widehat{\Psi}) \subset \mathrm{supp} \widehat{\Psi} \subset B_{\infty}(-\nu, \varepsilon/2). \qquad (142) $$

The translation operator $\mathcal{T}_{\nu}$ shifts the support with respect to $\nu$, which yields (16). □

## A.2 Proof of Proposition 1

*Proof.* Using the 2D Plancherel formula, we compute

$$ \begin{aligned} \| \mathcal{T}_h F_0 - F_0 \|_{L^2}^2 &= \frac{1}{4\pi^2} \left\| \widehat{\mathcal{T}_h F_0} - \widehat{F_0} \right\|_{L^2}^2 \\ &= \frac{1}{4\pi^2} \iint_{B_{\infty}(\varepsilon/2)} |\widehat{F_0}(\xi)|^2 |e^{-i\langle h, \xi \rangle} - 1|^2 d^2\xi \\ &= \frac{1}{4\pi^2} \iint_{B_{\infty}(\varepsilon/2)} |\widehat{F_0}(\xi)|^2 (2 - 2 \cos \langle h, \xi \rangle) d^2\xi \\ &\le \frac{1}{4\pi^2} \iint_{B_{\infty}(\varepsilon/2)} |\widehat{F_0}(\xi)|^2 |\langle h, \xi \rangle|^2 d^2\xi, \end{aligned} $$

because $\mathrm{cost} \ge 1 - \frac{t^2}{2}$. Note that the integral is computed on a compact domain because, according to Lemma 1, $\mathrm{supp} \widehat{F_0} \subset B_{\infty}(\varepsilon/2)$. Next, we use the Cauchy-Schwarz inequality to compute:

$$ \begin{aligned} \forall \xi \in B_{\infty}(\varepsilon/2), & \quad |\langle h, \xi \rangle| \le \|h\|_1 \cdot \|\xi\|_{\infty} \\ & \le \frac{\varepsilon}{2} \|h\|_1. \end{aligned} $$

Therefore,

$$ \| \mathcal{T}_h F_0 - F_0 \|_{L^2}^2 \le \frac{\varepsilon}{4} \|h\|_1^2 \|F_0\|_{L^2}^2, \qquad (143) $$

which yields the result. □
According to (3.92), the conditional expected value of $\tilde{Q}_X^2$ remains the same whatever the outcome of $\tilde{S}_X^2$. Thus, the law of total expectation states that

$$ (3.93) \qquad \mathbb{E}[\tilde{Q}_X^2] = \mathbb{E}[\mathbb{E}[\tilde{Q}_X^2 | \tilde{S}_X^2]] = \gamma_q(m\theta)^2. $$

Since we have assumed Hypothesis 3.4, we can apply Proposition 3.5. Using the definition of $\tilde{P}_X$ (3.58) and $\tilde{Q}_X$ (3.61), we get

$$ (3.94) \qquad \tilde{P}_X \le \tilde{Q}_X + \beta_q(m\kappa). $$

Then,

$$ (3.95) \qquad \mathbb{E}[\tilde{P}_X^2] \le \mathbb{E}[\tilde{Q}_X^2] + 2\beta_q(m\kappa)\mathbb{E}[\tilde{Q}_X] + \beta_q(m\kappa)^2. $$

According to Jensen's inequality,

$$ (3.96) \qquad \mathbb{E}[\tilde{Q}_X] \le \sqrt{\mathbb{E}[\tilde{Q}_X^2]} = \gamma_q(m\theta). $$

Thus,

$$ (3.97) \qquad \mathbb{E}[\tilde{P}_X^2] \le \gamma_q(m\theta)^2 + 2\beta_q(m\kappa)\gamma_q(m\theta) + \beta_q(m\kappa)^2, $$

which yields (3.87).

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

Let us analyze the bound obtained in (3.87). The first term, $\beta_q(m\kappa)$, accounts for the localized property of the convolution filter $W$. This term decreases linearly with the product $m\kappa$. In the limit case where $\kappa = 0$ (infinite, nonlocal filter), we get $\beta_q(m\kappa) = 0$. Note that a smaller subsampling factor $m$ allows for a larger bandwidth $\kappa$. Besides, $\beta_q(m\kappa)$ increases linearly with the size of the max pooling grid, which is characterized by $q$. The second term, $\gamma_q(m\theta)$, accounts for the discrete nature of the max pooling grid. It strongly depends on the characteristic frequency $\theta$, as illustrated in Figure 5. According to (3.82), this term has a polynomial decay when $q$ increases. However, increasing the size of the max pooling grid also results in increasing the term $\beta_q(m\kappa)$, as explained above. Therefore, a tradeoff must be found to get an optimal bound.

**4. Shift Invariance of $\mathbb{R}Max$ Outputs.** In this section, we present the main theoretical claim of this paper. Based on the previous results, we provide a probabilistic measure of shift invariance for $\mathbb{R}Max$ operators. First, we consider the following lemma.

**Lemma 4.1.** If Hypotheses 3.7 and 3.8 are satisfied, then they are also true with $X \leftarrow \mathcal{T}_u X$, for any $u \in \mathbb{R}^2$.

*Proof.* First, we show that, for any $x \in \mathbb{R}^2$,

$$ (4.1) \qquad \mathcal{T}_{uX}(x) = \mathcal{T}_{su}\mathcal{T}_X(x); $$

$$ (4.2) \qquad \mathcal{Z}_{uX}(x) = \mathcal{T}_{su}\mathcal{Z}_X(x). $$
### A.3 Proof of Lemma 2

*Proof.* Since $F \in \mathcal{V}^{(s)}$, the two-dimensional version of Shannon's sampling theorem (Mallat, 2009, Theorem 3.11, p. 81) yields

$$F = \sum_{n \in \mathbb{Z}^2} X[n] \Phi_n^{(s)}, \quad \text{and} \quad \hat{F} = \sum_{n \in \mathbb{Z}^2} X[n] \widehat{\Phi}_n^{(s)}. \qquad (144)$$

Moreover, using (19), we can show that, for any $\xi \in B_{\infty}(\pi/s)$,

$$\widehat{\Phi}_n^{(s)}(\xi) = \widehat{\Phi}^{(s)}(\xi) e^{-i\langle s\xi, n \rangle} = s e^{-i\langle s\xi, n \rangle}. \qquad (145)$$

Therefore, plugging (145) into (144) proves (22).

Then, by combining (22) with the Plancherel formula, we get

$$\begin{aligned} \|F\|_{L^2}^2 &= \frac{1}{4\pi^2} \|\hat{F}\|_{L^2}^2 \\ &= \frac{1}{4\pi^2} \iint_{B_{\infty}(\pi/s)} |\hat{F}(\xi)|^2 d^2\xi \\ &= \frac{1}{4\pi^2} \iint_{B_{\infty}(\pi/s)} |s \hat{X}(s\xi)|^2 d^2\xi. \end{aligned}$$

The integral is performed on $B_{\infty}(\pi/s)$ because $F \in \mathcal{V}^{(s)}$. Then, by applying the change of variable $\xi' \leftarrow s\xi$, we get

$$\begin{aligned} \|F\|_{L^2}^2 &= \frac{1}{4\pi^2} \iint_{B_{\infty}(\pi)} |\hat{X}(\xi')|^2 d^2\xi' \\ &= \frac{1}{4\pi^2} \|\hat{X}\|_{L^2}^2 = \|X\|_2^2, \end{aligned}$$

hence (23), which concludes the proof. □

### A.4 Proof of Proposition 2

*Proof.* First, $F_X$ and $\Psi_W$ are well defined because $X \in l_2^2(\mathbb{Z}^2)$ and $W \in l_2^2(\mathbb{Z}^2)$. By construction, $F_X$ and $\Psi_W \in \mathcal{V}^{(s)}$. Therefore, according to Shannon's sampling theorem (Mallat, 2009, Theorem 3.11, p. 81),

$$F_X := s \sum_{n \in \mathbb{Z}^2} F_X(sn) \Phi_n^{(s)} \quad \text{and} \quad \Psi_W := s \sum_{n \in \mathbb{Z}^2} \Psi_W(sn) \Phi_n^{(s)}. \qquad (146)$$

By uniqueness of decompositions in an orthonormal basis, we get (26). Moreover, using (22) in Lemma 2, we get, for any $\xi \in B_{\infty}(\pi/s)$,

$$\hat{\Psi}_W(\xi) = s \hat{W}(s\xi). \qquad (147)$$

Since $\hat{\Psi}_W(\xi) = 0$ outside $B_{\infty}(\pi/s)$, (147) is true for any $\xi \in \mathbb{R}^2$. Therefore, by hypothesis on W,

$$\mathrm{supp} \hat{\Psi}_W \subset B_{\infty}(\theta/s, \kappa/(2s)), \qquad (148)$$

which yields (25).

We now prove (27). For $n \in \mathbb{Z}^2$, we compute:

$$\begin{aligned} (F_X * \bar{\Psi}_W)(msn) &= \iint_{\mathbb{R}^2} F_X(msn - x) \bar{\Psi}_W(x) d^2x \\ &= \iint_{\mathbb{R}^2} \sum_{p \in \mathbb{Z}^2} X[p] \Phi_p^{(s)}(msn - x) \bar{\Psi}_W(x) d^2x \\ &= \sum_{p \in \mathbb{Z}^2} X[p] \iint_{\mathbb{R}^2} \Phi_p^{(s)}(msn - x) \bar{\Psi}_W(x) d^2x. \end{aligned}$$