\section{Introduction}
Recent research has shown that model performances can drop drastically when natural test images are altered \cite{szegedy2013intriguing,hendrycks2021natural,hendrycks2021many}.
One of these scenarios is when common image corruptions are added to the images. These corruptions include noise attributed to weather conditions, noisy environments, blurring effects and digital artifacts \cite{hendrycks2019benchmarking}. Another case where models can fail is adversarial examples where a malicious party can craft imperceptible perturbations to images to influence the model's prediction \cite{szegedy2013intriguing,carlini2017towards,papernot2018cleverhans,croce2019minimally}. While research to build models more robust to these two situations started independently, recent studies have started to draw a connection between them. Interestingly, models trained to be robust against adversarial examples have shown mixed results in common corruptions: improving accuracy for some corruptions types while doing poorly for others.

Though studies have shown a link between data augmentation strategies and robustness against corruptions with different frequency components \cite{yin2019fourier}, there is no study on how \emph{direct} changes to a model's Fourier profile would affect its robustness. Here, we aim to directly alter the Fourier profile of a model to study its direct effect on both the model's adversarial and corruption robustness. To achieve this, we investigate the model's Jacobian which represents a visual map of pixel importance in a particular input image \cite{smilkov2017smoothgrad}. Intuitively, a model would be relying on low (high)-frequency features when its Jacobians have a large ratio of low (high)-frequency components. While observing the Fourier spectra of natural images such as SVHN, CIFAR-10 and CIFAR-100, and the Jacobians of standard-trained models (Figure~\ref{fig:Fourier Spectra imgs and models}), these images have a much larger component of low-frequency features than their standard-trained models. This mismatch of frequency profiles motivates us to train models to bias towards low-frequency features through its Jacobians and study its effect on robustness.

To quantify the frequency profile of a model, we propose a frequency bias term that computes a scalar value from a Fourier spectrum of 2-D inputs such as an image or Jacobian to improve frequency evaluation beyond the visual inspection of Fourier spectra. Through this differentiable frequency bias term, we can use Jacobian frequency regularization (JaFR) to explicitly train a model to bias more heavily on low- or high-frequency features. Through our experiments, we find that a more direct change in a model's frequency profile towards low-frequency regions to match the frequency profile of the training data can boost clean accuracy, adversarial robustness and common corruptions in certain settings while trading off performance against low-frequency noise. Conversely, regularizing the model towards high-frequency regions can boost performance against low-frequency noise while sacrificing accuracy under high-frequency noise and adversarial examples. Our main goal here is not to claim state-of-the-art robustness but to elucidate a more direct link between models' frequency characteristics and robustness against noise.
All in all, the core contributions of this paper are:
\begin{itemize}
    \item We propose a frequency bias term to measure the frequency bias of a Fourier spectrum.
    \item We show that biasing the Jacobians of models towards low or high frequency have implications on model robustness against adversarial robustness and an array of corruptions.
    \item To achieve this, we propose Jacobian frequency regularization (JaFR) to train model's Jacobians to have a larger or smaller weightage of low-frequency components.
    \item We conduct experiments on SVHN, CIFAR-10, CIFAR-100 and TinyImageNet to show how low-frequency bias in Jacobians can improve robustness against adversarial and high-frequency corruptions, albeit with tradeoffs in performance for low-frequency corruptions.
\end{itemize}

\begin{figure}[!htbp]
    \centering
    \includegraphics[width=0.65\linewidth]{images/Fourier_Spectra_imgs_and_models.pdf}
    \caption{Fourier spectra of image datasets and the Jacobians of models trained on them, showing a mismatch between the frequency profiles of the images and models that are trained on them. Arrows show the direction of increasing frequency. Clean images contain large ratios of low-frequency components, shown by the intensity on the top left of the spectra (top row) while their corresponding standard-trained models rely heavily on high-frequency features (bottom row).}
    \label{fig:Fourier Spectra imgs and models}
\end{figure}

\section{Background and Related Work}

\textbf{Adversarial Robustness} \label{sec:adversarial robustness}
Adversarial robustness measures how well a model is resistant to attacks by malicious actors. In such attacks, imperceptible perturbations could be crafted to form adversarial examples with the aim to control the prediction of neural networks \cite{szegedy2013intriguing}. This threat could undermine the deployment of deep learning models in mission-critical applications.
In a classification task, a model ($f$) parameterized by $\theta$ takes an input $\xinput$ to predict the probabilities for $k$ classes, i.e., $f(\xinput ; \theta): \xinput \mapsto \mathbb{R}^k$. In supervised setting of empirical risk minimization (ERM), given training samples $(\xinput, \mathbf{y}) \sim D$, the model's parameters are trained to minmize the standard cross-entropy loss:
\begin{equation}
 \mathcal{L}(\xinput, \mathbf{y}) = \mathbb{E}_{(\xinput, \mathbf{y}) \sim D} \left[ - \mathbf{y}^\top \log f (\xinput) \right]
\end{equation}
where $\mathbf{y} \in \mathbb{R}^k$ is the one-hot label for the input $\xinput$.
Though ERM can train models that display high accuracy on a holdout test set, their performance degrades under adversarial examples. Given an adversarial perturbation of magnitude $\varepsilon$, we say a model is robust against this attack if

\begin{multline}
    \argmax_{i \in C} f_i(\xinput ; \theta) = \argmax_{i \in C} f_i(\xinput + \delta ; \theta) ~, \\~~~ \forall \delta \in B_p (\varepsilon) = {\delta : \| \delta \|_p \leq \varepsilon}
\end{multline}

where $p = \infty$ in this paper and is the most widely studied scenario.

One of the most effective approaches to train models robust against adversarial examples is adversarial training (AT) \cite{goodfellow2014explaining,madry2017towards,zhang2019defense,qin2019adversarial,andriushchenko2020understanding,wu2020adversarial}. More details of related work in adversarial training are deferred to the Appendix.

\if 0
\textbf{Adversarial Training} One of the most effective approaches to train models robust against adversarial examples is adversarial training (AT) \cite{goodfellow2014explaining}. The intuition behind AT is to match the training distribution with the adversarial example distribution. This is achieved by crafting adversarial examples in each training iteration and using them as training samples to minimize the following loss:

\begin{equation} \label{eq:AT objective}
\mathcal{L}(\xinput, \mathbf{y}) = \mathbb{E}_{(\xinput, \mathbf{y}) \sim D} \left[ \max_{\delta \in B(\varepsilon)} \mathcal{L}(\xinput + \delta, \mathbf{y}) \right]
\end{equation}
where the inner maximization, $\max_{\delta \in B(\varepsilon)} \mathcal{L}(\xinput + \delta, \mathbf{y})$, is performed with gradient descent. The earliest form of AT, fast gradient sign method (FGSM) \cite{goodfellow2014explaining}, uses one single gradient to create adversarial training samples but models trained this way have shown to be vulnerable to stronger adversarial examples that are crafted with more gradient steps. The more recent projected gradient descent (PGD) adversarial training \cite{madry2017towards} and its variants \cite{zhang2019defense,qin2019adversarial,andriushchenko2020understanding,wu2020adversarial} are one of the stronger defenses which performs the following gradient step iteratively:
\begin{equation} \label{eq:pgd step}
    \delta \gets \mathrm{Proj} \left[ \delta - \eta ~ \text{sign} \left( \nabla_{\delta} \mathcal{L}(\xinput + \delta, \mathbf{y}) \right) \right ]
\end{equation}
where $\mathrm{Proj}(\xinput) = \argmin_{\zeta \in B(\varepsilon)} \| \xinput - \zeta \|$. Different from the AT-based techniques, JaFR does not rely on adversarial examples for training and is orthogonal to these AT-based techniques.
\fi

\textbf{Corruption Robustness} \label{sec:corruption robustness}
In contrast to adversarial robustness, corruption robustness \cite{hendrycks2019benchmarking} entails studying the performance of models when input data are corrupted with common-occurring noise, not necessarily those created by a malicious actor to control models' prediction. CIFAR-10-C and CIFAR-100-C are two benchmark datasets that study 19 corruption types, each with 5 levels of severity. These 19 corruption types can be grouped under the `Noise', `Blur', `Weather', `Digital' categories. There is a line of work that seek to improve performances under these common corruptions by mostly improving the training data augmentation \cite{geirhos2018imagenet,cubuk2018autoaugment,hendrycks2019augmix,lopes2019improving,hendrycks2020many,kapoor2020fourier,rusak2020simple,vasconcelos2020effective}. Some other works assemble expert models whose performance are finetuned for subsets of the corruptions \cite{lee2020compounding,saikia2021improving} to boost overall performance. Different from these work, our paper aims to study the effect of frequency bias on corruption robustness of different types that corrupts features of varying frequencies, rather than to propose a new way to better resist these corruptions.

\textbf{Link between Frequency and Robustness} \label{sec:fourier perspectives on robustness}
There is a line of work that seeks to understand how robust models respond to corruption and adversarial perturbations of various frequency profiles \cite{yin2019fourier,sharma2019effectiveness,ortiz2020hold,tsuzuku2019structural,vasconcelos2021impact}. Details of these works are in the Appendix.
In contrast to these prior works, our work here takes the frequency analysis in a different direction by studying models through the Fourier spectrum of their Jacobians rather than their test or training data. More concretely, with the original training data, we train models and bias the frequency profile of the model's Jacobians towards low-frequency regions to see its effect on model robustness.

\if 0
\textbf{Link between Frequency and Robustness} \label{sec:fourier perspectives on robustness}
Studies have shown that data augmentation techniques have improved robustness against certain corruptions but degrade performances on others \cite{yin2019fourier}. \cite{yin2019fourier} carried out frequency analysis which shows a connection between the Fourier spectrum of a particular corruption type and whether its performance improves or degrades. More specifically, they found that Gaussian data augmentation and adversarial training result in models that rely heavily on low-frequency features in images. These models are hence more resistant against high-frequency noise but more susceptible to low-frequency corruptions. The authors analyzed the Fourier profile of trained models by measuring performance when noise is created by components bounded within a region of the Fourier spectrum or when training data have certain frequency components removed. \cite{yin2019fourier} use a qualitative inspection of the Fourier spectrum to show which corruption types have a relative ratio of low or high-frequency components while we propose a frequency bias term to quantify it. In other works that explored the link between frequency and robustness, \cite{sharma2019effectiveness} showed that models trained against adversarial examples are equally vulnerable as standard-trained models when faced with adversarial perturbations that are constrained to the low-frequency regions, confirming the findings from \cite{yin2019fourier} that adversarially trained models' robustness is based on invariance to high-frequency signals. \cite{ortiz2020hold} showed that trained models can learn to be invariant to low or high-frequency features depending on how discriminative these features are to the classification label (i.e., how much they change the label). \cite{tsuzuku2019structural} found that universal adversarial perturbations that fool a range of CNN classifiers are a combination of a few Fourier basis functions.
Our work here takes the Fourier analysis in a different direction by studying models through the Fourier spectrum of their Jacobians rather than their test or training data. More concretely, with the original training data, we train models and bias the frequency profile of the model's Jacobians towards low-frequency regions to see its effect on model robustness.
\fi

\textbf{Jacobians of Robust Models} \label{sec:jacobians of robust models}
The Jacobian,
\begin{equation}
    J \coloneqq \nabla_{\xinput} \mathcal{L(\xinput, \mathbf{y})}
\end{equation}
defines how the model's prediction changes with an infinitesimally small change to the input $\xinput$. For image classification, Jacobians can be loosely interpreted as a map of which pixels affect the model's prediction the most and, hence, give an illustration of important regions in an input image \cite{smilkov2017smoothgrad,adebayo2018sanity,etmann2019connection,ilyas2019adversarial}. There is a line of work that seeks to improve the adversarial robustness of models by matching the Jacobians to a target distribution \cite{chan2019jacobian,chan2020thinks} or by constraining their magnitude \cite{ross2018improving,jakubovitz2018improving}. Rather than aiming to improve the adversarial robustness, the core aim of our paper here is to investigate the relationship between the Fourier profile of models and robustness against corruptions. Moreover, the regularizing effect of JaFR has a different mechanism that acts directly on the Fourier spectrum of the Jacobians. A detailed comparison is in the Appendix.

\if 0
\textbf{Jacobians of Robust Models} \label{sec:jacobians of robust models}
The Jacobian,
\begin{equation}
    J \coloneqq \nabla_{\xinput} \mathcal{L(\xinput, \mathbf{y})}
