QuantLib Usermeeting 2014

ECB calls for lighter treatment of high-quality Asset Backed Securities. How QuantLib might help?

by Michael von den Driesch and Matthias Groncki, 2014

In this notebook we show, how to

  • use our extensions to the QuantLib to setup an ABS Portfolio and run a multithreaded Monte-Carlo simulation of the asset side.

  • model the liability side and write a waterfall in Python.

  • calculate the npv of the tranches.

  • visualise the results of the simulation.

We will setup the following simplified ABS:

  • The asset side consits of a large number of SME loans and all loans have a maturity date in 8 years from now and pay a monthly amount of 550 EUR. The amount will be splitted into 500 EUR amortizing and 50 EUR interest payment.

  • The liability side consits of 3 tranches. The first loss pieces absorbs the first 20 % of the losses, followed by a mezzanine tranche which absorbs the next 20% of losses and the most senior tranch absorbs the remaining 60 %.

  • All interest payments which exceeds the interest on the mezzanine and senior tranche will be used to cover losses from defaults, in particiular the excess flow will be used to redeem the higer seniorty notes.

Remark: This extension of the QuantLib is not part of the official QuantLib. And it's still work in progress. This script is just a draft to demonstrate the functionality of the underlying C++ Library. Before any 'productive' use one should design a more generic framework / toolbox for modelling the liability side.

Disclaimer: THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.

In [1]: 
import os
import numpy as np
import pandas as pd
import datetime as dt
import matplotlib.pyplot as plt
# Import QuantLib 
import IKB.QL.QuantLib as ql
%matplotlib inline

Some helper functions for the visualisation

In [2]: 
def plot_expected_cashflows(amortizing, interest, rocoveries, defaults):
    # prepare data vectors
    nominal = np.max(defaults)+np.sum(amortizing)
    nominal_vector = np.repeat(nominal, len(amortizing))
    nominal_vector = nominal_vector - defaults
    repayment_cum = []
    for i in range(len(amortizing)):
        if i ==0:
            repayment_cum.append(amortizing[i])
        else:
            repayment_cum.append(repayment_cum[i-1] + amortizing[i])
    nominal_vector = [x-y for x, y in zip(nominal_vector, repayment_cum)]
    #Plot
    fig1 = plt.figure(1, figsize=(20,4), dpi=300)
    ind = np.arange(len(amortizing))
    width = 0.5
    p1 = plt.subplot2grid((1,4), (0,0), colspan=3)

    p1 = plt.bar(ind, defaults, width, color='r')
    p2 = plt.bar(ind, nominal_vector, width, color='y', bottom = defaults)
    temp_vc = [x+y for x, y in zip(defaults, nominal_vector)]
    p3 = plt.bar(ind, repayment_cum, width, color='b', bottom = (temp_vc))
    # plot notional of the first loss piece and the mezzanine tranche
    p5 = plt.axhline(y = (91000000), xmin=0, xmax=1, c='black', linewidth=1,zorder=1, hold=None, ls='dashed')
    p6 = plt.axhline(y = (91000000 + 91000000), xmin=0, xmax=1, c='black', linewidth=1,zorder=1, hold=None, ls='dashed')
        
    plt.title('Cashflow evolution', fontsize=14, color='black') 
    plt.tick_params(axis = 'x', which = 'both', bottom = 'off')
    plt.legend((p1[0], p2[0], p3[0]),('Defaults', 'Nominal', 'Amortizing'), loc=6)
    plt.grid(True)
    plt.show()

