Pith. sign in

REVIEW 3 major objections 4 minor 12 references

Joint modeling of low and high extremes using a multivariate extended generalized Pareto distribution

T0 review · 3 major / 4 minor · reviewed 2026-08-05 · deepseek-v4-flash

Pith's one-line read The paper proposes a multivariate extended generalized Pareto distribution that models the entire range of a random vector—low, moderate, and high—while keeping joint lower and upper tails compliant with extreme-value theory and governed by

desk verdict New threshold-free multivariate eGPD construction is a real contribution, but the application's neural credible intervals run outside the training prior and need fixing. read the letter →

arxiv 2509.05982 v1 pith:VSDCWEEP submitted 2025-09-07 stat.ME

classification stat.ME MSC 62G32
keywords multivariateextremesextendedgeneralizedParetodistributionloweranduppertaildependenceregularvariationthreshold-freeextremevaluemodelingneuralBayesestimationsimulation-basedinferencerainfall
verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

This paper tries to establish that a multivariate eGPD can describe the full range of a random vector—its bulk, low extremes, and high extremes—without a threshold selection step, while remaining faithful to extreme-value theory in both tails. The construction multiplies a univariate eGPD radial variable R by a dynamically weighted sum of two simplex-valued vectors: L, which controls lower-tail dependence, and U, which controls upper-tail dependence. As R becomes large the weight shifts to U, and as R becomes small it shifts to L, so one generative mechanism interpolates smoothly through the middle of the distribution. If the claim is right, environmental risk studies of phenomena such as droughts and floods could use a single fitted model instead of gluing together separate bulk and tail models, and fast simulation makes inference practical via pretrained neural networks. Applied to daily rainfall in the Netherlands, the model reproduces observed dependence from low to high quantiles and matches classical bivariate GP fits in both tails.

What carries the argument

The object that carries the argument is the stochastic construction Y = R([1 − ω{F_R(R)}]L + ω{F_R(R)}U). R is a univariate eGPD with lower-tail shape κ and upper-tail shape ξ; L and U are positive random vectors on the L1 simplex, independent of R; ω is a cumulative distribution function on [0,1] acting as a smooth switch between L and U as the radius R moves from 0 to ∞. The identity ||Y||_1 = R makes the radius exactly the univariate eGPD, and the tail limits ω→1 and ω→0 transfer the angular dependence to U and L respectively. This is what breaks the dependence structure into two separately controllable pieces while preserving a smooth bulk.

What would settle it

Simulate from the model with known parameters and a weight function that approaches 1 slowly, such as ω(u) = 1 − (1 − u)^β with small β, then compute the multivariate angular measure or extremal coefficient at very high radial thresholds; if the observed tail dependence differs systematically from U's angular distribution, the claimed reduction of the upper tail to RU has failed. Analytically, check whether E[{R(1 − ω(F_R(R)))}^α L_j^α] is negligible relative to E[{R ω(F_R(R))}^α U_j^α] for α = 1/ξ; if the ratio does not vanish as the tail threshold increases, the separation does not hold.

Watch

Extended reading notes

Core claim

The central construction is Y = R([1 − ω{F_R(R)}]L + ω{F_R(R)}U), with R a univariate eGPD, L and U independent simplex-valued random vectors, and ω a cumulative distribution function on [0,1]. Because ||Y||_1 = R, the radius is exactly the univariate eGPD. As the radius tends to infinity, ω approaches 1 and Y behaves like RU; as it tends to zero, ω approaches 0 and Y behaves like RL. Using Breiman's lemma and its multivariate extension, the paper shows the upper joint tail is regularly varying with index 1/ξ and limiting angular measure determined by U, while the lower joint tail has index κ (or, in a special case, the index of L^{-1}) with angular measure determined by L. Margins are heavy

Load-bearing premise

The paper asserts that in the far upper tail the term [1 − ω{F_R(R)}]L becomes negligible compared with ω{F_R(R)}U, so Y behaves like RU, and symmetrically at low values Y behaves like RL; if the weight function does not reach its limits fast enough relative to R's tail, that separation of tail dependence into U and L may fail.

Editorial extensions

