McKinney Chapter 4 - Practice - Sec 04

FINA 6333 for Spring 2025

Author

Richard Herron 

import numpy as np
%precision 4
'%.4f'

Announcements

  1. Keep signing up for project teams on Canvas > People > Projects
  2. Claim a team, then ask to increase its limit if you want more than four teammates
  3. 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][3]
-1.9133
my_arr[2, 3]
-1.9133

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 * operator is elementwise multiplication with NumPy arrays.

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

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(my_arr, 0.5)
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.

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 * 1.000_001)
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
30.1 μs ± 6.69 μ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
22.1 ms ± 7.74 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!

  1. Slicing and brodacsting
  2. np.where()
  3. 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, tol=1e-6, step=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} iterations')
calc_irr(c)
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')