Example: "Custom" Radial Basis Function Approximation

Data is generated from the function
e^\frac{y}{2.0}\sin{x}\cos\frac{y}{2.0}
where a number of (x,y ) pairs make up a set of randomly chosen points. Random noise is added to the values, a "custom" Polyharmonic Spline radial basis function is specified
\varphi(r)=\left\{ \begin{array}{ll}r^{k}&k=1,3,5,\ldots\\r^{k}\ln{r}&k=2,4,6,\ldots\end{array} \right.
and a radial basis approximation of the noisy data is computed. The radial basis fit is then compared to the original function at another set of randomly chosen points. Both the average normalized error and the maximum normalized error are computed and printed.

In this example, the order of the Polyharmonic Spline, k =2. The function is sampled at 200 random points and the error is computed at 10000 random points.


import com.imsl.math.*;
import com.imsl.stat.Random;

public class RadialBasisEx2 {

    // The function to approximate
    static double fcn(double x[]) {
        return Math.exp((x[1]) / 2.0) * Math.sin(x[0]) - Math.cos((x[1]) / 2.0);
    }

    static public class PolyHarmonicSpline implements RadialBasis.Function {

        private int order = 3;
        private boolean isEven = false;

        public PolyHarmonicSpline(int order) {
            this.isEven = order % 2 == 0;
            this.order = order;
        }

        public double f(double x) {
            if (this.isEven) {
                return Math.pow(x, order) * Math.log(x);
            }
            return Math.pow(x, order);
        }

        public double g(double x) {
            if (order == 1) {
                return 1;
            }
            if (this.isEven) {
                return order * Math.pow(x, order - 1) * Math.log(x) + Math.pow(x, order - 1);
            }
            return order * Math.pow(x, order - 1);
        }

        public int getOrder() {
            return order;
        }

        public boolean isEvenOrder() {
            return isEven;
        }
    }

    public static void main(String args[]) {
        int nDim = 2;

        // Sample, with noise, the function at 100 randomly choosen points
        int nData = 200;
        double xData[][] = new double[nData][nDim];
        double fData[] = new double[nData];
        Random rand = new Random(123457);
        rand.setMultiplier(16807);
        double noise[] = new double[nData * nDim];
        for (int k = 0; k < nData; k++) {
            for (int i = 0; i < nDim; i++) {
                noise[k * 2 + i] = 1.0d - 2.0d * rand.nextDouble();
                xData[k][i] = 3 * noise[k * 2 + i];
            }
            // noisy sample
            fData[k] = fcn(xData[k]) + noise[k * 2] / 10;
        }

        // Compute the radial basis approximation using 100 centers
        int nCenters = 100;
        RadialBasis rb = new RadialBasis(nDim, nCenters);
        rb.setRadialFunction(new PolyHarmonicSpline(2));
        rb.update(xData, fData);

        // Compute the error at a randomly selected set of points
        int nTest = 10000;
        double maxError = 0.0;
        double aveError = 0.0;
        double maxMagnitude = 0.0;
        double x[][] = new double[nTest][nDim];
        noise = new double[nTest * nDim];

        for (int i = 0; i < nTest; i++) {
            for (int j = 0; j < nDim; j++) {
                noise[i * 2 + j] = 1.0d - 2.0d * (double) rand.nextDouble();
                x[i][j] = 3 * noise[i * 2 + j];
            }
            double error = Math.abs(fcn(x[i]) - rb.value(x[i]));
            maxMagnitude = Math.max(Math.abs(fcn(x[i])), maxMagnitude);
            aveError += error;
            maxError = Math.max(error, maxError);
        }
        aveError /= nTest;

        System.out.println("Average normalized error is " + aveError / maxMagnitude);
        System.out.println("Maximum normalized error is " + maxError / maxMagnitude);
        System.out.println("Using even order equation: " +
                ((PolyHarmonicSpline) rb.getRadialFunction()).isEvenOrder());
    }
}

Output

Average normalized error is 0.018055872051037286
Maximum normalized error is 0.25756531087544127
Using even order equation: true
Link to Java source.