NumPy Vectorization: Writing Quant Code That Runs 100x Faster

NumPy is Python’s foundation for numerical computing. Its core is the ndarray—a multi-dimensional array with C-level performance for element-wise operations. For quant finance, where you’re processing millions of price datapoints, NumPy vectorization is essential.

Why Loops Are Slow in Python

Python loops are interpreted. Every iteration has Python overhead. A loop over 1 million elements takes seconds. NumPy operations are compiled in C and operate on entire arrays at once.

Python loop:

import time
import numpy as np

n = 1_000_000
a = np.random.randn(n)

start = time.time()
result = []
for x in a:
    result.append(x * 2 + 1)
loop_time = time.time() - start

start = time.time()
vectorized = a * 2 + 1
vector_time = time.time() - start

print(f"Loop: {loop_time:.3f}s, Vectorized: {vector_time:.3f}s")
# Typical: Loop: 0.4s, Vectorized: 0.003s → 100x faster

The Key Principle: Array Operations

Instead of looping over elements, operate on entire arrays:

# Instead of this:
returns = []
for i in range(1, len(prices)):
    returns.append((prices[i] - prices[i-1]) / prices[i-1])

# Do this:
returns = (prices[1:] - prices[:-1]) / prices[:-1]

# Or this (pandas):
returns = prices.pct_change()

Both produce identical results, but the NumPy version is 50-100x faster.

Broadcasting: Operations Across Different Shapes

NumPy broadcasts arrays of different shapes together automatically:

# prices: shape (1000, 500) — 1000 days, 500 assets
# weights: shape (500,) — 1 weight per asset

# Multiply each day's prices by the weights
portfolio_value = prices @ weights  # (1000,) — sum of price*weight per day

# Or equivalently:
portfolio_value = (prices * weights).sum(axis=1)

Finance-Specific Examples

Computing Returns Matrix

def compute_returns_matrix(prices):
    """
    prices: DataFrame of shape (n_days, n_assets)
    returns: DataFrame of shape (n_days-1, n_assets)
    """
    return (prices.iloc[1:].values - prices.iloc[:-1].values) / prices.iloc[:-1].values

Portfolio Variance

def portfolio_variance(weights, cov_matrix):
    """
    Compute portfolio variance: w' * Cov * w

    weights: shape (n_assets,)
    cov_matrix: shape (n_assets, n_assets)
    """
    return weights @ cov_matrix @ weights

# Or using numpy:
var = np.dot(weights, np.dot(cov_matrix, weights))

Rolling IC (Information Coefficient)

def rolling_ic_vectorized(factor_values, returns, window=252):
    """
    Vectorized rolling IC computation.

    factor_values: shape (n_days, n_assets)
    returns: shape (n_days, n_assets)
    """
    # Align: factor at t predicts return t+1
    fv = factor_values[:-1]  # t = 0 to T-2
    ret = returns[1:]          # t+1 = 1 to T-1

    # Rolling IC using expanding window
    n = len(fv)
    ic = np.full(n, np.nan)

    for i in range(window, n):
        fv_window = fv[i-window:i].flatten()
        ret_window = ret[i-window:i].flatten()
        ic[i] = np.corrcoef(fv_window, ret_window)[0, 1]

    return ic

NumPy Random and Financial Simulation

# Generate correlated asset returns (Cholesky decomposition)
def simulate_correlated_returns(n_days, n_assets, mean_returns, cov_matrix):
    """
    Simulate returns from a multivariate normal distribution.
    """
    L = np.linalg.cholesky(cov_matrix)
    Z = np.random.standard_normal((n_days, n_assets))
    correlated = Z @ L.T + mean_returns
    return correlated

# Geometric Brownian Motion for asset prices
def gbm_prices(S0, mu, sigma, dt, n_paths):
    """
    Simulate stock prices using GBM.

    S0: initial price
    mu: drift (annual)
    sigma: volatility (annual)
    dt: time step (years)
    n_paths: number of simulation paths
    """
    n_steps = int(1 / dt)
    Z = np.random.standard_normal((n_steps, n_paths))
    drift = (mu - 0.5 * sigma**2) * dt
    diffusion = sigma * np.sqrt(dt) * Z
    log_returns = drift + diffusion
    prices = S0 * np.exp(np.cumsum(log_returns, axis=0))
    return prices

Common Pitfalls

Shape mismatches. When broadcasting fails, you get cryptic errors. Always print .shape when debugging.

In-place modification. NumPy arrays can be modified in-place unexpectedly. Use .copy() when you need to preserve data.

Memory for large arrays. A (10000, 5000) float64 array is 400 MB. For large datasets, use float32 or process in chunks.

References

• VanderPlas, Python Data Science Handbook, Ch. 2. O’Reilly, 2016. ISBN: 978-1-491-91205-8.
• Press, Numerical Recipes, Ch. 2. Cambridge, 2007. ISBN: 978-0-521-88068-8.
• McKinney, Python for Data Analysis, Ch. 4. O’Reilly, 2023. ISBN: 978-1-492-07231-9.