Genetic Programming

1             add
2           /     \
3         mul      sin
4        /   \       \
5       x    2.0      x

  1 # Genetic Programming Core: A Tiny, Dependency-Free Tree-GP Engine
  2 #
  3 # Provides the machinery every tree-based GP example in this chapter needs:
  4 # random tree creation ("full" and "grow" methods), subtree crossover,
  5 # subtree and point mutation, tournament selection, and a generational
  6 # evolve() loop with elitism. The engine knows nothing about what the trees
  7 # *mean* -- each example supplies its own function set, terminal sampler, and
  8 # fitness function, so the same core drives symbolic regression, Boolean
  9 # logic synthesis, program evolution, and L-system rule evolution.
 10 #
 11 # References:
 12 #   GP overview:    https://en.wikipedia.org/wiki/Genetic_programming
 13 #   Koza, J. (1992) "Genetic Programming: On the Programming of Computers
 14 #                   by Means of Natural Selection", MIT Press (full/grow
 15 #                   initialization and the standard benchmark problems)
 16 #   Bloat:          https://en.wikipedia.org/wiki/Genetic_programming#Bloat
 17 
 18 
 19 from __future__ import annotations
 20 
 21 import random
 22 from collections.abc import Callable, Mapping, Sequence
 23 from dataclasses import dataclass, field
 24 
 25 # A terminal sampler is a zero-argument callable that returns one terminal
 26 # symbol (variable name, literal constant, character-class tag, ...).
 27 # Examples use it to draw Ephemeral Random Constants (ERCs) on demand.
 28 TerminalSampler = Callable[[], str]
 29 
 30 
 31 @dataclass
 32 class Tree:
 33     """A GP expression tree: a name plus an ordered list of child subtrees.
 34 
 35     Leaves have no children; their `name` carries the terminal symbol
 36     (for example "x", "0.42", "digit"). Internal nodes name a function
 37     whose arity must agree with len(args).
 38     """
 39 
 40     name: str
 41     args: list[Tree] = field(default_factory=list)
 42 
 43     def copy(self) -> Tree:
 44         """Return a deep copy of this subtree."""
 45         return Tree(self.name, [child.copy() for child in self.args])
 46 
 47     @property
 48     def size(self) -> int:
 49         """Total number of nodes (functions + terminals)."""
 50         return 1 + sum(child.size for child in self.args)
 51 
 52     @property
 53     def depth(self) -> int:
 54         """Longest root-to-leaf path, counted in nodes."""
 55         return 1 + max((child.depth for child in self.args), default=0)
 56 
 57     @property
 58     def is_leaf(self) -> bool:
 59         return not self.args
 60 
 61     def __str__(self) -> str:
 62         """Prefix notation, e.g. (add x (mul 0.5 y))."""
 63         if self.is_leaf:
 64             return self.name
 65         inner = " ".join(str(arg) for arg in self.args)
 66         return f"({self.name} {inner})"
 67 
 68 
 69 # A mutator takes a tree and returns a mutated copy of it.
 70 Mutator = Callable[[Tree], Tree]
 71 
 72 # A fitness function maps a tree to a real number; evolve() MINIMIZES it.
 73 FitnessFunction = Callable[[Tree], float]
 74 
 75 
 76 def sequence_sampler(rng: random.Random, symbols: Sequence[str]) -> TerminalSampler:
 77     """Wrap a fixed terminal list into a TerminalSampler (uniform choice)."""
 78     pool = list(symbols)
 79 
 80     def sample() -> str:
 81         return rng.choice(pool)
 82 
 83     return sample
 84 
 85 
 86 def all_nodes(tree: Tree) -> list[tuple[Tree, Tree | None, int]]:
 87     """Listing of (node, parent, index-in-parent) triples, root first.
 88 
 89     Nodes are returned *by reference* so callers can splice subtrees in
 90     place. The root's parent is None and its index is -1. Every operator
 91     needs only each node together with its parent and index, so the exact
 92     visit order (depth-first, siblings right-to-left) is a nicety, not a
 93     contract.
 94     """
 95     out: list[tuple[Tree, Tree | None, int]] = [(tree, None, -1)]
 96     stack: list[Tree] = [tree]
 97     while stack:
 98         node = stack.pop()
 99         # Appending while iterating backwards lists each parent's children