If this is right

  • A single fitted model replaces the usual three-piece procedure (lower-tail GP, bulk, upper-tail GP) and avoids the threshold selection step in multivariate settings, where such choices are especially fragile.
  • Joint exceedance probabilities in both tails are EVT-compliant: high-level joint probabilities follow a regular variation limit whose angular measure is fixed by U, and low-level joint probabilities follow the analogous limit set by L.
  • Because both limiting angular measures are controlled by free parameters in L and U, practitioners can model drought-type joint lows and flood-type joint highs with different strengths of dependence in one model.
  • Fast simulation makes simulation-based inference practical; the paper's neural estimators provide point estimates or posteriors orders of magnitude faster than the hybrid likelihood-moment benchmark.
  • For rainfall, fits from the model reproduce observed χ-measures and quantile-quantile behavior from the lowest to highest values, matching classical bivariate GP fits in the tails while also covering the bulk.

Reading between the lines

Editorial extensions of the paper, not claims the author makes directly.

  • Because L and U can be any positive simplex-valued vectors, the construction offers a general template for full-range dependence modeling; a natural next step is to choose L or U so that one or both tails exhibit asymptotic independence, a case the paper explicitly leaves open.
  • The radial-angular separation is only asymptotic: for moderate values of R, both L and U contribute to dependence. If a user needs exact control at intermediate quantiles, the weight function ω would need to be chosen to shape that transition, and inference on θ_ω may be the weakest link (the paper's simulations show lower posterior contraction for θ_ω).
  • The paper notes that L and U contaminate marginal scales; when sites are not exchangeable, the component-specific scale extension it sketches would be needed. A testable extension is to fit that extended model and check whether the estimated L/U angular measures shift materially.
Share X Bluesky LinkedIn Reddit HN

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Referee Report

3 major / 4 minor

Summary. The paper proposes a new multivariate extended generalized Pareto distribution (eGPD) for jointly modeling low, moderate, and high values, with the stated aim of avoiding threshold selection while remaining compliant with multivariate extreme-value theory in both tails. The model is defined through the stochastic representation Y = R([1-ω{F_R(R)}]L + ω{F_R(R)}U), where R is a univariate eGPD radial variable, L and U are independent simplex-valued angular vectors, and ω is a distribution function on [0,1] that transitions from L to U as R moves from small to large values. The authors claim that the upper joint tail is governed by U, the lower joint tail by L, and that the margins behave like univariate eGPDs. Inference is performed using amortized neural Bayes estimators and neural posterior estimators, after also developing a hybrid likelihood--moment estimator. The methodology is illustrated on bivariate Dutch daily rainfall data from three stations. The central contribution is a generative full-range multivariate model with separately controlled lower- and upper-tail dependence.

Significance. If the mathematical claims hold, this is a useful contribution to the growing literature on bulk--tail models: it provides a simple, fast-to-simulate multivariate construction that avoids threshold selection and separates lower- and upper-tail dependence parameters. The authors are to be credited for making the code publicly available, for using amortized simulation-based inference, and for providing a fairly extensive simulation study with calibration diagnostics. The paper also honestly lists several limitations of the modeling framework, including the fact that L and U affect margins and that asymptotic independence is not covered. The main value is conceptual: a clean stochastic representation whose upper and lower extremal dependence can be tuned independently. However, as detailed below, the current manuscript contains a definitional ambiguity in the weight function, a missing proof of the central tail-equivalence assertion, and, most importantly, an application whose reported neural posterior estimates lie outside the training prior support and therefore cannot be taken at face value.

major comments (3)
  1. [Table 2 and §4.1] The NPE was trained with θ_U ∼ Unif(0.1,20) and θ_ω ∈ (0,0.5), but Table 2 reports for Ammerzoden–Zaltbommel a posterior median θ_U = 21.043 with 95% CI (19.358, 22.980), and for Giersbergen–Zaltbommel a θ_ω CI (0.417, 0.504) whose upper endpoint exceeds the parameter-space bound 0.5. Under the stated prior, the exact posterior has zero mass outside (0.1,20) and (0,0.5), respectively. These results therefore indicate that the normalizing flow is extrapolating outside its training support, so the reported posterior intervals and the conclusions drawn from them are not currently valid. The authors should either widen the priors (and retrain) or explicitly constrain the flow to respect the parameter bounds, and then reassess the application.
  2. [§4.1] The weight function is not well-defined as written. The text defines ω(u) as a Beta(3,3) CDF evaluated at (u−θ_ω)/(1−2θ_ω) with θ_ω ∈ (0,0.5). For u∈[0,1], this argument takes values below 0 (e.g., at u=0, it equals −θ_ω/(1−2θ_ω)) and above 1 (e.g., at u=1, it equals (1−θ_ω)/(1−2θ_ω) > 1). Since the Beta CDF is only defined on [0,1], the model (5) is not fully specified. The authors need to clarify whether the argument is clipped, renormalized, or replaced by a different transition function, and verify that the claimed limits ω(u)→0 as u→0 and ω(u)→1 as u→1 still hold.
  3. [§2.3.2] The central tail-separation result is asserted rather than proved. The text states that the upper tail of Y can be studied 'equivalently' through RU 'since ω(u)↑1 as u↑1', and similarly for the lower tail. This is a load-bearing step: it is needed to conclude that the joint upper (lower) tail of Y is fully governed by U (L). The convergence ω{F_R(R)}→1 in probability as R→∞ is plausible and likely sufficient, but the paper should supply a precise argument, e.g., a lemma showing that the remainder R(1−ω{F_R(R)})L (and the analogous lower-tail term) is asymptotically negligible in the relevant tail probability. Without such a proof, the main mathematical claim rests on an unverified heuristic.
minor comments (4)
  1. [§4.1] The variance-stabilizing transformation is written as 'h(x) = sign(x) log(1+|x|)−1', which is ambiguous. Please clarify the intended formula.
  2. [§2.3.2] Equation (6) and (7): it would be helpful to state explicitly that the moment formulas are only valid when the relevant moments of R exist, and to define ω_R before using it in (6).
  3. [§4.2] The SBC and z-score diagnostics are computed on test parameters drawn from the training prior. This is appropriate for checking calibration within the prior, but it does not validate behavior outside the prior support; a brief sentence acknowledging this limitation would be useful.
  4. [§5] The plots in Figure 5 and Appendix A are informative, but the lower-tail χ-measure appears to show systematic discrepancy at small levels; the caption mentions discretization of rainfall, but the paper would benefit from a more quantitative discussion of the magnitude of this lack-of-fit.

Circularity Check

0 steps flagged · score 0.0 of 10

No significant circularity: the model's tail properties follow from its explicit stochastic construction and are validated on independent simulations; the NPE prior-support issue is a validity concern, not circularity.

full rationale

The derivation chain is self-contained. The tail properties in Section 2.3.2 are consequences of the explicit stochastic representation (5), the exact identity ||Y||_1 = R in (8), and the assumed limits of the weight function ω; no fitted parameter is renamed as a prediction. The upper-tail equivalence with RU is asserted without a full proof, but it is a standard regular-variation/Slutsky consequence and does not assume the conclusion, so it is a mathematical gap rather than circularity. The neural estimators are trained on simulated data and validated on a separate test set (Section 4.2) with simulation-based calibration diagnostics; the application uses the trained NPE to produce posterior estimates, not predictions forced by a fitted input. The only caveat is non-circular: Table 2 reports credible intervals for θ_U and θ_ω that partially lie outside the training prior support (θ_U ~ Unif(0.1,20), θ_ω ~ Unif(0,0.5)), which undermines the validity of those NPE estimates but does not constitute circularity. Minor self-citations (Naveau et al. 2016; Sainsbury-Dale et al. 2024) provide background methods and are not load-bearing for the new model's mathematical construction.

Assumptions & free parameters 6 free parameters · 4 assumptions · 1 invented entities

The paper's central claims rest on six estimable parameters, two latent simplex vectors, and the unproved tail-equivalence step. The model itself postulates a specific generative mechanism; the EVT compliance of the tails is then derived, not assumed.

free parameters (6)
  • kappa (lower tail shape of R) = 1.116 (Ammerzoden-Zaltbommel)
    Estimated by NPE; controls lower tail index and interacts with the lower angular vector L.
  • sigma (scale of R) = 1.381 (Ammerzoden-Zaltbommel)
    Overall scale parameter of the eGPD radial variable; estimated from data.
  • xi (upper tail shape of R) = 0.193 (Ammerzoden-Zaltbommel)
    Controls the heavy upper tail index of all margins; estimated by NPE.
  • theta_L (lower dependence parameter) = 4.082 (Ammerzoden-Zaltbommel)
    Shape parameter of the Beta distribution for L; controls lower tail dependence.
  • theta_U (upper dependence parameter) = 21.043 (Ammerzoden-Zaltbommel)
    Shape parameter of the Beta distribution for U; controls upper tail dependence. Posterior median lies outside the training prior support.
  • theta_omega (bulk transition parameter) = 0.237 (Ammerzoden-Zaltbommel)
    Controls the transition region between lower and upper tail regimes in the weight function.
assumptions (4)
  • domain assumption R, L, and U are mutually independent in (5).
    Stated in Section 2.2; the tail characterizations via Breiman's lemma require it.
  • domain assumption Asymptotic dependence in both upper and lower tails.
    Stated in Section 1 and as a limitation in Section 6; the model is only EVT-compliant for asymptotically dependent tails.
  • ad hoc to paper Tail equivalence: Y behaves as RU in the upper tail and as RL in the lower tail.
    Asserted in Section 2.3.2 ('since omega(u) up arrow 1 as u up arrow 1') without a regular-variation proof of the remainder's negligibility.
  • standard math Breiman's lemma and the multivariate extension of Fougeres and Mercadier (2012).
    Used to derive the tail limit measures in Section 2.3.2.
invented entities (1)
  • Latent angular vectors L and U (on the L1 simplex)
    purpose: Control lower and upper tail dependence separately; U enters the upper angular limit and L the lower angular limit.
    They are unobserved latent variables internal to the model; no independent falsifiable handle is provided, though their distributions affect observable dependence.

how reviews work

0 comments
Cite this review

Pith. "Pith review of Joint modeling of low and high extremes using a multivariate extended generalized Pareto distribution." pith.science (2026). https://pith.science/paper/VSDCWEEP

@misc{pith2026250905982,
  author       = {Pith},
  title        = {Pith review of: Joint modeling of low and high extremes using a multivariate extended generalized Pareto distribution},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/VSDCWEEP}},
  note         = {Machine review of arXiv:2509.05982}
}
read the original abstract

