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.

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

Fixing propto tests in probability test framework #1764

Merged
merged 22 commits into from
Jul 22, 2020

Conversation

bbbales2
Copy link
Member

@bbbales2 bbbales2 commented Mar 5, 2020

This solves issue #1763

Tests

This is a bugfix to a function in the test framework. I did not add any tests for this. Do we have tests for the probability test framework?

Release Notes

Fixed propto tests in probability test framework (and fixed bugs this revealed in neg_binomial_lpdf and pareto_type_2_lpdf).

Checklist

@mcol
Copy link
Contributor

mcol commented Mar 6, 2020

This makes the prop tests for neg_binomial_lpdf, but I see that normal_lpdf fails as well, and perhaps other distributions too. Looking at the normal is much easier, as the distribution is much simpler and there are no numerical cutoffs. Right now I don't see anything wrong there, which may indicate that there's still something wrong in the distribution test, but before reaching that conclusion could you have a look at that too?

@mcol
Copy link
Contributor

mcol commented Mar 6, 2020

This is not an exhaustive list, but also cauchy_lpdf, chi_square_lpdf, exponential and gamma_lpdf fail, while bernoulli_lpmf, binomial_lpmf, poisson_lpmf and std_normal_lpdf pass.

@bbbales2
Copy link
Member Author

bbbales2 commented Mar 6, 2020

Good catch. I think the check here: https://github.com/stan-dev/math/blob/develop/test/prob/test_fixture_distr.hpp#L254

Needs changed from:

EXPECT_TRUE(reference_logprob_false - logprob_false
            == reference_logprob_true - logprob_true)

to:

EXPECT_NEAR(value_of(reference_logprob_false - reference_logprob_true),
            value_of(logprob_false - logprob_true), 1e-12)

This just checks to make sure the two numbers are within 1e-12 of each other (absolute tolerance).

With this in place only negative binomial fails.

@bob-carpenter
Copy link
Contributor

bob-carpenter commented Mar 6, 2020 via email

@mcol
Copy link
Contributor

mcol commented Mar 6, 2020

With this in place only negative binomial fails.

Great! I went back as far as 2.18.1 and the neg_binomial bug was already there, so it's nothing new. I couldn't go further back as I had some compilation errors I didn't think it was worth trying to fix (probably the makefiles for 2.17 were not as clever as what we are used to now).

@mcol mcol linked an issue Mar 6, 2020 that may be closed by this pull request
@bbbales2
Copy link
Member Author

bbbales2 commented Mar 6, 2020

I suspect the negative binomial is something @martinmodrak might know about. He got super deep on these functions a while back. I asked a question over here: #1497 (comment).

@martinmodrak
Copy link
Contributor

Yes, I believe the problem might be fixed with the new code, which, unfortunately, got a bit stuck, will try to move it forward a bit.

@bbbales2
Copy link
Member Author

bbbales2 commented May 2, 2020

@martinmodrak I removed the poisson stuff from neg_binomial in this branch and then updated the numerics a bit. I'm not sure I did it right though. Mind taking a look when you get the chance (no hurry)?

For n = 13 alpha = 1e11 and beta = 1e10, the version of negative_binomial in this branch is giving me:

-2.6185576442208003

R gives me:

-2.618557685656014211

with this code:

n = 13
alpha = 1e11
beta = 1e10
lp = dnbinom(n, mu = alpha / beta, size = alpha, log = TRUE)
sprintf("%.18f", lp)

Develop math was giving:

-2.618380836440565

So we improved but I'm not sure it's fixed yet.

@bbbales2
Copy link
Member Author

@martinmodrak with these values:

n = 13
alpha = 1.0e11
beta = 1.0e10

Can you run the code with lots of precision in Mathematica and tell me the result:

LogGamma[n + alpha] - LogGamma[alpha] - LogGamma[n + 1] +  alpha * Log[beta / (beta + 1)] - n * Log[beta + 1]

I'm trying to figure out how off the version of neg_binomial I implemented here is.

@martinmodrak
Copy link
Contributor

