Cupid/HyperbolicTan.m at master · xingzhe996/Cupid · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
classdef HyperbolicTan < dContinuous
    % Hyperbolic Tangent distribution with parameter Scale>0.
    %  See Strasburger (2002, Perception & Psychophysics.}

    properties(SetAccess = protected)
        Scale
    end

    methods (Static)

        function Reals = ParmsToReals(Parms,~)
            Reals = NumTrans.GT2Real(eps,Parms(1));
        end

        function Parms = RealsToParms(Reals,~)
            Parms = NumTrans.Real2GT(eps,Reals(1));
        end

    end

    methods

        function obj=HyperbolicTan(varargin)
            obj=obj@dContinuous('HyperbolicTan');
            obj.NDistParms = 1;
            obj.ParmTypes = 'r';
            obj.DefaultParmCodes = 'r';
            switch nargin
                case 0
                case 1
                    ResetParms(obj,varargin{1});
                otherwise
                    ME = MException('HyperbolicTan:Constructor', ...
                        'HyperbolicTan constructor needs 0 or 1 arguments.');
                    throw(ME);
            end
        end

        function []=ResetParms(obj,newparmvalues)
            ClearBeforeResetParmsC(obj);
            obj.Scale = newparmvalues(1);
            ReInit(obj);
        end

        function PerturbParms(obj,ParmCodes)
            % Perturb parameter values prior to estimation attempts.
            newScale = ifelse(ParmCodes(1)=='f', obj.Scale, 0.9*obj.Scale);
            obj.ResetParms(newScale);
        end

        function []=ReInit(obj)
            assert(obj.Scale>0,'HyperbolicTan Scale parameter must be > 0.');
            obj.Initialized = true;
            obj.LowerBound = InverseCDF(obj,obj.CDFNearlyZero);
            obj.UpperBound = InverseCDF(obj,obj.CDFNearlyOne);
            if (obj.NameBuilding)
                BuildMyName(obj);
            end
        end

        function thispdf=PDF(obj,X)
            [thispdf, InBounds, Done] = MaybeSplinePDF(obj,X);
            if Done
                return;
            end
            ExpPlus = exp(obj.Scale*X(InBounds));
            ExpMinus = 1 ./ ExpPlus;
            thispdf(InBounds) = 4 * obj.Scale ./ (ExpPlus + ExpMinus).^2;
            % for i=1:numel(X)
            %     if (X(i)>=obj.LowerBound) && (X(i)<=obj.UpperBound)
            %         ExpPlus = exp(obj.Scale*X(i));
            %         ExpMinus = 1/ExpPlus;
            %         thispdf(i) = 4 * obj.Scale / (ExpPlus + ExpMinus)^2;
            %     end
            % end
        end

        function thiscdf=CDF(obj,X)
            [thiscdf, InBounds, Done] = MaybeSplineCDF(obj,X);
            if Done
                return;
            end
            thiscdf(X>=obj.UpperBound) = 1;
            ExpPlus = exp(obj.Scale*X(InBounds));
            ExpMinus = 1 ./ ExpPlus;
            thiscdf(InBounds) = (ExpPlus-ExpMinus) ./ (ExpPlus+ExpMinus);
            % for i=1:numel(X)
            %     if X(i)<=obj.LowerBound
            %     elseif X(i)>=obj.UpperBound
            %         thiscdf(i) = 1;
            %     else
            %         ExpPlus = exp(obj.Scale*X(i));
            %         ExpMinus = 1/ExpPlus;
            %         thiscdf(i) = (ExpPlus - ExpMinus) / (ExpPlus+ExpMinus);
            %     end
            % end
        end

        function thisval=InverseCDF(obj,P)  % NWJEFF: Vectorize InverseCDFs. InBounds?
            [thisval, InBounds, Done] = MaybeSplineInvCDF(obj,P);
            if Done
                return;
            end
            OMP = 1 - P;
            S = sqrt(OMP.*(1+P)) ./ OMP;
            thisval = log(S) / obj.Scale;
        end

        function thisval=Mean(obj)
            if ~obj.Initialized
                error(UninitializedError(obj));
            end
            thisval = log(2) / obj.Scale;
        end

        function thisval=Variance(obj)
            if ~obj.Initialized
                error(UninitializedError(obj));
            end
            EXSqr = (pi/obj.Scale)^2 / 12.0;
            thisval = EXSqr - Mean(obj)^2;
        end

        function thisval=Median(obj)
            if ~obj.Initialized
                error(UninitializedError(obj));
            end
            thisval = 0.5493 / obj.Scale;
        end

    end  % methods

end  % class HyperbolicTan