Pith. sign in

REVIEW 4 major objections 8 minor 31 references

Fast Pricing of Energy Derivatives with Mean-reverting Jump-diffusion Processes

T0 review · 4 major / 8 minor · reviewed 2026-08-14 · deepseek-v4-flash

Challenge this review Re-run · record.json Download PDF Read on arXiv ↗ Review JSON

Pith's one-line read Two random draws—a Polya count and an Erlang size—replace Poisson jump-time simulation in energy price models, cutting Monte Carlo pricing time by up to thirty times while keeping the simulation exact.

desk verdict The Polya-Erlang simulation algorithm is exact and fast, but the paper's risk-neutral drift h(t) has sign errors that break forward-curve consistency; the pricing sections need correction before the headline speed claims can be trusted. read the letter →

arxiv 1908.03137 v2 pith:6IHCB7GU submitted 2019-08-08 q-fin.CP math.PR

Nicola Cufaro Petroni , Piergiacomo Sabino This is my paper · ORCID
classification q-fin.CPmath.PR MSC 91G2091G6065C0560J75
keywords Gamma-OUprocessexactsimulationjump-diffusionenergyderivativesAsianoptionsgasstorageswingMonteCarlopricing
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

The paper claims that the slow part of Monte Carlo pricing for mean-reverting jump-diffusion energy models—simulating the Gamma-OU jump component—can be replaced by an exact two-draw procedure: sample a Polya (negative binomial) random variable, then an Erlang random variable with rate set by the jump parameters. This avoids simulating Poisson jump times and individual jump sizes altogether. On Asian option pricing the authors report that path simulation becomes at least thirty times faster, and for gas storages and swing options the total Monte Carlo pricing time falls by about forty percent even though the stochastic-optimization step dominates. The paper's point is that such speedups make real-time pricing of energy derivatives feasible without changing the modeled price distribution.

What carries the argument

The load-bearing object is the $a$-remainder of a self-decomposable gamma law: for $a=e^{-kt}$, the law of the Gamma-OU innovation is an infinite Polya-weighted mixture of Erlang laws, equivalently a single Erlang draw with a Polya-distributed random index. This identity converts a compound Poisson simulation—which requires generating jump times, sorting them, and summing discounted jump sizes—into two parameterized random draws. The same machinery, with $a^2$ replacing $a$, handles the symmetric bilateral Gamma-OU process used for Laplace jumps.

What would settle it

Simulate a Gamma-OU process with Algorithm 3 for a fixed parameter set and compare the empirical characteristic function (or a large-sample histogram) against the closed form $\phi(u,t)=((\beta-iue^{-kt})/(\beta-iu))^{\lambda/k}$ and against the classical jump-time simulation at very high path counts; a statistically significant mismatch would refute the claimed exactness.

Watch

Extended reading notes

Core claim

On its own terms, the paper's central claim is that the increment of a Gamma-OU process with mean reversion $k$, jump intensity $\lambda$, and exponential jump rate $\beta$ has the same law as an Erlang variable $E_b(\beta/a)$ whose index $b$ is a Polya random variable $B(\alpha,1-a)$, where $\alpha=\lambda/k$ and $a=e^{-kt}$. This is the $a$-remainder representation of the self-decomposable gamma law, and it appears as characteristic function $\phi(u,t)=((\beta-iue^{-kt})/(\beta-iu))^{\lambda/k}$. The paper builds a simulation algorithm on this identity, and extends it to bilateral Gamma-OU processes (differences of two Gamma-OU processes) including the symmetric Laplace-jump case, where two independent Erlang draws are subtracted. The numerical experiments compare this procedure with the classical jump-time simulation and a recently proposed randomized-rate simulation, and find equal prices and root-mean-square errors with far lower CPU times.

Load-bearing premise

Algorithm 3 is exact only if the characteristic-function and Polya-Erlang mixture identities from the authors' earlier paper [13] are correct; the present paper cites that result rather than deriving or verifying it, so a flaw there would invalidate the exactness claim.

Editorial extensions