@bbbales2 Unfortunately, Mathematica (at least the free version available in Wolfram Cloud) refuses to compute this (allowed computation time exceeded).... To test for those big values of alpha and beta, the best I could do was to find analytical solutions for n=0 and n=1 (which then can be computed).

@andrjohns
Copy link
Collaborator

The raspberry pi has Mathematica freely available, I tested with those values (may have done it wrong, not at all familiar with mathematica) to 20 decimals:

In[1]:= n = 13                                                                          

Out[1]= 13

In[2]:= alpha = 100000000000                                                            

Out[2]= 100000000000

In[3]:= beta  = 10000000000                                                             

Out[3]= 10000000000
                                                                       
In[4]:= N[LogGamma[n + alpha] - LogGamma[alpha] - LogGamma[n + 1] +  alpha * Log[beta /  (beta + 1)] - n * Log[beta + 1],{\[Infinity],20}]                                       

Out[4]= -2.6185576442208289933

@bbbales2
Copy link
Member Author

@martinmodrak @andrjohns thanks!

What we have in this branch is: -2.6185576442208003
What Mathematica gives: -2.6185576442208289933

I'll look into the n = 0 / n = 1 thing. I just wanted a spot check that things weren't totally awry. The big fix was dropping the Poisson like @martinmodrak did previously so I think we're probably better off already.

@bbbales2
Copy link
Member Author

I added an n = 0 and n = 1 comparison against long double calculations to the tests but I only compared with EXPECT_FLOAT_EQ so it's not really testing for much precision.

The numbers are (printed out to 17 digits):

n = 0:
implementation: -9.9999999995
long double:    -9.9999999995
n = 1
implementation: -7.6974149066059496
long double:    -7.6974149066159543

@rok-cesnovar
Copy link
Member

Are we targeting this for the release? Given that its a test it probably doesnt matter?

@bbbales2
Copy link
Member Author

@rok-cesnovar it's bugfixes, so I wasn't worried about the feature freeze. If it passes it's ready for review though.

@rok-cesnovar
Copy link
Member

Oh right. Yeah its a bug fix definitely.

@stan-buildbot
Copy link
Contributor


Name Old Result New Result Ratio Performance change( 1 - new / old )
gp_pois_regr/gp_pois_regr.stan 4.15 4.13 1.0 0.48% faster
low_dim_corr_gauss/low_dim_corr_gauss.stan 0.02 0.02 1.01 0.54% faster
eight_schools/eight_schools.stan 0.09 0.09 1.03 2.72% faster
gp_regr/gp_regr.stan 0.19 0.19 1.01 0.77% faster
irt_2pl/irt_2pl.stan 5.32 5.3 1.0 0.23% faster
performance.compilation 87.01 85.88 1.01 1.3% faster
low_dim_gauss_mix_collapse/low_dim_gauss_mix_collapse.stan 8.52 8.53 1.0 -0.03% slower
pkpd/one_comp_mm_elim_abs.stan 26.64 27.88 0.96 -4.67% slower
sir/sir.stan 112.39 116.09 0.97 -3.3% slower
gp_regr/gen_gp_data.stan 0.05 0.05 1.01 0.8% faster
low_dim_gauss_mix/low_dim_gauss_mix.stan 3.29 3.34 0.98 -1.73% slower
pkpd/sim_one_comp_mm_elim_abs.stan 0.38 0.38 1.0 -0.17% slower
arK/arK.stan 1.82 1.84 0.99 -0.87% slower
arma/arma.stan 0.68 0.64 1.07 6.7% faster
garch/garch.stan 0.52 0.53 1.0 -0.42% slower
Mean result: 1.00219296502

Jenkins Console Log
Blue Ocean
Commit hash: 283a2b0


Machine information ProductName: Mac OS X ProductVersion: 10.11.6 BuildVersion: 15G22010

CPU:
Intel(R) Xeon(R) CPU E5-1680 v2 @ 3.00GHz

G++:
Configured with: --prefix=/Applications/Xcode.app/Contents/Developer/usr --with-gxx-include-dir=/usr/include/c++/4.2.1
Apple LLVM version 7.0.2 (clang-700.1.81)
Target: x86_64-apple-darwin15.6.0
Thread model: posix

