update on getting started
This commit is contained in:
mhjensen
Vendored
BIN
Binary file not shown.
@@ -0,0 +1,330 @@
|
||||
# import necessary packages
|
||||
import numpy as np
|
||||
import matplotlib.pyplot as plt
|
||||
from sklearn import datasets
|
||||
|
||||
|
||||
# ensure the same random numbers appear every time
|
||||
np.random.seed(0)
|
||||
|
||||
# display images in notebook
|
||||
plt.rcParams['figure.figsize'] = (12,12)
|
||||
|
||||
|
||||
# download MNIST dataset
|
||||
digits = datasets.load_digits()
|
||||
|
||||
# define inputs and labels
|
||||
inputs = digits.images
|
||||
labels = digits.target
|
||||
|
||||
# RGB images have a depth of 3
|
||||
# our images are grayscale so they should have a depth of 1
|
||||
inputs = inputs[:,:,:,np.newaxis]
|
||||
|
||||
print("inputs = (n_inputs, pixel_width, pixel_height, depth) = " + str(inputs.shape))
|
||||
print("labels = (n_inputs) = " + str(labels.shape))
|
||||
|
||||
|
||||
# choose some random images to display
|
||||
n_inputs = len(inputs)
|
||||
indices = np.arange(n_inputs)
|
||||
random_indices = np.random.choice(indices, size=5)
|
||||
|
||||
for i, image in enumerate(digits.images[random_indices]):
|
||||
plt.subplot(1, 5, i+1)
|
||||
plt.axis('off')
|
||||
plt.imshow(image, cmap=plt.cm.gray_r, interpolation='nearest')
|
||||
plt.title("Label: %d" % digits.target[random_indices[i]])
|
||||
plt.show()
|
||||
|
||||
from keras.utils import to_categorical
|
||||
from sklearn.model_selection import train_test_split
|
||||
|
||||
# representation of labels
|
||||
labels = to_categorical(labels)
|
||||
|
||||
# split into train and test data
|
||||
# one-liner from scikit-learn library
|
||||
train_size = 0.8
|
||||
test_size = 1 - train_size
|
||||
X_train, X_test, Y_train, Y_test = train_test_split(inputs, labels, train_size=train_size,
|
||||
test_size=test_size)
|
||||
|
||||
import tensorflow as tf
|
||||
|
||||
class ConvolutionalNeuralNetworkTensorflow:
|
||||
def __init__(
|
||||
self,
|
||||
X_train,
|
||||
Y_train,
|
||||
X_test,
|
||||
Y_test,
|
||||
n_filters=10,
|
||||
n_neurons_connected=50,
|
||||
n_categories=10,
|
||||
receptive_field=3,
|
||||
stride=1,
|
||||
padding=1,
|
||||
epochs=10,
|
||||
batch_size=100,
|
||||
eta=0.1,
|
||||
lmbd=0.0,
|
||||
):
|
||||
|
||||
self.global_step = tf.Variable(0, dtype=tf.int32, trainable=False, name='global_step')
|
||||
|
||||
self.X_train = X_train
|
||||
self.Y_train = Y_train
|
||||
self.X_test = X_test
|
||||
self.Y_test = Y_test
|
||||
|
||||
self.n_inputs, self.input_width, self.input_height, self.depth = X_train.shape
|
||||
|
||||
self.n_filters = n_filters
|
||||
self.n_downsampled = int(self.input_width*self.input_height*n_filters / 4)
|
||||
self.n_neurons_connected = n_neurons_connected
|
||||
self.n_categories = n_categories
|
||||
|
||||
self.receptive_field = receptive_field
|
||||
self.stride = stride
|
||||
self.strides = [stride, stride, stride, stride]
|
||||
self.padding = padding
|
||||
|
||||
self.epochs = epochs
|
||||
self.batch_size = batch_size
|
||||
self.iterations = self.n_inputs // self.batch_size
|
||||
self.eta = eta
|
||||
self.lmbd = lmbd
|
||||
|
||||
self.create_placeholders()
|
||||
self.create_CNN()
|
||||
self.create_loss()
|
||||
self.create_optimiser()
|
||||
self.create_accuracy()
|
||||
|
||||
def create_placeholders(self):
|
||||
with tf.name_scope('data'):
|
||||
self.X = tf.placeholder(tf.float32, shape=(None, self.input_width, self.input_height, self.depth), name='X_data')
|
||||
self.Y = tf.placeholder(tf.float32, shape=(None, self.n_categories), name='Y_data')
|
||||
|
||||
def create_CNN(self):
|
||||
with tf.name_scope('CNN'):
|
||||
|
||||
# Convolutional layer
|
||||
self.W_conv = self.weight_variable([self.receptive_field, self.receptive_field, self.depth, self.n_filters], name='conv', dtype=tf.float32)
|
||||
b_conv = self.weight_variable([self.n_filters], name='conv', dtype=tf.float32)
|
||||
z_conv = tf.nn.conv2d(self.X, self.W_conv, self.strides, padding='SAME', name='conv') + b_conv
|
||||
a_conv = tf.nn.relu(z_conv)
|
||||
|
||||
# 2x2 max pooling
|
||||
a_pool = tf.nn.max_pool(a_conv, [1, 2, 2, 1], [1, 2, 2, 1], padding='SAME', name='pool')
|
||||
|
||||
# Fully connected layer
|
||||
a_pool_flat = tf.reshape(a_pool, [-1, self.n_downsampled])
|
||||
self.W_fc = self.weight_variable([self.n_downsampled, self.n_neurons_connected], name='fc', dtype=tf.float32)
|
||||
b_fc = self.bias_variable([self.n_neurons_connected], name='fc', dtype=tf.float32)
|
||||
a_fc = tf.nn.relu(tf.matmul(a_pool_flat, self.W_fc) + b_fc)
|
||||
|
||||
# Output layer
|
||||
self.W_out = self.weight_variable([self.n_neurons_connected, self.n_categories], name='out', dtype=tf.float32)
|
||||
b_out = self.bias_variable([self.n_categories], name='out', dtype=tf.float32)
|
||||
self.z_out = tf.matmul(a_fc, self.W_out) + b_out
|
||||
|
||||
def create_loss(self):
|
||||
with tf.name_scope('loss'):
|
||||
softmax_loss = tf.reduce_mean(tf.nn.softmax_cross_entropy_with_logits_v2(labels=self.Y, logits=self.z_out))
|
||||
|
||||
regularizer_loss_conv = tf.nn.l2_loss(self.W_conv)
|
||||
regularizer_loss_fc = tf.nn.l2_loss(self.W_fc)
|
||||
regularizer_loss_out = tf.nn.l2_loss(self.W_out)
|
||||
regularizer_loss = self.lmbd*(regularizer_loss_conv + regularizer_loss_fc + regularizer_loss_out)
|
||||
|
||||
self.loss = softmax_loss + regularizer_loss
|
||||
|
||||
def create_accuracy(self):
|
||||
with tf.name_scope('accuracy'):
|
||||
probabilities = tf.nn.softmax(self.z_out)
|
||||
predictions = tf.argmax(probabilities, 1)
|
||||
labels = tf.argmax(self.Y, 1)
|
||||
|
||||
correct_predictions = tf.equal(predictions, labels)
|
||||
correct_predictions = tf.cast(correct_predictions, tf.float32)
|
||||
self.accuracy = tf.reduce_mean(correct_predictions)
|
||||
|
||||
def create_optimiser(self):
|
||||
with tf.name_scope('optimizer'):
|
||||
self.optimizer = tf.train.GradientDescentOptimizer(learning_rate=self.eta).minimize(self.loss, global_step=self.global_step)
|
||||
|
||||
def weight_variable(self, shape, name='', dtype=tf.float32):
|
||||
initial = tf.truncated_normal(shape, stddev=0.1)
|
||||
return tf.Variable(initial, name=name, dtype=dtype)
|
||||
|
||||
def bias_variable(self, shape, name='', dtype=tf.float32):
|
||||
initial = tf.constant(0.1, shape=shape)
|
||||
return tf.Variable(initial, name=name, dtype=dtype)
|
||||
|
||||
def fit(self):
|
||||
data_indices = np.arange(self.n_inputs)
|
||||
|
||||
with tf.Session() as sess:
|
||||
sess.run(tf.global_variables_initializer())
|
||||
for i in range(self.epochs):
|
||||
for j in range(self.iterations):
|
||||
chosen_datapoints = np.random.choice(data_indices, size=self.batch_size, replace=False)
|
||||
batch_X, batch_Y = self.X_train[chosen_datapoints], self.Y_train[chosen_datapoints]
|
||||
|
||||
sess.run([CNN.loss, CNN.optimizer],
|
||||
feed_dict={CNN.X: batch_X,
|
||||
CNN.Y: batch_Y})
|
||||
accuracy = sess.run(CNN.accuracy,
|
||||
feed_dict={CNN.X: batch_X,
|
||||
CNN.Y: batch_Y})
|
||||
step = sess.run(CNN.global_step)
|
||||
|
||||
self.train_loss, self.train_accuracy = sess.run([CNN.loss, CNN.accuracy],
|
||||
feed_dict={CNN.X: self.X_train,
|
||||
CNN.Y: self.Y_train})
|
||||
|
||||
self.test_loss, self.test_accuracy = sess.run([CNN.loss, CNN.accuracy],
|
||||
feed_dict={CNN.X: self.X_test,
|
||||
CNN.Y: self.Y_test})
|
||||
|
||||
epochs = 100
|
||||
batch_size = 100
|
||||
n_filters = 10
|
||||
n_neurons_connected = 50
|
||||
n_categories = 10
|
||||
|
||||
eta_vals = np.logspace(-5, 1, 7)
|
||||
lmbd_vals = np.logspace(-5, 1, 7)
|
||||
CNN_tf = np.zeros((len(eta_vals), len(lmbd_vals)), dtype=object)
|
||||
|
||||
for i, eta in enumerate(eta_vals):
|
||||
for j, lmbd in enumerate(lmbd_vals):
|
||||
CNN = ConvolutionalNeuralNetworkTensorflow(X_train, Y_train, X_test, Y_test,
|
||||
n_filters=n_filters, n_neurons_connected=n_neurons_connected,
|
||||
n_categories=n_categories, epochs=epochs, batch_size=batch_size,
|
||||
eta=eta, lmbd=lmbd)
|
||||
CNN.fit()
|
||||
|
||||
print("Learning rate = ", eta)
|
||||
print("Lambda = ", lmbd)
|
||||
print("Test accuracy: %.3f" % CNN.test_accuracy)
|
||||
print()
|
||||
|
||||
CNN_tf[i][j] = CNN
|
||||
|
||||
# visual representation of grid search
|
||||
# uses seaborn heatmap, could probably do this in matplotlib
|
||||
import seaborn as sns
|
||||
|
||||
sns.set()
|
||||
|
||||
train_accuracy = np.zeros((len(eta_vals), len(lmbd_vals)))
|
||||
test_accuracy = np.zeros((len(eta_vals), len(lmbd_vals)))
|
||||
|
||||
for i in range(len(eta_vals)):
|
||||
for j in range(len(lmbd_vals)):
|
||||
CNN = CNN_tf[i][j]
|
||||
|
||||
train_accuracy[i][j] = CNN.train_accuracy
|
||||
test_accuracy[i][j] = CNN.test_accuracy
|
||||
|
||||
|
||||
fig, ax = plt.subplots(figsize = (10, 10))
|
||||
sns.heatmap(train_accuracy, annot=True, ax=ax, cmap="viridis")
|
||||
ax.set_title("Training Accuracy")
|
||||
ax.set_ylabel("$\eta$")
|
||||
ax.set_xlabel("$\lambda$")
|
||||
plt.show()
|
||||
|
||||
fig, ax = plt.subplots(figsize = (10, 10))
|
||||
sns.heatmap(test_accuracy, annot=True, ax=ax, cmap="viridis")
|
||||
ax.set_title("Test Accuracy")
|
||||
ax.set_ylabel("$\eta$")
|
||||
ax.set_xlabel("$\lambda$")
|
||||
plt.show()
|
||||
|
||||
from keras.models import Sequential
|
||||
from keras.layers.convolutional import Conv2D
|
||||
from keras.layers.convolutional import MaxPooling2D
|
||||
from keras.layers import Flatten
|
||||
from keras.layers import Dense
|
||||
from keras.regularizers import l2
|
||||
from keras.optimizers import SGD
|
||||
|
||||
def create_convolutional_neural_network_keras(input_shape, receptive_field,
|
||||
n_filters, n_neurons_connected, n_categories,
|
||||
eta, lmbd):
|
||||
model = Sequential()
|
||||
model.add(Conv2D(n_filters, (receptive_field, receptive_field), input_shape=input_shape, padding='same',
|
||||
activation='relu', kernel_regularizer=l2(lmbd)))
|
||||
model.add(MaxPooling2D(pool_size=(2, 2)))
|
||||
model.add(Flatten())
|
||||
model.add(Dense(n_neurons_connected, activation='relu', kernel_regularizer=l2(lmbd)))
|
||||
model.add(Dense(n_categories, activation='softmax', kernel_regularizer=l2(lmbd)))
|
||||
|
||||
sgd = SGD(lr=eta)
|
||||
model.compile(loss='categorical_crossentropy', optimizer=sgd, metrics=['accuracy'])
|
||||
|
||||
return model
|
||||
|
||||
epochs = 100
|
||||
batch_size = 100
|
||||
input_shape = X_train.shape[1:4]
|
||||
receptive_field = 3
|
||||
n_filters = 10
|
||||
n_neurons_connected = 50
|
||||
n_categories = 10
|
||||
|
||||
eta_vals = np.logspace(-5, 1, 7)
|
||||
lmbd_vals = np.logspace(-5, 1, 7)
|
||||
|
||||
CNN_keras = np.zeros((len(eta_vals), len(lmbd_vals)), dtype=object)
|
||||
|
||||
for i, eta in enumerate(eta_vals):
|
||||
for j, lmbd in enumerate(lmbd_vals):
|
||||
CNN = create_convolutional_neural_network_keras(input_shape, receptive_field,
|
||||
n_filters, n_neurons_connected, n_categories,
|
||||
eta, lmbd)
|
||||
CNN.fit(X_train, Y_train, epochs=epochs, batch_size=batch_size, verbose=0)
|
||||
scores = CNN.evaluate(X_test, Y_test)
|
||||
|
||||
CNN_keras[i][j] = CNN
|
||||
|
||||
print("Learning rate = ", eta)
|
||||
print("Lambda = ", lmbd)
|
||||
print("Test accuracy: %.3f" % scores[1])
|
||||
print()
|
||||
|
||||
# visual representation of grid search
|
||||
# uses seaborn heatmap, could probably do this in matplotlib
|
||||
import seaborn as sns
|
||||
|
||||
sns.set()
|
||||
|
||||
train_accuracy = np.zeros((len(eta_vals), len(lmbd_vals)))
|
||||
test_accuracy = np.zeros((len(eta_vals), len(lmbd_vals)))
|
||||
|
||||
for i in range(len(eta_vals)):
|
||||
for j in range(len(lmbd_vals)):
|
||||
CNN = CNN_keras[i][j]
|
||||
|
||||
train_accuracy[i][j] = CNN.evaluate(X_train, Y_train)[1]
|
||||
test_accuracy[i][j] = CNN.evaluate(X_test, Y_test)[1]
|
||||
|
||||
|
||||
fig, ax = plt.subplots(figsize = (10, 10))
|
||||
sns.heatmap(train_accuracy, annot=True, ax=ax, cmap="viridis")
|
||||
ax.set_title("Training Accuracy")
|
||||
ax.set_ylabel("$\eta$")
|
||||
ax.set_xlabel("$\lambda$")
|
||||
plt.show()
|
||||
|
||||
fig, ax = plt.subplots(figsize = (10, 10))
|
||||
sns.heatmap(test_accuracy, annot=True, ax=ax, cmap="viridis")
|
||||
ax.set_title("Test Accuracy")
|
||||
ax.set_ylabel("$\eta$")
|
||||
ax.set_xlabel("$\lambda$")
|
||||
plt.show()
|
||||
@@ -0,0 +1,34 @@
|
||||
#The covariance matrix and its eigenvalues the hard way
|
||||
from random import random, seed
|
||||
import numpy as np
|
||||
|
||||
def covariance(x, y, n):
|
||||
sum = 0.0
|
||||
mean_x = np.mean(x)
|
||||
mean_y = np.mean(y)
|
||||
for i in range(0, n):
|
||||
sum += (x[(i)]-mean_x)*(y[i]-mean_y)
|
||||
return sum/n
|
||||
|
||||
n = 10
|
||||
|
||||
x = np.random.normal(size=n)
|
||||
y = np.random.normal(size=n)
|
||||
z = x*x*x+y*y 0.5*np.random.normal(size=n)
|
||||
covxx = covariance(x,x,n)
|
||||
covxy = covariance(x,y,n)
|
||||
covxz = covariance(x,z,n)
|
||||
covyy = covariance(y,y,n)
|
||||
covyz = covariance(y,z,n)
|
||||
covzz = covariance(z,z,n)
|
||||
|
||||
SigmaCov = np.array([ [covxx, covxy, covxz], [covxy, covyy, covyz], [covxz, covyz, covzz]])
|
||||
print(SigmaCov)
|
||||
|
||||
EigValues, EigVectors = np.linalg.eig(SigmaCov)
|
||||
# sort eigenvectors and eigenvalues
|
||||
permute = EigValues.argsort()
|
||||
EigValues = EigValues[permute]
|
||||
EigVectors = EigVectors[:,permute]
|
||||
print(EigValues)
|
||||
print(EigVectors)
|
||||
@@ -0,0 +1,38 @@
|
||||
# Simulation of financial transations with or without saving/taxation on transaction
|
||||
# If lambda =0.0, no saving/taxation
|
||||
# See Patriarca et al http://www.sciencedirect.com/science/article/pii/S0378437104004327
|
||||
#!/usr/bin/env python
|
||||
import numpy as np
|
||||
import matplotlib.mlab as mlab
|
||||
import matplotlib.pyplot as plt
|
||||
import random
|
||||
|
||||
# initialize the rng with a seed
|
||||
random.seed()
|
||||
# Hard coding of input parameters
|
||||
Agents = 500
|
||||
MCcounts = 1000
|
||||
Transactions = 100000
|
||||
startMoney = 1.0
|
||||
Lambda = 0.5
|
||||
FinancialAgents = startMoney*np.ones(Agents)
|
||||
for i in range (1, MCcounts, 1):
|
||||
for j in range (1, Transactions, 1):
|
||||
agent_i = int(Agents*random.random())
|
||||
agent_j = int(Agents*random.random())
|
||||
epsilon = random.random()
|
||||
if agent_i != agent_j:
|
||||
m1 = Lambda*FinancialAgents[agent_i] + (1-Lambda)*epsilon*(FinancialAgents[agent_i] + FinancialAgents[agent_j])
|
||||
m2 = Lambda*FinancialAgents[agent_j] + (1-Lambda)*(1-epsilon)*(FinancialAgents[agent_i] + FinancialAgents[agent_j])
|
||||
FinancialAgents[agent_i] = m1
|
||||
FinancialAgents[agent_j] = m2
|
||||
|
||||
# the histogram of the data
|
||||
n, bins, patches = plt.hist(FinancialAgents, 50, facecolor='green')
|
||||
|
||||
plt.xlabel('$x$')
|
||||
plt.ylabel('Distribution of wealth')
|
||||
plt.title(r'Money')
|
||||
plt.axis([0, 10, 0, 100])
|
||||
plt.grid(True)
|
||||
plt.show()
|
||||
@@ -0,0 +1,45 @@
|
||||
# Program to test the Metropolis algorithm with one particle at given temp in
|
||||
# one dimension
|
||||
#!/usr/bin/env python
|
||||
import numpy as np
|
||||
import matplotlib.mlab as mlab
|
||||
import matplotlib.pyplot as plt
|
||||
import random
|
||||
from math import sqrt, exp, log
|
||||
# initialize the rng with a seed
|
||||
random.seed()
|
||||
# Hard coding of input parameters
|
||||
MCcycles = 100000
|
||||
Temperature = 2.0
|
||||
beta = 1./Temperature
|
||||
InitialVelocity = -2.0
|
||||
CurrentVelocity = InitialVelocity
|
||||
Energy = 0.5*InitialVelocity*InitialVelocity
|
||||
VelocityRange = 10*sqrt(Temperature)
|
||||
VelocityStep = 2*VelocityRange/10.
|
||||
AverageEnergy = Energy
|
||||
AverageEnergy2 = Energy*Energy
|
||||
VelocityValues = np.zeros(MCcycles)
|
||||
# The Monte Carlo sampling with Metropolis starts here
|
||||
for i in range (1, MCcycles, 1):
|
||||
TrialVelocity = CurrentVelocity + (2.0*random.random() - 1.0)*VelocityStep
|
||||
EnergyChange = 0.5*(TrialVelocity*TrialVelocity -CurrentVelocity*CurrentVelocity);
|
||||
if random.random() <= exp(-beta*EnergyChange):
|
||||
CurrentVelocity = TrialVelocity
|
||||
Energy += EnergyChange
|
||||
VelocityValues[i] = CurrentVelocity
|
||||
AverageEnergy += Energy
|
||||
AverageEnergy2 += Energy*Energy
|
||||
#Final averages
|
||||
AverageEnergy = AverageEnergy/MCcycles
|
||||
AverageEnergy2 = AverageEnergy2/MCcycles
|
||||
Variance = AverageEnergy2 - AverageEnergy*AverageEnergy
|
||||
print(AverageEnergy, Variance)
|
||||
n, bins, patches = plt.hist(VelocityValues, 400, facecolor='green')
|
||||
|
||||
plt.xlabel('$v$')
|
||||
plt.ylabel('Velocity distribution P(v)')
|
||||
plt.title(r'Velocity histogram at $k_BT=2$')
|
||||
plt.axis([-5, 5, 0, 600])
|
||||
plt.grid(True)
|
||||
plt.show()
|
||||
@@ -0,0 +1,19 @@
|
||||
# ResamplingAnalysisScripts
|
||||
|
||||
## Sample Scripts for data Analysis
|
||||
So far this is a simple python script (should be made parallel...) to perform resampling of a data set. Methods used are __Bootstrapping__, __Jackknife__ and __Blocking__.
|
||||
|
||||
## Usage
|
||||
Simply run `python analysis.py FILENAME.xxx [NLINES]`
|
||||
|
||||
Where `FILENAME` is expected to have a 3 charachter extension `NLINES` (optional) is the number of lines in the file to read and process (default is the whole file, but it gets very slow above 2-3 hundred thousand entries)
|
||||
|
||||
Ouput is located into the `FILENAME/` folder.
|
||||
|
||||
If more than 10⁵ lines are specified the autocorrelation function won't be computed, as it would take too long.
|
||||
|
||||
The `gaussian.dat` dataset has been generated with numpy, as a proof of concept. It represents a normally distributed set of 5x10⁵ elements with `std = 0.05`. One will notice that the estimate on the error of the central value is greatly improved by all resampling methods.
|
||||
|
||||
`energy.dat` is an autocorrelated data set, with autocorrelation time of roughly 200. It is useful to see the use of blocking on this dataset as a convenient method to estimate the autocorrelation time (compare the elapsed time on the different methods).
|
||||
|
||||
In the `plaquette.dat` file there is a small data set (just 1000 samples) and it shows the strenght of using resampling methods to better estimate the error on the central value as opposed to the standard deviation.
|
||||
Executable
+9
@@ -0,0 +1,9 @@
|
||||
8 2
|
||||
0.001 -2.89017 0.00073621
|
||||
0.002 -2.88946 0.00052732
|
||||
0.005 -2.89067 0.00055038
|
||||
0.010 -2.89091 0.00040973
|
||||
0.015 -2.89084 0.00034278
|
||||
0.02 -2.89086 0.00029315
|
||||
0.025 -2.89059 0.00034278
|
||||
0.03 -2.89077 0.00025017
|
||||
Vendored
BIN
Binary file not shown.
Reference in New Issue
Block a user