[go: up one dir, main page]
1/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2
3/*
4Copyright (C) 2007 Ferdinando Ametrano
5Copyright (C) 2006 François du Vignaud
6
7This file is part of QuantLib, a free-software/open-source library
8for financial quantitative analysts and developers - http://quantlib.org/
9
10QuantLib is free software: you can redistribute it and/or modify it
11under the terms of the QuantLib license. You should have received a
12copy of the license along with this program; if not, please email
13<quantlib-dev@lists.sf.net>. The license is also available online at
14<http://quantlib.org/license.shtml>.
15
16This program is distributed in the hope that it will be useful, but WITHOUT
17ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18FOR A PARTICULAR PURPOSE. See the license for more details.
19*/
20
21#include "swapforwardmappings.hpp"
22#include "utilities.hpp"
23#include <ql/models/marketmodels/swapforwardmappings.hpp>
24#include <ql/models/marketmodels/correlations/timehomogeneousforwardcorrelation.hpp>
25#include <ql/models/marketmodels/curvestates/lmmcurvestate.hpp>
26#include <ql/models/marketmodels/evolutiondescription.hpp>
27#include <ql/models/marketmodels/evolvers/lognormalfwdratepc.hpp>
28#include <ql/models/marketmodels/models/flatvol.hpp>
29#include <ql/models/marketmodels/correlations/expcorrelations.hpp>
30#include <ql/models/marketmodels/browniangenerators/sobolbrowniangenerator.hpp>
31#include <ql/models/marketmodels/products/multistep/multistepcoterminalswaptions.hpp>
32#include <ql/models/marketmodels/accountingengine.hpp>
33#include <ql/models/marketmodels/models/cotswaptofwdadapter.hpp>
34#include <ql/models/marketmodels/curvestates/coterminalswapcurvestate.hpp>
35#include <ql/time/schedule.hpp>
36#include <ql/time/daycounters/simpledaycounter.hpp>
37#include <ql/math/statistics/sequencestatistics.hpp>
38#include <ql/pricingengines/blackcalculator.hpp>
39
40#include <ql/models/marketmodels/products/multistep/multistepswaption.hpp>
41
42#if defined(BOOST_MSVC)
43#include <float.h>
44//namespace { unsigned int u = _controlfp(_EM_INEXACT, _MCW_EM); }
45#endif
46
47using namespace QuantLib;
48using namespace boost::unit_test_framework;
49
50using std::fabs;
51using std::sqrt;
52
53namespace {
54
55 class MarketModelData{
56 public:
57 MarketModelData();
58 const std::vector<Time>& rateTimes(){return rateTimes_;}
59 const std::vector<Rate>& forwards(){return forwards_;}
60 const std::vector<Volatility>& volatilities(){return volatilities_;}
61 const std::vector<Rate>& displacements(){return displacements_;}
62 const std::vector<DiscountFactor>& discountFactors(){return discountFactors_;}
63 Size nbRates() const { return nbRates_; }
64
65 private:
66 std::vector<Time> rateTimes_, accruals_;
67 std::vector<Rate> forwards_;
68 std::vector<Spread> displacements_;
69 std::vector<Volatility> volatilities_;
70 std::vector<DiscountFactor> discountFactors_;
71 Size nbRates_;
72 };
73
74 MarketModelData::MarketModelData(){
75 // Times
76 Calendar calendar = NullCalendar();
77 Date todaysDate = Settings::instance().evaluationDate();
78 Date endDate = todaysDate + 9*Years; // change back
79 Schedule dates(todaysDate, endDate, Period(Semiannual),
80 calendar, Following, Following, DateGeneration::Backward, false);
81 nbRates_ = dates.size()-2;
82 rateTimes_ = std::vector<Time>(nbRates_+1);
83 //paymentTimes_ = std::vector<Time>(rateTimes_.size()-1);
84 accruals_ = std::vector<Time>(nbRates_);
85 DayCounter dayCounter = SimpleDayCounter();
86 for (Size i=1; i<nbRates_+2; ++i)
87 rateTimes_[i-1] = dayCounter.yearFraction(d1: todaysDate, d2: dates[i]);
88
89 displacements_ = std::vector<Rate>(nbRates_, .0);
90
91 forwards_ = std::vector<Rate>(nbRates_);
92 discountFactors_ = std::vector<Rate>(nbRates_+1);
93 discountFactors_[0] = 1.0; // .95; fdv1-> WHY ???????
94 for (Size i=0; i<nbRates_; ++i){
95 forwards_[i] = 0.03 + 0.0010*i;
96 accruals_[i] = rateTimes_[i+1] - rateTimes_[i];
97 discountFactors_[i+1] = discountFactors_[i]
98 /(1+forwards_[i]*accruals_[i]);
99 }
100 Volatility mktVols[] = {0.15541283,
101 0.18719678,
102 0.20890740,
103 0.22318179,
104 0.23212717,
105 0.23731450,
106 0.23988649,
107 0.24066384,
108 0.24023111,
109 0.23900189,
110 0.23726699,
111 0.23522952,
112 0.23303022,
113 0.23076564,
114 0.22850101,
115 0.22627951,
116 0.22412881,
117 0.22206569,
118 0.22009939
119 /*
120 0.2,
121 0.2,
122 0.2,
123 0.2,
124 0.2,
125 0.2,
126 0.2,
127 0.2,
128 0.2,
129 0.2,
130 0.2,
131 0.2,
132 0.2,
133 0.2,
134 0.2,
135 0.2,
136 0.2,
137 0.2,
138 0.2,
139 0.2
140 */
141
142 };
143 volatilities_ = std::vector<Volatility>(nbRates_);
144 for (Size i = 0; i < volatilities_.size(); ++i)
145 volatilities_[i] = mktVols[i];//.0;
146 }
147
148 ext::shared_ptr<SequenceStatisticsInc>
149 simulate(const std::vector<Real>& todaysDiscounts,
150 const ext::shared_ptr<MarketModelEvolver>& evolver,
151 const MarketModelMultiProduct& product) {
152 Size paths_;
153#ifdef _DEBUG
154 paths_ = 127;// //
155#else
156 paths_ = 32767; //262144-1; // //; // 2^15-1
157#endif
158
159 Size initialNumeraire = evolver->numeraires().front();
160 Real initialNumeraireValue = todaysDiscounts[initialNumeraire];
161
162 AccountingEngine engine(evolver, product, initialNumeraireValue);
163 ext::shared_ptr<SequenceStatisticsInc> stats(new
164 SequenceStatisticsInc(product.numberOfProducts()));
165 engine.multiplePathValues(stats&: *stats, numberOfPaths: paths_);
166 return stats;
167 }
168
169 MultiStepCoterminalSwaptions makeMultiStepCoterminalSwaptions(
170 const std::vector<Time>& rateTimes, Real strike ){
171 std::vector<Time> paymentTimes(rateTimes.begin(), rateTimes.end()-1);
172 std::vector<ext::shared_ptr<StrikedTypePayoff> > payoffs(paymentTimes.size());
173 for (auto& payoff : payoffs) {
174 payoff = ext::shared_ptr<StrikedTypePayoff>(
175 new PlainVanillaPayoff(Option::Call, strike));
176 }
177 return MultiStepCoterminalSwaptions (rateTimes,
178 paymentTimes, payoffs);
179
180 }
181
182}
183
184
185void SwapForwardMappingsTest::testForwardSwapJacobians()
186{
187 {
188 BOOST_TEST_MESSAGE("Testing forward-rate coinitial-swap Jacobian...");
189 MarketModelData marketData;
190 const std::vector<Time>& rateTimes = marketData.rateTimes();
191 const std::vector<Rate>& forwards = marketData.forwards();
192 const Size nbRates = marketData.nbRates();
193 LMMCurveState lmmCurveState(rateTimes);
194 lmmCurveState.setOnForwardRates(fwdRates: forwards);
195
196 Real bumpSize = 1e-8;
197
198 std::vector<Rate> bumpedForwards(forwards);
199
200 Matrix coinitialJacobian(nbRates,nbRates);
201
202 for (Size i=0; i < nbRates; ++i)
203 for (Size j=0; j < nbRates; ++j)
204 {
205 bumpedForwards = forwards;
206 bumpedForwards[j]+= bumpSize;
207 lmmCurveState.setOnForwardRates(fwdRates: bumpedForwards);
208 Real upRate = lmmCurveState.cmSwapRate(i: 0,spanningForwards: i+1);
209 bumpedForwards[j]-= 2.0*bumpSize;
210 lmmCurveState.setOnForwardRates(fwdRates: bumpedForwards);
211 Real downRate = lmmCurveState.cmSwapRate(i: 0,spanningForwards: i+1);
212 Real deriv = (upRate-downRate)/(2.0*bumpSize);
213 coinitialJacobian[i][j] = deriv;
214
215 }
216
217 Matrix modelJacobian(SwapForwardMappings::coinitialSwapForwardJacobian(cs: lmmCurveState));
218
219 Real errorTolerance = 1e-5;
220
221
222 for (Size i=0; i < nbRates; ++i)
223 for (Size j=0; j < nbRates; ++j)
224 if( fabs(x: modelJacobian[i][j]-coinitialJacobian[i][j]) > errorTolerance)
225 {
226 BOOST_TEST_MESSAGE("rate " << i
227 << ", sensitivity " << j
228 << ", formula value " << modelJacobian[i][j]
229 << " bumping value " << coinitialJacobian[i][j]
230 << "\n");
231
232 BOOST_ERROR("test failed");
233 }
234 }
235
236 {
237
238 BOOST_TEST_MESSAGE("Testing forward-rate constant-maturity swap Jacobian...");
239 MarketModelData marketData;
240 const std::vector<Time>& rateTimes = marketData.rateTimes();
241 const std::vector<Rate>& forwards = marketData.forwards();
242 const Size nbRates = marketData.nbRates();
243 LMMCurveState lmmCurveState(rateTimes);
244 lmmCurveState.setOnForwardRates(fwdRates: forwards);
245
246 Real bumpSize = 1e-8;
247
248 for( Size spanningForwards = 1; spanningForwards < nbRates; ++spanningForwards)
249 {
250
251 std::vector<Rate> bumpedForwards(forwards);
252
253 Matrix cmsJacobian(nbRates,nbRates);
254
255 for (Size i=0; i < nbRates; ++i)
256 for (Size j=0; j < nbRates; ++j)
257 {
258 bumpedForwards = forwards;
259 bumpedForwards[j]+= bumpSize;
260 lmmCurveState.setOnForwardRates(fwdRates: bumpedForwards);
261 Real upRate = lmmCurveState.cmSwapRate(i,spanningForwards);
262 bumpedForwards[j]-= 2.0*bumpSize;
263 lmmCurveState.setOnForwardRates(fwdRates: bumpedForwards);
264 Real downRate = lmmCurveState.cmSwapRate(i,spanningForwards);
265 Real deriv = (upRate-downRate)/(2.0*bumpSize);
266 cmsJacobian[i][j] = deriv;
267
268 }
269
270 Matrix modelJacobian(SwapForwardMappings::cmSwapForwardJacobian(cs: lmmCurveState, spanningForwards));
271
272 Real errorTolerance = 1e-5;
273
274
275 for (Size i=0; i < nbRates; ++i)
276 for (Size j=0; j < nbRates; ++j)
277 if( fabs(x: modelJacobian[i][j]-cmsJacobian[i][j]) > errorTolerance)
278 {
279 BOOST_TEST_MESSAGE(
280 "rate " << i
281 << ", sensitivity " << j
282 << ", formula value " << modelJacobian[i][j]
283 << " bumping value " << cmsJacobian[i][j]
284 << "\n");
285
286 BOOST_ERROR("test failed");
287
288 }
289 }
290
291 }
292}
293
294
295void SwapForwardMappingsTest::testForwardCoterminalMappings() {
296
297 BOOST_TEST_MESSAGE("Testing forward-rate coterminal-swap mappings...");
298 MarketModelData marketData;
299 const std::vector<Time>& rateTimes = marketData.rateTimes();
300 const std::vector<Rate>& forwards = marketData.forwards();
301 const Size nbRates = marketData.nbRates();
302 LMMCurveState lmmCurveState(rateTimes);
303 lmmCurveState.setOnForwardRates(fwdRates: forwards);
304
305 const Real longTermCorr=0.5;
306 const Real beta = .2;
307 Real strike = .03;
308 MultiStepCoterminalSwaptions product
309 = makeMultiStepCoterminalSwaptions(rateTimes, strike);
310
311 const EvolutionDescription& evolution = product.evolution();
312 const Size numberOfFactors = nbRates;
313 Spread displacement = marketData.displacements().front();
314 Matrix jacobian =
315 SwapForwardMappings::coterminalSwapZedMatrix(
316 cs: lmmCurveState, displacement);
317
318 Matrix correlations = exponentialCorrelations(rateTimes: evolution.rateTimes(),
319 longTermCorr,
320 beta);
321 ext::shared_ptr<PiecewiseConstantCorrelation> corr(new
322 TimeHomogeneousForwardCorrelation(correlations,
323 rateTimes));
324 ext::shared_ptr<MarketModel> smmMarketModel(new
325 FlatVol(marketData.volatilities(),
326 corr,
327 evolution,
328 numberOfFactors,
329 lmmCurveState.coterminalSwapRates(),
330 marketData.displacements()));
331
332 ext::shared_ptr<MarketModel>
333 lmmMarketModel(new CotSwapToFwdAdapter(smmMarketModel));
334
335 SobolBrownianGeneratorFactory generatorFactory(SobolBrownianGenerator::Diagonal);
336 std::vector<Size> numeraires(nbRates,
337 nbRates);
338 ext::shared_ptr<MarketModelEvolver> evolver(new LogNormalFwdRatePc
339 (lmmMarketModel, generatorFactory, numeraires));
340
341 ext::shared_ptr<SequenceStatisticsInc> stats =
342 simulate(todaysDiscounts: marketData.discountFactors(), evolver, product);
343 std::vector<Real> results = stats->mean();
344 std::vector<Real> errors = stats->errorEstimate();
345
346 const std::vector<DiscountFactor>& todaysDiscounts = marketData.discountFactors();
347 const std::vector<Rate>& todaysCoterminalSwapRates = lmmCurveState.coterminalSwapRates();
348 for (Size i=0; i<nbRates; ++i) {
349 const Matrix& cotSwapsCovariance = smmMarketModel->totalCovariance(endIndex: i);
350 //Matrix cotSwapsCovariance= jacobian * forwardsCovariance * transpose(jacobian);
351 //Time expiry = rateTimes[i];
352 ext::shared_ptr<PlainVanillaPayoff> payoff(
353 new PlainVanillaPayoff(Option::Call, strike+displacement));
354 //const std::vector<Time>& taus = lmmCurveState.rateTaus();
355 Real expectedSwaption = BlackCalculator(payoff,
356 todaysCoterminalSwapRates[i]+displacement,
357 std::sqrt(x: cotSwapsCovariance[i][i]),
358 lmmCurveState.coterminalSwapAnnuity(numeraire: i,i) *
359 todaysDiscounts[i]).value();
360 if (fabs(x: expectedSwaption-results[i]) > 0.0001)
361 BOOST_ERROR(
362 "expected\t" << expectedSwaption <<
363 "\tLMM\t" << results[i]
364 << "\tstdev:\t" << errors[i] <<
365 "\t" <<std::fabs(results[i]- expectedSwaption)/errors[i]);
366 }
367}
368
369void SwapForwardMappingsTest::testSwaptionImpliedVolatility()
370{
371
372 BOOST_TEST_MESSAGE("Testing implied swaption vol in LMM using HW approximation...");
373 MarketModelData marketData;
374 const std::vector<Time>& rateTimes = marketData.rateTimes();
375 const std::vector<Rate>& forwards = marketData.forwards();
376 const Size nbRates = marketData.nbRates();
377 LMMCurveState lmmCurveState(rateTimes);
378 lmmCurveState.setOnForwardRates(fwdRates: forwards);
379
380 const Real longTermCorr=0.5;
381 const Real beta = .2;
382 Real strike = .03;
383
384 for (Size startIndex = 1; startIndex+2 < nbRates; startIndex = startIndex+5)
385 {
386
387 Size endIndex = nbRates-2;
388
389 ext::shared_ptr<StrikedTypePayoff> payoff(new
390 PlainVanillaPayoff(Option::Call, strike));
391 MultiStepSwaption product(rateTimes, startIndex, endIndex,payoff );
392
393 const EvolutionDescription& evolution = product.evolution();
394 const Size numberOfFactors = nbRates;
395 Spread displacement = marketData.displacements().front();
396 Matrix jacobian =
397 SwapForwardMappings::coterminalSwapZedMatrix(
398 cs: lmmCurveState, displacement);
399
400 Matrix correlations = exponentialCorrelations(rateTimes: evolution.rateTimes(),
401 longTermCorr,
402 beta);
403 ext::shared_ptr<PiecewiseConstantCorrelation> corr(new
404 TimeHomogeneousForwardCorrelation(correlations,
405 rateTimes));
406 ext::shared_ptr<MarketModel> lmmMarketModel(new
407 FlatVol(marketData.volatilities(),
408 corr,
409 evolution,
410 numberOfFactors,
411 lmmCurveState.forwardRates(),
412 marketData.displacements()));
413
414
415 SobolBrownianGeneratorFactory generatorFactory(SobolBrownianGenerator::Diagonal);
416 std::vector<Size> numeraires(nbRates,
417 nbRates);
418 ext::shared_ptr<MarketModelEvolver> evolver(new LogNormalFwdRatePc
419 (lmmMarketModel, generatorFactory, numeraires));
420
421 ext::shared_ptr<SequenceStatisticsInc> stats =
422 simulate(todaysDiscounts: marketData.discountFactors(), evolver, product);
423 std::vector<Real> results = stats->mean();
424 std::vector<Real> errors = stats->errorEstimate();
425
426
427 Real estimatedImpliedVol = SwapForwardMappings::swaptionImpliedVolatility(volStructure: *lmmMarketModel,startIndex,endIndex);
428
429 Real swapRate = lmmCurveState.cmSwapRate(i: startIndex,spanningForwards: endIndex-startIndex);
430 Real swapAnnuity = lmmCurveState.cmSwapAnnuity(numeraire: startIndex,i: startIndex,spanningForwards: endIndex-startIndex)*marketData.discountFactors()[startIndex];
431
432 ext::shared_ptr<PlainVanillaPayoff> payoffDis( new PlainVanillaPayoff(Option::Call, strike+displacement));
433
434 Real expectedSwaption = BlackCalculator(payoffDis,
435 swapRate+displacement, estimatedImpliedVol *sqrt(x: rateTimes[startIndex]),
436 swapAnnuity).value();
437
438 Real error = expectedSwaption - results[0];
439 Real errorInSds = error/errors[0];
440 if (fabs(x: errorInSds) > 3.5 )
441 BOOST_ERROR(
442 "expected\t" << expectedSwaption <<
443 "\tLMM\t" << results[0]
444 << "\tstdev:\t" << errors[0] <<
445 "\t" <<errorInSds);
446 }
447
448}
449
450
451
452test_suite* SwapForwardMappingsTest::suite() {
453 auto* suite = BOOST_TEST_SUITE("swap-forward mappings tests");
454
455 suite->add(QUANTLIB_TEST_CASE(
456 &SwapForwardMappingsTest::testSwaptionImpliedVolatility));
457
458 suite->add(QUANTLIB_TEST_CASE(
459 &SwapForwardMappingsTest::testForwardSwapJacobians));
460 // suite->add(QUANTLIB_TEST_CASE(
461 // &SwapForwardMappingsTest::testForwardCoterminalMappings));
462 return suite;
463}
464
465

source code of quantlib/test-suite/swapforwardmappings.cpp