\end{equation}
defines how the model's prediction changes with an infinitesimally small change to the input $\xinput$. For image classification, Jacobians can be loosely interpreted as a map of which pixels affect the model's prediction the most and, hence, give an illustration of important regions in an input image \cite{smilkov2017smoothgrad,adebayo2018sanity,ilyas2019adversarial}. Previous studies have observed that adversarially robustness models trained with AT display Jacobians that are more salient than those from non-robust standard-trained models \cite{tsipras2018robustness}. Since then, there a line of work that studies the theoretical link between the saliency of the Jacobians and robustness \cite{etmann2019connection} and exploits this link to improve robustness by regularizing for Jacobians' saliency \cite{chan2019jacobian,chan2020thinks}. By using generative adversarial networks (GANs) to train a model's Jacobians to fit the distribution of either the input images \cite{chan2019jacobian} or of a robust teacher model \cite{chan2020thinks}, adversarial robustness can be improved. In contrast, the regularizing effect of JaFR acts directly on the Fourier spectrum of the Jacobians rather than fitting them to a target distribution through GANs. Furthermore, our work here also studies the effect of JaFR on common corruptions, on top of adversarial robustness.

Other work that regularizes the input gradients to boost adversarial robustness includes using double backpropagation \cite{drucker1991double} to minimize the input gradients' Frobenius norm \cite{ross2018improving,jakubovitz2018improving}. Those approaches aim to constrain the effect that changes at individual pixels have on the classifier's output but not the frequency profile of neural networks like our method. Rather than aiming to improve the adversarial robustness, the core aim of our paper here is to investigate the relationship between the Fourier profile of models and robustness against corruptions.
\fi

\section{Jacobian Frequency Bias}
\textbf{Motivation}
As mentioned in \S~\ref{sec:fourier perspectives on robustness}, it has been shown that there is a link between Fourier profile of input training data and robustness against corruptions with different frequency components. However, there is still no study on how changes to a neural network's Fourier profile would affect its robustness. Since the Jacobian of a model represents a visual map of pixel importance  \cite{smilkov2017smoothgrad}, it offers a medium for us to regularize the Fourier profile of the model. Intuitively, when a model's Jacobians are concentrated with low-frequency components, it places more importance on low-frequency features. Conversely, the model relies more on high-frequency features if its Jacobians have a relatively large proportion of high-frequency components. Here, we aim to study how changing the Fourier spectrum of the Jacobians would affect its robustness.

When analyzing the Fourier profile of Jacobians of adversarially robust models (Table~\ref{tab:cifar10 model fourier profile}), we see that it resembles the profile of the training data much more than the non-robust standard trained models. This raises the question of what would happen to model robustness if we directly train neural networks to have a low-frequency profile, similar to what we see in the images from SVHN, CIFAR-10 and CIFAR-100 (see Figure~\ref{fig:Fourier Spectra imgs and models}). This motivates a metric to quantitatively measure and control the Fourier profile of the neural network to more directly study the effect of Fourier profile on robustness. In the next sections, we propose the Jacobian frequency bias to achieve this goal and discuss how we can use it to train a neural network.

\subsection{Jacobian Frequency Bias}
Here, we present the measure of Jacobian frequency bias with an example of single-channel image classification where the input images are denoted as $\xinput \in \mathbb{R}^{hw}$. Given a training dataset ($\mathcal{D}_{\text{train}}$) where each training sample consists of an input image $\xinput$ and one-hot label vector of $k$ classes as $\mathbf{y} \in \mathbb{R}^k$, we can express $f_{\text{cls}} (\xinput) \in \mathbb{R}^k$ as the prediction of the classifier ($f_{\text{cls}}$), parameterized by $\theta$. Then, the classification cross entropy loss ($\mathcal{L}_{\text{cls}}$) is:
\begin{equation} \label{eq:cls xent}
    \mathcal{L}_{\text{cls}} = - \mathbf{y}^\top \log f_{\text{cls}} (\xinput)
\end{equation}

where $f_{\text{cls}}$ is the classifier model. The Jacobian matrix $J \in \mathbb{R}^{hw}$ of the model's classication loss value with respect to the input layer can be computed through backpropagation:
\begin{equation}
    J(\xinput) \coloneqq \nabla_{\xinput} \mathcal{L}_{\text{cls}} = \begin{bmatrix}
  \frac{\partial \mathcal{L}_{\text{cls}}}{\partial {\xinput_{1, 1}}} & \cdots & \frac{\partial \mathcal{L}_{\text{cls}}}{\partial {\xinput_{w,1}}} \\
  \vdots  & \ddots & \vdots  \\
  \frac{\partial \mathcal{L}_{\text{cls}}}{\partial {\xinput_{1,h}}} & \cdots & \frac{\partial \mathcal{L}_{\text{cls}}}{\partial {\xinput_{w,h}}}
 \end{bmatrix}
\end{equation}

We can then compute the Jacobian's Fourier spectrum to retrieve a frequency profile of the model. Since the input images are made up of discrete pixels, we can extract this information by applying a discrete Fourier transform ($\mathcal{F}$) to the Jacobian to get a map ($M$) of its frequency components' magnitude:
\begin{equation}
    M_{i,j} = | \mathcal{F}(J)[i,j] |
\end{equation}
In our experiments where the input images have 3 RGB channels, we compute the Fourier map for each other channel separately ($M_r$, $M_g$, $M_b$) and take the mean across these channels, i.e., $M_{i,j} = \frac{M_r[i,j] + M_g[i,j] + M_b[i,j]}{3} $.

Next, we propose to compute a scalar bias term ($\mathcal{B}_{\text{low}}$) from the 2-D map ($M$) to measure the relative bias (or ratio) of low-frequency components with respect to high-frequency components. One criterion for $\mathcal{B}_{\text{low}}$ is that the contribution of the frequency magnitude to this term should monotonically decrease as the frequency increase, i.e., larger high-frequency magnitudes results in lower $\mathcal{B}_{\text{low}}$ values. To satisfy this, we monotonically decrease the exponent value on the frequency magnitude as its frequency increases. For a 1-D scenario where $l$ is the dimension of the Fourier spectrum, we can express this bias term $\mathcal{B}_{\text{low}}$ as:
\begin{equation}
     \mathcal{B}_{\text{low}} = \Pi_i (M_i)^{\alpha_i},~~~~~ \alpha_i < \alpha_j,~~~~~ \forall i, j \in [1, l], ~~~~~  i < j
\end{equation}

where $M_1$ and $M_l$ are the magnitudes of the lowest and highest frequency components respectively.
To ensure that $\mathcal{B}_{\text{low}}$ measures the relative ratio between the frequencies rather than absolute values of the frequency components, we use the following constraint on the $\alpha$ so that it is independent of the sum of the components' magnitudes:
\begin{equation}
     \alpha_i = - \alpha_{(l-i+1)}, ~~~~~ \forall i \in [1, l]
\end{equation}

In all our experiments, we use an array of values whose values are evenly spaced with distance $k$ from $\alpha_1$ and $\alpha_l$ for the $\alpha$ values, i.e,
\begin{equation}
     \alpha_{i+1} = \alpha_i + k, ~~ \forall i \in [1, l-1]
\end{equation}
and use $\alpha_1 = 1, \alpha_l = -1$.
When generalizing the bias term to two axes, we can compute the sum of all bias terms along each row and column of the 2-D Fourier spectrum to give:
\begin{equation}
    \mathcal{B}_{\text{low}} = \left[ \sum_j \left[ \Pi_i (M_{i,j})^{\alpha_i} \right] +  \sum_i \left[ \Pi_j (M_{i,j})^{\alpha_j}\right] \right]
\end{equation}
In the next section, we discuss how $\mathcal{B}_{\text{low}}$ can be used to regularize neural networks' Jacobians to bias them towards low-frequency components.

\begin{figure}[!htbp]
    \centering
    \includegraphics[width=0.8\linewidth]{images/objectives_figure.pdf}
    \caption{Training architecture of Jacobian frequency regularization (JaFR). JaFR trains the model's Jacobian to bias towards low-frequency components by shifting its Fourier spectrum's intensity towards the low-frequency regions, leftwards for the horizontal and upwards for the vertical axis.}
    \label{fig:JaFR}
\end{figure}

\subsection{Jacobian Frequency Regularization (JaFR)} \label{sec:JaFR training}
To recall, the aim here is to control the Fourier spectrum of a model's Jacobians to study how its Fourier profile can affect model robustness. To achieve this, we propose Jacobian Frequency Regularization (JaFR) to bias the Fourier spectrum of a model towards the low-frequency region. Figure~\ref{fig:JaFR} shows a summary of how our proposed Jacobian Frequency Regularization (JaFR) trains a model to alter its frequency profile.
Since the operations involved in computing $\mathcal{B}_{\text{low}}$ are differentiable, we can simply incorporate it in the following loss function with the aim to maximize the low-frequency bias of a model:

\begin{multline}
    \mathcal{L}_{\text{freq}}  = - \log \mathcal{B}_{\text{low}} = \\ - \log \left[ \sum_j \left[ \Pi_i (M_{i,j})^{\alpha} \right] +  \sum_i \left[ \Pi_j (M_{i,j})^{\beta}\right] \right]
\end{multline}

Combining with the classification loss in Equation~\ref{eq:cls xent},  we can optimize through stochastic gradient descent to approximate the optimal parameters for the classifier $f_{\text{cls}}$ as follows,
\begin{equation}
\theta^* = \argmin_{\theta} (\mathcal{L}_{cls} + \lambda_{freq} \mathcal{L}_{freq})
\end{equation}
where $\lambda_{freq}$ determines the weight of JaFR in the model's training. A positive $\lambda_{freq}$ value regularizes the model to bias towards low-frequency features while a negative value (indicated as JaFR(-)) conversely steers the model towards high-frequency features. The $\mathcal{L}_{freq}$ term can be computed during each standard training iteration with an additional backpropagation step. Algorithm~\ref{algo:JaFR Training} summarizes the training of a classifier with JaFR.
In the next sections, we present the results on the effect of JaFR on model robustness under adversarial examples and corruption noises.

\begin{algorithm}
\footnotesize
 \caption{Jacobian Frequency Regularization Training}
 \label{algo:JaFR Training}

\textbf{Input:} Train data $\mathcal{D}_{\text{train}}$, learning rate $\gamma$
\begin{algorithmic}[1]
\FOR{ each training iteration }
\STATE Sample $(\xinput, \mathbf{y}) \sim \mathcal{D}_{\text{train}}$
\STATE $\mathcal{L}_{\text{cls}} \gets - \mathbf{y}^\top \log f_{\text{cls}} (\xinput)$ \hfill\textit{(1) Compute classification cross-entropy loss}
\STATE $J \gets \nabla_{\xinput} \mathcal{L}_{\text{cls}}$ \hfill\textit{(2) Compute Jacobian matrix}
\STATE $M_{i,j} \gets | \mathcal{F}(J)[i,j] | $ \hfill\textit{(3) Compute frequency magnitudes}
\STATE $\mathcal{L}_{\text{freq}} \gets - \log \left[ \sum_j \left[ \Pi_i (M_{i,j})^{\alpha} \right] +  \sum_i \left[ \Pi_j (M_{i,j})^{\beta}\right] \right] $ \hfill\textit{(4) Compute frequency bias}
\STATE $\theta \gets \theta - \gamma ~ \nabla_{\theta} (\mathcal{L}_{\text{cls}} + \lambda_{\text{freq}} \mathcal{L}_{\text{freq}})  $ \hfill\textit{(5) Update the classifier $f_{\text{cls}}$ to minimize $\mathcal{L}_{cls}$ and $\mathcal{L}_{\text{freq}}$}
\ENDFOR
\end{algorithmic}
\end{algorithm}

