A martingale correction method for Monte Carlo simulated interest rates

(1)

A martingale correction method for Monte Carlo

simulated interest rates

Ward van Santen

s2062313

University of Groningen, the Netherlands Faculty of Economics and Business

(2)

Master thesis Econometrics, Operations Research & Actuarial

Studies

Supervisor University of Groningen:

dr. D. Ronchetti

Supervisors Ortec Finance:

Ir. S.N. Singor

(3)

Abstract

(4)

1 Introduction

A typical client of Ortec Finance is an insurance company with a lot of products on the liability side of their balance sheet that are hard to value. These products are often insurance contracts with an uncertain pay-off scheme, that are typically characterized by embedded options and guarantees contractually promised to the policy holder. Certain rights may induce a profit for the policyholders, but never a loss. One can for example think of a minimum return guarantee on the policy-holders’ premiums, or the right to cancel a life insurance contract in exchange for its cash value.

In the valuation of these complex insurance products it is market practice to use so-called risk-neutral scenarios. Monte Carlo simulations are performed in order to generate those different paths.

However, in order for this simulation model to provide a reliable outcome, some properties have to hold. A crucial requirement is that the martingale property is replicated by the model. Generally, this property is not satisfied due to statistical noise in the simulation. In order to reduce the noise and improve the martingale property, the number of scenarios can be increased. However, this increases the computation time, which is highly undesirable.

There exist some methods that improve the martingale property, but these are mainly studied in the context of equity or short-term interest rates. Therefore, the main goal of this paper is to propose a method that effectively improves the martingale property for interest rate projections with a long horizon.

Before entering into more detail, the next section provides the reader with some theoretical background. In this framework the notions of bonds, swaps, and

risk-neutrality are explained. Additionally, the literature on martingale improving

(7)

1.1 Theoretical background

1.1.1 Bonds

First of all, the notion of a bond should be clarified. A bond is basically an agree-ment in which the holder lends an amount of money to the issuer, who is obliged to pay this back at some pre-specified moment in time, the maturity date. This contract contains terms about the interest payments, or coupons, that the issuer owes the holder. The way the coupon payments are made, distinguishes two kinds of bonds: coupon-bearing and zero coupon bonds. In case of a coupon-bearing bond the issuer has the obligation to pay the holder the principal, or its par value or face value, at maturity plus one or more periodically interest payments, most commonly payed semiannually or annually. On the contrary, a holder of a zero coupon bond only receives a single payment, the principal, at maturity.

According to Hull (2008), the theoretical price of a bond is calculated as the present value of all cash flows that are to be received by the holder of the bond. The discount factors are usually based on different interest rates, which are known as zero rates. An illustration will clarify the concepts.

Table 1 presents the zero rates for some given maturities. Furthermore,

sup-pose that the principal isAC1,000 and coupons are paid at a rate of 4% per annum

semiannually. Hence, for a coupon-bearing bond the first coupon payment yields

AC20 and is discounted at 0.9% for six months, the second coupon yields AC20 as

well and is discounted at 1.2% for one year, and so on. The theoretical price for this bond is thus:

20e−0.009·0.5+ 20e−0.012·1+ 20e−0.015·1.5+ 1020e−0.018·2 =AC1,043.16.

Moreover, the theoretical price for a zero coupon bond with a maturity of two years is:

1000e−0.018·2=AC964.64.

Note that continuously compounding is used in this example as well as in the remainder of this paper. For the difference between annually and continuously compounding one can consult Hull (2008).

Table 1: Continuously compounded zero rates.

Maturity (years) Zero rate (%)

0.5 0.9

1 1.2

1.5 1.5

2 1.8

(8)

notation for a zero coupon bond is P (t, T ) and it is defined as

P (t, T ) = N e−r(T −t), (1.1)

where N is its payoff, r is the zero rate, t is the moment at which the bond is valuated, and T is its maturity.

1.1.2 Swaps

According to Hull (2008) a swap is an agreement between two parties to exchange one or more cash flows at pre-specified dates in the future. The contract defines how these cash flows are calculated. The most common swap is a ’plain vanilla’ interest rate swap. Here both parties agree to exchange an amount of cash. Ac-cording to the contract one party is obliged to pay an amount based on a fixed, predetermined interest rate. The counterparty pays an amount based on a floating (uncertain) rate. The resulting cash flows are determined by the underlying no-tional. Hence, one party pays fixed cash flows at some future dates, while receiving variable cash flows at the same moments in time. Obviously, for the other party it is the other way around. The floating rate is often chosen to be the London Interbank Offered Rate (LIBOR).

Analogously to ’normal’ plain vanilla swaps there also exist forward starting inter-est rate swaps. As the notion sugginter-ests, these are swaps that start at some point in the future and have contractual specifications similar to regular swaps. Forward starting swaps play a key role in this study. However, before entering into more detail about these forward starting swaps, regular swaps are treated first.

There are basically two cases in which a party wants to take a position in a plain vanilla swap: one in which it wants to transform a floating obligation into a fixed and one in which it changes a fixed obligation into a floating. An illustration fol-lows in order to clarify both positions.

Suppose there is a party A that has an outstanding loan of AC1,000 for which

it has to pay LIBOR + 0.2% per year. Since this party wants to be sure of its obligations for the next few years, it takes a position in a swap in which it pays a fixed interest rate of 2% per year and where it receives LIBOR. Hence, party A, the so-called fixed-rate payer, has now locked in a fixed payment of 2.2% per

year over its AC1,000 loan and has protected itself against unexpected rises in the

LIBOR.

Suppose there is also a party B that has fixed obligations of 2.1% on a principal

of AC1,000. For some reason, this party expects the LIBOR to decrease in the

next few years. Therefore, it wants to speculate on this by interchanging its fixed obligation by a floating one. Hence, it takes a position as counterparty in the swap described above. As a result, party B, the floating-rate payer, now has a floating

obligation of LIBOR + 0.1% per year on its AC1,000 principal. It is important to

(9)

flows that the two parties owe each other. A schematic representation of this swap is given in Figure 1. Party A Party B LIBOR + 0.2% 2% LIBOR 2.1% Figure 1

In this paper the valuation of interest rate swaps is highly relevant. At this point it is therefore important to devote some attention to this topic. Hull (2008) discusses two different approaches for the valuation of interest rate swaps: one that regards the swap as the difference between two bonds and one that sees it as a portfolio of so-called forward rate agreements (FRAs). In this study the principle of the latter is used.

Forward rates are interest rates that can be locked in today for investments in future time periods. They are implied by current zero rates, the so-called term structure of discount factors. According to Hull (2008), an FRA is an agreement which states that a specific interest rate applies to either borrowing or lending a certain principal during some specified future time period. The underlying as-sumption is that both borrowing and lending would normally be done at LIBOR. Following Brigo and Mercurio (2006), an FRA is a contract that involves three time instants: the current time t, the expiry T > t (the reset date), and the

maturity S > T (the payment date). The holder of the contract receives an

interest-payment for the period [T, S). At the maturity S the payment based on the fixed rate K is exchanged for a payment based on a floating rate L(T, S), which is the spot rate that resets in T and matures at S. Likewise, the contract allows the holder to lock-in the interest rate for the period [T, S) at some desired level K.

Hence, at maturity S the holder receives N τ (T, S)K, where N is the notional and τ (T, S) is the year fraction between T and S. Discounting back to the current time t boils down to

N τ (T, S)KP (t, S),

where P (t, S) is the price of a zero coupon bond starting at t and maturing at S, which serves as discount factor. At the same time, the holder pays

N τ (T, S)L(T, S)P (t, S),

where, following Brigo and Mercurio (2006), L(T, S) is given by

(10)

