surface_lib/models/svi/
svi_model.rs

1// src/models/svi/svi_model.rs
2
3//! Stochastic Volatility Inspired (SVI) model implementation
4//!
5//! The SVI model provides a parametric form for the implied volatility surface
6//! that is both flexible and enforces no-arbitrage conditions. The total variance
7//! w(k,t) is given by:
8//!
9//! w(k) = a + b * (ρ(k-m) + sqrt((k-m)² + σ²))
10//!
11//! where k is log-moneyness, and the parameters are:
12//! - a: vertical shift (controls ATM level)
13//! - b: slope factor (controls overall level)
14//! - ρ: asymmetry parameter (skew, -1 < ρ < 1)
15//! - m: horizontal shift (ATM location)
16//! - σ: curvature parameter (controls smile curvature)
17
18use crate::models::traits::SurfaceModel;
19use anyhow::{anyhow, Result};
20use serde::{Deserialize, Serialize};
21
22/// Parameters for the SVI (Stochastic Volatility Inspired) model for a single maturity.
23///
24/// The SVI model provides a parametric form for implied volatility that ensures
25/// no calendar arbitrage and can be configured to avoid butterfly arbitrage.
26#[derive(Debug, Clone, PartialEq, Serialize, Deserialize)]
27pub struct SVIParams {
28    /// Time to maturity (years)
29    pub t: f64,
30    /// Vertical shift parameter (controls ATM variance level)
31    pub a: f64,
32    /// Slope factor (controls overall variance level)
33    pub b: f64,
34    /// Asymmetry parameter (skew, must be in (-1, 1))
35    pub rho: f64,
36    /// Horizontal shift parameter (ATM location)
37    pub m: f64,
38    /// Curvature parameter (controls smile curvature, must be > 0)
39    pub sigma: f64,
40}
41
42/// Helper function to validate SVI parameters for mathematical and no-arbitrage constraints.
43fn validate_svi_params(t: f64, a: f64, b: f64, rho: f64, m: f64, sigma: f64) -> Result<()> {
44    if t <= 0.0 || !t.is_finite() {
45        return Err(anyhow!(
46            "SVIParams validation: time to expiry (t={}) must be > 0 and finite",
47            t
48        ));
49    }
50    // Allow negative values for 'a' (vertical shift) provided they are finite. The
51    // no-arbitrage constraint below will ensure that the combination of parameters
52    // still yields non-negative total variance across strikes.
53    if !a.is_finite() {
54        return Err(anyhow!(
55            "SVIParams validation: parameter a (a={}) must be finite",
56            a
57        ));
58    }
59    if b <= 0.0 || !b.is_finite() {
60        return Err(anyhow!(
61            "SVIParams validation: parameter b (b={}) must be > 0 and finite",
62            b
63        ));
64    }
65    if rho <= -1.0 || rho >= 1.0 || !rho.is_finite() {
66        return Err(anyhow!(
67            "SVIParams validation: parameter rho (rho={}) must be in (-1, 1) and finite",
68            rho
69        ));
70    }
71    if !m.is_finite() {
72        return Err(anyhow!(
73            "SVIParams validation: parameter m (m={}) must be finite",
74            m
75        ));
76    }
77    if sigma <= 0.0 || !sigma.is_finite() {
78        return Err(anyhow!(
79            "SVIParams validation: parameter sigma (sigma={}) must be > 0 and finite",
80            sigma
81        ));
82    }
83
84    // Additional no-arbitrage constraint: a + b*sigma*sqrt(1-rho^2) >= 0
85    // This ensures the total variance is non-negative at the wings
86    let min_variance = a + b * sigma * (1.0 - rho * rho).sqrt();
87    if min_variance < 0.0 {
88        return Err(anyhow!(
89            "SVIParams validation: no-arbitrage constraint violated. a + b*sigma*sqrt(1-rho^2) = {} < 0", 
90            min_variance
91        ));
92    }
93
94    Ok(())
95}
96
97impl SVIParams {
98    /// Creates new SVI parameters with validation.
99    pub fn new(t: f64, a: f64, b: f64, rho: f64, m: f64, sigma: f64) -> Result<Self> {
100        validate_svi_params(t, a, b, rho, m, sigma)?;
101
102        Ok(Self {
103            t,
104            a,
105            b,
106            rho,
107            m,
108            sigma,
109        })
110    }
111
112    /// Validates the current parameter set.
113    pub fn validate(&self) -> Result<()> {
114        validate_svi_params(self.t, self.a, self.b, self.rho, self.m, self.sigma)
115    }
116}
117
118/// Represents the SVI volatility model for a single maturity slice.
119/// Contains the parameters and implements the SVI calculation logic.
120#[derive(Debug, Clone, PartialEq)]
121pub struct SVISlice {
122    pub params: SVIParams,
123}
124
125impl SVISlice {
126    /// Creates a new SVISlice instance from validated SVIParams.
127    pub fn new(params: SVIParams) -> Self {
128        // Validation should happen when creating SVIParams
129        Self { params }
130    }
131
132    /// Calculates total variance w(k) using the SVI formula:
133    /// w(k) = a + b * (ρ(k-m) + sqrt((k-m)² + σ²))
134    pub fn total_variance_at_k(&self, k: f64) -> f64 {
135        let k_minus_m = k - self.params.m;
136        let sqrt_term = (k_minus_m * k_minus_m + self.params.sigma * self.params.sigma).sqrt();
137
138        self.params.a + self.params.b * (self.params.rho * k_minus_m + sqrt_term)
139    }
140
141    /// Calculates implied volatility σ(k) from total variance.
142    pub fn implied_vol(&self, k: f64) -> f64 {
143        let total_var = self.total_variance_at_k(k);
144        if total_var <= 0.0 {
145            return 1e-6; // Small positive number to avoid issues
146        }
147        (total_var / self.params.t).sqrt()
148    }
149}
150
151// Define the 5-minute tolerance in years as a constant (matching Wing implementation)
152const FIVE_MINUTES_IN_YEARS: f64 = 5.0 / (60.0 * 24.0 * 365.0); // approx 9.51e-6
153
154// Implement the SurfaceModel trait for SVISlice
155impl SurfaceModel for SVISlice {
156    type Parameters = SVIParams;
157
158    fn parameters(&self) -> &Self::Parameters {
159        &self.params
160    }
161
162    /// Validates the internal parameters using the shared helper function.
163    fn validate_params(&self) -> Result<()> {
164        validate_svi_params(
165            self.params.t,
166            self.params.a,
167            self.params.b,
168            self.params.rho,
169            self.params.m,
170            self.params.sigma,
171        )
172    }
173
174    /// Calculates the model's implied total variance.
175    /// **Requires `t` to be within ~5 minutes of the slice's `params.t`.**
176    fn total_variance(&self, k: f64, t: f64) -> Result<f64> {
177        // Check if the provided time `t` is close enough to the slice's time `self.params.t`
178        if (t - self.params.t).abs() > FIVE_MINUTES_IN_YEARS {
179            return Err(anyhow!(
180                "SVISlice time mismatch: requested t={} is too far from slice t={}. Tolerance: {:.3e} years (~5 min)", 
181                t, self.params.t, FIVE_MINUTES_IN_YEARS
182            ));
183        }
184        if !k.is_finite() {
185            return Err(anyhow!("Log-moneyness k must be finite (k={})", k));
186        }
187
188        let total_var = self.total_variance_at_k(k);
189
190        if !total_var.is_finite() || total_var < 0.0 {
191            return Err(anyhow!(
192                "Calculated total variance is invalid: {} for k={}, t={}",
193                total_var,
194                k,
195                self.params.t
196            ));
197        }
198        Ok(total_var)
199    }
200
201    /// Checks for calendar spread arbitrage. Returns Ok(()) as it's not applicable for a single slice.
202    fn check_calendar_arbitrage(&self, _k: f64, _t1: f64, _t2: f64) -> Result<()> {
203        Ok(())
204    }
205
206    /// Checks for butterfly spread arbitrage violations at `k` and `t`.
207    /// Uses Gatheral's g(k) condition: g(k) = (1 - k*w'/(2*w))² - (w')²/4 * (1/w + 1/4) + w''/2 >= 0
208    /// **Requires `t` to be within ~5 minutes of the slice's `params.t`.**
209    fn check_butterfly_arbitrage_at_k(&self, k: f64, t: f64) -> Result<()> {
210        const EPSILON: f64 = 1e-5;
211        let tolerance = 1e-9; // Tolerance for g_k check
212
213        // Check if the provided time `t` is close enough to the slice's time `self.params.t`
214        if (t - self.params.t).abs() > FIVE_MINUTES_IN_YEARS {
215            return Err(anyhow!(
216                "SVISlice time mismatch for butterfly check: requested t={} is too far from slice t={}. Tolerance: {:.3e} years (~5 min)", 
217                t, self.params.t, FIVE_MINUTES_IN_YEARS
218            ));
219        }
220        if !k.is_finite() {
221            return Err(anyhow!(
222                "Butterfly check failed: k must be finite (k={})",
223                k
224            ));
225        }
226
227        // Use slice's exact time for consistency
228        let slice_t = self.params.t;
229        let w = self.total_variance(k, slice_t)?;
230        let w_p = self.total_variance(k - EPSILON, slice_t)?;
231        let w_n = self.total_variance(k + EPSILON, slice_t)?;
232
233        if w <= tolerance {
234            return Ok(()); // No arbitrage if variance is near zero
235        }
236
237        // Calculate first and second derivatives using finite differences
238        let w_k = (w_n - w_p) / (2.0 * EPSILON); // First derivative
239        let w_kk = (w_n - 2.0 * w + w_p) / (EPSILON * EPSILON); // Second derivative
240
241        // Gatheral's g(k) condition
242        let term1 = 1.0 - k * w_k / (2.0 * w);
243        let g_k = term1 * term1 - (w_k * w_k / 4.0) * (1.0 / w + 0.25) + w_kk / 2.0;
244
245        if g_k < -tolerance {
246            Err(anyhow!(
247                "Butterfly arbitrage detected at k={:.6}, t={:.4}. g(k) = {:.6e} < 0",
248                k,
249                t,
250                g_k
251            ))
252        } else {
253            Ok(())
254        }
255    }
256}
257
258/// Interpolates SVIParams across maturities using linear interpolation.
259/// Returns SVIParams for the requested time `t`.
260pub fn interpolate_svi_params(slices: &[(f64, SVIParams)], t: f64) -> SVIParams {
261    if slices.is_empty() {
262        panic!("Cannot interpolate SVI parameters with empty slice list");
263    }
264
265    // Clamp t to the range of provided slice times
266    let t_clamped = t.clamp(slices[0].0, slices.last().unwrap().0);
267
268    // If t matches an existing slice time exactly, return that slice's params
269    if let Some((_, params)) = slices
270        .iter()
271        .find(|(slice_t, _)| (*slice_t - t_clamped).abs() < 1e-9)
272    {
273        return params.clone();
274    }
275
276    // Find the bounding slices
277    let idx = slices.partition_point(|(slice_t, _)| *slice_t < t_clamped);
278    if idx == 0 {
279        return slices[0].1.clone();
280    }
281    if idx >= slices.len() {
282        return slices.last().unwrap().1.clone();
283    }
284
285    let (t0, params0) = &slices[idx - 1];
286    let (t1, params1) = &slices[idx];
287
288    // Linear interpolation weight
289    let weight1 = (t_clamped - t0) / (t1 - t0);
290    let weight0 = 1.0 - weight1;
291
292    // Interpolate each parameter
293    let a_interp = weight0 * params0.a + weight1 * params1.a;
294    let b_interp = weight0 * params0.b + weight1 * params1.b;
295    let rho_interp = weight0 * params0.rho + weight1 * params1.rho;
296    let m_interp = weight0 * params0.m + weight1 * params1.m;
297    let sigma_interp = weight0 * params0.sigma + weight1 * params1.sigma;
298
299    // Attempt to create new SVIParams with interpolated values
300    // Clamp values to ensure they remain valid
301    SVIParams::new(
302        t_clamped,
303        a_interp,                        // Allow negative values now
304        b_interp.max(1e-6),              // Ensure b > 0
305        rho_interp.clamp(-0.999, 0.999), // Ensure rho in (-1, 1)
306        m_interp,                        // m can be any finite value
307        sigma_interp.max(1e-6),          // Ensure sigma > 0
308    )
309    .unwrap_or_else(|err| {
310        // Fallback if interpolation results in invalid params
311        eprintln!(
312            "Warning: SVI interpolation failed for t={}: {}. Falling back to nearest slice.",
313            t_clamped, err
314        );
315        if (t_clamped - t0) < (t1 - t_clamped) {
316            params0.clone()
317        } else {
318            params1.clone()
319        }
320    })
321}
322
323/// Represents the full SVI volatility surface across multiple maturities.
324#[derive(Debug, Clone)]
325pub struct SVIModel {
326    // Slices sorted by time t
327    slices: Vec<(f64, SVIParams)>,
328    // Configurable tolerance for calendar arbitrage checks
329    calendar_arbitrage_tolerance: f64,
330}
331
332impl SVIModel {
333    /// Creates a new SVIModel (surface) from a vector of (time, params) tuples.
334    /// Sorts the slices by time and performs initial validation.
335    pub fn new(
336        mut slices: Vec<(f64, SVIParams)>,
337        calendar_arbitrage_tolerance: f64,
338    ) -> Result<Self> {
339        if slices.is_empty() {
340            return Err(anyhow!("SVIModel requires at least one slice"));
341        }
342
343        // Sort slices by time t
344        slices.sort_by(|a, b| a.0.partial_cmp(&b.0).unwrap_or(std::cmp::Ordering::Equal));
345
346        // Ensure times are unique and increasing
347        for i in 0..(slices.len() - 1) {
348            if (slices[i + 1].0 - slices[i].0).abs() < 1e-9 {
349                return Err(anyhow!("Duplicate slice time detected: {}", slices[i].0));
350            }
351        }
352
353        let model = Self {
354            slices,
355            calendar_arbitrage_tolerance,
356        };
357
358        // Perform initial validation of the surface
359        model.validate_params()?;
360        Ok(model)
361    }
362
363    /// Interpolates SVI parameters for a given time `t`.
364    fn interpolate_params(&self, t: f64) -> SVIParams {
365        interpolate_svi_params(&self.slices, t)
366    }
367}
368
369// Implement SurfaceModel for the SVIModel (Surface)
370impl SurfaceModel for SVIModel {
371    type Parameters = Vec<(f64, SVIParams)>;
372
373    fn parameters(&self) -> &Self::Parameters {
374        &self.slices
375    }
376
377    /// Validates parameters and checks calendar arbitrage across slices.
378    fn validate_params(&self) -> Result<()> {
379        if self.slices.is_empty() {
380            return Ok(());
381        }
382
383        // 1. Validate parameters for each individual slice
384        for (t, params) in &self.slices {
385            let slice = SVISlice::new(params.clone());
386            slice
387                .validate_params()
388                .map_err(|e| anyhow!("Invalid parameters for slice at t={}: {}", t, e))?;
389        }
390
391        // 2. Check for calendar arbitrage between adjacent slices
392        // We sample a few k points and ensure total variance is non-decreasing
393        if self.slices.len() > 1 {
394            let k_samples = vec![-0.5, -0.2, 0.0, 0.2, 0.5]; // Sample points
395
396            for i in 0..(self.slices.len() - 1) {
397                let (t1, params1) = &self.slices[i];
398                let (t2, params2) = &self.slices[i + 1];
399
400                let slice1 = SVISlice::new(params1.clone());
401                let slice2 = SVISlice::new(params2.clone());
402
403                for &k in &k_samples {
404                    let w1 = slice1.total_variance_at_k(k);
405                    let w2 = slice2.total_variance_at_k(k);
406
407                    if w2 < w1 - self.calendar_arbitrage_tolerance {
408                        eprintln!("Warning: Calendar arbitrage detected between t1={:.4} and t2={:.4} at k={:.3}. w1={:.6}, w2={:.6}", 
409                                t1, t2, k, w1, w2);
410                        // Note: We print a warning but don't return an error to allow flexibility
411                    }
412                }
413            }
414        }
415
416        Ok(())
417    }
418
419    /// Calculates total variance by interpolating parameters for time `t`.
420    fn total_variance(&self, k: f64, t: f64) -> Result<f64> {
421        let mut interpolated_params = self.interpolate_params(t);
422        // Set the time `t` on the interpolated parameters to the requested time
423        interpolated_params.t = t;
424
425        let temp_slice = SVISlice::new(interpolated_params);
426        temp_slice.total_variance(k, t)
427    }
428
429    /// Checks calendar arbitrage between two times at a given k.
430    fn check_calendar_arbitrage(&self, k: f64, t1: f64, t2: f64) -> Result<()> {
431        if t1 >= t2 {
432            return Err(anyhow!(
433                "Calendar check requires t1 < t2, got t1={}, t2={}",
434                t1,
435                t2
436            ));
437        }
438
439        let w1 = self.total_variance(k, t1)?;
440        let w2 = self.total_variance(k, t2)?;
441
442        if w2 < w1 - self.calendar_arbitrage_tolerance {
443            Err(anyhow!(
444                "Calendar arbitrage detected at k={:.6}: w(t1={:.4})={:.6} > w(t2={:.4})={:.6}",
445                k,
446                t1,
447                w1,
448                t2,
449                w2
450            ))
451        } else {
452            Ok(())
453        }
454    }
455
456    /// Checks butterfly arbitrage at a specific k and t by creating a temporary slice.
457    fn check_butterfly_arbitrage_at_k(&self, k: f64, t: f64) -> Result<()> {
458        let mut interpolated_params = self.interpolate_params(t);
459        interpolated_params.t = t;
460
461        let temp_slice = SVISlice::new(interpolated_params);
462        temp_slice.check_butterfly_arbitrage_at_k(k, t)
463    }
464}
465
466#[cfg(test)]
467mod tests {
468    use super::*;
469
470    // Helper function to create valid test parameters
471    fn create_test_svi_params() -> SVIParams {
472        SVIParams::new(
473            0.25, // t = 3 months
474            0.04, // a = 4% base variance level
475            0.2,  // b = slope factor
476            -0.3, // rho = negative skew
477            0.0,  // m = ATM at log-moneyness 0
478            0.2,  // sigma = curvature
479        )
480        .unwrap()
481    }
482
483    #[test]
484    fn test_svi_params_validation() {
485        // Valid parameters should work
486        let valid_params = SVIParams::new(0.25, 0.04, 0.2, -0.3, 0.0, 0.2);
487        assert!(valid_params.is_ok());
488
489        // Test invalid parameters
490        assert!(SVIParams::new(-0.1, 0.04, 0.2, -0.3, 0.0, 0.2).is_err()); // negative t
491        assert!(SVIParams::new(0.25, 0.04, -0.1, -0.3, 0.0, 0.2).is_err()); // negative b
492        assert!(SVIParams::new(0.25, 0.04, 0.2, -1.0, 0.0, 0.2).is_err()); // rho <= -1
493        assert!(SVIParams::new(0.25, 0.04, 0.2, 1.0, 0.0, 0.2).is_err()); // rho >= 1
494        assert!(SVIParams::new(0.25, 0.04, 0.2, -0.3, 0.0, -0.1).is_err()); // negative sigma
495    }
496
497    #[test]
498    fn test_svi_total_variance_calculation() {
499        let params = create_test_svi_params();
500        let slice = SVISlice::new(params.clone());
501
502        // Test ATM (k=0)
503        let w_atm = slice.total_variance_at_k(0.0);
504        let expected_atm = params.a + params.b * params.sigma; // rho*0 + sqrt(0^2 + sigma^2) = sigma
505        assert!((w_atm - expected_atm).abs() < 1e-10);
506
507        // Test positive k
508        let k_pos = 0.2;
509        let w_pos = slice.total_variance_at_k(k_pos);
510        let expected_pos = params.a
511            + params.b
512                * (params.rho * k_pos + (k_pos * k_pos + params.sigma * params.sigma).sqrt());
513        assert!((w_pos - expected_pos).abs() < 1e-10);
514    }
515
516    #[test]
517    fn test_svi_implied_volatility() {
518        let params = create_test_svi_params();
519        let slice = SVISlice::new(params);
520
521        // Test that implied volatility is positive and reasonable
522        let iv_atm = slice.implied_vol(0.0);
523        assert!(iv_atm > 0.0);
524        assert!(iv_atm < 10.0); // Sanity check
525
526        let iv_otm_call = slice.implied_vol(0.3);
527        assert!(iv_otm_call > 0.0);
528
529        let iv_otm_put = slice.implied_vol(-0.3);
530        assert!(iv_otm_put > 0.0);
531
532        // With negative rho, we expect some skew (put vol > call vol for same |k|)
533        assert!(iv_otm_put > iv_otm_call);
534    }
535}