Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Optimise von_mises_lpdf gradient calc #3010

Merged
merged 2 commits into from
Jan 31, 2024
Merged

Optimise von_mises_lpdf gradient calc #3010

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
1dd2dca
Use bessel identity to improve grad calc
andrjohns Jan 20, 2024
7024fae
Update von mises test to use AD framework
andrjohns Jan 20, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion stan/math/prim/prob/von_mises_lpdf.hpp
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -81,7 +81,7 @@ return_type_t<T_y, T_loc, T_scale> von_mises_lpdf(T_y const& y, T_loc const& mu,
if (!is_constant_all<T_scale>::value) {
edge<2>(ops_partials).partials_
= cos_mu_minus_y
- modified_bessel_first_kind(-1, kappa_val)
- modified_bessel_first_kind(1, kappa_val)
/ modified_bessel_first_kind(0, kappa_val);
}

Expand Down
122 changes: 21 additions & 101 deletions test/unit/math/mix/prob/von_mises_test.cpp
View file
Open in desktop
Original file line number Diff line number Diff line change
@@ -1,103 +1,23 @@
#include <stan/math/mix.hpp>
#include <test/unit/math/rev/fun/util.hpp>
#include <gtest/gtest.h>
#include <vector>

std::vector<double> test_von_mises_lpdf(double y, double mu, double kappa) {
using stan::math::var;
using stan::math::von_mises_lpdf;

var y_var = y;
var mu_var = mu;
var kappa_var = kappa;

std::vector<var> x;
x.push_back(y_var);
x.push_back(mu_var);
x.push_back(kappa_var);

var logp = von_mises_lpdf<false>(y_var, mu_var, kappa_var);
std::vector<double> grad;
logp.grad(x, grad);
return grad;
}

TEST(ProbAgradDistributionsVonMises, derivatives) {
using stan::math::fvar;
using stan::math::von_mises_lpdf;

std::vector<double> grad = test_von_mises_lpdf(0, 1, 0);

fvar<double> lp = von_mises_lpdf<false>(0, 1, fvar<double>(0, 1));
EXPECT_FLOAT_EQ(grad[2], lp.d());

fvar<double> kappa1(0, 1);
EXPECT_FLOAT_EQ(von_mises_lpdf(0, 1, kappa1).val_, -1.8378770664093453390819);
EXPECT_FLOAT_EQ(von_mises_lpdf(0, 1, kappa1).d_, 0.54030230586813976501);

fvar<fvar<double>> kappa2(0);
EXPECT_NO_THROW(von_mises_lpdf(0, 1, kappa2));
}

TEST(ProbAgradDistributionsVonMises, FvarVar_1stDeriv) {
using stan::math::fvar;
using stan::math::var;
using stan::math::von_mises_lpdf;

fvar<var> y_(0, 1);
double mu(1);
double kappa(0);

fvar<var> logp = von_mises_lpdf(y_, mu, kappa);

std::vector<stan::math::var> y{y_.val_};
std::vector<double> g;
logp.val_.grad(y, g);
EXPECT_FLOAT_EQ(0, g[0]);
}

TEST(ProbAgradDistributionsVonMises, FvarVar_2ndDeriv1) {
using stan::math::fvar;
using stan::math::var;
using stan::math::von_mises_lpdf;

double y_(0);
fvar<var> mu(1, 1);
double kappa(0);
fvar<var> logp = von_mises_lpdf(y_, mu, kappa);

std::vector<stan::math::var> y{mu.val_};
std::vector<double> g;
logp.d_.grad(y, g);
EXPECT_FLOAT_EQ(0, g[0]);
}

TEST(ProbAgradDistributionsVonMises, FvarVar_2ndDeriv2) {
using stan::math::fvar;
using stan::math::var;
using stan::math::von_mises_lpdf;

double y_(0);
double mu(1);
fvar<var> kappa(0, 1);
fvar<var> logp = von_mises_lpdf(y_, mu, kappa);

std::vector<stan::math::var> y{kappa.val_};
std::vector<double> g;
logp.d_.grad(y, g);
// Fails: g[0] is -nan
// EXPECT_FLOAT_EQ(0, g[0]);
}

// This test once failed sanitizer checks -- nothing explicitly tested in the
// code itself
TEST(ProbAgradDistributionsVonMises, sanitizer_error_fixed) {
using stan::math::var;
double y = boost::math::constants::third_pi<double>();
double mu = boost::math::constants::sixth_pi<double>();
std::vector<var> kappa = {0.5};

auto lp = stan::math::von_mises_lpdf(y, mu, kappa);

lp.grad();
#include <test/unit/math/test_ad.hpp>

TEST(mathMixScalFun, von_mises_lpdf) {
auto f = [](const auto& y, const auto& mu, const auto& kappa) {
return stan::math::von_mises_lpdf(y, mu, kappa);
};

Eigen::VectorXd y(2);
y << 1.0, 2.0;
Eigen::VectorXd mu(2);
mu << 1.0, 0.5;
Eigen::VectorXd kappa(2);
kappa << 1.0, 0.5;

stan::test::expect_ad(f, y[0], mu[0], kappa[0]);
stan::test::expect_ad(f, y[0], mu, kappa);
stan::test::expect_ad(f, y, mu[0], kappa);
stan::test::expect_ad(f, y, mu, kappa[0]);
stan::test::expect_ad(f, y[0], mu[0], kappa);
stan::test::expect_ad(f, y, mu[0], kappa[0]);
stan::test::expect_ad(f, y[0], mu, kappa[0]);
}