source: trunk/src/RegressionLocal.cc @ 245

Last change on this file since 245 was 245, checked in by Peter, 17 years ago

fixed major bug

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 2.5 KB
Line 
1// $Id: RegressionLocal.cc 245 2005-02-22 15:21:18Z peter $
2
3#include "RegressionLocal.h"
4
5#include "Regression.h"
6#include "RegressionKernel.h"
7#include "RegressionLinear.h"
8#include "vector.h"
9
10#include "algorithm"
11
12namespace theplu {
13namespace statistics {
14
15  RegressionLocal::RegressionLocal(Regression& r, 
16                                   RegressionKernel& k)
17    : kernel_(&k), regressor_(&r)
18  {
19  }
20
21  void RegressionLocal::fit(const double f, const u_int step_size)
22  {
23    sort(data_.begin(), data_.end());
24
25    // number of points on each side
26    size_t nof_points = static_cast<size_t>(f*data_.size());
27    double width = 0;
28    int min_index = 0;
29    for (size_t i=0; i<data_.size(); i+=step_size) {
30      // determining boundaries for window/kernel
31     
32      // extreme left case
33      if (data_[i].first <
34          (data_[nof_points-1].first+data_[nof_points-1].first)/2 ){
35        min_index=0;
36        width = data_[nof_points-1].first-data_[i].first;
37      }
38      // extreme right case
39      else if (data_[i].first > (data_[data_.size()-1].first+
40                            data_[data_.size()-nof_points].first)/2 ){
41        min_index = data_.size() - nof_points;
42        width = data_[i].first - data_[data_.size()-nof_points].first;
43      }
44      else {
45        // stepping forward until x is in the middle
46        while(data_[min_index+nof_points-1].first-data_[i].first
47              < data_[i].first-data_[min_index].first)
48          min_index++;
49        width = data_[min_index+nof_points-1].first-data_[i].first;
50        // Peter, should check if x gets closer to centrum if we step
51        // back one step
52      }
53
54
55      // copying data
56      // Peter, too much copying. Move the copying outside loop and
57      // use a vector view (when available in gslapi::vector class).
58      gslapi::vector x(nof_points);
59      gslapi::vector y(nof_points);
60      gslapi::vector w(nof_points);
61      for (size_t j=0; j<x.size(); j++)
62        x(j)=data_[min_index+j].first;
63      for (size_t j=0; j<y.size(); j++){
64        y(j)=data_[min_index+j].second;
65      }
66      // calculating weights
67      for (size_t j=0; j<y.size(); j++)
68        w(j) = kernel_->weight( (x(j)-(data_[min_index].first+
69                                       data_[min_index+nof_points-1].first) /2)
70                                /width );
71     
72      // fitting the regressor locally
73      regressor_->fit(x,y,w);
74      x_.push_back(data_[i].first);
75      y_.push_back(0.0);
76      y_err_.push_back(0.0);
77      regressor_->predict(x_[x_.size()-1], y_[y_.size()-1], 
78                          y_err_[y_err_.size()-1]);
79    }
80  }
81
82  std::ostream& RegressionLocal::print(std::ostream& s) const
83  {
84    for (size_t i=0; i<x_.size(); i++) {
85      regressor_->print(s);
86      s << std::endl;
87    }
88    return s;
89  }
90
91
92}} // of namespace cpptools and namespace theplu
Note: See TracBrowser for help on using the repository browser.