---
title: "COMET Flows: Towards Generative Modeling of Multivariate Extremes and Tail Dependence"
authors:
  - "Andrew McDonald"
  - "Pang-Ning Tan"
  - "Lifeng Luo"
doi: "10.24963/ijcai.2022/462"
url: "https://ijcai.org/proceedings/2022/462"
pdf: "https://ijcai.org/proceedings/2022/0462.pdf"
venue: "Proceedings of the Thirty-First International Joint Conference on Artificial Intelligence"
year: 2022
pages: "3328-3334"
track: "Main Track"
keywords:
  - "Machine Learning: Applications"
  - "Multidisciplinary Topics and Applications: Physical Science"
  - "Data Mining: Mining Spatial and/or Temporal Data"
  - "Machine Learning: Probabilistic Machine Learning"
  - "Machine Learning: Unsupervised Learning"
---

# COMET Flows: Towards Generative Modeling of Multivariate Extremes and Tail Dependence

Andrew McDonald, Pang-Ning Tan, Lifeng Luo

*Proceedings of the Thirty-First International Joint Conference on Artificial Intelligence*, pages 3328-3334. DOI: [10.24963/ijcai.2022/462](https://doi.org/10.24963/ijcai.2022/462)

## Abstract

Normalizing flows—a popular class of deep generative models—often fail to represent extreme phenomena observed in real-world processes. In particular, existing normalizing flow architectures struggle to model multivariate extremes, characterized by heavy-tailed marginal distributions and asymmetric tail dependence among variables. In light of this shortcoming, we propose COMET (COpula Multivariate ExTreme) Flows, which decompose the process of modeling a joint distribution into two parts: (i) modeling its marginal distributions, and (ii) modeling its copula distribution. COMET Flows capture heavy-tailed marginal distributions by combining a parametric tail belief at extreme quantiles of the marginals with an empirical kernel density function at mid-quantiles. In addition, COMET Flows capture asymmetric tail dependence among multivariate extremes by viewing such dependence as inducing a low-dimensional manifold structure in feature space. Experimental results on both synthetic and real-world datasets demonstrate the effectiveness of COMET flows in capturing both heavy-tailed marginals and asymmetric tail dependence compared to other state-of-the-art baseline architectures. All code is available at https://github.com/andrewmcdonald27/COMETFlows.

## Full text

## 1 Introduction

Modeling multivariate extremes involves two main challenges at odds with the usual assumptions of independence and Gaussianity in learning theory: (i) the challenge of modeling heavy-tailed marginal distributions [Resnick, 2007], and (ii) the challenge of modeling tail dependence within the joint distribution [Joe, 2014]. In this work, we propose COMET (COpula Multivariate ExTreme) Flows-a novel deep generative modeling architecture combining Extreme Value Theory and Copula Theory with normalizing flowsin an effort to more accurately model multivariate extremes exhibiting tail dependence.

Deep generative models estimate an unknown probability distribution from an input data for the purposes of density estimation and sampling. Popularized by advances in generative adversarial networks [Goodfellow et al., 2014], variational autoencoders [Kingma and Welling, 2014], autoregressive networks and normalizing flows, deep generative models offer a promising path forward as limits to labeled data constrains progress in supervised learning. Among deep generative models, normalizing flows set the gold standard of tractability, allowing exact sampling and exact computation of log-density by modeling a transformation between data vectors x ∼ p(x) and latent vectors z ∼ p(z). Given sufficient data and architectural complexity, normalizing flows are universal density approximators [Kobyzev et al., 2020].

Figure 1: Vanilla normalizing flows such as RealNVP [Dinh et al., 2017] struggle to model heavy-tailed, tail-dependent multivariate probability distributions. COMET Flows succeed.

![Image](https://ijcai.org/proceedings/2022/fulltext/462/image_000000_b68382e3caa9424cb2f803d95901048ee647c2019bf1e5549ad14908aa8c201f.png)

Yet, normalizing flows still struggle to address challenges (i) and (ii) when modeling multivariate extremes (Figure 1). As discussed at length in [Wiese et al., 2019; Jaini et al., 2020], vanilla normalizing flows mapping heavy-tailed data distributions to light-tailed (e.g., Gaussian) latent distributions cannot be Lipschitz-bounded and hence become numerically unstable in practice. Likewise, viewing tail dependence among two or more dimensions as inducing a lowdimensional manifold structure in a region of feature space, we may reasonably expect vanilla normalizing flows to perform poorly in that region of space [Cunningham et al., 2020; Brehmer and Cranmer, 2020; Kim et al., 2020; Horvat and Pfister, 2021]. In this work, we address challenge (i) by employing a copula transformation to decouple heavy-tailed marginal distributions pi(xi) from the joint dependence structure in p(x), and address challenge (ii) by using a conditional perturbation of inputs during training to capture tail dependence in a stable manner.

The value of addressing challenges (i) and (ii) is apparent in a number of real-world applications, including climate science and meteorology (Figure 2). Extreme weather has become more frequent and severe in recent years due to climate change. Of particular concern is the projected increase in compound extreme events [Zscheischler et al., 2018], in which simultaneous multivariate extremes pose a greater risk to life and property than comparable univariate extremes spread over time. Accurate deep generative models of such compound extremes allow one to ask how likely an event is to occur, and may even be used to simulate extremes and their effects under various climate change scenarios [Puchko et al., 2019; Ayala et al., 2020; Klemmer et al., 2021].

Our main contributions are as follows:

1. We demonstrate that vanilla normalizing flows parameterized by Lipschitz-bounded neural networks are unable to capture heavy-tailedness and tail dependence in multivariate data;
2. We present COMET Flows, a novel architecture designed to capture asymmetric heavy-tailedness and tail dependence in multivariate data; and,
3. We empirically show that COMET Flows present a promising approach for modeling multivariate extremes using synthetic, climate, and benchmark data.

## 2 Preliminaries

Normalizing Flows Normalizing flows [Kobyzev et al., 2020; Papamakarios et al., 2021] offer a systematic approach to modeling multivariate probability distributions with deep learning. In particular, normalizing flows facilitate tractable sampling and density estimation under a d-dimensional probability distribution p(x) by parameterizing the transformation of each data vector x ∈ Rd sampled from p(x) to a latent representation z ∈ Rd following some simple probability distribution p(z) using an invertible neural network fθ : Rd → Rd . The network is trained to maximize the log-likelihood of observed data {xi}Ni=1 , which may be computed in closed form using the change-in-variables formula from calculus as a sum of the log-likelihood of the latent vector z = fθ(x) and the log-Jacobian determinant log |dfθ(x)/dx| of fθ:

<!-- formula-not-decoded -->

Following training, the normalizing flow fθ may be used for density estimation by evaluating p(x) under (1). Assuming fθ is invertible, we may sample from p(x) by sampling z ∼ p(z) and considering f -1 θ (z).

Architecturally, normalizing flows present two challenges: (i) they must be invertible in order to allow sampling in the generative (z → x) direction, and (ii) they must allow for efficient computation of the Jacobian determinant |df(x)/dx| in order to be trained and to allow density estimation in the normalizing (x → z) direction. The model we propose in this work satisfies both (i) and (ii) with a coupling-layer based flow backbone [Dinh et al., 2017].

Extreme Value Theory Extreme Value Theory [Coles, 2001; Beirlant et al., 2005] characterizes the behavior of extremes in precise statistical terms, and may be used to model the realizations of a random variable greater than some threshold. Indeed, under relatively mild assumptions [Coles 2001, Theorem 4.1] one may show that the excesses of a univariate random variable x ∼ p(x) over some threshold µ follows the generalized Pareto (GP) distribution with probability density function:

<!-- formula-not-decoded -->

Here, σ &gt; 0 is a scale parameter, ξ ∈ R is a tail heaviness parameter, and p(x) has support on x ∈ [µ,∞) for ξ ≥ 0; alternatively, on x ∈ [µ, µ - σ/ξ] for ξ &lt; 0.

Note that the GP distribution may be used to model left excesses (i.e., values under a threshold µ) as right excesses by transforming x → -x and µ → -µ. In this work, we utilize the GP distribution to model heavy-tailed marginals of the random vector x ∼ p(x) in both left and right tails.

Copulas Copulas [Nelsen, 2006; Joe, 2014] describe the dependence structure between uniform random variables, and in particular may be used to model the d-dimensional dependence between the components of a multivariate probability distribution after transforming the marginals to be uniformly distributed. Through the lens of Sklar's Theorem [Nelsen 2006, Theorem 2.3.3], a copula C(u) : [0, 1]d → [0, 1] is uniquely defined by a multivariate cumulative distribution function (CDF) F(x) : Rd → [0, 1] having univariate marginal CDFs Fi(xi) : R → [0, 1] through the relationship

<!-- formula-not-decoded -->

If F is continuous with well-defined marginal inverse CDFs F -1 i (xi), then (3) may be equivalently expressed as

<!-- formula-not-decoded -->

Taking partial derivatives of (3) with respect to each component gives us a similar expression relating p(x) : Rd → R to the copula density function c(u) : [0, 1]d → R which is often more convenient in practice when we are interested in density estimation:

<!-- formula-not-decoded -->

<!-- formula-not-decoded -->

In this work, we leverage copulas to model heavy-tailed marginals separately from the dependence structure of the components of a multivariate probability distribution.

Tail Dependence Tail dependence [Joe, 2014] characterizes the degree to which the behavior of two components of a multivariate probability distribution p(x) coincide at high (u → 1-) and low (u → 0+) quantiles. More formally, let p(x) be a multivariate probability distribution with joint CDF F andmarginal CDFs Fi, i = 1, . . . , d. Then the coefficients of upper and lower tail dependence λUij , λLij describing the tail relationship of components i and j are given by

Figure 2: Multivariate extremes in climate science and meteorology-such as TXx, TXn, prcptot, and r10mm in the CLIMDEX [Sillmann et al., 2013] dataset-often exhibit heavy-tails and tail dependence that vanilla normalizing flows struggle to model.

![Image](https://ijcai.org/proceedings/2022/fulltext/462/image_000001_13d334e504900ef50fa333f99ee26c8abde80cd07dd8517b032288942bf98acb.png)

<!-- formula-not-decoded -->

<!-- formula-not-decoded -->

As conditional probabilities, both λUij , λLij take values in [0, 1]. If λUij &gt; 0, we say that components i and j of p(x) exhibit upper tail dependence; if λLij &gt; 0, we say that components i and j of p(x) exhibit lower tail dependence.

In the extreme case of tail dependence between components xi and xj, it is possible that F -1 i (u) = xi = xj = F -1 j (u) for u &lt; a ≈ 0 or u &gt; b ≈ 1, such that there exists a local manifold structure of dimension Rd-1 embedded in Rd in the ij-th tail of p(x). More generally, if there exists such co-linearity between k components of x in the upper or lower tail of p(x), the support of p(x) may be a manifold of dimension Rd-k with measure zero in Rd . Such manifolds are difficult to model using normalizing flows given the requirement that the learned mapping fθ : Rd → Rd be invertible-in fact, if fθ successfully maps a manifold of dimension Rd-k ⊂ Rd to a latent space of Rd , it necessarily has an exploding Jacobian determinant. As a result, the problem of manifold learning with normalizing flows remains an active area of research [Cunningham et al., 2020; Brehmer and Cranmer, 2020; Kim et al., 2020; Horvat and Pfister, 2021]-yet to the best of the authors' knowledge, this work is the first to consider the connection between manifold learning, normalizing flows, and tail dependencies in multivariate probability distributions.

## 3 RelatedWork

Extremes and Normalizing Flows [Wiese et al., 2019; Jaini et al., 2020] prove the inability of normalizing flows to capture heavy-tailed marginal distributions, showing that any mapping of heavy-tailed distributions to light-tailed (e.g., Gaussian) distributions cannot be Lipschitz-bounded. To address this issue, [Wiese et al., 2019] propose Copula &amp; Marginal Flows (CM Flows) which disentangle heavy-tailed marginals from the joint distribution by employing a hybrid deep-parametric marginal flow mapping inputs xi ∈ R to uniform space ui ∈ [0, 1], followed by a RealNVP-inspired [Dinh et al., 2017] copula flow, mapping the (non-uniform) copula distribution to a latent uniform distribution. However, the authors omit the implementation of marginal flows in experiments, consider only synthetic bivariate data on [0, 1]d , and suggest the use of Deep Dense Sigmoidal Flows [Huang et al., 2018] for constructing marginal flows despite the fact these are not analytically invertible and hence do not permit easily tractable sampling.

Alternatively, [Jaini et al., 2020] propose Tail Adaptive Flows (TAFs) which replace the latent isotropic Gaussian p(z) with a latent isotropic Student's t distribution having learnable degrees-of-freedom parameter ν &gt; 1, arguing that mapping heavy tails of the input distribution p(x) to heavy tails of the latent distribution will maintain Lipschitzboundedness. However, the degrees of freedom parameter ν and the log-likelihood of the data under the model (which depend on one another) must be optimized simultaneously. Furthermore, by using an isotropic Student's t distribution in latent space, Tail Adaptive Flows implicitly assume symmetry in tail heaviness and hence are unable to capture asymmetric tail heaviness across dimensions.

Copulas and Normalizing Flows [Bossemeyer, 2020] develops a generalization of CM Flows to higher dimensions using vine copulas [Czado, 2019]. Orthogonally, [Kamthe et al., 2021] present Copula Flows motivated by the task of synthetic data generation, and utilize Neural Spine Flows [Durkan et al., 2019] to design a generative model able to handle both discrete and continuous data. In [Laszkiewicz et al., 2021], the authors propose to optimize a normalizing flow against a copula distribution in latent space instead of a standard isotropic Gaussian. Other relevant works at the intersection of generative modeling and copulas include [Ling et al., 2020], which considers deep estimation of Archimedean copulas, and [Tagasovska et al., 2019], in which vine copulas are fit upon the latent space of an autoencoder to construct a generative model.

Manifolds and Normalizing Flows A novel contribution of our work is to view tail dependence in multivariate probability distributions as inducing a local manifold structure in feature space. Recent works on manifold learning in normalizing flows entails density estimation on a d-dimensional manifold embedded in D &gt; d dimensional space. Noisy Injective Flows presented in [Cunningham et al., 2020] utilize stochastic inversion for sampling. M-Flows are developed in [Brehmer and Cranmer, 2020] to carry out simultaneous manifold learning and density estimation. SoftFlow [Kim et al., 2020] estimates a conditional distribution of data given a perturbation noise level as input to the flow. This idea is advanced by [Horvat and Pfister, 2021] who derive theoretical conditions under which it is possible to recover the underlying manifold using a noising inflation-deflation approach.

Figure 3: Architecture of COMET Flows, a normalizing flow architecture designed to capture asymmetric and heavy-tailed marginals in fm, along with asymmetric and manifold-like tail dependence structures in fc. The composition fc ◦ fm allows for tractable log-density estimation, while the composition of analytical inverses f -1 m ◦ f -1 c allows for tractable sampling.

![Image](https://ijcai.org/proceedings/2022/fulltext/462/image_000002_2f4ce677664ed764ba64e16e38c26a0813d1e31e07e47bdf750a66aebe701cc3.png)

## 4 Proposed Framework: COMET Flows

Vanilla normalizing flows struggle to capture (i) the heavytailed marginals and (ii) asymmetric tail dependencies in multivariate probability distributions. To address these weaknesses, we propose COMET (COpula Multivariate ExTreme) Flows, a normalizing flow architecture that decouples marginal distributions from the joint dependence structure within a multivariate probability distribution using a copula, models heavy-tailed marginals using a combination of kernel density and generalized Pareto (GP) distribution functions, and captures tail dependence by considering cases of λUij &gt; 0, λLij &gt; 0 as inducing a low-dimensional manifold in feature space.

COMET Flows has two main architectural components: a marginal transform fm : Rd → [0, 1]d mapping each marginal distribution pi(xi) to a uniform distribution such that the dependence structure among dimensions is expressed as a copula, and a copula transform fc : [0, 1]d → Rd mapping this copula distribution to a tractable latent distribution. Thus, COMET Flows define the normalizing flow (x → z) by the composition fc ◦ fm, and the generative flow (z → x) by the composition f -1 m ◦ f -1 c . We describe the details of each component in the following subsections, and present a schematic diagram of COMET Flows in Figure 3.

Marginal Transform fm In order to capture heavy-tailed marginal distributions pi(xi), we construct the marginal transform fm as a semi-parametric mixture in which the center of each marginal distribution is modeled using a 1dimensional Gaussian kernel density function, and the tails are modeled using a parametrically-fitted GP distribution motivated by the threshold exceedance paradigm of EVT.

More formally, let ai, bi ∈ [0, 1] be the lower and upper quantiles beyond which to apply parametric GP tails in dimension i, and let αi, βi ∈ R be the values in data space corresponding to ai, bi obtained through an inverse empirical CDF transform. Then fm : Rd → [0, 1]d is defined to be a decoupled mapping acting on each dimension i = 1, . . . , d of its input independently of other dimensions j ̸= i, with i-th component fm,i given by

<!-- formula-not-decoded -->

Here, Ai and Bi are the rescaled CDFs of parametric GP distributions fit to the data in the left and right tails of pi(xi), respectively, and Ki is the rescaled CDF of a 1D kernel density estimate fit to points in the center of pi(xi):

<!-- formula-not-decoded -->

<!-- formula-not-decoded -->

<!-- formula-not-decoded -->

Here, ˜ Ai, ˜Ki, and ˜Bi refer to the original CDFs fit to points xi &lt; αi, xi ∈ [αi, βi], and xi &gt; βi, respectively. By construction, the inverse marginal transform f -1 m,i : [0, 1] → R is well-defined by the rescaled inverse CDFs (A-1 i , K-1 i , B-1 i ) on (ui &lt; ai, ui ∈ [ai, bi], ui &gt; bi).

Note that the only learnable portion of fm,i are the parameters of the GP distributions used to model each tail, which can be fit using maximum likelihood estimation (e.g., using scipy.stats.genpareto.fit). This is by design, and vastly simplifies the learning dynamics of the end-to-end fc ◦ fm model: the GP parameters of each fm,i may be fit and fixed before fc is trained. This comes in contrast to the approaches of [Wiese et al., 2019; Bossemeyer, 2020], who allow each fm,i to be learnable flows but struggle to successfully train the coupled fc ◦ fm model.

Copula Transform fc Given a copula distribution in the unit hypercube [0, 1]d , we construct fc to be a normalizing flow with codomain Rd . To guarantee that the generative direction (z → u) maps vectors into the unit hypercube, we apply a logit transform as the first layer of fc such that the last layer in the generative direction comprises a sigmoid transform. As the logit transform has no learnable parameters, we use fc to refer only to the learnable transformation applied after the logit transform in what follows.

The learnable transformation fc is implemented as a coupling layer-based RealNVP [Dinh et al., 2017], repeatedly alternating dimensions 1:k and k + 1:d and applying

<!-- formula-not-decoded -->

<!-- formula-not-decoded -->

in each layer, where s, t : Rd → Rk , 1 &lt; k &lt; d, and ⊙ denotes elementwise multiplication. With Ik representing the k × k identity matrix and ∗ being any arbitrary (d - k) × k matrix, the Jacobian of each coupling layer takes the form

<!-- formula-not-decoded -->

and has a readily-computable log determinant ∑di=1 s(x1:d)i. Because this layer is invertible with tractable log-Jacobian determinant regardless of the form of functions s and t, we implement s and t as concatsquash [Onken and Ruthotto, 2020] networks conditioned on a contextual input σ:

<!-- formula-not-decoded -->

<!-- formula-not-decoded -->

To justify this construction, we recall that modeling multivariate extremes presents two main challenges: capturing (i) heavy-tailed marginal distribution and (ii) asymmetric tail dependencies. The marginal transform fm presented previously addresses (i), leaving the copula transform fc to address (ii). We accomplish this by viewing tail dependencies as inducing a manifold structure within subspaces of the unit hypercube, and draw inspiration from SoftFlow [Kim et al., 2020].

More precisely, we interpret the conditioning input σ to the concatsquash networks s, t to be a noise level, and perturb inputs x to s and t by x + ϵ where ϵ ∼ Nd(0, σId) for a noise level chosen uniformly at random σ ∼ U[0, σmax]. This noise level σ is passed to the flow as a conditioning contextual input through networks s and t throughout training, allowing the network to learn the conditional distribution p(x|σ) while sidestepping the exploding gradients introduced by the manifold structure of tail dependencies. Yet, this approach still allows the network to recover such tail dependencies at inference time by setting σ = 0.

## 5 Experimental Evaluation

## 5.1 Baselines

RealNVP A coupling layer-based RealNVP [Dinh et al., 2017] implementation optimized against an isotropic Gaussian latent distribution.

SoftFlow A modified SoftFlow [Kim et al., 2020] implementation in which the FFJORD [Grathwohl et al., 2019] backbone is replaced with a RealNVP backbone.

Tail Adaptive Flow (TAF) A modified TAF [Jaini et al., 2020] implementation composed of a RealNVP flow backbone with an isotropic Student's t latent distribution, parameterized by degrees-of-freedom parameter ν. We consider ν to be a fixed hyperparameter as no algorithm is presented to learn ν in the original work.

Copula Marginal Flow (CMF) A modified CMF [Wiese et al., 2019; Bossemeyer, 2020] implementation in which the unspecified tail distribution and non-analytically-invertible Deep Dense Sigmoidal Flows [Huang et al., 2018] used for fm are replaced with an analytically-invertible Gaussian KDE and generalized Pareto tail mixture model for fm.

Figure 4: Samples from models trained on the 8D synthetic dataset. Note that RealNVP and SoftFlow blur the manifold structure of tail dependence and exhibit light tails; TAFs generate noisy samples; and CMFs capture structure well but exhibit light tails in comparison to COMET Flows. Labels denote marginal (xi, xj), 1 ≤ i, j ≤ 8.

![Image](https://ijcai.org/proceedings/2022/fulltext/462/image_000003_a6a9fd808f8bce3a4a4e5324e12b94d7065589de23eb0e3512f9d38b8a65d823.png)

## 5.2 Data

Synthetic Data We design a synthetic dataset exhibiting heavy tails and asymmetric, heterogeneous tail dependence by sampling from a uniform distribution U on [0, 1]8 with probability 0.95, and from independent univariate GP marginal distributions X parameterized by (µ, σ, ξ) = (0, 1, 1) with probability 0.05. These samples are combined such that x1, x2 exhibit a heavy upper tail with tail dependence; x3, x4 exhibit a heavy lower tail with tail dependence; x5, x6 exhibit heavy upper and lower tails, each with tail dependence; and x7, x8 are perfectly collinear with heavy upper and lower tails. We generate 200,000 vectors in x ∈ R8 for training, 25,000 vectors for validation, and 25,000 for a heldout test set. Key bivariate marginals of this distribution are shown in the rightmost column of Figure 4.

Climate Data To evaluate the practical applicability of COMET Flows, we apply them to a subset of CLIMDEX 1 extreme climate indices [Sillmann et al., 2013]. Motivated by the projected increase in extreme temperatures and rainfall under climate change, we extract data from a model run under the RCP 8.5 scenario of continued high greenhouse gas emissions, and aim to model the joint distribution of such extremes over North America. We consider the region encompassing 25◦N-75◦N and 170◦W-50◦W and simulation data from the years 2091-2100, in which the effects of climate change on extremes will be most pronounced. We choose four indices of extreme temperature (TXx, TNx, TXn, TNn), along with four indices of extreme precipitation (prcptot, r1mm, r10mm, and r20mm), then add latitude and longitude to form a 10dimensional dataset characterizing joint extremes and their spatial distribution. Indices are reported at every 0.94◦ latitude and 1.25◦ longitude, giving a total of 97 × 53 = 5,141 vectors x ∈ R10 per year. We use data from 2091-2098 for training, 2099 for validation, and 2100 for testing. A subset of this data was presented in Figure 2.

[1 https://www.climdex.org/](https://www.climdex.org/)

Table 1: Comparison of baselines and proposed model performance, measured by average negative log-likelhood on a held-out test set. All models are trained until validation loss fails to improve for two consecutive epochs. Standard deviations in NLL over three training runs for each model were all &lt; ±0.1 and are omitted for space. Columns HT and TD denote whether each model is designed to handle heavy tails and tail dependence, respectively.

| Model - Name   | Model - Copula                         | Model - Latent                                     | Model - HT   | Model - TD   | NLL - Synthetic       | NLL - CLIMDEX        | NLL - POWER            | NLL - GAS               |
|----------------|----------------------------------------|----------------------------------------------------|--------------|--------------|-----------------------|----------------------|------------------------|-------------------------|
| RealNVP        | ✗                                      | Gaussian                                           | ✗            | ✗            | 7.483                 | 11.036               | 5.754                  | 7.425                   |
| SoftFlow       | ✗                                      | Gaussian                                           | ✗            | ✓            | 7.563                 | 10.765               | 5.790                  | 7.484                   |
| TAF            | ✗                                      | Student (ν = 16) Student (ν = 32) Student (ν = 64) | ✓            | ✗            | 7.046 7.169 7.483     | 10.679 10.606 10.568 | 5.500 5.615 5.690      | 6.790 7.099 7.252       |
| CMF            | (0.01, 0.99) (0.05, 0.95) (0.10, 0.90) | Gaussian                                           | ✓            | ✗            | -11.500 -9.816 -8.175 | -5.847 -3.072 -1.482 | -11.492 -10.518 -9.793 | -16.349 -15.354 -13.424 |
| COMET          | (0.01, 0.99) (0.05, 0.95) (0.10, 0.90) | Gaussian                                           | ✓            | ✓            | -9.592 -7.716 -6.023  | -6.729 -3.723 -2.300 | -10.955 -10.387 -9.643 | -15.944 -15.051 -13.116 |

Benchmark Data To contextualize the performance of COMET Flows, we consider the POWER and GAS benchmark datasets originally presented in [Papamakarios et al., 2017]. We follow the standard preprocessing steps provided in the MAF GitHub 2 repository. Note that we do not expect our baseline results to match state-of-the-art negative loglikelihood values reported in [Kobyzev et al., 2020] due to the simplified and adjusted architectures we implement.

## 5.3 Results

Synthetic Data As shown in Figure 4 and Table 1, the copulabased CM and COMET models dominate in both the normalizing (x → z) and generative (z → x) directions on synthetic data due to their functional form: by construction, they are naturally able to account for the asymetrically heavy-tailed margins of training data where vanilla normalizing flows cannot, and where TAF are forced to assume symmetry in tail heaviness. In particular, COMET's ability to capture the lowdimensional manifold structure in the upper tail region of feature space is highlighted by its qualitative performance in the generative direction (Figure 4), and its ability to avoid the instabilities and exploding gradients associated with manifold structures in feature space by learning perturbed manifolds conditional on noise is apparent by its quantitative performance in the normalizing direction (Figure 1). COMET's edge over CMF is most evident in its ability to capture tail dependence, though both significantly outperform TAF, SoftFlow, and RealNVP.

Climate Data COMET's advantage in the normalizing direction is more evident in this experiment. Key bivariate marginals in Figure 2 appear to exhibit local manifold structures (e.g., between tails of TXn and TNn), explaining COMET's superior performance. This observation suggests that explicitly accounting for such manifold structures in the joint tails of multivariate probability distributions is profitable in practice, particularly when redundancy across dimensions is likely to exist as an artifact of the latent physical processes generating the observed data.

[2 https://github.com/gpapamak/maf](https://github.com/gpapamak/maf)

Benchmark Data COMET is not as performant as CMF on the benchmark data in the normalizing direction, likely due to the absence of heavy-tails and tail dependence in POWER and GAS. Nevertheless, COMET outperforms TAF, SoftFlow, and RealNVP, suggesting that it remains competitive even on datasets it is not explicitly designed to handle.

## 6 Conclusion

Modeling multivariate extremes involves two fundamental challenges: (i) modeling heavy-tailed marginal distributions, and (ii) modeling tail dependence within the joint distribution. To address these challenges, we propose COMET Flows, a novel normalizing flow architecture which draws inspiration from Extreme Value Theory and Copula Theory, and which views tail dependence as inducing a local manifold structure in feature space. We demonstrate the efficacy of COMET Flows through experiments on synthetic, climate, and benchmark datasets, observing competitive performance in sampling and density estimation.

## Acknowledgements

This research is partially supported by the National Science Foundation under grant IIS-2006633. Any use of trade, firm, or product names is for descriptive purposes only and does not imply endorsement by the U.S. Government. The authors thank Asadullah Hill Galib and Tyler Wilson for their input in discussions regarding the paper.

## References

1. [Ayala et al., 2020] Alexis Ayala, Christopher Drazic, Brian Hutchinson, Ben Kravitz, and Claudia Tebaldi. Loosely Conditioned Emulation of Global Climate Models With Generative Adversarial Networks. In NeurIPS Work- shop on Tackling Climate Change with Machine Learning, 2020.
2. [Beirlant et al., 2005] Jan Beirlant, Yuri Goegebeur, Jozef Teugels, Johan Segers, Daniel De Waal, and Chris Ferro. Statistics of Extremes: Theory and Applications. Springer, 2005.
3. [Bossemeyer, 2020] Leonie Bossemeyer. CM Flows -Cop- ula Density Estimation with Normalizing Flows. Master's thesis, University of Edinburgh, 2020.
4. [Brehmer and Cranmer, 2020] Johann Brehmer and Kyle Cranmer. Flows for simultaneous manifold learning and density estimation. In NeurIPS, 2020.
5. [Coles, 2001] Stuart Coles. An Introduction to Statistical Modeling of Extreme Values. Springer, 2001.
6. [Cunningham et al., 2020] Edmond Cunningham, Renos Zabounidis, Abhinav Agrawal, Ina Fiterau, and Daniel Sheldon. Normalizing Flows Across Dimensions. In ICML Workshop on INNF, 2020.
7. [Czado, 2019] Claudia Czado. Analyzing Dependent Data with Vine Copulas. Springer, 2019.
8. [Dinh et al., 2017] Laurent Dinh, Jascha Sohl-Dickstein, and Samy Bengio. Density estimation using Real NVP. In ICLR, 2017.
9. [Durkan et al., 2019] Conor Durkan, Artur Bekasov, Iain Murray, and George Papamakarios. Neural spline flows. In NeurIPS, 2019.
10. [Goodfellow et al., 2014] Ian Goodfellow, Jean Pouget- Abadie, Mehdi Mirza, Bing Xu, David Warde-Farley, Sherjil Ozair, Aaron Courville, and Yoshua Bengio. Gen- erative adversarial networks. In NeurIPS, 2014.
11. [Grathwohl et al., 2019] Will Grathwohl, Ricky T.Q. Chen, Jesse Bettencourt, Ilya Sutskever, and David Duvenaud. FFJORD: Free-form continuous dynamics for scalable re- versible generative models. In ICLR, 2019.
12. [Horvat and Pfister, 2021] Christian Horvat and Jean-Pascal Pfister. Density estimation on low-dimensional manifolds: an inflation-deflation approach. arXiv:2105.12152, 2021.
13. [Huang et al., 2018] Chin-Wei Huang, David Krueger, Alexandre Lacoste, and Aaron Courville. Neural autoregressive flows. In ICML, 2018.
14. [Jaini et al., 2020] Priyank Jaini, Ivan Kobyzev, Yaoliang Yu, and Marcus A. Brubaker. Tails of lipschitz triangu- lar flows. In ICML, 2020.
15. [Joe, 2014] Harry Joe. Dependence Modeling with Copulas. CRC Press, 2014.
16. [Kamthe et al., 2021] Sanket Kamthe, Samuel Assefa, and Marc Deisenroth. Copula Flows for Synthetic Data Gen- eration. arXiv:2101.00598, 2021.
17. [Kim et al., 2020] Hyeongju Kim, Hyeonseung Lee, Woo Hyun Kang, Joun Yeop Lee, and Nam Soo Kim. SoftFlow: Probabilistic Framework for Normalizing Flow on Manifolds. In NeurIPS, 2020.
18. [Kingma and Welling, 2014] Diederik P. Kingma and Max Welling. Auto-encoding variational bayes. In ICLR, 2014.
19. [Klemmer et al., 2021] Konstantin Klemmer, Sudipan Saha, Matthias Kahl, Tianlin Xu, and Xiao Xiang Zhu. Gen- erative modeling of spatio-temporal weather patterns with extreme event conditioning. In ICLR Workshop on Model- ing Oceans and Climate Change, 2021.
20. [Kobyzev et al., 2020] Ivan Kobyzev, Simon Prince, and Marcus Brubaker. Normalizing Flows: An Introduction and Review of Current Methods. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2020.
21. [Laszkiewicz et al., 2021] Mike Laszkiewicz, Johannes Led- erer, and Asja Fischer. Copula-Based Normalizing Flows. In ICML Workshop on INNF, 2021.
22. [Ling et al., 2020] Chun Kai Ling, Fei Fang, and J. Zico Kolter. Deep Archimedean Copulas. In NeurIPS, 2020.
23. [Nelsen, 2006] Roger B. Nelsen. An Introduction to Copu- las. Springer, 2006.
24. [Onken and Ruthotto, 2020] Derek Onken and Lars Ruthotto. Discretize-optimize vs. optimize-discretize for time-series regression and continuous normalizing flows. arXiv:2005.13420, 2020.
25. [Papamakarios et al., 2017] George Papamakarios, Theo Pavlakou, and Iain Murray. Masked autoregressive flow for density estimation. In NeurIPS, 2017.
26. [Papamakarios et al., 2021] George Papamakarios, Eric Nal- isnick, Danilo Jimenez Rezende, Shakir Mohamed, and Balaji Lakshminarayanan. Normalizing Flows for Prob- abilistic Modeling and Inference. Journal of Machine Learning Research, 2021.
27. [Puchko et al., 2019] Alexandra Puchko, Robert Link, Brian Hutchinson, Ben Kravitz, and Abigail Snyder. Deep- ClimGAN: A High-Resolution Climate Data Generator. In NeurIPS Workshop on Tackling Climate Change with Ma- chine Learning, 2019.
28. [Resnick, 2007] Sidney I. Resnick. Heavy Tail Phenomena: Probabilistic and Statistical Modeling. Springer, 2007.
29. [Sillmann et al., 2013] Jana Sillmann, Viatcheslav V. Kharin, Francis W. Zwiers, Xuebin Zhang, and David Bronaugh. Climate extremes indices in the CMIP5 multimodel ensemble: Part 2. Future climate projections. Journal of Geophysical Research: Atmospheres, 2013.
30. [Tagasovska et al., 2019] Natasa Tagasovska, Damien Ack- erer, and Thibault Vatter. Copulas as high-dimensional generative models: Vine copula autoencoders. In NeurIPS, 2019.
31. [Wiese et al., 2019] Magnus Wiese, Robert Knobloch, and Ralf Korn. Copula & Marginal Flows: Disentangling the Marginal from its Joint. arXiv:1907.03361, 2019.
32. [Zscheischler et al., 2018] Jakob Zscheischler, Seth Westra, Bart JJM Van Den Hurk, Sonia I Seneviratne, Philip J Ward, Andy Pitman, Amir AghaKouchak, David N Bresch, Michael Leonard, Thomas Wahl, et al. Future cli- mate risk from compound events. Nature Climate Change, 2018.
