yat/statistics/AveragerPair.cc

Code
Comments
Other
Rev Date Author Line
138 20 Aug 04 peter 1 // $Id$
138 20 Aug 04 peter 2
675 10 Oct 06 jari 3 /*
2119 12 Dec 09 peter 4   Copyright (C) 2004 Jari Häkkinen, Peter Johansson
1746 23 Jan 09 peter 5   Copyright (C) 2005 Peter Johansson
4359 23 Aug 23 peter 6   Copyright (C) 2006 Jari Häkkinen, Peter Johansson
4359 23 Aug 23 peter 7   Copyright (C) 2007 Peter Johansson
4359 23 Aug 23 peter 8   Copyright (C) 2008 Jari Häkkinen, Peter Johansson
2919 19 Dec 12 peter 9   Copyright (C) 2011, 2012 Peter Johansson
138 20 Aug 04 peter 10
1437 25 Aug 08 peter 11   This file is part of the yat library, http://dev.thep.lu.se/yat
675 10 Oct 06 jari 12
675 10 Oct 06 jari 13   The yat library is free software; you can redistribute it and/or
675 10 Oct 06 jari 14   modify it under the terms of the GNU General Public License as
1486 09 Sep 08 jari 15   published by the Free Software Foundation; either version 3 of the
675 10 Oct 06 jari 16   License, or (at your option) any later version.
675 10 Oct 06 jari 17
675 10 Oct 06 jari 18   The yat library is distributed in the hope that it will be useful,
675 10 Oct 06 jari 19   but WITHOUT ANY WARRANTY; without even the implied warranty of
675 10 Oct 06 jari 20   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
675 10 Oct 06 jari 21   General Public License for more details.
675 10 Oct 06 jari 22
675 10 Oct 06 jari 23   You should have received a copy of the GNU General Public License
1487 10 Sep 08 jari 24   along with yat. If not, see <http://www.gnu.org/licenses/>.
675 10 Oct 06 jari 25 */
675 10 Oct 06 jari 26
2881 18 Nov 12 peter 27 #include <config.h>
2881 18 Nov 12 peter 28
680 11 Oct 06 jari 29 #include "AveragerPair.h"
680 11 Oct 06 jari 30 #include "Averager.h"
675 10 Oct 06 jari 31
2565 26 Sep 11 peter 32 #include <limits>
138 20 Aug 04 peter 33 #include <utility>
138 20 Aug 04 peter 34
138 20 Aug 04 peter 35 namespace theplu {
680 11 Oct 06 jari 36 namespace yat {
683 11 Oct 06 jari 37 namespace statistics {
703 18 Dec 06 jari 38
703 18 Dec 06 jari 39   AveragerPair::AveragerPair(void)
2565 26 Sep 11 peter 40     : xy_centered_(0.0)
703 18 Dec 06 jari 41   {
703 18 Dec 06 jari 42   }
703 18 Dec 06 jari 43
2565 26 Sep 11 peter 44
703 18 Dec 06 jari 45   AveragerPair::AveragerPair(const AveragerPair& a)
2565 26 Sep 11 peter 46     : x_(a.x_), y_(a.y_), xy_centered_(a.xy_centered_)
703 18 Dec 06 jari 47   {
703 18 Dec 06 jari 48   }
703 18 Dec 06 jari 49
2565 26 Sep 11 peter 50
1295 12 May 08 jari 51   void AveragerPair::add(const double x, const double y, const long n)
703 18 Dec 06 jari 52   {
2565 26 Sep 11 peter 53     xy_add(x, y, 0, n);
2644 11 Nov 11 peter 54     x_.add(x,n);
2644 11 Nov 11 peter 55     y_.add(y,n);
703 18 Dec 06 jari 56   }
703 18 Dec 06 jari 57
2565 26 Sep 11 peter 58
718 26 Dec 06 jari 59   double AveragerPair::ccc(void) const
718 26 Dec 06 jari 60   {
1122 22 Feb 08 peter 61     return ((2*covariance()) /
1122 22 Feb 08 peter 62             ((x_.variance()+y_.variance()) +
1122 22 Feb 08 peter 63              (x_.mean()-y_.mean())*(x_.mean()-y_.mean())));
718 26 Dec 06 jari 64   }
718 26 Dec 06 jari 65
2565 26 Sep 11 peter 66
718 26 Dec 06 jari 67   double AveragerPair::correlation(void) const
2644 11 Nov 11 peter 68   {
2644 11 Nov 11 peter 69     return sum_xy_centered() /
2566 26 Sep 11 peter 70       std::sqrt(x_.sum_xx_centered() * y_.sum_xx_centered());
718 26 Dec 06 jari 71   }
718 26 Dec 06 jari 72
2565 26 Sep 11 peter 73
718 26 Dec 06 jari 74   double AveragerPair::covariance(void) const
718 26 Dec 06 jari 75   {
2565 26 Sep 11 peter 76     return xy_centered_ / n();
718 26 Dec 06 jari 77   }
718 26 Dec 06 jari 78
2565 26 Sep 11 peter 79
718 26 Dec 06 jari 80   double AveragerPair::mean_xy(void) const
718 26 Dec 06 jari 81   {
2565 26 Sep 11 peter 82     if (!n())
2565 26 Sep 11 peter 83       return std::numeric_limits<double>::quiet_NaN();
2565 26 Sep 11 peter 84     return sum_xy()/n();
718 26 Dec 06 jari 85   }
718 26 Dec 06 jari 86
2565 26 Sep 11 peter 87
718 26 Dec 06 jari 88   double AveragerPair::msd(void) const
718 26 Dec 06 jari 89   {
772 25 Feb 07 peter 90     return sum_squared_deviation()/n();
718 26 Dec 06 jari 91   }
718 26 Dec 06 jari 92
2565 26 Sep 11 peter 93
1295 12 May 08 jari 94   long AveragerPair::n(void) const
718 26 Dec 06 jari 95   {
718 26 Dec 06 jari 96     return x_.n();
718 26 Dec 06 jari 97   }
718 26 Dec 06 jari 98
2565 26 Sep 11 peter 99
703 18 Dec 06 jari 100   void AveragerPair::reset(void)
703 18 Dec 06 jari 101   {
2644 11 Nov 11 peter 102     x_.reset();
2644 11 Nov 11 peter 103     y_.reset();
2565 26 Sep 11 peter 104     xy_centered_ = 0.0;
703 18 Dec 06 jari 105   }
703 18 Dec 06 jari 106
2565 26 Sep 11 peter 107
703 18 Dec 06 jari 108   const AveragerPair& AveragerPair::operator=(const AveragerPair& a)
703 18 Dec 06 jari 109   {
2565 26 Sep 11 peter 110     x_=a.x_; y_=a.y_; xy_centered_=a.xy_centered_;
703 18 Dec 06 jari 111     return *this;
703 18 Dec 06 jari 112   }
703 18 Dec 06 jari 113
2565 26 Sep 11 peter 114
772 25 Feb 07 peter 115   double AveragerPair::sum_squared_deviation(void) const
772 25 Feb 07 peter 116   {
772 25 Feb 07 peter 117     return x_averager().sum_xx()+y_averager().sum_xx()-2*sum_xy() ;
772 25 Feb 07 peter 118   }
772 25 Feb 07 peter 119
2565 26 Sep 11 peter 120
718 26 Dec 06 jari 121   double AveragerPair::sum_xy(void) const
718 26 Dec 06 jari 122   {
2565 26 Sep 11 peter 123     return xy_centered_ + x_.sum_x()*y_.mean();;
718 26 Dec 06 jari 124   }
718 26 Dec 06 jari 125
2565 26 Sep 11 peter 126
718 26 Dec 06 jari 127   double AveragerPair::sum_xy_centered(void) const
718 26 Dec 06 jari 128   {
2565 26 Sep 11 peter 129     return xy_centered_;
718 26 Dec 06 jari 130   }
718 26 Dec 06 jari 131
2565 26 Sep 11 peter 132
2565 26 Sep 11 peter 133   void AveragerPair::xy_add(double mx, double my, double xy_centered, long n)
2565 26 Sep 11 peter 134   {
2565 26 Sep 11 peter 135     if (!n)
2565 26 Sep 11 peter 136       return;
2565 26 Sep 11 peter 137     /*
2644 11 Nov 11 peter 138       For derivation of update formula refer to
2565 26 Sep 11 peter 139       http://www.janinebennett.org/index_files/ParallelStatisticsAlgorithms.pdf
2565 26 Sep 11 peter 140     */
2565 26 Sep 11 peter 141     xy_centered_ += xy_centered;
2565 26 Sep 11 peter 142     if (this->n())
2565 26 Sep 11 peter 143       xy_centered_ += static_cast<double>(n*this->n()) / (n+this->n()) *
2565 26 Sep 11 peter 144         (mx-x_.mean()) * (my-y_.mean());
2565 26 Sep 11 peter 145
2565 26 Sep 11 peter 146   }
2565 26 Sep 11 peter 147
2565 26 Sep 11 peter 148
718 26 Dec 06 jari 149   const Averager& AveragerPair::x_averager(void) const
718 26 Dec 06 jari 150   {
718 26 Dec 06 jari 151     return x_;
718 26 Dec 06 jari 152   }
718 26 Dec 06 jari 153
2565 26 Sep 11 peter 154
718 26 Dec 06 jari 155   const Averager& AveragerPair::y_averager(void) const
718 26 Dec 06 jari 156   {
718 26 Dec 06 jari 157     return y_;
718 26 Dec 06 jari 158   }
718 26 Dec 06 jari 159
2565 26 Sep 11 peter 160
138 20 Aug 04 peter 161   const AveragerPair& AveragerPair::operator+=(const AveragerPair& a)
138 20 Aug 04 peter 162   {
2565 26 Sep 11 peter 163     xy_add(a.x_averager().mean(), a.y_averager().mean(), a.xy_centered_,a.n());
138 20 Aug 04 peter 164     x_+=a.x_averager();
138 20 Aug 04 peter 165     y_+=a.y_averager();
138 20 Aug 04 peter 166     return *this;
138 20 Aug 04 peter 167   }
138 20 Aug 04 peter 168
680 11 Oct 06 jari 169 }}} // of namespace statistics, yat, and theplu

Generated on Sat Oct 14 02:23:36 2023 (UTC) by svndigest 0.11.2pre (dev build 1672)