pip install numpy scipy matplotlib pandas

import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import math
import logging
from collections import deque
import time
import random

def create_market_state(
    S0: float = 100.0,
    q0: int = 0,
    T: float = 1.0,
    dt: float = 0.01,
    gamma: float = 0.1,
    kappa: float = 1.5,
    sigma: float = 0.2,
    lambda_0: float = 1.0
) -> dict:
    """
    Initializes the market environment and the market maker's initial state.

    Parameters
    ----------
    S0 : float, optional
        Initial mid-price of the asset, defaults to 100.0.
    q0 : int, optional
        Initial inventory of the market maker, defaults to 0.
    T : float, optional
        Total time horizon for the simulation, defaults to 1.0.
    dt : float, optional
        Time step increment for the simulation, defaults to 0.01.
    gamma : float, optional
        Risk aversion parameter, defaults to 0.1.
    kappa : float, optional
        Order book liquidity parameter, related to the sensitivity of order arrival rates to quote deviations, defaults to 1.5.
    sigma : float, optional
        Volatility of the underlying asset price, defaults to 0.2.
    lambda_0 : float, optional
        Baseline market order arrival rate, defaults to 1.0.

    Returns
    -------
    dict
        A dictionary representing the initial market state.
    """
    logging.basicConfig(level=logging.INFO, format='%(asctime)s - %(levelname)s - %(message)s')

    state = {
        'current_time': 0.0,
        'mid_price': S0,
        'inventory': q0,
        'total_time': T,
        'time_step': dt,
        'risk_aversion': gamma,
        'kappa': kappa,
        'volatility': sigma,
        'lambda_0': lambda_0,
        'cash': 0.0, # Initial cash, updated with trades
        'trades': [], # To store trade history
        'inventory_history': [q0],
        'mid_price_history': [S0],
        'reservation_price_history': [],
        'bid_quote_history': [],
        'ask_quote_history': [],
        'spread_history': [],
        'time_history': [0.0]
    }
    logging.info(f"Market state initialized: S0={S0}, q0={q0}, T={T}, dt={dt}")
    return state

def calculate_reservation_price(state: dict) -> dict:
    """
    Calculates the market maker's reservation price based on the Avellaneda-Stoikov model.

    Parameters
    ----------
    state : dict
        Current market state dictionary.

    Returns
    -------
    dict
        Updated state with the calculated reservation price.
    """
    S = state['mid_price']
    q = state['inventory']
    gamma = state['risk_aversion']
    sigma = state['volatility']
    T = state['total_time']
    t = state['current_time']

    # Avellaneda-Stoikov reservation price formula
    reservation_price = S - q * gamma * sigma**2 * (T - t)
    state['reservation_price'] = reservation_price
    # Removed: state['reservation_price_history'].append(reservation_price) -> now handled in run_simulation_step
    logging.debug(f"Time: {t:.2f}, Mid Price: {S:.2f}, Inventory: {q}, Reservation Price: {reservation_price:.2f}")
    return state

def calculate_optimal_spread(state: dict) -> dict:
    """
    Calculates the optimal spread for the market maker's quotes.

    Parameters
    ----------
    state : dict
        Current market state dictionary.

    Returns
    -------
    dict
        Updated state with the calculated optimal spread.
    """
    gamma = state['risk_aversion']
    sigma = state['volatility']
    kappa = state['kappa']
    T = state['total_time']
    t = state['current_time']

    # Avoid division by zero or log of zero if T - t is very small
    time_to_maturity = max(T - t, 1e-9) # Use a small epsilon instead of 0

    # Avellaneda-Stoikov optimal *total* spread formula
    optimal_spread = gamma * sigma**2 * time_to_maturity + (2 / kappa) * np.log(1 + kappa / gamma)

    # Ensure spread is non-negative (can be negative if inventory is heavily skewed)
    optimal_spread = max(0.0, optimal_spread)

    state['optimal_spread'] = optimal_spread
    # Removed: state['spread_history'].append(optimal_spread) -> now handled in run_simulation_step
    logging.debug(f"Time: {t:.2f}, Optimal Spread: {optimal_spread:.4f}")
    return state

