added matlab code for ELM

learning_folder/elmmodel/data/Ax.mat

@@ -0,0 +1,8 @@
+6
+6
+0.1666666666666667 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 
+0.0000000000000000 0.5000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 
+0.0000000000000000 0.0000000000000000 0.2348327403806639 0.0000000000000000 0.0000000000000000 0.0000000000000000 
+0.0000000000000000 0.0000000000000000 0.0000000000000000 0.2777777777777778 0.0000000000000000 0.0000000000000000 
+0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.3571428571428572 0.0000000000000000 
+0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.2136752136752137 

learning_folder/elmmodel/data/Ay.mat

@@ -0,0 +1,8 @@
+6
+6
+14.1999999999999993 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 
+0.0000000000000000 82.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 
+0.0000000000000000 0.0000000000000000 61.8000000000000043 0.0000000000000000 0.0000000000000000 0.0000000000000000 
+0.0000000000000000 0.0000000000000000 0.0000000000000000 68.4000000000000057 0.0000000000000000 0.0000000000000000 
+0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 12.1999999999999993 0.0000000000000000 
+0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 0.0000000000000000 60.3999999999999915 

learning_folder/elmmodel/data/bx.vec

@@ -0,0 +1,2 @@
+6
+0.0000000000000000 -2.1000000000000001 -2.1230699999999998 -3.1400000000000001 -1.3999999999999999 -1.5700000000000003 

learning_folder/elmmodel/data/by.vec

@@ -0,0 +1,2 @@
+6
+-7.9999999999999991 -58.2000000000000028 -37.6000000000000014 -27.0000000000000000 -6.2000000000000002 -32.3999999999999986 

learning_folder/elmmodel/info

@@ -0,0 +1 @@
+elm+datascale

learning_folder/learner/algorithms/ELM.m

