Making a start on Hilbert Spaces and matrix calculations using CPU

Mike's Notes

Working notes on finding a way for Pipi to use genetic algorithms with Hilbert Spaces and matrix calculations. A desktop exercise with Google at the moment. Computational experiments on CPU to come. Some generated Python and CFML sample code below to start playing with. 

It looks promising; it might help address limitations of other methods for using fitness landscapes. The output curves do look better.

This might void the need to use a Cerebras wafer for these fitness landscape calculations.

The next step is to watch free maths lectures on YouTube, then do a lot of drawings and process diagrams. Then play with code and see what happens. 😎

The maths formulas were in LaTeX format and were converted to non-LaTeX format using the free online LatexConverter.

Note

Nothing has been verified or tested. 😎

Resources

References

  • Reference

Repository

  • Home > Ajabbi Research > Library >
  • Home > Handbook > 

Last Updated

03/09/2026

Making a start on Hilbert Spaces and matrix calculations using CPU

By: Mike Peters & Google Search AI Mode
On a Sandy Beach: 30/08/2026

Mike Peters: Mike is the inventor and architect of Pipi and the founder of Ajabbi.

Google Search AI Mode: Free and surprisingly good.

...

Hilbert Spaces

Hilbert spaces use matrix calculations extensively, especially when the space is finite-dimensional or when dealing with linear operators. Here is how they connect:

  1. Finite-Dimensional Spaces
    • Direct Equivalence: An (n)-dimensional Hilbert space is isomorphic to (ℂⁿ) (complex numbers) or (ℝⁿ) (real numbers).
    • Vectors: Elements in the space are represented as column vectors.
    • Operators: Linear operators acting on the space are represented as (n × n) matrices.
    • Inner Product: The inner product is calculated using matrix multiplication (specifically, the conjugate transpose of a vector multiplied by another vector).
  2. Infinite-Dimensional Spaces
    • Infinite Matrices: In infinite-dimensional Hilbert spaces (like L²), the space of square-integrable functions, operators can be represented as infinite matrices.
    • Basis Dependence: By choosing an orthonormal basis, any linear operator can be written as a matrix with infinitely many rows and columns.
  3. Quantum Mechanics (Matrix Mechanics)
    • State Vectors: Quantum states are vectors in a Hilbert space.
    • Observables: Physical properties (like position, momentum, or spin) are Hermitian operators.
    • Heisenberg's Formulation: Werner Heisenberg's formulation of quantum mechanics is literally called "matrix mechanics" because it uses matrix calculations to solve Hilbert space problems.

In machine learning, Hilbert spaces—specifically Reproducible Kernel Hilbert Spaces (RKHS)—rely heavily on matrix calculations to make complex, non-linear data linearly separable.

Instead of explicitly mapping data points to an infinite-dimensional space, ML algorithms use the Kernel Trick to compute everything via an N × N matrix, where N is the number of data samples.

The Gram (Kernel) Matrix

The core matrix calculation in RKHS is the Gram Matrix (or Kernel Matrix), denoted as K.

  • Formula: (Kᵢⱼ = k(xᵢ, xⱼ) = ⟨ Φ(xᵢ), Φ(xⱼ) ⟩_{H})
  • Meaning: Each entry i,j is the inner product of data points (xᵢ) and (xⱼ) in the Hilbert space (H).
  • Properties: It is always a symmetric and positive semi-definite (PSD) matrix.

Core Matrix Operations in RKHS

Machine learning algorithms manipulate this kernel matrix using standard linear algebra:

  • Matrix Inversion: Used in Kernel Ridge Regression to solve for weights: (Kv = λ v).
  • Eigenvalue Decomposition: Used in Kernel PCA to find non-linear principal components by solving (Kv = λ v).
  • Quadratic Programming: Used in Support Vector Machines (SVMs) to maximise margins using matrix-vector products involving (K).

Computational Bottlenecks