def visualise_notes(notes):
    NoteNames = ['Senior Tranche','Mezzanine Tranche','First Loss Piece']
    fig = plt.figure(1, figsize=(16,15), dpi=300)
    index_note = np.arange(len(notes[0].PrinPaid))
    width = 0.5
    note_nominal   = []
    note_prinpaid  = []
    note_pdl       = []
    note_pdl_paid  = []
    note_intpaid   = []
    note_prinpaid_cum = []
    for i, note in enumerate(notes):
        note_nominal   = notes[i].PrinEnd
        note_prinpaid  = notes[i].PrinPaid
        note_pdl       = notes[i].PDL
        note_intpaid   = notes[i].IntPaid
        note_pdl_paid  = notes[i].PDLPaid
        note_prinpaid_cum = [notes[i].PrinStart[0]-(x+y) for x, y in zip(note_nominal, note_prinpaid)]
        note_nominal_with_pdl = [(x-y) for x, y in zip(note_nominal,note_pdl)]
        p1 = plt.subplot2grid((5,4),(i,0), colspan=3)
        p1 = plt.bar(index_note, note_pdl, width, color='r')
        p3 = plt.bar(index_note, note_nominal_with_pdl, width, color='y', bottom=note_pdl)
        tempvector = [x+y for x, y in zip(note_pdl, note_nominal_with_pdl)]
        p4 = plt.bar(index_note, note_prinpaid, width, color='b', bottom=tempvector)
        tempvector = [x+y for x, y in zip(note_prinpaid, tempvector)]
        p6 = plt.bar(index_note, note_prinpaid_cum, width, color='c', bottom=tempvector)
        plt.title(NoteNames[i]+' Notional evolution', fontsize=14, color='black') 
        plt.tick_params(axis = 'x', which = 'left', bottom = 'off')

        p29 = plt.twinx()
        p29 = plt.plot(notes[i].IntPaid, color='b', linewidth=2)
        plt.ylabel('Interest')
        plt.legend((p1[0], p3[0], p4[0]),('PDL', 'Nominal', 'Period Repay'), loc=6)
        plt.grid(True)
    plt.grid(True)
    plt.tight_layout()
    plt.show()
    
def plot_note_hist(sim_results, caption="PV"):
    NoteNames = ['Senior Tranche','Mezzanine Tranche','First Loss Piece']
    fig = plt.figure(1, figsize=(10,9), dpi=300)
    for i, npvs in enumerate(sim_results):
        idx = 0 if i < 2 else 1
        idy = i % 2
        p2 = plt.subplot2grid((2,2),(idx,idy))
        p2 = plt.hist(npvs, 30, normed=1, facecolor='g', alpha=0.75)
        mean_pv = np.mean(npvs)
        p7 = plt.axvline(mean_pv, color='r', linestyle='dashed', linewidth=2.3)
        plt.xlabel(caption)
        plt.ylabel('Frequency')
        plt.title(NoteNames[i])
        plt.grid(True)
    plt.tight_layout()
    plt.show()

def plot_note_boxplot(sim_results, caption="PV"):
    NoteNames = ['Senior Tranche','Mezzanine Tranche','First Loss Piece']
    fig = plt.figure(1, figsize=(10,9), dpi=300)
    for i, npvs in enumerate(sim_results):
        idx = 0 if i < 2 else 1
        idy = i % 2
        p2 = plt.subplot2grid((2,2),(idx,idy))
        p2 = plt.boxplot(npvs,0,'gD',1, 0.75)
        plt.xlabel('PV')
        plt.ylabel('')
        plt.title(NoteNames[i])
        plt.grid(True)
    plt.tight_layout()
    plt.show()
In [3]: 
def plot_default_distribution(defaultStock):
    # distribution of defaults
    defaultsMax = []
    for defaultvector in defaultStock:
        defaultsMax.append(np.max(defaultvector))
    mean = np.mean(defaultsMax)
    stdev = np.std(defaultsMax)
    print ('Expected Lifetime Defaults:                '+ str(round(mean,2)))
    print ('Standard Deviation of Lifetime Defaults:   '+ str(round(stdev,2)))
    p22 = plt.hist(defaultsMax, 50, normed=True, facecolor='g', alpha=0.75)   
    #p22 = plt.legend((mean, stdev))
    p22 = plt.xlabel('Lifetime Defaults')
    p22 = plt.ylabel('Frequency')
    p22 = plt.title('Lifetime Defaults Distribution')
    p22 = plt.grid(True)
    plt.show()

Asset-Side Modeling

Setting up ABS Portfolio

Assumptions:

  • We assume that our portfolio consists of 5000 obligors and each obligors has two loans.

  • Each contract pays a fixed amount of 550 EUR per month for the next 91 months.

  • Notional: 455,000,000.00 EUR

  • Global correlation between the assets and the market factor

  • All Obligor have the same rating

Setup the relevant dates and the collectionGrid