If this is right

  • For Asian options, the path-generation CPU time with Algorithm 3 is roughly one thirtieth of the time needed by the jump-time algorithms at the same number of simulations, with indistinguishable prices and root-mean-square errors.
  • For gas storages and swing options, total Least-Squares Monte Carlo pricing time drops by about 40 percent, because the jump-path simulation, previously a dominant cost, becomes a small fraction of the stochastic-optimization step.
  • At 500,000 simulated paths, the Asian-option example falls from about 48 minutes to about 2 minutes of CPU time in the paper's experiments.
  • Time-dependent jump intensities are handled by setting a different Polya parameter at each time step, preserving exactness under stepwise-constant intensities.
  • The same two-draw construction covers the three market models studied: exponential jumps, asymmetric bilateral exponential jumps, and symmetric Laplace jumps.

Reading between the lines

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

  • Because Algorithm 3 bypasses jump-time sampling entirely, its advantage should widen as jump intensity or time-step size increases: the number of Erlang draws is fixed while the number of Poisson arrivals to be sorted grows. A benchmark across intensities would make this explicit.
  • The $a$-remainder mixture suggests a route toward backward (bridge) simulation of Gamma-OU paths, which the paper lists as future work; a Gamma-OU bridge would further speed Least-Squares Monte Carlo by removing the sequential forward path constraint.
  • The same self-decomposability argument may extend the two-draw construction to other OU processes driven by self-decomposable laws (for example tempered-stable or variance-gamma type), provided the $a$-remainder is known in closed form; the paper does not test this generalization.
  • Since the innovation is now a function of a single Polya index and one Erlang draw, common-random-number or quasi-Monte Carlo variance reduction might be applied more directly than with jump-time algorithms; the paper does not explore this.
Share X Bluesky LinkedIn Reddit HN

Signed reviews

No signed human review yet.

Editorial analysis

A structured set of objections, weighed in public.

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

Referee Report

4 major / 8 minor

Summary. This paper proposes exact and fast simulation methods for Gamma-OU and bilateral Gamma-OU processes, based on representing the innovation as a Polya-mixed Erlang variable, and applies them to Monte Carlo pricing of Asian options, gas storages, and swing options under three jump-diffusion market models. The authors compare Algorithm 3 (and Algorithm 6 for the symmetric bilateral case) with the standard Cont-Tankov arrival-by-arrival algorithm and with the Qu-Dassios-Zhao method, reporting large CPU-time gains.

Significance. The underlying simulation identity is mathematically sound: the characteristic function in Eq. (20) is correct for the Gamma-OU process, and the mixture representation (18)-(19) justifies Algorithm 3. I verified that Algorithm 6 also produces the correct symmetric bilateral Gamma-OU law. If the pricing drift and timing claims are corrected, the method is a genuinely useful contribution because it eliminates Poisson arrival simulation and sorting, and the speed advantage in path generation is real. The numerical pricing results are plausible and the convergence of prices across algorithms supports the exactness claim. The paper would be strengthened by making the proofs of (20) self-contained rather than citing the authors' own companion paper [13].

