[go: up one dir, main page]
1/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2
3/*
4 Copyright (C) 2006 Ferdinando Ametrano
5 Copyright (C) 2006 Marco Bianchetti
6 Copyright (C) 2006 Cristina Duminuco
7 Copyright (C) 2006 StatPro Italia srl
8
9 This file is part of QuantLib, a free-software/open-source library
10 for financial quantitative analysts and developers - http://quantlib.org/
11
12 QuantLib is free software: you can redistribute it and/or modify it
13 under the terms of the QuantLib license. You should have received a
14 copy of the license along with this program; if not, please email
15 <quantlib-dev@lists.sf.net>. The license is also available online at
16 <http://quantlib.org/license.shtml>.
17
18 This program is distributed in the hope that it will be useful, but WITHOUT
19 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20 FOR A PARTICULAR PURPOSE. See the license for more details.
21*/
22
23#include "marketmodel_cms.hpp"
24#include "utilities.hpp"
25#include <ql/models/marketmodels/correlations/timehomogeneousforwardcorrelation.hpp>
26#include <ql/models/marketmodels/curvestates/lmmcurvestate.hpp>
27#include <ql/models/marketmodels/curvestates/cmswapcurvestate.hpp>
28#include <ql/models/marketmodels/evolvers/lognormalcmswapratepc.hpp>
29#include <ql/legacy/libormarketmodels/lmlinexpcorrmodel.hpp>
30#include <ql/legacy/libormarketmodels/lmextlinexpvolmodel.hpp>
31#include <ql/models/marketmodels/models/flatvol.hpp>
32#include <ql/models/marketmodels/models/abcdvol.hpp>
33#include <ql/models/marketmodels/correlations/expcorrelations.hpp>
34#include <ql/models/marketmodels/accountingengine.hpp>
35#include <ql/models/marketmodels/products/multistep/multistepcoterminalswaptions.hpp>
36#include <ql/models/marketmodels/products/multistep/multistepcoterminalswaps.hpp>
37#include <ql/models/marketmodels/products/multiproductcomposite.hpp>
38#include <ql/models/marketmodels/browniangenerators/sobolbrowniangenerator.hpp>
39#include <ql/models/marketmodels/swapforwardmappings.hpp>
40#include <ql/time/schedule.hpp>
41#include <ql/time/calendars/nullcalendar.hpp>
42#include <ql/time/daycounters/simpledaycounter.hpp>
43#include <ql/pricingengines/blackformula.hpp>
44#include <ql/pricingengines/blackcalculator.hpp>
45#include <ql/utilities/dataformatters.hpp>
46#include <ql/math/statistics/sequencestatistics.hpp>
47#include <ql/math/statistics/convergencestatistics.hpp>
48#include <iostream>
49#include <sstream>
50
51using namespace QuantLib;
52using namespace boost::unit_test_framework;
53
54namespace market_model_cms_test {
55
56 Date todaysDate, startDate, endDate;
57 std::vector<Time> rateTimes;
58 std::vector<Real> accruals;
59 Calendar calendar;
60 DayCounter dayCounter;
61 std::vector<Rate> todaysForwards, todaysCMSwapRates;
62 std::vector<Real> cMSwapAnnuity;
63 Spread displacement;
64 std::vector<DiscountFactor> todaysDiscounts;
65 std::vector<Volatility> volatilities, blackVols;
66 Real a, b, c, d;
67 Real longTermCorrelation, beta;
68 Size measureOffset_;
69 unsigned long seed_;
70 Size paths_, trainingPaths_;
71 bool printReport_ = false;
72 Size spanningForwards;
73
74 void setup() {
75
76 // Times
77 calendar = NullCalendar();
78 todaysDate = Settings::instance().evaluationDate();
79 //startDate = todaysDate + 5*Years;
80 endDate = todaysDate + 10*Years;
81 Schedule dates(todaysDate, endDate, Period(Semiannual),
82 calendar, Following, Following,
83 DateGeneration::Backward, false);
84 rateTimes = std::vector<Time>(dates.size()-1);
85
86 accruals = std::vector<Real>(rateTimes.size()-1);
87 dayCounter = SimpleDayCounter();
88 for (Size i=1; i<dates.size(); ++i)
89 rateTimes[i-1] = dayCounter.yearFraction(d1: todaysDate, d2: dates[i]);
90
91 for (Size i=1; i<rateTimes.size(); ++i)
92 accruals[i-1] = rateTimes[i] - rateTimes[i-1];
93
94 // Rates & displacement
95 todaysForwards = std::vector<Rate>(accruals.size());
96 displacement = 0.02;
97 for (Size i=0; i<todaysForwards.size(); ++i)
98 todaysForwards[i] = 0.03 + 0.0010*i;
99 LMMCurveState curveState_lmm(rateTimes);
100 curveState_lmm.setOnForwardRates(fwdRates: todaysForwards);
101 // until ConstantMaturitySwap is ready
102 spanningForwards = todaysForwards.size();
103 todaysCMSwapRates = curveState_lmm.cmSwapRates(spanningForwards);
104
105 // Discounts
106 todaysDiscounts = std::vector<DiscountFactor>(rateTimes.size());
107 todaysDiscounts[0] = 0.95;
108 for (Size i=1; i<rateTimes.size(); ++i)
109 todaysDiscounts[i] = todaysDiscounts[i-1] /
110 (1.0+todaysForwards[i-1]*accruals[i-1]);
111
112 // Swaption Volatilities
113 Volatility mktVols[] = {0.15541283,
114 0.18719678,
115 0.20890740,
116 0.22318179,
117 0.23212717,
118 0.23731450,
119 0.23988649,
120 0.24066384,
121 0.24023111,
122 0.23900189,
123 0.23726699,
124 0.23522952,
125 0.23303022,
126 0.23076564,
127 0.22850101,
128 0.22627951,
129 0.22412881,
130 0.22206569,
131 0.22009939
132 };
133 a = -0.0597;
134 b = 0.1677;
135 c = 0.5403;
136 d = 0.1710;
137 volatilities = std::vector<Volatility>(todaysCMSwapRates.size());
138 blackVols = std::vector<Volatility>(todaysCMSwapRates.size());
139 for (Size i=0; i<todaysCMSwapRates.size(); i++) {
140 volatilities[i] = todaysCMSwapRates[i]*mktVols[i]/
141 (todaysCMSwapRates[i]+displacement);
142 blackVols[i]= mktVols[i];
143 }
144
145 // Cap/Floor Correlation
146 longTermCorrelation = 0.5;
147 beta = 0.2;
148 measureOffset_ = 5;
149
150 // Monte Carlo
151 seed_ = 42;
152
153 #ifdef _DEBUG
154 paths_ = 127;
155 trainingPaths_ = 31;
156 #else
157 paths_ = 32767; //262144-1; //; // 2^15-1
158 trainingPaths_ = 8191; // 2^13-1
159 #endif
160 }
161
162 ext::shared_ptr<SequenceStatisticsInc>
163 simulate(const ext::shared_ptr<MarketModelEvolver>& evolver,
164 const MarketModelMultiProduct& product) {
165 Size initialNumeraire = evolver->numeraires().front();
166 Real initialNumeraireValue = todaysDiscounts[initialNumeraire];
167
168 AccountingEngine engine(evolver, product, initialNumeraireValue);
169 ext::shared_ptr<SequenceStatisticsInc> stats(
170 new SequenceStatisticsInc(product.numberOfProducts()));
171 engine.multiplePathValues(stats&: *stats, numberOfPaths: paths_);
172 return stats;
173 }
174
175
176 enum MarketModelType { ExponentialCorrelationFlatVolatility,
177 ExponentialCorrelationAbcdVolatility/*,
178 CalibratedMM*/
179 };
180
181 std::string marketModelTypeToString(MarketModelType type) {
182 switch (type) {
183 case ExponentialCorrelationFlatVolatility:
184 return "Exp. Corr. Flat Vol.";
185 case ExponentialCorrelationAbcdVolatility:
186 return "Exp. Corr. Abcd Vol.";
187 //case CalibratedMM:
188 // return "CalibratedMarketModel";
189 default:
190 QL_FAIL("unknown MarketModelEvolver type");
191 }
192 }
193
194 ext::shared_ptr<MarketModel> makeMarketModel(
195 const EvolutionDescription& evolution,
196 Size numberOfFactors,
197 MarketModelType marketModelType,
198 Spread rateBump = 0.0,
199 Volatility volBump = 0.0) {
200
201 std::vector<Time> fixingTimes(evolution.rateTimes());
202 fixingTimes.pop_back();
203 ext::shared_ptr<LmVolatilityModel> volModel(new
204 LmExtLinearExponentialVolModel(fixingTimes, 0.5, 0.6, 0.1, 0.1));
205 ext::shared_ptr<LmCorrelationModel> corrModel(new
206 LmLinearExponentialCorrelationModel(evolution.numberOfRates(),
207 longTermCorrelation, beta));
208 std::vector<Rate> bumpedRates(todaysCMSwapRates.size());
209 LMMCurveState curveState_lmm(rateTimes);
210 curveState_lmm.setOnForwardRates(fwdRates: todaysForwards);
211 std::vector<Rate> usedRates =
212 curveState_lmm.cmSwapRates(spanningForwards);
213 std::transform(first: usedRates.begin(), last: usedRates.end(),
214 result: bumpedRates.begin(),
215 unary_op: [=](Rate r){ return r + rateBump; });
216
217 std::vector<Volatility> bumpedVols(volatilities.size());
218 std::transform(first: volatilities.begin(), last: volatilities.end(),
219 result: bumpedVols.begin(),
220 unary_op: [=](Volatility v){ return v + volBump; });
221 Matrix correlations = exponentialCorrelations(rateTimes: evolution.rateTimes(),
222 longTermCorr: longTermCorrelation,
223 beta);
224 ext::shared_ptr<PiecewiseConstantCorrelation> corr(new
225 TimeHomogeneousForwardCorrelation(correlations,
226 evolution.rateTimes()));
227 switch (marketModelType) {
228 case ExponentialCorrelationFlatVolatility:
229 return ext::shared_ptr<MarketModel>(new
230 FlatVol(bumpedVols,
231 corr,
232 evolution,
233 numberOfFactors,
234 bumpedRates,
235 std::vector<Spread>(bumpedRates.size(),
236 displacement)));
237 case ExponentialCorrelationAbcdVolatility:
238 return ext::shared_ptr<MarketModel>(new
239 AbcdVol(0.0,0.0,1.0,1.0,
240 bumpedVols,
241 corr,
242 evolution,
243 numberOfFactors,
244 bumpedRates,
245 std::vector<Spread>(bumpedRates.size(),
246 displacement)));
247 //case CalibratedMM:
248 // return ext::shared_ptr<MarketModel>(new
249 // CalibratedMarketModel(volModel, corrModel,
250 // evolution,
251 // numberOfFactors,
252 // bumpedForwards,
253 // displacement));
254 default:
255 QL_FAIL("unknown MarketModel type");
256 }
257 }
258
259 enum MeasureType { ProductSuggested, Terminal,
260 MoneyMarket, MoneyMarketPlus };
261
262 std::string measureTypeToString(MeasureType type) {
263 switch (type) {
264 case ProductSuggested:
265 return "ProductSuggested measure";
266 case Terminal:
267 return "Terminal measure";
268 case MoneyMarket:
269 return "Money Market measure";
270 case MoneyMarketPlus:
271 return "Money Market Plus measure";
272 default:
273 QL_FAIL("unknown measure type");
274 }
275 }
276
277 std::vector<Size> makeMeasure(const MarketModelMultiProduct& product,
278 MeasureType measureType) {
279 std::vector<Size> result;
280 const EvolutionDescription& evolution(product.evolution());
281 switch (measureType) {
282 case ProductSuggested:
283 result = product.suggestedNumeraires();
284 break;
285 case Terminal:
286 result = terminalMeasure(evolution);
287 if (!isInTerminalMeasure(evolution, numeraires: result)) {
288 BOOST_ERROR("failure in verifying Terminal measure:\n"
289 << to_stream(result));
290 }
291 break;
292 case MoneyMarket:
293 result = moneyMarketMeasure(evolution);
294 if (!isInMoneyMarketMeasure(evolution, numeraires: result)) {
295 BOOST_ERROR("failure in verifying MoneyMarket measure:\n"
296 << to_stream(result));
297 }
298 break;
299 case MoneyMarketPlus:
300 result = moneyMarketPlusMeasure(evolution, offset: measureOffset_);
301 if (!isInMoneyMarketPlusMeasure(evolution, numeraires: result,
302 offset: measureOffset_)) {
303 BOOST_ERROR("failure in verifying MoneyMarketPlus(" <<
304 measureOffset_ << ") measure:\n" <<
305 to_stream(result));
306 }
307 break;
308 default:
309 QL_FAIL("unknown measure type");
310 }
311 checkCompatibility(evolution, numeraires: result);
312 if (printReport_) {
313 BOOST_TEST_MESSAGE(" " << measureTypeToString(measureType) << ": "
314 << to_stream(result));
315 }
316 return result;
317 }
318
319 enum EvolverType { Ipc, Pc, NormalPc};
320
321 std::string evolverTypeToString(EvolverType type) {
322 switch (type) {
323 case Ipc:
324 return "iterative predictor corrector";
325 case Pc:
326 return "predictor corrector";
327 case NormalPc:
328 return "predictor corrector for normal case";
329 default:
330 QL_FAIL("unknown MarketModelEvolver type");
331 }
332 }
333
334 ext::shared_ptr<MarketModelEvolver> makeMarketModelEvolver(
335 const ext::shared_ptr<MarketModel>& marketModel,
336 const std::vector<Size>& numeraires,
337 const BrownianGeneratorFactory& generatorFactory,
338 EvolverType evolverType,
339 Size initialStep = 0) {
340 switch (evolverType) {
341 case Pc:
342 return ext::shared_ptr<MarketModelEvolver>(new
343 LogNormalCmSwapRatePc(spanningForwards,
344 marketModel, generatorFactory,
345 numeraires,
346 initialStep));
347 default:
348 QL_FAIL("unknown ConstantMaturitySwapMarketModelEvolver type");
349 }
350 }
351
352
353 void
354 checkCMSAndSwaptions(const SequenceStatisticsInc& stats,
355 const Rate fixedRate,
356 const std::vector<ext::shared_ptr<StrikedTypePayoff> >& displacedPayoff,
357 const ext::shared_ptr<MarketModel>&, // marketModel,
358 const std::string& config) {
359 std::vector<Real> results = stats.mean();
360 std::vector<Real> errors = stats.errorEstimate();
361 std::vector<Real> discrepancies(todaysForwards.size());
362
363 Size N = todaysForwards.size();
364 // check Swaps
365 Real maxError = QL_MIN_REAL;
366 LMMCurveState curveState_lmm(rateTimes);
367 curveState_lmm.setOnForwardRates(fwdRates: todaysForwards);
368
369 std::vector<Real> expectedNPVs(todaysCMSwapRates.size());
370 Real errorThreshold = 0.5;
371 for (Size i=0; i<N; ++i) {
372 Real expectedNPV = curveState_lmm.cmSwapAnnuity(numeraire: i, i, spanningForwards)
373 * (todaysCMSwapRates[i]-fixedRate) * todaysDiscounts[i];
374 expectedNPVs[i] = expectedNPV;
375 discrepancies[i] = (results[i]-expectedNPVs[i])/errors[i];
376 maxError = std::max(a: std::fabs(x: discrepancies[i]), b: maxError);
377 }
378
379 if (maxError > errorThreshold) {
380 BOOST_TEST_MESSAGE(config);
381 for (Size i=0; i<N; ++i) {
382 BOOST_TEST_MESSAGE(io::ordinal(i+1) << " CMS NPV: "
383 << io::rate(results[i])
384 << " +- " << io::rate(errors[i])
385 << "; expected: " << io::rate(expectedNPVs[i])
386 << "; discrepancy/error = "
387 << discrepancies[N-1-i]
388 << " standard errors");
389 }
390 BOOST_ERROR("test failed");
391 }
392
393 // check Swaptions
394 maxError = 0;
395
396 std::vector<Rate> expectedSwaptions(N);
397 for (Size i=0; i<N; ++i) {
398 Real expectedSwaption =
399 BlackCalculator(displacedPayoff[i],
400 todaysCMSwapRates[i]+displacement,
401 volatilities[i]*std::sqrt(x: rateTimes[i]),
402 curveState_lmm.cmSwapAnnuity(numeraire: i,i, spanningForwards)
403 * todaysDiscounts[i]).value();
404 expectedSwaptions[i] = expectedSwaption;
405 discrepancies[i] = (results[N+i]-expectedSwaptions[i])/errors[N+i];
406 maxError = std::max(a: std::fabs(x: discrepancies[i]), b: maxError);
407 }
408 errorThreshold = 2.0;
409
410 if (maxError > errorThreshold) {
411 BOOST_TEST_MESSAGE(config);
412 for (Size i=1; i<=N; ++i) {
413 BOOST_TEST_MESSAGE(io::ordinal(i) << " Swaption: "
414 << io::rate(results[2*N-i])
415 << " +- " << io::rate(errors[2*N-i])
416 << "; expected: " << io::rate(expectedSwaptions[N-i])
417 << "; discrepancy/error = "
418 << io::percent(discrepancies[N-i])
419 << " standard errors");
420 }
421 BOOST_ERROR("test failed");
422 }
423 }
424
425}
426
427
428void MarketModelCmsTest::testMultiStepCmSwapsAndSwaptions() {
429
430 BOOST_TEST_MESSAGE("Testing exact repricing of "
431 "multi-step constant maturity swaps and swaptions "
432 "in a lognormal constant maturity swap market model...");
433
434 using namespace market_model_cms_test;
435
436 setup();
437
438 Real fixedRate = 0.04;
439
440 // swaps
441 std::vector<Time> swapPaymentTimes(rateTimes.begin()+1, rateTimes.end());
442 // until ConstantMaturitySwap is ready
443 MultiStepCoterminalSwaps swaps(rateTimes, accruals, accruals,
444 swapPaymentTimes,
445 fixedRate);
446 // swaptions
447 std::vector<Time> swaptionPaymentTimes(rateTimes.begin(), rateTimes.end()-1);
448 std::vector<ext::shared_ptr<StrikedTypePayoff> >
449 displacedPayoff(todaysForwards.size()), undisplacedPayoff(todaysForwards.size());
450 for (Size i=0; i<undisplacedPayoff.size(); ++i) {
451 displacedPayoff[i] = ext::shared_ptr<StrikedTypePayoff>(new
452 PlainVanillaPayoff(Option::Call, fixedRate+displacement));
453
454 undisplacedPayoff[i] = ext::shared_ptr<StrikedTypePayoff>(new
455 PlainVanillaPayoff(Option::Call, fixedRate));
456 }
457
458 // until ConstantMaturitySwap is ready
459 MultiStepCoterminalSwaptions swaptions(rateTimes,
460 swaptionPaymentTimes,
461 undisplacedPayoff);
462 MultiProductComposite product;
463 product.add(swaps);
464 product.add(swaptions);
465 product.finalize();
466
467 EvolutionDescription evolution = product.evolution();
468
469 MarketModelType marketModels[] = {// CalibratedMM,
470 ExponentialCorrelationFlatVolatility,
471 ExponentialCorrelationAbcdVolatility };
472
473 for (auto& j : marketModels) {
474
475 Size testedFactors[] = { /*4, 8,*/ todaysForwards.size()};
476 for (unsigned long factors : testedFactors) {
477 // Composite's ProductSuggested is the Terminal one
478 MeasureType measures[] = { // ProductSuggested,
479 Terminal,
480 // MoneyMarketPlus,
481 MoneyMarket};
482 for (auto& measure : measures) {
483 std::vector<Size> numeraires = makeMeasure(product, measureType: measure);
484
485 ext::shared_ptr<MarketModel> marketModel = makeMarketModel(evolution, numberOfFactors: factors, marketModelType: j);
486
487 EvolverType evolvers[] = { Pc/*, Ipc*/ };
488
489 ext::shared_ptr<MarketModelEvolver> evolver;
490 Size stop = isInTerminalMeasure(evolution, numeraires) ? 0 : 1;
491 for (Size i=0; i<LENGTH(evolvers)-stop; i++) {
492 for (Size n=0; n<1; n++) {
493 //MTBrownianGeneratorFactory generatorFactory(seed_);
494 SobolBrownianGeneratorFactory generatorFactory(
495 SobolBrownianGenerator::Diagonal, seed_);
496 evolver = makeMarketModelEvolver(marketModel,
497 numeraires,
498 generatorFactory,
499 evolverType: evolvers[i]);
500 std::ostringstream config;
501 config << marketModelTypeToString(type: j) << ", " << factors
502 << (factors > 1 ?
503 (factors == todaysForwards.size() ? " (full) factors, " :
504 " factors, ") :
505 " factor,")
506 << measureTypeToString(type: measure) << ", "
507 << evolverTypeToString(type: evolvers[i]) << ", "
508 << "MT BGF";
509 if (printReport_)
510 BOOST_TEST_MESSAGE(" " << config.str());
511
512 ext::shared_ptr<SequenceStatisticsInc> stats = simulate(evolver, product);
513 checkCMSAndSwaptions(stats: *stats, fixedRate,
514 displacedPayoff, marketModel,config: config.str());
515
516 }
517 }
518 }
519 }
520 }
521}
522
523
524
525// --- Call the desired tests
526test_suite* MarketModelCmsTest::suite(SpeedLevel speed) {
527 auto* suite = BOOST_TEST_SUITE("CMS Market-model tests");
528
529 if (speed == Slow) {
530 suite->add(QUANTLIB_TEST_CASE(
531 &MarketModelCmsTest::testMultiStepCmSwapsAndSwaptions));
532 }
533
534 return suite;
535}
536

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