บทนำ: ทำไม Black-Scholes ถึงไม่เพียงพอสำหรับคริปโตออปชัน

ในโลกของการเงินแบบดั้งเดิม Black-Scholes Model ถือเป็นหนึ่งในโมเดลที่ประสบความสำเร็จมากที่สุดในประวัติศาสตร์ แต่เมื่อนำมาใช้กับคริปโตออปชัน โมเดลนี้กลับมีจุดอ่อนสำคัญหลายประการ โดยเฉพาะเรื่อง **"Volatility Smile"** และ **"Volatility Skew"** ที่ทำให้ราคาที่คำนวณได้คลาดเคลื่อนจากความเป็นจริงอย่างมาก บทความนี้จะพาคุณเข้าใจพื้นฐาน Black-Scholes ที่ปรับปรุงแล้วสำหรับคริปโต พร้อมวิธีจัดการ Volatility Smile ที่ทำให้ราคาออปชันในตลาดคริปโตมีรูปแบบที่ผิดเพี้ยนจากทฤษฎีแบบ "รอยยิ้ม" อย่างเห็นได้ชัด จากประสบการณ์การสร้างระบบเทรดออปชันของตัวเองมาหลายปี ผมพบว่าการใช้ Black-Scholes แบบดิบๆ มักให้ผลลัพธ์ที่ผิดพลาดมากสำหรับคริปโต เนื่องจากความผันผวนที่สูงและการกระจุกตัวของสภาพคล่องในบาง Strike Price

พื้นฐาน Black-Scholes Model สำหรับคริปโตออปชัน

Black-Scholes Model แบบดั้งเดิมมีสมมติฐานหลักว่า ราคาสินทรัพย์เคลื่อนที่แบบ Geometric Brownian Motion และ Volatility คงที่ตลอดอายุสัญญา สำหรับ Call Option สูตรพื้นฐานคือ:
import numpy as np
from scipy.stats import norm

def black_scholes_basic(S, K, T, r, sigma, option_type='call'):
    """
    Black-Scholes Model พื้นฐานสำหรับคริปโตออปชัน
    
    Parameters:
    - S: ราคาปัจจุบันของสินทรัพย์ (USD)
    - K: Strike Price (USD)
    - T: เวลาที่เหลือจนหมดอายุ (ปี)
    - r: อัตราดอกเบี้ยปราศจากความเสี่ยง (ต่อปี)
    - sigma: ความผันผวน (ต่อปี)
    - option_type: 'call' หรือ 'put'
    
    Returns:
    - price: ราคาออปชันตามทฤษฎี
    """
    d1 = (np.log(S / K) + (r + 0.5 * sigma**2) * T) / (sigma * np.sqrt(T))
    d2 = d1 - sigma * np.sqrt(T)
    
    if option_type == 'call':
        price = S * norm.cdf(d1) - K * np.exp(-r * T) * norm.cdf(d2)
    else:
        price = K * np.exp(-r * T) * norm.cdf(-d2) - S * norm.cdf(-d1)
    
    return price

ตัวอย่าง: Bitcoin Call Option

ราคา BTC ปัจจุบัน $67,000, Strike $70,000, หมดอายุ 30 วัน

Volatility 80% ต่อปี, อัตราดอกเบี้ย 5%

S = 67000 K = 70000 T = 30 / 365 r = 0.05 sigma = 0.80 price = black_scholes_basic(S, K, T, r, sigma, 'call') print(f"ราคา Call Option ตาม Black-Scholes: ${price:.2f}")

ผลลัพธ์จะแตกต่างจากราคาตลาดจริงมาก เนื่องจากไม่คำนึงถึง Volatility Smile

ปัญหา Volatility Smile ในตลาดคริปโต

Volatility Smile เกิดจากปรากฏการณ์ที่ Options ที่มี Strike Price ต่างกันมี Implied Volatility (IV) ที่ต่างกัน แม้จะมีความน่าจะเป็นที่จะเคลื่อนที่ไปถึงราคานั้นเท่ากันก็ตาม ในตลาดคริปโต Volatility Smile มีลักษณะเฉพาะ: **ลักษณะเฉพาะของ Volatility Smile ในคริปโต:** - **Downside Volatility Premium**: ความผันผวงสำหรับ Put Options มักสูงกว่า Call Options มาก (Skew ลง) - **Short-Dated Extremes**: Options ระยะสั้นมี Smile ที่รุนแรงกว่าระยะยาว - **Jump Risk Premium**: ราคาคริปโตมีโอกาสล้มเหลวหรือพุ่งสูงผิดปกติได้ง่าย - **Liquidity Skew**: Strike Price ที่ไกลจาก Spot มีสภาพคล่องต่ำ ทำให้ IV บิดเบี้ยว
import numpy as np
import matplotlib.pyplot as plt

