Wednesday, November 2, 2016

Simple numpy and Cython example

In [1]:
from IPython.core.display import display, HTML
display(HTML("<style>.container { width:50% !important; }</style>"))
In [2]:
%load_ext cython
In [4]:
# Pass compile options to %%cython
# -a means annotate the source code to show sections
# which have Python interaction
In [6]:
%%cython -a
cimport numpy as np


cdef double _dot(double* a, double* b, int n):
    cdef int i
    cdef double sum


    sum = 0
    for i in range(n):
        sum += a[i] * b[i]
    
    return sum


cpdef double cy_dot(np.ndarray a, np.ndarray b):
    cdef int n = a.shape[0]

    cdef double* a_ptr = <double*> a.data
    cdef double* b_ptr = <double*> b.data

    return _dot(a_ptr, b_ptr, n)
Out[6]:
Cython: _cython_magic_3e017d4f9980d56ce613de0b5ec2e832.pyx

Generated by Cython 0.23.4

Yellow lines hint at Python interaction.
Click on a line that starts with a "+" to see the C code that Cython generated for it.

+01: cimport numpy as np
  __pyx_t_1 = PyDict_New(); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
  __Pyx_GOTREF(__pyx_t_1);
  if (PyDict_SetItem(__pyx_d, __pyx_n_s_test, __pyx_t_1) < 0) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 1; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
  __Pyx_DECREF(__pyx_t_1); __pyx_t_1 = 0;
 02: 
 03: 
+04: cdef double _dot(double* a, double* b, int n):
static double __pyx_f_46_cython_magic_3e017d4f9980d56ce613de0b5ec2e832__dot(double *__pyx_v_a, double *__pyx_v_b, int __pyx_v_n) {
  int __pyx_v_i;
  double __pyx_v_sum;
  double __pyx_r;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("_dot", 0);
/* … */
  /* function exit code */
  __pyx_L0:;
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}
 05:     cdef int i
 06:     cdef double sum
 07: 
 08: 
+09:     sum = 0
  __pyx_v_sum = 0.0;
+10:     for i in range(n):
  __pyx_t_1 = __pyx_v_n;
  for (__pyx_t_2 = 0; __pyx_t_2 < __pyx_t_1; __pyx_t_2+=1) {
    __pyx_v_i = __pyx_t_2;
+11:         sum += a[i] * b[i]
    __pyx_v_sum = (__pyx_v_sum + ((__pyx_v_a[__pyx_v_i]) * (__pyx_v_b[__pyx_v_i])));
  }
 12: 
+13:     return sum
  __pyx_r = __pyx_v_sum;
  goto __pyx_L0;
 14: 
 15: 
