Nguyen-5
The problem
Symbolic regression: find a formula for y in x from the points of
y = sin(x²) cos(x) − 1,
at 20 values of x drawn at random from [−1, 1]. It's the fifth of the twelve Nguyen problems (Uy et al. 2011), the first with sines and cosines. The formula is found when it's recovered exactly: an error at the level of rounding on the 20 training points and on 100 test points from [−1, 1] that the search never sees.
The problem's name, target, sampling and function set are as McDermott et al. (2012) restate them, not yet checked against the paper. Its points come from a fixed seed of genoxide's random stream; another library's differ.
The Python version builds the same islands with gx.gp and evaluates the trees in Rust
(problem.regression()): it prints the same output and writes the same trace.
What makes it hard
The function set, Koza's (addition, subtraction, multiplication, protected division, sine, cosine, the exponential and a protected logarithm, and x), has no constants. The shape sin(x²) cos(x) is a tree of 6 nodes, but the − 1 has to be built too, from x / x (a protected division that is 1 everywhere) or cos(x − x), and a search on the tree's own error has to find both at once. Without the − 1, the best trees fit the curve shifted by one, and their error stays large.
Representation
A tree of a genetic program (gp::Tree), of the primitives of gp::regression::problems::Nguyen5:
add, sub, mul, the protected division pdiv, sin, cos, exp, the protected logarithm
plog and x. gp::Gp sets Koza's limits and initialization: depth at most 17, at most 1024
nodes, and ramped half-and-half with depths 2 to 6.
The fitness is the root mean squared error on the 20 points after linear scaling
(gp::regression::Regression, with scaling on, its default): the error of a + b × tree, with a and
b fitted to the points by least squares for each tree (Keijzer 2003, Improving symbolic regression
with interval arithmetic and linear scaling, EuroGP 2003: 70-82). The search then only needs the
shape, sin(x²) cos(x); the fitted a = −1 and b = 1 give the rest. Without scaling, the same search
recovered the formula in 8 of 20 runs I tried; with it, in 94 of 100. A tree with a value that isn't
finite is invalid.
Algorithm
The same as the Nguyen-1 page and Koza's quartic: eight genetic algorithms of 500 trees each, as islands in a ring that pass their two best trees on every 10 generations, each with tournaments of 7, subtree crossover at a rate of 0.9 and subtree mutation at a rate of 0.1. The run stops at an RMSE of 10⁻¹⁰ times the standard deviation of the 20 values, or after 200 generations.
Output
The first two lines give the problem and the setting. Then come how the run stopped, after how many
generations and evaluations, the RMSE of the best tree on the training and the test points after
its scaling, and the tree with its scaling, as a + b * (tree) (Regression::display).
Good results
The optimum is the formula itself, recovered exactly. The run of output.txt found it after 2
generations and 10,670 evaluations, as
−1 + 1 × (mul(cos(x), sin(mul(mul(x, x), pdiv(x, x))))),
that is cos(x) · sin(x · x · 1) − 1, in 11 nodes, with an RMSE of 0 on the training and the test points: the scaling's intercept and slope are exactly −1 and 1.
Over 100 runs (the islands' seeds 100 s + 0 to 7 for s = 1 to 100, the first being output.txt's,
with the same data), 94 recovered the formula, after 2 generations and 10,600 evaluations in the
median, and 7 generations at most. The other 6 ended after 200 generations on near misses of 166 to
283 nodes, with test errors from 3 × 10⁻⁵ to 0.2.
Known optimum: sin(x²) cos(x) − 1, exactly (an RMSE of 0, up to rounding)
Source: examples/nguyen_5
Interactive run: tachsin.gr/projects/genoxide/examples/nguyen-5
cargo run --release --example nguyen_5
//! Nguyen-5: find the formula sin(x²) cos(x) − 1 from 20 points of it, by genetic programming.
//!
//! The fifth of Nguyen's twelve symbolic regression problems (Uy et al. 2011): trees of Koza's
//! functions (+, −, ×, protected division, sin, cos, exp and a protected logarithm) and the
//! variable x, evolved by a genetic algorithm with subtree crossover and mutation, fitted to the
//! root mean squared error on the points after linear scaling, which supplies the − 1. The run stops at exact recovery: an error at the level of
//! rounding, on the 20 training points and on 100 test points in [−1, 1], with the expression
//! printed.
//!
//! With `GENOXIDE_TRACE=<file>`, it also writes a trace of its run for the plot on the example's
//! page, with `trace.rs`.
//!
//! ```text
//! cargo run --release --example nguyen_5
//! ```
mod trace;
use genoxide::gp::regression::problems::{Nguyen5, RegressionProblem};
use genoxide::gp::{Gp, SubtreeCrossover, SubtreeMutation, Tree};
use genoxide::prelude::*;
// the islands, their trees, and the generations between migrations
const ISLANDS: u64 = 8;
const POPULATION: usize = 500;
const INTERVAL: u64 = 10;
fn main() -> Result<()> {
// Nguyen's function set and sampling: 20 training points uniform in [-1, 1], from a fixed
// seed, and 100 test points from another; the RMSE after linear scaling: a + b * tree, with
// a and b fitted by least squares, since the set has no constant for the - 1
let problem = Nguyen5::new();
let regression = problem.regression().clone();
let dataset = problem.dataset();
let (training, test) = (dataset.training(), dataset.test().expect("a test sample"));
// exact recovery: an error of at most 1e-10 of the values' standard deviation
let tolerance = 1e-10 * training.deviation();
println!("Nguyen-5: sin(x^2) cos(x) - 1 from 20 points in [-1, 1]");
println!(
"{ISLANDS} islands of {POPULATION} trees, until the RMSE is at most {tolerance:.2e}\n"
);
// Koza's limits and initialization: depth 17, ramped half-and-half of depths 2 to 6
let gp = Gp::builder(problem.primitives().clone()).build()?;
let set = gp.primitives().clone();
let islands = (0..ISLANDS)
.map(|island| {
let seed = 100 + island;
// Koza's even division among the depths and methods, without duplicates
let initial =
gp.ramped_half_and_half(POPULATION, &mut StreamRng::seed_from_u64(seed))?;
Ga::builder(gp.clone())
.population_size(POPULATION)
.initial_genomes(initial)
.select(Tournament::new(7)?)
.crossover(SubtreeCrossover::new())
.mutate(SubtreeMutation::new())
.crossover_rate(0.9)
.mutation_rate(0.1)
.minimize()
.seed(seed)
.build()
})
.collect::<Result<Vec<_>>>()?;
let islands = Islands::builder(islands)
.interval(INTERVAL)
.migrants(2)
.build()?;
// with GENOXIDE_TRACE=<file>, a trace of the run for the plot on the example's page
let mut trace = trace::Trace::from_env("nguyen_5", &problem, ®ression, -1.0..=1.0);
let outcome = Engine::new(islands, regression.clone())
.stop_when(Stop::target(tolerance).or(Stop::generations(200)))
.on_generation(|snapshot| trace.record(snapshot))
.run()?;
let best: &Tree = outcome.best_genome();
let rmse = |sample| regression.error(best, sample).unwrap_or(f64::NAN);
println!(
"{:?} after {} generations and {} evaluations",
outcome.stop_reason(),
outcome.generations(),
outcome.evaluations()
);
println!("RMSE on the 20 training points: {:.2e}", rmse(training));
println!("RMSE on the 100 test points: {:.2e}", rmse(test));
println!(
"\nthe expression, {} nodes of depth {}, and its scaling:\n{}",
best.len(),
best.depth(&set),
regression.display(best)
);
trace.write();
Ok(())
}
python examples/nguyen_5/main.py
"""Nguyen-5: find the formula sin(x^2) cos(x) - 1 from 20 points of it, by genetic
programming.
The fifth of Nguyen's twelve symbolic regression problems (Uy et al. 2011): trees of Koza's
functions (+, -, x, protected division, sin, cos, exp and a protected logarithm) and the
variable x, evolved by a genetic algorithm with subtree crossover and mutation, fitted to the
root mean squared error on the points after linear scaling, which supplies the - 1. The run
stops at exact recovery: an error at the level of rounding, on the 20 training points and on 100
test points in [-1, 1], with the expression printed.
The trees and their errors are computed in Rust (``gx.gp``), with genoxide's portable math, so
the run is the Rust example's, to the bit, on every platform.
With ``GENOXIDE_TRACE=<file>``, it also writes a trace of its run for the plot on the example's
page, with trace.py.
python examples/nguyen_5/main.py
"""
import genoxide as gx
from trace import Trace
# the islands, their trees, and the generations between migrations
ISLANDS = 8
POPULATION = 500
INTERVAL = 10
def scientific(value):
"""``value`` with 2 decimals and an exponent, as Rust's ``{:.2e}`` writes it: 8.76e-11,
0.00e0."""
if value != value:
return "NaN"
mantissa, exponent = f"{value:.2e}".split("e")
return f"{mantissa}e{int(exponent)}"
# Nguyen's function set and sampling: 20 training points uniform in [-1, 1], from a fixed seed,
# and 100 test points from another; the RMSE after linear scaling: a + b * tree, with a and b
# fitted by least squares, since the set has no constant for the - 1
problem = gx.gp.regression.problems.Nguyen5()
regression = problem.regression()
dataset = problem.dataset()
training, test = dataset.training, dataset.test
# exact recovery: an error of at most 1e-10 of the values' standard deviation
tolerance = 1e-10 * training.deviation()
print("Nguyen-5: sin(x^2) cos(x) - 1 from 20 points in [-1, 1]")
print(f"{ISLANDS} islands of {POPULATION} trees, until the RMSE is at most {scientific(tolerance)}")
print()
# Koza's limits and initialization: depth 17, ramped half-and-half of depths 2 to 6
gp = gx.gp.Gp(problem.primitives())
islands = []
for island in range(ISLANDS):
seed = 100 + island
islands.append(
gx.Ga(
gp,
population_size=POPULATION,
# Koza's even division among the depths and methods, without duplicates
initial_genomes=gp.ramped_half_and_half(POPULATION, seed),
select=gx.Tournament(7),
crossover=gx.gp.SubtreeCrossover(),
mutation=gx.gp.SubtreeMutation(),
crossover_rate=0.9,
mutation_rate=0.1,
objective="minimize",
seed=seed,
)
)
islands = gx.Islands(islands, interval=INTERVAL, migrants=2)
# with GENOXIDE_TRACE=<file>, a trace of the run for the plot on the example's page
trace = Trace("nguyen_5", problem, regression, (-1.0, 1.0))
result = islands.run(
regression, target=tolerance, generations=200, on_generation=trace.on_generation
)
best = result.best_genome
def rmse(sample):
error = regression.error(best, sample)
return float("nan") if error is None else error
print(
f"{result.stop_reason.capitalize()} after {result.generations} generations and "
f"{result.evaluations} evaluations"
)
print(f"RMSE on the 20 training points: {scientific(rmse(training))}")
print(f"RMSE on the 100 test points: {scientific(rmse(test))}")
print()
print(f"the expression, {len(best)} nodes of depth {best.depth}, and its scaling:")
print(regression.display(best))
trace.write()
What it prints, from a seeded run:
Nguyen-5: sin(x^2) cos(x) - 1 from 20 points in [-1, 1]
8 islands of 500 trees, until the RMSE is at most 1.66e-11
Target after 2 generations and 10670 evaluations
RMSE on the 20 training points: 0.00e0
RMSE on the 100 test points: 0.00e0
the expression, 11 nodes of depth 4, and its scaling:
-1 + 1 * (mul(cos(x), sin(mul(mul(x, x), pdiv(x, x)))))