Because (K) scales with the number of data points, scaling Hilbert space methods requires specialised matrix approximations:

  • Nyström Approximation: Factorises a low-rank sub-matrix to approximate the massive N × N matrix.
  • Random Fourier Features: Avoids the exact matrix calculation entirely by mapping data into a randomised, lower-dimensional space.

Graph image

Genetic Algorithms

Combining genetic algorithms (GAs) with Reproducing Kernel Hilbert Spaces (RKHS) is entirely possible, highly effective, and actively used in machine learning research.

In this hybrid framework, the fitness landscape is modelled as a continuous function belonging to an RKHS. This approach addresses a major limitation of traditional GAs: evaluating an individual's true fitness can be computationally expensive (e.g., running a physics simulation or a wet-lab experiment). The RKHS serves as an intelligent surrogate model.

Graph image

Core Mathematical Engines At Work

To make this framework highly effective, three core components work in tandem:

  • The Representer Theorem: This is the mathematical backbone of your RKHS surrogate. It guarantees that your continuous fitness approximation (̂{f}(x)) can be written as a finite linear combination of kernel evaluations centred at your sampled points:

(ᶠ(x)=∑ ᵢ₌₁ᴺα ᵢK(x,xᵢ))

  • This scales the search space completely independently of its true dimensionality, reducing prediction cost to a simple (O(N)) vector dot product.
  • Informed Exploration vs. Exploitation: Instead of letting the GA evaluate individuals blindly on the surrogate, you can use the RKHS variance (uncertainty) to construct an Acquisition Function (like Expected Improvement or Upper Confidence Bound). The GA then maximises this acquisition function rather than the raw estimated fitness, forcing the algorithm to intelligently search unmapped areas of the landscape.
  • Gram Matrix Regularisation: As new true data points are evaluated, they enter the Gram Matrix. To prevent numerical instability or overfitting as (N) grows, a small regularisation ridge ((λ I)) is maintained, smoothly adjusting the landscape's rigidity without rebuilding the framework from scratch.

How the Combination Works

[GA Population] ---> [Evaluate on RKHS Surrogate] ---> [Select & Crossover] ^ | |__________________ [Update RKHS with True Data] _________|


  1. The RKHS as the Fitness Landscape: You treat the unknown fitness landscape (f(x)) as a smooth function within an RKHS. By evaluating a small set of initial points, you use Kernel Ridge Regression or Kriging (Gaussian Process Regression) to build a continuous, analytical approximation of the landscape.
  2. Genetic Search on the Surrogate: The GA searches this approximated RKHS landscape. Because evaluating the RKHS kernel matrix is mathematically cheap (O(N)) for a new prediction, the GA can evolve through thousands of generations in seconds.
  3. Adaptive Sampling (Bayesian Optimisation): The best candidates found by the GA are evaluated using the true, expensive fitness function. These new data points are added back into the Gram matrix, updating the RKHS landscape to make it more accurate where it matters most.

Key Advantages of This Hybrid Approach

  • Global Search in Infinite Dimensions: Traditional GAs struggle with high-dimensional continuous spaces. By mapping the search space into an RKHS via a kernel (like the Gaussian RBF kernel), the GA implicitly searches a highly expressive, infinite-dimensional feature space.
  • Guaranteed Smoothness: Functions in an RKHS are mathematically smooth and bounded [6]. This eliminates extreme, discontinuous noise in the fitness landscape, allowing the GA to find gradient-like paths toward global optima without getting trapped in microscopic local minima.
  • Exact Kernel Alignment: You can use a GA to optimise the hyperparameters of the RKHS kernel itself (e.g., finding the optimal length-scale (σ ) of an RBF kernel) to perfectly align the geometry of the Hilbert space with your actual problem landscape.

Real-World Use Cases

  • Protein Design & Biochemistry: Designing new molecules where the "fitness" is binding affinity. The RKHS models the chemical landscape based on known proteins, and the GA mutates amino acid sequences to find optimal binders.
  • Hyperparameter Optimisation: Optimising complex neural network architectures where the evaluation landscape is highly non-linear and expensive to compute.
  • Robotics Control: Evolving locomotion parameters for robotics where physical testing is slow or dangerous.

