AN EFFICIENT LATTICE ALGORITHM FOR THE LIBOR MARKET MODEL

The LIBOR Market Model has become one of the most popular models for pricing interest rate products. It is commonly believed that Monte-Carlo simulation is the only viable method available for the LIB

AN EFFICIENT LATTICE ALGORITHM FOR THE LIBOR MARKET MODEL

Tim Xiao

Journal of Derivatives, 19 (1) 25-40, Fall 2011

ABSTRACT

The LIBOR Market Model has become one of the most popular models for pricing interest rate products. It is commonly believed that Monte-Carlo simulation is the only viable method available for the LIBOR Market Model. In this article, however, we propose a lattice approach to price interest rate products within the LIBOR Market Model by introducing a shifted forward measure and several novel fast drift approximation methods. This model should achieve the best performance without losing much accuracy. Moreover, the calibration is almost automatic and it is simple and easy to implement. Adding this model to the valuation toolkit is actually quite useful; especially for risk management or in the case there is a need for a quick turnaround.

Key Words: LIBOR Market Model, lattice model, tree model, shifted forward measure, drift approximation, risk management, calibration, callable exotics, callable bond, callable capped floater swap, callable inverse floater swap, callable range accrual swap.

The LIBOR Market Model (LMM) is an interest rate model based on evolving LIBOR market forward rates under a risk-neutral forward probability measure. In contrast to models that evolve the instantaneous short rates (e.g., Hull-White, Black-Karasinski models) or instantaneous forward rates (e.g., Heath-Jarrow-Morton (HJM) model), which are not directly observable in the market, the objects modeled using the LMM are market observable quantities. The explicit modeling of market forward rates allows for a natural formula for interest rate option volatility that is consistent with the market practice of using the formula of Black for caps. It is generally considered to have more desirable theoretical calibration properties than short rate or instantaneous forward rate models.

In general, it is believed that Monte Carlo simulation is the only viable numerical method available for the LMM (see Piterbarg [2003]). The Monte Carlo simulation is computationally expensive, slowly converging, and notoriously difficult to use for calculating sensitivities and hedges. Another notable weakness is its inability to determine how far the solution is from optimality in any given problem.

In this paper, we propose a lattice approach within the LMM. The model has similar accuracy to the current pricing models in the market, but is much faster. Some other merits of the model are that calibration is almost automatic and the approach is less complex and easier to implement than other current approaches.

We introduce a shifted forward measure that uses a variable substitution to shift the center of a forward rate distribution to zero. This ensures that the distribution is symmetric and can be represented by a relatively small number of discrete points. The shift transformation is the key to achieve high accuracy in relatively few discrete finite nodes. In addition, we present several fast and novel drift approximation approaches. Other concepts used in the model are probability distribution structure exploitation, numerical integration and the long jump technique (we only position nodes at times when decisions need to be made).

This model is actually quite useful for risk management because normally full-revaluations of an entire portfolio under hundreds of thousands of different future scenarios are required for a short time window (see FinPricing (2011)). Without an efficient algorithm, one cannot properly capture and manage the risk exposed by the portfolio.

The rest of this paper is organized as follows: The LMM is discussed in Section I. In Section II, the lattice model is elaborated. The calibration is presented in Section III. The numerical implementation is detailed in Section IV, which will enhance the reader’s understanding of the model and its practical implementation. The conclusions are provided in Section V.

1. LIBOR MARKET MODEL

Let (,,,) be a filtered probability space satisfying the usual conditions, where denotes a sample space, denotes a -algebra, denotes a probability measure, and denotes a filtration. Consider an increasing maturity structure from which expiry-maturity pairs of dates (,) for a family of spanning forward rates are taken. For any time , we define a right-continuous mapping function by . The simply compounded forward rate reset at t for forward period (,) is defined by

(1)

where denotes the time t price of a zero-coupon bond maturing at time T and is the accrual factor or day count fraction for period (,).

Inverting this relationship (1), we can express a zero coupon bond price in terms of forward rates as:

(2)

LIBOR Market Model Dynamics

Consider a zero coupon bond numeraire whose maturity coincides with the maturity of the forward rate. The measure associated with is called forward measure. Terminal measure is a forward measure where the maturity of the bond numeraire matches the terminal date .

For brevity, we discuss the one-factor LMM only. The one-factor LMM (Brace et al. [1997]) under forward measure can be expressed as

If , (3a)

If , (3b)

If , (3c)

where is a Brownian motion.

There is no requirement for what kind of instantaneous volatility structure should be chosen during the life of the caplet. All that is required is (see Hull-White [2000]):

(4)

