Repository: EC-KitY/EC-KitY (language: python)

Issue:
---TYPED GP---
* Added TreeNode class with two types: FunctionNode and TerminalNode and type field.
* Added type checks to the relevant places (adding node, replacing subtree - crossover and mutation).
* Added range of tries to each random addition/replacement until it would match by type - in order to avoid infinite loops and waste of running time. The ranges are proportional to the items we rand from or to.
* Added range for creating a tree to avoid situations of building a tree which we aren't able to satisfy type wise.
* Added utils to get the function's parameter types and return value type (they need to be implicitly written in the function definition for this to work).
* Added TYPED variable for the examples for enable choosing running with types or without.
* Added the option to define types of functions and terminals at the algorithm setting.
* Updated the imports in relevant files to call the typed version: import grow_typed.py, full_typed.py and tree_typed_individual_node.py.


Discussion / hints:
Added requested functionality, unittests and validated all system tests pass. 

For the BNF parser, parameters and grammer conventions should be decided - doesn't block this pr.
The parser would create new functions and terminals according to recursive search in the received BNF in a text file and the user will be able to use them to create a tree.

Fix patch:
diff --git a/.github/workflows/python-package.yml b/.github/workflows/python-package.yml
index 948efc1c..3837f49d 100644
--- a/.github/workflows/python-package.yml
+++ b/.github/workflows/python-package.yml
@@ -5,14 +5,13 @@ name: Python package
 
 on:
   push:
-    branches: [ main, develop ]
+    branches: [main, develop]
   pull_request:
-    branches: [ main, develop ]
+    branches: [main, develop]
   workflow_dispatch:
 
 jobs:
   build:
-
     runs-on: ubuntu-latest
     strategy:
       fail-fast: false
@@ -20,23 +19,24 @@ jobs:
         python-version: ["3.8", "3.9", "3.10", "3.11"]
 
     steps:
-    - uses: actions/checkout@v4.1.5
-    - name: Set up Python ${{ matrix.python-version }}
-      uses: actions/setup-python@v5.1.0
-      with:
-        python-version: ${{ matrix.python-version }}
-    - name: Install dependencies
-      run: |
-        python -m pip install --upgrade pip
-        python -m pip install flake8 pytest matplotlib
-        python -m pip install scikit-learn>=1.1 --use-pep517
-        if [ -f requirements.txt ]; then pip install -r requirements.txt; fi
-    - name: Lint with flake8
-      run: |
-        # stop the build if there are Python syntax errors or undefined names
-        flake8 . --count --select=E9,F63,F7,F82 --show-source --statistics
-        # exit-zero treats all errors as warnings. The GitHub editor is 127 chars wide
-        flake8 . --count --exit-zero --max-complexity=10 --max-line-length=127 --statistics
-    - name: Test with pytest
-      run: |
-        pytest
+      - uses: actions/checkout@v4.1.5
+      - name: Set up Python ${{ matrix.python-version }}
+        uses: actions/setup-python@v5.1.0
+        with:
+          python-version: ${{ matrix.python-version }}
+      - name: Install dependencies
+        run: |
+          python -m pip install --upgrade pip
+          python -m pip install flake8 pytest coverage matplotlib
+          python -m pip install scikit-learn>=1.1 --use-pep517
+          if [ -f requirements.txt ]; then pip install -r requirements.txt; fi
+      - name: Lint with flake8
+        run: |
+          # stop the build if there are Python syntax errors or undefined names
+          flake8 . --count --select=E9,F63,F7,F82 --show-source --statistics
+          # exit-zero treats all errors as warnings. The GitHub editor is 127 chars wide
+          flake8 . --count --exit-zero --max-complexity=10 --max-line-length=127 --statistics
+      - name: Test with pytest and coverage
+        run: |
+          coverage run -m pytest
+          coverage report --fail-under=85
diff --git a/.github/workflows/run-examples.yml b/.github/workflows/run-examples.yml
index d7bd0df6..262029e5 100644
--- a/.github/workflows/run-examples.yml
+++ b/.github/workflows/run-examples.yml
@@ -35,7 +35,10 @@ jobs:
         find examples -type f -name "*.py" -exec sed -E -i 's/max_generation\s*=\s*[0-9]+/max_generation=3/g' {} \;
         find examples -type f -name "*.py" -exec sed -E -i 's/population_size\s*=\s*[0-9]+/population_size=6/g' {} \;
 
-    - name: Find and run Python files in parallel
+    - name: Find and run Python files sequentially
       run: |
