QuantLib: a free/open-source library for quantitative finance
fully annotated source code - version 1.34
Loading...
Searching...
No Matches
capletcoterminalswaptioncalibration.cpp
Go to the documentation of this file.
1/* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
2
3/*
4 Copyright (C) 2007 Ferdinando Ametrano
5 Copyright (C) 2007 Mark Joshi
6
7 This file is part of QuantLib, a free-software/open-source library
8 for financial quantitative analysts and developers - http://quantlib.org/
9
10 QuantLib is free software: you can redistribute it and/or modify it
11 under the terms of the QuantLib license. You should have received a
12 copy 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
16 This program is distributed in the hope that it will be useful, but WITHOUT
17 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
18 FOR A PARTICULAR PURPOSE. See the license for more details.
19*/
20
25
26namespace QuantLib {
27
29 const EvolutionDescription& evolution,
30 const ext::shared_ptr<PiecewiseConstantCorrelation>& corr,
31 const std::vector<ext::shared_ptr<
33 displacedSwapVariances,
34 const std::vector<Volatility>& mktCapletVols,
35 const ext::shared_ptr<CurveState>& cs,
36 Spread displacement,
37 const std::vector<Real>& alpha,
38 bool lowestRoot,
39 bool useFullApprox)
40 : CTSMMCapletCalibration(evolution, corr, displacedSwapVariances,
41 mktCapletVols, cs, displacement),
42 alpha_(alpha), lowestRoot_(lowestRoot),
43 useFullApprox_(useFullApprox) {
44
46 "mismatch between number of rates (" << numberOfRates_ <<
47 ") and alpha (" << alpha.size() << ")");
48
49 }
50
52 const EvolutionDescription& evolution,
54 const std::vector<ext::shared_ptr<
56 displacedSwapVariances,
57 const std::vector<Volatility>& capletVols,
58 const CurveState& cs,
59 Spread displacement,
60
61 const std::vector<Real>& alpha,
62 bool lowestRoot,
63 bool useFullAprox,
64
65 Size numberOfFactors,
66 //Size maxIterations,
67 //Real tolerance,
68
69 std::vector<Matrix>& swapCovariancePseudoRoots) {
70
72 displacedSwapVariances, capletVols, cs);
73
74 Size numberOfSteps = evolution.numberOfSteps();
75 Size numberOfRates = evolution.numberOfRates();
76 const std::vector<Time>& rateTimes = evolution.rateTimes();
77
78 QL_REQUIRE(numberOfFactors<=numberOfRates,
79 "number of factors (" << numberOfFactors <<
80 ") cannot be greater than numberOfRates (" <<
81 numberOfRates << ")");
82 QL_REQUIRE(numberOfFactors>0,
83 "number of factors (" << numberOfFactors <<
84 ") must be greater than zero");
85
86 Natural failures = 0;
87 Real extraMultiplier = useFullAprox ? 1.0 : 0.0;
88
89 // factor reduction
90 std::vector<Matrix> corrPseudo(corr.times().size());
91 for (Size i=0; i<corrPseudo.size(); ++i)
92 corrPseudo[i] = rankReducedSqrt(corr.correlation(i),
93 numberOfFactors, 1.0,
95
96 Matrix zedMatrix =
98 Matrix invertedZedMatrix = inverse(zedMatrix);
99
100
101
102 // do alpha part
103 // first modify variances to take account of alpha
104 // then rescale so total variance is unchanged
105 Matrix swapTimeInhomogeneousVariances(numberOfSteps, numberOfRates, 0.0);
106 std::vector<Real> originalVariances(numberOfRates, 0.0);
107 std::vector<Real> modifiedVariances(numberOfRates, 0.0);
108 const std::vector<Time>& evolutionTimes = evolution.evolutionTimes();
109 for (Size i=0; i<numberOfSteps; ++i) {
110 Real s = (i==0 ? 0.0 : evolutionTimes[i-1]);
111 for (Size j=i; j<numberOfRates; ++j) {
112 const std::vector<Real>& var = displacedSwapVariances[j]->variances();
113 originalVariances[j]+=var[i];
114 swapTimeInhomogeneousVariances[i][j] = var[i]/
115 ((1.0+alpha[j]*s)*(1.0+alpha[j]*s));
116 modifiedVariances[j]+=swapTimeInhomogeneousVariances[i][j];
117 }
118 }
119
120 for (Size i=0; i<numberOfSteps; ++i)
121 for (Size j=i; j<numberOfRates; ++j)
122 swapTimeInhomogeneousVariances[i][j] *= originalVariances[j]/
123 modifiedVariances[j];
124
125
126 // compute swap covariances for caplet approximation formula
127 // without taking into account A and B
128 std::vector<Matrix> CovarianceSwapPseudos(numberOfSteps);
129 std::vector<Matrix> CovarianceSwapCovs(numberOfSteps); // this is total cov
130 std::vector<Matrix> CovarianceSwapMarginalCovs(numberOfSteps); // this is cov for one step
131
132 for (Size i=0; i<numberOfSteps; ++i) {
133 CovarianceSwapPseudos[i] = corrPseudo[i];
134 for (Size j=0; j<numberOfRates; ++j)
135 for (Size k=0; k < CovarianceSwapPseudos[i].columns(); ++k)
136 CovarianceSwapPseudos[i][j][k] *=
137 std::sqrt(swapTimeInhomogeneousVariances[i][j]);
138
139 CovarianceSwapMarginalCovs[i] = CovarianceSwapPseudos[i] *
140 transpose(CovarianceSwapPseudos[i]);
141
142 CovarianceSwapCovs[i] = CovarianceSwapMarginalCovs[i];
143 if (i>0)
144 CovarianceSwapCovs[i]+= CovarianceSwapCovs[i-1];
145 }
146
147 // compute partial variances and covariances which will take A and B coefficients
148 //const std::vector<Time>& taus = evolution.rateTaus();
149 std::vector<Real> totVariance(numberOfRates, 0.0);
150 std::vector<Real> almostTotVariance(numberOfRates, 0.0);
151 std::vector<Real> almostTotCovariance(numberOfRates, 0.0);
152 std::vector<Real> leftCovariance(numberOfRates, 0.0);
153 for (Size i=0; i<numberOfRates; ++i) {
154 for (Size jj=0; jj<=i; ++jj)
155 totVariance[i] += displacedSwapVariances[i]->variances()[jj];
156 Integer j;
157 for (j=0; j<=static_cast<Integer>(i)-1; ++j)
158 almostTotVariance[i] += swapTimeInhomogeneousVariances[j][i];
159 for (j=0; j<=static_cast<Integer>(i)-2; ++j) {
160 const Matrix& thisPseudo = corrPseudo[j];
161 Real correlation = 0.0;
162 for (Size k=0; k<numberOfFactors; ++k)
163 correlation += thisPseudo[i-1][k]*thisPseudo[i][k];
164 almostTotCovariance[i] += correlation *
165 std::sqrt(swapTimeInhomogeneousVariances[j][i] *
166 swapTimeInhomogeneousVariances[j][i-1]);
167 }
168 if (i>0) {
169 const Matrix& thisPseudo = corrPseudo[j];
170 Real correlation = 0.0;
171 for (Size k=0; k<numberOfFactors; ++k)
172 correlation += thisPseudo[i-1][k]*thisPseudo[i][k];
173 leftCovariance[i] = correlation *
174 std::sqrt(swapTimeInhomogeneousVariances[j][i] *
175 swapTimeInhomogeneousVariances[j][i-1]);
176 }
177 }
178
179
180 // multiplier up to rate reset previous time
181 // the first element is not used
182 std::vector<Real> a(numberOfSteps, 1.0);
183 // multiplier afterward
184 std::vector<Real> b(numberOfSteps);
185
186 b[0]=displacedSwapVariances[0]->variances()[0]/swapTimeInhomogeneousVariances[0][0];
187 // main loop where work is done
188 for (Size i=1; i<numberOfSteps; ++i) {
189 Integer j=0;
190 // up date variances to take account of last a and b computed
191 for (; j <= static_cast<Integer>(i)-2; j++)
192 swapTimeInhomogeneousVariances[j][i-1]*= a[i-1]*a[i-1];
193 swapTimeInhomogeneousVariances[j][i-1]*= b[i-1]*b[i-1];
194
195 Real w0 = invertedZedMatrix[i-1][i-1];
196 Real w1 = -invertedZedMatrix[i-1][i];
197 Real v1t1 = capletVols[i-1]*capletVols[i-1]*rateTimes[i-1];
198
199 // now compute contribution from lower right corner
200 Real extraConstantPart =0.0;
201 // part of caplet approximation formula coming from later rates
202 for (Size k = i+1; k < numberOfSteps; ++k)
203 for (Size l = i+1; l < numberOfSteps; ++l)
204 extraConstantPart+=invertedZedMatrix[i-1][k] *
205 CovarianceSwapCovs[i-1][k][l] *
206 invertedZedMatrix[i-1][l];
207
208 // now compute contribution from top row excluding first two columns and its transpose
209 // we divide into two parts as one part is multiplied by a[i-1] and the other by b[i-1]
210 // a lot could be precomputed when we need to optimize
211 for (Size k = i+1; k < numberOfSteps; ++k)
212 {
213 if (i > 1)
214 {
215 extraConstantPart+=invertedZedMatrix[i-1][i-1] *
216 CovarianceSwapCovs[i-2][i-1][k] *
217 invertedZedMatrix[i-1][k]*a[i-1];
218
219 extraConstantPart+=invertedZedMatrix[i-1][k] *
220 CovarianceSwapCovs[i-2][k][i-1] *
221 invertedZedMatrix[i-1][i-1]*a[i-1];
222 }
223
224 extraConstantPart+=invertedZedMatrix[i-1][i-1] *
225 CovarianceSwapMarginalCovs[i-1][i-1][k] *
226 invertedZedMatrix[i-1][k]*b[i-1];
227
228 extraConstantPart+=invertedZedMatrix[i-1][k] *
229 CovarianceSwapCovs[i-1][k][i-1] *
230 invertedZedMatrix[i-1][i-1]*b[i-1];
231
232 }
233
234 // we also get an extra linear part, this corresponds to row i, and columns j>i+1, and transpose
235 Real extraLinearPart=0.0;
236 for (Size k = i+1; k < numberOfSteps; ++k)
237 {
238 extraLinearPart+=invertedZedMatrix[i-1][k] *
239 CovarianceSwapCovs[i-1][k][i] *
240 invertedZedMatrix[i-1][i];
241
242 extraLinearPart+=invertedZedMatrix[i-1][i] *
243 CovarianceSwapCovs[i-1][i][k] *
244 invertedZedMatrix[i-1][k];
245 }
246
247
248
249 Real constantPart = w0*w0*totVariance[i-1] +
250 extraConstantPart*extraMultiplier-v1t1;
251 Real linearPart = -2*w0*w1*(a[i-1]*almostTotCovariance[i] +
252 b[i-1]*leftCovariance[i]) +extraLinearPart*extraMultiplier ;
253 Real quadraticPart = w1*w1*almostTotVariance[i];
254 Real disc = linearPart*linearPart-4.0*constantPart*quadraticPart;
255
256 Real root, minimum = -linearPart/(2.0*quadraticPart);
257 bool rightUsed = false;
258 if (disc <0.0) {
259 ++failures;
260 // pick up the minimum vol for the caplet
261 root = minimum;
262 } else if (lowestRoot) {
263 root = (-linearPart-std::sqrt(disc))/(2.0*quadraticPart);
264 } else {
265 if (minimum>1.0)
266 root = (-linearPart-std::sqrt(disc))/(2.0*quadraticPart);
267 else {
268 rightUsed = true;
269 root = (-linearPart+std::sqrt(disc))/(2.0*quadraticPart);
270 }
271 }
272
273 Real varianceFound = root*root*almostTotVariance[i];
274 Real varianceToFind = totVariance[i]-varianceFound;
275 Real mult = varianceToFind/swapTimeInhomogeneousVariances[i][i];
276 if (mult<=0.0 && rightUsed) {
277 root = (-linearPart-std::sqrt(disc))/(2.0*quadraticPart);
278 varianceFound = root*root*almostTotVariance[i];
279 varianceToFind = totVariance[i]-varianceFound;
280 mult = varianceToFind/swapTimeInhomogeneousVariances[i][i];
281 }
282
283 if (mult<0.0) // no solution...
284 {
285 ++failures;
286 a[i]=root;
287 b[i]=0.0;
288 }
289 else
290 {
291 a[i]=root;
292 b[i]=std::sqrt(mult);
293 }
294
295 QL_ENSURE(root>=0.0,
296 "negative root -- it should have not happened");
297
298 }
299
300 {
301 Integer i = numberOfSteps;
302 Integer j=0;
303 for (; j <= static_cast<Integer>(i)-2; j++)
304 swapTimeInhomogeneousVariances[j][i-1]*= a[i-1]*a[i-1];
305 swapTimeInhomogeneousVariances[j][i-1]*= b[i-1]*b[i-1];
306 }
307
308 // compute the results
309 swapCovariancePseudoRoots.resize(numberOfSteps);
310 for (Size k=0; k<numberOfSteps; ++k) {
311 swapCovariancePseudoRoots[k] = corrPseudo[k];
312 for (Size j=0; j<numberOfRates; ++j) {
313 Real coeff = std::sqrt(swapTimeInhomogeneousVariances[k][j]);
314 for (Size i=0; i<numberOfFactors; ++i)
315 swapCovariancePseudoRoots[k][j][i]*=coeff;
316 }
317 QL_ENSURE(swapCovariancePseudoRoots[k].rows()==numberOfRates,
318 "step " << k
319 << " abcd vol wrong number of rows: "
320 << swapCovariancePseudoRoots[k].rows()
321 << " instead of " << numberOfRates);
322 QL_ENSURE(swapCovariancePseudoRoots[k].columns()==numberOfFactors,
323 "step " << k
324 << " abcd vol wrong number of columns: "
325 << swapCovariancePseudoRoots[k].columns()
326 << " instead of " << numberOfFactors);
327 }
328
329 return failures;
330 }
331
333 Natural numberOfFactors,
334 Natural ,
335 Real ) {
336
338 *corr_,
340 // not mktCapletVols_ but...
342 *cs_,
344
345 alpha_,
348
349 numberOfFactors,
350
352 }
353
354}
ext::shared_ptr< PiecewiseConstantCorrelation > corr_
ext::shared_ptr< CurveState > cs_
static void performChecks(const EvolutionDescription &evolution, const PiecewiseConstantCorrelation &corr, const std::vector< ext::shared_ptr< PiecewiseConstantVariance > > &displacedSwapVariances, const std::vector< Volatility > &mktCapletVols, const CurveState &cs)
std::vector< ext::shared_ptr< PiecewiseConstantVariance > > displacedSwapVariances_
std::vector< Volatility > usedCapletVols_
std::vector< Matrix > swapCovariancePseudoRoots_
Natural calibrationImpl_(Natural numberOfFactors, Natural, Real) override
static Natural calibrationFunction(const EvolutionDescription &evolution, const PiecewiseConstantCorrelation &corr, const std::vector< ext::shared_ptr< PiecewiseConstantVariance > > &displacedSwapVariances, const std::vector< Volatility > &capletVols, const CurveState &cs, Spread displacement, const std::vector< Real > &alpha, bool lowestRoot, bool useFullApprox, Size numberOfFactors, std::vector< Matrix > &swapCovariancePseudoRoots)
CTSMMCapletOriginalCalibration(const EvolutionDescription &evolution, const ext::shared_ptr< PiecewiseConstantCorrelation > &corr, const std::vector< ext::shared_ptr< PiecewiseConstantVariance > > &displacedSwapVariances, const std::vector< Volatility > &capletVols, const ext::shared_ptr< CurveState > &cs, Spread displacement, const std::vector< Real > &alpha, bool lowestRoot, bool useFullApprox)
Curve state for market-model simulations
Definition: curvestate.hpp:41
Market-model evolution description.
const std::vector< Time > & rateTimes() const
const std::vector< Time > & evolutionTimes() const
Matrix used in linear algebra.
Definition: matrix.hpp:41
virtual const Matrix & correlation(Size i) const
virtual const std::vector< Time > & times() const =0
static Matrix coterminalSwapZedMatrix(const CurveState &cs, Spread displacement)
#define QL_ENSURE(condition, message)
throw an error if the given post-condition is not verified
Definition: errors.hpp:130
#define QL_REQUIRE(condition, message)
throw an error if the given pre-condition is not verified
Definition: errors.hpp:117
ext::function< Real(Real)> b
QL_REAL Real
real number
Definition: types.hpp:50
unsigned QL_INTEGER Natural
positive integer
Definition: types.hpp:43
QL_INTEGER Integer
integer number
Definition: types.hpp:35
Real Spread
spreads on interest rates
Definition: types.hpp:74
std::size_t Size
size of a container
Definition: types.hpp:58
Definition: any.hpp:35
Matrix rankReducedSqrt(const Matrix &matrix, Size maxRank, Real componentRetainedPercentage, SalvagingAlgorithm::Type sa)
Definition: pseudosqrt.cpp:427
Matrix transpose(const Matrix &m)
Definition: matrix.hpp:700
Matrix inverse(const Matrix &m)
Definition: matrix.cpp:44
pseudo square root of a real symmetric matrix
Real alpha
Definition: sabr.cpp:200
Utility functions for mapping between swap rate and forward rate.