Jpp test-rotations-old
the software that should make you happy
Loading...
Searching...
No Matches
JMEstimator.hh
Go to the documentation of this file.
1#ifndef __JFIT__JMESTIMATOR__
2#define __JFIT__JMESTIMATOR__
3
4#include <cmath>
5
6#include "JLang/JException.hh"
7
8/**
9 * \file
10 * Maximum likelihood estimator (M-estimators).
11 * \author mdejong
12 */
13namespace JFIT {}
14namespace JPP { using namespace JFIT; }
15
16namespace JFIT {
17
19
20 /**
21 * Interface for maximum likelihood estimator (M-estimator).
22 */
23 struct JMEstimator {
24 /**
25 * Virtual destructor.
26 */
27 virtual ~JMEstimator()
28 {}
29
30
31 /**
32 * Get maximum likelihood estimate.
33 *
34 * \param z deviation
35 * \return likelihood
36 */
37 virtual double getRho(const double z) const = 0;
38
39
40 /**
41 * Get derivative of maximum likelihood estimate.
42 *
43 * \param z deviation
44 * \return derivative
45 */
46 virtual double getPsi(const double z) const = 0;
47 };
48
49
50 /**
51 * Normal M-estimator.
52 *
53 * This estimator is based on a Gaussian PDF and produces the standard chi2.
54 */
56 public JMEstimator
57 {
58 virtual double getRho(const double z) const { return 0.5*z*z; }
59 virtual double getPsi(const double z) const { return z; }
60 };
61
62
63 /**
64 * Lorentzian M-estimator.
65 *
66 * This estimator prouces a logarithmic dependence for large deviations.
67 */
69 public JMEstimator
70 {
71 virtual double getRho(const double z) const { return log (1.0 + 0.5*z*z); }
72 virtual double getPsi(const double z) const { return z / (1.0 + 0.5*z*z); }
73 };
74
75
76 /**
77 * Linear M-estimator.
78 *
79 * This estimator produces a linear dependence for large deviations.
80 */
82 public JMEstimator
83 {
84 virtual double getRho(const double z) const { return sqrt(1.0 + 0.5*z*z) - 1.0; }
85 virtual double getPsi(const double z) const { return 0.5 * z / sqrt(1.0 + 0.5*z*z); }
86 };
87
88
89 /**
90 * Null M-estimator.
91 *
92 * This is not an estimator at all, but an object that just returns whatever it is given.\n
93 * It is introduced so that the user can directly access the likelihood calculated by JRegressor<JEnergy>.
94 */
96 public JMEstimator
97 {
98 virtual double getRho(const double z) const { return z; }
99 virtual double getPsi(const double z) const { return 1.0; }
100 };
101
102
103 /**
104 * Tukey's biweight M-estimator.
105 *
106 * This estimator produces a redescending dependence for large deviations.
107 */
109 public JMEstimator
110 {
111 /**
112 * Constructor.
113 *
114 * \param k standard deviation
115 */
116 JMEstimatorTukey(const double k) :
117 k(k)
118 {}
119
120 virtual double getRho(const double z) const override
121 {
122 const double w = 0.5 * k*k / 3.0;
123
124 if (fabs(z) < k) {
125
126 const double u = z/k;
127 const double v = 1.0 - u*u;
128
129 return w * (1.0 - v*v*v);
130 }
131
132 return w;
133 }
134
135 virtual double getPsi(const double z) const override
136 {
137 if (fabs(z) < k) {
138
139 const double u = z/k;
140 const double v = 1.0 - u*u;
141
142 return z * v*v;
143 }
144
145 return 0.0;
146 }
147
148 double k;
149 };
150
151
152 /**
153 * Normal M-estimator with background.
154 */
156 public JMEstimator
157 {
158 /**
159 * Constructor.
160 *
161 * \param p background probability
162 */
164 p(p)
165 {}
166
167 virtual double getRho(const double z) const override
168 {
169 const double w = exp(-0.5*z*z);
170
171 return -log(w + p);
172 }
173
174 virtual double getPsi(const double z) const override
175 {
176 const double w = exp(-0.5*z*z);
177
178 return z * w / (w + p);
179 }
180
181 double p;
182 };
183
184
185 /**
186 * Definition of the various M-Estimators available to use.
187 */
189 EM_NORMAL = 0,
190 EM_LORENTZIAN = 1,
191 EM_LINEAR = 2,
192 EM_NULL = 3,
195 };
196
197
198 /**
199 * Get M-Estimator.
200 *
201 * Note that for M-estimators with an additional parameter, defaults are used.
202 *
203 * \param type type
204 * \return pointer to newly created M-Estimator (may throw exception)
205 */
206 inline JMEstimator* getMEstimator(const int type)
207 {
208 switch (type) {
209
210 case EM_NORMAL:
211 return new JMEstimatorNormal();
212
213 case EM_LORENTZIAN:
214 return new JMEstimatorLorentzian();
215
216 case EM_LINEAR:
217 return new JMEstimatorLinear();
218
219 case EM_NULL:
220 return new JMEstimatorNull();
221
222 case EM_TUKEY:
223 return new JMEstimatorTukey(5.0);
224
226 return new JMEstimatorNormalWithBackground(1.0e-5);
227
228 default:
229 THROW(JPointerException, "Invalid M-Estimator type " << type);
230 }
231 }
232}
233
234#endif
Exceptions.
#define THROW(JException_t, A)
Marco for throwing exception with std::ostream compatible message.
Exception for accessing an invalid pointer.
Auxiliary classes and methods for linear and iterative data regression.
Definition JEnergy.hh:15
JMEstimator * getMEstimator(const int type)
Get M-Estimator.
JMEstimator_t
Definition of the various M-Estimators available to use.
@ EM_LORENTZIAN
@ EM_NORMALWITHBACKGROUND
@ EM_NORMAL
@ EM_LINEAR
This name space includes all other name spaces (except KM3NETDAQ, KM3NET and ANTARES).
Linear M-estimator.
virtual double getRho(const double z) const
Get maximum likelihood estimate.
virtual double getPsi(const double z) const
Get derivative of maximum likelihood estimate.
Lorentzian M-estimator.
virtual double getPsi(const double z) const
Get derivative of maximum likelihood estimate.
virtual double getRho(const double z) const
Get maximum likelihood estimate.
Normal M-estimator with background.
virtual double getRho(const double z) const override
Get maximum likelihood estimate.
virtual double getPsi(const double z) const override
Get derivative of maximum likelihood estimate.
JMEstimatorNormalWithBackground(const double p)
Constructor.
Normal M-estimator.
virtual double getPsi(const double z) const
Get derivative of maximum likelihood estimate.
virtual double getRho(const double z) const
Get maximum likelihood estimate.
Null M-estimator.
virtual double getRho(const double z) const
Get maximum likelihood estimate.
virtual double getPsi(const double z) const
Get derivative of maximum likelihood estimate.
Tukey's biweight M-estimator.
virtual double getRho(const double z) const override
Get maximum likelihood estimate.
virtual double getPsi(const double z) const override
Get derivative of maximum likelihood estimate.
JMEstimatorTukey(const double k)
Constructor.
Interface for maximum likelihood estimator (M-estimator).
virtual double getRho(const double z) const =0
Get maximum likelihood estimate.
virtual double getPsi(const double z) const =0
Get derivative of maximum likelihood estimate.
virtual ~JMEstimator()
Virtual destructor.