where denotes the market Black caplet volatility and denotes the strike. Given this equation, it is obviously not possible to uniquely pin down the instantaneous volatility function. In fact, this specification allows an infinite number of choices. People often assume that a forward rate has a piecewise constant instantaneous volatility. Here we choose the forward rate has constant instantaneous volatility regardless of t (see Brigo-Mercurio [2006]).

Shifted Forward Measure

The is a Martingale or driftless under its own measure . The solution to equation (3b) can be expressed as

(5)

where is the current (spot) forward rate. Under the volatility assumption described above, equation (5) can be further expressed as

(6)

Alternatively, we can reach the same Martingale conclusion by directly deriving the expectation of the forward rate (6); that is

(7)

where , are both Brownian motions with a normal distribution (0, t) at time t, is the expectation conditional on the , and the variable substitution used for derivation is

(8)

This variable substitution that ensures that the distribution is centered on zero and symmetry is the key to achieve high accuracy when we express the LMM in discrete finite form and use numerical integration to calculate the expectation. As a matter of fact, without this linear transformation, a lattice method in the LMM either does not exist or introduces too much error for longer maturities.

After applying this variable substitution (8), equation (6) can be expressed as

(9)

Since the LMM models the complete forward curve directly, it is essential to bring everything under a common measure. The terminal measure is a good choice for this purpose, although this is by no means the only choice. The forward rate dynamic under terminal measure is given by

(10)

The solution to equation (10) can be expressed as

(11a)

where the drift is given by

(11b)

where is the drift term.

Applying (8) to (11a), we have the forward rate dynamic under the shifted terminal measure as

(12)

Drift Approximation

Under terminal measure, the drifts of forward rate dynamics are state-dependent, which gives rise to sufficiently complicated non-lognormal distributions. This means that an explicit analytic solution to the forward rate stochastic differential equations cannot be obtained. Therefore, most work on the topic has focused on ways to approximate the drift, which is the fundamental trickiness in implementing the Market Model.

Our model works backwards recursively from forward rate N down to forward rate k. The N-th forward rate without drift can be determined exactly. By the time it takes to calculate the k-th forward rate , all forward rates from to at time t are already known. Therefore, the drift calculation (11b) is to estimate the integrals containing forward rate dynamics , for j=k+1,…,N, with known beginning and end points given by and . For completeness, we list all possible solutions below.

Frozen Drift (FD). Replace the random forward rates in the drift by their deterministic initial values, i.e.,

(13)

Arithmetic Average of the Forward Rates (AAFR). Apply the midpoint rule (rectangle rule) to the random forward rates in the drift, i.e.,

(14)

Arithmetic Average of the Drift Terms (AADT). Apply the midpoint rule to the random drift terms, i.e.,

(15)

Geometric Average of the Forward Rates (GAFR). Replace the random forward rates in the drift by their geometric averages, i.e.,

(16)

Geometric Average of the Drift Terms (GADT). Replace the random drift terms by their geometric averages, i.e.,

(17)

Conditional Expectation of the Forward Rate (CEFR). In addition to the two endpoints, we can further enhance our estimate based on the dynamics of the forward rates. The forward rate follows the dynamic (9) (The drift term is ignored). We can derive the expectation of the forward rate conditional on the two endpoints and replace the random forward rate in the drift by the conditional expectation of the forward rate.

Proposition 1. Assume the forward rate follows the dynamic (9), with the two known endpoints given by and . Based on the conditional expectation of the forward rate , the drift of can be expressed as

(18a)

where the conditional expectation of the forward rate is given by

(18b)

Proof. See Appendix A.

Conditional Expectation of the Drift Term (CEDT). Similarly, we can calculate the conditional expectation of the drift term and replace the random drift term by the conditional expectation.

Proposition 2. Assume the forward rate follows the dynamic (9), with the two known endpoints given by and . Based on the conditional expectation of the drift term , the drift of can be expressed as

(19a)

where the conditional expectation of the drift term is given by

(19b)

(19c)

(19d)

Proof. See Appendix A.

The accuracy and performance of these drift approximation methods are discussed in section IV.

1. THE LATTICE PROCEDURE IN THE LMM

The “lattice” is the generic term for any graph we build for the pricing of financial products. Each lattice is a layered graph that attempts to transform a continuous-time and continuous-space underlying process into a discrete-time and discrete-space process, where the nodes at each level represent the possible values of the underlying process in that period.

There are two primary types of lattices for pricing financial products: tree lattices and grid lattices (or rectangular lattices or Markov chain lattices). The tree lattices, e.g., traditional binomial tree, assume that the underlying process has two possible outcomes at each stage. In contrast with the binomial tree lattice, the grid lattices (see Amin [1993], Gandhi-Hunt [1997], Martzoukos-Trigeorgis [2002], Hagan [2005], and Das [2011]) shown in Exhibit 1, which permit the underlying process to change by multiple states, are built in a rectangular finite difference grid (not to be confused with finite difference numerical methods for solving partial differential equations). The grid lattices are more realistic and convenient for the implementation of a Markov chain solution.