In [4]: 
eval_date = ql.Date(28,11,2014)
maturity = ql.Date(30,6,2022)
replenish_end_date = ql.Date(28,11,2014)
replenish_end_time = ql.ActualActual().yearFraction(eval_date,replenish_end_date)
ql.Settings.instance().setEvaluationDate(eval_date)
EvaluationDate = dt.date(eval_date.year(), eval_date.month(), eval_date.dayOfMonth())
maturityABS = dt.date(maturity.year(), maturity.month(), maturity.dayOfMonth())
In [5]: 
grid = ql.Schedule(eval_date, maturity, ql.Period("1M"), ql.NullCalendar(), ql.Unadjusted, ql.Unadjusted, ql.DateGeneration.Forward, False)
actact = ql.ActualActual()
collectionGrid = [actact.yearFraction(eval_date, x) for x in grid]

** Construction of the PD curves**

In [6]: 
# Setup dummy default curves
defaultCurveRepo = ql.AbsRatingManager()
for rating in range(1,8):
    dates = []
    cum_pds = []
    for tenor in range(1,11):
        date = ql.TARGET().advance(eval_date, tenor, ql.Years, ql.Unadjusted)
        pd_ = 1-np.exp(-rating/360*tenor)
        dates.append(date)
        cum_pds.append(pd_)
    defaultCurveRepo.appendDefaultCurve(rating, dates.copy(), cum_pds.copy())

Setup of a portfolio

In [7]: 
# Setup dummy Portfolio
def setup_dummy_portfolio(correlation, recovery, recovery_delay, replenish_end, number_entities, lease_per_entity, num_periods):
    portfolio = ql.AbsDefaultablePortfolio()
    for entity in range(number_entities):
        entity_name = 'FirmNr '+str(entity)
        rating = 7
        entity = ql.AbsDefaultableEntity(entity_name, rating, correlation)
        for lease in range(lease_per_entity):
            interest_rate = 0.04
            leg = []
            for i in range(0,num_periods):
                leg.append(ql.SimpleCashFlow(50.0, ql.Date(1,12,2014)+ql.Period("%iM"%i)))
                leg.append(ql.AmortizingPayment(500.0, ql.Date(1,12,2014)+ql.Period("%iM"%i)))
            lease = ql.Lease(eval_date, ql.TARGET(), 0.0, leg)
            underlying = ql.AbsUnderlying(lease, ql.Date(), recovery, recovery_delay)
            entity.push_back(underlying)
        portfolio.push_back(entity)
    return portfolio
In [8]: 
entity_num = 5000
lease_num = 2
recovery = 0
correlation = 0
rec_delay = 0
portfolio = setup_dummy_portfolio(correlation, recovery, rec_delay, replenish_end_date, entity_num, lease_num, len(collectionGrid))

Monte Carlo Simulation

Define the result vectors

In [9]: 
interests = ql.DoubleMatrix()
amortizing = ql.DoubleMatrix()
recoveries = ql.DoubleMatrix()
defaultStock = ql.DoubleMatrix()
defaultStockinReplenishment = ql.DoubleMatrix()
defaultTimes = ql.Double3DMatrix()

Set a seed for each thread and create the Monte-Carlo-Engine

In [10]: 
# Setup simulatin engine
num_simulations = 10000
num_threads = 8
seeeds = [i*100000000 + 1 for i in range(0,num_threads)]
mce = ql.AbsMonteCarloEngine(defaultCurveRepo, portfolio, num_threads , seeeds, replenish_end_time)

Simulation of the assets

In [11]: 
mce.performSimulationBThreads(collectionGrid, num_simulations, 8, amortizing, interests, recoveries, defaultStock, defaultStockiRneplenishment)

Calculate the expected cashflows

In [12]: 
amortizing_avg = np.mean(amortizing, axis=0)
interests_avg = np.mean(interests, axis=0)
recoveries_avg = np.mean(recoveries, axis=0)
defaults_avg = np.mean(defaultStock, axis=0)
In [13]: 
mc_sim = {}
mc_sim_avg = {}
mc_sim[correlation] = [np.array(amortizing), np.array(interests), np.array(recoveries), np.array(defaultStock) ]
mc_sim_avg[correlation] = [amortizing_avg, interests_avg, recoveries_avg, defaults_avg ]

Visualistation of the simulated cashflows

In [14]: 
plot_expected_cashflows(amortizing_avg, interests_avg, recoveries_avg, defaults_avg)
In [15]: 
plot_default_distribution(defaultStock)

Expected Lifetime Defaults:                31693677.8
Standard Deviation of Lifetime Defaults:   1319842.32

** What will happen if we increase the correlation?**

