DiscreteHedging.cpp

This is an example of using the QuantLib Monte Carlo framework.

00001 /* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */
00002 
00020 /*  This example computes profit and loss of a discrete interval hedging
00021     strategy and compares with the results of Derman & Kamal's (Goldman Sachs
00022     Equity Derivatives Research) Research Note: "When You Cannot Hedge
00023     Continuously: The Corrections to Black-Scholes"
00024     http://www.ederman.com/emanuelderman/GSQSpapers/when_you_cannot_hedge.pdf
00025 
00026     Suppose an option hedger sells an European option and receives the
00027     Black-Scholes value as the options premium.
00028     Then he follows a Black-Scholes hedging strategy, rehedging at discrete,
00029     evenly spaced time intervals as the underlying stock changes. At
00030     expiration, the hedger delivers the option payoff to the option holder,
00031     and unwinds the hedge. We are interested in understanding the final
00032     profit or loss of this strategy.
00033 
00034     If the hedger had followed the exact Black-Scholes replication strategy,
00035     re-hedging continuously as the underlying stock evolved towards its final
00036     value at expiration, then, no matter what path the stock took, the final
00037     P&L would be exactly zero. When the replication strategy deviates from
00038     the exact Black-Scholes method, the final P&L may deviate from zero. This
00039     deviation is called the replication error. When the hedger rebalances at
00040     discrete rather than continuous intervals, the hedge is imperfect and the
00041     replication is inexact. The more often hedging occurs, the smaller the
00042     replication error.
00043 
00044     We examine the range of possibilities, computing the replication error.
00045 */
00046 
00047 // the only header you need to use QuantLib
00048 #include <ql/quantlib.hpp>
00049 #include <iostream>
00050 
00051 using namespace QuantLib;
00052 
00053 #if defined(QL_ENABLE_SESSIONS)
00054 namespace QuantLib {
00055 
00056     Integer sessionId() { return 0; }
00057 
00058 }
00059 #endif
00060 
00061 
00062 /* The ReplicationError class carries out Monte Carlo simulations to evaluate
00063    the outcome (the replication error) of the discrete hedging strategy over
00064    different, randomly generated scenarios of future stock price evolution.
00065 */
00066 class ReplicationError
00067 {
00068 public:
00069     ReplicationError(Option::Type type,
00070                      Time maturity,
00071                      Real strike,
00072                      Real s0,
00073                      Volatility sigma,
00074                      Rate r)
00075     : maturity_(maturity), payoff_(type, strike), s0_(s0),
00076       sigma_(sigma), r_(r) {
00077 
00078         // value of the option
00079         DiscountFactor rDiscount = std::exp(-r_*maturity_);
00080         DiscountFactor qDiscount = 1.0;
00081         Real forward = s0_*qDiscount/rDiscount;
00082         Real variance = sigma_*sigma_*maturity_;
00083         boost::shared_ptr<StrikedTypePayoff> payoff(
00084                                              new PlainVanillaPayoff(payoff_));
00085         BlackFormula black(forward,rDiscount,variance,payoff);
00086         std::cout << "Option value: " << black.value() << std::endl;
00087 
00088         // store option's vega, since Derman and Kamal's formula needs it
00089         vega_ = black.vega(maturity_);
00090 
00091         std::cout << std::endl;
00092         std::cout <<
00093             "        |        | P&L  \t|  P&L    | Derman&Kamal | P&L"
00094             "      \t| P&L" << std::endl;
00095 
00096         std::cout <<
00097             "samples | trades | Mean \t| Std Dev | Formula      |"
00098             " skewness \t| kurt." << std::endl;
00099 
00100         std::cout << "---------------------------------"
00101             "----------------------------------------------" << std::endl;
00102     }
00103 
00104     // the actual replication error computation
00105     void compute(Size nTimeSteps, Size nSamples);
00106 private:
00107     Time maturity_;
00108     PlainVanillaPayoff payoff_;
00109     Real s0_;
00110     Volatility sigma_;
00111     Rate r_;
00112     Real vega_;
00113 };
00114 
00115 // The key for the MonteCarlo simulation is to have a PathPricer that
00116 // implements a value(const Path& path) method.
00117 // This method prices the portfolio for each Path of the random variable
00118 class ReplicationPathPricer : public PathPricer<Path> {
00119   public:
00120     // real constructor
00121     ReplicationPathPricer(Option::Type type,
00122                           Real underlying,
00123                           Real strike,
00124                           Rate r,
00125                           Time maturity,
00126                           Volatility sigma)
00127     : type_(type), underlying_(underlying),
00128       strike_(strike), r_(r), maturity_(maturity), sigma_(sigma) {
00129         QL_REQUIRE(strike_ > 0.0, "strike must be positive");
00130         QL_REQUIRE(underlying_ > 0.0, "underlying must be positive");
00131         QL_REQUIRE(r_ >= 0.0,
00132                    "risk free rate (r) must be positive or zero");
00133         QL_REQUIRE(maturity_ > 0.0, "maturity must be positive");
00134         QL_REQUIRE(sigma_ >= 0.0,
00135                    "volatility (sigma) must be positive or zero");
00136 
00137     }
00138     // The value() method encapsulates the pricing code
00139     Real operator()(const Path& path) const;
00140 
00141   private:
00142     Option::Type type_;
00143     Real underlying_, strike_;
00144     Rate r_;
00145     Time maturity_;
00146     Volatility sigma_;
00147 };
00148 
00149 
00150 // Compute Replication Error as in the Derman and Kamal's research note
00151 int main(int, char* [])
00152 {
00153     try {
00154         QL_IO_INIT
00155 
00156         Time maturity = 1.0/12.0;   // 1 month
00157         Real strike = 100;
00158         Real underlying = 100;
00159         Volatility volatility = 0.20; // 20%
00160         Rate riskFreeRate = 0.05; // 5%
00161         ReplicationError rp(Option::Call, maturity, strike, underlying,
00162                 volatility, riskFreeRate);
00163 
00164         Size scenarios = 50000;
00165         Size hedgesNum;
00166 
00167         hedgesNum = 21;
00168         rp.compute(hedgesNum, scenarios);
00169 
00170         hedgesNum = 84;
00171         rp.compute(hedgesNum, scenarios);
00172 
00173         return 0;
00174     } catch (std::exception& e) {
00175         std::cout << e.what() << std::endl;
00176         return 1;
00177     } catch (...) {
00178         std::cout << "unknown error" << std::endl;
00179         return 1;
00180     }
00181 }
00182 
00183 
00184 /* The actual computation of the Profit&Loss for each single path.
00185 
00186    In each scenario N rehedging trades spaced evenly in time over
00187    the life of the option are carried out, using the Black-Scholes
00188    hedge ratio.
00189 */
00190 Real ReplicationPathPricer::operator()(const Path& path) const {
00191 
00192     Size n = path.length()-1;
00193     QL_REQUIRE(n>0, "the path cannot be empty");
00194 
00195     // discrete hedging interval
00196     Time dt = maturity_/n;
00197 
00198     // For simplicity, we assume the stock pays no dividends.
00199     Rate stockDividendYield = 0.0;
00200 
00201     // let's start
00202     Time t = 0;
00203 
00204     // stock value at t=0
00205     Real stock = underlying_;
00206 
00207     // money account at t=0
00208     Real money_account = 0.0;
00209 
00210     /************************/
00211     /*** the initial deal ***/
00212     /************************/
00213     // option fair price (Black-Scholes) at t=0
00214     DiscountFactor rDiscount = std::exp(-r_*maturity_);
00215     DiscountFactor qDiscount = std::exp(-stockDividendYield*maturity_);
00216     Real forward = stock*qDiscount/rDiscount;
00217     Real variance = sigma_*sigma_*maturity_;
00218     boost::shared_ptr<StrikedTypePayoff> payoff(
00219                                        new PlainVanillaPayoff(type_,strike_));
00220     BlackFormula black(forward,rDiscount,variance,payoff);
00221     // sell the option, cash in its premium
00222     money_account += black.value();
00223     // compute delta
00224     Real delta = black.delta(stock);
00225     // delta-hedge the option buying stock
00226     Real stockAmount = delta;
00227     money_account -= stockAmount*stock;
00228 
00229     /**********************************/
00230     /*** hedging during option life ***/
00231     /**********************************/
00232     for (Size step = 0; step < n-1; step++){
00233 
00234         // time flows
00235         t += dt;
00236 
00237         // accruing on the money account
00238         money_account *= std::exp( r_*dt );
00239 
00240         // stock growth:
00241         stock = path[step+1];
00242 
00243         // recalculate option value at the current stock value,
00244         // and the current time to maturity
00245         rDiscount = std::exp(-r_*(maturity_-t));
00246         qDiscount = std::exp(-stockDividendYield*(maturity_-t));
00247         forward = stock*qDiscount/rDiscount;
00248         variance = sigma_*sigma_*(maturity_-t);
00249         BlackFormula black(forward,rDiscount,variance,payoff);
00250 
00251         // recalculate delta
00252         delta = black.delta(stock);
00253 
00254         // re-hedging
00255         money_account -= (delta - stockAmount)*stock;
00256         stockAmount = delta;
00257     }
00258 
00259     /*************************/
00260     /*** option expiration ***/
00261     /*************************/
00262     // last accrual on my money account
00263     money_account *= std::exp( r_*dt );
00264     // last stock growth
00265     stock = path[n];
00266 
00267     // the hedger delivers the option payoff to the option holder
00268     Real optionPayoff = PlainVanillaPayoff(type_, strike_)(stock);
00269     money_account -= optionPayoff;
00270 
00271     // and unwinds the hedge selling his stock position
00272     money_account += stockAmount*stock;
00273 
00274     // final Profit&Loss
00275     return money_account;
00276 }
00277 
00278 
00279 // The computation over nSamples paths of the P&L distribution
00280 void ReplicationError::compute(Size nTimeSteps, Size nSamples)
00281 {
00282     QL_REQUIRE(nTimeSteps>0, "the number of steps must be > 0");
00283 
00284     // hedging interval
00285     // Time tau = maturity_ / nTimeSteps;
00286 
00287     /* Black-Scholes framework: the underlying stock price evolves
00288        lognormally with a fixed known volatility that stays constant
00289        throughout time.
00290     */
00291     Date today = Date::todaysDate();
00292     DayCounter dayCount = Actual365Fixed();
00293     Handle<Quote> stateVariable(
00294                           boost::shared_ptr<Quote>(new SimpleQuote(s0_)));
00295     Handle<YieldTermStructure> riskFreeRate(
00296                           boost::shared_ptr<YieldTermStructure>(
00297                                       new FlatForward(today, r_, dayCount)));
00298     Handle<YieldTermStructure> dividendYield(
00299                           boost::shared_ptr<YieldTermStructure>(
00300                                       new FlatForward(today, 0.0, dayCount)));
00301     Handle<BlackVolTermStructure> volatility(
00302                           boost::shared_ptr<BlackVolTermStructure>(
00303                                new BlackConstantVol(today, sigma_, dayCount)));
00304     boost::shared_ptr<StochasticProcess1D> diffusion(
00305                          new BlackScholesProcess(stateVariable, dividendYield,
00306                                                  riskFreeRate, volatility));
00307 
00308     // Black Scholes equation rules the path generator:
00309     // at each step the log of the stock
00310     // will have drift and sigma^2 variance
00311     PseudoRandom::rsg_type rsg =
00312         PseudoRandom::make_sequence_generator(nTimeSteps, 0);
00313 
00314     bool brownianBridge = false;
00315 
00316     typedef SingleVariate<PseudoRandom>::path_generator_type generator_type;
00317     boost::shared_ptr<generator_type> myPathGenerator(new
00318         generator_type(diffusion, maturity_, nTimeSteps,
00319                        rsg, brownianBridge));
00320 
00321     // The replication strategy's Profit&Loss is computed for each path
00322     // of the stock. The path pricer knows how to price a path using its
00323     // value() method
00324     boost::shared_ptr<PathPricer<Path> > myPathPricer(new
00325         ReplicationPathPricer(payoff_.optionType(), s0_,
00326             payoff_.strike(), r_, maturity_, sigma_));
00327 
00328     // a statistics accumulator for the path-dependant Profit&Loss values
00329     Statistics statisticsAccumulator;
00330 
00331     // The OneFactorMontecarloModel generates paths using myPathGenerator
00332     // each path is priced using myPathPricer
00333     // prices will be accumulated into statisticsAccumulator
00334     OneFactorMonteCarloOption MCSimulation(myPathGenerator,
00335                                            myPathPricer,
00336                                            statisticsAccumulator,
00337                                            false);
00338 
00339     // the model simulates nSamples paths
00340     MCSimulation.addSamples(nSamples);
00341 
00342     // the sampleAccumulator method of OneFactorMonteCarloOption_old
00343     // gives access to all the methods of statisticsAccumulator
00344     Real PLMean  = MCSimulation.sampleAccumulator().mean();
00345     Real PLStDev = MCSimulation.sampleAccumulator().standardDeviation();
00346     Real PLSkew  = MCSimulation.sampleAccumulator().skewness();
00347     Real PLKurt  = MCSimulation.sampleAccumulator().kurtosis();
00348 
00349     // Derman and Kamal's formula
00350     Real theorStD = std::sqrt(M_PI/4/nTimeSteps)*vega_*sigma_;
00351 
00352 
00353     std::cout << std::fixed
00354               << nSamples << "\t| "
00355               << nTimeSteps << "\t | "
00356               << std::setprecision(3) << PLMean << " \t| "
00357               << std::setprecision(2) << PLStDev << " \t  | "
00358               << std::setprecision(2) << theorStD << " \t | "
00359               << std::setprecision(2) << PLSkew << " \t| "
00360               << std::setprecision(2) << PLKurt << std::endl;
00361 }

QuantLib.org
QuantLib
Hosted by
SourceForge.net Logo
Documentation generated by
doxygen