+16: cpdef double cy_dot(np.ndarray a, np.ndarray b):
static PyObject *__pyx_pw_46_cython_magic_3e017d4f9980d56ce613de0b5ec2e832_1cy_dot(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/
static double __pyx_f_46_cython_magic_3e017d4f9980d56ce613de0b5ec2e832_cy_dot(PyArrayObject *__pyx_v_a, PyArrayObject *__pyx_v_b, CYTHON_UNUSED int __pyx_skip_dispatch) {
  int __pyx_v_n;
  double *__pyx_v_a_ptr;
  double *__pyx_v_b_ptr;
  double __pyx_r;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("cy_dot", 0);
/* … */
  /* function exit code */
  __pyx_L0:;
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}

/* Python wrapper */
static PyObject *__pyx_pw_46_cython_magic_3e017d4f9980d56ce613de0b5ec2e832_1cy_dot(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds); /*proto*/
static PyObject *__pyx_pw_46_cython_magic_3e017d4f9980d56ce613de0b5ec2e832_1cy_dot(PyObject *__pyx_self, PyObject *__pyx_args, PyObject *__pyx_kwds) {
  PyArrayObject *__pyx_v_a = 0;
  PyArrayObject *__pyx_v_b = 0;
  PyObject *__pyx_r = 0;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("cy_dot (wrapper)", 0);
  {
    static PyObject **__pyx_pyargnames[] = {&__pyx_n_s_a,&__pyx_n_s_b,0};
    PyObject* values[2] = {0,0};
    if (unlikely(__pyx_kwds)) {
      Py_ssize_t kw_args;
      const Py_ssize_t pos_args = PyTuple_GET_SIZE(__pyx_args);
      switch (pos_args) {
        case  2: values[1] = PyTuple_GET_ITEM(__pyx_args, 1);
        case  1: values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
        case  0: break;
        default: goto __pyx_L5_argtuple_error;
      }
      kw_args = PyDict_Size(__pyx_kwds);
      switch (pos_args) {
        case  0:
        if (likely((values[0] = PyDict_GetItem(__pyx_kwds, __pyx_n_s_a)) != 0)) kw_args--;
        else goto __pyx_L5_argtuple_error;
        case  1:
        if (likely((values[1] = PyDict_GetItem(__pyx_kwds, __pyx_n_s_b)) != 0)) kw_args--;
        else {
          __Pyx_RaiseArgtupleInvalid("cy_dot", 1, 2, 2, 1); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
        }
      }
      if (unlikely(kw_args > 0)) {
        if (unlikely(__Pyx_ParseOptionalKeywords(__pyx_kwds, __pyx_pyargnames, 0, values, pos_args, "cy_dot") < 0)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
      }
    } else if (PyTuple_GET_SIZE(__pyx_args) != 2) {
      goto __pyx_L5_argtuple_error;
    } else {
      values[0] = PyTuple_GET_ITEM(__pyx_args, 0);
      values[1] = PyTuple_GET_ITEM(__pyx_args, 1);
    }
    __pyx_v_a = ((PyArrayObject *)values[0]);
    __pyx_v_b = ((PyArrayObject *)values[1]);
  }
  goto __pyx_L4_argument_unpacking_done;
  __pyx_L5_argtuple_error:;
  __Pyx_RaiseArgtupleInvalid("cy_dot", 1, 2, 2, PyTuple_GET_SIZE(__pyx_args)); {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L3_error;}
  __pyx_L3_error:;
  __Pyx_AddTraceback("_cython_magic_3e017d4f9980d56ce613de0b5ec2e832.cy_dot", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __Pyx_RefNannyFinishContext();
  return NULL;
  __pyx_L4_argument_unpacking_done:;
  if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_a), __pyx_ptype_5numpy_ndarray, 1, "a", 0))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
  if (unlikely(!__Pyx_ArgTypeTest(((PyObject *)__pyx_v_b), __pyx_ptype_5numpy_ndarray, 1, "b", 0))) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
  __pyx_r = __pyx_pf_46_cython_magic_3e017d4f9980d56ce613de0b5ec2e832_cy_dot(__pyx_self, __pyx_v_a, __pyx_v_b);
  int __pyx_lineno = 0;
  const char *__pyx_filename = NULL;
  int __pyx_clineno = 0;

  /* function exit code */
  goto __pyx_L0;
  __pyx_L1_error:;
  __pyx_r = NULL;
  __pyx_L0:;
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}

static PyObject *__pyx_pf_46_cython_magic_3e017d4f9980d56ce613de0b5ec2e832_cy_dot(CYTHON_UNUSED PyObject *__pyx_self, PyArrayObject *__pyx_v_a, PyArrayObject *__pyx_v_b) {
  PyObject *__pyx_r = NULL;
  __Pyx_RefNannyDeclarations
  __Pyx_RefNannySetupContext("cy_dot", 0);
  __Pyx_XDECREF(__pyx_r);
  __pyx_t_1 = PyFloat_FromDouble(__pyx_f_46_cython_magic_3e017d4f9980d56ce613de0b5ec2e832_cy_dot(__pyx_v_a, __pyx_v_b, 0)); if (unlikely(!__pyx_t_1)) {__pyx_filename = __pyx_f[0]; __pyx_lineno = 16; __pyx_clineno = __LINE__; goto __pyx_L1_error;}
  __Pyx_GOTREF(__pyx_t_1);
  __pyx_r = __pyx_t_1;
  __pyx_t_1 = 0;
  goto __pyx_L0;

  /* function exit code */
  __pyx_L1_error:;
  __Pyx_XDECREF(__pyx_t_1);
  __Pyx_AddTraceback("_cython_magic_3e017d4f9980d56ce613de0b5ec2e832.cy_dot", __pyx_clineno, __pyx_lineno, __pyx_filename);
  __pyx_r = NULL;
  __pyx_L0:;
  __Pyx_XGIVEREF(__pyx_r);
  __Pyx_RefNannyFinishContext();
  return __pyx_r;
}
+17:     cdef int n = a.shape[0]
  __pyx_v_n = (__pyx_v_a->dimensions[0]);
 18: 
