classdef IsolationForest < anomaly.AnomalyDetector
    %IsolationForest Isolation forest.
    %   IsolationForest is a model for anomaly detection.
    %
    %   An object of this class cannot be created by calling the constructor.
    %   Use IFOREST to create an IsolationForest object by fitting the model
    %   to the data.
    %
    %   IsolationForest properties:
    %       CategoricalPredictors       - Indices of categorical predictors.
    %       ContaminationFraction       - Fraction of anomalies in the training data.
    %       NumLearners                 - Number of trees in isolation forest.
    %       NumObservationsPerLearner   - Number of observations used for each tree in isolation forest.
    %       PredictorNames              - Names of predictors used for this model.
    %       ScoreThreshold              - Threshold on anomaly score used to identify anomalies.
    %
    %   IsolationForest methods:
    %       isanomaly                   - Identify anomalies in data.

    %   Copyright 2021-2022 The MathWorks, Inc.

    properties(GetAccess=public,SetAccess=private)
        %NumLearners Number of trees in isolation forest.
        %   The NumLearners property is a positive integer scalar set to the
        %   number of trees in the forest.
        NumLearners

        %NumObservationsPerLearner Number of observations used for each tree in isolation forest.
        %   The NumObservationsPerLearner property is a positive integer scalar
        %   set to the number of observations used to train each isolation tree.
        NumObservationsPerLearner
    end

    properties(GetAccess=private,SetAccess=private)
        % Array of IsolationTree objects.
        Trees

        % 1-by-D logical array for D predictor. This array stores indices
        % of predictors that were not used in training an isolation forest.
        IsBadVariable
    end

    methods
        function [isAnomaly,scores] = isanomaly(obj,X,varargin)
            %ISANOMALY Identify anomalies in data.
            %   ISANOM = ISANOMALY(FOREST,TBL) returns logical array ISANOM with N
            %   elements for N observations in table TBL. TBL contains the predictor
            %   variables. FOREST is an object of class IsolationForest trained on a
            %   table of predictors.
            %
            %   ISANOM = ISANOMALY(FOREST,X) accepts X as a matrix of predictors with
            %   one row per observation and one column per predictor for model FOREST
            %   trained on a predictor matrix. You must use this syntax if you passed a
            %   matrix to IFOREST.
            %
            %   [ISANOM,SCORES] = ISANOMALY(...) also returns numeric array SCORES of
            %   anomaly scores with N elements using any of the input arguments in the
            %   previous syntaxes. The values in this array are between 0 and 1. Value
            %   close to 0 indicates a normal observation, and value close to 1
            %   indicates an anomaly.
            %
            %   ISANOM = ISANOMALY(...,Name=Value) specifies optional parameter
            %   name-value arguments:
            %       ScoreThreshold      - Numeric scalar between 0 and 1. Observations
            %                             with scores above this threshold are
            %                             identified as anomalies. Default: value of
            %                             the ScoreThreshold property
            %       UseParallel         - Logical flag to specify if computation should
            %                             be carried out using a parallel pool.
            %                             Default: false
            %
            %   See also IsolationForest, iforest, ScoreThreshold.

            [isAnomaly,scores] = obj.isanomalyHelper(X,varargin{:});
        end
    end

    methods(Access=protected) % Constructor and score methods
        function obj = IsolationForest(modelParams,dataSummary,trees)

            obj = obj@anomaly.AnomalyDetector(dataSummary.CategoricalPredictors,...
                dataSummary.PredictorNames,dataSummary.WasInputTable,...
                dataSummary.VariableRange,modelParams.ContaminationFraction);

            obj.NumLearners               = modelParams.NumLearners;
            obj.NumObservationsPerLearner = modelParams.NumObservationsPerLearner;
            obj.Trees                     = trees;
            obj.IsBadVariable             = dataSummary.IsBadVariable;
        end

        function s = score(obj,X,varargin)

            [varargin{:}] = convertStringsToChars(varargin{:});

            args = {'useparallel'};
            defs = {        false};
            [usepar] = internal.stats.parseArgs(args,defs,varargin{:});
            validateattributes(usepar,{'logical'},{'scalar'},'isanomaly','UseParallel')

            N = size(X,1);

            trees = obj.Trees;
            numlearners = obj.NumLearners;
            numdraw = double(obj.NumObservationsPerLearner);

            % Estimate depth
            ncats = cellfun(@numel,obj.VariableRange);
            if usepar
                dep = zeros(N,numlearners,'like',X);
                parfor j=1:numlearners
                    dep(:,j) = depth(trees(j),X,ncats);
                end
                dep = sum(dep,2);
            else
                dep = zeros(N,1,'like',X);
                for j=1:numlearners
                    dep = dep + depth(trees(j),X,ncats);
                end
            end
            dep = dep./double(numlearners);

            % Outlier score
            s = pow2(-dep/anomaly.iforest.cfactor(numdraw));
        end
    end

    methods(Hidden)
        function str = toStruct(obj)
            % Convert to a struct for codegen.

            str = toStruct@anomaly.AnomalyDetector(obj);
            str.NumLearners = obj.NumLearners;
            str.NumObservationsPerLearner = obj.NumObservationsPerLearner;
            str.IsBadVariable = obj.IsBadVariable;
            
            trees = num2cell(obj.Trees);
            % collect tree parameters for homogenization
            treeParams = cellfun(@treeCodeGenParameters,trees);
            maxtreeDepth = max(cellfun(@treeDepth,trees));
            % maxtreeDepth is always a power of 2, so this assertion should
            % always pass.
            assert(mod(log2(maxtreeDepth),1)==0,'maxtreeDepth must be a power of 2.');
            treeDepth = char("uint"+maxtreeDepth);
            
            fnames = fieldnames(treeParams);
            maxValues = struct;
            for ii = 1:numel(fnames)
                maxValues.(fnames{ii}) = max([treeParams.(fnames{ii})]);
            end

            str.Trees = cellfun(@(x)toStruct(x,maxValues,treeDepth),trees,'UniformOutput',false);
            
            str.FromStructFcn = 'IsolationForest.fromStruct';
            str.DepthType = treeDepth;

            str.TreeFromStructFcn = str.Trees{1}.FromStructFcn;

        end	
    end
	
    methods(Static,Hidden)
        function obj = fromStruct(str)
            % Reconstruct object from codegen struct

            [dataSummary,modelParams] = anomaly.AnomalyDetector.parametersFromStruct(str);


            trees = cellfun(@anomaly.iforest.IsolationTree.fromStruct,str.Trees);
            threshold = str.ScoreThreshold;
            numtrees = str.NumLearners;
            numobspertree = str.NumObservationsPerLearner;
            badvars = str.IsBadVariable;

            dataSummary.IsBadVariable = badvars;
            modelParams.NumLearners = numtrees;
            modelParams.NumObservationsPerLearner = numobspertree;

            obj = IsolationForest(modelParams,dataSummary,trees);

            obj = setScoreThreshold(obj,threshold);
        end

        function name = matlabCodegenRedirect(~)
            name = 'anomaly.coder.IsolationForest';
        end
    end	

    methods (Access=protected) % Display and validation methods

        function propList = propsForDisplay(obj)
            propList = propsForDisplay@anomaly.AnomalyDetector(obj); % get common properties
            propList.('NumLearners')                = obj.NumLearners;
            propList.('NumObservationsPerLearner')  = obj.NumObservationsPerLearner;
        end

        function X = validateXTest(obj,X)

            X = validateXTest@anomaly.AnomalyDetector(obj,X);

            isbadvar = obj.IsBadVariable;

            if any(isbadvar)
                X(:,isbadvar) = [];
            end
        end
    end

    methods (Static=true,Hidden=true)
        function [mdl,isAnomaly,scores] = fitandscore(X, varargin)

            Dorig = size(X,2);

            % prepareData may change the number of predictors after
            % removing constant predictors
            [X,dataSummary,extraArgs] = ...
                IsolationForest.prepareData(X, varargin{:});

            % Validate iforest parameters
            modelParams = IsolationForest.validateModelParameters(extraArgs{:});

            % Generate samples for trees
            N = size(X,1);
            numdraw = min(N,modelParams.NumObservationsPerLearner);
            modelParams.NumObservationsPerLearner = numdraw;

            % Generate trees
            iscat = false(1,Dorig);
            iscat(dataSummary.CategoricalPredictors) = true;
            iscat(dataSummary.IsBadVariable) = [];

            numlearners = modelParams.NumLearners;
            trees = repmat(anomaly.iforest.IsolationTree,numlearners,1);
            cls = anomaly.utils.inttype(numdraw);
            hasnan = any(isnan(X(:)));
            for j=1:numlearners
                trees(j) = anomaly.iforest.IsolationTree.fit(X,iscat,cls,numdraw,hasnan);
            end

            mdl = IsolationForest(modelParams,dataSummary,trees);

            % Compute scores
            % UseParallel is not stored as a property of mdl
            scores = score(mdl,X,'useparallel',modelParams.UseParallel);

            % Computation of the score threshold and isAnomaly
            [mdl,isAnomaly] = anomaliesInXTrain(mdl,scores);
        end
    end

    methods (Access=protected,Static=true)
        function validateScoreThreshold(thresh)
            validateattributes(thresh,{'single','double'},{'scalar','real','>=',0,'<=',1},'isanomaly','ScoreThreshold')
        end

        function [X,dataSummary,extraArgs] = prepareData(X, varargin)
            % prepareData Prepare data for creating IsolationForest

            % prepareData from base class does not modify number of
            % observation and predictors.
            [X,dataSummary,extraArgs1] = anomaly.AnomalyDetector.prepareData(X,varargin{:},'OrdinalIsCategorical',false);

            internal.stats.checkSupportedNumeric('X',X);

            % Take out ordinal predictors to catch invalid inputs later
            args = {'ordinalpredictors'};
            defs = {                [] };
            [~,~,extraArgs] = internal.stats.parseArgs(args,defs,extraArgs1{:});

            Dorig = size(X,2);
            iscat = false(1,Dorig);
            iscat(dataSummary.CategoricalPredictors) = true;

            % Number of predictors may have changed here
            [X,badvars]    = IsolationForest.removeConstantPredictors(X);
            iscat(badvars) = [];
            dataSummary.IsBadVariable = badvars;

            % Warn about categorical variables with more than 64 levels
            for d=1:size(X,2)
                if iscat(d)
                    if numel(levels(classreg.learning.internal.ClassLabel(X(:,d)))) > 64
                        warning(message('stats:anomaly:iforest:TooManyCategoricalLevels'))
                        break
                    end
                end
            end
        end

        function [X,badvars] = removeConstantPredictors(X)
            % Remove constant predictors

            badvars = min(X)==max(X) | all(isnan(X),1);
            if any(badvars)
                X(:,badvars) = [];
            end
        end

        function modelParams = validateModelParameters(varargin)

            [modelParams,varargin] = anomaly.AnomalyDetector.validateModelParameters('iforest',varargin{:});

            % 1.1 Set parameter defaults.
            dfltNumLearners                 = 100;
            dfltNumObservationsPerLearner   = 256;
            dfltUseParallel                 = false;

            % 1.2 Parse optional name/value pairs.
            paramNames = {'numlearners', 'numobservationsperlearner', 'useparallel'};
            paramDflts = {dfltNumLearners, dfltNumObservationsPerLearner, dfltUseParallel};

            [numlearners,numobslearner,usepar] = internal.stats.parseArgs(paramNames,paramDflts,varargin{:});
            validateattributes(numlearners,{'numeric'},{'scalar','real','positive','integer'},'iforest','NumLearners')
            validateattributes(numobslearner,{'numeric'},{'scalar','real','>',2,'integer'},'iforest','NumObservationsPerLearner')
            usepar = internal.stats.parseOnOff(usepar,'UseParallel');

            % create modelParams structure
            modelParams.NumLearners               = numlearners;
            modelParams.NumObservationsPerLearner = numobslearner;
            modelParams.UseParallel               = usepar;
        end
    end
end

% LocalWords:  ISANOM ordinalpredictors maxtreeDepth