Implementation Challenges

  • The Matrix Bottleneck Returns: As your GA discovers more points and updates the RKHS landscape, your Gram matrix (K) grows. Inverting this matrix to update the landscape scales at O(N³). You will need to use Nyström approximations or sparse online kernel methods if your total evaluation points exceed 10,000.
  • Exploration vs. Exploitation: GAs are aggressive exploiters. If your RKHS surrogate model is inaccurate in an unexplored region, the GA might get stuck optimising a "fake" peak. You must include an uncertainty metric (like the variance in a Gaussian Process) to force the GA to explore unknown areas.

Kernel Ridge Regression or Kriging (Gaussian Process Regression)

Kernel Ridge Regression (KRR) and Kriging (Gaussian Process Regression, or GPR) are both kernel-based methods for learning non-linear functions from data, but they differ fundamentally in their underlying statistical framework, objective functions, and outputs.

Key Differences

  • Core Approach: KRR minimises a regularised mean-squared error loss function in a Reproducing Kernel Hilbert Space (RKHS). GPR uses a probabilistic (Bayesian) approach, defining a Gaussian process prior over functions and updating it with a likelihood function based on observed data.
  • Uncertainty Estimation: KRR outputs point predictions only. GPR (and Kriging) naturally quantifies uncertainty, providing full posterior distributions, variance estimates, and confidence intervals.
  • Hyperparameter Optimisation: KRR typically optimises kernel parameters using grid search with cross-validation on a loss function. GPR optimises hyperparameters via gradient ascent on the marginal likelihood.
  • Terminology & Origin: Kriging originated in geostatistics (mining) to find the Best Linear Unbiased Predictor (BLUP), whereas GPR stems from machine learning and stochastic processes. Mathematically, standard Kriging is equivalent to GPR under matching covariance and prior assumptions

When to Use Which

  • Use Kernel Ridge Regression when: You need a fast, deterministic, regularised non-linear regressor and do not require predictive uncertainty bounds.
  • Use Kriging / Gaussian Process Regression when: You need confidence intervals for predictions, want to optimise hyperparameters via marginal likelihood, or are modelling spatial/geostatistical data.

Python Code (not tested)

Here is a clean, modular Python template using scikit-learn to stitch a Genetic Algorithm loop to a Gaussian Process (RKHS surrogate) framework.

import numpy as np
from sklearn.gaussian_process import GaussianProcessRegressor
from sklearn.gaussian_process.kernels import Matern

# --- 1. CONFIGURATION & SIMULATION BACKEND ---
BOUNDS = np.array([[-5.0, 5.0], [-5.0, 5.0]])  # 2D search space boundaries
N_DIM = BOUNDS.shape[0]
POP_SIZE = 20
GENERATIONS = 10
SURROGATE_MAX_ITER = 5   # Active learning/Bayesian optimization loops

def true_expensive_fitness(x):
    """
    Represents your expensive simulation or physical experiment.
    Expects a 1D array of shape (N_DIM,). Returns a scalar.
    """
    # Example: Multi-modal Ackley function (minimization converted to maximization)
    return -1 * (20 * np.exp(-0.2 * np.sqrt(0.5 * np.sum(x**2))) + 
                 np.exp(0.5 * np.sum(np.cos(2 * np.pi * x))) - 20 - np.e)

# --- 2. SURROGATE MODEL (RKHS / GAUSSIAN PROCESS) ---
# Using Matérn kernel with explicit noise handling (alpha) as regularization
kernel = Matern(nu=2.5)
surrogate = GaussianProcessRegressor(kernel=kernel, alpha=1e-6, normalize_y=True, random_state=42)

# --- 3. GENETIC ALGORITHM CORE ENGINE (OPERATING ON SURROGATE) ---
def init_population(pop_size, bounds):
    return np.random.uniform(bounds[:, 0], bounds[:, 1], size=(pop_size, len(bounds)))

