From b83de157c5be81e0e1ee7df463574ce6cff9c944 2010-02-22 21:25:40
From: Brian Granger <ellisonbg@gmail.com>
Date: 2010-02-22 21:25:40
Subject: [PATCH] Improved Monte-Carlo options pricing example.
---
diff --git a/docs/examples/kernel/mcdriver.py b/docs/examples/kernel/mcdriver.py
index 6f493ee..fe5f6b7 100644
--- a/docs/examples/kernel/mcdriver.py
+++ b/docs/examples/kernel/mcdriver.py
@@ -1,71 +1,148 @@
#!/usr/bin/env python
"""Run a Monte-Carlo options pricer in parallel."""
+#-----------------------------------------------------------------------------
+# Imports
+#-----------------------------------------------------------------------------
+
+import sys
+import time
from IPython.kernel import client
import numpy as np
from mcpricer import price_options
+from matplotlib import pyplot as plt
+
+#-----------------------------------------------------------------------------
+# Setup parameters for the run
+#-----------------------------------------------------------------------------
+
+def ask_question(text, the_type, default):
+ s = '%s [%r]: ' % (text, the_type(default))
+ result = raw_input(s)
+ if result:
+ return the_type(result)
+ else:
+ return the_type(default)
+
+cluster_profile = ask_question("Cluster profile", str, "default")
+price = ask_question("Initial price", float, 100.0)
+rate = ask_question("Interest rate", float, 0.05)
+days = ask_question("Days to expiration", int, 260)
+paths = ask_question("Number of MC paths", int, 10000)
+n_strikes = ask_question("Number of strike values", int, 5)
+min_strike = ask_question("Min strike price", float, 90.0)
+max_strike = ask_question("Max strike price", float, 110.0)
+n_sigmas = ask_question("Number of volatility values", int, 5)
+min_sigma = ask_question("Min volatility", float, 0.1)
+max_sigma = ask_question("Max volatility", float, 0.4)
+
+strike_vals = np.linspace(min_strike, max_strike, n_strikes)
+sigma_vals = np.linspace(min_sigma, max_sigma, n_sigmas)
+
+#-----------------------------------------------------------------------------
+# Setup for parallel calculation
+#-----------------------------------------------------------------------------
# The MultiEngineClient is used to setup the calculation and works with all
# engine.
-mec = client.MultiEngineClient(profile='mycluster')
+mec = client.MultiEngineClient(profile=cluster_profile)
# The TaskClient is an interface to the engines that provides dynamic load
# balancing at the expense of not knowing which engine will execute the code.
-tc = client.TaskClient(profile='mycluster')
+tc = client.TaskClient(profile=cluster_profile)
# Initialize the common code on the engines. This Python module has the
# price_options function that prices the options.
mec.run('mcpricer.py')
-# Define the function that will make up our tasks. We basically want to
-# call the price_options function with all but two arguments (K, sigma)
-# fixed.
-def my_prices(K, sigma):
- S = 100.0
- r = 0.05
- days = 260
- paths = 100000
- return price_options(S, K, sigma, r, days, paths)
-
-# Create arrays of strike prices and volatilities
-nK = 10
-nsigma = 10
-K_vals = np.linspace(90.0, 100.0, nK)
-sigma_vals = np.linspace(0.1, 0.4, nsigma)
-
-# Submit tasks to the TaskClient for each (K, sigma) pair as a MapTask.
-# The MapTask simply applies a function (my_prices) to the arguments:
-# my_prices(K, sigma) and returns the result.
+#-----------------------------------------------------------------------------
+# Perform parallel calculation
+#-----------------------------------------------------------------------------
+
+print "Running parallel calculation over strike prices and volatilities..."
+print "Strike prices: ", strike_vals
+print "Volatilities: ", sigma_vals
+sys.stdout.flush()
+
+# Submit tasks to the TaskClient for each (strike, sigma) pair as a MapTask.
+t1 = time.time()
taskids = []
-for K in K_vals:
+for strike in strike_vals:
for sigma in sigma_vals:
- t = client.MapTask(my_prices, args=(K, sigma))
+ t = client.MapTask(
+ price_options,
+ args=(price, strike, sigma, rate, days, paths)
+ )
taskids.append(tc.run(t))
print "Submitted tasks: ", len(taskids)
+sys.stdout.flush()
# Block until all tasks are completed.
tc.barrier(taskids)
+t2 = time.time()
+t = t2-t1
+
+print "Parallel calculation completed, time = %s s" % t
+print "Collecting results..."
# Get the results using TaskClient.get_task_result.
results = [tc.get_task_result(tid) for tid in taskids]
# Assemble the result into a structured NumPy array.
-prices = np.empty(nK*nsigma,
+prices = np.empty(n_strikes*n_sigmas,
dtype=[('ecall',float),('eput',float),('acall',float),('aput',float)]
)
+
for i, price_tuple in enumerate(results):
prices[i] = price_tuple
-prices.shape = (nK, nsigma)
-K_vals, sigma_vals = np.meshgrid(K_vals, sigma_vals)
+
+prices.shape = (n_strikes, n_sigmas)
+strike_mesh, sigma_mesh = np.meshgrid(strike_vals, sigma_vals)
-def plot_options(sigma_vals, K_vals, prices):
+print "Results are available: strike_mesh, sigma_mesh, prices"
+print "To plot results type 'plot_options(sigma_mesh, strike_mesh, prices)'"
+
+#-----------------------------------------------------------------------------
+# Utilities
+#-----------------------------------------------------------------------------
+
+def plot_options(sigma_mesh, strike_mesh, prices):
"""
- Make a contour plot of the option price in (sigma, K) space.
+ Make a contour plot of the option price in (sigma, strike) space.
"""
- from matplotlib import pyplot as plt
- plt.contourf(sigma_vals, K_vals, prices)
+ plt.figure(1)
+
+ plt.subplot(221)
+ plt.contourf(sigma_mesh, strike_mesh, prices['ecall'])
+ plt.axis('tight')
plt.colorbar()
- plt.title("Option Price")
+ plt.title('European Call')
+ plt.ylabel("Strike Price")
+
+ plt.subplot(222)
+ plt.contourf(sigma_mesh, strike_mesh, prices['acall'])
+ plt.axis('tight')
+ plt.colorbar()
+ plt.title("Asian Call")
+
+ plt.subplot(223)
+ plt.contourf(sigma_mesh, strike_mesh, prices['eput'])
+ plt.axis('tight')
+ plt.colorbar()
+ plt.title("European Put")
plt.xlabel("Volatility")
plt.ylabel("Strike Price")
+
+ plt.subplot(224)
+ plt.contourf(sigma_mesh, strike_mesh, prices['aput'])
+ plt.axis('tight')
+ plt.colorbar()
+ plt.title("Asian Put")
+ plt.xlabel("Volatility")
+
+
+
+
+
+