This article presents a grid lattice model for the LMM. To illustrate the lattice algorithm, we use a callable exotic as an example. Callable exotics are a class of interest rate derivatives that have Bermudan style provisions that allow for early exercise into various underlying interest rate products. In general, a callable exotic can be decomposed into an underlying instrument and an embedded Bermudan option.

We will simplify some of the definitions of the universe of instruments we will be dealing with for brevity. Assume the payoff of a generic underlying instrument is a stream of payments for i=1,…,N, where is the structured coupon. The callable exotic is a Bermudan style option to enter the underlying instrument on any of a sequence of notification dates . For any notification date , we define a right-continuous mapping function by . If the option is exercised at t, the reduced price of the underlying instrument, from the structured coupon payer’s perspective, is given by

(20)

where the ratio is usually called the reduced value of the underlying instrument or the reduced exercise value or the reduced intrinsic value.

Lattice approaches are ideal for pricing early exercise products, given their “backward-in-time” nature. Bermudan pricing is usually done by building a lattice to carry out a dynamic programming calculation via backward induction and is standard. The lattice model described below also uses backward induction but exploits the Gaussian structure to gain extra efficiencies.

First we need to create the lattice. The random process we are going to model in the lattice is the LMM (12). Unlike traditional trees, we only position nodes at the determination dates (the payment and exercise dates). At each determination date, the continuous-time stochastic equation (12) shall be discretized into a discrete-time scheme. Such discretized schemes basically convert the Brownian motion into discrete variables. There is no restriction on discretization schemes. At any determination date t, for instance, we discretize the Brownian motion to be equally spaced as a grid of nodes , for i = 1,…,. The number of nodes and the space between nodes at each determination date can vary depending on the length of time and the accuracy requirement. The nodes should cover a certain number of standard deviations of the Gaussian distribution to guarantee a certain level of accuracy. We have the discrete form of the forward rate as

(21)

The zero-coupon bond (2) can be expressed in discrete form as

(22)

We now have expressions for the forward rate (21) and discount bond (22), conditional on being in the state at time t, and from these we can perform valuation for the underlying instrument.

At the maturity date, the value of the underlying instrument is equal to the payoff, i.e.,

(23)

The underlying state process in the LMM (11) is a Brownian motion. The transition probability density from state (, ) to state (, ) is given by

(24)

Applying the variable substitution (8), equation (24) can be expressed as

(25)

Equation (20) can be further expressed as a conditional value on any state (, ) as:

(26)

This is a convolution integral. Some fast integration algorithms, e.g., Cubic Spline Integration, Fast Fourier Transform (FFT), etc., can be used for optimization. We use the Trapezoidal Rule Integration in this paper for ease of illustration.

Incomplete information handling. Convolution is widely used in Electrical Engineering, particularly in signal processing. The important part is that the far left and far right parts of the output are based on incomplete information. Any models that try to compute the transition values using integration will be inaccurate if this problem is not solved, especially for longer maturities and multiple exercise dates. Our solution is to extend the input nodes by padding the far end values on each side and only take the original range of the output nodes.

Next, we determine the option values in each final notification node. On the last exercise date, if we have not already exercised, the reduced option value in any state is given by

(27)

Then, we conduct the backward induction process that is performed by iteratively rolling back a series of long jumps from the final exercise date across notification dates and exercise opportunities until we reach the valuation date. Assume that in the previous rollback step , we calculated the reduced option value: . Now, we go to . The reduced option value at is

(28a)

where the reduced continuation value is given by

(28b)

We repeat the rollback procedure and eventually work our way through the first exercise date. Then the present value of the Bermudan option is found by a final integration given by

(29)

The present value or the price of the callable exotic from the coupon payer’s perspective is:

(30)

This framework can be used to price any interest rate products in the LMM setting and can be easily extended to the Swap Market Model (SMM).

1. Calibration

First, if we choose the LMM as the central model, we need to price interest rate derivatives that depend on either or both of cap and swaption markets. Second, we will undoubtedly use various swaptions to hedge a callable exotic. It is a reasonable expectation that the calibrated model we intend to use to price our exotic, will at least correctly price the market instruments that we intend to hedge with. Therefore, in an exotic derivative pricing situation, recovery of both cap and swaption markets might be desired.

The calibration of the LMM to caplet prices is quite straightforward. However, it is very difficult, if not impossible, to perfectly recover both cap and swaption markets. Fortunately for the LMM, there also exist extremely accurate approximate formulas for swaptions implied volatility, e.g., Rebonato's formula.

