Merge pull request #191 from jakobj/enh/mul-add-parameters-3.0

Allow multiple parameters for nodes and add example with parametrized additive node

.travis.yml

@@ -40,6 +40,7 @@ script:
        python examples/example_evo_regression.py || exit 1 ;
        python examples/example_differential_evo_regression.py || exit 1;
        python examples/example_caching.py || exit 1;
+       python examples/example_parametrized_nodes.py || exit 1;
     fi
 after_success:
   - coveralls

cgp/cartesian_graph.py

@@ -320,7 +320,7 @@ def to_torch(self) -> "torch.nn.Module":
             node.format_output_str_torch(self)
 
         active_nodes_by_hidden_column_idx = self._determine_active_nodes()
-        all_parameter_str = []
+        all_parameter_str: List[List[str]] = []
         for hidden_column_idx in sorted(active_nodes_by_hidden_column_idx):
             for node in active_nodes_by_hidden_column_idx[hidden_column_idx]:
                 node.format_output_str_torch(self)
@@ -337,8 +337,9 @@ def __init__(self):
         super().__init__()
 
 """
-        for s in all_parameter_str:
-            class_str += "        " + s
+        for parameter_str in all_parameter_str:
+            for s in parameter_str:
+                class_str += "        " + s + "\n"
 
         func_str = f"""\
 

cgp/node.py

@@ -172,7 +172,7 @@ class OperatorNode(Node):
     _def_torch_output: str
     _def_sympy_output: str
 
-    _parameter_str: str
+    _parameter_str: List[str]
 
     def __init_subclass__(cls: Type["OperatorNode"]) -> None:
         super().__init_subclass__()
@@ -208,7 +208,7 @@ def initial_value(cls, parameter_name: str) -> float:
         return cls._initial_values["<" + parameter_prefix + ">"]()
 
     @property
-    def parameter_str(self) -> str:
+    def parameter_str(self) -> List[str]:
         return self._parameter_str
 
     def __call__(self, x: List[float], graph: "CartesianGraph") -> None:
@@ -253,10 +253,10 @@ def format_output_str_numpy(self, graph: "CartesianGraph") -> None:
 
     def format_output_str_torch(self, graph: "CartesianGraph") -> None:
         if not hasattr(self, "_def_torch_output"):
-            self.format_output_str(graph)
+            output_str = self._format_output_str(self._def_output, graph)
         else:
             output_str = self._format_output_str(self._def_torch_output, graph)
-            self._output_str = self._replace_parameter_names_in_output_str_with_members(output_str)
+        self._output_str = self._replace_parameter_names_in_output_str_with_members(output_str)
 
     def _replace_parameter_names_in_output_str_with_members(self, output_str: str) -> str:
         g = re.findall("<([a-z]+[0-9]+)>", output_str)
@@ -273,11 +273,11 @@ def format_output_str_sympy(self, graph: "CartesianGraph") -> None:
             self._output_str = self._format_output_str(self._def_sympy_output, graph)
 
     def format_parameter_str(self) -> None:
-        parameter_str_list = []
+        parameter_str = []
         for parameter_name in self._parameter_names:
             parameter_name_with_idx = parameter_name[1:-1] + str(self._idx)