In [16]: 
entity_num = 5000
lease_num = 2
recovery = 0
correlation = 0.5
rec_delay = 0
portfolio = setup_dummy_portfolio(correlation, recovery, rec_delay, replenish_end_date, entity_num, lease_num, len(collectionGrid))
mce = ql.AbsMonteCarloEngine(defaultCurveRepo, portfolio, num_threads , seeeds, replenish_end_time)
mce.performSimulationBThreads(collectionGrid, num_simulations, 8, amortizing, interests, recoveries, defaultStock, defaultStockinReplenishment)
amortizing_avg = np.mean(amortizing, axis=0)
interests_avg = np.mean(interests, axis=0)
recoveries_avg = np.mean(recoveries, axis=0)
defaults_avg = np.mean(defaultStock, axis=0)
mc_sim[correlation] = [np.array(amortizing), np.array(interests), np.array(recoveries), np.array(defaultStock) ]
mc_sim_avg[correlation] = [amortizing_avg, interests_avg, recoveries_avg, defaults_avg ]
plot_default_distribution(defaultStock)

Expected Lifetime Defaults:                32695699.1
Standard Deviation of Lifetime Defaults:   52304189.23

And if we increase the correlation even further?

In [17]: 
entity_num = 5000
lease_num = 2
recovery = 0
correlation = 1
rec_delay = 0
portfolio = setup_dummy_portfolio(correlation, recovery, rec_delay, replenish_end_date, entity_num, lease_num, len(collectionGrid))
mce = ql.AbsMonteCarloEngine(defaultCurveRepo, portfolio, num_threads , seeeds, replenish_end_time)
mce.performSimulationBThreads(collectionGrid, num_simulations, 8, amortizing, interests, recoveries, defaultStock, defaultStockinReplenishment)
amortizing_avg = np.mean(amortizing, axis=0)
interests_avg = np.mean(interests, axis=0)
recoveries_avg = np.mean(recoveries, axis=0)
defaults_avg = np.mean(defaultStock, axis=0)
mc_sim[correlation] = [np.array(amortizing), np.array(interests), np.array(recoveries), np.array(defaultStock) ]
mc_sim_avg[correlation] = [amortizing_avg, interests_avg, recoveries_avg, defaults_avg ]
plot_default_distribution(defaultStock)

Expected Lifetime Defaults:                33187500.0
Standard Deviation of Lifetime Defaults:   96005454.76
In [18]: 
plot_expected_cashflows(amortizing_avg, interests_avg, recoveries_avg, defaults_avg)

Plot a simulation run without any defaults

In [19]: 
sim = 11
plot_expected_cashflows(amortizing[sim], interests[sim], recoveries[sim], defaultStock[sim])

Plot a simulation run in which all entities defaults at same time

In [20]: 
sim = 1
plot_expected_cashflows(amortizing[sim], interests[sim], recoveries[sim], defaultStock[sim])
In [21]: 
sim = 12
plot_expected_cashflows(amortizing[sim], interests[sim], recoveries[sim], defaultStock[sim])

Liability-Side Modeling

- Totals Liability-Side                     [455.000.000] 
- Senior Tranche                            [273.000.000]  .... EURIBOR 6M + 50bps
- Mezzanine Tranche                         [ 91.000.000]  .... EURIBOR 6M + 100bps
- First Loss Piece                          [ 91.000.000]

Load marketdata

Setup some flat yield curves

In [22]: 
discountCurve = ql.FlatForward(eval_date, 0.003, ql.Actual365Fixed())
fwd_EUR_6_EURIBOR = ql.FlatForward(eval_date, 0.003, ql.Actual365Fixed())

Define some helper classes for storeing the results of the liability side simulation

In [23]: 
class clNote:
    """ 
    Container representing a trance
    """
    def __init__(self, IntRate, PrinStart):
        self.IntRate         = IntRate
        self.PrinStart       = np.zeros(NumPeriods)
        self.PrinPayable     = np.zeros(NumPeriods)
        self.PrinPaid        = np.zeros(NumPeriods)
        self.PrinEnd         = np.zeros(NumPeriods)
        self.IntPayable      = np.zeros(NumPeriods)
        self.IntPaid         = np.zeros(NumPeriods)
        self.Cashflow        = np.zeros(NumPeriods)
        self.PDL             = np.zeros(NumPeriods)
        self.PDLPayable      = np.zeros(NumPeriods)
        self.PDLPaid         = np.zeros(NumPeriods)
        self.excessPDL       = np.zeros(NumPeriods)
        self.Unpaid          = np.zeros(NumPeriods)
        self.excessAmort     = np.zeros(NumPeriods)
        self.Prepays         = np.zeros(NumPeriods)
        self.PrinStart[0]    = PrinStart
        
    def NPV(self, rates_data):
        return np.sum(self.Cashflow * rates_data.DF)
    
    def defaults(self):
        return self.PDL[-1]
    
