Automatic Differentiation
 
Loading...
Searching...
No Matches
neg_binomial_2_cdf.hpp
Go to the documentation of this file.
1#ifndef STAN_MATH_PRIM_PROB_NEG_BINOMIAL_2_CDF_HPP
2#define STAN_MATH_PRIM_PROB_NEG_BINOMIAL_2_CDF_HPP
3
18#include <limits>
19
20namespace stan {
21namespace math {
22
23template <typename T_n, typename T_location, typename T_precision>
25 const T_n& n, const T_location& mu, const T_precision& phi) {
27 using T_n_ref = ref_type_t<T_n>;
28 using T_mu_ref = ref_type_t<T_location>;
29 using T_phi_ref = ref_type_t<T_precision>;
30 static constexpr const char* function = "neg_binomial_2_cdf";
31 check_consistent_sizes(function, "Random variable", n, "Location parameter",
32 mu, "Precision Parameter", phi);
33
34 T_n_ref n_ref = n;
35 T_mu_ref mu_ref = mu;
36 T_phi_ref phi_ref = phi;
37
38 check_positive_finite(function, "Location parameter", mu_ref);
39 check_positive_finite(function, "Precision parameter", phi_ref);
40 check_not_nan(function, "Random variable", n_ref);
41
42 if (size_zero(n, mu, phi)) {
43 return 1.0;
44 }
45
46 T_partials_return P(1.0);
47 auto ops_partials = make_partials_propagator(mu_ref, phi_ref);
48
49 scalar_seq_view<T_n_ref> n_vec(n_ref);
50 scalar_seq_view<T_mu_ref> mu_vec(mu_ref);
51 scalar_seq_view<T_phi_ref> phi_vec(phi_ref);
52 size_t size_phi = stan::math::size(phi);
53 size_t size_n_phi = max_size(n, phi);
54 size_t max_size_seq_view = max_size(n, mu, phi);
55
56 // Explicit return for extreme values
57 // The gradients are technically ill-defined, but treated as zero
58 for (size_t i = 0; i < stan::math::size(n); i++) {
59 if (n_vec.val(i) < 0) {
60 return ops_partials.build(0.0);
61 }
62 }
63
65 T_precision>
66 digamma_phi_vec(size_phi);
68 T_precision>
69 digamma_sum_vec(size_n_phi);
70
72 for (size_t i = 0; i < size_phi; i++) {
73 digamma_phi_vec[i] = digamma(phi_vec.val(i));
74 }
75 for (size_t i = 0; i < size_n_phi; i++) {
76 const T_partials_return n_dbl = n_vec.val(i);
77 const T_partials_return phi_dbl = phi_vec.val(i);
78 digamma_sum_vec[i] = digamma(n_dbl + phi_dbl + 1);
79 }
80 }
81
82 for (size_t i = 0; i < max_size_seq_view; i++) {
83 // Explicit results for extreme values
84 // The gradients are technically ill-defined, but treated as zero
85 if (n_vec.val(i) == std::numeric_limits<int>::max()) {
86 return ops_partials.build(1.0);
87 }
88
89 const T_partials_return n_dbl_p1 = n_vec.val(i) + 1;
90 const T_partials_return mu_dbl = mu_vec.val(i);
91 const T_partials_return phi_dbl = phi_vec.val(i);
92 const T_partials_return inv_mu_plus_phi = inv(mu_dbl + phi_dbl);
93 const T_partials_return p_dbl = phi_dbl * inv_mu_plus_phi;
94 const T_partials_return d_dbl = square(inv_mu_plus_phi);
95 const T_partials_return P_i = inc_beta(phi_dbl, n_dbl_p1, p_dbl);
96 const T_partials_return inc_beta_ddz_i
98 ? 0
99 : inc_beta_ddz(phi_dbl, n_dbl_p1, p_dbl) * d_dbl / P_i;
100
101 P *= P_i;
102
104 partials<0>(ops_partials)[i] -= inc_beta_ddz_i * phi_dbl;
105 }
106
108 partials<1>(ops_partials)[i]
109 += inc_beta_dda(phi_dbl, n_dbl_p1, p_dbl, digamma_phi_vec[i],
110 digamma_sum_vec[i])
111 / P_i
112 + inc_beta_ddz_i * mu_dbl;
113 }
114 }
115
117 for (size_t i = 0; i < stan::math::size(mu); ++i) {
118 partials<0>(ops_partials)[i] *= P;
119 }
120 }
121
123 for (size_t i = 0; i < size_phi; ++i) {
124 partials<1>(ops_partials)[i] *= P;
125 }
126 }
127
128 return ops_partials.build(P);
129}
130
131} // namespace math
132} // namespace stan
133#endif
VectorBuilder allocates type T1 values to be used as intermediate values.
scalar_seq_view provides a uniform sequence-like wrapper around either a scalar or a sequence of scal...
typename return_type< Ts... >::type return_type_t
Convenience type for the return type of the specified template parameters.
int64_t size(const T &m)
Returns the size (number of the elements) of a matrix_cl or var_value<matrix_cl<T>>.
Definition size.hpp:19
T inc_beta_ddz(T a, T b, T z)
Returns the partial derivative of the regularized incomplete beta function, I_{z}(a,...
bool size_zero(const T &x)
Returns 1 if input is of length 0, returns 0 otherwise.
Definition size_zero.hpp:19
void check_consistent_sizes(const char *)
Trivial no input case, this function is a no-op.
T inc_beta_dda(T a, T b, T z, T digamma_a, T digamma_ab)
Returns the partial derivative of the regularized incomplete beta function, I_{z}(a,...
fvar< T > inc_beta(const fvar< T > &a, const fvar< T > &b, const fvar< T > &x)
Definition inc_beta.hpp:19
void check_not_nan(const char *function, const char *name, const T_y &y)
Check if y is not NaN.
return_type_t< T_location, T_precision > neg_binomial_2_cdf(const T_n &n, const T_location &mu, const T_precision &phi)
int64_t max_size(const T1 &x1, const Ts &... xs)
Calculate the size of the largest input.
Definition max_size.hpp:20
fvar< T > inv(const fvar< T > &x)
Definition inv.hpp:12
auto make_partials_propagator(Ops &&... ops)
Construct an partials_propagator.
void check_positive_finite(const char *function, const char *name, const T_y &y)
Check if y is positive and finite.
fvar< T > digamma(const fvar< T > &x)
Return the derivative of the log gamma function at the specified argument.
Definition digamma.hpp:23
fvar< T > square(const fvar< T > &x)
Definition square.hpp:12
typename ref_type_if< true, T >::type ref_type_t
Definition ref_type.hpp:55
typename partials_return_type< Args... >::type partials_return_t
The lgamma implementation in stan-math is based on either the reentrant safe lgamma_r implementation ...
Extends std::true_type when instantiated with zero or more template parameters, all of which extend t...