100         # right-to-left, and pushes the leftmost child last so the walk
101         # itself continues left-to-right.
102         for index in range(len(node.args) - 1, -1, -1):
103             child = node.args[index]
104             out.append((child, node, index))
105             stack.append(child)
106     return out
107 
108 
109 def random_tree(
110     rng: random.Random,
111     functions: Mapping[str, int],
112     terminals: TerminalSampler,
113     max_depth: int,
114     min_depth: int = 1,
115     method: str = "half",
116 ) -> Tree:
117     """Create one random valid tree.
118 
119     Parameters:
120         rng        -- seeded random.Random for reproducibility
121         functions  -- map of function name -> arity (nonterminals)
122         terminals  -- sampler producing terminal symbols
123         max_depth  -- hard cap on root-to-leaf path length (in nodes)
124         min_depth  -- trees must reach at least this depth
125         method     -- "full" (perfect trees), "grow" (mixed shapes) or
126                       "half" (coin flip between the two, Koza's default)
127 
128     The "grow" method stops a branch with probability 0.4 once the minimum
129     depth is reached, so branches end at different depths and the tree comes
130     out irregular -- the shape GP is known for.
131     """
132     if method not in ("full", "grow", "half"):
133         raise ValueError(f"unknown initialization method {method!r}")
134     if max_depth < 1:
135         raise ValueError(f"max_depth must be at least 1, got {max_depth}")
136     if min_depth > max_depth:
137         raise ValueError(f"min_depth {min_depth} exceeds max_depth {max_depth}")
138     if method == "half":
139         method = "full" if rng.random() < 0.5 else "grow"
140 
141     names = list(functions)
142 
143     def build(depth_left: int, depth: int) -> Tree:
144         can_branch = depth_left > 1 and bool(names)
145         if not can_branch:
146             return Tree(terminals())
147         must_branch = depth < min_depth
148         if not must_branch and method == "grow" and rng.random() < 0.4:
149             return Tree(terminals())
150         name = rng.choice(names)
151         arity = functions[name]
152         return Tree(name, [build(depth_left - 1, depth + 1) for _ in range(arity)])
153 
154     return build(max_depth, 1)
155 
156 
157 def crossover(
158     rng: random.Random,
159     parent_a: Tree,
160     parent_b: Tree,
161     max_depth: int | None = None,
162 ) -> tuple[Tree, Tree]:
163     """Subtree crossover: swap one random subtree between two clones.
164 
165     Both offspring are returned; the parents are never modified. If
166     `max_depth` is given, a swap that would make either offspring deeper
167     than the cap is skipped (the clones come back unchanged), so bloat
168     stays bounded without disabling crossover.
169     """
170     child_a, child_b = parent_a.copy(), parent_b.copy()
171 
172     spots_a = all_nodes(child_a)
173     spots_b = all_nodes(child_b)
174     node_a, parent_na, index_a = rng.choice(spots_a)
175     node_b, parent_nb, index_b = rng.choice(spots_b)
176 
177     if max_depth is not None:
178         # Refuse the swap if either side would exceed the depth cap.
179         depth_a_rest = parent_a.depth - node_a.depth  # context above node_a
180         depth_b_rest = parent_b.depth - node_b.depth
181         new_depth_a = depth_a_rest + node_b.depth
182         new_depth_b = depth_b_rest + node_a.depth
183         if max(new_depth_a, new_depth_b) > max_depth:
184             return child_a, child_b
185 
186     # Splicing at the root means the whole other subtree becomes the child.
187     if parent_na is None:
188         child_a = node_b.copy()
189     else:
190         parent_na.args[index_a] = node_b.copy()
191     if parent_nb is None:
192         child_b = node_a.copy()
193     else:
194         parent_nb.args[index_b] = node_a.copy()
195     return child_a, child_b
196 
197 
198 def mutate_subtree(
199     rng: random.Random,
200     tree: Tree,
201     functions: Mapping[str, int],
202     terminals: TerminalSampler,
203     max_depth: int = 4,
204 ) -> Tree:
205     """Subtree mutation: replace one random node's subtree with fresh randoms.
206 
207     The replacement subtree respects `max_depth`, but the whole tree can
208     still end up deeper than the original where the graft happened.
209     """
210     child = tree.copy()
211     node, _parent, _index = rng.choice(all_nodes(child))
212     new = random_tree(rng, functions, terminals, max_depth=max_depth, method="grow")
213     # Overwrite in place: reusing the node object keeps the tree connected.
214     node.name, node.args = new.name, new.args
215     return child
216 
217 
218 def mutate_point(
219     rng: random.Random,
220     tree: Tree,
221     functions: Mapping[str, int],
222     terminals: TerminalSampler,
223 ) -> Tree:
224     """Point mutation: swap a function for another of the SAME arity, or a
225     terminal for another terminal. Structure (shape) is preserved."""
226     child = tree.copy()
227     node, _parent, _index = rng.choice(all_nodes(child))
228     if node.is_leaf:
229         node.name = terminals()
230     else:
231         same_arity = [n for n, a in functions.items() if a == len(node.args)]
232         if same_arity:
233             node.name = rng.choice(same_arity)
234     return child
235 
236 
237 def mutate_hoist(rng: random.Random, tree: Tree) -> Tree:
238     """Hoist mutation: replace a node by one of its own descendants, pruning
239     the surrounding context. A cheap anti-bloat operator."""
240     child = tree.copy()
241     node, _parent, _index = rng.choice(all_nodes(child))
242     pick, _p, _i = rng.choice(all_nodes(node))
243     node.name = pick.name
244     node.args = [arg.copy() for arg in pick.args]
245     return child
246 
247 
248 def tournament(
249     rng: random.Random,
250     fitnesses: Sequence[float],
251     population: Sequence[Tree],
252     k: int = 3,
253 ) -> Tree:
254     """Tournament selection (minimization). Draw k DISTINCT individuals at
255     random, return the one with the LOWEST fitness. Indices refer into the
256     parallel `fitnesses` and `population` sequences."""
257     n = len(population)
258     if n == 0:
259         raise ValueError("cannot select from an empty population")
260     if k < 1:
261         raise ValueError(f"tournament size must be at least 1, got {k}")
262     contenders = rng.sample(range(n), min(k, n))
263     return population[min(contenders, key=fitnesses.__getitem__)]
264 
265 
266 @dataclass
267 class EvolutionResult:
268     """Outcome of one evolve() run.
269 
270     `history[i]` is the best-so-far fitness reported at the start of
271     generation `i`, and the final entry repeats `best_fitness`.
272     """
273 
274     best: Tree
275     best_fitness: float
276     generations: int
277     solved: bool
278     history: list[float] = field(default_factory=list)
279 
280 
281 def evolve(
282     rng: random.Random,
283     population: list[Tree],
284     fitness: FitnessFunction,
285     mutate: Mutator,
286     *,
287     generations: int,
288     crossover_rate: float = 0.8,
289     mutation_rate: float = 0.15,
290     tournament_size: int = 7,
291     elitism: int = 1,
292     max_depth: int = 17,
293     target: float | None = None,
294     reporter: Callable[[int, float, Tree], None] | None = None,
295 ) -> EvolutionResult:
296     """Generational GP loop (minimizes `fitness`).
297 
298     Each generation: rank the population, optionally stop at `target`,
299     carry the best `elitism` trees unchanged, then fill the rest of the
300     next generation by tournament-selected crossover, mutating each child
301     with probability `mutation_rate`.
302 
303     Parameters:
304         rng             -- seeded random source
305         population      -- initial list of Trees (mutated in no way)
306         fitness         -- objective to minimize
307         mutate          -- the mutator chosen by the example
308         generations     -- budget of generations
309         crossover_rate  -- probability of crossover per child (else clone)
310         mutation_rate   -- probability of mutating each crossed/cloned child
311         tournament_size -- k for selection tournaments
312         elitism         -- number of best trees copied unchanged each gen;
313                            an elitism at or above the population size leaves
314                            no room for children, so nothing evolves
315         max_depth       -- depth cap enforced on crossover; mutators apply
316                            their own (relative) limits, so a mutated child
317                            can still exceed this value
318         target          -- stop early once best fitness <= target
319         reporter        -- called as reporter(generation, best_fitness, best_tree)
320 
321     Returns:
322         EvolutionResult with the best tree ever seen. `generations` counts
323         the reproduction rounds actually completed, so an early stop at
324         generation `g` reports `g`, not `g + 1`.
325     """
326     if not population:
327         raise ValueError("population must not be empty")
328     if generations < 0:
329         raise ValueError("generations must not be negative")
330     if tournament_size < 1:
331         raise ValueError("tournament_size must be at least 1")
332     if elitism < 0:
333         raise ValueError("elitism must not be negative")
334     n = len(population)
335     fitnesses = [fitness(tree) for tree in population]
336     best_index = min(range(n), key=fitnesses.__getitem__)
337     best_overall = (fitnesses[best_index], population[best_index].copy())
338     result = EvolutionResult(
339         best=best_overall[1],
340         best_fitness=best_overall[0],
341         generations=0,
342         solved=False,
343     )
344 
345     completed = 0
346     for generation in range(generations):
347         result.history.append(best_overall[0])
348         if reporter is not None:
349             reporter(generation, best_overall[0], best_overall[1])
350         if target is not None and best_overall[0] <= target:
351             result.solved = True
352             break
353 
354         order = sorted(range(n), key=fitnesses.__getitem__)
355         children: list[Tree] = [population[i].copy() for i in order[:elitism]]
356         while len(children) < n:
357             parent_a = tournament(rng, fitnesses, population, tournament_size)
358             parent_b = tournament(rng, fitnesses, population, tournament_size)
359             if rng.random() < crossover_rate:
360                 child, _sibling = crossover(rng, parent_a, parent_b, max_depth=max_depth)
361             else:
362                 child = parent_a.copy()
363             if rng.random() < mutation_rate:
364                 child = mutate(child)
365             children.append(child)
366 
367         population = children
368         fitnesses = [fitness(tree) for tree in population]
369         gen_best_index = min(range(n), key=fitnesses.__getitem__)
370         if fitnesses[gen_best_index] < best_overall[0]:
371             best_overall = (
372                 fitnesses[gen_best_index],
373                 population[gen_best_index].copy(),
374             )
375         completed = generation + 1
376 
377     result.generations = completed
378     result.best = best_overall[1]
379     result.best_fitness = best_overall[0]
380     result.history.append(result.best_fitness)
381     return result