+19:     cdef double* a_ptr = <double*> a.data
  __pyx_v_a_ptr = ((double *)__pyx_v_a->data);
+20:     cdef double* b_ptr = <double*> b.data
  __pyx_v_b_ptr = ((double *)__pyx_v_b->data);
 21: 
+22:     return _dot(a_ptr, b_ptr, n)
  __pyx_r = __pyx_f_46_cython_magic_3e017d4f9980d56ce613de0b5ec2e832__dot(__pyx_v_a_ptr, __pyx_v_b_ptr, __pyx_v_n);
  goto __pyx_L0;
In [5]:
import numpy as np
a = np.array([-2., 0., 5.])
b = np.array([1., 1., 1.])
print a.dot(b)
3.0
In [8]:
print cy_dot(a, b)
3.0

Sunday, September 4, 2016

Connecting to your computer from somewhere else

1. The easy way: Tmate.io
2. The harder way, ssh reverse tunneling
https://juntx.wordpress.com/2014/07/28/use-amazon-ec2-and-ssh-reverse-tunneling-to-connect-computers-behind-firewall-or-nat/
http://www.acooke.org/cute/ReverseRem0.html

https://help.ubuntu.com/community/SSH/OpenSSH/PortForwarding
https://thepcspy.com/read/making-ssh-secure/
http://www.howtogeek.com/75007/stupid-geek-tricks-use-your-ssh-config-file-to-create-aliases-for-hosts/

Start remote
From target computer
ssh into ec2
ssh -R 2222:localhost:22 ec2-ip/host-name

ssh into ec2
from ec2 run ssh -p 2222 me@localhost

Wednesday, August 24, 2016

Deep Learning Frameworks The Easy Way (NVIDIA-Docker)

Install drivers and cuda.

