| 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 | |
| 20 | #include "gsr.hpp" |
| 21 | #include "utilities.hpp" |
| 22 | #include <ql/processes/gsrprocess.hpp> |
| 23 | #include <ql/models/shortrate/onefactormodels/gsr.hpp> |
| 24 | #include <ql/instruments/nonstandardswap.hpp> |
| 25 | #include <ql/instruments/nonstandardswaption.hpp> |
| 26 | #include <ql/pricingengines/swaption/gaussian1dswaptionengine.hpp> |
| 27 | #include <ql/pricingengines/swaption/gaussian1djamshidianswaptionengine.hpp> |
| 28 | #include <ql/pricingengines/swaption/gaussian1dnonstandardswaptionengine.hpp> |
| 29 | #include <ql/indexes/swap/euriborswap.hpp> |
| 30 | #include <ql/termstructures/yield/flatforward.hpp> |
| 31 | #include <ql/time/calendars/target.hpp> |
| 32 | #include <ql/processes/hullwhiteprocess.hpp> |
| 33 | #include <ql/models/shortrate/onefactormodels/hullwhite.hpp> |
| 34 | #include <ql/models/shortrate/calibrationhelpers/swaptionhelper.hpp> |
| 35 | #include <ql/quotes/simplequote.hpp> |
| 36 | #include <ql/pricingengines/swaption/jamshidianswaptionengine.hpp> |
| 37 | #include <ql/time/daycounters/actual360.hpp> |
| 38 | #include <ql/time/daycounters/thirty360.hpp> |
| 39 | #include <ql/indexes/ibor/euribor.hpp> |
| 40 | #include <ql/termstructures/volatility/swaption/swaptionconstantvol.hpp> |
| 41 | #include <ql/instruments/makevanillaswap.hpp> |
| 42 | #include <ql/math/optimization/levenbergmarquardt.hpp> |
| 43 | |
| 44 | using namespace QuantLib; |
| 45 | using boost::unit_test_framework::test_suite; |
| 46 | |
| 47 | using std::fabs; |
| 48 | |
| 49 | void GsrTest::testGsrProcess() { |
| 50 | |
| 51 | BOOST_TEST_MESSAGE("Testing GSR process..." ); |
| 52 | |
| 53 | Date refDate = Settings::instance().evaluationDate(); |
| 54 | |
| 55 | // constant reversion, constant volatility, test conditional expectation and |
| 56 | // variance against |
| 57 | // existing HullWhiteForwardProcess |
| 58 | // technically we test two representations of the same constant reversion |
| 59 | // and volatility structure, |
| 60 | // namely with and without step dates |
| 61 | |
| 62 | Real tol = 1E-8; |
| 63 | |
| 64 | Real reversion = 0.01; |
| 65 | Real modelvol = 0.01; |
| 66 | |
| 67 | Handle<YieldTermStructure> yts0(ext::shared_ptr<YieldTermStructure>( |
| 68 | new FlatForward(0, TARGET(), 0.00, Actual365Fixed()))); |
| 69 | |
| 70 | std::vector<Date> stepDates0; |
| 71 | std::vector<Real> vols0(1, modelvol); |
| 72 | std::vector<Real> reversions0(1, reversion); |
| 73 | |
| 74 | std::vector<Date> stepDates1; |
| 75 | for (Size i = 1; i < 60; i++) |
| 76 | stepDates1.push_back(x: refDate + (i * 6 * Months)); |
| 77 | std::vector<Real> vols1(stepDates1.size() + 1, modelvol); |
| 78 | std::vector<Real> reversions1(stepDates1.size() + 1, reversion); |
| 79 | |
| 80 | Real T = 10.0; |
| 81 | do { |
| 82 | |
| 83 | ext::shared_ptr<Gsr> model( |
| 84 | new Gsr(yts0, stepDates0, vols0, reversions0, T)); |
| 85 | ext::shared_ptr<StochasticProcess1D> gsrProcess = |
| 86 | model->stateProcess(); |
| 87 | ext::shared_ptr<Gsr> model2( |
| 88 | new Gsr(yts0, stepDates1, vols1, reversions1, T)); |
| 89 | ext::shared_ptr<StochasticProcess1D> gsrProcess2 = |
| 90 | model2->stateProcess(); |
| 91 | |
| 92 | ext::shared_ptr<HullWhiteForwardProcess> hwProcess( |
| 93 | new HullWhiteForwardProcess(yts0, reversion, modelvol)); |
| 94 | hwProcess->setForwardMeasureTime(T); |
| 95 | |
| 96 | Real w, t, xw, hwVal, gsrVal, gsr2Val; |
| 97 | |
| 98 | t = 0.5; |
| 99 | do { |
| 100 | w = 0.0; |
| 101 | do { |
| 102 | xw = -0.1; |
| 103 | do { |
| 104 | hwVal = hwProcess->expectation(t0: w, x0: xw, dt: t - w); |
| 105 | gsrVal = gsrProcess->expectation(t0: w, x0: xw, dt: t - w); |
| 106 | gsr2Val = gsrProcess2->expectation(t0: w, x0: xw, dt: t - w); |
| 107 | if (fabs(x: hwVal - gsrVal) > tol) |
| 108 | BOOST_ERROR( |
| 109 | "Expectation E^{T=" |
| 110 | << T << "}(x(" << t << ") | x(" << w << ") = " << xw |
| 111 | << " is different in HullWhiteProcess(" << hwVal |
| 112 | << ") and GsrProcess (" << gsrVal << ")" ); |
| 113 | if (fabs(x: hwVal - gsr2Val) > tol) |
| 114 | BOOST_ERROR( |
| 115 | "Expectation E^{T=" |
| 116 | << T << "}(x(" << t << ") | x(" << w << ") = " << xw |
| 117 | << " is different in HullWhiteProcess(" << hwVal |
| 118 | << ") and GsrProcess2 (" << gsr2Val << ")" ); |
| 119 | |
| 120 | hwVal = hwProcess->variance(t0: w, x0: xw, dt: t - w); |
| 121 | gsrVal = gsrProcess->variance(t0: w, x0: xw, dt: t - w); |
| 122 | gsr2Val = gsrProcess2->variance(t0: w, x0: xw, dt: t - w); |
| 123 | if (fabs(x: hwVal - gsrVal) > tol) |
| 124 | BOOST_ERROR("Variance V((x(" |
| 125 | << t << ") | x(" << w << ") = " << xw |
| 126 | << " is different in HullWhiteProcess(" |
| 127 | << hwVal << ") and GsrProcess (" << gsrVal |
| 128 | << ")" ); |
| 129 | if (fabs(x: hwVal - gsr2Val) > tol) |
| 130 | BOOST_ERROR("Variance V((x(" |
| 131 | << t << ") | x(" << w << ") = " << xw |
| 132 | << " is different in HullWhiteProcess(" |
| 133 | << hwVal << ") and GsrProcess2 (" << gsr2Val |
| 134 | << ")" ); |
| 135 | xw += 0.01; |
| 136 | } while (xw <= 0.1); |
| 137 | w += t / 5.0; |
| 138 | } while (w <= t - 0.1); |
| 139 | t += T / 20.0; |
| 140 | } while (t <= T - 0.1); |
| 141 | T += 10.0; |
| 142 | } while (T <= 30.0); |
| 143 | |
| 144 | // time dependent reversion and volatility (test cases to be added) |
| 145 | |
| 146 | Array times(2); |
| 147 | Array vols(3); |
| 148 | Array reversions(3); |
| 149 | |
| 150 | times[0] = 1.0; |
| 151 | times[1] = 2.0; |
| 152 | vols[0] = 0.2; |
| 153 | vols[1] = 0.3; |
| 154 | vols[2] = 0.4; |
| 155 | reversions[0] = 0.50; |
| 156 | reversions[1] = 0.80; |
| 157 | reversions[2] = 1.30; |
| 158 | |
| 159 | GsrProcess p(times, vols, reversions); |
| 160 | p.setForwardMeasureTime(10.0); |
| 161 | |
| 162 | // add more test cases here ... |
| 163 | } |
| 164 | |
| 165 | void GsrTest::testGsrModel() { |
| 166 | |
| 167 | BOOST_TEST_MESSAGE("Testing GSR model..." ); |
| 168 | |
| 169 | Date refDate = Settings::instance().evaluationDate(); |
| 170 | |
| 171 | Real modelvol = 0.01; |
| 172 | Real reversion = 0.01; |
| 173 | |
| 174 | std::vector<Date> stepDates; // no step dates |
| 175 | std::vector<Real> vols(1, modelvol); |
| 176 | std::vector<Real> reversions(1, reversion); |
| 177 | |
| 178 | std::vector<Date> stepDates1; // artificial step dates (should yield the |
| 179 | // same result) |
| 180 | for (Size i = 1; i < 60; i++) |
| 181 | stepDates1.push_back(x: refDate + (i * 6 * Months)); |
| 182 | std::vector<Real> vols1(stepDates1.size() + 1, modelvol); |
| 183 | std::vector<Real> reversions1(stepDates1.size() + 1, reversion); |
| 184 | |
| 185 | Handle<YieldTermStructure> yts(ext::shared_ptr<YieldTermStructure>( |
| 186 | new FlatForward(0, TARGET(), 0.03, Actual365Fixed()))); |
| 187 | ext::shared_ptr<Gsr> model( |
| 188 | new Gsr(yts, stepDates, vols, reversions, 50.0)); |
| 189 | ext::shared_ptr<Gsr> model2( |
| 190 | new Gsr(yts, stepDates1, vols1, reversions1, 50.0)); |
| 191 | ext::shared_ptr<HullWhite> hw(new HullWhite(yts, reversion, modelvol)); |
| 192 | |
| 193 | // test zerobond prices against existing HullWhite model |
| 194 | // technically we test two representations of the same constant reversion |
| 195 | // and volatility structure, |
| 196 | // namely with and without step dates |
| 197 | |
| 198 | Real tol0 = 1E-8; |
| 199 | |
| 200 | Real w, t, xw; |
| 201 | |
| 202 | w = 0.1; |
| 203 | do { |
| 204 | t = w + 0.1; |
| 205 | do { |
| 206 | xw = -0.10; |
| 207 | do { |
| 208 | Real yw = |
| 209 | (xw - model->stateProcess()->expectation(t0: 0.0, x0: 0.0, dt: w)) / |
| 210 | model->stateProcess()->stdDeviation(t0: 0.0, x0: 0.0, dt: w); |
| 211 | Real rw = xw + 0.03; // instantaneous forward is 0.03 |
| 212 | Real gsrVal = model->zerobond(T: t, t: w, y: yw); |
| 213 | Real gsr2Val = model2->zerobond(T: t, t: w, y: yw); |
| 214 | Real hwVal = hw->discountBond(now: w, maturity: t, rate: rw); |
| 215 | if (fabs(x: gsrVal - hwVal) > tol0) |
| 216 | BOOST_ERROR("Zerobond P(" |
| 217 | << w << "," << t << " | x=" << xw << " / y=" |
| 218 | << yw << ") is different in HullWhite (" |
| 219 | << hwVal << ") and Gsr (" << gsrVal << ")" ); |
| 220 | if (fabs(x: gsr2Val - hwVal) > tol0) |
| 221 | BOOST_ERROR("Zerobond P(" |
| 222 | << w << "," << t << " | x=" << xw << " / y=" |
| 223 | << yw << ") is different in HullWhite (" |
| 224 | << hwVal << ") and Gsr2 (" << gsr2Val << ")" ); |
| 225 | xw += 0.01; |
| 226 | } while (xw <= 0.10); |
| 227 | t += 2.5; |
| 228 | } while (t <= 50.0); |
| 229 | w += 5.0; |
| 230 | } while (w <= 50.0); |
| 231 | |
| 232 | // test standard, nonstandard and jamshidian engine against existing Hull |
| 233 | // White Jamshidian engine |
| 234 | |
| 235 | Date expiry = TARGET().advance(date: refDate, period: 5 * Years); |
| 236 | Period tenor = 10 * Years; |
| 237 | ext::shared_ptr<SwapIndex> swpIdx(new EuriborSwapIsdaFixA(tenor, yts)); |
| 238 | Real forward = swpIdx->fixing(fixingDate: expiry); |
| 239 | |
| 240 | ext::shared_ptr<VanillaSwap> underlying = swpIdx->underlyingSwap(fixingDate: expiry); |
| 241 | ext::shared_ptr<VanillaSwap> underlyingFixed = |
| 242 | MakeVanillaSwap(10 * Years, swpIdx->iborIndex(), forward) |
| 243 | .withEffectiveDate(swpIdx->valueDate(fixingDate: expiry)) |
| 244 | .withFixedLegCalendar(cal: swpIdx->fixingCalendar()) |
| 245 | .withFixedLegDayCount(dc: swpIdx->dayCounter()) |
| 246 | .withFixedLegTenor(t: swpIdx->fixedLegTenor()) |
| 247 | .withFixedLegConvention(bdc: swpIdx->fixedLegConvention()) |
| 248 | .withFixedLegTerminationDateConvention( |
| 249 | bdc: swpIdx->fixedLegConvention()); |
| 250 | ext::shared_ptr<Exercise> exercise(new EuropeanExercise(expiry)); |
| 251 | ext::shared_ptr<Swaption> stdswaption( |
| 252 | new Swaption(underlyingFixed, exercise)); |
| 253 | ext::shared_ptr<NonstandardSwaption> nonstdswaption( |
| 254 | new NonstandardSwaption(*stdswaption)); |
| 255 | |
| 256 | stdswaption->setPricingEngine(ext::shared_ptr<PricingEngine>( |
| 257 | new JamshidianSwaptionEngine(hw, yts))); |
| 258 | Real HwJamNpv = stdswaption->NPV(); |
| 259 | |
| 260 | nonstdswaption->setPricingEngine(ext::shared_ptr<PricingEngine>( |
| 261 | new Gaussian1dNonstandardSwaptionEngine(model, 64, 7.0, true, false))); |
| 262 | stdswaption->setPricingEngine(ext::shared_ptr<PricingEngine>( |
| 263 | new Gaussian1dSwaptionEngine(model, 64, 7.0, true, false))); |
| 264 | Real GsrNonStdNpv = nonstdswaption->NPV(); |
| 265 | Real GsrStdNpv = stdswaption->NPV(); |
| 266 | stdswaption->setPricingEngine(ext::shared_ptr<PricingEngine>( |
| 267 | new Gaussian1dJamshidianSwaptionEngine(model))); |
| 268 | Real GsrJamNpv = stdswaption->NPV(); |
| 269 | |
| 270 | if (fabs(x: HwJamNpv - GsrNonStdNpv) > 0.00005) |
| 271 | BOOST_ERROR( |
| 272 | "Jamshidian HW NPV (" |
| 273 | << HwJamNpv |
| 274 | << ") deviates from Gaussian1dNonstandardSwaptionEngine NPV (" |
| 275 | << GsrNonStdNpv << ")" ); |
| 276 | if (fabs(x: HwJamNpv - GsrStdNpv) > 0.00005) |
| 277 | BOOST_ERROR("Jamshidian HW NPV (" |
| 278 | << HwJamNpv |
| 279 | << ") deviates from Gaussian1dSwaptionEngine NPV (" |
| 280 | << GsrStdNpv << ")" ); |
| 281 | if (fabs(x: HwJamNpv - GsrJamNpv) > 0.00005) |
| 282 | BOOST_ERROR("Jamshidian HW NPV (" |
| 283 | << HwJamNpv |
| 284 | << ") deviates from Gaussian1dJamshidianEngine NPV (" |
| 285 | << GsrJamNpv << ")" ); |
| 286 | } |
| 287 | |
| 288 | test_suite *GsrTest::suite() { |
| 289 | auto* suite = BOOST_TEST_SUITE("GSR model tests" ); |
| 290 | suite->add(QUANTLIB_TEST_CASE(&GsrTest::testGsrProcess)); |
| 291 | suite->add(QUANTLIB_TEST_CASE(&GsrTest::testGsrModel)); |
| 292 | return suite; |
| 293 | } |
| 294 | |