# This file shows how to implement a single hidden layer neural network for
# performing binary classification on the GPU using cudamat.

from __future__ import division
import pdb
import time
import numpy as np
import cudamat as cm
from cudamat import learn as cl
import util

# initialize CUDA
cm.cublas_init()

# load data
util.load('mnist49.dat', globals())

# Put training data onto the GPU.
dat_train = dat_train/255.
dat_train = dat_train - (np.mean(dat_train, 1)+10**-8)[:, np.newaxis]
dev_train = cm.CUDAMatrix(dat_train)
dev_lbl = cm.CUDAMatrix(lbl_train)

# training parameters
epsilon = 0.01
momentum = 0.9

num_epochs = 30
batch_size = 128
num_batches = dat_train.shape[1]//batch_size

# model parameters
dim_in = dat_train.shape[0]
dim_out = 1
num_hid = 1024

# initialize weights
w_w1 = cm.CUDAMatrix(dim_in ** -0.5 * np.random.randn(dim_in, num_hid))
w_b1 = cm.CUDAMatrix(np.zeros((num_hid, 1)))
w_w2 = cm.CUDAMatrix(num_hid ** -0.5 * np.random.randn(num_hid, dim_out))
w_b2 = cm.CUDAMatrix(np.zeros((dim_out, 1)))

# initialize weight update matrices
wu_w1 = cm.empty(w_w1.shape).assign(0)
wu_b1 = cm.empty(w_b1.shape).assign(0)
wu_w2 = cm.empty(w_w2.shape).assign(0)
wu_b2 = cm.empty(w_b2.shape).assign(0)

# initialize temporary storage
h = cm.empty((num_hid, batch_size))
out = cm.empty((dim_out, batch_size))
delta = cm.empty((num_hid, batch_size))

# Train neural network.
start_time = time.time()
for epoch in range(num_epochs):
    print("Epoch %i" % (epoch + 1))
    err = []

    for batch in range(num_batches):
        # get current minibatch
        inp = dev_train.slice(batch*batch_size,(batch + 1)*batch_size)
        target = dev_lbl.slice(batch*batch_size,(batch + 1)*batch_size)

        # forward pass
        cm.dot(w_w1.T, inp, target = h)

        h.add_col_vec(w_b1)
        h.apply_sigmoid()

        cm.dot(w_w2.T, h, target = out)

        out.add_col_vec(w_b2)
        out.apply_sigmoid()

        # back prop errors
        out.subtract(target) # compute error

        # gradients for w_w2 and w_b2
        wu_w2.add_dot(h, out.T, beta = momentum)
        wu_b2.add_sums(out, axis = 1, beta = momentum)

        # compute delta
        cm.dot(w_w2, out, target = delta)

        # delta = delta * h * (1 - h)
        cl.mult_by_sigmoid_deriv(delta, h)

        # gradients for w_w1 and w_b1
        wu_w1.add_dot(inp, delta.T, beta = momentum)
        wu_b1.add_sums(delta, axis = 1, beta = momentum)

        # update weights
        w_w1.subtract_mult(wu_w1, epsilon/batch_size)
        w_b1.subtract_mult(wu_b1, epsilon/batch_size)
        w_w2.subtract_mult(wu_w2, epsilon/batch_size)
        w_b2.subtract_mult(wu_b2, epsilon/batch_size)

        # calculate error on current minibatch 
        err.append(np.abs(out.asarray())>0.5)

    print("Training misclassification rate: %f" % np.mean(err))
    print("Time: %f" % (time.time() - start_time))

# Evaluate neural network on test data.

# Load test data onto the GPU.
dat_test = dat_test/255.
dat_test = dat_test - np.mean(dat_test, 1)[:, np.newaxis]
dev_test = cm.CUDAMatrix(dat_test)
dev_lbl = cm.CUDAMatrix(lbl_test)

# Initalize temporary storage.
h = cm.empty((num_hid, dat_test.shape[1]))
out = cm.empty((dim_out, dat_test.shape[1]))

# forward pass
cm.dot(w_w1.T, dev_test, target = h)

h.add_col_vec(w_b1)
h.apply_sigmoid()

cm.dot(w_w2.T, h, target = out)

out.add_col_vec(w_b2)
out.apply_sigmoid()

# compute error
out.subtract(dev_lbl)

print("Testing misclassification rate: %f" % np.mean(np.abs(out.asarray())>0.5))

cm.cublas_shutdown()

