classdef (CaseInsensitiveProperties, TruncatedProperties) idSupportVectorMachine <...
      idnlfun & nlident.internal.NonlinearMap & nlident.internal.UnstructuredFcn & ... 
      idoptions.MappingFcnEstimationOptionsInterface
   %idSupportVectorMachine Support Vector Machine function for use in nonlinear system identification.
   %
   %  Requires Statistics and Machine Learning Toolbox(TM)
   %
   % idSupportVectorMachine is a nonlinear mapping function that consists of a
   % nonlinear function, which is a Support Vector Machine.
   %
   %                  --------------------                 
   % Inputs (X) ---> | Nonlinear Function | ----> Output (y)  
   %                  --------------------           
   %
   % Mathematically, a Support Vector Machine is an R^m->R^1 function y = S(X)
   %  where:
   %     X is m-by-1 column vector of inputs  
   %     S(X) = \sum{n=1}^{N} alpha_n * G(x_n,X) + b
   %
   %     x_n is a support vector of the trained model
   %     alpha_n is the weight associated with the corresponding support vector
   %     G is the kernel function operating on the input X
   %     b is the offset of the trained model  
   %     N is the number of support vectors of the trained model
   %
   %  Usage:
   %   S = idSupportVectorMachine
   %   creates an @idSupportVectorMachine object, S, that represents a multi-input,
   %   single-output nonlinear function. This object functions as a type of
   %   mapping function for the Nonlinear ARX models; it maps the model's
   %   regressors to its output. It facilitates the use of a support vector machine
   %   of Statistics and Machine Learning Toolbox(TM) for system
   %   identification.
   %
   %   S = idSupportVectorMachine(KernelName)
   %   creates a @idSupportVectorMachine object, S, that uses the specified
   %   Kernel. Specify KernelName as one of 'Gaussian','rbf','Linear',
   %   'Polynomial'.
   %   Default: 'Gaussian'
   %
   %   S = idSupportVectorMachine(KernelName,EpsilonMargin)
   %   creates a @idSupportVectorMachine object, S, that uses the specified
   %   Kernel. EpsilonMargin is a nonnegative scalar which indicates half the
   %   width of the epsilon-insensitive band. 
   %
   % Properties of a @idSupportVectorMachine object:
   %          Inputs: Inputs of the function - their names, mean, and range.
   %                  The mean and range are determined during estimation.
   %                  The names are assigned when the @idSupportVectorMachine object is
   %                  assigned as a mapping function in the @idnlarx object.
   %         Outputs: Output of the function - its name, mean, and range.  
   %  KernelFunction: Name of the kernel function being employed
   %                  by the Support Vector Machine to compute the Gram matrix.
   %                  It can be one of the choices as outlined in FITRSVM function.
   %                  Default: 'Gaussian'
   %     KernelScale: A positive scalar that is used to divide all elements 
   %                  of the predictor matrix. Default: 'auto'
   %   EpsilonMargin: A nonnegative scalar that indicates half the width of
   %                  the epsilon-insensitive band. Default: 'auto'
   % PolynomialOrder: A positive integer that specifies the order of the
   %                  kernel function when the latter is set to
   %                  'Polynomial'. Default: 3
   %    KernelOffset: A nonnegative scalar that is added to each element of
   %                  the Gram matrix. Default: 0
   %   BoxConstraint: A positive scalar which is the upper bound for the
   %                  absolute value of the Alpha coefficients. Default: 1
   %
   %      Parameters: Parameters of trained Support Vector Machine model.
   %                  * Alpha: Weight associated with each support
   %                           vector of a trained Support Vector Machine
   %                           Model. Use this property to provide initial
   %                           estimates for the Alpha coefficients. Alpha
   %                           can be set before estimation but cannot be
   %                           modified afterwards. Default: []
   %      * LinearCoefficient: Read-only property indicating
   %                           the linear coefficients for the primal
   %                           optimization problem. This is determined
   %                           while estimating the Support Vector Machine
   %                           model. This property is non-empty only when
   %                           model is estimated using a 'Linear' kernel
   %                           function.
   %                   * Bias: Read-only property indicating the
   %                           offset term in the Support Vector
   %                           Machine model.
   %        * IsSupportVector: Read-only property indicating
   %                           whether an observation is a support vector.
   %
   %  Free: Option to update the parameters of idSupportVectorMachine during 
   %        nonlinear ARX model estimation, specified as true or false. When  
   %        free is true, the estimation process updates the idSupportVectorMachine
   %        model when it estimates the idnlarx model that contains it. 
   %        When Free is false, the idSupportVectorMachine model is fixed 
   %        during estimation. Setting Free to false is useful when you are
   %        using a previously estimated idSupportVectorMachine model as
   %        a mapping function for nlarx.
   %
   % EstimationOptions: Optimization options used during estimation. The
   %     supported options are: 
   %     * OutlierFraction: Numeric scalar in the interval [0, 1)
   %                        indicating the expected proportion of outliers
   %                        in training data.
   %           * CacheSize: Positive scalar indicating the cache size in
   %                        megabytes that the software reserves in memory
   %                        for training the model. Default: 1000
   %              * Solver: Solver used for parameter estimation. Specify
   %                        as one of:
   %                        'SMO' : Sequential Minimal Optimization
   %                        'ISDA': Iterative Single Data Algorithm
   %
   %   Use idSupportVectorMachine as the value of the "OutputFcn" property of IDNLARX
   %   models. Use the EVALUATE function to compute the output
   %   of the function to a given vector of inputs.
   %
   %  Example:
   %    load iddata1.mat z1;
   %    % Prepare the output function
   %    s = idSupportVectorMachine;
   %    s.KernelScale = 5;
   %    % Estimate the Nonlinear ARX model
   %    sys = nlarx(z1,[4 4 1],s);
   %    compare(z1,sys)
   %
   %   In this example, Fcn = sys.OutputFcn is an idSupportVectorMachine with
   %   default settings. The estimated model sys has an idSupportVectorMachine
   %   object for the output function. Fcn implements a mapping from
   %   regressor vector to scalar reals. The value of the function is
   %   computed by evaluate(Fcn,x), where x is the regressor matrix. The
   %   command COMPARE is used to compare the output of the estimated
   %   nonlinear ARX model with the idSupportVectorMachine nonlinearity with the
   %   measured data.
   %
   %   Type "Fcn = idSupportVectorMachine; Fcn.EstimationOptions" to view the
   %   idSupportVectorMachine estimation options.
   %
   % See also idnlarx, idnlhw, nlarx, nlhw, evaluate, idSigmoidNetwork,
   % idTreePartition, idCustomNetwork, idWaveletNetwork, idLinear,
   % idGaussianProcess, idTreeEnsemble.

   % Copyright 2021-2023 The MathWorks, Inc.

   properties (Constant, Hidden)
      NAME = getString(message('Ident:idnlfun:strSVMFcn'));
      HasLinearFcn = false;
      KernelFunctionList = {'Gaussian','rbf','Linear','Polynomial'}
   end

   properties (Access = protected)
      prvKernelFunction  char {mustBeMember(prvKernelFunction,{'Gaussian','rbf',...
         'Linear','Polynomial'})} = 'Gaussian'
      prvEpsilonMargin (1,1) {mustBeNumeric, mustBeNonnegative, mustBeFinite, mustBeReal} = 0.1
   end

   properties (Dependent, AbortSet)
      KernelFunction char {mustBeMember(KernelFunction,{'Gaussian', 'rbf', 'Linear', 'Polynomial'})}      
      EpsilonMargin
   end

   properties (AbortSet)
      PolynomialOrder (1,1) {mustBeInteger, mustBePositive, mustBeFinite, mustBeReal} = 3
      KernelOffset (1,1) {mustBeNumeric, mustBeNonnegative, mustBeFinite, mustBeReal} = 0
      BoxConstraint (1,1) {mustBeNumeric, mustBeNonnegative, mustBeFinite, mustBeReal} = 1
      KernelScale = 'auto'
   end   

   methods (Hidden)
      [yhat, jcb, dy_x] = getJacobian(nlobj, x, doParJac)
      this = soinitialize(this, yvec, regmat, varargin);
      yhat = soevaluate(this, regmat, varargin);
   end

   methods
      function this  = idSupportVectorMachine(KernelName,EpsilonMargin)
         % Constructor

         narginchk(0, 2);
         this.EstimationOptions = idoptions.idSupportVectorRegressionOptions;
         this.Parameters = nlident.internal.SVMParameters;

         if nargin == 0
            KernelName = 'Gaussian';
         end
         this.KernelFunction = ltipack.matchKey(KernelName,this.KernelFunctionList);

         if nargin == 2
            this.EpsilonMargin = EpsilonMargin;
         end
      end
   end

   methods (Hidden)
      function str = getInfoString(this)
         %getInfoString  1-line display string
         str = getString(message('Ident:idnlfun:strSVMDisplay',this.KernelFunction));
      end

      function this = initializeIOs(this,In,Out,varargin)
         % initialize names
         this = initializeIOs@idnlfun(this,In,Out,varargin{:});
      end

      function k = getInputIndex(this)
         % Indices of fcn inputs that are used relative to a global list.
         k = 1:getNumInputs(this);
      end

      function this = soreset(this)
         this.prvNonlinearFcn = [];
      end

      function this = renameInputs(this, Old, New)
         this = renameInputs@nlident.internal.NonlinearMap(this, Old, New);
      end

      function this = renameOutputs(this, Old, New)
         this = renameOutputs@nlident.internal.NonlinearMap(this, Old, New);
      end

      function NL = getNonlinearFcn(sys)
         % app/loadobj use
         NL = getNonlinearFcn@idnlfun(sys);
      end

      function sys = setNonlinearFcn(sys,NL)
         % app/loadobj use
         sys = setNonlinearFcn@idnlfun(sys,NL);
      end
   end

   % Block for dependent properties
   methods
      % Name
      function value = get.KernelFunction(this)
         value = this.prvKernelFunction;
      end

      function this = set.KernelFunction(this,value)
         this.prvKernelFunction = value;
      end

      % Kernel Scale
      function this = set.KernelScale(this,value)
         if isa(value,'char') && strcmpi(value,'auto')
            this.KernelScale = 'auto';
         else
            validateattributes(value,{'numeric'},{'positive','finite'});
            this.KernelScale = value;
         end       
      end

      % Epsilon Margin
      function value = get.EpsilonMargin(this)
         if this.prvEpsilonMargin == 0.1 && ~isa(this.prvNonlinearFcn,'RegressionSVM')
            value = 'auto';
         else
            value = this.prvEpsilonMargin;
         end
      end

      function this = set.EpsilonMargin(this,value)
         if isequal(value,'auto')
            value = 0.1;
         end
         this.prvEpsilonMargin = value;
      end
   end

   methods (Access = protected)
      function this = setParameters_(this, value)
         if isa(value,'nlident.internal.SVMParameters')
            this.Parameters_ = value;
         end
      end

      function status = soisinitialized(this)
         status = soisinitialized@idnlfun(this) && ...
            isa(this.prvNonlinearFcn,'RegressionSVM');
      end

      function status = soisdifferentiable(varargin)
         % Is the nonlinear mapping function differentiable?
         % FALSE here since Jacobian is not implemented yet.
         status = false;
      end

      function this = soSetFree(this,Free)
         this.Free = any(Free);
      end

      function [th, Free, MinMax] = sogetParameterVector(this)
         %sogetParameterVector returns the parameter vector of a single object.
         %
         %  th = sogetParameterVector(nlobj)
         %  [th, free_flag, min_max_bounds] = sogetParameterVector(nlobj)
         th = getParameterVector(this.Parameters);
         if nargout>1
            Free = this.Free & true(size(th));
            n = length(th);
            MinMax = [-Inf(n,1),Inf(n,1)];
         end
      end

      function [this, S] = sosetParameterVector(this, th, varargin)
         %sosetParameterVector Set the parameters of a single object.
         % this = sosetParameterVector(this, th) sets all the parameters,
         % fixed or free.
         % this = sosetParameterVector(this, th, 'free') sets only the free
         % parameters; th must contain only the free parameters.
         % S is a structure of all parameters, fixed or free.

         % no-op (no settable parameters yet)
         S = [];
      end

      function n = sonparams(this, varargin)
         %sonparams Number of parameters.
         % n = sonparams(obj): total number of parameters
         % n = sonparams(obj,'free'): number of free parameters only
         if isa(this.prvNonlinearFcn,'RegressionSVM')
            n = length(sogetParameterVector(this));
         else
            n = NaN;
         end
      end
      
      function PV = getPVFromNonlinearFunction(this)
         % Create PV pairs for fitrsvm
         PV = {'KernelFunction',this.KernelFunction};
         PV = [PV, {'KernelScale',this.KernelScale}];
         

         if strcmpi(this.KernelFunction,'Polynomial')
            PV = [PV, {'PolynomialOrder',this.PolynomialOrder}];
         end

         if this.KernelOffset ~= 0
            PV = [PV, {'KernelOffset',this.KernelOffset}];
         end

         if this.BoxConstraint ~= 1
            PV = [PV, {'BoxConstraint',this.BoxConstraint}];
         end

         if ~strcmpi(this.EpsilonMargin,'auto')
            PV = [PV, {'Epsilon',this.prvEpsilonMargin}];
         end

         PV = getPVFromNonlinearFunctionParameters(...
            this.Parameters,PV);
      end

      function this = updateNonlinearFunctionPostEstimation(this)
         % Update NonlinearFunction using estimated model
         mdl = this.prvNonlinearFcn;
         this.EpsilonMargin = mdl.Epsilon;
         this.KernelScale = mdl.Impl.KernelParameters.Scale;
         this.KernelOffset = mdl.Impl.KernelParameters.Offset;

         % Update Alpha, Beta, Bias, and IsSupportVector
         this.Parameters = updateSupportVectorParametersPostEstimation(...
            this.Parameters,mdl);
      end  

      function boo = soisLinearInFreeParameters(~)
         boo = false; % not used 
      end

      function [g, np0] = sogetParGroupsForSparseEstimation(this, varargin)
         [g, np0] = sogetParGroupsForSparseEstimation@nlident.internal.NonlinearMap(this, varargin{:});
      end

      function this = soremoveInputs(this, ~)
         % remove inputs of indicated names
         % no -op
      end

      function this = soundoProjection(this)
         % no-op
      end
   end

   methods (Hidden)
      function this = setpar(this, varargin)
         % no-op
      end

      function [Props,AsgnVals] = pnames(varargin)
         % Default implementation.
         Props = '';
         AsgnVals = {''};
      end

      function mdl = getRegressionSVMModel(this)
         mdl = this.prvNonlinearFcn;
      end

      function sys = reload(sys,s,ver)
         % common elements of loadobj operations required R2022a onwards.
         % no-op, since object introduced in R2022a
         sys = reload@nlident.internal.NonlinearMap(sys,s,ver);
         sys = reload@nlident.internal.UnstructuredFcn(sys,s,ver);
         sys = reload@idoptions.MappingFcnEstimationOptionsInterface(sys,s);
         sys = reload@idnlfun(sys,s,ver);
      end 
   end

   methods (Static, Hidden)
      function sys = loadobj(s)
         v = nlident.internal.FcnIO.utGetVersion(s);
         if isstruct(s)
            updatewarn
            Warn = ctrlMsgUtils.SuspendWarnings; %#ok<NASGU>
            sys = idSupportVectorMachine(s.prvKernelFunction, s.prvEpsilonMargin);
            sys.PolynomialOrder = s.PolynomialOrder;
            sys.KernelOffset = s.KernelOffset;
            sys.BoxConstraint = s.BoxConstraint;
            sys = reload(sys,s,v);
         else
            sys = s;
         end
      end
   end
end