In most risk assessment studies, it is important to accurately capture the entire distribution of the multivariate random vector of interest from low to high values. For example, in climate sciences, low precipitation events may lead to droughts, while heavy rainfall may generate large floods, and both of these extreme scenarios can have major impacts on the safety of people and infrastructure, as well as agricultural or other economic sectors. In the univariate case, the extended generalized Pareto distribution (eGPD) was specifically developed to accurately model low, moderate, and high precipitation intensities, while bypassing the threshold selection procedure usually conducted in extreme-value analyses. In this work, we extend this approach to the multivariate case. The proposed multivariate eGPD has the following appealing properties: (1) its marginal distributions behave like univariate eGPDs; (2) its lower and upper joint tails comply with multivariate extreme-value theory, with key parameters separately controlling dependence in each joint tail; and (3) the model allows for fast simulation and is thus amenable to simulation-based inference. We propose estimating model parameters by leveraging modern neural approaches, where a neural network, once trained, can provide point estimates, credible intervals, or full posterior approximations in a fraction of a second. Our new methodology is illustrated by application to daily rainfall times series data from the Netherlands. The proposed model is shown to provide satisfactory marginal and dependence fits from low to high quantiles.

Figures

Figures reproduced from arXiv: 2509.05982 by the authors.

Figure 1
Figure 1. Sample of size n = 2000 from the multivariate eGPD model (5) with the parametric spec￾ification specified in Section 4.1, with θ = (κ, σ, ξ, θL, θU , θω) ⊤ set to (3, 1, 0.05, 10, 20, 0.25)⊤ (top), (0.3, 1, 0.05, 10, 0.5, 0.25)⊤ (middle row), and (3, 1, 0.2, 4, 0.5, 0.25)⊤ (bottom), plotted on the original marginal scale (left), or rank-transformed to the Unif(0, 1) scale (middle column) with a corresponding bivaria… view at source ↗
Figure 2
Figure 2. Recovery plots for the likelihood-moment-based estimator (top), NBEs (middle), and posterior [PITH_FULL_IMAGE:figures/full_fig_p024_2.png] view at source ↗
Figure 3
Figure 3. Diagnostic plots assessing the quality of the approximate posterior distributions obtained from the [PITH_FULL_IMAGE:figures/full_fig_p025_3.png] view at source ↗
Figures from the paper (3 more)
Figure 4
Figure 4. Figure 4: Pairwise scatterplots of fall–winter daily scaled precipitation for the three station pairs: Ammer [PITH_FULL_IMAGE:figures/full_fig_p026_4.png]
Figure 5
Figure 5. Figure 5: Diagnostics for the Ammerzoden–Zaltbommel fit. Top: Real (red) and simulated (blue) scatterplot based on the estimated parameters (left), and the χ-dependence measures for the upper (middle) and lower (right) tails, computed based on observed (red) and simulated (blue)…
Figure 6
Figure 6. Figure 6: QQ-plots for marginal (left and middle) and joint (right; summarized by the sum across stations) distributions, comparing classical bivariate GP models fitted to data from the upper tail (top) or lower tail (bottom), versus our multivariate eGPD model fitted to the ent…