We introduced a parameter and set where denotes the market Black caplet volatility. One can choose different for different . For simplicity we describe one situation here. By choosing , we have perfectly calibrated the LMM to the caplet prices in the market. However, our goal is to select a to minimize the sum of the squared differences of the volatilities derived from the market and the volatilities implied by our model for both caps and swaptions combined.

In the optimization, we use Rebonato’s formula for an efficient expression of the model swaption volatilities, given by

(31a)

where =1 under one-factor LMM. The swap rate is given by

(31b)

(31c)

Assume the calibration containing caplets and G swaptions. The error minimization is given by

(32)

where denotes the market Black swaption volatility. The optimization can be found at a stationary point where the first derivative is zero; that is,

(33)

In terms of forward volatilities, we use the time-homogeneity assumption of the volatility structure, where a forward volatility for an option is the same or close to the spot volatility of the option with the same time to expiry. The time-homogeneous volatility structure can avoid non-stationary behavior.

In the LMM, forward swap rates are generally not lognormal. Such deviation from the lognormal paradigm however turns out to be extremely small. Rebonato [1999] shows that the pricing errors of swaptions caused by the lognormal approximation are well within the market bid/ask spread. For most short maturity interest rate products, we can use the lattice model without calibration (33). However, for longer maturity or deeply in the money (ITM) or out of the money (OTM) exotics we may need to use the calibration and even some specific skew/smile adjustment techniques to achieve high accuracy.

1. NUMERICAL IMPLEMENTATION

In this section, we will elaborate on more details of the implementation. We will start with a simple callable bond for the purpose of an easy illustration and then move on to some typical callable exotics, e.g., callable capped floater swap and callable range accrual swap. The reader should be able to implement and replicate the model after reading this section.

Callable Bond

A callable bond is a bond with an option that allows the issuer to retain the privilege of redeeming the bond at some points before the bond reaches the maturity date. For ease of illustration, we choose a very simple callable bond with a one-year maturity, a quarterly payment frequency, a $100 principal amount (A), and a 4% annual coupon rate (the quarterly coupon ). The call dates are 6 months, 9 months, and 12 months. The call price (H) is 100% of the principal. The bond spread () is 0.002. Let the valuation date be 0. A detailed description of the callable bond and current (spot) market data is shown in Exhibit 2.

For a short-term maturity callable bond, our lattice model can reach high accuracy even without calibration (33) and incomplete information handling. Therefore, we set and . The valuation procedure for a callable bond consists of 4 steps:

Step 1: Create the lattice. Based on the long jump technique, we position nodes only at the determination (payment/exercise) dates. The number of nodes and the space between nodes at each determination date may vary depending on the length of time and the accuracy requirement. To simplify the illustration, we choose the same settings across the lattice, with a grid space (space between nodes) , and a number of nodes S=7. It covers standard deviations for a standard normal distribution. The nodes are equally spaced and symmetric, as shown in Exhibit 3.

Step 2: Find the option value at each final node. At the final maturity date , the payoff of the callable bond in any state is given by

(34)

where A denotes the principal amount, C denotes the bond coupon, and H denotes the call price. The option values at the maturity are equal to the payoffs as shown in Exhibit 3.

Step 3: Find the option value at earlier nodes. Let us go to the penultimate notification date . The option value in any state is given by

(35)

Equation (35) can be further expressed in the form of reduced value as

(36a)

where denotes the reduced continuation value in state at given by

(36b)

where denotes the bond spread. Similarly we can compute the reduced callable bond values at . All intermediate reduced values are shown in Exhibit 3.

Step 4: Compute the final integration. The final integral at valuation date 0 is calculated as

(37)

Moreover, we need to add the present value of the coupon at into the final price. The final callable bond value is given by

(38)

The pseudo-code is supplied in Appendix B for the implementation program. The convergence results shown in Exhibit 4 indicate what occurs for a given grid space when we increase the number of nodes S. The speed of convergence is very fast, ensuring that a small number of grids are sufficient. All calculations are converged to 100.7518. One sanity check is that the callable bond price should be close to the straight bond price if the call prices become very high. Both of them are computed as 103.3536.

Callable capped floater swap

A callable capped floater swap has two legs: a regular floating leg and a structured coupon leg. The structured coupon rate of the j-th period () is given by

(39)

where is the notional amount, is the rate cap, is the rate floor, is the spread and is the scale factor. For > 0, it is called a callable capped floater swap. For < 0, it is called a callable inverse floater swap.

We choose a real middle life trade with more than 10 years remaining in its lifetime. The floating leg has a quarterly payment frequency with step-down notionals and step-up spreads. The structured coupon leg has a semi-annually payment frequency with varying notionals, spreads, scales, rate caps, and rate floors. The call schedule is semi-annual.

Callable range accrual swap

A callable range accrual swap has two legs: a regular floating leg and a structured coupon leg. The structured coupon rate of the j-th period () is given by