def generate_volatility_smile(strikes, spot_price, atm_vol, skew_factor=0.3):
    """
    สร้าง Volatility Smile สำหรับคริปโตออปชัน
    
    Parameters:
    - strikes: รายการ Strike Prices
    - spot_price: ราคา Spot ปัจจุบัน
    - atm_vol: ATM Implied Volatility
    - skew_factor: ความเข้มข้นของ Skew
    
    Returns:
    - implied_vols: IV สำหรับแต่ละ Strike
    """
    implied_vols = []
    moneyness = np.log(strikes / spot_price)
    
    for m in moneyness:
        # SABR-inspired volatility smile formula
        # ใช้ quadratic smile แทน flat volatility
        smile_component = skew_factor * m + 0.15 * m**2
        iv = atm_vol * (1 + smile_component)
        implied_vols.append(iv)
    
    return np.array(implied_vols)

ตัวอย่าง: ETH Options หมดอายุ 7 วัน

spot_eth = 3500 strikes = np.array([3000, 3200, 3400, 3500, 3700, 3900, 4100]) atm_vol = 1.20 # 120% annualized ivs = generate_volatility_smile(strikes, spot_eth, atm_vol, skew_factor=0.4) print("=== Volatility Smile สำหรับ ETH Options ===") print(f"{'Strike':<10} {'Moneyness':<12} {'IV (%)':<10} {'สถานะ':<10}") print("-" * 45) for k, iv in zip(strikes, ivs): moneyness = (k / spot_eth - 1) * 100 if k < spot_eth: status = "ITM Put" elif k > spot_eth: status = "OTM Put" else: status = "ATM" print(f"${k:<9} {moneyness:>+8.1f}% {iv*100:>6.1f}% {status}")

วาดกราฟ Volatility Smile

plt.figure(figsize=(12, 6)) plt.plot(moneyness := (strikes / spot_eth - 1) * 100, ivs * 100, 'b-o', linewidth=2) plt.axhline(y=atm_vol * 100, color='r', linestyle='--', label=f'Flat Vol = {atm_vol*100:.0f}%') plt.xlabel('Moneyness (%)') plt.ylabel('Implied Volatility (%)') plt.title('Volatility Smile - ETH Options (7DTE)') plt.legend() plt.grid(True, alpha=0.3) plt.savefig('volatility_smile.png', dpi=150) plt.show() print("\n📊 สังเกตได้ว่า OTM Puts (moneyness ติดลบ) มี IV สูงกว่าอย่างมีนัยสำคัญ")

การปรับ Black-Scholes ด้วย Local Volatility Model

วิธีหนึ่งในการแก้ปัญหา Volatility Smile คือการใช้ Local Volatility Model ที่ให้ Volatility เปลี่ยนตามราคาและเวลา ทำให้สะท้อน Volatility Surface ได้แม่นยำกว่า
import numpy as np
from scipy.interpolate import RectBivariateSpline, interp1d
from scipy.stats import norm

