Home  · Classes  · Annotated Classes  · Modules  · Members  · Namespaces  · Related Pages
LinearRegression.h
Go to the documentation of this file.
1 // --------------------------------------------------------------------------
2 // OpenMS -- Open-Source Mass Spectrometry
3 // --------------------------------------------------------------------------
4 // Copyright The OpenMS Team -- Eberhard Karls University Tuebingen,
5 // ETH Zurich, and Freie Universitaet Berlin 2002-2013.
6 //
7 // This software is released under a three-clause BSD license:
8 // * Redistributions of source code must retain the above copyright
9 // notice, this list of conditions and the following disclaimer.
10 // * Redistributions in binary form must reproduce the above copyright
11 // notice, this list of conditions and the following disclaimer in the
12 // documentation and/or other materials provided with the distribution.
13 // * Neither the name of any author or any participating institution
14 // may be used to endorse or promote products derived from this software
15 // without specific prior written permission.
16 // For a full list of authors, refer to the file AUTHORS.
17 // --------------------------------------------------------------------------
18 // THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
19 // AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
20 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
21 // ARE DISCLAIMED. IN NO EVENT SHALL ANY OF THE AUTHORS OR THE CONTRIBUTING
22 // INSTITUTIONS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
23 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
24 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS;
25 // OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY,
26 // WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR
27 // OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF
28 // ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
29 //
30 // --------------------------------------------------------------------------
31 // $Maintainer: Clemens Groepl $
32 // $Authors: $
33 // --------------------------------------------------------------------------
34 
35 #ifndef OPENMS_MATH_STATISTICS_LINEARREGRESSION_H
36 #define OPENMS_MATH_STATISTICS_LINEARREGRESSION_H
37 
38 #include <OpenMS/CONCEPT/Types.h>
40 
41 #include <iostream>
42 #include <vector>
43 #include <gsl/gsl_fit.h>
44 #include <gsl/gsl_statistics.h>
45 #include <gsl/gsl_cdf.h>
46 
47 namespace OpenMS
48 {
49  namespace Math
50  {
69  class OPENMS_DLLAPI LinearRegression
70  {
71 public:
72 
75  intercept_(0),
76  slope_(0),
77  x_intercept_(0),
78  lower_(0),
79  upper_(0),
80  t_star_(0),
81  r_squared_(0),
82  stand_dev_residuals_(0),
83  mean_residuals_(0),
84  stand_error_slope_(0),
85  chi_squared_(0),
86  rsd_(0)
87  {
88  }
89 
92  {
93  }
94 
111  template <typename Iterator>
112  void computeRegression(double confidence_interval_P, Iterator x_begin, Iterator x_end, Iterator y_begin);
113 
129  template <typename Iterator>
130  void computeRegressionNoIntercept(double confidence_interval_P, Iterator x_begin, Iterator x_end, Iterator y_begin);
131 
148  template <typename Iterator>
149  void computeRegressionWeighted(double confidence_interval_P, Iterator x_begin, Iterator x_end, Iterator y_begin, Iterator w_begin);
150 
151 
153  DoubleReal getIntercept() const;
155  DoubleReal getSlope() const;
157  DoubleReal getXIntercept() const;
159  DoubleReal getLower() const;
161  DoubleReal getUpper() const;
163  DoubleReal getTValue() const;
165  DoubleReal getRSquared() const;
167  DoubleReal getStandDevRes() const;
169  DoubleReal getMeanRes() const;
171  DoubleReal getStandErrSlope() const;
173  DoubleReal getChiSquared() const;
175  DoubleReal getRSD() const;
176 
177 protected:
178 
180  double intercept_;
182  double slope_;
184  double x_intercept_;
186  double lower_;
188  double upper_;
190  double t_star_;
192  double r_squared_;
200  double chi_squared_;
202  double rsd_;
203 
204 
206  void computeGoodness_(double * X, double * Y, int N, double confidence_interval_P);
207 
209  template <typename Iterator>
210  void iteratorRange2Arrays_(Iterator x_begin, Iterator x_end, Iterator y_begin, double * x_array, double * y_array);
211 
213  template <typename Iterator>
214  void iteratorRange3Arrays_(Iterator x_begin, Iterator x_end, Iterator y_begin, Iterator w_begin, double * x_array, double * y_array, double * w_array);
215 
216 private:
217 
219  LinearRegression(const LinearRegression & arg);
221  LinearRegression & operator=(const LinearRegression & arg);
222  };
223 
224  template <typename Iterator>
225  void LinearRegression::computeRegression(double confidence_interval_P, Iterator x_begin, Iterator x_end, Iterator y_begin)
226  {
227  int N = int(distance(x_begin, x_end));
228 
229  double * X = new double[N];
230  double * Y = new double[N];
231  iteratorRange2Arrays_(x_begin, x_end, y_begin, X, Y);
232 
233  double cov00, cov01, cov11;
234 
235  // Compute the unweighted linear fit.
236  // Get the intercept and the slope of the regression Y_hat=intercept_+slope_*X
237  // and the value of Chi squared, the covariances of the intercept and the slope
238  int error = gsl_fit_linear(X, 1, Y, 1, N, &intercept_, &slope_, &cov00, &cov01, &cov11, &chi_squared_);
239 
240  if (!error)
241  {
242  computeGoodness_(X, Y, N, confidence_interval_P);
243  }
244 
245  delete[] X;
246  delete[] Y;
247 
248  if (error)
249  {
250  throw Exception::UnableToFit(__FILE__, __LINE__, __PRETTY_FUNCTION__, "UnableToFit-LinearRegression", "Could not fit a linear model to the data");
251  }
252  }
253 
254  template <typename Iterator>
255  void LinearRegression::computeRegressionNoIntercept(double confidence_interval_P, Iterator x_begin, Iterator x_end, Iterator y_begin)
256  {
257  int N = int(distance(x_begin, x_end));
258 
259  double * X = new double[N];
260  double * Y = new double[N];
261  iteratorRange2Arrays_(x_begin, x_end, y_begin, X, Y);
262 
263  double cov;
264 
265  // Compute the linear fit.
266  // Get the intercept and the slope of the regression Y_hat=intercept_+slope_*X
267  // and the value of Chi squared, the covariances of the intercept and the slope
268  int error = gsl_fit_mul(X, 1, Y, 1, N, &slope_, &cov, &chi_squared_);
269  intercept_ = 0.0;
270 
271  if (!error)
272  {
273  computeGoodness_(X, Y, N, confidence_interval_P);
274  }
275 
276  delete[] X;
277  delete[] Y;
278 
279  if (error)
280  {
281  throw Exception::UnableToFit(__FILE__, __LINE__, __PRETTY_FUNCTION__, "UnableToFit-LinearRegression", "Could not fit a linear model to the data");
282  }
283  }
284 
285  template <typename Iterator>
286  void LinearRegression::computeRegressionWeighted(double confidence_interval_P, Iterator x_begin, Iterator x_end, Iterator y_begin, Iterator w_begin)
287  {
288  int N = int(distance(x_begin, x_end));
289 
290  double * X = new double[N];
291  double * Y = new double[N];
292  double * W = new double[N];
293  iteratorRange3Arrays_(x_begin, x_end, y_begin, w_begin, X, Y, W);
294 
295  double cov00, cov01, cov11;
296 
297  // Compute the weighted linear fit.
298  // Get the intercept and the slope of the regression Y_hat=intercept_+slope_*X
299  // and the value of Chi squared, the covariances of the intercept and the slope
300  int error = gsl_fit_wlinear(X, 1, W, 1, Y, 1, N, &intercept_, &slope_, &cov00, &cov01, &cov11, &chi_squared_);
301 
302  if (!error)
303  {
304  computeGoodness_(X, Y, N, confidence_interval_P);
305  }
306 
307  delete[] X;
308  delete[] Y;
309  delete[] W;
310 
311  if (error)
312  {
313  throw Exception::UnableToFit(__FILE__, __LINE__, __PRETTY_FUNCTION__, "UnableToFit-LinearRegression", "Could not fit a linear model to the data");
314  }
315  }
316 
317  template <typename Iterator>
318  void LinearRegression::iteratorRange2Arrays_(Iterator x_begin, Iterator x_end, Iterator y_begin, double * x_array, double * y_array)
319  {
320  int i = 0;
321  while (x_begin < x_end)
322  {
323  x_array[i] = *x_begin;
324  y_array[i] = *y_begin;
325  ++x_begin;
326  ++y_begin;
327  ++i;
328  }
329  }
330 
331  template <typename Iterator>
332  void LinearRegression::iteratorRange3Arrays_(Iterator x_begin, Iterator x_end, Iterator y_begin, Iterator w_begin, double * x_array, double * y_array, double * w_array)
333  {
334  int i = 0;
335  while (x_begin < x_end)
336  {
337  x_array[i] = *x_begin;
338  y_array[i] = *y_begin;
339  w_array[i] = *w_begin;
340  ++x_begin;
341  ++y_begin;
342  ++w_begin;
343  ++i;
344  }
345  }
346 
347  } // namespace Math
348 } // namespace OpenMS
349 
350 
351 #endif
virtual ~LinearRegression()
Destructor.
Definition: LinearRegression.h:91
LinearRegression()
Constructor.
Definition: LinearRegression.h:74
double rsd_
the relative standard deviation
Definition: LinearRegression.h:202
double upper_
The upper bound of the confidence intervall.
Definition: LinearRegression.h:188
double chi_squared_
The value of the Chi Squared statistic.
Definition: LinearRegression.h:200
void computeGoodness_(double *X, double *Y, int N, double confidence_interval_P)
Computes the goodness of the fitted regression line.
void computeRegressionWeighted(double confidence_interval_P, Iterator x_begin, Iterator x_end, Iterator y_begin, Iterator w_begin)
This function computes the best-fit linear regression coefficients of the model for the weighted da...
Definition: LinearRegression.h:286
This class offers functions to perform least-squares fits to a straight line model, .
Definition: LinearRegression.h:69
double x_intercept_
The intercept of the fitted line with the x-axis.
Definition: LinearRegression.h:184
double t_star_
The value of the t-statistic.
Definition: LinearRegression.h:190
double stand_error_slope_
The standard error of the slope.
Definition: LinearRegression.h:198
double mean_residuals_
Mean of residuals.
Definition: LinearRegression.h:196
double slope_
The slope of the fitted line.
Definition: LinearRegression.h:182
void iteratorRange3Arrays_(Iterator x_begin, Iterator x_end, Iterator y_begin, Iterator w_begin, double *x_array, double *y_array, double *w_array)
Copy the distance(x_begin,x_end) elements starting at x_begin, y_begin and w_begin into the arrays x_...
Definition: LinearRegression.h:332
Exception used if an error occurred while fitting a model to a given dataset.
Definition: Exception.h:662
void computeRegressionNoIntercept(double confidence_interval_P, Iterator x_begin, Iterator x_end, Iterator y_begin)
This function computes the best-fit linear regression coefficient of the model for the dataset ...
Definition: LinearRegression.h:255
double stand_dev_residuals_
The standard deviation of the residuals.
Definition: LinearRegression.h:194
void iteratorRange2Arrays_(Iterator x_begin, Iterator x_end, Iterator y_begin, double *x_array, double *y_array)
Copies the distance(x_begin,x_end) elements starting at x_begin and y_begin into the arrays x_array a...
Definition: LinearRegression.h:318
double lower_
The lower bound of the confidence intervall.
Definition: LinearRegression.h:186
double r_squared_
The squared correlation coefficient (Pearson)
Definition: LinearRegression.h:192
double intercept_
The intercept of the fitted line with the y-axis.
Definition: LinearRegression.h:180
void computeRegression(double confidence_interval_P, Iterator x_begin, Iterator x_end, Iterator y_begin)
This function computes the best-fit linear regression coefficients of the model for the dataset ...
Definition: LinearRegression.h:225

OpenMS / TOPP release 1.11.1 Documentation generated on Thu Nov 14 2013 11:19:16 using doxygen 1.8.5