(40a)

where

(40b)

where R is the fixed rate, and are the accrual range of the j-th period, is the LIBOR rate, is the range accrual index term, is the total number of the business days in the j-th period.

We choose a real 10 years maturity trade. The floating leg has a quarterly payment frequency and the structured coupon leg has a semi-annually payment frequency with varying accrual ranges. It starts with the first call opportunity being in 3 years from inception, and then every year until the last possibility being 9 years from inception. The range accrual index term is 6 months.

The lattice implementation procedure for a callable capped floater swap or a callable range accrual swap is quite similar to the one for a callable bond except the valuation for the underlying instrument.

The convergence diagrams of pricing calculations are shown in Exhibits 5 and 6. Each curve in the diagrams represents the convergence behavior for a given grid space as nodes are increased. All of the lattice results are well converged. If the grid space is smaller, the algorithm has better convergence accuracy but a slower convergence rate, and vice verse.

We benchmarked our model under different drift approximation methods with several standard market approaches, e.g., the regression-based Monte Carlo in the full LMM and the HJM trinomial tree. The model comparisons for the accuracy and speed are shown in Exhibits 7 and 8. With regards to accuracy, as expected, the FD performs very badly. AAFR and GAFR do a little better but errors go in different directions. The same conclusions can be drawn for AADT and GADT. Both CEFR and CEDT are the best. In terms of CPU times, FD, AAFR, AADT, GAFR and GADT are the same. But CEFR and CEDT are slower, especially in the callable range accrual swap case.

1. CONCLUSION

In this paper, we proposed a lattice model in the LMM to price interest rate products. Conclusions can be drawn, supported by the previous sections. First, the model is quite stable. The fast convergence behavior requires fewer discretization nodes. Second, this model has almost equivalent accuracy to the current pricing models in the market. Third, the implementation of the model is relatively easy. The calibration is very simple and straightforward. Finally, the performance of the model is probably the best among all known approaches at the time of writing.

We use the following techniques in our model: shifted forward measure, drift approximation, probability distribution structure exploitation, long jump, numerical integration, incomplete information handling, and calibration. Combining these techniques, the model achieves sufficient accuracy in relatively few time steps and discrete nodes, which makes it a very efficient method.

For ease of illustration, we present the lattice model based on the Trapezoidal Rule integration. A better but slightly more complicated solution is to spline the payoff functions. The cubic spline of the option payoffs can achieve higher accuracy, especially for Greeks calculations, and higher speed. Although cubic spline takes some time, the lattice will require much fewer nodes (23 ~ 28 nodes are good enough) and can perform a much faster integration. In general, the spline method can provide a speedup factor around 3 ~ 5 times.

We have implemented the lattice model to price a variety of interest rate exotics. The algorithm can always achieve a fast convergence rate. The accuracy, however, is a bit trickier, depending on many factors: drift approximation approaches, numerical integration schemes, volatility selections, and calibration, etc. Some work, such as calibration, is more of an art than a science.

REFERENCE

Amin, K. “Jump diffusion option valuation in discrete time.” Journal of Finance, Vol. 48, No. 5 (1993), pp. 1833-1863.

Brace, A., D. Gatarek, and M. Musiela. “The market model of interest rate dynamics.” Mathematical Finance, Vol. 7, No. 4 (1997), pp. 127-155.

Brigo, D., and F. Mercurio. “Interest Rate Models – Theory and Practice with Smiles, Inflation and Credit.” Second Edition, Springer Finance, 2006.

Das, S. “Random lattices for option pricing problems in finance.” Journal of Investment Management, Vol. 9, No.2 (2011), pp. 134-152.

FinPricing, Capital market solution, https://finpricing.com/faq/fxCurve.html

Gandhi, S. and P. Hunt. “Numerical option pricing using conditioned diffusions,” Mathematics of Derivative Securities, Cambridge University Press, Cambridge, 1997.

Hagan, P. “Accrual swaps and range notes.” Bloomberg Technical Report, 2005.

https://github.com/timxiao1203/lmm/raw/master/Lattice-in-LMM-4.pdf

Hull. J., and A. White. “Forward rate volatilities, swap rate volatilities and the implementation of the Libor Market Model.” Journal of Fixed Income, Vol. 10, No. 2 (2000), 46-62.

Martzoukos, H., and L. Trigeorgis. “Real (investment) options with multiple sources of rare events.” European Journal of Operational Research, 136 (2002), 696-706.

Piterbarg, V. “A Practitioner’s guide to pricing and hedging callable LIBOR exotics in LIBOR Market Models.” SSRN Working paper, 2003.

Rebonato, R. “Calibrating the BGM model.” RISK, March (1999), 74-79.

APPENDIX A:

Proof of Proposition 1. We rewrite (9) as

(A1)

