Stan Math Library  2.20.0
reverse mode automatic differentiation
ibeta.hpp
Go to the documentation of this file.
1 #ifndef STAN_MATH_REV_SCAL_FUN_IBETA_HPP
2 #define STAN_MATH_REV_SCAL_FUN_IBETA_HPP
3 
4 #include <stan/math/rev/meta.hpp>
5 #include <stan/math/rev/core.hpp>
9 
10 namespace stan {
11 namespace math {
12 
13 namespace internal {
19 inline double ibeta_hypergeometric_helper(double a, double b, double z,
20  double precision = 1e-8,
21  double max_steps = 1000) {
22  double val = 0;
23  double diff = 1;
24  double k = 0;
25  double a_2 = a * a;
26  double bprod = 1;
27  while (std::abs(diff) > precision && ++k < max_steps) {
28  val += diff;
29  bprod *= b + k - 1.0;
30  diff = a_2 * std::pow(a + k, -2) * bprod * std::pow(z, k) / tgamma(k + 1);
31  }
32  return val;
33 }
34 
35 class ibeta_vvv_vari : public op_vvv_vari {
36  public:
37  ibeta_vvv_vari(vari* avi, vari* bvi, vari* xvi)
38  : op_vvv_vari(ibeta(avi->val_, bvi->val_, xvi->val_), avi, bvi, xvi) {}
39  void chain() {
40  double a = avi_->val_;
41  double b = bvi_->val_;
42  double c = cvi_->val_;
43 
45  using std::log;
46  using std::pow;
47  using std::sin;
48  avi_->adj_ += adj_ * (log(c) - digamma(a) + digamma(a + b)) * val_
49  - tgamma(a) * tgamma(a + b) / tgamma(b) * pow(c, a)
50  / tgamma(1 + a) / tgamma(1 + a)
51  * ibeta_hypergeometric_helper(a, 1 - b, c);
52  bvi_->adj_ += adj_
53  * (tgamma(b) * tgamma(a + b) / tgamma(a) * pow(1 - c, b)
54  * ibeta_hypergeometric_helper(b, 1 - a, 1 - c)
55  / tgamma(b + 1) / tgamma(b + 1)
56  + ibeta(b, a, 1 - c)
57  * (digamma(b) - digamma(a + b) - log(1 - c)));
58  cvi_->adj_ += adj_ * boost::math::ibeta_derivative(a, b, c);
59  }
60 };
61 class ibeta_vvd_vari : public op_vvd_vari {
62  public:
63  ibeta_vvd_vari(vari* avi, vari* bvi, double x)
64  : op_vvd_vari(ibeta(avi->val_, bvi->val_, x), avi, bvi, x) {}
65  void chain() {
66  double a = avi_->val_;
67  double b = bvi_->val_;
68  double c = cd_;
69 
71  using std::log;
72  using std::pow;
73  using std::sin;
74  avi_->adj_ += adj_ * (log(c) - digamma(a) + digamma(a + b)) * val_
75  - tgamma(a) * tgamma(a + b) / tgamma(b) * pow(c, a)
76  / tgamma(1 + a) / tgamma(1 + a)
77  * ibeta_hypergeometric_helper(a, 1 - b, c);
78  bvi_->adj_ += adj_
79  * (tgamma(b) * tgamma(a + b) / tgamma(a) * pow(1 - c, b)
80  * ibeta_hypergeometric_helper(b, 1 - a, 1 - c)
81  / tgamma(b + 1) / tgamma(b + 1)
82  + ibeta(b, a, 1 - c)
83  * (digamma(b) - digamma(a + b) - log(1 - c)));
84  }
85 };
86 class ibeta_vdv_vari : public op_vdv_vari {
87  public:
88  ibeta_vdv_vari(vari* avi, double b, vari* xvi)
89  : op_vdv_vari(ibeta(avi->val_, b, xvi->val_), avi, b, xvi) {}
90  void chain() {
91  double a = avi_->val_;
92  double b = bd_;
93  double c = cvi_->val_;
94 
97  using boost::math::ibeta;
98  using std::log;
99  using std::pow;
100  using std::sin;
101  avi_->adj_ += adj_ * (log(c) - digamma(a) + digamma(a + b)) * val_
102  - tgamma(a) * tgamma(a + b) / tgamma(b) * pow(c, a)
103  / tgamma(1 + a) / tgamma(1 + a)
104  * ibeta_hypergeometric_helper(a, 1 - b, c);
105  cvi_->adj_ += adj_ * boost::math::ibeta_derivative(a, b, c);
106  }
107 };
108 class ibeta_vdd_vari : public op_vdd_vari {
109  public:
110  ibeta_vdd_vari(vari* avi, double b, double x)
111  : op_vdd_vari(ibeta(avi->val_, b, x), avi, b, x) {}
112  void chain() {
113  double a = avi_->val_;
114  double b = bd_;
115  double c = cd_;
116 
118  using boost::math::digamma;
119  using boost::math::ibeta;
120  using std::log;
121  using std::pow;
122  using std::sin;
123  avi_->adj_ += adj_ * (log(c) - digamma(a) + digamma(a + b)) * val_
124  - tgamma(a) * tgamma(a + b) / tgamma(b) * pow(c, a)
125  / tgamma(1 + a) / tgamma(1 + a)
126  * ibeta_hypergeometric_helper(a, 1 - b, c);
127  }
128 };
129 class ibeta_dvv_vari : public op_dvv_vari {
130  public:
131  ibeta_dvv_vari(double a, vari* bvi, vari* xvi)
132  : op_dvv_vari(ibeta(a, bvi->val_, xvi->val_), a, bvi, xvi) {}
133  void chain() {
134  double a = ad_;
135  double b = bvi_->val_;
136  double c = cvi_->val_;
137 
139  using boost::math::digamma;
140  using boost::math::ibeta;
141  using std::log;
142  using std::pow;
143  using std::sin;
144  bvi_->adj_ += adj_
145  * (tgamma(b) * tgamma(a + b) / tgamma(a) * pow(1 - c, b)
146  * ibeta_hypergeometric_helper(b, 1 - a, 1 - c)
147  / tgamma(b + 1) / tgamma(b + 1)
148  + ibeta(b, a, 1 - c)
149  * (digamma(b) - digamma(a + b) - log(1 - c)));
150  cvi_->adj_ += adj_ * boost::math::ibeta_derivative(a, b, c);
151  }
152 };
153 class ibeta_dvd_vari : public op_dvd_vari {
154  public:
155  ibeta_dvd_vari(double a, vari* bvi, double x)
156  : op_dvd_vari(ibeta(a, bvi->val_, x), a, bvi, x) {}
157  void chain() {
158  double a = ad_;
159  double b = bvi_->val_;
160  double c = cd_;
161 
163  using boost::math::digamma;
164  using boost::math::ibeta;
165  using std::log;
166  using std::pow;
167  using std::sin;
168  bvi_->adj_ += adj_
169  * (tgamma(b) * tgamma(a + b) / tgamma(a) * pow(1 - c, b)
170  * ibeta_hypergeometric_helper(b, 1 - a, 1 - c)
171  / tgamma(b + 1) / tgamma(b + 1)
172  + ibeta(b, a, 1 - c)
173  * (digamma(b) - digamma(a + b) - log(1 - c)));
174  }
175 };
176 class ibeta_ddv_vari : public op_ddv_vari {
177  public:
178  ibeta_ddv_vari(double a, double b, vari* xvi)
179  : op_ddv_vari(ibeta(a, b, xvi->val_), a, b, xvi) {}
180  void chain() {
181  double a = ad_;
182  double b = bd_;
183  double c = cvi_->val_;
184 
185  cvi_->adj_ += adj_ * boost::math::ibeta_derivative(a, b, c);
186  }
187 };
188 } // namespace internal
189 
208 inline var ibeta(const var& a, const var& b, const var& x) {
209  return var(new internal::ibeta_vvv_vari(a.vi_, b.vi_, x.vi_));
210 }
211 
212 } // namespace math
213 } // namespace stan
214 #endif
ibeta_ddv_vari(double a, double b, vari *xvi)
Definition: ibeta.hpp:178
fvar< T > abs(const fvar< T > &x)
Definition: abs.hpp:14
ibeta_dvv_vari(double a, vari *bvi, vari *xvi)
Definition: ibeta.hpp:131
ibeta_vvd_vari(vari *avi, vari *bvi, double x)
Definition: ibeta.hpp:63
ibeta_vdv_vari(vari *avi, double b, vari *xvi)
Definition: ibeta.hpp:88
void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
Definition: ibeta.hpp:180
fvar< T > log(const fvar< T > &x)
Definition: log.hpp:12
The variable implementation base class.
Definition: vari.hpp:30
Independent (input) and dependent (output) variables for gradients.
Definition: var.hpp:33
void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
Definition: ibeta.hpp:39
friend class var
Definition: vari.hpp:32
void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
Definition: ibeta.hpp:133
const double val_
The value of this variable.
Definition: vari.hpp:38
void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
Definition: ibeta.hpp:112
void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
Definition: ibeta.hpp:157
fvar< T > sin(const fvar< T > &x)
Definition: sin.hpp:11
ibeta_vvv_vari(vari *avi, vari *bvi, vari *xvi)
Definition: ibeta.hpp:37
vari * vi_
Pointer to the implementation of this variable.
Definition: var.hpp:45
ibeta_vdd_vari(vari *avi, double b, double x)
Definition: ibeta.hpp:110
double ibeta_hypergeometric_helper(double a, double b, double z, double precision=1e-8, double max_steps=1000)
Calculates the generalized hypergeometric 3F2(a, a, b; a + 1, a + 1; z).
Definition: ibeta.hpp:19
double ibeta(double a, double b, double x)
The normalized incomplete beta function of a, b, and x.
Definition: ibeta.hpp:24
double e()
Return the base of the natural logarithm.
Definition: constants.hpp:87
ibeta_dvd_vari(double a, vari *bvi, double x)
Definition: ibeta.hpp:155
void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
Definition: ibeta.hpp:90
var ibeta(const var &a, const var &b, const var &x)
The normalized incomplete beta function of a, b, and x.
Definition: ibeta.hpp:208
void chain()
Apply the chain rule to this variable based on the variables on which it depends. ...
Definition: ibeta.hpp:65
double pi()
Return the value of pi.
Definition: constants.hpp:80
double adj_
The adjoint of this variable, which is the partial derivative of this variable with respect to the ro...
Definition: vari.hpp:44
fvar< T > pow(const fvar< T > &x1, const fvar< T > &x2)
Definition: pow.hpp:16
fvar< T > tgamma(const fvar< T > &x)
Return the result of applying the gamma function to the specified argument.
Definition: tgamma.hpp:21
fvar< T > digamma(const fvar< T > &x)
Return the derivative of the log gamma function at the specified argument.
Definition: digamma.hpp:23

     [ Stan Home Page ] © 2011–2018, Stan Development Team.