@@ -0,0 +1,150 @@
+classdef ELM < Learner
+    %ELM Extreme Learning Machine
+    %   Single layer feed-forward network
+
+    properties
+        hidDim = 10;    %number of neurons in the hidden layer
+        wInp = [];      %weights from input layer to hidden layer
+        wOut = [];      %weights from hidden layer to output layer
+        a = [];         %slope parameters of activation functions
+        b = [];         %bias parameters of activation functions
+
+        reg = 1e-6;     %regularization parameter for ridge regression
+        BIP = 1;        %flag for using Batch Intrinsic Plasticity (BIP)
+        mu = 0.2;       %desired mean activity parameter for BIP
+
+        batchSize = 5000; %number of samples used for a minibatch
+
+        confest=[];
+    end
+
+    methods
+        function l = ELM(inpDim, outDim, spec)
+            l = l@Learner(inpDim, outDim, spec);
+        end
+
+        function init(l, X)
+            l.a = ones(1,l.hidDim);
+            l.b = 2 * rand(1,l.hidDim) - ones(1,l.hidDim);
+            l.wInp = 2 * rand(l.hidDim,l.inpDim) - ones(l.hidDim,l.inpDim);
+            l.wOut = 2 * rand(l.hidDim,l.outDim) - ones(l.hidDim,l.outDim);
+
+            if nargin > 1 && l.BIP
+                if iscell(X)
+                    X = l.normalizeIO(cell2mat(X));
+                else
+                    X = l.normalizeIO(X);
+                end
+                l.bip(X);
+            end
+        end
+
+        function train(l, X, Y)
+            if iscell(X)
+                [X Y] = l.normalizeIO(cell2mat(X), cell2mat(Y));
+            else
+                [X Y] = l.normalizeIO(X, Y);
+            end
+            if size(X,1) < l.batchSize
+                hs = l.calcHiddenStates(X);
+                l.wOut = (hs'*hs + l.reg * eye(l.hidDim))\(hs' * Y);
+            else
+                HTH = zeros(l.hidDim, l.hidDim);
+                HTY = zeros(l.hidDim, l.outDim);                
+                numBatches = floor(size(X,1)/l.batchSize);
+                rest = mod(size(X,1),l.batchSize);                
+                for b=1:numBatches
+                    H = l.calcHiddenStates(X((b-1)*l.batchSize+1:b*l.batchSize,:));
+                    HTH = HTH + H' * H;
+                    HTY = HTY + H' * Y((b-1)*l.batchSize+1:b*l.batchSize,:);
+                end
+                if rest > 0
+                    H = l.calcHiddenStates(X(numBatches*l.batchSize+1:end,:));
+                    HTH = HTH + H' * H;
+                    HTY = HTY + H' * Y(numBatches*l.batchSize+1:end,:);
+                end                
+                HTH = HTH + diag(repmat(l.reg, 1, l.hidDim));
+                l.wOut = HTH\HTY;
+            end
+        end
+
+        function [Y] = apply(l, X)
+            conf=[];
+            if iscell(X)
+                Y = cell(length(X),1);
+                for i=1:length(X)
+                    Y{i} = l.apply(X{i});
+                end
+            else
+                X = range2norm(X, l.inpRange, l.inpOffset);
+                if size(X,1) < l.batchSize
+                    H = l.calcHiddenStates(X);
+                    Y = H * l.wOut;
+
+
+                else
+                    Y = zeros(size(X,1), l.outDim);
+                    numBatches = floor(size(X,1)/l.batchSize);
+                    rest = mod(size(X,1),l.batchSize);
+                    for b=1:numBatches
+                        H = l.calcHiddenStates(X((b-1)*l.batchSize+1:b*l.batchSize,:));
+                        Y((b-1)*l.batchSize+1:b*l.batchSize,:) = H * l.wOut;
+                    end
+                    if rest > 0
+                        H = l.calcHiddenStates(X(numBatches*l.batchSize+1:end,:));
+                        Y(numBatches*l.batchSize+1:end,:) = H * l.wOut;
+                    end
+                end
+                Y = norm2range(Y, l.outRange, l.outOffset);
+                l.out = Y(end,:);
+            end
+        end
+
+
+        function [H G] = calcHiddenStates(l, X)
+            numSamples = size(X,1);
+
+            %TODO this is extremely costly, make mex call for fermi fct
+            atemp = repmat(l.a, numSamples, 1);
+            btemp = repmat(l.b, numSamples, 1);
+            %
+
+            G = X * l.wInp';
+
+            H = 1./(1+exp(-atemp .* (G - btemp)));
+        end
+
+
+
+        %batch intrinsic plasticity
+        function bip(l, X)
+            numSamples = size(X,1);
+
+            PD = ProbDistUnivParam('exponential', l.mu);
+            G = X * l.wInp';
+
+            for hn = 1:l.hidDim
+                targets = random(PD, 1, numSamples);
+
+                hightars = length([targets(targets>=1),targets(targets<=0)]);
+                while hightars > 0
+                    further_targets = random(PD, 1, hightars);
+                    targets = [targets(targets<1&targets>0),further_targets];
+                    hightars = length([targets(targets>=1),targets(targets<=0)]);
+                end
+                targets = sort(targets);
+
+                s = sort(G(:,hn));
+                Phi = [s,ones(numSamples,1)];
+                targetsinv = -log(1./targets - 1); %apply inverse activation function
+
+                w = pinv(Phi)*targetsinv';                
+
+                l.a(hn) = w(1);
+                l.b(hn) = w(2);
+            end
+        end
+    end
+
+end
+

learning_folder/learner/algorithms/Learner.m

@@ -0,0 +1,70 @@
+classdef Learner < handle
+    %LEARNER Interface class for function approximator models.
+    %   Detailed explanation goes here
+
+    properties
+        inpDim = 0; %dimension of input
+        outDim = 0; %dimension of output
+        out = [];
+        outvar=NaN;
+
+        inpOffset = [];
+        inpRange = [];
+        outOffset = [];
+        outRange = [];
+    end
+
+    methods
+        function l = Learner(inpDim, outDim, spec)
+            l.inpDim = inpDim;
+            l.outDim = outDim;
+            l.out = zeros(1,l.outDim);
+
+            if nargin > 2
+                for s=1:size(spec,1)
+                    if ~strcmp(spec{s,1}, 'class')
+                        try
+                            l.(spec{s,1}) = spec{s,2};
+                        catch me
+                            disp(me.message);
+                        end
+                    end
+                end
+            end
+        end
+
+        function [X Y] = normalizeIO(l, X, Y)
+            [X l.inpOffset l.inpRange] = normalize(X);
+            if exist('Y', 'var')
+                [Y l.outOffset l.outRange] = normalize(Y);
+            end
+        end
+
+        function init(l, X)
+        end
+
+        % expect X and Y to be either matrices or cell-arrays of matrices
+        function train(l, X, Y)
+        end
+
+        function Y = apply(l, X)
+
+        end
+
+        function new = copy(this)
+            new = feval(class(this));
+            p = properties(this);
+            for i = 1:length(p)
+                if( isa(this.(p{i}), 'handle'))
+                    new.(p{i}) = this.(p{i}).copy();
+                elseif isa(this.(p{i}),'cell')
+                    new.(p{i}) = deepCopyCellArray(numel(this.(p{i})));
+                else
+                    new.(p{i}) = this.(p{i});
+                end
+            end
+        end
+    end
+
+end
+

