Multivariate normal implementation

[MikeOMa]

Implementation of multivariate normal distribution just for regression with LogScore.

It does not use the scipy implementation primarily as I found scipy was much slower.

Got the idea from a paper I mentioned in the doc string of distns.multivariate_normal. I recreated one of the experiments from the paper in experiments/multivariate_normal.py and results seem very similar.

Here's a sanity check I did against a tensorflow implementation for a 3 dimensional multivariate normal, I can share my calculations for the Fisher Information soon if needed, they are currently very messy. mvn_comparison.pdf

Deleted the old implementation

I did override the old MultivariateNormal class. I salvaged a few bits (the .fit method) but the rest is all new. I could not see any score implementations so I assume this was legacy code.

changes to ngboost.py I had to add a self.multi_output attribute which defaults to False if the Dist does not have a multi_output attribute. Also, changed __set_state__ and __get_state__ to allow pickling as the MultivariateNormal implementation is a function factory like k_categorical.

P.S. This is my first pull request so I am fully expecting something to go wrong.

EDIT: Here's the derivations might have typos let me know if its' not clear MVN_DERIV.pdf

[MikeOMa]

Can't fully understand why it failed

(test_36) $ python -m poetry run pre-commit run --hook-stage manual --all-files
Trim Trailing Whitespace.................................................Passed
Fix End of Files.........................................................Passed
Check Yaml...............................................................Passed
Check JSON...............................................................Passed
Check for added large files..............................................Passed
Check for case conflicts.................................................Passed
Check for merge conflicts................................................Passed
Debug Statements (Python)................................................Passed
black....................................................................Passed
isort....................................................................Passed
Flake8...................................................................Passed
pylint on nboost*........................................................Passed
pylint on tests*.........................................................Passed

I ran this on my PC and on the second run all the tests passed. First run I got

isort....................................................................Failed
- hook id: isort
- files were modified by this hook

....

...lib/python3.8/site-packages/setuptools/distutils_patch.py:25: UserWarning: Distutils was imported before Setuptools. This usage is discouraged and may exhibit undesirable behaviors or errors. Please use Setuptools' objects directly or at least import Setuptools first.

Alongside a pile of file changes

do I just push the changed files and try again?

[alejandroschuler]

@MikeOMa yeah I think if you push again it will be fine

[alejandroschuler]

Nicely done @MikeOMa. Could you add tests for d_score and metric that numerically test against finite differences and the default metric implemented in ngboost.scores.LogScore, respectively? I want to add that kind of testing for all the distributions eventually so this would be a great example to have for future PRs too. lmk if that makes sense.

[ryan-wolbeck]

@MikeOMa Thanks for this Mike!

************* Module tests.test_distns
tests/test_distns.py:146:0: R0914: Too many local variables (17/15) (too-many-locals)

-----------------------------------
Your code has been rated at 9.92/10

------------------------------------
Your code has been rated at 10.00/10

make: *** [Makefile:12: lint] Error 1
Error: Process completed with exit code 2.

Can you triage the linter failures here? Thanks!

[ryan-wolbeck]

Can't fully understand why it failed

(test_36) $ python -m poetry run pre-commit run --hook-stage manual --all-files
Trim Trailing Whitespace.................................................Passed
Fix End of Files.........................................................Passed
Check Yaml...............................................................Passed
Check JSON...............................................................Passed
Check for added large files..............................................Passed
Check for case conflicts.................................................Passed
Check for merge conflicts................................................Passed
Debug Statements (Python)................................................Passed
black....................................................................Passed
isort....................................................................Passed
Flake8...................................................................Passed
pylint on nboost*........................................................Passed
pylint on tests*.........................................................Passed

I ran this on my PC and on the second run all the tests passed. First run I got

isort....................................................................Failed
- hook id: isort
- files were modified by this hook

....

...lib/python3.8/site-packages/setuptools/distutils_patch.py:25: UserWarning: Distutils was imported before Setuptools. This usage is discouraged and may exhibit undesirable behaviors or errors. Please use Setuptools' objects directly or at least import Setuptools first.

Alongside a pile of file changes

do I just push the changed files and try again?

Typically the linter will fix these after they are run the first time, here it's just looking for proper ordering of the imports

