Merge branch 'develop' of https://github.com/QMCSoftware/QMCSoftware …

…into develop

* 'develop' of https://github.com/QMCSoftware/QMCSoftware:
  asian put option. "AsianCall" --> "AsianOption"

demos/importance_sampling.ipynb

@@ -275,7 +275,7 @@
      "data": {
       "text/plain": [
        "Solution: 1.7863         \n",
-       "AsianCall (Integrand Object)\n",
+       "AsianOption (Integrand Object)\n",
        "    volatility      2^(-1)\n",
        "    start_price     30\n",
        "    strike_price    35\n",
@@ -314,7 +314,7 @@
    "source": [
     "distribution = Sobol(dimension)\n",
     "measure = BrownianMotion(distribution)\n",
-    "integrand = AsianCall(measure)\n",
+    "integrand = AsianOption(measure)\n",
     "solution1,data1 = CubQMCSobolG(integrand, abs_tol).integrate()\n",
     "data1"
    ]
@@ -357,7 +357,7 @@
      "data": {
       "text/plain": [
        "Solution: 1.7851         \n",
-       "AsianCall (Integrand Object)\n",
+       "AsianOption (Integrand Object)\n",
        "    volatility      2^(-1)\n",
        "    start_price     30\n",
        "    strike_price    35\n",
@@ -397,7 +397,7 @@
     "drift = 1\n",
     "distribution = Sobol(dimension)\n",
     "measure = BrownianMotion(distribution,drift=drift)\n",
-    "integrand = AsianCall(measure)\n",
+    "integrand = AsianOption(measure)\n",
     "solution2,data2 = CubQMCSobolG(integrand, abs_tol).integrate()\n",
     "data2"
    ]

demos/integration_examples.ipynb

@@ -110,7 +110,7 @@
      "output_type": "stream",
      "text": [
       "Solution: 6.2708         \n",
-      "AsianCall (Integrand Object)\n",
+      "AsianOption (Integrand Object)\n",
       "    volatility      2^(-1)\n",
       "    start_price     30\n",
       "    strike_price    25\n",
@@ -147,7 +147,7 @@
    "source": [
     "distribution = IIDStdUniform(dimension=16, seed=7)\n",
     "measure = BrownianMotion(distribution)\n",
-    "integrand = AsianCall(\n",
+    "integrand = AsianOption(\n",
     "    measure = measure,\n",
     "    volatility = 0.5,\n",
     "    start_price = 30,\n",
@@ -187,7 +187,7 @@
      "output_type": "stream",
      "text": [
       "Solution: 6.2717         \n",
-      "AsianCall (Integrand Object)\n",
+      "AsianOption (Integrand Object)\n",
       "    volatility      2^(-1)\n",
       "    start_price     30\n",
       "    strike_price    25\n",
@@ -224,7 +224,7 @@
    "source": [
     "distribution = IIDStdGaussian(seed=7)\n",
     "measure = BrownianMotion(distribution)\n",
-    "integrand = AsianCall(measure,\n",
+    "integrand = AsianOption(measure,\n",
     "        volatility = 0.5,\n",
     "        start_price = 30,\n",
     "        strike_price = 25,\n",

qmcpy/integrand/__init__.py

@@ -1,4 +1,4 @@
-from .asian_call import AsianCall
+from .asian_option import AsianOption
 from .european_option import EuropeanOption
 from .keister import Keister
 from .linear import Linear

qmcpy/integrand/asian_call.py → qmcpy/integrand/asian_option.py

@@ -2,16 +2,16 @@
 from ._integrand import Integrand
 from ..true_measure import BrownianMotion
 from ..util import ParameterError
-from numpy import array, exp, log, maximum, repeat
+from numpy import *
 
 
-class AsianCall(Integrand):
+class AsianOption(Integrand):
     """
     >>> dd = Sobol(4,seed=7)
     >>> m = BrownianMotion(dd)
-    >>> ac = AsianCall(m)
+    >>> ac = AsianOption(m)
     >>> ac
-    AsianCall (Integrand Object)
+    AsianOption (Integrand Object)
         volatility      2^(-1)
         start_price     30
         strike_price    35
@@ -27,9 +27,9 @@ class AsianCall(Integrand):
     >>> dd2 = Sobol(seed=7)
     >>> m2 = BrownianMotion(dd2,drift=1)
     >>> level_dims = [2,4,8]
-    >>> ac2 = AsianCall(m2,multi_level_dimensions=level_dims)
+    >>> ac2 = AsianOption(m2,multi_level_dimensions=level_dims)
     >>> ac2
-    AsianCall (Integrand Object)
+    AsianOption (Integrand Object)
         volatility      2^(-1)
         start_price     30
         strike_price    35
@@ -51,7 +51,7 @@ class AsianCall(Integrand):
                   'interest_rate','mean_type', 'dimensions', 'dim_fracs']
 
     def __init__(self, measure, volatility=0.5, start_price=30., strike_price=35.,\
-                 interest_rate=0., mean_type='arithmetic', multi_level_dimensions=None):
+                 interest_rate=0., call_put='call', mean_type='arithmetic', multi_level_dimensions=None):
         """
         Args:
             measure (TrueMeasure): A BrownianMotion TrueMeasure object
@@ -64,12 +64,15 @@ def __init__(self, measure, volatility=0.5, start_price=30., strike_price=35.,\
                 Leave as None for single-level problems
         """
         if not isinstance(measure,BrownianMotion):
-            raise ParameterError('AsianCall measure must be a BrownianMotion instance')
+            raise ParameterError('AsianOption measure must be a BrownianMotion instance')
         self.measure = measure
         self.volatility = float(volatility)
         self.start_price = float(start_price)
         self.strike_price = float(strike_price)
         self.interest_rate = float(interest_rate)
+        self.call_put = call_put.lower()
+        if self.call_put not in ['call','put']:
+            raise ParameterError("call_put must be either 'call' or 'put'")
         self.mean_type = mean_type.lower()
         if self.mean_type not in ['arithmetic', 'geometric']:
             raise ParameterError("mean_type must either 'arithmetic' or 'geometric'")
@@ -85,7 +88,7 @@ def __init__(self, measure, volatility=0.5, start_price=30., strike_price=35.,\
             self.dim_fracs = [0.]
             self.leveltype = 'single'
         self.exercise_time = self.measure.time_vector[-1]
-        super(AsianCall,self).__init__()        
+        super(AsianOption,self).__init__()        
 
     def _get_discounted_payoffs(self, stock_path, dimension):
         """
@@ -108,9 +111,12 @@ def _get_discounted_payoffs(self, stock_path, dimension):
                        log(stock_path[:, :-1]).sum(1) +
                        log(stock_path[:, -1]) / 2.) /
                       float(dimension))
-        y = maximum(avg - self.strike_price, 0.) * \
-            exp(-self.interest_rate * self.exercise_time)
-        return y
+        if self.call_put == 'call':
+            y_raw = maximum(avg - self.strike_price, 0)
+        else: # put
+            y_raw = maximum(self.strike_price - avg, 0)
+        y_adj = y_raw * exp(-self.interest_rate * self.exercise_time)
+        return y_adj
 
     def g(self, x, l=0):
         """ See abstract method. """
@@ -119,6 +125,8 @@ def g(self, x, l=0):
         s_fine = self.start_price * exp(
             (self.interest_rate - self.volatility ** 2 / 2.) *
             self.measure.time_vector + self.volatility * x)
+        for xx,yy in zip(*where(s_fine<0)): # if stock becomes <=0, 0 out rest of path
+            s_fine[xx,yy:] = 0
         y = self._get_discounted_payoffs(s_fine, dimension)
         if dim_frac > 0:
             s_course = s_fine[:, int(dim_frac - 1):: int(dim_frac)]

qmcpy/integrand/european_option.py

@@ -2,7 +2,7 @@
 from ..true_measure import BrownianMotion
 from ..discrete_distribution import Sobol
 from ..util import ParameterError
-from numpy import exp, maximum, log, sqrt
+from numpy import *
 from scipy.stats import norm 
 
 
@@ -51,14 +51,15 @@ def __init__(self, measure, volatility=0.5, start_price=30, strike_price=35,
 
     def g(self, x):
         """ See abstract method. """
-        s_last = self.start_price * exp(
+        s = self.start_price * exp(
             (self.interest_rate - self.volatility ** 2 / 2) *
-            self.measure.time_vector[-1] + self.volatility * x[:,-1])
+            self.measure.time_vector + self.volatility * x)
+        for xx,yy in zip(*where(s<0)): # if stock becomes <=0, 0 out rest of path
+            s[xx,yy:] = 0
         if self.call_put == 'call':
-            y_raw = maximum(s_last - self.strike_price, 0)
+            y_raw = maximum(s[:,-1] - self.strike_price, 0)
         else: # put
-            s_last = maximum(s_last,0) # avoid negative stock values
-            y_raw = maximum(self.strike_price - s_last, 0)
+            y_raw = maximum(self.strike_price - s[:,-1], 0)
         y_adj = y_raw * exp(-self.interest_rate * self.exercise_time)
         return y_adj
 

qmcpy/integrand/ml_call_options.py

@@ -52,7 +52,7 @@ def __init__(self, measure, option='european', volatility=.2,
             t_final (float): exercise time
         """
         if not (isinstance(measure,Gaussian) and (measure.mu==0).all() and (measure.sigma==eye(measure.d)).all()):
-            raise ParameterError('AsianCall measure must be a Gaussian instance with mean 0 and variance 1')
+            raise ParameterError('AsianOption measure must be a Gaussian instance with mean 0 and variance 1')
         options = ['european','asian']
         self.option = option.lower()
         if self.option not in options:

qmcpy/stopping_criterion/cub_mc_clt.py

@@ -2,7 +2,7 @@
 from ..accumulate_data import MeanVarData
 from ..discrete_distribution import IIDStdGaussian
 from ..true_measure import Gaussian, BrownianMotion
-from ..integrand import Keister, AsianCall
+from ..integrand import Keister, AsianOption
 from ..util import MaxSamplesWarning
 from numpy import array, ceil, floor, maximum
 from scipy.stats import norm
@@ -45,7 +45,7 @@ class CubMCCLT(StoppingCriterion):
         error_hat       0.050
         confid_int      [1.785 1.885]
         time_integrate  ...
-    >>> ac = AsianCall(
+    >>> ac = AsianOption(
     ...     measure = BrownianMotion(IIDStdGaussian()),
     ...     multi_level_dimensions = [2,4,8])
     >>> sc = CubMCCLT(ac,abs_tol=.05)

test/fasttests/test_integrands.py

@@ -8,21 +8,21 @@
 else: import unittest2 as unittest
 
 
-class TestAsianCall(unittest.TestCase):
-    """ Unit tests for AsianCall Integrand. """
+class TestAsianOption(unittest.TestCase):
+    """ Unit tests for AsianOption Integrand. """
 
     def test_f(self):
         distribution = Sobol(dimension=2)
         measure = BrownianMotion(distribution)
-        integrand = AsianCall(measure)
+        integrand = AsianOption(measure)
         samples = integrand.measure.distribution.gen_samples(4)
         y = integrand.f(samples).squeeze()
         self.assertTrue(y.shape==(4,))        
 
     def test__dim_at_level(self):
         distribution = Sobol(dimension=4)
         measure = BrownianMotion(distribution)
-        integrand = AsianCall(measure, multi_level_dimensions=[4,8])
+        integrand = AsianOption(measure, multi_level_dimensions=[4,8])
         self.assertTrue(integrand._dim_at_level(0)==4)
         self.assertTrue(integrand._dim_at_level(1)==8)
 

test/longtests/test_integrate.py

@@ -30,7 +30,7 @@ def test_asian_option_multi_level(self):
         abs_tol =.05
         distribution = IIDStdGaussian()
         measure = BrownianMotion(distribution)
-        integrand = AsianCall(measure,multi_level_dimensions=[4,16,64])
+        integrand = AsianOption(measure,multi_level_dimensions=[4,16,64])
         solution,data = CubMCCLT(integrand, abs_tol).integrate()
         true_value = 1.7845
         self.assertTrue(abs(solution - true_value) < abs_tol)

workouts/integration_examples/asian_option_multi_level.py

@@ -20,7 +20,7 @@ def asian_option_multi_level(
     # CubMCCLT
     distribution = IIDStdGaussian(seed=7)
     measure = BrownianMotion(distribution)
-    integrand = AsianCall(
+    integrand = AsianOption(
         measure = measure,
         volatility = volatility,
         start_price = start_price,

workouts/integration_examples/asian_option_single_level.py

@@ -12,6 +12,7 @@ def asian_option_single_level(
     start_price = 30,
     strike_price = 25,
     interest_rate = .01,
+    call_put = 'call',
     mean_type = 'geometric',
     abs_tol = .1):
 
@@ -21,35 +22,35 @@ def asian_option_single_level(
     # CubMCCLT
     distribution = IIDStdGaussian(dimension, seed=7)
     measure = BrownianMotion(distribution)
-    integrand = AsianCall(measure, volatility, start_price, strike_price, interest_rate, mean_type)
+    integrand = AsianOption(measure, volatility, start_price, strike_price, interest_rate, call_put, mean_type)
     solution,data = CubMCCLT(integrand,abs_tol=abs_tol).integrate()
     print('%s%s'%(data,bar))
 
     # CubQMCCLT
     distribution = Lattice(dimension, randomize=True, seed=7, backend='MPS')
     measure = BrownianMotion(distribution)
-    integrand = AsianCall(measure, volatility, start_price, strike_price, interest_rate, mean_type)
+    integrand = AsianOption(measure, volatility, start_price, strike_price, interest_rate, call_put, mean_type)
     solution,data = CubQMCCLT(integrand,abs_tol=abs_tol).integrate()
     print('%s%s'%(data,bar))
 
     # CubMCG
     distribution = IIDStdGaussian(dimension, seed=7)
     measure = BrownianMotion(distribution)
-    integrand = AsianCall(measure, volatility, start_price, strike_price, interest_rate, mean_type)
+    integrand = AsianOption(measure, volatility, start_price, strike_price, interest_rate, call_put, mean_type)
     solution,data = CubMCG(integrand,abs_tol=abs_tol).integrate()
     print('%s%s'%(data,bar))
 
     # CubQMCLatticeG
     distribution = Lattice(dimension=dimension, randomize=True, seed=7, backend='GAIL')
     measure = BrownianMotion(distribution)
-    integrand = AsianCall(measure, volatility, start_price, strike_price, interest_rate, mean_type)
+    integrand = AsianOption(measure, volatility, start_price, strike_price, interest_rate, call_put, mean_type)
     solution,data = CubQMCLatticeG(integrand,abs_tol=abs_tol).integrate()
     print('%s%s'%(data,bar))
 
     # CubQMCSobolG
     distribution = Sobol(dimension=dimension, randomize=True, seed=7, backend='QRNG')
     measure = BrownianMotion(distribution)
-    integrand = AsianCall(measure, volatility, start_price, strike_price, interest_rate, mean_type)
+    integrand = AsianOption(measure, volatility, start_price, strike_price, interest_rate, call_put, mean_type)
     solution,data = CubQMCSobolG(integrand,abs_tol=abs_tol).integrate()
     print('%s%s'%(data,bar))
 

workouts/mc_vs_qmc/integrations_asian_call.py

@@ -5,56 +5,56 @@
 def cubmcclt_iidstduniform(dimension, abs_tol, drift):
     distribution = IIDStdUniform(dimension,seed=7)
     measure = BrownianMotion(distribution, drift=drift)
-    integrand = AsianCall(measure)
+    integrand = AsianOption(measure)
     solution,data = CubMCCLT(integrand, abs_tol).integrate()
     return data
 
 def cubmcclt_iidstdgaussian(dimension, abs_tol, drift):
     distribution = IIDStdGaussian(dimension,seed=7)
     measure = BrownianMotion(distribution, drift=drift)
-    integrand = AsianCall(measure)
+    integrand = AsianOption(measure)
     solution,data = CubMCCLT(integrand, abs_tol).integrate()
     return data
 
 def cubmcg_iidstduniform(dimension, abs_tol, drift):
     distribution = IIDStdUniform(dimension,seed=7)
     measure = BrownianMotion(distribution, drift=drift)
-    integrand = AsianCall(measure)
+    integrand = AsianOption(measure)
     solution,data = CubMCG(integrand, abs_tol).integrate()
     return data
 
 def cubmcg_iidstdgaussian(dimension, abs_tol, drift):
     distribution = IIDStdGaussian(dimension,seed=7)
     measure = BrownianMotion(distribution, drift=drift)
-    integrand = AsianCall(measure)
+    integrand = AsianOption(measure)
     solution,data = CubMCG(integrand, abs_tol).integrate()
     return data
 
 def cubqmcclt_lattice(dimension, abs_tol, drift):
     distribution = Lattice(dimension, seed=7, backend="MPS")
     measure = BrownianMotion(distribution, drift=drift)
-    integrand = AsianCall(measure)
+    integrand = AsianOption(measure)
     solution,data = CubQMCCLT(integrand, abs_tol).integrate()
     return data
 
 def cubqmcclt_sobol(dimension, abs_tol, drift):
     distribution = Sobol(dimension, seed=7, backend="QRNG")
     measure = BrownianMotion(distribution, drift=drift)
-    integrand = AsianCall(measure)
+    integrand = AsianOption(measure)
     solution,data = CubQMCCLT(integrand, abs_tol).integrate()
     return data
 
 def cubqmclatticeg(dimension, abs_tol, drift):
     distribution = Lattice(dimension, seed=7, backend="GAIL")
     measure = BrownianMotion(distribution, drift=drift)
-    integrand = AsianCall(measure)
+    integrand = AsianOption(measure)
     solution,data = CubQMCLatticeG(integrand, abs_tol,).integrate()
     return data
 
 def cubqmcsobolg(dimension, abs_tol, drift):
     distribution = Sobol(dimension, seed=7, backend="QRNG")
     measure = BrownianMotion(distribution, drift=drift)
-    integrand = AsianCall(measure)
+    integrand = AsianOption(measure)
     solution,data = CubQMCSobolG(integrand, abs_tol).integrate()
     return data