def get_optimal_quotes(state: dict) -> dict:
    """
    Determines the optimal bid and ask quotes based on the reservation price and optimal spread.

    Parameters
    ----------
    state : dict
        Current market state dictionary, must contain 'reservation_price' and 'optimal_spread'.

    Returns
    -------
    dict
        Updated state with the calculated optimal bid and ask quotes.
    """
    reservation_price = state['reservation_price']
    optimal_spread = state['optimal_spread']

    optimal_bid = reservation_price - optimal_spread / 2
    optimal_ask = reservation_price + optimal_spread / 2

    state['optimal_bid'] = optimal_bid
    state['optimal_ask'] = optimal_ask
    # Removed: state['bid_quote_history'].append(optimal_bid) -> now handled in run_simulation_step
    # Removed: state['ask_quote_history'].append(optimal_ask) -> now handled in run_simulation_step

    logging.debug(f"Time: {state['current_time']:.2f}, Optimal Bid: {optimal_bid:.2f}, Optimal Ask: {optimal_ask:.2f}")
    return state

import random
import numpy as np # Adding np import for robustness, though it's typically available globally

def simulate_order_arrival(state: dict) -> tuple[str, float | None]:
    """
    Simulates the arrival of market orders (buy or sell).

    Parameters
    ----------
    state : dict
        Current market state dictionary.

    Returns
    -------
    tuple
        A tuple (trade_type, trade_price). trade_type is 'buy', 'sell', or 'none'.
        trade_price is the price at which the trade would occur, or None if no trade.
    """
    S = state['mid_price']
    optimal_bid = state['optimal_bid']
    optimal_ask = state['optimal_ask']
    lambda_0 = state['lambda_0']
    kappa = state['kappa']
    dt = state['time_step']

    # Calculate depths (distances from mid-price)
    delta_b = S - optimal_bid
    delta_a = optimal_ask - S

    # Calculate arrival rates for buy and sell orders
    lambda_buy = lambda_0 * np.exp(-kappa * delta_a) # Market buy hits our ask
    lambda_sell = lambda_0 * np.exp(-kappa * delta_b) # Market sell hits our bid

    # Simulate arrivals using Poisson probabilities P = lambda * dt
    # Add random jitter to simulate real-world arrival times / retry delays if this were a retry mechanism
    # For simple Poisson simulation, we compare random number to probability.
    jitter = random.uniform(-0.01, 0.01) # Small jitter for event probability

    # Simulate market buy order (hits our ask)
    if random.random() < (lambda_buy * dt + jitter):
        logging.debug(f"Time: {state['current_time']:.2f}, Market buy order arrived at Ask: {optimal_ask:.2f}")
        return 'buy', optimal_ask

    # Simulate market sell order (hits our bid)
    if random.random() < (lambda_sell * dt + jitter):
        logging.debug(f"Time: {state['current_time']:.2f}, Market sell order arrived at Bid: {optimal_bid:.2f}")
        return 'sell', optimal_bid

    return 'none', None

def execute_trade(state: dict, trade_type: str, trade_price: float) -> dict:
    """
    Updates the market maker's state (inventory and cash) after a trade.

    Parameters
    ----------
    state : dict
        Current market state dictionary.
    trade_type : str
        The type of trade ('buy' from market maker's perspective means market maker sells, 'sell' means market maker buys).
    trade_price : float
        The price at which the trade occurred.

    Returns
    -------
    dict
        Updated state dictionary after the trade.
    """
    if trade_type == 'buy': # Market order buys from us (we sell from inventory)
        state['inventory'] -= 1
        state['cash'] += trade_price
        logging.info(f"Time: {state['current_time']:.2f}, Sold 1 unit at {trade_price:.2f}. New Inventory: {state['inventory']}, New Cash: {state['cash']:.2f}")
    elif trade_type == 'sell': # Market order sells to us (we buy into inventory)
        state['inventory'] += 1
        state['cash'] -= trade_price
        logging.info(f"Time: {state['current_time']:.2f}, Bought 1 unit at {trade_price:.2f}. New Inventory: {state['inventory']}, New Cash: {state['cash']:.2f}")

    state['trades'].append({
        'time': state['current_time'],
        'type': trade_type,
        'price': trade_price,
        'inventory_after': state['inventory'],
        'cash_after': state['cash']
    })
    state['inventory_history'].append(state['inventory'])
    return state