1   x=0.000  y=-0.012794
2   x=0.128  y=+0.511572
3   x=0.256  y=+0.889584
4   x=0.385  y=+1.208017

  1 # Example 1 -- Symbolic Regression: Recovering a Damped Oscillator
  2 #
  3 # GP evolves computer programs (expression trees) that fit noisy numeric
  4 # data. Unlike curve fitting with a FIXED model shape, GP searches over
  5 # model STRUCTURES: the tree decides which operators appear, where, and
  6 # with which constants. The result is a human-readable formula.
  7 #
  8 # This demo hides the true generating function -- a damped sine wave
  9 #     y(x) = 2 * exp(-x/3) * sin(2x) + noise
 10 # sampled at 40 points on [0, 5] -- and asks GP to approximate it. The
 11 # function set below deliberately has NO `exp` primitive, so the evolved
 12 # formula cannot be an exact copy of the truth; GP must find an equally
 13 # good structural approximation. That is typical: symbolic regression
 14 # rewards "right shape", not "textbook answer".
 15 #
 16 # Two classic GP knobs are on display:
 17 #   * Parsimony pressure -- fitness = error + lambda * tree size, which
 18 #     discourages bloat (the runaway growth of useless subtrees).
 19 #   * Protected operators -- division returns 1.0 when the denominator is
 20 #     near zero, so every random tree evaluates to a finite number.
 21 #
 22 # References:
 23 #   Symbolic regression: https://en.wikipedia.org/wiki/Symbolic_regression
 24 #   Bloat control:       Luke & Panait (2006), "A Comparison of Bloat
 25 #                        Control Methods for Genetic Programming",
 26 #                        https://doi.org/10.1162/evco.2006.14.3.309
 27 #   GP field guide:      https://www.gp-field-guide.org.uk/
 28 
 29 from __future__ import annotations
 30 
 31 import math
 32 import random
 33 from collections.abc import Callable, Sequence
 34 from typing import Final
 35 
 36 from gp_core import (
 37     TerminalSampler,
 38     Tree,
 39     evolve,
 40     mutate_point,
 41     mutate_subtree,
 42     random_tree,
 43 )
 44 
 45 SEED: Final = 7  # chosen so this run converges quickly and repeatably
 46 GENERATIONS: Final = 40
 47 
 48 # ---------------------------------------------------------------------------
 49 # The hidden ground truth and the data we are allowed to see.
 50 # ---------------------------------------------------------------------------
 51 
 52 
 53 def target(x: float) -> float:
 54     """The secret: a damped oscillation GP must approximate."""
 55     return 2.0 * math.exp(-x / 3.0) * math.sin(2.0 * x)
 56 
 57 
 58 NOISE_SIGMA: Final = 0.05
 59 XS: Final[list[float]] = [5.0 * i / 39 for i in range(40)]
 60 
 61 
 62 def make_dataset(rng: random.Random) -> list[tuple[float, float]]:
 63     """Sample the target at XS with Gaussian measurement noise."""
 64     return [(x, target(x) + rng.gauss(0.0, NOISE_SIGMA)) for x in XS]
 65 
 66 
 67 # ---------------------------------------------------------------------------
 68 # Function set: name -> (arity, evaluation). The raw operators are partial:
 69 # sin overflows the domain at huge inputs, mul can overflow, div can divide
 70 # by zero. Division is protected here and evaluate() maps every non-finite
 71 # result to 0.0, so every random tree still yields a finite float.
 72 # ---------------------------------------------------------------------------
 73 
 74 ArityAndEval = tuple[int, Callable[[Sequence[float]], float]]
 75 
 76 FUNCTIONS: Final[dict[str, ArityAndEval]] = {
 77     "add": (2, lambda a: a[0] + a[1]),
 78     "sub": (2, lambda a: a[0] - a[1]),
 79     "mul": (2, lambda a: a[0] * a[1]),
 80     "div": (  # protected: |denominator| < 1e-10 -> 1.0
 81         2,
 82         lambda a: a[0] / a[1] if abs(a[1]) > 1e-10 else 1.0,
 83     ),
 84     "sin": (1, lambda a: math.sin(a[0])),
 85     "cos": (1, lambda a: math.cos(a[0])),
 86     "neg": (1, lambda a: -a[0]),
 87     "abs": (1, lambda a: abs(a[0])),
 88 }
 89 
 90 FUNCTION_ARITIES: Final[dict[str, int]] = {
 91     name: arity for name, (arity, _eval) in FUNCTIONS.items()
 92 }
 93 
 94 
 95 def terminal_sampler(rng: random.Random) -> TerminalSampler:
 96     """Terminals = the input variable x plus Ephemeral Random Constants.
 97 
 98     An ERC is re-drawn every time it is needed, so no two trees ever share
 99     the same constant set: GP re-invents its own numeric literals every
100     generation and selection tunes them statistically.
101     """
102 
103     def sample() -> str:
104         if rng.random() < 0.25:
105             return "x"
106         return f"{rng.uniform(-3.0, 3.0):.4f}"
107 
108     return sample
109 
110 
111 # ---------------------------------------------------------------------------
112 # Evaluating a tree on one input, and the parsimony-pressured fitness.
113 # ---------------------------------------------------------------------------
114 
115 
116 def evaluate(tree: Tree, x: float) -> float:
117     """Recursive post-order evaluation. Leaves are 'x' or a float literal."""
118     if tree.is_leaf:
119         return x if tree.name == "x" else float(tree.name)
120     arity, op = FUNCTIONS[tree.name]
121     if len(tree.args) != arity:  # defensive: corrupted tree
122         return 0.0
123     values = [evaluate(child, x) for child in tree.args]
124     result = op(values)
125     return result if math.isfinite(result) else 0.0
126 
127 
128 PARSIMONY: Final = 5e-4  # fitness penalty per tree node
129 
130 
131 class RegressionFitness:
132     """Mean squared error plus a size penalty (parsimony pressure)."""
133 
134     def __init__(self, data: list[tuple[float, float]]) -> None:
135         self.data = data
136 
137     def __call__(self, tree: Tree) -> float:
138         squared = sum((evaluate(tree, x) - y) ** 2 for x, y in self.data)
139         return squared / len(self.data) + PARSIMONY * tree.size
140 
141 
142 # ---------------------------------------------------------------------------
143 # Pretty printing: prefix tree -> readable infix formula.
144 # ---------------------------------------------------------------------------
145 
146 
147 def render(tree: Tree) -> str:
148     """Infix for binaries, call-style for unaries, plain for leaves."""
149     if tree.is_leaf:
150         return tree.name
151     (arity, _op) = FUNCTIONS[tree.name]
152     kids = [render(child) for child in tree.args]
153     if arity == 2:
154         symbol = {"add": " + ", "sub": " - ", "mul": " * ", "div": " / "}.get(tree.name)
155         if symbol is not None:
156             return f"({kids[0]}{symbol}{kids[1]})"
157         return f"{tree.name}({kids[0]}, {kids[1]})"
158     return f"{tree.name}({kids[0]})"
159 
160 
161 # ---------------------------------------------------------------------------
162 # The run.
163 # ---------------------------------------------------------------------------
164 
165 
166 def build_population(rng: random.Random, size: int, max_depth: int) -> list[Tree]:
167     terminals = terminal_sampler(rng)
168     return [
169         random_tree(
170             rng,
171             FUNCTION_ARITIES,
172             terminals,
173             max_depth=max_depth,
174             min_depth=2,
175             method="half",
176         )
177         for _ in range(size)
178     ]
179 
180 
181 def main() -> None:
182     rng = random.Random(SEED)
183     data = make_dataset(rng)
184     scoring = RegressionFitness(data)
185 
186     population = build_population(rng, size=400, max_depth=8)
187     terminals = terminal_sampler(rng)
188 
189     def mutate(tree: Tree) -> Tree:
190         roll = rng.random()
191         if roll < 0.55:
192             return mutate_subtree(rng, tree, FUNCTION_ARITIES, terminals, max_depth=4)
193         return mutate_point(rng, tree, FUNCTION_ARITIES, terminals)
194 
195     def report(generation: int, penalty_fitness: float, best: Tree) -> None:
196         # The run plateaus early, so print checkpoints only.
197         if generation % 5 and generation != GENERATIONS - 1:
198             return
199         mse = sum((evaluate(best, x) - y) ** 2 for x, y in data) / len(data)
200         print(f"gen {generation:3d}  best MSE {mse:8.5f}  size {best.size:3d}  {render(best)[:60]}")
201 
202     print(f"# Symbolic regression: fit {len(data)} noisy samples of a damped oscillator")
203     print(f"# hidden truth: y = 2*exp(-x/3)*sin(2x), noise sd = {NOISE_SIGMA}\n")
204 
205     result = evolve(
206         rng,
207         population,
208         scoring,
209         mutate,
210         generations=GENERATIONS,
211         crossover_rate=0.8,
212         mutation_rate=0.15,
213         tournament_size=7,
214         elitism=2,
215         max_depth=12,
216         reporter=report,
217     )
218 
219     mse = sum((evaluate(result.best, x) - y) ** 2 for x, y in data) / len(data)
220     signal_power = sum(y * y for _, y in data) / len(data)
221     print(f"\nBest after {result.generations} generations:")
222     print(f"  formula   y = {render(result.best)}")
223     print(f"  nodes     {result.best.size}, depth {result.best.depth}")
224     print(
225         f"  MSE {mse:.5f} vs signal power {signal_power:.5f}  "
226         f"(explains {100 * (1 - mse / signal_power):.1f}% of variance)"
227     )
228 
229     print("\n     x      y*      y_hat   error")
230     for x in XS[::4]:
231         y_hat = evaluate(result.best, x)
232         print(f"  {x:5.2f}  {target(x):6.3f}  {y_hat:6.3f}  {y_hat - target(x):+6.3f}")
233 
234 
235 if __name__ == "__main__":
236     main()

