From 3f71bae3d531f953928fb6d06cbdcd839bcfe7c3 Mon Sep 17 00:00:00 2001
From: Michael Hush <MichaelRHush@gmail.com>
Date: Wed, 5 Oct 2016 15:39:34 +1100
Subject: [PATCH 1/5] Adding differential evolution

Started modifying the differential evolution code from scipy for use in
M-LOOP.
---
 mloop/learners.py | 499 +++++++++++++++++++++++++++++++++++++++++++++++++++++-
 1 file changed, 498 insertions(+), 1 deletion(-)

diff --git a/mloop/learners.py b/mloop/learners.py
index 03fdc28..7cfdaf9 100644
--- a/mloop/learners.py
+++ b/mloop/learners.py
@@ -1178,7 +1178,504 @@ def find_local_minima(self):
         
         self.has_local_minima = True
         self.log.info('Search completed')
-        
+
     
+class DifferentialEvolutionLearner(Learner, threading.Thread):
+
+    """
+    Adaption of the differential evolution algorithm in scipy. 
     
+    Args:
+        params_out_queue (queue): Queue for parameters sent to controller.
+        costs_in_queue (queue): Queue for costs for gaussian process. This must be tuple
+        end_event (event): Event to trigger end of learner.
+        
+    Keyword Args:
+        evolution_strategy (Optional [string]): the differential evolution strategy to use, options are 'best1', 'best1', 'rand1' and 'rand2'. The default is 'rand2'.
+        population_size (Optional [int]): multiplier proportional to the number of parameters in a generation. The generation population is set to population_size * parameter_num. Default 15.
+        mutation_scale (Optional [tuple]): The mutation scale when picking new points. Otherwise known as differential weight. When provided as a tuple (min,max) a mutation constant is picked randomly in the interval. Default (0.5,1.0).
+        crossover_probability (Optional [float]): The recombination constand or crossover probability, the probability a new points will be added to the population.
+        restart_tolerance (Optional [float]): when the current population have a spread less than the initial tolerance, namely stdev(curr_pop) < restart_tolerance stdev(init_pop), it is likely the population is now in a minima, and so the search is started again.
+        trust_region (Optional [float or array]): The trust region defines the maximum distance the learner will travel from the current best set of parameters. If None, the learner will search everywhere. If a float, this number must be between 0 and 1 and defines maximum distance the learner will venture as a percentage of the boundaries. If it is an array, it must have the same size as the number of parameters and the numbers define the maximum absolute distance that can be moved along each direction. 
+        
+    Attributes:
+        has_trust_region (bool): Whether the learner has a trust region. 
+        
+    """
+
+    # Dispatch of mutation strategy method (binomial or exponential).
+    
+    def __init__(self, 
+                 strategy='rand2', 
+                 population_size=15,
+                 restart_tolerance=0.01, 
+                 mutation=(0.5, 1), 
+                 recombination=0.7, 
+                 **kwargs):
+        
+        super(NelderMeadLearner,self).__init__(**kwargs)
+        
+        if strategy == 'best1':
+            self.mutation_func = self._best1
+        elif strategy == 'best2':
+            self.mutation_func = self._best2
+        elif strategy == 'rand1':
+            self.mutation_func = self._rand1
+        elif strategy == 'rand2':
+            self.mutation_func = self._rand2
+        else:
+            self.log.error("Please select a valid mutation strategy")
+            raise ValueError
+        
+        self.strategy = strategy
+        self.restart_tolerance = restart_tolerance
+
+        # Mutation constant should be in [0, 2). If specified as a sequence
+        # then dithering is performed.
+        self.scale = mutation
+        if (not np.all(np.isfinite(mutation)) or
+                np.any(np.array(mutation) >= 2) or
+                np.any(np.array(mutation) < 0)):
+            raise ValueError('The mutation constant must be a float in '
+                             'U[0, 2), or specified as a tuple(min, max)'
+                             ' where min < max and min, max are in U[0, 2).')
+
+        self.dither = None
+        if hasattr(mutation, '__iter__') and len(mutation) > 1:
+            self.dither = [mutation[0], mutation[1]]
+            self.dither.sort()
+
+        self.cross_over_probability = recombination
+
+        self.func = func
+        self.args = args
+
+        # convert tuple of lower and upper bounds to limits
+        # [(low_0, high_0), ..., (low_n, high_n]
+        #     -> [[low_0, ..., low_n], [high_0, ..., high_n]]
+        self.limits = np.array(bounds, dtype='float').T
+        if (np.size(self.limits, 0) != 2 or not
+                np.all(np.isfinite(self.limits))):
+            raise ValueError('bounds should be a sequence containing '
+                             'real valued (min, max) pairs for each value'
+                             ' in x')
+
+        if maxiter is None:  # the default used to be None
+            maxiter = 1000
+        self.maxiter = maxiter
+        if maxfun is None:  # the default used to be None
+            maxfun = np.inf
+        self.maxfun = maxfun
+
+        # population is scaled to between [0, 1].
+        # We have to scale between parameter <-> population
+        # save these arguments for _scale_parameter and
+        # _unscale_parameter. This is an optimization
+        self.__scale_arg1 = 0.5 * (self.limits[0] + self.limits[1])
+        self.__scale_arg2 = np.fabs(self.limits[0] - self.limits[1])
+
+        self.parameter_count = np.size(self.limits, 1)
+
+        self.random_number_generator = _make_random_gen(seed)
+
+        # default population initialization is a latin hypercube design, but
+        # there are other population initializations possible.
+        self.num_population_members = popsize * self.parameter_count
+
+        self.population_shape = (self.num_population_members,
+                                 self.parameter_count)
+
+        self._nfev = 0
+        if init == 'latinhypercube':
+            self.init_population_lhs()
+        elif init == 'random':
+            self.init_population_random()
+        else:
+            raise ValueError("The population initialization method must be one"
+                             "of 'latinhypercube' or 'random'")
+
+        self.disp = disp
+
+    def init_population_lhs(self):
+        """
+        Initializes the population with Latin Hypercube Sampling.
+        Latin Hypercube Sampling ensures that each parameter is uniformly
+        sampled over its range.
+        """
+        rng = self.random_number_generator
+
+        # Each parameter range needs to be sampled uniformly. The scaled
+        # parameter range ([0, 1)) needs to be split into
+        # `self.num_population_members` segments, each of which has the following
+        # size:
+        segsize = 1.0 / self.num_population_members
+
+        # Within each segment we sample from a uniform random distribution.
+        # We need to do this sampling for each parameter.
+        samples = (segsize * rng.random_sample(self.population_shape)
+
+        # Offset each segment to cover the entire parameter range [0, 1)
+                   + np.linspace(0., 1., self.num_population_members,
+                                 endpoint=False)[:, np.newaxis])
+
+        # Create an array for population of candidate solutions.
+        self.population = np.zeros_like(samples)
+
+        # Initialize population of candidate solutions by permutation of the
+        # random samples.
+        for j in range(self.parameter_count):
+            order = rng.permutation(range(self.num_population_members))
+            self.population[:, j] = samples[order, j]
+
+        # reset population energies
+        self.population_energies = (np.ones(self.num_population_members) *
+                                    np.inf)
+
+        # reset number of function evaluations counter
+        self._nfev = 0
+
+    def init_population_random(self):
+        """
+        Initialises the population at random.  This type of initialization
+        can possess clustering, Latin Hypercube sampling is generally better.
+        """
+        rng = self.random_number_generator
+        self.population = rng.random_sample(self.population_shape)
+
+        # reset population energies
+        self.population_energies = (np.ones(self.num_population_members) *
+                                    np.inf)
+
+        # reset number of function evaluations counter
+        self._nfev = 0
+
+    @property
+    def x(self):
+        """
+        The best solution from the solver
+
+        Returns
+        -------
+        x : ndarray
+            The best solution from the solver.
+        """
+        return self._scale_parameters(self.population[0])
+
+    @property
+    def convergence(self):
+        """
+        The standard deviation of the population energies divided by their
+        mean.
+        """
+        return (np.std(self.population_energies) /
+                np.abs(np.mean(self.population_energies) + _MACHEPS))
+
+    def solve(self):
+        """
+        Runs the DifferentialEvolutionSolver.
+
+        Returns
+        -------
+        res : OptimizeResult
+            The optimization result represented as a ``OptimizeResult`` object.
+            Important attributes are: ``x`` the solution array, ``success`` a
+            Boolean flag indicating if the optimizer exited successfully and
+            ``message`` which describes the cause of the termination. See
+            `OptimizeResult` for a description of other attributes.  If `polish`
+            was employed, and a lower minimum was obtained by the polishing,
+            then OptimizeResult also contains the ``jac`` attribute.
+        """
+        nit, warning_flag = 0, False
+        status_message = _status_message['success']
+
+        # The population may have just been initialized (all entries are
+        # np.inf). If it has you have to calculate the initial energies.
+        # Although this is also done in the evolve generator it's possible
+        # that someone can set maxiter=0, at which point we still want the
+        # initial energies to be calculated (the following loop isn't run).
+        if np.all(np.isinf(self.population_energies)):
+            self._calculate_population_energies()
+
+        # do the optimisation.
+        for nit in range(1, self.maxiter + 1):
+            # evolve the population by a generation
+            try:
+                next(self)
+            except StopIteration:
+                warning_flag = True
+                status_message = _status_message['maxfev']
+                break
+
+            if self.disp:
+                print("differential_evolution step %d: f(x)= %g"
+                      % (nit,
+                         self.population_energies[0]))
+
+            # stop when the fractional s.d. of the population is less than tol
+            # of the mean energy
+            convergence = self.convergence
+
+            if (self.callback and
+                    self.callback(self._scale_parameters(self.population[0]),
+                                  convergence=self.tol / convergence) is True):
+
+                warning_flag = True
+                status_message = ('callback function requested stop early '
+                                  'by returning True')
+                break
+
+            if convergence < self.tol or warning_flag:
+                break
+
+        else:
+            status_message = _status_message['maxiter']
+            warning_flag = True
+
+        DE_result = OptimizeResult(
+            x=self.x,
+            fun=self.population_energies[0],
+            nfev=self._nfev,
+            nit=nit,
+            message=status_message,
+            success=(warning_flag is not True))
+
+        if self.polish:
+            result = minimize(self.func,
+                              np.copy(DE_result.x),
+                              method='L-BFGS-B',
+                              bounds=self.limits.T,
+                              args=self.args)
+
+            self._nfev += result.nfev
+            DE_result.nfev = self._nfev
+
+            if result.fun < DE_result.fun:
+                DE_result.fun = result.fun
+                DE_result.x = result.x
+                DE_result.jac = result.jac
+                # to keep internal state consistent
+                self.population_energies[0] = result.fun
+                self.population[0] = self._unscale_parameters(result.x)
+
+        return DE_result
+
+    def _calculate_population_energies(self):
+        """
+        Calculate the energies of all the population members at the same time.
+        Puts the best member in first place. Useful if the population has just
+        been initialised.
+        """
+        for index, candidate in enumerate(self.population):
+            if self._nfev > self.maxfun:
+                break
+
+            parameters = self._scale_parameters(candidate)
+            self.population_energies[index] = self.func(parameters,
+                                                        *self.args)
+            self._nfev += 1
+
+        minval = np.argmin(self.population_energies)
+
+        # put the lowest energy into the best solution position.
+        lowest_energy = self.population_energies[minval]
+        self.population_energies[minval] = self.population_energies[0]
+        self.population_energies[0] = lowest_energy
+
+        self.population[[0, minval], :] = self.population[[minval, 0], :]
+
+    def __iter__(self):
+        return self
+
+    def __next__(self):
+        """
+        Evolve the population by a single generation
+
+        Returns
+        -------
+        x : ndarray
+            The best solution from the solver.
+        fun : float
+            Value of objective function obtained from the best solution.
+        """
+        # the population may have just been initialized (all entries are
+        # np.inf). If it has you have to calculate the initial energies
+        if np.all(np.isinf(self.population_energies)):
+            self._calculate_population_energies()
+
+        if self.dither is not None:
+            self.scale = (self.random_number_generator.rand()
+                          * (self.dither[1] - self.dither[0]) + self.dither[0])
+
+        for candidate in range(self.num_population_members):
+            if self._nfev > self.maxfun:
+                raise StopIteration
+
+            # create a trial solution
+            trial = self._mutate(candidate)
+
+            # ensuring that it's in the range [0, 1)
+            self._ensure_constraint(trial)
+
+            # scale from [0, 1) to the actual parameter value
+            parameters = self._scale_parameters(trial)
+
+            # determine the energy of the objective function
+            energy = self.func(parameters, *self.args)
+            self._nfev += 1
+
+            # if the energy of the trial candidate is lower than the
+            # original population member then replace it
+            if energy < self.population_energies[candidate]:
+                self.population[candidate] = trial
+                self.population_energies[candidate] = energy
+
+                # if the trial candidate also has a lower energy than the
+                # best solution then replace that as well
+                if energy < self.population_energies[0]:
+                    self.population_energies[0] = energy
+                    self.population[0] = trial
+
+        return self.x, self.population_energies[0]
+
+    def next(self):
+        """
+        Evolve the population by a single generation
+
+        Returns
+        -------
+        x : ndarray
+            The best solution from the solver.
+        fun : float
+            Value of objective function obtained from the best solution.
+        """
+        # next() is required for compatibility with Python2.7.
+        return self.__next__()
+
+    def _scale_parameters(self, trial):
+        """
+        scale from a number between 0 and 1 to parameters.
+        """
+        return self.__scale_arg1 + (trial - 0.5) * self.__scale_arg2
+
+    def _unscale_parameters(self, parameters):
+        """
+        scale from parameters to a number between 0 and 1.
+        """
+        return (parameters - self.__scale_arg1) / self.__scale_arg2 + 0.5
+
+    def _ensure_constraint(self, trial):
+        """
+        make sure the parameters lie between the limits
+        """
+        for index, param in enumerate(trial):
+            if param > 1 or param < 0:
+                trial[index] = self.random_number_generator.rand()
+
+    def _mutate(self, candidate):
+        """
+        create a trial vector based on a mutation strategy
+        """
+        trial = np.copy(self.population[candidate])
+
+        rng = self.random_number_generator
+
+        fill_point = rng.randint(0, self.parameter_count)
+
+        if (self.strategy == 'randtobest1exp' or
+                self.strategy == 'randtobest1bin'):
+            bprime = self.mutation_func(candidate,
+                                        self._select_samples(candidate, 5))
+        else:
+            bprime = self.mutation_func(self._select_samples(candidate, 5))
+
+        if self.strategy in self._binomial:
+            crossovers = rng.rand(self.parameter_count)
+            crossovers = crossovers < self.cross_over_probability
+            # the last one is always from the bprime vector for binomial
+            # If you fill in modulo with a loop you have to set the last one to
+            # true. If you don't use a loop then you can have any random entry
+            # be True.
+            crossovers[fill_point] = True
+            trial = np.where(crossovers, bprime, trial)
+            return trial
+
+        elif self.strategy in self._exponential:
+            i = 0
+            while (i < self.parameter_count and
+                   rng.rand() < self.cross_over_probability):
+
+                trial[fill_point] = bprime[fill_point]
+                fill_point = (fill_point + 1) % self.parameter_count
+                i += 1
+
+            return trial
+
+    def _best1(self, samples):
+        """
+        best1bin, best1exp
+        """
+        r0, r1 = samples[:2]
+        return (self.population[0] + self.scale *
+                (self.population[r0] - self.population[r1]))
+
+    def _rand1(self, samples):
+        """
+        rand1bin, rand1exp
+        """
+        r0, r1, r2 = samples[:3]
+        return (self.population[r0] + self.scale *
+                (self.population[r1] - self.population[r2]))
+
+    def _best2(self, samples):
+        """
+        best2bin, best2exp
+        """
+        r0, r1, r2, r3 = samples[:4]
+        bprime = (self.population[0] + self.scale *
+                  (self.population[r0] + self.population[r1] -
+                   self.population[r2] - self.population[r3]))
+
+        return bprime
+
+    def _rand2(self, samples):
+        """
+        rand2bin, rand2exp
+        """
+        r0, r1, r2, r3, r4 = samples
+        bprime = (self.population[r0] + self.scale *
+                  (self.population[r1] + self.population[r2] -
+                   self.population[r3] - self.population[r4]))
+
+        return bprime
+
+    def _select_samples(self, candidate, number_samples):
+        """
+        obtain random integers from range(self.num_population_members),
+        without replacement.  You can't have the original candidate either.
+        """
+        idxs = list(range(self.num_population_members))
+        idxs.remove(candidate)
+        self.random_number_generator.shuffle(idxs)
+        idxs = idxs[:number_samples]
+        return idxs
+
+
+def _make_random_gen(seed):
+    """Turn seed into a np.random.RandomState instance
+
+    If seed is None, return the RandomState singleton used by np.random.
+    If seed is an int, return a new RandomState instance seeded with seed.
+    If seed is already a RandomState instance, return it.
+    Otherwise raise ValueError.
+    """
+    if seed is None or seed is np.random:
+        return np.random.mtrand._rand
+    if isinstance(seed, (numbers.Integral, np.integer)):
+        return np.random.RandomState(seed)
+    if isinstance(seed, np.random.RandomState):
+        return seed
+    raise ValueError('%r cannot be used to seed a numpy.random.RandomState'
+                     ' instance' % seed)
+
+
 