Clang:
Apple LLVM version 7.0.2 (clang-700.1.81)
Target: x86_64-apple-darwin15.6.0
Thread model: posix

@bbbales2
Copy link
Member Author

@andrjohns @SteveBronder this is ready to be reviewed.

There's three things going on here:

  1. A fix to the original issue: propto tests not testing the correct things in probability unit tests #1763

  2. Changed how neg_binomial worked to fix a problem this caught (avoids the poisson approximation and uses some new numerics copying from: Fixing negative binomial phi cutoff #1497)

  3. A fix a problem with templating in pareto_type_2 this caught

@andrjohns
Copy link
Collaborator

I can review this one tonight

Copy link
Collaborator

@andrjohns andrjohns left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Mostly comments on the ordering of some of the calcs, and some basic q's around the tests. Otherwise looks great!

I didn't check the any of the numerics, since I assumed the approach was reviewed in Martin's original PR (plus the tests are passing)

stan/math/prim/prob/neg_binomial_lpmf.hpp Show resolved Hide resolved
stan/math/prim/prob/neg_binomial_lpmf.hpp Show resolved Hide resolved
stan/math/prim/prob/neg_binomial_lpmf.hpp Outdated Show resolved Hide resolved
stan/math/prim/prob/neg_binomial_lpmf.hpp Show resolved Hide resolved
stan/math/prim/prob/pareto_type_2_lpdf.hpp Outdated Show resolved Hide resolved
test/prob/test_fixture_distr.hpp Show resolved Hide resolved
test/unit/math/mix/prob/neg_binomial_test.cpp Outdated Show resolved Hide resolved
@stan-buildbot
Copy link
Contributor


Name Old Result New Result Ratio Performance change( 1 - new / old )
gp_pois_regr/gp_pois_regr.stan 4.15 4.15 1.0 0.0% slower
low_dim_corr_gauss/low_dim_corr_gauss.stan 0.02 0.02 1.02 1.6% faster
eight_schools/eight_schools.stan 0.09 0.09 1.0 -0.12% slower
gp_regr/gp_regr.stan 0.2 0.2 1.0 -0.12% slower
irt_2pl/irt_2pl.stan 5.37 5.32 1.01 0.87% faster
performance.compilation 88.54 85.76 1.03 3.14% faster
low_dim_gauss_mix_collapse/low_dim_gauss_mix_collapse.stan 8.58 8.52 1.01 0.73% faster
pkpd/one_comp_mm_elim_abs.stan 27.53 27.1 1.02 1.59% faster
sir/sir.stan 115.06 115.35 1.0 -0.25% slower
gp_regr/gen_gp_data.stan 0.05 0.05 1.01 0.78% faster
low_dim_gauss_mix/low_dim_gauss_mix.stan 3.32 3.29 1.01 1.02% faster
pkpd/sim_one_comp_mm_elim_abs.stan 0.41 0.42 0.98 -2.41% slower
arK/arK.stan 1.84 1.82 1.01 0.9% faster
arma/arma.stan 0.63 0.61 1.03 3.37% faster
garch/garch.stan 0.53 0.52 1.01 1.18% faster
Mean result: 1.00843475832

Jenkins Console Log
Blue Ocean
Commit hash: 55cf763


Machine information ProductName: Mac OS X ProductVersion: 10.11.6 BuildVersion: 15G22010

CPU:
Intel(R) Xeon(R) CPU E5-1680 v2 @ 3.00GHz

G++:
Configured with: --prefix=/Applications/Xcode.app/Contents/Developer/usr --with-gxx-include-dir=/usr/include/c++/4.2.1
Apple LLVM version 7.0.2 (clang-700.1.81)
Target: x86_64-apple-darwin15.6.0
Thread model: posix

Clang:
Apple LLVM version 7.0.2 (clang-700.1.81)
Target: x86_64-apple-darwin15.6.0
Thread model: posix

Copy link
Collaborator

@andrjohns andrjohns left a comment

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

All looks good to me, merge when ready

@bbbales2 bbbales2 merged commit bf2d9a3 into develop Jul 22, 2020
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

propto tests not testing the correct things in probability unit tests
7 participants