-
Notifications
You must be signed in to change notification settings - Fork 7
/
Copy pathExWaldMSM.m
124 lines (106 loc) · 5.1 KB
/
ExWaldMSM.m
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
classdef ExWaldMSM < ExWald % Old version superceded by simple convolution which seems much faster
% ExWaldMSM(WaldMean, WaldSD, ExpMean): 3-parameter Ex-Wald distribution, the sum of a 2-parameter Wald, specified in terms of its mean & sd,
% plus an exponential, specified in terms of its mean. The Wald's sigma is 1.
properties(SetAccess = protected)
WaldMean, WaldSD, ExpMean
end
properties(SetAccess = public)
MinWaldMean, MinWaldSD, MinExpMean
StartParmsMLECandidateProportions
end
methods (Static)
function parms = MomsToParms(WaldMean,WaldVar,ExpMean)
% Return values of 3 distribution parameters yielding specified
% mean and variance of normal and mean of exponential.
% Used with ExHelpStartParmsMLE
parms = zeros(3,1);
parms(1) = WaldMean;
parms(2) = sqrt(WaldVar);
parms(3) = ExpMean;
end
end % methods (Static)
methods
function obj=ExWaldMSM(varargin)
obj=obj@ExWald;
obj.FamilyName = 'ExWaldMSM';
obj.NDistParms = 3;
obj.ParmNames{1} = 'WaldMean';
obj.ParmNames{2} = 'WaldSD';
obj.ParmNames{3} = 'ExpMean';
obj.ParmTypes = 'rrr';
obj.DefaultParmCodes = 'rrr';
obj.MinWaldMean = 0.001; % These lower limits are to prevent numerical errors
obj.MinWaldSD = 0.001; % in parameter searching. .01 worked for all.
obj.MinExpMean = 0.001;
obj.WaldBasis.InverseCDFTol = 0.0005; % Tolerance on how close the CDF value must be to the target.
obj.WaldBasish.InverseCDFTol = 0.0005; % Tolerance on how close the CDF value must be to the target.
% This would improve speed but I don't want to hard-code it in; user may call it after making the dist.
% obj.SearchOptions.TolFun = 1e-4;
% obj.SearchOptions.TolX = 1e-4;
NSteps = 10;
obj.StartParmsMLECandidateProportions = ( (1:NSteps) - 0.5) / NSteps;
obj.StartParmsMLECandidateProportions = [0.001 obj.StartParmsMLECandidateProportions 0.999];
switch nargin
case 0
case 3
ResetParms(obj,[varargin{:}]);
otherwise
ME = MException('ExWaldMSM:Constructor', ...
'ExWaldMSM constructor needs 0 or 3 arguments.');
throw(ME);
end
end
function []=ResetParms(obj,newparmvalues)
ClearBeforeResetParmsC(obj);
obj.WaldMean = newparmvalues(1);
obj.WaldSD = newparmvalues(2);
obj.ExpMean = newparmvalues(3);
if obj.WaldMean<=obj.MinWaldMean
error(['ExWaldMSM WaldMean must be > ' num2str(obj.MinWaldMean) '.']);
end
if obj.WaldSD<obj.MinWaldSD
error(['ExWaldMSM WaldSD must be > ' num2str(obj.MinWaldSD) '.']);
end
if obj.ExpMean<obj.MinExpMean
error(['ExWaldMSM ExpMean must be > ' num2str(obj.MinExpMean) '.']);
end
PassA = sqrt(obj.WaldMean)*obj.WaldMean/obj.WaldSD;
PassMu = PassA / obj.WaldMean;
PassRate = 1/obj.ExpMean;
HoldNameBuilding = obj.NameBuilding;
obj.NameBuilding = false;
obj.NDistParms = 4; % Temporarily for ExWald processing.
ResetParms@ExWald(obj,[PassMu 1 PassA PassRate]);
obj.NDistParms = 3; % Return to the true value for ExWaldMSM.
obj.NameBuilding = HoldNameBuilding;
ReInit(obj);
end
function PerturbParms(obj,ParmCodes)
% Perturb parameter values a little bit, e.g., prior to estimation attempts for testing.
newWaldMean = ifelse(ParmCodes(1)=='f',obj.WaldMean,1.05*obj.WaldMean);
newWaldSD = ifelse(ParmCodes(2)=='f',obj.WaldSD, 1.05*obj.WaldSD);
newExpMean = ifelse(ParmCodes(3)=='f',obj.ExpMean, 1.05*obj.ExpMean);
obj.ResetParms([newWaldMean newWaldSD newExpMean]);
end
function []=ReInit(obj)
HoldNameBuilding = obj.NameBuilding;
obj.NameBuilding = false;
ReInit@ExWald(obj);
obj.NameBuilding = HoldNameBuilding;
if (obj.NameBuilding)
BuildMyName(obj);
end
end
function Reals = ParmsToReals(obj,Parms,~)
Reals = [NumTrans.GT2Real(obj.MinWaldMean,Parms(1)) NumTrans.GT2Real(obj.MinWaldSD,Parms(2)) NumTrans.GT2Real(obj.MinExpMean,Parms(3))];
end
function Parms = RealsToParms(obj,Reals,~)
Parms = [NumTrans.Real2GT(obj.MinWaldMean,Reals(1)) NumTrans.Real2GT(obj.MinWaldSD,Reals(2)) NumTrans.Real2GT(obj.MinExpMean,Reals(3))];
end
function parms = StartParmsMLE(obj,X)
HoldParms = obj.ParmValues;
parms = ExHelpStartParmsMLE(obj,X);
obj.ResetParms(HoldParms);
end
end % methods
end % class ExWaldMSM