In the general Brownian Bridge case when the Wiener process has =a and =b, the distribution of at time is normal given by

(A2)

In our case: , , a=0, b=, , thus (A2) can be expressed as

(A3)

Let . According to the linear transformation rule, is a normal given by

(A4)

Let . By definition, is a lognormal given by . According to the characterizations of the lognormal distribution, the mean and variance of are

(A5a)

(A5b)

We have the conditional expectation of the forward rate as

(A6)

Proof of Proposition 2. Let where is defined above. According to the linear transformation rule, is a lognormal given by . The mean and variance of are

(A7a)

(A7b)

On the other hand, according to the characterizations of the lognormal distribution, the mean and variance of are

(A8a)

(A8b)

Solving the equation (A8a) and (A8b), we get

(A9a)

(A9b)

We know the first negative moment of the lognormal is and have the conditional expectation of the drift term as

(A10)

where , are given by (A7a) and (A7b).

APPENDIX B:

The following pseudo-code (C++) demonstrates how to implement the model to price a callable bond. For the purpose of an easy illustration, we choose the same settings (the number of nodes and the grid space) across the lattice and use the Trapezoidal Rule for numerical integration.

// 2*numNodes = 2*mNumNodes = the number of nodes (S); gap = mGap = the grid space (Phi)

double priceCallableBond (BondTrade* bd, CallableBond* cb, int numNodes, double gap) {

double pv;

cb->fillLattice();

// The last exercise

CallSchedule& cs = bd->callSch[numCallSch-1];

if (cs.term == bd->cFlow[numCashFlow-1].endDate) // The last exercise is at maturity

for (int i= -numNodes; i <= numNodes; i++)

cs.reducedValue[i+numNodes] = min (cs.callPrice,

bd->cFlow[numCashFlow-1].reducedPayoff[i+numNodes]);

else { // The last exercise is before maturity

for (int i= -numNodes; i <= numNodes; i++) {

pv = 0;

for (int j = bd->numCF-1; (bd->cFlow[j].endDate >= cs.term) && (j >= 0); j--) {

CashFlow& cf = bd->cFlow[j];

(cf.endDate == cs.term) ? pv += cf.reducedPayoff[i+numNodes]

: pv += exp(-bondSpread*(cf.endDate-cs.term)) * cb->integral(i,

cs.vol, cf.vol, cf.endDate, cs.term, cf.reducedPayoff);

}

cs.reducedValue[i+numNodes] = min (cs.callPrice/cs.df[i+numNodes], pv);

}

}

if (numCallSch > 1) { // The remaining exercises

for (int i = numCallSch - 2; i>=0; i--) {

CallSchedule& cs = bd->callSch[i];

CallSchedule& preCs = bd->callSch[i+1];

for (int j = -numNodes; j <= numNodes; j++) {

pv = exp(-bondSpread * (preCs.term - cs.term))

* cb->integral (j, cs.vol, preCs.vol, preCs.term, cs.term, preCs.reducedValue);

for (int k=bd->numCF-1; k >= 0; k--) // Count intermediate coupons

if ((bd->cFlow[k].endDate < preCs.term) && (bd->cFlow[k].endDate >= cs.term))

pv += bd->cFlow[k].reducedPayoff[j+numNodes]

* exp (-bondSpread*(bd->cFlow[k].endDate - cs.term));

cs.reducedValue[j+numNodes] = min (cs.callPrice/cs.df[j+numNodes], pv);

}

}

}

// The final integral

CallSchedule& preCs = bd->callSch[0];

pv = cb->integral (0, 0, preCs.vol, preCs.term, 0, preCs.reducedValue) *exp(-bondSpread*(preCs.term));

pv *= bd->cFlow[bd->numCF-1].endDf; // endDf: discount factor from 0 to the end date

for (int k=bd->numCF-1; k >= 0; k--) // Count intermediate coupons

if ((bd->cFlow[k].endDate < preCs.term))

pv += bd->cFlow[k].coupon * bd->cFlow[k].endDf * exp(-bondSpread * bd->cFlow[k].endDate);

return pv;

}

void CallableBond::fillLattice() {

for (int i = mTrade->numCF-1; i>=0; i--) {

CashFlow& cf = mTrade->cFlow[i];

if (cf.endDate < mTrade->callSch[0].term) break;

for (int j = -mNumNodes; j <= mNumNodes; j++)

fillNode(i, j, cf.startDate, mDrift);

}

}

