QuantLib: a free/open-source library for quantitative finance
fully annotated source code - version 1.34
Loading...
Searching...
No Matches
creditriskplus.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) 2013 Peter Caspers
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
21#include <map>
22#include <utility>
23
24using std::sqrt;
25
26namespace QuantLib {
27
28 CreditRiskPlus::CreditRiskPlus(std::vector<Real> exposure,
29 std::vector<Real> defaultProbability,
30 std::vector<Size> sector,
31 std::vector<Real> relativeDefaultVariance,
32 Matrix correlation,
33 const Real unit)
34 : exposure_(std::move(exposure)), pd_(std::move(defaultProbability)),
35 sector_(std::move(sector)), relativeDefaultVariance_(std::move(relativeDefaultVariance)),
36 correlation_(std::move(correlation)), unit_(unit) {
37
38 m_ = exposure_.size();
39
40 QL_REQUIRE(m_ > 0, "no exposures given");
41 QL_REQUIRE(m_ == pd_.size(), "number of exposures ("
42 << m_
43 << ") must be equal to number of pds ("
44 << pd_.size() << ")");
45 QL_REQUIRE(m_ == sector_.size(),
46 "number of exposures ("
47 << m_
48 << ") must be equal to number of exposure sectors ("
49 << sector_.size() << ")");
50
53 "correlation matrix (" << n_ << "," << correlation_.columns()
54 << ") must be a square matrix");
55
57 "number of relative default variances ("
58 << relativeDefaultVariance_.size() << ")"
59 << " must be equal to number of sectors (" << n_ << ")");
60
61 exposureSum_ = 0.0;
62 el_ = 0.0;
63 el2_ = 0.0;
64 for (Size i = 0; i < m_; ++i) {
65 QL_REQUIRE(exposure_[i] >= 0.0, "exposure #"
66 << i << " is negative ("
67 << exposure_[i] << ")");
68 QL_REQUIRE(pd_[i] > 0.0, "pd #" << i << " is negative (" << pd_[i]
69 << ")");
70 QL_REQUIRE(sector_[i] < n_, "sector #" << i << " (" << sector_[i]
71 << ") is out of range 0..."
72 << (n_ - 1));
74 el_ += pd_[i] * exposure_[i];
75 el2_ += pd_[i] * exposure_[i]*exposure_[i];
76 }
77
78 QL_REQUIRE(unit_ > 0.0, "loss unit (" << unit_ << ") must be positive");
79
80 compute();
81 }
82
84
85 Size i = 0;
86 Real sum = loss_[0];
87 while(i < upperIndex_-1 && sum < p) {
88 ++i;
89 sum += loss_[i];
90 }
91
92 if(loss_[0] >= p)
93 return 0.0;
94
95 Real p1 = sum - loss_[i];
96 Real p2 = sum >= p ? sum : 1.0;
97 Real l1 = (i - 1) * unit_;
98 Real l2 = i * unit_;
99
100 return l1 + (p - p1) / (p2 - p1) * (l2 - l1);
101 }
102
104
105 std::vector<Real> sectorPdSum_, sectorSpecTerms_;
106
107 sectorPdSum_ = std::vector<Real>(n_, 0.0);
108 sectorExposure_ = std::vector<Real>(n_, 0.0);
109 sectorEl_ = std::vector<Real>(n_, 0.0);
110 sectorSpecTerms_ = std::vector<Real>(n_, 0.0);
111 sectorUl_ = std::vector<Real>(n_, 0.0);
112 marginalLoss_ = std::vector<Real>(m_, 0.0);
113
114 std::vector<Real> pdAdj(m_, 0.0);
115
116 // compute exposure bands
117
118 unsigned long maxNu_ = 0;
119 upperIndex_ = 0;
120
121 // map of nuC_ to expected loss
122 std::map<unsigned long, Real, std::less<> > epsNuC_;
123
124 std::map<unsigned long, Real, std::less<> >::iterator iter;
125
126 for (Size k = 0; k < m_; ++k) {
127 auto exUnit = (unsigned long)(std::floor(0.5 + exposure_[k] / unit_)); // round
128 if (exposure_[k] > 0 && exUnit == 0)
129 exUnit = 1; // but avoid zero exposure
130 if (exUnit > maxNu_)
131 maxNu_ = exUnit;
132 pdAdj[k] = exposure_[k] > 0.0
133 ? exposure_[k] * pd_[k] / (exUnit * unit_)
134 : Real(0.0); // adjusted pd
135 Real el = exUnit * pdAdj[k];
136 if (exUnit > 0) {
137 iter = epsNuC_.find(exUnit);
138 if (iter == epsNuC_.end()) {
139 epsNuC_.insert(std::pair<unsigned long, Real>(exUnit, el));
140 } else {
141 (*iter).second += el;
142 }
143 upperIndex_ += exUnit;
144 }
145 }
146
147 // compute per sector figures
148
149 Real pdSum_ = 0;
150 for (Size k = 0; k < m_; ++k) {
151 pdSum_ += pdAdj[k];
152 sectorPdSum_[sector_[k]] += pd_[k];
154 sectorEl_[sector_[k]] += exposure_[k] * pd_[k];
155 }
156
157 for (Size i = 0; i < n_; ++i) {
158
159 // precompute sector specific terms (formula 15 in [1])
160
161 sectorSpecTerms_[i] += relativeDefaultVariance_[i] * sectorEl_[i];
162 for (Size j = 0; j < n_; ++j) {
163 if (j != i) {
164 sectorSpecTerms_[i] +=
165 correlation_[i][j] *
166 std::sqrt(relativeDefaultVariance_[i] *
168 sectorEl_[j];
169 }
170 }
171 }
172
173 // compute synthetic standard deviation (formula 12 in [1])
174
175 ul_ = 0.0;
176 for (Size i = 0; i < n_; ++i) {
177 sectorUl_[i] =
179 ul_ += sectorUl_[i];
180 for (Size j = 0; j < n_; ++j) {
181 if (j != i) {
182 ul_ += correlation_[i][j] *
183 std::sqrt(relativeDefaultVariance_[i] *
185 sectorEl_[i] * sectorEl_[j];
186 }
187 }
188 }
189
190 Real matchUl_ = ul_; // formula 13 in [1], rhs
191 for (Size k = 0; k < m_; ++k) {
192 Real tmp = pd_[k] * exposure_[k] * exposure_[k];
193 sectorUl_[sector_[k]] += tmp;
194 ul_ += tmp;
195 }
196 ul_ = std::sqrt(ul_);
197 for (Size i = 0; i < n_; ++i)
198 sectorUl_[i] = std::sqrt(sectorUl_[i]);
199
200 // compute risk contributions (formula 15 in [1])
201
202 for (Size k = 0; k < m_; ++k) {
203 marginalLoss_[k] = pd_[k] * exposure_[k] / ul_ *
204 (sectorSpecTerms_[sector_[k]] + exposure_[k]);
205 }
206
207 // compute sigmaC_ and deduced figures
208
209 Real sigmaC_ = pdSum_ * sqrt(matchUl_ / (el_ * el_));
210 Real alphaC_ = pdSum_ * pdSum_ / (sigmaC_ * sigmaC_);
211 Real betaC_ = sigmaC_ * sigmaC_ / pdSum_;
212 Real pC_ = betaC_ / (1.0 + betaC_);
213
214 // compute loss distribution
215
216 loss_.clear();
217 loss_.push_back(std::pow(1.0 - pC_, alphaC_)); // A(0)
218
219 Real res;
220 for (unsigned long n = 0; n < upperIndex_ - 1; ++n) { // compute A(n+1)
221 // recursively
222 res = 0.0;
223 for (unsigned long j = 0;
224 j <= std::min<unsigned long>(maxNu_ - 1, n); ++j) {
225 iter = epsNuC_.find(j + 1);
226 if (iter != epsNuC_.end()) {
227 res += (*iter).second * loss_[n - j] * alphaC_;
228 if (j <= n - 1)
229 res += (*iter).second / ((Real)(j + 1)) *
230 ((Real)(n - j)) * loss_[n - j];
231 }
232 }
233 loss_.push_back(res * pC_ / (pdSum_ * ((Real)(n + 1))));
234 }
235 }
236}
const Real correlation_
std::vector< Real > loss_
const std::vector< Real > exposure_
CreditRiskPlus(std::vector< Real > exposure, std::vector< Real > defaultProbability, std::vector< Size > sector, std::vector< Real > relativeDefaultVariance, Matrix correlation, Real unit)
std::vector< Real > sectorEl_
std::vector< Real > sectorExposure_
const std::vector< Real > relativeDefaultVariance_
const std::vector< Real > pd_
const std::vector< Size > sector_
std::vector< Real > marginalLoss_
std::vector< Real > sectorUl_
Matrix used in linear algebra.
Definition: matrix.hpp:41
Size rows() const
Definition: matrix.hpp:504
Size columns() const
Definition: matrix.hpp:508
Extended CreditRisk+ Model.
#define QL_REQUIRE(condition, message)
throw an error if the given pre-condition is not verified
Definition: errors.hpp:117
QL_REAL Real
real number
Definition: types.hpp:50
std::size_t Size
size of a container
Definition: types.hpp:58
Definition: any.hpp:35
STL namespace.
Real pd_