Create a mortgage amortization table in Python

Home ownership is considered one piece of the American Dream by many. Today we'll look at another seemingly integral part of the American Dream: A home mortgage.

Like it or not, most people don't buy their houses with cash, instead having to take out a loan through a bank. A key part of that loan is the amortization schedule, which is a schedule of payments over the life of the loan, showing what portion of each payment is going towards the balance of the loan, and what portion is paying interest on the loan.

For this I created two functions: amortization_table prints out the amortization schedule, while pmt calculates the monthly payment amount.

Note that due to the nature of the beast, getting the numbers to be exactly right could be nigh unto impossible, due to rounding and such. But the idea's the important thing!

Let's put the following in mortgage.py:

# Mortgage amortization

from decimal import *

def amortization_table(principal, rate, term):
    ''' Prints the amortization table for a loan.

    Prints the amortization table for a loan given
    the principal, the interest rate (as an APR), and
    the term (in months).'''

    payment = pmt(principal, rate, term)
    begBal = principal

    # Print headers
    print 'Pmt no'.rjust(6), ' ', 'Beg. bal.'.ljust(13), ' ',
    print 'Payment'.ljust(9), ' ', 'Principal'.ljust(9), ' ',
    print 'Interest'.ljust(9), ' ', 'End. bal.'.ljust(13)
    print ''.rjust(6, '-'), ' ', ''.ljust(13, '-'), ' ',
    print ''.rjust(9, '-'), ' ', ''.ljust(9, '-'), ' ',
    print ''.rjust(9, '-'), ' ', ''.ljust(13, '-'), ' '
    # Print data
    for num in range(1, term + 1):
        
        interest = round(begBal * (rate / (12 * 100.0)), 2)
        applied = round(payment - interest, 2)
        endBal = round(begBal - applied, 2)
        
        print str(num).center(6), ' ',
        print '{0:,.2f}'.format(begBal).rjust(13), ' ',
        print '{0:,.2f}'.format(payment).rjust(9), ' ',
        print '{0:,.2f}'.format(applied).rjust(9), ' ',
        print '{0:,.2f}'.format(interest).rjust(9), ' ',
        print '{0:,.2f}'.format(endBal).rjust(13)

        begBal = endBal
    
def pmt(principal, rate, term):
    '''Calculates the payment on a loan.

    Returns the payment amount on a loan given
    the principal, the interest rate (as an APR),
    and the term (in months).'''
    
    ratePerTwelve = rate / (12 * 100.0)
    
    result = principal * (ratePerTwelve / (1 - (1 + ratePerTwelve) ** (-term)))

    # Convert to decimal and round off to two decimal
    # places.
    result = Decimal(result)
    result = round(result, 2)
    return result

What happens if we run this? To use the amortization_table function, we call it with the loan amount, the interest rate (as an APR), and the term of the loan in months. Let's try a loan for $150,000 at 4% for fifteen years (180 months):

>>> import mortgage
>>> mortgage.amortization_table(150000, 4, 180)
Pmt no   Beg. bal.       Payment     Principal   Interest    End. bal.    
------   -------------   ---------   ---------   ---------   -------------  
  1         150,000.00    1,109.53      609.53      500.00      149,390.47
  2         149,390.47    1,109.53      611.56      497.97      148,778.91
  3         148,778.91    1,109.53      613.60      495.93      148,165.31
  4         148,165.31    1,109.53      615.65      493.88      147,549.66
  5         147,549.66    1,109.53      617.70      491.83      146,931.96
...output deleted...
 175          6,580.75    1,109.53    1,087.59       21.94        5,493.16
 176          5,493.16    1,109.53    1,091.22       18.31        4,401.94
 177          4,401.94    1,109.53    1,094.86       14.67        3,307.08
 178          3,307.08    1,109.53    1,098.51       11.02        2,208.57
 179          2,208.57    1,109.53    1,102.17        7.36        1,106.40
 180          1,106.40    1,109.53    1,105.84        3.69            0.56

And if we check the handy-dandy calculator at BankRate.com, we find that our numbers are quite nearly the same!

But why aren't they the same? Well, I think it's because I'm switching back and forth between Decimal and float. I need to learn how to use Decimal exclusively, as it is more accurate...

What does this mean? We'll probably revisit this later.

Calculating orbital velocity with Python

Earlier we calculated escape velocity for a rocket using Python, and today, we are going to look at calculating orbital velocity using Python. For a background on some of the terms I use, you may want to refer to the above linked article.

Orbital velocity is what the Space
Shuttle needs to reach to get to orbit.

Orbital velocity is very similar to escape velocity, except that it is instead what is necessary to achieve orbit, rather than totally escape a planet's gravity well. That is, if we're trying to get into orbit from the ground (like the Space Shuttle), we calculate orbital velocity. If we're trying to get to Mars (i.e. escape Earth's gravity well), we calculate escape velocity (as well as Mars' escape velocity in that case).

The formula for orbital velocity is as follows:

Vo ≈ sqrt[ GM / r ]

where

Vo = Orbital velocity (m/s)
G = Gravitational constant
M = Mass of planet (kg)
r = Distance from planet's center of gravity the rocket is launching at (m)

We can easily write a Python function to solve for this. In fact, it's almost identical to the one we wrote for calculating escape velocity:

from decimal import *

_GRAVITATIONAL_CONSTANT = Decimal('6.67384E-11')

def orbital_velocity(mass, radius):
    '''Calculates the delta-V orbital velocity.

    Returns delta-V in meters per second for the orbital
    velocity for a given object given its mass in kilograms
    and its radius in meters.'''

    if mass == Decimal('0') or radius == Decimal('0'):
        # User gave a mass or radius which is unrealistic,
        # abort.
        return None
    
    result = Decimal((_GRAVITATIONAL_CONSTANT * mass) / radius)
    result = result.sqrt()

    return result

Let's run this using our figures for Earth (mean radius of 6,371 km, mass of 5.9736E24 kg), assuming the above code is in sci.py:

>>> from decimal import *
>>> import sci
>>> sci.orbital_velocity(Decimal('5.9736E24'), Decimal('6371000'))
Decimal('7910.467706232366495042659538')

So, 7.9 km/s is approximately the orbital velocity according to our function. And, if we check the Wikipedia article for orbital velocity, it mentions that for low Earth orbit, the orbital velocity ranges from 6.9 km/s to 7.8 km/s! That's what I call close enough for government work.

Present and future values of a sum with Python

Not only am I interested in science and computer stuff, but I'm also interested in the worlds of business and finance. That's why today I'll look at computing the present and future values of a lump sum using Python.

Present and future values revolve around the idea of time-value of money: the idea that I'd rather have money now than later. Because of this, when someone says they will give me $1 a year from now, I have to discount that dollar to the present, or calculate its present value, to see how much it's worth right now to me.

Likewise, when calculating the future value of a sum, I'm seeing how much it's worth in the future.

Future value is useful when dealing with compound interest. For example, if I put $1,000 in a savings account right now that's paying 0.8% interest, how much money will be in the account in ten years?

While we could use a financial calculator (like an HP 12C financial calculator [affiliate link]), or one of the plethora of online calculators to solve this, I think it'd be more fun to write a Python function to calculate it for us!

import decimal

def fv(pv, rate, n):
    '''Calculates future value of present value.

    Returns the future value of a sum given the
    present value, a rate, and the number of periods.'''

    result = Decimal(pv * (1 + rate) ** n)

    return result

Note that I used some abbreviations for these variables: fv is the future value, pv is the present value, rate is for the interest rate per period (as a decimal), and n is the number of periods.

The function is very simple, and all it's doing is implementing the equation for the future value of a lump sum: one plus the interest rate, all raised to the number of periods, times the present value.

If we try this out in Python, we get the following (assuming you had the above code in finance.py):

>>> import finance
>>> from decimal import *
>>> finance.fv(Decimal(1000), Decimal(0.008), Decimal(10))
Decimal('1082.942308472838656347119069')

$1,083 isn't a lot of money, so we'll probably want to consider putting our money elsewhere, but the function is the important thing.

Now let's write the function for calculating the present value given a future value:

import decimal

def pv(fv, rate, n):
    '''Calculates present value of future value.

    Returns the present value of a sum given the
    future value, a rate, and the number of periods.'''

    result = Decimal(fv * (1 + rate) ** -n)

    return result

You may notice that this function is basically the inverse of the above function. Instead of raising to the number of periods, we raise to the number of periods as a negative number (thus working backwards), and we multiply it all by the future value as opposed to the present value.

Does it work? Let's try the inverse of the problem we just did. If we have $1,082.94 ten years from now, after compounding at 0.8% annual interest, how much did we start with? Note that due to rounding there may be a slight error.

>>> import finance
>>> from decimal import *
>>> finance.pv(Decimal(1082.94), Decimal(0.008), Decimal(10))
Decimal('999.9978683325782541661697781')

As I rounded the future value to two decimal places when entering it into the present value function, that may have introduced a small amount of error.

These functions only allow us to calculate the value of a lump sum. Later we'll create functions that allow us to calculate the value of a sum plus some stream of payments.

Calculating escape velocity with Python

As a kid I was fascinated by rockets, the ideals of space exploration, and the promise of a brighter tomorrow. While it looks like our dreams are on hold for now with the evisceration of NASA, we can still have some fun by running through some of the necessary calculations ourselves using Python, a free scripting language.

For this, I'm assuming you are already familiar with Python. If not, you will probably want to look at the Python 2.7.x tutorial.

Today we'll be calculating escape velocity. Escape velocity is the speed which is necessary to escape the influence of gravity from a planet or other heavenly body. That is, if a rocket wants to launch from Earth and then leave Earth's orbit, it has to have enough fuel to reach escape velocity.

The equation for escape velocity is as follows:

Ve = sqrt[ (2 * G * M) / R ]

Ve = Escape velocity (m/s)
G = Gravitational constant (N)
M = Mass of the planet/heavenly body (kg)
R = Distance from the planet/heavenly body's center of gravity we're interested in solving for (m)

Note that R is a little funky. If we're launching at sea level, we would just use the radius of the planet as R, but if we're higher up, like on a mountain, the escape velocity will be less.

How can we implement this in Python? To increase re-usability, I decided to make a function to solve for escape velocity, given a planet's mass and radius:

from decimal import *

GRAVITATIONAL_CONSTANT = Decimal('6.67384E-11')

def escape_velocity(mass, radius):
    '''Calculates the delta-V escape velocity.

    Returns delta-V in meters per second for the escape
    velocity for a given object given its mass in kilograms
    and its radius in meters.'''

    if mass == Decimal('0') or radius == Decimal('0'):
        # User gave a mass or radius which is unrealistic,
        # abort.
        return None
    
    result = Decimal((2 * GRAVITATIONAL_CONSTANT * mass) / radius)
    result = result.sqrt()

    return result

Let's work through this piece by piece.

I decided to use the decimal module, as I wanted greater precision than the built-in float was going to give me. Binary floating point numbers can't represent decimals exactly, so I decided that since we're dealing with such large and small numbers here, we had better use a unit type that's designed for that type of work.

On the next line, I define the gravitational constant, using the Decimal type.

Now we get to the actual function. For best results, the user should pass in Decimals for both mass and radius.

Python then checks to see if either mass or radius are zero, which is highly unrealistic, so it aborts. Otherwise it continues.

The function then calculates the result, and returns it as a Decimal.

How does this look in practice? Let's say we were in the Python interpreter, and the previous code was in sci.py in the local directory. We'll calculate the escape velocity for Earth, which according to Wikipedia has a mass of 5.9736E24 kg, a mean radius of 6,371 km (or 6,371,000 m), and an escape velocity of 11.186 km/s (or 11,186 m/s):

>>> import sci
>>> from decimal import *
>>> sci.escape_velocity(Decimal('5.9736E24'), Decimal('6371000'))
Decimal('11187.09071486820095554468233')

As we can see, it is calculating an escape velocity of 11,187 m/s. Variation could be introduced depending on what radius we use for Earth; I used mean radius as opposed to the equatorial radius or polar radius.

Let's try another celestial object, the Moon:

>>> import sci
>>> from decimal import *
>>> sci.escape_velocity(Decimal('7.3477E22'), Decimal('1737100'))
Decimal('2376.108079493511717080823539')

This is based off a mass of 7.3477E22 kg, and a mean radius of 1,737.10 km. Wikipedia gives an answer of 2.38 km/s (or 2,380 m/s), which is again close to our calculation.