-            parameter_str_list.append(
+            parameter_str.append(
                 f"self._{parameter_name_with_idx} = torch.nn.Parameter("
-                + f"torch.DoubleTensor([<{parameter_name_with_idx}>]))\n"
+                + f"torch.DoubleTensor([<{parameter_name_with_idx}>]))"
             )
-        self._parameter_str = "\n".join(parameter_str_list)
+        self._parameter_str = parameter_str

examples/example_parametrized_nodes.py

@@ -0,0 +1,172 @@
+"""
+Example for evolutionary regression with parametrized nodes
+===========================================================
+
+Example demonstrating the use of Cartesian genetic programming for a
+regression task that requires fine tuning of constants in parametrized
+nodes. This is achieved by introducing a new node, "ParametrizedAdd"
+which produces a scaled and shifted version of the sum of its inputs.
+
+"""
+
+import functools
+import math
+import matplotlib.pyplot as plt
+import numpy as np
+import scipy.constants
+import torch
+
+import cgp
+
+
+# %%
+# We first define a new node that adds two inputs then scales and
+# finally shifts the result. The scale ("w") and shift factors ("b")
+# are parameters that are adapted by local search. We need to define
+# the arity of the node, callables for the initial values for the
+# parameters and the operation of the node as a string. In this string
+# parameters are enclosed in angle brackets, inputs are denoted by "x_i"
+# with i representing their corresponding index.
+
+
+class ParametrizedAdd(cgp.OperatorNode):
+    """A node that adds its two inputs.
+
+    The result of addition is scaled by w and shifted by b. Both these
+    parameters can be adapted via local search are passed on from
+    parents to their offspring.
+
+    """
+
+    _arity = 2
+    _initial_values = {"<w>": lambda: 1.0, "<b>": lambda: 0.0}
+    _def_output = "<w> * (x_0 + x_1) + <b>"
+
+
+# %%
+# We define a target function which contains numerical constants that
+# are not available as constants for the search and need to be found
+# by local search on parameterized nodes.
+
+
+def f_target(x):
+    return math.pi * (x[:, 0] + x[:, 1]) + math.e
+
+
+# %%
+# Then we define a differentiable(!) inner objective function for the
+# evolution. This function accepts a torch class as a parameter. It
+# returns the mean-squared error between the output of the forward
+# pass of this class and the target function evaluated on a set of
+# random points. This inner objective is then used by actual objective
+# function to determine the fitness of the individual.
+
+
+def inner_objective(f, seed):
+    torch.manual_seed(seed)
+    batch_size = 500
+    x = torch.DoubleTensor(batch_size, 2).uniform_(-5, 5)
+    y = f(x)
+    return torch.nn.MSELoss()(f_target(x), y[:, 0])
+
+
+def objective(individual, seed):
+    if individual.fitness is not None:
+        return individual
+
+    f = individual.to_torch()
+    loss = inner_objective(f, seed)
+    individual.fitness = -loss.item()
+
+    return individual
+
+
+# %%
+# Next, we define the parameters for the population, the genome of
+# individuals, the evolutionary algorithm, and the local search. Note
+# that we add the custom node defined above as a primitive.
+
+population_params = {"n_parents": 1, "mutation_rate": 0.04, "seed": 818821}
+
+genome_params = {
+    "n_inputs": 2,
+    "n_outputs": 1,
+    "n_columns": 5,
+    "n_rows": 1,
+    "levels_back": None,
+    "primitives": (ParametrizedAdd, cgp.Add, cgp.Sub, cgp.Mul),
+}
+
+ea_params = {"n_offsprings": 4, "tournament_size": 1, "n_processes": 2}
+
+evolve_params = {"max_generations": 500, "min_fitness": 0.0}
+
+local_search_params = {"lr": 1e-3, "gradient_steps": 9}
+
+# %%
+# We then create a Population instance and instantiate the local search
+# and evolutionary algorithm.
+
+pop = cgp.Population(**population_params, genome_params=genome_params)
+
+local_search = functools.partial(
+    cgp.local_search.gradient_based,
+    objective=functools.partial(inner_objective, seed=population_params["seed"]),
+    **local_search_params,
+)
+
+ea = cgp.ea.MuPlusLambda(**ea_params, local_search=local_search)
+
+
+# %%
+# We define a recording callback closure for bookkeeping of the progression of the evolution.
+
+history = {}
+history["fitness_champion"] = []
+history["expr_champion"] = []
+
+
+def recording_callback(pop):
+    history["fitness_champion"].append(pop.champion.fitness)
+    history["expr_champion"].append(pop.champion.to_sympy())
+
+
+# %%
+# We fix the seed for the objective function to make sure results are
+# comparable across individuals and, finally, we call the `evolve`
+# method to perform the evolutionary search.
+
+obj = functools.partial(objective, seed=population_params["seed"])
+
+cgp.evolve(
+    pop, obj, ea, **evolve_params, print_progress=True, callback=recording_callback,
+)
+
+# %%
+# After the evolutionary search has ended, we print the expression
+# with the highest fitness and plot the progression of the search.
+
+print(f"Final expression {pop.champion.to_sympy()[0]} with fitness {pop.champion.fitness}")
+
+print("Best performing expression per generation (for fitness increase > 0.5):")
+old_fitness = -np.inf
+for i, (fitness, expr) in enumerate(zip(history["fitness_champion"], history["expr_champion"])):
+    delta_fitness = fitness - old_fitness
+    if delta_fitness > 0.5:
+        print(f"{i:3d}: {fitness}, {expr}")
+        old_fitness = fitness
+print(f"{i:3d}: {fitness}, {expr}")
+
+width = 9.0
+
+fig = plt.figure(figsize=(width, width / scipy.constants.golden))
+
+ax_fitness = fig.add_subplot(111)
+ax_fitness.set_xlabel("Generation")
+ax_fitness.set_ylabel("Fitness")
+ax_fitness.set_yscale("symlog")
+
+ax_fitness.axhline(0.0, color="k")
+ax_fitness.plot(history["fitness_champion"], lw=2)
+
+plt.savefig("example_parametrized_nodes.pdf", dpi=300)

test/test_node.py

@@ -405,6 +405,42 @@ class CustomParameter(cgp.Parameter):
     assert y != pytest.approx(1.0)
 
 
+def test_multiple_parameters_per_node():
+
+    p = 3.1415
+    q = 2.7128
+
+    class DoubleParameter(cgp.OperatorNode):
+        _arity = 0
+        _initial_values = {"<p>": lambda: p, "<q>": lambda: q}
+        _def_output = "<p> + <q>"
+        _def_numpy_output = "np.ones(x.shape[0]) * (<p> + <q>)"
+        _def_torch_output = "torch.ones(1).expand(x.shape[0]) * (<p> + <q>)"
+
+    genome_params = {
+        "n_inputs": 1,
+        "n_outputs": 1,
+        "n_columns": 1,
+        "n_rows": 1,
+        "levels_back": None,
+    }
+    primitives = (DoubleParameter,)
+    genome = cgp.Genome(**genome_params, primitives=primitives)
+    # f(x) = p + q
+    genome.dna = [
+        ID_INPUT_NODE,
+        ID_NON_CODING_GENE,
+        0,
+        0,
+        ID_OUTPUT_NODE,
+        1,
+    ]
+    f = cgp.CartesianGraph(genome).to_func()
+    y = f([0.0])[0]
+
+    assert y == pytest.approx(p + q)
+
+
 def test_raise_broken_def_output():
     with pytest.raises(SyntaxError):