class Accounts:
    """
    Container collects all cashflows from a simulation
    """
    def __init__(self, Interests, Principals, Recoveries, Defaults):
        self.Interest      = np.array(Interests)
        self.Principal     = np.array(Principals)
        self.Recoveries    = np.array(Recoveries)
        self.cumDefaults   = np.array(Defaults)
        self.ADA           = self.Interest + self.Recoveries
        self.APDA          = self.Principal
        self.prdDefaults   = np.zeros(self.Interest.shape)

class RatesData:
    """
    Interest rate container, contains all needed fixings and
    discount factors
    """
    def __init__(self):
        self.SOPeriod      = np.zeros(NumPeriods)
        self.EOPeriod      = np.zeros(NumPeriods)
        self.euribor_1m    = np.zeros(NumPeriods)
        self.DF            = np.zeros(NumPeriods)
        
class Fee:
    """
    Cashflow container for payable fees
    """
    def __init__(self):
        self.FeePayable    = np.zeros(NumPeriods)
        self.FeePaid       = np.zeros(NumPeriods)
        
class Results:
    """
    Result container contains all npvs and 
    """
    def __init__(self, rates_data):
        self.senTranche  = []
        self.mezzTranche  = []
        self.flp  = []
        self.rates_data = rates_data
    
    def get_expected_cashflows(self, i):
        return (self.senTranche[i].Cashflow,
                self.mezzTranche[i].Cashflow,
                self.flp[i].Cashflow
                )
    
    def get_note(self, i):
        return (self.senTranche[i],
                self.mezzTranche[i],
                self.flp[i]
                )
    
    def get_senior_npvs(self):
        return [x.NPV(self.rates_data) for x in self.senTranche]
    
    def get_mezz_npvs(self):
        return [x.NPV(self.rates_data) for x in self.mezzTranche]
    
    def get_flp_npvs(self):
        return [x.NPV(self.rates_data) for x in self.flp]
    
    
    def get_senior_defaults(self):
        return [x.defaults() for x in self.senTranche]
    
    def get_mezz_defaults(self):
        return [x.defaults() for x in self.mezzTranche]
    
    def get_flp_defaults(self):
        return [x.defaults() for x in self.flp]
    
    
    def defaults(self):
        return [self.get_senior_defaults(),
                self.get_mezz_defaults(),
                self.get_flp_defaults()]
    
    def npvs(self):
        return [self.get_senior_npvs(),
                self.get_mezz_npvs(),
                self.get_flp_npvs()]
    
    def get_expected_value(self, senority):
        if senority == "sen":
            return np.mean(self.get_senior_npvs())
        elif senority == "mezz":
            return np.mean(self.get_mezz_npvs())
        elif senority == "flp":
            return np.mean(self.get_flp_npvs())

Load relevant fixing and discount factors

In [24]: 
def qlDate(dat):
    return ql.Date(dat.day, dat.month, dat.year)


NumPeriods = len(collectionGrid)-1
offset = pd.core.datetools.MonthEnd()
offsetStart = pd.core.datetools.MonthBegin()

NoIssuerEventOfDefault_start = [True  for x in range(0,NumPeriods)]
        
#Defining dataframe with cash accounts
columnNames2 = ['SOPeriod', 'EOPeriod', 'StartDateQl', 'EndDateQl','Interest',
                'Principal','Recoveries', 'ADA','APDA','PDL_cum','Defaults','euribor_1m', 'DF']
