First upload of code

This commit is contained in:
R.J. Renken 2019-04-04 08:35:34 +02:00
commit 41d09392ca
2 changed files with 142 additions and 0 deletions

120
CREClass.m Normal file
View File

@ -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

22
TestCREClass.m Normal file
View File

@ -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