From 5e089a2c449964405817c30ee53233a44a14b889 Mon Sep 17 00:00:00 2001
From: Michael Hush <MichaelRHush@gmail.com>
Date: Thu, 6 Oct 2016 13:59:45 +1100
Subject: [PATCH 2/5] Further refining the differential evolution

Continuing the implementation.
---
 mloop/learners.py | 169 ++++++++++++++----------------------------------------
 1 file changed, 42 insertions(+), 127 deletions(-)

diff --git a/mloop/learners.py b/mloop/learners.py
index 7cfdaf9..028c30b 100644
--- a/mloop/learners.py
+++ b/mloop/learners.py
@@ -259,6 +259,7 @@ class RandomLearner(Learner, threading.Thread):
     Keyword Args:
         min_boundary (Optional [array]): If set to None, overrides default learner values and sets it to a set of value 0. Default None.
         max_boundary (Optional [array]): If set to None overides default learner values and sets it to an array of value 1. Default None.
+        first_params (Optional [array]): The first parameters to test. If None will just randomly sample the initial condition.
         trust_region (Optional [float or array]): The trust region defines the maximum distance the learner will travel from the current best set of parameters. If None, the learner will search everywhere. If a float, this number must be between 0 and 1 and defines maximum distance the learner will venture as a percentage of the boundaries. If it is an array, it must have the same size as the number of parameters and the numbers define the maximum absolute distance that can be moved along each direction.   
     '''
     
@@ -1181,8 +1182,7 @@ def find_local_minima(self):
 
     
 class DifferentialEvolutionLearner(Learner, threading.Thread):
-
-    """
+    '''
     Adaption of the differential evolution algorithm in scipy. 
     
     Args:
@@ -1191,30 +1191,44 @@ class DifferentialEvolutionLearner(Learner, threading.Thread):
         end_event (event): Event to trigger end of learner.
         
     Keyword Args:
+        first_params (Optional [array]): The first parameters to test. If None will just randomly sample the initial condition. Default None.
+        trust_region (Optional [float or array]): The trust region defines the maximum distance the learner will travel from the current best set of parameters. If None, the learner will search everywhere. If a float, this number must be between 0 and 1 and defines maximum distance the learner will venture as a percentage of the boundaries. If it is an array, it must have the same size as the number of parameters and the numbers define the maximum absolute distance that can be moved along each direction. 
         evolution_strategy (Optional [string]): the differential evolution strategy to use, options are 'best1', 'best1', 'rand1' and 'rand2'. The default is 'rand2'.
         population_size (Optional [int]): multiplier proportional to the number of parameters in a generation. The generation population is set to population_size * parameter_num. Default 15.
         mutation_scale (Optional [tuple]): The mutation scale when picking new points. Otherwise known as differential weight. When provided as a tuple (min,max) a mutation constant is picked randomly in the interval. Default (0.5,1.0).
-        crossover_probability (Optional [float]): The recombination constand or crossover probability, the probability a new points will be added to the population.
+        cross_over_probability (Optional [float]): The recombination constand or crossover probability, the probability a new points will be added to the population.
         restart_tolerance (Optional [float]): when the current population have a spread less than the initial tolerance, namely stdev(curr_pop) < restart_tolerance stdev(init_pop), it is likely the population is now in a minima, and so the search is started again.
-        trust_region (Optional [float or array]): The trust region defines the maximum distance the learner will travel from the current best set of parameters. If None, the learner will search everywhere. If a float, this number must be between 0 and 1 and defines maximum distance the learner will venture as a percentage of the boundaries. If it is an array, it must have the same size as the number of parameters and the numbers define the maximum absolute distance that can be moved along each direction. 
         
     Attributes:
         has_trust_region (bool): Whether the learner has a trust region. 
         
