Showing posts with label python. Show all posts
Showing posts with label python. Show all posts

Saturday, May 16, 2020

Running python cgi scripts on the Raspberry Pi nginx

Basically, the setup of python plugin cgi is here at https://www.takaitra.com/running-python-cgi-scripts-on-the-raspberry-pi/
and the enhanced functions of uwsgi-cgi documentation here https://uwsgi-docs.readthedocs.io/en/latest/CGI.html
Except the followings:

# Build and install the uwsgi with the cgi plugin
wget https://projects.unbit.it/downloads/uwsgi-latest.tar.gz
tar zxvf uwsgi-latest.tar.gz 
cd uwsgi-2.0.18
# compile as cgi plugin
make PROFILE=cgi
sudo cp uwsgi /usr/local/bin/



# Create the file /etc/uwsgi.ini
plugins = cgi
# change to unix sock
socket = /tmp/uwsgi.sock
#socket = 127.0.0.1:9000
module = pyindex
cgi = /var/www/html/cgi-bin
#cgi = /usr/share/nginx/www
cgi-allowed-ext = .py
cgi-helper = .py=python
logger=file:/tmp/uwsgi-error.log
uid = www-data
gid = www-data



# Add a location to the /etc/nginx/sites-available/default
location ~ \.py$ {
  # uwsgi_pass 127.0.0.1:9000;
  # change to unix sock
  uwsgi_pass unix:/tmp/uwsgi.sock;
  include uwsgi_params;
  uwsgi_modifier1 9;
}


Test this python script to show the temperature of Raspberry Pi in a html web page

/var/www/html/cgi-bin/temp.py    Select all
#!/usr/bin/env python import os # Return CPU temperature as a character string def getCPUtemperature(): res = os.popen('vcgencmd measure_temp').readline() return(res.replace("temp=","").replace("'C\n","")) #We have to print a valid HTTP header first so the browser will know how to decode the data print "Content-type: text/html\n\n" temp1=getCPUtemperature() print temp1