major comments (4)
  1. [§4.1–§4.3 and Eq. (24)] The printed risk-neutral drift h(t) is not martingale-consistent. From Eq. (20) with the decomposition Y=Y1-Y2, the characteristic function of the negative-jump component Y2 at u=i is ((β2+e^{-kt})/(β2+1))^{λ2/k}; therefore the second term in the §4.1 formula must be -λ2/k log((β2+e^{-kt})/(β2+1)), not the printed -λ2/k log((β2-e^{-kt})/(β2-1)). The same error appears in §4.2. In §4.3, Eq. (24) should read ((β²+u²e^{-2kt})/(β²+u²))^{λ/(2k)}; the printed minus signs make the expression not a characteristic function and are inconsistent with Eqs. (25)-(26), and the diffusion term in the §4.3 h(t_m) is printed with the wrong sign. As a result, a reader implementing the closed forms as written will not reproduce E[S(T)]=F(0,T); for the Case 1 parameters the bias in the forward is on the order of 5-10% at T=1.
  2. [§4.3, time-dependent intensity] The formula h(t_m) = ... - λ_m/(2k) log((β²-e^{-2kt_m})/(β²-1)) uses only the current intensity λ_m together with the total maturity t_m. For a stepwise time-dependent Poisson intensity, the deterministic adjustment is the sum over past increments of -log φ_increment(-i); it cannot be written with a single λ_m and total t_m unless λ is constant. Please provide the correct cumulative formula or explain the approximation and its impact on forward-curve consistency.
  3. [Tables 2, 5, and 8] The tables do not support the speed headlines. In Table 2, the CPU ratios Alg2/Alg3 range from about 21 to 28 and Alg4/Alg3 from about 26 to 34; the abstract's 'at least thirty times faster' and the conclusion's 'factor larger than forty' are not reproducible, since the largest ratio in the table is about 34. For storages and swings, the worst-case overall CPU ratios in Tables 5 and 8 are about 1.3-1.4 (for example Alg2/Alg3 at N=10^4 in Table 5), corresponding to roughly 25-30% time savings, not '40% faster in the worst case.' The direction of the advantage is clear, but the reported magnitudes should be recomputed or stated precisely.
  4. [Algorithm 1] Algorithm 1 as printed cannot be implemented as written. Line 2 uses Δt_m and φ_jumps(u,t_m) instead of the total time t_m and the argument -i, and line 8 omits the forward F(0,t_m); if taken literally, the simulated spot does not satisfy S(t)=F(0,t)e^{h(t)+H(t)}. Please rewrite the pseudocode so that h(t_m) is the closed-form total drift (with the corrected signs) and S(t_m) includes F(0,t_m).
minor comments (8)
  1. [Table 1] Table 1 appears misaligned: as printed, the header has nine entries and the data row has nine numbers, which would set p=50 and λ1=0.6; please reformat so p, λ1, λ2, β1, and β2 are unambiguous.
  2. [Algorithm 5] Algorithm 5 uses E1(λJ) for the jump-size exponentials; these should be E1(β) or a clearly defined parameter, not λJ.
  3. [Algorithm 7] Algorithm 7 has inconsistent notation: line 3 defines β^{(r)}_n but line 4 uses i=1,...,n and β^{(d)}_i; use N and consistent subscripts throughout.
  4. [Tables 2, 5, and 8] CPU times in Tables 2, 5, and 8 are single measurements; report averages over several runs with standard deviations or interquartile ranges.
  5. [Table 2 caption] The caption of Table 2 says RMSEs are 'divided by the squared root of the number of simulations'; the numbers look like standard errors, so please clarify the definition.
  6. [§3.1, reference [13]] The short derivation of Eq. (20) and the mixture (18)-(19) is delegated to the companion paper [13]; including a self-contained proof would avoid relying on an unpublished reference.
  7. [§4.1] The sentence 'almost one hour' for N=5×10^5 is only true for Algorithm 4 (3469 s); Algorithm 2 takes about 48 minutes, so the wording should be adjusted.
  8. [Figures 3b and 5b] Figures 3b and 5b do not say whether the plotted ratios are for LSMC or PATH times; the text discusses both, and the 'worst-case' statements are not visible in the figures.

Circularity Check

0 steps flagged · score 0.0 of 10

No significant circularity: the exact-simulation identity is a parameter-free mathematical result, and the speed advantages are empirical measurements.

full rationale

The paper's central methodological step is the representation of a Gamma-OU increment as a Polya-Erlang mixture, Eqs. (18)-(20), and the corresponding Algorithms 3 and 6. Although these formulas are imported from the authors' own companion paper [13], they are parameter-free mathematical identities with explicitly stated assumptions (alpha = lambda/k, a = e^{-kt}); they do not encode the prices, CPU times, or any fitted quantity of the present paper. The exactness of Algorithm 3 is therefore not a fitted input renamed as a prediction. The speed comparisons in Tables 2, 5, and 8 are empirical measurements of wall-clock time against independent baselines (Cont-Tankov Algorithm 6.2 and the Qu-Dassios-Zhao conditional-compound-Poisson representation), and the prices agree across the algorithms within Monte Carlo error. The h(t) formulas in Sections 4.1-4.3 are derived from the martingale condition (3); possible sign errors in the printed closed forms for downward jumps would be correctness defects, not circularity, because h(t) is not a fitted parameter used to produce the alleged prediction. Likewise, the overstatement that the speed-up reaches 'thirty times' and 'larger than forty' is a numerical reporting issue, not a definitional reduction. No circular step was found.

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