absDateRange = pd.date_range(EvaluationDate, maturityABS, freq = 'M')
dfAccounts = pd.DataFrame(0, index = absDateRange, columns = columnNames2)
dfAccounts.EOPeriod = absDateRange
dfAccounts.SOPeriod = dfAccounts.EOPeriod.map(lambda x: offsetStart.rollback(dfAccounts.EOPeriod[x] - pd.DateOffset(months=0)))
dfAccounts.EOPeriod = dfAccounts.EOPeriod.map(lambda x: x.strftime('%Y.%m.%d'))
dfAccounts.SOPeriod = dfAccounts.SOPeriod.map(lambda x: x.to_pydatetime())
dfAccounts.SOPeriod = dfAccounts.SOPeriod.map(lambda x: x.strftime('%Y.%m.%d'))
dfAccounts['EndDateQl'] = dfAccounts.EOPeriod.map(lambda x: qlDate(dt.datetime.strptime(str(x),'%Y.%m.%d')))
dfAccounts['StartDateQl'] = dfAccounts.SOPeriod.map(lambda x: qlDate(dt.datetime.strptime(str(x),'%Y.%m.%d')))
dfAccounts.index = range(0,len(absDateRange))

#forward rate + discount factors
def AssignFwdRate(dfAcc):
    fwdRate = 0.
    fwdRate = fwd_EUR_6_EURIBOR.forwardRate(dfAcc['StartDateQl'],dfAcc['EndDateQl'],ql.Actual360(),ql.Continuous).rate()
    return fwdRate

dfAccounts['euribor_1m'] = dfAccounts[1:].apply(AssignFwdRate, axis=1)
dfAccounts.euribor_1m.loc[0] = 0.005
dfAccounts['DF'] = dfAccounts.EndDateQl.map(discountCurve.discount)

#decomposing the dataframe into the lists :)

clRatesData = RatesData()
clRatesData.SOPeriod      = dfAccounts.SOPeriod.values
clRatesData.EOPeriod      = dfAccounts.EOPeriod.values
clRatesData.euribor_1m    = dfAccounts.euribor_1m.values
clRatesData.DF            = dfAccounts.DF.values

print ('done '+ str(dt.datetime.now()))

done 2014-12-09 16:34:03.426979

Initiate Cashflow Waterfall