def crossover(parent1, parent2):
    # Blend crossover (BLX-alpha)
    alpha = 0.5
    gamma = (1 + 2 * alpha) * np.random.random(size=N_DIM) - alpha
    return parent1 + gamma * (parent2 - parent1)

def mutate(individual, bounds, rate=0.2, scale=0.5):
    for i in range(len(individual)):
        if np.random.random() < rate:
            individual[i] += np.random.normal(0, scale)
            individual[i] = np.clip(individual[i], bounds[i, 0], bounds[i, 1])
    return individual

def run_ga_on_surrogate(model, bounds, pop_size, generations):
    """
    Standard GA that queries the cheap RKHS surrogate model instead of the true fitness.
    """
    population = init_population(pop_size, bounds)
    
    for _ in range(generations):
        # O(N) evaluation using the trained kernel model
        # predict() returns (mean, std). We maximize the predicted mean fitness.
        fitnesses, _ = model.predict(population, return_std=True)
        
        # Rank-based selection
        idx = np.argsort(fitnesses)[::-1]
        population = population[idx]
        
        # Breed next generation
        next_gen = list(population[:2])  # Elitist preservation
        while len(next_gen) < pop_size:
            p1, p2 = population[np.random.randint(0, 5)], population[np.random.randint(0, 5)]
            child = crossover(p1, p2)
            child = mutate(child, bounds)
            next_gen.append(child)
            
        population = np.array(next_gen)
        
    # Return best candidate found by GA in this surrogate landscape
    best_idx = np.argmax(model.predict(population, return_std=False))
    return population[best_idx]

# --- 4. THE ACTIVE LEARNING HYBRID LOOP ---
# Step A: Initialize with a small Latin Hypercube or random experimental design
X_train = init_population(pop_size=10, bounds=BOUNDS)
y_train = np.array([true_expensive_fitness(x) for x in X_train])

print("Starting Hybrid GA-RKHS Optimization Engine...\n")

for loop in range(SURROGATE_MAX_ITER):
    print(f"--- Iteration {loop + 1}/{SURROGATE_MAX_ITER} ---")
    print(f"Gram Matrix Size: {X_train.shape[0]} points")
    
    # Step B: Fit/Update RKHS surrogate model
    surrogate.fit(X_train, y_train)
    
    # Step C: Let the GA aggressively explore the surrogate landscape
    # This costs thousands of operations but takes milliseconds
    proposed_candidate = run_ga_on_surrogate(surrogate, BOUNDS, POP_SIZE, GENERATIONS)
    
    # Step D: Adaptive Sampling (Evaluate the candidate on the true simulation)
    true_score = true_expensive_fitness(proposed_candidate)
    print(f"GA Proposed Vector: {np.round(proposed_candidate, 4)}")
    print(f"True Fitness Score: {true_score:.4f}")
    
    # Step E: Augment data (Updates the Gram matrix for the next iteration)
    X_train = np.vstack([X_train, proposed_candidate])
    y_train = np.append(y_train, true_score)
    print(f"Current Best Score Found: {np.max(y_train):.4f}\n")

print("Optimization Complete.")

Enhancements 

  • Regularisation (Ridge): The alpha=1e-6 in GaussianProcessRegressor serves as the λI identity-matrix offset. This prevents numerical inversion crashes if the GA proposes points structurally close to existing data.
  • Normalisation Engine: normalize_y=True transforms unstable fitness functions to a standard normal distribution mapping, stabilising kernel hyperparameter tuning during .fit().
  • Scalability Hand-off: run_ga_on_surrogate can run for massive populations or generations; it never interacts with true_expensive_fitness, protecting your heavy computation budget.

CFML Code (not tested)

<cfscript>
/**
 * CFML Hybrid GA-RKHS Optimization Engine
 */

