Skip to content

Commit

Permalink
added brownian bridge for pricing american options
Browse files Browse the repository at this point in the history
  • Loading branch information
anthonymakarewicz committed Aug 16, 2024
1 parent afa746e commit e02dffc
Show file tree
Hide file tree
Showing 3 changed files with 159 additions and 58 deletions.
3 changes: 2 additions & 1 deletion include/solver/monte_carlo/mc_american.h
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,8 @@ namespace OptionPricer {
const unsigned int& steps);

double calculate_price(const unsigned long& N) const override;
double standardPrice(const unsigned long& N, const double &dt, const double &discountFactor) const ;
double brownianBridgePrice(const unsigned long& N, const double &dt, const double &discountFactor) const ;

private:
std::shared_ptr<BasisFunctionStrategy> basisFunctionStrategy_;
Expand All @@ -28,5 +30,4 @@ namespace OptionPricer {

}


#endif //MC_AMERICAN_H
6 changes: 4 additions & 2 deletions main.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,8 @@
#include "random/distribution/standard_normal_distribution.h"
#include "random/number_generator/sobol_quasi_random_generator.h"
#include <Eigen/Dense>
#include <gtest/internal/gtest-param-util.h>
#include <random/number_generator/random_generator.h>
#include <solver/monte_carlo/basis_function/chebyshev.h>
#include <solver/monte_carlo/basis_function/laguerre.h>
#include <solver/monte_carlo/regression/lasso.h>
Expand All @@ -32,7 +34,7 @@ int main() {
std::string ticker = "AAPL";
double T = 1.0;
double K = 100.0;
double S = 90.0;
double S = 100.0;
double sigma = 0.15;
double r = 0.03;
int dim = 150;
Expand All @@ -47,7 +49,7 @@ int main() {
params.setParameter("K", K);

auto normal = std::make_shared<StandardNormalDistribution>();
auto generator = std::make_shared<SobolGenerator>(normal, dim);
auto generator = std::make_shared<RNGenerator>(normal);
auto brownianMotion = std::make_shared<GeometricBrownianMotionModel>(ticker, marketData);
auto laguerre = std::make_shared<LaguerreBasisFunction>(4);
auto regression = std::make_shared<LeastSquaresRegression>();
Expand Down
208 changes: 153 additions & 55 deletions src/solver/monte_carlo/mc_american.cpp
Original file line number Diff line number Diff line change
@@ -1,4 +1,7 @@
#include "solver/monte_carlo/mc_american.h"

#include <random/number_generator/random_generator.h>

#include "Eigen/Dense"

namespace OptionPricer {
Expand All @@ -16,76 +19,171 @@ namespace OptionPricer {
regressionStrategy_(std::move(regressionStrategy)) {}

double AmericanMCPricer::calculate_price(const unsigned long &N) const {
/* This algorithm uses the Least Square Monte Carlo (LSMC) method as described in
* Longstaff and Schwartz (2001) for pricing American options.
*/
const double dt = option_->getT() / static_cast<double>(steps_);
Eigen::MatrixXd paths(N, steps_ + 1);
const double dt = option_->getT() / static_cast<double>(steps_);
const double discountFactor = exp(-marketData_->getR() * dt);

// Step 1: Simulate paths
for (int i = 0; i < N; ++i) {
paths(i, 0) = marketData_->getStockData(option_->getID())->getPrice();
for (int j = 1; j <= steps_; ++j) {
double z = generator_->generate(j - 1);
paths(i, j) = paths(i, j - 1) * stockModel_->simulateStepPrice(dt, z);
}
if (std::dynamic_pointer_cast<RNGenerator>(generator_)) {
return brownianBridgePrice(N, dt, discountFactor);
}

return standardPrice(N, dt, discountFactor);
}

double AmericanMCPricer::standardPrice(const unsigned long &N,
const double &dt,
const double &discountFactor) const {
/* This algorithm uses the Least Square Monte Carlo (LSMC) method as described in
* Longstaff and Schwartz (2001) for pricing American options.
*/
const double S = marketData_->getStockData(option_->getID())->getPrice();
Eigen::MatrixXd paths(N, steps_+1);
paths.col(0).setConstant(S);

// Step 1: Simulate paths
for (int i = 0; i < N; ++i) {
for (int j = 1; j <= steps_; ++j) {
double z = generator_->generate(j - 1);
paths(i, j) = paths(i, j - 1) * stockModel_->simulateStepPrice(dt, z);
}
}

// Step 2: Initialize americanPrices with the payoff at maturity
Eigen::VectorXd americanPrices = paths.col(steps_).unaryExpr([this](const double& path){
return option_->payoff(path);
});

// Step 3: Perform backward induction
for (int j = static_cast<int>(steps_) - 1; j >= 1; --j) {
std::vector<double> inTheMoneyPaths, discountedCashFlows;
std::vector<bool> isInTheMoney;

// Step 2: Initialize americanPrices with the payoff at maturity
Eigen::VectorXd americanPrices(N);
for (int i = 0; i < N; ++i) {
americanPrices(i) = option_->payoff(paths(i, steps_));
}
const double exerciseValue = option_->payoff(paths(i, j));

// Step 3: Perform backward induction
for (int step = static_cast<int>(steps_) - 1; step >= 0; --step) {
std::vector<double> inTheMoneyPaths, discountedCashFlows;
std::vector<bool> isInTheMoney;
// Collect In-The-Money paths
if (exerciseValue > 0.0) {
inTheMoneyPaths.push_back(paths(i, j)); //
discountedCashFlows.push_back(americanPrices(i) * discountFactor);
isInTheMoney.push_back(true);
} else {
isInTheMoney.push_back(false);
}
}

// Step 4: Regression to estimate continuation values with basis functions
if (!inTheMoneyPaths.empty()) {
Eigen::VectorXd vInTheMoneyPaths = Eigen::Map<Eigen::VectorXd>(inTheMoneyPaths.data(),
inTheMoneyPaths.size());
Eigen::VectorXd vDiscountedCashFlows = Eigen::Map<Eigen::VectorXd>(discountedCashFlows.data(),
discountedCashFlows.size());

const Eigen::MatrixXd inTheMoneyBases = basisFunctionStrategy_->generate(vInTheMoneyPaths);
const Eigen::VectorXd coeffs = regressionStrategy_->solve(inTheMoneyBases, vDiscountedCashFlows);

// Step 5: Determine whether to exercise or continue
int idxITM = 0;
for (int i = 0; i < N; ++i) {
const double exerciseValue = option_->payoff(paths(i, step));
if (exerciseValue > 0.0) {
inTheMoneyPaths.push_back(paths(i, step)); //
discountedCashFlows.push_back(americanPrices(i) * exp(-marketData_->getR() * dt));
isInTheMoney.push_back(true);
if (isInTheMoney[i]) {
const double exerciseValue = option_->payoff(paths(i, j));
const double continuationValue = inTheMoneyBases.row(idxITM) * coeffs;
idxITM++;
if (exerciseValue > continuationValue) {
americanPrices(i) = exerciseValue;
} else {
americanPrices(i) *= discountFactor;
}
} else {
isInTheMoney.push_back(false);
americanPrices(i) *= discountFactor;
}
}
}
}

// Step 6: Compute the option price as the expected discounted payoffs
return americanPrices.mean() * discountFactor;
}

double AmericanMCPricer::brownianBridgePrice(const unsigned long &N,
const double &dt,
const double &discountFactor) const{
/* This algorithm is an extension of the Longstaff and Schwartz (2001)
* which introduce a brownian bridge to store only the simulated paths for 1 time step
* instead thus reducing memory complexity from O(N x m) to O(N) where N is the number of
* simulations and m the number of steps.
*/
const double S = marketData_->getStockData(option_->getID())->getPrice();
const double sigma = marketData_->getStockData(option_->getID())->getSigma();

// Step 1: Simulate paths at maturity
Eigen::VectorXd paths(N);
for (int i = 0; i < N; ++i) {
const double z = generator_->generate(0); // For each path generate a steps_ dimensional point
paths(i) = stockModel_->simulatePriceAtMaturity(option_->getT(), z);
}

// Step 2: Initialize americanPrices with the payoff at maturity
Eigen::VectorXd americanPrices = paths.unaryExpr([this](const double& path) {
return option_->payoff(path);
});

// Step 3: Perform backward induction
for (int j = static_cast<int>(steps_) - 1; j > 0; --j) {
std::vector<double> inTheMoneyPaths, discountedCashFlows;
std::vector<bool> isInTheMoney;

// Generate stock prices backward using Brownian Bridge
for (int i = 0; i < N; ++i) {
const double t_j = static_cast<double>(j) * dt; // Time at current step
const double t_jp1 = static_cast<double>(j+1) * dt; // Time at the next step
const double X_next = log(paths(i) / S);
const double z = generator_->generate(static_cast<int>(steps_) - j);

// Use the Brownian Bridge interpolation formula
const double X_curr = (t_j/t_jp1) * X_next + sqrt(t_j*dt/t_jp1) * sigma * z;
paths(i) = S * exp(X_curr);
const double exerciseValue = option_->payoff(paths(i));

// Collect In-The-Money paths
if (exerciseValue > 0.0) {
inTheMoneyPaths.push_back(paths(i));
discountedCashFlows.push_back(americanPrices(i) * discountFactor);
isInTheMoney.push_back(true);
} else {
isInTheMoney.push_back(false);
}
}

// Step 4: Regression to estimate continuation values with basis functions
if (!inTheMoneyPaths.empty()) {
Eigen::VectorXd vInTheMoneyPaths = Eigen::Map<Eigen::VectorXd>(inTheMoneyPaths.data(),
inTheMoneyPaths.size());
Eigen::VectorXd vDiscountedCashFlows = Eigen::Map<Eigen::VectorXd>(discountedCashFlows.data(),
discountedCashFlows.size());

// Step 4: Regression to estimate continuation values with basis functions
if (!inTheMoneyPaths.empty()) {
Eigen::VectorXd vInTheMoneyPaths = Eigen::Map<Eigen::VectorXd>(inTheMoneyPaths.data(),
inTheMoneyPaths.size());
Eigen::VectorXd vDiscountedCashFlows = Eigen::Map<Eigen::VectorXd>(discountedCashFlows.data(),
discountedCashFlows.size());

const Eigen::MatrixXd inTheMoneyBases = basisFunctionStrategy_->generate(vInTheMoneyPaths);
const Eigen::VectorXd coeffs = regressionStrategy_->solve(inTheMoneyBases, vDiscountedCashFlows);

// Generate basis functions for all paths at the current step
//const Eigen::MatrixXd allPathsBases = basisFunctionStrategy_->generate(paths.col(step));
int idxITM = 0;
// Step 5: Determine whether to exercise or continue
for (int i = 0; i < N; ++i) {
if (isInTheMoney[i]) {
const double exerciseValue = option_->payoff(paths(i, step));
const double continuationValue = inTheMoneyBases.row(idxITM) * coeffs;
idxITM++;
if (exerciseValue > continuationValue) {
americanPrices(i) = exerciseValue;
} else {
americanPrices(i) *= exp(-marketData_->getR() * dt);
}
const Eigen::MatrixXd inTheMoneyBases = basisFunctionStrategy_->generate(vInTheMoneyPaths);
const Eigen::VectorXd coeffs = regressionStrategy_->solve(inTheMoneyBases, vDiscountedCashFlows);

// Step 5: Determine whether to exercise or continue
int idxITM = 0;
for (int i = 0; i < N; ++i) {
if (isInTheMoney[i]) {
const double exerciseValue = option_->payoff(paths(i));
const double continuationValue = coeffs.dot(inTheMoneyBases.row(idxITM));
++idxITM;
if (exerciseValue > continuationValue) {
americanPrices(i) = exerciseValue;
} else {
americanPrices(i) *= exp(-marketData_->getR() * dt);
americanPrices(i) *= discountFactor;
}
} else {
americanPrices(i) *= discountFactor;
}
}
}

// Step 6: Compute the option price as the expected discounted payoffs
return (americanPrices.sum() / static_cast<double>(N));
}

}
// Step 6: Compute the option price as the expected discounted payoffs
return americanPrices.mean() * discountFactor;
}

}

0 comments on commit e02dffc

Please sign in to comment.