ProbabilitySelector.java
001 /*
002  * Java Genetic Algorithm Library (jenetics-2.0.2).
003  * Copyright (c) 2007-2014 Franz Wilhelmstötter
004  *
005  * Licensed under the Apache License, Version 2.0 (the "License");
006  * you may not use this file except in compliance with the License.
007  * You may obtain a copy of the License at
008  *
009  *      http://www.apache.org/licenses/LICENSE-2.0
010  *
011  * Unless required by applicable law or agreed to in writing, software
012  * distributed under the License is distributed on an "AS IS" BASIS,
013  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
014  * See the License for the specific language governing permissions and
015  * limitations under the License.
016  *
017  * Author:
018  *    Franz Wilhelmstötter (franz.wilhelmstoetter@gmx.at)
019  */
020 package org.jenetics;
021 
022 import static java.lang.Math.abs;
023 import static java.lang.String.format;
024 import static java.util.Objects.requireNonNull;
025 import static org.jenetics.util.math.pow;
026 import static org.jenetics.util.math.statistics.sum;
027 import static org.jenetics.util.math.ulpDistance;
028 
029 import java.util.Random;
030 
031 import org.jenetics.util.Factory;
032 import org.jenetics.util.Function;
033 import org.jenetics.util.RandomRegistry;
034 
035 
036 /**
037  * Probability selectors are a variation of fitness proportional selectors and
038  * selects individuals from a given population based on it's selection
039  * probability <i>P(i)</i>.
040  <p>
041  <img src="doc-files/FitnessProportionalSelection.svg" width="400" alt="Selection">
042  <p>
043  * Fitness proportional selection works as shown in the figure above. The
044  * runtime complexity of the implemented probability selectors is
045  <i>O(n+</i>log<i>(n))</i> instead of <i>O(n<sup>2</sup>)</i> as for the naive
046  * approach: <i>A binary (index) search is performed on the summed probability
047  * array.</i>
048  *
049  @author <a href="mailto:franz.wilhelmstoetter@gmx.at">Franz Wilhelmstötter</a>
050  @since 1.0
051  @version 2.0 — <em>$Date: 2014-08-27 $</em>
052  */
053 public abstract class ProbabilitySelector<
054     extends Gene<?, G>,
055     extends Comparable<? super C>
056 >
057     implements Selector<G, C>
058 {
059     private static final long MAX_ULP_DISTANCE = pow(1010);
060 
061     private final Function<double[]double[]> _revert;
062 
063     /**
064      * Create a new {@code ProbabilitySelector} with the given {@code sorting}
065      * flag. <em>This flag must set to {@code true} if the selector
066      * implementation is sorting the population in the
067      {@link #probabilities(Population, int)} method.</em>
068      *
069      @param sorted {@code true} if the implementation is sorting the
070      *        population when calculating the selection probabilities,
071      *        {@code false} otherwise.
072      */
073     protected ProbabilitySelector(final boolean sorted) {
074         _revert = sorted ?
075             new Function<double[]double[]>() {
076                 @Override public double[] apply(final double[] values) {
077                     return revert(values);
078                 }
079             :
080             new Function<double[]double[]>() {
081                 @Override public double[] apply(final double[] values) {
082                     return sortAndRevert(values);
083                 }
084             };
085     }
086 
087     /**
088      * Create a new selector with {@code sorted = false}.
089      */
090     protected ProbabilitySelector() {
091         this(false);
092     }
093 
094     @Override
095     public Population<G, C> select(
096         final Population<G, C> population,
097         final int count,
098         final Optimize opt
099     ) {
100         requireNonNull(population, "Population");
101         requireNonNull(opt, "Optimization");
102         if (count < 0) {
103             throw new IllegalArgumentException(format(
104                 "Selection count must be greater or equal then zero, but was %s.",
105                 count
106             ));
107         }
108 
109         final Population<G, C> selection = new Population<>(count);
110 
111         if (count > 0) {
112             final double[] probabilities = probabilities(population, count, opt);
113             assert (population.size() == probabilities.length:
114                 "Population size and probability length are not equal.";
115             assert (sum2one(probabilities)) "Probabilities doesn't sum to one.";
116 
117             incremental(probabilities);
118             final Factory<Phenotype<G, C>> factory = factory(
119                 population, probabilities, RandomRegistry.getRandom()
120             );
121 
122             selection.fill(factory, count);
123             assert (count == selection.size());
124         }
125 
126         return selection;
127     }
128 
129     private static <
130         extends Gene<?, G>,
131         extends Comparable<? super C>
132     >
133     Factory<Phenotype<G, C>> factory(
134         final Population<G, C> population,
135         final double[] probabilities,
136         final Random random
137     ) {
138         return new Factory<Phenotype<G, C>>() {
139             @Override
140             public Phenotype<G, C> newInstance() {
141                 return select(population, probabilities, random);
142             }
143         };
144     }
145 
146     private static <
147         extends Gene<?, G>,
148         extends Comparable<? super C>
149     >
150     Phenotype<G, C> select(
151         final Population<G, C> population,
152         final double[] probabilities,
153         final Random random
154     ) {
155         final double value = random.nextDouble();
156         return population.get(indexOf(probabilities, value));
157     }
158 
159     /**
160      * This method takes the probabilities from the
161      {@link #probabilities(Population, int)} method and inverts it if needed.
162      *
163      @param population The population.
164      @param count The number of phenotypes to select.
165      @param opt Determines whether the individuals with higher fitness values
166      *        or lower fitness values must be selected. This parameter
167      *        determines whether the GA maximizes or minimizes the fitness
168      *        function.
169      @return Probability array.
170      */
171     protected final double[] probabilities(
172         final Population<G, C> population,
173         final int count,
174         final Optimize opt
175     ) {
176         return requireNonNull(opt== Optimize.MINIMUM ?
177             _revert.apply(probabilities(population, count)) :
178             probabilities(population, count);
179     }
180 
181     // Package private for testing.
182     static double[] sortAndRevert(final double[] probabilities) {
183         final int N = probabilities.length;
184         final int[] indexes = sort(probabilities);
185         final double[] result = new double[N];
186 
187         for (int i = 0; i < N; ++i) {
188             result[indexes[N - i - 1]] = probabilities[indexes[i]];
189         }
190 
191         return result;
192     }
193 
194     private static final int INSERTION_SORT_THRESHOLD = 75;
195 
196     // Package private for testing.
197     static int[] sort(final double[] values) {
198         return values.length < INSERTION_SORT_THRESHOLD ?
199             insertionSort(values:
200             quickSort(values);
201     }
202 
203     private static int[] indexes(final int length) {
204         final int[] indexes = new int[length];
205         for (int i = 0; i < indexes.length; ++i) {
206             indexes[i= i;
207         }
208         return indexes;
209     }
210 
211     // Package private for testing.
212     static int[] quickSort(final double[] array) {
213         final int[] indexes = indexes(array.length);
214         quickSort(array, indexes, 0, array.length - 1);
215         return indexes;
216     }
217 
218     private static void quickSort(
219         final double[] array,
220         final int[] indexes,
221         final int left, final int right
222     ) {
223         if (right > left) {
224             final int j = partition(array, indexes, left, right);
225             quickSort(array, indexes, left, j - 1);
226             quickSort(array, indexes, j + 1, right);
227         }
228     }
229 
230     private static int partition(
231         final double[] array, final int[] indexes,
232         final int left, final int right
233     ) {
234         final double pivot = array[indexes[left]];
235         int i = left;
236         int j = right + 1;
237 
238         while (true) {
239             do ++i; while (i < right && array[indexes[i]] < pivot);
240             do --j; while (j > left && array[indexes[j]] > pivot);
241             if (j <= ibreak;
242             swap(indexes, i, j);
243         }
244         swap(indexes, left, j);
245 
246         return j;
247     }
248 
249     private static void swap(final int[] indexes, final int i, final int j) {
250         final int temp = indexes[i];
251         indexes[i= indexes[j];
252         indexes[j= temp;
253     }
254 
255     // Package private for testing.
256     static int[] insertionSort(final double[] array) {
257         final int[] indexes = indexes(array.length);
258 
259         for (int sz = array.length, i = 1; i < sz; ++i) {
260             int j = i;
261             while (j > 0) {
262                 if (array[indexes[j - 1]] > array[indexes[j]]) {
263                     swap(indexes, j - 1, j);
264                 else {
265                     break;
266                 }
267                 --j;
268             }
269         }
270 
271         return indexes;
272     }
273 
274     /**
275      <p>
276      * Return an Probability array, which corresponds to the given Population.
277      * The probability array and the population must have the same size. The
278      * population is not sorted. If a subclass needs a sorted population, the
279      * subclass is responsible to sort the population.
280      </p>
281      * The implementer always assumes that higher fitness values are better. The
282      * base class inverts the probabilities ({@code p = 1.0 - p }) if the GA is
283      * supposed to minimize the fitness function.
284      *
285      @param population The <em>unsorted</em> population.
286      @param count The number of phenotypes to select. <i>This parameter is not
287      *        needed for most implementations.</i>
288      @return Probability array. The returned probability array must have the
289      *         length {@code population.size()} and <strong>must</strong> sum to
290      *         one. The returned value is checked with
291      *         {@code assert(Math.abs(math.sum(probabilities) - 1.0) < 0.0001)}
292      *         in the base class.
293      */
294     protected abstract double[] probabilities(
295         final Population<G, C> population,
296         final int count
297     );
298 
299     /**
300      * Check if the given probabilities sum to one.
301      *
302      @param probabilities the probabilities to check.
303      @return {@code true} if the sum of the probabilities are within the error
304      *         range, {@code false} otherwise.
305      */
306     static boolean sum2one(final double[] probabilities) {
307         final double sum = sum(probabilities);
308         return abs(ulpDistance(sum, 1.0)) < MAX_ULP_DISTANCE;
309     }
310 
311     /**
312      * Perform a binary-search on the summed probability array.
313      */
314     static int indexOf(final double[] incremental, final double v) {
315         int imin = 0;
316         int imax = incremental.length;
317 
318         while (imax > imin) {
319             int imid = (imin + imax>>> 1;
320 
321             if (imid == 0) {
322                 return imid;
323             else if (incremental[imid>= v && incremental[imid - 1< v) {
324                 return imid;
325             else if (incremental[imid<= v) {
326                 imin = imid + 1;
327             else if (incremental[imid> v) {
328                 imax = imid;
329             }
330         }
331 
332         return incremental.length - 1;
333     }
334 
335     /**
336      * In-place summation of the probability array.
337      */
338     static double[] incremental(final double[] values) {
339         for (int i = 1; i < values.length; ++i) {
340             values[i= values[i - 1+ values[i];
341         }
342         return values;
343     }
344 
345     public static double[] revert(final double[] array) {
346         for (int i = 0, j = array.length - 1; i < j; ++i, --j) {
347             swap(array, i, j);
348         }
349 
350         return array;
351     }
352 
353     public static void swap(final double[] array, final int i, final int j) {
354         final double temp = array[i];
355         array[i= array[j];
356         array[j= temp;
357     }
358 
359 }