\begin{table*}[!htbp]
    \centering
    \footnotesize
        \begin{tabular}{ l|cc|cccc|c }
         \hline
         ~ & Standard & JaFR & FGSM AT & FGSM AT & FGSM AT & PGD AT & $\mathcal{B}_{\text{low}}$ \\
         ~ & ~ & ~ & + JaFR(-) & ~ & + JaFR & ~ & ~ \\
        \hline
        Clean & 93.11$\pm$0.20 & 93.13$\pm$0.11 & 83.16 $\pm$1.01 & 84.80$\pm$1.37 & 79.94$\pm$0.22 & 79.31$\pm$0.23 & - \\
        \hline
        mCE & 100.00 & 95.84 & 127.51 & 124.32 & 134.37 & 135.71 & - \\
        \hline
        Fog & 86.23$\pm$0.45 & \textbf{86.60$\pm$0.28}  &   \textbf{67.15$\pm$3.81}   & 63.00$\pm$4.03  & 56.60$\pm$0.72 & 56.36$\pm$1.25 & 12.85\\
        Saturate & \textbf{89.27$\pm$0.16}  & 89.07$\pm$0.28 &  77.55$\pm$0.76   & \textbf{79.54$\pm$1.07} & 76.61$\pm$0.22 & 76.29$\pm$0.21 & 12.28\\
        Contrast & \textbf{75.10$\pm$1.08} & 73.45$\pm$0.95 &   \textbf{50.93$\pm$3.10}   & 44.06$\pm$3.97 & 41.96$\pm$0.61 & 41.19$\pm$0.88 & 12.14\\
        Bright & \textbf{91.56$\pm$0.15} & 91.45$\pm$0.18 &  81.02$\pm$0.88  & \textbf{83.04$\pm$1.20} & 76.95$\pm$0.42 & 76.01$\pm$0.58 & 12.01\\
        Snow & 79.54$\pm$0.71  & \textbf{81.06$\pm$0.15} &  75.05$\pm$1.52   & \textbf{77.83$\pm$1.43} & 74.28$\pm$0.27 & 73.73$\pm$0.21 & 11.53\\
        Frost & 75.41$\pm$0.63 & \textbf{78.97$\pm$0.56} &   75.82$\pm$1.64  & \textbf{77.70$\pm$1.44} & 70.50$\pm$0.54 & 69.06$\pm$1.18 & 10.61\\
        Motion & \textbf{75.67$\pm$1.35} & 75.13$\pm$1.21 &  67.33$\pm$2.01  & 65.15$\pm$3.47 & \textbf{72.65$\pm$0.41}  & 72.15$\pm$0.78 & 10.61\\
        Zoom & 76.42$\pm$1.65 & \textbf{76.48$\pm$1.20}  &   71.39$\pm$1.87  & 69.37$\pm$3.91 & \textbf{75.38$\pm$0.44} & 74.78$\pm$0.67 & 10.37\\
        Elastic & 81.62$\pm$0.60 & \textbf{81.98$\pm$0.55} &  73.88$\pm$1.45   & 73.90$\pm$2.27 & \textbf{74.86$\pm$0.37} & 74.25$\pm$0.58 & 10.05\\
        Pixel & \textbf{72.90$\pm$0.70} & 72.45$\pm$0.89 &   79.78$\pm$0.72   & \textbf{81.20$\pm$0.55} & 78.15$\pm$0.14 & 77.47$\pm$0.33 & 8.2\\
        Gauss. B & 73.38$\pm$1.58 & \textbf{73.53$\pm$1.16}  &  70.56$\pm$1.51   & 68.49$\pm$3.93 & \textbf{74.39$\pm$0.38} & 73.71$\pm$0.63 & 7.82\\
        Defocus & 81.57$\pm$0.83 & \textbf{81.76$\pm$0.71} &  74.22$\pm$1.37  & 73.42$\pm$3.01 & \textbf{76.15$\pm$0.29} & 75.46$\pm$0.54 & 7.49\\
        Glass & 50.18$\pm$1.75 & \textbf{55.12$\pm$0.58} &  66.28$\pm$2.89   & 70.18$\pm$1.07 & \textbf{73.75$\pm$0.28} & 73.33$\pm$0.49 & 7.01\\
        Spatter & 80.77$\pm$0.11 & \textbf{82.29$\pm$0.45} &  76.46$\pm$1.68  & \textbf{78.49$\pm$2.04} & 75.83$\pm$0.32 & 75.29$\pm$0.32 & 6.77\\
        JPEG & 77.76$\pm$0.56 & \textbf{79.56$\pm$0.60} &  80.38$\pm$0.77   & \textbf{81.47$\pm$1.54} & 78.19$\pm$0.18 & 77.58$\pm$0.24 & 6.51\\
        Speckle & 62.79$\pm$2.44 & \textbf{66.74$\pm$0.51} &  68.86$\pm$4.83  & 74.07$\pm$4.08 & \textbf{77.36$\pm$0.21} & 76.22$\pm$0.54 & 3.76\\
        Shot & 59.63$\pm$2.75 & \textbf{64.32$\pm$0.76} &  68.32$\pm$5.06   & 73.73$\pm$4.19 & \textbf{77.35$\pm$0.24} & 76.19$\pm$0.42 & 3.68\\
        Impulse  & 56.50$\pm$1.64 & \textbf{58.75$\pm$0.59} &   57.73$\pm$6.09  & 64.51$\pm$8.91 & 71.55$\pm$0.54 & \textbf{71.56$\pm$0.74} & 3.63\\
        Gauss. N & 48.70$\pm$3.31 & \textbf{54.34$\pm$1.00} &  62.91$\pm$6.74   & 69.40$\pm$5.27 & \textbf{76.59$\pm$0.33} & 75.15$\pm$0.49 & 3.62\\
        \hline
        \end{tabular}
    \caption{Accuracy values ($\uparrow$ better) and mCE ($\downarrow$ better) for different models under CIFAR-10 corruptions. The corruption types are arranged with descending order of low-frequency bias $\mathcal{B}_{\text{low}}$.}
\label{tab:cifar10 corruption}
\end{table*}

\if 0
\begin{table*}[!htbp]
    \centering

        \begin{tabular}{ l|cc|ccc|c }
         \hline
         ~ & Standard & JaFR & FGSM AT & FGSM AT & PGD AT & $\mathcal{B}_{\text{low}}$ \\
         ~ & ~ & ~ & ~ & + JaFR & ~ & ~ \\
        \hline
        Clean & 93.11$\pm$0.20 & 93.13$\pm$0.11 & 84.80$\pm$1.37 & 79.94$\pm$0.22 & 79.31$\pm$0.23 & - \\
        \hline
        mCE & 100.00 & 95.84 & 124.32 & 134.37 & 135.71 & - \\
        \hline
        Fog & 86.23$\pm$0.45 & \textbf{86.60$\pm$0.28}  & \textbf{63.00$\pm$4.03}  & 56.60$\pm$0.72 & 56.36$\pm$1.25 & 12.85\\
        Saturate & \textbf{89.27$\pm$0.16}  & 89.07$\pm$0.28 & \textbf{79.54$\pm$1.07} & 76.61$\pm$0.22 & 76.29$\pm$0.21 & 12.28\\
        Contrast & \textbf{75.10$\pm$1.08} & 73.45$\pm$0.95 & \textbf{44.06$\pm$3.97} & 41.96$\pm$0.61 & 41.19$\pm$0.88 & 12.14\\
        Bright & \textbf{91.56$\pm$0.15} & 91.45$\pm$0.18 & \textbf{83.04$\pm$1.20} & 76.95$\pm$0.42 & 76.01$\pm$0.58 & 12.01\\
        Snow & 79.54$\pm$0.71  & \textbf{81.06$\pm$0.15} & \textbf{77.83$\pm$1.43} & 74.28$\pm$0.27 & 73.73$\pm$0.21 & 11.53\\
        Frost & 75.41$\pm$0.63 & \textbf{78.97$\pm$0.56} & \textbf{77.70$\pm$1.44} & 70.50$\pm$0.54 & 69.06$\pm$1.18 & 10.61\\
        Motion & \textbf{75.67$\pm$1.35} & 75.13$\pm$1.21 & 65.15$\pm$3.47 & \textbf{72.65$\pm$0.41}  & 72.15$\pm$0.78 & 10.61\\
        Zoom & 76.42$\pm$1.65 & \textbf{76.48$\pm$1.20}  & 69.37$\pm$3.91 & \textbf{75.38$\pm$0.44} & 74.78$\pm$0.67 & 10.37\\
        Elastic & 81.62$\pm$0.60 & \textbf{81.98$\pm$0.55} & 73.90$\pm$2.27 & \textbf{74.86$\pm$0.37} & 74.25$\pm$0.58 & 10.05\\
        Pixel & \textbf{72.90$\pm$0.70} & 72.45$\pm$0.89 & \textbf{81.20$\pm$0.55} & 78.15$\pm$0.14 & 77.47$\pm$0.33 & 8.2\\
        Gauss. B & 73.38$\pm$1.58 & \textbf{73.53$\pm$1.16}  & 68.49$\pm$3.93 & \textbf{74.39$\pm$0.38} & 73.71$\pm$0.63 & 7.82\\
        Defocus & 81.57$\pm$0.83 & \textbf{81.76$\pm$0.71} & 73.42$\pm$3.01 & \textbf{76.15$\pm$0.29} & 75.46$\pm$0.54 & 7.49\\
        Glass & 50.18$\pm$1.75 & \textbf{55.12$\pm$0.58} & 70.18$\pm$1.07 & \textbf{73.75$\pm$0.28} & 73.33$\pm$0.49 & 7.01\\
        Spatter & 80.77$\pm$0.11 & \textbf{82.29$\pm$0.45} & 78.49$\pm$2.04 & \textbf{75.83$\pm$0.32} & 75.29$\pm$0.32 & 6.77\\
        JPEG & 77.76$\pm$0.56 & \textbf{79.56$\pm$0.60} & \textbf{81.47$\pm$1.54} & 78.19$\pm$0.18 & 77.58$\pm$0.24 & 6.51\\
        Speckle & 62.79$\pm$2.44 & \textbf{66.74$\pm$0.51} & 74.07$\pm$4.08 & \textbf{77.36$\pm$0.21} & 76.22$\pm$0.54 & 3.76\\
        Shot & 59.63$\pm$2.75 & \textbf{64.32$\pm$0.76} & 73.73$\pm$4.19 & \textbf{77.35$\pm$0.24} & 76.19$\pm$0.42 & 3.68\\
        Impulse  & 56.50$\pm$1.64 & \textbf{58.75$\pm$0.59} & 64.51$\pm$8.91 & 71.55$\pm$0.54 & \textbf{71.56$\pm$0.74} & 3.63\\
        Gauss. N & 48.70$\pm$3.31 & \textbf{54.34$\pm$1.00} & 69.40$\pm$5.27 & \textbf{76.59$\pm$0.33} & 75.15$\pm$0.49 & 3.62\\
        \hline
        \end{tabular}
    \caption{Accuracy values ($\uparrow$ better) and mCE ($\downarrow$ better) for different models under CIFAR-10 corruptions. The corruption types are arranged with descending order of low-frequency bias $\mathcal{B}_{\text{low}}$.}
\label{tab:cifar10 corruption}
\end{table*}
\fi

\section{Experiments}
We conduct experiments across 4 image datasets (SVHN, CIFAR-10 and CIFAR-100, TinyImageNet) \cite{netzer2011reading,krizhevsky2009learning} to study the effect of JaFR on model robustness against adversarial examples and common image corruptions.
We largely follow the training setting as \cite{andriushchenko2020understanding} where the PreAct ResNet-18 architecture \cite{he2016identity} is used for all models. To benchmark JaFR's effect on adversarial robustness, we compare with a relatively weak AT baseline (FGSM AT) and 2 strong AT baselines (FGSM AT + GradAlign \cite{andriushchenko2020understanding} and PGD AT \cite{madry2017towards}). All AT models use $\epsilon=\frac{8}{255}$. The dataset-specific parameters are further detailed in the following sections. Experiments are run on Nvidia V100 GPUs.

\textbf{SVHN}
All models were trained for 15 epochs, with a cyclic learning rate schedule and a max learning rate of 0.05. We use $\lambda_{freq}=0.001$ and $\lambda_{freq}=0.05$ for the JaFR and FGSM AT + JaFR model respectively. FGSM + GradAlign model uses the same hyperparameters in \cite{andriushchenko2020understanding} for training.

\textbf{CIFAR-10}
The training setup largely follows the SVHN experiments where the PreAct ResNet-18 architecture is used for all models. All models were trained for 30 epochs, with a cyclic learning schedule. The standard-trained, JaFR and FGSM models use a max learning rate of 0.2 while the FGSM + JaFR, FGSM + GradAlign and PGD AT models used a max learning rate of 0.3. Both JaFR and FGSM AT + JaFR models use $\lambda_{freq}=0.001$.

\textbf{CIFAR-100}
The training setup is similar to that in CIFAR-10, except for the learning rates and $\lambda_{freq}=0.001$ values. The max learning rates for all models are 0.3 except for the JaFR (0.1) and standard-trained (0.2) models.

\textbf{TinyImageNet}
The training setup is similar to that in CIFAR-100, except for the $\lambda_{freq}=0.002$ value.

