Implemented HMM Archetype (#44)

* Fixed Swap Nodes in mutate_rows per issue #41

Swapped based on state, row[0], then sorted.  As per our conversation in issue #41 thread.

* Updated for fixed swap in mutation

* Moved HMMParams class and implemented tests.

* Implemented PSO for HMM.

* Implemented PSO for HMM.

* Deleted CSVs from bin.

* Changed mutation_rate to mutation_probability in test.

* Added tests/docstrings to HMM, and fixed bugs in HMM.

* Fixed typo

* Changed test strategies from hard-to-pickle strategies to basic strategies.

* Made read_vector, create_vector_bounds, and vector_to_instance more instance-specific

* Made vector changes to FSM too.

* Updated test_hmm for vector_to_instance change.

* Updated test_gambler for vector_to_instance change.

* Fixed test_vector_to_instance in test_fsm.

* Made minor formatting fixes.

* Put conversion back in.

* Cleaned up docstrings.

Author: gaffney2010

bin/hmm_evolve.py

@@ -2,11 +2,11 @@
 Hidden Markov Model Evolver
 
 Usage:
-    fsm_evolve.py [-h] [--generations GENERATIONS] [--population POPULATION]
+    hmm_evolve.py [-h] [--generations GENERATIONS] [--population POPULATION]
     [--mu MUTATION_RATE] [--bottleneck BOTTLENECK] [--processes PROCESSORS]
     [--output OUTPUT_FILE] [--objective OBJECTIVE] [--repetitions REPETITIONS]
     [--turns TURNS] [--noise NOISE] [--nmoran NMORAN]
-    [--states NUM_STATES]
+    [--states NUM_STATES] [--algorithm ALGORITHM]
 
 Options:
     -h --help                   Show help
@@ -22,190 +22,14 @@
     --noise NOISE               Match noise [default: 0.00]
     --nmoran NMORAN             Moran Population Size, if Moran objective [default: 4]
     --states NUM_STATES         Number of FSM states [default: 5]
+    --algorithm ALGORITHM       Which algorithm to use (EA for evolutionary algorithm or PS for
+                                particle swarm algorithm) [default: EA]
 """
 
-#####
-# This is a potential candidate for PSO optimization, which will require
-# combining the matrices.
-#####
-
-import random
-from random import randrange, choice
-
 from docopt import docopt
-import numpy as np
-
-from axelrod import Action
-from axelrod.strategies.hmm import HMMPlayer
-from axelrod_dojo import Params, Population, prepare_objective
-
-C, D = Action.C, Action.D
-
-
-def copy_lists(rows):
-    new_rows = list(map(list, rows))
-    return new_rows
-
-def random_vector(size):
-    """Create a random vector of values in [0, 1] that sums to 1."""
-    vector = []
-    s = 1
-    for _ in range(size - 1):
-        r = s * random.random()
-        vector.append(r)
-        s -= r
-    vector.append(s)
-    return vector
-
-def normalize_vector(vec):
-    s = sum(vec)
-    vec = [v / s for v in vec]
-    return vec
-
-def mutate_row(row, mutation_rate):
-    randoms = np.random.random(len(row))
-    for i in range(len(row)):
-        if randoms[i] < mutation_rate:
-            ep = random.uniform(-1, 1) / 4
-            row[i] += ep
-            if row[i] < 0:
-                row[i] = 0
-            if row[i] > 1:
-                row[i] = 1
-    return row
-
-
-class HMMParams(Params):
-
-    def __init__(self, num_states, mutation_rate=None, transitions_C=None,
-                 transitions_D=None, emission_probabilities=None,
-                 initial_state=0, initial_action=C):
-        self.PlayerClass = HMMPlayer
-        self.num_states = num_states
-        if mutation_rate is None:
-            self.mutation_rate = 10 / (num_states ** 2)
-        else:
-            self.mutation_rate = mutation_rate
-        if transitions_C is None:
-            self.randomize()
-        else:
-            # Make sure to copy the lists
-            self.transitions_C = copy_lists(transitions_C)
-            self.transitions_D = copy_lists(transitions_D)
-            self.emission_probabilities = list(emission_probabilities)
-            self.initial_state = initial_state
-            self.initial_action = initial_action
-
 
-    def player(self):
-        player = self.PlayerClass(self.transitions_C, self.transitions_D,
-                                  self.emission_probabilities,
-                                  self.initial_state, self.initial_action)
-        return player
-
-    def copy(self):
-        return HMMParams(self.num_states, self.mutation_rate,
-                         self.transitions_C, self.transitions_D,
-                         self.emission_probabilities,
-                         self.initial_state, self.initial_action)
-
-    @staticmethod
-    def random_params(num_states):
-        t_C = []
-        t_D = []
-        for _ in range(num_states):
-            t_C.append(random_vector(num_states))
-            t_D.append(random_vector(num_states))
-        initial_state = randrange(num_states)
-        # initial_action = choice([C, D])
-        initial_action = C
-        return t_C, t_D, initial_state, initial_action
-
-    def randomize(self):
-        t_C, t_D, initial_state, initial_action = self.random_params(self.num_states)
-        self.emission_probabilities = [random.random() for _ in range(self.num_states)]
-        self.transitions_C = t_C
-        self.transitions_D = t_D
-        self.initial_state = initial_state
-        self.initial_action = initial_action
-
-    @staticmethod
-    def mutate_rows(rows, mutation_rate):
-        for i, row in enumerate(rows):
-            row = mutate_row(row, mutation_rate)
-            rows[i] = normalize_vector(row)
-        return rows
-
-    def mutate(self):
-        self.transitions_C = self.mutate_rows(
-            self.transitions_C, self.mutation_rate)
-        self.transitions_D = self.mutate_rows(
-            self.transitions_D, self.mutation_rate)
-        self.emission_probabilities = mutate_row(
-            self.emission_probabilities, self.mutation_rate)
-        if random.random() < self.mutation_rate / 10:
-            self.initial_action = self.initial_action.flip()
-        if random.random() < self.mutation_rate / (10 * self.num_states):
-            self.initial_state = randrange(self.num_states)
-        # Change node size?
-
-    @staticmethod
-    def crossover_rows(rows1, rows2):
-        num_states = len(rows1)
-        crosspoint = randrange(num_states)
-        new_rows = copy_lists(rows1[:crosspoint])
-        new_rows += copy_lists(rows2[crosspoint:])
-        return new_rows
-
-    @staticmethod
-    def crossover_weights(w1, w2):
-        crosspoint = random.randrange(len(w1))
-        new_weights = list(w1[:crosspoint]) + list(w2[crosspoint:])
-        return new_weights
-
-    def crossover(self, other):
-        # Assuming that the number of states is the same
-        t_C = self.crossover_rows(self.transitions_C, other.transitions_C)
-        t_D = self.crossover_rows(self.transitions_D, other.transitions_D)
-        emissions = self.crossover_weights(
-            self.emission_probabilities, other.emission_probabilities)
-        return HMMParams(self.num_states, self.mutation_rate,
-                         t_C, t_D, emissions,
-                         self.initial_state, self.initial_action)
-
-    @staticmethod
-    def repr_rows(rows):
-        ss = []
-        for row in rows:
-            ss.append("_".join(list(map(str, row))))
-        return "|".join(ss)
-
-    def __repr__(self):
-        return "{}:{}:{}:{}:{}".format(
-            self.initial_state,
-            self.initial_action,
-            self.repr_rows(self.transitions_C),
-            self.repr_rows(self.transitions_D),
-            self.repr_rows([self.emission_probabilities])
-        )
-
-    @classmethod
-    def parse_repr(cls, s):
-        def parse_matrix(sm):
-            rows = []
-            lines = sm.split('|')
-            for line in lines:
-                row = line.split('_')
-                row = list(map(float, row))
-                rows.append(row)
-            return row
-        lines = s.split(':')
-        initial_state = int(lines[0])
-        initial_action = lines[1]
-        t_C = parse_matrix(lines[2])
-        t_D = parse_matrix(lines[3])
-        ps = parse_matrix(lines[4])
-        return cls(t_C, t_D, ps, initial_state, initial_action)
+from axelrod_dojo import HMMParams, Population, prepare_objective
+from axelrod_dojo.algorithms.particle_swarm_optimization import PSO
 
 
 if __name__ == '__main__':
@@ -227,13 +51,33 @@ def parse_matrix(sm):
     noise = float(arguments['--noise'])
     nmoran = int(arguments['--nmoran'])
 
-    # FSM
+    # HMM
     num_states = int(arguments['--states'])
-    param_kwargs = {"num_states": num_states}
-
-    objective = prepare_objective(name, turns, noise, repetitions, nmoran)
-    population = Population(HMMParams, param_kwargs, population, objective,
-                            output_filename, bottleneck,
-                            mutation_probability,
-                            processes=processes)
-    population.run(generations)
+    params_kwargs = {"num_states": num_states}
+
+    if arguments['--algorithm'] == "PS":
+        objective = prepare_objective(name, turns, noise, repetitions, nmoran)
+        pso = PSO(HMMParams, params_kwargs, objective=objective,
+                  population=population, generations=generations,
+                  size=num_states)
+
+        xopt_helper, fopt = pso.swarm()
+        xopt = HMMParams(num_states=num_states)
+        xopt.read_vector(xopt_helper, num_states)
+    else:
+        objective = prepare_objective(name, turns, noise, repetitions, nmoran)
+        population = Population(HMMParams, params_kwargs, population, objective,
+                                output_filename, bottleneck, mutation_probability,
+                                processes=processes)
+        population.run(generations)
+        
+        # Get the best member of the population to output.
+        scores = population.score_all()
+        record, record_holder = 0, -1
+        for i, s in enumerate(scores):
+            if s >= record:
+                record = s
+                record_holder = i
+        xopt, fopt = population.population[record_holder], record
+    
+    print("Best Score: {} {}".format(fopt, xopt))

bin/pso_evolve.py

@@ -54,13 +54,13 @@
     plays = int(arguments['--plays'])
     op_plays = int(arguments['--op_plays'])
     op_start_plays = int(arguments['--op_start_plays'])
-    param_kwargs = {"plays": plays,
+    params_kwargs = {"plays": plays,
                     "op_plays": op_plays,
                     "op_state_plays": op_start_plays}
 
     objective = prepare_objective(name, turns, noise, repetitions, nmoran)
 
-    pso = PSO(GamblerParams, param_kwargs, objective=objective,
+    pso = PSO(GamblerParams, params_kwargs, objective=objective,
               population=population, generations=generations)
 
     xopt, fopt = pso.swarm()

docs/background/genetic_algorithm.rst

@@ -51,9 +51,44 @@ The crossover and mutation are implemented in the following way:
 - Crossover: this is done by taking a randomly selected number of target
   state/actions
   pairs from one individual and the rest from the other.
-- Mutation: given a mutation probability :math:`delta` each target state/action
+- Mutation: given a mutation probability :math:`\delta` each target state/action
   has a probability :math:`\delta` of being randomly changed to one of the other
   states or actions. Furthermore the **initial** action has a probability of
   being swapped of :math:`\delta\times 10^{-1}` and the **initial** state has a
   probability of being changed to another random state of :math:`\delta \times
   10^{-1} \times N` (where :math:`N` is the number of states).
+
+Hidden Markov models
+---------------------
+
+A hidden Markov model is made up of the following:
+
+- a mapping from a state/action pair to a probability of defect or cooperation.
+- a cooperation transition matrix, the probability of transitioning to each
+  state, given current state and an opponent cooperation.
+- a defection transition matrix, the probability of transitioning to each
+  state, given current state and an opponent defection.
+- an initial state/action pair.
+
+(See [Harper2017]_ for more details.)
+
+The crossover and mutation are implemented in the following way:
+
+- Crossover: this is done by taking a randomly selected number of rows from
+  one cooperation transition matrix and the rest from the other to form a target
+  cooperation transition matrix; then a different number of randomly selected
+  rows from one defection transition matrix and the rest from the other; and
+  then a randomly select number of entries from one state/part -> probability
+  mapping and the rest from the other.
+- Mutation: given a mutation probability :math:`delta` each cell of both
+  transition matrices and the state/part -> probability mapping have probability
+  :math:`delta` of being increased by :math:`varepsilon`, where
+  :math:`varepsilon` is randomly drawn uniformly from :math:`[-0.25, 0.25]`
+  (A negative number would decrease.)  Then the transition matrices and mapping
+  are adjusted so that no cell is outside :math:`[0, 1]` and the transition
+  matrices are normalized so that each row adds to 1. Furthermore the
+  **initial** action has a probability of being swapped of
+  :math:`\delta\times 10^{-1}` and the **initial** state has a probability of
+  being changed to another random state of
+  :math:`\delta \times 10^{-1} \times N` (where :math:`N` is the number of
+  states).

src/axelrod_dojo/__init__.py

@@ -1,7 +1,9 @@
 from .version import __version__
 from .archetypes.fsm import FSMParams
+from .archetypes.hmm import HMMParams
 from .archetypes.gambler import GamblerParams
 from .algorithms.genetic_algorithm import Population
+from .algorithms.particle_swarm_optimization import PSO
 from .utils import (prepare_objective,
                     load_params,
                     Params,

src/axelrod_dojo/algorithms/particle_swarm_optimization.py

@@ -41,7 +41,7 @@ def swarm(self):
 
         def objective_function(vector):
             params.receive_vector(vector=vector)
-            instance_generation_function = 'vector_to_instance'
+            instance_generation_function = 'player'
 
             return - score_params(params=params, objective=self.objective,
                                   opponents_information=self.opponents_information,

src/axelrod_dojo/archetypes/fsm.py

@@ -140,34 +140,30 @@ def parse_repr(cls, s):
         return cls(num_states, rows, initial_state, initial_action)
 
     def receive_vector(self, vector):
-        """Receives a vector and creates an instance attribute called
-        vector."""
-        self.vector = vector
-
-    def vector_to_instance(self):
-        """Turns the attribute vector in to a FSM player instance.
+        """
+        Read a serialized vector into the set of FSM parameters (less initial
+        state).  Then assign those FSM parameters to this class instance.
 
         The vector has three parts. The first is used to define the next state 
         (for each of the player's states - for each opponents action).
 
         The second part is the player's next moves (for each state - for 
         each opponent's actions).
 
-        Finally, a probability to determine the player's first move."""
-
-        num_states = int((len(self.vector) - 1) / 4)
-        state_scale = self.vector[:num_states * 2]
-        next_states = [int(s * (num_states - 1)) for s in state_scale]
-        actions = self.vector[num_states * 2: -1]
-        starting_move = C if round(self.vector[-1]) == 0 else D
+        Finally, a probability to determine the player's first move.
+        """
+        state_scale = vector[:self.num_states * 2]
+        next_states = [int(s * (self.num_states - 1)) for s in state_scale]
+        actions = vector[self.num_states * 2: -1]
+        
+        self.initial_action = C if round(vector[-1]) == 0 else D
+        self.initial_state = 1
 
-        fsm = []
+        self.rows = []
         for i, (initial_state, action) in enumerate(
-                itertools.product(range(num_states), [C, D])):
+                itertools.product(range(self.num_states), [C, D])):
             next_action = C if round(actions[i]) == 0 else D
-            fsm.append([initial_state, action, next_states[i], next_action])
-
-        return FSMPlayer(fsm, initial_action=starting_move)
+            self.rows.append([initial_state, action, next_states[i], next_action])
 
     def create_vector_bounds(self):
         """Creates the bounds for the decision variables."""