class LocalVolatilityModel:
    """
    Local Volatility Model สำหรับคริปโตออปชัน
    ใช้ Dupire's Formula ในการคำนวณ Local Volatility
    """
    
    def __init__(self, strikes, maturities, implied_vols):
        """
        Parameters:
        - strikes: รายการ Strike Prices
        - maturities: รายการ Maturities (ปี)
        - implied_vols: 2D array ของ IV surface [maturity, strike]
        """
        self.strikes = strikes
        self.maturities = maturities
        self.iv_surface = implied_vols
        
        # สร้าง interpolator สำหรับ IV surface
        self.iv_interpolator = RectBivariateSpline(
            maturities, strikes, implied_vols, 
            kx=2, ky=2
        )
    
    def get_local_vol(self, S, T):
        """
        คำนวณ Local Volatility ที่ราคา S และเวลา T
        ใช้ Dupire's Formula
        
        dupire_vol² = [∂C/∂T + (r-q)K∂C/∂K + qC] / [½K²∂²C/∂K²]
        """
        K = S  # ใช้ราคาปัจจุบันเป็น Strike โดยประมาณ
        T = max(T, 1e-6)  # ป้องกัน T=0
        
        # หา IV ณ เวลานี้
        iv = self.iv_interpolator(T, K, grid=False)
        iv = max(iv, 1e-6)
        
        # คำนวณ Local Volatility (simplified)
        # ใน implementation จริงต้องใช้ numerical differentiation
        local_vol = iv * (1 + 0.3 * (1 - np.exp(-2*T)))  # Mean reversion effect
        
        return local_vol
    
    def price_option(self, S, K, T, r, q, option_type='call'):
        """
        กำหนดราคาออปชันด้วย Local Volatility
        ใช้ Monte Carlo Simulation
        """
        dt = T / 100  # 100 time steps
        n_paths = 10000
        n_steps = int(T / dt)
        
        # Simulate paths with local vol
        Z = np.random.standard_normal((n_paths, n_steps))
        
        paths = np.zeros((n_paths, n_steps + 1))
        paths[:, 0] = S
        
        for i in range(n_steps):
            t = i * dt
            vol = self.get_local_vol(paths[:, i], t + dt)
            paths[:, i + 1] = paths[:, i] * np.exp(
                (r - q - 0.5 * vol**2) * dt + vol * np.sqrt(dt) * Z[:, i]
            )
        
        # Calculate payoff
        final_prices = paths[:, -1]
        if option_type == 'call':
            payoffs = np.maximum(final_prices - K, 0)
        else:
            payoffs = np.maximum(K - final_prices, 0)
        
        # Discount to present
        price = np.exp(-r * T) * np.mean(payoffs)
        
        return price, np.std(payoffs) / np.sqrt(n_paths)  # price, std error


ตัวอย่างการใช้งาน

strikes = np.array([3000, 3500, 4000, 4500, 5000]) maturities = np.array([7, 14, 30, 60]) / 365

สร้าง IV surface (ตัวอย่าง)