The central simulation algorithm is built on the self-decomposability of the gamma law and the a-remainder representation, sourced primarily from the authors' own prior publication [13]. The numerical experiments use hand-picked market parameters but no fitted values, and no new entities are introduced.

free parameters (1)
  • Market model parameters = Tables 1, 3, and 6
    The spot price parameters (mean reversion rates, jump intensities, jump size rates) are chosen by hand as realistic values, not calibrated to market data. The reported speedups are demonstrated only for these settings.
assumptions (4)
  • standard math Gamma law is self-decomposable, so for any a in (0,1) the decomposition X = aY + Z_a holds
    Standard result (Sato), used in Section 3.1 to derive the a-remainder mixture representation.
  • standard math The chf of a Gamma-OU process at time t is ((β - i u e^{-kt})/(β - i u))^{λ/k}
    Given in eq. (20) and follows from the definition of the OU process and self-decomposability; also derived in the authors' companion paper [13].
  • domain assumption Spot price can be decomposed as S(t) = F(0,t) exp(h(t)+X(t)+Y(t)) with independent factors and h(t) = -log φ_H(-i,t)
    Market modeling assumption using risk-neutral arguments from Hambly et al. (eqs. (1)-(4)).
  • domain assumption Time-dependent jump intensity is stepwise constant on the simulation grid
    Assumed in Section 3.2 to handle non-homogeneous Poisson processes; allows using different λ_m per time step.

how reviews work

0 comments
Cite this review

Pith. "Pith review of Fast Pricing of Energy Derivatives with Mean-reverting Jump-diffusion Processes." pith.science (2026). https://pith.science/paper/6IHCB7GU

@misc{pith2026190803137,
  author       = {Pith},
  title        = {Pith review of: Fast Pricing of Energy Derivatives with Mean-reverting Jump-diffusion Processes},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/6IHCB7GU}},
  note         = {Machine review of arXiv:1908.03137}
}
read the original abstract

Most energy and commodity markets exhibit mean-reversion and occasional distinctive price spikes, which results in demand for derivative products which protect the holder against high prices. To this end, in this paper we present exact and fast methodologies for the simulation of the spot price dynamics modeled as the exponential of the sum of an Ornstein-Uhlenbeck and an independent pure jump process, where the latter one is driven by a compound Poisson process with (bilateral) exponentially distributed jumps. These methodologies are finally applied to the pricing of Asian options, gas storages and swings under different combinations of jump-diffusion market models, and the apparent computational advantages of the proposed procedures are emphasized.

Figures

Figures reproduced from arXiv: 1908.03137 by the authors.

Figure 1
Figure 1. Asian options. 0 100 200 300 400 18 20 22 24 26 28 30 Price Trajectory: Market Model Case 1 (a) Sample Trajectory. 103 104 105 106 Number of Simulations 101 CPU time ratio Ratios of CPU Times Alg 2 / Alg 3 Alg 4 / Alg 3 Alg 4 / Alg 2 (b) Ratio of CPU times. pricing of the Asian option above, is accomplished in almost two minutes whereas, it takes almost one hour with the other alternatives. Figure 1b clearly shows t… view at source ↗
Figure 2
Figure 2. Gas Storages. 0 100 200 300 400 time 0 5 10 15 20 25 C(t) Feasible Volumes Storage (a) Feasible Volumes of the fast churn storage 0 100 200 300 400 16 18 20 22 24 26 28 30 Price Trajectory: Market Model Case 2 (b) Sample Trajectory [PITH_FULL_IMAGE:figures/full_fig_p015_2.png] view at source ↗
Figure 3
Figure 3. Gas Storage Results. 103 104 105 Price 281 282 283 284 285 286 287 288 Number of Simulations Gas Storage Value Alg 2 / Alg 3 Alg 4 / Alg 3 Alg 4 / Alg 2 (a) Values Gas Storage. 103 104 105 Number Of simulations 1 1.2 1.4 1.6 1.8 2 2.2 Ratio simulation times Ratios of CPU Times Alg 2 / Alg 3 Alg 4 / Alg 3 Alg 4 / Alg 2 (b) Ratio of CPU Times. 15 [PITH_FULL_IMAGE:figures/full_fig_p015_3.png] view at source ↗
Figures from the paper (2 more)
Figure 4
Figure 4. Figure 4: Market Model. 0 100 200 300 400 30 35 40 45 50 55 60 65 Time-dependent Intensity: Market Model Case 3 (a) Time-dependent Intensity 0 100 200 300 400 18 20 22 24 26 28 Price Trajectory: Market Model Case 3 (b) Sample Trajectory [PITH_FULL_IMAGE:figures/full_fig_p018_4.png]
Figure 5
Figure 5. Figure 5: Swings. 0 100 200 300 400 0 20 40 60 80 100 120 Feasible Volumes Swings (a) Feasible Volumes of a 120-120 Swing. 103 104 105 Number Of simulations 1 1.2 1.4 1.6 1.8 2 Ratio simulation times Ratios of CPU Times Alg 2 / Alg 3 Alg 4 / Alg 3 Alg 4 / Alg 2 (b) Swing Values …