def update_mid_price(state: dict) -> dict:
    """
    Simulates the movement of the underlying asset's mid-price using Geometric Brownian Motion.

    Parameters
    ----------
    state : dict
        Current market state dictionary.

    Returns
    -------
    dict
        Updated state with the new mid-price.
    """
    S = state['mid_price']
    sigma = state['volatility']
    dt = state['time_step']

    # Geometric Brownian Motion update (assuming zero drift for short-term MM)
    dW = np.random.normal(0, np.sqrt(dt))
    new_S = S + sigma * S * dW # For simplicity, omitting drift mu*S*dt term

    # Ensure price is non-negative
    state['mid_price'] = max(1e-6, new_S)
    state['mid_price_history'].append(state['mid_price'])
    logging.debug(f"Time: {state['current_time']:.2f}, New Mid Price: {state['mid_price']:.2f}")
    return state

def run_simulation_step(state: dict) -> dict:
    """
    Executes a single time step of the Avellaneda-Stoikov market making simulation.

    Parameters
    ----------
    state : dict
        Current market state dictionary.

    Returns
    -------
    dict
        Updated state dictionary after one simulation step.
    """
    state['current_time'] += state['time_step']
    state['time_history'].append(state['current_time'])

    # 1. Update mid-price
    state = update_mid_price(state)

    # 2. Calculate reservation price
    state = calculate_reservation_price(state)
    # Ensure reservation_price history is updated for every step
    state['reservation_price_history'].append(state['reservation_price'])

    # 3. Calculate optimal spread
    state = calculate_optimal_spread(state)
    # Ensure spread history is updated for every step
    state['spread_history'].append(state['optimal_spread'])

    # 4. Get optimal quotes
    state = get_optimal_quotes(state)
    # Ensure bid/ask quote histories are updated for every step
    state['bid_quote_history'].append(state['optimal_bid'])
    state['ask_quote_history'].append(state['optimal_ask'])

    # 5. Simulate order arrival
    trade_type, trade_price = simulate_order_arrival(state)

    # 6. Execute trade if any
    if trade_type != 'none':
        state = execute_trade(state, trade_type, trade_price)

    # Ensure inventory history is updated for every step even if no trade occurred
    if len(state['inventory_history']) < len(state['time_history']):
        state['inventory_history'].append(state['inventory'])

    logging.debug(f"--- End of Step {state['current_time'] / state['time_step']:.0f} --- ")
    return state

def calculate_pnl(state: dict) -> float:
    """
    Calculates the market maker's total Profit and Loss (PnL).

    PnL = Final Cash + Final Inventory Value (at current mid-price).

    Parameters
    ----------
    state : dict
        The final market state dictionary.

    Returns
    -------
    float
        The calculated total PnL.
    """
    final_cash = state['cash']
    final_inventory = state['inventory']
    final_mid_price = state['mid_price']

    final_pnl = final_cash + final_inventory * final_mid_price
    logging.info(f"Final PnL: {final_pnl:.2f}")
    return final_pnl

import logging
import numpy as np
import pandas as pd

# Reset logging level to INFO for cleaner output during demonstration
logging.getLogger().setLevel(logging.INFO)

# Initialize market state
initial_state = create_market_state(
    S0=100.0,
    q0=0,
    T=1.0,
    dt=0.01,
    gamma=0.1,  # Risk aversion
    kappa=1.5,  # Order book liquidity parameter
    sigma=0.2,  # Volatility
    lambda_0=1.0 # Baseline order arrival rate
)