As a consequence, the total value of an FRA at time t is F RA(t, T, S, τ (T, S), N, K) = N P (t, S)τ (T, S) [K − L(T, S)] = N P (t, S) τ (T, S)K − 1 P (T, S)+ 1 = N [P (t, S)τ (T, S)K − P (t, T ) + P (t, S)] = N P (t, S)τ (T, S) [K − Fs(t; T, S)] , where Fs(t; T, S) = 1 τ (T, S) P (t, T ) P (t, S) − 1

is the simply-compounded forward interest rate that serves as an estimate of the future spot rate L(T, S).

According to Brigo and Mercurio (2006), the interest rate swap is a generalization of the FRA. A forward starting interest rate swap is a contract that exchanges payments between a fixed and a floating leg, starting from some future time

in-stant. Hence, the contract contains pre-specified dates Tα+1, ..., Tβ at which the

fixed leg pays out N · τi· K and the floating leg pays out an amount corresponding

to the actual interest rate L(Ti−1, Ti). In case the fixed leg is paid and the floating

leg is received the forward starting swap is termed a payer forward swap (PFS). Its discounted payoff at time t is given by

β

X

i=α+1

D(t, Ti)N τi L(Ti−1, Ti) − K,

where D(t, Ti) is the stochastic discount factor that discounts for the period [t, Ti].

Conversely, the discounted payoff for a receiver forward swap (RFS) is expressed as

β

X

i=α+1

D(t, Ti)N τi K − L(Ti−1, Ti).

By seeing this last contract as a portfolio of FRAs, an RFS can be valued by:

RF S(t, T, τ, N, K) = β X i=α+1 F RA(t, Ti−1, Ti, τi, N, K) = N β X i=α+1 τiP (t, Ti) K − Fs(t; Ti−1, Ti).

Assuming the above interest rate swap to be fair, or ’at par’, Brigo and Mercurio (2006) derive the following definition of a forward swap rate:

Sα,β(t) =

P (t, Tα) − P (t, Tβ)

Pβ

i=α+1τiP (t, Ti)

(11)

where Tα is the year at which the swap starts and Tβ the maturity of the swap.

The value of an RFS, which is used as the price of a swap in the remainder of this paper, can thus be restated as follows:

RF S(t, T, τ, N ) = N Sα,β(0) − Sα,β(t) β X i=α+1 τiP (t, Ti) = N P (0, Tα) − P (0, Tβ) Pβ i=α+1τiP (0, Ti) − P (t, Tα) − P (t, Tβ) Pβ i=α+1τiP (t, Ti) ! _β X i=α+1 τiP (t, Ti). (1.3) Note that the definition of the value of a forward starting swap can serve as the

value of a swap starting today by initiating Tα = 0.

Analogously to zero coupon bonds there also exist zero coupon swaps. In con-trast to the plain vanilla swaps described above, zero coupon swaps only involve one single payment at the maturity date. Therefore, the value of a forward starting zero coupon swap is given by:

RF ZCS(t, T, τ, N ) = N P (0, Tα) − P (0, Tβ) τ P (0, Tβ) − P (t, Tα) − P (t, Tβ) τ P (t, Tβ) ! τ P (t, Tβ). (1.4)

Note that the valuation method used above and in the remainder of this study is termed DCF (Discounted Cash Flow). After the financial crisis so-called OIS (Overnight Indexed swap) discounting, or dual curve discounting, has emerged. Siliadin (2013) describes this as the practice of using one interest rate curve to project the cash flows of an interest rate swap and another curve to discount them. Since Ortec Finance does not use dual curve discounting and an implemen-tation of this method would highly complicate this study, a formal treatment of this concept is outside the scope of this paper. The reader can consult Mercurio (2009) for a full treatment.

1.1.3 Risk-neutrality

(12)

Simply put, a risk-neutral investor is neither risk-seeking nor risk-averse and thus indifferent to risk when making a decision about an investment. This means that he only looks at the potential gains of an investment. For example, a risk-neutral

investor is indifferent between an investment with a 100% certain gain of AC1,000,

a 50% chance of gettingAC2,000 (and AC0 otherwise), and a 10% probability of

re-ceivingAC10,000 (andAC0 otherwise), since the expected return of all three options

is AC1,000.

The concepts of risk-neutral probabilities and a risk-neutral measure are best ex-plained by use of an example. However, before introducing this illustration, the notion of an option should be treated.

An option is an example of a (financial) derivative, which is the collective term for financial products that derive their value from an underlying entity. There are lots of variations on options, but for now a general introduction will do. The two most commonly traded types of options are call and put options. A call option gives the holder the right, but not the obligation to buy the underlying asset for a fixed price, the strike price. Similarly, a put option gives the holder the right, but not the obligation to sell the underlying asset for a fixed price. Both options contain a contractual expiration date, the maturity. Furthermore, there is a distinction between American options, that can be exercised at any time before maturity, and European options, that can only be exercised at maturity. In case of a European call option the holder makes a profit when the price of the underlying asset is higher than the strike price, because the holder can then buy the asset for the strike price and immediately sell it for a higher price in the market. Similarly, an investor who holds a European put option makes a profit when the actual price is lower than the strike price.

A so-called binomial tree, as presented in Figure 2, is used to clarify an (overly simplified) option pricing example. Suppose an investor has a call option on a

stock S with a strike price equal to AC110. Furthermore, assume that the risk-free

interest rate, which is the theoretical rate of return of an investment with zero risk, is 1.2%. In this one-step case, the price of S can either go up or down. As

the figure indicates, the initial stock price S0 is AC100 and it either increases to

AC120 after an upward movement or decreases to AC80 after an downward

move-ment. Hence, the value of the option isAC10 after the stock price moves up, or AC0

when it moves down.

Now define p as the probability of an upward move and hence 1 − p as the prob-ability of a downward move. Following Hull (2008), one can now argue that the expected return on the stock in a risk-neutral world must be the risk-free rate of 1.2%. Hence, assuming a one-year period, the following equation should be satisfied in the risk-neutral world:

(13)

Solving this equation yields p = 0.47. Hence, the expected value of the call option is now:

0.47 · 10 + 0.53 · 0 =AC4.70.

Since this is the value in one year, this should be discounted with the risk-free rate to obtain the value today:

4.70e−0.012·1=AC4.65.

Note that p and 1 − p are risk-neutral probabilities.

S0 = 100 S1 = 80 c = 0 S1 = 120 c = 10 1 − p p Figure 2

Again, the above situation is an extremely simplified example of risk-neutrality and only serves as intuitive support. Nevertheless, it enables the opportunity to introduce some more general statements concerning risk neutrality.

First of all, risk-neutral probabilities can be seen as probabilities of future states which are adjusted for risk. These probabilities can for example be used to calcu-late asset values. The benefit of risk-neutral pricing is that once the corresponding probabilities are computed, these can be used to price any asset based on its ex-pected payoff.

Furthermore, it is important to point out two of the key assumptions when working in a risk-neutral world: the absence of arbitrage and the complete market assump-tion. Roughly speaking, the absence of arbitrage is equivalent to the impossibility of making a risk-less profit. This implies that there cannot be two portfolios with the same payoff at some future date, but a different price today. According to Brigo and Mercurio (2006), a market is complete if and only if every contingent claim is attainable. This basically means that any asset’s payoff structure can be replicated by a portfolio of other assets.

(14)

these two probability measures are equivalent, denoted by Q ∼ P, if they agree on the events with probability zero and the events with probability one.

It is now time to introduce the central concept of this paper: a martingale.

Follow-ing Bj¨ork (2009), a given discrete-time process {Xn: n = 0, 1, ...} is a discrete-time

martingale with respect to the conditioning set Fn, denoting the information

avail-able at time n, if

E[Xm|Fn] = Xn, for all n ≤ m. (1.5)

Moreover, Bj¨ork (2009) presents an economically more easy interpretable theorem:

assuming that there exists a risk-free asset S and the risk-free rate is r, then the market is free of arbitrage if and only if there exists a measure Q ∼ P such that

S0 =

1

1 + rE

Q_[S

1], (1.6)

