source:trunk/lib/statistics/Linear.cc@430

Last change on this file since 430 was 430, checked in by Peter, 18 years ago

changed interface of Regression::Local

• Property svn:eol-style set to `native`
• Property svn:keywords set to `Author Date Id Revision`
File size: 1.4 KB
Line
1// \$Id: Linear.cc 430 2005-12-08 22:53:08Z peter \$
2
3#include <c++_tools/statistics/Linear.h>
4
5#include <c++_tools/statistics/AveragerPair.h>
6#include <c++_tools/gslapi/vector.h>
7
8#include <gsl/gsl_fit.h>
9
10
11namespace theplu {
12namespace statistics {
13namespace regression {
14
15
16  void Linear::fit(const gslapi::vector& x, const gslapi::vector& y)
17  {
18    ap_.reset();
19    for (size_t i=0; i<x.size(); i++)
21
22    alpha_ = ap_.y_averager().mean();
23    beta_ = ap_.covariance() / ap_.x_averager().variance();
24
25    // estimating the noise level, i.e. the conditional variance of y
26    // given x, Var(y|x).
27    double Q = (ap_.y_averager().sum_xsqr_centered() - ap_.sum_xy_centered() *
28                ap_.sum_xy_centered()/ap_.x_averager().sum_xsqr_centered() ) ;
29    s2_ = Q/(x.size()-2);
30    r2_= 1;
31    alpha_var_ = s2_ / x.size();
32    beta_var_ = s2_ / ap_.x_averager().sum_xsqr_centered();
33    m_x_ = ap_.x_averager().mean();
34  }
35
36  void Linear::predict(const double x, double& y, double& y_err) const
37  {
38    y = alpha_ + beta_ * (x-m_x_);
39    y_err = sqrt( alpha_var_+beta_var_*(x-m_x_)*(x-m_x_)+s2_ );
40  }
41
43  {
44    s << "# column 1: x\n"
45      << "# column 2: y\n"
46      << "# column 3: y_err\n"
47      << "# column 4: alpha\n"
48      << "# column 5: alpha_err\n"
49      << "# column 6: beta\n"
50      << "# column 7: beta_err\n"
51      << "# column 8: s2 (var(y|x))\n"
52      << "# column 9: r2 (coefficient of determination)";
53    return s;
54  }
55
56
57}}} // of namespaces regression, statisitcs and thep
Note: See TracBrowser for help on using the repository browser.