market_state = initial_state.copy()

# Calculate and append initial (t=0) values for history lists
# This ensures all history lists start with a consistent value at t=0
market_state = calculate_reservation_price(market_state)
market_state['reservation_price_history'].append(market_state['reservation_price'])

market_state = calculate_optimal_spread(market_state)
market_state['spread_history'].append(market_state['optimal_spread'])

market_state = get_optimal_quotes(market_state)
market_state['bid_quote_history'].append(market_state['optimal_bid'])
market_state['ask_quote_history'].append(market_state['optimal_ask'])

# Run simulation
num_steps = int(market_state['total_time'] / market_state['time_step'])
print(f"Running simulation for {num_steps} steps...")

for _ in range(num_steps):
    market_state = run_simulation_step(market_state)

print("Simulation finished.")

# Calculate final PnL
final_pnl = calculate_pnl(market_state)
print(f"Market Maker's Final PnL: {final_pnl:.2f}")

# Convert history lists to a DataFrame for easier plotting and analysis
simulation_df = pd.DataFrame({
    'Time': market_state['time_history'],
    'Mid_Price': market_state['mid_price_history'],
    'Inventory': market_state['inventory_history'],
    'Reservation_Price': market_state['reservation_price_history'],
    'Bid_Quote': market_state['bid_quote_history'],
    'Ask_Quote': market_state['ask_quote_history'],
    'Optimal_Spread': market_state['spread_history']
})

import matplotlib.pyplot as plt

# Set a larger figure size for better readability
plt.figure(figsize=(15, 10))

# Plot 1: Mid-Price, Reservation Price, Bid and Ask Quotes
plt.subplot(3, 1, 1)
plt.plot(simulation_df['Time'], simulation_df['Mid_Price'], label='Mid Price', color='blue')
plt.plot(simulation_df['Time'], simulation_df['Reservation_Price'], label='Reservation Price', linestyle='--', color='purple')
plt.plot(simulation_df['Time'], simulation_df['Bid_Quote'], label='Optimal Bid Quote', color='green')
plt.plot(simulation_df['Time'], simulation_df['Ask_Quote'], label='Optimal Ask Quote', color='red')
plt.title('Mid Price, Reservation Price, and Optimal Quotes Over Time')
plt.xlabel('Time')
plt.ylabel('Price')
plt.legend()
plt.grid(True)

# Plot 2: Inventory Level
plt.subplot(3, 1, 2)
plt.plot(simulation_df['Time'], simulation_df['Inventory'], label='Inventory Level', color='orange')
plt.title('Market Maker Inventory Over Time')
plt.xlabel('Time')
plt.ylabel('Inventory')
plt.legend()
plt.grid(True)

# Plot 3: Optimal Spread
plt.subplot(3, 1, 3)
plt.plot(simulation_df['Time'], simulation_df['Optimal_Spread'], label='Optimal Spread', color='brown')
plt.title('Optimal Spread Over Time')
plt.ylabel('Spread')
plt.legend()
plt.grid(True)

plt.tight_layout()
plt.show()

summary_stats = simulation_df[['Mid_Price', 'Inventory', 'Optimal_Spread', 'Reservation_Price']].describe()
print("\nSummary Statistics of Simulation:\n")
print(summary_stats)

import pandas as pd

trades_df = pd.DataFrame(market_state['trades'])
if not trades_df.empty:
    print("\nTrade Activity:\n")
    print(trades_df.head())
    print(f"\nTotal Trades: {len(trades_df)}")
    buy_trades = trades_df[trades_df['type'] == 'sell'] # Market Maker buys
    sell_trades = trades_df[trades_df['type'] == 'buy'] # Market Maker sells
    print(f"Market Maker Bought: {len(buy_trades)} times")
    print(f"Market Maker Sold: {len(sell_trades)} times")
else:
    print("No trades occurred during the simulation.")