void CallableBond::fillNode(int cI, int nI, double vT, DriftAppx flag) {

int numCF = mTrade->numCF;

double avgF, expon, fwdt, drift = 0;

CashFlow& fl = mTrade->cFlow[cI];

if (cI == numCF-1) { // At maturity

fl.df[nI + mNumNodes] = 1.0;

fl.reducedPayoff[nI + mNumNodes] = fl.notional + fl.coupon;

}

else if (fl.startDate <= 0) // Starting before valuation date)

fl.reducedPayoff[nI + mNumNodes] = fl.coupon * fl.endDf / mTrade->cFlow[numCF-1].endDf;

else {

fl.df[nI + mNumNodes] = 1.0;

for (int i = numCF - 1; i > cI; i--) {

CashFlow& cf = mTrade->cFlow[i];

expon = (cf.vol * cf.vol * vT / 2) + cf.vol * nI * mGap;

fwdt = cf.fwd0 * exp(-drift + expon);

switch (flag) { // The other cases are similar to either AAFR or CEFR

case AAFR: // Arithemic Average Fwd Rate

avgF = 0.5 * (cf.fwd0 + fwdt);

drift += vT * fl.vol * cf.vol * cf.delta * avgF / (1 + cf.delta * avgF);

break;

case CEFR: // Conditional Expectation of Fwd Rate

drift += fl.vol * cf.vol * integralFwd(cf.fwd0, fwdt, 0, vT, cf.vol, cf.delta);

break;

default:

break;

}

fl.df[nI + mNumNodes] /= (1 + fwdt * cf.delta); // df: discount factor maturing at maturity

}

fl.reducedPayoff[nI + mNumNodes] = fl.coupon / fl.df[nI + mNumNodes];

}

}

// Gauss-Legendre integration for drift

const double xArray[] = {0, 0.1488743389, 0.4333953941, 0.6794095682, 0.8650633666, 0.9739065285};

const double wArray[] = {0, 0.2955242247, 0.2692667193, 0.2190863625, 0.1494513491, 0.0666713443};

double CallableBond::integralFwd(double F0, double Ft, double a, double b, double vol, double delta) {

double xm = 0.5 * (b + a);

double xr = 0.5 * (b - a);

double ss = 0, dx = 0;

for (int j = 1; j <= 5; j++) {

dx = xr * xArray[j];

ss += wArray[j] * (expectFwd(F0, Ft, (xm + dx), b, vol, delta)

+ expectFwd(F0, Ft, (xm - dx), b, vol, delta));

}

return ss * xr;

}

double CallableBond::expectFwd(double F0, double Ft, double s, double t, double vol, double delta) {

double mean = F0 * pow ((Ft / F0), (s / t)) * exp(0.5 * vol * vol * s * (t - s) / t);

return delta * mean / (1 + delta * mean);-

}

// Trapezoidal Rule Integration

double CallableBond::integral (int curPos, double curVol, double preVol, double preTerm,

double curTerm, double* value){

double diffPos, tmpV, sum = 0;

for (int k = -mNumNodes; k <= mNumNodes; k++) {

diffPos = k*mGap - curPos*mGap + preVol * preTerm - curVol * curTerm;

tmpV = value[k+mNumNodes] * exp (-diffPos * diffPos/(2*(preTerm - curTerm)));

((k == -mNumNodes) || (k == mNumNodes)) ? sum += 0.5 * tmpV : sum += tmpV;

}

return sum * mGap / sqrt(2 * PI * (preTerm - curTerm));

}

EXHIBIT 1. The Grid/Rectangular Lattice

This exhibit defines the state space for the underlying process over the first two discrete time periods. The starting state at valuation date 0 is the single root of the lattice. At each date the underlying process is discretized into a number of vertical nodes/states indexed by j. The value denotes the underlying process in state j at date . The node , for instance, can evolve to any discrete state in with certain transition probabilities. For a Brownian motion, the transition probability can be easily determined by (25).

EXHIBIT 2: The Callable Bond and Associated Spot Market Data

The callable bond has a one-year maturity, a $100 principal, a quarterly payment frequency, and a 4% annual coupon rate. Delta = (end date – start date)/365 (day count: ACT/365). The discount bond matures at the end date . The call dates are 6, 9, and 12 months.

Cash flow index	1	2	3	4
Start date (days)	0	92	181	273
End date (days)	92 ()	181 ()	273 ()	365 ()
Delta (years)	0.252055	0.243836	0.252055	0.252055
Payoff ($)	1	1	1	101
Call Schedule (days)	-	181	273	365
Discount bond	0.999313	0.998557	0.997293	0.995667
Black Volatility	-	0.337631	0.344218	0.350878

EXHIBIT 3: The LMM Lattice Structure of the Callable Bond.

The callable bond is defined in Exhibit 2. denotes the reduced value of the callable bond at any node (i, j). denotes the coupon at . is the value calculated by the final integration. is the final callable bond value that is equal to plus the present value of . The grid space is and the number of nodes is . This lattice has 3 steps and 7 nodes.

0

EXHIBIT 4: The Convergence Results for the Callable Bond.