where S0 and S1 are the asset prices at t = 0 and t = 1, respectively.

Equation (1.6) is also called the risk-neutral pricing formula. Its economic inter-pretation is that the asset prices today are obtained as the discounted expected value of the assets prices tomorrow. Note that while one should discount with the risk-free rate, in real-life ’the’ risk-free rate is generally hard to find.

Once an asset satisfies the risk-neutral pricing formula, such an asset is said to possess the martingale property. As indicated in the introduction, in Monte Carlo simulations this property is often not perfectly replicated, which is highly undesir-able. The literature offers several methods that improve this martingale property. The most relevant are discussed in the next section.

1.1.4 Methods on improving the martingale property

Before entering the realm of variance reduction techniques and correction methods, the reader should be familiarized with the principle of Monte Carlo simulation. This algorithm is generally used to solve problems that are analytically unsolvable. More specifically, it is employed as a tool for the discretization of continuous time processes, which is exactly what it is used for in this study.

(15)

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 Time in years

Simulated paths (N = 100): Standardized stock price (S/B)

S/B

Figure 3

In case of an underlying risk-neutral model, one can test whether the generated val-ues satisfy the martingale property. Using Monte Carlo simulation, this property mostly does not hold due to statistical noise. In general, this property improves as the number of scenarios increase. However, higher numbers of scenarios are in practice undesirable, since this increases the simulation time. There are several methods in the literature that are able to improve the martingale property. Variance reduction techniques are a powerful tool for increasing the accuracy of simulations. A well-known and commonly-used technique is antithetic sampling. This method can be applied when Gaussian variates drive a Monte Carlo

sim-ulation. According to J¨ackel (2001), one can make use of the fact that for any

drawn path its mirror image has equal probability. Without going into too much detail, in practice this means that when drawing random numbers from a standard normal distribution one can constantly take the ’counterpart’ of a drawn number, while neither affecting the mean nor the variance of the underlying distribution.

Hence, if an evaluation driven by a Gaussian variate vector draw zi is represented

by vi = v(zi), one also uses ˜vi = v(−zi). Likewise, the underlying mean of the

random number is then always exactly equal to zero, which improves the accuracy of the simulation and therefore the martingale property. Note that antithetic sam-pling is a rather basic tool which is widely used in Monte Carlo studies.

(16)

1. Simulate asset prices ˆSi(tj) for i = 1, 2, ..., N with N the number of paths

and j = 1, 2, ..., M , where M is the number of time steps. Note that ˆSi(tj)

denotes the asset price without EMS adjustment.

2. Take the simulated return from tj−1 to tj:

ˆ Si(tj)

ˆ Si(tj−1)

.

3. Create a temporary asset price at time tj: Zi(tj, n) = Si∗(tj−1, n)

ˆ Si(tj)

ˆ Si(tj−1)

,

where S_i∗(tj−1, n) represents the EMS-corrected asset price of the previous

time step.

4. Calculate the discounted sample average Z0(tj, n) =

1

Ne

−r·tjPN

i=1Zi(tj, n).

5. Calculate the EMS-corrected asset price at tj by Si∗(tj, n) = S0

Zi(tj, n)

Z0(tj, n)

.

6. Move on to tj+1 and repeat the steps above.

The authors show that the EMS method contributes to a substantial error re-duction in the price estimates of several derivatives, such as the earlier discussed European call options.

In the papers of Van Haastrecht et al. (2009) and Andersen (2007) a correc-tion method is proposed that enforces the martingale property on asset prices.

The former authors design martingale corrections for the Sch¨obel-Zhu model and

some extensions, whereas Andersen (2007) introduces martingale corrections for different schemes applied to the Heston stochastic volatility model. Despite the fact that both papers study different models, their martingale corrections are sim-ilar. They both make use of, how Van Haastrecht et al. (2009) introduces it, the exponentially affine in expectation property. The idea of the martingale correction of both papers is basically as follows.

The expected value of some Monte Carlo simulated asset ˆX generally does not

satisfy the martingale property. Hence,

M 6= E[eA· ˆX(t+∆)| ˆX(t)],

where M is the corresponding martingale value and A can have any real-valued form. Due to this convenient form it is easy to adjust A such that the above equation holds. Thus,

M = E[eA∗· ˆX(t+∆)| ˆX(t)],

where A∗ includes the martingale correction.

(17)

The studies of Van Haastrecht et al. (2009) and Andersen (2007) are inspired by the paper of Glasserman and Zhao (2000). They developed a similar method which enforces martingale properties on simulated LIBOR and swap rates in order to keep the discretized model arbitrage-free and thus risk-neutral. The LIBOR and swap rate models used in this study are examples of so-called forward rate models. These are models that are typically used for ’short’ time horizons, say up to a maximum of 30 years. The numerical results suggest that their proposed method generally gives better results than more standard discretizations.

1.2 Contribution and research question

In the previous section several methods on improving the martingale property in a Monte Carlo simulation framework were discussed. General variance reduction techniques such as antithetic sampling are able to improve martingale properties for a wide range of models. However, their power is limited. The methods pro-posed by Glasserman and Zhao (2000), Andersen (2007), and Van Haastrecht et al. (2009) are effective tools in improving martingale properties. Their techniques are explicitly developed for either equity or forward rate models.

Hence, the literature offers some relevant martingale correction methods. However, there are no papers that specifically study martingale improvements for short rate models. A crucial feature that distinguishes forward rate and short rate models is their horizon. The former is more suitable for short-term projections, whereas the latter is typically used for longer horizons. In times where interest rates have become increasingly relevant, it seems surprising that there is such a gap in the literature.

Therefore, the goal of this paper is to propose a martingale correction method for a widely-used short rate model: the two-factor Hull-White model. The tech-nique should have the flexibility of being applicable to a wide range of interest rate dependent products. As Glasserman and Zhao (2000), Andersen (2007), and Van Haastrecht et al. (2009) made promising contributions to the literature by developing martingale correction methods for equity and forward rate models, the method proposed in this paper is inspired by their studies.

The research question that serves as guidance in this paper is as follows:

Is it possible to derive a martingale correction method for the two-factor Hull-White short rate model that improves martingale properties of Monte Carlo simulated interest rates?

2 Two-factor Hull-White model

(18)

because short rate models are more suitable for projecting interest rates for a long horizon. This is relevant for insurance companies, who are typically interested in the long term. In this context, one should think of a long horizon as a horizon of 100 years or even longer. Moreover, Hull and White (1994) argue that the two-factor Hull-White model provides a richer pattern of term structure movements and volatility structures than its one-factor variant.

The literature offers multiple alternative short rate models. However, the two-factor Hull-White model possesses some favorable features. Firstly, it is possible due to the underlying Gaussian distributions to derive explicit formulas for a num-ber of financial instruments, such as bonds and interest rate derivatives. More specifically, the Hull-White model offers an analytical formula for the price of a zero coupon bond. As the remainder of this paper shows, this proves to be a highly convenient feature. Finally, the different model parameters provide a large degree of flexibility and give insight into the dynamic behavior of the term structure. In their paper Hull and White (1994) propose the two-factor Hull-White model, which can be represented by the following processes:

dr(t) = [θ(t) + u(t) − ar(t)]dt + σ1dz1(t) (2.1a)

du(t) = −bu(t)dt + σ2dz2(t), (2.1b)

where dt is an infinitesimal time difference, a and b are mean reversion parameters,

σ1 and σ2are volatility parameters, dz1 and dz2 are so-called Wiener processes

cor-related with ρ, u(0) = 0, and r(0) is based on the underlying term structure as will be discussed in section 5.1. Furthermore, θ(t) is a deterministic parametric function of time and is chosen to make the model consistent with the initial term structure. In this study the two-factor Hull-White model is taken as the model for the short rate under the risk-neutral probability measure.

The correlated Wiener processes dz1 and dz2 are determined by applying the

Cholesky decomposition to uncorrelated Wiener processes dW1 and dW2. This