\if 0
\begin{table*}[!htbp]
    \centering
    \footnotesize

        \begin{tabular}{ l|cccccc }
         \hline
         Corruption & fog & saturate & contrast & brightness & snow & frost \\
         \hline
         $\mathcal{B}_{\text{low}}$ & 12.85 & 12.28 & 12.14 & 12.01 & 11.53 & 10.61  \\
         $\mathbb{E} \left[| \mathcal{F}(C(\xinput) - \xinput) | \right]$ &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/fog.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/saturate.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/contrast.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/brightness.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/snow.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/frost.png}} \\
         \vspace{0.01cm} \\
         $ C(\xinput) $ &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/fog.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/saturate.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/contrast.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/brightness.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/snow.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/frost.png}} \\
        \end{tabular}

        \vspace{1em}

        \begin{tabular}{ l|cccccc }
         \hline
         Corruption &  spatter & jpeg & speckle & shot & impulse & gauss noise \\
         \hline
         $\mathcal{B}_{\text{low}}$ & 6.77 & 6.51 & 3.76 & 3.68 & 3.63 & 3.62 \\
         $\mathbb{E} \left[| \mathcal{F}(C(\xinput) - \xinput) | \right]$ &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/spatter.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/jpeg_compression.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/speckle_noise.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/shot_noise.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/impulse_noise.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_fourier_map/gaussian_noise.png}} \\
         \vspace{0.01cm} \\
         $ C(\xinput) $ &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/spatter.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/jpeg_compression.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/speckle_noise.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/shot_noise.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/impulse_noise.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_c_sample_img/gaussian_noise.png}} \\
        \end{tabular}

    \caption{CIFAR-10 corruptions arranged according to their frequency bias, i.e., `fog' with largest ratio of low-frequency components indicated by its bigger $\mathcal{B}_{\text{low}}$ and frequency spectrum that shows concentration of intensity at the top-left.}

\label{tab:cifar10 corruption stats}
\end{table*}
\fi

\if 0
\begin{table*}[!htbp]
    \centering

        \begin{tabular}{ l|cc|ccc|c }
         \hline
         ~ & Standard & JaFR & FGSM AT & FGSM AT & PGD AT & $\mathcal{B}_{\text{low}}$ \\
         ~ & ~ & ~ & ~ & + JaFR & ~ & ~ \\
        \hline
        Clean & 72.27$\pm$0.31 & 74.59$\pm$0.37 & 52.08$\pm$2.85 & 50.49$\pm$0.21 & 49.57$\pm$0.44 & - \\
        \hline
        mCE & 100.00 & 96.29 & 120.03 & 114.70 & 116.50 & - \\
        \hline
        Fog & 58.58$\pm$0.54 & \textbf{61.97$\pm$0.52} & \textbf{31.83$\pm$3.36} & 27.39$\pm$0.45  & 26.69$\pm$0.78 & 12.91\\
        Saturate & 59.48$\pm$0.26 & \textbf{61.11$\pm$0.24} & \textbf{39.98$\pm$2.12} & 38.93$\pm$0.09 &  38.53$\pm$0.45 & 12.33\\
        Contrast & 48.88$\pm$0.72 & \textbf{51.78$\pm$0.66} & \textbf{22.27$\pm$2.84} & 20.22$\pm$0.41 & 19.38$\pm$0.74 & 12.15\\
        Bright & 67.69$\pm$0.20 & \textbf{70.19$\pm$0.34} & \textbf{48.63$\pm$2.58} & 44.87$\pm$0.31 & 44.28$\pm$0.49 & 11.98\\
        Snow & 50.46$\pm$0.28 & \textbf{54.50$\pm$0.49} &  41.96$\pm$2.81 & \textbf{42.44$\pm$0.49} &  42.25$\pm$0.78 & 11.56\\
        Frost & 45.82$\pm$0.95 & \textbf{52.18$\pm$0.47} & \textbf{42.06$\pm$2.45} & 37.22$\pm$0.65 & 36.54$\pm$1.07 & 10.77\\
        Motion & 51.22$\pm$0.56 & \textbf{52.22$\pm$0.96} & 36.73$\pm$4.48 & \textbf{43.24$\pm$0.36} & 42.16$\pm$0.63 & 10.63\\
        Zoom & 50.44$\pm$1.18 & \textbf{51.25$\pm$1.54} & 39.13$\pm$4.59 & \textbf{45.40$\pm$0.25} & 44.19$\pm$0.71 & 10.34\\
        Elastic & 56.66$\pm$0.31 & \textbf{58.45$\pm$0.48} & 42.17$\pm$3.44 & \textbf{44.88$\pm$0.24} & 43.93$\pm$0.47 & 10.03\\
        Pixel & 48.19$\pm$1.08 & \textbf{50.41$\pm$0.49} & 48.82$\pm$2.31 & \textbf{48.89$\pm$0.24} & 47.86$\pm$0.48 & 8.11\\
        Gauss. B & 48.64$\pm$1.22 & \textbf{48.78$\pm$1.49} & 38.42$\pm$4.46 & \textbf{44.43$\pm$0.24 }& 43.31$\pm$0.59 & 7.66\\
        Defocus & 56.52$\pm$0.89 & \textbf{57.44$\pm$1.22} & 42.12$\pm$3.96 & \textbf{46.15$\pm$0.21} & 44.99$\pm$0.49 & 7.38\\
        Glass & 21.32$\pm$0.37 & \textbf{24.79$\pm$0.46} & 38.01$\pm$3.00 & \textbf{45.24$\pm$0.24} & 44.33$\pm$0.40 & 7.02\\
        Spatter & 54.94$\pm$0.56 & \textbf{57.21$\pm$0.88} & 44.91$\pm$2.89 & \textbf{45.70$\pm$0.30} & 44.91$\pm$0.39 & 6.76\\
        JPEG & \textbf{52.04$\pm$0.43} & 50.40$\pm$0.65 & 48.33$\pm$2.79 & \textbf{48.50$\pm$0.33} & 47.48$\pm$0.43 & 6.6\\
        Speckle & 32.23$\pm$0.79 & \textbf{34.56$\pm$0.58} & 34.58$\pm$4.34 & \textbf{46.26$\pm$0.43} & 44.46$\pm$0.68 & 3.86\\
        Shot & 31.49$\pm$0.84 & \textbf{33.46$\pm$0.61} & 33.95$\pm$4.76 & \textbf{47.09$\pm$0.37} & 45.43$\pm$0.57 & 3.74\\
        Gauss. N& 24.21$\pm$0.78 & \textbf{25.14$\pm$0.56} & 28.47$\pm$6.33 & \textbf{46.32$\pm$0.44} & 44.57$\pm$0.49 & 3.68\\
        Impulse & \textbf{26.52$\pm$1.54} & 24.20$\pm$2.07 & 27.31$\pm$6.27 & \textbf{37.40$\pm$0.87} & 37.10$\pm$0.94 & 3.64\\
        \hline
        \end{tabular}
    \caption{Accuracy values ($\uparrow$ better) and mCE ($\downarrow$ better) for different models under CIFAR-100 corruptions. The corruption types are arranged with descending order of low-frequency bias $\mathcal{B}_{\text{low}}$.}
\label{tab:cifar100 corruption}
\end{table*}
\fi

\subsection{Common Corruption Robustness}
CIFAR-10-C and CIFAR-100-C are common corruption benchmarks where the CIFAR-10/CIFAR-100 test sets are corrupted with an array of 19 corruption types with 5 levels of severity. The accuracy values are reported as the average over the 5 severities of each corruption. Table~\ref{tab:cifar10 corruption stats} (in Appendix) shows top-6 corruptions that have the highest ratio of low- and high-frequency components as measured by $\mathcal{B}_{\text{low}}$.

\textbf{JaFR improves corruption robustness for standard-trained models.} From Table~\ref{tab:cifar10 corruption} and \ref{tab:cifar100 corruption} (Appendix), we observe that the JaFR model achieves the best overall corruption robustness for both CIFAR-10 and -100 images, as indicated by its lowest relative mean corruption errors (mCE). The JaFR outperforms the second-best model (standard-trained) in 14 out of 19 corruptions for CIFAR-10-C and 17 out of 19 corruptions for CIFAR-100-C test samples. This is an encouraging sign for direct frequency regularization of neural networks as an approach for better corruption robustness.

When comparing the spectra of the standard and JaFR model in Table~\ref{tab:cifar10 model fourier profile}, we observe that the intensities shift from the highest-frequency (bottom-right) region to lower frequency regions (left and top) with JaFR. This marks a lower reliance on high-frequency features in the model, which can confer better robustness against high-frequency corruptions. Indeed, with JaFR, we see the biggest improvements for high-frequency corruptions (e.g., Glass, Shot and Gaussian noise). Moreover, we can see that the spectra of our JaFR-only model, unlike the adversarially trained models (e.g., PGD AT), do not concentrate the intensities on the lowest-frequency region (top-left) which may result in an over-reliance on low-frequency features and make the model more susceptible to low-frequency corruptions. During our experiments, larger $\lambda_{freq}$ values for JaFR do concentrate the intensities even more in the low-frequency regions but result in poorer overall performance. This balance in reliance on both high and low-frequency features may explain the improved overall performance of JaFR against corruption.

\begin{table*}[!htbp]
    \centering
    \footnotesize
    \setlength{\tabcolsep}{0.4em}
        \begin{tabular}{ l|cccccccc }
         \hline
         Model & Standard & JaFR & FGSM AT & FGSM AT & FGSM AT & FGSM AT & PGD AT & Original \\
         ~ & ~ & ~ & ~ & + JaFR(-) & + JaFR & + GradAlign & ~ & Image \\
         \hline
         $\mathbb{E} \left[| \mathcal{F}(J) | \right]$ &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/standard.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/jafr.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/fgsm.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/fgsm_n_jafrminus.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/fgsm_n_jafr.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/fgsm_n_gradalign.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/pgd.png}} & - \\
         \vspace{0.02cm} \\
         $ J $ &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/standard.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/jafr.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/fgsm.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/fgsm_n_jafrminus.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/fgsm_n_jafr.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/fgsm_n_gradalign.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/pgd.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/original.png}} \\
         \hline
         $\mathcal{B}_{\text{low}}$ & 1.56 & 16.71 & 2.68 & 0.326 & 15.41 & 7.88 & 7.49 & -\\
        \end{tabular}

    \caption{CIFAR-10 models' low-frequency bias ($\mathcal{B}_{\text{low}}$), frequency spectra and Jacobians.}
\label{tab:cifar10 model fourier profile}
\end{table*}

\if 0
\begin{table*}[!htbp]
    \centering
    \footnotesize
    \setlength{\tabcolsep}{0.4em}
        \begin{tabular}{ l|ccccccc }
         \hline
         Model & Standard & JaFR & FGSM AT & FGSM AT & FGSM AT & PGD AT & Original \\
         ~ & ~ & ~ & ~ & + JaFR & + GradAlign & ~ & Image \\
         \hline
         $\mathbb{E} \left[| \mathcal{F}(J) | \right]$ &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/standard.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/jafr.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/fgsm.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/fgsm_n_jafr.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/fgsm_n_gradalign.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_fourier_map/pgd.png}} & - \\
         \vspace{0.02cm} \\
         $ J $ &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/standard.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/jafr.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/fgsm.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/fgsm_n_jafr.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/fgsm_n_gradalign.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/pgd.png}} &
         \parbox[c]{4em}{\includegraphics[width=4em]{images/c10_model_ig_sample10/original.png}} \\
         \hline
         $\mathcal{B}_{\text{low}}$ & 1.56 & 16.71 & 2.68 & 15.41 & 7.88 & 7.49 & -\\
        \end{tabular}

    \caption{CIFAR-10 models' low-frequency bias ($\mathcal{B}_{\text{low}}$), frequency spectra and Jacobians.}
\label{tab:cifar10 model fourier profile}
\end{table*}
\fi

\textbf{Tradeoff between robustness against low- and high-frequency corruption exists.} For both CIFAR-10-C (Table~\ref{tab:cifar10 corruption}) and CIFAR-100-C (Table~\ref{tab:cifar100 corruption}), we observe that when JaFR is added to FGSM AT, the performance against low-frequency corruptions (e.g., Fog \& Contrast) drops when the accuracy against high-frequency corruptions (e.g., Impulse \& Gaussian noise) improves. We also see such behavior when JaFR is added to the standard trained model in CIFAR-10-C. This aligns with previous studies that show such a tradeoff between corruptions with different frequency profiles.

To study the opposite effect of JaFR, a JaFR(-) model variant is added where the $\lambda_{freq}$ term has a negative value to investigate the effect of biasing the model towards high-frequency components.
We can see such tradeoff emerging with an opposite effect when JaFR(-) is combined with FGSM AT to bias the model towards high frequency features, where the performance against high-frequency corruptions (e.g., Impulse \& Gaussian noise) drops when the accuracy against low-frequency corruptions (e.g., Fog \& Contrast) improves. This suggests that direct frequency regularization of models is a promising way to tailoring robustness against a set of corruptions that is more relevant for their applications.

\textbf{FGSM AT + JaFR outperforms PGD AT in almost all corruptions.} From our experiments, apart from one case of Impulse corruption in CIFAR-10-C, the FGSM AT + JaFR outperforms the PGD AT model in all the other 37 corruptions scenarios. We speculate that the PGD AT model might have overfitted to resist adversarial perturbations, making them excessively reliant on a small subset of low-frequency features that can be easily disrupted by the common corruptions.

\if 0
\textbf{CIFAR-10-C}
From Table~\ref{tab:cifar10 corruption}, we can see that JaFR outperforms the standard-trained model in 12 out of 19 corruptions. There is a performance drop when JaFR is applied for the three high-frequency corruptions  (brightness, contrast and fog). This shows that training models to focus on low-frequency features makes them more susceptible to low-frequency noises.

FGSM AT + JaFR outperforms FGSM for 10 out of 19 corruption types. These performance gains are concentrated on corruption types with high-frequency noises (e.g., gaussian blur and motion blur) while the performance drops from FGSM to FGSM + JaFR are the highest for the low-frequency corruption types such as brightness, contrast and fog. This corroborates previous findings that training models to be invariant to low-frequency features results in poorer performances in low-frequency corruptions
[to check!]The FGSM AT + JaFR model outperforms the PGD AT model in almost all of the corruption types in CIFAR-10, except for fog corruption. We speculate that the PGD AT model might have overfitted to adversarial perturbations, making them excessively invariant to other important high-frequency features.
\fi

\if 0
\textbf{CIFAR-100-C}
The JaFR model outperforms the standard trained model on 18 out of 19 corruption types, except for the `jpeg compression' corruption, as shown in Table~\ref{tab:cifar100 corruption}. The low-frequency corruption types (e.g., brightness, contrast and fog corruption), which have seen performance drop with JaFR in CIFAR-10, have shown performance boost with JaFR. However, the improvement for these types are the smallest among the corruption types which suggest that the detrimental effect of JaFR is offset up the increase in clean accuracy, explaining the net performance gain we see in them.

We observe that FGSM AT + JaFR outperforms FGSM for 12 out of 19 corruption types. For the other 7 corruption types (contrast, fog, frost, pixelate, saturate, snow and spatter) where FGSM AT + JaFR shows lower performance, we see similar performance drop in the CIFAR-10-C results (Table~\ref{tab:cifar10 corruption}). A commonality amongst these 7 corruptions is their low-frequency nature, indicating JaFR trades off robustness against low-frequency noise for gains under high-frequency corruption.

The FGSM AT + JaFR model outperforms the PGD model in almost all of the corruption types here, except for frost corruption. Similar to CIFAR-10-C, this could be explained by the overfitting of the PGD AT model to be in resistant to high-frequency adversarial perturbation, causing it to not rely on other important high-frequency features in the input images.
\fi

\subsection{Adversarial Robustness}
We evaluate the adversarial robustness of models with FGSM \cite{goodfellow2014explaining} and PGD \cite{madry2017towards} attacks. PGD uses 50 gradient iterations and 10 restarts with a step size of $\alpha_{adv}=\epsilon/4$.

\textbf{Biasing towards low frequency can boost adversarial robustness in weakly adversarially trained models.}
For all the four image datasets, when JaFR is combined with FGSM, there is a significant boost in both clean and adversarial accuracy values from the FGSM AT model for all attack types (see Table~\ref{tab:svhn}, \ref{tab:cifar10}, \ref{tab:cifar100} and \ref{tab:tinyimagenet}). This observation is similar to what is observed for a recent defense, GradAlign, which also requires combining with FGSM to improve the adversarial robustness of the model. The need of using FGSM training samples to see an improvement in JaFR's adversarial robustness indicates that strong adversarial examples are more than just high-frequency noise and are able to find regions of high error in the uneven loss landscape of a non-adversarially trained model \cite{liu2020loss}. For the CIFAR-10 and -100, the improvement of adding JaFR to FGSM is large enough to achieve a robustness level that is competitive with the strong PGD AT baseline. In Table C, we observe that JaFR scales well to a larger dataset, with a smaller amount of computation than PGD AT. Conversely, combining with JaFR(-) to bias the FGSM AT model towards high frequency with a negative $\lambda_{freq}$ term was observed in our experiments to worsen its adversarial robustness against PGD attacks. Experiments (Table~\ref{tab:more attacks} in Appendix) on more advanced attacks such as CW \cite{carlini2017towards} and AutoAttack \cite{croce2020reliable} show that the robustness gain from JaFR can resist even stronger attacks and is not due to gradient masking. Moreover, from Table~\ref{tab:tinyimagenet}, we observe that JaFR scales well, with a smaller amount of computation than the PGD AT.

\begin{table}[!htbp]
    \centering

    \footnotesize
        \begin{tabular}{ l|ccc }
         \hline

         Model & Clean & FGSM & PGD \\
         \hline
        Standard & 96.62$\pm$0.05 &  21.86$\pm$0.90 & 0.17$\pm$0.02 \\
        JaFR &  96.58$\pm$0.07  &  21.94$\pm$1.16 & 0.14$\pm$0.03  \\
        FGSM AT &  92.33$\pm$0.20  &  89.16$\pm$5.87 & 0.52$\pm$0.90   \\

        FGSM + JaFR & 87.22$\pm$4.63 & 59.4$\pm$12.75 & 19.64$\pm$16.26   \\
        FGSM + GA & 92.51$\pm$0.26 & 59.52$\pm$0.26 & 46.63$\pm$0.15 \\
        PGD AT & 92.13$\pm$0.51 & 62.38$\pm$0.86 & 56.79$\pm$0.29 \\
         \hline
        \end{tabular}
    \caption{SVHN accuracy (\%) on clean and adversarial test samples.}
\label{tab:svhn}
\end{table}

\if 0
\begin{table}[!htbp]
    \centering

    \footnotesize
        \begin{tabular}{ l|ccc }
         \hline

         Model & Clean & FGSM & PGD \\
         \hline
        Standard & 96.62$\pm$0.05 &  21.86$\pm$0.90 & 0.17$\pm$0.02 \\
        JaFR &  96.58$\pm$0.07  &  21.94$\pm$1.16 & 0.14$\pm$0.03  \\
        FGSM AT &  92.33$\pm$0.20  &  89.16$\pm$5.87 & 0.52$\pm$0.90   \\
        FGSM + JaFR & 87.22$\pm$4.63 & 59.4$\pm$12.75 & 19.64$\pm$16.26   \\
        FGSM + GA & 92.51$\pm$0.26 & 59.52$\pm$0.26 & 46.63$\pm$0.15 \\
        PGD AT & 92.13$\pm$0.51 & 62.38$\pm$0.86 & 56.79$\pm$0.29 \\
         \hline
        \end{tabular}
    \caption{SVHN accuracy (\%) on clean and adversarial test samples.}
\label{tab:svhn}
\end{table}
\fi

\begin{table}[!htbp]
    \centering
    \footnotesize
        \begin{tabular}{ l|ccc }
         \hline

         Model & Clean & FGSM & PGD \\
         \hline
        Standard & 93.11$\pm$0.20 & 16.04$\pm$0.85 & 0 \\
        JaFR & 93.13$\pm$0.11 & 15.86$\pm$0.88 & 0 \\
        FGSM AT & 84.80$\pm$1.37 & 85.14$\pm$4.01 & 0.01$\pm$0.02 \\

        FGSM + JaFR & 79.94$\pm$0.22 & 53.12$\pm$0.25 & 46.32$\pm$0.15 \\
        FGSM + GA & 80.07$\pm$0.21 & 54.26$\pm$0.55 & 46.97$\pm$0.15 \\
        PGD AT & 79.31$\pm$0.23 & 54.00$\pm$0.72 & 49.63$\pm$0.20 \\
         \hline
        \end{tabular}
    \caption{CIFAR-10 accuracy (\%) on clean and adversarial test samples.}
\label{tab:cifar10}
\end{table}

\if 0
\begin{table}[!htbp]
    \centering
    \footnotesize
        \begin{tabular}{ l|ccc }
         \hline

         Model & Clean & FGSM & PGD \\
         \hline
        Standard & 93.11$\pm$0.20 & 16.04$\pm$0.85 & 0 \\
        JaFR & 93.13$\pm$0.11 & 15.86$\pm$0.88 & 0 \\
        FGSM AT & 84.80$\pm$1.37 & 85.14$\pm$4.01 & 0.01$\pm$0.02 \\
        FGSM + JaFR & 79.94$\pm$0.22 & 53.12$\pm$0.25 & 46.32$\pm$0.15 \\
        FGSM + GA & 80.07$\pm$0.21 & 54.26$\pm$0.55 & 46.97$\pm$0.15 \\
        PGD AT & 79.31$\pm$0.23 & 54.00$\pm$0.72 & 49.63$\pm$0.20 \\
         \hline
        \end{tabular}
    \caption{CIFAR-10 accuracy (\%) on clean and adversarial test samples.}
\label{tab:cifar10}
\end{table}
\fi

\if 0
\textbf{JaFR can improve clean accuracy when performance is low.}
For CIFAR-100, the JaFR model outperforms the standard-trained model by a larger margin as compared to what is observed in the SVHN and CIFAR-10 experiments, as shown in Table~\ref{tab:cifar100}. This suggests that when the model performance is relatively low (72.27\% in CIFAR-100 vs 91.11\% in CIFAR-10), the model can achieve better performance when it relies on (low-frequency) components that more closely match the frequency profile of input images.
\fi

\if 0
\textbf{Biasing towards high frequency is detrimental to adversarial robustness.}
Combining with JaFR(-) to bias the FGSM AT model towards high frequency with a negative $\lambda_{freq}$ term worsens its adversarial robustness, in contrast to using a positive $\lambda_{freq}$ term. We can see this from the drop in accuracy values for FGSM + JaFR(-) from the FGSM AT model, in Table~\ref{tab:svhn}, \ref{tab:cifar10} and \ref{tab:cifar100}, under the FGSM and PGD adversarial examples. This aligns with previous findings \cite{yin2019fourier} that dependence on high-frequency features results in higher susceptibility to adversarial examples.
\fi

\textbf{JaFR removes high-frequency components in models' Jacobian.} From Table~\ref{tab:cifar10 model fourier profile}, we can see that models with JaFR have the highest low-frequency bias ($\mathcal{B}_{\text{low}}$) values, indicating that the training regularization is successful in increasing the low-frequency components of the model's Jacobians. Furthermore, from the Jacobian's Fourier spectra ($\mathbb{E} \left[| \mathcal{F}(J) | \right]$), intensities of the spectra shift from the high-frequency regions (right and bottom) towards the left and top regions of the spectra which indicate the low-frequency components along the horizontal and vertical axis respectively. In contrast, using a negative $\lambda_{freq}$ term (in FGSM AT + JaFR(-)) concentrates the spectra towards highest-frequency (bottom-right) region and drastically lowers the ($\mathcal{B}_{\text{low}}$) value which indicates a strong reliance in high-frequency features.

\begin{table}[!htbp]
    \centering
    \footnotesize
        \begin{tabular}{ l|ccc }
         \hline

         Model & Clean & FGSM & PGD \\
         \hline
        Standard & 72.27$\pm$0.31 & 8.04$\pm$0.84 & 0.13$\pm$0.03 \\
        JaFR & 74.59$\pm$0.37 & 8.60$\pm$0.81 & 0.08$\pm$0.04 \\
        FGSM AT & 52.08$\pm$2.85 & 31.43$\pm$2.15 & 0.04$\pm$0.04\\

        FGSM + JaFR & 50.49$\pm$0.21 & 26.96$\pm$0.36 & 23.42$\pm$0.27 \\
        FGSM + GA & 50.68$\pm$0.28 & 27.44$\pm$0.68 & 23.93$\pm$0.19  \\
        PGD AT & 49.57$\pm$0.44 & 27.68$\pm$0.77 & 25.24$\pm$0.25 \\
         \hline
        \end{tabular}
    \caption{CIFAR-100 accuracy (\%) on clean and adversarial test samples.}
\label{tab:cifar100}
\end{table}

\begin{table}[!htbp]
    \centering

    \footnotesize
        \begin{tabular}{ l|ccc|c }
         \hline

         Model & Clean & FGSM & PGD & Complex. (s/iter) \\
         \hline
        Standard & 63.59 & 2.23 & 0.05 & 0.1575 \\
        JaFR & 63.44 & 1.87 & 0.04 & 0.6015  \\
        FGSM AT & 45.89 & 17.54 & 0 & 0.4245 \\

        FGSM + JaFR &  47.43 & 20.3 & 18.86 & 0.6027  \\
        FGSM + GA & 40.69 & 19.86 & 17.81 & 0.6003  \\
        PGD AT & 40 & 17.54 & 19.91 & 1.34 \\
         \hline
        \end{tabular}
    \caption{Model accuracy and computational complexity on the TinyImageNet dataset.}
\label{tab:tinyimagenet}
\end{table}

\if 0
\begin{table}[!htbp]
    \centering
    \footnotesize
        \begin{tabular}{ l|ccc }
         \hline

         Model & Clean & FGSM & PGD \\
         \hline
        Standard & 72.27$\pm$0.31 & 8.04$\pm$0.84 & 0.13$\pm$0.03 \\
        JaFR & 74.59$\pm$0.37 & 8.60$\pm$0.81 & 0.08$\pm$0.04 \\
        FGSM AT & 52.08$\pm$2.85 & 46.08$\pm$18.03 & 0.04$\pm$0.04\\
        FGSM + JaFR & 50.49$\pm$0.21 & 26.96$\pm$0.36 & 23.42$\pm$0.27 \\
        FGSM + GA & 50.68$\pm$0.28 & 27.44$\pm$0.68 & 23.93$\pm$0.19  \\
        PGD AT & 49.57$\pm$0.44 & 27.68$\pm$0.77 & 25.24$\pm$0.25 \\
         \hline
        \end{tabular}
    \caption{CIFAR-100 accuracy (\%) on clean and adversarial test samples.}
\label{tab:cifar100}
\end{table}
\fi

\textbf{JaFR improves saliency of Jacobians.} When JaFR is added to the FGSM AT model, we observe that its Jacobians ($J$) become more salient (see Table~\ref{tab:cifar10 model fourier profile}). Together with the boost in adversarial robustness from the FGSM AT model, this corroborates previous results that saliency of Jacobian is correlated with adversarial robustness \cite{etmann2019connection}.

\if 0
\textbf{SVHN}
With only JaFR, we observe no significant change in both clean and adversarial accuracy values from the standard-trained model. In contrast, when JaFR is combined with FGSM, there is a significant boost in both clean and adversarial accuracy values from the FGSM AT model. This observation is similar to what is observed for a recent defense, GradAlign, which also requires combining with FGSM to improve the adversarial robustness of the model. The frequency components are weighted more on the low-frequency region [new figure: show fourier map] after adding JaFR as compared to the FGSM model, showing that the JaFR is achieving what it set out to do: shifting the components towards to low-frequency region.
\fi

\if 0
\textbf{CIFAR-10}
From Table~\ref{tab:cifar10}, we observe that the model slightly outperforms the standard trained model in clean accuracy when JaFR is added. However, JaFR does not offer any significant benefit for adversarial robustness on top of the standard trained model, similar to what is observed for SVHN. In contrast, when JaFR is combined with FGSM, there is a significant boost in both clean and adversarial accuracy values from the FGSM AT model, achieving robustness levels that are competitive with the strong PGD AT baseline. As compared to the FGSM AT + GradAlign, the performance of clean and adversarial accuracy values are closer to each other than in SVHN, suggesting that the robustness improvement from low-frequency bias may vary across image datasets.
\fi

\if 0
\textbf{CIFAR-100}
For CIFAR-100, the JaFR model outperforms the standard-trained model by a larger margin as compared to what is observed in the CIFAR-10 experiments, as shown in Table~\ref{tab:cifar100}. This suggests that when the model performance is relatively low (71.9\% here vs 92.9\% in CIFAR-10), the model can achieve better performance when it relies on (low-frequency) components that more closely match the (frequency components of) input images. Similar to CIFAR-10 results, by JaFR alone, there is no gain in adversarial robustness over the standard-trained model. However, when FGSM AT is present, we observe a significant boost to the model's performance under adversarial examples, to levels that are competitive to the PGD AT model.
\fi

\section{Conclusion}
Model robustness is growing more important as deep learning models are gaining wider adoption. Here, we delve further into the link between the frequency characteristic of a model and its robustness by proposing Jacobian frequency bias. Through this term, we can control the distribution of high- and low-frequency components in the model's Jacobian and find that it can affect both corruption and adversarial robustness. We hope that our findings here will open an avenue for future work to explore other frequency-focused approaches to improve model robustness.

\bibliographystyle{named}


\section{Introduction}

\subsection{Asset return distributions}
There is a long history of researchers observing that
the tails of asset returns are not well modeled by a Gaussian distribution,
going back to the thesis of Fama in 1965 \cite{Fama},
who observed that
while somewhat symmetric, the tails of the return distribution were much
heavier than those of a Gaussian distribution. Additionally, asset returns are
skewed, violating normality \cite{neuberger-skewness}.
There is also long history of researchers proposing alternative distributions
to model asset returns, including elliptical distributions \cite{bingham-2001},
and Gaussian mixtures (GMs) \cite{ball-1983,akgiray-1987}, the focus of this paper.

GMs can in principle approximate any continuous
distribution; for asset returns, it has been observed that good
approximations can be obtained with just a handful of mixture
components \cite{kon-1984}.
We can interpret the components of a GM return distribution as market regimes,
with a latent variable that represents the active regime, and
a return distribution that is Gaussian, given the regime.
Many authors have observed that the
correlations among asset returns can change during different market regimes,
for example, increased correlations during bear markets
\cite{bear-correlation,ang2004regimes,ang2015regimes}.  A GM can model such regime-dependent
correlation structures. Another desirable atribute of a GM is that it can model
skewness in return distributions, for which many authors have argued real investors
exhibit preferences \cite{arditti1967risk,scott1980direction}.

GM return models arise in a hidden Markov Gaussian model \cite{ryden-1998} of returns,
which models the regimes as Markovian, and the returns as
Gaussian, given the regime. 1/
In such a model the means and covariances
of the Gaussian components corresponding to the regimes are fixed, but the
component weights change in each period \cite{gupta-2012}; but in each period, the
asset return distribution, conditioned on the past returns, is GM, so the methods
in this paper can be applied.

\subsection{Mean-variance versus expected utility}

\paragraph{Mean-variance portfolio construction.}
In mean-variance portfolio construction, pioneered by Markowitz \cite{Markowitz1952},
portfolio constuction is viewed as an optimization problem with
two main objectives: the mean or expected return of the portfolio, and the risk, taken
to be the variance of the portfolio return.  These objectives are
combined into a risk-adjusted return using a positive weight parameter,
interpreted as setting the level of risk aversion.
Mean-variance portfolio construction can be carried out
analytically, when there are very simple constraints,
or numerically, with realistic constraints,
by solving a convex optimization problem such as a
quadratic program \cite{grinold-1999,cvxbook,boyd-2017,osqp}.
With current convex optimization methods, mean-variance construction
can be done reliably
and quickly for up to thousands of assets, and many more
when a factor risk model is used.
These optimization problems can be solved in well under one second,
allowing back-tests (what-if simulations, based on real or
simulated data) to be carried out quickly \cite{boyd-2017,cvxpygen,mathworks}.

One obvious criticism of mean-variance portfolio construction
is that the (quadratic) objective function penalizes returns that are well
above the mean (a desirable event) just as much as returns
that are well below the mean (an undesirable event) \cite{hanoch1970efficient}.
Another is that it only uses the first two moments of the
return distribution, and so cannot take into account skewed or
fat-tailed distributions.  Nevertheless it is very widely used
in practice. There is work in analyzing portfolios under higher order Taylor
approximations of utilities \cite{jondeau2006optimal}, but these are not used in practice due to
both semantic and computational complexity.

\paragraph{Expected utility portfolio construction.}
In work that predates mean-variance portfolio construction,
Von Neumann and Morgenstern \cite{utility} introduced the notion of utility to
model decision making with uncertain outcomes. A utility function specifies a
value indexing an agent's preference for each specific outcome; their
theory posits that the agent makes a choice so as to maximize her expected utility.

Portfolio construction by expected utility maximization also frames the
problem as an optimization problem.  The trader specifies a utility
function that is concave and increasing, and the objective (to be
maximized) is the expected utility under the return distribution.
This formulation avoids the awkward situation in mean-variance portfolio
construction where high portfolio returns are
considered bad.  Expected utility maxmization better captures the asymmetry in downside
versus upside risks than mean-variance optimization.  Since
the return distribution is arbitrary, expected utility can directly handle
return distributions with skew or fat tails.

Expected utility maximization, like mean-variance optimization,
leads to a convex optimization problem,
more specifically, a stochastic optimization problem \cite{shapiro}.
Almost all expected utility methods for constructing portfolios
work with samples of the asset returns.
This can be considered an advantage, since it means that such methods
can work with any return distribution from which we can sample returns.
The disadvantage is that sample-based optimization, while tractable,
can be slow compared to mean-variance methods, and scales poorly with problem size.

There are several related portfolio construction methods that rely on
return samples and stochastic convex optimization.
One is based on conditional value at risk (CVaR) \cite{Rockafellar-CVaR,CVaR-gen}.
A more recently proposed method uses entropic value at risk
(EVaR) \cite{EVaR-theory,Cajas}, which we address in \S\ref{s-EVaR}.
Both of these are coherent measures of risk \cite{Artzner,rockafellar-coherent},
and result in convex stochastic optimization problems.

Sample based stochastic optimization methods are used in practice \cite{grinold1999mean},
but far less often than methods based on mean-variance optimization,
which do not involve samples.  This is part because solving sample based stochastic
convex optimization problems is tractable, but far more involved than solving
the convex optimization problems that do not involve return samples,
{\it e.g.}, mean-variance optimization or the methods proposed in this paper.

\paragraph{Comparison.}
These two main approaches, mean-variance optimization and expected utility maximization,
are not as different as they might seem.
Levy and Markowitz \cite{levy-mark-utility} show that maximizing a second order
Taylor approximation of a utility function is equivalent to mean-variance
optimization.  So very roughly speaking, mean-variance optimization is the
second order approximation of expected utility optimization.

When the returns are Gaussian, and we use an exponential utility,
mean-variance and expected utility optimization are
not merely close, but exactly the same.
It appears that Merton was the first to note this connection
\cite{Merton1969}, but his observation does not seem to be
mentioned often after that (see also \S\ref{s-special}).

\subsection{This paper}
In this paper we consider a GM model for asset returns,
and maximize expected utility with a generic utility function,
the exponential utility \cite{saha-1993}.
We refer to this type of portfolio construction as EGM,
for exponential utility with Gaussian mixture returns.
We show that the EGM portfolio construction problem can be solved
exactly as a convex optimization problem, without the need for
any samples from the distribution or other approximations.
The EGM portfolio construction problem is not only convex,
but is easily specified in just a few lines
of code in domain-specific languages (DSLs) for convex optimization
such as CVXPY \cite{diamond-2016}, CVX \cite{grant2014cvx},
or CVXR \cite{cvxr2020}.
Thus EGM combines the efficiency, reliability, and scalability
of mean-variance optimization with the ability of expected utility
maximization to handle non-Gaussian returns and the asymmetry
in our preferences.
When the GM has only one component, our return model is Gaussian,
and EGM reduces to mean-variance optimization.  Thus we can think of EGM as an
extension of mean-variance optimization, or as a special case
of expected utility where the
problem can be solved exactly, without any return samples.

We also show that EGM is closely related to portfolio construction methods
based on the entropic value at risk (EVaR).  With GM return model,
we show that EVaR portfolio construction problem leads to a convex
optimization problem that, like EGM, does not involve sample based
stochastic optimization.

\subsection{Previous and related work}
\paragraph{Portfolio construction with Gaussian mixture returns.}
In \cite{buckley-gmm},
Buckley et al.\ consider a two-component Gaussian mixture
of tranquil and distressed regimes, and analyze several objectives, including
Markowitz, Sharpe ratio, exponential utility, and lower partial moments.
In \S3.3.6 of their paper
they derive the closed form expression for expected
exponential utility under Gaussian mixture returns,
but do not observe that maximizing exponential utility leads to a convex problem.
Studying single period portfolios consisting of a risk free asset and a risky
asset, Prigent and Kaffel analyze optimal portfolios under arbitrary utility
functions, and show on historical data that GM return models lead to
significantly different portfolios than those arising from a Gaussian
return model \cite{kaffel-2014}.
\paragraph{EVaR portfolio construction.} In recent work \cite{Cajas},
Cajas develops a disciplined convex
(DCP) formulation of EVaR, with return samples,
which allows it to be used as
either the objective or as a constraint in
portfolio optimization problems specified using DSLs for convex optimization
such as CVXPY.
Since a return distribution that takes on a finite set of values
({\it e.g.}, the empirical distribution of samples) is a special case
of GM, we can consider EGM (with EVaR) as a generalization of
Cajas' formulation.

\subsection{Outline}
We describe the GM return model is \S\ref{s-gmm-ret},
and in \S\ref{s-port-opt} we show that portfolio optimization with
exponential utility is a convex optimization problem.
In \S\ref{s-EVaR} we show that the closely related function EVaR is also convex,
so minimizing it, or adding a limit on it as a constraint,
results in a convex optimization problem.

\section{Gaussian mixture return model}\label{s-gmm-ret}

\subsection{Asset return distribution}
We let $r\in {\mbox{\bf R}}^n$ denote the return of $n$ assets over some specific period.
We model $r$ as having a GM distribution with $k$ components,
\[
r\sim \textbf{GM}(\{\mu_i,\Sigma_i,\pi_i\}_{i=1}^k),
\]
where $\pi_i\in {\mbox{\bf R}}$ are the (positive) component probabilities,
$\mu_i \in {\mbox{\bf R}}^n$ are the component means, and
$\Sigma_i \in {\mbox{\bf R}}^{n \times n}$ are the (symmetric positive definite)
component covariance matrices.

The GM return distribution includes two interesting special cases.  When there is
only one component, it reduces to Gaussian, with $r \sim \mathcal N(\mu_1, \Sigma_1)$.
Another special case arises when $\Sigma_1 = \cdots = \Sigma_k=0$.  Here
$r$ takes on only the values $\mu_1, \ldots, \mu_k$, with probabilities
$\pi_1, \ldots, \pi_k$. We refer to this as a finite values return distribution.

\subsection{Portfolio return distribution}
Let $w\in{\mbox{\bf R}}^{n}$ denote the weights in an investment portfolio, with
$\mathbf 1^Tw =1$, where $\mathbf 1$ is the vector with all entries one. For $w_i\geq 0$,
$w_i$ is the fraction of the total portfolio value invested in asset $i$; for $w_i<0$,
$-w_i$ is the fraction of total portfolio value that is held in a short position in
asset $i$.
The portfolio return is $R=w^Tr$.  This scalar random variable is also GM with
component probabilities $\pi_i$, and means and variances
\[
\nu_i = w^T \mu_i, \quad \sigma_i^2 = w^T\Sigma_i w,  \quad i=1, \ldots, k.
\]

We observe that various quantities associated with
the portfolio return $R$ can be evaluated analytically,
without the need for Monte Carlo or other sampling methods.
For example its cumulative distribution
function (CDF) is given by
\begin{equation}\label{e-gmm-cdf}
\Phi_R(w,a) = \sum_{i=1}^k \pi_i \Phi \left(\frac{a-\nu_i}{\sigma_i}\right) =
\sum_{i=1}^k \pi_i \Phi \left(\frac{a-w^T\mu_i}{(w^T \Sigma_i w)^{1/2}}\right),
\end{equation}
where $\Phi$ is the CDF of a standard Gaussian.

\subsection{Moment and cumulant generating functions}
Two other quantities we will encounter later are the moment generating function
\begin{equation}\label{e-mgf}
M(w,t) = \mathop{\bf E{}}\exp (t R) =
\sum_{i=1}^k\pi_i \exp \left( t \nu_i+\frac{t^2}{2}\sigma_i^2\right) =
\sum_{i=1}^k\pi_i \exp \left( t \mu_i^Tw+\frac{t^2}{2}w^T\Sigma_i w \right),
\end{equation}
where we use $\mathop{\bf E{}} \exp Z = \exp (\mu+\sigma^2/2)$ for $Z \sim \mathcal N(\mu,\sigma^2)$,
and the cumulant generating function
\begin{eqnarray}
\nonumber
K(w,t) &=& \log \mathop{\bf E{}}\exp (tR)\\
\nonumber &=& \log \left(
\sum_{i=1}^k\pi_i \exp \left( t\nu_i+\frac{t^2}{2}\sigma_i^2\right)\right)\\
&=&
\log \left(
\sum_{i=1}^k\pi_i \exp \left( t\mu_i^Tw+\frac{t^2}{2}w^T\Sigma_i w \right)
\right).
\label{e-cgf}
\end{eqnarray}
We observe for future use the identity
\begin{equation}\label{e-ident}
K(w,t)=K(tw,1),
\end{equation}
{\it i.e.}, the parameter $t$ simply multiplies the argument $w$.

\section{Portfolio optimization with exponential utility}\label{s-port-opt}

\subsection{Expected exponential utility}
Our objective is to choose $w$ to maximize the expected exponential utility
$\mathop{\bf E{}} U_\gamma(R)$, where
\[
U_\gamma(a) = 1-\exp(-\gamma a),
\]
with $\gamma>0$ the risk aversion parameter.
Using \eqref{e-mgf}, we can express this as
\[
\mathop{\bf E{}} U_\gamma(R) = 1-\mathop{\bf E{}} \exp(-\gamma w^Tr)
= 1-  M(w,-\gamma).
\]
It follows that we can maximize $\mathop{\bf E{}} U_\gamma(R)$ by minimizing
the moment generating function $M(w,-\gamma)$,
or equivalently the cumulant generating function
\begin{equation}\label{e-obj}
K(w,-\gamma) = \log \left( \sum_{i=1}^k\exp \left( \log\pi_i -\gamma\mu_i^Tw+
\frac{\gamma^2}{2}w^T\Sigma_i w \right) \right).
\end{equation}

\paragraph{Convexity.}
The function $K(w,-\gamma)$ is a convex function of $w$.
To see this, we note that for each $i$,
$\log \pi_i- \gamma \mu_i^Tw + \frac{\gamma^2}{2} x^T \Sigma_i w$ is a convex
quadratic function of $w$, and therefore convex.
The function $K(w,-\gamma)$ is the log-sum-exp function (also called the soft-max
function),
\begin{equation}\label{e-log-sum-exp}
S(u) = \log\left(\sum_{i=1}^k\exp u_i\right),
\end{equation}
of these arguments.
The log-sum-exp function is convex and increasing in all arguments, so
the composition $f_\gamma$ is convex \cite[\S3.1.5]{cvxbook}.

\subsection{EGM portfolio construction}
Our portfolio construction optimization problem has the form
\begin{equation}\label{e-prob}
\begin{array}{ll}
\mbox{minimize} &  K(w,-\gamma)\\
\mbox{subject to} & \mathbf 1^Tw =1, \quad w \in \mathcal W,
\end{array}
\end{equation}
where $\mathcal W$ is a convex set of portfolio constraints.
This is evidently a convex optimization problem.
One implication is that we can efficiently solve this problem globally using
a variety of methods.

\paragraph{DSL specification.}
The EGM problem \eqref{e-prob} is not just convex.
It is readily specified in domain-specific languages (DSLs) for convex
optimization, since all such systems include the log-sum-exp function,
and all such systems can handle the convex function composition rules
that we used to establish convexity of $K(w,-\gamma)$ in $w$.
No special methods (or gradient or other derivatives) are needed;
the function $K(w,-\gamma)$ can be specified in a DSL by just typing it
in as is.
As a simple example, CVXPY code for specifying
the EGM construction problem \eqref{e-prob} and solving it,
with a long only portfolio ({\it i.e.}, $w \geq 0$), is given below.
(This code snippet is also available at the repository
\url{https://github.com/cvxgrp/exp_util_gm_portfolio_opt}.)

\bigskip

\begin{lstlisting}
import cvxpy as cvx

def K(w):
  u = cvx.vstack([cvx.log(pi[i])
    - gamma * mus[i] @ w
    + (gamma**2/2) * cvx.quad_form(w, Sigmas[i])])
  return cvx.log_sum_exp(u)

w = cvx.Variable(n)
objective = cvx.Minimize(K(w))
constraints = [ w >= 0, cvx.sum(w) == 1 ]
egm_prob = cvx.Problem(objective, constraints)
egm_prob.solve()
w.value
\end{lstlisting}

\bigskip

Here it is assumed that \verb|n|, \verb|pi|, and \verb|gamma| are
constants corresponding to $n$, $\pi$ and $\gamma$,
and \verb|mus| and \verb|Sigmas| are lists
of the $\mu_i$ and $\Sigma_i$, respectively.
In lines 3--7 the objective $K(w,-\gamma)$ is formed,
and in lines 9--12 the EGM optimization problem is formed.
The problem is solved in line 13, which populates \verb|w.value| with optimal
weights.
In this simple example our portfolio constraint set $\mathcal W$ is simple.
One of the advantages of using a DSL is that more complex constraints
can be added by just appending them to the list of constraints defined in line~11.

\paragraph{Soft-max interpretation.}
We can give an interpretation of the objective $K(w,-\gamma)$ in \eqref{e-obj}
in terms of the soft-max function,
which can be thought of as a smooth approximation to the max, since it
satisfies
\begin{equation}\label{e-soft-max}
\max_i u_i  \leq S(u) \leq \max_i u_i  + \log k.
\end{equation}

The objective \eqref{e-obj} can be expressed as
\[
K(w,-\gamma) = S(u), \qquad
u_i = \log \pi_i + \gamma \left( -\mu_i^Tw+ \frac{\gamma}{2}w^T\Sigma_i w\right),
\quad i=1, \ldots, k.
\]
We recognize $-\mu_i^Tw + \frac{\gamma}{2} w^T \Sigma_i w$ as the negative
risk-adjusted return of the portfolio under the $i$th Gaussian component.
Thus $K(w,-\gamma)$ is the soft-max of these negative risk-adjusted returns,
offset by the terms $\log \pi_i$, and scaled by $\gamma$.
Roughly speaking, our objective is an approximation of the maximum of
the negative risk-adjusted returns under the component distributions.

From \eqref{e-soft-max} we have
\begin{eqnarray}
\label{e-K-lower-bnd}
K(w,-\gamma) &\geq&
\max_{i=1,\ldots, k} \left( \log \pi_i-\gamma \mu_i^Tw+ \frac{\gamma^2}{2}w^T\Sigma_i w \right),\\
\label{e-K-upper-bnd}
K(w,-\gamma) &\leq& \log k +
\max_{i=1,\ldots,k}
\left( \log \pi_i-\gamma \mu_i^Tw+ \frac{\gamma^2}{2}w^T\Sigma_i w \right).
\end{eqnarray}

\subsection{Special cases}\label{s-special}
\paragraph{Gaussian returns.}
When $k=1$ our GM return distribution reduces to Gaussian,
and the problem \eqref{e-prob} reduces to the standard Markowitz problem
\cite{Markowitz1952,mark-book}
\[
\begin{array}{ll}
\mbox{maximize} &  \mu_1^Tw - \frac{\gamma }{2} w^T\Sigma_1 w\\
\mbox{subject to} & \mathbf 1^Tw =1, \quad w \in \mathcal W.
\end{array}
\]

\paragraph{Finite values returns.}
When $\Sigma_i=0$, so $r$ takes on only the values $\mu_1, \ldots, \mu_k$,
the problem \eqref{e-prob} can be expressed as
\begin{equation}\label{e-prob-finite-vals}
\begin{array}{ll}
\mbox{minimize} & \log\left( \sum_{i=1}^k \pi_i \exp (- \gamma \mu_i^Tw) \right) \\
\mbox{subject to} & \mathbf 1^Tw =1, \quad w \in \mathcal W.
\end{array}
\end{equation}

\subsection{Simple example}
To illustrate the difference between EGM and mean-variance portfolios,
we construct a very simple example for which
the two portfolios can be analytically found.
We take a finite values distribution with $n=2$ assets and $k=2$ components,
with $\Sigma_1 = \Sigma_2=0$,
\[
\mu_1 = \left[\begin{array}{r} -1 \\ 0 \end{array}\right], \qquad
\mu_2 = \left[\begin{array}{r} 1 \\ 0 \end{array}\right].
\]
Thus $r=(-1,0)$ with probability $\pi_1$ and $r=(1,0)$ with probability $\pi_2$.
The first asset is risky, and the second is riskless, with zero return.
We take $\mathcal W= {\mbox{\bf R}}^2$, so the only constraint on
the portfolio weight is $w_1+w_2=1$.

\paragraph{Markowitz portfolio.}
The mean and covariance of $r$ are
\[
\mu = \left[\begin{array}{c} 1-2\pi_1 \\ 0 \end{array}\right], \qquad
\Sigma = \left[ \begin{array}{cc} 4\pi_1(1-\pi_1) & 0 \\ 0 & 0 \end{array}\right].
\]
The Markowitz optimal portfolio is
\[
w^\text{M}_1 = \frac{1-2\pi_1}{4 \gamma \pi_1(1-\pi_1)},
\]
with $w^\text{M}_2 = 1-w^\text{M}_1$.

\paragraph{EGM portfolio.}
The EGM portfolio minimizes
\[
\pi_1 \exp (\gamma w_1) + (1-\pi_1)\exp (-\gamma w_1),
\]
so the EGM portfolio is
\[
w^\text{E}_1 = \frac{\log(1/\pi_1 -1)}{2\gamma},
\]
with $w^\text{E}_2 = 1-w^\text{E}_1$.

\paragraph{Comparison.}
The two portfolios are the same for $\pi_1=1/2$,
with $w^\text{M} = w^\text{E} =(0,1)$.
They are not too far from each other for other values of $\pi_1$ and $\gamma$,
but can differ substantially for others.
For example with $\pi_1 =0.05$ and $\gamma=1$,
the Markowitz and EGM portfolios are
\[
w^\text{M}=(4.74,-3.74), \qquad
w^\text{E}=(1.47,-.47).
\]
The value at risk 5\% is 4.74 for the Markowitz porfolio compared
to 1.47 for the EGM portfolio.

\subsection{High and low risk aversion limits}
\paragraph{High risk aversion limit.}
Here we consider the case where $\gamma \to \infty$.
Dividing \eqref{e-K-lower-bnd} and
\eqref{e-K-upper-bnd} by $\gamma$, we find that
\[
\frac{K(w, -\gamma)}{\gamma} =
\max_{i=1, \ldots, k} \left(-\mu_i^T w + \frac{\gamma}{2} w^T
\Sigma_i w \right) + O(1/\gamma).
\]
So for large risk aversion parameter $\gamma$,
the EGM portfolio construction problem \eqref{e-prob} is approximately
\[
\begin{array}{ll}
\mbox{minimize} &  \max_{i=1, \ldots, k} \left( -\mu_i^T w + \frac{\gamma}{2}
w^T \Sigma_i w \right)\\
\mbox{subject to} & \mathbf 1^Tw =1, \quad w \in \mathcal W.
\end{array}
\]
Thus in the limit of high risk aversion,
the EGM portfolio minimizes the maximum of the risk adjusted returns
under each of the components, regardless of the $\pi_i$.
This is similar to solving a minimax Markowitz problem, where we
use as a risk model the maximum risk over a set of covariance matrices
(see \cite[\S4.2, p.\;30]{boyd-2017}).

\paragraph{Low risk aversion limit.}
Here we consider the case where $\gamma \to 0$.
We start with the well known expansion
\[
\frac{1}{\gamma} \log \left(\sum_{i=1}^k \pi_i \exp \gamma z_i \right)
= \sum_{i=1}^k \pi_i z_i + \frac{\gamma}{2} \left(
\sum_{i=1}^k \pi_i z_i^2 - \left( \sum_{i=1}^k \pi_i z_i\right)^2 \right)
+ O(\gamma^2),
\]
for any $z_i$. (We recognize the first term on the righthand side
as the mean of $z$, and the
second as $\gamma/2$ times the variance of $z$, when $z$ is a random variable
taking values $z_1, \ldots, z_k$ with probabilities $\pi_1, \ldots, \pi_k$.)
Substituting $z_i = -\mu_i^T w + \frac{\gamma}{2} w^T \Sigma_i w$ we obtain
\[
\frac{K(w,-\gamma)}{\gamma} = -\mu^T w + \frac{\gamma}{2} w^T \Sigma w + O(\gamma^2),
\]
where $\mu$ and $\Sigma$ are the mean and covariance of $r$,
\[
\mu = \mathop{\bf E{}} r = \sum_{i=1}^k \pi_i \mu_i, \qquad
\Sigma = \mathop{\bf E{}} rr^T -(\mathop{\bf E{}} r)(\mathop{\bf E{}} r)^T =
\sum_{i=1}^k \pi_i \left( \Sigma_i + (\mu_i - \mu)(\mu_i-\mu)^T \right).
\]
This has a very nice interpretation: in the limit of small risk aversion,
EGM reduces to Markowitz, using the mean and covariance of the return.

\section{Portfolio optimization with entropic value at risk} \label{s-EVaR}
\subsection{Entropic value at risk}
The traditional measure of downside risk
is the value at risk (VaR) with probability $\alpha$, which is the $(1-\alpha)$
quantile of the negative return $-R$,
\[
\text{VaR}_\alpha(R) = -\inf\{x\in{\mbox{\bf R}}\mid \mathop{\bf Prob}(R\leq x)>\alpha\}.
\]
(We are typically interested in values such as $\alpha=0.05$ or $\alpha = 0.01$.)
For example if the value at risk of a portfolio with probability $5\%$ is $15\%$,
the probability of a loss exceeding $15\%$ ({\it i.e.}, $R\leq -0.15$) is $5\%$.
Value at risk is interpretable and widely used,
but it is not a coherent risk measure \cite{Rockafellar-CVaR}.
For example, VaR is not sub-additive, so the sum of two
portfolios can have a higher VaR than the sum of the component VaRs.
Several coherent risk measures have been proposed, including the
conditional value at risk $\text{CVaR}_\alpha$ \cite{Rockafellar-CVaR}
and entropic value at risk $\text{EVaR}_\alpha$ \cite{EVaR-theory}.

The entropic value at risk $\text{EVaR}_\alpha$ is the tightest
Chernoff upper bound on $\text{VaR}_\alpha$, which can be expressed in
terms of the cumulant generating function as
\[
\text{EVaR}_{\alpha}(R) = \inf_{\lambda>0} \frac{K(w,-\lambda) -\log \alpha}{\lambda}
\geq \text{VaR}_\alpha(R)
\]
(It is also an upper bound on $\text{CVaR}_\alpha$.)
Cajas \cite{Cajas} and Shen et al.~\cite{shen2022minimizing}
describe convex optimization problems involving EVaR with the expectation replaced
with its sample approximation.

\paragraph{Minimum EVaR portfolio.}
To minimize $\text{EVaR}_\alpha(R)$, we solve the
optimization problem
\begin{equation}\label{e-EVaR-prob}
\begin{array}{ll}
\mbox{minimize} &  \frac{K(w,-\lambda) -\log \alpha}{\lambda} \\
\mbox{subject to} & \mathbf 1^Tw =1, \quad w \in \mathcal W, \quad \lambda>0,
\end{array}
\end{equation}
with variables $w\in {\mbox{\bf R}}^n$ and $\lambda \in {\mbox{\bf R}}$.

There is a close connection between this problem and the
exponential utility maximization problem \eqref{e-prob}.
Suppose that $w^\star$ and $\lambda^\star$ are optimal for
\eqref{e-EVaR-prob}.  Then $w^\star$ is also optimal for the
exponential utility problem \eqref{e-prob}, with risk aversion
parameter $\lambda^\star$.
Thus we can think of minimizing $\text{EVaR}_\alpha(R)$ as
simply choosing a value of the risk aversion parameter in EGM.
(This choice of parameter depends on $\alpha$.)
We will refer to portfolio construction using \eqref{e-EVaR-prob} also as
EGM, since any such portfolio is optimal for EGM with some
value of risk aversion, and also, conveniently, entropic and
exponential both start with E.

\subsection{Convex formulation}
The problem \eqref{e-EVaR-prob} is not a convex optimization problem
since the objective is not jointly convex in $w$ and $\lambda$.
But a change of variable can give us an equivalent convex problem.
Instead of using the variable $\lambda$, we use
the new variable $\delta = 1/\lambda$, and the problem \eqref{e-EVaR-prob}
becomes
\begin{equation}\label{e-EVaR-prob-cvx}
\begin{array}{ll}
\mbox{minimize} &  \delta K(w/\delta,-1) - \delta \log \alpha \\
\mbox{subject to} & \mathbf 1^Tw =1, \quad w \in \mathcal W, \quad \delta>0,
\end{array}
\end{equation}
with variables $w\in {\mbox{\bf R}}^n$ and $\delta \in {\mbox{\bf R}}$.
(We use the identity \eqref{e-ident} above.)
This objective is jointly convex in the variables $w$ and $\delta$, since it
is the perspective function of $K(w,-1)$, which is convex in $w$
\cite[\S 3.2.6]{cvxbook}, so \eqref{e-EVaR-prob-cvx} is a convex
optimization problem, which is readily solved.
(The constraint $\delta>0$ is actually not needed, since the perspective
function is defined to be $+\infty$ if $\delta \leq 0$.)

Unfortunately, current DSLs for convex optimization do not automate
the creation of the perspective of a function, so the
problem \eqref{e-EVaR-prob-cvx} cannot simply be typed in; we must form
the the perspective function by hand, as outlined below
in appendix~\ref{a-graph-form}.

There are also simple methods that can be used to solve it, with a modest loss
in efficiency, that are immediately compatible with DSLs.
One method is alternating optimization, where we alternate between fixing
$\delta$ and optimizing over $w$ (easy with current DSLs),
and fixing $w$ and optimizing over $\delta$ (minimization of
a scalar convex function, which can be done by many simple methods).
To start we can replace $K$ with the lower bound
\eqref{e-K-lower-bnd} (or the upper bound \eqref{e-K-upper-bnd}), to
obtain the approximate problem
\begin{equation}\label{e-EVaR-prob-approx}
\begin{array}{ll}
\mbox{minimize} & \max_i \left( \delta \log (\pi_i/\alpha) - \mu_i^T w +
\frac{w^T\Sigma_i w}{2 \delta} \right) \\
\mbox{subject to} & \mathbf 1^Tw =1, \quad w \in \mathcal W, \quad \delta>0.
\end{array}
\end{equation}
This problem is convex, and also immediately representable in DSLs using
the quadratic-over-linear function for the last term in the objective.
(Here too the constraint $\delta>0$ is redundant, since the
quadratic-over-linear function is defined as $+\infty$ if the
denominator is not positive.)

\subsection{Special cases}
\paragraph{Gaussian returns.}
When $k=1$, our GM return distribution is Gaussian and we have
\begin{equation}\label{e-EVaR-gaussian}
\delta K(w/\delta, -1)-\delta \log \alpha =
-\delta \log \alpha - \mu_1^Tw+ \frac{w^T\Sigma_1w}{2 \delta},
\end{equation}
with variables $w$ and $\gamma$.  This objective is readily minimized using DSLs,
using the quadratic-over-linear function for the last term.

The value of $\delta$ that minimizes this, with fixed $w$, is
\[
\delta = \left( \frac{w^T\Sigma_1 w}{-2 \log \alpha}\right)^{1/2}.
\]
Thus we see that for Gaussian returns, the portfolio that minimizes
$\text{EVaR}_\alpha$ is in fact Markowitz, with the specific
choice of risk aversion parameter
\[
\gamma = \left( \frac{-2 \log \alpha}{w^T\Sigma_1 w} \right)^{1/2}.
\]
(This depends on $w$, so to find it we must solve the convex problem
with objective \eqref{e-EVaR-gaussian}.)
We see that as $\alpha$ decreases, the associated risk aversion increases, which makes sense.

Substituting the optimal value of $\delta$ into the objective \eqref{e-EVaR-gaussian},
we find that the objective is
\[
- \mu_1^T w + (-2 \log \alpha)^{1/2} \left( w^T\Sigma_1 w \right)^{1/2},
\]
plus a constant.
Thus we maximize a risk adjusted return, using the standard
deviation instead of the traditional variance as risk,
and the very specific risk aversion constant $(-2\log \alpha)^{1/2}$.

\paragraph{Finite values returns.}
When $\Sigma_i=0$, so $r$ takes on only the values $\mu_1, \ldots, \mu_k$,
the problem \eqref{e-EVaR-prob-cvx} can be expressed as
\begin{equation}\label{e-EVaR-prob-finite-vals}
\begin{array}{ll}
\mbox{minimize} & \delta \log\left( \sum_{i=1}^k (\pi_i/\alpha)
\exp (-\mu_i^Tw/\delta) \right) \\
\mbox{subject to} & \mathbf 1^Tw =1, \quad w \in \mathcal W.
\end{array}
\end{equation}

\section{Conclusions}

In this paper we have shown that two specific portfolio construction problems,
maximizing expected exponential utility and minimizing entropic value at risk,
with a Gaussian mixture return model, can be
formulated as convex optimization problems, and exactly solved
with no need for return samples or Monte Carlo approximations.
The resulting problems
are not much harder to solve than a mean-variance problem, but have the advantage
of directly handling return distributions with substantial skews or large tails.

Our focus in this paper is on the formulation of the these portfolio
construction problems as tractable
convex optimization problems that do not need return samples.
In a future paper we will report on practical portfolio construction
using these methods.

\clearpage