FiMa I - Homework Assignment No. 5

FiMa1_Hom5_2.ipynb

{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The following exercises will guide you towards implementations of pricing tools for \n",
    "European options in the CRR-model.\n",
    "\n",
    "Concerning the return of those exercises, please send an email with your jupyter (ipynb) file to either\n",
    "\n",
    "          [email protected]   or   [email protected]\n",
    "          \n",
    "until Wednesday, 20 November 2019, 12:15 p.m. Indicate your group members in the email."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Setting** We consider the Cox-Ross-Rubinstein model with time horizon $T>0$.\n",
    "Using the notation from the lecture, we assume\n",
    "$$\n",
    "\t-1 < d < r < u.\n",
    "$$\n",
    "to make the model arbitrage-free and complete. We want to tackle the following tasks:\n",
    "\n",
    "**Setting for Exercise 5.2.1**\n",
    "    \n",
    "Let $T=3$. Calculate the *Price*\n",
    "$$\n",
    "    P^{call}(0):=S^0(0)\\mathbb{E}^{\\mathbb{Q}}\\left[\\frac{(S^1(T)-K)_+}{(1+r)^T}\\right]\n",
    "$$\n",
    "at time $t=0$ for the plain-vanilla call option\n",
    "for parameters \n",
    "$$\n",
    "    S^0(0)=1,\\,\\, S^1(0)=120,\\,\\, d=-0.5,\\,\\, u=0.5,\\,\\,  r=0.1,\\,\\,  K=120\\,\n",
    "$$\n",
    "with unique martingale measure $\\mathbb{Q}$.\n",
    "\n",
    "**Setting for  Exercise 5.2.2**\n",
    "\n",
    "Now we consider a so called *up-and-out-Call* with *knock-out* price $B$ ($B=200$).\n",
    "This claim is a European Call on the asset $S^1$ with maturity $T$ and strike $K$ if\n",
    "the asset price never exceeds the knock-out price $B$ within the time period $[0,T]$.\n",
    "If, however, the asset price is larger than $B$ at any time point $t\\in \\{0,\\ldots,T\\}$,\n",
    "the holder gets nothing.\n",
    "\n",
    "a) Give a mathematical description of the discounted payoff of the\n",
    "up-and-out call-option (On your written homework). \n",
    "\n",
    "b) Calculate the discounted price $C(0)$.\n",
    " Use the same parameters as in 5.2.1.\n",
    "\n",
    "     \n",
    "**The goal of the programming exercise is to write Python code which does this computations for general $T$.**"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Toy example\n",
    "\n",
    "For your convinience, let us give the following python code that calculates (not particularly efficiently) $C(0)$ from Exercise 5.2.1 for $T=1$.  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#In this cell we use the \"Code\"-environment to code something in python. \n",
    "#More precisley we just load all necessary packages.\n",
    "import numpy as np \n",
    "#This package contains a lot of methods to implement mathoperations\n",
    "#in a fast way\n",
    "import matplotlib.pyplot as plt \n",
    "#This package is used for generating plots\n",
    "from scipy.stats import bernoulli \n",
    "#This package is fast for simulating bernoulli random variables\n",
    "from scipy.special import comb\n",
    "#This package allows to compute n over k fast\n",
    "import time \n",
    "#Used to stop the compution time"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#input variables:\n",
    "S0 = ???\n",
    "a = ???\n",
    "b = ???\n",
    "r = ???\n",
    "B = ???\n",
    "K = ???\n",
    "T = 1\n",
    "#Unique equivalent martingale measure:\n",
    "p_star = ???"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Calculating the price via the backwards formula from 5.1 iv)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Backward formula of the expected value for T=1 only!!!\n",
    "tic_1 = time.time()\n",
    "pi_1 = p_star*max(S0*(1+b)-K,0)/(1+r)+(1-p_star)*max(S0*(1+a)-K,0)/(1+r)\n",
    "toc_1 = time.time()\n",
    "print( \"The calculated value is %.2f.\" %(pi_1))\n",
    "print(\"The time needed was %.5f seconds.\" %(toc_1-tic_1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Calculating the price via Monte Carlo"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#The following function provides independent samples of the call at time T\n",
    "def sample_call_value_one(a,b,S0,K,r,M):\n",
    "    random_values = bernoulli.rvs(size=(M,1),p=p_star)\n",
    "    up_moves = random_values*(1+b)\n",
    "    down_moves = -1*(random_values-1)*(1+a)\n",
    "    moves = up_moves+down_moves\n",
    "    return ((np.maximum(S0*np.prod(moves, axis=1)-K,0))/(1+r)**1)  "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Value of pi via Monte Carlo\n",
    "n=10000\n",
    "mc_values = 50\n",
    "tic_2 = time.time()\n",
    "x_values = []\n",
    "pi_M = []\n",
    "pi_M_RMSE = []\n",
    "for i in range(1,mc_values,1):\n",
    "    M=n*i\n",
    "    if i == mc_values-1:\n",
    "        tic_2 = time.time()    \n",
    "    x_values.append([M])\n",
    "    values_call = sample_call_value_one(a,b,S0,K,r,M)\n",
    "    pi_M.append([np.mean(values_call)])\n",
    "    pi_M_RMSE.append([np.sqrt(np.var(values_call,ddof=1))/np.sqrt(n*i)])\n",
    "toc_2 = time.time()\n",
    "print( \"The last calculated value is %.2f.\" %(float(pi_M[len(pi_M)-1][0])))\n",
    "print(\"The time needed for the last run was %.5f seconds.\" %(toc_2-tic_2))  \n",
    "#The following will plot the calculated value for pi\n",
    "plt.subplot(2, 1, 1)\n",
    "plt.plot(x_values, pi_M, 'o')\n",
    "plt.plot(x_values, [pi_1]*len(x_values), '-')\n",
    "plt.plot(x_values, pi_1-1.96*np.array(pi_M_RMSE), '-')\n",
    "plt.plot(x_values, pi_1+1.96*np.array(pi_M_RMSE), '-')\n",
    "plt.title('i*n')\n",
    "plt.ylabel('Expected Value')\n",
    "#The following will plot the root mean square error\n",
    "plt.subplot(2, 1, 2)\n",
    "plt.plot(x_values, pi_M_RMSE, '-')\n",
    "plt.xlabel('i*n')\n",
    "plt.ylabel('Root Mean Square Error')\n",
    "plt.show()    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Exercise 5.2.1 (4 points)\n",
    "\n",
    "Calculate the searched Price of Exercise 5.2.1 above for $T=3$ (rather than $T=1$ as assumed \n",
    "in the code above) via jupyther in three ways:\n",
    "\n",
    "    1. Use the backwards algorithm.\n",
    "    \n",
    "    2. Via a combinatorial argument.\n",
    "\n",
    "    3. Use Monte-Carlo Simulation, stop the time and plot your results as in the toy example.\n",
    "    \n",
    "    4. Why did we not ask you to implement explicit integration over $\\Omega=\\{-1,1\\}^T$ here? Compare in general the numerical complexity, when T=50, for all of the above methods.\n",
    "    \n",
    "    \n",
    "   \n",
    "Next you should calculate the call option price in the **$N$-step CRR model**\n",
    "($N \\in \\mathbb{N}$). Specifically, we want to investigate the\n",
    "convergence of the computed prices to the prices obtained by the\n",
    "corresponding limiting *Black-Scholes pricing formula* as the\n",
    "number of time steps $N$ becomes large. \n",
    "\n",
    "The $N$-step CRR model with time step $T/N$ for $T = 1$, riskless interest rate $r=0$ (for simplicity) and risky asset\n",
    "$S^{(N)}$ whose dynamics are described by\n",
    "\\begin{equation*}\n",
    "S^{(N)}_0 := (S^{(N)})^1(0):=s:= 100 > 0, \\quad S^{(N)}_k:=(S^{(N)})^1(k):= s \\prod_{i=1}^k (1 + R^{(N)}_i) \\quad (k=1,\\ldots,N),\n",
    "\\end{equation*}\n",
    "where $1+R^{(N)}_k$, $k=1,\\ldots,N$, takes the values\n",
    "\\begin{equation*}\n",
    "d_N := e^{-0.3 \\sqrt{T/N}} \\quad \\text{and} \\quad u_N :=\n",
    "e^{0.3 \\sqrt{T/N}}.\n",
    "\\end{equation*}\n",
    "Furthermore the martingale\n",
    "measure $\\mathbb{Q}$ with respect to the filtration\n",
    "$(\\sigma(R^{(N)}_1,\\ldots, R^{(N)}_n))_{n=1,\\ldots,N}$ is given by\n",
    "\\begin{equation*}\n",
    "  q_N :=\\mathbb{Q}[R^{(N)}_k = b_N] = \\frac{1-d_N}{u_N - d_N}, \\quad\n",
    "  1-q_N =\\mathbb{Q}[R^{(N)}_k = a_N] = \\frac{u_N - 1}{u_N - d_N} \\quad (k=1,\\ldots,N),\n",
    "\\end{equation*}\n",
    "where $R^{(N)}_k$ $(k=1,2,\\ldots)$ are i.i.d. under $\\mathbb{Q}$. Note that\n",
    "the (time discrete) price process of the risky asset\n",
    "$(S_k^{(N)})_{k=0,1\\ldots,N}$ can be discretized *in space* via\n",
    "the formula\n",
    "\\begin{equation}\n",
    "  S^{i,N}_k := s d_N^{k-i} u_N^{i}, \\quad 0 \\leq i \\leq k, \\quad 0 \\label{disc}\n",
    "  \\leq k \\leq N,\n",
    "\\end{equation}\n",
    "where $i$ denotes the number of upward movements until time step\n",
    "$k$.    \n",
    "\n",
    "5. Implement the backwards recursion algorithm from 5.1 iv) which allows you to compute the arbitrage-free price $P^{put,(N)}(0)$ at time t = 0 of a European put option with strike K = 90 in the above N-step CRR model with step sizes $$N = 50,60,70,...,1980,1990,2000.$$ Plot the deviation of the computed prices $(P^{put,(N)}(0))_{N =50,60,...,2000}$ from the price obtained by using the corresponding limiting Black-Scholes pricing formula for put options. What do you observe?\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "\n",
    "For your convenience, see below the code lines that demand your attention."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "T = ???"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Backward formula for arbitrary $T$ (1.-4.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#backward formula of the expected value for general T\n",
    "tic_1 = time.time()\n",
    "#start the time\n",
    "asset_val = ???\n",
    "#Vector which contains all possible outcomes of the asset at T\n",
    "call_val = ???\n",
    "#Vector which contains all possible outcomes of the call at T\n",
    "#Using the recursion formula for pricing in the CRR model:\n",
    "for t in np.arange(T,0,-1):\n",
    "    ???\n",
    "pi_1 = ???\n",
    "toc_1 = time.time()\n",
    "#Ending the timer\n",
    "print( \"The calculated value is %.4f.\" %(pi_1))\n",
    "print(\"The time needed was %.4f seconds.\" %(toc_1-tic_1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Monte Carlo for arbitrary $T$ (1.-4.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#The following function provides independent samples of the call at time T\n",
    "def sample_call_value???\n",
    "    ???"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Value of pi via Monte Carlo\n",
    "n=10000\n",
    "mc_values = 50\n",
    "tic_2 = time.time()\n",
    "x_values = []\n",
    "pi_M = []\n",
    "pi_M_RMSE = []\n",
    "for i in range(1,mc_values,1):\n",
    "    M=n*i\n",
    "    if i == mc_values-1:\n",
    "        tic_2 = time.time()    \n",
    "    x_values.append([M])\n",
    "    values_call = ???\n",
    "    pi_M.append([np.mean(values_call)])\n",
    "    pi_M_RMSE.append([np.sqrt(np.var(values_call,ddof=1))/np.sqrt(n*i)])\n",
    "toc_2 = time.time()\n",
    "print( \"The last calculated value is %.2f.\" %(float(pi_M[len(pi_M)-1][0])))\n",
    "print(\"The time needed for the last run was %.5f seconds.\" %(toc_2-tic_2))  \n",
    "#The following will plot the calculated value for pi\n",
    "plt.subplot(2, 1, 1)\n",
    "plt.plot(x_values, pi_M, 'o')\n",
    "plt.plot(x_values, [pi_1]*len(x_values), '-')\n",
    "plt.plot(x_values, pi_1-1.96*np.array(pi_M_RMSE), '-')\n",
    "plt.plot(x_values, pi_1+1.96*np.array(pi_M_RMSE), '-')\n",
    "plt.title('i*n')\n",
    "plt.ylabel('Expected Value')\n",
    "#The following will plot the root mean square error\n",
    "plt.subplot(2, 1, 2)\n",
    "plt.plot(x_values, pi_M_RMSE, '-')\n",
    "plt.xlabel('i*n')\n",
    "plt.ylabel('Root Mean Square Error')\n",
    "plt.show()    "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Pricing in the N-step CRR model (5.)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "##########    Setting                                        ##########\n",
    "\n",
    "S0 = ???\n",
    "K = ???\n",
    "sigma = ???\n",
    "T = ???\n",
    "\n",
    "steps = np.arange(50,2001,10)\n",
    "\n",
    "\n",
    "##########    Price of the Put Option via CRR model          ##########\n",
    "\n",
    "#Starting the timer\n",
    "tic = time.time()\n",
    "\n",
    "#Calculating the quantities of the CRR model\n",
    "delta = ???\n",
    "u = ???\n",
    "d = ???\n",
    "p = ???\n",
    "q = ???\n",
    "n = ???\n",
    "\n",
    "Put_value = [];\n",
    "\n",
    "for it in np.arange(0,n,1):\n",
    "    #Vector which contains all possible outcomes of the asset\n",
    "    SN = ???\n",
    "    #Vector which contains all possible outputs of the put\n",
    "    P = ???\n",
    "    #Using the recursion formula for pricing in the CRR model\n",
    "    for t in np.arange(steps[it],0,-1):\n",
    "        ???\n",
    "      \n",
    "    Put_value.append(P[0])       \n",
    "\n",
    "#Ending the timer and printing the result\n",
    "toc = time.time()\n",
    "comp_time = toc - tic\n",
    "print (\"The computational time was %f seconds and the last approximated value is %f.\"\n",
    "       %(comp_time,Put_value[n-1]))\n",
    "\n",
    "##########    Price of the Put Option via formula           ##########\n",
    "\n",
    "#Starting the timer\n",
    "tic_2 = time.time()\n",
    "\n",
    "#Incredients of the Black sholes formula\n",
    "d1 = ???\n",
    "d2 = ???\n",
    "\n",
    "#Black sholes formula for the Put\n",
    "BS_val_Put = ???\n",
    "\n",
    "#End timer for the formula calculation\n",
    "toc_2 = time.time()\n",
    "comp_time_2 = toc_2 - tic_2\n",
    "\n",
    "##########    Comparison outputs    ########\n",
    "\n",
    "print (\"The value of the put computed with the formula is %f and it took %f seconds.\" \n",
    "       %(BS_val_Put,comp_time_2))\n",
    "\n",
    "print (\"The error between those two is %f and the time difference was %f seconds.\" \n",
    "       %(np.abs(Put_value[n-1]-BS_val_Put),comp_time-comp_time_2))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "%matplotlib inline\n",
    "\n",
    "plt.plot(steps, Put_value, steps, [BS_val_Put]*len(Put_value), 'b-')\n",
    "plt.ylabel('Put value')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Exercise 5.2.2 (4 points)\n",
    "\n",
    "Proceed similarly to Exercise 5.2.1 the price $P^{k.o.}(0)$ of the knock-out option from exercise 5.2.2 above for $T=3$\n",
    "and discuss your findings also for $T=2$ and $T=30$. \n",
    "\n",
    "Note that you will have to specify a computationaly reasonable backwards scheme and a suitable Monte Carlo sampling routine."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "T=???\n",
    "B=???"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Backward formula for arbitrary $T$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#backward formula of the expected value for general T\n",
    "tic_1 = time.time()\n",
    "#start the time\n",
    "???\n",
    "for t in np.arange(T-1,-1,-1):\n",
    "    ???\n",
    "pi_1 = ???\n",
    "toc_1 = time.time()\n",
    "#Ending the timer\n",
    "print( \"The calculated value is %.4f.\" %(pi_1))\n",
    "print(\"The time needed was %.4f seconds.\" %(toc_1-tic_1))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Monte Carlo for arbitrary $T$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#The following function provides independent samples of the knock-out option at time T\n",
    "def sample_knockoutcall_value???\n",
    "    ???"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "#Value of pi via Monte Carlo\n",
    "n=10000\n",
    "mc_values = 50\n",
    "tic_2 = time.time()\n",
    "x_values = []\n",
    "pi_M = []\n",
    "pi_M_RMSE = []\n",
    "for i in range(1,mc_values,1):\n",
    "    M=n*i\n",
    "    if i == mc_values-1:\n",
    "        tic_2 = time.time()    \n",
    "    x_values.append([M])\n",
    "    values_call = ???\n",
    "    pi_M.append([np.mean(values_call)])\n",
    "    pi_M_RMSE.append([np.sqrt(np.var(values_call,ddof=1))/np.sqrt(n*i)])\n",
    "toc_2 = time.time()\n",
    "print( \"The last calculated value is %.2f.\" %(float(pi_M[len(pi_M)-1][0])))\n",
    "print(\"The time needed for the last run was %.5f seconds.\" %(toc_2-tic_2))  \n",
    "#The following will plot the calculated value for pi\n",
    "plt.subplot(2, 1, 1)\n",
    "plt.plot(x_values, pi_M, 'o')\n",
    "plt.plot(x_values, [pi_1]*len(x_values), '-')\n",
    "plt.plot(x_values, pi_1-1.96*np.array(pi_M_RMSE), '-')\n",
    "plt.plot(x_values, pi_1+1.96*np.array(pi_M_RMSE), '-')\n",
    "plt.title('i*n')\n",
    "plt.ylabel('Expected Value')\n",
    "#The following will plot the root mean square error\n",
    "plt.subplot(2, 1, 2)\n",
    "plt.plot(x_values, pi_M_RMSE, '-')\n",
    "plt.xlabel('i*n')\n",
    "plt.ylabel('Root Mean Square Error')\n",
    "plt.show()    "
   ]
  }
 ],
 "metadata": {
  "anaconda-cloud": {},
  "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.7.5"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 1
}