gives: dz1 = dW1 (2.2a) dz2 = ρdW1+ p 1 − ρ2_dW 2, (2.2b) where dW1, dW2 ∼ N (0, 1).

The parametric function θ(t) is calculated by

θ(t) = Ft(0, t) + aF (0, t) + φt(0, t) + aφ(0, t), (2.3)

(19)

is thus today’s observed forward rate. There will be elaborated more on this in Chapter 5. Furthermore, φ(t, T ) is given by

φ(t, T ) = 1 2σ 2 1B(t, T ) 2 +1 2σ 2 2C(t, T ) 2 + ρσ1σ2B(t, T )C(t, T ), (2.4)

where B(t, T ) and C(t, T ) are defined as

B(t, T ) = 1 a1 − e −a(T −t) (2.5a) C(t, T ) = 1 a(a − b)e −a(T −t)₋ 1 b(a − b)e −b(T −t)₊ 1 ab. (2.5b)

The two-factor Hull-White model has the attractive feature that it has an analytic formula for the price of a zero coupon bond. According to Hull (2014), the value at time t for a zero coupon bond with a payoff of 1 at time T is

P (t, T ) = A(t, T )e−B(t,T )r(t)−C(t,T )u(t), (2.6)

where r(t) and u(t) are obtained from the Hull-White model, B(t, T ) and C(t, T ) are as defined by (2.5a) and (2.5b), respectively, and A(t, T ) is given by (A.1) in Appendix A.

The reader can consult Hull and White (1994) for a derivation of the analytical formula of the zero coupon bond price and the corresponding functions.

3 Methodology

3.1 Simulation framework

The two-factor Hull-White model is a continuous-time stochastic progress. The literature offers multiple discretization schemes that transform a continuous-time into a discrete-time processes. In this paper the most basic method, an Euler scheme, is implemented. The study of other discretization schemes is outside the scope of this paper.

Following Rouah (2011) on the Euler scheme, the discrete-time approximations of (2.1a) and (2.1b) have the following form:

r(k + 1) = r(k) + [θ(k) + u(k) − ar(k)]∆t + σ1

√

∆tdz1(k) (3.1a)

u(k + 1) = u(t) − bu(k)∆t + σ2

√

∆tdz2(k), (3.1b)

(20)

6.3 more attention will be devoted to this choice. Note that θ(k) is obtained using (2.3).

The random numbers are drawn in Matlab using the randn function, which gen-erates pseudo-random numbers that are standard normally distributed. Pseudo-random numbers are generated by an algorithm that produces a sequence whose properties approximate those of ’real’ random numbers. In this simulation frame-work the Mersenne Twister pseudo-random number generator is used, which is the default algorithm in Matlab. More specifically, Matlab uses the MT19937 version,

which is able to generate random sequences with a period of 219937−1. This version

is widely used and the default algorithm for many software systems.

In order to speed up the process the simulation is vectorized. This means that the random numbers are drawn as vectors, which enables the possibility of simulating N paths simultaneously. Hence, the vectorized discrete-time approximation has the following form:

r(k + 1) = r(k) + [θ(k) + u(k) − ar(k)]∆t + σ1

√

∆tdz1(k)

u(k + 1) = u(k) − bu(k)∆t + σ2

√

∆tdz2(k),

where the bold symbols correspond to N -dimensional vectors. Note that vector θ(k) contains N identical scalars obtained for step k.

Furthermore, it is important to note that vectors dz1 and dz2 are obtained using

antithetic sampling. This is easy to implement in Matlab when the simulation is vectorized. First of all, using the function randn, one draws a random matrix Wh ∈ RN/2×2 with N the number of scenarios. Note that an N/2 × 2 random

matrix is drawn since there are two stochastic processes that drive the Hull-White model. In the second step, one combines this matrix with its mirror image to obtain

W = Wh

−Wh

,

such that W ∈ RN ×2_{contains the standard normally distributed random numbers}

for a single discrete time step. Denote the first and second column by W1 and

W2, respectively. Hence, one has W1, W2 ∈ RN, whose means are now exactly

equal to zero. Using (2.2a) and (2.2b), dz1 and dz2 are thus obtained as:

dz1 = W1

dz2 = ρW1+

p

1 − ρ2_W

2.

(21)

3.2 Martingale tests

Interest rates are an important part of Ortec Finance’s risk-neutral scenario mod-els. Monte Carlo simulations are used to model the interest rates and this pro-cedure typically does not deliver martingales. Therefore, it is required that the martingale property is at least replicated ’well enough’. This means that the mar-tingale property has to hold within certain bounds.

Therefore, several martingale tests are designed in order to carry out tests on the simulated interest rates. The martingale tests not only serve as a crucial tool in validating the risk-neutrality of scenarios, but they also enable one to detect possi-ble bugs in the code and study the effect of both stochastic noise and discretization. In the remainder of this section the framework of martingale tests is discussed. First of all, a martingale test for zero coupon bonds is constructed. Secondly, martingale tests for forward starting swaps with different maturities are designed. In order to give the reader some more grip on the methodology, the theory is ac-companied with some visualizations of the martingale tests. Note that the purpose of this section is not to elaborate on how these figures and tables are obtained. A detailed treatment of how such results are acquired, is provided in the chapters that follow.

3.2.1 Martingale test for zero coupon bonds

In this test the simulated zero rate is compared with the zero rate corresponding to the underlying term structure. Due to stochastic noise in Monte Carlo simu-lations, these curves typically do not perfectly match. In the remainder of this section the procedure for obtaining the zero rates from the two-factor Hull-White model is discussed.

First of all, the value at t = 0 of a zero coupon bond with a payoff of 1 and maturity T is determined. Analogous to (1.1), this is given by:

π0(T ) = E[D(0, T )|F0], (3.3)

where

D(0, T ) = e−R0Tr(s) ds (3.4)

is the stochastic discount factor and F0 is the conditioning set at t = 0 as used in

equation (1.5). In discretized form this boils down to:

(22)

is the value of the discretized discount factor for the period [0, T ] in scenario i. Here N denotes the number of scenarios, k = T /∆t is the discretization step corresponding to maturity T , and r(i, j) is the short rate at step j in scenario i. Furthermore, note that the conditional expectation in equation (3.3) is estimated by just a sample average, without including any conditioning information set at t = 0.

Hence, π0(k) is calculated for each step k and used to estimate the model implied

zero rate as follows:

rmodel(T ) = −

ln(π0(k))

T (3.7)

for k = T /∆t and T = ∆t, 2∆t, ..., H years. Note that with slight abuse of no-tation, the dependence of k on T will be suppressed in the remainder of this paper. Finally, threshold bounds are defined in order to test whether the martingale property holds within some desirable range. The bounds are constantly linearly increasing over time, such that less deviation is allowed for in an earlier stage. As a benchmark the bounds are set at 0.01% for t = 0 and at 0.05% for the time horizon. Hence, in that case a deviation of 0.05% from the underlying term struc-ture is allowed for. This is graphically represented by Figure 4. Note that one can modify the threshold bounds to any desirable level.

0 10 20 30 40 50 60 70 80 90 −0.4 −0.3 −0.2 −0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 Maturity (years)

Benchmark: Zero rates and threshold bounds, July 2016 (10000 scenarios)

Interest rate (%)

Observed zero rate Simulated zero rate Thresholds

Figure 4

(23)

goal is to test whether the fit of the zero rate curve improves as the number of scenarios increases. Therefore, a consistent measure has to be used in order to properly test the performance. Since the width of the confidence interval decreases as the number of scenarios increases, confidence bounds are obviously no consis-tent tool for measuring the performance of a Monte Carlo simulated zero curve with respect to the underlying term structure.

The martingale test is evaluated by determining the number of simulated points of the zero curve that lie between the thresholds and by calculating the corresponding percentage over all simulated points. For example, in the above martingale test the simulated curve is between the thresholds for 99.37% of the time.