Discussion (0). Sign in to comment.

Reference graph

Works this paper leans on

12 extracted references · 6 canonical work pages

  1. [2]

    ArXiv preprint arXiv:2503.06255

    Gaussian mixture copulas for flexible dependence modelling in the body and tails of joint distributions. ArXiv preprint arXiv:2503.06255. André, L., Wadsworth, J., O’Hagan, A.,

  2. [4]

    Implications of Artificial Intelligence on Health Data Privacy and Confidentiality

    Fast likelihood-free parameter estimation for Lévy processes. ArXiv preprint arXiv.2501.01639. Davison, A.C., Huser, R.,

  3. [6]

    ArXiv preprint arXiv.2404.10899

    A variational neural Bayes framework for inference on intractable posterior distributions. ArXiv preprint arXiv.2404.10899. Mackay, E., Jonathan, P.,

  4. [7]

    ArXiv preprint arXiv:2310.12711

    Modelling multivariate extremes through angular-radial decomposition of the density function. ArXiv preprint arXiv:2310.12711. McFadden, D.,

  5. [8]

    ArXiv:2506.18508

    Theoretical guarantees for neural estimators in parametric statistics. ArXiv:2506.18508. Sainsbury-Dale, M.,

  6. [9]

    R package version 0.1.2,https://CRAN.R-project.org/package= NeuralEstimators

    NeuralEstimators: Likelihood-Free Parameter Estimation us- ing Neural Networks. R package version 0.1.2,https://CRAN.R-project.org/package= NeuralEstimators. Sainsbury-Dale, M., Richards, J., Zammit-Mangion, A., Huser, R., 2025a. Neural Bayes es- timators for irregular spatial data using graph neural networks. Journal of Computational and Graphical Statis...

  7. [10]

    ArXiv preprint arXiv:2412.07957

    Spatial scale-aware tail dependence modeling for high-dimensional spatial extremes. ArXiv preprint arXiv:2412.07957. Smith, R.L.,

  8. [2018]

    ArXiv preprint arXiv:1804.06788

    Validating Bayesian inference algorithms with simulation-based calibration. ArXiv preprint arXiv:1804.06788. Tencaliec, P., Favre, A.C., Naveau, P., Prieur, C., Nicolet, G.,