iv_surface = np.array([ [1.40, 1.20, 1.00, 1.20, 1.50], # 7DTE [1.30, 1.15, 1.00, 1.10, 1.35], # 14DTE [1.20, 1.10, 1.00, 1.05, 1.25], # 30DTE [1.10, 1.05, 1.00, 1.02, 1.15], # 60DTE ]) lv_model = LocalVolatilityModel(strikes, maturities, iv_surface)

กำหนดราคา ETH Call Option

S, K, T, r, q = 3500, 4000, 14/365, 0.05, 0.0 price, se = lv_model.price_option(S, K, T, r, q, 'call') print(f"ETH Call Option (K=${K}, 14DTE)") print(f"ราคาจาก Local Vol Model: ${price:.2f} (±${se:.2f})") print(f"Local Vol ณ ราคาปัจจุบัน: {lv_model.get_local_vol(S, T)*100:.1f}%")

การปรับแต่ง Volatility Skew ด้วย SABR Model

SABR Model (Stochastic Alpha Beta Rho) เป็นอีกหนึ่งวิธีที่นิยมใช้ในการอธิบาย Volatility Smile โดยเฉพาะในตลาดที่มี Skew รุนแรงอย่างคริปโต
import numpy as np
from scipy.optimize import minimize
from scipy.stats import norm

class SABRVolatilityModel:
    """
    SABR Model สำหรับการปรับ Volatility Smile
    
    Parameters:
    - alpha: Initial volatility (ATM vol scale)
    - beta: CEV exponent (0 < beta <= 1)
    - rho: Correlation between asset and vol
    - nu: Volatility of volatility
    """
    
    def __init__(self, alpha=0.5, beta=0.5, rho=-0.3, nu=0.5):
        self.alpha = alpha
        self.beta = beta
        self.rho = rho
        self.nu = nu
    
    def sabr_vol(self, F, K, T, alpha, beta, rho, nu):
        """
        คำนวณ SABR implied volatility ด้วย Hagan's formula
        
        Parameters:
        - F: Forward price
        - K: Strike price
        - T: Time to maturity
        """
        if abs(F - K) < 1e-10:  # ATM case
            FK_mid = F
            term1 = alpha / (FK_mid ** (1 - beta))
            term2 = 1 + ((1 - beta)**2 / 24 * alpha**2 / (FK_mid ** (2*(1-beta)))
                        + 0.25 * rho * beta * nu * alpha / (FK_mid ** (1 - beta))
                        + (2 - 3*rho**2) / 24 * nu**2) * T
            return term1 * term2
        
        FK = F * K
        FK_mid = np.sqrt(FK)
        log_FK = np.log(F / K)
        
        z = (nu / alpha) * FK_mid ** (1 - beta) * log_FK
        x_z = np.log((np.sqrt(1 - 2*rho*z + z**2) + z - rho) / (1 - rho))
        
        denom = 1 + (1 - beta)**2 / 24 * log_FK**2 + (1 - beta)**4 / 1920 * log_FK**4
        
        term1 = alpha / (FK ** ((1 - beta) / 2) * denom)
        term2 = z / x_z
        term3 = 1 + ((1 - beta)**2 / 24 * alpha**2 / (FK ** (1 - beta))
                    + 0.25 * rho * beta * nu * alpha / (FK ** ((1 - beta)/2))
                    + (2 - 3*rho**2) / 24 * nu**2) * T
        
        return term1 * term2 * term3
    
    def fit_to_market(self, strikes, market_ivs, forward, T):
        """
        Fit SABR parameters ให้ตรงกับ Market IV
        """
        def objective(params):
            alpha, beta, rho, nu = params
            if alpha <= 0 or nu <= 0 or beta <= 0 or beta > 1 or abs(rho) >= 1:
                return 1e10
            
            total_error = 0
            for K, mkt_iv in zip(strikes, market_ivs):
                try:
                    sabr_iv = self.sabr_vol(forward, K, T, alpha, beta, rho, nu)
                    total_error += (sabr_iv - mkt_iv) ** 2
                except:
                    total_error += 1e10
            return total_error
        
        # Initial guess
        x0 = [0.3, 0.5, -0.3, 0.4]
        bounds = [(0.01, 2), (0.01, 0.99), (-0.99, 0.99), (0.01, 3)]
        
        result = minimize(objective, x0, method='L-BFGS-B', bounds=bounds)
        self.alpha, self.beta, self.rho, self.nu = result.x
        
        return result.x
    
    def get_vol(self, F, K, T):
        """คืนค่า IV ที่คำนวณจาก fitted SABR"""
        return self.sabr_vol(F, K, T, self.alpha, self.beta, self.rho, self.nu)


ตัวอย่าง: BTC Options Smile

forward_btc = 67000 strikes = np.array([60000, 63000, 66000, 67000, 70000, 73000, 76000]) market_ivs = np.array([0.95, 0.85, 0.75, 0.70, 0.78, 0.88, 1.00]) # IV จากตลาด T = 30 / 365 sabr = SABRVolatilityModel() fitted_params = sabr.fit_to_market(strikes, market_ivs, forward_btc, T) print("=== SABR Fitted Parameters ===") print(f"Alpha (Initial Vol): {fitted_params[0]:.4f}") print(f"Beta (CEV exponent): {fitted_params[1]:.4f}") print(f"Rho (Correlation): {fitted_params[2]:.4f}") print(f"Nu (Vol of Vol): {fitted_params[3]:.4f}")

เปรียบเทียบ Market IV vs SABR IV

print("\n=== Market IV vs SABR IV ===") print(f"{'Strike':<10} {'Market IV':<12} {'SABR IV':<12} {'Error':<10}") print("-" * 45) total_error = 0 for K, mkt_iv in zip(strikes, market_ivs): sabr_iv = sabr.get_vol(forward_btc, K, T) error = sabr_iv - mkt_iv total_error += abs(error) print(f"${K:<9} {mkt_iv*100:>8.1f}% {sabr_iv*100:>8.1f}% {error*100:>+6.2f}%") print(f"\nAverage Absolute Error: {total_error/len(strikes)*100:.2f}%")

ข้อผิดพลาดที่พบบ่อยและวิธีแก้ไข

1. การใช้ Flat Volatility โดยไม่ปรับ Volatility Smile

ปัญหา: นักพัฒนาหลายคนใช้ Historical Volatility เดียวสำหรับทุก Strike และทุก Maturity ทำให้ราคาออปชันที่ OTM คลาดเคลื่อนอย่างมาก
# ❌ วิธีผิด: ใช้ Flat Volatility
def price_with_flat_vol(S, K, T, r, sigma, option_type='call'):
    d1 = (np.log(S / K) + (r + 0.5 * sigma**2) * T) / (sigma * np.sqrt(T))
    d2 = d1 - sigma * np.sqrt(T)
    if option_type == 'call':
        return S * norm.cdf(d1) - K * np.exp(-r * T) * norm.cdf(d2)
    else:
        return K * np.exp(-r * T) * norm.cdf(-d2) - S * norm.cdf(-d1)

✅ วิธีถูก: ใช้ Strike-Specific Volatility

def price_with_volatility_smile(S, K, T, r, iv_function, option_type='call'): """ iv_function: function(K) ที่คืนค่า IV สำหรับ Strike นั้นๆ สามารถใช้ SABR, Local Vol, หรือ Interpolated IV surface """ sigma = iv_function(K) # ใช้ IV เฉพาะของ Strike นี้ d1 = (np.log(S / K) + (r + 0.5 * sigma**2) * T) / (sigma * np.sqrt(T)) d2 = d1 - sigma * np.sqrt(T) if option_type == 'call': return S * norm.cdf(d1) - K * np.exp(-r * T) * norm.cdf(d2) else: return K * np.exp(-r * T) * norm.cdf(-d2) - S * norm.cdf(-d1)

ตัวอย่างการใช้งาน

from scipy.interpolate import interp1d strikes = np.array([60000, 65000, 70000, 75000, 80000]) ivs = np.array([0.95, 0.82, 0.70, 0.82, 0.95])

สร้าง interpolation function สำหรับ IV surface

iv_interpolator = interp1d(strikes, ivs, kind='quadratic', fill_value='extrapolate')

ราคาออปชัน OTM ด้วยวิธีที่ถูกต้อง

K_otm = 80000 # OTM Call price_correct = price_with_volatility_smile(67000, K_otm, 30/365, 0.05, iv_interpolator) price_wrong = price_with_flat_vol(67000, K_otm, 30/365, 0.05, 0.70) print(f"ราคาที่ถูกต้อง (ใช้ Smile): ${price_correct:.2f}") print(f"ราคาที่ผิดพลาด (Flat Vol): ${price_wrong:.2f}") print(f"ความคลาดเคลื่อน: ${abs(price_correct - price_wrong):.2f} ({abs(price_correct - price_wrong)/price_correct*100:.1f}%)")

2. การละเลย Jump Risk ในคริปโต

ปัญหา: Black-Scholes และ SABR ไม่ได้คำนึงถึง Extreme Jumps ที่เกิดขึ้นบ่อยในคริปโต เช่น การล้มของ FTX หรือ Elon Musk Tweet
# ❌ วิธีผิด: ใช้โมเดลที่ไม่มี Jump
def merton_jump_penalty(F, K, T, base_vol, jump_intensity=2, jump_size=0.3):
    """
    เพิ่ม Jump Risk Premium ให้กับ Base Volatility
    
    Parameters:
    - jump_intensity: lambda (จำนวน Jump ต่อปี)
    - jump_size: ขนาด Jump เฉลี่ย (เป็น % ของราคา)
    
    สูตร: σ²_adj = σ² + λ * (σ_j² + μ_j²)
    โดย σ_j และ μ_j คือ vol และ mean ของ log-jump
    """
    # Log-normal jump parameters
    mu_j = np.log(1 + jump_size) - 0.5 * (jump_size ** 2) / 2
    sigma_j = jump_size * 0.5  # Jump volatility
    
    # Variance from jumps
    jump_variance = jump_intensity * (sigma_j**2 + mu_j**2)
    
    # Adjusted total variance
    total_variance = base_vol**2 + jump_variance
    adjusted_vol = np.sqrt(total_variance)
    
    return adjusted_vol

✅ วิธีถูก: ใช้ Merton Jump-Diffusion Model

def price_with_jump_diffusion(S, K, T, r, sigma_base, jump_lambda=2, jump_mean=0.3, jump_vol=0.15): """ กำหนดราคาออปชันด้วย Merton Jump-Diffusion Model สูตร: C = Σ(p^n) * BS_Call(S*e^(n*μ_j + n²σ_j²/2), K, T_n, r, σ) โดย p^n = e^(-λT) * (λT)^n / n! (Poisson probability) """ n_steps = 50 jump_price = 0 for n in range(30): # ใช้ 30 jumps ก็เพียงพอ # Probability of n jumps poisson_prob = np.exp(-jump_lambda * T) * (jump_lambda * T)**n / np.math.factorial(n) # Jump impact on price mu_j = np.log(1 + jump_mean) - 0.5 * jump_vol**2 jump_adjustment = n * mu_j + 0.5 * n * jump_vol**2 # Effective forward price F_adj = S * np.exp(jump_adjustment + (r - 0.5 * sigma_base**2) * T) # BS price with adjusted forward d1 = (np.log(F_adj / K) + 0.5 * sigma_base**2 * T) / (sigma_base * np.sqrt(T)) d2 = d1 - sigma_base * np.sqrt(T) call_price = F_adj * norm.cdf(d1) - K * np.exp(-r * T) * norm.cdf(d2) jump_price += poisson_prob * call_price return jump_price

ตั