From 41d09392ca78ae2a36b93fddeab5a58b818f3b9c Mon Sep 17 00:00:00 2001
From: Remco Renken <r.j.renken@umcg.nl>
Date: Thu, 4 Apr 2019 08:35:34 +0200
Subject: [PATCH] First upload of code

---
 CREClass.m     | 120 +++++++++++++++++++++++++++++++++++++++++++++++++
 TestCREClass.m |  22 +++++++++
 2 files changed, 142 insertions(+)
 create mode 100644 CREClass.m
 create mode 100644 TestCREClass.m

diff --git a/CREClass.m b/CREClass.m
new file mode 100644
index 0000000..eeb698d
--- /dev/null
+++ b/CREClass.m
@@ -0,0 +1,120 @@
+classdef CREClass < handle
+    properties
+        Data; % the  data set to get the CRE off
+        nBin; %number of bins for the reference
+        Monitor=true; %set to true to get intermediat results/figures
+        UseHistProxy=true; % set this to get a proxy of the survival function using the hist function
+    end
+    properties %to be made read only
+        CRE
+        RB
+        P_RB
+    end
+    properties (Hidden)
+        WaitBarHandle;
+    end
+    properties (Constant,Hidden)
+        alpha=1; %needed for the scaled cummulative redual entropy as defined in zografoss. Will be implemented in the future. (I hope)
+        beta=1; %they are currently obsolete.
+        %        FeedBackOptions={'CRE','RB','All'}
+    end
+    %% constructor
+    methods
+        function obj = CREClass(varargin)
+            obj.Data=[];
+            obj.nBin=[];
+            %obj.FeedBack={'CRE'};
+            if nargin==0;return;end % if no input arguments are given, just create an empty instance
+            if ischar(varargin{1}) %assume number value pairs.
+                if nargin<=2; error('please use number value pairs' ); end
+                for k=2:2:nargin
+                    switch lower(varargin{k-1})
+                        case {'ref','data'}
+                            obj.Data=varargin{k};
+                        case {'nbin','bin'}
+                            obj.nBin=varargin{k};
+                        otherwise
+                            error('invalid variable name');
+                    end
+                end
+            else
+                if nargin>=1
+                    obj.Data=varargin{1};
+                end
+                if nargin==2
+                    obj.nBin=varargin{2};
+                end
+            end
+        end
+    end
+    %% set functions
+    methods
+        function set.Data(obj,val)
+            if isnumeric(val)
+                obj.Data=val;
+            else
+                error('please enter numeric array/matrix as reference');
+            end
+        end
+        function set.nBin(obj,val)
+            if isnumeric(val)
+                obj.nBin=val;
+            else
+                error('please enter numeric array or value as bin definition');
+            end
+        end
+    end
+    %% core function(s)
+    methods
+        
+        function Calc(obj)
+            if isempty(obj.Data)
+                warning('nothing to do');
+                return;
+            end
+            if obj.UseHistProxy
+                if isempty(obj.nBin)
+                    obj.nBin=floor(numel(obj.Data)./10);
+                end
+                [N,bin]=hist(obj.Data(:),obj.nBin);
+                N=transpose(N);
+                P=N./sum(N);
+                CP=cumsum(P);
+                FC=1-CP;
+                FC(FC<0)=0;
+                dl=ones(size(FC))./numel(FC);
+            else
+                q=sort(unique(obj.Data(:)));
+                FC=zeros(size(q));
+                dl=zeros(size(q));
+                for k=1:numel(q)
+                    FC(k)=sum(obj.Data(:)>=q(k));
+                    if k==1
+                        dl(k)=(q(k+1)-q(k))/2; %should I divide by 2?
+                    elseif k==numel(q)
+                        dl(k)=(q(k)-q(k-1))/2;
+                    else
+                        dl(k)=(q(k+1)-q(k))/2+(q(k)-q(k-1))/2;
+                    end
+                end
+                FC=FC./numel(obj.Data);
+            end
+            LogFC=log(FC);
+            LogFC(~isfinite(LogFC))=0;%the log of 0 is -inf. however in Fc.*logFc it should end up as 0. to avoid conflicts removing the -inf
+            if any(isnan(FC));keyboard;end
+            obj.CRE=-transpose(dl)*(FC.*LogFC)./(transpose(dl)*FC); %CRE zografos 
+            %% get the RB
+            dl(FC>0.5)=0; %set the weight for all in the histogram to the left side (i.e. <p50) to 0
+            CRE_Med_Inf=-transpose(dl)*(FC.*LogFC)./(transpose(dl)*FC); %CRE zografos 
+            ind=find(exp(-CRE_Med_Inf)>=FC,1,'first');
+            if obj.UseHistProxy
+                obj.P_RB=FC(ind);
+                obj.RB=bin(ind);
+            else
+                SortData=sort(obj.Data);
+                obj.RB=SortData(ind);
+                obj.P_RB=FC(ind);
+            end
+        end
+    end
+end
\ No newline at end of file
diff --git a/TestCREClass.m b/TestCREClass.m
new file mode 100644
index 0000000..91c4a72
--- /dev/null
+++ b/TestCREClass.m
@@ -0,0 +1,22 @@
+S=CREClass;
+S.Data=randn(100000,1);
+S.nBin=(1000);
+%% use histogram
+S.UseHistProxy=true;
+S.Calc;
+S
+%% Use my aproximation
+S.UseHistProxy=false;
+S.Calc;
+S
+std(S.Data)
+%% multiply S
+S.Data=100*S.Data;
+S.UseHistProxy=true;
+S.Calc
+S
+std(S.Data)
+%% ones again the "new" way
+S.UseHistProxy=false;
+S.Calc
+S
\ No newline at end of file