Skip to content

Commit 2ebea30

Browse files
pcaspersjenkins
pcaspers
authored and
jenkins
committed
QPR-12119 implement regressor model
1 parent faf6752 commit 2ebea30

2 files changed

Lines changed: 86 additions & 45 deletions

File tree

‎QuantExt/qle/pricingengines/mcmultilegbaseengine.cpp‎

Lines changed: 78 additions & 42 deletions
Original file line numberDiff line numberDiff line change
@@ -46,12 +46,13 @@ McMultiLegBaseEngine::McMultiLegBaseEngine(
4646
const Size calibrationSeed, const Size pricingSeed, const Size polynomOrder,
4747
const LsmBasisSystem::PolynomialType polynomType, const SobolBrownianGenerator::Ordering ordering,
4848
SobolRsg::DirectionIntegers directionIntegers, const std::vector<Handle<YieldTermStructure>>& discountCurves,
49-
const std::vector<Date>& simulationDates, const std::vector<Size>& externalModelIndices, const bool minimalObsDate)
49+
const std::vector<Date>& simulationDates, const std::vector<Size>& externalModelIndices, const bool minimalObsDate,
50+
const RegressorModel regressorModel)
5051
: model_(model), calibrationPathGenerator_(calibrationPathGenerator), pricingPathGenerator_(pricingPathGenerator),
5152
calibrationSamples_(calibrationSamples), pricingSamples_(pricingSamples), calibrationSeed_(calibrationSeed),
5253
pricingSeed_(pricingSeed), polynomOrder_(polynomOrder), polynomType_(polynomType), ordering_(ordering),
5354
directionIntegers_(directionIntegers), discountCurves_(discountCurves), simulationDates_(simulationDates),
54-
externalModelIndices_(externalModelIndices), minimalObsDate_(minimalObsDate) {
55+
externalModelIndices_(externalModelIndices), minimalObsDate_(minimalObsDate), regressorModel_(regressorModel) {
5556

5657
if (discountCurves_.empty())
5758
discountCurves_.resize(model_->components(CrossAssetModel::AssetType::IR));
@@ -828,8 +829,8 @@ void McMultiLegBaseEngine::calculate() const {
828829

829830
if (exercise_ != nullptr) {
830831
regModelUndExInto[counter] = RegressionModel(
831-
*t, cashflowInfo, [&cfStatus](std::size_t i) { return cfStatus[i] == CfStatus::done; },
832-
model_->stateProcess()->size());
832+
*t, cashflowInfo, [&cfStatus](std::size_t i) { return cfStatus[i] == CfStatus::done; }, **model_,
833+
regressorModel_);
833834
regModelUndExInto[counter].train(polynomOrder_, polynomType_, pathValueUndExInto, pathValuesRef,
834835
simulationTimes);
835836
}
@@ -838,8 +839,8 @@ void McMultiLegBaseEngine::calculate() const {
838839
auto exerciseValue = regModelUndExInto[counter].apply(model_->stateProcess()->initialValues(),
839840
pathValuesRef, simulationTimes);
840841
regModelContinuationValue[counter] = RegressionModel(
841-
*t, cashflowInfo, [&cfStatus](std::size_t i) { return cfStatus[i] == CfStatus::done; },
842-
model_->stateProcess()->size());
842+
*t, cashflowInfo, [&cfStatus](std::size_t i) { return cfStatus[i] == CfStatus::done; }, **model_,
843+
regressorModel_);
843844
regModelContinuationValue[counter].train(polynomOrder_, polynomType_, pathValueOption, pathValuesRef,
844845
simulationTimes,
845846
exerciseValue > RandomVariable(calibrationSamples_, 0));
@@ -849,23 +850,23 @@ void McMultiLegBaseEngine::calculate() const {
849850
exerciseValue > RandomVariable(calibrationSamples_, 0),
850851
pathValueUndExInto, pathValueOption);
851852
regModelOption[counter] = RegressionModel(
852-
*t, cashflowInfo, [&cfStatus](std::size_t i) { return cfStatus[i] == CfStatus::done; },
853-
model_->stateProcess()->size());
853+
*t, cashflowInfo, [&cfStatus](std::size_t i) { return cfStatus[i] == CfStatus::done; }, **model_,
854+
regressorModel_);
854855
regModelOption[counter].train(polynomOrder_, polynomType_, pathValueOption, pathValuesRef, simulationTimes);
855856
}
856857

