McKinney Chapter 4 - Practice - Sec 03

FINA 6333 for Spring 2025

Author

Richard Herron

import numpy as np

%precision 4

'%.4f'

Announcements

Keep signing up for project teams on Canvas > People > Projects
Claim a team, then ask to increase its limit if you want more than four teammates
Keep adding and voting for students choice topics

Five-Minute Review

NumPy is the foundation of scientific computing in Python! NumPy is also the foundation of pandas, so all these tools in this notebook are relevant to pandas, too!

Here are some random data for us to quickly review NumPy with. The np.random.seed() function makes our random number draws repeatable.

np.random.seed(42)
my_arr = np.random.randn(3, 5)
my_arr

array([[ 0.4967, -0.1383,  0.6477,  1.523 , -0.2342],
       [-0.2341,  1.5792,  0.7674, -0.4695,  0.5426],
       [-0.4634, -0.4657,  0.242 , -1.9133, -1.7249]])

Slicing works the same as with lists and tuples. As well, we can replace chained slices, like [i][j][k], with a Matlab notation, like [i. j. k]

my_arr[2][4]

-1.7249

my_arr[2, 4]

-1.7249

What if we want the first two rows?

my_arr[:2]

array([[ 0.4967, -0.1383,  0.6477,  1.523 , -0.2342],
       [-0.2341,  1.5792,  0.7674, -0.4695,  0.5426]])

What if we want the first two columns?

my_arr[:, :2]

array([[ 0.4967, -0.1383],
       [-0.2341,  1.5792],
       [-0.4634, -0.4657]])

What about linear algebra in NumPy?

A = np.array([[1, 2, 3], [4, 5, 6]])
A

array([[1, 2, 3],
       [4, 5, 6]])

B = np.array([[1, 2, 3], [4, 5, 6]])
B

array([[1, 2, 3],
       [4, 5, 6]])

We can use the .dot() method to chain matrix multiplication.

A.dot(B.T)

array([[14, 32],
       [32, 77]])

The @ operator is matrix multiplication with NumPy arrays.

A @ B.T

array([[14, 32],
       [32, 77]])

The * symbol is elementwise multiplication.

A * B

array([[ 1,  4,  9],
       [16, 25, 36]])

There is another way, too! The np.matmul() function is similar to the .dot() method!

np.matmul(A, B.T)

array([[14, 32],
       [32, 77]])

The np.where() and np.select() functions are the NumPy equivalents of Excel’s if().

Say we want to replace values in my_arr greater than 0.5 with 0.5.

my_arr

array([[ 0.4967, -0.1383,  0.6477,  1.523 , -0.2342],
       [-0.2341,  1.5792,  0.7674, -0.4695,  0.5426],
       [-0.4634, -0.4657,  0.242 , -1.9133, -1.7249]])

np.where(my_arr>0.5, 0.5, my_arr)

array([[ 0.4967, -0.1383,  0.5   ,  0.5   , -0.2342],
       [-0.2341,  0.5   ,  0.5   , -0.4695,  0.5   ],
       [-0.4634, -0.4657,  0.242 , -1.9133, -1.7249]])

Or, we could use np.minimum() instead.

np.minimum(0.5, my_arr)

array([[ 0.4967, -0.1383,  0.5   ,  0.5   , -0.2342],
       [-0.2341,  0.5   ,  0.5   , -0.4695,  0.5   ],
       [-0.4634, -0.4657,  0.242 , -1.9133, -1.7249]])

Now, say we want to replace values in my_arr greater than 0 with 0 and greater than 0.5 with 0.5.

np.where(
    my_arr>0.5, # first condition is most restrictive
    0.5, # if first condition True
    np.where( # otherwise
        my_arr>0, # second condition is less restrictive
        0, # if second condition True
        my_arr # otherwise
    )
)

array([[ 0.    , -0.1383,  0.5   ,  0.5   , -0.2342],
       [-0.2341,  0.5   ,  0.5   , -0.4695,  0.5   ],
       [-0.4634, -0.4657,  0.    , -1.9133, -1.7249]])

We have a better way to test more than one condition! np.select()!

np.select(
    condlist=[my_arr>0.5, my_arr>0],
    choicelist=[0.5, 0],
    default=my_arr
)

array([[ 0.    , -0.1383,  0.5   ,  0.5   , -0.2342],
       [-0.2341,  0.5   ,  0.5   , -0.4695,  0.5   ],
       [-0.4634, -0.4657,  0.    , -1.9133, -1.7249]])

Practice

Create a 1-dimensional array `a1` that counts from 0 to 24 by 1.

np.array(range(25))

array([ 0,  1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11, 12, 13, 14, 15, 16,
       17, 18, 19, 20, 21, 22, 23, 24])

a1 = np.arange(25)
a1