1 uv sync
2 make regression    # or: uv run python example_1_symbolic_regression.py

 1 # Symbolic regression: fit 40 noisy samples of a damped oscillator
 2 # hidden truth: y = 2*exp(-x/3)*sin(2x), noise sd = 0.05
 3 
 4 gen   0  best MSE  0.14600  size   7  neg((sin((x + x)) * -0.6569))
 5 gen   5  best MSE  0.08975  size   4  sin((x + x))
 6 gen  10  best MSE  0.08975  size   4  sin((x + x))
 7 gen  15  best MSE  0.08975  size   4  sin((x + x))
 8 gen  20  best MSE  0.08975  size   4  sin((x + x))
 9 gen  25  best MSE  0.08975  size   4  sin((x + x))
10 gen  30  best MSE  0.08975  size   4  sin((x + x))
11 gen  35  best MSE  0.08975  size   4  sin((x + x))
12 gen  39  best MSE  0.08975  size   4  sin((x + x))
13 
14 Best after 40 generations:
15   formula   y = sin((x + x))
16   nodes     4, depth 3
17   MSE 0.08975 vs signal power 0.56206  (explains 84.0% of variance)
18 
19      x      y*      y_hat   error
20    0.00   0.000   0.000  +0.000
21    0.51   1.441   0.855  -0.586
22    1.03   1.260   0.887  -0.373
23    1.54   0.077   0.065  -0.013
24    2.05  -0.827  -0.820  +0.008
25    2.56  -0.778  -0.915  -0.136
26    3.08  -0.092  -0.129  -0.036
27    3.59   0.472   0.781  +0.309
28    4.10   0.478   0.939  +0.461
29    4.62   0.083   0.193  +0.110

  1 # Example 2 -- Boolean Logic Synthesis with Tree GP
  2 #
  3 # Genetic programming on the canonical GP problem class: evolve a PROGRAM
  4 # (here, a combinational logic circuit) whose INPUT/OUTPUT behavior matches
  5 # a specification given as a truth table. Fitness counts wrong rows, so the
  6 # landscape is a staircase of partially-correct circuits -- exactly the kind
  7 # of "behavior, not structure" scoring that distinguishes GP from ordinary
  8 # parameter optimization.
  9 #
 10 # The specification is a machine-safety interlock:
 11 #
 12 #   ENABLE = (both operator keys pressed) AND (no fault)
 13 #            AND (at least one of the two redundant sensors agrees)
 14 #
 15 # Inputs: K1, K2 = operator keys; F = fault flag (active high);
 16 #         S1, S2 = redundant sensors. The interlock needs only ONE of the
 17 #         redundant sensors to agree, so the pair is free to disagree on
 18 #         any row and ENABLE still fires.
 19 # Hand-writing the gate netlist for this is error-prone; GP rediscovers a
 20 # correct netlist from 32 example rows in seconds.
 21 #
 22 # This is also a clean illustration of two GP subtleties:
 23 #   * Ephemeral Boolean constants (TRUE/FALSE leaves) -- without them, GP
 24 #     struggles badly on parity-like targets.
 25 #   * Behavior-only scoring: the evolved netlist in our run is NOT the
 26 #     textbook formula. It is a structurally different circuit that
 27 #     computes the same function, and the demo verifies it row by row.
 28 #     GP optimizes what a program DOES, never how it LOOKS.
 29 #
 30 # References:
 31 #   Koza, J. (1992) "Genetic Programming: On the Programming of Computers
 32 #                   by Means of Natural Selection", MIT Press (Boolean
 33 #                   problems)
 34 #   Parity difficulty in GP: Langdon & Poli (1998), "Why Building Blocks
 35 #                   Don't Work on Parity Problems", CSRP-98-17:
 36 #                   http://web4.cs.ucl.ac.uk/staff/W.Langdon/csrp-98-17/eq.html
 37 #   Logic synthesis:           https://en.wikipedia.org/wiki/Logic_synthesis
 38 
 39 from __future__ import annotations
 40 
 41 import itertools
 42 import random
 43 from collections.abc import Callable, Mapping, Sequence
 44 from typing import Final
 45 
 46 from gp_core import (
 47     TerminalSampler,
 48     Tree,
 49     evolve,
 50     mutate_hoist,
 51     mutate_point,
 52     mutate_subtree,
 53     random_tree,
 54 )
 55 
 56 SEED: Final = 12
 57 
 58 INPUTS: Final[tuple[str, ...]] = ("K1", "K2", "F", "S1", "S2")
 59 
 60 
 61 def specification(row: Sequence[bool]) -> bool:
 62     """The interlock we want GP to rediscover (row order = INPUTS)."""
 63     k1, k2, fault, s1, s2 = row
 64     return (k1 and k2) and not fault and (s1 or s2)
 65 
 66 
 67 TRUTH_TABLE: Final[list[tuple[tuple[bool, ...], bool]]] = [
 68     (row, specification(row)) for row in itertools.product((False, True), repeat=len(INPUTS))
 69 ]
 70 
 71 # The same cases keyed by input name, built once so fitness() does not
 72 # rebuild a dict for every candidate. Evaluating a tree is the hot path.
 73 CASES: Final[list[tuple[Mapping[str, bool], bool]]] = [
 74     (dict(zip(INPUTS, row, strict=True)), expected) for row, expected in TRUTH_TABLE
 75 ]
 76 
 77 # Gate-level function set. NOT is unary; AND/OR binary. (NAND-only synthesis
 78 # is possible too -- try it as an exercise!)
 79 GATES: Final[dict[str, tuple[int, Callable[[Sequence[bool]], bool]]]] = {
 80     "AND": (2, lambda a: a[0] and a[1]),
 81     "OR": (2, lambda a: a[0] or a[1]),
 82     "NOT": (1, lambda a: not a[0]),
 83 }
 84 
 85 FUNCTION_ARITIES: Final[dict[str, int]] = {name: arity for name, (arity, _eval) in GATES.items()}
 86 
 87 
 88 def terminals(rng: random.Random) -> TerminalSampler:
 89     """Inputs plus 20% chance of an Ephemeral Random Boolean constant."""
 90     pool = list(INPUTS)
 91 
 92     def sample() -> str:
 93         if rng.random() < 0.2:
 94             return rng.choice(("TRUE", "FALSE"))
 95         return rng.choice(pool)
 96 
 97     return sample
 98 
 99 