In order to evaluate the performance of the martingale test more thoroughly, some numerical results corresponding to Figure 4 are presented in Table 2. As it turns out later, it also serves as means for comparing this benchmark model with the corrected model.

Table 5 shows the difference and the absolute difference between the observed zero rates and the simulated zero rates. These differences are obtained for each time step in the simulation. The mean, minimum, and maximum are calculated over the whole time period and serve as clarification of the deviation from the observed zero curve.

Table 2: Hull-White benchmark

Performance ZCB test: difference observed and simulated zero rates (in 1/100 %).

Diff. Abs. diff.

Mean −5.45 · 10−5 _{1.56 · 10}−4

Minimum −5.02 · 10−4 5.04 · 10−9

Maximum 1.44 · 10−4 5.02 · 10−4

July 2016, 10,000 scenarios

In order to give the reader some intuition about how these differences evolve with respect to the maturity, a graphical representation is given by Figure 14 in Ap-pendix B.

3.2.2 Martingale tests for forward starting swaps

(24)

Using (1.3), the value at t = T of a d-year contract that starts at T is given by ψd(T, T ) = ST,d(0) − ST ,d(T ) d X j=1 P (T, T + j), (3.8)

where ST ,d(0) denotes the d-year contract rate for maturity T , ST ,d(T ) the actual

d-year swap rate for maturity T , and P (T, T + j) is the value of a zero coupon bond that starts at t = T and matures j year(s) later, for some integer j ≥ 1. Using (1.2) as guidance, the d-year contract rate for maturity T is defined as

ST ,d(0) =

P (0, T ) − P (0, T + d)

Pd

j=1P (0, T + j)

(3.9)

and the actual d-year swap rate for maturity T is given by

ST,d(T ) =

1 − P (T, T + d)

Pd

j=1P (T, T + j)

. (3.10)

Here the zero coupon bond price P (t, T ) is calculated by the analytical formula in (2.6) for T > t > 0. Furthermore, note that P (0, T ) is simply given by

P (0, T ) = e−robs(T )·T_, _(3.11)

where robs(T ) is the zero rate observed in the market for maturity T .

Note that ψd(T, T ) denotes the value of a contract at t = T . Hence, in order

to obtain the value today, ψd(0, T ), one has to discount ψd(T, T ) with the

stochas-tic discount factor D(0, T ), as given by (3.4).

Finally, the value of the contract at t = 0 is estimated by calculating the mean of the discounted prices over all scenarios:

ˆ ψd(0, T ) = 1 N N X i=1 ψd,i(T, T ) · Di,T, (3.12)

where ˆψd(0, T ) is the estimated value of the contract today, ψd,i(T, T ) is the value

at time T for scenario i, and Di,T is the discretized stochastic discount factor at

time T for scenario i, as given by (3.6).

A forward swap that is bought at par has a theoretical value at t = 0 that is equal

to zero. Therefore, the martingale test concerns testing whether ˆψd(0, T ) ≈ 0.

(25)

Figure 5 shows simulated swap prices for swaps with maturities of 1, 5, 10, 20, 30, 40, and 50 years. In the remainder of this paper this is the set that is used for the martingale tests. These maturities are selected for two reasons. First of all, this choice enables the opportunity to observe both short and longer term swaps. Secondly, it should give the reader a sufficient feeling on how to adjust and apply this method to any desirable preferences. Consequences of possible other choices are discussed in the conclusion.

Note that the forward starting swaps are only simulated for T = 1, 2, ..., H years. This choice is justified by the fact that the simulation of swaps is highly computa-tional intensive. Furthermore, the set {1, 2, ..., H} gives an adequate representation of the martingale behavior.

0 10 20 30 40 50 60 70 80 90 −0.1 −0.08 −0.06 −0.04 −0.02 0 0.02 0.04 0.06 0.08 0.1

Benchmark: Martingale test forward starting swaps, July 2016 (10000 scenarios)

Maturity (years)

Price forward swap with notional 1

FSS 1 year FSS 5 years FSS 10 years FSS 20 years FSS 30 years FSS 40 years FSS 50 years Figure 5

(26)

4 Correction methods

In this chapter several correction methods are discussed. Inspired by Glasserman and Zhao (2000), Andersen (2007), and Van Haastrecht et.al (2009), short rate r and stochastic variable u are corrected in such a way that the martingale properties for the simulated zero coupon bonds and the forward starting swaps improve.

4.1 Correcting short rate r based on ZCB martingale tests

In this section a correction for the short rate r is derived. In order to achieve an

improvement of the martingale properties, a scalar δZCB is determined for each

time step and applied to every scenario path. This δZCB corrects the short rate r

in such a way that the martingale property for zero coupon bonds (ZCB) perfectly holds. Furthermore, it is chosen to be an additive instead of a multiplicative fac-tor since this does not disturb the variance of the underlying normal distribution. This is crucial with regards to the use of the analytical formula for the price of a zero coupon bond.

In order to find δZCB that corrects the short rate r such that the

correspond-ing zero rate satisfies the martcorrespond-ingale property, the simulated zero coupon bond

price, π0(T ), is compared with the price observed in the market, π0,obs(T ). Note

that both π0(T ) and π0,obs(T ) denote today’s value of a zero coupon bond maturing

at T . Following (3.3) and (3.4), the analytic form of the value of a zero coupon bond is: π0(T ) = E h e−R0Tr(s)ds|F 0 i .

Moreover, according to (3.5) and (3.6), its discrete form is given by:

π0(k) = 1 N N X i=1 e−Pkj=1r(i,j)·∆t_,

where k = T /∆t is the discretization step corresponding to maturity T , r(i, j) the short rate in scenario i at step j, and N the number of scenarios.

Assume that the adjusted price is a martingale. Hence, ˜π0(k − 1) = π0,obs(k − 1).

The corrected price ˜π0(k − 1) can then be defined as

˜ π0(k − 1) = 1 N N X i=1 e−Pk−1j=1˜r(i,j)·∆t_,

where ˜r(i, j) is the corrected short rate for scenario i and step j. Hence, the

simulated price π0(k), uncorrected at step k, is given by:

(27)

Now the correction δZCB(T ), corresponding to the k-th step, is introduced as follows: π0,obs(k) = ˜π0(k) = 1 N N X i=1 e−Pk−1j=1r(i,j)·∆t−(r(i,k)+δ˜ ZCB(T ))·∆t = e−δZCB(T )·∆t 1 N N X i=1 e−Pk−1j=1r(i,j)·∆t−r(i,k)·∆t˜ _.

From here it is straightforward to derive δZCB(T ):

δZCB(T ) = − 1 ∆tln π0,obs(k) + 1 ∆t ln 1 N N X i=1 e−Pk−1j=1r(i,j)·∆t−r(i,k)·∆t˜ ! = − 1 ∆tln e −robs(T )·k·∆t + 1 ∆t ln 1 N N X i=1 e−Pk−1j=1r(i,j)·∆t−r(i,k)·∆t˜ ! = robs(T ) · k + 1 ∆t ln 1 N N X i=1 e−Pk−1j=1˜r(i,j)·∆t−r(i,k)·∆t ! ,

where robs(T ) is the observed zero rate corresponding to the k-th step.

Hence, the corrected short rate for scenario i at step k is given by: ˜ r(i, k) = r(i, k) + δZCB(T ), where δZCB(T ) = robs(T ) · k + 1 ∆t ln 1 N N X i=1 e−Pk−1j=1˜r(i,j)·∆t−r(i,k)·∆t ! ! (4.1) for T = ∆t, 2∆t, ..., H and k = T /∆t.

Consequently, ˜r(k) is an N -dimensional vector containing the corrected short rates

˜

r(i, k). This vector is used to estimate the model implied zero rate ˜rmodel(T ), as

defined by (3.7), which satisfies the martingale property. The corrected ˜rmodel is

(28)

