Black-Scholes Option Pricing Explained: Black-Scholes Model with Code, Charts & Intuition

Black-Scholes Option Pricing, Black-Scholes Model, Options Greeks explained with Python

🧠 Why Black-Scholes Option Pricing matters

If you trade options, you are not just trading direction.

You are trading:

• Probability
• Time
• Uncertainty

The Black-Scholes model converts all of this into:

Price + Greeks → actual P&L behavior

🎯 What you will build

From your notebook:

• Real data using yfinance
• A contract → model → output pipeline
• Visual understanding of:
• Delta
• Gamma
• Theta
• And finally: P&L intuition

🧠 Step 1: Black-Scholes Engine

import numpy as np
from scipy.stats import norm

def black_scholes_call(S, K, T, r, sigma):
    d1 = (np.log(S / K) + (r + 0.5 * sigma**2) * T) / (sigma * np.sqrt(T))
    d2 = d1 - sigma * np.sqrt(T)

    call_price = S * norm.cdf(d1) - K * np.exp(-r * T) * norm.cdf(d2)

    delta = norm.cdf(d1)
    gamma = norm.pdf(d1) / (S * sigma * np.sqrt(T))
    theta = (-S * norm.pdf(d1) * sigma / (2 * np.sqrt(T))
             - r * K * np.exp(-r * T) * norm.cdf(d2)) / 365
    vega = (S * norm.pdf(d1) * np.sqrt(T)) / 100

    return {
        "price": call_price,
        "delta": delta,
        "gamma": gamma,
        "theta": theta,
        "vega": vega
    }

Step 2: Get Market Data

import yfinance as yf

ticker = yf.Ticker("AAPL")
spot = ticker.history(period="1d")["Close"].iloc[-1]

expiry = ticker.options[0]
calls = ticker.option_chain(expiry).calls

👉 This gives you the real option chain

🧠 Step 3: Parse Contracts

def parse_contract_symbol(symbol):
    underlying = symbol[:-15]
    expiry_str = symbol[-15:-9]
    option_type = symbol[-9]
    strike_str = symbol[-8:]

    expiry = pd.to_datetime(expiry_str, format="%y%m%d")
    strike = int(strike_str) / 1000

    return pd.Series([underlying, expiry, option_type, strike])

Step 4: Build Selection Layer (df1)

Contracts = calls[[
    'underlying',
    'parsed_strike',
    'expiry',
    'option_type',
    'contract_size',
    'lastPrice',
    'impliedVolatility',
    'inTheMoney'
]].rename(columns={
    'parsed_strike': 'strike',
    'lastPrice': 'market_price'
})

👉 This is your decision layer (what to trade)

🧠 Step 5: Build Model Inputs (df2)

def build_bs_inputs(row, spot, r=0.05):
    delta = row['expiry'] - datetime.now()
    T = delta.total_seconds() / (365 * 24 * 60 * 60)

    if T <= 0:
        T = 1 / 365

    sigma = row['impliedVolatility']
    if pd.isna(sigma) or sigma < 0.01:
        sigma = 0.2

    return pd.DataFrame([{
        "Spot Price": spot,
        "Strike Price": row['strike'],
        "Time to Expiry": T,
        "Volatility": sigma,
        "Interest Rate": r
    }])

Step 6: Run Model

res = black_scholes_call(
    S=df2["Spot Price"].iloc[0],
    K=df2["Strike Price"].iloc[0],
    T=df2["Time to Expiry"].iloc[0],
    r=df2["Interest Rate"].iloc[0],
    sigma=df2["Volatility"].iloc[0]
)

📊 Output Interpretation

From your notebook:

Price  ≈ 23.86
Delta  = 1.00
Gamma  ≈ 0
Theta  ≈ -0.032
Vega   ≈ 0

🧠 Meaning

• Delta = 1 → behaves like stock
• Gamma ≈ 0 → no convexity
• Vega ≈ 0 → no uncertainty
• Theta small → mostly intrinsic

👉 This is a deep ITM option

📈 Step 7: Delta vs Stock Price

S_range = np.linspace(200, 300, 100)

deltas = [
    black_scholes_call(
        S,
        df2["Strike Price"].iloc[0],
        df2["Time to Expiry"].iloc[0],
        df2["Interest Rate"].iloc[0],
        df2["Volatility"].iloc[0]
    )['delta']
    for S in S_range
]

plt.plot(S_range, deltas)

Interpretation

• Delta moves 0 → 1
• ATM = steep region
• ITM/OTM = flat

👉 Transition from lottery → convex → stock

📈 Step 8: Gamma vs Stock Price

gammas = [
    black_scholes_call(
        S,
        df2["Strike Price"].iloc[0],
        df2["Time to Expiry"].iloc[0],
        df2["Interest Rate"].iloc[0],
        df2["Volatility"].iloc[0]
    )['gamma']
    for S in S_range
]

Interpretation

• Peak at ATM
• Low elsewhere

👉 Maximum convexity happens at the strike

📉 Step 9: Theta vs Time

T_range = np.linspace(0.001, 1, 100)

thetas = [
    black_scholes_call(
        df2["Spot Price"].iloc[0],
        df2["Strike Price"].iloc[0],
        T,
        df2["Interest Rate"].iloc[0],
        df2["Volatility"].iloc[0]
    )['theta']
    for T in T_range
]

plt.plot(T_range, thetas)

Interpretation

• Theta becomes more negative near expiry
• Decay accelerates
• Longer time = slower decay

👉 Time is not linear

Step 10: Theta vs Time (ATM vs ITM vs OTM)

🟦 ATM (Blue)

• Massive decay near expiry
• Highest uncertainty collapse
• 👉 Most dangerous zone

🟧 ITM (Orange)

• Stable decay
• Mostly intrinsic
• 👉 Stock-like behavior

🟩 OTM (Green)

• Small decay early
• Sudden drop near expiry
• 👉 Lottery behavior

📊 Step 11: Theta + Gamma Overlay (MOST IMPORTANT)

🟦 ATM (Battlefield)

• Gamma spikes
• Theta crashes

👉 Move fast → big gain
👉 No move → fast loss

🟧 ITM (Stable)

• Low gamma
• Mild theta

👉 Linear exposure

🟩 OTM (Lottery)

• Low gamma
• Late theta

👉 Nothing → then collapse

🔥 THE BIG INSIGHT

Gamma (reward) and Theta (cost) peak together

🎯 What charts prove

• ATM = high risk, high reward
• ITM = stable, low convexity
• OTM = low probability, asymmetric

⚠️ Critical trading insight

The most attractive trades (high gamma)
are also the ones bleeding the fastest (high theta)

🧠 Final Mental Model

Options = Probability × Time × Volatility

One-line takeaway

You are always paying Theta to own Gamma

Black-Scholes Option Pricing Explained: Black-Scholes Model with Code, Charts & Intuition was originally published in DataDrivenInvestor on Medium, where people are continuing the conversation by highlighting and responding to this story.