100 def evaluate(tree: Tree, row: Mapping[str, bool]) -> bool:
101     """Interpret the tree as a gate netlist over one input row."""
102     if tree.is_leaf:
103         if tree.name == "TRUE":
104             return True
105         if tree.name == "FALSE":
106             return False
107         return row[tree.name]
108     arity, op = GATES[tree.name]
109     if len(tree.args) != arity:
110         return False
111     return op([evaluate(child, row) for child in tree.args])
112 
113 
114 def fitness(tree: Tree) -> float:
115     """Fraction of truth-table rows the circuit gets WRONG (0.0 = perfect)."""
116     wrong = 0
117     for row, expected in CASES:
118         if evaluate(tree, row) != expected:
119             wrong += 1
120     return float(wrong) / len(CASES)
121 
122 
123 def render(tree: Tree) -> str:
124     """Human-readable gate expression, e.g. (NOT (K1 AND F))."""
125     if tree.is_leaf:
126         return tree.name
127     arity, _op = GATES[tree.name]
128     kids = [render(child) for child in tree.args]
129     if arity == 2:
130         return f"({kids[0]} {tree.name} {kids[1]})"
131     return f"({tree.name} {kids[0]})"
132 
133 
134 def build_population(rng: random.Random, size: int, max_depth: int) -> list[Tree]:
135     sampler = terminals(rng)
136     return [
137         random_tree(
138             rng,
139             FUNCTION_ARITIES,
140             sampler,
141             max_depth=max_depth,
142             min_depth=2,
143             method="grow",
144         )
145         for _ in range(size)
146     ]
147 
148 
149 def main() -> None:
150     rng = random.Random(SEED)
151     population = build_population(rng, size=250, max_depth=8)
152     sampler = terminals(rng)
153 
154     def mutate(tree: Tree) -> Tree:
155         roll = rng.random()
156         if roll < 0.45:
157             return mutate_subtree(rng, tree, FUNCTION_ARITIES, sampler, max_depth=4)
158         if roll < 0.9:
159             return mutate_point(rng, tree, FUNCTION_ARITIES, sampler)
160         return mutate_hoist(rng, tree)
161 
162     def report(generation: int, best_fitness: float, best: Tree) -> None:
163         print(
164             f"gen {generation:3d}  wrong rows {round(best_fitness * len(CASES)):2d}/32  "
165             f"size {best.size:3d}  {render(best)[:58]}"
166         )
167 
168     print("# Evolving a 5-input safety interlock from its 32-row truth table")
169     print("# inputs: K1, K2 (keys), F (fault), S1, S2 (sensors)\n")
170 
171     result = evolve(
172         rng,
173         population,
174         fitness,
175         mutate,
176         generations=60,
177         crossover_rate=0.85,
178         mutation_rate=0.12,
179         tournament_size=5,
180         elitism=1,
181         max_depth=12,
182         target=0.0,
183         reporter=report,
184     )
185 
186     print(
187         f"\nPerfect circuit found: {result.solved} "
188         f"(generation {result.generations}, {result.best.size} gates, depth {result.best.depth})"
189     )
190     print(f"  ENABLE = {render(result.best)}")
191 
192     # Double-check: re-evaluate EVERY row against the discovered netlist.
193     failures = [row for row, expected in CASES if evaluate(result.best, row) != expected]
194     print(f"  verification: {len(CASES) - len(failures)}/{len(CASES)} rows correct")
195 
196     print("\n  truth table (1 = ENABLE):")
197     header = "  " + " ".join(f"{name:>2}" for name in INPUTS) + " | spec | evolved"
198     print(header)
199     print("  " + "-" * (len(header) - 2))
200     for row, expected in CASES:
201         got = evaluate(result.best, row)
202         mark = "  " if got == expected else " <-- MISMATCH"
203         bits = " ".join(f"{int(row[name]):>2}" for name in INPUTS)
204         print(f"  {bits} |  {int(expected)}  |   {int(got)}{mark}")
205 
206 
207 if __name__ == "__main__":
208     main()