array([ 0,  1,  2,  3,  4,  5,  6,  7,  8,  9, 10, 11, 12, 13, 14, 15, 16,
       17, 18, 19, 20, 21, 22, 23, 24])

Create a 1-dimentional array `a2` that counts from 0 to 24 by 3.

a2 = np.arange(0, 25, 3)
a2

array([ 0,  3,  6,  9, 12, 15, 18, 21, 24])

Create a 1-dimentional array `a3` that counts from 0 to 100 by multiples of 3 or 5.

a3 = np.array([i for i in range(101) if (i%3==0) | (i%5==0)])
a3

array([  0,   3,   5,   6,   9,  10,  12,  15,  18,  20,  21,  24,  25,
        27,  30,  33,  35,  36,  39,  40,  42,  45,  48,  50,  51,  54,
        55,  57,  60,  63,  65,  66,  69,  70,  72,  75,  78,  80,  81,
        84,  85,  87,  90,  93,  95,  96,  99, 100])

zero_to_100 = np.arange(101)
a3_alt = zero_to_100[(zero_to_100 % 3 == 0) | (zero_to_100 % 5 == 0)]
a3_alt

array([  0,   3,   5,   6,   9,  10,  12,  15,  18,  20,  21,  24,  25,
        27,  30,  33,  35,  36,  39,  40,  42,  45,  48,  50,  51,  54,
        55,  57,  60,  63,  65,  66,  69,  70,  72,  75,  78,  80,  81,
        84,  85,  87,  90,  93,  95,  96,  99, 100])

(a3 == a3_alt).all()

np.True_

np.allclose(a3, a3_alt)

True

Create a 1-dimensional array `a4` that contains the squares of the even integers through 100,000.

%timeit np.arange(0, 100_001, 2) ** 2
a4 = np.arange(0, 100_001, 2) ** 2
a4

43.3 μs ± 10.3 μs per loop (mean ± std. dev. of 7 runs, 10,000 loops each)

array([          0,           4,          16, ...,  9999200016,
        9999600004, 10000000000], shape=(50001,))

%timeit a4_alt = np.array([i**2 for i in range(100_001) if i%2==0])
a4_alt = np.array([i**2 for i in range(100_001) if i%2==0])
a4_alt

29.5 ms ± 2.69 ms per loop (mean ± std. dev. of 7 runs, 100 loops each)

array([          0,           4,          16, ...,  9999200016,
        9999600004, 10000000000], shape=(50001,))

(a4 == a4_alt).all()

np.True_

np.allclose(a4, a4_alt*1.0000001)

True

Write a function `calc_pv()` that mimic Excel’s `PV` function.

Excel’s present value function is: =PV(rate, nper, pmt, [fv], [type])

The present value of an annuity payment is: \(PV_{pmt} = \frac{pmt}{rate} \times \left(1 - \frac{1}{(1+rate)^{nper}} \right)\)

The present value of a lump sum is: \(PV_{fv} = \frac{fv}{(1+rate)^{nper}}\)

def calc_pv(rate, nper, pmt=None, fv=None, type=0):
    # or we could set the pmt and fv defaults to 0
    if pmt is None:
        pmt = 0
    if fv is None:
        fv = 0
    if type not in [0, 1]:
        raise ValueError(f'type must be 0 or 1; you provided {type}')

    pv_pmt = (pmt / rate) * (1 - (1 + rate)**(-nper))
    pv_fv = fv * (1 + rate)**(-nper)
    pv = pv_pmt + pv_fv
    
    if type == 1:
        pv *= 1 + rate

    return -pv

calc_pv(rate=0.1, nper=10, pmt=10, fv=100)

-100.0000

calc_pv(rate=0.1, nper=10, pmt=10, fv=100, type=1)

-110.0000

Write a function `calc_fv()` that mimic Excel’s `FV` function.

Excel’s future value function is: =FV(rate, nper, pmt, [pv], [type])

def calc_fv(rate, nper, pmt=None, pv=None, type='END'):
    if pmt is None:
        pmt = 0
    if pv is None:
        pv = 0
    if type not in ['BGN', 'END']:
        raise ValueError('type must be BGN or END')
    
    fv_pmt = (pmt / rate) * ((1 + rate)**nper - 1)
    fv_pv = pv * (1 + rate)**nper
    fv = fv_pmt + fv_pv
    
    if type == 'BGN':
        fv *= 1 + rate
        
    return -fv

calc_fv(rate=0.1, nper=10, pmt=10, pv=-100)

100.0000

The rule of 72!

calc_fv(rate=0.072, nper=10, pmt=0, pv=-100)

200.4231

Replace the negative values in `data` with -1 and positive values with +1.

np.random.seed(42)
data = np.random.randn(4, 4)

data