-    """
+    '''
 
     # Dispatch of mutation strategy method (binomial or exponential).
     
     def __init__(self, 
-                 strategy='rand2', 
+                 first_params = None,
+                 trust_region = None,
+                 evolution_strategy='rand2', 
                  population_size=15,
+                 mutation_scale=(0.5, 1), 
+                 cross_over_probability=0.7, 
                  restart_tolerance=0.01, 
-                 mutation=(0.5, 1), 
-                 recombination=0.7, 
                  **kwargs):
         
         super(NelderMeadLearner,self).__init__(**kwargs)
         
+        if first_params is None:
+            self.first_params = None
+        else:
+            self.first_params = np.array(first_params, dtype=float)
+            if not self.check_num_params(self.first_params):
+                self.log.error('first_params has the wrong number of parameters:' + repr(self.first_params))
+                raise ValueError
+            if not self.check_in_boundary(self.first_params):
+                self.log.error('first_params is not in the boundary:' + repr(self.first_params))
+                raise ValueError
+        
         if strategy == 'best1':
             self.mutation_func = self._best1
         elif strategy == 'best2':
@@ -1224,7 +1238,7 @@ def __init__(self,
         elif strategy == 'rand2':
             self.mutation_func = self._rand2
         else:
-            self.log.error("Please select a valid mutation strategy")
+            self.log.error('Please select a valid mutation strategy')
             raise ValueError
         
         self.strategy = strategy
@@ -1232,41 +1246,18 @@ def __init__(self,
 
         # Mutation constant should be in [0, 2). If specified as a sequence
         # then dithering is performed.
-        self.scale = mutation
-        if (not np.all(np.isfinite(mutation)) or
-                np.any(np.array(mutation) >= 2) or
-                np.any(np.array(mutation) < 0)):
-            raise ValueError('The mutation constant must be a float in '
-                             'U[0, 2), or specified as a tuple(min, max)'
-                             ' where min < max and min, max are in U[0, 2).')
-
-        self.dither = None
-        if hasattr(mutation, '__iter__') and len(mutation) > 1:
-            self.dither = [mutation[0], mutation[1]]
-            self.dither.sort()
-
-        self.cross_over_probability = recombination
-
-        self.func = func
-        self.args = args
-
-        # convert tuple of lower and upper bounds to limits
-        # [(low_0, high_0), ..., (low_n, high_n]
-        #     -> [[low_0, ..., low_n], [high_0, ..., high_n]]
-        self.limits = np.array(bounds, dtype='float').T
-        if (np.size(self.limits, 0) != 2 or not
-                np.all(np.isfinite(self.limits))):
-            raise ValueError('bounds should be a sequence containing '
-                             'real valued (min, max) pairs for each value'
-                             ' in x')
-
-        if maxiter is None:  # the default used to be None
-            maxiter = 1000
-        self.maxiter = maxiter
-        if maxfun is None:  # the default used to be None
-            maxfun = np.inf
-        self.maxfun = maxfun
-
+        if len(mutation_scale) == 2 and (np.any(np.array(mutation_scale) <= 2) or np.any(np.array(mutation_scale) > 0)):
+            self.mutation_scale = mutation_scale
+        else:
+            self.log.error('Mutation scale must be a tuple with (min,max) between 0 and 2. mutation_scale:' + repr(mutation_scale))
+            raise ValueError
+        
+        if cross_over_probability <= 1 and cross_over_probability >= 0:
+            self.cross_over_probability = cross_over_probability
+        else:
+            self.log.error('Cross over probability must be between 0 and 1. cross_over_probability:' + repr(cross_over_probability))
+        
+        
         # population is scaled to between [0, 1].
         # We have to scale between parameter <-> population
         # save these arguments for _scale_parameter and
@@ -1276,8 +1267,6 @@ def __init__(self,
 
         self.parameter_count = np.size(self.limits, 1)
 
-        self.random_number_generator = _make_random_gen(seed)
-
         # default population initialization is a latin hypercube design, but
         # there are other population initializations possible.
         self.num_population_members = popsize * self.parameter_count
@@ -1286,59 +1275,15 @@ def __init__(self,
                                  self.parameter_count)
 
         self._nfev = 0
-        if init == 'latinhypercube':
-            self.init_population_lhs()
-        elif init == 'random':
-            self.init_population_random()
-        else:
-            raise ValueError("The population initialization method must be one"
-                             "of 'latinhypercube' or 'random'")
+        
+        self.sample_random_population()
 
         self.disp = disp
 
-    def init_population_lhs(self):
-        """
-        Initializes the population with Latin Hypercube Sampling.
-        Latin Hypercube Sampling ensures that each parameter is uniformly
-        sampled over its range.
-        """
-        rng = self.random_number_generator
-
-        # Each parameter range needs to be sampled uniformly. The scaled
-        # parameter range ([0, 1)) needs to be split into
-        # `self.num_population_members` segments, each of which has the following
-        # size:
-        segsize = 1.0 / self.num_population_members
-
-        # Within each segment we sample from a uniform random distribution.
-        # We need to do this sampling for each parameter.
-        samples = (segsize * rng.random_sample(self.population_shape)
-
-        # Offset each segment to cover the entire parameter range [0, 1)
-                   + np.linspace(0., 1., self.num_population_members,
-                                 endpoint=False)[:, np.newaxis])
-
-        # Create an array for population of candidate solutions.
-        self.population = np.zeros_like(samples)
-
-        # Initialize population of candidate solutions by permutation of the
-        # random samples.
-        for j in range(self.parameter_count):
-            order = rng.permutation(range(self.num_population_members))
-            self.population[:, j] = samples[order, j]
-
-        # reset population energies
-        self.population_energies = (np.ones(self.num_population_members) *
-                                    np.inf)
-
-        # reset number of function evaluations counter
-        self._nfev = 0
-
-    def init_population_random(self):
-        """
-        Initialises the population at random.  This type of initialization
-        can possess clustering, Latin Hypercube sampling is generally better.
-        """
+    def sample_random_population(self):
+        '''
+        Sample a new random set of variables
+        '''
         rng = self.random_number_generator
         self.population = rng.random_sample(self.population_shape)
 
@@ -1374,16 +1319,6 @@ def solve(self):
         """
         Runs the DifferentialEvolutionSolver.
 
-        Returns
-        -------
-        res : OptimizeResult
-            The optimization result represented as a ``OptimizeResult`` object.
-            Important attributes are: ``x`` the solution array, ``success`` a
-            Boolean flag indicating if the optimizer exited successfully and
-            ``message`` which describes the cause of the termination. See
-            `OptimizeResult` for a description of other attributes.  If `polish`
-            was employed, and a lower minimum was obtained by the polishing,
-            then OptimizeResult also contains the ``jac`` attribute.
         """
         nit, warning_flag = 0, False
         status_message = _status_message['success']
@@ -1395,7 +1330,7 @@ def solve(self):
         # initial energies to be calculated (the following loop isn't run).
         if np.all(np.isinf(self.population_energies)):
             self._calculate_population_energies()
-
+        
         # do the optimisation.
         for nit in range(1, self.maxiter + 1):
             # evolve the population by a generation
@@ -1659,23 +1594,3 @@ def _select_samples(self, candidate, number_samples):
         idxs = idxs[:number_samples]
         return idxs
 
-
-def _make_random_gen(seed):
-    """Turn seed into a np.random.RandomState instance
-
-    If seed is None, return the RandomState singleton used by np.random.
-    If seed is an int, return a new RandomState instance seeded with seed.
-    If seed is already a RandomState instance, return it.
-    Otherwise raise ValueError.
-    """
-    if seed is None or seed is np.random:
-        return np.random.mtrand._rand
-    if isinstance(seed, (numbers.Integral, np.integer)):
-        return np.random.RandomState(seed)
-    if isinstance(seed, np.random.RandomState):
-        return seed
-    raise ValueError('%r cannot be used to seed a numpy.random.RandomState'
-                     ' instance' % seed)
-
-
-

From 7788a8bd08f716deb5e1eae2a9ad4c0a1c8531df Mon Sep 17 00:00:00 2001
From: Michael Hush <MichaelRHush@gmail.com>
Date: Mon, 10 Oct 2016 18:03:29 +1100
Subject: [PATCH 3/5] First version of DEComplete

First version of differential evolution complete. Next testing will be
investigated
---
 mloop/learners.py | 495 +++++++++++++++++++-----------------------------------
 1 file changed, 174 insertions(+), 321 deletions(-)

diff --git a/mloop/learners.py b/mloop/learners.py
index 028c30b..d4eeafd 100644
--- a/mloop/learners.py
+++ b/mloop/learners.py
@@ -8,6 +8,7 @@
 
 import threading
 import numpy as np
+import random
 import numpy.random as nr
 import scipy.optimize as so
 import logging
@@ -318,7 +319,6 @@ def run(self):
         self._shut_down()
         self.log.debug('Ended Random Learner')
     
-    
 class NelderMeadLearner(Learner, threading.Thread):
     '''
     Nelder-Mead learner. Executes the Nelder-Mead learner algorithm and stores the needed simplex to estimate the next points.
@@ -1193,7 +1193,7 @@ class DifferentialEvolutionLearner(Learner, threading.Thread):
     Keyword Args:
         first_params (Optional [array]): The first parameters to test. If None will just randomly sample the initial condition. Default None.
         trust_region (Optional [float or array]): The trust region defines the maximum distance the learner will travel from the current best set of parameters. If None, the learner will search everywhere. If a float, this number must be between 0 and 1 and defines maximum distance the learner will venture as a percentage of the boundaries. If it is an array, it must have the same size as the number of parameters and the numbers define the maximum absolute distance that can be moved along each direction. 
-        evolution_strategy (Optional [string]): the differential evolution strategy to use, options are 'best1', 'best1', 'rand1' and 'rand2'. The default is 'rand2'.
+        evolution_strategy (Optional [string]): the differential evolution strategy to use, options are 'best1', 'best1', 'rand1' and 'rand2'. The default is 'best2'.
         population_size (Optional [int]): multiplier proportional to the number of parameters in a generation. The generation population is set to population_size * parameter_num. Default 15.
         mutation_scale (Optional [tuple]): The mutation scale when picking new points. Otherwise known as differential weight. When provided as a tuple (min,max) a mutation constant is picked randomly in the interval. Default (0.5,1.0).
         cross_over_probability (Optional [float]): The recombination constand or crossover probability, the probability a new points will be added to the population.
@@ -1229,6 +1229,8 @@ def __init__(self,
                 self.log.error('first_params is not in the boundary:' + repr(self.first_params))
                 raise ValueError
         
+        self._set_trust_region(trust_region)
+        
         if strategy == 'best1':
             self.mutation_func = self._best1
         elif strategy == 'best2':
@@ -1244,8 +1246,6 @@ def __init__(self,
         self.strategy = strategy
         self.restart_tolerance = restart_tolerance
 
-        # Mutation constant should be in [0, 2). If specified as a sequence
-        # then dithering is performed.
         if len(mutation_scale) == 2 and (np.any(np.array(mutation_scale) <= 2) or np.any(np.array(mutation_scale) > 0)):
             self.mutation_scale = mutation_scale
         else:
@@ -1257,340 +1257,193 @@ def __init__(self,
         else:
             self.log.error('Cross over probability must be between 0 and 1. cross_over_probability:' + repr(cross_over_probability))
         
+        if population_size >= 5:
+            self.population_size = population_size
+        else:
+            self.log.error('Population size must be greater or equal to 5:' + repr(population_size))
         
-        # population is scaled to between [0, 1].
-        # We have to scale between parameter <-> population
-        # save these arguments for _scale_parameter and
-        # _unscale_parameter. This is an optimization
-        self.__scale_arg1 = 0.5 * (self.limits[0] + self.limits[1])
-        self.__scale_arg2 = np.fabs(self.limits[0] - self.limits[1])
-
-        self.parameter_count = np.size(self.limits, 1)
-
-        # default population initialization is a latin hypercube design, but
-        # there are other population initializations possible.
-        self.num_population_members = popsize * self.parameter_count
-
-        self.population_shape = (self.num_population_members,
-                                 self.parameter_count)
-
-        self._nfev = 0
+        self.num_population_members = self.population_size * self.num_params
         
-        self.sample_random_population()
-
-        self.disp = disp
-
-    def sample_random_population(self):
+        self.first_sample = True
+        
+        self.params_generations = []
+        self.costs_generations = []
+        
+        self.min_index = 0
+        self.init_std = 0
+        self.curr_std = 0
+    
+    def run(self):
+        '''
+        Runs the Differential Evolution Learner.
+        '''
+        try:
+        
+            generate_population(self)
+        
+            while not self.end_event.is_set():
+                
+                self.next_generation()
+                
+                if self.curr_std < self.restart_tolerance * self.init_std:
+                    self.generate_population()
+            
+        except LearnerInterupt:
+            return
+        
+    def generate_population(self):
         '''
         Sample a new random set of variables
         '''
-        rng = self.random_number_generator
-        self.population = rng.random_sample(self.population_shape)
-
-        # reset population energies
-        self.population_energies = (np.ones(self.num_population_members) *
-                                    np.inf)
-
-        # reset number of function evaluations counter
-        self._nfev = 0
-
-    @property
-    def x(self):
-        """
-        The best solution from the solver
-
-        Returns
-        -------
-        x : ndarray
-            The best solution from the solver.
-        """
-        return self._scale_parameters(self.population[0])
-
-    @property
-    def convergence(self):
-        """
-        The standard deviation of the population energies divided by their
-        mean.
-        """
-        return (np.std(self.population_energies) /
-                np.abs(np.mean(self.population_energies) + _MACHEPS))
-
-    def solve(self):
-        """
-        Runs the DifferentialEvolutionSolver.
-
-        """
-        nit, warning_flag = 0, False
-        status_message = _status_message['success']
-
-        # The population may have just been initialized (all entries are
-        # np.inf). If it has you have to calculate the initial energies.
-        # Although this is also done in the evolve generator it's possible
-        # that someone can set maxiter=0, at which point we still want the
-        # initial energies to be calculated (the following loop isn't run).
-        if np.all(np.isinf(self.population_energies)):
-            self._calculate_population_energies()
-        
-        # do the optimisation.
-        for nit in range(1, self.maxiter + 1):
-            # evolve the population by a generation
-            try:
-                next(self)
-            except StopIteration:
-                warning_flag = True
-                status_message = _status_message['maxfev']
-                break
-
-            if self.disp:
-                print("differential_evolution step %d: f(x)= %g"
-                      % (nit,
-                         self.population_energies[0]))
-
-            # stop when the fractional s.d. of the population is less than tol
-            # of the mean energy
-            convergence = self.convergence
-
-            if (self.callback and
-                    self.callback(self._scale_parameters(self.population[0]),
-                                  convergence=self.tol / convergence) is True):
-
-                warning_flag = True
-                status_message = ('callback function requested stop early '
-                                  'by returning True')
-                break
-
-            if convergence < self.tol or warning_flag:
-                break
-
+        
+        self.population = []
+        self.population_costs = []
+        self.min_index = 0
+        
+        if self.first_params is not None and self.first_sample:
+            curr_params = self.first_params
         else:
-            status_message = _status_message['maxiter']
-            warning_flag = True
-
-        DE_result = OptimizeResult(
-            x=self.x,
-            fun=self.population_energies[0],
-            nfev=self._nfev,
-            nit=nit,
-            message=status_message,
-            success=(warning_flag is not True))
-
-        if self.polish:
-            result = minimize(self.func,
-                              np.copy(DE_result.x),
-                              method='L-BFGS-B',
-                              bounds=self.limits.T,
-                              args=self.args)
-
-            self._nfev += result.nfev
-            DE_result.nfev = self._nfev
-
-            if result.fun < DE_result.fun:
-                DE_result.fun = result.fun
-                DE_result.x = result.x
-                DE_result.jac = result.jac
-                # to keep internal state consistent
-                self.population_energies[0] = result.fun
-                self.population[0] = self._unscale_parameters(result.x)
-
-        return DE_result
-
-    def _calculate_population_energies(self):
-        """
-        Calculate the energies of all the population members at the same time.
-        Puts the best member in first place. Useful if the population has just
-        been initialised.
-        """
-        for index, candidate in enumerate(self.population):
-            if self._nfev > self.maxfun:
-                break
-
-            parameters = self._scale_parameters(candidate)
-            self.population_energies[index] = self.func(parameters,
-                                                        *self.args)
-            self._nfev += 1
-
-        minval = np.argmin(self.population_energies)
-
-        # put the lowest energy into the best solution position.
-        lowest_energy = self.population_energies[minval]
-        self.population_energies[minval] = self.population_energies[0]
-        self.population_energies[0] = lowest_energy
-
-        self.population[[0, minval], :] = self.population[[minval, 0], :]
-
-    def __iter__(self):
-        return self
-
-    def __next__(self):
-        """
-        Evolve the population by a single generation
-
-        Returns
-        -------
-        x : ndarray
-            The best solution from the solver.
-        fun : float
-            Value of objective function obtained from the best solution.
-        """
-        # the population may have just been initialized (all entries are
-        # np.inf). If it has you have to calculate the initial energies
-        if np.all(np.isinf(self.population_energies)):
-            self._calculate_population_energies()
-
-        if self.dither is not None:
-            self.scale = (self.random_number_generator.rand()
-                          * (self.dither[1] - self.dither[0]) + self.dither[0])
-
-        for candidate in range(self.num_population_members):
-            if self._nfev > self.maxfun:
-                raise StopIteration
-
-            # create a trial solution
-            trial = self._mutate(candidate)
-
-            # ensuring that it's in the range [0, 1)
-            self._ensure_constraint(trial)
-
-            # scale from [0, 1) to the actual parameter value
-            parameters = self._scale_parameters(trial)
-
-            # determine the energy of the objective function
-            energy = self.func(parameters, *self.args)
-            self._nfev += 1
-
-            # if the energy of the trial candidate is lower than the
-            # original population member then replace it
-            if energy < self.population_energies[candidate]:
-                self.population[candidate] = trial
-                self.population_energies[candidate] = energy
-
-                # if the trial candidate also has a lower energy than the
-                # best solution then replace that as well
-                if energy < self.population_energies[0]:
-                    self.population_energies[0] = energy
-                    self.population[0] = trial
-
-        return self.x, self.population_energies[0]
-
-    def next(self):
-        """
+            curr_params = self.min_boundary + nr.rand(self.num_params) * self.diff_boundary
+        
+        try:
+            curr_cost = self.put_params_and_get_cost(curr_params)
+        except LearnerInterrupt:
+            self.log.info('DELearner ended during first sample of population.')
+            raise
+        
+        self.population.append(curr_params)
+        self.population_costs.append(curr_cost)
+        
+        for index in range(1, self.num_population_members):
+            
+            if self.has_trust_region:
+                temp_min = np.maximum(self.min_boundary,self.population[self.min_index] - self.trust_region)
+                temp_max = np.minimum(self.max_boundary,self.population[self.min_index] + self.trust_region)
+                curr_params = temp_min + nr.rand(self.num_params) * (temp_max - temp_min)
+            else:
+                curr_params = self.min_boundary + nr.rand(self.num_params) * self.diff_boundary
+            
+            try:
+                curr_cost = self.put_params_and_get_cost(curr_params)
+            except LearnerInterrupt:
+                self.log.info('DELearner ended during initial sample of population.')
+                raise
+            
+            self.population.append(curr_params)
+            self.population_costs.append(curr_cost)
+            
+            if curr_cost < self.population_costs[self.min_index]:
+                self.min_index = index
+        
+        self.population = np.array(self.population)
+        self.population_costs = np.array(self.population_costs)
+        
+        self.init_std = self.std(self.population_costs)
+        self.curr_std = self.init_std
+        
+        self.params_generations.append(np.copy(self.population))
+        self.costs_generations.append(np.copy(self.population_costs))
+    
+    def next_generation(self):
+        '''
         Evolve the population by a single generation
+        '''
+        
+        self.curr_scale = nr.uniform(self.mutation_scale[0], self.mutation_scale[1])
+        
+        for index in range(self.num_population_members):
+            
+            curr_params = self.mutate(index)
 