The callable bond is defined in Exhibit 2. and drift approximation is AADT. Each curve represents the convergence behavior for a given grid space (phi) as nodes are added. All calculations are converged to 100.7518.

.

EXHIBIT 5: The Convergence Results for the Callable Capped Floater Swap

The callable capped floater swap has more than 10 years remaining in its lifetime. The floating leg has a quarterly payment frequency. The structural leg has a semi-annually payment frequency. The call schedule is semi-annual.=1 and drift approximation is CEDT. Each curve represents the convergence behavior for a given grid space (phi) as nodes (N) are added.

EXHIBIT 6: The Convergence Results for the Callable Range Accrual Swap

The callable range accrual swap has 10 years maturity. The floating leg has a quarterly payment frequency. The structural leg has a semi-annually payment frequency. There are 7 call opportunities. =1 and drift approximation is CEDT. Each curve represents the convergence behavior for a given grid space (phi) as nodes are added.

EXHIBIT 7: The Benchmark Results for the Callable Capped Floater Swap

This exhibit presents the results for model comparison. We benchmark the lattice model under different drift approximation methods with several standard market approaches, e.g., the regression-based Monte Carlo in the full LMM and the HJM trinomial tree, for both accuracy and speed. The trade is the same as the one in Exhibit 5. The grid space is =1/8 and the number of nodes is S=200. PC denotes Predictor-Corrector. The column ‘Dif from MC’ = 1 – (current row price) / (price of MC in LMM). All computational times are denoted in seconds on a computer with a 2.33 GHz Duo Core CPU.

Model		Drift	Steps n	Calls	Nodes/Paths	Price	Err from MC	Run time
MC in LMM	-	PC	40	20	1 million	4,546,863.3	0	290.32
HJM tri-tree	-	-	1979	20	2n+1	4,602,136.3	1.22%	15.01
Our Model	1	FD	40	20	200	4,822,728.4	6.07%	0.32
1	AAFR	40	20	200	4,637,263.2	1.99%	0.32
1	AADT	40	20	200	4,637,718.1	2.00%	0.32
1	GAFR	40	20	200	4,698,215.6	3.33%	0.32
1	GADT	40	20	200	4,698,441.3	3.33%	0.32
1	CEFR	40	20	200	4,665,210.3	2.60%	0.38
1	CEDT	40	20	200	4,665,552.4	2.61%	0.39
0.99	FD	40	20	200	4,708,768.9	3.56%	0.32
0.99	AAFR	40	20	200	4,504,989.2	-0.92%	0.32
0.99	AADT	40	20	200	4,505,426.3	-0.91%	0.32
0.99	GAFR	40	20	200	4,609,779.5	1.38%	0.32
0.99	GADT	40	20	200	4,609,996.6	1.39%	0.32
0.99	CEFR	40	20	200	4,563,689.2	0.37%	0.38
0.99	CEDT	40	20	200	4,563,730.9	0.37%	0.39

EXHIBIT 8: The Benchmark Results for the Callable Range Accrual Swap

This exhibit presents the results for model comparison. We benchmark the lattice model under different drift approximation methods with several standard market approaches, e.g., the regression-based Monte Carlo in the full LMM and the HJM trinomial tree, for both accuracy and speed. The trade is the same as the one in Exhibit 6. The grid space is =1/8 and the number of nodes is S=200. The column ‘Dif from MC’ = 1 – (current row price) / (price of MC in LMM). All computational times are denoted in seconds on a computer with a 2.33 GHz Duo Core CPU.

Model		Drift	Steps n	Calls	Nodes/Paths	Price	Dif from MC	Run time
MC in LMM	-	Euler	1801	7	1 million	585793.2	0.00%	2372.21
HJM tri-tree	-	-	1801	7	2n+1	582167.8	-0.62%	15.62
Our Model	1	FD	1801	7	200	648365.4	10.68%	0.21
1	AAFR	1801	7	200	602482.2	2.85%	0.21
1	AADT	1801	7	200	602742.1	2.89%	0.21
1	GAFR	1801	7	200	616318.6	5.21%	0.21
1	GADT	1801	7	200	616425.3	5.23%	0.21
1	CEFR	1801	7	200	598253.3	2.13%	2.21
1	CEDT	1801	7	200	598372.4	2.15%	2.35
0.99	FD	1801	7	200	609373.9	4.03%	0.21
0.99	AAFR	1801	7	200	579337.2	-1.10%	0.21
0.99	AADT	1801	7	200	579386.3	-1.09%	0.21
0.99	GAFR	1801	7	200	591981.5	1.06%	0.21
0.99	GADT	1801	7	200	591917.6	1.05%	0.21
0.99	CEFR	1801	7	200	588918.9	0.53%	2.21
0.99	CEDT	1801	7	200	588935.7	0.54%	2.35