[((296, 312), 'cudamat.cublas_init', 'cm.cublas_init', ([], {}), '()\n', (310, 312), True, 'import cudamat as cm\n'), ((506, 530), 'cudamat.CUDAMatrix', 'cm.CUDAMatrix', (['dat_train'], {}), '(dat_train)\n', (519, 530), True, 'import cudamat as cm\n'), ((541, 565), 'cudamat.CUDAMatrix', 'cm.CUDAMatrix', (['lbl_train'], {}), '(lbl_train)\n', (554, 565), True, 'import cudamat as cm\n'), ((1260, 1291), 'cudamat.empty', 'cm.empty', (['(num_hid, batch_size)'], {}), '((num_hid, batch_size))\n', (1268, 1291), True, 'import cudamat as cm\n'), ((1298, 1329), 'cudamat.empty', 'cm.empty', (['(dim_out, batch_size)'], {}), '((dim_out, batch_size))\n', (1306, 1329), True, 'import cudamat as cm\n'), ((1338, 1369), 'cudamat.empty', 'cm.empty', (['(num_hid, batch_size)'], {}), '((num_hid, batch_size))\n', (1346, 1369), True, 'import cudamat as cm\n'), ((1408, 1419), 'time.time', 'time.time', ([], {}), '()\n', (1417, 1419), False, 'import time\n'), ((3067, 3090), 'cudamat.CUDAMatrix', 'cm.CUDAMatrix', (['dat_test'], {}), '(dat_test)\n', (3080, 3090), True, 'import cudamat as cm\n'), ((3101, 3124), 'cudamat.CUDAMatrix', 'cm.CUDAMatrix', (['lbl_test'], {}), '(lbl_test)\n', (3114, 3124), True, 'import cudamat as cm\n'), ((3161, 3199), 'cudamat.empty', 'cm.empty', (['(num_hid, dat_test.shape[1])'], {}), '((num_hid, dat_test.shape[1]))\n', (3169, 3199), True, 'import cudamat as cm\n'), ((3206, 3244), 'cudamat.empty', 'cm.empty', (['(dim_out, dat_test.shape[1])'], {}), '((dim_out, dat_test.shape[1]))\n', (3214, 3244), True, 'import cudamat as cm\n'), ((3261, 3295), 'cudamat.dot', 'cm.dot', (['w_w1.T', 'dev_test'], {'target': 'h'}), '(w_w1.T, dev_test, target=h)\n', (3267, 3295), True, 'import cudamat as cm\n'), ((3338, 3367), 'cudamat.dot', 'cm.dot', (['w_w2.T', 'h'], {'target': 'out'}), '(w_w2.T, h, target=out)\n', (3344, 3367), True, 'import cudamat as cm\n'), ((3535, 3555), 'cudamat.cublas_shutdown', 'cm.cublas_shutdown', ([], {}), '()\n', (3553, 3555), True, 'import cudamat as cm\n'), ((888, 910), 'numpy.zeros', 'np.zeros', (['(num_hid, 1)'], {}), '((num_hid, 1))\n', (896, 910), True, 'import numpy as np\n'), ((1007, 1029), 'numpy.zeros', 'np.zeros', (['(dim_out, 1)'], {}), '((dim_out, 1))\n', (1015, 1029), True, 'import numpy as np\n'), ((833, 865), 'numpy.random.randn', 'np.random.randn', (['dim_in', 'num_hid'], {}), '(dim_in, num_hid)\n', (848, 865), True, 'import numpy as np\n'), ((951, 984), 'numpy.random.randn', 'np.random.randn', (['num_hid', 'dim_out'], {}), '(num_hid, dim_out)\n', (966, 984), True, 'import numpy as np\n'), ((1076, 1096), 'cudamat.empty', 'cm.empty', (['w_w1.shape'], {}), '(w_w1.shape)\n', (1084, 1096), True, 'import cudamat as cm\n'), ((1115, 1135), 'cudamat.empty', 'cm.empty', (['w_b1.shape'], {}), '(w_b1.shape)\n', (1123, 1135), True, 'import cudamat as cm\n'), ((1154, 1174), 'cudamat.empty', 'cm.empty', (['w_w2.shape'], {}), '(w_w2.shape)\n', (1162, 1174), True, 'import cudamat as cm\n'), ((1193, 1213), 'cudamat.empty', 'cm.empty', (['w_b2.shape'], {}), '(w_b2.shape)\n', (1201, 1213), True, 'import cudamat as cm\n'), ((1746, 1775), 'cudamat.dot', 'cm.dot', (['w_w1.T', 'inp'], {'target': 'h'}), '(w_w1.T, inp, target=h)\n', (1752, 1775), True, 'import cudamat as cm\n'), ((1842, 1871), 'cudamat.dot', 'cm.dot', (['w_w2.T', 'h'], {'target': 'out'}), '(w_w2.T, h, target=out)\n', (1848, 1871), True, 'import cudamat as cm\n'), ((2182, 2213), 'cudamat.dot', 'cm.dot', (['w_w2', 'out'], {'target': 'delta'}), '(w_w2, out, target=delta)\n', (2188, 2213), True, 'import cudamat as cm\n'), ((2263, 2297), 'cudamat.learn.mult_by_sigmoid_deriv', 'cl.mult_by_sigmoid_deriv', (['delta', 'h'], {}), '(delta, h)\n', (2287, 2297), True, 'from cudamat import learn as cl\n'), ((3020, 3040), 'numpy.mean', 'np.mean', (['dat_test', '(1)'], {}), '(dat_test, 1)\n', (3027, 3040), True, 'import numpy as np\n'), ((449, 470), 'numpy.mean', 'np.mean', (['dat_train', '(1)'], {}), '(dat_train, 1)\n', (456, 470), True, 'import numpy as np\n'), ((2835, 2847), 'numpy.mean', 'np.mean', (['err'], {}), '(err)\n', (2842, 2847), True, 'import numpy as np\n'), ((2873, 2884), 'time.time', 'time.time', ([], {}), '()\n', (2882, 2884), False, 'import time\n')]