// --- 1. CONFIGURATION & CONFIG MATRIX ---
variables.bounds = [ {min: -5.0, max: 5.0}, {min: -5.0, max: 5.0} ]; // 2D Search Space
variables.popSize = 20;
variables.generations = 10;
variables.surrogateMaxIter = 5;

/**
 * Represents your expensive simulation or backend physical experiment.
 * Minimisation problem converted to Maximisation.
 */
public numeric function trueExpensiveFitness(required array x) {
    // Standard multi-modal benchmark proxy
    var sumSq = 0;
    var sumCos = 0;
    var n = arrayLen(arguments.x);
    
    for (var i = 1; i <= n; i++) {
        sumSq += arguments.x[i] ^ 2;
        sumCos += cos(2 * pi() * arguments.x[i]);
    }
    
    var term1 = -20.0 * exp(-0.2 * sqrt(0.5 * sumSq));
    var term2 = -exp(0.5 * sumCos);
    return -1 * (term1 + term2 - 20 - exp(1));
}

// --- 2. SURROGATE ENGINE (Continuous Analytical Proxy) ---
/**
 * Approximates f(x) using existing Gram dataset via standard RBF/IDW Kernel.
 */
public numeric function predictSurrogateFitness(required array x, required array xTrain, required array yTrain) {
    var totalWeight = 0;
    var weightedSum = 0;
    var p = 2; // Power parameter for spatial continuity
    var regularizationRidge = 1e-6; // Prevents division by zero on exact matches
    
    for (var i = 1; i <= arrayLen(arguments.xTrain); i++) {
        var distSq = 0;
        for (var d = 1; d <= arrayLen(arguments.x); d++) {
            distSq += (arguments.x[d] - arguments.xTrain[i][d]) ^ 2;
        }
        var dist = sqrt(distSq) + regularizationRidge;
        var weight = 1.0 / (dist ^ p);
        
        totalWeight += weight;
        weightedSum += weight * arguments.yTrain[i];
    }
    
    return weightedSum / totalWeight;
}

// --- 3. GENETIC ALGORITHM CORE ENGINE (OPERATING ON CHEAP SURROGATE) ---
public array function initPopulation(required numeric popSize, required array bounds) {
    var pop = [];
    for (var i = 1; i <= arguments.popSize; i++) {
        var ind = [];
        for (var d = 1; d <= arrayLen(arguments.bounds); d++) {
            arrayAppend(ind, rand() * (arguments.bounds[d].max - arguments.bounds[d].min) + arguments.bounds[d].min);
        }
        arrayAppend(pop, ind);
    }
    return pop;
}

public array function crossover(required array p1, required array p2) {
    var child = [];
    var alpha = 0.5; // Blend Crossover parameter
    for (var d = 1; d <= arrayLen(arguments.p1); d++) {
        var gamma = (1 + 2 * alpha) * rand() - alpha;
        arrayAppend(child, arguments.p1[d] + gamma * (arguments.p2[d] - arguments.p1[d]));
    }
    return child;
}

public array function mutate(required array individual, required array bounds, numeric rate=0.2, numeric scale=0.5) {
    var mutated = duplicate(arguments.individual);
    for (var d = 1; d <= arrayLen(mutated); d++) {
        if (rand() < arguments.rate) {
            // Box-Muller transform for normal distribution mutation
            var u1 = rand(); var u2 = rand();
            if(u1 == 0) u1 = 0.0001;
            var normalRandom = sqrt(-2.0 * log(u1)) * cos(2.0 * pi() * u2);
            
            mutated[d] += normalRandom * arguments.scale;
            // Clip boundaries
            if (mutated[d] < arguments.bounds[d].min) mutated[d] = arguments.bounds[d].min;
            if (mutated[d] > arguments.bounds[d].max) mutated[d] = arguments.bounds[d].max;
        }
    }
    return mutated;
}