Discussion (0). Continue with ORCID to comment.

Reference graph

Works this paper leans on

31 extracted references · 31 canonical work pages

  1. [13]

    Cufaro Petroni and P

    N. Cufaro Petroni and P. Sabino. Gamma Related OrnsteinUhlen beck Processes and their Simulation. available at: https://arxiv.org/abs/2003.0881 0, 2020

    work page arXiv 2003

  2. [1]

    Bardou, S

    O. Bardou, S. Bouthemy, and G. Pag´ es. Optimal Quantization f or the Pricing of Swing Options. Applied Mathematical Finance , 16(2):183–217, 2009

    work page 2009

  3. [2]

    Barndorff-Nielsen and N

    O.E. Barndorff-Nielsen and N. Shephard. Non-Gaussian Ornstein -Uhlenbeck- based Models and Some of Their Uses in Financial Economics. Journal of the Royal Statistical Society: Series B , 63(2):167–241, 2001

    work page 2001

  4. [3]

    Ben-Ameur, M

    H. Ben-Ameur, M. Breton, L. Karoui, and P. L’Ecuyer. A Dynamic Program- ming Approach for Pricing Options Embedded in Bonds. Journal of Economic Dynamics and Control , 31(7):2212–2233, July 2007

    work page 2007

  5. [4]

    Benth and A

    F.E. Benth and A. Pircalabu. A non-gaussian Ornstein-Uhlenbeck Model for Pricing Wind Power Futures. Applied Mathematical Finance , 25(1), 2018

    work page 2018

  6. [5]

    Benth and J

    F.E. Benth and J. Kallsen T. Meyer-Brandis. A non-Gaussian Orns tein- Uhlenbeck Process for Electricity Spot Price Modeling and Derivative s Pricing. Applied Mathematical Finance , 14(2):153–169, 2007

    work page 2007

  7. [6]

    D. P. Bertsekas. Dynamic Programming and Optimal Control, Volume I . Athena Scientific, Belmont, Mass., third edition, 2005

    work page 2005

  8. [7]

    Boogert and C

    A. Boogert and C. de Jong. Gas Storage Valuation Using a Monte C arlo Method. Journal of Derivatives , 15:81–91, 2008

    work page 2008