-        NUM_PROCS=$(nproc)
-        find examples -type f -name "*.py" | parallel --halt now,fail=1 -j $NUM_PROCS 'echo "Running {}"; python {}'
+        find examples -type f -name "*.py" | while read -r file; do
+        echo "Running $file"
+        cat "$file"
+        python "$file" || exit 1
+        done
diff --git a/README.md b/README.md
index ff022bcf..f02bcb40 100644
--- a/README.md
+++ b/README.md
@@ -44,11 +44,23 @@ You can run an EA with just 3 lines of code. The problem being solved herein is
 
 Additional information on this problem can be found in the [Symbolic Regression Tutorial](https://github.com/EC-KitY/EC-KitY/wiki/Tutorial:-Symbolic-Regression).
 ```python
-from eckity.algorithms.simple_evolution import SimpleEvolution
-from eckity.subpopulation import Subpopulation
-from examples.treegp.non_sklearn_mode.symbolic_regression.sym_reg_evaluator import SymbolicRegressionEvaluator
-
-algo = SimpleEvolution(Subpopulation(SymbolicRegressionEvaluator()))
+from eckity import Subpopulation
+from eckity.algorithms import SimpleEvolution
+from eckity.base.untyped_functions import f_add, f_sub, f_mul, f_div
+from eckity.creators import FullCreator
+from eckity.genetic_operators import SubtreeCrossover, SubtreeMutation
+from examples.treegp.basic_mode.symbolic_regression import SymbolicRegressionEvaluator
+
+algo = SimpleEvolution(
+    Subpopulation(
+        SymbolicRegressionEvaluator(),
+        creator=FullCreator(
+            terminal_set=['x', 'y', 'z'],
+            function_set=[f_add, f_sub, f_mul, f_div]
+        ),
+        operators_sequence=[SubtreeCrossover(), SubtreeMutation()]
+    )
+)
 algo.evolve()
 print(f'algo.execute(x=2,y=3,z=4): {algo.execute(x=2, y=3, z=4)}')
 ```
diff --git a/eckity/__init__.py b/eckity/__init__.py
index d351cc5d..5a44cafb 100644
--- a/eckity/__init__.py
+++ b/eckity/__init__.py
@@ -2,7 +2,7 @@
 EC-KitY: Evolutionary Computation Tool Kit in Python
 """
 
-__version__ = "0.3.4"
+__version__ = "0.4.0"
 
 from .before_after_publisher import BeforeAfterPublisher
 from .individual import Individual
diff --git a/eckity/algorithms/algorithm.py b/eckity/algorithms/algorithm.py
index 4cdfdc0e..1bd3361f 100644
--- a/eckity/algorithms/algorithm.py
+++ b/eckity/algorithms/algorithm.py
@@ -213,7 +213,7 @@ def initialize(self) -> None:
         Initialize seed, Executor and relevant operators
         """
         self.set_random_seed(self.random_seed)
-        logger.info("random seed = %f", self.random_seed)
+        logger.info("random seed = %d", self.random_seed)
         self.population_evaluator.set_executor(self.executor)
 
         for field in self.__dict__.values():
diff --git a/eckity/algorithms/simple_evolution.py b/eckity/algorithms/simple_evolution.py
index 409ce140..fd509494 100644
--- a/eckity/algorithms/simple_evolution.py
+++ b/eckity/algorithms/simple_evolution.py
@@ -85,7 +85,7 @@ def __init__(
         self,
         population,
         statistics=None,
-        breeder: SimpleBreeder = None,
+        breeder: SimpleBreeder = SimpleBreeder(),
         population_evaluator: SimplePopulationEvaluator = None,
         max_generation=500,
         events=None,
diff --git a/eckity/base/__init__.py b/eckity/base/__init__.py
index 244b0d3b..e69de29b 100644
--- a/eckity/base/__init__.py
+++ b/eckity/base/__init__.py
@@ -1,1 +0,0 @@
-# turns current folder into module
diff --git a/eckity/base/typed_functions.py b/eckity/base/typed_functions.py
new file mode 100644
index 00000000..3ce0c5d6
--- /dev/null
+++ b/eckity/base/typed_functions.py
@@ -0,0 +1,199 @@
+"""
+This module implements functions used in the function nodes of a GP tree.
+Note: all functions work on numpy arrays.
+"""
+
+from typing import Any
+import numpy as np
+from .types import t_argmax
+
+
+def add2floats(x: float, y: float) -> float:
+    """x+y"""
+    return np.add(x, y)
+
+
+def sub2floats(x: float, y: float) -> float:
+    """x-y"""
+    return np.subtract(x, y)
+
+
+def mul2floats(x: float, y: float) -> float:
+    """x*y"""
+    return np.multiply(x, y)
+
+
+def div2floats(x: float, y: float) -> float:
+    """protected division: if abs(y) > 0.001 return x/y else return 0"""
+    with np.errstate(divide="ignore", invalid="ignore"):
+        return np.where(np.abs(y) > 0.001, np.divide(x, y), 0.0)
+
+
+def sqrt_float(x: float) -> float:
+    """protected square root: sqrt(abs(x))"""
+    return np.sqrt(np.absolute(x))
+
+
+def log_float(x: float) -> float:
+    """protected log: if abs(x) > 0.001 return log(abs(x)) else return 0"""
+    with np.errstate(divide="ignore", invalid="ignore"):
+        return np.where(np.abs(x) > 0.001, np.log(np.abs(x)), 0.0)
+
+
+def abs_float(x: float) -> float:
+    """absolute value of x"""
+    return np.absolute(x)
+
+
+def neg_float(x: float) -> float:
+    """negative of x"""
+    return np.negative(x)
+
+
+def inv_float(x: float) -> float:
+    """protected inverse: if abs(x) > 0.001 return 1/x else return 0"""
+    return div2floats(1, x)
+
+
+def max2floats(x: float, y: float) -> float:
+    """maximum(x,y)"""
+    return np.maximum(x, y)
+
+
+def min2floats(x: float, y: float) -> float:
+    """minimum(x,y)"""
+    return np.minimum(x, y)
+
+
+def sin_float(x: float) -> float:
+    """sin(x)"""
+    return np.sin(x)
+
+
+def cos_float(x: float) -> float:
+    """cos(x)"""
+    return np.cos(x)
+
+
+def tan_float(x: float) -> float:
+    """tan(x)"""
+    return np.tan(x)
+
+
+def iflte0_floats(x: float, y: float, z: float) -> float:
+    """if x <= 0 return y else return z"""
+    return np.where(x <= 0, y, z)
+
+
+def ifgt0_floats(x: float, y: float, z: float) -> float:
+    """if x > 0 return y else return z"""
+    return np.where(x > 0, y, z)
+
+
+def iflte_floats(x: float, y: float, z: float, w: float) -> float:
+    """if x <= y return z else return w"""
+    return np.where(x <= y, z, w)
+
+
+def ifgt_floats(x: float, y: float, z: float, w: float) -> float:
+    """if x > y return z else return w"""
+    return np.where(x > y, z, w)
+
+
+def and2floats(x: float, y: float) -> float:
+    """x and y"""
+    return np.logical_and(x, y)
+
+
+def or2floats(x: float, y: float) -> float:
+    """x or y"""
+    return np.logical_or(x, y)
+
+
+def not2floats(x: float) -> float:
+    """not x"""
+    return np.logical_not(x).astype(float)
+
+
+def if_then_else(test: bool, dit: Any, dif: Any) -> float:
+    """if test return dit else return dif"""
+    return np.where(test, dit, dif)
+
+
+def argmax2floats(x0: float, x1: float) -> t_argmax:
+    return np.argmax([x0, x1], axis=0)
+
+
+def and2ints(x: int, y: int) -> int:
+    """x and y"""
+    return np.logical_and(x, y)
+
+
+def or2ints(x: int, y: int) -> int:
+    """x or y"""
+    return np.logical_or(x, y)
+
+
+def not2ints(x: int) -> int:
+    """not x"""
+    return np.logical_not(x).astype(int)
+
+
+def if_then_else3ints(test: int, dit: int, dif: int) -> int:
+    """if test return dit else return dif"""
+    return np.where(test, dit, dif)
+
+
+def and2bools(x: bool, y: bool) -> bool:
+    """x and y"""
+    return np.bitwise_and(x, y)
+
+
+def or2bools(x: bool, y: bool) -> bool:
+    """x or y"""
+    return np.bitwise_or(x, y)
+
+
+def not2bools(x: bool) -> bool:
+    """not x"""
+    return np.logical_not(x)
+
+
+def if_then_else3bools(test: bool, dit: bool, dif: bool) -> bool:
+    """if test return dit else return dif"""
+    return np.where(test, dit, dif)
+
+
+__all__ = [
+    "add2floats",
+    "sub2floats",
+    "mul2floats",
+    "div2floats",
+    "sqrt_float",
+    "log_float",
+    "abs_float",
+    "neg_float",
+    "inv_float",
+    "max2floats",
+    "min2floats",
+    "sin_float",
+    "cos_float",
+    "tan_float",
+    "iflte0_floats",
+    "ifgt0_floats",
+    "iflte_floats",
+    "ifgt_floats",
+    "and2floats",
+    "or2floats",
+    "not2floats",
+    "if_then_else",
+    "argmax2floats",
+    "and2ints",
+    "or2ints",
+    "not2ints",
+    "if_then_else3ints",
+    "and2bools",
+    "or2bools",
+    "not2bools",
+    "if_then_else3bools",
+]
diff --git a/eckity/base/types.py b/eckity/base/types.py
new file mode 100644
index 00000000..7bcd0168
--- /dev/null
+++ b/eckity/base/types.py
@@ -0,0 +1,1 @@
+t_argmax = type("argmax", (int,), {})
diff --git a/eckity/genetic_encodings/gp/tree/functions.py b/eckity/base/untyped_functions.py
similarity index 79%
rename from eckity/genetic_encodings/gp/tree/functions.py
rename to eckity/base/untyped_functions.py
index 1bb6879d..c94c5d92 100644
--- a/eckity/genetic_encodings/gp/tree/functions.py
+++ b/eckity/base/untyped_functions.py
@@ -23,8 +23,8 @@ def f_mul(x, y):
 
 def f_div(x, y):
     """protected division: if abs(y) > 0.001 return x/y else return 0"""
-    with np.errstate(divide='ignore', invalid='ignore'):
-        return np.where(np.abs(y) > 0.001, np.divide(x, y), 0.)
+    with np.errstate(divide="ignore", invalid="ignore"):
+        return np.where(np.abs(y) > 0.001, np.divide(x, y), 0.0)
 
 
 def f_sqrt(x):
@@ -34,8 +34,8 @@ def f_sqrt(x):
 
 def f_log(x):
     """protected log: if abs(x) > 0.001 return log(abs(x)) else return 0"""
-    with np.errstate(divide='ignore', invalid='ignore'):
-        return np.where(np.abs(x) > 0.001, np.log(np.abs(x)), 0.)
+    with np.errstate(divide="ignore", invalid="ignore"):
+        return np.where(np.abs(x) > 0.001, np.log(np.abs(x)), 0.0)
 
 
 def f_abs(x):
@@ -50,8 +50,6 @@ def f_neg(x):
 
 def f_inv(x):
     """protected inverse: if abs(x) > 0.001 return 1/x else return 0"""
-    # with np.errstate(divide='ignore', invalid='ignore'):
-    #     return np.where(np.abs(x) > 0.001, 1. / x, 0.)
     return f_div(1, x)
 
 
@@ -120,6 +118,27 @@ def f_if_then_else(test, dit, dif):
     return np.where(test, dit, dif)
 
 
-
-full_function_set = [f_add, f_sub, f_mul, f_div, f_sqrt, f_log, f_abs, f_neg, f_inv, f_max, f_min, f_sin, f_cos, f_tan,
-                     f_iflte0, f_ifgt0, f_iflte, f_ifgt, f_and, f_or, f_not, f_if_then_else]
+__all__ = [
+    "f_add",
+    "f_sub",
+    "f_mul",
+    "f_div",
+    "f_sqrt",
+    "f_log",
+    "f_abs",
+    "f_neg",
+    "f_inv",
+    "f_max",
+    "f_min",
+    "f_sin",
+    "f_cos",
+    "f_tan",
+    "f_iflte0",
+    "f_ifgt0",
+    "f_iflte",
+    "f_ifgt",
+    "f_and",
+    "f_or",
+    "f_not",
+    "f_if_then_else",
+]
diff --git a/eckity/base/utils.py b/eckity/base/utils.py
index 2a4facf9..624f52ab 100644
--- a/eckity/base/utils.py
+++ b/eckity/base/utils.py
@@ -3,18 +3,19 @@
 """
 
 from inspect import getfullargspec
+from typing import Callable
 
 
-def arity(func):
-	"""
-	Parameters
-	----------
-	func : function
-		A function.
+def arity(func: Callable) -> int:
+    """
+    Parameters
+    ----------
+    func : function
+            A function.
 
-	Returns
-	-------
-	arity : int
-		The function's arity.
-	"""
-	return len(getfullargspec(func)[0])
+    Returns
+    -------
+    arity : int
+            The function's arity.
+    """
+    return len(getfullargspec(func)[0])
diff --git a/eckity/creators/__init__.py b/eckity/creators/__init__.py
index cb87c921..4ff69fb1 100644
--- a/eckity/creators/__init__.py
+++ b/eckity/creators/__init__.py
@@ -10,4 +10,4 @@
 from .gp_creators.tree_creator import GPTreeCreator
 from .gp_creators.grow import GrowCreator
 from .gp_creators.full import FullCreator
-from .gp_creators.ramped_hh import RampedHalfAndHalfCreator
+from .gp_creators.half import HalfCreator
diff --git a/eckity/creators/creator.py b/eckity/creators/creator.py
index 184dd133..d30d02fd 100644
--- a/eckity/creators/creator.py
+++ b/eckity/creators/creator.py
@@ -1,5 +1,7 @@
 from abc import abstractmethod
+from typing import Any, Dict, List
 
+from eckity import Individual
 from eckity.event_based_operator import Operator
 from eckity.fitness.simple_fitness import SimpleFitness
 
@@ -11,14 +13,17 @@ def __init__(self, events=None, fitness_type=SimpleFitness):
         self.fitness_type = fitness_type
 
     @abstractmethod
-    def create_individuals(self, n_individuals, higher_is_better):
+    def create_individuals(
+        self, n_individuals: int, higher_is_better: bool
+    ) -> List[Individual]:
         pass
 
-    def apply_operator(self, payload):
+    def apply_operator(self, payload: Any) -> Any:
         return self.create_individuals(payload)
 
-    def event_name_to_data(self, event_name):
-        if event_name == "after_operator":
-            return {"created_individuals": self.created_individuals}
-        else:
-            return {}
+    def event_name_to_data(self, event_name: str) -> Dict[str, Any]:
+        return (
+            {"created_individuals": self.created_individuals}
+            if event_name == "after_operator"
+            else {}
+        )
diff --git a/eckity/creators/gp_creators/full.py b/eckity/creators/gp_creators/full.py
index 6890338a..ff600664 100644
--- a/eckity/creators/gp_creators/full.py
+++ b/eckity/creators/gp_creators/full.py
@@ -1,18 +1,26 @@
+from typing import Any, Callable, Dict, List, Optional, Tuple, Union
+
 from overrides import overrides
 
 from eckity.creators.gp_creators.tree_creator import GPTreeCreator
+from eckity.genetic_encodings.gp import (
+    TerminalNode,
+    FunctionNode,
+)
+from eckity.genetic_encodings.gp.tree.utils import get_func_types
 
 
 class FullCreator(GPTreeCreator):
     def __init__(
         self,
-        init_depth=None,
-        function_set=None,
-        terminal_set=None,
-        erc_range=None,
-        bloat_weight=0.1,
-        events=None,
-        update_parents=False,
+        init_depth: Tuple[int, int] = None,
+        function_set: List[Callable] = None,
+        terminal_set: Union[Dict[Any, type], List[Any]] = None,
+        erc_range: Union[Tuple[int, int], Tuple[float, float]] = None,
+        bloat_weight: float = 0.0,
+        events: List[str] = None,
+        root_type: Optional[type] = None,
+        update_parents: bool = False,
     ):
         """
         Tree creator using the full method
@@ -28,9 +36,6 @@ def __init__(
         terminal_set : list
                 List of terminals used in the GP-tree leaves. The default is None.
 
-        erc_range : (float, float)
-                Range of values for ephemeral random constant (ERC). The default is None.
-
         bloat_weight : float
                 Bloat control weight to punish large trees. Bigger values make a bigger punish.
 
@@ -41,58 +46,57 @@ def __init__(
             init_depth=init_depth,
             function_set=function_set,
             terminal_set=terminal_set,
-            erc_range=erc_range,
             bloat_weight=bloat_weight,
+            erc_range=erc_range,
             events=events,
+            root_type=root_type,
             update_parents=update_parents,
         )
 
     @overrides
-    def create_tree(self, tree_ind, max_depth=5):
-        """
-        Create a random tree using the full method, and assign it to the given individual.
-
-        Parameters
-        ----------
-        tree_ind: Tree
-                An individual of GP Tree representation with an initially empty tree
-
-        max_depth: int
-                Maximum depth of tree. The default is 5.
-
-        Returns
-        -------
-        None.
-        """
-        self._create_tree(tree_ind, max_depth, 0)
-
-    def _create_tree(self, tree_ind, max_depth=5, depth=0):
+    def create_tree(
+        self,
+        tree,
+        random_function: Callable[[type], Optional[FunctionNode]],
+        random_terminal: Callable[[type], Optional[TerminalNode]],
+        depth: int = 0,
+        node_type: Optional[type] = None,
+    ) -> None:
         """
-        Recursively create a random tree using the full method.
+        Recursively create a random tree using the full method
 
         Parameters
         ----------
-        max_depth : int
-                Maximum depth of tree. The default is 5.
-
-        depth : int, optional
-                Current depth in recursive process. The default is 0.
+        depth: int
+                Current depth in recursive process.
 
         Returns
         -------
         None.
 
         """
-        is_func = False
+        max_depth = self.init_depth[1]
 
         if depth >= max_depth:
-            node = tree_ind.random_terminal()
-        else:
-            node = tree_ind.random_function()
-            is_func = True
+            node = random_terminal(node_type)
+            self._assert_node_created(node)
 
-        tree_ind.add_tree(node)
-
-        if is_func:
-            for i in range(tree_ind.arity[node]):
-                self._create_tree(tree_ind, max_depth, depth=depth + 1)
+            # add the new node to the tree of the given individual
+            tree.append(node)
+        else:
+            node = random_function(node_type)
+            self._assert_node_created(node)
+
+            # add the new node to the tree of the given individual
+            tree.append(node)
+
+            # recursively add argument nodes to the tree
+            func_types = get_func_types(node.function)[:-1]
+            for t in func_types:
+                self.create_tree(
+                    tree,
+                    random_function,
+                    random_terminal,
+                    depth=depth + 1,
+                    node_type=t,
+                )
diff --git a/eckity/creators/gp_creators/grow.py b/eckity/creators/gp_creators/grow.py
index 67512f80..65b71cbc 100644
--- a/eckity/creators/gp_creators/grow.py
+++ b/eckity/creators/gp_creators/grow.py
@@ -1,104 +1,118 @@
-from overrides import overrides
 import random
+from typing import Any, Callable, Dict, List, Optional, Tuple, Union
+
+from overrides import overrides
 
 from eckity.creators.gp_creators.tree_creator import GPTreeCreator
+from eckity.genetic_encodings.gp import (
+    FunctionNode,
+    TerminalNode,
+    TreeNode,
+)
 
 
 class GrowCreator(GPTreeCreator):
-	def __init__(self,
-				 init_depth=None,
-				 function_set=None,
-				 terminal_set=None,
-				 erc_range=None,
-				 bloat_weight=0.1,
-				 events=None,
-				 update_parents=False):
-		"""
-		Tree creator using the grow method
-
-		Parameters
-		----------
-		init_depth : (int, int)
-		Min and max depths of initial random trees. The default is None.
-
-		function_set : list
-			List of functions used as internal nodes in the GP-tree. The default is None.
-
-		terminal_set : list
-			List of terminals used in the GP-tree leaves. The default is None.
-
-		erc_range : (float, float)
-			Range of values for ephemeral random constant (ERC). The default is None.
-
-		bloat_weight : float
-			Bloat control weight to punish large trees. Bigger values make a bigger punish.
-
-		events : list
-			List of events related to this class
-		"""
-		super().__init__(
-			init_depth=init_depth,
-			function_set=function_set,
-			terminal_set=terminal_set,
-			erc_range=erc_range,
-			bloat_weight=bloat_weight,
-			events=events,
-			update_parents=update_parents
-		)
-
-	@overrides
-	def create_tree(self, tree_ind, max_depth=5):
-		"""
-		Create a random tree using the grow method, and assign it to the given individual.
-
-		Parameters
-		----------
-		tree_ind: Tree
-			An individual of GP Tree representation with an initially empty tree
-
-		max_depth: int
-			Maximum depth of tree. The default is 5.
-
-		Returns
-		-------
-		None.
-		"""
-		self._create_tree(tree_ind, max_depth, 0)
-
-	def _create_tree(self, tree_ind, max_depth=5, depth=0):
-		"""
-		Recursively create a random tree using the grow method
-
-		Parameters
-		----------
-		max_depth : int
-			Maximum depth of tree. The default is 2.
-		depth : int, optional
-			Current depth in recursive process. The default is 0.
-
-		Returns
-		-------
-		None.
-
-		"""
-		is_func = False
-
-		if depth < self.init_depth[0]:
-			node = tree_ind.random_function()
-			is_func = True
-		elif depth >= max_depth:
-			node = tree_ind.random_terminal()
-		else:  # intermediate depth, grow
-			if random.random() > 0.5:
-				node = tree_ind.random_function()
-				is_func = True
-			else:
-				node = tree_ind.random_terminal()
-
-		# add the new node to the tree of the given individual
-		tree_ind.add_tree(node)
-
-		if is_func:
-			# recursively add arguments to the function node, according to its arity
-			for i in range(tree_ind.arity[node]):
-				self._create_tree(tree_ind, max_depth, depth=depth + 1)
+    def __init__(
+        self,
+        init_depth: Tuple[int, int] = None,
+        function_set: List[Callable] = None,
+        terminal_set: Union[Dict[Any, type], List[Any]] = None,
+        bloat_weight: float = 0.0,
+        erc_range: Union[Tuple[int, int], Tuple[float, float]] = None,
+        p_prune: float = 0.5,
+        events: List[str] = None,
+        root_type: Optional[type] = None,
+        update_parents: bool = False
+    ):
+        """
+        Tree creator using the grow method
+
+        Parameters
+        ----------
+        init_depth : (int, int)
+        Min and max depths of initial random trees. The default is None.
+
+        function_set : list, default=None
+                List of functions used as internal nodes in the GP-tree.
+
+        terminal_set : list, default=None.
+                List of terminals used in the GP-tree leaves.
+
+        bloat_weight : float, default=0.0
+                Bloat control weight to punish large trees.
+
+        p_prune : float, default=0.5
+                Probability of pruning the tree at each level.
+
+        events : list
+                List of events related to this class
+        """
+        super().__init__(
+            init_depth=init_depth,
+            function_set=function_set,
+            terminal_set=terminal_set,
+            bloat_weight=bloat_weight,
+            events=events,
+            root_type=root_type,
+            update_parents=update_parents,
+            erc_range=erc_range,
+        )
+        self.p_prune = p_prune
+
+    @overrides
+    def create_tree(
+        self,
+        tree: List[TreeNode],
+        random_function: Callable[[type], Optional[FunctionNode]],
+        random_terminal: Callable[[type], Optional[TerminalNode]],
+        depth: int = 0,
+        node_type: Optional[type] = None,
+    ) -> None:
+        """
+        Generate a list of TreeNodes in-place.
+
+        Parameters
+        ----------
+        tree : List[TreeNode]
+            List of tree nodes representing a tree of a Tree individual.
+        random_function : Callable[[type], Optional[FunctionNode]]
+            Random FunctionNode generator.
+        random_terminal : Callable[[type], Optional[TerminalNode]]
+            Random TerminalNode generator.
+        depth : int, optional
+            current depth of the tree, by default 0
+        node_type : Optional[type], optional
+            _description_, by default None
+        """
+        if tree is None:
+            tree = []
+
+        min_depth, max_depth = self.init_depth
+
+        is_func = False
+        if depth < min_depth:
+            node = random_function(node_type)
+
+            is_func = True
+        elif depth >= max_depth:
+            node = random_terminal(node_type)
+        else:  # intermediate depth, grow
+            if random.random() < self.p_prune:
+                node = random_terminal(node_type)
+            else:
+                node = random_function(node_type)
+                is_func = True
+
+        # add the new node to the tree of the given individual
+        tree.append(node)
+
+        self._assert_node_created(node)
+
+        if is_func:
+            self._add_children(
+                tree=tree,
+                fn_node=node,
+                random_function=random_function,
+                random_terminal=random_terminal,
+                depth=depth,
+            )
diff --git a/eckity/creators/gp_creators/half.py b/eckity/creators/gp_creators/half.py
new file mode 100644
index 00000000..f24ed9a5
--- /dev/null
+++ b/eckity/creators/gp_creators/half.py
@@ -0,0 +1,141 @@
+from typing import Any, Callable, Dict, List, Optional, Tuple, Union
+
+from overrides import override
+
+from eckity.creators.gp_creators.full import FullCreator
+from eckity.creators.gp_creators.grow import GrowCreator
+from eckity.creators.gp_creators.tree_creator import GPTreeCreator
+from eckity.genetic_encodings.gp.tree.tree_individual import Tree
+
+
+class HalfCreator(GPTreeCreator):
+    def __init__(
+        self,
+        grow_creator: GrowCreator = None,
+        full_creator: FullCreator = None,
+        init_depth: Tuple[int, int] = None,
+        function_set: List[Callable] = None,
+        terminal_set: Union[Dict[Any, type], List[Any]] = None,
+        bloat_weight: float = 0.0,
+        erc_range: Union[Tuple[int, int], Tuple[float, float]] = None,
+        events: List[str] = None,
+        root_type: Optional[type] = None,
+    ):
+        """
+        Tree creator that creates trees using the Ramped Half and Half method
+
+        Parameters
+        ----------
+        grow_creator: GrowCreator
+                a tree creator that creates trees using the grow method
+
+        full_creator: FullCreator
+                a tree creator that creates trees using the full method
+
+        init_depth : (int, int)
+        Min and max depths of initial random trees. The default is None.
+
+        function_set : list
+                List of functions used as internal nodes in the GP tree. The default is None.
+
+        terminal_set : list or dict
+                List of terminals used in the GP-tree leaves. The default is None.
+
+        bloat_weight : float
+                Bloat control weight to punish large trees. Bigger values make a bigger punish.
+
+        events : list
+                List of events related to this class
+        """
+        super().__init__(
+            init_depth=init_depth,
+            function_set=function_set,
+            terminal_set=terminal_set,
+            bloat_weight=bloat_weight,
+            erc_range=erc_range,
+            events=events,
+            root_type=root_type,
+        )
+
+        # assign default creators
+        if grow_creator is None:
+            grow_creator = GrowCreator(
+                init_depth=self.init_depth,
+                function_set=self.function_set,
+                terminal_set=self.terminal_set,
+                events=self.events,
+                root_type=root_type,
+                erc_range=self.erc_range,
+            )
+        if full_creator is None:
+            full_creator = FullCreator(
+                init_depth=self.init_depth,
+                function_set=self.function_set,
+                terminal_set=self.terminal_set,
+                events=self.events,
+                root_type=root_type,
+                erc_range=self.erc_range,
+            )
+
+        self.grow_creator = grow_creator
+        self.full_creator = full_creator
+
+    @override
+    def create_individuals(
+        self, n_individuals: int, higher_is_better: bool
+    ) -> List[Tree]:
+        """
+        Initialize the subpopulation individuals using ramped half-and-half method.
+
+        Parameters
+        ----------
+        n_individuals: int
+                number of individuals to create
+
+        higher_is_better: bool
+                determines if the fitness of the created individuals should be minimized or maximized
+
+        Returns
+        -------
+
+        """
+        min_depth, max_depth = self.init_depth
+
+        # if pop size is 100 and we want depths 2,3,4,5,6 then group_size is 10:
+        # 10 'grow' with depth 2, 10 'full' with depth 2, 10 'grow' with depth 3, 10 'full' with depth 3, etc.
+        group_size = int(n_individuals / (2 * (max_depth + 1 - min_depth)))
+
+        if group_size == 0:
+            raise ValueError(
+                "Incompatible population size and init_depth. "
+                "Population size must be "
+                f"at least {2 * (max_depth + 1 - min_depth)}."
+            )
+
+        individuals = []
+
+        for depth in range(min_depth, max_depth + 1):
+            # first create `group_size` individuals using grow method
+            self.grow_creator.init_depth = (min_depth, depth)
+            grown_inds = self.grow_creator.create_individuals(
+                group_size, higher_is_better
+            )
+
+            # then create `group_size` individuals using full method
+            self.full_creator.init_depth = (min_depth, depth)
+            full_inds = self.full_creator.create_individuals(
+                group_size, higher_is_better
+            )
+
+            individuals.extend(grown_inds + full_inds)
+
+        # might need to add a few because 'group_size' may have been a float that was truncated
+        n_missing = n_individuals - len(individuals)
+        if n_missing > 0:
+            missing_inds = self.full_creator.create_individuals(
+                n_missing, higher_is_better
+            )
+            individuals.extend(missing_inds)
+
+        self.created_individuals = individuals
+        return individuals
diff --git a/eckity/creators/gp_creators/ramped_hh.py b/eckity/creators/gp_creators/ramped_hh.py
deleted file mode 100644
index 0937d439..00000000
--- a/eckity/creators/gp_creators/ramped_hh.py
+++ /dev/null
@@ -1,122 +0,0 @@
-from eckity.creators.gp_creators.full import FullCreator
-from eckity.creators.gp_creators.grow import GrowCreator
-from eckity.creators.gp_creators.tree_creator import GPTreeCreator
-from eckity.genetic_encodings.gp.tree.tree_individual import Tree
-from eckity.fitness.gp_fitness import GPFitness
-
-
-class RampedHalfAndHalfCreator(GPTreeCreator):
-	def __init__(self,
-				 grow_creator=None,
-				 full_creator=None,
-				 init_depth=None,
-				 function_set=None,
-				 terminal_set=None,
-				 erc_range=None,
-				 bloat_weight=0.1,
-				 events=None,
-				 update_parents=False):
-		"""
-		Tree creator that creates trees using the Ramped Half and Half method
-
-		Parameters
-		----------
-		grow_creator: GrowCreator
-			a tree creator that creates trees using the grow method
-
-		full_creator: FullCreator
-			a tree creator that creates trees using the full method
-
-		init_depth : (int, int)
-		Min and max depths of initial random trees. The default is None.
-
-		function_set : list
-			List of functions used as internal nodes in the GP tree. The default is None.
-
-		terminal_set : list
-			List of terminals used in the GP-tree leaves. The default is None.
-
-		erc_range : (float, float)
-			Range of values for ephemeral random constant (ERC). The default is None.
-
-		bloat_weight : float
-			Bloat control weight to punish large trees. Bigger values make a bigger punish.
-
-		events : list
-			List of events related to this class
-		"""
-		super().__init__(init_depth=init_depth,
-						 function_set=function_set,
-						 terminal_set=terminal_set,
-						 erc_range=erc_range,
-						 bloat_weight=bloat_weight,
-						 events=events,
-						 update_parents=update_parents)
-
-		# assign default creators
-		if grow_creator is None:
-			grow_creator = GrowCreator(init_depth=self.init_depth, function_set=self.function_set,
-									   terminal_set=self.terminal_set, erc_range=self.erc_range, events=self.events)
-		if full_creator is None:
-			full_creator = FullCreator(init_depth=self.init_depth, function_set=self.function_set,
-									   terminal_set=self.terminal_set, erc_range=self.erc_range, events=self.events)
-
-		self.grow_creator, self.full_creator = grow_creator, full_creator
-
-		self.init_method = None  # current creator in use (either grow or full)
-
-	def create_individuals(self, n_individuals, higher_is_better):
-		"""
-		Initialize the subpopulation individuals using ramped half-and-half method.
-
-		Parameters
-		----------
-		n_individuals: int
-			number of individuals to create
-
-		higher_is_better: bool
-			determines if the fitness of the created individuals should be minimized or maximized
-
-		Returns
-		-------
-
-		"""
-
-		min_depth, max_depth = self.init_depth[0], self.init_depth[1]
-
-		# if pop size is 100 and we want depths 2,3,4,5,6 then group_size is 10:
-		# 10 'grow' with depth 2, 10 'full' with depth 2, 10 'grow' with depth 3, 10 'full' with depth 3, etc.
-		group_size = int(n_individuals / (max_depth + 1 - min_depth) / 2)
-
-		individuals = []
-
-		for depth in range(min_depth, max_depth + 1):
-			for i in range(group_size):
-				# as explained above, first create (group_size) individuals using grow method
-				self.init_method = self.grow_creator
-				self._create_individuals(individuals, depth, higher_is_better)
-
-				# then create (group_size) individuals using full method
-				self.init_method = self.full_creator
-				self._create_individuals(individuals, depth, higher_is_better)
-
-		# might need to add a few because 'group_size' may have been a float that was truncated
-		self.init_method = self.full_creator
-		for i in range(n_individuals - len(individuals)):
-			self._create_individuals(individuals, max_depth, higher_is_better)
-
-		# TODO we don't need this event since creators have before/after operator events
-		# self.publish("after_creation")
-
-		self.created_individuals = individuals
-		return individuals
-
-	def _create_individuals(self, individuals, max_depth, higher_is_better):
-		t = Tree(init_depth=self.init_depth, function_set=self.function_set,
-				 terminal_set=self.terminal_set, erc_range=self.erc_range,
-				 fitness=GPFitness(bloat_weight=self.bloat_weight, higher_is_better=higher_is_better))
-		self.create_tree(t, max_depth=max_depth)
-		individuals.append(t)
-
-	def create_tree(self, tree_ind, max_depth):
-		self.init_method.create_tree(tree_ind, max_depth)
diff --git a/eckity/creators/gp_creators/tree_creator.py b/eckity/creators/gp_creators/tree_creator.py
index 5e160344..c2acdcef 100644
--- a/eckity/creators/gp_creators/tree_creator.py
+++ b/eckity/creators/gp_creators/tree_creator.py
@@ -1,23 +1,32 @@
+from typing import Any, Callable, Dict, List, Tuple, Union, Optional
+
+from overrides import override
+
+from eckity.base.utils import arity
 from eckity.creators.creator import Creator
-from eckity.fitness.simple_fitness import SimpleFitness
-from eckity.genetic_encodings.gp.tree.functions import f_add, f_sub, f_mul, f_div
 from eckity.fitness.gp_fitness import GPFitness
-from eckity.genetic_encodings.gp.tree.tree_individual import Tree
-
-from abc import abstractmethod
+from eckity.fitness.simple_fitness import SimpleFitness
+from eckity.genetic_encodings.gp import (
+    FunctionNode,
+    TerminalNode,
+    Tree,
+    TreeNode,
+)
+from eckity.genetic_encodings.gp.tree.utils import get_func_types
 
 
 class GPTreeCreator(Creator):
     def __init__(
         self,
-        init_depth=None,
-        function_set=None,
-        terminal_set=None,
-        erc_range=None,
-        fitness_type=SimpleFitness,
-        bloat_weight=0.1,
-        events=None,
-        update_parents=False,
+        init_depth: Tuple[int, int] = None,
+        function_set: List[Callable] = None,
+        terminal_set: Union[Dict[Any, type], List[Any]] = None,
+        erc_range: Union[Tuple[int, int], Tuple[float, float]] = None,
+        fitness_type: type = SimpleFitness,
+        bloat_weight: float = 0.0,
+        events: List[str] = None,
+        root_type: Optional[type] = None,
+        update_parents: bool = False
     ):
         if events is None:
             events = ["after_creation"]
@@ -27,38 +36,138 @@ def __init__(
             init_depth = (2, 4)
 
         if function_set is None:
-            function_set = [f_add, f_sub, f_mul, f_div]
+            raise ValueError("function_set must be provided")
 
         if terminal_set is None:
-            terminal_set = ["x", "y", "z", 0, 1, -1]
+            raise ValueError("terminal_set must be provided")
 
         self.init_depth = init_depth
         self.function_set = function_set
         self.terminal_set = terminal_set
-        self.erc_range = erc_range
         self.bloat_weight = bloat_weight
+        self.root_type = root_type
         self.update_parents = update_parents
+        self.erc_range = erc_range
 
-    def create_individuals(self, n_individuals, higher_is_better):
+    @override
+    def create_individuals(
+        self, n_individuals: int, higher_is_better: bool
+    ) -> List[Tree]:
         individuals = [
             Tree(
                 function_set=self.function_set,
                 terminal_set=self.terminal_set,
-                erc_range=self.erc_range,
                 fitness=GPFitness(
                     bloat_weight=self.bloat_weight,
-                    higher_is_better=higher_is_better
+                    higher_is_better=higher_is_better,
                 ),
-                init_depth=self.init_depth,
+                root_type=self.root_type,
                 update_parents=self.update_parents,
+                erc_range=self.erc_range,
             )
             for _ in range(n_individuals)
         ]
         for ind in individuals:
-            self.create_tree(ind, self.init_depth[1])
+            self.create_tree(
+                ind.tree,
+                ind.random_function,
+                ind.random_terminal,
+                node_type=self.root_type,
+            )
         self.created_individuals = individuals
         return individuals
 
-    @abstractmethod
-    def create_tree(self, tree_ind, max_depth):
+    def create_tree(
+        self,
+        tree: List[TreeNode],
+        random_function: Callable[[type], Optional[FunctionNode]],
+        random_terminal: Callable[[type], Optional[TerminalNode]],
+        depth: int = 0,
+        node_type: Optional[type] = None,
+    ) -> None:
+        """
+        Recursively build the tree representation
+        of an existing Tree individual.
+        
+        (This method is not abstract as it is not required in HalfCreator)
+
+        Parameters
+        ----------
+        tree_ind : Tree
+            Individual to create the tree representation for
+        """
         pass
+
+    def _add_children(
+        self,
+        tree: List[TreeNode],
+        fn_node: FunctionNode,
+        random_function: Callable[[type], Optional[FunctionNode]],
+        random_terminal: Callable[[type], Optional[TerminalNode]],
+        depth: int,
+    ) -> None:
+        """
+        Recursively add children to a function node.
+
+        Parameters
+        ----------
+        tree : List[TreeNode]
+            Total representation of the tree.
+        fn_node : FunctionNode
+            Function node, parent of the children to be created.
+        random_function : Callable[[type], Optional[FunctionNode]]
+            Random FunctionNode generator.
+        random_terminal : Callable[[type], Optional[TerminalNode]]
+            Random TerminalNode generator.
+        depth : int
+            current depth of the tree
+        """
+        func_types = get_func_types(fn_node.function)
+        for i in range(arity(fn_node.function)):
+            self.create_tree(
+                tree,
+                random_function,
+                random_terminal,
+                depth=depth + 1,
+                node_type=func_types[i],
+            )
+
+    def _assert_node_created(
+        self,
+        node: Optional[TreeNode],
+        node_type: Optional[type] = None
+    ) -> None:
+        """
+        Assert that a TreeNode was created successfully.
+
+        Parameters
+        ----------
+        node : Optional[TreeNode]
+            Generated tree node.
+        node_type : Optional[type], optional
+            Generated tree node type, by default None
+
+        Raises
+        ------
+        ValueError
+            Raised if a node was not generated successfully.
+        """
+        if node is None:
+            # optionally add type info for typed case
+            type_info = (
+                f"with type {node_type}" if node_type is not None else ""
+            )
+
+            # optionally add erc_range info if it was defined
+            erc_info = (
+                f"erc range ({self.erc_range}) ,"
+                if self.erc_range is not None
+                else ""
+            )
+
+            raise ValueError(
+                f"Could not generate node {type_info} "
+                f"due to a mismatch in terminals ({self.terminal_set}), "
+                f"functions ({self.function_set}), {erc_info} "
+                f"or init_depth ({self.init_depth}) configuration."
+            )
diff --git a/eckity/event_based_operator.py b/eckity/event_based_operator.py
index 1a4a6d93..fed6f5c5 100644
--- a/eckity/event_based_operator.py
+++ b/eckity/event_based_operator.py
@@ -2,7 +2,8 @@
 This module implements the Operator class
 """
 
-from abc import abstractmethod, ABC
+from abc import ABC, abstractmethod
+from typing import Any
 
 from eckity.before_after_publisher import BeforeAfterPublisher
 
@@ -34,7 +35,24 @@ def act(self, payload=None):
         -------
         the return value of the operator implemented in the sub-class
         """
-        return self.act_and_publish_before_after(lambda: self.apply_operator(payload))
+        self._assert_arity(payload)
+        return self.act_and_publish_before_after(
+            lambda: self.apply_operator(payload)
+        )
+
+    def _assert_arity(self, payload: Any) -> None:
+        if not isinstance(payload, list):
+            if self.arity != 1:
+                raise ValueError(
+                    f"Received payload of type {type(payload)}, "
+                    f"expected arity of 1 but got {self.arity}."
+                )
+
+        elif (n_individuals := len(payload)) != self.arity:
+            raise ValueError(
+                f"Expected {self.arity} individuals,"
+                f"but received {n_individuals}."
+            )
 
     def get_operator_arity(self):
         """
diff --git a/eckity/fitness/gp_fitness.py b/eckity/fitness/gp_fitness.py
index 9687d6eb..cfc50e9a 100644
--- a/eckity/fitness/gp_fitness.py
+++ b/eckity/fitness/gp_fitness.py
@@ -65,6 +65,10 @@ def get_augmented_fitness(self, individual):
         """
         score = self.get_pure_fitness()
 
+        if self.bloat_weight == 0:
+            # no bloat
+            return score
+
         # subtract bloat value from the fitness score if it should be maximized,
         # otherwise add bloat value to fitness score
         return (
diff --git a/eckity/genetic_encodings/gp/__init__.py b/eckity/genetic_encodings/gp/__init__.py
index 40dd6f6e..8c7c485b 100644
--- a/eckity/genetic_encodings/gp/__init__.py
+++ b/eckity/genetic_encodings/gp/__init__.py
@@ -1,1 +1,2 @@
 from .tree.tree_individual import Tree
+from .tree.tree_node import TreeNode, FunctionNode, TerminalNode
diff --git a/eckity/genetic_encodings/gp/tree/tree_individual.py b/eckity/genetic_encodings/gp/tree/tree_individual.py
index b501acc6..58fa8d44 100644
--- a/eckity/genetic_encodings/gp/tree/tree_individual.py
+++ b/eckity/genetic_encodings/gp/tree/tree_individual.py
@@ -3,65 +3,100 @@
 """
 
 import logging
-import numpy as np
+import random
 from numbers import Number
-from random import randint, uniform, random
+from typing import Any, Callable, Dict, List, Optional, Tuple, Union
 
-from eckity.base.utils import arity
-from eckity.genetic_encodings.gp.tree.utils import _generate_args
+import numpy as np
 
+from eckity.fitness import Fitness, GPFitness
+from eckity.genetic_encodings.gp.tree.tree_node import (
+    FunctionNode,
+    TerminalNode,
+    TreeNode,
+)
 from eckity.individual import Individual
-from eckity.genetic_encodings.gp.tree.functions import f_add, f_sub, f_mul, f_div
+
+from .utils import generate_args, get_func_types, get_return_type
 
 logger = logging.getLogger(__name__)
 
+
 class Tree(Individual):
     """
     A tree optimized for genetic programming operations.
     It is represented by a list of nodes in depth-first order.
     There are two types of nodes: functions and terminals.
 
-    (tree is not meant as a stand-alone -- parameters are supplied through the call from the Tree Creators)
-
-    Parameters
-    ----------
-    init_depth : (int, int)
-        Min and max depths of initial random trees. The default is None.
-        
-    function_set : list
-        List of functions used as internal nodes in the GP tree. The default is None.
-        
-    terminal_set : list
-        List of terminals used in the GP-tree leaves. The default is None.
-        
-    erc_range : (float, float)
-        Range of values for ephemeral random constant (ERC). The default is None.
+    (tree is not meant as a stand-alone,
+    parameters are supplied through the call from the Tree Creators)
     """
-    def __init__(self,
-                 fitness,
-                 function_set=None,
-                 terminal_set=None,
-                 erc_range=None,
-                 init_depth=(1, 2),
-                 update_parents=False):
+
+    def __init__(
+        self,
+        fitness: Fitness = GPFitness(),
+        function_set: List[Callable] = None,
+        terminal_set: Union[Dict[Any, type], List[Any]] = None,
+        tree: List[TreeNode] = None,
+        erc_range: Optional[
+            Union[Tuple[float, float], Tuple[int, int]]
+        ] = None,
+        root_type: Optional[type] = None,
+        update_parents: bool = False
+    ):
+        """
+        GP Tree Individual.
+
+        Parameters
+        ----------
+        fitness : Fitness
+            Manages fitness state, by default GPFitness
+        function_set: List[Callable], default=None
+            List of functions used as internal nodes in the GP tree.
+        terminal_set : Union[Dict[Any, type], List[Any]], optional
+            Mapping of terminal nodes and their types.
+            In the untyped case, all types are None.
+            Lists are treated as untyped, and will be assigned None.
+        tree : List[TreeNode], optional
+            Actual tree representation, by default None
+        erc_range : tuple of float or int, optional
+            Range of Ephemeral random constant values, by default None
+        root_type: type, optional
+            Root node type, by default None
+
+        Raises
+        ------
+        ValueError
+            If typed function is used with untyped terminals.
+        """
         super().__init__(fitness, update_parents=update_parents)
-        if function_set is None:
-            function_set = [f_add, f_sub, f_mul, f_div]
 
-        if terminal_set is None:
-            terminal_set = ['x', 'y', 'z', 0, 1, -1]
+        self.erc_range = erc_range
+
+        function_set, terminal_set = self._handle_input_types(
+            function_set, terminal_set, root_type
+        )
 
         self.function_set = function_set
         self.terminal_set = terminal_set
-        self.n_terminals = len(terminal_set)
-        self.arity = dict([(func, arity(func)) for func in self.function_set])
-        self.vars = [t for t in terminal_set if not isinstance(t, Number)]
-        self.erc_range = erc_range
-        self.n_functions = len(self.function_set)
-        self.init_depth = init_depth
-        self.tree = []
 
-    def size(self):
+        # actual tree representation
+        if tree is None:
+            tree = []
+        self.tree = tree
+
+        # this is the type of the execution result of the program (tree)
+        self.root_type = root_type
+
+    @property
+    def root(self) -> TreeNode:
+        return self.tree[0]
+
+    @property
+    def erc_type(self) -> Optional[Union[int, float]]:
+        return type(self.erc_range[0]) if self.erc_range else None
+
+    def size(self) -> int:
         """
         Compute size of tree.
 
@@ -72,28 +107,43 @@ def size(self):
         """
         return len(self.tree)
 
-    def add_tree(self, node):
-        self.tree.append(node)
-
-    def empty_tree(self):
-        self.tree = []
+    def add_tree(self, node: TreeNode) -> None:
+        """Add a node to the tree following the defined type constrains"""
+        if (
+            self.size() == 0 and node.node_type == self.root_type
+        ) or self._should_add([0], node):
+            self.tree.append(node)
+        else:
+            raise ValueError(f"Could not add node {node} to tree {self.tree}")
 
-    def _depth(self, pos, depth):
-        """Recursively compute depth 
-           (pos is a size-1 list so as to pass "by reference" on successive recursive calls)."""
+    def _should_add(self, pos: int, new_node: TreeNode) -> bool:
+        """
+        Recursively find the parent function of the new node and check it
+        matches by type to the expected parameter. (pos is a size-1 list
+        so as to pass "by reference" on successive recursive calls).
+        """
+        if (
+            pos[0] == self.size()
+        ):  # reached the place the new node should be placed
+            return True
 
         node = self.tree[pos[0]]
-
-        depths = []
-        if node in self.function_set:
-            for i in range(self.arity[node]):
+        res = None
+        if isinstance(node, FunctionNode):
+            func_types = get_func_types(node.function)
+            for i in range(node.n_args):
                 pos[0] += 1
-                depths.append(1 + self._depth(pos, depth + 1))
-            return max(depths)
-        else:
-            return 0
+                res = self._should_add(pos, node)
+                if res:
+                    return func_types[i] == new_node.node_type
+                elif res is not None:
+                    return res
+        return res
 
-    def depth(self):
+    def empty_tree(self) -> None:
+        self.tree = []
+
+    def depth(self) -> int:
         """
         Compute depth of tree (maximal path length to a leaf).
 
@@ -102,56 +152,92 @@ def depth(self):
         int
             tree depth.
         """
+        return self._depth([0], d=0)
 
-        return self._depth([0], 0)
-
-    def random_function(self):
-        """select a random function"""
-        return self.function_set[randint(0, self.n_functions - 1)]
-
-    def random_terminal(self):
-        """Select a random terminal or create an ERC terminal"""
-        if self.erc_range is None:
-            node = self.terminal_set[randint(0, self.n_terminals - 1)]
-        else:
-            if random() > 0.5:
-                node = self.terminal_set[randint(0, self.n_terminals - 1)]
-            else:
-                node = round(uniform(*self.erc_range), 4)
-        return node
-
-    def _execute(self, pos, **kwargs):
-        """Recursively execute the tree by traversing it in a depth-first order 
-           (pos is a size-1 list so as to pass "by reference" on successive recursive calls)."""
+    def _depth(self, pos, d) -> int:
+        """Recursively compute depth
+        (pos is a size-1 list so as to pass "by reference"
+        on successive recursive calls).
+        """
 
         node = self.tree[pos[0]]
 
-        if node in self.function_set:
-            arglist = []
-            for i in range(self.arity[node]):
+        depths = []
+        if isinstance(node, FunctionNode):
+            for _ in range(node.n_args):
                 pos[0] += 1
-                res = self._execute(pos, **kwargs)
-                arglist.append(res)
-            return node(*arglist)
-        else:  # terminal
-            if isinstance(node, Number):  # terminal is a constant
-                return node
-            else:  # terminal is a variable, return its value
-                return kwargs[node]
+                depths.append(1 + self._depth(pos, d + 1))
+            return max(depths)
+        else:
+            return 0
 
-    def execute(self, *args, **kwargs):
+    def random_function(
+        self,
+        node_type: Optional[type] = None
+    ) -> Optional[FunctionNode]:
+        functions_types = {
+            func: get_return_type(func) for func in self.function_set
+        }
+        relevant_functions = [
+            func
+            for func in self.function_set
+            if functions_types[func] == node_type
+        ]
+
+        # Return None in case there are no functions of the given type
+        if not relevant_functions:
+            return None
+
+        func = random.choice(relevant_functions)
+        return FunctionNode(func)
+
+    def random_terminal(
+        self,
+        node_type: Optional[type] = None
+    ) -> Optional[TerminalNode]:
+        """Select a random terminal, including constants from ERC range"""
+        relevant_terminals = [
+            term
+            for term, term_type in self.terminal_set.items()
+            if term_type == node_type
+        ]
+
+        if self.erc_range is not None and (
+            node_type is None or issubclass(node_type, Number)
+        ):
+            relevant_terminals.append(
+                random.uniform(*self.erc_range)
+                if type(self.erc_range[0]) is float
+                else random.randint(*self.erc_range)
+            )
+
+        # Return None in case there are no terminals of the given type
+        if not relevant_terminals:
+            return None
+
+        terminal = random.choice(relevant_terminals)
+
+        # erc terminal will be typeless in untyped case,
+        # and int/float in typed case
+        default_type = None if node_type is None else type(terminal)
+        node_type = self.terminal_set.get(terminal, default_type)
+
+        return TerminalNode(terminal, node_type=node_type)
+
+    def execute(self, *args, **kwargs) -> object:
         """
-        Execute the program (tree). 
+        Execute the program (tree).
         Input is a numpy array or keyword arguments (but not both).
 
         Parameters
         ----------
         args : arguments
             A numpy array.
-        
+
         kwargs : keyword arguments
-            Input to program, including every variable in the terminal set as a keyword argument.
-            For example, if `terminal_set=['x', 'y', 'z', 0, 1, -1]` 
+            Input to program, including every variable
+            in the terminal set as a keyword argument.
+            For example, if `terminal_set=['x', 'y', 'z']`
             then call `execute(x=..., y=..., z=...)`.
 
         Returns
@@ -160,114 +246,295 @@ def execute(self, *args, **kwargs):
             Result of tree execution.
         """
 
+        if not self.tree:
+            raise ValueError("Tree is empty, cannot execute.")
+
         reshape = False
-        if args != ():  # numpy array -- convert to kwargs
+        if args:  # numpy array -- convert to kwargs
             try:
                 X = args[0]
-                kwargs = _generate_args(X)
+                kwargs = generate_args(X)
                 reshape = True
             except Exception:
-                raise ValueError(f'Bad argument to tree.execute, must be numpy array or kwargs: {args}')
+                raise ValueError(
+                    f"Bad argument to tree.execute, "
+                    f"must be np array or kwargs: {args}"
+                )
 
         kw = list(kwargs.keys())
 
-        bad_vars = [item for item in kw if item not in self.vars]
+        bad_vars = [item for item in kw if item not in self.terminal_set]
         if len(bad_vars) > 0:
-            raise ValueError(f'tree.execute received variable arguments not in terminal set: {bad_vars}')
+            raise ValueError(
+                f"received variable arguments not in terminal set: {bad_vars}"
+            )
 
-        missing_vars = [item for item in self.vars if item not in kw]
+        missing_vars = [item for item in self.terminal_set if item not in kw]
         if len(missing_vars) > 0:
             raise ValueError(
-                f'Some variable terminals were not passed to tree.execute as keyword arguments: {missing_vars}')
+                f"Missing variable terminals as execute kwargs: {missing_vars}"
+            )
 
         res = self._execute([0], **kwargs)
+
         if reshape and (isinstance(res, Number) or res.shape == np.shape(0)):
             # sometimes a tree degenrates to a scalar value
             res = np.full_like(X[:, 0], res)
         return res
 
-    def random_subtree(self):
-        # select a random node index from this individual's tree
-        rnd_i = randint(0, self.size() - 1)
-        # find the subtree's end
-        end_i = self._find_subtree_end([rnd_i])
-        # now we have a random subtree from this individual
-        return self.tree[rnd_i:end_i + 1]
-
-    def _find_subtree_end(self, pos):
-        """find index of final node of subtree that starts at `pos` 
-          (pos is a size-1 list so as to pass "by reference" on successive recursive calls)."""
+    def _execute(self, pos, **kwargs):
+        """
+        Recursively execute the tree by traversing it in a depth-first order
+        (pos is a size-1 list so as to pass "by reference"
+        on successive recursive calls).
+        """
 
         node = self.tree[pos[0]]
 
-        if node in self.function_set:
-            for i in range(self.arity[node]):
+        if isinstance(node, FunctionNode):
+            arglist = []
+            for _ in range(node.n_args):
                 pos[0] += 1
-                self._find_subtree_end(pos)
+                res = self._execute(pos, **kwargs)
+                arglist.append(res)
+            return node.function(*arglist)
+        else:  # terminal
+            if node.value in kwargs:
+                # terminal is a variable, return its value
+                return kwargs[node.value]
+            else:  # terminal is a constant
+                return node.value
 
-        return pos[0]
+    def filter_tree(self, filter_func: Callable) -> None:
+        return [node for node in self.tree if filter_func(node)]
 
-    def replace_subtree(self, subtree):
+    def get_random_numeric_node(self) -> Optional[TerminalNode]:
         """
-        Replace the subtree starting at `index` with `subtree`
+        Get a random leaf node of the tree.
 
         Parameters
         ----------
-        subtree - new subtree to replace the some existing subtree in this individual's tree
+        node_type : type, default=None
+            Type of node to return.
 
         Returns
         -------
-        None
+        TreeNode
+            Random leaf node.
         """
-
-        index = randint(0, self.size() - 1)  # select a random node (index)
-        end_i = self._find_subtree_end([index])
-        if isinstance(self.tree[end_i], list):
-            logger.debug(self.tree[end_i], list)
-        left_part = self.tree[:index]
-        right_part = self.tree[(end_i + 1):]
-        self.tree = left_part + subtree + right_part
-
-    def _node_label(self, node):
+        erc_nodes = self.filter_tree(
+            lambda node: isinstance(node, TerminalNode)
+            and isinstance(node.value, Number)
+        )
+        return random.choice(erc_nodes) if erc_nodes else None
+
+    def random_subtree(
+        self,
+        node_type: Optional[type] = None
+    ) -> Optional[List[TreeNode]]:
+        relevant_nodes = self.filter_tree(
+            lambda node: node.node_type is None  # untyped case
+            or (
+                node_type is None and node != self.root
+            )  # typed case with first invocation
+            or (
+                node.node_type == node_type and node != self.root
+            )  # typed case with subsequent invocations
+        )
+        if not relevant_nodes:
+            return None
+        subtree_root = random.choice(relevant_nodes)
+
+        return self._get_subtree_by_root(subtree_root)
+
+    def _get_subtree_by_root(self, subtree_root: TreeNode) -> List[TreeNode]:
+        start_i = self.tree.index(subtree_root)
+
+        if start_i == 0:
+            return self.tree
+
+        end_i = self._find_subtree_end([start_i])
+        return self.tree[start_i: end_i + 1]
+
+    def replace_subtree(
+        self, old_subtree: List[TreeNode], new_subtree: List[TreeNode]
+    ) -> None:
         """
-        return a string label for the node
+        Replace the subtree starting at `index` with `subtree`
 
         Parameters
         ----------
-        node - some node in the tree's function/terminal set
+        subtree - new subtree to replace the existing subtree in the tree
 
         Returns
         -------
-        node name - either a terminal (x0, x1,...) or a function (f_add, f_or, ...)
+        None
         """
-        return node.__name__ if node in self.function_set else str(node)
+        start_i = self.tree.index(old_subtree[0])
+        end_i = start_i + len(old_subtree)
+        left_part = self.tree[:start_i]
+        right_part = self.tree[end_i:]
+        self.tree = left_part + new_subtree + right_part
 
-    def _str_rec(self, prefix, pos, result, use_python_syntax):
-        """Recursively produce a simple textual printout of the tree 
-        (pos is a size-1 list so as to pass "by reference" on successive recursive calls)."""
+    def _find_subtree_end(self, pos):
+        """find index of final node of subtree that starts at `pos`
+        (pos is a size-1 list so as to pass "by reference" on successive recursive calls).
+        """
 
         node = self.tree[pos[0]]
-        if node in self.function_set:
-            result.append(f'{prefix}{self._node_label(node)}{"(" if use_python_syntax else ""}\n')
-            for i in range(self.arity[node]):
+
+        if (
+            isinstance(node, FunctionNode)
+            and node.function in self.function_set
+        ):
+            for _ in range(node.n_args):
                 pos[0] += 1
-                self._str_rec(prefix + "   ", pos, result, use_python_syntax)
-                if use_python_syntax:
-                    result.append(',')
-                if use_python_syntax or i < self.arity[node] - 1:
-                    result.append('\n')
-            if use_python_syntax:
-                result.append(prefix + ')')
-        else:  # terminal
-            result.append(f'{prefix}{self._node_label(node)}')
+                self._find_subtree_end(pos)
 
-    def __str__(self, use_python_syntax=False):
-        if use_python_syntax:
-            result = [f"def func_{self.id}({', '.join(self.terminal_set)}):\n   return "]
-        else:
-            result = []
-        self._str_rec("   " if use_python_syntax else "", [0], result, use_python_syntax)
-        return ''.join(result)
+        return pos[0]
+
+    def _handle_input_types(
+        self,
+        function_set: List[Callable],
+        terminal_set: Union[Dict[Any, type], List[str]],
+        root_type: type
+    ):
+        if function_set is None:
+            raise ValueError("Function set must be provided.")
+
+        if terminal_set is None:
+            raise ValueError("Terminal set must be provided.")
+
+        for t in function_set:
+            if not isinstance(t, Callable):
+                raise ValueError(
+                    f"Functions must be Callble, but {t} is of type {type(t)}"
+                )
+
+        if not isinstance(terminal_set, (list, dict)):
+            raise ValueError(
+                "Terminal set must be a list or a dictionary, "
+                f"got {type(terminal_set)}."
+            )
+
+        # untyped case - convert to dictionary of Nones.
+        if isinstance(terminal_set, list):
+            # check if any function has type hints
+            if any(f.__annotations__ for f in function_set):
+                raise ValueError(
+                    "Detected typed function with untyped terminal set. \
+                    Please provide a dictionary with types for terminals."
+                )
+
+            terminal_set = {t: None for t in terminal_set}
+            return function_set, terminal_set
+
+        # typed case - check every value is a type
+        for v in terminal_set.values():
+            if not isinstance(v, type):
+                raise ValueError(
+                    "Values in terminal set dictionary must be types, "
+                    f"but {v} is of type {type(v)}."
+                )
+
+        function_return_types = {
+            get_return_type(f) for f in function_set
+        }
+        if root_type not in function_return_types:
+            raise ValueError(
+                f"Detected a mismatch between root_type ({root_type}) "
+                f"and function set ({function_set}).\n"
+                f"Root type must be the return type of at least one function."
+            )
+        if self.erc_type is not None and self.erc_type not in function_return_types:
+            raise ValueError(
+                f"Detected a mismatch between ERC type ({self.erc_type}) "
+                f"and function set ({function_set}).\n"
+                f"ERC range should not be defined if there are no numeric functions."
+            )
+
+        # check terminals and functions type intersection
+        function_arg_types = {
+            t for func in function_set for t in get_func_types(func)[:-1]
+        }
+        terminal_types = set(terminal_set.values())
+
+        if self.erc_type:
+            terminal_types.add(self.erc_type)
+
+        if function_arg_types != terminal_types:
+            raise ValueError(
+                f"Function args type hints ({function_arg_types}) "
+                f"must match terminal types ({terminal_types})."
+            )
+
+        return function_set, terminal_set
+
+    def __str__(self) -> str:
+        """
+        Return a simple textual representation of the tree.
+
+        Returns
+        -------
+        str
+            String representation of the tree.
+
+        Examples
+        --------
+        untyped case:
+        >>> ut = Tree(terminal_set=[x, y, z], function_set=[f_add])
+        >>> ut.tree = [FunctionNode(f_add), TerminalNode(1), TerminalNode(2)]
+        >>> tree.show()
+        def func_0(x, y, z):
+            return f_add(x, y)
+
+        typed case:
+        >>> tt = Tree(terminal_set={x: float, y: float},
+        ...           function_set=[typed_add])
+        >>> tt.tree = [FunctionNode(typed_add),
+        ...            TerminalNode(1.0, float),
+        ...            TerminalNode(2.0, float)]
+        >>> str(tree)
+        def func_0(x: float, y: float) -> float:
+            return typed_add(1.0, 2.0)
+        """
+        args = (
+            list(self.terminal_set.keys())
+            if None in self.terminal_set.values()
+            else [f"{k}: {v.__name__}" for k, v in self.terminal_set.items()]
+        )
+        ret_type_str = (
+            ""
+            if None in self.terminal_set.values()
+            else f" -> {self.root.node_type.__name__}"
+        )
+
+        result = [
+            f"def func_{str(self.id)}({', '.join(args)}){ret_type_str}:\n  return "
+        ]
+        self._str_rec("", [0], result)
+        return "".join(result)
+
+    def _str_rec(self, prefix: str, pos: List[int], result: str) -> str:
+        """Recursively produce a simple textual printout of the tree
+        (pos is a size-1 list so as to pass "by reference" on successive
+         recursive calls).
+        """
+
+        node = self.tree[pos[0]]
+        if isinstance(node, FunctionNode):
+            result.append(f"{prefix}{str(node)}(\n")
+            prefix = prefix + "  " if pos[0] == 0 else prefix
+            for i in range(node.n_args):
+                pos[0] += 1
+                self._str_rec(prefix + "  ", pos, result)
+                if i < node.n_args - 1:
+                    result.append(",\n")
+                elif 0 < i < node.n_args:
+                    result.append("\n")
+            result.append(prefix + ")")
+        else:  # terminal
+            result.append(f"{prefix}{str(node)}")
 
     def show(self):
         """
@@ -277,6 +544,7 @@ def show(self):
         -------
         None.
         """
-        logger.info('\n' + str(self))
+        logger.info("\n" + str(self))
 
-# end class tree
+    def __repr__(self):
+        return str(self)
diff --git a/eckity/genetic_encodings/gp/tree/tree_node.py b/eckity/genetic_encodings/gp/tree/tree_node.py
new file mode 100644
index 00000000..6cc07f18
--- /dev/null
+++ b/eckity/genetic_encodings/gp/tree/tree_node.py
@@ -0,0 +1,102 @@
+from abc import ABC, abstractmethod
+from typing import Any, Callable, Optional
+
+from overrides import override
+
+from eckity.base.utils import arity
+
+from .utils import get_func_types
+
+
+class TreeNode(ABC):
+    """
+    GP Node
+
+    Attributes
+    ----------
+    node_type : type
+        node type
+    """
+
+    def __init__(self, node_type: Optional[type] = None) -> None:
+        self.node_type = node_type
+
+    @abstractmethod
+    def __str__(self) -> str:
+        pass
+
+    def __eq__(self, other):
+        return (
+            isinstance(other, TreeNode) and self.node_type is other.node_type
+        )
+
+    def __repr__(self):
+        return str(self)
+
+
+class FunctionNode(TreeNode):
+    def __init__(
+        self,
+        function: Callable,
+    ) -> None:
+        # infer the return type of the function
+        func_types = get_func_types(function)
+        return_type = func_types[-1] if func_types else None
+        self.n_args = arity(function)
+
+        if 0 < len(func_types) < self.n_args + 1:
+            raise ValueError(
+                f"Function {function.__name__} has missing type hints."
+                f"Please provide type hints for all arguments and return type."
+            )
+
+        super().__init__(return_type)
+
+        self.function = function
+
+    @override
+    def __eq__(self, other: object) -> bool:
+        """
+        Compare two FunctionNodes for equality.
+        Function nodes are equal if they have the same function.
+
+        Parameters
+        ----------
+        other : object
+            The object to compare to
+
+        Returns
+        -------
+        bool
+            True if the two FunctionNodes are equal, False otherwise
+        """
+        return (
+            super().__eq__(other)
+            and isinstance(other, FunctionNode)
+            and self.function == other.function
+        )
+
+    @override
+    def __str__(self) -> str:
+        return self.function.__name__
+
+
+class TerminalNode(TreeNode):
+    def __init__(
+        self,
+        value: Any,
+        node_type: Optional[type] = None
+    ) -> None:
+        super().__init__(node_type)
+        self.value = value
+
+    @override
+    def __eq__(self, other):
+        return (
+            super().__eq__(other)
+            and isinstance(other, TerminalNode)
+            and self.value is other.value
+        )
+
+    def __str__(self):
+        return str(self.value)
diff --git a/eckity/genetic_encodings/gp/tree/utils.py b/eckity/genetic_encodings/gp/tree/utils.py
index 165b8950..8e979650 100644
--- a/eckity/genetic_encodings/gp/tree/utils.py
+++ b/eckity/genetic_encodings/gp/tree/utils.py
@@ -1,11 +1,19 @@
 """
-This module implements some utility functions. 
+This module implements some utility functions.
 """
 
+from typing import Callable, Dict, List, Union, get_type_hints
 
-def create_terminal_set(X):
+import numpy as np
+
+from eckity.base.utils import arity
+
+
+def create_terminal_set(
+    X: np.ndarray, typed=False
+) -> Union[List[str], Dict[str, type]]:
     """
-    Create a terminal set from a 2D-shaped numpy array. 
+    Create a terminal set from a 2D-shaped numpy array.
 
     Example: \n
         X = array([[  4,   7,  -7, -10],  \n
@@ -25,15 +33,19 @@ def create_terminal_set(X):
     Terminal set (a list).
 
     """
-    
-    return [f'x{i}' for i in range(X.shape[1])]
 
+    features = [f"x{i}" for i in range(X.shape[1])]
+    if not typed:
+        return features
+    # convert numpy dtypes to python types
+    return {x: type(X[0, i].item()) for i, x in enumerate(features)}
 
-def _generate_args(X):
+
+def generate_args(X: np.ndarray) -> Dict[str, np.ndarray]:
     """
     Generate keyword arguments from a 2d array for passing to GPTree.execute.
-    
-    Example: \n 
+
+    Example: \n
         X = array([[  4,   7,  -7, -10],  \n
                    [  7,  -3,   3,  -8],  \n
                    [  8,  -5,  -3,  -1]])
@@ -42,16 +54,54 @@ def _generate_args(X):
         {'x0': array([4, 7, 8]), \n
          'x1': array([ 7, -3, -5]), \n
          'x2': array([-7,  3, -3]), \n
-         'x3': array([-10,  -8,  -1])}     
+         'x3': array([-10,  -8,  -1])}
 
     Returns
     -------
     kwargs : dict
-        A keyword dictionary that includes a value for every variable x_i in the terminal set. 
+        A keyword dictionary that includes a value for every variable x_i in the terminal set.
 
     """
-    
-    kwargs = dict()
-    for i in range(X.shape[1]):
-        kwargs[f'x{i}'] = X[:, i]
-    return kwargs
+
+    return {f"x{i}": X[:, i] for i in range(X.shape[1])}
+
+
+def get_func_types(f: Callable) -> List[type]:
+    """
+    Return list of function types in the following format:
+    [type_arg_1, type_arg_2, ..., type_arg_n, return_type]
+
+    Parameters
+    ----------
+    f : Callable
+        function (builtin or user-defined)
+
+    Returns
+    -------
+    List[type]
+        List of function types, sorted by argument order
+        with the return type as the last element.
+        For untyped functions, None is used.
+
+    Examples
+    --------
+    >>> def f(x: int, y: float) -> float:
+    ...     return x + y
+    >>> get_func_types(f)
+    [int, float, float]
+
+    >>> def f(x, y):
+    ...     return x + y
+    >>> get_func_types(f)
+    [None, None, None]
+    """
+    params_types: Dict = get_type_hints(f)
+    type_list = list(params_types.values())
+    if not type_list:
+        # If we don't have type hints, assign None
+        type_list = [None] * (arity(f) + 1)
+    return type_list
+
+
+def get_return_type(func: Callable) -> type:
+    return get_type_hints(func).get("return", None)
diff --git a/eckity/genetic_operators/crossovers/subtree_crossover.py b/eckity/genetic_operators/crossovers/subtree_crossover.py
index 15fc99c5..8af06599 100644
--- a/eckity/genetic_operators/crossovers/subtree_crossover.py
+++ b/eckity/genetic_operators/crossovers/subtree_crossover.py
@@ -1,47 +1,105 @@
-from eckity.genetic_operators.genetic_operator import GeneticOperator
 
+from typing import Any, List, Tuple, Optional
 
-class SubtreeCrossover(GeneticOperator):
-    def __init__(self, probability=1, arity=2, events=None):
+from overrides import override
+
+from eckity.genetic_encodings.gp import Tree, TreeNode
+from eckity.genetic_operators.failable_operator import FailableOperator
+
+
+class SubtreeCrossover(FailableOperator):
+    def __init__(self, probability=1.0, arity=2, events=None, attempts=1):
+        super().__init__(
+            probability=probability,
+            arity=arity,
+            events=events,
+            attempts=attempts,
+        )
         self.individuals = None
         self.applied_individuals = None
-        super().__init__(probability=probability, arity=arity, events=events)
 
-    def apply(self, individuals):
+    @override
+    def attempt_operator(
+        self, payload: Any, attempt_num: int
+    ) -> Tuple[bool, Any]:
         """
-        Perform subtree crossover between this tree and `other` tree:
-            1. Select random node from `other` tree
-            2. Get subtree rooted at selected node
-            1. Select a random node in this tree
-            2. Place `other` selected subtree at this node, replacing current subtree
+        Perform subtree crossover between a list of trees in a cyclic manner.
+        Meaning, the second individual will have a subtree from the first,
+        and the first individual will have a subtree from the last individual.
 
         Parameters
         ----------
-        individuals
-        select_func: callable
-        Selection method used to receive additional individuals to perform crossover on
+        payload: List[Individual]
+            List of Trees to perform crossover on
 
         individual: Tree
         tree individual to perform crossover on
 
         Returns
         -------
-        a new, modified individual
+        List
+            List of individuals after crossover (modified in-place)
         """
+        individuals = payload
 
-        assert len(individuals) == self.arity, f'Expected individuals list of size {self.arity}, got {len(individuals)}'
+        if len(individuals) != self.arity:
+            raise ValueError(
+                f"Expected individuals of size {self.arity}, "
+                f"got {len(individuals)}."
+            )
 
         self.individuals = individuals
 
-        # select a random subtree from each individual's tree
-        subtrees = [ind.random_subtree() for ind in individuals]
-
-        # assign the next individual's subtree to the current individual's tree
-        for i in range(len(individuals) - 1):
-            individuals[i].replace_subtree(subtrees[i+1])
+        subtrees: Optional[List[List[TreeNode]]] = self._pick_subtrees(
+            individuals
+        )
 
-        # to complete the crossover circle, assign the first subtree to the last individual
-        individuals[-1].replace_subtree(subtrees[0])
+        if subtrees is None:
+            return False, individuals
 
+        self._swap_subtrees(individuals, subtrees)
         self.applied_individuals = individuals
-        return individuals
+
+        return True, individuals
+
+    @staticmethod
+    def _pick_subtrees(
+        individuals: List[Tree],
+    ) -> Optional[List[List[TreeNode]]]:
+        # select a random subtree from first individual tree
+        first_subtree: List[TreeNode] = individuals[0].random_subtree()
+
+        if first_subtree is None:
+            # failed attempt
+            return None
+
+        m_type: type = first_subtree[0].node_type
+
+        # now select a random subtree from the rest of the individuals
+        # with regards to the type of the first subtree
+        rest_subtrees = [ind.random_subtree(m_type) for ind in individuals[1:]]
+
+        # fails if any subtree doesn't contain a node with of type `m_type`
+        if None in rest_subtrees:
+            return None
+
+        subtrees = [first_subtree] + rest_subtrees
+        return subtrees
+
+    @staticmethod
+    def _swap_subtrees(
+        individuals: List[TreeNode], subtrees: List[List[TreeNode]]
+    ) -> None:
+        """
+        Replace subtrees for all individuals in a cyclic manner
+        For n subtrees (st_1, st_2, ..., st_n):
+        st_n receives the subtree of st_n-1
+        st_n-1 receives the subtree of st_n-2
+        ...
+        st_2 receives the subtree of st_1
+        st_1 receives the subtree of st_n
+        """
+        for i in range(len(individuals) - 1, -1, -1):
+            individuals[i].replace_subtree(
+                old_subtree=subtrees[i], new_subtree=subtrees[i - 1]
+            )
diff --git a/eckity/genetic_operators/failable_operator.py b/eckity/genetic_operators/failable_operator.py
index 84bb0698..cb9d5f0e 100644
--- a/eckity/genetic_operators/failable_operator.py
+++ b/eckity/genetic_operators/failable_operator.py
@@ -1,5 +1,5 @@
 from abc import abstractmethod
-from typing import List, Tuple
+from typing import Any, List, Tuple
 
 from eckity.genetic_operators.genetic_operator import GeneticOperator
 
@@ -70,8 +70,8 @@ def apply(self, payload: object) -> object:
 
     @abstractmethod
     def attempt_operator(
-        self, payload: object, attempt_num: int
-    ) -> Tuple[bool, object]:
+        self, payload: Any, attempt_num: int
+    ) -> Tuple[bool, Any]:
         """
         A single attempt of the operator
 
@@ -90,12 +90,11 @@ def attempt_operator(
         """
         pass
 
-    @abstractmethod
-    def on_fail(self, payload: object) -> object:
+    def on_fail(self, payload: Any) -> Any:
         """
         What to do when all operator attempts failed
-
         This method is called once all operator attempts have failed
+        By default, this method does nothing.
 
         Parameters
         -------
@@ -105,7 +104,7 @@ def on_fail(self, payload: object) -> object:
 
         Returns
         -------
-        object
-            result value
+        Any
+            failure handling mechanism, does nothing by default
         """
-        pass
+        return payload
diff --git a/eckity/genetic_operators/mutations/erc_mutation.py b/eckity/genetic_operators/mutations/erc_mutation.py
index addd3616..083bdd17 100644
--- a/eckity/genetic_operators/mutations/erc_mutation.py
+++ b/eckity/genetic_operators/mutations/erc_mutation.py
@@ -1,14 +1,29 @@
-from random import choice, gauss
+from typing import Any, List, Tuple
 
-from eckity.event_based_operator import Operator
-from eckity.genetic_operators.genetic_operator import GeneticOperator
+from overrides import override
+import random
 
+from eckity.genetic_encodings.gp import Tree
+from eckity.genetic_operators import FailableOperator
 
-class ERCMutation(GeneticOperator):
-    def __init__(self, probability=1, arity=1, events=None):
+
+class ERCMutation(FailableOperator):
+    def __init__(
+        self,
+        probability=1.0,
+        arity=1,
+        mu=0,
+        sigma=1,
+        events=None,
+    ):
         super().__init__(probability=probability, arity=arity, events=events)
+        self.mu = mu
+        self.sigma = sigma
 
-    def apply(self, individuals):
+    @override
+    def attempt_operator(
+        self, payload: Any, attempt_num: int
+    ) -> Tuple[bool, Any]:
         """
         Perform ephemeral random constant (ERC) mutation: select an ERC node at random
         and add Gaussian noise to it.
@@ -17,11 +32,16 @@ def apply(self, individuals):
         -------
         None.
         """
-        for j in range(len(individuals)):
-            erc_indexes = [i for i, node in enumerate(individuals[j].tree) if
-                           isinstance(node, (int, float)) and node not in individuals[j].terminal_set]
-            if len(erc_indexes) > 0:
-                m_point = choice(erc_indexes)
-                individuals[j].tree[m_point] = round(individuals[j].tree[m_point] + gauss(mu=0, sigma=1), 4)
+        individuals: List[Tree] = payload
+
+        leaves = [ind.get_random_numeric_node() for ind in individuals]
+
+        if None in leaves:
+            return False, individuals
+
+        mu, sigma = self.mu, self.sigma
+        for terminal in leaves:
+            terminal.value = terminal.value + random.gauss(mu, sigma)
+
         self.applied_individuals = individuals
-        return individuals
+        return True, individuals
diff --git a/eckity/genetic_operators/mutations/subtree_mutation.py b/eckity/genetic_operators/mutations/subtree_mutation.py
index db1dde22..2811ad67 100644
--- a/eckity/genetic_operators/mutations/subtree_mutation.py
+++ b/eckity/genetic_operators/mutations/subtree_mutation.py
@@ -1,41 +1,80 @@
-from random import randint
+from typing import Any, List, Optional, Tuple
+
+from overrides import override
 
 from eckity.creators.gp_creators.grow import GrowCreator
-from eckity.genetic_operators.genetic_operator import GeneticOperator
+from eckity.genetic_encodings.gp import Tree, TreeNode
+from eckity.genetic_operators import FailableOperator
 
 
-class SubtreeMutation(GeneticOperator):
-    def __init__(self, probability=1, arity=1, init_depth=None, events=None):
-        super().__init__(probability=probability, arity=arity, events=events)
+class SubtreeMutation(FailableOperator):
+    def __init__(
+        self,
+        arity=1,
+        probability: float = 1.0,
+        init_depth: Tuple[int, int] = (2, 4),
+        events=None,
+        attempts=1,
+    ):
+        super().__init__(
+            probability=probability, arity=1, events=events, attempts=attempts
+        )
         self.init_depth = init_depth
+        self.tree_creator = None
 
-    def apply(self, individuals):
+    @override
+    def attempt_operator(
+        self, payload: Any, attempt_num: int
+    ) -> Tuple[bool, Any]:
         """
-        Perform subtree mutation: select a subtree at random to be replaced by a new, randomly generated subtree.
+        Perform subtree mutation: select a subtree at random
+        to be replaced by a new, randomly generated subtree.
 
         Returns
         -------
-        None.
+        Tuple[bool, Any]
+            A tuple containing a boolean indicating whether the operator was
+            successful and a list of the individuals.
         """
+        individuals: List[Tree] = payload
 
-        for ind in individuals:
-            init_depth = (
-                (ind.init_depth[0], randint(ind.init_depth[0], ind.init_depth[1]))
-                if self.init_depth is None
-                else self.init_depth
-            )
-            tree_creator = GrowCreator(
-                init_depth=init_depth,
-                function_set=ind.function_set,
-                terminal_set=ind.terminal_set,
-                erc_range=ind.erc_range,
-            )
+        old_subtrees: List[Optional[List[TreeNode]]] = [
+            ind.random_subtree() for ind in individuals
+        ]
 
-            # TODO refactor dummy individual creation, only the tree should be generated
-            subtree_individual = tree_creator.create_individuals(
-                1, ind.fitness.higher_is_better
-            )[0]
-            ind.replace_subtree(subtree_individual.tree)
+        if None in old_subtrees:
+            return False, individuals
+
+        self._swap_subtrees(individuals, old_subtrees)
 
         self.applied_individuals = individuals
-        return individuals
+        return True, individuals
+
+    def _swap_subtrees(
+        self, individuals: List[Tree], old_subtrees: List[List[TreeNode]]
+    ) -> None:
+
+        if self.tree_creator is None:
+            self.tree_creator = GrowCreator(
+                init_depth=self.init_depth,
+                function_set=individuals[0].function_set,
+                terminal_set=individuals[0].terminal_set,
+            )
+
+        for ind, old_subtree in zip(individuals, old_subtrees):
+            # generate a random tree with the same root type
+            # of the old subtree to not cause type errors
+
+            new_subtree = []
+            self.tree_creator.create_tree(
+                new_subtree,
+                ind.random_function,
+                ind.random_terminal,
+                depth=0,
+                node_type=old_subtree[0].node_type,
+            )
+
+            # replace the old subtree with the newly generated one
+            ind.replace_subtree(
+                old_subtree=old_subtree, new_subtree=new_subtree
+            )
diff --git a/eckity/sklearn_compatible/classification_evaluator.py b/eckity/sklearn_compatible/classification_evaluator.py
index 4e40a32a..3135911a 100644
--- a/eckity/sklearn_compatible/classification_evaluator.py
+++ b/eckity/sklearn_compatible/classification_evaluator.py
@@ -2,12 +2,16 @@
 This module implements the fitness evaluation class, which delivers the fitness function.
 You will need to implement such a class to work with your own problem and fitness function.
 """
+
 import numpy as np
 from sklearn.metrics import accuracy_score
 
-from eckity.evaluators.simple_individual_evaluator import SimpleIndividualEvaluator
+from eckity.evaluators.simple_individual_evaluator import (
+    SimpleIndividualEvaluator,
+)
+from scipy.special import expit as sigmoid
 
-CLASSIFICATION_THRESHOLD = 0
+CLF_METHODS = ("sigmoid", "argmax", "softmax")
 
 
 class ClassificationEvaluator(SimpleIndividualEvaluator):
@@ -16,11 +20,25 @@ class ClassificationEvaluator(SimpleIndividualEvaluator):
     All simple classes assume only one sub-population.
     """
 
-    def __init__(self, X=None, y=None, metric=accuracy_score):
+    def __init__(
+        self,
+        X=None,
+        y=None,
+        metric=accuracy_score,
+        n_classes=2,
+        clf_method=CLF_METHODS[0],
+    ):
         super().__init__()
         self.X = X
         self.y = y
         self.metric = metric
+        self.n_classes = n_classes
+
+        if clf_method not in CLF_METHODS:
+            raise ValueError(
+                "classification method must be one of", CLF_METHODS
+            )
+        self.clf_method = clf_method
 
     def set_context(self, context):
         """
@@ -60,4 +78,31 @@ def evaluate_individual(self, individual):
         return self.metric(y_true=self.y, y_pred=y_pred)
 
     def classify_individual(self, individual):
-        return np.where(individual.execute(self.X) > CLASSIFICATION_THRESHOLD, 1, 0)
+        clf_method_to_function = {
+            "sigmoid": self._clf_sigmoid,
+            "argmax": self._clf_argmax,
+            "softmax": self._clf_softmax,
+        }
+        selected_func = clf_method_to_function[self.clf_method]
+        return selected_func(individual)
+
+    def _clf_sigmoid(self, individual):
+        # normalize execute results between 0 and 1
+        probs = sigmoid(individual.execute(self.X))
+        # Create thresholds: 1/N, 2/N, ..., (N-1)/N
+        thresholds = np.linspace(0, 1, self.n_classes + 1)[1:-1]
+        return np.digitize(probs, thresholds)
+
+    def _clf_argmax(self, individual):
+        return self._clf_root_func(individual, "argmax")
+
+    def _clf_softmax(self, individual):
+        return self._clf_root_func(individual, "softmax")
+
+    def _clf_root_func(self, individual, method):
+        # assumes individual is a GP tree with argmax function in depth 1
+        if method not in individual.root.function.__name__:
+            raise ValueError(
+                f"Individual must have {method} function in depth 0 to classify."
+            )
+        return individual.execute(self.X)
diff --git a/eckity/sklearn_compatible/sk_classifier.py b/eckity/sklearn_compatible/sk_classifier.py
index 14330f7b..b65f3eda 100644
--- a/eckity/sklearn_compatible/sk_classifier.py
+++ b/eckity/sklearn_compatible/sk_classifier.py
@@ -2,7 +2,9 @@
 from sklearn.utils.validation import check_is_fitted
 
 from eckity.sklearn_compatible.sklearn_wrapper import SklearnWrapper
-from eckity.sklearn_compatible.classification_evaluator import ClassificationEvaluator
+from eckity.sklearn_compatible.classification_evaluator import (
+    ClassificationEvaluator,
+)
 
 
 class SKClassifier(SklearnWrapper, ClassifierMixin):
@@ -25,15 +27,17 @@ def predict(self, X):
         # Check is fit had been called
         check_is_fitted(self)
 
-        clf_eval: ClassificationEvaluator = self.algorithm.get_individual_evaluator()
+        clf_eval: ClassificationEvaluator = (
+            self.algorithm.get_individual_evaluator()
+        )
 
-        # y doesn't matter since we only need execute result and evolution has already finished
+        # ignore y since we only need execute result and evolution is finished
         clf_eval.set_context((X, None))
 
         return clf_eval.classify_individual(self.algorithm.best_of_run_)
 
     def predict_proba(self, X):
-        raise NotImplementedError('not implemented yet')
+        raise NotImplementedError("not implemented yet")
 
     def predict_log_proba(self, X):
-        raise NotImplementedError('not implemented yet')
+        raise NotImplementedError("not implemented yet")
diff --git a/eckity/subpopulation.py b/eckity/subpopulation.py
index 61c48382..a7d2d1ef 100644
--- a/eckity/subpopulation.py
+++ b/eckity/subpopulation.py
@@ -4,10 +4,7 @@
 import numpy as np
 
 from eckity.creators.creator import Creator
-from eckity.creators.gp_creators.full import FullCreator
-from eckity.genetic_operators import SubtreeCrossover, TournamentSelection
-from eckity.genetic_operators.mutations.erc_mutation import ERCMutation
-from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
+from eckity.genetic_operators import TournamentSelection
 
 logger = logging.getLogger(__name__)
 
@@ -79,13 +76,12 @@ def __init__(
 
         # verify valid creators and creation probability inputs
         if creators is None:
-            full_cr = FullCreator()
-            creators = [full_cr]
+            raise ValueError("Must specify at least one creator")
         elif isinstance(creators, Creator):
             creators = [creators]
         elif isinstance(creators, list):
             if len(creators) == 0:
-                raise ValueError("Creators list cannot be empty")
+                raise ValueError("Creators cannot be empty list")
             for creator in creators:
                 if not isinstance(creator, Creator):
                     raise ValueError(
@@ -93,8 +89,8 @@ def __init__(
                     )
         else:
             raise ValueError(
-                "Creators must be either a Creator or a list of Creators\n "
-                f"Received creators with unexpected type of {type(creators)}",
+                "Creators must be either a Creator or a list of Creators. "
+                f"Got unexpected type of {type(creators)}",
             )
 
         if pcr is None:
@@ -112,11 +108,7 @@ def __init__(
 
         # set default args
         if operators_sequence is None:
-            operators_sequence = [
-                SubtreeCrossover(arity=2, probability=0.9),
-                SubtreeMutation(arity=1, probability=0.7),
-                ERCMutation(arity=1, probability=0.1),
-            ]
+            raise ValueError("Must specify at least one operator")
         if selection_methods is None:
             selection_methods = [
                 (
@@ -139,8 +131,8 @@ def __init__(
 
         if self.n_elite == 0 and elitism_rate > 0.0:
             logger.warning(
-                "Detected elitism_rate > 0 but 0 elites. \
-                    Try increasing elitism_rate."
+                "Detected elitism_rate > 0 but 0 elites. "
+                "Try increasing elitism_rate."
             )
 
         self.individuals = individuals
diff --git a/examples/treegp/non_sklearn_mode/__init__.py b/examples/treegp/basic_mode/__init__.py
similarity index 100%
rename from examples/treegp/non_sklearn_mode/__init__.py
rename to examples/treegp/basic_mode/__init__.py
diff --git a/examples/treegp/basic_mode/multiplexer.py b/examples/treegp/basic_mode/multiplexer.py
new file mode 100644
index 00000000..b55c1216
--- /dev/null
+++ b/examples/treegp/basic_mode/multiplexer.py
@@ -0,0 +1,397 @@
+from time import time
+
+
+from eckity.algorithms.simple_evolution import SimpleEvolution
+from eckity.creators.gp_creators.full import FullCreator
+from eckity.base.untyped_functions import (
+    f_and,
+    f_or,
+    f_not,
+    f_if_then_else,
+)
+from eckity.genetic_operators.crossovers.subtree_crossover import (
+    SubtreeCrossover,
+)
+from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
+from eckity.genetic_operators.selections.tournament_selection import (
+    TournamentSelection,
+)
+from eckity.statistics.best_average_worst_statistics import (
+    BestAverageWorstStatistics,
+)
+from eckity.subpopulation import Subpopulation
+from eckity.termination_checkers.threshold_from_target_termination_checker import (
+    ThresholdFromTargetTerminationChecker,
+)
+
+from itertools import product
+from numbers import Number
+
+import pandas as pd
+import numpy as np
+
+from eckity.evaluators.simple_individual_evaluator import (
+    SimpleIndividualEvaluator,
+)
+
+NUM_SELECT_ENTRIES = 3
+NUM_INPUT_ENTRIES = 2 ** NUM_SELECT_ENTRIES
+NUM_COLUMNS = NUM_SELECT_ENTRIES + NUM_INPUT_ENTRIES
+NUM_ROWS = 2 ** NUM_COLUMNS
+
+
+def _target_func(s0, s1, s2, d0, d1, d2, d3, d4, d5, d6, d7):
+    """
+    Truth table of a 8-3 mux gate
+
+    Returns the value of d_i, where i is the decimal value of the binary
+    number gained from joining s0, s1 and s2 digits (see examples below).
+
+    Parameters
+    ----------
+    s0-s2: int
+        select values for the mux gate
+
+    d0-d7: int
+        input values for the mux gate
+
+    Returns
+    -------
+    int
+        0 (False) or 1 (True), depends on the values of the given parameters
+
+    Examples
+    -------
+    _target_func(s0=0, s1=0, s2=0, d0=1, ...) = 1 (value of input d0)
+    _target_func(s0=0, s1=0, s2=1, d0=1, d1=0, ...) = 0 (value of input d1)
+    """
+    return (
+            ((not s0) and (not s1) and (not s2) and d0)
+            or ((not s0) and (not s1) and s2 and d1)
+            or ((not s0) and s1 and (not s2) and d2)
+            or ((not s0) and s1 and s2 and d3)
+            or (s0 and (not s1) and (not s2) and d4)
+            or (s0 and (not s1) and s2 and d5)
+            or (s0 and s1 and (not s2) and d6)
+            or (s0 and s1 and s2 and d7)
+    )
+
+
+class MuxEvaluator(SimpleIndividualEvaluator):
+    """
+    Evaluator class for the Multiplexer problem, responsible of defining
+    a fitness evaluation method and evaluating it.
+
+    Attributes
+    -------
+    inputs: pd.DataFrame
+        Input columns representing all possible combinations
+        of all select and input values.
+
+    output: pandas Series
+        All possible output values. Values depend on the matching rows
+        from the inputs data frame.
+    """
+
+    def __init__(self):
+        super().__init__()
+
+        # construct a truth table of all combinations of ones and zeros
+        values = [
+            list(x) + [_target_func(*x)]
+            for x in product([0, 1], repeat=_target_func.__code__.co_argcount)
+        ]
+        truth_tbl = pd.DataFrame(
+            values,
+            columns=(list(_target_func.__code__.co_varnames) + ["output"]),
+        )
+
+        # split dataframe to input columns and an output column
+        self.inputs = truth_tbl.iloc[:, :-1]
+        self.output = truth_tbl["output"]
+
+    def evaluate_individual(self, individual):
+        """
+        Compute the fitness value of a given individual.
+
+        Fitness evaluation is done calculating the accuracy between the tree
+        execution result and the optimal result (multiplexer truth table).
+
+        Parameters
+        ----------
+        individual: Tree
+            The individual to compute the fitness value for.
+
+        Returns
+        -------
+        float
+            The evaluated fitness value of the given individual.
+            The value ranges from 0 (worst case) to 1 (best case).
+        """
+
+        # select entries columns
+        select_entries = {
+            f"s{i}": self.inputs[f"s{i}"] for i in range(NUM_SELECT_ENTRIES)
+        }
+        # input entries columns
+        input_entries = {
+            f"d{i}": self.inputs[f"d{i}"] for i in range(NUM_INPUT_ENTRIES)
+        }
+
+        exec_res = individual.execute(**select_entries, **input_entries)
+
+        # sometimes execute will return a single scalar (in constant trees)
+        if isinstance(exec_res, Number) or exec_res.shape == ():
+            exec_res = np.full((NUM_ROWS,), exec_res)
+
+        # Return normalized number of correct results
+        return np.mean(exec_res == self.output)
+
+
+def main():
+    """
+    The goal in the Multiplexer (Mux) problem to create a GP tree that approximates a Multiplexer logical gate.
+
+    Expected run time: less than a minute (on 2 cores, 2.5 GHz CPU)
+    Example output is provided at the end of the file.
+
+    References
+    ----------
+    DEAP Multiplexer Example: https://deap.readthedocs.io/en/master/examples/gp_multiplexer.html
+    """
+
+    start_time = time()
+
+    # The terminal set of the tree will contain the mux inputs (d0-d7 in a 8-3 mux gate),
+    # 3 select lines (s0-s2 in a 8-3 mux gate) and the constants 0 and 1
+    untyped_select_terminals = [f"s{i}" for i in range(NUM_SELECT_ENTRIES)]
+    untyped_input_terminals = [f"d{i}" for i in range(NUM_INPUT_ENTRIES)]
+    terminal_set = untyped_select_terminals + untyped_input_terminals
+
+    # Logical functions: and, or, not and if-then-else
+    function_set = [f_and, f_or, f_not, f_if_then_else]
+    # Initialize SimpleEvolution instance
+    algo = SimpleEvolution(
+        Subpopulation(
+            creators=FullCreator(
+                init_depth=(2, 4),
+                terminal_set=terminal_set,
+                function_set=function_set,
+                bloat_weight=0.00001,
+                erc_range=(0, 1)
+            ),
+            population_size=40,
+            # user-defined fitness evaluation method
+            evaluator=MuxEvaluator(),
+            # this is a maximization problem (fitness is accuracy), so higher fitness is better
+            higher_is_better=True,
+            elitism_rate=0.0,
+            # genetic operators sequence to be applied in each generation
+            operators_sequence=[
+                SubtreeCrossover(probability=0.8, arity=2),
+                SubtreeMutation(probability=0.1, arity=1),
+            ],
+            selection_methods=[
+                # (selection method, selection probability) tuple
+                (
+                    TournamentSelection(
+                        tournament_size=7, higher_is_better=True
+                    ),
+                    1,
+                )
+            ],
+        ),
+        max_workers=1,
+        max_generation=40,
+        termination_checker=ThresholdFromTargetTerminationChecker(
+            optimal=1, threshold=0.01
+        ),
+        statistics=BestAverageWorstStatistics(),
+        random_seed=10,
+    )
+
+    # evolve the generated initial population
+    algo.evolve()
+
+    # execute the best individual after the evolution process ends
+    exec1 = algo.execute(
+        s0=0, s1=0, s2=1, d0=0, d1=0, d2=1, d3=1, d4=1, d5=0, d6=0, d7=1
+    )
+    exec3 = algo.execute(
+        s0=0, s1=1, s2=1, d0=0, d1=0, d2=1, d3=1, d4=1, d5=0, d6=0, d7=1
+    )
+    exec7 = algo.execute(
+        s0=1, s1=1, s2=1, d0=0, d1=0, d2=1, d3=1, d4=1, d5=0, d6=0, d7=1
+    )
+    print("execute(s0=0, s1=1, s2=1, d1=0): expected = 0, actual =", exec1)
+    print("execute(s0=1, s1=1, s2=1, d3=1): expected = 1, actual =", exec3)
+    print("execute(s0=1, s1=1, s2=1, d7=1): expected = 1, actual =", exec7)
+
+    print("best pure fitness:", algo.best_of_run_.get_pure_fitness())
+
+    print("total time:", time() - start_time)
+
+
+if __name__ == "__main__":
+    main()
+
+"""
+Example output (s - select inputs, d - data inputs):
+f_if_then_else
+    f_and
+        f_not
+            f_if_then_else
+            f_or
+                f_and
+                    f_not
+                        s1
+                    f_if_then_else
+                        f_not
+                        f_or
+                            s0
+                            f_or
+                                f_and
+                                    s0
+                                    d3
+                                f_and
+                                    d0
+                                    s1
+                        f_and
+                        d1
+                        s2
+                        f_and
+                        f_or
+                            d4
+                            f_and
+                                f_not
+                                    s1
+                                f_if_then_else
+                                    f_not
+                                    f_or
+                                        s0
+                                        d5
+                                    f_and
+                                    d5
+                                    f_if_then_else
+                                        f_not
+                                            d2
+                                        f_not
+                                            d4
+                                        f_if_then_else
+                                            d5
+                                            d7
+                                            d6
+                                    f_and
+                                    d1
+                                    f_if_then_else
+                                        d1
+                                        s2
+                                        d0
+                        f_and
+                            d5
+                            d5
+                f_and
+                    f_if_then_else
+                        f_if_then_else
+                        s2
+                        s2
+                        s2
+                        d7
+                        d2
+                    s1
+            d7
+            f_if_then_else
+                f_and
+                    f_if_then_else
+                        d2
+                        f_if_then_else
+                        d3
+                        d7
+                        d7
+                        f_if_then_else
+                        d0
+                        s1
+                        d7
+                    s1
+                s1
+                d7
+        d7
+    s1
+    f_or
+        f_and
+            f_not
+            s1
+            f_if_then_else
+            f_not
+                f_or
+                    s0
+                    f_or
+                        f_and
+                        d4
+                        d0
+                        f_and
+                        f_if_then_else
+                            d2
+                            f_if_then_else
+                                d3
+                                f_and
+                                    f_or
+                                    f_and
+                                        d2
+                                        f_and
+                                            1
+                                            d7
+                                    d1
+                                    f_if_then_else
+                                    d5
+                                    f_or
+                                        d0
+                                        f_and
+                                            d0
+                                            d4
+                                    f_if_then_else
+                                        f_if_then_else
+                                            d0
+                                            d0
+                                            d7
+                                        f_or
+                                            d2
+                                            d7
+                                        s0
+                                d7
+                            f_if_then_else
+                                d0
+                                s1
+                                d7
+                        s1
+            f_and
+                d1
+                f_if_then_else
+                    s2
+                    s2
+                    d0
+            f_and
+                f_or
+                    d4
+                    f_and
+                        f_not
+                        s1
+                        f_if_then_else
+                        s2
+                        s2
+                        s2
+                f_and
+                    d5
+                    d5
+        f_and
+            f_if_then_else
+            f_if_then_else
+                f_and
+                    d5
+                    d5
+                s2
+                s2
+            d7
+            d2
+            s1
+"""
diff --git a/examples/treegp/non_sklearn_mode/symbolic_regression/sym_reg_main.py b/examples/treegp/basic_mode/symbolic_regression.py
similarity index 51%
rename from examples/treegp/non_sklearn_mode/symbolic_regression/sym_reg_main.py
rename to examples/treegp/basic_mode/symbolic_regression.py
index 1942a2d0..2531f033 100644
--- a/examples/treegp/non_sklearn_mode/symbolic_regression/sym_reg_main.py
+++ b/examples/treegp/basic_mode/symbolic_regression.py
@@ -3,34 +3,82 @@
 This is a non-sklearn setting so we use `evolve` and `execute`.
 """
 
+import numpy as np
+import pandas as pd
+
 from eckity.algorithms.simple_evolution import SimpleEvolution
-from eckity.breeders.simple_breeder import SimpleBreeder
-from eckity.creators.gp_creators.ramped_hh import RampedHalfAndHalfCreator
-from eckity.genetic_encodings.gp.tree.functions import (
-    f_add,
-    f_mul,
-    f_sub,
-    f_div,
-    f_sqrt,
-    f_log,
-    f_abs,
-    f_max,
-    f_min,
-    f_inv,
-    f_neg,
+from eckity.base.untyped_functions import *
+from eckity.creators.gp_creators.half import HalfCreator
+from eckity.evaluators.simple_individual_evaluator import (
+    SimpleIndividualEvaluator,
+)
+from eckity.genetic_operators.crossovers.subtree_crossover import (
+    SubtreeCrossover,
 )
-from eckity.genetic_operators.crossovers.subtree_crossover import SubtreeCrossover
 from eckity.genetic_operators.mutations.erc_mutation import ERCMutation
 from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
-from eckity.genetic_operators.selections.tournament_selection import TournamentSelection
-from eckity.statistics.best_average_worst_statistics import BestAverageWorstStatistics
+from eckity.genetic_operators.selections.tournament_selection import (
+    TournamentSelection,
+)
+from eckity.statistics.best_average_worst_statistics import (
+    BestAverageWorstStatistics,
+)
 from eckity.subpopulation import Subpopulation
 from eckity.termination_checkers.threshold_from_target_termination_checker import (
     ThresholdFromTargetTerminationChecker,
 )
-from examples.treegp.non_sklearn_mode.symbolic_regression.sym_reg_evaluator import (
-    SymbolicRegressionEvaluator,
-)
+
+
+class SymbolicRegressionEvaluator(SimpleIndividualEvaluator):
+    """
+    This module implements the Evaluator class, which delivers the fitness function.
+    You will need to implement such a class to work with your own problem and fitness function.
+    """
+
+    def __init__(self):
+        super().__init__()
+
+        data = np.random.uniform(-100, 100, size=(200, 3))
+        self.df = pd.DataFrame(data, columns=["x", "y", "z"])
+        self.df["target"] = self.target_func(
+            self.df["x"], self.df["y"], self.df["z"]
+        )
+
+    @staticmethod
+    def target_func(x, y, z):
+        """
+        True regression function, the individuals
+        Parameters
+        ----------
+        x, y, z: float
+            Values to the parameters of the function.
+
+        Returns
+        -------
+        float
+            The result of target function activation.
+        """
+        return x + 2 * y + 3 * z
+
+    def evaluate_individual(self, individual):
+        """
+        Parameters
+        ----------
+        individual : Tree
+            An individual program tree in the gp population, whose fitness needs to be computed.
+            Makes use of GPTree.execute, which runs the program.
+            Calling `gptree.execute` must use keyword arguments that match the terminal-set variables.
+            For example, if the terminal set includes `x` and `y` then the call is `gptree.execute(x=..., y=...)`.
+
+        Returns
+        -------
+        float
+            fitness value
+        """
+        x, y, z = self.df["x"], self.df["y"], self.df["z"]
+        return np.mean(
+            np.abs(individual.execute(x=x, y=y, z=z) - self.df["target"])
+        )
 
 
 def main():
@@ -68,16 +116,17 @@ def main():
         f_neg,
     ]
 
-    # terminal set, consisted of variables and constants
-    terminal_set = ["x", "y", "z", 0, 1, -1]
+    # terminal set, consisted of variables
+    terminal_set = ["x", "y", "z"]
 
     # Initialize the evolutionary algorithm
     algo = SimpleEvolution(
         Subpopulation(
-            creators=RampedHalfAndHalfCreator(
+            creators=HalfCreator(
                 init_depth=(2, 4),
                 terminal_set=terminal_set,
                 function_set=function_set,
+                erc_range=(-1.0, 1.0),
                 bloat_weight=0.0001,
             ),
             population_size=200,
@@ -88,19 +137,22 @@ def main():
             elitism_rate=0.05,
             # genetic operators sequence to be applied in each generation
             operators_sequence=[
-                SubtreeCrossover(probability=0.9, arity=2),
-                SubtreeMutation(probability=0.2, arity=1),
-                ERCMutation(probability=0.05, arity=1),
+                SubtreeCrossover(probability=0.9),
+                SubtreeMutation(probability=0.2),
+                ERCMutation(probability=0.05),
             ],
             selection_methods=[
                 # (selection method, selection probability) tuple
-                (TournamentSelection(tournament_size=4, higher_is_better=False), 1)
+                (
+                    TournamentSelection(
+                        tournament_size=4
+                    ),
+                    1,
+                )
             ],
         ),
-        breeder=SimpleBreeder(),
         max_workers=4,
-        max_generation=500,
-        # random_seed=0,
+        max_generation=40,
         termination_checker=ThresholdFromTargetTerminationChecker(
             optimal=0, threshold=0.001
         ),
diff --git a/examples/treegp/basic_mode/typed_multiplexer.py b/examples/treegp/basic_mode/typed_multiplexer.py
new file mode 100644
index 00000000..3d82dd7a
--- /dev/null
+++ b/examples/treegp/basic_mode/typed_multiplexer.py
@@ -0,0 +1,393 @@
+from time import time
+from eckity.algorithms.simple_evolution import SimpleEvolution
+from eckity.base.typed_functions import and2bools, or2bools, not2bools, if_then_else3bools
+from eckity.creators.gp_creators.full import FullCreator
+from eckity.genetic_operators.crossovers.subtree_crossover import (
+    SubtreeCrossover,
+)
+from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
+from eckity.genetic_operators.selections.tournament_selection import (
+    TournamentSelection,
+)
+from eckity.statistics.best_average_worst_statistics import (
+    BestAverageWorstStatistics,
+)
+from eckity.subpopulation import Subpopulation
+from eckity.termination_checkers.threshold_from_target_termination_checker import (
+    ThresholdFromTargetTerminationChecker,
+)
+
+from itertools import product
+from numbers import Number
+
+import pandas as pd
+import numpy as np
+
+from eckity.evaluators.simple_individual_evaluator import (
+    SimpleIndividualEvaluator,
+)
+
+NUM_SELECT_ENTRIES = 3
+NUM_INPUT_ENTRIES = 2 ** NUM_SELECT_ENTRIES
+NUM_COLUMNS = NUM_SELECT_ENTRIES + NUM_INPUT_ENTRIES
+NUM_ROWS = 2 ** NUM_COLUMNS
+
+
+def _target_func(s0, s1, s2, d0, d1, d2, d3, d4, d5, d6, d7):
+    """
+    Truth table of a 8-3 mux gate
+
+    Returns the value of d_i, where i is the decimal value of the binary
+    number gained from joining s0, s1 and s2 digits (see examples below).
+
+    Parameters
+    ----------
+    s0-s2: bool
+        select values for the mux gate
+
+    d0-d7: bool
+        input values for the mux gate
+
+    Returns
+    -------
+    bool
+        False or True, depends on the values of the given parameters
+
+    Examples
+    -------
+    _target_func(s0=False, s1=False, s2=False, d0=True, ...) = True (value of input d0)
+    _target_func(s0=False, s1=False, s2=True, d0=True, d1=False, ...) = False (value of input d1)
+    """
+    return (
+            ((not s0) and (not s1) and (not s2) and d0)
+            or ((not s0) and (not s1) and s2 and d1)
+            or ((not s0) and s1 and (not s2) and d2)
+            or ((not s0) and s1 and s2 and d3)
+            or (s0 and (not s1) and (not s2) and d4)
+            or (s0 and (not s1) and s2 and d5)
+            or (s0 and s1 and (not s2) and d6)
+            or (s0 and s1 and s2 and d7)
+    )
+
+
+class MuxEvaluator(SimpleIndividualEvaluator):
+    """
+    Evaluator class for the Multiplexer problem, responsible of defining
+    a fitness evaluation method and evaluating it.
+
+    Attributes
+    -------
+    inputs: pd.DataFrame
+        Input columns representing all possible combinations
+        of all select and input values.
+
+    output: pandas Series
+        All possible output values. Values depend on the matching rows
+        from the inputs data frame.
+    """
+
+    def __init__(self):
+        super().__init__()
+
+        # construct a truth table of all combinations of ones and zeros
+        values = [
+            list(map(bool, x)) + [_target_func(*map(bool, x))]
+            for x in product([False, True], repeat=_target_func.__code__.co_argcount)
+        ]
+        truth_tbl = pd.DataFrame(
+            values,
+            columns=(list(_target_func.__code__.co_varnames) + ["output"]),
+        )
+
+        # split dataframe to input columns and an output column
+        self.inputs = truth_tbl.iloc[:, :-1]
+        self.output = truth_tbl["output"]
+
+    def evaluate_individual(self, individual):
+        """
+        Compute the fitness value of a given individual.
+
+        Fitness evaluation is done calculating the accuracy between the tree
+        execution result and the optimal result (multiplexer truth table).
+
+        Parameters
+        ----------
+        individual: Tree
+            The individual to compute the fitness value for.
+
+        Returns
+        -------
+        float
+            The evaluated fitness value of the given individual.
+            The value ranges from 0 (worst case) to 1 (best case).
+        """
+
+        # select entries columns
+        select_entries = {
+            f"s{i}": self.inputs[f"s{i}"] for i in range(NUM_SELECT_ENTRIES)
+        }
+        # input entries columns
+        input_entries = {
+            f"d{i}": self.inputs[f"d{i}"] for i in range(NUM_INPUT_ENTRIES)
+        }
+
+        exec_res = individual.execute(**select_entries, **input_entries)
+
+        if isinstance(exec_res, Number) or exec_res.shape == ():
+            exec_res = np.full((NUM_ROWS,), bool(exec_res))
+
+        # Normalize number of correct results
+        return np.mean(exec_res == self.output)
+
+
+def main():
+    """
+    The goal in the Multiplexer (Mux) problem to create a GP tree that approximates a Multiplexer logical gate.
+
+    Expected run time: less than a minute (on 2 cores, 2.5 GHz CPU)
+    Example output is provided at the end of the file.
+
+    References
+    ----------
+    DEAP Multiplexer Example: https://deap.readthedocs.io/en/master/examples/gp_multiplexer.html
+    """
+
+    start_time = time()
+
+    # The terminal set of the tree will contain the mux inputs (d0-d7 in a 8-3 mux gate),
+    # 3 select lines (s0-s2 in a 8-3 mux gate) and the constants 0 and 1
+    # Terminal set for boolean inputs
+    typed_select_terminals = {f"s{i}": bool for i in range(NUM_SELECT_ENTRIES)}
+    typed_input_terminals = {f"d{i}": bool for i in range(NUM_INPUT_ENTRIES)}
+    terminal_set = {**typed_select_terminals, **typed_input_terminals}
+
+    # Logical functions: and, or, not and if-then-else
+    function_set = [and2bools, or2bools, not2bools, if_then_else3bools]
+
+    # Initialize SimpleEvolution instance
+    algo = SimpleEvolution(
+        Subpopulation(
+            creators=FullCreator(
+                init_depth=(2, 4),
+                terminal_set=terminal_set,
+                function_set=function_set,
+                bloat_weight=0.00001,
+                root_type=bool,
+            ),
+            population_size=40,
+            # user-defined fitness evaluation method
+            evaluator=MuxEvaluator(),
+            # this is a maximization problem (fitness is accuracy), so higher fitness is better
+            higher_is_better=True,
+            # genetic operators sequence to be applied in each generation
+            operators_sequence=[
+                SubtreeCrossover(probability=0.8, arity=2),
+                SubtreeMutation(probability=0.1, arity=1),
+            ],
+            selection_methods=[
+                # (selection method, selection probability) tuple
+                (
+                    TournamentSelection(
+                        tournament_size=7, higher_is_better=True
+                    ),
+                    1,
+                )
+            ],
+        ),
+        max_workers=1,
+        max_generation=40,
+        termination_checker=ThresholdFromTargetTerminationChecker(
+            optimal=1, threshold=0.01
+        ),
+        statistics=BestAverageWorstStatistics(),
+        random_seed=10,
+    )
+
+    # evolve the generated initial population
+    algo.evolve()
+
+    # execute the best individual after the evolution process ends
+    exec1 = algo.execute(
+        s0=False, s1=False, s2=True, d0=False, d1=False, d2=True, d3=True,
+        d4=True, d5=False, d6=False, d7=True
+    )
+    exec3 = algo.execute(
+        s0=False, s1=True, s2=True, d0=False, d1=False, d2=True, d3=True,
+        d4=True, d5=False, d6=False, d7=True
+    )
+    exec7 = algo.execute(
+        s0=True, s1=True, s2=True, d0=False, d1=False, d2=True, d3=True,
+        d4=True, d5=False, d6=False, d7=True
+    )
+    print("execute(s0=False, s1=True, s2=True, d1=False): expected = False, actual =", exec1)
+    print("execute(s0=True, s1=True, s2=True, d3=True): expected = True, actual =", exec3)
+    print("execute(s0=True, s1=True, s2=True, d7=True): expected = True, actual =", exec7)
+
+    print("best pure fitness:", algo.best_of_run_.get_pure_fitness())
+
+    print("total time:", time() - start_time)
+
+
+if __name__ == "__main__":
+    main()
+
+"""
+Example output (s - select inputs, d - data inputs):
+f_if_then_else
+    f_and
+        f_not
+            f_if_then_else
+            f_or
+                f_and
+                    f_not
+                        s1
+                    f_if_then_else
+                        f_not
+                        f_or
+                            s0
+                            f_or
+                                f_and
+                                    s0
+                                    d3
+                                f_and
+                                    d0
+                                    s1
+                        f_and
+                        d1
+                        s2
+                        f_and
+                        f_or
+                            d4
+                            f_and
+                                f_not
+                                    s1
+                                f_if_then_else
+                                    f_not
+                                    f_or
+                                        s0
+                                        d5
+                                    f_and
+                                    d5
+                                    f_if_then_else
+                                        f_not
+                                            d2
+                                        f_not
+                                            d4
+                                        f_if_then_else
+                                            d5
+                                            d7
+                                            d6
+                                    f_and
+                                    d1
+                                    f_if_then_else
+                                        d1
+                                        s2
+                                        d0
+                        f_and
+                            d5
+                            d5
+                f_and
+                    f_if_then_else
+                        f_if_then_else
+                        s2
+                        s2
+                        s2
+                        d7
+                        d2
+                    s1
+            d7
+            f_if_then_else
+                f_and
+                    f_if_then_else
+                        d2
+                        f_if_then_else
+                        d3
+                        d7
+                        d7
+                        f_if_then_else
+                        d0
+                        s1
+                        d7
+                    s1
+                s1
+                d7
+        d7
+    s1
+    f_or
+        f_and
+            f_not
+            s1
+            f_if_then_else
+            f_not
+                f_or
+                    s0
+                    f_or
+                        f_and
+                        d4
+                        d0
+                        f_and
+                        f_if_then_else
+                            d2
+                            f_if_then_else
+                                d3
+                                f_and
+                                    f_or
+                                    f_and
+                                        d2
+                                        f_and
+                                            1
+                                            d7
+                                    d1
+                                    f_if_then_else
+                                    d5
+                                    f_or
+                                        d0
+                                        f_and
+                                            d0
+                                            d4
+                                    f_if_then_else
+                                        f_if_then_else
+                                            d0
+                                            d0
+                                            d7
+                                        f_or
+                                            d2
+                                            d7
+                                        s0
+                                d7
+                            f_if_then_else
+                                d0
+                                s1
+                                d7
+                        s1
+            f_and
+                d1
+                f_if_then_else
+                    s2
+                    s2
+                    d0
+            f_and
+                f_or
+                    d4
+                    f_and
+                        f_not
+                        s1
+                        f_if_then_else
+                        s2
+                        s2
+                        s2
+                f_and
+                    d5
+                    d5
+        f_and
+            f_if_then_else
+            f_if_then_else
+                f_and
+                    d5
+                    d5
+                s2
+                s2
+            d7
+            d2
+            s1
+"""
diff --git a/examples/treegp/basic_mode/typed_symbolic_regression.py b/examples/treegp/basic_mode/typed_symbolic_regression.py
new file mode 100644
index 00000000..2f031f7b
--- /dev/null
+++ b/examples/treegp/basic_mode/typed_symbolic_regression.py
@@ -0,0 +1,172 @@
+"""
+A simple example optimizing a three-variable function.
+This is a non-sklearn setting so we use `evolve` and `execute`.
+"""
+
+import numpy as np
+import pandas as pd
+
+from eckity.algorithms.simple_evolution import SimpleEvolution
+from eckity.base.typed_functions import *
+from eckity.creators.gp_creators.half import HalfCreator
+from eckity.evaluators.simple_individual_evaluator import (
+    SimpleIndividualEvaluator,
+)
+from eckity.genetic_operators.crossovers.subtree_crossover import (
+    SubtreeCrossover,
+)
+from eckity.genetic_operators.mutations.erc_mutation import ERCMutation
+from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
+from eckity.genetic_operators.selections.tournament_selection import (
+    TournamentSelection,
+)
+from eckity.statistics.best_average_worst_statistics import (
+    BestAverageWorstStatistics,
+)
+from eckity.subpopulation import Subpopulation
+from eckity.termination_checkers.threshold_from_target_termination_checker import (
+    ThresholdFromTargetTerminationChecker,
+)
+
+
+class SymbolicRegressionEvaluator(SimpleIndividualEvaluator):
+    """
+    This module implements the Evaluator class, which delivers the fitness function.
+    You will need to implement such a class to work with your own problem and fitness function.
+    """
+
+    def __init__(self):
+        super().__init__()
+
+        data = np.random.uniform(-100, 100, size=(200, 3))
+        self.df = pd.DataFrame(data, columns=["x", "y", "z"])
+        self.df["target"] = self.target_func(
+            self.df["x"], self.df["y"], self.df["z"]
+        )
+
+    @staticmethod
+    def target_func(x, y, z):
+        """
+        True regression function, the individuals
+        Parameters
+        ----------
+        x, y, z: float
+            Values to the parameters of the function.
+
+        Returns
+        -------
+        float
+            The result of target function activation.
+        """
+        return x + 2 * y + 3 * z
+
+    def evaluate_individual(self, individual):
+        """
+        Parameters
+        ----------
+        individual : Tree
+            An individual program tree in the gp population, whose fitness needs to be computed.
+            Makes use of GPTree.execute, which runs the program.
+            Calling `gptree.execute` must use keyword arguments that match the terminal-set variables.
+            For example, if the terminal set includes `x` and `y` then the call is `gptree.execute(x=..., y=...)`.
+
+        Returns
+        -------
+        float
+            fitness value
+        """
+        x, y, z = self.df["x"], self.df["y"], self.df["z"]
+        return np.mean(
+            np.abs(individual.execute(x=x, y=y, z=z) - self.df["target"])
+        )
+
+
+def main():
+    """
+    Solve the Symbolic Regression problem by approximating the regression target function.
+
+    Expected runtime: less than a minute (on 2 cores, 2.5 GHz CPU)
+    Example of an optimal evolved tree:
+    f_add
+       f_add
+          f_add
+             f_add
+                y
+                z
+             x
+          z
+       f_add
+          y
+          z
+    """
+
+    # each node of the GP tree is either a terminal or a function
+    # function nodes, each has two children (which are its operands)
+    function_set = [
+        add2floats,
+        mul2floats,
+        sub2floats,
+        div2floats,
+        sqrt_float,
+        log_float,
+        abs_float,
+        max2floats,
+        min2floats,
+        inv_float,
+        neg_float,
+    ]
+
+    # terminal set, consisted of variables
+    terminal_set = {"x": float, "y": float, "z": float}
+
+    # Initialize the evolutionary algorithm
+    algo = SimpleEvolution(
+        Subpopulation(
+            creators=HalfCreator(
+                init_depth=(2, 4),
+                terminal_set=terminal_set,
+                function_set=function_set,
+                erc_range=(-1.0, 1.0),
+                bloat_weight=0.0001,
+                root_type=float,
+            ),
+            population_size=200,
+            # user-defined fitness evaluation method
+            evaluator=SymbolicRegressionEvaluator(),
+            # minimization problem (fitness is MAE), so higher fitness is worse
+            higher_is_better=False,
+            elitism_rate=0.05,
+            # genetic operators sequence to be applied in each generation
+            operators_sequence=[
+                SubtreeCrossover(probability=0.9),
+                SubtreeMutation(probability=0.2),
+                ERCMutation(probability=0.05),
+            ],
+            selection_methods=[
+                # (selection method, selection probability) tuple
+                (
+                    TournamentSelection(
+                        tournament_size=4
+                    ),
+                    1,
+                )
+            ],
+        ),
+        max_workers=4,
+        max_generation=40,
+        termination_checker=ThresholdFromTargetTerminationChecker(
+            optimal=0, threshold=0.001
+        ),
+        statistics=BestAverageWorstStatistics(),
+    )
+
+    # evolve the generated initial population
+    algo.evolve()
+
+    # execute the best individual after the evolution process ends, by assigning numeric values to the variable
+    # terminals in the tree
+    print(f"algo.execute(x=2,y=3,z=4): {algo.execute(x=2, y=3, z=4)}")
+
+
+if __name__ == "__main__":
+    main()
diff --git a/examples/treegp/non_sklearn_mode/multiplexer/mux_evaluator.py b/examples/treegp/non_sklearn_mode/multiplexer/mux_evaluator.py
deleted file mode 100644
index 3c0ea00e..00000000
--- a/examples/treegp/non_sklearn_mode/multiplexer/mux_evaluator.py
+++ /dev/null
@@ -1,120 +0,0 @@
-from itertools import product
-from numbers import Number
-
-import pandas as pd
-import numpy as np
-
-from eckity.evaluators.simple_individual_evaluator import SimpleIndividualEvaluator
-
-NUM_SELECT_ENTRIES = 3
-NUM_INPUT_ENTRIES = 2**NUM_SELECT_ENTRIES
-NUM_COLUMNS = NUM_SELECT_ENTRIES + NUM_INPUT_ENTRIES
-NUM_ROWS = 2**NUM_COLUMNS
-
-
-def _target_func(s0, s1, s2, d0, d1, d2, d3, d4, d5, d6, d7):
-    """
-    Truth table of a 8-3 mux gate
-
-    Returns the value of d_i, where i is the decimal value of the binary
-    number gained from joining s0, s1 and s2 digits (see examples below).
-
-    Parameters
-    ----------
-    s0-s2: int
-        select values for the mux gate
-
-    d0-d7: int
-        input values for the mux gate
-
-    Returns
-    -------
-    int
-        0 (False) or 1 (True), depends on the values of the given parameters
-
-    Examples
-    -------
-    _target_func(s0=0, s1=0, s2=0, d0=1, ...) = 1 (value of input d0)
-    _target_func(s0=0, s1=0, s2=1, d0=1, d1=0, ...) = 0 (value of input d1)
-    """
-    return (
-        ((not s0) and (not s1) and (not s2) and d0)
-        or ((not s0) and (not s1) and s2 and d1)
-        or ((not s0) and s1 and (not s2) and d2)
-        or ((not s0) and s1 and s2 and d3)
-        or (s0 and (not s1) and (not s2) and d4)
-        or (s0 and (not s1) and s2 and d5)
-        or (s0 and s1 and (not s2) and d6)
-        or (s0 and s1 and s2 and d7)
-    )
-
-
-class MuxEvaluator(SimpleIndividualEvaluator):
-    """
-    Evaluator class for the Multiplexer problem, responsible of defining
-    a fitness evaluation method and evaluating it.
-
-    Attributes
-    -------
-    inputs: pd.DataFrame
-        Input columns representing all possible combinations
-        of all select and input values.
-
-    output: pandas Series
-        All possible output values. Values depend on the matching rows
-        from the inputs data frame.
-    """
-
-    def __init__(self):
-        super().__init__()
-
-        # construct a truth table of all combinations of ones and zeros
-        values = [
-            list(x) + [_target_func(*x)]
-            for x in product([0, 1], repeat=_target_func.__code__.co_argcount)
-        ]
-        truth_tbl = pd.DataFrame(
-            values,
-            columns=(list(_target_func.__code__.co_varnames) + ["output"])
-        )
-
-        # split dataframe to input columns and an output column
-        self.inputs = truth_tbl.iloc[:, :-1]
-        self.output = truth_tbl["output"]
-
-    def evaluate_individual(self, individual):
-        """
-        Compute the fitness value of a given individual.
-
-        Fitness evaluation is done calculating the accuracy between the tree
-        execution result and the optimal result (multiplexer truth table).
-
-        Parameters
-        ----------
-        individual: Tree
-            The individual to compute the fitness value for.
-
-        Returns
-        -------
-        float
-            The evaluated fitness value of the given individual.
-            The value ranges from 0 (worst case) to 1 (best case).
-        """
-
-        # select entries columns
-        select_entries = {
-            f"s{i}": self.inputs[f"s{i}"] for i in range(NUM_SELECT_ENTRIES)
-        }
-        # input entries columns
-        input_entries = {
-            f"d{i}": self.inputs[f"d{i}"] for i in range(NUM_INPUT_ENTRIES)
-        }
-
-        exec_res = individual.execute(**select_entries, **input_entries)
-
-        # sometimes execute will return a single scalar (in constant trees)
-        if isinstance(exec_res, Number) or exec_res.shape == ():
-            exec_res = np.full((NUM_ROWS,), exec_res)
-
-        # Return normalized number of correct results
-        return np.mean(exec_res == self.output)
diff --git a/examples/treegp/non_sklearn_mode/multiplexer/mux_main.py b/examples/treegp/non_sklearn_mode/multiplexer/mux_main.py
deleted file mode 100644
index 9d93ca4a..00000000
--- a/examples/treegp/non_sklearn_mode/multiplexer/mux_main.py
+++ /dev/null
@@ -1,262 +0,0 @@
-from time import time
-
-from eckity.algorithms.simple_evolution import SimpleEvolution
-from eckity.breeders.simple_breeder import SimpleBreeder
-from eckity.creators.gp_creators.full import FullCreator
-from eckity.genetic_encodings.gp.tree.functions import (
-    f_and,
-    f_or,
-    f_not,
-    f_if_then_else,
-)
-from eckity.genetic_operators.crossovers.subtree_crossover import SubtreeCrossover
-from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
-from eckity.genetic_operators.selections.tournament_selection import TournamentSelection
-from eckity.statistics.best_average_worst_statistics import BestAverageWorstStatistics
-from eckity.subpopulation import Subpopulation
-from eckity.termination_checkers.threshold_from_target_termination_checker import (
-    ThresholdFromTargetTerminationChecker,
-)
-from examples.treegp.non_sklearn_mode.multiplexer.mux_evaluator import (
-    MuxEvaluator,
-    NUM_SELECT_ENTRIES,
-    NUM_INPUT_ENTRIES,
-)
-
-
-def main():
-    """
-    The goal in the Multiplexer (Mux) problem to create a GP tree that approximates a Multiplexer logical gate.
-
-    Expected run time: less than a minute (on 2 cores, 2.5 GHz CPU)
-    Example output (s - select inputs, d - data inputs):
-    f_if_then_else
-       f_and
-          f_not
-             f_if_then_else
-                f_or
-                   f_and
-                      f_not
-                         s1
-                      f_if_then_else
-                         f_not
-                            f_or
-                               s0
-                               f_or
-                                  f_and
-                                     s0
-                                     d3
-                                  f_and
-                                     d0
-                                     s1
-                         f_and
-                            d1
-                            s2
-                         f_and
-                            f_or
-                               d4
-                               f_and
-                                  f_not
-                                     s1
-                                  f_if_then_else
-                                     f_not
-                                        f_or
-                                           s0
-                                           d5
-                                     f_and
-                                        d5
-                                        f_if_then_else
-                                           f_not
-                                              d2
-                                           f_not
-                                              d4
-                                           f_if_then_else
-                                              d5
-                                              d7
-                                              d6
-                                     f_and
-                                        d1
-                                        f_if_then_else
-                                           d1
-                                           s2
-                                           d0
-                            f_and
-                               d5
-                               d5
-                   f_and
-                      f_if_then_else
-                         f_if_then_else
-                            s2
-                            s2
-                            s2
-                         d7
-                         d2
-                      s1
-                d7
-                f_if_then_else
-                   f_and
-                      f_if_then_else
-                         d2
-                         f_if_then_else
-                            d3
-                            d7
-                            d7
-                         f_if_then_else
-                            d0
-                            s1
-                            d7
-                      s1
-                   s1
-                   d7
-          d7
-       s1
-       f_or
-          f_and
-             f_not
-                s1
-             f_if_then_else
-                f_not
-                   f_or
-                      s0
-                      f_or
-                         f_and
-                            d4
-                            d0
-                         f_and
-                            f_if_then_else
-                               d2
-                               f_if_then_else
-                                  d3
-                                  f_and
-                                     f_or
-                                        f_and
-                                           d2
-                                           f_and
-                                              1
-                                              d7
-                                        d1
-                                     f_if_then_else
-                                        d5
-                                        f_or
-                                           d0
-                                           f_and
-                                              d0
-                                              d4
-                                        f_if_then_else
-                                           f_if_then_else
-                                              d0
-                                              d0
-                                              d7
-                                           f_or
-                                              d2
-                                              d7
-                                           s0
-                                  d7
-                               f_if_then_else
-                                  d0
-                                  s1
-                                  d7
-                            s1
-                f_and
-                   d1
-                   f_if_then_else
-                      s2
-                      s2
-                      d0
-                f_and
-                   f_or
-                      d4
-                      f_and
-                         f_not
-                            s1
-                         f_if_then_else
-                            s2
-                            s2
-                            s2
-                   f_and
-                      d5
-                      d5
-          f_and
-             f_if_then_else
-                f_if_then_else
-                   f_and
-                      d5
-                      d5
-                   s2
-                   s2
-                d7
-                d2
-             s1
-
-    References
-    ----------
-    DEAP Multiplexer Example: https://deap.readthedocs.io/en/master/examples/gp_multiplexer.html
-    """
-    start_time = time()
-
-    # The terminal set of the tree will contain the mux inputs (d0-d7 in a 8-3 mux gate),
-    # 3 select lines (s0-s2 in a 8-3 mux gate) and the constants 0 and 1
-    select_terminals = [f"s{i}" for i in range(NUM_SELECT_ENTRIES)]
-    input_terminals = [f"d{i}" for i in range(NUM_INPUT_ENTRIES)]
-    terminal_set = select_terminals + input_terminals + [0, 1]
-
-    # Logical functions: and, or, not and if-then-else
-    function_set = [f_and, f_or, f_not, f_if_then_else]
-
-    # Initialize SimpleEvolution instance
-    algo = SimpleEvolution(
-        Subpopulation(
-            creators=FullCreator(
-                init_depth=(2, 4),
-                terminal_set=terminal_set,
-                function_set=function_set,
-                bloat_weight=0.00001,
-            ),
-            population_size=40,
-            # user-defined fitness evaluation method
-            evaluator=MuxEvaluator(),
-            # this is a maximization problem (fitness is accuracy), so higher fitness is better
-            higher_is_better=True,
-            elitism_rate=0.0,
-            # genetic operators sequence to be applied in each generation
-            operators_sequence=[
-                SubtreeCrossover(probability=0.8, arity=2),
-                SubtreeMutation(probability=0.1, arity=1),
-            ],
-            selection_methods=[
-                # (selection method, selection probability) tuple
-                (TournamentSelection(tournament_size=7, higher_is_better=True), 1)
-            ],
-        ),
-        breeder=SimpleBreeder(),
-        max_workers=1,
-        max_generation=40,
-        termination_checker=ThresholdFromTargetTerminationChecker(
-            optimal=1, threshold=0.01
-        ),
-        statistics=BestAverageWorstStatistics(),
-        random_seed=10,
-    )
-
-    # evolve the generated initial population
-    algo.evolve()
-
-    # execute the best individual after the evolution process ends
-    exec1 = algo.execute(
-        s0=0, s1=0, s2=1, d0=0, d1=0, d2=1, d3=1, d4=1, d5=0, d6=0, d7=1
-    )
-    exec3 = algo.execute(
-        s0=0, s1=1, s2=1, d0=0, d1=0, d2=1, d3=1, d4=1, d5=0, d6=0, d7=1
-    )
-    exec7 = algo.execute(
-        s0=1, s1=1, s2=1, d0=0, d1=0, d2=1, d3=1, d4=1, d5=0, d6=0, d7=1
-    )
-    print("execute(s0=0, s1=1, s2=1, d1=0): expected value = 0, actual value =", exec1)
-    print("execute(s0=0, s1=0, s2=1, d3=1): expected value = 1, actual value =", exec3)
-    print("execute(s0=1, s1=1, s2=1, d7=1): expected value = 1, actual value =", exec7)
-
-    print("total time:", time() - start_time)
-
-
-if __name__ == "__main__":
-    main()
diff --git a/examples/treegp/non_sklearn_mode/symbolic_regression/__init__.py b/examples/treegp/non_sklearn_mode/symbolic_regression/__init__.py
deleted file mode 100644
index e69de29b..00000000
diff --git a/examples/treegp/non_sklearn_mode/symbolic_regression/sym_reg_evaluator.py b/examples/treegp/non_sklearn_mode/symbolic_regression/sym_reg_evaluator.py
deleted file mode 100644
index cde0a1d0..00000000
--- a/examples/treegp/non_sklearn_mode/symbolic_regression/sym_reg_evaluator.py
+++ /dev/null
@@ -1,58 +0,0 @@
-"""
-This module implements the fitness class, which delivers the fitness function.
-You will need to implement such a class to work with your own problem and fitness function.
-"""
-
-import numpy as np
-import pandas as pd
-
-from eckity.evaluators.simple_individual_evaluator import SimpleIndividualEvaluator
-
-
-def _target_func(x, y, z):
-    """
-    True regression function, the individuals
-    Parameters
-    ----------
-    x, y, z: float
-        Values to the parameters of the function.
-
-    Returns
-    -------
-    float
-        The result of target function activation.
-    """
-    return x + 2 * y + 3 * z
-
-
-class SymbolicRegressionEvaluator(SimpleIndividualEvaluator):
-    """
-    Compute the fitness of an individual.
-    """
-
-    def __init__(self):
-        super().__init__()
-
-        # np.random.seed(0)
-
-        data = np.random.uniform(-100, 100, size=(200, 3))
-        self.df = pd.DataFrame(data, columns=['x', 'y', 'z'])
-        self.df['target'] = _target_func(self.df['x'], self.df['y'], self.df['z'])
-
-    def evaluate_individual(self, individual):
-        """
-        Parameters
-        ----------
-        individual : Tree
-            An individual program tree in the gp population, whose fitness needs to be computed.
-            Makes use of GPTree.execute, which runs the program.
-            Calling `gptree.execute` must use keyword arguments that match the terminal-set variables.
-            For example, if the terminal set includes `x` and `y` then the call is `gptree.execute(x=..., y=...)`.
-
-        Returns
-        -------
-        float
-            fitness value
-        """
-        x, y, z = self.df['x'], self.df['y'], self.df['z']
-        return np.mean(np.abs(individual.execute(x=x, y=y, z=z) - self.df['target']))
diff --git a/examples/treegp/sklearn_mode/breast_cancer.py b/examples/treegp/sklearn_mode/breast_cancer.py
new file mode 100644
index 00000000..c2625a02
--- /dev/null
+++ b/examples/treegp/sklearn_mode/breast_cancer.py
@@ -0,0 +1,115 @@
+"""
+Solving a sklearn_mode problem created through scikit-learn's `load breast cancer`.
+This is an sklearn setting so we use `fit` and `predict`.
+"""
+
+from time import time
+
+from sklearn.datasets import load_breast_cancer
+from sklearn.metrics import accuracy_score
+from sklearn.model_selection import train_test_split
+
+from eckity import Subpopulation
+from eckity.algorithms import SimpleEvolution
+from eckity.base.untyped_functions import (
+    f_add,
+    f_div,
+    f_mul,
+    f_sub,
+)
+from eckity.genetic_encodings.gp.tree.utils import create_terminal_set
+from eckity.genetic_operators import (
+    SubtreeCrossover,
+    SubtreeMutation,
+    TournamentSelection,
+)
+from eckity.sklearn_compatible import ClassificationEvaluator, SKClassifier
+from eckity.statistics import BestAverageWorstSizeTreeStatistics
+from eckity.termination_checkers import ThresholdFromTargetTerminationChecker
+from eckity.creators import HalfCreator
+
+
+def main():
+    """
+    Evolve a GP Tree that classifies breast cancer,
+    using a dataset from sklearn
+    """
+    start_time = time()
+
+    # load the brest cancer dataset from sklearn
+    X, y = load_breast_cancer(return_X_y=True)
+
+    # Automatically generate a terminal set.
+    # Since there are 30 features, set terminal_set to: ['x0', 'x1', ..., 'x29']
+    terminal_set = create_terminal_set(X, typed=False)
+
+    # Define function set
+    function_set = [
+        f_add,
+        f_sub,
+        f_mul,
+        f_div,
+    ]
+
+    # Initialize SimpleEvolution instance
+    algo = SimpleEvolution(
+        Subpopulation(
+            creators=HalfCreator(
+                init_depth=(2, 4),
+                terminal_set=terminal_set,
+                function_set=function_set,
+                bloat_weight=0.0001,
+            ),
+            population_size=100,
+            evaluator=ClassificationEvaluator(metric=accuracy_score),
+            # maximization problem (fitness is accuracy), so higher fitness is better
+            higher_is_better=True,
+            elitism_rate=0.05,
+            # genetic operators sequence to be applied in each generation
+            operators_sequence=[
+                SubtreeCrossover(probability=0.9, arity=2),
+                SubtreeMutation(probability=0.2, arity=1),
+            ],
+            selection_methods=[
+                # (selection method, selection probability) tuple
+                (
+                    TournamentSelection(
+                        tournament_size=4, higher_is_better=True
+                    ),
+                    1,
+                )
+            ],
+        ),
+        max_workers=1,
+        max_generation=100,
+        # optimal fitness is 1, evolution ("training") process will be finished when best fitness <= threshold
+        termination_checker=ThresholdFromTargetTerminationChecker(
+            optimal=1, threshold=0.03
+        ),
+        statistics=BestAverageWorstSizeTreeStatistics(),
+    )
+    # wrap the basic evolutionary algorithm with a sklearn-compatible classifier
+    classifier = SKClassifier(algo)
+
+    # split brest cancer dataset to train and test set
+    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2)
+
+    # fit the model (perform evolution process)
+    classifier.fit(X_train, y_train)
+
+    print("best individual:\n", algo.best_of_run_.root)
+
+    # check training set results
+    print(
+        f"\nbest pure fitness over training set: {algo.best_of_run_.get_pure_fitness()}"
+    )
+
+    # check test set results by computing the accuracy score between the prediction result and the test set result
+    test_score = accuracy_score(y_test, classifier.predict(X_test))
+    print(f"test score: {test_score}")
+
+    print(f"Total runtime: {time() - start_time} seconds.")
+
+
+if __name__ == "__main__":
+    main()
diff --git a/examples/treegp/sklearn_mode/breast_cancer/__init__.py b/examples/treegp/sklearn_mode/breast_cancer/__init__.py
deleted file mode 100644
index e69de29b..00000000
diff --git a/examples/treegp/sklearn_mode/breast_cancer/breast_cancer_main.py b/examples/treegp/sklearn_mode/breast_cancer/breast_cancer_main.py
deleted file mode 100644
index 9237a58b..00000000
--- a/examples/treegp/sklearn_mode/breast_cancer/breast_cancer_main.py
+++ /dev/null
@@ -1,177 +0,0 @@
-"""
-Solving a sklearn_mode problem created through scikit-learn's `load breast cancer`.
-This is an sklearn setting so we use `fit` and `predict`.
-"""
-
-from time import time
-from sklearn.datasets import load_breast_cancer
-from sklearn.metrics import accuracy_score
-from sklearn.model_selection import train_test_split
-
-from eckity.algorithms.simple_evolution import SimpleEvolution
-from eckity.sklearn_compatible.sk_classifier import SKClassifier
-from eckity.breeders.simple_breeder import SimpleBreeder
-from eckity.creators.gp_creators.ramped_hh import RampedHalfAndHalfCreator
-from eckity.genetic_encodings.gp.tree.functions import (
-    f_add,
-    f_mul,
-    f_sub,
-    f_div,
-    f_neg,
-    f_sqrt,
-    f_log,
-    f_abs,
-    f_inv,
-    f_max,
-    f_min,
-)
-from eckity.genetic_encodings.gp.tree.utils import create_terminal_set
-from eckity.genetic_operators.crossovers.subtree_crossover import SubtreeCrossover
-from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
-from eckity.genetic_operators.selections.tournament_selection import TournamentSelection
-from eckity.statistics.best_avg_worst_size_tree_statistics import (
-    BestAverageWorstSizeTreeStatistics,
-)
-from eckity.subpopulation import Subpopulation
-from eckity.termination_checkers.threshold_from_target_termination_checker import (
-    ThresholdFromTargetTerminationChecker,
-)
-
-# Adding your own functions
-from eckity.sklearn_compatible.classification_evaluator import ClassificationEvaluator
-
-
-def main():
-    """
-    Evolve a GP Tree that classifies breast cancer, using a dataset from sklearn
-
-    Expected running time: ~6 minutes (on 2 cores, 2.5 GHz CPU)
-    Example output (with over 92% accuracy on test set):
-    f_div
-       f_log
-          x10
-       f_min
-          f_div
-             f_log
-                f_sqrt
-                   f_sqrt
-                      x7
-             f_max
-                x19
-                x23
-          f_div
-             f_log
-                f_mul
-                   f_add
-                      f_abs
-                         f_neg
-                            f_sqrt
-                               x19
-                      f_sub
-                         f_inv
-                            x23
-                         x26
-                   f_log
-                      x4
-             f_max
-                f_add
-                   f_abs
-                      f_max
-                         f_sub
-                            f_inv
-                               x23
-                            x26
-                         x23
-                   f_add
-                      x12
-                      f_sqrt
-                         f_sqrt
-                            f_log
-                               f_sqrt
-                                  f_sqrt
-                                     x7
-                f_inv
-                   f_mul
-                      x19
-                      x25
-    """
-    start_time = time()
-
-    # load the brest cancer dataset from sklearn
-    X, y = load_breast_cancer(return_X_y=True)
-
-    # Automatically generate a terminal set.
-    # Since there are 5 features, set terminal_set to: ['x0', 'x1', 'x2', ..., 'x9']
-    terminal_set = create_terminal_set(X)
-
-    # Define function set
-    function_set = [
-        f_add,
-        f_mul,
-        f_sub,
-        f_div,
-        f_sqrt,
-        f_log,
-        f_abs,
-        f_neg,
-        f_inv,
-        f_max,
-        f_min,
-    ]
-
-    # Initialize SimpleEvolution instance
-    algo = SimpleEvolution(
-        Subpopulation(
-            creators=RampedHalfAndHalfCreator(
-                init_depth=(2, 4),
-                terminal_set=terminal_set,
-                function_set=function_set,
-                bloat_weight=0.0001,
-            ),
-            population_size=1000,
-            evaluator=ClassificationEvaluator(),
-            # maximization problem (fitness is accuracy), so higher fitness is better
-            higher_is_better=True,
-            elitism_rate=0.05,
-            # genetic operators sequence to be applied in each generation
-            operators_sequence=[
-                SubtreeCrossover(probability=0.9, arity=2),
-                SubtreeMutation(probability=0.2, arity=1),
-            ],
-            selection_methods=[
-                # (selection method, selection probability) tuple
-                (TournamentSelection(tournament_size=4, higher_is_better=True), 1)
-            ],
-        ),
-        breeder=SimpleBreeder(),
-        max_workers=1,
-        max_generation=1000,
-        # optimal fitness is 1, evolution ("training") process will be finished when best fitness <= threshold
-        termination_checker=ThresholdFromTargetTerminationChecker(
-            optimal=1, threshold=0.03
-        ),
-        statistics=BestAverageWorstSizeTreeStatistics(),
-    )
-    # wrap the basic evolutionary algorithm with a sklearn-compatible classifier
-    classifier = SKClassifier(algo)
-
-    # split brest cancer dataset to train and test set
-    X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2)
-
-    # fit the model (perform evolution process)
-    classifier.fit(X_train, y_train)
-
-    # check training set results
-    print(
-        f"\nbest pure fitness over training set: {algo.best_of_run_.get_pure_fitness()}"
-    )
-
-    # check test set results by computing the accuracy score between the prediction result and the test set result
-    test_score = accuracy_score(y_test, classifier.predict(X_test))
-    print(f"test score: {test_score}")
-
-    print(f"Total runtime: {time() - start_time} seconds.")
-
-
-if __name__ == "__main__":
-    main()
diff --git a/examples/treegp/sklearn_mode/breast_cancer/classification_gridsearch.py b/examples/treegp/sklearn_mode/classification_gridsearch.py
similarity index 68%
rename from examples/treegp/sklearn_mode/breast_cancer/classification_gridsearch.py
rename to examples/treegp/sklearn_mode/classification_gridsearch.py
index 98eda2b8..ad04950f 100644
--- a/examples/treegp/sklearn_mode/breast_cancer/classification_gridsearch.py
+++ b/examples/treegp/sklearn_mode/classification_gridsearch.py
@@ -1,72 +1,91 @@
-from sklearn.datasets import load_breast_cancer
+import numpy as np
+from sklearn.datasets import make_classification
 from sklearn.model_selection import GridSearchCV
 
 from eckity.algorithms.simple_evolution import SimpleEvolution
-from eckity.sklearn_compatible.sk_classifier import SKClassifier
-from eckity.breeders.simple_breeder import SimpleBreeder
-from eckity.creators.gp_creators.ramped_hh import RampedHalfAndHalfCreator
-from eckity.genetic_encodings.gp.tree.functions import (
+from eckity.base.untyped_functions import (
+    f_abs,
     f_add,
-    f_mul,
-    f_sub,
     f_div,
-    f_neg,
-    f_sqrt,
-    f_log,
-    f_abs,
     f_inv,
+    f_log,
     f_max,
     f_min,
+    f_mul,
+    f_neg,
+    f_sqrt,
+    f_sub,
 )
+from eckity.breeders.simple_breeder import SimpleBreeder
+from eckity.creators.gp_creators.half import HalfCreator
 from eckity.genetic_encodings.gp.tree.utils import create_terminal_set
-from eckity.genetic_operators.crossovers.subtree_crossover import SubtreeCrossover
+from eckity.genetic_operators.crossovers.subtree_crossover import (
+    SubtreeCrossover,
+)
 from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
-from eckity.genetic_operators.selections.tournament_selection import TournamentSelection
-from eckity.statistics.best_average_worst_statistics import BestAverageWorstStatistics
+from eckity.genetic_operators.selections.tournament_selection import (
+    TournamentSelection,
+)
+from eckity.sklearn_compatible.classification_evaluator import (
+    ClassificationEvaluator,
+)
+from eckity.sklearn_compatible.sk_classifier import SKClassifier
+from eckity.statistics.best_average_worst_statistics import (
+    BestAverageWorstStatistics,
+)
 from eckity.subpopulation import Subpopulation
 from eckity.termination_checkers.threshold_from_target_termination_checker import (
     ThresholdFromTargetTerminationChecker,
 )
-from eckity.sklearn_compatible.classification_evaluator import ClassificationEvaluator
+
+# Adding your own types and functions
+
+t_argmax = type("argmax", (int,), {})
+
+
+def argmax3(x0: float, x1: float, x2: float) -> t_argmax:
+    return np.argmax([x0, x1, x2])
 
 
 def main():
     """
-    Demonstrate sklearn compatibility through use of grid search, solving the breast cancer problem.
+    Demonstrate sklearn compatibility through use of grid search, solving a classification problem.
     """
 
     # load the brest cancer dataset from sklearn
-    X, y = load_breast_cancer(return_X_y=True)
+    X, y = make_classification(
+        n_samples=100, n_features=10, n_classes=3, n_informative=5
+    )
 
     # Automatically generate a terminal set.
-    # Since there are 5 features, set terminal_set to: ['x0', 'x1', 'x2', ..., 'x9']
-    terminal_set = create_terminal_set(X)
+    # Since there are 10 features, set terminal_set to: ['x0', 'x1', ..., 'x9']
+    terminal_set = create_terminal_set(X, typed=False)
 
     # Define function set
     function_set = [
+        f_abs,
         f_add,
-        f_mul,
-        f_sub,
         f_div,
-        f_sqrt,
-        f_log,
-        f_abs,
-        f_neg,
         f_inv,
+        f_log,
         f_max,
         f_min,
+        f_mul,
+        f_neg,
+        f_sqrt,
+        f_sub,
     ]
 
     # Initialize SimpleEvolution instance
     algo = SimpleEvolution(
         Subpopulation(
-            creators=RampedHalfAndHalfCreator(
+            creators=HalfCreator(
                 init_depth=(2, 4),
                 terminal_set=terminal_set,
                 function_set=function_set,
                 bloat_weight=0.0001,
             ),
-            population_size=1000,
+            population_size=10,
             evaluator=ClassificationEvaluator(),
             higher_is_better=True,
             elitism_rate=0.05,
@@ -75,12 +94,16 @@ def main():
                 SubtreeMutation(probability=0.2, arity=1),
             ],
             selection_methods=[
-                (TournamentSelection(tournament_size=4, higher_is_better=True), 1)
+                (
+                    TournamentSelection(
+                        tournament_size=4, higher_is_better=True
+                    ),
+                    1,
+                )
             ],
         ),
-        breeder=SimpleBreeder(),
         max_workers=4,
-        max_generation=1000,
+        max_generation=10,
         termination_checker=ThresholdFromTargetTerminationChecker(
             optimal=1, threshold=0.03
         ),
diff --git a/examples/treegp/sklearn_mode/symbolic_regression/regression_gridsearch.py b/examples/treegp/sklearn_mode/regression_gridsearch.py
similarity index 85%
rename from examples/treegp/sklearn_mode/symbolic_regression/regression_gridsearch.py
rename to examples/treegp/sklearn_mode/regression_gridsearch.py
index d08ebe80..8c34d17d 100644
--- a/examples/treegp/sklearn_mode/symbolic_regression/regression_gridsearch.py
+++ b/examples/treegp/sklearn_mode/regression_gridsearch.py
@@ -1,19 +1,25 @@
 from time import time
 
+
 from sklearn.datasets import make_regression
 from sklearn.model_selection import GridSearchCV
 
 from eckity.algorithms.simple_evolution import SimpleEvolution
 from eckity.sklearn_compatible.sk_regressor import SKRegressor
-from eckity.breeders.simple_breeder import SimpleBreeder
-from eckity.creators.gp_creators.ramped_hh import RampedHalfAndHalfCreator
-from eckity.genetic_encodings.gp.tree.functions import f_add, f_mul, f_sub
+from eckity.creators.gp_creators.half import HalfCreator
+from eckity.base.untyped_functions import f_add, f_mul, f_sub
 from eckity.genetic_encodings.gp.tree.utils import create_terminal_set
-from eckity.genetic_operators.crossovers.subtree_crossover import SubtreeCrossover
+from eckity.genetic_operators.crossovers.subtree_crossover import (
+    SubtreeCrossover,
+)
 from eckity.genetic_operators.mutations.erc_mutation import ERCMutation
 from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
-from eckity.genetic_operators.selections.tournament_selection import TournamentSelection
-from eckity.statistics.best_average_worst_statistics import BestAverageWorstStatistics
+from eckity.genetic_operators.selections.tournament_selection import (
+    TournamentSelection,
+)
+from eckity.statistics.best_average_worst_statistics import (
+    BestAverageWorstStatistics,
+)
 from eckity.subpopulation import Subpopulation
 from eckity.termination_checkers.threshold_from_target_termination_checker import (
     ThresholdFromTargetTerminationChecker,
@@ -39,7 +45,7 @@ def main():
 
     # Automatically generate a terminal set.
     # Since there are 5 features, set terminal_set to: ['x0', 'x1', 'x2', 'x3', 'x4']
-    terminal_set = create_terminal_set(X)
+    terminal_set = create_terminal_set(X, typed=False)
 
     # Set function set to binary addition, binary multiplication and binary subtraction
     function_set = [f_add, f_mul, f_sub]
@@ -47,12 +53,12 @@ def main():
     # Initialize Simple Evolutionary Algorithm instance
     algo = SimpleEvolution(
         Subpopulation(
-            creators=RampedHalfAndHalfCreator(
+            creators=HalfCreator(
                 init_depth=(2, 4),
                 terminal_set=terminal_set,
                 function_set=function_set,
-                erc_range=(-100, 100),
                 bloat_weight=0.0001,
+                erc_range=(-100.0, 100.0),
             ),
             population_size=1000,
             # user-defined fitness evaluation method
@@ -68,15 +74,19 @@ def main():
             ],
             selection_methods=[
                 # (selection method, selection probability) tuple
-                (TournamentSelection(tournament_size=4, higher_is_better=False), 1)
+                (
+                    TournamentSelection(
+                        tournament_size=4, higher_is_better=False
+                    ),
+                    1,
+                )
             ],
         ),
-        breeder=SimpleBreeder(),
         max_workers=1,
         max_generation=1000,
         # optimal fitness is 0, evolution ("training") process will be finished when best fitness <= threshold
         termination_checker=ThresholdFromTargetTerminationChecker(
-            optimal=0, threshold=0.01
+            threshold=0.01
         ),
         statistics=BestAverageWorstStatistics(),
     )
diff --git a/examples/treegp/sklearn_mode/symbolic_regression/regression_pipeline.py b/examples/treegp/sklearn_mode/regression_pipeline.py
similarity index 86%
rename from examples/treegp/sklearn_mode/symbolic_regression/regression_pipeline.py
rename to examples/treegp/sklearn_mode/regression_pipeline.py
index 7b378358..82b2e1a3 100644
--- a/examples/treegp/sklearn_mode/symbolic_regression/regression_pipeline.py
+++ b/examples/treegp/sklearn_mode/regression_pipeline.py
@@ -1,5 +1,6 @@
 from time import time
 
+
 from sklearn.decomposition import PCA
 from sklearn.datasets import make_regression
 from sklearn.model_selection import train_test_split
@@ -9,15 +10,20 @@
 
 from eckity.algorithms.simple_evolution import SimpleEvolution
 from eckity.sklearn_compatible.sk_regressor import SKRegressor
-from eckity.breeders.simple_breeder import SimpleBreeder
-from eckity.creators.gp_creators.ramped_hh import RampedHalfAndHalfCreator
-from eckity.genetic_encodings.gp.tree.functions import f_add, f_mul, f_sub
+from eckity.creators.gp_creators.half import HalfCreator
+from eckity.base.untyped_functions import f_add, f_mul, f_sub
 from eckity.genetic_encodings.gp.tree.utils import create_terminal_set
-from eckity.genetic_operators.crossovers.subtree_crossover import SubtreeCrossover
+from eckity.genetic_operators.crossovers.subtree_crossover import (
+    SubtreeCrossover,
+)
 from eckity.genetic_operators.mutations.erc_mutation import ERCMutation
 from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
-from eckity.genetic_operators.selections.tournament_selection import TournamentSelection
-from eckity.statistics.best_average_worst_statistics import BestAverageWorstStatistics
+from eckity.genetic_operators.selections.tournament_selection import (
+    TournamentSelection,
+)
+from eckity.statistics.best_average_worst_statistics import (
+    BestAverageWorstStatistics,
+)
 from eckity.subpopulation import Subpopulation
 from eckity.termination_checkers.threshold_from_target_termination_checker import (
     ThresholdFromTargetTerminationChecker,
@@ -42,7 +48,7 @@ def main():
 
     # Automatically generate a terminal set.
     # Since there are 5 features, set terminal_set to: ['x0', 'x1', 'x2', 'x3', 'x4']
-    terminal_set = create_terminal_set(X)
+    terminal_set = create_terminal_set(X, typed=False)
 
     # Set function set to binary addition, binary multiplication and binary subtraction
     function_set = [f_add, f_mul, f_sub]
@@ -50,12 +56,12 @@ def main():
     # Initialize Simple Evolutionary Algorithm instance
     algo = SimpleEvolution(
         Subpopulation(
-            creators=RampedHalfAndHalfCreator(
+            creators=HalfCreator(
                 init_depth=(2, 4),
                 terminal_set=terminal_set,
                 function_set=function_set,
-                erc_range=(-100, 100),
                 bloat_weight=0.0001,
+                erc_range=(-100.0, 100.0),
             ),
             population_size=1000,
             # user-defined fitness evaluation method
@@ -71,15 +77,19 @@ def main():
             ],
             selection_methods=[
                 # (selection method, selection probability) tuple
-                (TournamentSelection(tournament_size=4, higher_is_better=False), 1)
+                (
+                    TournamentSelection(
+                        tournament_size=4, higher_is_better=False
+                    ),
+                    1,
+                )
             ],
         ),
-        breeder=SimpleBreeder(),
         max_workers=1,
         max_generation=1000,
         # optimal fitness is 0, evolution ("training") process will be finished when best fitness <= threshold
         termination_checker=ThresholdFromTargetTerminationChecker(
-            optimal=0, threshold=0.01
+            threshold=0.01
         ),
         statistics=BestAverageWorstStatistics(),
     )
diff --git a/examples/treegp/sklearn_mode/root_function_creator.py b/examples/treegp/sklearn_mode/root_function_creator.py
new file mode 100644
index 00000000..e9f43d9e
--- /dev/null
+++ b/examples/treegp/sklearn_mode/root_function_creator.py
@@ -0,0 +1,50 @@
+from typing import Callable, List, Optional
+
+from overrides import override
+
+from eckity.creators import FullCreator
+from eckity.genetic_encodings.gp import FunctionNode, TerminalNode, TreeNode
+
+
+class RootFunctionCreator(FullCreator):
+    """
+    Creator for trees with an immutable root function
+    (e.g. argmax/softmax for classification)
+    """
+
+    def __init__(
+        self,
+        root_function: Callable = None,
+        init_depth=None,
+        function_set=None,
+        terminal_set=None,
+        bloat_weight=0.1,
+        events=None,
+    ) -> None:
+        super().__init__(
+            init_depth=init_depth,
+            function_set=function_set,
+            terminal_set=terminal_set,
+            bloat_weight=bloat_weight,
+            events=events,
+        )
+        self.root_function = root_function
+
+    @override
+    def create_tree(
+        self,
+        tree: List[TreeNode],
+        random_function: Callable[[type], Optional[FunctionNode]],
+        random_terminal: Callable[[type], Optional[TerminalNode]],
+        depth: int = 0,
+        node_type: Optional[type] = None,
+    ) -> None:
+        root = FunctionNode(self.root_function)
+
+        self._add_children(
+            tree=tree,
+            fn_node=root,
+            random_function=random_function,
+            random_terminal=random_terminal,
+            depth=0,
+        )
diff --git a/examples/treegp/sklearn_mode/symbolic_regression/sklearn_sym_reg_main.py b/examples/treegp/sklearn_mode/symbolic_regression.py
similarity index 62%
rename from examples/treegp/sklearn_mode/symbolic_regression/sklearn_sym_reg_main.py
rename to examples/treegp/sklearn_mode/symbolic_regression.py
index 6f961c36..408d1d52 100644
--- a/examples/treegp/sklearn_mode/symbolic_regression/sklearn_sym_reg_main.py
+++ b/examples/treegp/sklearn_mode/symbolic_regression.py
@@ -11,16 +11,21 @@
 from sklearn.metrics import mean_absolute_error
 
 from eckity.algorithms.simple_evolution import SimpleEvolution
-from eckity.breeders.simple_breeder import SimpleBreeder
-from eckity.creators.gp_creators.ramped_hh import RampedHalfAndHalfCreator
-from eckity.genetic_encodings.gp.tree.functions import f_add, f_mul, f_sub
+from eckity.creators.gp_creators.half import HalfCreator
+from eckity.base.untyped_functions import f_add, f_mul, f_sub
 from eckity.genetic_encodings.gp.tree.utils import create_terminal_set
-from eckity.genetic_operators.crossovers.subtree_crossover import SubtreeCrossover
+from eckity.genetic_operators.crossovers.subtree_crossover import (
+    SubtreeCrossover,
+)
 from eckity.genetic_operators.mutations.erc_mutation import ERCMutation
 from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
-from eckity.genetic_operators.selections.tournament_selection import TournamentSelection
+from eckity.genetic_operators.selections.tournament_selection import (
+    TournamentSelection,
+)
 from eckity.sklearn_compatible.sk_regressor import SKRegressor
-from eckity.statistics.best_average_worst_statistics import BestAverageWorstStatistics
+from eckity.statistics.best_average_worst_statistics import (
+    BestAverageWorstStatistics,
+)
 from eckity.subpopulation import Subpopulation
 from eckity.termination_checkers.threshold_from_target_termination_checker import (
     ThresholdFromTargetTerminationChecker,
@@ -41,68 +46,6 @@ def main():
     """
     Solve a regression problem imported from sklearn `make_regression` function, using GP Trees.
     Expected run time: ~25 minutes (on 2 cores, 2.5 GHz CPU)
-    Example output (with an error of 0.09 on test set):
-    f_sub
-       f_mul
-          x3
-          57.9788
-       f_sub
-          f_mul
-             f_add
-                x1
-                x0
-             -69.8759
-          f_add
-             f_sub
-                f_add
-                   f_add
-                      f_add
-                         f_add
-                            f_add
-                               f_add
-                                  x1
-                                  x2
-                               x2
-                            f_add
-                               f_add
-                                  x2
-                                  x2
-                               f_sub
-                                  x1
-                                  x0
-                         x1
-                      f_add
-                         f_add
-                            f_add
-                               f_sub
-                                  f_sub
-                                     x1
-                                     f_mul
-                                        x1
-                                        -1.0477
-                                  x1
-                               x2
-                            x4
-                         x1
-                   x1
-                f_mul
-                   -33.4406
-                   f_add
-                      x4
-                      x4
-             f_sub
-                f_mul
-                   x1
-                   f_sub
-                      x2
-                      x2
-                f_mul
-                   f_add
-                      x2
-                      x0
-                   f_sub
-                      -19.0099
-                      -1.0477
     """
     start_time = time()
 
@@ -111,7 +54,7 @@ def main():
 
     # Automatically generate a terminal set.
     # Since there are 5 features, set terminal_set to: ['x0', 'x1', 'x2', 'x3', 'x4']
-    terminal_set = create_terminal_set(X)
+    terminal_set = create_terminal_set(X, typed=False)
 
     # Set function set to binary addition, binary multiplication and binary subtraction
     function_set = [f_add, f_mul, f_sub]
@@ -119,12 +62,12 @@ def main():
     # Initialize Simple Evolutionary Algorithm instance
     algo = SimpleEvolution(
         Subpopulation(
-            creators=RampedHalfAndHalfCreator(
+            creators=HalfCreator(
                 init_depth=(2, 4),
                 terminal_set=terminal_set,
                 function_set=function_set,
-                erc_range=(-100, 100),
                 bloat_weight=0.0001,
+                erc_range=(-100.0, 100.0),
             ),
             population_size=1000,
             # user-defined fitness evaluation method
@@ -140,10 +83,14 @@ def main():
             ],
             selection_methods=[
                 # (selection method, selection probability) tuple
-                (TournamentSelection(tournament_size=4, higher_is_better=False), 1)
+                (
+                    TournamentSelection(
+                        tournament_size=4, higher_is_better=False
+                    ),
+                    1,
+                )
             ],
         ),
-        breeder=SimpleBreeder(),
         max_workers=1,
         max_generation=1000,
         # optimal fitness is 0, evolution ("training") process will be finished when best fitness <= threshold
diff --git a/examples/treegp/sklearn_mode/symbolic_regression/__init__.py b/examples/treegp/sklearn_mode/symbolic_regression/__init__.py
deleted file mode 100644
index e69de29b..00000000
diff --git a/examples/treegp/sklearn_mode/typed_classification_gridsearch.py b/examples/treegp/sklearn_mode/typed_classification_gridsearch.py
new file mode 100644
index 00000000..5d1972b2
--- /dev/null
+++ b/examples/treegp/sklearn_mode/typed_classification_gridsearch.py
@@ -0,0 +1,128 @@
+import numpy as np
+from sklearn.datasets import make_classification
+from sklearn.model_selection import GridSearchCV
+
+from eckity.algorithms.simple_evolution import SimpleEvolution
+from eckity.base.typed_functions import (
+    abs_float,
+    add2floats,
+    div2floats,
+    inv_float,
+    log_float,
+    max2floats,
+    min2floats,
+    mul2floats,
+    neg_float,
+    sqrt_float,
+    sub2floats,
+)
+from eckity.creators.gp_creators.half import HalfCreator
+from eckity.genetic_encodings.gp.tree.utils import create_terminal_set
+from eckity.genetic_operators.crossovers.subtree_crossover import (
+    SubtreeCrossover,
+)
+from eckity.genetic_operators.mutations.subtree_mutation import SubtreeMutation
+from eckity.genetic_operators.selections.tournament_selection import (
+    TournamentSelection,
+)
+from eckity.sklearn_compatible.classification_evaluator import (
+    ClassificationEvaluator,
+)
+from eckity.sklearn_compatible.sk_classifier import SKClassifier
+from eckity.statistics.best_average_worst_statistics import (
+    BestAverageWorstStatistics,
+)
+from eckity.subpopulation import Subpopulation
+from eckity.termination_checkers.threshold_from_target_termination_checker import (
+    ThresholdFromTargetTerminationChecker,
+)
+
+# Adding your own types and functions
+
+t_argmax = type("argmax", (int,), {})
+
+
+def argmax3(x0: float, x1: float, x2: float) -> t_argmax:
+    return np.argmax([x0, x1, x2])
+
+
+def main():
+    """
+    Demonstrate sklearn compatibility through use of grid search, solving a classification problem.
+    """
+
+    # load the brest cancer dataset from sklearn
+    X, y = make_classification(
+        n_samples=100, n_features=10, n_classes=3, n_informative=5
+    )
+
+    # Automatically generate a terminal set.
+    # Since there are 10 features, set terminal_set to: ['x0', 'x1', ..., 'x9']
+    terminal_set = create_terminal_set(X, typed=True)
+
+    # Define function set
+    function_set = [
+        abs_float,
+        add2floats,
+        div2floats,
+        inv_float,
+        log_float,
+        max2floats,
+        min2floats,
+        mul2floats,
+        neg_float,
+        sqrt_float,
+        sub2floats,
+    ]
+
+    # Initialize SimpleEvolution instance
+    algo = SimpleEvolution(
+        Subpopulation(
+            creators=HalfCreator(
+                init_depth=(2, 4),
+                terminal_set=terminal_set,
+                function_set=function_set,
+                bloat_weight=0.0001,
+                root_type=float,
+            ),
+            population_size=10,
+            evaluator=ClassificationEvaluator(),
+            higher_is_better=True,
+            elitism_rate=0.05,
+            operators_sequence=[
+                SubtreeCrossover(probability=0.9, arity=2),
+                SubtreeMutation(probability=0.2, arity=1),
+            ],
+            selection_methods=[
+                (
+                    TournamentSelection(
+                        tournament_size=4, higher_is_better=True
+                    ),
+                    1,
+                )
+            ],
+        ),
+        max_workers=4,
+        max_generation=10,
+        termination_checker=ThresholdFromTargetTerminationChecker(
+            optimal=1, threshold=0.03
+        ),
+        statistics=BestAverageWorstStatistics(),
+    )
+    classifier = SKClassifier(algo)
+
+    print("Showcasing GridSearchCV...")
+
+    # Grid search parameters.
+    # The Grid Search model will fit the classifier with each parameter value
+    parameters = {"max_generation": [10, 20, 30]}
+
+    model = GridSearchCV(classifier, parameters)
+    model.fit(X, y)
+
+    print(f"best params: {model.best_params_}")
+    print(f"best score:  {model.best_score_}")
+
+
+if __name__ == "__main__":
+    main()
diff --git a/setup.py b/setup.py
index 46ffe595..f40c98f1 100644
--- a/setup.py
+++ b/setup.py
@@ -6,6 +6,9 @@
 with open("README.md", "r", encoding="utf-8") as fh:
     long_description = fh.read()
 
+with open("requirements.txt", 'r', encoding="utf-8") as fh:
+    requirements = [line.strip() for line in fh]
+
 setup(
     name='eckity',
     version=VERSION,
@@ -20,5 +23,5 @@
     },
     license='GNU GPLv3',
     packages=find_packages(),
-    install_requires=['numpy>=1.14.6', 'overrides>=6.1.0', 'pandas>=0.25.0'],
+    install_requires=requirements,
 )


Test patch:
diff --git a/eckity/algorithms/tests/test_simple_evolution.py b/eckity/algorithms/tests/test_simple_evolution.py
index 2c3966ea..83553bc6 100644
--- a/eckity/algorithms/tests/test_simple_evolution.py
+++ b/eckity/algorithms/tests/test_simple_evolution.py
@@ -1,5 +1,9 @@
 from eckity import Subpopulation
 from eckity.evaluators import SimpleIndividualEvaluator
+from eckity.creators import FullCreator
+from eckity.genetic_operators.mutations.identity_transformation import (
+    IdentityTransformation,
+)
 
 from ..simple_evolution import SimpleEvolution
 
@@ -11,5 +15,13 @@ def evaluate_individual(self, ind):
 
 def test_initialize():
     # creates a Tree-GP experiment by default
-    algo = SimpleEvolution(Subpopulation(DummyIndividualEvaluator()))
+    algo = SimpleEvolution(
+        Subpopulation(
+            DummyIndividualEvaluator(),
+            creators=FullCreator(
+                function_set=[lambda x: x], terminal_set=["x"]
+            ),
+            operators_sequence=[IdentityTransformation()],
+        )
+    )
     algo.initialize()
diff --git a/eckity/base/tests/test_base_utils.py b/eckity/base/tests/test_base_utils.py
new file mode 100644
index 00000000..df821dd9
--- /dev/null
+++ b/eckity/base/tests/test_base_utils.py
@@ -0,0 +1,18 @@
+import pytest
+
+from eckity.base.utils import arity
+from eckity.base.typed_functions import add2floats, or2floats
+from eckity.base.untyped_functions import f_mul, f_not
+
+
+@pytest.mark.parametrize(
+    "func, expected",
+    [
+        (add2floats, 2),
+        (or2floats, 2),
+        (f_mul, 2),
+        (f_not, 1),
+    ]
+)
+def test_arity(func, expected):
+    assert arity(func) == expected
diff --git a/eckity/base/tests/test_typed_functions.py b/eckity/base/tests/test_typed_functions.py
new file mode 100644
index 00000000..faf42071
--- /dev/null
+++ b/eckity/base/tests/test_typed_functions.py
@@ -0,0 +1,130 @@
+import numpy as np
+import pytest
+
+from eckity.base.typed_functions import *
+
+
+def test_add2floats():
+    assert add2floats(1.0, 2.0) == 3.0
+    assert add2floats(-1.0, -2.0) == -3.0
+    assert add2floats(0.0, 0.0) == 0.0
+
+
+def test_sub2floats():
+    assert sub2floats(5.0, 3.0) == 2.0
+    assert sub2floats(-1.0, -1.0) == 0.0
+    assert sub2floats(0.0, 1.0) == -1.0
+
+
+def test_mul2floats():
+    assert mul2floats(3.0, 2.0) == 6.0
+    assert mul2floats(-2.0, -2.0) == 4.0
+    assert mul2floats(0.0, 100.0) == 0.0
+
+
+def test_div2floats():
+    assert div2floats(6.0, 2.0) == 3.0
+    assert div2floats(1.0, 0.0) == 0.0
+    assert div2floats(0.0, 1.0) == 0.0
+
+
+def test_sqrt_float():
+    assert sqrt_float(4.0) == 2.0
+    assert sqrt_float(-4.0) == 2.0
+    assert sqrt_float(0.0) == 0.0
+
+
+def test_log_float():
+    assert log_float(1.0) == 0.0
+    assert log_float(0.0) == 0.0
+    assert log_float(-1.0) == 0.0
+
+
+def test_abs_float():
+    assert abs_float(-3.0) == 3.0
+    assert abs_float(3.0) == 3.0
+    assert abs_float(0.0) == 0.0
+
+
+def test_neg_float():
+    assert neg_float(5.0) == -5.0
+    assert neg_float(-5.0) == 5.0
+    assert neg_float(0.0) == 0.0
+
+
+def test_inv_float():
+    assert inv_float(2.0) == 0.5
+    assert inv_float(0.0) == 0.0
+    assert inv_float(-2.0) == -0.5
+
+
+def test_max2floats():
+    assert max2floats(3.0, 5.0) == 5.0
+    assert max2floats(-1.0, -2.0) == -1.0
+    assert max2floats(0.0, 0.0) == 0.0
+
+
+def test_min2floats():
+    assert min2floats(3.0, 5.0) == 3.0
+    assert min2floats(-1.0, -2.0) == -2.0
+    assert min2floats(0.0, 0.0) == 0.0
+
+
+def test_sin_float():
+    assert np.isclose(sin_float(np.pi / 2), 1.0)
+    assert np.isclose(sin_float(0.0), 0.0)
+    assert np.isclose(sin_float(-np.pi / 2), -1.0)
+
+
+def test_cos_float():
+    assert np.isclose(cos_float(0.0), 1.0)
+    assert np.isclose(cos_float(np.pi), -1.0)
+    assert np.isclose(cos_float(np.pi / 2), 0.0)
+
+
+def test_tan_float():
+    assert np.isclose(tan_float(0.0), 0.0)
+    assert np.isclose(tan_float(np.pi / 4), 1.0)
+    assert np.isclose(tan_float(-np.pi / 4), -1.0)
+
+
+def test_iflte0_floats():
+    assert iflte0_floats(0.0, 1.0, 2.0) == 1.0
+    assert iflte0_floats(-1.0, 1.0, 2.0) == 1.0
+    assert iflte0_floats(1.0, 1.0, 2.0) == 2.0
+
+
+def test_ifgt0_floats():
+    assert ifgt0_floats(1.0, 1.0, 2.0) == 1.0
+    assert ifgt0_floats(0.0, 1.0, 2.0) == 2.0
+    assert ifgt0_floats(-1.0, 1.0, 2.0) == 2.0
+
+
+def test_iflte_floats():
+    assert iflte_floats(1.0, 2.0, 3.0, 4.0) == 3.0
+    assert iflte_floats(2.0, 2.0, 3.0, 4.0) == 3.0
+    assert iflte_floats(3.0, 2.0, 3.0, 4.0) == 4.0
+
+
+def test_ifgt_floats():
+    assert ifgt_floats(3.0, 2.0, 3.0, 4.0) == 3.0
+    assert ifgt_floats(2.0, 2.0, 3.0, 4.0) == 4.0
+    assert ifgt_floats(1.0, 2.0, 3.0, 4.0) == 4.0
+
+
+def test_and2floats():
+    assert and2floats(1.0, 1.0) == 1.0
+    assert and2floats(0.0, 1.0) == 0.0
+    assert and2floats(0.0, 0.0) == 0.0
+
+
+def test_or2floats():
+    assert or2floats(1.0, 0.0) == 1.0
+    assert or2floats(0.0, 1.0) == 1.0
+    assert or2floats(0.0, 0.0) == 0.0
+
+
+def test_not2floats():
+    assert not2floats(0.0) == 1.0
+    assert not2floats(1.0) == 0.0
+    assert not2floats(-1.0) == 0.0
diff --git a/eckity/base/tests/test_untyped_functions.py b/eckity/base/tests/test_untyped_functions.py
new file mode 100644
index 00000000..67c50c0e
--- /dev/null
+++ b/eckity/base/tests/test_untyped_functions.py
@@ -0,0 +1,143 @@
+import numpy as np
+
+from eckity.base.untyped_functions import *
+
+
+def test_f_add():
+    assert np.array_equal(f_add(1, 2), 3)
+    assert np.array_equal(
+        f_add(np.array([1, 2]), np.array([3, 4])), np.array([4, 6])
+    )
+
+
+def test_f_sub():
+    assert np.array_equal(f_sub(5, 3), 2)
+    assert np.array_equal(
+        f_sub(np.array([10, 20]), np.array([5, 5])), np.array([5, 15])
+    )
+
+
+def test_f_mul():
+    assert np.array_equal(f_mul(3, 4), 12)
+    assert np.array_equal(
+        f_mul(np.array([2, 3]), np.array([4, 5])), np.array([8, 15])
+    )
+
+
+def test_f_div():
+    assert np.array_equal(f_div(6, 3), 2)
+    assert np.array_equal(
+        f_div(np.array([1, 2]), np.array([0.001, 0.002])),
+        np.array([0, 1000]),
+    )
+    assert np.array_equal(f_div(1, 0), 0)
+    assert np.array_equal(f_div(1, 0.0001), 0)
+
+
+def test_f_sqrt():
+    assert np.array_equal(f_sqrt(4), 2)
+    assert np.array_equal(f_sqrt(-4), 2)
+    assert np.array_equal(f_sqrt(np.array([4, -9])), np.array([2, 3]))
+
+
+def test_f_log():
+    assert np.array_equal(f_log(1), 0)
+    assert np.array_equal(f_log(0), 0)
+    assert np.array_equal(
+        f_log(np.array([0.001, 10])), np.array([0, np.log(10)])
+    )
+
+
+def test_f_abs():
+    assert np.array_equal(f_abs(-5), 5)
+    assert np.array_equal(f_abs(np.array([-1, 2, -3])), np.array([1, 2, 3]))
+
+
+def test_f_neg():
+    assert np.array_equal(f_neg(5), -5)
+    assert np.array_equal(f_neg(np.array([1, -2, 3])), np.array([-1, 2, -3]))
+
+
+def test_f_inv():
+    assert np.array_equal(f_inv(2), 0.5)
+    assert np.array_equal(f_inv(0), 0)
+    assert np.array_equal(f_inv(np.array([2, 0.0001])), np.array([0.5, 0]))
+
+
+def test_f_max():
+    assert np.array_equal(f_max(1, 2), 2)
+    assert np.array_equal(
+        f_max(np.array([1, 3]), np.array([2, 2])), np.array([2, 3])
+    )
+
+
+def test_f_min():
+    assert np.array_equal(f_min(1, 2), 1)
+    assert np.array_equal(
+        f_min(np.array([1, 3]), np.array([2, 2])), np.array([1, 2])
+    )
+
+
+def test_f_sin():
+    assert np.allclose(f_sin(0), 0)
+    assert np.allclose(f_sin(np.pi / 2), 1)
+
+
+def test_f_cos():
+    assert np.allclose(f_cos(0), 1)
+    assert np.allclose(f_cos(np.pi), -1)
+
+
+def test_f_tan():
+    assert np.allclose(f_tan(0), 0)
+    assert np.allclose(f_tan(np.pi / 4), 1)
+
+
+def test_f_iflte0():
+    assert np.array_equal(f_iflte0(-1, 2, 3), 2)
+    assert np.array_equal(f_iflte0(1, 2, 3), 3)
+
+
+def test_f_ifgt0():
+    assert np.array_equal(f_ifgt0(1, 2, 3), 2)
+    assert np.array_equal(f_ifgt0(-1, 2, 3), 3)
+
+
+def test_f_iflte():
+    assert np.array_equal(f_iflte(1, 2, 3, 4), 3)
+    assert np.array_equal(f_iflte(3, 2, 3, 4), 4)
+
+
+def test_f_ifgt():
+    assert np.array_equal(f_ifgt(3, 2, 3, 4), 3)
+    assert np.array_equal(f_ifgt(1, 2, 3, 4), 4)
+
+
+def test_f_and():
+    assert np.array_equal(f_and(1, 0), 0)
+    assert np.array_equal(
+        f_and(np.array([1, 0]), np.array([0, 1])), np.array([0, 0])
+    )
+
+
+def test_f_or():
+    assert np.array_equal(f_or(1, 0), 1)
+    assert np.array_equal(
+        f_or(np.array([1, 0]), np.array([0, 1])), np.array([1, 1])
+    )
+
+
+def test_f_not():
+    assert np.array_equal(f_not(1), 0)
+    assert np.array_equal(f_not(np.array([1, 0])), np.array([0, 1]))
+
+
+def test_f_if_then_else():
+    assert np.array_equal(f_if_then_else(True, 1, 0), 1)
+    assert np.array_equal(f_if_then_else(False, 1, 0), 0)
+    assert np.array_equal(
+        f_if_then_else(
+            np.array([True, False]), np.array([1, 2]), np.array([3, 4])
+        ),
+        np.array([1, 4]),
+    )
diff --git a/eckity/breeders/tests/test_simple_breeder.py b/eckity/breeders/tests/test_simple_breeder.py
index 2feac36c..adfa1755 100644
--- a/eckity/breeders/tests/test_simple_breeder.py
+++ b/eckity/breeders/tests/test_simple_breeder.py
@@ -4,6 +4,7 @@
 from eckity import Subpopulation
 from eckity.algorithms import SimpleEvolution
 from eckity.breeders import SimpleBreeder
+from eckity.creators import FullCreator
 from eckity.evaluators import SimpleIndividualEvaluator
 from eckity.genetic_operators import GeneticOperator
 
@@ -30,6 +31,10 @@ def test_incompatible_operator_arities():
     algo = SimpleEvolution(
         Subpopulation(
             DummyIndividualEvaluator(),
+            creators=FullCreator(
+                function_set=[lambda x: x],
+                terminal_set=['x']
+            ),
             population_size=10,
             operators_sequence=[DummyCrossover(arity=3)],
         ),
diff --git a/examples/treegp/non_sklearn_mode/multiplexer/__init__.py b/eckity/creators/gp_creators/tests/__init__.py
similarity index 100%
rename from examples/treegp/non_sklearn_mode/multiplexer/__init__.py
rename to eckity/creators/gp_creators/tests/__init__.py
diff --git a/eckity/creators/gp_creators/tests/test_full.py b/eckity/creators/gp_creators/tests/test_full.py
new file mode 100644
index 00000000..728c041c
--- /dev/null
+++ b/eckity/creators/gp_creators/tests/test_full.py
@@ -0,0 +1,27 @@
+from ....base.untyped_functions import f_add
+from ....fitness import SimpleFitness
+from ....genetic_encodings.gp import FunctionNode, TerminalNode, Tree
+from ..full import FullCreator
+
+
+def test_create_tree():
+    init_depth = (1, 1)
+    function_set = [f_add]
+    terminal_set = ["x"]
+
+    creator = FullCreator(init_depth, function_set, terminal_set)
+    tree_ind = Tree(
+        SimpleFitness(),
+        function_set=function_set,
+        terminal_set=terminal_set,
+        erc_range=None,
+    )
+    creator.create_tree(
+        tree_ind.tree, tree_ind.random_function, tree_ind.random_terminal
+    )
+
+    assert tree_ind.tree == [
+        FunctionNode(function_set[0]),
+        TerminalNode(terminal_set[0]),
+        TerminalNode(terminal_set[0]),
+    ]
diff --git a/eckity/creators/gp_creators/tests/test_grow.py b/eckity/creators/gp_creators/tests/test_grow.py
new file mode 100644
index 00000000..9bb21d4c
--- /dev/null
+++ b/eckity/creators/gp_creators/tests/test_grow.py
@@ -0,0 +1,60 @@
+from typing import Any, Callable, Dict, List, Union
+
+import pytest
+
+from eckity.base.typed_functions import sqrt_float
+from eckity.base.untyped_functions import f_log
+from eckity.creators import GrowCreator
+from eckity.genetic_encodings.gp import (
+    FunctionNode,
+    TerminalNode,
+    Tree,
+    TreeNode,
+)
+
+
+@pytest.mark.parametrize(
+    "function_set, terminal_set, expected",
+    [
+        (
+            [sqrt_float],
+            {1.0: float},
+            [FunctionNode(sqrt_float), TerminalNode(1.0, float)],
+        ),
+        (
+            [f_log],
+            [1.0],
+            [FunctionNode(f_log), TerminalNode(1.0)],
+        ),
+    ],
+)
+def test_add_children(
+    function_set: List[Callable],
+    terminal_set: Union[Dict[Any, type], List[Any]],
+    expected: List[TreeNode],
+):
+    grow_creator = GrowCreator(
+        init_depth=(1, 1),
+        function_set=function_set,
+        terminal_set=terminal_set,
+    )
+    
+    node_type = float if isinstance(terminal_set, dict) else None
+    tree_ind = Tree(
+        function_set=function_set,
+        terminal_set=terminal_set,
+        erc_range=None,
+        root_type=node_type,
+    )
+
+    node = tree_ind.random_function(node_type=node_type)
+    tree_ind.add_tree(node)
+
+    grow_creator._add_children(
+        tree_ind.tree,
+        node,
+        tree_ind.random_function,
+        tree_ind.random_terminal,
+        depth=0,
+    )
+    assert tree_ind.tree == expected
diff --git a/eckity/creators/gp_creators/tests/test_half.py b/eckity/creators/gp_creators/tests/test_half.py
new file mode 100644
index 00000000..a92cbab4
--- /dev/null
+++ b/eckity/creators/gp_creators/tests/test_half.py
@@ -0,0 +1,36 @@
+import pytest
+
+from eckity.base.untyped_functions import f_add
+from eckity.creators import HalfCreator
+
+
+def test_zero_group_size():
+    half_creator = HalfCreator(
+        init_depth=(2, 4), function_set=[f_add], terminal_set=["x"]
+    )
+    with pytest.raises(ValueError) as exc_info:
+        half_creator.create_individuals(4, True)
+
+    error_str = "Incompatible population size and init_depth."
+    assert error_str in str(exc_info.value)
+
+
+def test_create_individuals():
+    half_creator = HalfCreator(
+        init_depth=(3, 4), function_set=[f_add], terminal_set=["x"]
+    )
+    inds = half_creator.create_individuals(20, True)
+    assert len(inds) == 20
+
+    min_depth, max_depth = half_creator.init_depth
+    assert all(min_depth <= ind.depth() <= max_depth for ind in inds)
+
+    # assert at least 25% of the individuals are full trees of size 2^4
+    # and another 25% are full trees of size 2^5 - 1
+    sizes = [ind.size() for ind in inds]
+    max_depth = half_creator.init_depth[1]
+    assert sum(size == 2**max_depth - 1 for size in sizes) >= len(inds) / 4
+    assert (
+        sum(size == 2 ** (max_depth + 1) - 1 for size in sizes)
+        >= len(inds) / 4
+    )
diff --git a/eckity/genetic_encodings/gp/tree/tests/test_tree.py b/eckity/genetic_encodings/gp/tree/tests/test_tree.py
new file mode 100644
index 00000000..2625be71
--- /dev/null
+++ b/eckity/genetic_encodings/gp/tree/tests/test_tree.py
@@ -0,0 +1,442 @@
+from typing import List, Optional
+
+import pytest
+
+from eckity.base.typed_functions import (
+    add2floats,
+    div2floats,
+    mul2floats,
+    sub2floats,
+    or2bools,
+    and2bools,
+)
+from eckity.base.untyped_functions import f_add, f_div, f_mul, f_sub
+from eckity.fitness.gp_fitness import GPFitness
+from eckity.genetic_encodings.gp import (
+    FunctionNode,
+    TerminalNode,
+    Tree,
+    TreeNode,
+)
+
+
+class TestTree:
+    """
+    This class contains test cases for Tree and TreeNode classes
+    """
+
+    typed_functions = [add2floats, sub2floats, mul2floats, div2floats]
+    typed_terminals = {"x": float, "y": float}
+
+    untyped_functions = [f_add, f_sub, f_mul, f_div]
+    untyped_terminals = ["x", "y"]
+    typed_tree = Tree(
+        fitness=GPFitness(),
+        function_set=typed_functions,
+        terminal_set=typed_terminals,
+        root_type=float,
+    )
+    untyped_tree = Tree(
+        fitness=GPFitness(),
+        function_set=untyped_functions,
+        terminal_set=untyped_terminals,
+    )
+
+    @pytest.fixture
+    def setup(self) -> None:
+        """
+        Empties both trees before each test
+        Returns
+        -------
+        None
+        """
+        self.typed_tree.empty_tree()
+        self.untyped_tree.empty_tree()
+
+    def test_add_tree_root_typed(self, setup):
+        """
+        Test that add_child method adds child to the tree
+        """
+        typed_child = TerminalNode(1.0, float)
+        self.typed_tree.add_tree(typed_child)
+        assert self.typed_tree.size() == 1
+        assert self.typed_tree.root == typed_child
+        assert self.typed_tree.root.node_type is float
+
+    def test_add_tree_root_untyped(self, setup):
+        """
+        Test that add_child method adds child to the tree
+        """
+        untyped_child = TerminalNode(1)
+        self.untyped_tree.add_tree(untyped_child)
+        assert self.untyped_tree.root == untyped_child
+        assert self.untyped_tree.size() == 1
+        assert self.untyped_tree.root.node_type is None
+
+    def test_add_tree_inner_typed(self, setup):
+        """
+        Test that add_child method adds child to the tree
+        """
+        typed_parent = FunctionNode(add2floats)
+        children = [TerminalNode(float(i), float) for i in range(2)]
+
+        self.typed_tree.add_tree(typed_parent)
+        assert self.typed_tree.size() == 1
+        assert self.typed_tree.root == typed_parent
+        assert self.typed_tree.root.node_type is float
+
+        for i, child in enumerate(children):
+            self.typed_tree.add_tree(child)
+            assert self.typed_tree.size() == i + 2
+            assert self.typed_tree.tree[i + 1] == child
+            assert self.typed_tree.tree[i + 1].node_type is float
+
+    def test_add_tree_inner_untyped(self, setup):
+        """
+        Test that add_tree method adds child to the tree
+        """
+        untyped_child = TerminalNode(1)
+        self.untyped_tree.add_tree(untyped_child)
+        assert self.untyped_tree.root == untyped_child
+        assert self.untyped_tree.size() == 1
+        assert self.untyped_tree.root.node_type is None
+
+    @pytest.mark.parametrize(
+        "root, child",
+        [
+            (
+                FunctionNode(add2floats),
+                TerminalNode("1", str),
+            ),
+            (
+                FunctionNode(add2floats),
+                TerminalNode(1, None),
+            ),
+            (
+                FunctionNode(f_add),
+                TerminalNode(1, int),
+            ),
+        ],
+    )
+    def test_add_tree_bad_type(self, root, child):
+        tree = Tree(
+            erc_range=None,
+            function_set=[add2floats],
+            terminal_set={"x": float},
+            root_type=float,
+        )
+        tree.add_tree(root)
+
+        with pytest.raises(ValueError) as e:
+            tree.add_tree(child)
+        assert "Could not add node" in str(e.value)
+
+    def test_replace_subtree(self, setup):
+        """
+        Test that replace_subtree replaces subtree with another subtree
+        """
+
+        self.typed_tree.tree = list(
+            [
+                FunctionNode(add2floats),
+                TerminalNode(1.0, float),
+                TerminalNode(2.0, float),
+            ],
+        )
+
+        new_subtree = list(
+            [
+                FunctionNode(sub2floats),
+                TerminalNode(3.0, float),
+                TerminalNode(4.0, float),
+            ],
+        )
+
+        old_subtree = [self.typed_tree.tree[1]]
+
+        self.typed_tree.replace_subtree(old_subtree, new_subtree)
+        assert self.typed_tree.tree[1 : 1 + len(new_subtree)] == new_subtree
+        assert self.typed_tree.size() == 5
+
+    @pytest.mark.parametrize(
+        "typed, tree, expected",
+        [
+            (True, [TerminalNode(1, int)], 0),
+            (False, [TerminalNode(1)], 0),
+            (
+                True,
+                [
+                    FunctionNode(add2floats),
+                    TerminalNode(1.0, float),
+                    TerminalNode(2.0, float),
+                ],
+                1,
+            ),
+            (
+                False,
+                [
+                    FunctionNode(f_add),
+                    TerminalNode(1),
+                    TerminalNode(2),
+                ],
+                1,
+            ),
+            (
+                True,
+                [
+                    FunctionNode(add2floats),
+                    FunctionNode(add2floats),
+                    TerminalNode(1.0, float),
+                    TerminalNode(2.0, float),
+                    TerminalNode(2.0, float),
+                ],
+                2,
+            ),
+            (
+                False,
+                [
+                    FunctionNode(f_add),
+                    FunctionNode(f_add),
+                    TerminalNode(1),
+                    TerminalNode(2),
+                    TerminalNode(2),
+                ],
+                2,
+            ),
+        ],
+    )
+    def test_depth(
+        self, setup, typed: bool, tree: List[TerminalNode], expected: bool
+    ):
+        tree_ind = self.typed_tree if typed else self.untyped_tree
+        tree_ind.tree = tree
+        assert tree_ind.depth() == expected
+
+    @pytest.mark.parametrize(
+        "typed, tree, expected",
+        [
+            (
+                True,
+                [TerminalNode(1, int)],
+                f"def func_{typed_tree.id}(x: float, y: float) -> int:\n  return 1",
+            ),
+            (
+                False,
+                [TerminalNode(1)],
+                f"def func_{untyped_tree.id}(x, y):\n  return 1",
+            ),
+            (
+                True,
+                [
+                    FunctionNode(add2floats),
+                    TerminalNode(1.0, float),
+                    TerminalNode(2.0, float),
+                ],
+                f"def func_{typed_tree.id}(x: float, y: float) -> float:\n  return add2floats(\n    1.0,\n    2.0\n  )"
+            ),
+            (
+                False,
+                [
+                    FunctionNode(f_add),
+                    TerminalNode(1),
+                    TerminalNode(2),
+                ],
+                f"def func_{untyped_tree.id}(x, y):\n  return f_add(\n    1,\n    2\n  )"
+            ),
+            (
+                True,
+                [
+                    FunctionNode(add2floats),
+                    TerminalNode(1.0, float),
+                    FunctionNode(sub2floats),
+                    TerminalNode(2.0, float),
+                    TerminalNode(3.0, float),
+                ],
+                f"def func_{typed_tree.id}(x: float, y: float) -> float:\n  return add2floats(\n    1.0,\n    sub2floats(\n      2.0,\n      3.0\n    )\n  )"
+            ),
+            (
+                False,
+                [
+                    FunctionNode(f_add),
+                    TerminalNode(1),
+                    FunctionNode(f_sub),
+                    TerminalNode(2),
+                    TerminalNode(3),
+                ],
+                f"def func_{untyped_tree.id}(x, y):\n  return f_add(\n    1,\n    f_sub(\n      2,\n      3\n    )\n  )"
+            ),
+        ],
+    )
+    def test_str(self, setup, typed, tree, expected):
+        """
+        Test that str method returns a string representation of the tree
+        """
+        tree_ind = self.typed_tree if typed else self.untyped_tree
+        tree_ind.tree = tree
+
+        tree_str = str(tree_ind)
+        assert tree_str == expected
+
+    @pytest.mark.parametrize(
+        "typed, tree, root_indices, expected_results",
+        [
+            (
+                True,
+                [
+                    FunctionNode(add2floats),
+                    TerminalNode(1.0, float),
+                    FunctionNode(sub2floats),
+                    TerminalNode(2.0, float),
+                    TerminalNode(3.0, float),
+                ],
+                [0, 1, 2, 3, 4],
+                [
+                    [
+                        FunctionNode(add2floats),
+                        TerminalNode(1.0, float),
+                        FunctionNode(sub2floats),
+                        TerminalNode(2.0, float),
+                        TerminalNode(3.0, float),
+                    ],
+                    [TerminalNode(1.0, float)],
+                    [
+                        FunctionNode(sub2floats),
+                        TerminalNode(2.0, float),
+                        TerminalNode(3.0, float),
+                    ],
+                    [TerminalNode(2.0, float)],
+                    [TerminalNode(3.0, float)],
+                ],
+            ),
+            (
+                False,
+                [
+                    FunctionNode(f_add),
+                    TerminalNode(1),
+                    FunctionNode(f_sub),
+                    TerminalNode(2),
+                    TerminalNode(3),
+                ],
+                [0, 1, 2, 3, 4],
+                [
+                    [
+                        FunctionNode(f_add),
+                        TerminalNode(1),
+                        FunctionNode(f_sub),
+                        TerminalNode(2),
+                        TerminalNode(3),
+                    ],
+                    [TerminalNode(1)],
+                    [
+                        FunctionNode(f_sub),
+                        TerminalNode(2),
+                        TerminalNode(3),
+                    ],
+                    [TerminalNode(2)],
+                    [TerminalNode(3)],
+                ],
+            ),
+        ],
+    )
+    def test_get_subtree_by_root(
+        self,
+        setup,
+        typed,
+        tree,
+        root_indices: List[int],
+        expected_results: List[List[TreeNode]],
+    ):
+        tree_ind = self.typed_tree if typed else self.untyped_tree
+        tree_ind.tree = tree
+
+        for root_idx, expected in zip(root_indices, expected_results):
+            subtree_root = tree[root_idx]
+            assert tree_ind._get_subtree_by_root(subtree_root) == expected
+
+    @pytest.mark.parametrize(
+        "typed",
+        [True, False]
+    )
+    def test_random_function(
+        self,
+        typed: bool,
+    ):
+        tree = self.typed_tree if typed else self.untyped_tree
+        function_set = (
+            self.typed_functions
+            if typed
+            else self.untyped_functions
+        )
+        node_type = float if typed else None
+
+        for i in range(10):
+            function_node = tree.random_function(node_type)
+            assert function_node.function in function_set
+
+        class MyType:
+            pass
+
+        assert tree.random_function(MyType) is None
+
+    @pytest.mark.parametrize(
+        "typed",
+        [True, False]
+    )
+    def test_random_terminal(
+        self,
+        typed: bool,
+    ):
+        tree = self.typed_tree if typed else self.untyped_tree
+        terminal_set = (
+            self.typed_terminals
+            if typed
+            else self.untyped_terminals
+        )
+        node_type = float if typed else None
+
+        for i in range(10):
+            terminal = tree.random_terminal(node_type)
+            assert terminal.value in terminal_set
+
+        class MyType:
+            pass
+
+        assert tree.random_terminal(MyType) is None
+
+    @pytest.mark.parametrize(
+        "kwargs, error_message_substring",
+        [
+            # supplying int terminal with boolean functions only
+            (
+                {
+                    "function_set": [and2bools, or2bools],
+                    "terminal_set": {1: int, False: bool},
+                    "root_type": bool
+                },
+                "must match terminal types"
+            ),
+            # supplying root_type int with no functions that return int
+            (
+                {
+                    "function_set": [and2bools, or2bools],
+                    "terminal_set": {True: bool, False: bool},
+                    "root_type": int
+                },
+                "Root type must be the return type of at least one function."
+            ),
+            # supplying erc_range with no numeric functions
+            (
+                {
+                    "erc_range": (1, 2),
+                    "function_set": [and2bools, or2bools],
+                    "terminal_set": {True: bool, False: bool},
+                    "root_type": bool
+                },
+                "ERC range should not be defined if there are no numeric functions."
+            ),
+        ],
+    )
+    def test_init_bad_args(self, kwargs, error_message_substring):
+        with pytest.raises(ValueError) as e:
+            Tree(**kwargs)
+        assert error_message_substring in str(e.value)
diff --git a/eckity/genetic_encodings/gp/tree/tests/test_tree_node.py b/eckity/genetic_encodings/gp/tree/tests/test_tree_node.py
new file mode 100644
index 00000000..c30b55d7
--- /dev/null
+++ b/eckity/genetic_encodings/gp/tree/tests/test_tree_node.py
@@ -0,0 +1,61 @@
+import pytest
+
+from eckity.base.typed_functions import add2floats
+from eckity.base.untyped_functions import f_add
+from eckity.genetic_encodings.gp import FunctionNode, TerminalNode
+from eckity.genetic_encodings.gp.tree.utils import get_func_types
+
+
+@pytest.mark.parametrize(
+    "function, expected_types",
+    [
+        (add2floats, [float, float, float]),
+        (f_add, [None, None, None]),
+    ],
+)
+def test_get_func_types(function, expected_types):
+    """
+    Test that get_func_types returns the correct function types
+    """
+    func_types = get_func_types(function)
+    assert func_types == expected_types
+
+
+def test_missing_type_hints():
+    """
+    define some bad functions
+    not using pytest decorator because lambda expressions have no type hints
+    """
+
+    # missing return type hint
+    def add_no_return_type(x: int, y: int):
+        return x + y
+
+    # missing x type hint
+    def add_no_x_type(x, y: int) -> int:
+        return x + y
+
+    # missing y type hint
+    def add_no_y_type(x: int, y) -> int:
+        return x + y
+
+    for func in [add_no_return_type, add_no_x_type, add_no_y_type]:
+        with pytest.raises(ValueError) as excinfo:
+            FunctionNode(func)
+        assert "missing type hints" in str(excinfo)
+
+
+@pytest.mark.parametrize(
+    "node1, node2, expected",
+    [
+        (TerminalNode(1, int), TerminalNode(1, int), True),
+        (TerminalNode(1, int), TerminalNode(True, bool), False),
+        (TerminalNode(0, int), TerminalNode(False, bool), False),
+        (TerminalNode(1, int), TerminalNode(2, int), False),
+        (TerminalNode(1), TerminalNode(1), True),
+        (TerminalNode(1), TerminalNode(True), False),
+    ],
+)
+def test_eq(node1, node2, expected):
+    assert (node1 == node2) == expected
+    assert (node2 == node1) == expected
diff --git a/eckity/genetic_encodings/gp/tree/tests/test_tree_utils.py b/eckity/genetic_encodings/gp/tree/tests/test_tree_utils.py
new file mode 100644
index 00000000..7512e179
--- /dev/null
+++ b/eckity/genetic_encodings/gp/tree/tests/test_tree_utils.py
@@ -0,0 +1,39 @@
+import numpy as np
+import pytest
+from eckity.genetic_encodings.gp.tree.utils import (
+    create_terminal_set,
+    generate_args,
+)
+
+
+@pytest.mark.parametrize(
+    "X, typed, expected",
+    [
+        (
+            np.array([[4, 7, -7, -10], [7, -3, 3, -8], [8, -5, -3, -1]]),
+            False,
+            ["x0", "x1", "x2", "x3"],
+        ),
+        (
+            np.array([[4, 7, -7, -10], [7, -3, 3, -8], [8, -5, -3, -1]]),
+            True,
+            {"x0": int, "x1": int, "x2": int, "x3": int},
+        ),
+    ],
+)
+def test_create_terminal_set(X, typed, expected):
+    assert create_terminal_set(X, typed) == expected
+
+
+def test_generate_args():
+    X = np.array([[4, 7, -7, -10], [7, -3, 3, -8], [8, -5, -3, -1]])
+    expected = {
+        "x0": np.array([4, 7, 8]),
+        "x1": np.array([7, -3, -5]),
+        "x2": np.array([-7, 3, -3]),
+        "x3": np.array([-10, -8, -1]),
+    }
+    args = generate_args(X)
+    assert args.keys() == expected.keys()
+    for key in args:
+        assert np.array_equal(args[key], expected[key])
diff --git a/eckity/genetic_operators/crossovers/tests/test_subtree_crossover.py b/eckity/genetic_operators/crossovers/tests/test_subtree_crossover.py
new file mode 100644
index 00000000..6c52196e
--- /dev/null
+++ b/eckity/genetic_operators/crossovers/tests/test_subtree_crossover.py
@@ -0,0 +1,229 @@
+from typing import List
+
+import pytest
+
+from eckity.fitness import GPFitness
+from eckity.genetic_encodings.gp import FunctionNode, TerminalNode, Tree, TreeNode
+from eckity.genetic_operators import SubtreeCrossover
+
+
+# Define custom functions
+def typed_optional_inc(n: int, inc: bool) -> int:
+    return n + 1 if inc else n
+
+
+def untyped_optional_inc(n, inc):
+    return n + 1 if inc else n
+
+
+@pytest.mark.parametrize(
+    "individuals, subtrees, expected",
+    [
+        (
+            [
+                Tree(
+                    fitness=GPFitness(),
+                    function_set=[typed_optional_inc],
+                    terminal_set={True: bool, False: bool},
+                    erc_range=(1, 3),
+                    tree=[
+                        FunctionNode(typed_optional_inc),
+                        TerminalNode(1, int),
+                        TerminalNode(True, bool),
+                    ],
+                    root_type=int,
+                ),
+                Tree(
+                    fitness=GPFitness(),
+                    function_set=[typed_optional_inc],
+                    terminal_set={True: bool, False: bool},
+                    erc_range=(1, 3),
+                    tree=[
+                        FunctionNode(typed_optional_inc),
+                        TerminalNode(1, int),
+                        TerminalNode(False, bool),
+                    ],
+                    root_type=int,
+                ),
+            ],
+            [
+                [
+                    FunctionNode(typed_optional_inc),
+                    TerminalNode(1, int),
+                    TerminalNode(True, bool),
+                ],
+                [
+                    FunctionNode(typed_optional_inc),
+                    TerminalNode(1, int),
+                    TerminalNode(False, bool),
+                ],
+            ],
+            [
+                [
+                    FunctionNode(typed_optional_inc),
+                    TerminalNode(1, int),
+                    TerminalNode(False, bool),
+                ],
+                [
+                    FunctionNode(typed_optional_inc),
+                    TerminalNode(1, int),
+                    TerminalNode(True, bool),
+                ],
+            ]
+        ),
+        (
+            [
+                Tree(
+                    fitness=GPFitness(),
+                    function_set=[typed_optional_inc],
+                    terminal_set={True: bool, False: bool},
+                    erc_range=(1, 3),
+                    tree=[
+                        FunctionNode(typed_optional_inc),
+                        TerminalNode(1, int),
+                        TerminalNode(True, bool),
+                    ],
+                    root_type=int,
+                ),
+                Tree(
+                    fitness=GPFitness(),
+                    function_set=[typed_optional_inc],
+                    terminal_set={True: bool, False: bool},
+                    erc_range=(1, 3),
+                    tree=[
+                        FunctionNode(typed_optional_inc),
+                        TerminalNode(1, int),
+                        TerminalNode(False, bool),
+                    ],
+                    root_type=int,
+                ),
+            ],
+            [
+                [
+                    TerminalNode(True, bool),
+                ],
+                [
+                    TerminalNode(False, bool),
+                ],
+            ],
+            [
+                [
+                    FunctionNode(typed_optional_inc),
+                    TerminalNode(1, int),
+                    TerminalNode(False, bool),
+                ],
+                [
+                    FunctionNode(typed_optional_inc),
+                    TerminalNode(1, int),
+                    TerminalNode(True, bool),
+                ],
+            ]
+        ),
+        ######################################################################
+        # Untyped case
+        ######################################################################
+        (
+            [
+                Tree(
+                    fitness=GPFitness(),
+                    function_set=[untyped_optional_inc],
+                    terminal_set=[True, False],
+                    erc_range=(1, 3),
+                    tree=[
+                        FunctionNode(untyped_optional_inc),
+                        TerminalNode(1),
+                        TerminalNode(True),
+                    ],
+                ),
+                Tree(
+                    fitness=GPFitness(),
+                    function_set=[untyped_optional_inc],
+                    terminal_set=[True, False],
+                    erc_range=(1, 3),
+                    tree=[
+                        FunctionNode(untyped_optional_inc),
+                        TerminalNode(1),
+                        TerminalNode(False),
+                    ],
+                ),
+            ],
+            [
+                [
+                    FunctionNode(untyped_optional_inc),
+                    TerminalNode(1),
+                    TerminalNode(True),
+                ],
+                [
+                    FunctionNode(untyped_optional_inc),
+                    TerminalNode(1),
+                    TerminalNode(False),
+                ],
+            ],
+            [
+                [
+                    FunctionNode(untyped_optional_inc),
+                    TerminalNode(1),
+                    TerminalNode(False),
+                ],
+                [
+                    FunctionNode(untyped_optional_inc),
+                    TerminalNode(1),
+                    TerminalNode(True),
+                ],
+            ]
+        ),
+        (
+            [
+                Tree(
+                    fitness=GPFitness(),
+                    function_set=[untyped_optional_inc],
+                    terminal_set=[True, False],
+                    erc_range=(1, 3),
+                    tree=[
+                        FunctionNode(untyped_optional_inc),
+                        TerminalNode(1),
+                        TerminalNode(True),
+                    ],
+                ),
+                Tree(
+                    fitness=GPFitness(),
+                    function_set=[untyped_optional_inc],
+                    terminal_set=[True, False],
+                    erc_range=(1, 3),
+                    tree=[
+                        FunctionNode(untyped_optional_inc),
+                        TerminalNode(1),
+                        TerminalNode(False),
+                    ],
+                ),
+            ],
+            [
+                [
+                    TerminalNode(True),
+                ],
+                [
+                    TerminalNode(False),
+                ],
+            ],
+            [
+                [
+                    FunctionNode(untyped_optional_inc),
+                    TerminalNode(1),
+                    TerminalNode(False),
+                ],
+                [
+                    FunctionNode(untyped_optional_inc),
+                    TerminalNode(1),
+                    TerminalNode(True),
+                ],
+            ]
+        ),
+    ],
+)
+def test_swap_subtrees(
+    individuals: List[Tree],
+    subtrees: List[List[TreeNode]],
+    expected: List[Tree],
+):
+    SubtreeCrossover._swap_subtrees(individuals, subtrees)
+    assert [ind.tree for ind in individuals] == expected
diff --git a/eckity/genetic_operators/crossovers/tests/test_vector_crossover.py b/eckity/genetic_operators/crossovers/tests/test_vector_crossover.py
new file mode 100644
index 00000000..fcaf3033
--- /dev/null
+++ b/eckity/genetic_operators/crossovers/tests/test_vector_crossover.py
@@ -0,0 +1,26 @@
+import random
+from eckity.genetic_encodings.ga.int_vector import IntVector
+from eckity.genetic_operators.crossovers.vector_k_point_crossover import (
+    VectorKPointsCrossover,
+)
+from eckity.fitness.simple_fitness import SimpleFitness
+
+
+def test_vector_two_point_crossover():
+    random.seed(0)
+    length = 4
+    v1 = IntVector(SimpleFitness(), length, bounds=(1, 10))
+    v1.set_vector(list(range(1, 5)))
+    v2 = IntVector(SimpleFitness(), length, bounds=(1, 10))
+    v2.set_vector(list(range(5, 9)))
+
+    # random sample will return [2, 3]
+    expected_v1 = [5, 2, 7, 8]
+    expected_v2 = [1, 6, 3, 4]
+
+    crossover = VectorKPointsCrossover(k=2)
+    crossover.apply_operator([v1, v2])
+    assert v1.vector == expected_v1
+    assert v2.vector == expected_v2
+    assert v1.applied_operators == ["VectorKPointsCrossover"]
+    assert v2.applied_operators == ["VectorKPointsCrossover"]
diff --git a/eckity/genetic_operators/mutations/tests/test_subtree_mutation.py b/eckity/genetic_operators/mutations/tests/test_subtree_mutation.py
new file mode 100644
index 00000000..06260f0c
--- /dev/null
+++ b/eckity/genetic_operators/mutations/tests/test_subtree_mutation.py
@@ -0,0 +1,55 @@
+import pytest
+
+from eckity.fitness import GPFitness
+from eckity.genetic_encodings.gp import FunctionNode, TerminalNode, Tree
+from eckity.genetic_operators import SubtreeMutation
+
+
+def equal2ints(a: int, b: int) -> bool:
+    return a == b
+
+
+def inc_int(x: int) -> int:
+    return x + 1
+
+
+def f_equal(a, b):
+    return a == b
+
+
+@pytest.mark.parametrize(
+    "tree",
+    [
+        Tree(
+            fitness=GPFitness(),
+            function_set=[equal2ints, inc_int],
+            terminal_set={"x": int, "y": int},
+            tree=[
+                FunctionNode(function=equal2ints),
+                TerminalNode(value=1, node_type=int),
+                TerminalNode(value=2, node_type=int),
+            ],
+            root_type=bool,
+        ),
+        Tree(
+            fitness=GPFitness(),
+            function_set=[f_equal],
+            terminal_set=["x", "y"],
+            tree=[
+                FunctionNode(function=equal2ints),
+                TerminalNode(value=1),
+                TerminalNode(value=2),
+            ],
+        ),
+    ],
+)
+def test_subtree_mutation_success(tree):
+    subtree_mutation = SubtreeMutation(probability=1.0)
+
+    tree_copy = tree.clone()
+
+    # Perform subtree crossover
+    subtree_mutation.apply([tree_copy])
+
+    assert tree.root == tree_copy.root
+    assert tree != tree_copy
diff --git a/eckity/genetic_operators/mutations/tests/test_mutations.py b/eckity/genetic_operators/mutations/tests/test_vector_mutation.py
similarity index 100%
rename from eckity/genetic_operators/mutations/tests/test_mutations.py
rename to eckity/genetic_operators/mutations/tests/test_vector_mutation.py
diff --git a/eckity/sklearn_compatible/tests/test_sklearn_compatibility.py b/eckity/sklearn_compatible/tests/test_sklearn_compatibility.py
new file mode 100644
index 00000000..cd5c51b0
--- /dev/null
+++ b/eckity/sklearn_compatible/tests/test_sklearn_compatibility.py
@@ -0,0 +1,37 @@
+"""
+Tests for EC-KitY are currently under active development.
+For now, this class is a placeholder for future tests to be added.
+"""
+
+import pytest
+from sklearn.exceptions import NotFittedError
+
+from eckity import Subpopulation
+from eckity.algorithms import SimpleEvolution
+from eckity.base.typed_functions import sqrt_float
+from eckity.creators import GrowCreator
+from eckity.genetic_operators import SubtreeCrossover, SubtreeMutation
+from eckity.sklearn_compatible import ClassificationEvaluator, SKClassifier
+
+
+class TestSklearnCompatibility:
+    algo = SimpleEvolution(
+        Subpopulation(
+            evaluator=ClassificationEvaluator(),
+            creators=GrowCreator(
+                init_depth=(1, 1),
+                function_set=[sqrt_float],
+                terminal_set={"x": float},
+            ),
+            operators_sequence=[
+                SubtreeCrossover(),
+                SubtreeMutation(),
+            ],
+        )
+    )
+    clf = SKClassifier(algo)
+
+    def test_predict_without_fit(self):
+        X = [1, 2, 3]
+        with pytest.raises(NotFittedError):
+            self.clf.predict(X)
diff --git a/tests/test_reproducibility.py b/tests/test_reproducibility.py
index 435a238b..e17fe77f 100644
--- a/tests/test_reproducibility.py
+++ b/tests/test_reproducibility.py
@@ -1,12 +1,15 @@
-import random
-import numpy as np
 import itertools
-
+import random
 from time import sleep
 
+import numpy as np
+
+from eckity import Subpopulation
 from eckity.algorithms import SimpleEvolution
+from eckity.base.untyped_functions import f_add, f_div, f_mul, f_sub
+from eckity.creators import FullCreator
 from eckity.evaluators import SimpleIndividualEvaluator
-from eckity import Subpopulation
+from eckity.genetic_operators import SubtreeCrossover, SubtreeMutation
 
 
 class RandomIndividualEvaluator(SimpleIndividualEvaluator):
@@ -16,9 +19,17 @@ def evaluate_individual(self, ind):
 
 class TestReproducibility:
     def test_reproducibility(self):
-        # creates a Tree-GP experiment by default
+        # creates a Tree-GP algorithm
         algo1 = SimpleEvolution(
-            Subpopulation(RandomIndividualEvaluator(), population_size=10),
+            Subpopulation(
+                RandomIndividualEvaluator(),
+                population_size=10,
+                creators=FullCreator(
+                    function_set=[f_add, f_sub, f_mul, f_div],
+                    terminal_set=["x", "y", "z"],
+                ),
+                operators_sequence=[SubtreeCrossover(), SubtreeMutation()],
+            ),
             random_seed=42,
             max_generation=10,
         )
@@ -27,7 +38,15 @@ def test_reproducibility(self):
 
         # creates a Tree-GP experiment by default
         algo2 = SimpleEvolution(
-            Subpopulation(RandomIndividualEvaluator(), population_size=10),
+            Subpopulation(
+                RandomIndividualEvaluator(),
+                population_size=10,
+                creators=FullCreator(
+                    function_set=[f_add, f_sub, f_mul, f_div],
+                    terminal_set=["x", "y", "z"],
+                ),
+                operators_sequence=[SubtreeCrossover(), SubtreeMutation()],
+            ),
             random_seed=42,
             max_generation=10,
         )
@@ -46,7 +65,14 @@ def test_random_seeds(self):
         seeds = []
         for i in range(n_reps):
             seed = SimpleEvolution(
-                Subpopulation(RandomIndividualEvaluator()),
+                Subpopulation(
+                    RandomIndividualEvaluator(),
+                    creators=FullCreator(
+                        function_set=[f_add, f_sub, f_mul, f_div],
+                        terminal_set=["x", "y", "z"],
+                    ),
+                    operators_sequence=[SubtreeCrossover()],
+                ),
             ).random_seed
             seeds.append(seed)
 
diff --git a/tests/test_sanity.py b/tests/test_sanity.py
new file mode 100644
index 00000000..e9f107f8
--- /dev/null
+++ b/tests/test_sanity.py
@@ -0,0 +1,50 @@
+"""
+Sanity tests for some basic examples.
+"""
+
+import pytest
+import subprocess
+
+
+def run_example(example_path: str):
+    result = subprocess.run(
+        ["python", "-m", example_path], capture_output=True
+    )
+    return result
+
+
+@pytest.mark.parametrize(
+    "module_path, n_reps, expected_fitness, higher_is_better",
+    [
+        ("examples.treegp.basic_mode.multiplexer", 3, 0.7, True),
+    ],
+)
+def test_example(
+    module_path: str,
+    n_reps: int,
+    expected_fitness: float,
+    higher_is_better: bool,
+):
+    fitness_scores = []
+    for _ in range(n_reps):
+        # run the example
+        result = subprocess.run(
+            ["python", "-m", module_path],
+            capture_output=True,
+        )
+        assert result.returncode == 0
+
+        # parse best pure fitness:
+        output = result.stdout.decode("utf-8")
+        lines = output.split("\n")
+        for line in lines:
+            if "best pure fitness:" in line:
+                best_pure_fitness = float(line.split(":")[-1].strip())
+                fitness_scores.append(best_pure_fitness)
+
+    mean_fitness = sum(fitness_scores) / n_reps
+
+    if higher_is_better:
+        assert mean_fitness >= expected_fitness
+    else:
+        assert mean_fitness <= expected_fitness
diff --git a/tests/test_sklearn_compatibility.py b/tests/test_sklearn_compatibility.py
deleted file mode 100644
index 0c87fecb..00000000
--- a/tests/test_sklearn_compatibility.py
+++ /dev/null
@@ -1,21 +0,0 @@
-"""
-Tests for EC-KitY are currently under active development.
-For now, this class is a placeholder for future tests to be added.
-"""
-
-import pytest
-from sklearn.exceptions import NotFittedError
-
-from eckity.algorithms import SimpleEvolution
-from eckity.sklearn_compatible import ClassificationEvaluator
-from eckity.sklearn_compatible import SKClassifier
-from eckity.subpopulation import Subpopulation
-
-
-class TestSklearnCompatibility:
-    def test_predict_without_fit(self):
-        algo = SimpleEvolution(Subpopulation(evaluator=ClassificationEvaluator()))
-        clf = SKClassifier(algo)
-        X = [1, 2, 3]
-        with pytest.raises(NotFittedError):
-            clf.predict(X)
diff --git a/tests/test_subpopulation.py b/tests/test_subpopulation.py
index 06dd2aba..3690150b 100644
--- a/tests/test_subpopulation.py
+++ b/tests/test_subpopulation.py
@@ -1,6 +1,8 @@
 import logging
 from eckity.subpopulation import Subpopulation
-from examples.treegp.non_sklearn_mode.symbolic_regression.sym_reg_evaluator import (
+from eckity.creators import FullCreator
+from eckity.genetic_operators import IdentityTransformation
+from examples.treegp.basic_mode.symbolic_regression import (
     SymbolicRegressionEvaluator,
 )
 
@@ -9,6 +11,11 @@ def test_zero_elites_warnning(caplog):
     with caplog.at_level(logging.WARNING):
         Subpopulation(
             SymbolicRegressionEvaluator(),
+            creators=FullCreator(
+                function_set=[lambda x: x],
+                terminal_set=['x'],
+            ),
+            operators_sequence=[IdentityTransformation()],
             population_size=10,
             elitism_rate=0.01,
         )