Package com.opengamma.strata.math.impl.interpolation

Class PenaltyMatrixGenerator

  • public abstract class PenaltyMatrixGenerator
extends Object
    The k^th order difference matrix will act on a vector to produce the k^th order difference series. Related to this is the k^th order finite difference differentiation matrix which acts on the values of a function evaluated at a (non-uniformly spaced) set of points, to give the estimate of the k^th order derivative at those points (for unit spaces points these are identical). One use of these matrices, is to provide matrices that penalise the gradient or curvature of data on a non-uniform grid. This should work for an arbitrary number of dimensions (for data that has been flattened to a vector).

    • Constructor Detail

      • PenaltyMatrixGenerator

        public PenaltyMatrixGenerator()

    • Method Detail

      • getDifferenceMatrix

        public static DoubleMatrix getDifferenceMatrix​(int m,
                                               int k)
        get the k^th order difference matrix, D, which acts on a vector, x, of length m to produce the k^th order difference vector. The first k rows of D are set to zero - because of this D_1 * D_1 != D_2

        For example, for x = {x1,x2,....xn}, D_1 * x = {0, x2-x1, x3-x2,...., xn - x(n-1)} and D_2 * x = {0,0, x3-2*x2+x1,....,xn-2*x(n-1)+x(n-2)} etc

        Parameters:
        m - Length of the vector
        k - Difference order. Require m > k
        Returns:
        The k^th order difference matrix

      • getPenaltyMatrix

        public static DoubleMatrix getPenaltyMatrix​(int m,
                                            int k)
        get the k^th order penalty matrix, P. This is defined as P = (D^T)*D , where D is the k^th order difference matrix (see getDifferenceMatrix), so the scalar amount (x^T)*P*x = |Dx|^2 is greater the more k^th order differences there are in the vector, x. This can then act as a penalty on x in some optimisation routine where x is the vector of (fit) parameters.
        Parameters:
        m - Length of the vector.
        k - Difference order. Require m > k
        Returns:
        The k^th order penalty matrix, P

      • getPenaltyMatrix

        public static DoubleMatrix getPenaltyMatrix​(int[] numElements,
                                            int k,
                                            int index)
        Assume a tensor has been flattened to a vector as {A_{0,0}, A_{0,1},...._A_{0,m}, A_{1,0}, A_{1,1},...._A_{1,m},...,A_{n,0}, A_{n,1},...._A_{n,m}} (see flattenMatrix(com.opengamma.strata.collect.array.DoubleMatrix)) that is, the last index changes most rapidly. This produces a penalty matrix that acts on a given set of indexes only

        To produce a penalty matrix that acts on multiple indexes, produce one for each set of indexes and add them together (scaling if necessary).

        Parameters:
        numElements - The range of each index. In the example above, this would be {n,m}
        k - Difference order. Require size[indices] > k
        index - Which set of indices does the matrix act on
        Returns:
        A penalty matrix

      • getPenaltyMatrix

        public static DoubleMatrix getPenaltyMatrix​(int[] numElements,
                                            int[] k,
                                            double[] lambda)
        Assume a tensor has been flattened to a vector as {A_{0,0}, A_{0,1},...._A_{0,m}, A_{1,0}, A_{1,1},...._A_{1,m},...,A_{n,0}, A_{n,1},...._A_{n,m}} (see flattenMatrix(com.opengamma.strata.collect.array.DoubleMatrix)) that is, the last index changes most rapidly. This produces the sum of penalty matrices (or order given by k) with each scaled by lambda.
        Parameters:
        numElements - The range of each index. In the example above, this would be {n,m}
        k - The difference order for each dimension
        lambda - The scaling for each dimension
        Returns:
        A penalty matrix

      • getDerivativeMatrix

        public static DoubleMatrix getDerivativeMatrix​(double[] x,
                                               int k,
                                               boolean includeEnds)
        Get the kth order finite difference derivative matrix, D_k(x), for a non-uniform set of points. For a (1D) set of points, x, (which are not necessarily uniformly spaced), and a set of values at those points, y (i.e. y=f(x)) the vector y_k = D_k(x) * y is the finite difference estimate of the k^th order derivative (d^k y/ dx^k) at the points, x.

        K = 0, trivially return the identity matrix; for k = 1 and 2, this is a three point estimate. K > 2 is not implemented

        Parameters:
        x - A non-uniform set of points
        k - The order Only first and second order are currently implemented
        includeEnds - If true the matrix includes an estimate for the derivative at the end points; if false the first and last rows of the matrix are empty
        Returns:
        The derivative matrix

      • getPenaltyMatrix

        public static DoubleMatrix getPenaltyMatrix​(double[] x,
                                            int k)
        get a k^th order penalty matrix,P, for a non-uniform grid, x. $P = D^T D$ where D is the kth order finite difference matrix (not including the ends). Given values y = f(x) at the grid points, y^T*P*y = |Dy|^2 is the sum of squares of the k^th order derivative estimates at the points.

        Note: In certain applications we may need an estimate of $\frac{1}{b-a}\int_a^b \left(\frac{d^k y}{dx^k}\right)^2 dx$ i.e. the RMS value of the k^th order derivative, between a & b - for a uniform x this forms a crude approximation.

        Parameters:
        x - A non-uniform set of points
        k - order The order Only first and second order are currently implemented
        Returns:
        The k^th order penalty matrix

      • getPenaltyMatrix

        public static DoubleMatrix getPenaltyMatrix​(double[][] x,
                                            int k,
                                            int index)
        Get a kth order penalty matrix for a non-uniform grid whose values have been flattened to a vector.
        Parameters:
        x - the grid positions in each dimension
        k - the (finite difference) order
        index - which index to act on
        Returns:
        a penalty matrix

      • getPenaltyMatrix

        public static DoubleMatrix getPenaltyMatrix​(double[][] x,
                                            int[] k,
                                            double[] lambda)
        Get a penalty for a non-uniform grid whose values have been flattened to a vector. This is the sum of penalty matrices that act on each index, scaled by an amount lambda
        Parameters:
        x - the grid positions in each dimension
        k - the (finite difference) order in each dimension
        lambda - the strength of the penalty in each dimension
        Returns:
        a penalty matrix

      • flattenMatrix

        public static DoubleArray flattenMatrix​(DoubleMatrix aMatrix)
        for a matrix {{A_{0,0}, A_{0,1},...._A_{0,m},{A_{1,0}, A_{1,1},...._A_{1,m},...,{A_{n,0}, A_{n,1},...._A_{n,m}} flattened to a vector {A_{0,0}, A_{0,1},...._A_{0,m}, A_{1,0}, A_{1,1},...._A_{1,m},...,A_{n,0}, A_{n,1},...._A_{n,m}}.
        Parameters:
        aMatrix - A matrix
        Returns:
        a the flattened matrix

      • getMatrixForFlattened

        public static DoubleMatrix getMatrixForFlattened​(int[] numElements,
                                                 DoubleMatrix m,
                                                 int index)
        Assume a tensor has been flattened to a vector as {A_{0,0}, A_{0,1},...._A_{0,m}, A_{1,0}, A_{1,1},...._A_{1,m},...,A_{n,0}, A_{n,1},...._A_{n,m}} (see flattenMatrix(com.opengamma.strata.collect.array.DoubleMatrix)) that is, the last index changes most rapidly. Given a matrix, M, that acts on the elements of one index only, i.e. $$y_{i, i_1, i_2, \dots,i_{k-1}, i_{k+1},\dots, i_n} = \sum_{i_k=0}^{N_k-1} M_{i,i_k} x_{i_1, i_2, \dots,i_k,\dots, i_n} $$ form the larger matrix that acts on the flattened vector.
        Parameters:
        numElements - The number of elements in each index. In the example above, this would be {n,m}
        m - the matrix M
        index - Which index does the matrix act on
        Returns:
        A (larger) matrix which acts on the flattened vector