In [25]: 
# Cashflow Waterfall
def CFWaterfall(interests, amortizing, recoveries, defaultStock):
    clSenior     = clNote(0.0050, 273000000.00)   
    clMezz     = clNote(0.01, 91000000.00)
    clFLP     = clNote(0.0 , 91000000.00) 
    clFee   = Fee()
    clAccounts = Accounts(interests, amortizing, recoveries, defaultStock)
    for i in range(0,len(clAccounts.Interest)-1):
        if i > 0: clAccounts.prdDefaults[i] = clAccounts.cumDefaults[i] - clAccounts.cumDefaults[i-1]
    iMax=len(interests)-1
    for i in range(0,len(clAccounts.ADA)):
        #loss allocation
        if i==0:
            clFLP.PDL[i] = min(clFLP.PrinStart[i], clAccounts.prdDefaults[i])
            clFLP.excessPDL[i] = -min(0, clFLP.PrinStart[i]-clAccounts.prdDefaults[i])
            clMezz.PDL[i] = -min(clMezz.PrinStart[i], clFLP.excessPDL[i])
            clMezz.excessPDL[i] = -min(0, clMezz.PrinStart[i]-clFLP.excessPDL[i])
            clSenior.PDL[i] = -min(clSenior.PrinStart[i], clMezz.excessPDL[i])
            clSenior.excessPDL[i] = -min(0, clSenior.PrinStart[i]-clMezz.excessPDL[i])
        else:
            clFLP.PDL[i] = clFLP.PDL[i-1] + min(clFLP.PrinStart[i]-clFLP.PDL[i-1],clAccounts.prdDefaults[i])
            clFLP.excessPDL[i] =  -min(0, clFLP.PrinStart[i]-clFLP.PDL[i-1]-clAccounts.prdDefaults[i])
            clMezz.PDL[i] = clMezz.PDL[i-1] + min(clMezz.PrinStart[i]-clMezz.PDL[i-1], clFLP.excessPDL[i])
            clMezz.excessPDL[i] = -min(0, clMezz.PrinStart[i]-clMezz.PDL[i-1]-clFLP.excessPDL[i])
            clSenior.PDL[i] = clSenior.PDL[i-1] + min(clSenior.PrinStart[i] - clSenior.PDL[i-1], clMezz.excessPDL[i])
            clSenior.excessPDL[i] = -min(0, clSenior.PrinStart[i]- clSenior.PDL[i-1]- clMezz.excessPDL[i])
        #fees
        clFee.FeePayable[i] = 100000# side.fee.calculateFees(clAccount, i)
        clFee.FeePaid[i] = max(0, min(clAccounts.ADA[i], clFee.FeePayable[i]))
        clAccounts.ADA[i] = clAccounts.ADA[i]-clFee.FeePaid[i]
        #interest for first tranche
        if i==0:
            clSenior.IntPayable[i] = clSenior.PrinStart[i]*(clRatesData.euribor_1m[i]+clSenior.IntRate)*1/12
        else:
            clSenior.IntPayable[i] = clSenior.PrinStart[i]*(clRatesData.euribor_1m[i]+clSenior.IntRate)*1/12 + clSenior.Unpaid[i-1]
        clSenior.IntPaid[i] = max(0, min(clSenior.IntPayable[i], clAccounts.ADA[i]))
        clSenior.Unpaid[i] = max(0, clSenior.IntPayable[i] - clAccounts.ADA[i])
        clSenior.Cashflow[i] = clSenior.Cashflow[i] + clSenior.IntPaid[i]
        clAccounts.ADA[i] = clAccounts.ADA[i] - clSenior.IntPaid[i]
        #interest for second tranche
        if i==0:
            clMezz.IntPayable[i] = clMezz.PrinStart[i]*(clRatesData.euribor_1m[i]+clMezz.IntRate)*1/12
        else:
            clMezz.IntPayable[i] = clMezz.PrinStart[i]*(clRatesData.euribor_1m[i]+clMezz.IntRate)*1/12 + clMezz.Unpaid[i-1]
        clMezz.IntPaid[i] = max(0, min(clMezz.IntPayable[i], clAccounts.ADA[i]))
        clMezz.Unpaid[i] = max(0, clMezz.IntPayable[i] - clAccounts.ADA[i])
        clMezz.Cashflow[i] = clMezz.Cashflow[i] + clMezz.IntPaid[i]
        clAccounts.ADA[i] = clAccounts.ADA[i] - clMezz.IntPaid[i]
        #redeemption of defaults: most senior tranche
        clSenior.PDLPaid[i] = min(clSenior.PDL[i], clAccounts.ADA[i])
        clAccounts.ADA[i] = clAccounts.ADA[i] - clSenior.PDLPaid[i]
        clAccounts.APDA[i] = clAccounts.APDA[i] + clSenior.PDLPaid[i]
        clSenior.PDL[i] = clSenior.PDL[i] - clSenior.PDLPaid[i]
        #redeemption of defaults: mezzanine tranche
        clMezz.PDLPaid[i] = min(clMezz.PDL[i], clAccounts.ADA[i])
        clAccounts.ADA[i] = clAccounts.ADA[i] - clMezz.PDLPaid[i]
        clAccounts.APDA[i] = clAccounts.APDA[i] + clMezz.PDLPaid[i]
        clMezz.PDL[i] = clMezz.PDL[i] - clMezz.PDLPaid[i]
        #redeemption of defaults: first-loss-piece
        clFLP.PDLPaid[i] = min(clFLP.PDL[i], clAccounts.ADA[i])
        clAccounts.ADA[i] = clAccounts.ADA[i] - clFLP.PDLPaid[i]
        clAccounts.APDA[i] = clAccounts.APDA[i] + clFLP.PDLPaid[i]
        clFLP.PDL[i] = clFLP.PDL[i] - clFLP.PDLPaid[i]
        #interest for the first-loss-piece
        clFLP.IntPayable[i] = clAccounts.ADA[i]
        clFLP.IntPaid[i] = max(0,min(clFLP.IntPayable[i], clAccounts.ADA[i]))
        clFLP.Cashflow[i] = clFLP.Cashflow[i] + clFLP.IntPaid[i]
        clAccounts.ADA[i] = clAccounts.ADA[i] - clFLP.IntPaid[i]
        #amortizing
        # First tranche
        clSenior.PrinPayable[i] = max(0, min(clSenior.PrinStart[i], clAccounts.APDA[i]))
        clSenior.PrinPaid[i] = min(clSenior.PrinPayable[i], clAccounts.APDA[i])
        clSenior.PrinEnd[i] = max(0, clSenior.PrinStart[i] - clSenior.PrinPaid[i])
        if i<iMax: clSenior.PrinStart[i+1] = clSenior.PrinEnd[i]
        clSenior.Cashflow[i] = clSenior.Cashflow[i] + clSenior.PrinPaid[i]
        clAccounts.APDA[i] = clAccounts.APDA[i] - clSenior.PrinPaid[i]
        # Second tranche
        clMezz.PrinPayable[i] = max(0, min(clMezz.PrinStart[i], clAccounts.APDA[i]))
        clMezz.PrinPaid[i] = min(clMezz.PrinPayable[i], clAccounts.APDA[i])
        clMezz.PrinEnd[i] = max(0, clMezz.PrinStart[i] - clMezz.PrinPaid[i])
        if i<iMax: clMezz.PrinStart[i+1] = clMezz.PrinEnd[i]
        clMezz.Cashflow[i] = clMezz.Cashflow[i] + clMezz.PrinPaid[i]
        clAccounts.APDA[i] = clAccounts.APDA[i] - clMezz.PrinPaid[i]
        # Excess flow for the first loss piece
        clFLP.PrinPayable[i] = max(0, min(clFLP.PrinStart[i], clAccounts.APDA[i]))
        clFLP.PrinPaid[i] = min(clFLP.PrinPayable[i], clAccounts.APDA[i])
        clFLP.PrinEnd[i] = max(0, clFLP.PrinStart[i] - clFLP.PrinPaid[i])
        if i<iMax: clFLP.PrinStart[i+1] = clFLP.PrinEnd[i]
        clFLP.Cashflow[i] = clFLP.Cashflow[i] + clFLP.PrinPaid[i]
        clAccounts.APDA[i] = clAccounts.APDA[i] - clFLP.PrinPaid[i]
    return clSenior, clMezz, clFLP

