yat/regression/NegativeBinomial.cc

Code
Comments
Other
Rev Date Author Line
3659 13 Jul 17 peter 1 // $Id$
3659 13 Jul 17 peter 2
3659 13 Jul 17 peter 3 /*
4207 26 Aug 22 peter 4   Copyright (C) 2017, 2022 Peter Johansson
3659 13 Jul 17 peter 5
3659 13 Jul 17 peter 6   This file is part of the yat library, http://dev.thep.lu.se/yat
3659 13 Jul 17 peter 7
3659 13 Jul 17 peter 8   The yat library is free software; you can redistribute it and/or
3659 13 Jul 17 peter 9   modify it under the terms of the GNU General Public License as
3659 13 Jul 17 peter 10   published by the Free Software Foundation; either version 3 of the
3659 13 Jul 17 peter 11   License, or (at your option) any later version.
3659 13 Jul 17 peter 12
3659 13 Jul 17 peter 13   The yat library is distributed in the hope that it will be useful,
3659 13 Jul 17 peter 14   but WITHOUT ANY WARRANTY; without even the implied warranty of
3659 13 Jul 17 peter 15   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
3659 13 Jul 17 peter 16   General Public License for more details.
3659 13 Jul 17 peter 17
3659 13 Jul 17 peter 18   You should have received a copy of the GNU General Public License
3659 13 Jul 17 peter 19   along with yat. If not, see <http://www.gnu.org/licenses/>.
3659 13 Jul 17 peter 20 */
3659 13 Jul 17 peter 21
3659 13 Jul 17 peter 22 #include <config.h>
3659 13 Jul 17 peter 23
3659 13 Jul 17 peter 24 #include "NegativeBinomial.h"
3659 13 Jul 17 peter 25
3663 20 Jul 17 peter 26 #include "detail/NegativeBinomialFitter.h"
3663 20 Jul 17 peter 27
3659 13 Jul 17 peter 28 #include <yat/utility/Matrix.h>
3659 13 Jul 17 peter 29 #include <yat/utility/Vector.h>
3659 13 Jul 17 peter 30 #include <yat/utility/VectorBase.h>
3659 13 Jul 17 peter 31 #include <yat/utility/VectorConstView.h>
3659 13 Jul 17 peter 32
3659 13 Jul 17 peter 33 #include <cassert>
3659 13 Jul 17 peter 34 #include <cmath>
3659 13 Jul 17 peter 35
3659 13 Jul 17 peter 36 namespace theplu {
3659 13 Jul 17 peter 37 namespace yat {
3659 13 Jul 17 peter 38 namespace regression {
3659 13 Jul 17 peter 39
3663 20 Jul 17 peter 40   double NegativeBinomial::alpha(void) const
3663 20 Jul 17 peter 41   {
3663 20 Jul 17 peter 42     return alpha_;
3663 20 Jul 17 peter 43   }
3663 20 Jul 17 peter 44
3663 20 Jul 17 peter 45
3663 20 Jul 17 peter 46   const utility::Matrix& NegativeBinomial::covariance(void) const
3663 20 Jul 17 peter 47   {
3663 20 Jul 17 peter 48     return covariance_;
3663 20 Jul 17 peter 49   }
3663 20 Jul 17 peter 50
3663 20 Jul 17 peter 51
4130 20 Jan 22 peter 52   void NegativeBinomial::fit(const utility::Matrix& x,
3659 13 Jul 17 peter 53                              const utility::VectorBase& y)
3659 13 Jul 17 peter 54   {
4130 20 Jan 22 peter 55     fit2(x, y);
4130 20 Jan 22 peter 56   }
4130 20 Jan 22 peter 57
4130 20 Jan 22 peter 58
4130 20 Jan 22 peter 59   void NegativeBinomial::fit2(const utility::MatrixBase& X,
4130 20 Jan 22 peter 60                              const utility::VectorBase& y)
4130 20 Jan 22 peter 61   {
3659 13 Jul 17 peter 62     assert(X.rows() == y.size());
3663 20 Jul 17 peter 63
3663 20 Jul 17 peter 64     utility::Matrix X2(X.rows(), X.columns()+1, 1.0);
3663 20 Jul 17 peter 65     for (size_t i=0; i<X.columns(); ++i)
3663 20 Jul 17 peter 66       X2.column_view(i+1) = X.column_const_view(i);
3663 20 Jul 17 peter 67
3663 20 Jul 17 peter 68     detail::NegativeBinomialFitter fitter(X2, y, beta_, covariance_, alpha_);
3663 20 Jul 17 peter 69     assert(beta_.size() == X.columns() + 1);
3659 13 Jul 17 peter 70   }
3659 13 Jul 17 peter 71
3659 13 Jul 17 peter 72
3659 13 Jul 17 peter 73   const utility::Vector& NegativeBinomial::fit_parameters(void) const
3659 13 Jul 17 peter 74   {
3659 13 Jul 17 peter 75     return beta_;
3659 13 Jul 17 peter 76   }
3659 13 Jul 17 peter 77
3659 13 Jul 17 peter 78
3659 13 Jul 17 peter 79   double NegativeBinomial::predict(const utility::VectorBase& x) const
3659 13 Jul 17 peter 80   {
3659 13 Jul 17 peter 81     assert(x.size() + 1 == beta_.size());
3659 13 Jul 17 peter 82     using utility::VectorConstView;
3659 13 Jul 17 peter 83     double logm = beta_(0) + VectorConstView(beta_, 1, x.size()) * x;
3659 13 Jul 17 peter 84     return std::exp(logm);
3659 13 Jul 17 peter 85   }
3659 13 Jul 17 peter 86
3659 13 Jul 17 peter 87 }}} // of namespaces regression, yat, and theplu