001/*
002 * Java Genetic Algorithm Library (jenetics-8.1.0).
003 * Copyright (c) 2007-2024 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@gmail.com)
019 */
020package io.jenetics;
021
022import static java.lang.Math.abs;
023import static java.lang.String.format;
024import static java.util.Objects.requireNonNull;
025import static io.jenetics.internal.math.Basics.pow;
026import static io.jenetics.internal.math.Basics.ulpDistance;
027
028import java.util.Comparator;
029import java.util.function.Function;
030
031import io.jenetics.internal.math.DoubleAdder;
032import io.jenetics.internal.util.Arrays;
033import io.jenetics.util.ISeq;
034import io.jenetics.util.MSeq;
035import io.jenetics.util.ProxySorter;
036import io.jenetics.util.RandomRegistry;
037import io.jenetics.util.Seq;
038
039/**
040 * Probability selectors are a variation of fitness proportional selectors and
041 * selects individuals from a given population based on its selection
042 * probability <i>P(i)</i>.
043 * <p>
044 * <img src="doc-files/FitnessProportionalSelection.svg" width="400" alt="Selection">
045 * <p>
046 * Fitness proportional selection works as shown in the figure above. The
047 * runtime complexity of the implemented probability selectors is
048 * <i>O(n+</i>log<i>(n))</i> instead of <i>O(n<sup>2</sup>)</i> as for the naive
049 * approach: <i>A binary (index) search is performed on the summed probability
050 * array.</i>
051 *
052 * @author <a href="mailto:franz.wilhelmstoetter@gmail.com">Franz Wilhelmstötter</a>
053 * @since 1.0
054 * @version 4.0
055 */
056public abstract class ProbabilitySelector<
057        G extends Gene<?, G>,
058        C extends Comparable<? super C>
059>
060        implements Selector<G, C>
061{
062        private static final int SERIAL_INDEX_THRESHOLD = 35;
063
064        private static final long MAX_ULP_DISTANCE = pow(10, 10);
065
066        protected final Comparator<Phenotype<G, C>> POPULATION_COMPARATOR = (a, b) ->
067                Optimize.MAXIMUM.<C>descending().compare(a.fitness(), b.fitness());
068
069        protected final boolean _sorted;
070        protected final Function<double[], double[]> _reverter;
071
072
073        /**
074         * Create a new {@code ProbabilitySelector} with the given {@code sorting}
075         * flag. <em>This flag must set to {@code true} if the selector
076         * implementation is sorting the population in the
077         * {@link #probabilities(Seq, int)} method.</em>
078         *
079         * @param sorted {@code true} if the implementation is sorting the
080         *        population when calculating the selection probabilities,
081         *        {@code false} otherwise.
082         */
083        protected ProbabilitySelector(final boolean sorted) {
084                _sorted = sorted;
085                _reverter = sorted ? Arrays::revert : ProbabilitySelector::sortAndRevert;
086        }
087
088        /**
089         * Create a new selector with {@code sorted = false}.
090         */
091        protected ProbabilitySelector() {
092                this(false);
093        }
094
095        @Override
096        public ISeq<Phenotype<G, C>> select(
097                final Seq<Phenotype<G, C>> population,
098                final int count,
099                final Optimize opt
100        ) {
101                requireNonNull(population, "Population");
102                requireNonNull(opt, "Optimization");
103                if (count < 0) {
104                        throw new IllegalArgumentException(format(
105                                "Selection count must be greater or equal then zero, but was %s.",
106                                count
107                        ));
108                }
109
110                final MSeq<Phenotype<G, C>> selection = MSeq
111                        .ofLength(population.isEmpty() ? 0 : count);
112
113                if (count > 0 && !population.isEmpty()) {
114                        final Seq<Phenotype<G, C>> pop = _sorted
115                                ? population.asISeq().copy().sort(POPULATION_COMPARATOR)
116                                : population;
117
118
119                        final double[] prob = probabilities(pop, count, opt);
120                        assert pop.size() == prob.length
121                                : "Population size and probability length are not equal.";
122
123                        checkAndCorrect(prob);
124                        assert sum2one(prob) : "Probabilities doesn't sum to one.";
125
126                        incremental(prob);
127
128                        final var random = RandomRegistry.random();
129                        selection.fill(() -> pop.get(indexOf(prob, random.nextDouble())));
130                }
131
132                return selection.toISeq();
133        }
134
135        /**
136         * This method takes the probabilities from the
137         * {@link #probabilities(Seq, int)} method and inverts it if needed.
138         *
139         * @param population The population.
140         * @param count The number of phenotypes to select.
141         * @param opt Determines whether the individuals with higher fitness values
142         *        or lower fitness values must be selected. This parameter
143         *        determines whether the GA maximizes or minimizes the fitness
144         *        function.
145         * @return Probability array.
146         */
147        protected final double[] probabilities(
148                final Seq<Phenotype<G, C>> population,
149                final int count,
150                final Optimize opt
151        ) {
152                return requireNonNull(opt) == Optimize.MINIMUM
153                        ? _reverter.apply(probabilities(population, count))
154                        : probabilities(population, count);
155        }
156
157        // Package private for testing.
158        static double[] sortAndRevert(final double[] array) {
159                final int[] indexes = ProxySorter.sort(array);
160
161                // Copy the elements in reversed order.
162                final double[] result = new double[array.length];
163                for (int i = 0; i < result.length; ++i) {
164                        result[indexes[i]] = array[indexes[result.length - 1 - i]];
165                }
166
167                return result;
168        }
169
170        /**
171         * <p>
172         * Return an Probability array, which corresponds to the given Population.
173         * The probability array and the population must have the same size. The
174         * population is not sorted. If a subclass needs a sorted population, the
175         * subclass is responsible to sort the population.
176         * </p>
177         * The implementer always assumes that higher fitness values are better. The
178         * base class inverts the probabilities, by reverting the returned
179         * probability array, if the GA is supposed to minimize the fitness function.
180         *
181         * @param population The <em>unsorted</em> population.
182         * @param count The number of phenotypes to select. <i>This parameter is not
183         *        needed for most implementations.</i>
184         * @return Probability array. The returned probability array must have the
185         *         length {@code population.size()} and <strong>must</strong> sum to
186         *         one. The returned value is checked with
187         *         {@code assert(Math.abs(math.sum(probabilities) - 1.0) < 0.0001)}
188         *         in the base class.
189         */
190        protected abstract double[] probabilities(
191                final Seq<Phenotype<G, C>> population,
192                final int count
193        );
194
195        /**
196         * Checks if the given probability values are finite. If not, all values are
197         * set to the same probability.
198         *
199         * @param probabilities the probabilities to check.
200         */
201        private static void checkAndCorrect(final double[] probabilities) {
202                boolean ok = true;
203                for (int i = probabilities.length; --i >= 0 && ok;) {
204                        ok = Double.isFinite(probabilities[i]);
205                }
206
207                if (!ok) {
208                        final double value = 1.0/probabilities.length;
209                        for (int i = probabilities.length; --i >= 0;) {
210                                probabilities[i] = value;
211                        }
212                }
213        }
214
215        /**
216         * Check if the given probabilities sum to one.
217         *
218         * @param probabilities the probabilities to check.
219         * @return {@code true} if the sum of the probabilities is within the error
220         *         range, {@code false} otherwise.
221         */
222        static boolean sum2one(final double[] probabilities) {
223                final double sum = probabilities.length > 0
224                        ? DoubleAdder.sum(probabilities)
225                        : 1.0;
226                return abs(ulpDistance(sum, 1.0)) < MAX_ULP_DISTANCE;
227        }
228
229        static boolean eq(final double a, final double b) {
230                return abs(ulpDistance(a, b)) < MAX_ULP_DISTANCE;
231        }
232
233        static int indexOf(final double[] incr, final double v) {
234                return incr.length <= SERIAL_INDEX_THRESHOLD
235                        ? indexOfSerial(incr, v)
236                        : indexOfBinary(incr, v);
237        }
238
239        /**
240         * Perform a binary-search on the summed probability array.
241         */
242        static int indexOfBinary(final double[] incr, final double v) {
243                int imin = 0;
244                int imax = incr.length;
245                int index = -1;
246
247                while (imax > imin && index == -1) {
248                        final int imid = (imin + imax) >>> 1;
249
250                        if (imid == 0 || (incr[imid] >= v && incr[imid - 1] < v)) {
251                                index = imid;
252                        } else if (incr[imid] <= v) {
253                                imin = imid + 1;
254                        } else if (incr[imid] > v) {
255                                imax = imid;
256                        }
257                }
258
259                return index;
260        }
261
262        /**
263         * Perform a serial-search on the summed probability array.
264         */
265        static int indexOfSerial(final double[] incr, final double v) {
266                int index = -1;
267                for (int i = 0; i < incr.length && index == -1; ++i) {
268                        if (incr[i] >= v) {
269                                index = i;
270                        }
271                }
272
273                return index;
274        }
275
276        /**
277         * In-place summation of the probability array.
278         */
279        static double[] incremental(final double[] values) {
280                final DoubleAdder adder = new DoubleAdder(values[0]);
281                for (int i = 1; i < values.length; ++i) {
282                        values[i] = adder.add(values[i]).doubleValue();
283                }
284                return values;
285        }
286
287}