QuantLib: a free/open-source library for quantitative finance
fully annotated source code - version 1.34
Loading...
Searching...
No Matches
saddlepointlossmodel.hpp
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) 2014 Jose Aparicio
5
6 This file is part of QuantLib, a free-software/open-source library
7 for financial quantitative analysts and developers - http://quantlib.org/
8
9 QuantLib is free software: you can redistribute it and/or modify it
10 under the terms of the QuantLib license. You should have received a
11 copy of the license along with this program; if not, please email
12 <quantlib-dev@lists.sf.net>. The license is also available online at
13 <http://quantlib.org/license.shtml>.
14
15 This program is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17 FOR A PARTICULAR PURPOSE. See the license for more details.
18*/
19
20#ifndef quantlib_saddle_point_lossmodel_hpp
21#define quantlib_saddle_point_lossmodel_hpp
22
23#include <ql/tuple.hpp>
29
30namespace QuantLib {
31
32 /*! \brief Saddle point portfolio credit default loss model.\par
33 Default Loss model implementing the Saddle point expansion
34 integrations on several default risk metrics. Codepence is dealt
35 through a latent model making the integrals conditional to the latent
36 model factor. Latent variables are integrated indirectly.\par
37 See:\par
38 <b>Taking to the saddle</b> by R.Martin, K.Thompson and C.Browne; RISK JUNE
39 2001; p.91\par
40 <b>The saddlepoint method and portfolio optionalities</b> R.Martin in Risk
41 December 2006\par
42 <b>VAR: who contributes and how much?</b> R.Martin, K.Thompson and
43 C.Browne RISK AUGUST 2001\par
44 <b>Shortfall: Who contributes and how much?</b> R. J. Martin, Credit Suisse
45 January 3, 2007 \par
46 <b>Don't Fall from the Saddle: the Importance of Higher Moments of Credit
47 Loss Distributions</b> J.Annaert, C.Garcia Joao Batista, J.Lamoot,
48 G.Lanine February 2006, Gent University\par
49 <b>Analytical techniques for synthetic CDOs and credit default risk
50 measures</b> A. Antonov, S. Mechkovy, and T. Misirpashaevz;
51 NumeriX May 23, 2005 \par
52 <b>Computation of VaR and VaR contribution in the Vasicek portfolio credit
53 loss model: a comparative study</b> X.Huang, C.W.Oosterlee, M.Mesters
54 Journal of Credit Risk (75-96) Volume 3/ Number 3, Fall 2007 \par
55 <b>Higher-order saddlepoint approximations in the Vasicek portfolio credit
56 loss model</b> X.Huang, C.W.Oosterlee, M.Mesters Journal of
57 Computational Finance (93-113) Volume 11/Number 1, Fall 2007 \par
58 While more expensive, a high order expansion is used here; see the paper by
59 Antonov et al for the terms retained.\par
60 For a discussion of an alternative to fix the error at low loss levels
61 (more relevant to pricing than risk metrics) see: \par
62 <b>The hybrid saddlepoint method for credit portfolios</b> by A.Owen,
63 A.McLeod and K.Thompson; in Risk, August 2009. This is not implemented here
64 though (yet?...)\par
65 For the more general context mathematical theory see: <b>Saddlepoint
66 approximations with applications</b> by R.W. Butler, Cambridge series in
67 statistical and probabilistic mathematics. 2007 \par
68 \todo Some portfolios show instabilities in the high order expansion terms.
69 \todo Methods here are calling and integrating using the unconditional
70 probabilities without inverting them first; quite a lot of calls to
71 the copula inversion can be avoided; this should improve performance.
72 \todo Revise the model for stability of the saddle point calculation. The
73 search for the point does not convege in extreme cases; e.g. very high
74 value of all the factors; factors for each variable not ordered from
75 high to low,...
76 */
77
78 /* The treatment of recovery wont work with random recoveries, they should
79 be passed to the conditional methods in the same way as the probabilities.
80 */
81
82 /*
83 TO DO:
84 -> Failing when the tranche upper loss limit goes over the max attainable
85 loss.
86
87 - With 15 quadrature points things are OK but 25 gives me -1#IND000 errors
88 (over region around the EL I think)
89 - Silly bug when calling some methods on todays date (zero time).
90 ProbDef = 0 there
91 - VaR <- tranched?????!
92 - ESF <- tranched????!!
93 - VaR split
94 - ESF split?
95
96 When introducing defaults; somewhere, (after an update?) there should be
97 a check that: copula_->basketSize() EQUALS remainingBasket_.size()
98 */
99 template<class CP>
101 public:
103 const ext::shared_ptr<ConstantLossLatentmodel<CP> >& m)
104 : copula_(m) { }
105 protected:
106 // ----------- Cumulants and derivatives auxiliary functions ---------
107
108 /*! Returns the cumulant generating function (zero-th order
109 expansion term) conditional to the mkt factor:
110 \f$ K = \sum_j ln(1-p_j + p_j e^{N_j \times lgd_j \times s}) \f$
111 */
113 const std::vector<Real>& invUncondProbs,
114 Real lossFraction,// saddle pt
115 const std::vector<Real>& mktFactor) const;
116 /*! Returns the first derivative of the cumulant generating function
117 (first order expansion term) conditional to the mkt factor:
118 \f$ K1 = \sum_j \frac{p_j \times N_j \times LGD_j \times
119 e^{N_j \times LGD_j \times s}} \
120 {1-p_j + p_j e^{N_j \times LGD_j \times s}} \f$
121 One of its properties is that its value at zero is the portfolio
122 expected loss (in fractional units). Its value at infinity is the
123 max attainable portfolio loss. To be understood conditional to the
124 market factor.
125 */
127 const std::vector<Real>& invUncondProbs,
128 Real saddle, // in fract loss units... humm not really
129 const std::vector<Real>& mktFactor) const;
130 /*! Returns the second derivative of the cumulant generating function
131 (first order expansion term) conditional to the mkt factor:
132 \f$ K2 = \sum_j \frac{p_j \times (N_j \times LGD_j)^2 \times
133 e^{N_j \times LGD_j \times s}}
134 {1-p_j + p_j e^{N_j \times LGD_j \times s}}
135 - (\frac{p_j \times N_j \times LGD_j \times e^{N_j \times
136 LGD_j \times s}}
137 {1-p_j + p_j e^{N_j \times LGD_j \times s}})^2 \f$
138 */
140 const std::vector<Real>& invUncondProbs,
141 Real saddle,
142 const std::vector<Real>& mktFactor) const;
144 const std::vector<Real>& invUncondProbs,
145 Real saddle,
146 const std::vector<Real>& mktFactor) const;
148 const std::vector<Real>& invUncondProbs,
149 Real saddle,
150 const std::vector<Real>& mktFactor) const ;
151 /*! Returns the cumulant and second to fourth derivatives together.
152 Included for optimization, most methods work on expansion of these
153 terms.
154 Alternatively use a local private buffer member? */
155 ext::tuple<Real, Real, Real, Real> CumGen0234DerivCond(
156 const std::vector<Real>& invUncondProbs,
157 Real saddle,
158 const std::vector<Real>& mktFactor) const;
159 ext::tuple<Real, Real> CumGen02DerivCond(
160 const std::vector<Real>& invUncondProbs,
161 Real saddle,
162 const std::vector<Real>& mktFactor) const;
163
164 /* Unconditional cumulants. Because this class integrates the various
165 statistics it provides in indirect mode they are never used.
166 Provided for completeness/extendability */
167 /*! Returns the cumulant generating function (zero-th order expansion
168 term) weighting the conditional value by the prob density of the
169 market factor, called by integrations */
170 Real CumulantGenerating(const Date& date, Real s) const;
171 Real CumGen1stDerivative(const Date& date, Real s) const;
172 Real CumGen2ndDerivative(const Date& date, Real s) const;
173 Real CumGen3rdDerivative(const Date& date, Real s) const;
174 Real CumGen4thDerivative(const Date& date, Real s) const;
175
176 // -------- Saddle point search functions ---------------------------
180 const std::vector<Real>& mktFactor_;
181 const std::vector<Real>& invUncondProbs_;
182 public:
183 //! The passed target is in fractional loss units
185 const Real target,
186 const std::vector<Real>& invUncondProbs,
187 const std::vector<Real>& mktFactor
188 )
189 : me_(me),
190 targetValue_(target),
191 mktFactor_(mktFactor),
192 invUncondProbs_(invUncondProbs)
193 {}
194 Real operator()(const Real x) const {
197 }
200 mktFactor_);
201 }
202 };
203
204 /*! Calculates the mkt-fct-conditional saddle point for the loss level
205 given and the probability passed.
206 The date is implicitly given through the probability. Performance
207 requires to pass the probabilities for that date. Otherwise once we
208 integrate this over the market factor we would be computing the same
209 probabilities over and over. While this works fine here some models
210 of the recovery rate might require the date.
211
212 The passed lossLevel is in total portfolio loss fractional units.
213
214 \todo Improve convergence speed (which is bad at the moment).See
215 discussion in several places; references above and The Oxford
216 Handbook of CD, sect 2.9
217 */
219 const std::vector<Real>& invUncondProbs,
220 Real lossLevel,
221 const std::vector<Real>& mktFactor,
222 Real accuracy = 1.0e-3,//1.e-4
223 Natural maxEvaluations = 50
224 ) const;
225
230 public:
232 const SaddlePointLossModel& me,
233 const Real target,
234 const Date& date)
235 : me_(me), targetValue_(1.-target), date_(date) {}
236 /*! The passed x is the _tranche_ loss fraction */
237 Real operator()(const Real x) const {
238 return me_.probOverLoss(date_, x) - targetValue_;
239 }
240 };
241 // Functionality, Provides various portfolio statistics---------------
242 public:
243 /*! Returns the loss amount at the requested date for which the
244 probability of lossing that amount or less is equal to the value passed.
245 */
246 Real percentile(const Date& d, Probability percentile) const override;
247
248 protected:
249 /*! Conditional (on the mkt factor) prob of a loss fraction of
250 the the tranched portfolio.
251
252 The trancheLossFract parameter is the fraction over the
253 tranche notional and must be in [0,1].
254 */
256 const std::vector<Real>& invUncondProbs,
257 Real trancheLossFract,
258 const std::vector<Real>& mktFactor) const;
260 const std::vector<Real>& invUncondProbs,
261 Real loss,
262 const std::vector<Real>& mktFactor) const;
263 public:
264 Probability probOverLoss(const Date& d, Real trancheLossFract) const override;
265
266 std::map<Real, Probability> lossDistribution(const Date& d) const override;
267
268 protected:
269 /*!
270 Probability of having losses in the portfolio due to default
271 events equal or larger than a given absolute loss value on a
272 given date conditional to the latent model factor.
273 The integral expression on the expansion is the first order
274 integration as presented in several references, see for instance;
275 equation 8 in R.Martin, K.Thompson, and C. Browne 's
276 'Taking to the Saddle', Risk Magazine, June 2001, page 91
277
278 The passed loss is in absolute value.
279 */
281 const std::vector<Real>& invUncondProbs,
282 Real loss,
283 const std::vector<Real>& mktFactor) const;
284 public:
285 Probability probOverPortfLoss(const Date& d, Real loss) const;
286 Real expectedTrancheLoss(const Date& d) const override;
287
288 protected:
289 /*!
290 Probability density of having losses in the total portfolio (untranched)
291 due to default events equal to a given value on a given date conditional
292 to the latent model factor.
293 Based on the integrals of the expected shortfall.
294 */
295 Probability probDensityCond(const std::vector<Real>& invUncondProbs,
296 Real loss, const std::vector<Real>& mktFactor) const;
297 public:
298 Probability probDensity(const Date& d, Real loss) const;
299 protected:
300 std::vector<Real> splitLossCond(
301 const std::vector<Real>& invUncondProbs,
302 Real loss, std::vector<Real> mktFactor) const;
304 const std::vector<Real>& invUncondProbs,
305 Real lossPerc, const std::vector<Real>& mktFactor) const;
307 const std::vector<Real>& invUncondProbs,
308 Real lossPerc, Probability percentile,
309 const std::vector<Real>& mktFactor) const;
310 std::vector<Real> expectedShortfallSplitCond(
311 const std::vector<Real>& invUncondProbs,
312 Real lossPerc, const std::vector<Real>& mktFactor) const;
313 public:
314 /*! Sensitivities of the individual names to a given portfolio loss
315 value due to defaults. It returns ratios to the total structure
316 notional, which aggregated add up to the requested loss value.
317 Notice then that it refers to the total portfolio, not the tranched
318 basket.
319 \todo Fix this.
320 \par
321 see equation 8 in <b>VAR: who contributes and how much?</b> by
322 R.Martin, K.Thompson, and C. Browne in Risk Magazine, August 2001
323
324 The passed loss is the loss amount level at which we want
325 to request the sensitivity. Equivalent to a percentile.
326 */
327 std::vector<Real> splitVaRLevel(const Date& date, Real loss) const override;
328 Real expectedShortfall(const Date& d, Probability percentile) const override;
329
330 protected:
332 const std::vector<Real>& invUncondProbs,
333 const std::vector<Real>& mktFactor) const;
335 const std::vector<Real>& invUncondProbs,
336 const std::vector<Real>& mktFactor) const;
337
338 void resetModel() override {
339 remainingNotionals_ = basket_->remainingNotionals();
340 remainingNotional_ = basket_->remainingNotional();
341 attachRatio_ = std::min(basket_->remainingAttachmentAmount()
342 / basket_->remainingNotional(), 1.);
343 detachRatio_ = std::min(basket_->remainingDetachmentAmount()
344 / basket_->remainingNotional(), 1.);
345 copula_->resetBasket(basket_.currentLink());
346 }
347
348 const ext::shared_ptr<ConstantLossLatentmodel<CP> > copula_;
349 // cached todays arguments values
351 mutable std::vector<Real> remainingNotionals_;
353 // remaining basket levels:
355 /*
356 // Just for testing the ESF direct integration, not for release,
357 // this is very inneficient:
358 class ESFIntegrator {
359 public:
360 ESFIntegrator(const SaddlePointLossModel& me,
361 const Date& date,
362 Real lossPercentileFract//,
363 //const std::vector<Real>& mktFactor
364 )
365 : me_(me), date_(date),lossPercentileFract_(lossPercentileFract)
366 {}
367
368
369 Real operator()(Real x) const {
370 return me_.densityTrancheLoss(date_, x + lossPercentileFract_)
371 * (x + lossPercentileFract_);
372 }
373
374 Real lossPercentileFract_;
375 Date date_;
376 // const std::vector<Real>& mktFactor_;
377 const SaddlePointLossModel& me_;
378 };
379 */
380 };
381
382
383 // -- Inlined integrations------------------------------------------------
384
385 // Unconditional Moments and derivatives. --------------------------------
386 template<class CP>
388 const Date& date, Real s) const
389 {
390 std::vector<Real> invUncondProbs =
391 basket_->remainingProbabilities(date);
392 for(Size i=0; i<invUncondProbs.size(); i++)
393 invUncondProbs[i] =
394 copula_->inverseCumulativeY(invUncondProbs[i], i);
395
396 return copula_->integratedExpectedValue(
397 [&](const std::vector<Real>& v1) {
398 return CumulantGeneratingCond(invUncondProbs, s, v1);
399 });
400 }
401
402 template<class CP>
404 const Date& date, Real s) const
405 {
406 std::vector<Real> invUncondProbs =
407 basket_->remainingProbabilities(date);
408 for(Size i=0; i<invUncondProbs.size(); i++)
409 invUncondProbs[i] =
410 copula_->inverseCumulativeY(invUncondProbs[i], i);
411
412 return copula_->integratedExpectedValue(
413 [&](const std::vector<Real>& v1) {
414 return CumGen1stDerivativeCond(invUncondProbs, s, v1);
415 });
416 }
417
418 template<class CP>
420 const Date& date, Real s) const
421 {
422 std::vector<Real> invUncondProbs =
423 basket_->remainingProbabilities(date);
424 for(Size i=0; i<invUncondProbs.size(); i++)
425 invUncondProbs[i] =
426 copula_->inverseCumulativeY(invUncondProbs[i], i);
427
428 return copula_->integratedExpectedValue(
429 [&](const std::vector<Real>& v1) {
430 return CumGen2ndDerivativeCond(invUncondProbs, s, v1);
431 });
432 }
433
434 template<class CP>
436 const Date& date, Real s) const
437 {
438 std::vector<Real> invUncondProbs =
439 basket_->remainingProbabilities(date);
440 for(Size i=0; i<invUncondProbs.size(); i++)
441 invUncondProbs[i] =
442 copula_->inverseCumulativeY(invUncondProbs[i], i);
443
444 return copula_->integratedExpectedValue(
445 [&](const std::vector<Real>& v1) {
446 return CumGen3rdDerivativeCond(invUncondProbs, s, v1);
447 });
448 }
449
450 template<class CP>
452 const Date& date, Real s) const
453 {
454 std::vector<Real> invUncondProbs =
455 basket_->remainingProbabilities(date);
456 for(Size i=0; i<invUncondProbs.size(); i++)
457 invUncondProbs[i] =
458 copula_->inverseCumulativeY(invUncondProbs[i], i);
459
460 return copula_->integratedExpectedValue(
461 [&](const std::vector<Real>& v1) {
462 return CumGen4thDerivativeCond(invUncondProbs, s, v1);
463 });
464 }
465
466 template<class CP>
468 const Date& d, Real trancheLossFract) const
469 {
470 // avoid computation:
471 if (trancheLossFract >=
472 // time dependent soon:
473 basket_->detachmentAmount()) return 0.;
474
475 std::vector<Real> invUncondProbs =
476 basket_->remainingProbabilities(d);
477 for(Size i=0; i<invUncondProbs.size(); i++)
478 invUncondProbs[i] =
479 copula_->inverseCumulativeY(invUncondProbs[i], i);
480
481 return copula_->integratedExpectedValue(
482 [&](const std::vector<Real>& v1) {
483 return probOverLossCond(invUncondProbs, trancheLossFract, v1);
484 });
485 }
486
487 template<class CP>
489 const Date& d, Real loss) const
490 {
491 std::vector<Probability> invUncondProbs =
492 basket_->remainingProbabilities(d);
493 for(Size i=0; i<invUncondProbs.size(); i++)
494 invUncondProbs[i] =
495 copula_->inverseCumulativeY(invUncondProbs[i], i);
496
497 return copula_->integratedExpectedValue(
498 [&](const std::vector<Real>& v1) {
499 return probOverLossPortfCond(invUncondProbs, loss, v1);
500 });
501 }
502
503 template<class CP>
505 const Date& d) const
506 {
507 std::vector<Real> invUncondProbs =
508 basket_->remainingProbabilities(d);
509 for(Size i=0; i<invUncondProbs.size(); i++)
510 invUncondProbs[i] =
511 copula_->inverseCumulativeY(invUncondProbs[i], i);
512
513 return copula_->integratedExpectedValue(
514 [&](const std::vector<Real>& v1) {
515 return conditionalExpectedTrancheLoss(invUncondProbs, v1);
516 });
517 }
518
519 template<class CP>
521 const Date& d, Real loss) const
522 {
523 std::vector<Real> invUncondProbs =
524 basket_->remainingProbabilities(d);
525 for(Size i=0; i<invUncondProbs.size(); i++)
526 invUncondProbs[i] =
527 copula_->inverseCumulativeY(invUncondProbs[i], i);
528
529 return copula_->integratedExpectedValue(
530 [&](const std::vector<Real>& v1) {
531 return probDensityCond(invUncondProbs, loss, v1);
532 });
533 }
534
535 template<class CP>
536 inline std::vector<Real> SaddlePointLossModel<CP>::splitVaRLevel(const Date& date, Real s) const
537 {
538 std::vector<Real> invUncondProbs =
539 basket_->remainingProbabilities(date);
540 for(Size i=0; i<invUncondProbs.size(); i++)
541 invUncondProbs[i] =
542 copula_->inverseCumulativeY(invUncondProbs[i], i);
543
544 return copula_->integratedExpectedValueV(
545 [&](const std::vector<Real>& v1) {
546 return splitLossCond(invUncondProbs, s, v1);
547 });
548 }
549
550
551
552
553
554
555
556 /* ------------------------------------------------------------------------
557 Conditional Moments and derivatives.
558
559 Notice that in all this methods the date dependence is implicitly
560 present in the unconditional probabilities. But, as in other LMs, it
561 is redundant and expensive to perform the call to the probabilities in
562 these methods since they are integrands.
563 ---------------------------------------------------------------------- */
564
565 template<class CP>
567 const std::vector<Real>& invUncondProbs,
568 Real lossFraction,
569 const std::vector<Real>& mktFactor) const
570 {
571 const Size nNames = remainingNotionals_.size();
572 Real sum = 0.;
573
574 for(Size iName=0; iName < nNames; iName++) {
575 Probability pBuffer =
576 copula_->conditionalDefaultProbabilityInvP(
577 invUncondProbs[iName], iName, mktFactor);
578 sum += std::log(1. - pBuffer +
579 pBuffer * std::exp(remainingNotionals_[iName] *
580 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
581 iName, mktFactor)) * lossFraction / remainingNotional_));
582 }
583 return sum;
584 }
585
586 template<class CP>
588 const std::vector<Real>& invUncondProbs,
589 Real saddle,
590 const std::vector<Real>& mktFactor) const
591 {
592 const Size nNames = remainingNotionals_.size();
593 Real sum = 0.;
594
595 for(Size iName=0; iName < nNames; iName++) {
596 Probability pBuffer =
597 copula_->conditionalDefaultProbabilityInvP(
598 invUncondProbs[iName], iName, mktFactor);
599 // loss in fractional units
600 Real lossInDef = remainingNotionals_[iName] *
601 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
602 iName, mktFactor)) / remainingNotional_;
603 Real midFactor = pBuffer * std::exp(lossInDef * saddle);
604 sum += lossInDef * midFactor / (1.-pBuffer + midFactor);
605 }
606 return sum;
607 }
608
609 template<class CP>
611 const std::vector<Real>& invUncondProbs,
612 Real saddle,
613 const std::vector<Real>& mktFactor) const
614 {
615 const Size nNames = remainingNotionals_.size();
616 Real sum = 0.;
617
618 for(Size iName=0; iName < nNames; iName++) {
619 Probability pBuffer =
620 copula_->conditionalDefaultProbabilityInvP(
621 invUncondProbs[iName], iName, mktFactor);
622 // loss in fractional units
623 Real lossInDef = remainingNotionals_[iName] *
624 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
625 iName, mktFactor)) / remainingNotional_;
626 Real midFactor = pBuffer * std::exp(lossInDef * saddle);
627 Real denominator = 1.-pBuffer + midFactor;
628 sum += lossInDef * lossInDef * midFactor / denominator -
629 std::pow(lossInDef * midFactor / denominator , 2.);
630 }
631 return sum;
632 }
633
634 template<class CP>
636 const std::vector<Real>& invUncondProbs,
637 Real saddle,
638 const std::vector<Real>& mktFactor) const
639 {
640 const Size nNames = remainingNotionals_.size();
641 Real sum = 0.;
642
643 for(Size iName=0; iName < nNames; iName++) {
644 Probability pBuffer =
645 copula_->conditionalDefaultProbabilityInvP(
646 invUncondProbs[iName], iName, mktFactor);
647 Real lossInDef = remainingNotionals_[iName] *
648 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
649 iName, mktFactor)) / remainingNotional_;
650
651 const Real midFactor = pBuffer * std::exp(lossInDef * saddle);
652 const Real denominator = 1.-pBuffer + midFactor;
653
654 const Real& suma0 = denominator;
655 const Real suma1 = lossInDef * midFactor;
656 const Real suma2 = lossInDef * suma1;
657 const Real suma3 = lossInDef * suma2;
658
659 sum += (suma3 + (2.*std::pow(suma1, 3.)/suma0 -
660 3.*suma1*suma2)/suma0)/suma0;
661 }
662 return sum;
663 }
664
665 template<class CP>
667 const std::vector<Real>& invUncondProbs,
668 Real saddle,
669 const std::vector<Real>& mktFactor) const
670 {
671 const Size nNames = remainingNotionals_.size();
672 Real sum = 0.;
673
674 for(Size iName=0; iName < nNames; iName++) {
675 Probability pBuffer =
676 copula_->conditionalDefaultProbabilityInvP(
677 invUncondProbs[iName], iName, mktFactor);
678 Real lossInDef = remainingNotionals_[iName] *
679 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
680 iName, mktFactor)) / remainingNotional_;
681
682 Real midFactor = pBuffer * std::exp(lossInDef * saddle);
683 Real denominator = 1.-pBuffer + midFactor;
684
685 const Real& suma0 = denominator;
686 const Real suma1 = lossInDef * midFactor;
687 const Real suma2 = lossInDef * suma1;
688 const Real suma3 = lossInDef * suma2;
689 const Real suma4 = lossInDef * suma3;
690
691 sum += (suma4 + (-4.*suma1*suma3 - 3.*suma2*suma2 +
692 (12.*suma1*suma1*suma2 -
693 6.*std::pow(suma1,4.)/suma0)/suma0)/suma0)/suma0;
694 }
695 return sum;
696 }
697
698 template<class CP>
699 ext::tuple<Real, Real, Real, Real> SaddlePointLossModel<CP>::CumGen0234DerivCond(
700 const std::vector<Real>& invUncondProbs,
701 Real saddle,
702 const std::vector<Real>& mktFactor) const
703 {
704 const Size nNames = remainingNotionals_.size();
705 Real deriv0 = 0.,
706 //deriv1 = 0.,
707 deriv2 = 0.,
708 deriv3 = 0.,
709 deriv4 = 0.;
710 for(Size iName=0; iName < nNames; iName++) {
711 Probability pBuffer =
712 copula_->conditionalDefaultProbabilityInvP(
713 invUncondProbs[iName], iName, mktFactor);
714 Real lossInDef = remainingNotionals_[iName] *
715 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
716 iName, mktFactor)) / remainingNotional_;
717
718 Real midFactor = pBuffer * std::exp(lossInDef * saddle);
719 Real denominator = 1.-pBuffer + midFactor;
720
721 const Real& suma0 = denominator;
722 const Real suma1 = lossInDef * midFactor;
723 const Real suma2 = lossInDef * suma1;
724 const Real suma3 = lossInDef * suma2;
725 const Real suma4 = lossInDef * suma3;
726
727 // To do: optimize these:
728 deriv0 += std::log(suma0);
729 //deriv1 += suma1 / suma0;
730 deriv2 += suma2 / suma0 - std::pow(suma1 / suma0 , 2.);
731 deriv3 += (suma3 + (2.*std::pow(suma1, 3.)/suma0 -
732 3.*suma1*suma2)/suma0)/suma0;
733 deriv4 += (suma4 + (-4.*suma1*suma3 - 3.*suma2*suma2 +
734 (12.*suma1*suma1*suma2 -
735 6.*std::pow(suma1,4.)/suma0)/suma0)/suma0)/suma0;
736 }
737 return {deriv0, deriv2, deriv3, deriv4};
738 }
739
740 template<class CP>
742 const std::vector<Real>& invUncondProbs,
743 Real saddle,
744 const std::vector<Real>& mktFactor) const
745 {
746 const Size nNames = remainingNotionals_.size();
747 Real deriv0 = 0.,
748 //deriv1 = 0.,
749 deriv2 = 0.;
750 for(Size iName=0; iName < nNames; iName++) {
751 Probability pBuffer =
752 copula_->conditionalDefaultProbabilityInvP(
753 invUncondProbs[iName], iName, mktFactor);
754 Real lossInDef = remainingNotionals_[iName] *
755 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
756 iName, mktFactor)) / remainingNotional_;
757
758 Real midFactor = pBuffer * std::exp(lossInDef * saddle);
759 Real denominator = 1.-pBuffer + midFactor;
760
761 const Real& suma0 = denominator;
762 const Real suma1 = lossInDef * midFactor;
763 const Real suma2 = lossInDef * suma1;
764
765 // To do: optimize these:
766 deriv0 += std::log(suma0);
767 //deriv1 += suma1 / suma0;
768 deriv2 += suma2 / suma0 - std::pow(suma1 / suma0 , 2.);
769 }
770 return {deriv0, deriv2};
771 }
772
773 // ----- Saddle point search ----------------------------------------------
774
775 template<class CP>
777 const std::vector<Real>& invUncondPs,
778 Real lossLevel, // in total portfolio loss fractional unit
779 const std::vector<Real>& mktFactor,
780 Real accuracy,
781 Natural maxEvaluations
782 ) const
783 {
784 // \to do:
785 // REQUIRE that loss level is below the max loss attainable in
786 // the portfolio, otherwise theres no solution...
787 SaddleObjectiveFunction f(*this, lossLevel, invUncondPs, mktFactor);
788
789 Size nNames = remainingNotionals_.size();
790 std::vector<Real> lgds;
791 for(Size iName=0; iName<nNames; iName++)
792 lgds.push_back(remainingNotionals_[iName] *
793 (1.-copula_->conditionalRecoveryInvP(invUncondPs[iName], iName,
794 mktFactor)) );
795
796 // computed limits:
797 // position of the name with the largest relative exposure loss (i.e.:
798 // largest: N_i LGD_i / N_{total})
799 Size iNamMax = std::distance(lgds.begin(),
800 std::max_element(lgds.begin(), lgds.end()) );
801 // gap to be considered zero at the negative side of the logistic
802 // inversion:
803 static const Real deltaMin = 1.e-5;
804 //
805 Probability pMaxName = copula_->conditionalDefaultProbabilityInvP(
806 invUncondPs[iNamMax], iNamMax, mktFactor);
807 // aproximates the saddle pt corresponding to this minimum; finds
808 // it by using only the smallest logistic term and thus this is
809 // smaller than the true value:
810 Real saddleMin = 1./(lgds[iNamMax]/remainingNotional_) *
811 std::log(deltaMin*(1.-pMaxName)/
812 (pMaxName*lgds[iNamMax]/remainingNotional_-pMaxName*deltaMin));
813 // and the associated minimum loss is approximately: (this is thence
814 // the minimum loss we can resolve/invert)
815 Real minLoss =
816 CumGen1stDerivativeCond(invUncondPs, saddleMin, mktFactor);
817
818 // If we are below the loss resolution it returns approximating
819 // by the minimum/maximum attainable point. Typically the functionals
820 // to integrate will have a low dependency on this point.
821 if(lossLevel < minLoss) return saddleMin;
822
823 Real saddleMax = 1./(lgds[iNamMax]/remainingNotional_) *
824 std::log((lgds[iNamMax]/remainingNotional_
825 -deltaMin)*(1.-pMaxName)/(pMaxName*deltaMin));
826 Real maxLoss =
827 CumGen1stDerivativeCond(invUncondPs, saddleMax, mktFactor);
828 if(lossLevel > maxLoss) return saddleMax;
829
830 Brent solverBrent;
831 Real guess = (saddleMin+saddleMax)/2.;
832 /*
833 (lossLevel -
834 CumGen1stDerivativeCond(invUncondPs, lossLevel, mktFactor))
835 /CumGen2ndDerivativeCond(invUncondPs, lossLevel, mktFactor);
836 if(guess > saddleMax) guess = (saddleMin+saddleMax)/2.;
837 */
838 solverBrent.setMaxEvaluations(maxEvaluations);
839 return solverBrent.solve(f, accuracy, guess, saddleMin, saddleMax);
840 }
841
842
843 // ----- Statistics -------------------------------------------------------
844
845
846 template<class CP>
848 Probability percentile) const
849 {
850 //this test should be in the calling basket...?
851 QL_REQUIRE(percentile >=0. && percentile <=1.,
852 "Incorrect percentile value.");
853
854 // still this does not tackle the situation where we have cumulated
855 // losses from previous defaults:
856 if(d <= Settings::instance().evaluationDate()) return 0.;
857
858 // Trivial cases when the percentile is outside the prob range
859 // associated to the tranche limits:
860 if(percentile <= 1.-probOverLoss(d, 0.)) return 0.;
861 if(percentile >= 1.-probOverLoss(d, 1.))
862 return basket_->remainingTrancheNotional();
863
864 SaddlePercObjFunction f(*this, percentile, d);
865 Brent solver;
866 solver.setMaxEvaluations(100);
867 Real minVal = QL_EPSILON;
868
869 Real maxVal = 1.-QL_EPSILON;
870 Real guess = 0.5;
871
872 Real solut = solver.solve(f, 1.e-4, guess, minVal, maxVal);
873 return basket_->remainingTrancheNotional() * solut;
874 }
875
876 template<class CP>
878 const std::vector<Real>& invUncondPs,
879 Real trancheLossFract,
880 const std::vector<Real>& mktFactor) const {
881 Real portfFract = attachRatio_ + trancheLossFract *
882 (detachRatio_-attachRatio_);// these are remaining ratios
883
884 // for non-equity losses add the probability jump at zero tranche
885 // losses (since this method returns prob of losing more or
886 // equal to)
887 ////////////////--- if(trancheLossFract <= QL_EPSILON) return 1.;
888 return
889 probOverLossPortfCond(invUncondPs,
890 //below; should substract realized loses. Use remaining amounts??
891 portfFract * basket_->basketNotional(),
892 mktFactor);
893 }
894
895 template<class CP>
896 std::map<Real, Probability> SaddlePointLossModel<CP>::lossDistribution(const Date& d) const {
897 std::map<Real, Probability> distrib;
898 static constexpr double numPts = 500.;
899 for(Real lossFraction=1./numPts; lossFraction<0.45;
900 lossFraction+= 1./numPts)
901 distrib.insert(std::make_pair<Real, Probability>(
902 lossFraction * remainingNotional_ ,
903 1.-probOverPortfLoss(d, lossFraction* remainingNotional_ )));
904 return distrib;
905 }
906
907 /* NOTICE THIS IS ON THE TOTAL PORTFOLIO ---- UNTRANCHED..............
908 Probability of having losses in the portfolio due to default
909 events equal or larger than a given absolute loss value on a
910 given date conditional to the latent model factor.
911 The integral expression on the expansion is the first order
912 integration as presented in several references, see for instance;
913 equation 8 in R.Martin, K.Thompson, and C. Browne 's
914 'Taking to the Saddle', Risk Magazine, June 2001, page 91
915
916 The loss is passed in absolute value.
917 */
918 template<class CP>
920 const std::vector<Real>& invUncondProbs,
921 Real loss,
922 const std::vector<Real>& mktFactor) const
923 {
924 /* This is taking in the unconditional probabilites non inverted. See if
925 the callers can be written taking the inversion already; if they are
926 performing it thats a perf hit. At least this can be seen to be true
927 for the recovery call (but rand rr are not intended to be used yet)
928 */
929 // return probOverLossPortfCond1stOrder(d, loss, mktFactor);
930 if (loss <= QL_EPSILON) return 1.;
931
932 Real relativeLoss = loss / remainingNotional_;
933 if (relativeLoss >= 1.-QL_EPSILON) return 0.;
934
935 const Size nNames = remainingNotionals_.size();
936
937 Real averageRecovery_ = 0.;
938 for(Size iName=0; iName < nNames; iName++)
939 averageRecovery_ += copula_->conditionalRecoveryInvP(
940 invUncondProbs[iName], iName, mktFactor);
941 averageRecovery_ = averageRecovery_ / nNames;
942
943 Real maxAttLossFract = 1.-averageRecovery_;
944 if(relativeLoss > maxAttLossFract) return 0.;
945
946 Real saddlePt = findSaddle(invUncondProbs,
947 relativeLoss, mktFactor);
948
949 ext::tuple<Real, Real, Real, Real> cumulants =
950 CumGen0234DerivCond(invUncondProbs,
951 saddlePt, mktFactor);
952 Real baseVal = ext::get<0>(cumulants);
953 Real secondVal = ext::get<1>(cumulants);
954 Real K3Saddle = ext::get<2>(cumulants);
955 Real K4Saddle = ext::get<3>(cumulants);
956
957 Real saddleTo2 = saddlePt * saddlePt;
958 Real saddleTo3 = saddleTo2 * saddlePt;
959 Real saddleTo4 = saddleTo3 * saddlePt;
960 Real saddleTo6 = saddleTo4 * saddleTo2;
961 Real K3SaddleTo2 = K3Saddle*K3Saddle;
962
963 if(saddlePt > 0.) { // <-> (loss > condEL)
964 Real exponent = baseVal - relativeLoss * saddlePt +
965 .5 * saddleTo2 * secondVal;
966 if( std::abs(exponent) > 700.) return 0.;
967 return
968 std::exp(exponent)
969 * CumulativeNormalDistribution()(-std::abs(saddlePt)*
970 std::sqrt(/*saddleTo2 **/secondVal))
971
972 // high order corrections:
973 * (1. - saddleTo3*K3Saddle/6. + saddleTo4*K4Saddle/24. +
974 saddleTo6*K3SaddleTo2/72.)
975 /*
976 // FIX ME: this term introduces at times numerical
977 // instabilty (shows up in percentile computation)
978 + (3.*secondVal*(1.-secondVal*saddleTo2)*
979 (saddlePt*K4Saddle-4.*K3Saddle)
980 - saddlePt*K3SaddleTo2*(3.-saddleTo2*secondVal +
981 saddleTo4*secondVal*secondVal))
982 / (72.*M_SQRTPI*M_SQRT_2*std::pow(secondVal, 5./2.) )
983 */
984 ;
985 }else if(saddlePt==0.){// <-> (loss == condEL)
986 return .5;
987 }else {// <->(loss < condEL)
988 Real exponent = baseVal - relativeLoss * saddlePt +
989 .5 * saddleTo2 * secondVal;
990 if( std::abs(exponent) > 700.) return 0.;
991 return
992 1.-
993 std::exp(exponent)
994 * CumulativeNormalDistribution()(-std::abs(saddlePt)
995 * std::sqrt(/*saddleTo2 **/secondVal))// static call?
996
997 // high order corrections:
998 * (1. - saddleTo3*K3Saddle/6. + saddleTo4*K4Saddle/24. +
999 saddleTo6*K3SaddleTo2/72.)
1000 /*
1001 + (3.*secondVal*(1.-secondVal*saddleTo2)*
1002 (saddlePt*K4Saddle-4.*K3Saddle)
1003 - saddlePt*K3SaddleTo2*(3.-saddleTo2*secondVal +
1004 saddleTo4*secondVal*secondVal))
1005 / (72.*M_SQRTPI*M_SQRT_2*std::pow(secondVal, 5./2.) )
1006 */
1007 ;
1008 }
1009 }
1010
1011 template<class CP>
1012 // cheaper; less terms retained; yet the cost lies in the saddle point calc
1014 const std::vector<Real>& invUncondPs,
1015 Real loss,
1016 const std::vector<Real>& mktFactor) const
1017 {
1018 if (loss <= QL_EPSILON) return 1.;
1019 const Size nNames = remainingNotionals_.size();
1020
1021 Real relativeLoss = loss / remainingNotional_;
1022 if(relativeLoss >= 1.-QL_EPSILON) return 0.;
1023
1024 // only true for constant recovery models......?
1025 Real averageRecovery_ = 0.;
1026 for(Size iName=0; iName < nNames; iName++)
1027 averageRecovery_ +=
1028 copula_->conditionalRecoveryInvP(invUncondPs[iName], iName,
1029 mktFactor);
1030 averageRecovery_ = averageRecovery_ / nNames;
1031
1032 Real maxAttLossFract = 1.-averageRecovery_;
1033 if(relativeLoss > maxAttLossFract) return 0.;
1034
1035 Real saddlePt = findSaddle(invUncondPs,
1036 relativeLoss, mktFactor);
1037
1038 ext::tuple<Real, Real> cumulants =
1039 CumGen02DerivCond(invUncondPs,
1040 saddlePt, mktFactor);
1041 Real baseVal = ext::get<0>(cumulants);
1042 Real secondVal = ext::get<1>(cumulants);
1043
1044 Real saddleTo2 = saddlePt * saddlePt;
1045
1046 if(saddlePt > 0.) { // <-> (loss > condEL)
1047 Real exponent = baseVal - relativeLoss * saddlePt +
1048 .5 * saddleTo2 * secondVal;
1049 if( std::abs(exponent) > 700.) return 0.;
1050 return
1051 // dangerous exponential; fix me
1052 std::exp(exponent)
1053 /* std::exp(baseVal - relativeLoss * saddlePt
1054 + .5 * saddleTo2 * secondVal)*/
1055 * CumulativeNormalDistribution()(-std::abs(saddlePt)*
1056 std::sqrt(/*saddleTo2 **/secondVal));
1057 }else if(saddlePt==0.){// <-> (loss == condEL)
1058 return .5;
1059 }else {// <->(loss < condEL)
1060 Real exponent = baseVal - relativeLoss * saddlePt +
1061 .5 * saddleTo2 * secondVal;
1062 if( std::abs(exponent) > 700.) return 0.;
1063
1064 return
1065 1.-
1066 /* std::exp(baseVal - relativeLoss * saddlePt
1067 + .5 * saddleTo2 * secondVal)*/
1068 std::exp(exponent)
1069 * CumulativeNormalDistribution()(-std::abs(saddlePt)*
1070 std::sqrt(/*saddleTo2 **/secondVal));
1071 }
1072 }
1073
1074
1075 /*! NOTICE THIS IS ON THE TOTAL PORTFOLIO ---- UNTRANCHED
1076 Probability density of having losses in the portfolio due to default
1077 events equal to a given value on a given date conditional to the w
1078 latent model factor.
1079 Based on the integrals of the expected shortfall. See......refernce.
1080 */
1081 template<class CP>
1083 const std::vector<Real>& invUncondPs,
1084 Real loss,
1085 const std::vector<Real>& mktFactor) const
1086 {
1087 if (loss <= QL_EPSILON) return 0.;
1088
1089 Real relativeLoss = loss / remainingNotional_;
1090 Real saddlePt = findSaddle(invUncondPs,
1091 relativeLoss, mktFactor);
1092
1093 ext::tuple<Real, Real, Real, Real> cumulants =
1094 CumGen0234DerivCond(invUncondPs,
1095 saddlePt, mktFactor);
1096 /// access them directly rather than through this copy
1097 Real K0Saddle = ext::get<0>(cumulants);
1098 Real K2Saddle = ext::get<1>(cumulants);
1099 Real K3Saddle = ext::get<2>(cumulants);
1100 Real K4Saddle = ext::get<3>(cumulants);
1101 /* see, for instance R.Martin "he saddle point method and portfolio
1102 optionalities." in Risk December 2006 p.93 */
1103 //\todo the exponentials below are dangerous and agressive, tame them.
1104 return
1105 (
1106 1.
1107 + K4Saddle
1108 /(8.*std::pow(K2Saddle, 2.))
1109 - 5.*std::pow(K3Saddle,2.)
1110 /(24.*std::pow(K2Saddle, 3.))
1111 ) * std::exp(K0Saddle - saddlePt * relativeLoss)
1112 / (std::sqrt(2. * M_PI * K2Saddle));
1113 }
1114
1115 /* NOTICE THIS IS ON THE TOTAL PORTFOLIO ---- UNTRANCHED..
1116 Sensitivities of the individual names to a given portfolio loss value
1117 due to defaults. It returns ratios to the total structure notional,
1118 which aggregated add up to the requested loss value.
1119
1120 see equation 8 in 'VAR: who contributes and how much?' by R.Martin,
1121 K.Thompson, and C. Browne in Risk Magazine, August 2001
1122
1123 The passed loss is the loss amount level at which we want to
1124 request the sensitivity. Equivalent to a percentile.
1125 */
1126 template<class CP>
1128 const std::vector<Real>& invUncondProbs,
1129 Real loss,
1130 std::vector<Real> mktFactor) const
1131 {
1132 const Size nNames = remainingNotionals_.size();
1133 std::vector<Real> condContrib(nNames, 0.);
1134 if (loss <= QL_EPSILON) return condContrib;
1135
1136 Real saddlePt = findSaddle(invUncondProbs, loss / remainingNotional_,
1137 mktFactor);
1138
1139 for(Size iName=0; iName<nNames; iName++) {
1140 Probability pBuffer =
1141 copula_->conditionalDefaultProbabilityInvP(
1142 invUncondProbs[iName], iName, mktFactor);
1143 Real lossInDef = remainingNotionals_[iName] *
1144 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
1145 iName, mktFactor));
1146 Real midFactor = pBuffer *
1147 std::exp(lossInDef * saddlePt/ remainingNotional_);
1148 Real denominator = 1.-pBuffer + midFactor;
1149
1150 condContrib[iName] = lossInDef * midFactor / denominator;
1151 }
1152 return condContrib;
1153 }
1154
1155 template<class CP>
1157 const std::vector<Real>& invUncondProbs,
1158 const std::vector<Real>& mktFactor) const
1159 {
1160 const Size nNames = remainingNotionals_.size();
1161 Real eloss = 0.;
1162 /// USE STL.....-------------------
1163 for(Size iName=0; iName < nNames; iName++) {
1164 Probability pBuffer = copula_->conditionalDefaultProbabilityInvP(
1165 invUncondProbs[iName], iName, mktFactor);
1166 eloss += pBuffer * remainingNotionals_[iName] *
1167 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
1168 iName, mktFactor));
1169 }
1170 return eloss;
1171 }
1172
1173 template<class CP>
1175 const std::vector<Real>& invUncondProbs,
1176 const std::vector<Real>& mktFactor) const
1177 {
1178 const Size nNames = remainingNotionals_.size();
1179 Real eloss = 0.;
1180 /// USE STL.....-------------------
1181 for(Size iName=0; iName < nNames; iName++) {
1182 Probability pBuffer = copula_->conditionalDefaultProbabilityInvP(
1183 invUncondProbs[iName], iName, mktFactor);
1184 eloss +=
1185 pBuffer * remainingNotionals_[iName] *
1186 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
1187 iName, mktFactor));
1188 }
1189 return std::min(
1190 std::max(eloss - attachRatio_ * remainingNotional_, 0.),
1191 (detachRatio_ - attachRatio_) * remainingNotional_);
1192 }
1193
1194 template<class CP>
1196 const std::vector<Real>& invUncondProbs,
1197 Real lossPerc, const std::vector<Real>& mktFactor) const
1198 {
1199 const Size nNames = remainingNotionals_.size();
1200 std::vector<Real> lgds;
1201 for(Size iName=0; iName<nNames; iName++)
1202 lgds.push_back(remainingNotionals_[iName] *
1203 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
1204 iName, mktFactor)));
1205 std::vector<Real> vola(nNames, 0.), mu(nNames, 0.);
1206 Real volaTot = 0., muTot = 0.;
1207 for(Size iName=0; iName < nNames; iName++) {
1208 Probability pBuffer = copula_->conditionalDefaultProbabilityInvP(
1209 invUncondProbs[iName], iName, mktFactor);
1210 mu[iName] = lgds[iName] * pBuffer / remainingNotionals_[iName];
1211 muTot += lgds[iName] * pBuffer;
1212 vola[iName] = lgds[iName] * lgds[iName] * pBuffer * (1.-pBuffer)
1213 / remainingNotionals_[iName];
1214 volaTot += lgds[iName] * lgds[iName] * pBuffer * (1.-pBuffer) ;
1215 }
1216 for (Size iName=0; iName < nNames; iName++)
1217 vola[iName] = vola[iName] / volaTot;
1218
1219 std::vector<Real> esfPartition(nNames, 0.);
1220 for(Size iName=0; iName < nNames; iName++) {
1221 Real uEdisp = (lossPerc-muTot)/std::sqrt(volaTot);
1222 esfPartition[iName] = mu[iName]
1223 * CumulativeNormalDistribution()(uEdisp) // static call?
1224 + vola[iName] * NormalDistribution()(uEdisp);
1225 }
1226 return esfPartition;
1227 }
1228
1229 template<class CP>
1231 const std::vector<Real>& invUncondProbs,
1232 Real lossPerc, // value
1233 Probability percentile,
1234 const std::vector<Real>& mktFactor) const
1235 {
1236 /* TO DO: this is too crude, a general expression valid for all
1237 situations is possible (with no extra cost as long as the loss limits
1238 are checked).
1239 */
1240 //tranche correction term:
1241 Real correctionTerm = 0.;
1242 Real probLOver = probOverLossPortfCond(invUncondProbs,
1243 basket_->detachmentAmount(), mktFactor);
1244 if(basket_->attachmentAmount() > QL_EPSILON) {
1245 if(lossPerc < basket_->attachmentAmount()) {
1246 correctionTerm = ( (basket_->detachmentAmount()
1247 - 2.*basket_->attachmentAmount())*
1248 probOverLossPortfCond(invUncondProbs, lossPerc,
1249 mktFactor)
1250 + basket_->attachmentAmount() * probLOver )/(1.-percentile);
1251 }else{
1252 correctionTerm = ( (percentile-1)*basket_->attachmentAmount()
1253 + basket_->detachmentAmount() * probLOver
1254 )/(1.-percentile);
1255 }
1256 }
1257
1258 return expectedShortfallFullPortfolioCond(invUncondProbs,
1259 std::max(lossPerc, basket_->attachmentAmount()), mktFactor)
1260 + expectedShortfallFullPortfolioCond(invUncondProbs,
1261 basket_->detachmentAmount(), mktFactor)
1262 - correctionTerm;
1263 }
1264
1265 template<class CP>
1267 const std::vector<Real>& invUncondProbs,
1268 Real lossPerc, // value
1269 const std::vector<Real>& mktFactor) const
1270 {
1271 /* This version is based on: Martin 2006 paper and on the expression
1272 in 'SaddlePoint approximation of expected shortfall for transformed
1273 means' S.A. Broda and M.S.Paolella , Amsterdam School of Economics
1274 discussion paper, available online.
1275 */
1276 Real lossPercRatio = lossPerc /remainingNotional_;
1277 Real elCond = 0.;
1278 const Size nNames = remainingNotionals_.size();
1279
1280 /// use stl algorthms
1281 for(Size iName=0; iName < nNames; iName++) {
1282 Probability pBuffer = copula_->conditionalDefaultProbabilityInvP(
1283 invUncondProbs[iName], iName, mktFactor);
1284 elCond += pBuffer * remainingNotionals_[iName] *
1285 (1.-copula_->conditionalRecoveryInvP(invUncondProbs[iName],
1286 iName, mktFactor));
1287 }
1288 Real saddlePt = findSaddle(invUncondProbs, lossPercRatio, mktFactor);
1289
1290 // Martin 2006:
1291 return
1292 elCond * probOverLossPortfCond(invUncondProbs, lossPerc, mktFactor)
1293 + (lossPerc - elCond) * probDensityCond(invUncondProbs, lossPerc,
1294 mktFactor) /saddlePt;
1295
1296 // calling the EL tranche
1297 // return elCond - expectedEquityLossCond(d, lossPercRatio, mktFactor);
1298
1299 /*
1300 // Broda and Paolella:
1301 Real elCondRatio = elCond / remainingNotional_;
1302
1303 ext::tuple<Real, Real, Real, Real> cumulants =
1304 CumGen0234DerivCond(uncondProbs,
1305 saddlePt, mktFactor);
1306 Real K0Saddle = ext::get<0>(cumulants);///USE THEM DIRECTLY
1307 Real K2Saddle = ext::get<1>(cumulants);
1308
1309 Real wq = std::sqrt(2. * saddlePt * lossPercRatio - 2. * K0Saddle);
1310 //std::sqrt(-2. * saddlePt * lossPerc + 2. * K0Saddle);????
1311 Real factor = 1.;
1312 if(saddlePt<0.) {
1313 wq = -wq;
1314 factor = -1.;
1315 }
1316
1317 Real numNames = static_cast<Real>(nNames);
1318
1319 Real term1 = CumulativeNormalDistribution()(wq)// * std::sqrt(numNames)
1320 * elCond ;
1321 Real term2 = .5 * M_2_SQRTPI * M_SQRT1_2 * (1./std::sqrt(numNames))
1322 * exp(-wq*wq * numNames/2.)*(elCond/wq -
1323 lossPerc/(saddlePt * std::sqrt(K2Saddle)));
1324 return term1 + term2;
1325 */
1326 }
1327
1328 template<class CP>
1330 Probability percProb) const
1331 {
1332 // assuming I have the tranched one.
1333 Real lossPerc = percentile(d, percProb);
1334
1335 // check the trivial case when the loss is over the detachment limit
1336 // to avoid computation:
1337 Real trancheAmount = basket_->trancheNotional() *
1338 (detachRatio_-attachRatio_);
1339 //assumed the amount includes the realized loses
1340 if(lossPerc >= trancheAmount) return trancheAmount;
1341 //SHOULD CHECK NOW THE OPPOSITE LIMIT ("zero" losses)....
1342 std::vector<Real> invUncondProbs =
1343 basket_->remainingProbabilities(d);
1344 for(Size i=0; i<invUncondProbs.size(); i++)
1345 invUncondProbs[i] =
1346 copula_->inverseCumulativeY(invUncondProbs[i], i);
1347
1348 // Integrate with the tranche or the portfolio according to the limits.
1349 return copula_->integratedExpectedValue(
1350 [&](const std::vector<Real>& v1) {
1351 return expectedShortfallFullPortfolioCond(invUncondProbs, lossPerc, v1);
1352 }) / (1.-percProb);
1353
1354 /* test:?
1355 return std::inner_product(integrESFPartition.begin(),
1356 integrESFPartition.end(), remainingNotionals_.begin(), 0.);
1357 */
1358
1359 }
1360
1361
1362
1363}
1364
1365#endif
basket of issuers and related notionals
Brent 1-D solver.
Brent 1-D solver
Definition: brent.hpp:37
Cumulative normal distribution function.
Concrete date class.
Definition: date.hpp:125
RelinkableHandle< Basket > basket_
Normal distribution function.
SaddleObjectiveFunction(const SaddlePointLossModel &me, const Real target, const std::vector< Real > &invUncondProbs, const std::vector< Real > &mktFactor)
The passed target is in fractional loss units.
SaddlePercObjFunction(const SaddlePointLossModel &me, const Real target, const Date &date)
Saddle point portfolio credit default loss model.
SaddlePointLossModel(const ext::shared_ptr< ConstantLossLatentmodel< CP > > &m)
Real CumulantGenerating(const Date &date, Real s) const
Probability probDensityCond(const std::vector< Real > &invUncondProbs, Real loss, const std::vector< Real > &mktFactor) const
Probability probOverLossCond(const std::vector< Real > &invUncondProbs, Real trancheLossFract, const std::vector< Real > &mktFactor) const
Real CumGen1stDerivative(const Date &date, Real s) const
Real CumGen3rdDerivative(const Date &date, Real s) const
Real conditionalExpectedLoss(const std::vector< Real > &invUncondProbs, const std::vector< Real > &mktFactor) const
ext::tuple< Real, Real > CumGen02DerivCond(const std::vector< Real > &invUncondProbs, Real saddle, const std::vector< Real > &mktFactor) const
Probability probOverPortfLoss(const Date &d, Real loss) const
Real expectedShortfallTrancheCond(const std::vector< Real > &invUncondProbs, Real lossPerc, Probability percentile, const std::vector< Real > &mktFactor) const
Real expectedShortfallFullPortfolioCond(const std::vector< Real > &invUncondProbs, Real lossPerc, const std::vector< Real > &mktFactor) const
Real CumulantGeneratingCond(const std::vector< Real > &invUncondProbs, Real lossFraction, const std::vector< Real > &mktFactor) const
Real CumGen1stDerivativeCond(const std::vector< Real > &invUncondProbs, Real saddle, const std::vector< Real > &mktFactor) const
Real CumGen3rdDerivativeCond(const std::vector< Real > &invUncondProbs, Real saddle, const std::vector< Real > &mktFactor) const
std::vector< Real > splitVaRLevel(const Date &date, Real loss) const override
Probability probOverLossPortfCond1stOrder(const std::vector< Real > &invUncondProbs, Real loss, const std::vector< Real > &mktFactor) const
Real CumGen4thDerivative(const Date &date, Real s) const
std::vector< Real > splitLossCond(const std::vector< Real > &invUncondProbs, Real loss, std::vector< Real > mktFactor) const
Real CumGen2ndDerivative(const Date &date, Real s) const
Probability probOverLossPortfCond(const std::vector< Real > &invUncondProbs, Real loss, const std::vector< Real > &mktFactor) const
void resetModel() override
Concrete models do now any updates/inits they need on basket reset.
Real percentile(const Date &d, Probability percentile) const override
Real expectedShortfall(const Date &d, Probability percentile) const override
Expected shortfall given a default loss percentile.
Real CumGen2ndDerivativeCond(const std::vector< Real > &invUncondProbs, Real saddle, const std::vector< Real > &mktFactor) const
const ext::shared_ptr< ConstantLossLatentmodel< CP > > copula_
Real expectedTrancheLoss(const Date &d) const override
ext::tuple< Real, Real, Real, Real > CumGen0234DerivCond(const std::vector< Real > &invUncondProbs, Real saddle, const std::vector< Real > &mktFactor) const
Probability probDensity(const Date &d, Real loss) const
Real findSaddle(const std::vector< Real > &invUncondProbs, Real lossLevel, const std::vector< Real > &mktFactor, Real accuracy=1.0e-3, Natural maxEvaluations=50) const
Real conditionalExpectedTrancheLoss(const std::vector< Real > &invUncondProbs, const std::vector< Real > &mktFactor) const
Probability probOverLoss(const Date &d, Real trancheLossFract) const override
std::map< Real, Probability > lossDistribution(const Date &d) const override
Full loss distribution.
Real CumGen4thDerivativeCond(const std::vector< Real > &invUncondProbs, Real saddle, const std::vector< Real > &mktFactor) const
std::vector< Real > expectedShortfallSplitCond(const std::vector< Real > &invUncondProbs, Real lossPerc, const std::vector< Real > &mktFactor) const
static Settings & instance()
access to the unique instance
Definition: singleton.hpp:104
void setMaxEvaluations(Size evaluations)
Definition: solver1d.hpp:238
Real solve(const F &f, Real accuracy, Real guess, Real step) const
Definition: solver1d.hpp:84
#define QL_REQUIRE(condition, message)
throw an error if the given pre-condition is not verified
Definition: errors.hpp:117
Date d
#define QL_EPSILON
Definition: qldefines.hpp:178
QL_REAL Real
real number
Definition: types.hpp:50
unsigned QL_INTEGER Natural
positive integer
Definition: types.hpp:43
Real Probability
probability
Definition: types.hpp:82
std::size_t Size
size of a container
Definition: types.hpp:58
#define M_PI
Definition: any.hpp:35
Newton 1-D solver.
Maps tuple to either the boost or std implementation.