0 10 20 30 40 50 60 70 80 90 −0.4 −0.3 −0.2 −0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 Maturity (years)

δ_ZCB−correction: Zero rates and threshold bounds, July 2016 (10 scenarios)

Interest rate (%)

Observed zero rate Simulated zero rate Thresholds

Figure 6

As the above graph indicates, the martingale test for zero coupon bonds improves greatly. Table 3, that corresponds to Figure 6, illustrates that more clearly. Com-paring these results to those presented in Table 2, shows the great improvement that has been achieved.

Table 3: δZCB-correction.

Performance ZCB test: difference observed and simulated zero rates (in 1/100 %).

Diff. Abs. diff.

Mean 1.34 · 10−18 1.34 · 10−17

Minimum −6.55 · 10−15 ₀

Maximum 1.17 · 10−14 1.17 · 10−14

July 2016, 10 scenarios

(29)

4.2 Correcting stochastic variable u based on one-year F SS

martingale tests

Analogous to the previous section, the idea is to derive a scalar δF SS1 for each time

step, which corrects the stochastic variable u in such a way that the martingale property for forward starting swaps (FSS) with a maturity of one year perfectly holds. Likewise, the martingale property for forward starting swaps with longer maturities are expected to improve as well.

The correction has the following form: ˜

P (T, T + 1) = e−δ(T )P (T, T + 1),

where ˜P (T, T + 1) denotes the martingale value for a zero coupon bond starting

at T and maturing at T + 1, P (T, T + 1) is its simulated value, and δ(T ) is defined by

δ(T ) = C(T, T + 1)δF SS1(T ),

where δF SS1(T ) is the correction for u(T ). Hence, using (2.6), this implies

˜

P (T, T + 1) = A(T, T + 1)e−B(T,T +1)r(T )−C(T,T +1)(u(T )+δF SS1(T ))_.

Following equations (3.8)-(3.12), for one-year forward starting swaps it holds that: ˆ ψ1(0, T ) = 1 N N X i=1 ST ,1(0) − ST,1i (T ) Pi(T, T + 1)Di,T = 1 N N X i=1 ST,1(0) − 1 Pi(T, T + 1) − 1 Pi(T, T + 1)Di,T = 1 N N X i=1 ST,1(0) + 1Pi(T, T + 1) − 1 Di,T,

where S_{T ,1}i (T ) denotes the one-year swap rate at T in scenario i. Moreover, note

that

Pi(T, T + 1) = A(T, T + 1)e−B(T,T +1)r(i,k)−C(T,T +1)u(i,k)

and

Di,T = e−

Pk

j=1r(i,j)·∆t

for k = T /∆t.

(30)

Hence, for ˜Pi(T, T + 1) the ’martingale value’, finding a suitable δF SS1(T ) requires solving: 1 N N X i=1 ST ,1(0) + 1 _˜ Pi(T, T + 1) − 1 Di,T = 1 N N X i=1

ST,1(0) + 1A(T, T + 1)e−B(T,T +1)r(i,k)−C(T,T +1)(u(i,k)+δF SS1(T ))− 1

Di,T = 1 N N X i=1

ST,1(0) + 1e−C(T,T +1)δF SS1(T )A(T, T + 1)e−B(T,T +1)r(i,k)−C(T,T +1)u(i,k)− 1

Di,T = 1 N N X i=1

ST,1(0) + 1e−C(T,T +1)δF SS1(T )Pi(T, T + 1)Di,T − Di,T

= 1 N N X i=1 ST,1(0) + 1e−C(T,T +1)δF SS1(T )Pi(T, T + 1)Di,T − 1 N N X i=1 Di,T = ST ,1(0) + 1e−C(T,T +1)δF SS1(T ) 1 N N X i=1 Pi(T, T + 1)Di,T − 1 N N X i=1 Di,T = 0.

Hence, the corrected u for scenario i at step k is given by: ˜ u(i, k) = u(i, k) + δF SS1(T ), where δF SS1(T ) = − 1 C(T, T + 1)ln    1 N PN i=1Di,T ST ,1(0) + 1 1 N PN i=1Pi(T, T + 1)Di,T    (4.2) for T = ∆t, 2∆t, ..., H.

Note that by substituting r(i, k) for ˜r(i, k) in Di,T and Pi(T, T + 1), one obtains

δF SS1(T ) that enforces the martingale property for one-year forward starting swaps

by correcting stochastic variable u, while allowing the short rate r to be corrected as well. Hence, it is possible to apply both corrections simultaneously such that two martingale properties are satisfied.

By correcting u with δF SS1, the martingale property for forward starting swaps

(31)

Table 4: δF SS1- and δZCB-correction.

Performance F SS1: price contract.

Price ( ˆψ1) Abs. price (| ˆψ1|)

Mean −9.62 · 10−18 _{1.04 · 10}−16

Minimum −2.94 · 10−16 _{6.94 · 10}−19

Maximum 3.09 · 10−16 3.09 · 10−16

July 2016, 10 scenarios

Moreover, the swaps with a longer maturity seem to improve as well. Figure 7 shows the resulting martingale tests.

0 10 20 30 40 50 60 70 80 90 −0.1 −0.08 −0.06 −0.04 −0.02 0 0.02 0.04 0.06 0.08 0.1 δ_ZCB− and δ_FSS 1

−correction: Martingale test forward starting swaps, July 2016 (10000 scenarios)

Maturity (years)

(32)

4.3 Correcting stochastic variable u based on d-year F SZCS

martingale tests

A similar derivation as in the previous section was attempted for longer maturities.

However, as is shown in Appendix C, an analytical derivation of δF SSd that

per-fectly corrects forward starting swaps with a maturity longer than one year turned

out unfruitful. Therefore, in this section a correction δF SZCSd is derived, which is

based on a martingale test for forward starting zero coupon swaps (FSZCS) with a maturity of d years.

The goal of the derivation of this correction is twofold. First of all, this should give a better understanding on how stochastic variable u affects interest rates with longer maturities. For this particular reason it is irrelevant whether it is market convention to use zero coupon swaps or not. Secondly, it is expected that the δF SZCSd-correction improves the martingale tests for forward starting swaps as well, because of the similar behavior of both swaps.

As discussed in section 1.1.2, zero coupon swaps only pay one lump sum at the swaps maturity. Using equation (1.4), the price of a forward starting zero coupon swap with a d-year maturity is thus given by:

ˆ ψZC,d(0, T ) = 1 N N X i=1 ST ,d(0) − ST,di (T ) Pi(T, T + d)Di,T = 1 N N X i=1 ST,d(0) − Pi(T, T ) − Pi(T, T + d) Pi(T, T + d) ! Pi(T, T + d)Di,T = 1 N N X i=1 ST,d(0) − 1 Pi(T, T + d) − 1 ! Pi(T, T + d)Di,T = 1 N N X i=1 ST,d(0) + 1Pi(T, T + d) − 1 Di,T,

where S_{T ,d}i (T ) denotes the d-year swap rate at T in scenario i and

ST ,d(0) =

P (0, T ) − P (0, T + d) P (0, T + d) with P (0, T ) as in (3.11).

Note the similarity with the derivation of δF SS1 in section 4.2. Following the

same procedure, it is thus straightforward to derive δF SZCSd. Hence, the corrected

u for scenario i at step k is given by: ˜

(33)

where δF SZCSd(T ) = − 1 C(T, T + d)ln    1 N PN i=1Di,T ST ,d(0) + 1 1 N PN i=1Pi(T, T + d)Di,T    (4.3) for T = ∆t, 2∆t, ..., H. Obviously, δF SZCSd(T ) = δF SS1(T ) if d = 1.

Imagine the case in which u is corrected with δF SZCS10. Then it is to be

ex-pected that the martingale property for forward starting swaps with a maturity of

ten years improves with respect to the case where u is corrected with δF SS1. Figure

8 indeed suggests an improvement for forward starting swaps with a maturity of ten years. However, note that the martingale property of the ten-year swap does