-        Returns
-        -------
-        x : ndarray
-            The best solution from the solver.
-        fun : float
-            Value of objective function obtained from the best solution.
-        """
-        # next() is required for compatibility with Python2.7.
-        return self.__next__()
-
-    def _scale_parameters(self, trial):
-        """
-        scale from a number between 0 and 1 to parameters.
-        """
-        return self.__scale_arg1 + (trial - 0.5) * self.__scale_arg2
-
-    def _unscale_parameters(self, parameters):
-        """
-        scale from parameters to a number between 0 and 1.
-        """
-        return (parameters - self.__scale_arg1) / self.__scale_arg2 + 0.5
-
-    def _ensure_constraint(self, trial):
-        """
-        make sure the parameters lie between the limits
-        """
-        for index, param in enumerate(trial):
-            if param > 1 or param < 0:
-                trial[index] = self.random_number_generator.rand()
-
-    def _mutate(self, candidate):
-        """
-        create a trial vector based on a mutation strategy
-        """
-        trial = np.copy(self.population[candidate])
-
-        rng = self.random_number_generator
-
-        fill_point = rng.randint(0, self.parameter_count)
+            try:
+                curr_cost = self.put_params_and_get_cost(curr_params)
+            except LearnerInterrupt:
+                self.log.info('DELearner ended during initial sample of population.')
+                raise
+            
+            if curr_cost < self.population_costs[index]:
+                self.population[index] = curr_params
+                self.population_costs[index] = curr_cost
+                
+                if curr_cost < self.population_costs[self.min_index]:
+                    self.min_index = index
+        
+        self.curr_std = self.std(self.population_costs)
+        
+        self.params_generations.append(np.copy(self.population))
+        self.costs_generations.append(np.copy(self.population_costs))
 
-        if (self.strategy == 'randtobest1exp' or
-                self.strategy == 'randtobest1bin'):
-            bprime = self.mutation_func(candidate,
-                                        self._select_samples(candidate, 5))
+    def mutate(self, index):
+        '''
+        Mutate the parameters at index.
+        
+        Args:
+            index (int): Index of the point to be mutated.
+        '''
+        
+        fill_point = nr.randint(0, self.parameter_count)
+        candidate_params = self.mutation_func(index)
+        crossovers = nr.rand(self.parameter_count) < self.cross_over_probability
+        crossovers[fill_point] = True
+        mutated_params = np.where(crossovers, candidate_params, self.population[index])
+        print(mutated_params)
+        
+        if self.has_trust_region:
+            temp_min = np.maximum(self.min_boundary,self.population[self.min_index] - self.trust_region)
+            temp_max = np.minimum(self.max_boundary,self.population[self.min_index] + self.trust_region)
+            rand_params = temp_min + nr.rand(self.num_params) * (temp_max - temp_min)
         else:
-            bprime = self.mutation_func(self._select_samples(candidate, 5))
-
-        if self.strategy in self._binomial:
-            crossovers = rng.rand(self.parameter_count)
-            crossovers = crossovers < self.cross_over_probability
-            # the last one is always from the bprime vector for binomial
-            # If you fill in modulo with a loop you have to set the last one to
-            # true. If you don't use a loop then you can have any random entry
-            # be True.
-            crossovers[fill_point] = True
-            trial = np.where(crossovers, bprime, trial)
-            return trial
-
-        elif self.strategy in self._exponential:
-            i = 0
-            while (i < self.parameter_count and
-                   rng.rand() < self.cross_over_probability):
+            rand_params = self.min_boundary + nr.rand(self.num_params) * self.diff_boundary
+        print(rand_params)
+        
+        projected_params = np.where(np.logical_or(mutated_params < self.min_boundary, mutated_params > self.max_boundary),rand_params,mutated_params)
+        print(projected_params)
+        
+        return projected_params
 
-                trial[fill_point] = bprime[fill_point]
-                fill_point = (fill_point + 1) % self.parameter_count
-                i += 1
+    def _best1(self, index):
+        '''
+        Use best parameters and two others to generate mutation.
+        
+        Args:
+            index (int): Index of member to mutate. 
+        '''
+        r0, r1 = random.sample(range(self.num_population_members).remove(index),2)
+        samples[:2]
+        return (self.population[self.min_index] + self.curr_scale *(self.population[r0] - self.population[r1]))
 
-            return trial
+    def _rand1(self, index):
+        '''
+        Use three random parameters to generate mutation.
+        
+        Args:
+            index (int): Index of member to mutate. 
+        '''
+        r0, r1, r2 = random.sample(range(self.num_population_members).remove(index),3)
+        return (self.population[r0] + self.curr_scale * (self.population[r1] - self.population[r2]))
 
-    def _best1(self, samples):
-        """
-        best1bin, best1exp
-        """
-        r0, r1 = samples[:2]
-        return (self.population[0] + self.scale *
-                (self.population[r0] - self.population[r1]))
+    def _best2(self, index):
+        '''
+        Use best parameters and four others to generate mutation.
+        
+        Args:
+            index (int): Index of member to mutate. 
+        '''
+        r0, r1, r2, r3 = random.sample(range(self.num_population_members).remove(index),4)
+        return self.population[self.min_index] + self.curr_scale * (self.population[r0] + self.population[r1] - self.population[r2] - self.population[r3])
 
-    def _rand1(self, samples):
-        """
-        rand1bin, rand1exp
-        """
-        r0, r1, r2 = samples[:3]
-        return (self.population[r0] + self.scale *
-                (self.population[r1] - self.population[r2]))
+    def _rand2(self, index):
+        '''
+        Use five random parameters to generate mutation.
+        
+        Args:
+            index (int): Index of member to mutate. 
+        '''
+        r0, r1, r2, r3, r4 = random.sample(range(self.num_population_members).remove(index),5)
+        return self.population[r0] + self.curr_scale * (self.population[r1] + self.population[r2] - self.population[r3] - self.population[r4])
 
-    def _best2(self, samples):
-        """
-        best2bin, best2exp
-        """
-        r0, r1, r2, r3 = samples[:4]
-        bprime = (self.population[0] + self.scale *
-                  (self.population[r0] + self.population[r1] -
-                   self.population[r2] - self.population[r3]))
 
-        return bprime
 
-    def _rand2(self, samples):
-        """
-        rand2bin, rand2exp
-        """
-        r0, r1, r2, r3, r4 = samples
-        bprime = (self.population[r0] + self.scale *
-                  (self.population[r1] + self.population[r2] -
-                   self.population[r3] - self.population[r4]))
 
-        return bprime
 
-    def _select_samples(self, candidate, number_samples):
-        """
-        obtain random integers from range(self.num_population_members),
-        without replacement.  You can't have the original candidate either.
-        """
-        idxs = list(range(self.num_population_members))
-        idxs.remove(candidate)
-        self.random_number_generator.shuffle(idxs)
-        idxs = idxs[:number_samples]
-        return idxs
 

From 5084b9063ca2b736524d6f477f913b214dd21cde Mon Sep 17 00:00:00 2001
From: Michael Hush <MichaelRHush@gmail.com>
Date: Tue, 11 Oct 2016 18:05:06 +1100
Subject: [PATCH 4/5] DE learner complete and added tests.

The differential evolution controller has been implemented and is now
set at the default trainer for gaussian process. Tests have been added
to the automated suite and there is some basic visualizations available.
Still have to debug a possible issue with GP andexamples logging and
extras.
---
 .../differential_evolution_complete_config.txt     |  19 +
 examples/differential_evolution_simple_config.txt  |  15 +
 mloop/controllers.py                               |  57 ++-
 mloop/learners.py                                  | 558 +++++++++++----------
 mloop/visualizations.py                            | 118 ++++-
 tests/test_examples.py                             |  16 +
 6 files changed, 511 insertions(+), 272 deletions(-)
 create mode 100644 examples/differential_evolution_complete_config.txt
 create mode 100644 examples/differential_evolution_simple_config.txt

diff --git a/examples/differential_evolution_complete_config.txt b/examples/differential_evolution_complete_config.txt
new file mode 100644
index 0000000..88a8547
--- /dev/null
+++ b/examples/differential_evolution_complete_config.txt
@@ -0,0 +1,19 @@
+#Differential Evolution Complete Options
+#---------------------------------------
+
+#General options
+max_num_runs = 500              #number of planned runs
+target_cost = 0.1               #cost to beat
+
+#Differential evolution controller options
+controller_type = 'differential_evolution'
+num_params = 2                  #number of parameters
+min_boundary = [-1.2,-2]        #minimum boundary
+max_boundary = [10.0,4]         #maximum boundary
+trust_region = [3.2,3.1]        #maximum move distance from best params
+first_params = None          	#first parameters to try if None a random set of parameters is chosen
+evolution_strategy='best2'		#evolution strategy can be 'best1', 'best2', 'rand1' and 'rand2'. Best uses the best point, rand uses a random one, the number indicates the number of directions added.
+population_size=10				#a multiplier for the population size of a generation
+mutation_scale=(0.4, 1.1)		#the minimum and maximum value for the mutation scale factor. Each generation is randomly selected from this. Each value must be between 0 and 2.
+cross_over_probability=0.8     	#the probability a parameter will be resampled during a mutation in a new generation
+restart_tolerance=0.02			#the fraction the standard deviation in the costs of the population must reduce from the initial sample, before the search is restarted.
\ No newline at end of file
diff --git a/examples/differential_evolution_simple_config.txt b/examples/differential_evolution_simple_config.txt
new file mode 100644
index 0000000..d4615a0
--- /dev/null
+++ b/examples/differential_evolution_simple_config.txt
@@ -0,0 +1,15 @@
+#Differential Evolution Basic Options
+#------------------------------------
+
+#General options
+max_num_runs = 500              #number of planned runs
+target_cost = 0.1               #cost to beat
+
+#Differential evolution controller options
+controller_type = 'differential_evolution'
+num_params = 1                  #number of parameters
+min_boundary = [-4.8]           #minimum boundary
+max_boundary = [10.0]           #maximum boundary
+trust_region = 0.6              #maximum % move distance from best params
+first_params = [5.3]            #first parameters to try
+
diff --git a/mloop/controllers.py b/mloop/controllers.py
index bb1ebae..98d3a23 100644
--- a/mloop/controllers.py
+++ b/mloop/controllers.py
@@ -11,8 +11,8 @@
 import logging
 import os
 
-controller_dict = {'random':1,'nelder_mead':2,'gaussian_process':3}
-number_of_controllers = 3
+controller_dict = {'random':1,'nelder_mead':2,'gaussian_process':3,'differential_evolution':4}
+number_of_controllers = 4
 default_controller_archive_filename = 'controller_archive'
 default_controller_archive_file_type = 'txt'
 
@@ -47,6 +47,8 @@ def create_controller(interface,
     controller_type = str(controller_type)
     if controller_type=='gaussian_process':
         controller = GaussianProcessController(interface, **controller_config_dict)
+    elif controller_type=='differential_evolution':
+        controller = DifferentialEvolutionController(interface, **controller_config_dict)
     elif controller_type=='nelder_mead':
         controller = NelderMeadController(interface, **controller_config_dict)
     elif controller_type=='random':
@@ -489,6 +491,37 @@ def _next_params(self):
         self.learner_costs_queue.put(cost)
         return self.learner_params_queue.get()
 
+class DifferentialEvolutionController(Controller):
+    '''
+    Controller for the differential evolution learner. 
+    
+    Args:
+        params_out_queue (queue): Queue for parameters to next be run by experiment.
+        costs_in_queue (queue): Queue for costs (and other details) that have been returned by experiment.
+        **kwargs (Optional [dict]): Dictionary of options to be passed to Controller parent class and differential evolution learner.
+    '''
+    def __init__(self, interface,
+                **kwargs):
+        super(DifferentialEvolutionController,self).__init__(interface, **kwargs)    
+        
+        self.learner = mll.DifferentialEvolutionLearner(start_datetime = self.start_datetime,
+                                                        **self.remaining_kwargs)
+        
+        self._update_controller_with_learner_attributes()
+        self.out_type.append('differential_evolution')
+    
+    def _next_params(self):
+        '''
+        Gets next parameters from differential evolution learner.
+        '''
+        if self.curr_bad:
+            cost = float('inf')
+        else:
+            cost = self.curr_cost       
+        self.learner_costs_queue.put(cost)
+        return self.learner_params_queue.get()
+
+
 
 
 class GaussianProcessController(Controller):
@@ -506,7 +539,7 @@ class GaussianProcessController(Controller):
     '''
     
     def __init__(self, interface, 
-                 training_type='random',
+                 training_type='differential_evolution',
                  num_training_runs=None,
                  no_delay=True,
                  num_params=None,
@@ -553,6 +586,18 @@ def __init__(self, interface,
                                                  learner_archive_filename='training_learner_archive',
                                                  learner_archive_file_type=learner_archive_file_type,
                                                  **self.remaining_kwargs)
+        
+        elif self.training_type == 'differential_evolution':
+            self.learner = mll.DifferentialEvolutionLearner(start_datetime=self.start_datetime,
+                                                            num_params=num_params,
+                                                            min_boundary=min_boundary,
+                                                            max_boundary=max_boundary,
+                                                            trust_region=trust_region,
+                                                            evolution_strategy='rand2',
+                                                            learner_archive_filename='training_learner_archive',
+                                                            learner_archive_file_type=learner_archive_file_type,
+                                                            **self.remaining_kwargs)    
+        
         else:
             self.log.error('Unknown training type provided to Gaussian process controller:' + repr(training_type))
         
@@ -601,12 +646,12 @@ def _next_params(self):
         '''
         Gets next parameters from training learner.
         '''
-        if self.training_type == 'nelder_mead':
+        if self.training_type == 'differential_evolution' or self.training_type == 'nelder_mead':
             #Copied from NelderMeadController
-            if self.curr_bad:
+            if self.last_training_bad:
                 cost = float('inf')
             else:
-                cost = self.curr_cost       
+                cost = self.last_training_cost       
             self.learner_costs_queue.put(cost)
             temp = self.learner_params_queue.get()
             
diff --git a/mloop/learners.py b/mloop/learners.py
index d4eeafd..af44655 100644
--- a/mloop/learners.py
+++ b/mloop/learners.py
@@ -8,6 +8,7 @@
 
 import threading
 import numpy as np
+import math
 import random
 import numpy.random as nr
 import scipy.optimize as so
@@ -549,14 +550,305 @@ def run(self):
         self._shut_down()
         self.log.info('Ended Nelder-Mead')
 
-def update_archive(self):
+    def update_archive(self):
         '''
         Update the archive.
         '''
-        self.archive_dict.update({'archive_type':'nelder_mead_learner',
-                                  'simplex_parameters':self.simplex_params,
+        self.archive_dict.update({'simplex_parameters':self.simplex_params,
                                   'simplex_costs':self.simplex_costs})
 
+class DifferentialEvolutionLearner(Learner, threading.Thread):
+    '''
+    Adaption of the differential evolution algorithm in scipy. 
+    
+    Args:
+        params_out_queue (queue): Queue for parameters sent to controller.
+        costs_in_queue (queue): Queue for costs for gaussian process. This must be tuple
+        end_event (event): Event to trigger end of learner.
+        
+    Keyword Args:
+        first_params (Optional [array]): The first parameters to test. If None will just randomly sample the initial condition. Default None.
+        trust_region (Optional [float or array]): The trust region defines the maximum distance the learner will travel from the current best set of parameters. If None, the learner will search everywhere. If a float, this number must be between 0 and 1 and defines maximum distance the learner will venture as a percentage of the boundaries. If it is an array, it must have the same size as the number of parameters and the numbers define the maximum absolute distance that can be moved along each direction. 
+        evolution_strategy (Optional [string]): the differential evolution strategy to use, options are 'best1', 'best1', 'rand1' and 'rand2'. The default is 'best2'.
+        population_size (Optional [int]): multiplier proportional to the number of parameters in a generation. The generation population is set to population_size * parameter_num. Default 15.
+        mutation_scale (Optional [tuple]): The mutation scale when picking new points. Otherwise known as differential weight. When provided as a tuple (min,max) a mutation constant is picked randomly in the interval. Default (0.5,1.0).
+        cross_over_probability (Optional [float]): The recombination constand or crossover probability, the probability a new points will be added to the population.
+        restart_tolerance (Optional [float]): when the current population have a spread less than the initial tolerance, namely stdev(curr_pop) < restart_tolerance stdev(init_pop), it is likely the population is now in a minima, and so the search is started again.
+        
+    Attributes:
+        has_trust_region (bool): Whether the learner has a trust region.
+        num_population_members (int): The number of parameters in a generation. 
+        params_generations (list): History of the parameters generations. A list of all the parameters in the population, for each generation created.
+        costs_generations (list): History of the costs generations. A list of all the costs in the population, for each generation created.
+        init_std (float): The initial standard deviation in costs of the population. Calucalted after sampling (or resampling) the initial population.
+        curr_std (float): The current standard devation in costs of the population. Calculated after sampling each generation.
+    '''
+
+    def __init__(self, 
+                 first_params = None,
+                 trust_region = None,
+                 evolution_strategy='best1', 
+                 population_size=15,
+                 mutation_scale=(0.5, 1), 
+                 cross_over_probability=0.7, 
+                 restart_tolerance=0.01, 
+                 **kwargs):
+        
+        super(DifferentialEvolutionLearner,self).__init__(**kwargs)
+        
+        if first_params is None:
+            self.first_params = float('nan')
+        else:
+            self.first_params = np.array(first_params, dtype=float)
+            if not self.check_num_params(self.first_params):
+                self.log.error('first_params has the wrong number of parameters:' + repr(self.first_params))
+                raise ValueError
+            if not self.check_in_boundary(self.first_params):
+                self.log.error('first_params is not in the boundary:' + repr(self.first_params))
+                raise ValueError
+        
+        self._set_trust_region(trust_region)
+        
+        if evolution_strategy == 'best1':
+            self.mutation_func = self._best1
+        elif evolution_strategy == 'best2':
+            self.mutation_func = self._best2
+        elif evolution_strategy == 'rand1':
+            self.mutation_func = self._rand1
+        elif evolution_strategy == 'rand2':
+            self.mutation_func = self._rand2
+        else:
+            self.log.error('Please select a valid mutation strategy')
+            raise ValueError
+        
+        self.evolution_strategy = evolution_strategy
+        self.restart_tolerance = restart_tolerance
+
+        if len(mutation_scale) == 2 and (np.any(np.array(mutation_scale) <= 2) or np.any(np.array(mutation_scale) > 0)):
+            self.mutation_scale = mutation_scale
+        else:
+            self.log.error('Mutation scale must be a tuple with (min,max) between 0 and 2. mutation_scale:' + repr(mutation_scale))
+            raise ValueError
+        
+        if cross_over_probability <= 1 and cross_over_probability >= 0:
+            self.cross_over_probability = cross_over_probability
+        else:
+            self.log.error('Cross over probability must be between 0 and 1. cross_over_probability:' + repr(cross_over_probability))
+        
+        if population_size >= 5:
+            self.population_size = population_size
+        else:
+            self.log.error('Population size must be greater or equal to 5:' + repr(population_size))
+        
+        self.num_population_members = self.population_size * self.num_params
+        
+        self.first_sample = True
+        
+        self.params_generations = []
+        self.costs_generations = []
+        self.generation_count = 0
+        
+        self.min_index = 0
+        self.init_std = 0
+        self.curr_std = 0
+        
+        self.archive_dict.update({'archive_type':'differential_evolution',
+                                  'evolution_strategy':self.evolution_strategy,
+                                  'mutation_scale':self.mutation_scale,
+                                  'cross_over_probability':self.cross_over_probability,
+                                  'population_size':self.population_size,
+                                  'num_population_members':self.num_population_members,
+                                  'restart_tolerance':self.restart_tolerance,
+                                  'first_params':self.first_params,
+                                  'has_trust_region':self.has_trust_region,
+                                  'trust_region':self.trust_region})
+        
+        
+    def run(self):
+        '''
+        Runs the Differential Evolution Learner.
+        '''
+        try:
+        
+            self.generate_population()
+        
+            while not self.end_event.is_set():
+                
+                self.next_generation()
+                
+                if self.curr_std < self.restart_tolerance * self.init_std:
+                    self.generate_population()
+            
+        except LearnerInterrupt:
+            return
+    
+    def save_generation(self):
+        '''
+        Save history of generations.
+        '''
+        self.params_generations.append(np.copy(self.population))
+        self.costs_generations.append(np.copy(self.population_costs))
+        self.generation_count += 1
+    
+    def generate_population(self):
+        '''
+        Sample a new random set of variables
+        '''
+        
+        self.population = []
+        self.population_costs = []
+        self.min_index = 0
+        
+        if (not math.isnan(self.first_params)) and self.first_sample:
+            curr_params = self.first_params
+            self.first_sample = False
+        else:
+            curr_params = self.min_boundary + nr.rand(self.num_params) * self.diff_boundary
+        
+        curr_cost = self.put_params_and_get_cost(curr_params)
+        
+        self.population.append(curr_params)
+        self.population_costs.append(curr_cost)
+        
+        for index in range(1, self.num_population_members):
+            
+            if self.has_trust_region:
+                temp_min = np.maximum(self.min_boundary,self.population[self.min_index] - self.trust_region)
+                temp_max = np.minimum(self.max_boundary,self.population[self.min_index] + self.trust_region)
+                curr_params = temp_min + nr.rand(self.num_params) * (temp_max - temp_min)
+            else:
+                curr_params = self.min_boundary + nr.rand(self.num_params) * self.diff_boundary
+
+            curr_cost = self.put_params_and_get_cost(curr_params)
+            
+            self.population.append(curr_params)
+            self.population_costs.append(curr_cost)
+            
+            if curr_cost < self.population_costs[self.min_index]:
+                self.min_index = index
+        
+        self.population = np.array(self.population)
+        self.population_costs = np.array(self.population_costs)
+        
+        self.init_std = np.std(self.population_costs)
+        self.curr_std = self.init_std
+        
+        self.save_generation()
+        
+    def next_generation(self):
+        '''
+        Evolve the population by a single generation
+        '''
+        
+        self.curr_scale = nr.uniform(self.mutation_scale[0], self.mutation_scale[1])
+        
+        for index in range(self.num_population_members):
+            
+            curr_params = self.mutate(index)
+
+            curr_cost = self.put_params_and_get_cost(curr_params)
+            
+            if curr_cost < self.population_costs[index]:
+                self.population[index] = curr_params
+                self.population_costs[index] = curr_cost
+                
+                if curr_cost < self.population_costs[self.min_index]:
+                    self.min_index = index
+        
+        self.curr_std = np.std(self.population_costs)
+        
+        self.save_generation()
+
+    def mutate(self, index):
+        '''
+        Mutate the parameters at index.
+        
+        Args:
+            index (int): Index of the point to be mutated.
+        '''
+        
+        fill_point = nr.randint(0, self.num_params)
+        candidate_params = self.mutation_func(index)
+        crossovers = nr.rand(self.num_params) < self.cross_over_probability
+        crossovers[fill_point] = True
+        mutated_params = np.where(crossovers, candidate_params, self.population[index])
+        
+        if self.has_trust_region:
+            temp_min = np.maximum(self.min_boundary,self.population[self.min_index] - self.trust_region)
+            temp_max = np.minimum(self.max_boundary,self.population[self.min_index] + self.trust_region)
+            rand_params = temp_min + nr.rand(self.num_params) * (temp_max - temp_min)
+        else:
+            rand_params = self.min_boundary + nr.rand(self.num_params) * self.diff_boundary
+        
+        projected_params = np.where(np.logical_or(mutated_params < self.min_boundary, mutated_params > self.max_boundary), rand_params, mutated_params)
+        
+        return projected_params
+
+    def _best1(self, index):
+        '''
+        Use best parameters and two others to generate mutation.
+        
+        Args:
+            index (int): Index of member to mutate. 
+        '''
+        r0, r1 = self.random_index_sample(index, 2)
+        return (self.population[self.min_index] + self.curr_scale *(self.population[r0] - self.population[r1]))
+
+    def _rand1(self, index):
+        '''
+        Use three random parameters to generate mutation.
+        
+        Args:
+            index (int): Index of member to mutate. 
+        '''
+        r0, r1, r2 = self.random_index_sample(index, 3)
+        return (self.population[r0] + self.curr_scale * (self.population[r1] - self.population[r2]))
+
+    def _best2(self, index):
+        '''
+        Use best parameters and four others to generate mutation.
+        
+        Args:
+            index (int): Index of member to mutate. 
+        '''
+        r0, r1, r2, r3 = self.random_index_sample(index, 4)
+        return self.population[self.min_index] + self.curr_scale * (self.population[r0] + self.population[r1] - self.population[r2] - self.population[r3])
+
+    def _rand2(self, index):
+        '''
+        Use five random parameters to generate mutation.
+        
+        Args:
+            index (int): Index of member to mutate. 
+        '''
+        r0, r1, r2, r3, r4 = self.random_index_sample(index, 5)
+        return self.population[r0] + self.curr_scale * (self.population[r1] + self.population[r2] - self.population[r3] - self.population[r4])
+
+    def random_index_sample(self, index, num_picks):
+        '''
+        Randomly select a num_picks of indexes, without index.
+        
+        Args:
+            index(int): The index that is not included
+            num_picks(int): The number of picks.
+        '''
+        rand_indexes = list(range(self.num_population_members))
+        rand_indexes.remove(index)
+        return random.sample(rand_indexes, num_picks)
+        
+    def update_archive(self):
+        '''
+        Update the archive.
+        '''
+        self.archive_dict.update({'params_generations':self.params_generations,
+                                  'costs_generations':self.costs_generations,
+                                  'population':self.population,
+                                  'population_costs':self.population_costs,
+                                  'init_std':self.init_std,
+                                  'curr_std':self.curr_std,
+                                  'generation_count':self.generation_count})
+
+
 
 class GaussianProcessLearner(Learner, mp.Process):
     '''
@@ -1181,266 +1473,6 @@ def find_local_minima(self):
         self.log.info('Search completed')
 
     
-class DifferentialEvolutionLearner(Learner, threading.Thread):
-    '''
-    Adaption of the differential evolution algorithm in scipy. 
-    
-    Args:
-        params_out_queue (queue): Queue for parameters sent to controller.
-        costs_in_queue (queue): Queue for costs for gaussian process. This must be tuple
-        end_event (event): Event to trigger end of learner.
-        
-    Keyword Args:
-        first_params (Optional [array]): The first parameters to test. If None will just randomly sample the initial condition. Default None.
-        trust_region (Optional [float or array]): The trust region defines the maximum distance the learner will travel from the current best set of parameters. If None, the learner will search everywhere. If a float, this number must be between 0 and 1 and defines maximum distance the learner will venture as a percentage of the boundaries. If it is an array, it must have the same size as the number of parameters and the numbers define the maximum absolute distance that can be moved along each direction. 
-        evolution_strategy (Optional [string]): the differential evolution strategy to use, options are 'best1', 'best1', 'rand1' and 'rand2'. The default is 'best2'.
-        population_size (Optional [int]): multiplier proportional to the number of parameters in a generation. The generation population is set to population_size * parameter_num. Default 15.
-        mutation_scale (Optional [tuple]): The mutation scale when picking new points. Otherwise known as differential weight. When provided as a tuple (min,max) a mutation constant is picked randomly in the interval. Default (0.5,1.0).
-        cross_over_probability (Optional [float]): The recombination constand or crossover probability, the probability a new points will be added to the population.
-        restart_tolerance (Optional [float]): when the current population have a spread less than the initial tolerance, namely stdev(curr_pop) < restart_tolerance stdev(init_pop), it is likely the population is now in a minima, and so the search is started again.
-        
-    Attributes:
-        has_trust_region (bool): Whether the learner has a trust region. 
-        
-    '''
-
-    # Dispatch of mutation strategy method (binomial or exponential).
-    
-    def __init__(self, 
-                 first_params = None,
-                 trust_region = None,
-                 evolution_strategy='rand2', 
-                 population_size=15,
-                 mutation_scale=(0.5, 1), 
-                 cross_over_probability=0.7, 
-                 restart_tolerance=0.01, 
-                 **kwargs):
-        
-        super(NelderMeadLearner,self).__init__(**kwargs)
-        
-        if first_params is None:
-            self.first_params = None
-        else:
-            self.first_params = np.array(first_params, dtype=float)
-            if not self.check_num_params(self.first_params):
-                self.log.error('first_params has the wrong number of parameters:' + repr(self.first_params))
-                raise ValueError
-            if not self.check_in_boundary(self.first_params):
-                self.log.error('first_params is not in the boundary:' + repr(self.first_params))
-                raise ValueError
-        
-        self._set_trust_region(trust_region)
-        
-        if strategy == 'best1':
-            self.mutation_func = self._best1
-        elif strategy == 'best2':
-            self.mutation_func = self._best2
-        elif strategy == 'rand1':
-            self.mutation_func = self._rand1
-        elif strategy == 'rand2':
-            self.mutation_func = self._rand2
-        else:
-            self.log.error('Please select a valid mutation strategy')
-            raise ValueError
-        
-        self.strategy = strategy
-        self.restart_tolerance = restart_tolerance
-
-        if len(mutation_scale) == 2 and (np.any(np.array(mutation_scale) <= 2) or np.any(np.array(mutation_scale) > 0)):
-            self.mutation_scale = mutation_scale
-        else:
-            self.log.error('Mutation scale must be a tuple with (min,max) between 0 and 2. mutation_scale:' + repr(mutation_scale))
-            raise ValueError
-        
-        if cross_over_probability <= 1 and cross_over_probability >= 0:
-            self.cross_over_probability = cross_over_probability
-        else:
-            self.log.error('Cross over probability must be between 0 and 1. cross_over_probability:' + repr(cross_over_probability))
-        
-        if population_size >= 5:
-            self.population_size = population_size
-        else:
-            self.log.error('Population size must be greater or equal to 5:' + repr(population_size))
-        
-        self.num_population_members = self.population_size * self.num_params
-        
-        self.first_sample = True
-        
-        self.params_generations = []
-        self.costs_generations = []
-        
-        self.min_index = 0
-        self.init_std = 0
-        self.curr_std = 0
-    
-    def run(self):
-        '''
-        Runs the Differential Evolution Learner.
-        '''
-        try:
-        
-            generate_population(self)
-        
-            while not self.end_event.is_set():
-                
-                self.next_generation()
-                
-                if self.curr_std < self.restart_tolerance * self.init_std:
-                    self.generate_population()
-            
-        except LearnerInterupt:
-            return
-        
-    def generate_population(self):
-        '''
-        Sample a new random set of variables
-        '''
-        
-        self.population = []
-        self.population_costs = []
-        self.min_index = 0
-        
-        if self.first_params is not None and self.first_sample:
-            curr_params = self.first_params
-        else:
-            curr_params = self.min_boundary + nr.rand(self.num_params) * self.diff_boundary
-        
-        try:
-            curr_cost = self.put_params_and_get_cost(curr_params)
-        except LearnerInterrupt:
-            self.log.info('DELearner ended during first sample of population.')
-            raise
-        
-        self.population.append(curr_params)
-        self.population_costs.append(curr_cost)
-        
-        for index in range(1, self.num_population_members):
-            
-            if self.has_trust_region:
-                temp_min = np.maximum(self.min_boundary,self.population[self.min_index] - self.trust_region)
-                temp_max = np.minimum(self.max_boundary,self.population[self.min_index] + self.trust_region)
-                curr_params = temp_min + nr.rand(self.num_params) * (temp_max - temp_min)
-            else:
-                curr_params = self.min_boundary + nr.rand(self.num_params) * self.diff_boundary
-            
-            try:
-                curr_cost = self.put_params_and_get_cost(curr_params)
-            except LearnerInterrupt:
-                self.log.info('DELearner ended during initial sample of population.')
-                raise
-            
-            self.population.append(curr_params)
-            self.population_costs.append(curr_cost)
-            
-            if curr_cost < self.population_costs[self.min_index]:
-                self.min_index = index
-        
-        self.population = np.array(self.population)
-        self.population_costs = np.array(self.population_costs)
-        
-        self.init_std = self.std(self.population_costs)
-        self.curr_std = self.init_std
-        
-        self.params_generations.append(np.copy(self.population))
-        self.costs_generations.append(np.copy(self.population_costs))
-    
-    def next_generation(self):
-        '''
-        Evolve the population by a single generation
-        '''
-        
-        self.curr_scale = nr.uniform(self.mutation_scale[0], self.mutation_scale[1])
-        
-        for index in range(self.num_population_members):
-            
-            curr_params = self.mutate(index)
-
-            try:
-                curr_cost = self.put_params_and_get_cost(curr_params)
-            except LearnerInterrupt:
-                self.log.info('DELearner ended during initial sample of population.')
-                raise
-            
-            if curr_cost < self.population_costs[index]:
-                self.population[index] = curr_params
-                self.population_costs[index] = curr_cost
-                
-                if curr_cost < self.population_costs[self.min_index]:
-                    self.min_index = index
-        
-        self.curr_std = self.std(self.population_costs)
-        
-        self.params_generations.append(np.copy(self.population))
-        self.costs_generations.append(np.copy(self.population_costs))
-
-    def mutate(self, index):
-        '''
-        Mutate the parameters at index.
-        
-        Args:
-            index (int): Index of the point to be mutated.
-        '''
-        
-        fill_point = nr.randint(0, self.parameter_count)
-        candidate_params = self.mutation_func(index)
-        crossovers = nr.rand(self.parameter_count) < self.cross_over_probability
-        crossovers[fill_point] = True
-        mutated_params = np.where(crossovers, candidate_params, self.population[index])
-        print(mutated_params)
-        
-        if self.has_trust_region:
-            temp_min = np.maximum(self.min_boundary,self.population[self.min_index] - self.trust_region)
-            temp_max = np.minimum(self.max_boundary,self.population[self.min_index] + self.trust_region)
-            rand_params = temp_min + nr.rand(self.num_params) * (temp_max - temp_min)
-        else:
-            rand_params = self.min_boundary + nr.rand(self.num_params) * self.diff_boundary
-        print(rand_params)
-        
-        projected_params = np.where(np.logical_or(mutated_params < self.min_boundary, mutated_params > self.max_boundary),rand_params,mutated_params)
-        print(projected_params)
-        
-        return projected_params
-
-    def _best1(self, index):
-        '''
-        Use best parameters and two others to generate mutation.
-        
-        Args:
-            index (int): Index of member to mutate. 
-        '''
-        r0, r1 = random.sample(range(self.num_population_members).remove(index),2)
-        samples[:2]
-        return (self.population[self.min_index] + self.curr_scale *(self.population[r0] - self.population[r1]))
-
-    def _rand1(self, index):
-        '''
-        Use three random parameters to generate mutation.
-        
-        Args:
-            index (int): Index of member to mutate. 
-        '''
-        r0, r1, r2 = random.sample(range(self.num_population_members).remove(index),3)
-        return (self.population[r0] + self.curr_scale * (self.population[r1] - self.population[r2]))
-
-    def _best2(self, index):
-        '''
-        Use best parameters and four others to generate mutation.
-        
-        Args:
-            index (int): Index of member to mutate. 
-        '''
-        r0, r1, r2, r3 = random.sample(range(self.num_population_members).remove(index),4)
-        return self.population[self.min_index] + self.curr_scale * (self.population[r0] + self.population[r1] - self.population[r2] - self.population[r3])
-
-    def _rand2(self, index):
-        '''
-        Use five random parameters to generate mutation.
-        
-        Args:
-            index (int): Index of member to mutate. 
-        '''
-        r0, r1, r2, r3, r4 = random.sample(range(self.num_population_members).remove(index),5)
-        return self.population[r0] + self.curr_scale * (self.population[r1] + self.population[r2] - self.population[r3] - self.population[r4])
 
 
 
diff --git a/mloop/visualizations.py b/mloop/visualizations.py
index cc9fe04..9f47743 100644
--- a/mloop/visualizations.py
+++ b/mloop/visualizations.py
@@ -11,12 +11,12 @@
 import logging
 import matplotlib.pyplot as plt
 import matplotlib as mpl
-from mloop.controllers import GaussianProcessController
 
 figure_counter = 0
 cmap = plt.get_cmap('hsv')
 run_label = 'Run number'
 cost_label = 'Cost'
+generation_label = 'Generation number'
 scale_param_label = 'Min (0) to max (1) parameters'
 param_label = 'Parameter'
 log_length_scale_label = 'Log of length scale'
@@ -38,12 +38,19 @@ def show_all_default_visualizations(controller, show_plots=True):
     log.debug('Creating controller visualizations.')
     create_contoller_visualizations(controller.total_archive_filename, 
                                     file_type=controller.controller_archive_file_type)
-    if isinstance(controller, GaussianProcessController):
+    
+    if isinstance(controller, mlc.DifferentialEvolutionController):
+        log.debug('Creating differential evolution visualizations.')
+        create_differential_evolution_learner_visualizations(controller.learner.total_archive_filename, 
+                                                             file_type=controller.learner.learner_archive_file_type)
+        
+    if isinstance(controller, mlc.GaussianProcessController):
         log.debug('Creating gaussian process visualizations.')
         plot_all_minima_vs_cost_flag = bool(controller.gp_learner.has_local_minima)
         create_gaussian_process_learner_visualizations(controller.gp_learner.total_archive_filename, 
                                                        file_type=controller.gp_learner.learner_archive_file_type,
                                                        plot_all_minima_vs_cost=plot_all_minima_vs_cost_flag)
+        
     log.info('Showing visualizations, close all to end MLOOP.')
     if show_plots:
         plt.show()
@@ -225,6 +232,111 @@ def plot_parameters_vs_cost(self):
             artists.append(plt.Line2D((0,1),(0,0), color=self.param_colors[ind],marker='o',linestyle=''))
         plt.legend(artists,[str(x) for x in range(1,self.num_params+1)], loc=legend_loc)
 
+def create_differential_evolution_learner_visualizations(filename,
+                                                         file_type='pkl',
+                                                         plot_params_vs_generations=True,
+                                                         plot_costs_vs_generations=True):
+    '''
+    Runs the plots from a differential evolution learner file.
+    
+    Args:
+        filename (Optional [string]): Filename for the differential evolution archive. Must provide datetime or filename. Default None.
+        
+    Keyword Args:
+        file_type (Optional [string]): File type 'pkl' pickle, 'mat' matlab or 'txt' text.
+        plot_params_generations (Optional [bool]): If True plot parameters vs generations, else do not. Default True. 
+        plot_costs_generations (Optional [bool]): If True plot costs vs generations, else do not. Default True. 
+    '''
+    visualization = DifferentialEvolutionVisualizer(filename, file_type=file_type)
+    if plot_params_vs_generations:
+        visualization.plot_params_vs_generations()
+    if plot_costs_vs_generations:
+        visualization.plot_costs_vs_generations()
+
+class DifferentialEvolutionVisualizer():
+    '''
+    DifferentialEvolutionVisualizer creates figures from a differential evolution archive. 
+    
+    Args:
+        filename (String): Filename of the DifferentialEvolutionVisualizer archive.
+    
+    Keyword Args:
+        file_type (String): Can be 'mat' for matlab, 'pkl' for pickle or 'txt' for text. Default 'pkl'. 
+    
+    '''
+    def __init__(self, filename, 
+                 file_type ='pkl', 
+                 **kwargs):
+        
+        self.log = logging.getLogger(__name__)
+        
+        self.filename = str(filename)
+        self.file_type = str(file_type)
+        if not mlu.check_file_type_supported(self.file_type):
+            self.log.error('GP training file type not supported' + repr(self.file_type))
+        learner_dict = mlu.get_dict_from_file(self.filename, self.file_type)
+        
+        if 'archive_type' in learner_dict and not (learner_dict['archive_type'] == 'differential_evolution'):
+            self.log.error('The archive appears to be the wrong type.' + repr(learner_dict['archive_type']))
+            raise ValueError
+        self.archive_type = learner_dict['archive_type']
+        
+        self.num_generations = int(learner_dict['generation_count'])
+        self.num_population_members = int(learner_dict['num_population_members'])
+        self.num_params = int(learner_dict['num_params'])
+        self.min_boundary = np.squeeze(np.array(learner_dict['min_boundary']))
+        self.max_boundary = np.squeeze(np.array(learner_dict['max_boundary']))
+        self.params_generations = np.array(learner_dict['params_generations'])
+        self.costs_generations = np.array(learner_dict['costs_generations'])
+          
+        self.finite_flag = True
+        self.param_scaler = lambda p: (p-self.min_boundary)/(self.max_boundary - self.min_boundary)
+        self.scaled_params_generations = np.array([[self.param_scaler(self.params_generations[inda,indb,:]) for indb in range(self.num_population_members)] for inda in range(self.num_generations)])
+        
+        self.gen_numbers = np.arange(1,self.num_generations+1)
+        self.param_colors = _color_list_from_num_of_params(self.num_params)
+        self.gen_plot = np.array([np.full(self.num_population_members, ind, dtype=int) for ind in self.gen_numbers]).flatten()
+        
+    def plot_costs_vs_generations(self):
+        '''
+        Create a plot of the costs versus run number.
+        '''
+        if self.costs_generations.size == 0:
+            self.log.warning('Unable to plot DE: costs vs generations as the initial generation did not complete.')
+            return
+        
+        global figure_counter, cost_label, generation_label
+        figure_counter += 1
+        plt.figure(figure_counter)
+        plt.plot(self.gen_plot,self.costs_generations.flatten(),marker='o',linestyle='',color='k')
+        plt.xlabel(generation_label)
+        plt.ylabel(cost_label)
+        plt.title('Differential evolution: Cost vs generation number.')
+        
+    def plot_params_vs_generations(self):
+        '''
+        Create a plot of the parameters versus run number.
+        '''
+        if self.params_generations.size == 0:
+            self.log.warning('Unable to plot DE: params vs generations as the initial generation did not complete.')
+            return
+        
+        global figure_counter, generation_label, scale_param_label, legend_loc
+        figure_counter += 1
+        plt.figure(figure_counter)
+        
+        for ind in range(self.num_params):
+            plt.plot(self.gen_plot,self.params_generations[:,:,ind].flatten(),marker='o',linestyle='',color=self.param_colors[ind])
+            plt.ylim((0,1))
+        plt.xlabel(generation_label)
+        plt.ylabel(scale_param_label)
+        
+        plt.title('Differential evolution: Params vs generation number.')
+        artists=[]
+        for ind in range(self.num_params):
+            artists.append(plt.Line2D((0,1),(0,0), color=self.param_colors[ind],marker='o',linestyle=''))
+        plt.legend(artists,[str(x) for x in range(1,self.num_params+1)],loc=legend_loc)
+        
 def create_gaussian_process_learner_visualizations(filename,
                                                    file_type='pkl',
                                                    plot_cross_sections=True,
@@ -234,7 +346,7 @@ def create_gaussian_process_learner_visualizations(filename,
     Runs the plots from a gaussian process learner file.
     
     Args:
-        filename (Optional [string]): Filename for the controller archive. Must provide datetime or filename. Default None.
+        filename (Optional [string]): Filename for the gaussian process archive. Must provide datetime or filename. Default None.
         
     Keyword Args:
         file_type (Optional [string]): File type 'pkl' pickle, 'mat' matlab or 'txt' text.
diff --git a/tests/test_examples.py b/tests/test_examples.py
index e9727e7..d56cbe6 100644
--- a/tests/test_examples.py
+++ b/tests/test_examples.py
@@ -33,6 +33,8 @@ def test_controller_config(self):
     def test_extras_config(self):
         controller = mll.launch_from_file(mlu.mloop_path+'/../examples/extras_config.txt', 
                                           num_params=1,
+                                          min_boundary = [-1.0],
+                                          max_boundary = [1.0],
                                           target_cost = 0.1,
                                           interface_type = 'test',
                                           no_delay = False, 
@@ -42,6 +44,8 @@ def test_extras_config(self):
     def test_logging_config(self):
         controller = mll.launch_from_file(mlu.mloop_path+'/../examples/logging_config.txt',
                                           num_params=1,
+                                          min_boundary = [-1.0],
+                                          max_boundary = [1.0],
                                           target_cost = 0.1,
                                           interface_type = 'test',
                                           no_delay = False, 
@@ -70,6 +74,18 @@ def test_nelder_mead_complete_config(self):
                                           **self.override_dict)
         self.asserts_for_cost_and_params(controller)
     
+    def test_differential_evolution_simple_config(self):
+        controller = mll.launch_from_file(mlu.mloop_path+'/../examples/differential_evolution_simple_config.txt',
+                                          interface_type = 'test',
+                                          **self.override_dict)
+        self.asserts_for_cost_and_params(controller)
+    
+    def test_differential_evolution_complete_config(self):
+        controller = mll.launch_from_file(mlu.mloop_path+'/../examples/differential_evolution_complete_config.txt',
+                                          interface_type = 'test',
+                                          **self.override_dict)
+        self.asserts_for_cost_and_params(controller)
+    
     def test_gaussian_process_simple_config(self):
         controller = mll.launch_from_file(mlu.mloop_path+'/../examples/gaussian_process_simple_config.txt',
                                           interface_type = 'test',

From 4c44e1f8c40947cfcaf50293ea326a1c590bb1ba Mon Sep 17 00:00:00 2001
From: Michael Hush <MichaelRHush@gmail.com>
Date: Wed, 12 Oct 2016 11:34:05 +1100
Subject: [PATCH 5/5] Completed Differential Evolution

Differential Evolution now been added to M-LOOP and is set to the
default trainer for the gaussian process. Tests and examples have been
added. The installation section of the documentation has also been
updated.
---
 docs/contributing.rst |  2 +-
 docs/install.rst      | 78 +++++++++++++++++++++++++++++++++++++++++++--------
 mloop/controllers.py  |  4 +--
 mloop/learners.py     |  3 +-
 4 files changed, 71 insertions(+), 16 deletions(-)

diff --git a/docs/contributing.rst b/docs/contributing.rst
index e7f98a2..984e6f0 100644
--- a/docs/contributing.rst
+++ b/docs/contributing.rst
@@ -5,7 +5,7 @@ Contributing
 
 If you use M-LOOP please consider contributing to the project. There are many quick and easy ways to help out.
 
-- If you use M-LOOP be sure to cite paper where it first used: `'Fast machine-learning online optimization of ultra-cold-atom experiments', Sci Rep 6, 25890 (2016) <http://www.nature.com/srep/>`_. 
+- If you use M-LOOP be sure to cite the paper where it first used: `'Fast machine-learning online optimization of ultra-cold-atom experiments', Sci Rep 6, 25890 (2016) <http://www.nature.com/srep/>`_. 
 - Star and watch the `M-LOOP github <https://github.com/michaelhush/M-LOOP/watchers>`_.
 - Make a suggestion on what features you would like added, or report an issue, on the `github <https://github.com/michaelhush/M-LOOP/watchers>`_ or by `email <mailto:MichaelRHush@gmail.com>`_. 
 - Contribute your own code to the `M-LOOP github <https://github.com/michaelhush/M-LOOP/watchers>`_, this could be the interface you designed, more options or a completely new solver.
diff --git a/docs/install.rst b/docs/install.rst
index 142a56e..f74647c 100644
--- a/docs/install.rst
+++ b/docs/install.rst
@@ -2,13 +2,15 @@
 
 Installation
 ============
-M-LOOP is available on PyPI and can be installed with your favorite package manager. However, we currently recommend you install from the source code to ensure you have the latest improvements and bug fixes. 
+M-LOOP is available on PyPI and can be installed with your favorite package manager. Simply search for 'M-LOOP' and install. For those new to python, we also provide a more comprehensive list of instructions on how to install below.
 
 The installation process involves three steps.
 
 1. Get a Python distribution with the standard scientific packages. We recommend installing :ref:`sec-anaconda`.
-2. Install the development version of :ref:`sec-M-LOOP`.
-3. :ref:`Test<sec-Testing>` your M-LOOP install.
+2. Install the latest release of :ref:`sec-M-LOOP`.
+3. (Optional) :ref:`Test<sec-Testing>` your M-LOOP install.
+
+If you are having any trouble with the installation you may need to check your the :ref:`package dependencies<sec-dependencies>` have been correctly installed. If you ares still having trouble, you can report an issue on the `Link github <https://github.com/michaelhush/M-LOOP/issues>`_. 
 
 .. _sec-anaconda:
 
@@ -20,13 +22,30 @@ https://www.continuum.io/downloads
 
 Follow the installation instructions they provide.
 
-M-LOOP is targeted at python 3.\* but also supports 2.7. Please use python 3.\* if you do not have a reason to use 2.7, see :ref:`sec-py3vpy2` for details.
+M-LOOP is targeted at python 3 but also supports 2. Please use python 3 if you do not have a reason to use 2, see :ref:`sec-py3vpy2` for details.
 
 .. _sec-m-loop:
 
 M-LOOP
 ------
-M-LOOP can be installed from the source code with three commands::
+
+You have two options when installing M-LOOP, you can get the last stable release using pip or you can install from source to get the latest features and bug fixes.
+
+Installing with pip
+^^^^^^^^^^^^^^^^^^^
+
+M-LOOP can be installed with pip with a single command::
+
+   pip install M-LOOP
+   
+If you are using linux or MacOS you may need admin privileges to run the command. To update M-LOOP to the latest version use::
+
+   pip install M-LOOP --upgrade
+
+Installing from source
+^^^^^^^^^^^^^^^^^^^^^^
+
+M-LOOP can be installed from the latest source code with three commands::
 
    git clone git://github.com/michaelhush/M-LOOP.git
    cd ./M-LOOP
@@ -42,19 +61,56 @@ in the M-LOOP directory.
 
 .. _sec-Testing:
 
-Test Installation
------------------
+Testing
+-------
 
-To test your M-LOOP installation use the command::
+If you have installed from source, to test you installation use the command::
 
    python setup.py test
    
 In the M-LOOP source code directory. The tests should take around five minutes to complete. If you find a error please consider :ref:`sec-contributing` to the project and report a bug on the `github <https://github.com/michaelhush/M-LOOP>`_.
 
+If you installed M-LOOP using pip, you will not need to test your installation. 
+
+.. _sec-dependencies:
+
+Dependencies
+------------
+M-LOOP requires the following packages to run correctly.
+
+============   =======
+Package        Version
+============   =======
+docutils       >=0.3
+matplotlib     >=1.5
+numpy          >=1.11
+pip            >=7.0  
+pytest         >=2.9
+setuptools     >=26   
+scikit-learn   >=0.18
+scipy          >=0.17 
+============   =======  
+
+These packages should be automatically installed by pip or the script setup.py when you install M-LOOP.
+
+However if you are using Anaconda some packages that are managed by the conda command may not be correctly updated, even if your installation passes all the tests. In this case you will have to update these packages yourself manually. You can check what packages you have installed and their version with the command::
+
+   conda list
+   
+To install a package that is missing, say for example pytest, use the command::
+
+   conda install pytest
+   
+To update a package to the latest version, say for example scikit-learn, use the command::
+
+   conda update scikit-learn
+
+Once you install and update all the required packages with conda M-LOOP should run correctly. 
+
 Documentation
 -------------
 
-If you would also like a local copy of the documentation enter the docs folder and use the command::
+The latest documentation will always be available here online. If you would also like a local copy of the documentation enter the docs folder and use the command::
 
    make html
    
@@ -65,6 +121,6 @@ Which will generate the documentation in docs/_build/html.
 Python 3 vs 2
 -------------
 
-M-LOOP is developed in python 3.\* and it gets the best performance in this environment. This is primarily because other packages that M-LOOP uses, like numpy, run fastest in python 3. The tests typically take about 20% longer to complete in python 2 than 3.
+M-LOOP is developed in python 3 and it gets the best performance in this environment. This is primarily because other packages that M-LOOP uses, like numpy, run fastest in python 3. The tests typically take about 20% longer to complete in python 2 than 3.
 
-If you have a specific reason to stay in a python 2.7 environment, you may use other packages which are not python 3 compatible, then you can still use M-LOOP without upgrading to 3.\*. However, if you do not have a specific reason to stay with python 2, it is highly recommended you use the latest python 3.\* package.
+If you have a specific reason to stay in a python 2 environment (you may use other packages which are not python 3 compatible) then you can still use M-LOOP without upgrading to 3. However, if you do not have a specific reason to stay with python 2, it is highly recommended you use the latest python 3 package.
diff --git a/mloop/controllers.py b/mloop/controllers.py
index 98d3a23..e4b6964 100644
--- a/mloop/controllers.py
+++ b/mloop/controllers.py
@@ -583,7 +583,7 @@ def __init__(self, interface,
                                                  num_params=num_params,
                                                  min_boundary=min_boundary,
                                                  max_boundary=max_boundary,
-                                                 learner_archive_filename='training_learner_archive',
+                                                 learner_archive_filename=None,
                                                  learner_archive_file_type=learner_archive_file_type,
                                                  **self.remaining_kwargs)
         
@@ -594,7 +594,7 @@ def __init__(self, interface,
                                                             max_boundary=max_boundary,
                                                             trust_region=trust_region,
                                                             evolution_strategy='rand2',
-                                                            learner_archive_filename='training_learner_archive',
+                                                            learner_archive_filename=None,
                                                             learner_archive_file_type=learner_archive_file_type,
                                                             **self.remaining_kwargs)    
         
diff --git a/mloop/learners.py b/mloop/learners.py
index af44655..14556f2 100644
--- a/mloop/learners.py
+++ b/mloop/learners.py
@@ -8,7 +8,6 @@
 
 import threading
 import numpy as np
-import math
 import random
 import numpy.random as nr
 import scipy.optimize as so
@@ -699,7 +698,7 @@ def generate_population(self):
         self.population_costs = []
         self.min_index = 0
         
-        if (not math.isnan(self.first_params)) and self.first_sample:
+        if np.all(np.isfinite(self.first_params)) and self.first_sample:
             curr_params = self.first_params
             self.first_sample = False
         else: