Conversation
Jenkins Console Log Machine informationNo LSB modules are available. Distributor ID: Ubuntu Description: Ubuntu 20.04.3 LTS Release: 20.04 Codename: focalCPU: G++: Clang: |
|
Rtools & ARM64 CI / rev tests (windows-11-arm) (pull_request) Seems a test on ode fails, but this might be unrelated to this PR? My mac can pass the rev tests. |
|
@lingium it passed on re-run, so I think it was just random noise |
SteveBronder
left a comment
SteveBronder
left a comment
There was a problem hiding this comment.
Two minor comments but imo looks good!
| scalar_seq_view<T_w> w_vec(w); | ||
| size_t size_w = stan::math::size(w); | ||
|
|
||
| VectorBuilder<true, int, T_alpha> output(size_w); | ||
| for (size_t n = 0; n < size_w; ++n) { | ||
| double p = stan::math::exp(-w_vec.val(n)); | ||
| double odds_ratio_p | ||
| = stan::math::exp(stan::math::log(p) - stan::math::log1m(p)); | ||
| output[n] = neg_binomial_rng(1.0, odds_ratio_p, rng) + 1; | ||
| } | ||
|
|
||
| return output.data(); |
There was a problem hiding this comment.
Since you use neg_binomial_rng I think we can use the vectorized version of everything here
| scalar_seq_view<T_w> w_vec(w); | |
| size_t size_w = stan::math::size(w); | |
| VectorBuilder<true, int, T_alpha> output(size_w); | |
| for (size_t n = 0; n < size_w; ++n) { | |
| double p = stan::math::exp(-w_vec.val(n)); | |
| double odds_ratio_p | |
| = stan::math::exp(stan::math::log(p) - stan::math::log1m(p)); | |
| output[n] = neg_binomial_rng(1.0, odds_ratio_p, rng) + 1; | |
| } | |
| return output.data(); | |
| auto p = stan::math::exp(-w); | |
| auto odds_ratio_p = stan::math::exp(stan::math::log(p) - stan::math::log1m(p)); | |
| return neg_binomial_rng(1.0, odds_ratio_p, rng) + 1; |
There was a problem hiding this comment.
Seems neg_binomial_rng(...) + 1; cannot pass compile, let's keep the loop version?
There was a problem hiding this comment.
How about like this?
/** \ingroup prob_dists
* Return a yule-simon random variate with the given shape parameter,
* using the given random number generator.
*
* alpha can be a scalar or a one-dimensional container.
*
* @tparam T_alpha A scalar, `std::vector` or Eigen vector
* @tparam RNG type of random number generator
*
* @param alpha (Sequence of) shape parameter(s)
* @param rng random number generator
* @return (Sequence of) yule-simon random variate(s)
* @throw std::domain_error if alpha is nonpositive
*/
template <typename T_alpha, typename RNG>
inline auto yule_simon_rng(T_alpha&& alpha, RNG& rng) {
static constexpr const char *function = "yule_simon_rng";
decltype(auto) alpha_ref = to_ref(std::forward<T_alpha>(alpha));
check_positive_finite(function, "Shape parameter", alpha_ref);
auto w = exponential_rng(std::forward<decltype(alpha_ref)>(alpha_ref), rng);
auto w_arr = as_array_or_scalar(w);
const auto p = stan::math::exp(-w_arr);
const auto odds_ratio_p
= stan::math::exp(stan::math::log(p) - stan::math::log1m(p));
if constexpr (is_stan_scalar_v<T_alpha>) {
return neg_binomial_rng(1.0, odds_ratio_p, rng) + 1;
} else {
return to_array_1d(as_array_or_scalar(neg_binomial_rng(1.0, std::move(odds_ratio_p), rng)) + 1);
}
}
Jenkins Console Log Machine informationNo LSB modules are available. Distributor ID: Ubuntu Description: Ubuntu 20.04.3 LTS Release: 20.04 Codename: focalCPU: G++: Clang: |
4 participants
Summary
With this PR the rng of Yule-Simon distribution are added.
See issue #3239
For$$Y\sim \text{Yule-Simon}(\alpha)$$ :
It has a Beta–Geometric mixture. To sample form$$Y$$ , we do:
Stan does not have Geometric yet so we use
neg_binomial_rng(1.0, odds_ratio_p).Tests
Test is written follow the guide.
Side Effects
No.
Release notes
yule_simon_rngis available if merged.Checklist
Copyright holder: Zhi Ling
The copyright holder is typically you or your assignee, such as a university or company. By submitting this pull request, the copyright holder is agreeing to the license the submitted work under the following licenses:
- Code: BSD 3-clause (https://opensource.org/licenses/BSD-3-Clause)
- Documentation: CC-BY 4.0 (https://creativecommons.org/licenses/by/4.0/)
the basic tests are passing
./runTests.py test/unit)make test-headers)make test-math-dependencies)make doxygen)make cpplint)the code is written in idiomatic C++ and changes are documented in the doxygen
the new changes are tested