Show all 12 references
  1. [2019]

    secondary

    Natural catastrophes and man-made disasters in 2018: “secondary” perils on the frontline. (Tech. Rep. No. 2). Säilynoja, T., Bürkner, P., Vehtari, A.,

  2. [2023]

    (Eds.), Contribution of Working Groups I, II, and III to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change

    Climate Change 2023: Synthesis Report, in: Core Writing Team, H.L., Romero, J. (Eds.), Contribution of Working Groups I, II, and III to the Sixth Assessment Report of the Intergovernmental Panel on Climate Change. IPCC, Geneva, Switzerland, pp. 35–115. 33 Kakampakou, L., Simps...

  3. [2024]

    Computational Statistics & Data Analysis 189, 107841

    Joint modelling of the body and tail of bivariate data. Computational Statistics & Data Analysis 189, 107841. Ardizzone, L., Kruse, J., Rother, C., Köthe, U., 2019a. Analyzing inverse problems with invertibleneuralnetworks, in: Proceedingsofthe7thInternationalConferenceonLearn...

  4. [2025]

    ArXiv preprint arXiv:2503.23156

    Neural Bayes inference for complex bivariate extremal dependence models. ArXiv preprint arXiv:2503.23156. André, L.M., Tawn, J.A.,

Pith tools

Reviewed August 5, 2026 · model on record in the stance chip above.