yat/regression/MultiDimensional.cc

Code
Comments
Other
Rev Date Author Line
386 14 Aug 05 jari 1 // $Id$
386 14 Aug 05 jari 2
675 10 Oct 06 jari 3 /*
2119 12 Dec 09 peter 4   Copyright (C) 2005 Jari Häkkinen
2119 12 Dec 09 peter 5   Copyright (C) 2006, 2007, 2008 Jari Häkkinen, Peter Johansson
4359 23 Aug 23 peter 6   Copyright (C) 2011, 2012, 2022, 2023 Peter Johansson
386 14 Aug 05 jari 7
1437 25 Aug 08 peter 8   This file is part of the yat library, http://dev.thep.lu.se/yat
675 10 Oct 06 jari 9
675 10 Oct 06 jari 10   The yat library is free software; you can redistribute it and/or
675 10 Oct 06 jari 11   modify it under the terms of the GNU General Public License as
1486 09 Sep 08 jari 12   published by the Free Software Foundation; either version 3 of the
675 10 Oct 06 jari 13   License, or (at your option) any later version.
675 10 Oct 06 jari 14
675 10 Oct 06 jari 15   The yat library is distributed in the hope that it will be useful,
675 10 Oct 06 jari 16   but WITHOUT ANY WARRANTY; without even the implied warranty of
675 10 Oct 06 jari 17   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
675 10 Oct 06 jari 18   General Public License for more details.
675 10 Oct 06 jari 19
675 10 Oct 06 jari 20   You should have received a copy of the GNU General Public License
1487 10 Sep 08 jari 21   along with yat. If not, see <http://www.gnu.org/licenses/>.
675 10 Oct 06 jari 22 */
675 10 Oct 06 jari 23
2881 18 Nov 12 peter 24 #include <config.h>
2881 18 Nov 12 peter 25
680 11 Oct 06 jari 26 #include "MultiDimensional.h"
750 17 Feb 07 jari 27 #include "yat/utility/Exception.h"
1121 22 Feb 08 peter 28 #include "yat/utility/Matrix.h"
1021 01 Feb 08 peter 29 #include "yat/utility/VectorBase.h"
1120 21 Feb 08 peter 30 #include "yat/utility/Vector.h"
675 10 Oct 06 jari 31
731 06 Jan 07 peter 32 #include <cassert>
731 06 Jan 07 peter 33
386 14 Aug 05 jari 34 namespace theplu {
680 11 Oct 06 jari 35 namespace yat {
386 14 Aug 05 jari 36 namespace regression {
386 14 Aug 05 jari 37
386 14 Aug 05 jari 38
703 18 Dec 06 jari 39   MultiDimensional::MultiDimensional(void)
703 18 Dec 06 jari 40     : chisquare_(0), work_(NULL)
703 18 Dec 06 jari 41   {
703 18 Dec 06 jari 42   }
703 18 Dec 06 jari 43
703 18 Dec 06 jari 44
703 18 Dec 06 jari 45   MultiDimensional::~MultiDimensional(void)
703 18 Dec 06 jari 46   {
3907 09 May 20 peter 47     gsl_multifit_linear_free(work_);
703 18 Dec 06 jari 48   }
703 18 Dec 06 jari 49
703 18 Dec 06 jari 50
1121 22 Feb 08 peter 51   const utility::Matrix& MultiDimensional::covariance(void) const
726 04 Jan 07 peter 52   {
726 04 Jan 07 peter 53     return covariance_;
726 04 Jan 07 peter 54   }
726 04 Jan 07 peter 55
726 04 Jan 07 peter 56
4125 14 Jan 22 peter 57   void MultiDimensional::fit(const utility::Matrix& x,
1021 01 Feb 08 peter 58                              const utility::VectorBase& y)
386 14 Aug 05 jari 59   {
4130 20 Jan 22 peter 60     fit2(x, y);
4130 20 Jan 22 peter 61   }
4130 20 Jan 22 peter 62
4130 20 Jan 22 peter 63
4130 20 Jan 22 peter 64   void MultiDimensional::fit2(const utility::MatrixBase& x,
4130 20 Jan 22 peter 65                               const utility::VectorBase& y)
4130 20 Jan 22 peter 66   {
731 06 Jan 07 peter 67     assert(x.rows()==y.size());
1098 18 Feb 08 jari 68     covariance_.resize(x.columns(),x.columns());
2534 31 Jul 11 peter 69     fit_parameters_.resize(x.columns());
3907 09 May 20 peter 70     gsl_multifit_linear_free(work_);
750 17 Feb 07 jari 71     if (!(work_=gsl_multifit_linear_alloc(x.rows(),fit_parameters_.size())))
750 17 Feb 07 jari 72       throw utility::GSL_error("MultiDimensional::fit failed to allocate memory");
750 17 Feb 07 jari 73
750 17 Feb 07 jari 74     int status = gsl_multifit_linear(x.gsl_matrix_p(), y.gsl_vector_p(),
750 17 Feb 07 jari 75                                      fit_parameters_.gsl_vector_p(),
750 17 Feb 07 jari 76                                      covariance_.gsl_matrix_p(), &chisquare_,
750 17 Feb 07 jari 77                                      work_);
750 17 Feb 07 jari 78     if (status)
750 17 Feb 07 jari 79       throw utility::GSL_error(std::string("MultiDimensional::fit",status));
750 17 Feb 07 jari 80
739 12 Jan 07 peter 81     s2_ = chisquare_/(x.rows()-x.columns());
386 14 Aug 05 jari 82   }
386 14 Aug 05 jari 83
2534 31 Jul 11 peter 84
1120 21 Feb 08 peter 85   const utility::Vector& MultiDimensional::fit_parameters(void) const
4200 19 Aug 22 peter 86   {
4200 19 Aug 22 peter 87     return fit_parameters_;
713 21 Dec 06 peter 88   }
386 14 Aug 05 jari 89
713 21 Dec 06 peter 90
713 21 Dec 06 peter 91   double MultiDimensional::chisq(void) const
713 21 Dec 06 peter 92   {
713 21 Dec 06 peter 93     return chisquare_;
713 21 Dec 06 peter 94   }
713 21 Dec 06 peter 95
713 21 Dec 06 peter 96
4200 19 Aug 22 peter 97   double MultiDimensional::predict(const utility::VectorBase& x) const
718 26 Dec 06 jari 98   {
731 06 Jan 07 peter 99     assert(x.size()==fit_parameters_.size());
718 26 Dec 06 jari 100     return fit_parameters_ * x;
718 26 Dec 06 jari 101   }
718 26 Dec 06 jari 102
718 26 Dec 06 jari 103
1021 01 Feb 08 peter 104   double MultiDimensional::prediction_error2(const utility::VectorBase& x) const
586 19 Jun 06 peter 105   {
739 12 Jan 07 peter 106     return standard_error2(x) + s2_;
586 19 Jun 06 peter 107   }
586 19 Jun 06 peter 108
586 19 Jun 06 peter 109
1021 01 Feb 08 peter 110   double MultiDimensional::standard_error2(const utility::VectorBase& x) const
586 19 Jun 06 peter 111   {
4291 03 Feb 23 peter 112     return x * covariance_ * x;
586 19 Jun 06 peter 113   }
586 19 Jun 06 peter 114
681 11 Oct 06 jari 115 }}} // of namespaces regression, yat, and theplu

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