array([[ 0.4967, -0.1383,  0.6477,  1.523 ],
       [-0.2342, -0.2341,  1.5792,  0.7674],
       [-0.4695,  0.5426, -0.4634, -0.4657],
       [ 0.242 , -1.9133, -1.7249, -0.5623]])

We have at least three good solutions!

Slicing and brodacsting
np.where()
np.select()

First, here is the slicing and broadcasting solution.

data[data < 0] = -1
data[data > 0] = +1

data # here "1." and "-1." indicate that these values are floats

array([[ 1., -1.,  1.,  1.],
       [-1., -1.,  1.,  1.],
       [-1.,  1., -1., -1.],
       [ 1., -1., -1., -1.]])

Second, here is the np.where() solution. NumPy’s np.where() function has the same logic as Excel’s if() function! I will recreate data so we start from the same point.

np.random.seed(42)
data = np.random.randn(4, 4)

np.where(data < 0, -1, np.where(data > 0, +1, data))

array([[ 1., -1.,  1.,  1.],
       [-1., -1.,  1.,  1.],
       [-1.,  1., -1., -1.],
       [ 1., -1., -1., -1.]])

The nested np.where() function calls can get confusing! An option is to insert white space!

np.random.seed(42)
data = np.random.randn(4, 4)

np.where(
    data < 0, # condition
    -1, # result if True
    np.where(data > 0, +1, data) # result if False
)

array([[ 1., -1.,  1.,  1.],
       [-1., -1.,  1.,  1.],
       [-1.,  1., -1., -1.],
       [ 1., -1., -1., -1.]])

Third, here is the np.select() solution. NumPy’s np.select() function lets us test many conditions! I will recreate data so we start from the same point.

np.random.seed(42)
data = np.random.randn(4, 4)

np.select(
    condlist=[data<0, data>0],
    choicelist=[-1, +1],
    default=data
)

array([[ 1., -1.,  1.,  1.],
       [-1., -1.,  1.,  1.],
       [-1.,  1., -1., -1.],
       [ 1., -1., -1., -1.]])

Write a function `calc_n()` that calculates the number of payments that generate `x`% of the present value of a perpetuity.

The present value of a growing perpetuity is \(PV = \frac{C_1}{r - g}\), and the present value of a growing annuity is \(PV = \frac{C_1}{r - g}\left[ 1 - \left( \frac{1 + g}{1 + r} \right)^t \right]\).

We can use the growing annuity and perpetuity formulas to show: \(x = \left[ 1 - \left( \frac{1 + g}{1 + r} \right)^t \right]\).

Then: \(1 - x = \left( \frac{1 + g}{1 + r} \right)^t\).

Finally: \(t = \frac{\log(1-x)}{\log\left(\frac{1 + g}{1 + r}\right)}\)

We do not need to accept an argument c1 because \(C_1\) cancels out!

def npmts(x, r, g):
    return np.log(1-x) / np.log((1 + g) / (1 + r))

npmts(0.5, 0.1, 0.05)

14.9000

Write a function that calculates the internal rate of return of a NumPy array of cash flows.

Let us pick a set of cash flows c where we know the internal rate of return! For the following c, the IRR is 10%.

c = np.array([-100, 110])
r = 0.1

First we need an function to calculate net present value (NPV) from c and r! The npv() function below uses NumPy arrays to calculate NPV as: \[NPV = \sum_{t=0}^T \frac{c_t}{(1+r)^t}\]

def calc_npv(r, c):
    t = np.arange(len(c))
    return np.sum(c / (1 + r)**t)

calc_npv(r=r, c=c)

-0.0000

def calc_irr(c, guess=0, step=1e-6, tol=1e-6, max_iter=1_000, verbose=False):
    irr = guess
    for i in range(max_iter):
        npv = calc_npv(r=irr, c=c)
        if abs(npv) < tol:
            if verbose:
                print(f'IRR: {irr:0.4f}\nNPV: {npv:0.4f}\nIterations: {i+1:d}')
            return irr
        deriv = (calc_npv(r=irr+step, c=c) - npv) / step
        irr -= npv / deriv

    raise ValueError(f'NPV did not converge to zero after {i+1:d} iterations')

calc_irr(c=c, verbose=True)

IRR: 0.1000
NPV: 0.0000
Iterations: 4

0.1000

calc_irr(c=np.array([-100, 10, 10, 10, 10, 110]), verbose=True)

IRR: 0.1000
NPV: 0.0000
Iterations: 5

0.1000

# calc_irr(c=np.array([-2000, 1000, -500]))

Write a function `calc_returns()` that accepts NumPy arrays of prices and dividends and returns a NumPy array of returns.

prices = np.array([100, 150, 100, 50, 100, 150, 100, 150])
dividends = np.array([1, 1, 1, 1, 2, 2, 2, 2])

We want to slice our arrays to “lag” or “shift” them! For example, we slice the prices array to calculate capital gains as follows.