Show all 31 references
  1. [8]

    Cartea and M

    A. Cartea and M. Figueroa. Pricing in Electricity Markets: a Mean R everting Jump Diffusion Model with Seasonality. Applied Mathematical Finance, No. 4, December 2005, 12(4):313–335, 2005

    2005

  2. [9]

    Cont and P

    R. Cont and P. Tankov. Financial Modelling with Jump Processes . Chapman and Hall, 2004

    2004

  3. [10]

    Cufaro Petroni

    N. Cufaro Petroni. Self-decomposability and Self-similarity: a Co ncise Primer. Physica A, Statistical Mechanics and its Applications , 387(7-9):1875–1894, 2008

    2008

  4. [11]

    Cufaro Petroni and P

    N. Cufaro Petroni and P. Sabino. Coupling Poisson Processes b y Self- decomposability. Mediterranean Journal of Mathematics , 14(2):69, 2017

    2017

  5. [12]

    Cufaro Petroni and P

    N. Cufaro Petroni and P. Sabino. Pricing exchange options with correlated jump diffusion processes. Quantitative Finance , pages 1–13, 2018. 20

    2018

  6. [14]

    Cummins, G

    M. Cummins, G. Kiely, and B. Murphy. Gas Storage Valuation unde r L´ evy Processes using Fast Fourier Transform. Journal of Energy Markets , 4:43–86, 2017

    2017

  7. [15]

    S. Deng. Stochastic Models of Energy Commodity Prices and The ir Applica- tions: Mean-reversion with Jumps and Spikes. Citeseer, 2000

    2000

  8. [16]

    Hambly, S

    B. Hambly, S. Howison, and T. Kluge. Information-Based Models for Finance and Insurance. Quantitative Finance , 9(8):937–949, 2009

    2009

  9. [17]

    M. Kjaer. Pricing of Swing Options in a Mean Reverting Model with J umps. Applied Mathematical Finance , 15(5-6):479–502, 2008

    2008

  10. [18]

    T. Kluge. Pricing Swing Options and other Electricity Derivatives. Technical report, University of Oxford, 2006. PhD Thesis, Available at http:/ /perso- math.univ-mlv.fr/users/bally.vlad/publications.html

    2006

  11. [19]

    S. G. Kou. A Jump-Diffusion Model for Option Pricing. Manage. Sci. , 48(8):1086–1101, August 2002

    2002

  12. [20]

    K¨ uchler and S

    U. K¨ uchler and S. Tappe. Bilateral Gamma Distributions and Pro cesses in Fi- nancial Mathematics. Stochastic Processes and their Applications , 118(2):261– 283, 2008

    2008

  13. [21]

    Lindskog and J

    F. Lindskog and J. McNeil. Common poisson shock models: applicat ions to insurance and credit risk modelling. ASTIN Bulletin , 33(2):209–238, 2003

    2003

  14. [22]

    F. A. Longstaff and E.S. Schwartz. Valuing American Options by S imulation: a Simple Least-Squares Approach. Review of Financial Studies , 14(1):113–147, 2001

    2001

  15. [23]

    Lucia and E.S

    J.J. Lucia and E.S. Schwartz. Electricity Prices and Power Deriva tives: Evi- dence from the Nordic Power Exchange. Review of Derivatives Research, 5(1):5– 50, Jan 2002

    2002

  16. [24]

    McKenzie

    E. McKenzie. Innovation Distridution for Gamma and Negative Bin omial Au- toregressions. Scandinavian Journal of Statistics: Theory and Applicatio ns, 14(1):79–85, 1987

    1987

  17. [25]

    Meyer-Brandis and P

    T. Meyer-Brandis and P. Tankov. Multi-factor Jump-diffusion M odels of Electricity Prices. International Journal of Theoretical and Applied Finance , 11(5):503–528, 2008. 21

    2008

  18. [26]

    Pellegrino and P

    T. Pellegrino and P. Sabino. Enhancing Least Squares Monte Car lo with Dif- fusion Bridges: an Application to Energy Facilities. Quantitative Finance , 15(5):761–772, 2015

    2015

  19. [27]

    Y. Qu, A. Dassios, and H. Zhao. Exact Simulation of Gamma-drive n Ornstein- Uhlenbeck Processes with Finite and Infinite Activity Jumps. Journal of the Operational Research Society, 0(0):1–14, 2019

    2019

  20. [28]

    P. Sabino. Forward or Backward Simulations? A Comparative Stu dy. Quanti- tative Finance , 2020. In press

    2020

  21. [29]

    K. Sato. L´ evy Processes and Infinitely Divisible Distributions. Cambridge U.P., Cambridge, 1999

    1999

  22. [30]

    Schoutens

    W. Schoutens. L´ evy Processes in Finance: Pricing Financial Derivatives. John Wiley and Sons Inc, 2003

    2003

  23. [31]

    Schwartz and J.E

    P. Schwartz and J.E. Smith. Short-term Variations and Long-t erm Dynamics in Commodity Prices. Management Science , 46(7):893–911, 2000. 22

    2000

Pith tools

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