Install Docker from the official site (NOT the Ubuntu repos)
sudo apt-get install docker.io DO NOT DO THIS
(What will happen is nvidia-docker will complain during installation that it can't find docker)
see https://github.com/NVIDIA/nvidia-docker/issues/134

Then install nvidia-docker
https://github.com/NVIDIA/nvidia-docker

Finally go to Docker Hub and find a suitable Docker image. e.g. "tensorflow/tensorflow:latest-gpu"

From what I understand these deep learning + gpu docker images are based on Nvidia's official CUDA images at https://hub.docker.com/r/nvidia/cuda

Notes:
Theano + GTX1060 needs the latest version of cuda, v8.0.
The latest binaries for Tensorflow + GTX1060 only work with CUDA v7.5

That means I need two different docker images, one for Tensorflow and one for Theano.

An alternative would be to build Tensorflow from source so that it can use v8.0 but I'm not experienced with that.

Friday, August 19, 2016

The Road to VR

  1. Buy Vive.
  2. Buy computer parts.
  3. Assemble computer.
  4. Install Windows 10.
    1. Flash USB with Windows 10 ISO.
  5. Install motherboard drivers.
    1. Install LAN driver first by copying LAN driver from CD to a flash drive then moving the driver files to the VR PC.
    2. Download and install the rest of the motherboard-related drivers.
      • USB3.0, Audio, Intel Rapid Storage Technology, etc. ... man there's a lot of crap
  6. Install GPU drivers.
  7. Install Steam.
  8. Setup Vive.
  9. FUN!

Sunday, August 14, 2016

"Neuronlike Adaptive Elements That Can Solve Difficult Learning Control Problems"

This is a Python implementation of the method presented in the paper "Neuronlike Adaptive Elements That Can Solve Difficult Learning Control Problems", published in 1983.

The original C implementation written by the authors can be found here.

The "ASE" and "ACE" are neurons which:

  1. have a set of weights which determins the output of the element
    • ASE's weight determine the outputed action
    • ACE's weights determine the corrected reward
  2. update weights according to a
    • reward input
    • exponentially-decaying eligibility trace

In the paper the raw observation vectors are first passed to a decoder which creates the $x$ seen below. $x$ is one-hot encoded.

So even though the below formulas are written in terms of vectors and dot products, anything involving $x$ is really just an index.

e.g. $\mathbf{w}^\intercal \mathbf{x}$

is equivalent to

"Let i be the index where there is a 1 in w. Return x[i]"

Decision function of the ASE:

$ y(t) = \text{sign}(\mathbf{w}^\intercal \mathbf{x} + \text{noise}(t))$

Update rules:

$\mathbf{w}(t+1) = \mathbf{w} + \alpha r(t) \mathbf{e}(t) $

$\mathbf{e}(t+1) = \delta \mathbf{e}(t) + (1-\delta) y(t)\mathbf{x}(t) $

ACE Output (an improved prediction):

$\hat{r}(t) = r(t) + \gamma p(t) - p(t-1)$

where $p$ is the prediction of the reward $p(t) = \mathbf{v}^\intercal \mathbf{x}$

Update rules:

$\mathbf{v}(t+1) = \mathbf{v}(t) + \beta{\hat{r}}\mathbf{\bar{x}}(t)$

$\mathbf{\bar{x}}(t+1) = \lambda\mathbf{\bar{x}}(t) + (1-\lambda) \mathbf{x}$

In [10]:
import numpy as np
In [11]:
import numpy as np
from sklearn.utils import check_random_state

class ASE(object):
    def __init__(self, n_input_dim, learning_rate, decay_rate, sigma=0.01, random_state=None):
        self.n_input_dim = n_input_dim
        self.learning_rate = learning_rate
        self.decay_rate = decay_rate
        self.sigma = sigma
        self.random_state = random_state
        
        self.reset_values()
    
    def reset_values(self):
        # Weights
        self.w = np.zeros(self.n_input_dim)

        # Eligibility
        self.e = np.zeros(self.n_input_dim)
        
        self.rs = check_random_state(self.random_state)
        
    def step(self, x, reward=None):
        '''
        Updates the weights and eligibility trace.
        
        Parameters
        ----------
        'x' : the state vector 
        'reward' : the reward value
        
        Returns
        -------
        Action : 0 or 1
        '''
        # Probabilistic action
        action = self.w.dot(x) + self.rs.randn()*self.sigma
        action = np.clip(action, -50, 50)
        action = int((1 + np.exp(-action)) ** (-1) > 0.5)

        if reward is not None:
            self.w += self.learning_rate * reward * self.e       

            self.e *= self.decay_rate
            self.e += (1. - self.decay_rate) * (action*2 - 1) * x
        
        return action
        
In [20]:
class ACE(object):

    def __init__(self, n_input_dim, learning_rate, decay_rate, discount_factor, random_state=None):
        self.n_input_dim = n_input_dim
        self.learning_rate = learning_rate
        self.decay_rate = decay_rate
        self.discount_factor = discount_factor
        self.random_state = random_state
        
        self.reset_values()
        
    def reset_values(self):
        self.rs = check_random_state(self.random_state)
        self.w = np.zeros(self.n_input_dim)
        self.trace = np.zeros(self.n_input_dim)
        self.prev_p = 0.
        
    def step(self, x, reward=None, done=False):
        if done:
            p = 0.
        else:
            p = self.w.dot(x)
        
        if reward is None:
            self.prev_p = p
            return None
        
        else:
            revised_reward = reward + self.discount_factor * p - self.prev_p
            
            self.w += self.learning_rate * revised_reward * self.trace
            self.trace *= self.decay_rate
            self.trace += (1. - self.decay_rate) * x
            
            self.prev_p = p
            return revised_reward
In [21]:
ONE_DEGREE = 1. * np.pi / 180
SIX_DEGREES = 6. * np.pi / 180
FIFTY_DEGREES = 50. * np.pi / 180

def get_box(observation):
    x, x_dot, theta, theta_dot = observation
    box=0

    bin_edges = [[-0.8, 0.8],
                 [-0.5, 0.5],
                 [-SIX_DEGREES, -ONE_DEGREE, 0, ONE_DEGREE, SIX_DEGREES],
                 [-FIFTY_DEGREES, FIFTY_DEGREES]]
    
    box = 0
    for s, edges in zip(observation, bin_edges):

        i = np.digitize([s], edges)[0]
        box = box * (len(edges)+1) + i
    
    vec = np.zeros(162)
    vec[box] = 1.
    
    return vec
In [22]:
import gym

N_TRIALS = 150
MAX_STEPS = 100000
TERMINATE_ON_MAX_STEPS = True
SEED = 12345


env = gym.make('CartPole-v0')
ase = ASE(n_input_dim=162, learning_rate=1000, decay_rate=0.9, random_state=SEED)
ace = ACE(n_input_dim=162, learning_rate=0.5, decay_rate=0.8, discount_factor=0.95, random_state=SEED)

for trial in range(1, N_TRIALS+1):
    obs = env.reset()
    reward = None
    done = False
    
    for t in range(1, MAX_STEPS+1):
        x = get_box(obs)
        revised_reward = ace.step(x, reward, done)
        #print reward, revised_reward
        action = ase.step(x, revised_reward)
        
        if done:
            break
        
        obs, _, done, _ = env.step(action)

        if done:
            reward = -1
        else:
            reward = 0
    
    if trial % 10 == 0 :
        print "Trial {}: survived {} steps".format(trial, t)
    if TERMINATE_ON_MAX_STEPS and t == MAX_STEPS:
        print "Terminated after {} trials. Successfully balanced pole for MAX_STEPS={}".format(trial, MAX_STEPS)
        break
INFO:gym.envs.registration:Making new env: CartPole-v0
[2016-08-14 18:48:41,473] Making new env: CartPole-v0
Trial 10: survived 23 steps
Trial 20: survived 60 steps
Trial 30: survived 268 steps
Trial 40: survived 815 steps
WARNING:gym.core:Observation '[ 2.40528917  0.38981101  0.01398455  0.26120958]' is not contained within observation space 'Box(4,)'.
[2016-08-14 18:48:44,502] Observation '[ 2.40528917  0.38981101  0.01398455  0.26120958]' is not contained within observation space 'Box(4,)'.
Trial 50: survived 354 steps
Trial 60: survived 1194 steps
Trial 70: survived 657 steps
Trial 80: survived 587 steps
Trial 90: survived 1967 steps
Terminated after 91 trials. Successfully balanced pole for MAX_STEPS=100000

This method succesfully solves Cart-Pole in under 100 iterations!

What if I didn't use any boxes? What if I just used the raw 4-element observation vector??

The formulas are all writen in terms of vectors and dot products so maybe I can skip the state-space, one-hot encoding and just pass the raw observation vector?

In [23]:
import gym

N_TRIALS = 1500
MAX_STEPS = 100000
TERMINATE_ON_MAX_STEPS = True
SEED = 12345


env = gym.make('CartPole-v0')
ase = ASE(n_input_dim=4, learning_rate=10, decay_rate=0.5, random_state=SEED)
ace = ACE(n_input_dim=4, learning_rate=0.01, decay_rate=0.5, discount_factor=0.05, random_state=SEED)

for trial in range(1,N_TRIALS+1):
    obs = env.reset()
    reward = None
    done = False
    
    for t in range(1, MAX_STEPS+1):
        x = obs
        revised_reward = ace.step(x, reward, done)
        #print reward, revised_reward
        action = ase.step(x, revised_reward)
        
        if done:
            break
        
        obs, _, done, _ = env.step(action)

        if done:
            reward = -1
        else:
            reward = 0
            
    if trial % 10 == 0 :
        print "Trial {}: survived {} steps".format(trial, t)
    if TERMINATE_ON_MAX_STEPS and t == MAX_STEPS:
        print "Terminated after {} trials. Successfully balanced pole for MAX_STEPS={}".format(trial, MAX_STEPS)
        break
INFO:gym.envs.registration:Making new env: CartPole-v0
[2016-08-14 18:49:36,056] Making new env: CartPole-v0
Trial 10: survived 64 steps
Trial 20: survived 47 steps
Trial 30: survived 11 steps
Trial 40: survived 32 steps
Trial 50: survived 23 steps
Trial 60: survived 53 steps
Trial 70: survived 33 steps
Trial 80: survived 74 steps
Trial 90: survived 25 steps
Trial 100: survived 25 steps
Trial 110: survived 86 steps
Trial 120: survived 25 steps
Trial 130: survived 31 steps
Trial 140: survived 17 steps
Trial 150: survived 61 steps
Trial 160: survived 54 steps
Trial 170: survived 22 steps
Trial 180: survived 81 steps
Trial 190: survived 17 steps
Trial 200: survived 42 steps
Trial 210: survived 32 steps
Trial 220: survived 13 steps
Trial 230: survived 50 steps
Trial 240: survived 23 steps
Trial 250: survived 18 steps
Trial 260: survived 55 steps
Trial 270: survived 11 steps
Trial 280: survived 49 steps
Trial 290: survived 18 steps
Trial 300: survived 67 steps
Trial 310: survived 43 steps
Trial 320: survived 31 steps
Trial 330: survived 52 steps
Trial 340: survived 47 steps
Trial 350: survived 33 steps
Trial 360: survived 34 steps
Trial 370: survived 30 steps
Trial 380: survived 29 steps
Trial 390: survived 107 steps
Trial 400: survived 11 steps
Trial 410: survived 63 steps
Trial 420: survived 58 steps
Trial 430: survived 84 steps
Trial 440: survived 42 steps
Trial 450: survived 39 steps
Trial 460: survived 32 steps
Trial 470: survived 50 steps
Trial 480: survived 59 steps
Trial 490: survived 33 steps
Trial 500: survived 11 steps
Trial 510: survived 37 steps
Trial 520: survived 11 steps
Trial 530: survived 21 steps
Trial 540: survived 22 steps
Trial 550: survived 87 steps
Trial 560: survived 67 steps
Trial 570: survived 67 steps
Trial 580: survived 44 steps
Trial 590: survived 64 steps
Trial 600: survived 26 steps
Trial 610: survived 11 steps
Trial 620: survived 15 steps
Trial 630: survived 48 steps
Trial 640: survived 39 steps
Trial 650: survived 52 steps
Trial 660: survived 31 steps
Trial 670: survived 72 steps
Trial 680: survived 30 steps
Trial 690: survived 10 steps
Trial 700: survived 49 steps
Trial 710: survived 29 steps
Trial 720: survived 40 steps
Trial 730: survived 44 steps
Trial 740: survived 42 steps
Trial 750: survived 36 steps
Trial 760: survived 39 steps
Trial 770: survived 31 steps
Trial 780: survived 28 steps
Trial 790: survived 43 steps
Trial 800: survived 50 steps
Trial 810: survived 44 steps
Trial 820: survived 75 steps
Trial 830: survived 29 steps
Trial 840: survived 45 steps
Trial 850: survived 31 steps
Trial 860: survived 64 steps
Trial 870: survived 11 steps
Trial 880: survived 29 steps
Trial 890: survived 25 steps
Trial 900: survived 31 steps
Trial 910: survived 31 steps
Trial 920: survived 23 steps
Trial 930: survived 42 steps
Trial 940: survived 33 steps
Trial 950: survived 41 steps
Trial 960: survived 25 steps
Trial 970: survived 45 steps
Trial 980: survived 27 steps
Trial 990: survived 45 steps
Trial 1000: survived 41 steps
Trial 1010: survived 42 steps
Trial 1020: survived 33 steps
Trial 1030: survived 52 steps
Trial 1040: survived 67 steps
Trial 1050: survived 26 steps
Trial 1060: survived 47 steps
Trial 1070: survived 40 steps
Trial 1080: survived 58 steps
Trial 1090: survived 31 steps
Trial 1100: survived 31 steps
Trial 1110: survived 31 steps
Trial 1120: survived 98 steps
Trial 1130: survived 30 steps
Trial 1140: survived 74 steps
Trial 1150: survived 59 steps
Trial 1160: survived 43 steps
Trial 1170: survived 38 steps
Trial 1180: survived 39 steps
Trial 1190: survived 51 steps
Trial 1200: survived 56 steps
Trial 1210: survived 50 steps
Trial 1220: survived 35 steps
Trial 1230: survived 43 steps
Trial 1240: survived 49 steps
Trial 1250: survived 50 steps
Trial 1260: survived 32 steps
Trial 1270: survived 27 steps
Trial 1280: survived 50 steps
Trial 1290: survived 28 steps
Trial 1300: survived 25 steps
Trial 1310: survived 50 steps
Trial 1320: survived 45 steps
Trial 1330: survived 26 steps
Trial 1340: survived 33 steps
Trial 1350: survived 41 steps
Trial 1360: survived 51 steps
Trial 1370: survived 11 steps
Trial 1380: survived 29 steps
Trial 1390: survived 35 steps
Trial 1400: survived 40 steps
Trial 1410: survived 34 steps
Trial 1420: survived 34 steps
Trial 1430: survived 41 steps
Trial 1440: survived 85 steps
Trial 1450: survived 45 steps
Trial 1460: survived 56 steps
Trial 1470: survived 40 steps
Trial 1480: survived 25 steps
Trial 1490: survived 20 steps
Trial 1500: survived 28 steps

Nope. It doesn't work. It seems this method is too strongly reliant on the input vector being a one-hot encoding of discrete states. This is not surprising since the paper was written as a successor to another method which had the name "Boxes".

What about a less sophisticated discretization of the space?

np.linspace

In [4]:
from itertools import product
import numpy as np

import gym

GRANULARITY = 8

env = gym.make('CartPole-v0')

# Split each dimension up into `GRANULARITY` number of buckets
# store these buckets in `linspaces`
bucket_edges = []
ranges = zip(env.observation_space.low, env.observation_space.high)
for low, high in ranges:
    if np.isneginf(low):
        low = -5
    if np.isinf(high):
        high = 5
    bucket_edges.append(np.linspace(low, high, num=GRANULARITY+1, endpoint=True)[1:-1])
INFO:gym.envs.registration:Making new env: CartPole-v0
[2016-08-14 18:35:24,993] Making new env: CartPole-v0
In [6]:
bucket_edges
Out[6]:
[array([-1.8, -1.2, -0.6,  0. ,  0.6,  1.2,  1.8]),
 array([-3.75, -2.5 , -1.25,  0.  ,  1.25,  2.5 ,  3.75]),
 array([-0.31415927, -0.20943951, -0.10471976,  0.        ,  0.10471976,
         0.20943951,  0.31415927]),
 array([-3.75, -2.5 , -1.25,  0.  ,  1.25,  2.5 ,  3.75])]
In [25]:
def get_box(observation):
    x, x_dot, theta, theta_dot = observation
    box=0

    bin_edges = bucket_edges
    
    box = 0
    for s, edges in zip(observation, bin_edges):

        i = np.digitize([s], edges)[0]
        box = box * (len(edges)+1) + i
    
    vec = np.zeros(GRANULARITY**observation.shape[0])
    vec[box] = 1.
    
    return vec
In [26]:
import gym

N_TRIALS = 1000
MAX_STEPS = 100000
TERMINATE_ON_MAX_STEPS = True
SEED = 12345


env = gym.make('CartPole-v0')
ase = ASE(n_input_dim=GRANULARITY**4, learning_rate=1000, decay_rate=0.9, random_state=SEED)
ace = ACE(n_input_dim=GRANULARITY**4, learning_rate=0.5, decay_rate=0.8, discount_factor=0.95, random_state=SEED)

for trial in range(1, N_TRIALS):
    obs = env.reset()
    reward = None
    done = False
    
    for t in range(1, MAX_STEPS+1):
        x = get_box(obs)
        revised_reward = ace.step(x, reward, done)
        #print reward, revised_reward
        action = ase.step(x, revised_reward)
        
        if done:
            break
        
        obs, _, done, _ = env.step(action)

        if done:
            reward = -1
        else:
            reward = 0
    
    if trial % 10 == 0 :
        print "Trial {}: survived {} steps".format(trial, t)
    if TERMINATE_ON_MAX_STEPS and t == MAX_STEPS:
        print "Terminated after {} trials. Successfully balanced pole for MAX_STEPS={}".format(trial, MAX_STEPS)
        break
INFO:gym.envs.registration:Making new env: CartPole-v0
[2016-08-14 18:50:37,159] Making new env: CartPole-v0
Trial 10: survived 10 steps
Trial 20: survived 72 steps
Trial 30: survived 133 steps
Trial 40: survived 129 steps
WARNING:gym.core:Observation '[ 2.4481542   2.44331325 -0.01609566 -0.74575892]' is not contained within observation space 'Box(4,)'.
[2016-08-14 18:50:37,834] Observation '[ 2.4481542   2.44331325 -0.01609566 -0.74575892]' is not contained within observation space 'Box(4,)'.
Trial 50: survived 116 steps
Trial 60: survived 285 steps
Trial 70: survived 11 steps
Trial 80: survived 165 steps
Trial 90: survived 219 steps
Trial 100: survived 191 steps
Trial 110: survived 254 steps
Trial 120: survived 169 steps
Trial 130: survived 196 steps
Trial 140: survived 102 steps
Trial 150: survived 211 steps
Trial 160: survived 152 steps
Trial 170: survived 129 steps
Trial 180: survived 174 steps
Trial 190: survived 299 steps
Trial 200: survived 130 steps
Trial 210: survived 165 steps
Trial 220: survived 184 steps
Trial 230: survived 184 steps
Trial 240: survived 242 steps
Trial 250: survived 419 steps
Trial 260: survived 151 steps
Trial 270: survived 166 steps
Trial 280: survived 372 steps
Trial 290: survived 319 steps
Trial 300: survived 155 steps
Trial 310: survived 165 steps
Trial 320: survived 234 steps
Trial 330: survived 362 steps
Trial 340: survived 143 steps
Trial 350: survived 409 steps
Trial 360: survived 773 steps
Trial 370: survived 405 steps
Trial 380: survived 386 steps
Trial 390: survived 320 steps
Trial 400: survived 393 steps
Trial 410: survived 124 steps
Trial 420: survived 128 steps
Trial 430: survived 359 steps
Trial 440: survived 155 steps
Trial 450: survived 979 steps
Trial 460: survived 229 steps
Trial 470: survived 370 steps
Trial 480: survived 361 steps
Trial 490: survived 1150 steps
Trial 500: survived 168 steps
Trial 510: survived 411 steps
Trial 520: survived 751 steps
Trial 530: survived 187 steps
Trial 540: survived 5854 steps
Terminated after 541 trials. Successfully balanced pole for MAX_STEPS=100000

Their method works with this very rough discretization but (unsurprisingly) it takes much longer to learn.

Saturday, August 6, 2016

Installing Theano on Windows

My build:
OS:Windows 10 x64
GPU: GTX 1060



  1. Install Visual Studio 2013 Community Edition
    • Community, because it's free
    • 2013 because 2015 is not supported by CUDA yet.
  2. Install CUDA 7.5 8.0 (Currenly only CUDA 8.0 RC is available. You have to register with NVIDIA in order to download it.) Installing and testing 8.0 is basically the same as installing 7.5.

    1. Choose Custom during the NVIDA CUDA install. Don't install drivers for the GPU, just install CUDA stuff.
      (I installed the wrong version of CUDA the first time. The GTX 1060 need CUDA 8.)

      My GPU drivers are newer than the ones packaged with this CUDA install. I can ignore this warning.
      Don't install the outdated GPU drivers.
    2. Test CUDA install. Compile one of the samples and run it.
      Open a CUDA sample project .sln, e.g. deviceQuery or bandwidthTest.

      Running the just-built device

  3. Install gcc via the TDM-GCC compiler suite.
    • Install the 64bit version.
    • TDM-GCC is a compiler suite for Windows which consists of components that are all individually available elsewhere. TDM-GCC include GCC and MinGW. TDM-GCC includes command-line tools only.
  4. Install WinPython, one of many scientific Python distributions. (WinPython-64bit-3.4.4.3Qt5.exe)
    • By design, WinPython unzips itself to a local directory but does not mess with the system environment variables.
  5. Create a startup shell script for Python.
  6. Create a link library for gcc? But there is already a .a file in the target location.
  7. Clone and install latest version of Theano. Optional, because Theano is included in WinPy.
    • pip install git+git://github.com/Theano/Theano.git
  8. Create a .theanorc (or .theanorc.txt) for GPU usage. Put it in the settings folder of the WinPython directory.

    .theanorc

    [global]
    device = gpu
    floatX = float32

    [nvcc]
    flags = -LC:\WinPython-64bit-3.4.4.3Qt5\python-3.4.4.amd64\libs
    compiler_binddir=C:\Program Files (x86)\Microsoft Visual Studio 12.0\VC\bin\amd64
  9. Create and run theano-test.py from the command line to see if it works.

    theano-test.py

    import numpy as np
    import time
    import theano
    A = np.random.rand(1000,10000).astype(theano.config.floatX)
    B = np.random.rand(10000,1000).astype(theano.config.floatX)
    np_start = time.time()
    AB = A.dot(B)
    np_end = time.time()
    X,Y = theano.tensor.matrices('XY')
    mf = theano.function([X,Y],X.dot(Y))
    t_start = time.time()
    tAB = mf(A,B)
    t_end = time.time()
    print("NP time: %f[s], theano time: %f[s] (times should be close when run on CPU!)" %(
                                               np_end-np_start, t_end-t_start))
    print("Result difference: %f" % (np.abs(AB-tAB).max(), ))
    Note: The first time it runs, there will be a lot of debug output as theano compiles the code. Subsequent runs are much quieter.