prices[1:] - prices[:-1]

array([ 50, -50, -50,  50,  50, -50,  50])

def calc_returns(p, d):
    return (p[1:] - p[:-1] + d[1:]) / p[:-1]

calc_returns(p=prices, d=dividends)

array([ 0.51  , -0.3267, -0.49  ,  1.04  ,  0.52  , -0.32  ,  0.52  ])

Rewrite the function `calc_returns()` as `calc_returns_2()` so it returns NumPy arrays of returns, capital gains yields, and dividend yields.

def calc_returns_2(p, d):
    cg = p[1:] / p[:-1] - 1
    dp = d[1:] / p[:-1]
    r = cg + dp
    return {'r': r, 'cg': cg, 'dp': dp}

calc_returns_2(p=prices, d=dividends)

{'r': array([ 0.51  , -0.3267, -0.49  ,  1.04  ,  0.52  , -0.32  ,  0.52  ]),
 'cg': array([ 0.5   , -0.3333, -0.5   ,  1.    ,  0.5   , -0.3333,  0.5   ]),
 'dp': array([0.01  , 0.0067, 0.01  , 0.04  , 0.02  , 0.0133, 0.02  ])}

Write a function `rescale()` to rescale and shift numbers so that they cover the range [0, 1]

Input: np.array([18.5, 17.0, 18.0, 19.0, 18.0])
Output: np.array([0.75, 0.0, 0.5, 1.0, 0.5])

x = np.array([18.5, 17.0, 18.0, 19.0, 18.0])
x

array([18.5, 17. , 18. , 19. , 18. ])

def rescale(x):
    return (x - x.min()) / (x.max() - x.min())

rescale(x=x)

array([0.75, 0.  , 0.5 , 1.  , 0.5 ])

Write a function `calc_portval()` that accepts a dictionary of prices and share holdings and returns the portfolio value

First convert your dictionary to a NumPy array with one column for prices and another for shares.

data = {
    "AAPL": (150.25, 10),  # (price, shares)
    "GOOGL": (2750.00, 2),
    "MSFT": (300.75, 5)
}

def calc_portval(data):
    x = np.array(list(data.values()))
    return x.prod(axis=1).sum()

calc_portval(data=data)

8506.2500

Write functions `calc_var()` and `calc_std()` that calculate variance and standard deviation.

NumPy’s .var() and .std() methods return population statistics (i.e., denominators of \(n\)). The pandas equivalents return sample statistics (denominators of \(n-1\)), which are more appropriate for financial data analysis where we have a sample instead of a population.

Your calc_var() and calc_std() functions should have a sample argument that is True by default so both functions return sample statistics by default.

def calc_var(x, sample=True):
    sq_err = (x - x.mean()) ** 2
    den = len(x)
    if sample:
        den -= 1
    return sq_err.sum() / den

We can re-use calc_var() in our calc_std() function.

def calc_std(x, sample=True):
    return calc_var(x=x, sample=sample) ** 0.5

np.random.seed(42)
arr = np.random.randn(1_000_000)
arr

array([ 0.4967, -0.1383,  0.6477, ..., -0.113 ,  1.4691,  0.4764],
      shape=(1000000,))

calc_var(arr)

1.0004

calc_std(arr)

1.0002

arr.var(ddof=1) == calc_var(arr)

np.True_

arr.std(ddof=1) == calc_std(arr)

np.True_

Write a function `calc_ret()` to convert quantitative returns to qualitative returns

Returns within one standard deviation of the mean are “Medium”. Returns less than one standard deviation below the mean are “Low”, and returns greater than one standard deviation above the mean are “High”.

np.random.seed(42)
returns = np.random.randn(10)

def calc_ret(r):
    mu = r.mean()
    sigma = r.std(ddof=1)
    return np.select(
        condlist=[
            r<(mu-sigma),
            r<=(mu+sigma),
            r>(mu+sigma)
        ],
        choicelist=[
            'Low',
            'Medium',
            'High'
        ],
        default=''
    )

calc_ret(returns)

array(['Medium', 'Medium', 'Medium', 'High', 'Medium', 'Medium', 'High',
       'Medium', 'Low', 'Medium'], dtype='<U6')

Announcements

Five-Minute Review

Practice

Create a 1-dimensional array a1 that counts from 0 to 24 by 1.

Create a 1-dimentional array a2 that counts from 0 to 24 by 3.

Create a 1-dimentional array a3 that counts from 0 to 100 by multiples of 3 or 5.

Create a 1-dimensional array a4 that contains the squares of the even integers through 100,000.

Write a function calc_pv() that mimic Excel’s PV function.

Write a function calc_fv() that mimic Excel’s FV function.

Replace the negative values in data with -1 and positive values with +1.

Write a function calc_n() that calculates the number of payments that generate x% of the present value of a perpetuity.