not perfectly hold, since δF SZCS10 is based on zero coupon swaps, while ’normal’

swaps are tested. Furthermore, note that the one-year swap does no longer per-fectly possess the martingale property. Nevertheless, Figure 8 suggests an overall

improvement with respect to the δF SS1-correction.

0 10 20 30 40 50 60 70 80 90 −0.1 −0.08 −0.06 −0.04 −0.02 0 0.02 0.04 0.06 0.08 0.1 δ_ZCB− and δ_FSZCS 10

Maturity (years)

An illustration of the effect of the δF SZCS50-correction is given in Figure 9. It

(34)

0 10 20 30 40 50 60 70 80 90 −0.1 −0.08 −0.06 −0.04 −0.02 0 0.02 0.04 0.06 0.08 0.1 δ_ZCB− and δ_FSZCS 50

Maturity (years)

5 Historical data

After setting up the framework, Ortec Finance provided the calibrated input pa-rameters for the two-factor Hull-White model for the period January 1999-July

2016. Hence, 211 sets of monthly calibrated parameters a, b, σ1, σ2, and ρ are

available. Typical values for a, b, σ1, σ2, and ρ are 0.6, 0.05, 0.006, 0.006, and -0.9,

respectively. Ortec Finance uses Monte Carlo calibration instead of an analytical calibration method, because this offers more flexibility when applying to other models.

In order to obtain the corresponding vector θ, zero rates observed in the market are used as underlying term structure. For each month in the period January 1999-July 2016 yearly data on zero rates are used, all with a time horizon of 89 years. These bootstrapped rates are obtained from Bloomberg.

5.1 Obtaining time-dependent Hull-White paramater θ

(35)

t = 1, 2, ..., 89 years, one needs to obtain the rates between those time steps. Spline interpolation is a useful tool for obtaining intermediate values. Ametrano and Bianchetti (2009) argue in favor of spline interpolation on log-discounts. Therefore, there will be proceeded as follows.

1. Set the observed zero rate for t = 0 equal to NaN, say z(0) = NaN, where z denotes the zero rate.

2. Calculate the log-discounts for the zero rates:

y(t) = ln e−z(t)t for t = 1, 2, ..., 89 years.

3. Set y(0) = 0.

4. Apply spline interpolation to y for time steps ∆t, 2∆t, ..., 89 years. 5. Convert y back to obtain the ’observed’ zero rates:

robs(t) =

−y(t)

t for t = ∆t, 2∆t, ..., 89 years.

Using these zero rates as the underlying term structure, θ(t) is obtained according equation (2.3). An implementation for the Nelson Siegel term structure can be found in Appendix D.

Note that the starting value r(0) in the process dr of the Hull-White model is chosen as r(0) = z(0.01).

5.2 Extrapolating the zero rate curve

At this point some attention has to be devoted to the extrapolation of the zero rate curve since this is required to calculate the prices of the forward starting swaps. In section 3.3.2 the pricing of forward starting swaps was discussed. Inspecting the d-year contract rate

ST ,d(0) = P (0, T ) − P (0, T + d) Pd j=1P (0, T + j) , where P (0, T ) = e−robs(T )·T

shows that the zero rate curve has to be extrapolated in order to be able to calculate the prices of forward starting swaps for T = 1, 2, ..., 89 years.

(36)

Hence, the last observed forward rate is used as the forward rate for the period extrapolated on. Note that the forward rate is defined as

F (0, t) = −∂ ∂ty(t) = −∂ ∂tln e −z(t)·t = ∂ ∂tz(t) · t

and the discount factor is

P (0, t) = ey(t)

= eln(e−z(t)·t)

= e−z(t)·t.

Hence, by constantly extrapolating the forward rate F (0, t), one can obtain z(t),

P (0, t), and consequently robs(t) for t = 89 + ∆t, 89 + 2∆t, ..., Hu, where Hu is

the preferred ultimate time horizon in years. An extrapolated zero rate curve for

Hu = 200 years is presented in Figure 10.

0 20 40 60 80 100 120 140 160 180 200 −0.4 −0.3 −0.2 −0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8

Extrapolation zero rate curve, July 2016

Maturity (years)

Interest rate (%)

Extrapolated zero rate curve Zero rate curve

(37)

6 Validation

In this section several tests are presented that were performed in order to serve as a justification of the simulation framework. First of all, by eliminating stochastic el-ements the performance of a deterministic model is evaluated. Secondly, the effect of discretization is studied. Finally, the influence of the differentiation methods is discussed.

6.1 Deterministic model

A useful tool in justifying the correctness of the model is eliminating stochasticity. In the risk-neutral two-factor Hull-White model this should imply that the mar-tingale properties perfectly hold.

In order to perform a proper test, all historical moments in the period January 1999-July 2016 are evaluated. Hence, the calibrated values for a and b are used and

σ1 = σ2 = 0. Consequently, it is obvious that the calibrated value for ρ becomes

redundant.

The performance of the martingale test for a specific historical moment is mea-sured by the maximum value it deviates from its martingale value over the 89-year simulation period. By inspecting the overall performance over this historical pe-riod, a clear picture of its deterministic behavior will be achieved.

The results of this test for the benchmark case are presented in Table 5. The per-formance of the zero coupon bond martingale test is measured by the maximum absolute difference from the zero curve (given in 1/100%) and the forward starting swap is measured by its maximum absolute price. Note that it was only necessary to generate a single scenario in the simulation since there was no stochasticity involved.

Table 5: Performance ZCB and FSS martingale tests in deterministic model (∆t = 0.01).

avg. maxh max. maxh min. maxh

ZCB 3.34 · 10−5 1.19 · 10−4 4.69 · 10−6 F SS1 3.96 · 10−5 9.09 · 10−5 7.13 · 10−6 F SS5 7.93 · 10−5 1.86 · 10−4 1.31 · 10−5 F SS10 8.21 · 10−5 1.92 · 10−4 1.34 · 10−5 F SS20 8.21 · 10−5 1.92 · 10−4 1.34 · 10−5 F SS30 8.22 · 10−5 1.92 · 10−4 1.34 · 10−5 F SS40 8.22 · 10−5 1.93 · 10−4 1.34 · 10−5 F SS50 8.23 · 10−5 1.93 · 10−4 1.34 · 10−5

In Table 5 maxh refers to a 211 × 1 vector containing the maximum values for

the period January 1999-July 2016 for a specific martingale test. Hence, there are

(38)

tests. The table presents the average, maximum, and minimum of each vector for

all martingale tests. Note that the unusual notation of vector maxh is used here

as well as in the remainder of this paper to emphasize its purpose.

Note that graphically the deterministic model would suggest that the martingale properties are well satisfied. However, Table 5 shows that there is a significant error that cannot be ignored. Two important factors that are potentially accountable for these errors are discussed in the next sections.

6.2 Differentiation effects

The method of differentiation plays a crucial role in the discretization approxi-mation, since the results are typically very sensitive to any differentiation errors. This sensitivity particularly has its roots in the estimation of θ(t), which is highly dependent on the derivatives. Since θ(t) is the key factor in the risk-neutral char-acter of the Hull-White model, the relevance of a correct approximation is obvious. First of all, attention is devoted to the different built-in functions that Matlab offers with respect to derivatives. Four different functions are tested: diff, gradi-ent, mkpp, and fnder. Table 8-10 in Appendix E present some results. It appears that the function diff performs worse than gradient. Moreover, it turns out that the functions fnder and mkpp deliver the same results as well. The function fnder is preferred because of its relative ease of use.

The tables indicate that the functions gradient and fnder perform quite similar. The function fnder performs slightly better with respect to the forward starting swaps, whereas the martingale property of the zero coupon bonds holds better for the gradient function. Since the difference of the latter is relatively bigger, the gradient function was preferred over fnder.