public array function runGaOnSurrogate(required array xTrain, required array yTrain, required array bounds, required numeric popSize, required numeric generations) {
    var population = initPopulation(arguments.popSize, arguments.bounds);
    
    for (var gen = 1; gen <= arguments.generations; gen++) {
        // Evaluate complete generation instantly on the analytical surrogate proxy
        var scoredPop = [];
        for (var i = 1; i <= arrayLen(population); i++) {
            arrayAppend(scoredPop, {
                vector: population[i],
                score: predictSurrogateFitness(population[i], arguments.xTrain, arguments.yTrain)
            });
        }
        
        // Rank-based Sort Descending
        arraySort(scoredPop, function(a, b) {
            return b.score > a.score ? 1 : (b.score < a.score ? -1 : 0);
        });
        
        // Rebuild next generation
        var nextGen = [ scoredPop[1].vector, scoredPop[2].vector ]; // Elitist strategy
        while (arrayLen(nextGen) < arguments.popSize) {
            var parent1 = scoredPop[randRange(1, 5)].vector;
            var parent2 = scoredPop[randRange(1, 5)].vector;
            var child = crossover(parent1, parent2);
            child = mutate(child, arguments.bounds);
            arrayAppend(nextGen, child);
        }
        
        // Update population pointer
        for (var k = 1; k <= arrayLen(nextGen); k++) {
            population[k] = nextGen[k];
        }
    }
    
    // Return best vector from final generation iteration
    return population[1];
}

// --- 4. EXECUTING THE HYBRID SEED LOOP ---
writeOutput("Starting CFML Hybrid GA-RKHS Optimization Engine...<br><br>");

// Step A: Initialize Gram Dataset
variables.X_train = initPopulation(10, variables.bounds);
variables.y_train = [];
for (var i = 1; i <= arrayLen(variables.X_train); i++) {
    arrayAppend(variables.y_train, trueExpensiveFitness(variables.X_train[i]));
}

// Step B: Active Learning Global Intercept Loop
for (var loop = 1; loop <= variables.surrogateMaxIter; loop++) {
    writeOutput("--- Iteration #loop#/#variables.surrogateMaxIter# ---<br>");
    writeOutput("Gram Matrix Dataset Size: " & arrayLen(variables.X_train) & " points<br>");
    
    // Step C: Let the GA scan the math surrogate space aggressively
    var proposedCandidate = runGaOnSurrogate(variables.X_train, variables.y_train, variables.bounds, variables.popSize, variables.generations);
    
    // Step D: Run costly physical simulation check on candidate
    var trueScore = trueExpensiveFitness(proposedCandidate);
    
    writeOutput("GA Proposed Vector: [" & arrayToList(proposedCandidate, ", ") & "]<br>");
    writeOutput("True Real-World Fitness Score: " & numberFormat(trueScore, "0.0000") & "<br>");
    
    // Step E: Augment data matrices
    arrayAppend(variables.X_train, proposedCandidate);
    arrayAppend(variables.y_train, trueScore);
    
    writeOutput("Current System Best Score: " & numberFormat(arrayMax(variables.y_train), "0.0000") & "<br><br>");
}

writeOutput("Optimization Execution Cycle Terminated.");
</cfscript>

Key Architectural Shifts in the CFML Variant

  • Script-Only Architecture: Written purely within <cfscript> blocks to remove legacy tag overhead and maximise processing throughput for vector arithmetic loops.
  • Native Box-Muller Generator: Because ColdFusion’s native rand() produces a flat uniform distribution, a standard Box-Muller transformation was written directly into the mutate() method to correctly inject normal (Gaussian) distribution mutations.
  • Surrogate Mapping Array Matrix: Replaces Python arrays with native array structs. The predictSurrogateFitness method models continuous field behaviours natively via distances weighted dynamically across the current known points.

Lec 50 Kernel Ridge Regression (KRR)

NPTEL - Indian Institute of Science, Bengaluru

YouTube: 09/2025 35:13
By ,
Tags , , , , , ,

No comments:

Post a Comment

Subscribe to: Post Comments (Atom)