MonteCarloPricing-Swaption.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": false
   },
   "source": [
    "#Swaption Pricing : Monte-Carlo Methods"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "*Copyright (c) 2015 Matthias Groncki*\n",
    "\n",
    "\n",
    "Redistribution and use in source and binary forms, with or without\n",
    "modification, are permitted provided that the following conditions are met:\n",
    "\n",
    "    - Redistributions of source code must retain the above copyright notice,\n",
    "    this list of conditions and the following disclaimer.\n",
    "\n",
    "THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS \"AS IS\"\n",
    "AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE\n",
    "IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE\n",
    "DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDERS OR CONTRIBUTORS BE LIABLE\n",
    "FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL\n",
    "DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR\n",
    "SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER\n",
    "CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY,\n",
    "OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE\n",
    "OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.\n",
    "\n",
    "This disclaimer is taken from the QuantLib license"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 103,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# import the used libraries\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "import QuantLib as ql\n",
    "%matplotlib inline"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 104,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Setting evaluation date\n",
    "today = ql.Date(7,4,2015)\n",
    "ql.Settings.instance().setEvaluationDate(today)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 105,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Setup Marketdata\n",
    "rate = ql.SimpleQuote(0.03)\n",
    "rate_handle = ql.QuoteHandle(rate)\n",
    "dc = ql.Actual365Fixed()\n",
    "yts = ql.FlatForward(today, rate_handle, dc)\n",
    "yts.enableExtrapolation()\n",
    "hyts = ql.RelinkableYieldTermStructureHandle(yts)\n",
    "t0_curve = ql.YieldTermStructureHandle(yts)\n",
    "euribor6m = ql.Euribor6M(hyts)\n",
    "cal = ql.TARGET()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 106,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Setup a dummy portfolio with two Swaps\n",
    "def makeSwap(start, maturity, nominal, fixedRate, index, typ=ql.VanillaSwap.Payer):\n",
    "    \"\"\"\n",
    "    creates a plain vanilla swap with fixedLegTenor 1Y\n",
    "    \n",
    "    parameter:\n",
    "        \n",
    "        start (ql.Date) : Start Date\n",
    "        \n",
    "        maturity (ql.Period) : SwapTenor\n",
    "        \n",
    "        nominal (float) : Nominal\n",
    "        \n",
    "        fixedRate (float) : rate paid on fixed leg\n",
    "        \n",
    "        index (ql.IborIndex) : Index\n",
    "        \n",
    "    return: tuple(ql.Swap, list<Dates>) Swap and all fixing dates\n",
    "    \n",
    "        \n",
    "    \"\"\"\n",
    "    end = ql.TARGET().advance(start, maturity)\n",
    "    fixedLegTenor = ql.Period(\"1y\")\n",
    "    fixedLegBDC = ql.ModifiedFollowing\n",
    "    fixedLegDC = ql.Thirty360(ql.Thirty360.BondBasis)\n",
    "    spread = 0.0\n",
    "    fixedSchedule = ql.Schedule(start,\n",
    "                                end, \n",
    "                                fixedLegTenor, \n",
    "                                index.fixingCalendar(), \n",
    "                                fixedLegBDC,\n",
    "                                fixedLegBDC, \n",
    "                                ql.DateGeneration.Backward,\n",
    "                                False)\n",
    "    floatSchedule = ql.Schedule(start,\n",
    "                                end,\n",
    "                                index.tenor(),\n",
    "                                index.fixingCalendar(),\n",
    "                                index.businessDayConvention(),\n",
    "                                index.businessDayConvention(),\n",
    "                                ql.DateGeneration.Backward,\n",
    "                                False)\n",
    "    swap = ql.VanillaSwap(typ, \n",
    "                          nominal,\n",
    "                          fixedSchedule,\n",
    "                          fixedRate,\n",
    "                          fixedLegDC,\n",
    "                          floatSchedule,\n",
    "                          index,\n",
    "                          spread,\n",
    "                          index.dayCounter())\n",
    "    return swap, [index.fixingDate(x) for x in floatSchedule][:-1]\n",
    "\n",
    "\n",
    "def makeSwaption(swap, callDates, settlement):\n",
    "    if len(callDates) == 1:\n",
    "        exercise = ql.EuropeanExercise(callDates[0])\n",
    "    else:\n",
    "        exercise = ql.BermudanExercise(callDates)\n",
    "    return ql.Swaption(swap, exercise, settlement)\n",
    "\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Lets start with an european swaption"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Create a forward starting plain vanilla swap using the helper function above and create a european swaption which allows us to to enter the swap in one year from today."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 107,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "settlementDate = today + ql.Period(\"1Y\")\n",
    "\n",
    "swaps = [makeSwap(settlementDate,\n",
    "                  ql.Period(\"5Y\"),\n",
    "                  1e6,\n",
    "                  0.03047096,\n",
    "                  euribor6m)\n",
    "        ]\n",
    "\n",
    "calldates = [euribor6m.fixingDate(settlementDate)]\n",
    " \n",
    "swaptions = [makeSwaption(swap, \n",
    "                          calldates, \n",
    "                          ql.Settlement.Physical) \n",
    "             for swap, fd in swaps]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Pricing of the underlying at time 0 using the QuantLib pricing engine"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 108,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Swap NPV at time 0: 0.0091\n"
     ]
    }
   ],
   "source": [
    "#%%timeit\n",
    "# Setup pricing engine and calculate the npv of the underlying swap\n",
    "engine = ql.DiscountingSwapEngine(hyts)\n",
    "for swap, fixingDates in swaps:\n",
    "    swap.setPricingEngine(engine)\n",
    "    print(\"Swap NPV at time 0: %.4f\" % swap.NPV())\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###Setup the Gaussian Shortrate model (a.k.a Hull White model)\n",
    "\n",
    "Don't worry about calibration, assume we know the calbriated model parameter"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 109,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "# Stochastic Process "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 1,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "ename": "NameError",
     "evalue": "name 'ql' is not defined",
     "output_type": "error",
     "traceback": [
      "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m",
      "\u001b[0;31mNameError\u001b[0m                                 Traceback (most recent call last)",
      "\u001b[0;32m<ipython-input-1-199b3d41ba3e>\u001b[0m in \u001b[0;36m<module>\u001b[0;34m()\u001b[0m\n\u001b[1;32m      1\u001b[0m \u001b[0;31m# Assume the model is already calibrated either historical or market implied\u001b[0m\u001b[0;34m\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m----> 2\u001b[0;31m volas = [ql.QuoteHandle(ql.SimpleQuote(0.0075)),\n\u001b[0m\u001b[1;32m      3\u001b[0m          ql.QuoteHandle(ql.SimpleQuote(0.0075))]\n\u001b[1;32m      4\u001b[0m \u001b[0mmeanRev\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0mql\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mQuoteHandle\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mql\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mSimpleQuote\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m0.002\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m      5\u001b[0m \u001b[0mmodel\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mql\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mGsr\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mt0_curve\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m[\u001b[0m\u001b[0mtoday\u001b[0m\u001b[0;34m+\u001b[0m\u001b[0;36m100\u001b[0m\u001b[0;34m]\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mvolas\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmeanRev\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;36m16.\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n",
      "\u001b[0;31mNameError\u001b[0m: name 'ql' is not defined"
     ]
    }
   ],
   "source": [
    "# Assume the model is already calibrated either historical or market implied\n",
    "volas = [ql.QuoteHandle(ql.SimpleQuote(0.0075)),\n",
    "         ql.QuoteHandle(ql.SimpleQuote(0.0075))]\n",
    "meanRev = [ql.QuoteHandle(ql.SimpleQuote(0.002))]\n",
    "model = ql.Gsr(t0_curve, [today+100], volas, meanRev, 16.)\n",
    "process = model.stateProcess()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###Calculate the swaption price using an integral pricing engine"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 111,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "swaptionEngine = ql.Gaussian1dSwaptionEngine(model)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 112,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Swaption NPV : 12898.40\n"
     ]
    }
   ],
   "source": [
    "for swaption in swaptions:\n",
    "    swaption.setPricingEngine(swaptionEngine)\n",
    "    print(\"Swaption NPV : %.2f\" % swaption.NPV())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "###Pricing with an Monte Carlo method"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "####Create a swap path pricer in Python \n",
    "\n",
    "Convert all Dates into times in years (using the same DayCounter as in the yieldTermStructure and store all fix cashflows in a numpy array."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 113,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "mcDC = yts.dayCounter()\n",
    "\n",
    "def timeFromReferenceFactory(daycounter, ref):\n",
    "    \"\"\"\n",
    "    returns a function, that calculate the time in years\n",
    "    from a the reference date *ref* to date *dat* \n",
    "    with respect to the given DayCountConvention *daycounter*\n",
    "    \n",
    "    Parameter:\n",
    "        dayCounter (ql.DayCounter)\n",
    "        ref (ql.Date)\n",
    "        \n",
    "    Return:\n",
    "    \n",
    "        f(np.array(ql.Date)) -> np.array(float)\n",
    "    \"\"\"\n",
    "    def impl(dat):\n",
    "        return daycounter.yearFraction(ref, dat)\n",
    "    return np.vectorize(impl)\n",
    "\n",
    "timeFromReference = timeFromReferenceFactory(mcDC, today)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 114,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "fixed_leg = swap.leg(0)\n",
    "n = len(fixed_leg)\n",
    "fixed_times = np.zeros(n)\n",
    "fixed_amounts = np.zeros(n)\n",
    "for i in range(n):\n",
    "    cf = fixed_leg[i]\n",
    "    fixed_times[i] = timeFromReference(cf.date())\n",
    "    fixed_amounts[i] = cf.amount()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 115,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "float_leg = swap.leg(1)\n",
    "n = len(float_leg)\n",
    "float_times = np.zeros(n)\n",
    "float_dcf = np.zeros(n)\n",
    "accrual_start_time = np.zeros(n)\n",
    "accrual_end_time = np.zeros(n)\n",
    "nominals = np.zeros(n)\n",
    "for i in range(n):\n",
    "    # convert base classiborstart_idx Cashflow to\n",
    "    # FloatingRateCoupon\n",
    "    cf = ql.as_floating_rate_coupon(float_leg[i])\n",
    "    iborIndex = cf.index()\n",
    "    value_date = cf.referencePeriodStart()\n",
    "    index_mat = cf.referencePeriodEnd()\n",
    "    # year fraction\n",
    "    float_dcf[i] = cf.accrualPeriod()\n",
    "    # calculate the start and end time\n",
    "    accrual_start_time[i] = timeFromReference(value_date)\n",
    "    accrual_end_time[i] = timeFromReference(index_mat)\n",
    "    # payment time\n",
    "    float_times[i] = timeFromReference(cf.date())\n",
    "    # nominals \n",
    "    nominals[i] = cf.nominal()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 116,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Store all times for which we need a discount factor in one array\n",
    "df_times = np.concatenate([fixed_times, \n",
    "                           accrual_start_time, \n",
    "                           accrual_end_time, \n",
    "                           float_times])\n",
    "df_times = np.unique(df_times)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 117,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Store indices of fix leg payment times in \n",
    "# the df_times array\n",
    "fix_idx = np.in1d(df_times, fixed_times, True)\n",
    "fix_idx = fix_idx.nonzero()\n",
    "# Indices of the floating leg payment times \n",
    "# in the df_times array\n",
    "float_idx = np.in1d(df_times, float_times, True)\n",
    "float_idx = float_idx.nonzero()\n",
    "# Indices of the accrual start and end time\n",
    "# in the df_times array\n",
    "accrual_start_idx = np.in1d(df_times, accrual_start_time, True)\n",
    "accrual_start_idx = accrual_start_idx.nonzero()\n",
    "accrual_end_idx = np.in1d(df_times, accrual_end_time, True)\n",
    "accrual_end_idx = accrual_end_idx.nonzero()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 118,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Convert call date to time\n",
    "callTimes = timeFromReference(calldates)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "####Test the t0 Swap Pricing"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 119,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "t = 0\n",
    "x_t = 0\n",
    "# Calculate all discount factors\n",
    "discount = np.vectorize(lambda T: model.zerobond(T, t, x_t))\n",
    "dfs = discount(df_times)\n",
    "# Calculate fixed leg npv\n",
    "fix_leg_npv = np.sum(fixed_amounts * dfs[fix_idx])\n",
    "# Estimate the index fixings\n",
    "index_fixings = (dfs[accrual_start_idx] / dfs[accrual_end_idx] - 1) \n",
    "index_fixings /= float_dcf\n",
    "# Calculate the floating leg npv\n",
    "float_leg_npv = np.sum(nominals * index_fixings * float_dcf * dfs[float_idx])\n",
    "npv = float_leg_npv - fix_leg_npv\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 120,
   "metadata": {
    "collapsed": false
   },
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "Swap NPV at time 0: 0.0091\n",
      "Error : 0.00000000\n"
     ]
    }
   ],
   "source": [
    "npv\n",
    "print(\"Swap NPV at time 0: %.4f\" % npv)\n",
    "print(\"Error : %.8f\" % (npv - swap.NPV()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "###Monte Carlo Simulation\n",
    "\n",
    "#### Generate time grid and paths"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 121,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Define evaluation grid\n",
    "date_grid = [today] + calldates\n",
    "\n",
    "date_grid = np.unique(np.sort(date_grid))\n",
    "time_grid = np.vectorize(lambda x: ql.Actual365Fixed().yearFraction(today, x))(date_grid)\n",
    "dt = time_grid[1:] - time_grid[:-1]\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 122,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# Random number generator\n",
    "seed = 1\n",
    "urng = ql.MersenneTwisterUniformRng(seed)\n",
    "usrg = ql.MersenneTwisterUniformRsg(len(time_grid)-1,urng)\n",
    "generator = ql.InvCumulativeMersenneTwisterGaussianRsg(usrg)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#%%timeit\n",
    "# Generate N paths\n",
    "M = 100000\n",
    "m = len(time_grid)\n",
    "x = np.zeros((M, m))\n",
    "y = np.zeros((M, m))\n",
    "numeraires = np.zeros((M,m))\n",
    "dfs = np.zeros((M, m, len(df_times)))\n",
    "\n",
    "for n in range(0,M):\n",
    "    numeraires[n, 0] = model.numeraire(0, 0)\n",
    "    \n",
    "for n in range(0,M):\n",
    "    dWs = generator.nextSequence().value()\n",
    "    for i in range(1, len(time_grid)):\n",
    "        t0 = time_grid[i-1]\n",
    "        t1 = time_grid[i]\n",
    "        e = process.expectation(t0, \n",
    "                                x[n,i-1], \n",
    "                                dt[i-1])\n",
    "        std = process.stdDeviation(t0,\n",
    "                                   x[n,i-1],\n",
    "                                   dt[i-1])\n",
    "        x[n,i] = e + dWs[i-1] * std \n",
    "        e_0_0 = process.expectation(0,0,t1)\n",
    "        std_0_0 = process.stdDeviation(0,0,t1)\n",
    "        y[n,i] = (x[n,i] - e_0_0) / std_0_0\n",
    "        numeraires[n ,i] = model.numeraire(t1, y[n, i])\n",
    "        dfs[n,i] = np.vectorize(lambda T : model.zerobond(T, t1, y[n,i]))(df_times)                            "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": 136,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "# plot the paths\n",
    "#for i in range(0, 1000):\n",
    "#    plt.plot(time_grid, x[i,:])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "#%%timeit\n",
    "# Calculate NPV of the underlying swap at expiry\n",
    "index_fixings = dfs[:,-1, accrual_start_idx][:,0,:] / dfs[:, -1, accrual_end_idx][:,0,:] - 1\n",
    "index_fixings /= float_dcf\n",
    "floatLeg_npv = np.sum(index_fixings * float_dcf * dfs[:,-1, float_idx][:,0,:] * nominals, \n",
    "                     axis = 1) \n",
    "fixedLeg_npv = np.sum(fixed_amounts * dfs[:, -1, fix_idx][:,0,:], axis=1)\n",
    "npv = (floatLeg_npv - fixedLeg_npv)\n",
    "# Apply payoff function \n",
    "npv[npv < 0] = 0\n",
    "# Deflate NPV\n",
    "npv = npv / numeraires[:,-1] \n",
    "npv = np.mean(npv) * numeraires[0,0]"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "print(\"MC NPV: %.4f\" % npv)\n",
    "print(\"Error : %.4f \" % (npv-swaption.NPV()))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.5.1"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}