857858
if (isXvaTime) {
858859
regModelUndDirty[counter] = RegressionModel(
859-
*t, cashflowInfo, [&cfStatus](std::size_t i) { return cfStatus[i] != CfStatus::open; },
860-
model_->stateProcess()->size());
860+
*t, cashflowInfo, [&cfStatus](std::size_t i) { return cfStatus[i] != CfStatus::open; }, **model_,
861+
regressorModel_);
861862
regModelUndDirty[counter].train(polynomOrder_, polynomType_, pathValueUndDirty, pathValuesRef,
862863
simulationTimes);
863864
}
864865

865866
if (exercise_ != nullptr) {
866867
regModelOption[counter] = RegressionModel(
867-
*t, cashflowInfo, [&cfStatus](std::size_t i) { return cfStatus[i] == CfStatus::done; },
868-
model_->stateProcess()->size());
868+
*t, cashflowInfo, [&cfStatus](std::size_t i) { return cfStatus[i] == CfStatus::done; }, **model_,
869+
regressorModel_);
869870
regModelOption[counter].train(polynomOrder_, polynomType_, pathValueOption, pathValuesRef, simulationTimes);
870871
}
871872

@@ -1048,29 +1049,43 @@ std::vector<QuantExt::RandomVariable> McMultiLegBaseEngine::MultiLegBaseAmcCalcu
10481049
McMultiLegBaseEngine::RegressionModel::RegressionModel(const Real observationTime,
10491050
const std::vector<CashflowInfo>& cashflowInfo,
10501051
const std::function<bool(std::size_t)>& cashflowRelevant,
1051-
const Size numberOfModelIndices)
1052-
: observationTime_(observationTime), numberOfModelIndices_(numberOfModelIndices) {
1052+
const CrossAssetModel& model,
1053+
const McMultiLegBaseEngine::RegressorModel regressorModel)
1054+
: observationTime_(observationTime) {
10531055

1054-
// determine (time, modelIndex) pairs to be included in the regressor
1056+
// we always include the full model state as of the observation time
1057+
1058+
for (Size m = 0; m < model.dimension(); ++m) {
1059+
regressorTimesModelIndices_.insert(std::make_pair(observationTime, m));
1060+
}
1061+
1062+
// for LaggedFX we add past fx states
1063+
1064+
if (regressorModel == McMultiLegBaseEngine::RegressorModel::LaggedFX) {
1065+
1066+
std::set<Size> modelFxIndices;
1067+
for (Size i = 1; i < model.components(CrossAssetModel::AssetType::IR); ++i) {
1068+
for (Size j = 0; j < model.stateVariables(CrossAssetModel::AssetType::FX, i - 1); ++j)
1069+
modelFxIndices.insert(model.pIdx(CrossAssetModel::AssetType::FX, i - 1, j));
1070+
}
10551071

1056-
if (false) { // switch off for the moment
10571072
for (Size i = 0; i < cashflowInfo.size(); ++i) {
10581073
if (!cashflowRelevant(i))
10591074
continue;
1075+
// add the cashflow simulation indices
10601076
for (Size j = 0; j < cashflowInfo[i].simulationTimes.size(); ++j) {
10611077
Real t = std::min(observationTime_, cashflowInfo[i].simulationTimes[j]);
1078+
// the simulation time might be zero, but then we want to skip the factors
1079+
if (QuantLib::close_enough(t, 0.0))
1080+
continue;
10621081
for (auto& m : cashflowInfo[i].modelIndices[j]) {
1063-
regressorTimesModelIndices_.insert(std::make_pair(t, m));
1082+
// add FX factors
1083+
if (modelFxIndices.find(m) != modelFxIndices.end())
1084+
regressorTimesModelIndices_.insert(std::make_pair(t, m));
10641085
}
10651086
}
10661087
}
10671088
}
1068-
1069-
// we always include the full model state as of the observation time
1070-
1071-
for (Size m = 0; m < numberOfModelIndices_; ++m) {
1072-
regressorTimesModelIndices_.insert(std::make_pair(observationTime, m));
1073-
}
10741089
}
10751090

10761091
void McMultiLegBaseEngine::RegressionModel::train(const Size polynomOrder,
@@ -1095,13 +1110,25 @@ void McMultiLegBaseEngine::RegressionModel::train(const Size polynomOrder,
10951110
regressor.push_back(paths[std::distance(pathTimes.begin(), pt)][modelIdx]);
10961111
}
10971112

1098-
// get the basis functions
1113+
if (!regressor.empty()) {
10991114

1100-
basisFns_ = RandomVariableLsmBasisSystem::multiPathBasisSystem(regressor.size(), polynomOrder, polynomType);
1115+
// get the basis functions
11011116

1102-
// compute the regression coefficients
1117+
basisFns_ = RandomVariableLsmBasisSystem::multiPathBasisSystem(regressor.size(), polynomOrder, polynomType);
11031118

1104-
regressionCoeffs_ = regressionCoefficients(regressand, regressor, basisFns_, filter);
1119+
// compute the regression coefficients
1120+
1121+
regressionCoeffs_ =
1122+
regressionCoefficients(regressand, regressor, basisFns_, filter, RandomVariableRegressionMethod::QR);
1123+
1124+
} else {
1125+
1126+
// an empty regressor is possible if there are no relevant cashflows, but then the regressand has to be zero too
1127+
1128+
QL_REQUIRE(close_enough_all(regressand, RandomVariable(regressand.size(), 0.0)),
1129+
"McMultiLegBaseEngine::RegressionModel::train(): internal error: regressand is not identically "
1130+
"zero, but no regressor was built.");
1131+
}
11051132

11061133
// update state of model
11071134

@@ -1117,31 +1144,39 @@ McMultiLegBaseEngine::RegressionModel::apply(const Array& initialState,
11171144

11181145
QL_REQUIRE(isTrained_, "McMultiLegBaseEngine::RegressionMdeol::apply(): internal error: model is not trained.");
11191146

1120-
// build initial state pointer
1147+
// determine sample size
11211148

11221149
QL_REQUIRE(!paths.empty() && !paths.front().empty(),
11231150
"McMultiLegBaseEngine::RegressionMdeol::apply(): paths are empty or have empty first component");
11241151
Size samples = paths.front().front()->size();
11251152

1126-
std::vector<RandomVariable> initialStateValues;
1127-
std::vector<const RandomVariable*> initialStatePointer;
1153+
// if we do not have regression coefficients, the regressand was zero
1154+
1155+
if (regressionCoeffs_.empty())
1156+
return RandomVariable(samples, 0.0);
1157+
1158+
// build initial state pointer
1159+
1160+
std::vector<RandomVariable> initialStateValues(initialState.size());
1161+
std::vector<const RandomVariable*> initialStatePointer(initialState.size());
11281162
for (Size j = 0; j < initialState.size(); ++j) {
1129-
initialStateValues.push_back(RandomVariable(samples, initialState[j]));
1130-
initialStatePointer.push_back(&initialStateValues.back());
1163+
initialStateValues[j] = RandomVariable(samples, initialState[j]);
1164+
initialStatePointer[j] = &initialStateValues[j];
11311165
}
11321166

11331167
// build the regressor
11341168

1135-
std::vector<const RandomVariable*> regressor;
1136-
std::vector<RandomVariable> tmp;
1169+
std::vector<const RandomVariable*> regressor(regressorTimesModelIndices_.size());
1170+
std::vector<RandomVariable> tmp(regressorTimesModelIndices_.size());
11371171

1172+
Size i = 0;
11381173
for (auto const& [t, modelIdx] : regressorTimesModelIndices_) {
11391174
auto pt = pathTimes.find(t);
11401175
if (pt != pathTimes.end()) {
11411176

11421177
// the time is a path time, no need to interpolate the path
11431178

1144-
regressor.push_back(paths[std::distance(pathTimes.begin(), pt)][modelIdx]);
1179+
regressor[i] = paths[std::distance(pathTimes.begin(), pt)][modelIdx];
11451180

11461181
} else {
11471182

@@ -1153,32 +1188,33 @@ McMultiLegBaseEngine::RegressionModel::apply(const Array& initialState,
11531188
// t is after last path time => flat extrapolation
11541189

11551190
if (t2 == pathTimes.end()) {
1156-
regressor.push_back(paths[pathTimes.size() - 1][modelIdx]);
1191+
regressor[i] = paths[pathTimes.size() - 1][modelIdx];
11571192
continue;
11581193
}
11591194

11601195
// t is before last path time
11611196

11621197
Real time2 = *t2;
1163-
std::vector<const RandomVariable*> s2 = paths[std::distance(pathTimes.begin(), t2)];
1198+
const RandomVariable* s2 = paths[std::distance(pathTimes.begin(), t2)][modelIdx];
11641199

11651200
Real time1;
1166-
std::vector<const RandomVariable*> s1;
1201+
const RandomVariable* s1;
11671202
if (t2 == pathTimes.begin()) {
11681203
time1 = 0.0;
1169-
s1 = initialStatePointer;
1204+
s1 = initialStatePointer[modelIdx];
11701205
} else {
11711206
time1 = *std::next(t2, -1);
1172-
s1 = paths[std::distance(pathTimes.begin(), std::next(t2, -1))];
1207+
s1 = paths[std::distance(pathTimes.begin(), std::next(t2, -1))][modelIdx];
11731208
}
11741209

11751210
// compute the interpolated state
11761211

11771212
RandomVariable alpha1(samples, (time2 - t) / (time2 - time1));
11781213
RandomVariable alpha2(samples, (t - time1) / (time2 - time1));
1179-
tmp.push_back(alpha1 * *s1[modelIdx] + alpha2 * *s2[modelIdx]);
1180-
regressor.push_back(&tmp.back());
1214+
tmp[i] = alpha1 * *s1 + alpha2 * *s2;
1215+
regressor[i] = &tmp[i];
11811216
}
1217+
++i;
11821218
}
11831219

11841220
// compute result and return it

‎QuantExt/qle/pricingengines/mcmultilegbaseengine.hpp‎

Lines changed: 8 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -52,6 +52,9 @@ struct McEngineStats : public QuantLib::Singleton<McEngineStats> {
5252
};
5353

5454
class McMultiLegBaseEngine {
55+
public:
56+
enum RegressorModel { Simple, LaggedFX };
57+
5558
protected:
5659
/*! The npv is computed in the model's base currency, discounting curves are taken from the model. simulationDates
5760
are additional simulation dates. The cross asset model here must be consistent with the multi path that is the
@@ -69,7 +72,8 @@ class McMultiLegBaseEngine {
6972
const SobolRsg::DirectionIntegers directionIntegers,
7073
const std::vector<Handle<YieldTermStructure>>& discountCurves = std::vector<Handle<YieldTermStructure>>(),
7174
const std::vector<Date>& simulationDates = std::vector<Date>(),
72-
const std::vector<Size>& externalModelIndices = std::vector<Size>(), const bool minimalObsDate = true);
75+
const std::vector<Size>& externalModelIndices = std::vector<Size>(), const bool minimalObsDate = true,
76+
const RegressorModel regressorModel = RegressorModel::Simple);
7377

7478
// run calibration and pricing (called from derived engines)
7579
void calculate() const;
@@ -96,6 +100,7 @@ class McMultiLegBaseEngine {
96100
std::vector<Date> simulationDates_;
97101
std::vector<Size> externalModelIndices_;
98102
bool minimalObsDate_;
103+
RegressorModel regressorModel_;
99104

100105
// the generated amc calculator
101106
mutable boost::shared_ptr<AmcCalculator> amcCalculator_;
@@ -122,7 +127,8 @@ class McMultiLegBaseEngine {
122127
public:
123128
RegressionModel() = default;
124129
RegressionModel(const Real observationTime, const std::vector<CashflowInfo>& cashflowInfo,
125-
const std::function<bool(std::size_t)>& cashflowRelevant, const Size numberOfModelIndices);
130+
const std::function<bool(std::size_t)>& cashflowRelevant, const CrossAssetModel& model,
131+
const RegressorModel regressorModel);
126132
// pathTimes must contain the observation time and the relevant cashflow simulation times
127133
void train(const Size polynomOrder, const LsmBasisSystem::PolynomialType polynomType,
128134
const RandomVariable& regressand, const std::vector<std::vector<const RandomVariable*>>& paths,
@@ -133,7 +139,6 @@ class McMultiLegBaseEngine {
133139

134140
private:
135141
Real observationTime_ = Null<Real>();
136-
Size numberOfModelIndices_ = Null<Size>();
137142
bool isTrained_ = false;
138143
std::set<std::pair<Real, Size>> regressorTimesModelIndices_;
139144
std::vector<std::function<RandomVariable(const std::vector<const RandomVariable*>&)>> basisFns_;

0 commit comments

Comments
 (0)