/var/www/html/temp.html    Select all
<html> <head> <title>Pi Temp</title> <script src="http://code.jquery.com/jquery-1.10.1.min.js"></script> </head> <body> <h1>Temp from Pi</h1> <script> $(document).ready(function () { var interval = 500; //number of milli seconds between each call var refresh = function() { $.ajax({ url: "temp.py", cache: false, success: function(html) { $('#pi-temp-here').html(html); setTimeout(function() { refresh(); }, interval); } }); }; refresh(); }); </script> <div id="pi-temp-here"></div> </body> </html>


// sed 's/<\([^>]*\)>/\<\1\>/g;'


Shell script    Select all
# Append video Group to www-date user sudo usermod -aG video www-data # reboot the Raspberry Pi and test the python cgi script sudo reboot http://127.0.0.1/temp.html




To install FastCGI for php in ngnix, please follow this guide -> https://getgrav.org/blog/raspberrypi-nginx-php7-dev

If you use buster, it will install the latest php-7.3, so change everything from 7.2 from this guide to 7.3 and the installation of packages will be
sudo apt-get update
sudo apt-get install php php-curl php-gd php-fpm php-cli php-opcache php-mbstring php-xml php-zip


#add in /etc/php/7.3/fpm/pool.d/www.conf
user = pi
group = pi


#reload web server and test
#check to ensure the /var/run/php/php7.3-fpm.sock file exists
sudo service nginx restart
sudo service php7.3-fpm restart




Saturday, April 4, 2020

How to Install TensorFlow with GPU Support on Windows 10 Notebook



Need Windows 10, 64 bits (Home Edition is OK) and NVIDIA Graphics Card with minimum Cuda capability of 3.0(2GB Graphic RAM, or more is better for larger dataset), i5/i7 Notebook with 8GB RAM (more is better) are fine. The latest gaming i7 laptop typically has 16GB RAM or more and NVIDIA Card with 6GB graphics.

This installation guide works for my i7 SAMSUNG notebook with GT 650M.

How to Install TensorFlow with GPU Support on Windows 10


After the installation, and go into tf-gpu python environment
conda activate tf-gpu
test run the followings
python tftest.py    Select all
#!/usr/bin/env python import tensorflow as tf with tf.compat.v1.Session() as sess: hello = tf.constant('Hello, TensorFlow! '+ tf.__version__) print (sess.run(hello).decode()) sess.close()



The following irislearn.py requires the additional installation of python modules
pip install matplotlib pandas sklearn
and download the iris.data from
https://archive.ics.uci.edu/ml/machine-learning-databases/iris/iris.data
or create the iris.data file from the content below
python irislearn.py    Select all
#!/usr/bin/env python # Check the versions of libraries # Python version import sys print('Python: {}'.format(sys.version)) # scipy import scipy print('scipy: {}'.format(scipy.__version__)) # numpy import numpy print('numpy: {}'.format(numpy.__version__)) # matplotlib import matplotlib print('matplotlib: {}'.format(matplotlib.__version__)) # pandas import pandas print('pandas: {}'.format(pandas.__version__)) # scikit-learn import sklearn print('sklearn: {}'.format(sklearn.__version__)) # Load libraries import pandas from pandas.plotting import scatter_matrix import matplotlib.pyplot as plt from sklearn import model_selection from sklearn.metrics import classification_report from sklearn.metrics import confusion_matrix from sklearn.metrics import accuracy_score from sklearn.linear_model import LogisticRegression from sklearn.tree import DecisionTreeClassifier from sklearn.neighbors import KNeighborsClassifier from sklearn.discriminant_analysis import LinearDiscriminantAnalysis from sklearn.naive_bayes import GaussianNB from sklearn.svm import SVC # Load dataset #url = "https://archive.ics.uci.edu/ml/machine-learning-databases/iris/iris.data" url = "iris.data" names = ['sepal-length', 'sepal-width', 'petal-length', 'petal-width', 'class'] dataset = pandas.read_csv(url, names=names) # shape print(dataset.shape) # head print(dataset.head(20)) # descriptions print(dataset.describe()) # class distribution print(dataset.groupby('class').size()) # box and whisker plots dataset.plot(kind='box', subplots=True, layout=(2,2), sharex=False, sharey=False) plt.suptitle("Box and Whisker Plots for inputs") plt.show() # histograms dataset.hist() plt.suptitle('Histograms for inputs') plt.show() # scatter plot matrix scatter_matrix(dataset) plt.suptitle('Scatter Plot Matrix for inputs') plt.show() # Split-out validation dataset array = dataset.values X = array[:,0:4] Y = array[:,4] validation_size = 0.20 seed = 7 X_train, X_validation, Y_train, Y_validation = model_selection.train_test_split(X, Y, test_size=validation_size, random_state=seed) # Test options and evaluation metric seed = 7 scoring = 'accuracy' # Spot Check Algorithms models = [] models.append(('LR', LogisticRegression(max_iter=10000))) models.append(('LDA', LinearDiscriminantAnalysis())) models.append(('KNN', KNeighborsClassifier())) models.append(('CART', DecisionTreeClassifier())) models.append(('NB', GaussianNB())) models.append(('SVM', SVC())) # evaluate each model in turn results = [] names = [] for name, model in models: kfold = model_selection.KFold(n_splits=10, random_state=seed, shuffle=True) cv_results = model_selection.cross_val_score(model, X_train, Y_train, cv=kfold, scoring=scoring) results.append(cv_results) names.append(name) msg = "%s: %f (%f)" % (name, cv_results.mean(), cv_results.std()) print(msg) # Compare Algorithms fig = plt.figure() fig.suptitle('Algorithm Comparison') ax = fig.add_subplot(111) plt.boxplot(results) ax.set_xticklabels(names) plt.show() # Make predictions on validation dataset knn = KNeighborsClassifier() knn.fit(X_train, Y_train) predictions = knn.predict(X_validation) print(accuracy_score(Y_validation, predictions)) print(confusion_matrix(Y_validation, predictions)) print(classification_report(Y_validation, predictions)) print("Make predictions on LogisticRegression Model") model = LogisticRegression(max_iter=10000) model.fit(X_train, Y_train) predictions = model.predict(X_validation) print(accuracy_score(Y_validation, predictions)) print(confusion_matrix(Y_validation, predictions)) print(classification_report(Y_validation, predictions)) # print prediction results on test data for i, prediction in enumerate(predictions): print ('Predicted: %s, Target: %s %s' % (prediction, Y_validation[i], '' if prediction==Y_validation[i] else '(WRONG!!!)'))



iris.data
iris.data    Select all 5.1,3.5,1.4,0.2,Iris-setosa 4.9,3.0,1.4,0.2,Iris-setosa 4.7,3.2,1.3,0.2,Iris-setosa 4.6,3.1,1.5,0.2,Iris-setosa 5.0,3.6,1.4,0.2,Iris-setosa 5.4,3.9,1.7,0.4,Iris-setosa 4.6,3.4,1.4,0.3,Iris-setosa 5.0,3.4,1.5,0.2,Iris-setosa 4.4,2.9,1.4,0.2,Iris-setosa 4.9,3.1,1.5,0.1,Iris-setosa 5.4,3.7,1.5,0.2,Iris-setosa 4.8,3.4,1.6,0.2,Iris-setosa 4.8,3.0,1.4,0.1,Iris-setosa 4.3,3.0,1.1,0.1,Iris-setosa 5.8,4.0,1.2,0.2,Iris-setosa 5.7,4.4,1.5,0.4,Iris-setosa 5.4,3.9,1.3,0.4,Iris-setosa 5.1,3.5,1.4,0.3,Iris-setosa 5.7,3.8,1.7,0.3,Iris-setosa 5.1,3.8,1.5,0.3,Iris-setosa 5.4,3.4,1.7,0.2,Iris-setosa 5.1,3.7,1.5,0.4,Iris-setosa 4.6,3.6,1.0,0.2,Iris-setosa 5.1,3.3,1.7,0.5,Iris-setosa 4.8,3.4,1.9,0.2,Iris-setosa 5.0,3.0,1.6,0.2,Iris-setosa 5.0,3.4,1.6,0.4,Iris-setosa 5.2,3.5,1.5,0.2,Iris-setosa 5.2,3.4,1.4,0.2,Iris-setosa 4.7,3.2,1.6,0.2,Iris-setosa 4.8,3.1,1.6,0.2,Iris-setosa 5.4,3.4,1.5,0.4,Iris-setosa 5.2,4.1,1.5,0.1,Iris-setosa 5.5,4.2,1.4,0.2,Iris-setosa 4.9,3.1,1.5,0.1,Iris-setosa 5.0,3.2,1.2,0.2,Iris-setosa 5.5,3.5,1.3,0.2,Iris-setosa 4.9,3.1,1.5,0.1,Iris-setosa 4.4,3.0,1.3,0.2,Iris-setosa 5.1,3.4,1.5,0.2,Iris-setosa 5.0,3.5,1.3,0.3,Iris-setosa 4.5,2.3,1.3,0.3,Iris-setosa 4.4,3.2,1.3,0.2,Iris-setosa 5.0,3.5,1.6,0.6,Iris-setosa 5.1,3.8,1.9,0.4,Iris-setosa 4.8,3.0,1.4,0.3,Iris-setosa 5.1,3.8,1.6,0.2,Iris-setosa 4.6,3.2,1.4,0.2,Iris-setosa 5.3,3.7,1.5,0.2,Iris-setosa 5.0,3.3,1.4,0.2,Iris-setosa 7.0,3.2,4.7,1.4,Iris-versicolor 6.4,3.2,4.5,1.5,Iris-versicolor 6.9,3.1,4.9,1.5,Iris-versicolor 5.5,2.3,4.0,1.3,Iris-versicolor 6.5,2.8,4.6,1.5,Iris-versicolor 5.7,2.8,4.5,1.3,Iris-versicolor 6.3,3.3,4.7,1.6,Iris-versicolor 4.9,2.4,3.3,1.0,Iris-versicolor 6.6,2.9,4.6,1.3,Iris-versicolor 5.2,2.7,3.9,1.4,Iris-versicolor 5.0,2.0,3.5,1.0,Iris-versicolor 5.9,3.0,4.2,1.5,Iris-versicolor 6.0,2.2,4.0,1.0,Iris-versicolor 6.1,2.9,4.7,1.4,Iris-versicolor 5.6,2.9,3.6,1.3,Iris-versicolor 6.7,3.1,4.4,1.4,Iris-versicolor 5.6,3.0,4.5,1.5,Iris-versicolor 5.8,2.7,4.1,1.0,Iris-versicolor 6.2,2.2,4.5,1.5,Iris-versicolor 5.6,2.5,3.9,1.1,Iris-versicolor 5.9,3.2,4.8,1.8,Iris-versicolor 6.1,2.8,4.0,1.3,Iris-versicolor 6.3,2.5,4.9,1.5,Iris-versicolor 6.1,2.8,4.7,1.2,Iris-versicolor 6.4,2.9,4.3,1.3,Iris-versicolor 6.6,3.0,4.4,1.4,Iris-versicolor 6.8,2.8,4.8,1.4,Iris-versicolor 6.7,3.0,5.0,1.7,Iris-versicolor 6.0,2.9,4.5,1.5,Iris-versicolor 5.7,2.6,3.5,1.0,Iris-versicolor 5.5,2.4,3.8,1.1,Iris-versicolor 5.5,2.4,3.7,1.0,Iris-versicolor 5.8,2.7,3.9,1.2,Iris-versicolor 6.0,2.7,5.1,1.6,Iris-versicolor 5.4,3.0,4.5,1.5,Iris-versicolor 6.0,3.4,4.5,1.6,Iris-versicolor 6.7,3.1,4.7,1.5,Iris-versicolor 6.3,2.3,4.4,1.3,Iris-versicolor 5.6,3.0,4.1,1.3,Iris-versicolor 5.5,2.5,4.0,1.3,Iris-versicolor 5.5,2.6,4.4,1.2,Iris-versicolor 6.1,3.0,4.6,1.4,Iris-versicolor 5.8,2.6,4.0,1.2,Iris-versicolor 5.0,2.3,3.3,1.0,Iris-versicolor 5.6,2.7,4.2,1.3,Iris-versicolor 5.7,3.0,4.2,1.2,Iris-versicolor 5.7,2.9,4.2,1.3,Iris-versicolor 6.2,2.9,4.3,1.3,Iris-versicolor 5.1,2.5,3.0,1.1,Iris-versicolor 5.7,2.8,4.1,1.3,Iris-versicolor 6.3,3.3,6.0,2.5,Iris-virginica 5.8,2.7,5.1,1.9,Iris-virginica 7.1,3.0,5.9,2.1,Iris-virginica 6.3,2.9,5.6,1.8,Iris-virginica 6.5,3.0,5.8,2.2,Iris-virginica 7.6,3.0,6.6,2.1,Iris-virginica 4.9,2.5,4.5,1.7,Iris-virginica 7.3,2.9,6.3,1.8,Iris-virginica 6.7,2.5,5.8,1.8,Iris-virginica 7.2,3.6,6.1,2.5,Iris-virginica 6.5,3.2,5.1,2.0,Iris-virginica 6.4,2.7,5.3,1.9,Iris-virginica 6.8,3.0,5.5,2.1,Iris-virginica 5.7,2.5,5.0,2.0,Iris-virginica 5.8,2.8,5.1,2.4,Iris-virginica 6.4,3.2,5.3,2.3,Iris-virginica 6.5,3.0,5.5,1.8,Iris-virginica 7.7,3.8,6.7,2.2,Iris-virginica 7.7,2.6,6.9,2.3,Iris-virginica 6.0,2.2,5.0,1.5,Iris-virginica 6.9,3.2,5.7,2.3,Iris-virginica 5.6,2.8,4.9,2.0,Iris-virginica 7.7,2.8,6.7,2.0,Iris-virginica 6.3,2.7,4.9,1.8,Iris-virginica 6.7,3.3,5.7,2.1,Iris-virginica 7.2,3.2,6.0,1.8,Iris-virginica 6.2,2.8,4.8,1.8,Iris-virginica 6.1,3.0,4.9,1.8,Iris-virginica 6.4,2.8,5.6,2.1,Iris-virginica 7.2,3.0,5.8,1.6,Iris-virginica 7.4,2.8,6.1,1.9,Iris-virginica 7.9,3.8,6.4,2.0,Iris-virginica 6.4,2.8,5.6,2.2,Iris-virginica 6.3,2.8,5.1,1.5,Iris-virginica 6.1,2.6,5.6,1.4,Iris-virginica 7.7,3.0,6.1,2.3,Iris-virginica 6.3,3.4,5.6,2.4,Iris-virginica 6.4,3.1,5.5,1.8,Iris-virginica 6.0,3.0,4.8,1.8,Iris-virginica 6.9,3.1,5.4,2.1,Iris-virginica 6.7,3.1,5.6,2.4,Iris-virginica 6.9,3.1,5.1,2.3,Iris-virginica 5.8,2.7,5.1,1.9,Iris-virginica 6.8,3.2,5.9,2.3,Iris-virginica 6.7,3.3,5.7,2.5,Iris-virginica 6.7,3.0,5.2,2.3,Iris-virginica 6.3,2.5,5.0,1.9,Iris-virginica 6.5,3.0,5.2,2.0,Iris-virginica 6.2,3.4,5.4,2.3,Iris-virginica 5.9,3.0,5.1,1.8,Iris-virginica



The following keraslearn.py requires the additional installation of python module as below
pip install keras
python keraslearn.py    Select all#!/usr/bin/env python from keras.models import Sequential from keras.layers import Dense import numpy import time # fix random seed for reproducibility numpy.random.seed(7) # load pima indians dataset #dataset = numpy.loadtxt("pima-indians-diabetes.csv", delimiter=",") dataset = numpy.loadtxt("pima-indians-diabetes.data", delimiter=",") # split into input (X) and output (Y) variables X = dataset[:,0:8] Y = dataset[:,8] # create model model = Sequential() model.add(Dense(12, input_dim=8, activation='relu')) model.add(Dense(1, activation='sigmoid')) # Compile model model.compile(loss='binary_crossentropy', optimizer='adam', metrics=['accuracy']) start_time=time.time() # Fit the model #model.fit(X, Y, epochs=150, batch_size=10) model.fit(X, Y, 10, 150) # parameters change to keras 1.2.2 # evaluate the model scores = model.evaluate(X, Y) print("\n%s: %.2f%%" % (model.metrics_names[1], scores[1]*100)) print("\nTraining took %.2f seconds\n" %(time.time()-start_time))



pima-indians-diabetes.data
pima-indians-diabetes.data    Select all 6,148,72,35,0,33.6,0.627,50,1 1,85,66,29,0,26.6,0.351,31,0 8,183,64,0,0,23.3,0.672,32,1 1,89,66,23,94,28.1,0.167,21,0 0,137,40,35,168,43.1,2.288,33,1 5,116,74,0,0,25.6,0.201,30,0 3,78,50,32,88,31.0,0.248,26,1 10,115,0,0,0,35.3,0.134,29,0 2,197,70,45,543,30.5,0.158,53,1 8,125,96,0,0,0.0,0.232,54,1 4,110,92,0,0,37.6,0.191,30,0 10,168,74,0,0,38.0,0.537,34,1 10,139,80,0,0,27.1,1.441,57,0 1,189,60,23,846,30.1,0.398,59,1 5,166,72,19,175,25.8,0.587,51,1 7,100,0,0,0,30.0,0.484,32,1 0,118,84,47,230,45.8,0.551,31,1 7,107,74,0,0,29.6,0.254,31,1 1,103,30,38,83,43.3,0.183,33,0 1,115,70,30,96,34.6,0.529,32,1 3,126,88,41,235,39.3,0.704,27,0 8,99,84,0,0,35.4,0.388,50,0 7,196,90,0,0,39.8,0.451,41,1 9,119,80,35,0,29.0,0.263,29,1 11,143,94,33,146,36.6,0.254,51,1 10,125,70,26,115,31.1,0.205,41,1 7,147,76,0,0,39.4,0.257,43,1 1,97,66,15,140,23.2,0.487,22,0 13,145,82,19,110,22.2,0.245,57,0 5,117,92,0,0,34.1,0.337,38,0 5,109,75,26,0,36.0,0.546,60,0 3,158,76,36,245,31.6,0.851,28,1 3,88,58,11,54,24.8,0.267,22,0 6,92,92,0,0,19.9,0.188,28,0 10,122,78,31,0,27.6,0.512,45,0 4,103,60,33,192,24.0,0.966,33,0 11,138,76,0,0,33.2,0.420,35,0 9,102,76,37,0,32.9,0.665,46,1 2,90,68,42,0,38.2,0.503,27,1 4,111,72,47,207,37.1,1.390,56,1 3,180,64,25,70,34.0,0.271,26,0 7,133,84,0,0,40.2,0.696,37,0 7,106,92,18,0,22.7,0.235,48,0 9,171,110,24,240,45.4,0.721,54,1 7,159,64,0,0,27.4,0.294,40,0 0,180,66,39,0,42.0,1.893,25,1 1,146,56,0,0,29.7,0.564,29,0 2,71,70,27,0,28.0,0.586,22,0 7,103,66,32,0,39.1,0.344,31,1 7,105,0,0,0,0.0,0.305,24,0 1,103,80,11,82,19.4,0.491,22,0 1,101,50,15,36,24.2,0.526,26,0 5,88,66,21,23,24.4,0.342,30,0 8,176,90,34,300,33.7,0.467,58,1 7,150,66,42,342,34.7,0.718,42,0 1,73,50,10,0,23.0,0.248,21,0 7,187,68,39,304,37.7,0.254,41,1 0,100,88,60,110,46.8,0.962,31,0 0,146,82,0,0,40.5,1.781,44,0 0,105,64,41,142,41.5,0.173,22,0 2,84,0,0,0,0.0,0.304,21,0 8,133,72,0,0,32.9,0.270,39,1 5,44,62,0,0,25.0,0.587,36,0 2,141,58,34,128,25.4,0.699,24,0 7,114,66,0,0,32.8,0.258,42,1 5,99,74,27,0,29.0,0.203,32,0 0,109,88,30,0,32.5,0.855,38,1 2,109,92,0,0,42.7,0.845,54,0 1,95,66,13,38,19.6,0.334,25,0 4,146,85,27,100,28.9,0.189,27,0 2,100,66,20,90,32.9,0.867,28,1 5,139,64,35,140,28.6,0.411,26,0 13,126,90,0,0,43.4,0.583,42,1 4,129,86,20,270,35.1,0.231,23,0 1,79,75,30,0,32.0,0.396,22,0 1,0,48,20,0,24.7,0.140,22,0 7,62,78,0,0,32.6,0.391,41,0 5,95,72,33,0,37.7,0.370,27,0 0,131,0,0,0,43.2,0.270,26,1 2,112,66,22,0,25.0,0.307,24,0 3,113,44,13,0,22.4,0.140,22,0 2,74,0,0,0,0.0,0.102,22,0 7,83,78,26,71,29.3,0.767,36,0 0,101,65,28,0,24.6,0.237,22,0 5,137,108,0,0,48.8,0.227,37,1 2,110,74,29,125,32.4,0.698,27,0 13,106,72,54,0,36.6,0.178,45,0 2,100,68,25,71,38.5,0.324,26,0 15,136,70,32,110,37.1,0.153,43,1 1,107,68,19,0,26.5,0.165,24,0 1,80,55,0,0,19.1,0.258,21,0 4,123,80,15,176,32.0,0.443,34,0 7,81,78,40,48,46.7,0.261,42,0 4,134,72,0,0,23.8,0.277,60,1 2,142,82,18,64,24.7,0.761,21,0 6,144,72,27,228,33.9,0.255,40,0 2,92,62,28,0,31.6,0.130,24,0 1,71,48,18,76,20.4,0.323,22,0 6,93,50,30,64,28.7,0.356,23,0 1,122,90,51,220,49.7,0.325,31,1 1,163,72,0,0,39.0,1.222,33,1 1,151,60,0,0,26.1,0.179,22,0 0,125,96,0,0,22.5,0.262,21,0 1,81,72,18,40,26.6,0.283,24,0 2,85,65,0,0,39.6,0.930,27,0 1,126,56,29,152,28.7,0.801,21,0 1,96,122,0,0,22.4,0.207,27,0 4,144,58,28,140,29.5,0.287,37,0 3,83,58,31,18,34.3,0.336,25,0 0,95,85,25,36,37.4,0.247,24,1 3,171,72,33,135,33.3,0.199,24,1 8,155,62,26,495,34.0,0.543,46,1 1,89,76,34,37,31.2,0.192,23,0 4,76,62,0,0,34.0,0.391,25,0 7,160,54,32,175,30.5,0.588,39,1 4,146,92,0,0,31.2,0.539,61,1 5,124,74,0,0,34.0,0.220,38,1 5,78,48,0,0,33.7,0.654,25,0 4,97,60,23,0,28.2,0.443,22,0 4,99,76,15,51,23.2,0.223,21,0 0,162,76,56,100,53.2,0.759,25,1 6,111,64,39,0,34.2,0.260,24,0 2,107,74,30,100,33.6,0.404,23,0 5,132,80,0,0,26.8,0.186,69,0 0,113,76,0,0,33.3,0.278,23,1 1,88,30,42,99,55.0,0.496,26,1 3,120,70,30,135,42.9,0.452,30,0 1,118,58,36,94,33.3,0.261,23,0 1,117,88,24,145,34.5,0.403,40,1 0,105,84,0,0,27.9,0.741,62,1 4,173,70,14,168,29.7,0.361,33,1 9,122,56,0,0,33.3,1.114,33,1 3,170,64,37,225,34.5,0.356,30,1 8,84,74,31,0,38.3,0.457,39,0 2,96,68,13,49,21.1,0.647,26,0 2,125,60,20,140,33.8,0.088,31,0 0,100,70,26,50,30.8,0.597,21,0 0,93,60,25,92,28.7,0.532,22,0 0,129,80,0,0,31.2,0.703,29,0 5,105,72,29,325,36.9,0.159,28,0 3,128,78,0,0,21.1,0.268,55,0 5,106,82,30,0,39.5,0.286,38,0 2,108,52,26,63,32.5,0.318,22,0 10,108,66,0,0,32.4,0.272,42,1 4,154,62,31,284,32.8,0.237,23,0 0,102,75,23,0,0.0,0.572,21,0 9,57,80,37,0,32.8,0.096,41,0 2,106,64,35,119,30.5,1.400,34,0 5,147,78,0,0,33.7,0.218,65,0 2,90,70,17,0,27.3,0.085,22,0 1,136,74,50,204,37.4,0.399,24,0 4,114,65,0,0,21.9,0.432,37,0 9,156,86,28,155,34.3,1.189,42,1 1,153,82,42,485,40.6,0.687,23,0 8,188,78,0,0,47.9,0.137,43,1 7,152,88,44,0,50.0,0.337,36,1 2,99,52,15,94,24.6,0.637,21,0 1,109,56,21,135,25.2,0.833,23,0 2,88,74,19,53,29.0,0.229,22,0 17,163,72,41,114,40.9,0.817,47,1 4,151,90,38,0,29.7,0.294,36,0 7,102,74,40,105,37.2,0.204,45,0 0,114,80,34,285,44.2,0.167,27,0 2,100,64,23,0,29.7,0.368,21,0 0,131,88,0,0,31.6,0.743,32,1 6,104,74,18,156,29.9,0.722,41,1 3,148,66,25,0,32.5,0.256,22,0 4,120,68,0,0,29.6,0.709,34,0 4,110,66,0,0,31.9,0.471,29,0 3,111,90,12,78,28.4,0.495,29,0 6,102,82,0,0,30.8,0.180,36,1 6,134,70,23,130,35.4,0.542,29,1 2,87,0,23,0,28.9,0.773,25,0 1,79,60,42,48,43.5,0.678,23,0 2,75,64,24,55,29.7,0.370,33,0 8,179,72,42,130,32.7,0.719,36,1 6,85,78,0,0,31.2,0.382,42,0 0,129,110,46,130,67.1,0.319,26,1 5,143,78,0,0,45.0,0.190,47,0 5,130,82,0,0,39.1,0.956,37,1 6,87,80,0,0,23.2,0.084,32,0 0,119,64,18,92,34.9,0.725,23,0 1,0,74,20,23,27.7,0.299,21,0 5,73,60,0,0,26.8,0.268,27,0 4,141,74,0,0,27.6,0.244,40,0 7,194,68,28,0,35.9,0.745,41,1 8,181,68,36,495,30.1,0.615,60,1 1,128,98,41,58,32.0,1.321,33,1 8,109,76,39,114,27.9,0.640,31,1 5,139,80,35,160,31.6,0.361,25,1 3,111,62,0,0,22.6,0.142,21,0 9,123,70,44,94,33.1,0.374,40,0 7,159,66,0,0,30.4,0.383,36,1 11,135,0,0,0,52.3,0.578,40,1 8,85,55,20,0,24.4,0.136,42,0 5,158,84,41,210,39.4,0.395,29,1 1,105,58,0,0,24.3,0.187,21,0 3,107,62,13,48,22.9,0.678,23,1 4,109,64,44,99,34.8,0.905,26,1 4,148,60,27,318,30.9,0.150,29,1 0,113,80,16,0,31.0,0.874,21,0 1,138,82,0,0,40.1,0.236,28,0 0,108,68,20,0,27.3,0.787,32,0 2,99,70,16,44,20.4,0.235,27,0 6,103,72,32,190,37.7,0.324,55,0 5,111,72,28,0,23.9,0.407,27,0 8,196,76,29,280,37.5,0.605,57,1 5,162,104,0,0,37.7,0.151,52,1 1,96,64,27,87,33.2,0.289,21,0 7,184,84,33,0,35.5,0.355,41,1 2,81,60,22,0,27.7,0.290,25,0 0,147,85,54,0,42.8,0.375,24,0 7,179,95,31,0,34.2,0.164,60,0 0,140,65,26,130,42.6,0.431,24,1 9,112,82,32,175,34.2,0.260,36,1 12,151,70,40,271,41.8,0.742,38,1 5,109,62,41,129,35.8,0.514,25,1 6,125,68,30,120,30.0,0.464,32,0 5,85,74,22,0,29.0,1.224,32,1 5,112,66,0,0,37.8,0.261,41,1 0,177,60,29,478,34.6,1.072,21,1 2,158,90,0,0,31.6,0.805,66,1 7,119,0,0,0,25.2,0.209,37,0 7,142,60,33,190,28.8,0.687,61,0 1,100,66,15,56,23.6,0.666,26,0 1,87,78,27,32,34.6,0.101,22,0 0,101,76,0,0,35.7,0.198,26,0 3,162,52,38,0,37.2,0.652,24,1 4,197,70,39,744,36.7,2.329,31,0 0,117,80,31,53,45.2,0.089,24,0 4,142,86,0,0,44.0,0.645,22,1 6,134,80,37,370,46.2,0.238,46,1 1,79,80,25,37,25.4,0.583,22,0 4,122,68,0,0,35.0,0.394,29,0 3,74,68,28,45,29.7,0.293,23,0 4,171,72,0,0,43.6,0.479,26,1 7,181,84,21,192,35.9,0.586,51,1 0,179,90,27,0,44.1,0.686,23,1 9,164,84,21,0,30.8,0.831,32,1 0,104,76,0,0,18.4,0.582,27,0 1,91,64,24,0,29.2,0.192,21,0 4,91,70,32,88,33.1,0.446,22,0 3,139,54,0,0,25.6,0.402,22,1 6,119,50,22,176,27.1,1.318,33,1 2,146,76,35,194,38.2,0.329,29,0 9,184,85,15,0,30.0,1.213,49,1 10,122,68,0,0,31.2,0.258,41,0 0,165,90,33,680,52.3,0.427,23,0 9,124,70,33,402,35.4,0.282,34,0 1,111,86,19,0,30.1,0.143,23,0 9,106,52,0,0,31.2,0.380,42,0 2,129,84,0,0,28.0,0.284,27,0 2,90,80,14,55,24.4,0.249,24,0 0,86,68,32,0,35.8,0.238,25,0 12,92,62,7,258,27.6,0.926,44,1 1,113,64,35,0,33.6,0.543,21,1 3,111,56,39,0,30.1,0.557,30,0 2,114,68,22,0,28.7,0.092,25,0 1,193,50,16,375,25.9,0.655,24,0 11,155,76,28,150,33.3,1.353,51,1 3,191,68,15,130,30.9,0.299,34,0 3,141,0,0,0,30.0,0.761,27,1 4,95,70,32,0,32.1,0.612,24,0 3,142,80,15,0,32.4,0.200,63,0 4,123,62,0,0,32.0,0.226,35,1 5,96,74,18,67,33.6,0.997,43,0 0,138,0,0,0,36.3,0.933,25,1 2,128,64,42,0,40.0,1.101,24,0 0,102,52,0,0,25.1,0.078,21,0 2,146,0,0,0,27.5,0.240,28,1 10,101,86,37,0,45.6,1.136,38,1 2,108,62,32,56,25.2,0.128,21,0 3,122,78,0,0,23.0,0.254,40,0 1,71,78,50,45,33.2,0.422,21,0 13,106,70,0,0,34.2,0.251,52,0 2,100,70,52,57,40.5,0.677,25,0 7,106,60,24,0,26.5,0.296,29,1 0,104,64,23,116,27.8,0.454,23,0 5,114,74,0,0,24.9,0.744,57,0 2,108,62,10,278,25.3,0.881,22,0 0,146,70,0,0,37.9,0.334,28,1 10,129,76,28,122,35.9,0.280,39,0 7,133,88,15,155,32.4,0.262,37,0 7,161,86,0,0,30.4,0.165,47,1 2,108,80,0,0,27.0,0.259,52,1 7,136,74,26,135,26.0,0.647,51,0 5,155,84,44,545,38.7,0.619,34,0 1,119,86,39,220,45.6,0.808,29,1 4,96,56,17,49,20.8,0.340,26,0 5,108,72,43,75,36.1,0.263,33,0 0,78,88,29,40,36.9,0.434,21,0 0,107,62,30,74,36.6,0.757,25,1 2,128,78,37,182,43.3,1.224,31,1 1,128,48,45,194,40.5,0.613,24,1 0,161,50,0,0,21.9,0.254,65,0 6,151,62,31,120,35.5,0.692,28,0 2,146,70,38,360,28.0,0.337,29,1 0,126,84,29,215,30.7,0.520,24,0 14,100,78,25,184,36.6,0.412,46,1 8,112,72,0,0,23.6,0.840,58,0 0,167,0,0,0,32.3,0.839,30,1 2,144,58,33,135,31.6,0.422,25,1 5,77,82,41,42,35.8,0.156,35,0 5,115,98,0,0,52.9,0.209,28,1 3,150,76,0,0,21.0,0.207,37,0 2,120,76,37,105,39.7,0.215,29,0 10,161,68,23,132,25.5,0.326,47,1 0,137,68,14,148,24.8,0.143,21,0 0,128,68,19,180,30.5,1.391,25,1 2,124,68,28,205,32.9,0.875,30,1 6,80,66,30,0,26.2,0.313,41,0 0,106,70,37,148,39.4,0.605,22,0 2,155,74,17,96,26.6,0.433,27,1 3,113,50,10,85,29.5,0.626,25,0 7,109,80,31,0,35.9,1.127,43,1 2,112,68,22,94,34.1,0.315,26,0 3,99,80,11,64,19.3,0.284,30,0 3,182,74,0,0,30.5,0.345,29,1 3,115,66,39,140,38.1,0.150,28,0 6,194,78,0,0,23.5,0.129,59,1 4,129,60,12,231,27.5,0.527,31,0 3,112,74,30,0,31.6,0.197,25,1 0,124,70,20,0,27.4,0.254,36,1 13,152,90,33,29,26.8,0.731,43,1 2,112,75,32,0,35.7,0.148,21,0 1,157,72,21,168,25.6,0.123,24,0 1,122,64,32,156,35.1,0.692,30,1 10,179,70,0,0,35.1,0.200,37,0 2,102,86,36,120,45.5,0.127,23,1 6,105,70,32,68,30.8,0.122,37,0 8,118,72,19,0,23.1,1.476,46,0 2,87,58,16,52,32.7,0.166,25,0 1,180,0,0,0,43.3,0.282,41,1 12,106,80,0,0,23.6,0.137,44,0 1,95,60,18,58,23.9,0.260,22,0 0,165,76,43,255,47.9,0.259,26,0 0,117,0,0,0,33.8,0.932,44,0 5,115,76,0,0,31.2,0.343,44,1 9,152,78,34,171,34.2,0.893,33,1 7,178,84,0,0,39.9,0.331,41,1 1,130,70,13,105,25.9,0.472,22,0 1,95,74,21,73,25.9,0.673,36,0 1,0,68,35,0,32.0,0.389,22,0 5,122,86,0,0,34.7,0.290,33,0 8,95,72,0,0,36.8,0.485,57,0 8,126,88,36,108,38.5,0.349,49,0 1,139,46,19,83,28.7,0.654,22,0 3,116,0,0,0,23.5,0.187,23,0 3,99,62,19,74,21.8,0.279,26,0 5,0,80,32,0,41.0,0.346,37,1 4,92,80,0,0,42.2,0.237,29,0 4,137,84,0,0,31.2,0.252,30,0 3,61,82,28,0,34.4,0.243,46,0 1,90,62,12,43,27.2,0.580,24,0 3,90,78,0,0,42.7,0.559,21,0 9,165,88,0,0,30.4,0.302,49,1 1,125,50,40,167,33.3,0.962,28,1 13,129,0,30,0,39.9,0.569,44,1 12,88,74,40,54,35.3,0.378,48,0 1,196,76,36,249,36.5,0.875,29,1 5,189,64,33,325,31.2,0.583,29,1 5,158,70,0,0,29.8,0.207,63,0 5,103,108,37,0,39.2,0.305,65,0 4,146,78,0,0,38.5,0.520,67,1 4,147,74,25,293,34.9,0.385,30,0 5,99,54,28,83,34.0,0.499,30,0 6,124,72,0,0,27.6,0.368,29,1 0,101,64,17,0,21.0,0.252,21,0 3,81,86,16,66,27.5,0.306,22,0 1,133,102,28,140,32.8,0.234,45,1 3,173,82,48,465,38.4,2.137,25,1 0,118,64,23,89,0.0,1.731,21,0 0,84,64,22,66,35.8,0.545,21,0 2,105,58,40,94,34.9,0.225,25,0 2,122,52,43,158,36.2,0.816,28,0 12,140,82,43,325,39.2,0.528,58,1 0,98,82,15,84,25.2,0.299,22,0 1,87,60,37,75,37.2,0.509,22,0 4,156,75,0,0,48.3,0.238,32,1 0,93,100,39,72,43.4,1.021,35,0 1,107,72,30,82,30.8,0.821,24,0 0,105,68,22,0,20.0,0.236,22,0 1,109,60,8,182,25.4,0.947,21,0 1,90,62,18,59,25.1,1.268,25,0 1,125,70,24,110,24.3,0.221,25,0 1,119,54,13,50,22.3,0.205,24,0 5,116,74,29,0,32.3,0.660,35,1 8,105,100,36,0,43.3,0.239,45,1 5,144,82,26,285,32.0,0.452,58,1 3,100,68,23,81,31.6,0.949,28,0 1,100,66,29,196,32.0,0.444,42,0 5,166,76,0,0,45.7,0.340,27,1 1,131,64,14,415,23.7,0.389,21,0 4,116,72,12,87,22.1,0.463,37,0 4,158,78,0,0,32.9,0.803,31,1 2,127,58,24,275,27.7,1.600,25,0 3,96,56,34,115,24.7,0.944,39,0 0,131,66,40,0,34.3,0.196,22,1 3,82,70,0,0,21.1,0.389,25,0 3,193,70,31,0,34.9,0.241,25,1 4,95,64,0,0,32.0,0.161,31,1 6,137,61,0,0,24.2,0.151,55,0 5,136,84,41,88,35.0,0.286,35,1 9,72,78,25,0,31.6,0.280,38,0 5,168,64,0,0,32.9,0.135,41,1 2,123,48,32,165,42.1,0.520,26,0 4,115,72,0,0,28.9,0.376,46,1 0,101,62,0,0,21.9,0.336,25,0 8,197,74,0,0,25.9,1.191,39,1 1,172,68,49,579,42.4,0.702,28,1 6,102,90,39,0,35.7,0.674,28,0 1,112,72,30,176,34.4,0.528,25,0 1,143,84,23,310,42.4,1.076,22,0 1,143,74,22,61,26.2,0.256,21,0 0,138,60,35,167,34.6,0.534,21,1 3,173,84,33,474,35.7,0.258,22,1 1,97,68,21,0,27.2,1.095,22,0 4,144,82,32,0,38.5,0.554,37,1 1,83,68,0,0,18.2,0.624,27,0 3,129,64,29,115,26.4,0.219,28,1 1,119,88,41,170,45.3,0.507,26,0 2,94,68,18,76,26.0,0.561,21,0 0,102,64,46,78,40.6,0.496,21,0 2,115,64,22,0,30.8,0.421,21,0 8,151,78,32,210,42.9,0.516,36,1 4,184,78,39,277,37.0,0.264,31,1 0,94,0,0,0,0.0,0.256,25,0 1,181,64,30,180,34.1,0.328,38,1 0,135,94,46,145,40.6,0.284,26,0 1,95,82,25,180,35.0,0.233,43,1 2,99,0,0,0,22.2,0.108,23,0 3,89,74,16,85,30.4,0.551,38,0 1,80,74,11,60,30.0,0.527,22,0 2,139,75,0,0,25.6,0.167,29,0 1,90,68,8,0,24.5,1.138,36,0 0,141,0,0,0,42.4,0.205,29,1 12,140,85,33,0,37.4,0.244,41,0 5,147,75,0,0,29.9,0.434,28,0 1,97,70,15,0,18.2,0.147,21,0 6,107,88,0,0,36.8,0.727,31,0 0,189,104,25,0,34.3,0.435,41,1 2,83,66,23,50,32.2,0.497,22,0 4,117,64,27,120,33.2,0.230,24,0 8,108,70,0,0,30.5,0.955,33,1 4,117,62,12,0,29.7,0.380,30,1 0,180,78,63,14,59.4,2.420,25,1 1,100,72,12,70,25.3,0.658,28,0 0,95,80,45,92,36.5,0.330,26,0 0,104,64,37,64,33.6,0.510,22,1 0,120,74,18,63,30.5,0.285,26,0 1,82,64,13,95,21.2,0.415,23,0 2,134,70,0,0,28.9,0.542,23,1 0,91,68,32,210,39.9,0.381,25,0 2,119,0,0,0,19.6,0.832,72,0 2,100,54,28,105,37.8,0.498,24,0 14,175,62,30,0,33.6,0.212,38,1 1,135,54,0,0,26.7,0.687,62,0 5,86,68,28,71,30.2,0.364,24,0 10,148,84,48,237,37.6,1.001,51,1 9,134,74,33,60,25.9,0.460,81,0 9,120,72,22,56,20.8,0.733,48,0 1,71,62,0,0,21.8,0.416,26,0 8,74,70,40,49,35.3,0.705,39,0 5,88,78,30,0,27.6,0.258,37,0 10,115,98,0,0,24.0,1.022,34,0 0,124,56,13,105,21.8,0.452,21,0 0,74,52,10,36,27.8,0.269,22,0 0,97,64,36,100,36.8,0.600,25,0 8,120,0,0,0,30.0,0.183,38,1 6,154,78,41,140,46.1,0.571,27,0 1,144,82,40,0,41.3,0.607,28,0 0,137,70,38,0,33.2,0.170,22,0 0,119,66,27,0,38.8,0.259,22,0 7,136,90,0,0,29.9,0.210,50,0 4,114,64,0,0,28.9,0.126,24,0 0,137,84,27,0,27.3,0.231,59,0 2,105,80,45,191,33.7,0.711,29,1 7,114,76,17,110,23.8,0.466,31,0 8,126,74,38,75,25.9,0.162,39,0 4,132,86,31,0,28.0,0.419,63,0 3,158,70,30,328,35.5,0.344,35,1 0,123,88,37,0,35.2,0.197,29,0 4,85,58,22,49,27.8,0.306,28,0 0,84,82,31,125,38.2,0.233,23,0 0,145,0,0,0,44.2,0.630,31,1 0,135,68,42,250,42.3,0.365,24,1 1,139,62,41,480,40.7,0.536,21,0 0,173,78,32,265,46.5,1.159,58,0 4,99,72,17,0,25.6,0.294,28,0 8,194,80,0,0,26.1,0.551,67,0 2,83,65,28,66,36.8,0.629,24,0 2,89,90,30,0,33.5,0.292,42,0 4,99,68,38,0,32.8,0.145,33,0 4,125,70,18,122,28.9,1.144,45,1 3,80,0,0,0,0.0,0.174,22,0 6,166,74,0,0,26.6,0.304,66,0 5,110,68,0,0,26.0,0.292,30,0 2,81,72,15,76,30.1,0.547,25,0 7,195,70,33,145,25.1,0.163,55,1 6,154,74,32,193,29.3,0.839,39,0 2,117,90,19,71,25.2,0.313,21,0 3,84,72,32,0,37.2,0.267,28,0 6,0,68,41,0,39.0,0.727,41,1 7,94,64,25,79,33.3,0.738,41,0 3,96,78,39,0,37.3,0.238,40,0 10,75,82,0,0,33.3,0.263,38,0 0,180,90,26,90,36.5,0.314,35,1 1,130,60,23,170,28.6,0.692,21,0 2,84,50,23,76,30.4,0.968,21,0 8,120,78,0,0,25.0,0.409,64,0 12,84,72,31,0,29.7,0.297,46,1 0,139,62,17,210,22.1,0.207,21,0 9,91,68,0,0,24.2,0.200,58,0 2,91,62,0,0,27.3,0.525,22,0 3,99,54,19,86,25.6,0.154,24,0 3,163,70,18,105,31.6,0.268,28,1 9,145,88,34,165,30.3,0.771,53,1 7,125,86,0,0,37.6,0.304,51,0 13,76,60,0,0,32.8,0.180,41,0 6,129,90,7,326,19.6,0.582,60,0 2,68,70,32,66,25.0,0.187,25,0 3,124,80,33,130,33.2,0.305,26,0 6,114,0,0,0,0.0,0.189,26,0 9,130,70,0,0,34.2,0.652,45,1 3,125,58,0,0,31.6,0.151,24,0 3,87,60,18,0,21.8,0.444,21,0 1,97,64,19,82,18.2,0.299,21,0 3,116,74,15,105,26.3,0.107,24,0 0,117,66,31,188,30.8,0.493,22,0 0,111,65,0,0,24.6,0.660,31,0 2,122,60,18,106,29.8,0.717,22,0 0,107,76,0,0,45.3,0.686,24,0 1,86,66,52,65,41.3,0.917,29,0 6,91,0,0,0,29.8,0.501,31,0 1,77,56,30,56,33.3,1.251,24,0 4,132,0,0,0,32.9,0.302,23,1 0,105,90,0,0,29.6,0.197,46,0 0,57,60,0,0,21.7,0.735,67,0 0,127,80,37,210,36.3,0.804,23,0 3,129,92,49,155,36.4,0.968,32,1 8,100,74,40,215,39.4,0.661,43,1 3,128,72,25,190,32.4,0.549,27,1 10,90,85,32,0,34.9,0.825,56,1 4,84,90,23,56,39.5,0.159,25,0 1,88,78,29,76,32.0,0.365,29,0 8,186,90,35,225,34.5,0.423,37,1 5,187,76,27,207,43.6,1.034,53,1 4,131,68,21,166,33.1,0.160,28,0 1,164,82,43,67,32.8,0.341,50,0 4,189,110,31,0,28.5,0.680,37,0 1,116,70,28,0,27.4,0.204,21,0 3,84,68,30,106,31.9,0.591,25,0 6,114,88,0,0,27.8,0.247,66,0 1,88,62,24,44,29.9,0.422,23,0 1,84,64,23,115,36.9,0.471,28,0 7,124,70,33,215,25.5,0.161,37,0 1,97,70,40,0,38.1,0.218,30,0 8,110,76,0,0,27.8,0.237,58,0 11,103,68,40,0,46.2,0.126,42,0 11,85,74,0,0,30.1,0.300,35,0 6,125,76,0,0,33.8,0.121,54,1 0,198,66,32,274,41.3,0.502,28,1 1,87,68,34,77,37.6,0.401,24,0 6,99,60,19,54,26.9,0.497,32,0 0,91,80,0,0,32.4,0.601,27,0 2,95,54,14,88,26.1,0.748,22,0 1,99,72,30,18,38.6,0.412,21,0 6,92,62,32,126,32.0,0.085,46,0 4,154,72,29,126,31.3,0.338,37,0 0,121,66,30,165,34.3,0.203,33,1 3,78,70,0,0,32.5,0.270,39,0 2,130,96,0,0,22.6,0.268,21,0 3,111,58,31,44,29.5,0.430,22,0 2,98,60,17,120,34.7,0.198,22,0 1,143,86,30,330,30.1,0.892,23,0 1,119,44,47,63,35.5,0.280,25,0 6,108,44,20,130,24.0,0.813,35,0 2,118,80,0,0,42.9,0.693,21,1 10,133,68,0,0,27.0,0.245,36,0 2,197,70,99,0,34.7,0.575,62,1 0,151,90,46,0,42.1,0.371,21,1 6,109,60,27,0,25.0,0.206,27,0 12,121,78,17,0,26.5,0.259,62,0 8,100,76,0,0,38.7,0.190,42,0 8,124,76,24,600,28.7,0.687,52,1 1,93,56,11,0,22.5,0.417,22,0 8,143,66,0,0,34.9,0.129,41,1 6,103,66,0,0,24.3,0.249,29,0 3,176,86,27,156,33.3,1.154,52,1 0,73,0,0,0,21.1,0.342,25,0 11,111,84,40,0,46.8,0.925,45,1 2,112,78,50,140,39.4,0.175,24,0 3,132,80,0,0,34.4,0.402,44,1 2,82,52,22,115,28.5,1.699,25,0 6,123,72,45,230,33.6,0.733,34,0 0,188,82,14,185,32.0,0.682,22,1 0,67,76,0,0,45.3,0.194,46,0 1,89,24,19,25,27.8,0.559,21,0 1,173,74,0,0,36.8,0.088,38,1 1,109,38,18,120,23.1,0.407,26,0 1,108,88,19,0,27.1,0.400,24,0 6,96,0,0,0,23.7,0.190,28,0 1,124,74,36,0,27.8,0.100,30,0 7,150,78,29,126,35.2,0.692,54,1 4,183,0,0,0,28.4,0.212,36,1 1,124,60,32,0,35.8,0.514,21,0 1,181,78,42,293,40.0,1.258,22,1 1,92,62,25,41,19.5,0.482,25,0 0,152,82,39,272,41.5,0.270,27,0 1,111,62,13,182,24.0,0.138,23,0 3,106,54,21,158,30.9,0.292,24,0 3,174,58,22,194,32.9,0.593,36,1 7,168,88,42,321,38.2,0.787,40,1 6,105,80,28,0,32.5,0.878,26,0 11,138,74,26,144,36.1,0.557,50,1 3,106,72,0,0,25.8,0.207,27,0 6,117,96,0,0,28.7,0.157,30,0 2,68,62,13,15,20.1,0.257,23,0 9,112,82,24,0,28.2,1.282,50,1 0,119,0,0,0,32.4,0.141,24,1 2,112,86,42,160,38.4,0.246,28,0 2,92,76,20,0,24.2,1.698,28,0 6,183,94,0,0,40.8,1.461,45,0 0,94,70,27,115,43.5,0.347,21,0 2,108,64,0,0,30.8,0.158,21,0 4,90,88,47,54,37.7,0.362,29,0 0,125,68,0,0,24.7,0.206,21,0 0,132,78,0,0,32.4,0.393,21,0 5,128,80,0,0,34.6,0.144,45,0 4,94,65,22,0,24.7,0.148,21,0 7,114,64,0,0,27.4,0.732,34,1 0,102,78,40,90,34.5,0.238,24,0 2,111,60,0,0,26.2,0.343,23,0 1,128,82,17,183,27.5,0.115,22,0 10,92,62,0,0,25.9,0.167,31,0 13,104,72,0,0,31.2,0.465,38,1 5,104,74,0,0,28.8,0.153,48,0 2,94,76,18,66,31.6,0.649,23,0 7,97,76,32,91,40.9,0.871,32,1 1,100,74,12,46,19.5,0.149,28,0 0,102,86,17,105,29.3,0.695,27,0 4,128,70,0,0,34.3,0.303,24,0 6,147,80,0,0,29.5,0.178,50,1 4,90,0,0,0,28.0,0.610,31,0 3,103,72,30,152,27.6,0.730,27,0 2,157,74,35,440,39.4,0.134,30,0 1,167,74,17,144,23.4,0.447,33,1 0,179,50,36,159,37.8,0.455,22,1 11,136,84,35,130,28.3,0.260,42,1 0,107,60,25,0,26.4,0.133,23,0 1,91,54,25,100,25.2,0.234,23,0 1,117,60,23,106,33.8,0.466,27,0 5,123,74,40,77,34.1,0.269,28,0 2,120,54,0,0,26.8,0.455,27,0 1,106,70,28,135,34.2,0.142,22,0 2,155,52,27,540,38.7,0.240,25,1 2,101,58,35,90,21.8,0.155,22,0 1,120,80,48,200,38.9,1.162,41,0 11,127,106,0,0,39.0,0.190,51,0 3,80,82,31,70,34.2,1.292,27,1 10,162,84,0,0,27.7,0.182,54,0 1,199,76,43,0,42.9,1.394,22,1 8,167,106,46,231,37.6,0.165,43,1 9,145,80,46,130,37.9,0.637,40,1 6,115,60,39,0,33.7,0.245,40,1 1,112,80,45,132,34.8,0.217,24,0 4,145,82,18,0,32.5,0.235,70,1 10,111,70,27,0,27.5,0.141,40,1 6,98,58,33,190,34.0,0.430,43,0 9,154,78,30,100,30.9,0.164,45,0 6,165,68,26,168,33.6,0.631,49,0 1,99,58,10,0,25.4,0.551,21,0 10,68,106,23,49,35.5,0.285,47,0 3,123,100,35,240,57.3,0.880,22,0 8,91,82,0,0,35.6,0.587,68,0 6,195,70,0,0,30.9,0.328,31,1 9,156,86,0,0,24.8,0.230,53,1 0,93,60,0,0,35.3,0.263,25,0 3,121,52,0,0,36.0,0.127,25,1 2,101,58,17,265,24.2,0.614,23,0 2,56,56,28,45,24.2,0.332,22,0 0,162,76,36,0,49.6,0.364,26,1 0,95,64,39,105,44.6,0.366,22,0 4,125,80,0,0,32.3,0.536,27,1 5,136,82,0,0,0.0,0.640,69,0 2,129,74,26,205,33.2,0.591,25,0 3,130,64,0,0,23.1,0.314,22,0 1,107,50,19,0,28.3,0.181,29,0 1,140,74,26,180,24.1,0.828,23,0 1,144,82,46,180,46.1,0.335,46,1 8,107,80,0,0,24.6,0.856,34,0 13,158,114,0,0,42.3,0.257,44,1 2,121,70,32,95,39.1,0.886,23,0 7,129,68,49,125,38.5,0.439,43,1 2,90,60,0,0,23.5,0.191,25,0 7,142,90,24,480,30.4,0.128,43,1 3,169,74,19,125,29.9,0.268,31,1 0,99,0,0,0,25.0,0.253,22,0 4,127,88,11,155,34.5,0.598,28,0 4,118,70,0,0,44.5,0.904,26,0 2,122,76,27,200,35.9,0.483,26,0 6,125,78,31,0,27.6,0.565,49,1 1,168,88,29,0,35.0,0.905,52,1 2,129,0,0,0,38.5,0.304,41,0 4,110,76,20,100,28.4,0.118,27,0 6,80,80,36,0,39.8,0.177,28,0 10,115,0,0,0,0.0,0.261,30,1 2,127,46,21,335,34.4,0.176,22,0 9,164,78,0,0,32.8,0.148,45,1 2,93,64,32,160,38.0,0.674,23,1 3,158,64,13,387,31.2,0.295,24,0 5,126,78,27,22,29.6,0.439,40,0 10,129,62,36,0,41.2,0.441,38,1 0,134,58,20,291,26.4,0.352,21,0 3,102,74,0,0,29.5,0.121,32,0 7,187,50,33,392,33.9,0.826,34,1 3,173,78,39,185,33.8,0.970,31,1 10,94,72,18,0,23.1,0.595,56,0 1,108,60,46,178,35.5,0.415,24,0 5,97,76,27,0,35.6,0.378,52,1 4,83,86,19,0,29.3,0.317,34,0 1,114,66,36,200,38.1,0.289,21,0 1,149,68,29,127,29.3,0.349,42,1 5,117,86,30,105,39.1,0.251,42,0 1,111,94,0,0,32.8,0.265,45,0 4,112,78,40,0,39.4,0.236,38,0 1,116,78,29,180,36.1,0.496,25,0 0,141,84,26,0,32.4,0.433,22,0 2,175,88,0,0,22.9,0.326,22,0 2,92,52,0,0,30.1,0.141,22,0 3,130,78,23,79,28.4,0.323,34,1 8,120,86,0,0,28.4,0.259,22,1 2,174,88,37,120,44.5,0.646,24,1 2,106,56,27,165,29.0,0.426,22,0 2,105,75,0,0,23.3,0.560,53,0 4,95,60,32,0,35.4,0.284,28,0 0,126,86,27,120,27.4,0.515,21,0 8,65,72,23,0,32.0,0.600,42,0 2,99,60,17,160,36.6,0.453,21,0 1,102,74,0,0,39.5,0.293,42,1 11,120,80,37,150,42.3,0.785,48,1 3,102,44,20,94,30.8,0.400,26,0 1,109,58,18,116,28.5,0.219,22,0 9,140,94,0,0,32.7,0.734,45,1 13,153,88,37,140,40.6,1.174,39,0 12,100,84,33,105,30.0,0.488,46,0 1,147,94,41,0,49.3,0.358,27,1 1,81,74,41,57,46.3,1.096,32,0 3,187,70,22,200,36.4,0.408,36,1 6,162,62,0,0,24.3,0.178,50,1 4,136,70,0,0,31.2,1.182,22,1 1,121,78,39,74,39.0,0.261,28,0 3,108,62,24,0,26.0,0.223,25,0 0,181,88,44,510,43.3,0.222,26,1 8,154,78,32,0,32.4,0.443,45,1 1,128,88,39,110,36.5,1.057,37,1 7,137,90,41,0,32.0,0.391,39,0 0,123,72,0,0,36.3,0.258,52,1 1,106,76,0,0,37.5,0.197,26,0 6,190,92,0,0,35.5,0.278,66,1 2,88,58,26,16,28.4,0.766,22,0 9,170,74,31,0,44.0,0.403,43,1 9,89,62,0,0,22.5,0.142,33,0 10,101,76,48,180,32.9,0.171,63,0 2,122,70,27,0,36.8,0.340,27,0 5,121,72,23,112,26.2,0.245,30,0 1,126,60,0,0,30.1,0.349,47,1 1,93,70,31,0,30.4,0.315,23,0



Friday, May 5, 2017

How to translate QuantLib C++ code to Python

Here is a demo of how QuantLib c++ code are translated to Python. This included the code for importing of csv file and construction of volatility surface and the timing of MCDiscreteArithmeticAPEngine. Slicing and manipulation of list/array is much easier in Python than that of C++ code. However, C++ is faster.

QuantLib C++ source code, AsianOption.cpp

AsianOption.cpp    Select all
// g++ -std=c++11 AsianOption.cpp -o AsianOption -lQuantLib #include <ql/quantlib.hpp> #include <boost/timer.hpp> #include <iostream> #include <iomanip> #include <fstream> #include <string> #include <boost/algorithm/string/split.hpp> #include <boost/algorithm/string/classification.hpp> #include <boost/lexical_cast.hpp> using namespace QuantLib; using namespace std; void GetStrikes(string &path, vector<Real> &strikes) { std::ifstream file(path); std::string line; std::vector<std::string> tokens; int linecount = 0; while (std::getline(file, line)) { std::stringstream stringStream(line); std::string content; int item = 0; if (linecount >= 1) while (std::getline(stringStream, content, ',')) { switch (item) { case 0: // strikes are on first column only strikes.push_back(boost::lexical_cast<double>(content)); break; default: break; } item++; } linecount++; } return; /* if hardcode csv data strikes.push_back(0.67263); strikes.push_back(0.71865); strikes.push_back(0.7487); strikes.push_back(0.77129); strikes.push_back(0.78984); strikes.push_back(0.80587); strikes.push_back(0.82034); strikes.push_back(0.83379); strikes.push_back(0.84658); strikes.push_back(0.85792); strikes.push_back(0.87354); strikes.push_back(0.89085); strikes.push_back(0.90904); strikes.push_back(0.9289); strikes.push_back(0.95157); strikes.push_back(0.97862); strikes.push_back(1.01337); strikes.push_back(1.06261); strikes.push_back(1.14631); */ } void GetExpiryDates(string &path, vector<Date> &expirations) { std::ifstream file(path); std::string line; std::vector<std::string> tokens; int linecount = 0; while (std::getline(file, line)) { std::stringstream stringStream(line); std::string content; int item = 0; if (linecount == 0) // expiration dates are on first row only while (std::getline(stringStream, content, ',')) { switch (item) { case 1: case 2: case 3: case 4: boost::algorithm::split(tokens, content, boost::algorithm::is_any_of("/")); expirations.push_back(Date(Day(boost::lexical_cast<int>(tokens.at(1))), Month(boost::lexical_cast<int>(tokens.at(0))), Year(boost::lexical_cast<int>(tokens.at(2))))); break; default: break; } item++; } linecount++; } return; /* if hardcode csv data expirations.push_back(Date(27, June, 2017)); expirations.push_back(Date(27, September, 2017)); expirations.push_back(Date(27, December, 2017)); expirations.push_back(Date(27, June, 2018)); */ } Matrix GetVolData(string &path, vector<Date> &expirations, vector<Real> &strikes) { // Matrix volMatrix(19, 4); Matrix volMatrix(strikes.size(), expirations.size()); std::ifstream file(path); std::string line; std::vector<std::string> tokens; int linecount = 0; while (std::getline(file, line)) { std::stringstream stringStream(line); std::string content; int item = 0; if (linecount >= 1) // vols are on second row onward while (std::getline(stringStream, content, ',')) { switch (item) { case 1: case 2: case 3: case 4: volMatrix[linecount-1][item-1] = boost::lexical_cast<double>(content); // vols are on second column onward break; default: break; } item++; } linecount++; } return volMatrix; /* if hardcode csv data //0.67263,0.144183920277296,0.139374695503699,0.135526204819277,0.12885 volMatrix[0][0] = 0.144183920277296; volMatrix[0][1] = 0.139374695503699; volMatrix[0][2] = 0.135526204819277; volMatrix[0][3] = 0.12885; //0.71865,0.133703802426343,0.129893056346044,0.126909006024096,0.12175 volMatrix[1][0] = 0.133703802426343; volMatrix[1][1] = 0.129893056346044; volMatrix[1][2] = 0.126909006024096; volMatrix[1][3] = 0.12175; //0.7487,0.126860526863085,0.123701764371087,0.121209416342412,0.11695 volMatrix[2][0] = 0.126860526863085; volMatrix[2][1] = 0.123701764371087; volMatrix[2][2] = 0.121209416342412; volMatrix[2][3] = 0.11695; // 0.77129,0.121720863192182,0.118881707209199,0.116979766476388,0.11381 volMatrix[3][0] = 0.121720863192182; volMatrix[3][1] = 0.118881707209199; volMatrix[3][2] = 0.116979766476388; volMatrix[3][3] = 0.11381; // 0.78984,0.117581136690647,0.115218428824572,0.113899219047619,0.11163 volMatrix[4][0] = 0.117581136690647; volMatrix[4][1] = 0.115218428824572; volMatrix[4][2] = 0.113899219047619; volMatrix[4][3] = 0.11163; // 0.80587,0.114363421052632,0.112523118729097,0.111637193240265,0.11019 volMatrix[5][0] = 0.114363421052632; volMatrix[5][1] = 0.112523118729097; volMatrix[5][2] = 0.111637193240265; volMatrix[5][3] = 0.11019; //0.82034,0.111728795180723,0.110489402985075,0.109987692307692,0.10921 volMatrix[6][0] = 0.111728795180723; volMatrix[6][1] = 0.110489402985075; volMatrix[6][2] = 0.109987692307692; volMatrix[6][3] = 0.10921; // 0.83379,0.109805703883495,0.109100413723512,0.108870460584588,0.1086 volMatrix[7][0] = 0.109805703883495; volMatrix[7][1] = 0.109100413723512; volMatrix[7][2] = 0.108870460584588; volMatrix[7][3] = 0.1086; // 0.84658,0.108581646586345,0.108250493273543,0.108197213114754,0.10829 volMatrix[8][0] = 0.108581646586345; volMatrix[8][1] = 0.108250493273543; volMatrix[8][2] = 0.108197213114754; volMatrix[8][3] = 0.10829; // 0.85792,0.108190964125561,0.107986172506739,0.10796631037213,0.10822 volMatrix[9][0] = 0.108190964125561; volMatrix[9][1] = 0.107986172506739; volMatrix[9][2] = 0.10796631037213; volMatrix[9][3] = 0.10822; // 0.87354,0.10859510460251,0.108310304612707,0.108232350773766,0.10849 volMatrix[10][0] = 0.10859510460251; volMatrix[10][1] = 0.108310304612707; volMatrix[10][2] = 0.108232350773766; volMatrix[10][3] = 0.10849; // 0.89085,0.110043016488846,0.109404567049808,0.109102906403941,0.10919 volMatrix[11][0] = 0.110043016488846; volMatrix[11][1] = 0.109404567049808; volMatrix[11][2] = 0.109102906403941; volMatrix[11][3] = 0.10919; // 0.90904,0.112447321958457,0.111343238289206,0.110615417475728,0.11036 volMatrix[12][0] = 0.112447321958457; volMatrix[12][1] = 0.111343238289206; volMatrix[12][2] = 0.110615417475728; volMatrix[12][3] = 0.11036; // 0.9289,0.115567066189624,0.113888152866242,0.112830993150685,0.11201 volMatrix[13][0] = 0.115567066189624; volMatrix[13][1] = 0.113888152866242; volMatrix[13][2] = 0.112830993150685; volMatrix[13][3] = 0.11201; // 0.95157,0.119454321849106,0.117151688909342,0.115569047072331,0.11433 volMatrix[14][0] = 0.119454321849106; volMatrix[14][1] = 0.117151688909342; volMatrix[14][2] = 0.115569047072331; volMatrix[14][3] = 0.11433; // 0.97862,0.123858310308183,0.121275916334661,0.119199029605263,0.11731 volMatrix[15][0] = 0.123858310308183; volMatrix[15][1] = 0.121275916334661; volMatrix[15][2] = 0.119199029605263; volMatrix[15][3] = 0.11731; // 1.01337,0.129434558979809,0.126231870274572,0.123929902439024,0.12145 volMatrix[16][0] = 0.129434558979809; volMatrix[16][1] = 0.126231870274572; volMatrix[16][2] = 0.123929902439024; volMatrix[16][3] = 0.12145; // 1.06261,0.137335982996812,0.133099606048548,0.12994278699187,0.12723 volMatrix[17][0] = 0.137335982996812; volMatrix[17][1] = 0.133099606048548; volMatrix[17][2] = 0.12994278699187; volMatrix[17][3] = 0.12723; // 1.14631,0.150767120085016,0.144773641066454,0.140163713821138,0.13547 volMatrix[18][0] = 0.150767120085016; volMatrix[18][1] = 0.144773641066454; volMatrix[18][2] = 0.140163713821138; volMatrix[18][3] = 0.13547; return volMatrix; */ } void asian() { // Calendar set up Calendar calendar = TARGET(); Date todaysDate(4, April, 2017); Settings::instance().evaluationDate() = todaysDate; DayCounter dayCounter = Actual360(); // Option parameters Asian FX Option::Type optionType(Option::Call); Average::Type averageType = Average::Arithmetic; Date maturity(4, April, 2018); Real strike = 0.74; Volatility volatility = 0.07053702474; Date obsStart(4, March, 2018); Real runningSum = 0; Size pastFixings = 0; vector<Date> fixingDates; for (Date incrementedDate = obsStart; incrementedDate <= maturity; incrementedDate += 1) { if (calendar.isBusinessDay(incrementedDate)) { fixingDates.push_back(incrementedDate); } } // Option parameters // European Exercise boost::shared_ptr<Exercise> europeanExercise( new EuropeanExercise(maturity)); // Payoff boost::shared_ptr<StrikedTypePayoff> payoffAsianOption( new PlainVanillaPayoff(Option::Type(optionType), strike)); // Model parameters Real underlying = 0.748571186; Spread dividendYield = 0.04125; Rate riskFreeRate = 0.0225377; // Market Data // Quote handling Handle<Quote> underlyingH( boost::shared_ptr<Quote>(new SimpleQuote(underlying))); // Yield term structure handling Handle<YieldTermStructure> flatTermStructure( boost::shared_ptr<YieldTermStructure>(new FlatForward(todaysDate, dividendYield, dayCounter))); // Dividend term structure handling Handle<YieldTermStructure> flatDividendTermStructure( boost::shared_ptr<YieldTermStructure>(new FlatForward(todaysDate, riskFreeRate, dayCounter))); // Volatility structure handling: constant volatility Handle<BlackVolTermStructure> flatVolTermStructure( boost::shared_ptr<BlackVolTermStructure>(new BlackConstantVol(todaysDate, calendar, volatility, dayCounter))); // Read csv file string path = "./VolMatrixA.csv"; vector<Real> strikes = {}; GetStrikes(path, strikes); vector<Date> expirations = {}; GetExpiryDates(path, expirations); cout << "strikes.size() " << strikes.size() << endl; cout << "expirations.size() " << expirations.size() << endl; // assert csv data BOOST_ASSERT_MSG(strikes.size() > 0, static_cast<std::stringstream&>(std::stringstream() << "No valid strikes.size() found! It is " << strikes.size()).str().c_str()); BOOST_ASSERT_MSG(expirations.size() > 0, static_cast<std::stringstream&>(std::stringstream() << "No valid expirations.size() found! It is " << expirations.size()).str().c_str()); Matrix volMatrix = GetVolData(path, expirations, strikes); // Volatility Surface BlackVarianceSurface volatilitySurface(Settings::instance().evaluationDate(), calendar, expirations, strikes, volMatrix, dayCounter); volatilitySurface.setInterpolation<Bicubic>(); volatilitySurface.enableExtrapolation(true); const boost::shared_ptr<BlackVarianceSurface> volatilitySurfaceH( new BlackVarianceSurface(volatilitySurface)); Handle<BlackVolTermStructure> volTermStructure(volatilitySurfaceH); // the BS equation behind boost::shared_ptr<BlackScholesMertonProcess> bsmProcess( new BlackScholesMertonProcess(underlyingH, flatDividendTermStructure, flatTermStructure, volTermStructure)); // Options DiscreteAveragingAsianOption discreteArithmeticAsianAverageOption( averageType, runningSum, pastFixings, fixingDates, payoffAsianOption, europeanExercise); // Outputting on the screen cout << "Option type = " << optionType << endl; cout << "Option maturity = " << maturity << endl; cout << "Underlying = " << underlying << endl; cout << "Strike = " << strike << endl; cout << "Risk-free interest rate = " << setprecision(4) << io::rate(riskFreeRate) << endl; cout << "Dividend yield = " << setprecision(4) << io::rate(dividendYield) << endl; cout << "Volatility = " << setprecision(4) << io::volatility(volatility) << endl; cout << "Time-length between successive fixings = weekly time step" << endl; cout << "Previous fixings = " << pastFixings << endl; cout << setprecision(10) << endl; boost::timer timer; // Pricing engine discreteArithmeticAsianAverageOption.setPricingEngine( boost::shared_ptr<PricingEngine>( MakeMCDiscreteArithmeticAPEngine<LowDiscrepancy>(bsmProcess) .withSamples(1500))); // Timer timer.restart(); try { cout << "Discrete ArithMC Price: " << discreteArithmeticAsianAverageOption.NPV() << endl; } catch (exception const& e) { cout << "Erreur: " << e.what() << endl; } cout << " in " << timer.elapsed() << " s" << endl; timer.restart(); } int main(int, char* []) { asian(); }


QuantLib Python source code, AsianOption.py

AsianOption.py    Select all
#!python2 #!/usr/bin/env python # AsianOption.py from QuantLib import * import csv import time # Calendar set up calendar = TARGET() todaysDate = Date(4, April, 2017) Settings.instance().evaluationDate = todaysDate dayCounter = Actual360() # Option parameters Asian FX optionType = Option.Call averageType = Average.Arithmetic maturity = Date(4, April, 2018) strike = 0.74 volatility = 0.07053702474 obsStart = Date(4, March, 2018) runningSum = 0 pastFixings = 0 fixingDates = [ Date(serial) for serial in range(obsStart.serialNumber(), maturity.serialNumber()) if calendar.isBusinessDay(Date(serial)) ] # Model parameters underlying = 0.748571186 dividendYield = 0.04125 riskFreeRate = 0.0225377 #settlementDate = todaysDate # Option parameters # European Exercise europeanExercise = EuropeanExercise(maturity) # Payoff payoffAsianOption = PlainVanillaPayoff(optionType, strike) # Market Data # Quote handling underlyingH = QuoteHandle(SimpleQuote(underlying)) # Yield term structure handling flatTermStructure = YieldTermStructureHandle(FlatForward(todaysDate, dividendYield, dayCounter)) # Dividend term structure handling flatDividendTermStructure = YieldTermStructureHandle(FlatForward(todaysDate, riskFreeRate, dayCounter)) # Volatility structure handling: constant volatility flatVolTermStructure = BlackVolTermStructureHandle(BlackConstantVol(Settings.instance().evaluationDate, calendar, volatility, dayCounter)) # Read csv file with open('VolMatrixA.csv', 'rb') as f: reader = csv.reader(f) csv_list = list(reader) expirations = [ Date(int(col.split("/")[1]), int(col.split("/")[0]), int(col.split("/")[2])) for col in csv_list[0][1:] ] # expirations are on first row[0] and for second column[1:] onward strikes = [ float(row[0]) for row in csv_list[1:] ] # strikes are for second row [1:] onward and on first column [0] """ # if hardcode csv data expirations = [Date(27, June, 2017), Date(27, September, 2017), Date(27, December, 2017), Date(27, June, 2018)] strikes = [0.67263, 0.71865, 0.7487, 0.77129, 0.78984, 0.80587, 0.82034, 0.83379, 0.84658, 0.85792, 0.87354, 0.89085, 0.90904, 0.9289, 0.95157, 0.97862, 1.01337, 1.06261, 1.14631] volMatrix = [ [0.144183920277296,0.139374695503699,0.135526204819277,0.12885], [0.133703802426343,0.129893056346044,0.126909006024096,0.12175], [0.126860526863085,0.123701764371087,0.121209416342412,0.11695], [0.121720863192182,0.118881707209199,0.116979766476388,0.11381], [0.117581136690647,0.115218428824572,0.113899219047619,0.11163], [0.114363421052632,0.112523118729097,0.111637193240265,0.11019], [0.111728795180723,0.110489402985075,0.109987692307692,0.10921], [0.109805703883495,0.109100413723512,0.108870460584588,0.1086], [0.108581646586345,0.108250493273543,0.108197213114754,0.10829], [0.108190964125561,0.107986172506739,0.10796631037213,0.10822], [0.10859510460251,0.108310304612707,0.108232350773766,0.10849], [0.110043016488846,0.109404567049808,0.109102906403941,0.10919], [0.112447321958457,0.111343238289206,0.110615417475728,0.11036], [0.115567066189624,0.113888152866242,0.112830993150685,0.11201], [0.119454321849106,0.117151688909342,0.115569047072331,0.11433], [0.123858310308183,0.121275916334661,0.119199029605263,0.11731], [0.129434558979809,0.126231870274572,0.123929902439024,0.12145], [0.137335982996812,0.133099606048548,0.12994278699187,0.12723], [0.150767120085016,0.144773641066454,0.140163713821138,0.13547] ] """ # assert csv data assert len(strikes) > 0, "No valid len(strikes) found ! It is " + str(len(strikes)) assert len(expirations) > 0, "No valid len(expirations) found ! It is " + str(len(expirations)) #volMatrix = Matrix(len(strikes), len(expirations)) volMatrix = [[float(y) for y in x[1:]] for x in csv_list[1:]] # vols are for second row [1:] and for second column [1:] onward print "len(strikes) ", len(strikes) print "len(expirations) ", len(expirations) # Volatility Surface volatilitySurface = BlackVarianceSurface(Settings.instance().evaluationDate, calendar, expirations, strikes, volMatrix, dayCounter) volatilitySurface.setInterpolation("Bicubic") volatilitySurface.enableExtrapolation() volTermStructure = BlackVolTermStructureHandle(volatilitySurface) # the BS equation behind bsmProcess = BlackScholesMertonProcess(underlyingH, flatDividendTermStructure, flatTermStructure, volTermStructure) # Options discreteArithmeticAsianAverageOption = DiscreteAveragingAsianOption(averageType, runningSum, pastFixings, fixingDates, payoffAsianOption, europeanExercise) # Outputting on the screen optionTypeKeyName = dict((v,k) for k, v in vars(Option).iteritems() if v == optionType) print "Option type = ", optionTypeKeyName[optionType] print "Option maturity = ", maturity print "Underlying = ", underlying print "Strike = ", strike print "Risk-free interest rate = ", '{0:.{prec}f}%'.format(riskFreeRate*100.00, prec=4) print "Dividend yield = ", '{0:.{prec}f}%'.format(dividendYield*100.00, prec=4) print "Volatility = ", '{0:.{prec}f}%'.format(volatility*100.00, prec=4) print "Time-length between successive fixings = weekly time step" print "Previous fixings = ", pastFixings print "" # Pricing engine engine = MCDiscreteArithmeticAPEngine(bsmProcess, "LowDiscrepancy", requiredSamples=1500) discreteArithmeticAsianAverageOption.setPricingEngine(engine) # Timer start = time.time() print "Discrete ArithMC Price: ", discreteArithmeticAsianAverageOption.NPV() print ' in ' + '{0:.2f}'.format(time.time() - start), 's\n'


VolMatrixA.csv file is the data source of the Volatility Surface

VolMatrixA.csv    Select all
Strike ,6/27/2017,9/27/2017,12/27/2017,6/27/2018 0.67263,0.144183920277296,0.139374695503699,0.135526204819277,0.12885 0.71865,0.133703802426343,0.129893056346044,0.126909006024096,0.12175 0.7487,0.126860526863085,0.123701764371087,0.121209416342412,0.11695 0.77129,0.121720863192182,0.118881707209199,0.116979766476388,0.11381 0.78984,0.117581136690647,0.115218428824572,0.113899219047619,0.11163 0.80587,0.114363421052632,0.112523118729097,0.111637193240265,0.11019 0.82034,0.111728795180723,0.110489402985075,0.109987692307692,0.10921 0.83379,0.109805703883495,0.109100413723512,0.108870460584588,0.1086 0.84658,0.108581646586345,0.108250493273543,0.108197213114754,0.10829 0.85792,0.108190964125561,0.107986172506739,0.10796631037213,0.10822 0.87354,0.10859510460251,0.108310304612707,0.108232350773766,0.10849 0.89085,0.110043016488846,0.109404567049808,0.109102906403941,0.10919 0.90904,0.112447321958457,0.111343238289206,0.110615417475728,0.11036 0.9289,0.115567066189624,0.113888152866242,0.112830993150685,0.11201 0.95157,0.119454321849106,0.117151688909342,0.115569047072331,0.11433 0.97862,0.123858310308183,0.121275916334661,0.119199029605263,0.11731 1.01337,0.129434558979809,0.126231870274572,0.123929902439024,0.12145 1.06261,0.137335982996812,0.133099606048548,0.12994278699187,0.12723 1.14631,0.150767120085016,0.144773641066454,0.140163713821138,0.13547


QuantLib Python source code, Gaussian1dModels.py (Gaussian1dModels.cpp in QuantLib Examples Folder)

Gaussian1dModels.py    Select all
#!python2 #!/usr/bin/env python #Gaussian1dModels.py import QuantLib as ql def printBasket(basket): print ("%-20s %-20s %-20s %-20s %-20s %-20s" % ("Expiry", "Maturity", "Nominal", "Rate", "MarketVol", "Pay/Rec")) print ("==================================================================================================================") for i in range(0, len(basket)): expiryDate = basket[i].swaptionExpiryDate() endDate = basket[i].swaptionMaturityDate() nominal = basket[i].swaptionNominal() vol = basket[i].volatility().value() rate = basket[i].swaptionStrike() #type = basket[i].swaption.type() print ("%-20s %-20s %-20f %-20f %-20f" % (str(expiryDate), str(endDate), nominal, rate, vol)) print("==================================================================================================================") def printModelCalibration(basket, volatility): print ("%-20s %-20s %-20s %-20s %-20s %-20s" % ("Expiry","Model sigma","ModelPrice","MarketPrice","Model impVol","Market impVol")) print ("=================================================================================================================") for i in range(0, len(basket)): expiryDate = basket[i].swaptionExpiryDate() modelValue = basket[i].modelValue() marketValue= basket[i].marketValue() impVol = basket[i].impliedVolatility(modelValue, 1e-6, 1000, 0.0, 2.0) vol = basket[i].volatility().value() print ("%-20s %-20f %-20f %-20f %-20f %-20f" % (str(expiryDate), volatility[i], modelValue, marketValue, impVol, vol)) print("==================================================================================================================") refDate = ql.Date(30, 4, 2014) # Date refDate(30, April, 2014); ql.Settings.instance().setEvaluationDate(refDate) # Settings::instance().evaluationDate() = refDate; forward6mQuote = ql.QuoteHandle(ql.SimpleQuote(0.025)) # Handle<Quote> forward6mQuote(boost::make_shared(0.025)); oisQuote = ql.QuoteHandle(ql.SimpleQuote(0.02)) # Handle<Quote> oisQuote(boost::make_shared(0.02)); volQuote = ql.QuoteHandle(ql.SimpleQuote(0.2)) # Handle<Quote> volQuote(boost::make_shared(0.2)); dc = ql.Actual365Fixed() yts6m = ql.FlatForward(refDate, forward6mQuote, dc) ytsOis= ql.FlatForward(refDate, oisQuote, dc) yts6m.enableExtrapolation() ytsOis.enableExtrapolation() hyts6m = ql.RelinkableYieldTermStructureHandle(yts6m) t0_curve = ql.YieldTermStructureHandle(yts6m) t0_Ois = ql.YieldTermStructureHandle(ytsOis) euribor6m = ql.Euribor6M(hyts6m) swaptionVol = ql.ConstantSwaptionVolatility(0, ql.TARGET(), ql.ModifiedFollowing, volQuote, ql.Actual365Fixed()) # Handle<SwaptionVolatilityStructure> swaptionVol(boost::make_shared(0, TARGET(), ModifiedFollowing, volQuote, Actual365Fixed())); effectiveDate = ql.TARGET().advance(refDate, ql.Period('2D')) # Date effectiveDate = TARGET().advance(refDate, 2 * Days); maturityDate = ql.TARGET().advance(effectiveDate, ql.Period('10Y')) # Date maturityDate = TARGET().advance(effectiveDate, 10 * Years); fixedSchedule = ql.Schedule(effectiveDate, maturityDate, ql.Period('1Y'), ql.TARGET(), ql.ModifiedFollowing, ql.ModifiedFollowing, ql.DateGeneration.Forward, False) # Schedule fixedSchedule(effectiveDate, maturityDate, 1 * Years, TARGET(), ModifiedFollowing, ModifiedFollowing, DateGeneration::Forward, false); floatSchedule = ql.Schedule(effectiveDate, maturityDate, ql.Period('6M'), ql.TARGET(), ql.ModifiedFollowing, ql.ModifiedFollowing, ql.DateGeneration.Forward, False) # Schedule floatingSchedule(effectiveDate, maturityDate, 6 * Months, TARGET(), ModifiedFollowing, ModifiedFollowing, DateGeneration::Forward, false); # Vector input for the NonstandardSwap obj fixedNominal = [1 for x in range(0,len(fixedSchedule)-1)] floatingNominal = [1 for x in range(0,len(floatSchedule)-1)] strike = [0.04 for x in range(0,len(fixedSchedule)-1)] spread = [0 for x in range(0,len(floatSchedule)-1)] gearing = [1 for x in range(0,len(floatSchedule)-1)] underlying = ql.NonstandardSwap(ql.VanillaSwap.Payer, fixedNominal, floatingNominal, fixedSchedule, strike, ql.Thirty360(), floatSchedule, euribor6m, gearing, spread, ql.Actual360(), False, False, ql.ModifiedFollowing) # boost::shared_ptr<NonstandardSwap> underlying = boost::make_shared<NonstandardSwap>(VanillaSwap(VanillaSwap::Payer, 1.0, fixedSchedule, strike, Thirty360(), floatingSchedule, euribor6m, 0.00, Actual360())); exerciseDates = [ql.TARGET().advance(x, -ql.Period('2D')) for x in fixedSchedule] exerciseDates = exerciseDates[1:-1] # std::vector<Date> exerciseDates; # for (Size i = 1; i < 10; ++i) # exerciseDates.push_back(TARGET().advance(fixedSchedule[i], -2 * Days)); exercise = ql.BermudanExercise(exerciseDates) # boost::shared_ptr<Exercise> exercise = boost::make_shared<BermudanExercise>(exerciseDates, false); swaption = ql.NonstandardSwaption(underlying,exercise,ql.Settlement.Physical) # boost::shared_ptr<NonstandardSwaption> swaption = boost::make_shared<NonstandardSwaption>(underlying, exercise); stepDates = exerciseDates[:-1] # std::vector<Date> stepDates(exerciseDates.begin(), exerciseDates.end() - 1); sigmas = [ql.QuoteHandle(ql.SimpleQuote(0.01)) for x in range(1, 10)] # std::vector<Real> sigmas(stepDates.size() + 1, 0.01); reversion = [ql.QuoteHandle(ql.SimpleQuote(0.01))] # Real reversion = 0.01; gsr = ql.Gsr(t0_curve, stepDates, sigmas, reversion) # boost::shared_ptr<Gsr> gsr = boost::make_shared<Gsr>(yts6m, stepDates, sigmas, reversion); swaptionEngine = ql.Gaussian1dSwaptionEngine(gsr, 64, 7.0, True, False, t0_Ois) # boost::shared_ptr<PricingEngine> swaptionEngine = boost::make_shared<Gaussian1dSwaptionEngine>(gsr, 64, 7.0, true, false, ytsOis); nonstandardSwaptionEngine = ql.Gaussian1dNonstandardSwaptionEngine(gsr, 64, 7.0, True, False, ql.QuoteHandle(ql.SimpleQuote(0)), t0_Ois) # boost::shared_ptr<PricingEngine> nonstandardSwaptionEngine = boost::make_shared<Gaussian1dNonstandardSwaptionEngine>(gsr, 64, 7.0, true, false, Handle<Quote>(), ytsOis); swaption.setPricingEngine(nonstandardSwaptionEngine) # swaption->setPricingEngine(nonstandardSwaptionEngine); swapBase = ql.EuriborSwapIsdaFixA(ql.Period('10Y'), t0_curve, t0_Ois) # boost::shared_ptr<SwapIndex> swapBase = boost::make_shared<EuriborSwapIsdaFixA>(10 * Years, yts6m, ytsOis); basket = swaption.calibrationBasket(swapBase, swaptionVol, 'Naive') # std::vector<boost::shared_ptr<CalibrationHelper> > basket = swaption->calibrationBasket(swapBase, *swaptionVol, BasketGeneratingEngine::Naive); for i in range(0, len(basket)): basket[i].setPricingEngine(swaptionEngine) # for (Size i = 0; i < basket.size(); ++i) basket[i]->setPricingEngine(swaptionEngine); method = ql.LevenbergMarquardt() # LevenbergMarquardt method; ec = ql.EndCriteria(1000, 10, 1e-8, 1e-8, 1e-8) # EndCriteria ec(1000, 10, 1E-8, 1E-8, 1E-8); gsr.calibrateVolatilitiesIterative(basket, method, ec) # gsr->calibrateVolatilitiesIterative(basket, method, ec); printBasket(basket) printModelCalibration(basket, gsr.volatility()) npv = swaption.NPV() print(npv)


Wednesday, March 29, 2017

How to install QuantLib Python for Windows 32 in offline installation

1. Download Python 2.7 Windows x86 MSI installer from

https://www.python.org/ftp/python/2.7.13/python-2.7.13.msi

2. Download required python packages and dependencies from

http://www.lfd.uci.edu/~gohlke/pythonlibs/#quantlib

QuantLib_Python‑1.9‑cp27‑cp27m‑win32.whl

Dependencies for matplotlib and others (win32)
six‑1.10.0‑py2.py3‑none‑any.whl
pyparsing‑2.2.0‑py2.py3‑none‑any.whl
packaging‑16.8‑py2.py3‑none‑any.whl
appdirs-1.4.3-py2.py3-none-any.whl
python_dateutil‑2.6.0‑py2.py3‑none‑any.whl
pytz‑2016.10‑py2.py3‑none‑any.whl
cycler‑0.10.0‑py2.py3‑none‑any.whl
setuptools‑34.3.3‑py2.py3‑none‑any.whl

numpy‑1.11.3+mkl‑cp27‑cp27m‑win32.whl
matplotlib‑1.5.3‑cp27‑cp27m‑win32.whl
xlrd-1.0.0-py2.py3-none-any.whl
pandas-0.19.2-cp27-cp27m-win32.whl
scipy‑0.19.0‑cp27‑cp27m‑win32.whl

3. Copy the above Installers and Packages to destination machine for offline installation

4. Install Python 2.7 and add Environment Variable PATH

SET PATH=C:\Python27;C:\Python27\Scripts;%PATH%

To be persistence, use sysdm.cpl to edit Advanced -> Environment Variables -> Path

or use py -2 to run python script under Windows OS without setting PATH
py -2 chap06.py

5. Install each and every the python packages above using pip

For example
pip install QuantLib_Python‑1.9‑cp27‑cp27m‑win32.whl

or if have python2 and python3 co-exist

py -2 -m pip install QuantLib_Python‑1.9‑cp27‑cp27m‑win32.whl

6. Microsoft Visual Studio not required

and no need to build your own QuantLib-Python library
However, a code editor like Microsoft VS Code is recommended.
https://code.visualstudio.com/download
(requirement : .NET Framework 4.5.2 for Windows 7)

Offline VS code extension for python can be downloaded from
https://donjayamanne.gallery.vsassets.io/_apis/public/gallery/publisher/donjayamanne/extension/python/0.6.0/assetbyname/Microsoft.VisualStudio.Services.VSIXPackage
see stack overflow discussion here http://stackoverflow.com/questions/37071388/how-to-install-vscode-extensions-offline

or NotePad ++
https://notepad-plus-plus.org/download/


7. Test

The code is borrowed from QuantLib Python Cookbook chapter 06 (requires QuantLib-Python, matplotlib)
chap06.py    Select all
#! python2 #!/usr/bin/env python # pylint: disable-msg=C0103 # Interest-rate Curves # Chapter 6 EONIA curve bootstrapping # Everything You Always Wanted to Know About Multiple Interest Rate Curve Bootstrapping but Were Afraid to Ask # In[1] import math # In[2] from QuantLib import * print ("\nOut[1]:") print ("QuantLib version", QuantLib.__version__) # In[3] # setup evaluationDate today = Date(11, December,2012) Settings.instance().evaluationDate = today # In[4] # setup DepositRateHelper for 0-2 days helpers = [ DepositRateHelper(QuoteHandle(SimpleQuote(rate/100)), Period(1,Days), fixingDays, TARGET(), Following, False, Actual360()) for rate, fixingDays in [(0.04, 0), (0.04, 1), (0.04, 2)] ] # DepositRateHelper (const Handle< Quote > &rate, # const Period &tenor, # Natural fixingDays, # const Calendar &calendar, # BusinessDayConvention convention, # bool endOfMonth, # const DayCounter &dayCounter) # DepositRateHelper (Rate rate, # const Period &tenor, # Natural fixingDays, # const Calendar &calendar, # BusinessDayConvention convention, # bool endOfMonth, # const DayCounter &dayCounter) # In[5] """ Eonia(const Handle< YieldTermStructure > &h=Handle< YieldTermStructure >()) Eonia (Euro Overnight Index Average) rate fixed by the ECB HKDHibor(const Period &tenor, const Handle< YieldTermStructure > &h=Handle< YieldTermStructure >()) """ eonia = Eonia() """ http://quant.stackexchange.com/questions/32345/quantlib-python-dual-curve-bootstrapping-example swap-rate helpers used to bootstrap the LIBOR curve can take a discount curve to use. In the old single-curve examples, a SwapRateHelper instance would be created as helper = SwapRateHelper(quoted_rate, tenor, calendar, fixedLegFrequency, fixedLegAdjustment, fixedLegDayCounter, Euribor6M()) and use the curve being bootstrapped for both forecast and discounting. To use dual-curve bootstrapping, instead, you'll have to build it as helper = SwapRateHelper(quoted_rate, tenor, calendar, fixedLegFrequency, fixedLegAdjustment, fixedLegDayCounter, Euribor6M(), QuoteHandle(), Period(0,Days), # needed as default value discountCurve) # the discountCurve argument would be a handle to the OIS curve that you bootstrapped previously SwapRateHelper (Rate rate, const Period &tenor, const Calendar &calendar, Frequency fixedFrequency, BusinessDayConvention fixedConvention, const DayCounter &fixedDayCount, const boost::shared_ptr< IborIndex > &iborIndex, const Handle< Quote < &spread=Handle< Quote >(), const Period &fwdStart=0 *Days, const Handle< YieldTermStructure > &discountingCurve=Handle< YieldTermStructure >()) In the above, the additional QuoteHandle() and Period(0,Days) arguments are, unfortunately, needed because the SWIG wrappers don't support keyword arguments for this constructor; and the discountCurve argument would be a handle to the OIS curve that you bootstrapped previously. When the swap-rate helpers are instantiated as above, they will use the LIBOR curve being bootstrapped for forecast and the OIS curve for discounting. """ # In[6] # Overnight Index Swap rate # setup OISRateHelper for 1,2,3 weeks and 1 month helpers += [ OISRateHelper(2, Period(*tenor), QuoteHandle(SimpleQuote(rate/100)), eonia) for rate, tenor in [(0.070, (1,Weeks)), (0.069, (2,Weeks)), (0.078, (3,Weeks)), (0.074, (1,Months))] ] # In[7] """ DatedOISRateHelper(const Date &startDate, const Date &endDate, const Handle< Quote > &fixedRate, const boost::shared_ptr< OvernightIndex > &overnightIndex) """ # setup DatedOISRateHelper helpers += [ DatedOISRateHelper(start_date, end_date, QuoteHandle(SimpleQuote(rate/100)), eonia) for rate, start_date, end_date in [(0.046, Date(16,January,2013), Date(13,February,2013)), (0.016, Date(13,February,2013), Date(13,March,2013)), (-0.007, Date(13,March,2013), Date(10,April,2013)), (-0.013, Date(10,April,2013), Date(8,May,2013)), (-0.014, Date(8,May,2013), Date(12,June,2013))] ] # In[8] """ Overnight Index Swap rate OISRateHelper(Natural settlementDays, const Period &tenor, const Handle< Quote > &fixedRate, const boost::shared_ptr< OvernightIndex > &overnightIndex) """ # setup OISRateHelper from 15 months to 30 years helpers += [ OISRateHelper(2, Period(*tenor), QuoteHandle(SimpleQuote(rate/100)), eonia) for rate, tenor in [(0.002, (15,Months)), (0.008, (18,Months)), (0.021, (21,Months)), (0.036, (2,Years)), (0.127, (3,Years)), (0.274, (4,Years)), (0.456, (5,Years)), (0.647, (6,Years)), (0.827, (7,Years)), (0.996, (8,Years)), (1.147, (9,Years)), (1.280, (10,Years)), (1.404, (11,Years)), (1.516, (12,Years)), (1.764, (15,Years)), (1.939, (20,Years)), (2.003, (25,Years)), (2.038, (30,Years))] ] # In[9] eonia_curve_c = PiecewiseLogCubicDiscount(0, TARGET(), helpers, Actual365Fixed()) """ # QuantLib-SWIG/SWIG/piecewiseyieldcurve.i %define export_piecewise_curve(Name, Base, Interpolator) export_piecewise_curve(PiecewiseLogCubicDiscount, Discount, MonotonicLogCubic); PiecewiseYieldCurve<Base,Interpolator>( settlementDays, # Integer settlementDays, calendar, # const Calendar& calendar, instruments, # const std::vector<boost::shared_ptr<RateHelper> >& instruments, dayCounter, # const DayCounter& dayCounter, jumps, # const std::vector<Handle<Quote> >& jumps=std::vector<Handle<Quote> >(), jumpDates, # const std::vector<Date>& jumpDates = std::vector<Date>(), accuracy, # Real accuracy = 1.0e-12, i # const Interpolator& i = Interpolator() ) """ eonia_curve_c.enableExtrapolation() # In[10] today = eonia_curve_c.referenceDate() end = today+Period(2,Years) dates = [ Date(serial) for serial in range(today.serialNumber(), end.serialNumber()+1) ] rates_c = [ eonia_curve_c.forwardRate(d, TARGET().advance(d, 1, Days), Actual360(), Simple).rate()*100 for d in dates ] # In[11] import matplotlib.pyplot as plt plt.title("Multiple Interest Rate Curve Bootstrapping") plt.plot(rates_c, '-') plt.ylabel('Rates') plt.show() # get spot rates spots = [] tenors = [] today = eonia_curve_c.referenceDate() end = today+Period(2,Years) dates = [ Date(serial) for serial in range(today.serialNumber(), end.serialNumber()+1) ] #for d in eonia_curve_c.dates(): # return boost::dynamic_pointer_cast<Name>(*self)->dates(); for d in dates: day_count = Actual360() yrs = day_count.yearFraction(today, d) compounding = Simple freq = Annual zero_rate = eonia_curve_c.zeroRate(yrs, compounding, freq) tenors.append(yrs) eq_rate = zero_rate.equivalentRate(day_count, compounding, freq, today, d).rate() spots.append(100*eq_rate) plt.title('Discount Curve') plt.plot(tenors[1::], spots[1::], linewidth=2.0) plt.xlabel('tenor (Y)') plt.ylabel('spot (%)') plt.show()





8. Test2, requires QuantLib, numpy, scipy, matplotlib

The code is borrowed from QuantLib Python Cookbook chapter 13
chap13.py    Select all
#! python2 #!/usr/bin/env python # pylint: disable-msg=C0103 # Interest-rate Models # Chapter 13 Thoughts on the Convergence of Hull-White Model Monte-Carlo Simulations # In[1] import QuantLib as ql print ("\nOut[1]:") print ("QuantLib version", ql.__version__) import matplotlib.pyplot as plt import numpy as np from scipy.integrate import simps, cumtrapz, romb # % matplotlib inline import math todays_date = ql.Date(15, 1, 2015) ql.Settings.instance().evaluationDate = todays_date # In[2] # <!-- collapse=True --> def get_path_generator(timestep, hw_process, length, low_discrepancy=False, brownian_bridge=True): """ Returns a path generator The `get_path_generator` function creates the a path generator. This function takes various inputs such as """ if low_discrepancy: usg = ql.UniformLowDiscrepancySequenceGenerator(timestep) rng = ql.GaussianLowDiscrepancySequenceGenerator(usg) seq = ql.GaussianSobolPathGenerator( hw_process, length, timestep, rng,brownian_bridge) else: usg = ql.UniformRandomSequenceGenerator(timestep, ql.UniformRandomGenerator()) rng = ql.GaussianRandomSequenceGenerator(usg) seq = ql.GaussianPathGenerator( hw_process, length, timestep, rng, brownian_bridge) return seq # In[3] # <!-- collapse=True --> def generate_paths(num_paths, timestep, seq): """ The `generate_paths` function uses the generic path generator produced by the `get_path_generator` function to return a tuple of the array of the points in the time grid and a matrix of the short rates generated." """ arr = np.zeros((num_paths, timestep+1)) for i in range(num_paths): sample_path = seq.next() path = sample_path.value() time = [path.time(j) for j in range(len(path))] value = [path[j] for j in range(len(path))] arr[i, :] = np.array(value) return np.array(time), arr # In[4] # <!-- collapse=True --> def generate_paths_zero_price(spot_curve_handle, a, sigma, timestep, length, num_paths, avg_grid_array, low_discrepancy=False, brownian_bridge=True): """ This function returns a tuple (T_array, F_array), where T_array is the array of points in the time grid, and F_array is the array of the average of zero prices observed from the simulation. The `generate_paths_zero_price` essentially is a wrapper around `generate_path_generator` and `generate_paths` taking all the required raw inputs. This function returns the average of zero prices from all the paths for different points in time. I wrote this out so that I can conveniently change all the required inputs and easily plot the results." """ hw_process = ql.HullWhiteProcess(spot_curve_handle, a, sigma) seq = get_path_generator( timestep, hw_process, length, low_discrepancy, brownian_bridge ) time, paths = generate_paths(num_paths, timestep, seq) avgs = [(time[j], (np.mean([math.exp(-simps(paths[i][0:j], time[0:j])) for i in range(num_paths)]))) for j in avg_grid_array ] return zip(*avgs) def generate_paths_discount_factors(spot_curve_handle, a, sigma, timestep, length, num_paths, avg_grid_array, low_discrepancy=False, brownian_bridge=True): """ This function returns a tuple (T_array, S_matrix), where T_array is the array of points in the time grid, and S_matrix is the matrix of the spot rates for each path in the different points in the time grid. """ hw_process = ql.HullWhiteProcess(spot_curve_handle, a, sigma) seq = get_path_generator( timestep, hw_process, length, low_discrepancy, brownian_bridge ) time, paths = generate_paths(num_paths, timestep, seq) arr = np.zeros((num_paths, len(avg_grid_array))) for i in range(num_paths): arr[i, :] = [np.exp(-simps(paths[i][0:j], time[0:j])) for j in avg_grid_array ] t_array = [time[j] for j in avg_grid_array] return t_array, arr def V(t,T, a, sigma): """ Variance of the integral of short rates, used below """ return sigma*sigma/a/a*(T-t + 2.0/a*math.exp(-a*(T-t)) - 1.0/(2.0*a)*math.exp(-2.0*a*(T-t)) - 3.0/(2.0*a) ) # In[5] # <!-- collapse=True --> # Here we vary sigma with fixed a and observe the error epsilon # define constants num_paths = 500 sigma_array = np.arange(0.01,0.1,0.03) a = 0.1 timestep = 180 length = 15 # in years forward_rate = 0.05 day_count = ql.Thirty360() avg_grid_array = np.arange(12, timestep+1, 12) # generate spot curve spot_curve = ql.FlatForward( todays_date, ql.QuoteHandle(ql.SimpleQuote(forward_rate)), day_count ) spot_curve_handle = ql.YieldTermStructureHandle(spot_curve) #initialize plots figure, axis = plt.subplots() plots = [] zero_price_theory = np.array([spot_curve.discount(j*float(length)/float(timestep)) for j in avg_grid_array]) for sigma in sigma_array: term, zero_price_empirical = generate_paths_zero_price( spot_curve_handle, a, sigma, timestep, length, num_paths, avg_grid_array ) plots += axis.plot( term, np.abs(zero_price_theory - np.array(zero_price_empirical)), lw=2, alpha=0.6, label="$\sigma=$"+str(sigma) ) # plot legend labels = [p.get_label() for p in plots] legend =axis.legend(plots,labels, loc=0)#, loc=0, bbox_to_anchor=(1.1,0.4)) axis.set_xlabel("T (years)", size=12) axis.set_ylabel("|$\epsilon(T)$|", size=12) axis.set_title("Out[5]:Discount Factor Error for $a=$%0.2f and Varying $\sigma$"%a, size=14) plt.show() # In[6] # <!-- collapse=True --> # Here we vary a with fixed sigma and observe the error epsilon # define constants num_paths = 500 sigma = 0.1 a_array = np.arange(0.1, 0.51, 0.1) timestep = 180 length = 15 # in years forward_rate = 0.05 day_count = ql.Thirty360() avg_grid_array = np.arange(12, timestep+1, 12) # generate spot curve spot_curve = ql.FlatForward( todays_date, ql.QuoteHandle(ql.SimpleQuote(forward_rate)), day_count ) spot_curve_handle = ql.YieldTermStructureHandle(spot_curve) #initialize plots figure, axis = plt.subplots() plots = [] zero_price_theory = np.array([spot_curve.discount(j*float(length)/float(timestep)) for j in avg_grid_array]) for a in a_array: term, zero_price_empirical = generate_paths_zero_price( spot_curve_handle, a, sigma, timestep, length, num_paths, avg_grid_array ) plots += axis.plot( term,np.abs(zero_price_theory - np.array(zero_price_empirical)), lw=2, alpha=0.6, label="a="+str(a) ) # plot legend labels = [p.get_label() for p in plots] legend =axis.legend(plots,labels, loc=0)#, loc=0, bbox_to_anchor=(1.1,0.4)) axis.set_xlabel("T (years)", size=12) axis.set_ylabel("|$\\epsilon(T)$|", size=12) axis.set_title("Out[6]:Discount Factor Error for $\sigma$=%0.2f and Varying $a$"%sigma, size=14) plt.show() # In[7] # <!-- collapse=True --> #define constants num_paths = 500 sigma = 0.02 a = 0.1 timestep = 180 length = 15 # in years forward_rate = 0.05 day_count = ql.Thirty360() avg_grid_array = np.arange(1, timestep+1, 12) # generate spot curve spot_curve = ql.FlatForward( todays_date, ql.QuoteHandle(ql.SimpleQuote(forward_rate)), day_count ) spot_curve_handle = ql.YieldTermStructureHandle(spot_curve) term, discount_factor_matrix = generate_paths_discount_factors( spot_curve_handle, a, sigma, timestep, length, num_paths, avg_grid_array ) vol = [np.var(discount_factor_matrix[:, i]) for i in range(len(term))] l1 = plt.plot(term, 100*np.sqrt(vol),"b", lw=2, alpha=0.6, label="Empirical") vol_theory = [100*np.sqrt(math.exp(V(0,T,a, sigma))-1.0) * spot_curve_handle.discount(T) for T in term] l2 = plt.plot(term, vol_theory,"r--", lw=2, alpha=0.6, label="Theory") plots = l1+l2 labels = [p.get_label() for p in plots] legend =plt.legend(plots,labels, loc=0) plt.xlabel("Time (Years)", size=12) plt.ylabel("$\sigma_D(0,T)$ (%)", size=12) plt.title("Out[7]:Standard Deviation of Discount Factors " "(a=%0.2f, $\sigma$=%0.2f)"%(a, sigma), size=14) plt.show()




9. Test3, requires QuantLib, numpy, panda, xlrd

The code is using QuantLib, numpy, panda and xlrd to read xls data
testreadxls.py    Select all
#! python2 #!/usr/bin/env python # pylint: disable-msg=C0103 # pylint: disable-msg=C0301 # testreadxls.py import os import pandas import numpy as np import QuantLib as ql xls_file = os.path.dirname(os.path.realpath(__file__)) + '/usd_market_data_2016-07-13.xls' #alternate way to mount google drive and read xlsx file data #import os #from google.colab import drive #drive.mount('/content/drive') #xls_file = '/content/drive/My Drive/usd_market_data_2016-07-13.xlsx' # if the file is in csv format #csv_file = '/content/drive/My Drive/usd_market_data_2016-07-13.csv' #csv = pandas.read_csv(csv_file) #print("\nCSV data\n") #print(csv) #print (csv.loc[:, 'dates']) #print (csv.iloc[0:, 1]) #discounts_csv = np.array(csv.iloc[0:,1].tolist()) #print(discounts_csv) xl = pandas.ExcelFile(xls_file) print (xl.sheet_names) #df = xl.parse("usdstd") #df = xl.parse(0) # alternate way to read first worksheet filename = xls_file sheetname = 'usdois' sheet = pandas.read_excel(filename, sheetname) # Read an Excel table into a pandas DataFrame #sheet = pandas.read_excel(filename, 0) # alternate way to read first worksheet instead of name print (sheet.columns) #tenors = np.array([Date.from_timestamp(d).t for d in sheet['dates']]) today = ql.Date(21,3,2016) act365 = ql.Actual365Fixed() dates = np.array([ql.Date(d.day, d.month, d.year) for d in sheet['dates']]) #dates = np.array([ql.Date(d.day, d.month, d.year) for d in df.iloc[0:,0]]) # use iloc to read first column, all rows print ("\nusd_market_data_2016-07-13.xls dates column") print (dates) tenors = np.array([act365.yearFraction(today,ql.Date(d.day, d.month, d.year)) for d in sheet['dates']]) print ("\nusd_market_data_2016-07-13.xls dates column convert to tenor") print (tenors) print ("\nusd_market_data_2016-07-13.xls discounts column") discounts = np.array(sheet['discounts'].tolist()) #discounts = np.array(df.iloc[0:,1].tolist()) # use iloc to read second column, all rows print (discounts)


Assume the usd_market_data_2016-07-13.xls sheet usdois (first worksheet) has the following data
https://mega.nz/#!m4Bh2TJD!KVZbvAn4D_sSR8-FYJtNASkiWvY2TWZ1LokJ-vcw-u4


10. Test3, requires QuantLib, pandas, matplotlib

The code demonstrates the download of Interest Rate xml market data from markit.com and Bootstrapping IR Curve in QuantLib
test_xml_download.py    Select all
#!python2 #!/usr/bin/env python # pylint: disable-msg=C0103 # pylint: disable-msg=C0301 # test_xml_download.py from QuantLib import * print ("\nQuantLib version", QuantLib.__version__) import urllib import zipfile import xml.etree.ElementTree as ET import pandas as pd import datetime as dt from matplotlib.dates import YearLocator, MonthLocator, DateFormatter from matplotlib.ticker import FuncFormatter import sys if sys.version_info[0] >= 3: from urllib.request import urlretrieve else: # Not Python 3 - today, it is most likely to be Python 2 # But note that this might need an update when Python 4 # might be around one day from urllib import urlretrieve def to_datetime(d): """ to_datetime """ return dt.datetime(d.year(), d.month(), d.dayOfMonth()) def format_rate(r, p=2): """ format_rate """ return '{0:.{prec}f}%'.format(r*100.00, prec=p) period_dict = {str(k)+'M': (k, Months) for k in range(1, 12, 1)} d2 = {str(k)+'Y': (k, Years) for k in range(1, 50, 1)} period_dict.update(d2) (ir_currency, ir_date) = ('USD', '20170331') url = 'https://www.markit.com/news/InterestRates_%s_%s.zip' % (ir_currency, ir_date) #filehandle, _ = urllib.urlretrieve(url) filehandle, _ = urlretrieve(url) zip_file_object = zipfile.ZipFile(filehandle, 'r') interest_rate_file = zip_file_object.open(zip_file_object.namelist()[1]) content = interest_rate_file.read() #print (content) root = ET.fromstring(content) print ("\nInfo on InterestRates_%s_%s" % (ir_currency, ir_date)) for a in [root.find(b) for b in ['./currency', './effectiveasof', './deposits/daycountconvention', './deposits/snaptime', './deposits/spotdate', './swaps/fixeddaycountconvention', './swaps/floatingdaycountconvention', './swaps/snaptime', './swaps/spotdate']]: print (a.tag, a.text) effectiveasof = root.find('./effectiveasof').text currency = root.find('./currency').text today = DateParser.parseFormatted(effectiveasof, '%Y-%m-%d') Settings.instance().evaluationDate = today print ("\nDeposit Rates") deposits_maturitydates = [a.text for a in root.findall('.//deposits/curvepoint/maturitydate')] deposits_parrates_text = [a.text for a in root.findall('.//deposits/curvepoint/parrate')] deposits_tenors = [a.text for a in root.findall('.//deposits/curvepoint/tenor')] deposits_periods = [period_dict[a] for a in deposits_tenors] print (pd.DataFrame(zip(deposits_tenors, deposits_maturitydates, map(float, deposits_parrates_text)), columns=('Tenor', 'Maturity', 'Parrate'))) maturity_dict = {} for t, m in zip(deposits_tenors, deposits_maturitydates): d = {t: m} maturity_dict.update(d) print ("\nSwap Rates") swaps_maturitydates = [a.text for a in root.findall('.//swaps/curvepoint/maturitydate')] swaps_parrates_text = [a.text for a in root.findall('.//swaps/curvepoint/parrate')] swaps_tenors = [a.text for a in root.findall('.//swaps/curvepoint/tenor')] swaps_periods = [period_dict[a] for a in swaps_tenors] print (pd.DataFrame(zip(swaps_tenors, swaps_maturitydates, map(float, swaps_parrates_text)), columns=('Tenor', 'Maturity', 'Parrate'))) for t, m in zip(swaps_tenors, swaps_maturitydates): d = {t: m} maturity_dict.update(d) #print (maturity_dict) euribor6m = Euribor6M() helpers = [DepositRateHelper(rate, Period(*tenor), 2, TARGET(), Following, False, Actual360()) for tenor, rate in zip(deposits_periods, map(float, deposits_parrates_text))] helpers += [SwapRateHelper(rate, Period(*tenor), TARGET(), Semiannual, Unadjusted, Thirty360(), euribor6m) for tenor, rate in zip(swaps_periods, map(float, swaps_parrates_text))] curve = PiecewiseLogCubicDiscount(2, TARGET(), helpers, Actual360()) curve.enableExtrapolation() spot = curve.referenceDate() #dates = [spot+Period(i,Months) for i in range(0, 30*12+1)] dates = [spot+Period(i, Years) for i in range(0, 30)] rates = [curve.forwardRate(d, euribor6m.maturityDate(d), Actual360(), Simple).rate() for d in dates] valid_dates = [d for d in dates if d >= spot] import matplotlib.pyplot as plt fig, ax = plt.subplots() fig.autofmt_xdate() ax.plot_date([to_datetime(d) for d in valid_dates], rates, '-') ax.set_xlim(to_datetime(min(dates)), to_datetime(max(dates))) ax.xaxis.set_major_locator(YearLocator(2, month=today.month(), day=today.dayOfMonth())) ax.xaxis.set_major_formatter(DateFormatter("%Y")) ax.xaxis.grid(True, 'major') ax.xaxis.grid(False, 'minor') ax.yaxis.set_major_formatter(FuncFormatter(lambda r, pos: format_rate(r))) plt.title("%s Interest Rate Curve Bootstrapping, effective %s" % (currency, effectiveasof)) plt.xlabel('Year') plt.ylabel('Rates') plt.show()






11. Compile QuantLib Python 1.10 for Mac

script.sh    Select all
# install Xcode and HomeBrew # ruby -e "$(curl -fsSL https://raw.githubusercontent.com/Homebrew/install/master/install)" # install boost brew install boost # Download QuantLib-1.10 cd $(HOME)/Downloads wget --no-check-certificate https://jaist.dl.sourceforge.net/project/quantlib/test/QuantLib-1.10.tar.gz tar xzvf QuantLib-1.10.tar.gz # Compile QuantLib 1.10 cd QuantLib-1.10 ./configure --prefix=/usr/local/ CXXFLAGS='-O2 -stdlib=libstdc++ -mmacosx-version-min=10.6' LDFLAGS='-stdlib=libstdc++ -mmacosx-version-min=10.6' make && sudo make install # test compile c++ g++ Bonds.cpp -std=c++11 -stdlib=libstdc++ -o Bonds -lQuantLib ./Bonds # Download QuantLib-SWIG 1.10 cd $(HOME)/Downloads wget --no-check-certificate https://jaist.dl.sourceforge.net/project/quantlib/test/QuantLib-SWIG-1.10.tar.gz tar xzvf QuantLib-SWIG-1.10.tar.gz # Compile QuantLib Python cd QuantLib-SWIG-1.10 ./configure CXXFLAGS='-O2 -stdlib=libstdc++ -mmacosx-version-min=10.6' make -C Python make -C Python check sudo make -C Python install # Package wheel file sudo -H python pip install setuptools wheel cd $(HOME)/Downloads/QuantLib-SWIG-1.10/Python # Modify setup.py and replace from distutils.core import setup, Extension # By try: from setuptools import setup, Extension except: from distutils.core import setup, Extension # python 2 wheel file python setup.py bdist_wheel ls dist/ # python 3 wheel file LDFLAGS="-arch x86_64 -bundle -flat_namespace -undefined suppress" CFLAGS="-fno-strict-aliasing -Wsign-compare -fno-common -static -arch x86_64 -Wno-shorten-64-to-32" python3 setup.py bdist_wheel ls dist/


QuantLib_Python-1.10-cp27-cp27m-macosx_10_10_intel.whl download

QuantLib_Python-1.10-cp36-cp36m-macosx_10_6_intel.whl download
For Python3, have to install libQuantLib to /usr/local/lib
cd /usr/local/lib
sudo tar xzvf ~/Download/libQuantLib.tgz
Download libQuantLib.tgz here https://mega.nz/#!rpJg0B6T!NlG3Ijgo0weNOi_CVOkCzpGCMRqHZkWSeF3S0sAR12o



12. For Windows 64-bit, Regeneration of python SWIG interface method
If you need to edit the SWIG interface file and add functionality in QuantLib Python, you need to regenerate.


script.cmd    Select all
# First, compile QuantLib using Visual Studio 2015. # Select 21 projects and right click select Property -> All Configurations Add C:\local\boost_1_59_0; to the Property Pages : VC++ Include directories to win32 platform and Add C:\local\boost_1_59_0\lib32-msvc-14.0; to the VC++ Libs directories to win32 platform # Select All Configurations Add C:\local\boost_1_59_0_64; to the Property Pages : VC++ Include directories to x64 platform and Add C:\local\boost_1_59_0_64\lib64-msvc-14.0; to the VC++ Libs directories to x64 platform # Add /wd4819 to Command Line : Additional Options to disable C4819 warnings when compiling Quantlib to both win32 and x64 platform # Build Release x64 Solution or Release win32 Solution. # setting of QuantLib Project directory, e.g. SET QL_DIR=C:\local\QuantLib-1.10\QuantLib-1.10 # or alternatively, edit the setup.py # QL_INSTALL_DIR = r'C:\local\QuantLib-1.10\QuantLib-1.10' # setting of compiled libraries location for QuantLib ('QuantLib-vc140-x64-mt.lib') and Boost, e.g. for x64 platform SET LIB=C:\local\boost_1_59_0_64\lib64-msvc-14.0 # setting of Boost Project directory, e.g. SET INCLUDE=C:\local\boost_1_59_0_64 # setting of SWIGWIN.exe PATH, e.g. SET PATH=C:\local\swigwin-3.0.12;%PATH% # setting of VISUAL STUDIO 2015 (Version 14) SET VS90COMNTOOLS=%VS140COMNTOOLS% # Edit SWIG interface file if any C:\local\QuantLib-1.10\QuantLib-SWIG-1.10\SWIG\*.i # e.g. edit piecewiseyieldcurve.i and add at the end # export_piecewise_curve(PiecewiseLogLinearDiscount,Discount,LogLinear); # Wrap the edited interface file again after edit py -2 setup.py wrap # compile using msvc # need to edit setup.py to add ,'/wd4819' to the extra_compile_args to disable warnings py -2 setup.py build --compiler=msvc # must use Administrator Command Prompt to install py -2 -m setup.py install --skip-build # generation of wheel file is same as above that is modify setup.py first # need to install setuptools and wheel py -2 -m pip install setuptools wheel py -2 setup.py test py -2 setup.py bdist_wheel dir dist


Visual Studio 2015 setup https://mega.nz/#!70ZwyARa!z4et3sKwguU16tEbCYTdwBC1VEkNJTTakryrDRBSKM8