Send generated cashflows through the waterfall

In [26]: 
result_notes = {}

for corr in [0, 0.5, 1]:
    result = Results(clRatesData)
    for i in range(0, len(amortizing)):
        sen,mez,flp = CFWaterfall(mc_sim[corr][1][i],mc_sim[corr][0][i],mc_sim[corr][2][i],mc_sim[corr][3][i])
        result.senTranche.append(sen)
        result.mezzTranche.append(mez)
        result.flp.append(flp)
    result_notes[corr] = result

avg_notes = {}
for corr in [0, 0.5, 1]:
    result = Results(clRatesData)
    sen, mez, flp = CFWaterfall(mc_sim_avg[corr][1],mc_sim_avg[corr][0],mc_sim_avg[corr][2],mc_sim_avg[corr][3])
    result.senTranche.append(sen)
    result.mezzTranche.append(mez)
    result.flp.append(flp)
    avg_notes[corr] = result
In [27]: 
for corr in [0, 0.5, 1]:
    print("Correlation = %f\n" % corr)
    print("NPV senior tranche:", result_notes[corr].get_expected_value("sen"))
    print("NPV mezzanine tranche:", result_notes[corr].get_expected_value("mezz"))
    print("NPV First Loss Piece:", result_notes[corr].get_expected_value("flp"))
    print()

Correlation = 0.000000

NPV senior tranche: 276218968.842
NPV mezzanine tranche: 95892710.762
NPV First Loss Piece: 79445640.304

Correlation = 0.500000

NPV senior tranche: 274590069.355
NPV mezzanine tranche: 91081661.4696
NPV First Loss Piece: 84808081.0592

Correlation = 1.000000

NPV senior tranche: 263724868.968
NPV mezzanine tranche: 86437979.3929
NPV First Loss Piece: 100425003.303

-------------------------------------------------------------------------------------------------------------------------

Plotting results liability side

Evolution of the notional and interest payment of the notes under zero correlation assumption

In [28]: 
visualise_notes(avg_notes[0].get_note(0))
In [29]: 
plot_note_hist(result_notes[0].npvs(), "NPV")
In [30]: 
plot_note_boxplot(result_notes[0].npvs(), "NPV")
In [31]: 
plot_note_hist(result_notes[0].defaults(), "Default")

Correlation = 0.5

In [32]: 
plot_note_hist(result_notes[0.5].npvs(), "NPV")
In [33]: 
plot_note_hist(result_notes[0.5].defaults(), "Exposure at Default")

Plot simulation run in which the defaults exceed the first loss piece

In [34]: 
visualise_notes(result_notes[0.5].get_note(12))

Correlation = 1

In [35]: 
plot_note_hist(result_notes[1].npvs(), "NPV")
In [36]: 
plot_note_hist(result_notes[1].defaults(), "Exposure at Default")
In [37]: 
visualise_notes(result_notes[1].get_note(11))
In [38]: 
visualise_notes(result_notes[1].get_note(1))
In [39]: 
visualise_notes(result_notes[1].get_note(12))