deeraw/options-pricing-api
0
1"""2Options Pricing API — Black-Scholes, Heston, Greeks, implied vol surface, arbitrage checks.3"""4from __future__ import annotations5 6import math7from typing import List, Literal, Optional8 9import numpy as np10from fastapi import FastAPI, HTTPException11from fastapi.middleware.cors import CORSMiddleware12from pydantic import BaseModel, Field13from scipy.optimize import brentq14from scipy.stats import norm15 16app = FastAPI(title="Options Pricing Lab API", version="1.0.0")17 18app.add_middleware(19 CORSMiddleware,20 allow_origins=["*"],21 allow_methods=["*"],22 allow_headers=["*"],23)24 25 26# ---------- Black-Scholes ----------27 28def _d1_d2(S: float, K: float, T: float, r: float, q: float, sigma: float):29 if T <= 0 or sigma <= 0:30 raise ValueError("T and sigma must be positive")31 vt = sigma * math.sqrt(T)32 d1 = (math.log(S / K) + (r - q + 0.5 * sigma * sigma) * T) / vt33 return d1, d1 - vt34 35 36def bs_price(S, K, T, r, q, sigma, option: str) -> float:37 d1, d2 = _d1_d2(S, K, T, r, q, sigma)38 if option == "call":39 return S * math.exp(-q * T) * norm.cdf(d1) - K * math.exp(-r * T) * norm.cdf(d2)40 return K * math.exp(-r * T) * norm.cdf(-d2) - S * math.exp(-q * T) * norm.cdf(-d1)41 42 43def bs_greeks(S, K, T, r, q, sigma, option: str) -> dict:44 d1, d2 = _d1_d2(S, K, T, r, q, sigma)45 pdf_d1 = norm.pdf(d1)46 disc_q = math.exp(-q * T)47 disc_r = math.exp(-r * T)48 sqrtT = math.sqrt(T)49 50 if option == "call":51 delta = disc_q * norm.cdf(d1)52 theta = (-S * disc_q * pdf_d1 * sigma / (2 * sqrtT)53 - r * K * disc_r * norm.cdf(d2)54 + q * S * disc_q * norm.cdf(d1))55 rho = K * T * disc_r * norm.cdf(d2)56 else:57 delta = -disc_q * norm.cdf(-d1)58 theta = (-S * disc_q * pdf_d1 * sigma / (2 * sqrtT)59 + r * K * disc_r * norm.cdf(-d2)60 - q * S * disc_q * norm.cdf(-d1))61 rho = -K * T * disc_r * norm.cdf(-d2)62 63 gamma = disc_q * pdf_d1 / (S * sigma * sqrtT)64 vega = S * disc_q * pdf_d1 * sqrtT65 return {66 "delta": delta,67 "gamma": gamma,68 "theta": theta / 365.0, # per calendar day69 "vega": vega / 100.0, # per 1 vol point70 "rho": rho / 100.0, # per 1% rate move71 }72 73 74def implied_vol(price, S, K, T, r, q, option: str) -> Optional[float]:75 intrinsic = max(0.0, (S * math.exp(-q * T) - K * math.exp(-r * T)) if option == "call"76 else (K * math.exp(-r * T) - S * math.exp(-q * T)))77 upper = S * math.exp(-q * T) if option == "call" else K * math.exp(-r * T)78 if price <= intrinsic + 1e-9 or price >= upper - 1e-9:79 return None80 try:81 return brentq(lambda v: bs_price(S, K, T, r, q, v, option) - price,82 1e-4, 5.0, xtol=1e-7, maxiter=200)83 except (ValueError, RuntimeError):84 return None85 86 87# ---------- Heston (Lewis 2001 / Carr-Madan style FFT-free integration) ----------88 89def heston_price(S, K, T, r, q, v0, kappa, theta, sigma_v, rho, option: str) -> float:90 """Heston European price via Gauss-Legendre quadrature on the Lewis integrand."""91 def char_func(phi):92 a = kappa * theta93 u = -0.594 b = kappa95 d = np.sqrt((rho * sigma_v * 1j * phi - b) ** 296 - sigma_v ** 2 * (2 * u * 1j * phi - phi ** 2))97 g = (b - rho * sigma_v * 1j * phi - d) / (b - rho * sigma_v * 1j * phi + d)98 C = ((r - q) * 1j * phi * T99 + (a / sigma_v ** 2) * ((b - rho * sigma_v * 1j * phi - d) * T100 - 2 * np.log((1 - g * np.exp(-d * T)) / (1 - g))))101 D = ((b - rho * sigma_v * 1j * phi - d) / sigma_v ** 2) \102 * ((1 - np.exp(-d * T)) / (1 - g * np.exp(-d * T)))103 return np.exp(C + D * v0 + 1j * phi * np.log(S))104 105 def integrand(phi, j):106 if j == 1:107 num = np.exp(-1j * phi * np.log(K)) * char_func(phi - 1j)108 den = 1j * phi * char_func(-1j)109 else:110 num = np.exp(-1j * phi * np.log(K)) * char_func(phi)111 den = 1j * phi112 return np.real(num / den)113 114 nodes, weights = np.polynomial.legendre.leggauss(96)115 a_lo, a_hi = 1e-5, 100.0116 x = 0.5 * (a_hi - a_lo) * nodes + 0.5 * (a_hi + a_lo)117 w = 0.5 * (a_hi - a_lo) * weights118 119 P1 = 0.5 + (1 / np.pi) * np.sum(w * np.array([integrand(xi, 1) for xi in x]))120 P2 = 0.5 + (1 / np.pi) * np.sum(w * np.array([integrand(xi, 2) for xi in x]))121 122 call = S * math.exp(-q * T) * P1 - K * math.exp(-r * T) * P2123 if option == "call":124 return float(max(call, 0.0))125 # put-call parity126 put = call - S * math.exp(-q * T) + K * math.exp(-r * T)127 return float(max(put, 0.0))128 129 130# ---------- Schemas ----------131 132class PriceReq(BaseModel):133 S: float = Field(..., gt=0, description="Spot")134 K: float = Field(..., gt=0, description="Strike")135 T: float = Field(..., gt=0, description="Years to expiry")136 r: float = Field(..., description="Risk-free rate")137 q: float = Field(0.0, description="Dividend yield")138 sigma: float = Field(..., gt=0, description="Vol")139 option: Literal["call", "put"] = "call"140 141 142class HestonReq(BaseModel):143 S: float = Field(..., gt=0)144 K: float = Field(..., gt=0)145 T: float = Field(..., gt=0)146 r: float = 0.05147 q: float = 0.0148 v0: float = Field(0.04, gt=0)149 kappa: float = Field(2.0, gt=0)150 theta: float = Field(0.04, gt=0)151 sigma_v: float = Field(0.5, gt=0)152 rho: float = Field(-0.7, ge=-1, le=1)153 option: Literal["call", "put"] = "call"154 155 156class IVReq(BaseModel):157 price: float = Field(..., gt=0)158 S: float; K: float; T: float; r: float159 q: float = 0.0160 option: Literal["call", "put"] = "call"161 162 163class SurfaceQuote(BaseModel):164 K: float; T: float; price: float165 option: Literal["call", "put"] = "call"166 167 168class SurfaceReq(BaseModel):169 S: float; r: float; q: float = 0.0170 quotes: List[SurfaceQuote]171 172 173# ---------- Routes ----------174 175@app.get("/health")176def health():177 return {"status": "ok", "service": "options-pricing-lab"}178 179 180@app.post("/price/black-scholes")181def price_bs(req: PriceReq):182 try:183 price = bs_price(req.S, req.K, req.T, req.r, req.q, req.sigma, req.option)184 greeks = bs_greeks(req.S, req.K, req.T, req.r, req.q, req.sigma, req.option)185 return {"model": "Black-Scholes", "price": price, "greeks": greeks}186 except Exception as e:187 raise HTTPException(400, str(e))188 189 190@app.post("/price/heston")191def price_heston(req: HestonReq):192 try:193 feller = 2 * req.kappa * req.theta - req.sigma_v ** 2194 price = heston_price(req.S, req.K, req.T, req.r, req.q,195 req.v0, req.kappa, req.theta, req.sigma_v, req.rho,196 req.option)197 return {198 "model": "Heston",199 "price": price,200 "feller_condition": feller,201 "feller_satisfied": feller > 0,202 }203 except Exception as e:204 raise HTTPException(400, str(e))205 206 207@app.post("/iv")208def iv(req: IVReq):209 sigma = implied_vol(req.price, req.S, req.K, req.T, req.r, req.q, req.option)210 if sigma is None:211 raise HTTPException(422, "Could not solve for IV (price outside no-arb bounds)")212 greeks = bs_greeks(req.S, req.K, req.T, req.r, req.q, sigma, req.option)213 return {"implied_vol": sigma, "greeks": greeks}214 215 216@app.post("/surface")217def surface(req: SurfaceReq):218 """219 Build implied vol surface from a list of (K, T, price) quotes.220 Returns a grid + arbitrage diagnostics:221 - vertical (butterfly) arb: convexity of C(K) at fixed T222 - calendar arb: total variance must be non-decreasing in T at fixed moneyness223 """224 points = []225 for q in req.quotes:226 sigma = implied_vol(q.price, req.S, q.K, q.T, req.r, req.q, q.option)227 if sigma is not None:228 points.append({229 "K": q.K, "T": q.T,230 "moneyness": q.K / req.S,231 "iv": sigma,232 "total_variance": sigma * sigma * q.T,233 })234 235 # group by tenor for butterfly check236 butterfly_violations: list[dict] = []237 by_T: dict[float, list[dict]] = {}238 for p in points:239 by_T.setdefault(p["T"], []).append(p)240 for T_, pts in by_T.items():241 pts_sorted = sorted(pts, key=lambda x: x["K"])242 # convexity check on call prices (re-price at IV)243 prices = [bs_price(req.S, p["K"], p["T"], req.r, req.q, p["iv"], "call")244 for p in pts_sorted]245 for i in range(1, len(pts_sorted) - 1):246 k_lo, k_mid, k_hi = (pts_sorted[i - 1]["K"],247 pts_sorted[i]["K"],248 pts_sorted[i + 1]["K"])249 w = (k_hi - k_mid) / (k_hi - k_lo)250 interp = w * prices[i - 1] + (1 - w) * prices[i + 1]251 if prices[i] > interp + 1e-6:252 butterfly_violations.append({253 "T": T_, "K": k_mid,254 "C_mid": prices[i], "C_interp": interp,255 })256 257 # calendar arbitrage check on total variance per moneyness bucket258 calendar_violations: list[dict] = []259 pts_by_money = sorted(points, key=lambda x: (round(x["moneyness"], 2), x["T"]))260 bucket: dict[float, list[dict]] = {}261 for p in pts_by_money:262 bucket.setdefault(round(p["moneyness"], 2), []).append(p)263 for m, ps in bucket.items():264 ps = sorted(ps, key=lambda x: x["T"])265 for i in range(1, len(ps)):266 if ps[i]["total_variance"] + 1e-8 < ps[i - 1]["total_variance"]:267 calendar_violations.append({268 "moneyness": m,269 "T_short": ps[i - 1]["T"], "T_long": ps[i]["T"],270 "var_short": ps[i - 1]["total_variance"],271 "var_long": ps[i]["total_variance"],272 })273 274 return {275 "points": points,276 "butterfly_violations": butterfly_violations,277 "calendar_violations": calendar_violations,278 "n_points": len(points),279 }280 281 282@app.post("/greeks/strip")283def greeks_strip(req: PriceReq):284 """Greeks across a strike strip — for plotting Δ/Γ/V profiles."""285 strikes = np.linspace(req.K * 0.6, req.K * 1.4, 41)286 rows = []287 for k in strikes:288 try:289 g = bs_greeks(req.S, float(k), req.T, req.r, req.q, req.sigma, req.option)290 p = bs_price(req.S, float(k), req.T, req.r, req.q, req.sigma, req.option)291 rows.append({"K": float(k), "price": p, **g})292 except Exception:293 continue294 return {"rows": rows}295