1 make logic    # or: uv run python example_2_boolean_logic.py

 1 # Evolving a 5-input safety interlock from its 32-row truth table
 2 # inputs: K1, K2 (keys), F (fault), S1, S2 (sensors)
 3 
 4 gen   0  wrong rows  3/32  size   4  (NOT (NOT (NOT TRUE)))
 5 gen   1  wrong rows  2/32  size  28  ((NOT ((F AND S1) OR ((F AND ((NOT K1) OR (S1 OR S2))) OR 
 6 gen   2  wrong rows  2/32  size  28  ((NOT ((F AND S1) OR ((F AND ((NOT K1) OR (S1 OR S2))) OR 
 7 gen   3  wrong rows  1/32  size  13  (K1 AND ((NOT (F OR (NOT (K2 AND (K1 AND K1))))) AND K1))
 8 gen   4  wrong rows  1/32  size  13  (K1 AND ((NOT (F OR (NOT (K2 AND (K1 AND K1))))) AND K1))
 9 gen   5  wrong rows  1/32  size  13  (K1 AND ((NOT (F OR (NOT (K2 AND (K1 AND K1))))) AND K1))
10 gen   6  wrong rows  1/32  size  13  (K1 AND ((NOT (F OR (NOT (K2 AND (K1 AND K1))))) AND K1))
11 gen   7  wrong rows  1/32  size  13  (K1 AND ((NOT (F OR (NOT (K2 AND (K1 AND K1))))) AND K1))
12 gen   8  wrong rows  0/32  size  18  (((NOT K1) OR (S1 OR S2)) AND ((NOT (F OR (NOT (K2 AND (K2
13 
14 Perfect circuit found: True (generation 8, 18 gates, depth 8)
15   ENABLE = (((NOT K1) OR (S1 OR S2)) AND ((NOT (F OR (NOT (K2 AND (K2 AND K1))))) AND K2))
16   verification: 32/32 rows correct
17 
18   truth table (1 = ENABLE):
19   K1 K2  F S1 S2 | spec | evolved
20   -------------------------------
21    0  0  0  0  0 |  0  |   0  
22    0  0  0  0  1 |  0  |   0  
23    0  0  0  1  0 |  0  |   0  
24    0  0  0  1  1 |  0  |   0  
25    0  0  1  0  0 |  0  |   0  
26    0  0  1  0  1 |  0  |   0  
27    0  0  1  1  0 |  0  |   0  
28    0  0  1  1  1 |  0  |   0  
29    0  1  0  0  0 |  0  |   0  
30    0  1  0  0  1 |  0  |   0  
31    0  1  0  1  0 |  0  |   0  
32    0  1  0  1  1 |  0  |   0  
33    0  1  1  0  0 |  0  |   0  
34    0  1  1  0  1 |  0  |   0  
35    0  1  1  1  0 |  0  |   0  
36    0  1  1  1  1 |  0  |   0  
37    1  0  0  0  0 |  0  |   0  
38    1  0  0  0  1 |  0  |   0  
39    1  0  0  1  0 |  0  |   0  
40    1  0  0  1  1 |  0  |   0  
41    1  0  1  0  0 |  0  |   0  
42    1  0  1  0  1 |  0  |   0  
43    1  0  1  1  0 |  0  |   0  
44    1  0  1  1  1 |  0  |   0  
45    1  1  0  0  0 |  0  |   0  
46    1  1  0  0  1 |  1  |   1  
47    1  1  0  1  0 |  1  |   1  
48    1  1  0  1  1 |  1  |   1  
49    1  1  1  0  0 |  0  |   0  
50    1  1  1  0  1 |  0  |   0  
51    1  1  1  1  0 |  0  |   0  
52    1  1  1  1  1 |  0  |   0  

1 positive: 09:15  23:59  00:01  17:42  08:05  12:00  01:30  20:47
2 negative: 24:00  29:00  09:5  9:15  0915  ab:cd  09:60  1:2:3  ''  09:1a  123:45

  1 # Example 3 -- Evolving Regular Expressions from String Examples
  2 #
  3 # GP where the "program" is a MATCHER instead of a computer: the evolved
  4 # tree is a pattern program that reads a string and returns accept/reject.
  5 # This flips the usual GP setup in an instructive way -- the phenotype is a
  6 # piece of text-processing logic, the fitness cases are example strings,
  7 # and the primitives are regex operators instead of arithmetic.
  8 #
  9 # Task: synthesize a validator for 24-hour clock times ("09:15", "23:59")
 10 # from POSITIVE and NEGATIVE examples only. Nobody tells GP the structure
 11 # of clock arithmetic; it has to discover that the hour's first digit is
 12 # limited to [0-2], that the minute tens digit is limited to [0-5], and so
 13 # on -- purely from the pain of accepting "24:00" or rejecting "09:15".
 14 #
 15 # Fitness is example-weighted error with a smoothing trick: an unmatched
 16 # positive is charged (1 - best_prefix_fraction) instead of a flat 1.0 --
 17 # how much of the string the pattern could match from the left. Accepting
 18 # a negative still costs a full 1.0, so the search cannot cheat by
 19 # accepting every string, while partial credit stops it from getting
 20 # trapped in plateaus like `(?:(?:[0-2]|:))+` (rejects everything, matches
 21 # almost nothing). A small size term breaks ties toward shorter, cleaner
 22 # patterns.
 23 #
 24 # The terminal alphabet is *pre-generalized*: it contains character-class
 25 # fragments like [0-3] and [4-9] alongside single literals, mimicking how
 26 # real GP-regex tools build patterns from library components.
 27 #
 28 # References:
 29 #   Bartoli, De Lorenzo, Medvet & Tarlao, "Playing Regex Golf with Genetic
 30 #   Programming" (GECCO 2014):
 31 #                   https://doi.org/10.1145/2576768.2598333
 32 #   Bartoli, De Lorenzo, Medvet & Tarlao, "On the Automatic Construction
 33 #   of Regular Expressions from Examples (GP vs. Humans 1-0)" (GECCO
 34 #   2016):          https://doi.org/10.1145/2908961.2930946
 35 #   Python re module:           https://docs.python.org/3/library/re.html
 36 
 37 from __future__ import annotations
 38 
 39 import functools
 40 import random
 41 import re
 42 from typing import Final
 43 
 44 from gp_core import (
 45     TerminalSampler,
 46     Tree,
 47     evolve,
 48     mutate_hoist,
 49     mutate_point,
 50     mutate_subtree,
 51     random_tree,
 52 )
 53 
 54 SEED: Final = 1  # chosen so this run converges quickly and repeatably
 55 
 56 POSITIVE: Final[list[str]] = [
 57     "09:15",
 58     "23:59",
 59     "00:01",
 60     "17:42",
 61     "08:05",
 62     "12:00",
 63     "01:30",
 64     "20:47",
 65 ]
 66 NEGATIVE: Final[list[str]] = [
 67     "24:00",  # hour out of range
 68     "29:00",  # hour out of range (also punishes [0-5]-for-hours)
 69     "09:5",  # single-digit minutes
 70     "9:15",  # single-digit hour
 71     "0915",  # missing separator
 72     "ab:cd",  # not digits
 73     "09:60",  # minute out of range
 74     "1:2:3",  # extra separator
 75     "",  # empty input
 76     "09:1a",  # non-digit tail
 77     "123:45",  # too many hour digits
 78 ]
 79 
 80 # (negative string, cost of accepting it). All equal here; the tuple form
 81 # leaves room to weight near-miss cases differently.
 82 NEGATIVE_WEIGHTS: Final[list[tuple[str, float]]] = [(s, 1.0) for s in NEGATIVE]
 83 
 84 # ---------------------------------------------------------------------------
 85 # Pattern-program function set. Every node compiles to a regex fragment;
 86 # unary ops wrap children in (?:...) so precedence is always safe.
 87 # ---------------------------------------------------------------------------
 88 
 89 PATTERN_OPS: Final[dict[str, int]] = {
 90     "cat": 2,  # concatenation:  a b
 91     "alt": 2,  # alternation:    a|b
 92     "opt": 1,  # optional:       a?
 93     "star": 1,  # Kleene star:    a*
 94     "plus": 1,  # one-or-more:    a+
 95 }
 96 
 97 # The quantifier operators, used by the anti-backtracking guard.
 98 QUANTIFIER_OPS: Final = frozenset({"opt", "star", "plus"})
 99 
100 # Terminal alphabet -> regex fragment. Pre-generalized classes do most of
101 # the "cleverness"; GP discovers WHERE they must appear. Deliberately NO
102 # "." wildcard: with prefix credit below, `.` would grant every positive
103 # full credit and flatten the whole fitness landscape into plateaus.
104 TERMINAL_FRAGMENTS: Final[dict[str, str]] = {
105     "D": r"\d",  # any digit
106     "h1": r"[0-1]",  # hour tens digit 0 or 1
107     "h2": r"[0-2]",  # ... alternative tens constraint
108     "h3": r"[0-3]",  # hour units digit when the tens digit is 2
109     "h9": r"[3-9]",  # hour units digit paired with a 0/1 tens digit
110     "m1": r"[0-5]",  # minute tens constraint
111     "m9": r"[6-9]",  # minute tens digit that makes the time invalid
112     "c": r":",  # the separator, as a literal
113     "d0": r"0",
114     "d1": r"1",
115     "d2": r"2",
116     "d5": r"5",
117     "d9": r"9",
118 }
119 
120 TERMINAL_POOL: Final[list[str]] = list(TERMINAL_FRAGMENTS)
121 
122 FUNCTION_ARITIES: Final[dict[str, int]] = dict(PATTERN_OPS)
123 
124 SIZE_PRESSURE: Final = 0.005  # fitness = example error + 0.5% per tree node
125 MAX_QUANTIFIERS: Final = 4  # anti-backtracking guard (see compile_tree)
126 MAX_SOURCE_LEN: Final = 300
127 
128 NEVER_MATCHES: Final = r"(?!a)a"  # impossible pattern: fails instantly
129 
130 
131 def terminals(rng: random.Random) -> TerminalSampler:
132     """Bias leaf sampling toward the meaningful classes, not bare literals."""
133     weights = {"D": 4, "h1": 4, "h2": 3, "h3": 4, "h9": 2, "m1": 4, "m9": 2, "c": 4}
134     pool: list[str] = []
135     for name in TERMINAL_POOL:
136         pool.extend([name] * weights.get(name, 1))
137 
138     def sample() -> str:
139         return rng.choice(pool)
140 
141     return sample
142 
143 
144 # ---------------------------------------------------------------------------
145 # Compile a pattern tree into an anchored regex string (and evaluate it).
146 # ---------------------------------------------------------------------------
147 
148 
149 def quantifier_stats(tree: Tree) -> tuple[int, bool]:
150     """Return (number of quantifier nodes, whether any quantifier nests).
151 
152     `nested` is True when a quantifier node has another quantifier anywhere
153     below it, the shape that makes the backtracking `re` engine exponential
154     even on 6-character strings: `(?:(?:(?:[0-5]*)*)+)*` can hang a CPU
155     forever. One post-order pass answers both questions.
156     """
157     total = 1 if tree.name in QUANTIFIER_OPS else 0
158     nested = False
159     for child in tree.args:
160         child_total, child_nested = quantifier_stats(child)
161         total += child_total
162         if child_nested or (tree.name in QUANTIFIER_OPS and child_total):
163             nested = True
164     return total, nested
165 
166 
167 def compile_tree(tree: Tree) -> str:
168     r"""Tree -> regex source string.
169 
170     Bloat guard: trees with more than `MAX_QUANTIFIERS` quantifier nodes,
171     NESTED quantifiers, or an absurdly long source compile to a
172     never-matching pattern. Without this, degenerate candidates like the
173     greedy `(?:\d|:)+` or `(?:(?:[0-5]*)*)+` make fitness evaluation take
174     exponential time on a six-character string -- a dramatic demonstration
175     that fitness evaluation must stay computable. The exact solution here
176     uses no quantifiers at all, so the guard only removes junk.
177 
178     The guard is a heuristic, not a proof: a single quantifier over an
179     ambiguous alternation (for example one digit class followed by another)
180     can still backtrack. It is enough because these fitness cases are at
181     most six characters.
182     """
183     count, nested = quantifier_stats(tree)
184     if count > MAX_QUANTIFIERS or nested:
185         return NEVER_MATCHES
186     source = _compile_node(tree)
187     return source if len(source) <= MAX_SOURCE_LEN else NEVER_MATCHES
188 
189 
190 def _compile_node(tree: Tree) -> str:
191     if tree.is_leaf:
192         return TERMINAL_FRAGMENTS.get(tree.name, re.escape(tree.name))
193     arity = PATTERN_OPS.get(tree.name)
194     if arity is None:  # defensive: unknown function
195         return NEVER_MATCHES
196     kids = [_compile_node(child) for child in tree.args]
197     if len(kids) != arity:  # defensive: malformed tree
198         return NEVER_MATCHES
199     match tree.name:
200         case "cat":
201             return kids[0] + kids[1]
202         case "alt":
203             return f"(?:{kids[0]}|{kids[1]})"
204         case "opt":
205             return f"(?:{kids[0]})?"
206         case "star":
207             return f"(?:{kids[0]})*"
208         case "plus":
209             return f"(?:{kids[0]})+"
210         case _:
211             raise AssertionError(f"unknown pattern op {tree.name!r}")
212 
213 
214 @functools.lru_cache(maxsize=8192)
215 def compile_cached(source: str) -> re.Pattern[str] | None:
216     """Compile anchored regex text once per distinct string; None if invalid."""
217     try:
218         return re.compile(f"^(?:{source})$")
219     except re.error:
220         return None
221 
222 
223 def prefix_credit(tree: Tree, text: str) -> float:
224     """Fraction of `text` the pattern can consume from the left (1.0 = match).
225 
226     Turns the staircase fitness into a ramp: a pattern that gets half of
227     "23:59" right is closer to the truth than one that cannot even start.
228     Note we must ask "can the WHOLE prefix k be matched?" (fullmatch on
229     text[:k]) -- a plain .match() succeeds as long as the pattern consumes
230     even one character, which would silently grant full credit.
231     """
232     if not text:
233         return 0.0
234     anchored = compile_cached(compile_tree(tree))
235     if anchored is None:
236         return 0.0
237     for k in range(len(text), 0, -1):
238         if anchored.fullmatch(text[:k]):
239             return k / len(text)
240     return 0.0
241 
242 
243 def matches(tree: Tree, text: str) -> bool:
244     pattern = compile_cached(compile_tree(tree))
245     return pattern is not None and pattern.fullmatch(text) is not None
246 
247 
248 def fitness(tree: Tree) -> float:
249     r"""Weighted example error + tiny pressure toward short patterns.
250 
251     Rejecting a positive costs (1 - prefix_credit); accepting a negative
252     costs its full weight -- see module docstring for why this mix keeps
253     both degenerate strategies (the greedy `(?:\d|:)+` and the
254     never-matching pattern) bad while still giving the search somewhere to
255     climb.
256     """
257     error = 0.0
258     for text in POSITIVE:
259         error += 1.0 - prefix_credit(tree, text)
260     for text, weight in NEGATIVE_WEIGHTS:
261         if matches(tree, text):
262             error += weight
263     return error + SIZE_PRESSURE * tree.size
264 
265 
266 def render(tree: Tree) -> str:
267     """Pattern source, kept short for progress lines."""
268     text = compile_tree(tree)
269     return text if len(text) <= 58 else text[:55] + "..."
270 
271 
272 def build_population(rng: random.Random, size: int, max_depth: int) -> list[Tree]:
273     sampler = terminals(rng)
274     return [
275         random_tree(
276             rng,
277             FUNCTION_ARITIES,
278             sampler,
279             max_depth=max_depth,
280             min_depth=3,
281             method="grow",
282         )
283         for _ in range(size)
284     ]
285 
286 
287 def main() -> None:
288     rng = random.Random(SEED)
289     population = build_population(rng, size=200, max_depth=9)
290     sampler = terminals(rng)
291 
292     def mutate(tree: Tree) -> Tree:
293         roll = rng.random()
294         if roll < 0.5:
295             return mutate_subtree(rng, tree, FUNCTION_ARITIES, sampler, max_depth=4)
296         if roll < 0.85:
297             return mutate_point(rng, tree, FUNCTION_ARITIES, sampler)
298         return mutate_hoist(rng, tree)
299 
300     def report(generation: int, best_fitness: float, best: Tree) -> None:
301         accepted = sum(1 for s in POSITIVE if matches(best, s))
302         rejected = sum(1 for s, _w in NEGATIVE_WEIGHTS if not matches(best, s))
303         print(
304             f"gen {generation:3d}  fit {best_fitness:6.3f}  "
305             f"pos {accepted}/{len(POSITIVE)}  neg {rejected}/{len(NEGATIVE)}  {render(best)[:40]}"
306         )
307 
308     print("# Evolving a regular expression for 24-hour clock times")
309     print(f"# {len(POSITIVE)} positive, {len(NEGATIVE)} negative example strings\n")
310 
311     result = evolve(
312         rng,
313         population,
314         fitness,
315         mutate,
316         generations=100,
317         crossover_rate=0.85,
318         mutation_rate=0.12,
319         tournament_size=5,
320         elitism=2,
321         max_depth=12,
322         target=SIZE_PRESSURE * 20,  # zero example error with <= 20 nodes
323         reporter=report,
324     )
325 
326     best_source = compile_tree(result.best)
327     print(f"\nEvolved pattern ({result.best.size} nodes, depth {result.best.depth}):")
328     print(f"  ^{best_source}$")
329     print(f"  all examples satisfied: {result.solved}")
330 
331     print("\n  input       expected  evolved")
332     for text in POSITIVE:
333         print(f"  {text:<10}  MATCH     {'MATCH' if matches(result.best, text) else 'no match'}")
334     for text in NEGATIVE:
335         shown = text if text else "<empty>"
336         verdict = "REJECTED" if not matches(result.best, text) else "accepted (!)"
337         print(f"  {shown:<10}  reject    {verdict}")
338 
339     # Show generalization: strings the GP system never saw.
340     print("\n  unseen probes:")
341     for text in ["05:59", "23:61", "22:00", "07:7", "11:11", "2:03"]:
342         verdict = "MATCH" if matches(result.best, text) else "reject"
343         truth = "valid" if re.fullmatch(r"([01]\d|2[0-3]):[0-5]\d", text) else "invalid"
344         print(f"  {text:<8} -> {verdict:<7} ({truth})")
345 
346 
347 if __name__ == "__main__":
348     main()

1 make regex    # or: uv run python example_3_regex_synthesis.py

 1 # Evolving a regular expression for 24-hour clock times
 2 # 8 positive, 11 negative example strings
 3 
 4 gen   0  fit  5.430  pos 0/8  neg 11/11  (?:9[3-9]|(?:[0-3])+)
 5 gen   1  fit  5.410  pos 0/8  neg 11/11  (?:[0-3])+
 6 gen   2  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
 7 gen   3  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
 8 gen   4  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
 9 gen   5  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
10 gen   6  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
11 gen   7  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
12 gen   8  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
13 gen   9  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
14 gen  10  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
15 gen  11  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
16 gen  12  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
17 gen  13  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
18 gen  14  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
19 gen  15  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
20 gen  16  fit  4.830  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:[0-3])+)
21 gen  17  fit  4.070  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:(?:5:)+|[0-1]\d:))
22 gen  18  fit  3.680  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:(?:[0-3][3-9]|(?:[0-3])
23 gen  19  fit  3.270  pos 0/8  neg 11/11  (?:[0-3][3-9]|(?:(?:5:)+|[0-5]\d:))
24 gen  20  fit  3.260  pos 0/8  neg 11/11  (?:(?:[0-3][3-9]|(?:[0-3])+)(?::)+|(?:[0
25 gen  21  fit  2.730  pos 0/8  neg 11/11  (?:[0-3](?:[0-3](?::)+|:)[3-9]|(?:(?:5:)
26 gen  22  fit  2.730  pos 0/8  neg 11/11  (?:[0-3](?:[0-3](?::)+|:)[3-9]|(?:(?:5:)
27 gen  23  fit  2.720  pos 0/8  neg 11/11  (?:[0-3](?:[0-3](?::)+|:)[3-9]|(?:(?:[0-
28 gen  24  fit  2.350  pos 0/8  neg 11/11  (?:(?:[0-3](?:[0-3](?::)+|[0-3])|(?:[0-3
29 gen  25  fit  1.755  pos 3/8  neg 11/11  (?:(?:[0-3](?:[0-3](?::)+|[0-3])|(?:[0-3
30 gen  26  fit  1.755  pos 3/8  neg 11/11  (?:(?:[0-3](?:[0-3](?::)+|[0-3])|(?:[0-3
31 gen  27  fit  1.755  pos 3/8  neg 11/11  (?:(?:[0-3](?:[0-3](?::)+|[0-3])|(?:[0-3
32 gen  28  fit  1.190  pos 8/8  neg 10/11  (?:[0-1]\d:(?:[0-5])+|(?:(?:(?:[0-3](?:[
33 gen  29  fit  1.170  pos 8/8  neg 10/11  (?:[0-1]\d:(?:[0-5])+|(?:(?:(?:[0-3](?:[
34 gen  30  fit  1.170  pos 8/8  neg 10/11  (?:[0-1]\d:(?:[0-5])+|(?:(?:(?:[0-3](?:[
35 gen  31  fit  1.170  pos 8/8  neg 10/11  (?:[0-1]\d:(?:[0-5])+|(?:(?:(?:[0-3](?:[
36 gen  32  fit  1.165  pos 8/8  neg 10/11  (?:[0-1]\d:(?:[0-5])+|(?:(?:(?:[0-3](?:[
37 gen  33  fit  1.135  pos 8/8  neg 10/11  (?:[0-1]\d:(?:[0-5])+|(?:(?:[0-3](?:[0-3
38 gen  34  fit  1.135  pos 8/8  neg 10/11  (?:[0-1]\d:(?:[0-5])+|(?:(?:[0-3](?:[0-3
39 gen  35  fit  1.125  pos 8/8  neg 10/11  (?:[0-1]\d:(?:[0-5])+|(?:(?:[0-3](?:[0-3
40 gen  36  fit  0.140  pos 8/8  neg 11/11  (?:[0-1]\d:(?:[0-5])+(?:[0-5])+|(?:(?:[0
41 gen  37  fit  0.140  pos 8/8  neg 11/11  (?:[0-1]\d:(?:[0-5])+(?:[0-5])+|(?:(?:[0
42 gen  38  fit  0.130  pos 8/8  neg 11/11  (?:[0-1]\d:(?:[0-5])+(?:[0-5])+|(?:[0-3]
43 gen  39  fit  0.130  pos 8/8  neg 11/11  (?:[0-1]\d:(?:[0-5])+(?:[0-5])+|(?:[0-3]
44 gen  40  fit  0.125  pos 8/8  neg 11/11  (?:[0-1]\d:[0-5](?:[0-5])+|(?:[0-3](?:[0
45 gen  41  fit  0.120  pos 8/8  neg 11/11  (?:[0-1]\d:[0-5](?:[0-5])+|(?:(?:[0-3]\d
46 gen  42  fit  0.120  pos 8/8  neg 11/11  (?:[0-1]\d:[0-5](?:[0-5])+|(?:(?:[0-3]\d
47 gen  43  fit  0.120  pos 8/8  neg 11/11  (?:[0-1]\d:[0-5](?:[0-5])+|(?:(?:[0-3]\d
48 gen  44  fit  0.105  pos 8/8  neg 11/11  (?:[0-1]\d:[0-5][0-5]|(?:[0-3]\d:|[0-5])
49 gen  45  fit  0.105  pos 8/8  neg 11/11  (?:[0-1]\d:[0-5][0-5]|(?:[0-3]\d:|[0-5])
50 gen  46  fit  0.105  pos 8/8  neg 11/11  (?:[0-1]\d:[0-5][0-5]|(?:[0-3]\d:|[0-5])
51 gen  47  fit  0.105  pos 8/8  neg 11/11  (?:[0-1]\d:[0-5][0-5]|(?:[0-3]\d:|[0-5])
52 gen  48  fit  0.095  pos 8/8  neg 11/11  (?:[0-1]\d:[0-5][0-5]|[0-3]\d:\d[3-9])
53 
54 Evolved pattern (19 nodes, depth 6):
55   ^(?:[0-1]\d:[0-5][0-5]|[0-3]\d:\d[3-9])$
56   all examples satisfied: True
57 
58   input       expected  evolved
59   09:15       MATCH     MATCH
60   23:59       MATCH     MATCH
61   00:01       MATCH     MATCH
62   17:42       MATCH     MATCH
63   08:05       MATCH     MATCH
64   12:00       MATCH     MATCH
65   01:30       MATCH     MATCH
66   20:47       MATCH     MATCH
67   24:00       reject    REJECTED
68   29:00       reject    REJECTED
69   09:5        reject    REJECTED
70   9:15        reject    REJECTED
71   0915        reject    REJECTED
72   ab:cd       reject    REJECTED
73   09:60       reject    REJECTED
74   1:2:3       reject    REJECTED
75   <empty>     reject    REJECTED
76   09:1a       reject    REJECTED
77   123:45      reject    REJECTED
78 
79   unseen probes:
80   05:59    -> MATCH   (valid)
81   23:61    -> reject  (invalid)
82   22:00    -> reject  (valid)
83   07:7     -> reject  (invalid)
84   11:11    -> MATCH   (valid)
85   2:03     -> reject  (invalid)