Secondly, a flat term structure was initiated in order to eliminate the effect of differentiation. Therefore, a constant zero rate of 1% was chosen as underlying term structure. Note that again a deterministic model is tested. Moreover, the same calibrated parameters a and b are used and θ(t) makes the model consistent with the underlying (flat) term structure.

(39)

Table 6: Flat term structure: 1% zero rate.

Performance ZCB and FSS martingale tests in deterministic model (∆t = 0.01).

avg. maxh max. maxh min. maxh

ZCB 1.10 · 10−14 1.10 · 10−14 1.10 · 10−14 F SS1 6.22 · 10−15 6.57 · 10−15 5.82 · 10−15 F SS5 1.27 · 10−14 1.51 · 10−14 1.04 · 10−14 F SS10 1.33 · 10−14 1.62 · 10−14 1.06 · 10−14 F SS20 1.34 · 10−14 1.63 · 10−14 1.06 · 10−14 F SS30 1.34 · 10−14 1.64 · 10−14 1.07 · 10−14 F SS40 1.35 · 10−14 1.64 · 10−14 1.08 · 10−14 F SS50 1.34 · 10−14 1.63 · 10−14 1.07 · 10−14

6.3 Discretization effects

As suggested above, it is likely that the discretization approximation is accountable for the errors in the martingale tests. Therefore, this section studies the effect of the size of the discrete time steps.

Since increasing the number of time steps non-linearly increases the amount of simulation time, it is almost impossible to perform a proper test over the whole historical period. As an illustration: simulating a single historical moment with ∆t = 1, ∆t = 0.1, ∆t = 0.01, ∆t = 0.001, and ∆t = 0.0001 costs about 2, 2, 3, 15, and 2600 seconds, respectively. Hence, one can imagine that performing the test using the current implementation for the whole historical period and for even more steps is too time-consuming. Nonetheless, Table 7 gives a clear picture of how the discretization process has its impact. It suggests that by increasing the number of time steps with a factor ten, the martingale property improves with that factor as

well. Note that the table only shows results for July 2016 and that maxh,jul16 is a

scalar analogous to an entry of vector maxh. Tests for other historical moments

gave the same result.

Table 7: Performance ZCB and FSS martingale tests in deterministic model. Effect number of time steps. July 2016.

(40)

Concluding, the choice for ∆t = 0.01 is justified from a time-consuming point of view. Furthermore, the impact of discretization is quite clear and therefore has to be kept in mind in the remainder of this paper. Finally, it has to be noted that an attempt was made to run a simulation for ∆t = 0.00001, but due to a memory issue this could not be completed.

7 Numerical results

This chapter is organized in three parts. The first section shows how the optimal

maturity d for δF SZCSd was found for the set of forward starting swap martingale

tests and how these optimal corrections vary through the period January 1999-July 2016. In the second section back-test results are presented in which the performance of the ’optimally corrected’ martingale tests are compared with the benchmark. The last section shows how the performance of the benchmark and the corrected model depend on the number of scenarios that is used.

Note once again that all results are obtained for the set of forward starting swaps with maturities of 1, 5, 10, 20, 30, 40, and 50 years. Moreover, the calibrated input parameters for the two-factor Hull-White model discussed in Chapter 5 are used to obtain the results.

7.1 Optimal corrections over time

As multiple δF SZCSd-corrections were available, it had to be decided which

correc-tion was the most effective. In determining a correccorrec-tion’s effectiveness two rules were initiated. Firstly, it was chosen that all martingale tests for forward starting swaps are labeled equally important. Secondly, it was decided to reduce the ex-tremes. Hence, one wants to minimize the worst performing martingale test.

By trying different δF SZCSd-corrections and eye-bolling the results, it became

im-mediately clear that the solution was not as simple as just choosing either δF SZCS1

or δF SZCS50. Therefore, a set of δF SZCSd-corrections was applied in order to get

grip on their functioning. Since it turned out that δF SZCSd-corrections for d > 50

did neither improve short-term nor long-term martingale tests with respect to δF SZCS50, it was chosen to study the set {δF SZCS1, δF SZCS2, ..., δF SZCS50}.

In order to meet both rules, a sum of squared errors (SSE)-like approach was per-formed to determine the optimal maturity on which a correction should be based. This approach simply calculates the integral value of each squared swap price. After that, the resulting squared integrals of the different swaps are added. This

procedure is repeated for every δF SZCSd. As a consequence, 50 different SSE values

(41)

Moreover, this method has the advantage of being a threshold-independent ap-proach, while improving the threshold-based test results. This was confirmed by

comparing the optimal maturity of δF SZCSd determined by both the SSE-like and

a ’highest minimum’ approach. The latter chooses the optimal d based on the δF SZCSd that implies the ’best worst performing’ martingale test. This can be

understood as follows. For each δF SZCSd there are seven martingale tests with a

score (proportion of time between thresholds). The minimum of these seven scores indicates the worst performing martingale test for that specific correction. Hence,

the ’highest minimum’ approach choses d corresponding to the δF SZCSd with the

highest minimum.

Figure 11 shows the optimal maturity d on which δF SZCSd is based for all historic

moments in the period January 1999-July 2016. It is interesting to note that there is a clear pattern observable that distinguishes the pre- and post-2008 periods.

0 5 10 15 20 25 30 35 40 45

Jan−99 Jan−00 Jan−01 Jan−02 Jan−03 Jan−04 Jan−05 Jan−06 Jan−07 Jan−08 Jan−09 Jan−10 Jan−11 Jan−12 Jan−13 Jan−14 Jan−15 Jan−16

Optimal maturity δ

FSZCS_d−correction

maturity d

Figure 11

The simulation procedure to determine these optimal d’s takes approximately eight hours for 1,000 scenarios and all the 211 historical dates. Keeping both computa-tional time and the robustness of the results in mind, 1,000 scenarios seem to be a fair choice.

7.2 Martingale tests over time

This section evaluates the martingale behavior over time for the benchmark and

the corrected model. Using (4.1) and (4.3), both δZCB and δF SZCSd (obviously,

based on its optimal maturity d) were implemented in the corrected model. The goal of this section is to study whether there exists a certain pattern in the mar-tingale performance over time and whether the corrected model incorporates any improvements.

(42)

were used. Furthermore, it should be stressed once again that all simulations are generated by the same sequence of random numbers. Moreover, the results in this section are based on the thresholds introduced in section 3.2.1 and 3.2.2. Hence, for the martingale test for zero coupon bonds the bounds are constantly linearly increasing over time. It was chosen to set the bounds at 0.01% at t = 0 and to let them increase up to 0.05% at the time horizon. In the case of the martingale test for forward starting swaps the thresholds were set at a constant level of 0.1% for all forward starting swaps.

Figure 12 presents the martingale tests over time. Graph (a) shows the

per-formance of the martingale tests for zero coupon bonds for both models. The Y-axis represents the proportion of time the simulated zero curve is between the two threshold bounds, whereas the X-axis stands for the historical date. The pro-portion is calculated by simply dividing the number of simulated points between the bounds by the total number of time steps.

Since δZCB causes the corrected model to perfectly enforce the martingale

prop-erty for zero coupon bonds, it is obvious that the simulated zero curve is strictly between the thresholds. This is obviously a great improvement with respect to the benchmark. The graph clearly indicates that the martingale property holds less well in the more volatile period starting in 2008. Note that the martingale test holds perfectly for the corrected model for any number of scenarios. The amount of 10,000 scenarios used here is thus rather high, but nicely envisions the martingale behavior over time.

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Martingale test zero coupon bonds

proportion between tresholds

benchmark δ−corrections (a) 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

Martingale test FSS 1

benchmark

δ−corrections

(43)

0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

benchmark δ−corrections (c) 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

benchmark δ−corrections (d) 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

benchmark δ−corrections (e) 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

benchmark δ−corrections (f) 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

benchmark δ−corrections (g) 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

benchmark

δ−corrections

(h) Figure 12

A martingale correction method for Monte Carlo simulated interest rates