[MikeOMa]

Thanks for the reviews. I understand it all and I'll get to it as soon as I can.

Nicely done @MikeOMa. Could you add tests for d_score and metric that numerically test against finite differences and the default metric implemented in ngboost.scores.LogScore, respectively? I want to add that kind of testing for all the distributions eventually so this would be a great example to have for future PRs too. lmk if that makes sense.

Metric part should be straightforward, I have tested the metric implementation against that in the past so just need to redo that code.

1. For the finite differences do you mean use something like scipy.optimize.approx_fprime?

2. Should I add this to test_distns.py or make a new test file?

[alejandroschuler]

@MikeOMa yes, exactly- scipy.optimize.approx_fprime would be perfect.

I hadn't thought it through but yes, absolutely- I guess if you add it to test_distns the tests will automatically run for all the distributions! The trick here will be to set the tolerances generously enough that moderate numerical differences don't trigger constant failures.

[MikeOMa]

Certainly looking for input on the new tests,

gradients

Went for a mean squared error between finite differences and the d_score implementation. This seems pretty accurate errors are typically around 1e-5.

Strangely Laplace CRPScore fails this test but passes the metric test. Absolutely no idea why and I haven't looked into it. The error was around 1e0 where everything else is around 1e-5 or 1e-6 I think the difference was in the 2nd parameter.

Metric

Went for |estimate - D.metric()| / |D.metric()| Where | is a frobenious norm. I also tried dividing by D.n_params**2 but this seemed to be even more seed dependent. I know this will pass if an incorrect calculation |D.metric()| is large relative to the estimate but I think that kind of error is easily detectable when model fitting as the scale from line search will always be tiny.

TFixedDFfails this test and every other version of this test that I have tried even with a very high number of samples

EDIT: Opted to add the test in a new file, test_distns was getting crowded.

Had trouble getting Exponential, LogNormal and T to work with the estimate_metric_error so I did not include them.

[alejandroschuler]

@MikeOMa great job with the new tests.

iirc I think @tonyduan did the laplace CRPS so maybe he could look into that. @angeldroth and @mzjp2 worked on the T so maybe they can check that out and see if there was any issue with the implementation. Any fixes to those should be their own PRs though.

What was the nature of the issue with exponential, etc.? I wouldn't worry about dealing with them in this PR (we can do it later) but I'm just curious to know.

[MikeOMa]

@MikeOMa great job with the new tests.

Thanks let me know if you want me to change anything

iirc I think @tonyduan did the laplace CRPS so maybe he could look into that. @angeldroth and @mzjp2 worked on the T so maybe they can check that out and see if there was any issue with the implementation. Any fixes to those should be their own PRs though.

I can't be sure! It might just be because of the metrics I was using to reject.

What was the nature of the issue with exponential, etc.? I wouldn't worry about dealing with them in this PR (we can do it later) but I'm just curious to know.

I didn't really understand how the uncensoring works. I did not try extremely hard though. Posted the error at the bottom! This is the implementation I tried:

Exponential = Exponential.uncensor(LogScore)
LogNormal = LogNormal.uncensor(LogScore)
TEST_METRIC: List[DistScore] = [
    (Exponential, LogScore),

I think if the following example worked then it could be added easily:

from ngboost.distns import LogNormal
from ngboost.manifold import manifold
import numpy as np 
from ngboost.scores import LogScore
uncensored_LN = LogNormal.uncensor(LogScore)
mf = manifold(LogScore, uncensored_LN)(np.array([5,5]).reshape(-1,1))
mf.sample(1)

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-15-1ab9cb776b8b> in <module>
      5 uncensored_LN = LogNormal.uncensor(LogScore)
      6 mf = manifold(LogScore, uncensored_LN)(np.array([5,5]).reshape(-1,1))
----> 7 mf.sample(1)

TypeError: 'NoneType' object is not callable

[MikeOMa]

I posted the derivations as an edit in the first post of this PR (last line)

[ryan-wolbeck]

Hey sorry I didn't get to this sooner but I think we'll want to do a version bump on this with the addition of the new distribution. I'll submit a quick pr for that

[alejandroschuler]

oh sure, my bad, should have waited to merge

[ryan-wolbeck]

@alejandroschuler no you are good :)