Working with Numpy arrays, getting an Error

Question

I'm trying to code functions, that compute financial greeks for B-S model.
I started with something like this:

def greeks_vanilla(S, K, r, q, t, T, sigma):
    import numpy as np
    from scipy.stats import norm

    tau = np.linspace(t, T, num = T)
    d1 = ( np.log(S/K) + (r - q + sigma**2 /2)*tau ) / (sigma * np.sqrt(tau))
    d2 = d1 - sigma*np.sqrt(tau)

    delta_call = np.exp(-q*tau) * norm.cdf(d1)
    delta_put = -np.exp(-q*tau) * norm.cdf(-d1)
    ...

    return {'d1':d1, 'd2': d2, 'delta_call': delta_call, 'delta_put': delta_put, ...}

'...' means, there are more greeks being computed, but it's not important here.
It was working fine, I had reasonable values, nice plots, etc.; however, my teacher told me, that he want to see those values versus not only time (tau on x-axis), but also versus S (S - stock price, on x-axis). In other words, I have to compute greeks for both tau and S changing.

I've tried following:

def greeks_vanilla(S, K, r, q, t, T, sigma):
    import numpy as np
    from scipy.stats import norm

    S = np.linspace(1, S, num = S)

    tau = np.linspace(t, T, num = T)
    d1 = ( np.log(S/K) + (r - q + sigma**2 /2)*tau ) / (sigma * np.sqrt(tau))
    d2 = d1 - sigma*np.sqrt(tau)

    delta_call = np.exp(-q*tau) * norm.cdf(d1)
    delta_put = -np.exp(-q*tau) * norm.cdf(-d1)
    ...

For both versions I initialized following params and run (a-variable):

S = 30.0 #Stock price
K = 50.0 #Strike price
r = 0.05 #Risk-free rate
q = 0.01 #also called delta, annual dividend yield
t = 1.0
T = 100.0
sigma = 0.15
a = greeks_vanilla(S, K, r, q, t, T, sigma)

It works fine in the first case, but when I want to vary S (second function), I'm getting following error:

     File "greeks.py", line 11, in greeks_vanilla
    d1 = ( np.log(S/K) + (r - q + sigma**2 /2)*tau ) / (sigma * np.sqrt(tau))
ValueError: operands could not be broadcast together with shapes (30,) (100,) 

I googled a bit, and it's look like there is something with Numpy library and it's data types (arrays). I'm not very skilled or experienced programmer (still learning), so I wasn't able to figure it by myself.
It looks like it works now only with S == 100 (S = T), but that's not desirable solution.
I've tried:

S = list(S)

But it only outputs:

File "greeks.py", line 11, in greeks_vanilla
    d1 = ( np.log(S/K) + (r - q + sigma**2 /2)*tau ) / (sigma * np.sqrt(tau))
TypeError: unsupported operand type(s) for /: 'list' and 'float'

Please, help me find working solution of this. I don't know if I should try to loop overs S (I tried and failed...), do some data type tricks, calculate it other way (how?) or anything else.

Answer 1

From what I can understand, you got this working for the 1D case, evaluating your function at each value in a np.linspace. Now you're trying to get it to work in 2D and are still using np.linspace.

Numpy wants to work on things element-wise. So when you give it two linspaces that are different sizes, it doesn't know how to map the elements from one linspace to elements in the other. But if you think about what you want, you should realize that you want to evaluate your function at all pair combinations (s, t) where s in S and t in tau, and all those pairs are points on a 2D grid. So the tool you want to be using is np.meshgrid.

Graphically, the way you've written it, you're trying to evaluate your expression only at the #s.

#...................
#...................
#...................
#...................
#...................
####################

and that's why numpy is complaining--it's saying I'm not sure which value in tau corresponds with what value of S. What I'd like to do is either match one value of tau with all of S or pair them up, but that only happens with len(tau)==len(S) as you saw. What you really want is to be evaluating your expression at all points on that grid, and that's what np.meshgrid will let you do.

Numpy's meshgrid can be tricky to understand how to use at first. If you feed it two vectors (e.g. your linspaces), it will return two 2D arrays that specify the coordinates at every grid point spanned by those two vectors. I imagine that you would probably like to do something like:

def greeks_vanilla(S, K, r, q, t, T, sigma):
    import numpy as np
    from scipy.stats import norm

    v1 = np.linspace(1, S, num = S)
    v2 = np.linspace(t, T, num = T)

    S, TAU = np.meshgrid(v1, v2)

    d1 = ( np.log(S/K) + (r - q + sigma**2 /2)*TAU ) / (sigma * np.sqrt(TAU))
    d2 = d1 - sigma*np.sqrt(TAU)

    delta_call = np.exp(-q*TAU) * norm.cdf(d1)
    delta_put = -np.exp(-q*TAU) * norm.cdf(-d1)
    ...

which doesn't raise the error anymore.