learning_folder/learner/metric/Metric.m

@@ -0,0 +1,22 @@
+classdef Metric < handle
+    %METRIC Summary of this class goes here
+    %   Detailed explanation goes here
+
+    properties
+    end
+
+    methods
+        function m = Metric()
+        end
+
+        function d = distance(m, a, b)
+            if size(b,1) < size(a,1)
+                d = sqrt(sum((a-repmat(b,size(a,1),1)).^2,2));
+            else
+                d = sqrt(sum((a-b).^2,2));
+            end
+        end
+    end
+
+end
+

learning_folder/learner/metric/mexSqEucDistance.cpp

@@ -0,0 +1,37 @@
+#include <math.h>
+#include <matrix.h>
+#include <mex.h>
+
+//input: vector a, vector b
+//output: vector sum((a-b).^2)
+
+void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
+    //dx = l.alpha * sum(repmat(l.h', 1, l.visDim) .* (l.center.C - repmat(x, l.hidDim, 1)), 1);
+
+    //get dimensions
+    int dim = mxGetDimensions(prhs[0])[1]; //row vector
+    int numSamples = mxGetDimensions(prhs[0])[0];
+    int numSamplesB = mxGetDimensions(prhs[1])[0];
+
+    //get inputs
+    double *a = mxGetPr(prhs[0]);
+    double *b = mxGetPr(prhs[1]);
+
+    //create output vector
+    plhs[0] = mxCreateDoubleMatrix(numSamples,1,mxREAL); //column vector
+    double *d = mxGetPr(plhs[0]);
+
+    if (numSamples == numSamplesB) {
+        for (int i=numSamples-1; i>=0; i--) {
+            for (int j=dim-1; j>=0; j--) {
+                d[i] += pow(a[j*numSamples+i] - b[j*numSamples+i], 2);
+            }
+        }
+    } else {
+        for (int i=numSamples-1; i>=0; i--) {
+            for (int j=dim-1; j>=0; j--) {
+                d[i] += pow(a[j*numSamples+i] - b[j], 2);
+            }
+        }
+    }
+}

learning_folder/learner/metric/mexWeightedSqEucDistance.cpp

@@ -0,0 +1,38 @@
+#include <math.h>
+#include <matrix.h>
+#include <mex.h>
+
+//input: vector a, vector b, vector w
+//output: vector sum(w .* (a-b).^2)
+
+void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
+    //dx = l.alpha * sum(repmat(l.h', 1, l.visDim) .* (l.center.C - repmat(x, l.hidDim, 1)), 1);
+
+    //get dimensions
+    int dim = mxGetDimensions(prhs[0])[1]; //row vector
+    int numSamples = mxGetDimensions(prhs[0])[0];
+    int numSamplesB = mxGetDimensions(prhs[1])[0];
+
+    //get inputs
+    double *a = mxGetPr(prhs[0]);
+    double *b = mxGetPr(prhs[1]);
+    double *w = mxGetPr(prhs[2]);
+
+    //create output vector
+    plhs[0] = mxCreateDoubleMatrix(numSamples,1,mxREAL); //column vector
+    double *d = mxGetPr(plhs[0]);
+
+    if (numSamples == numSamplesB) {
+        for (int i=numSamples-1; i>=0; i--) {
+            for (int j=dim-1; j>=0; j--) {
+                d[i] += w[j] * pow(a[j*numSamples+i] - b[j*numSamples+i], 2);
+            }
+        }
+    } else {
+        for (int i=numSamples-1; i>=0; i--) {
+            for (int j=dim-1; j>=0; j--) {
+                d[i] += w[j] * pow(a[j*numSamples+i] - b[j], 2);
+            }
+        }
+    }
+}