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

[Fix] decrease tolerance on forecast #1655

Merged
merged 3 commits into from Oct 22, 2018

Conversation

Projects
None yet
4 participants
@skoudoro
Copy link
Member

skoudoro commented Oct 20, 2018

The aim of this PR is to fix #1654 by decreasing the solver tolerance. At the moment, the current tolerance is 0.001, which I think is too high for getting any quality solutions.

  • Advantage: The solution is more reliable and stable, and hopefully it will fix the negative variables that violate our positivity constraints
  • Drawback: We have a tiny loss of speed.

As you can see here, cvxpy (and osqp) improved its solver. In fact, some unnecessary extra variables were introduced during the problem resolution. Below, the previous version verbose result:

with cvxpy v1.0.9
-----------------------------------------------------------------
           OSQP v0.4.1  -  Operator Splitting QP Solver
              (c) Bartolomeo Stellato,  Goran Banjac
        University of Oxford  -  Stanford University 2018
-----------------------------------------------------------------
problem:  variables n = 249, constraints m = 272
          nnz(P) + nnz(A) = 7247
settings: linear system solver = qdldl,
          eps_abs = 1.0e-03, eps_rel = 1.0e-03,
          eps_prim_inf = 1.0e-04, eps_dual_inf = 1.0e-04,
          rho = 1.00e-01 (adaptive),
          sigma = 1.00e-06, alpha = 1.60, max_iter = 4000
          check_termination: on (interval 25),
          scaling: on, scaled_termination: off
          warm start: on, polish: on

iter   objective    pri res    dua res    rho        time
   1   0.0000e+00   1.00e+00   1.11e+04   1.00e-01   1.61e-03s
 175   1.4997e-01   2.58e-04   9.77e-04   5.50e-01   7.11e-03s

status:               solved
solution polish:      unsuccessful
number of iterations: 175
optimal objective:    0.1500
run time:             7.85e-03s
optimal rho estimate: 1.53e-01

here, you can see that variables n = 249 is really strange because our input matrix, in this test, has a shape of (193, 28). It will make more sense to have 221 variables. This is what we get with the new cvxpy release as you can see below.

with cvxpy v1.0.10
-----------------------------------------------------------------
           OSQP v0.4.1  -  Operator Splitting QP Solver
              (c) Bartolomeo Stellato,  Goran Banjac
        University of Oxford  -  Stanford University 2018
-----------------------------------------------------------------
problem:  variables n = 221, constraints m = 244
          nnz(P) + nnz(A) = 7191
settings: linear system solver = qdldl,
          eps_abs = 1.0e-03, eps_rel = 1.0e-03,
          eps_prim_inf = 1.0e-04, eps_dual_inf = 1.0e-04,
          rho = 1.00e-01 (adaptive),
          sigma = 1.00e-06, alpha = 1.60, max_iter = 4000
          check_termination: on (interval 25),
          scaling: on, scaled_termination: off
          warm start: on, polish: on

iter   objective    pri res    dua res    rho        time
   1   0.0000e+00   1.00e+00   1.11e+04   1.00e-01   1.27e-03s
 200   1.3808e-01   4.50e-03   2.42e-05   1.00e-01   6.75e-03s
 250   1.4711e-01   1.49e-03   9.48e-04   5.71e-01   8.38e-03s

status:               solved
solution polish:      unsuccessful
number of iterations: 250
optimal objective:    0.1471
run time:             9.35e-03s
optimal rho estimate: 3.20e-01

To me, this indicates that the 2 versions of cvxpy are solving a different problem. Fixing the tolerance issue may be a good starting point for our test to stop failing. Does it make sense @maurozucchelli @arokem @Garyfallidis ?

@@ -479,7 +480,7 @@ def lb_forecast(sh_order):
diag_lb = np.zeros(n_c)
counter = 0
for l in range(0, sh_order + 1, 2):
for m in range(-l, l + 1):
for _ in range(-l, l + 1):

This comment has been minimized.

@arokem

arokem Oct 20, 2018

Member

Why not use the following?

Suggested change Beta
for _ in range(-l, l + 1):
for l in range(0, sh_order + 1, 2):
diag_lb[l//2] = (l * (l + 1)) ** 2

And get rid of the counter variable? Wouldn't that have the same effect?

This comment has been minimized.

@skoudoro

skoudoro Oct 20, 2018

Member

I do not think so, because in your case, we will not fill diag_lb, it misses some iteration or another trick

@arokem

This comment has been minimized.

Copy link
Member

arokem commented Oct 20, 2018

Yes. This looks like a good solution to me, and fixes the test failure on my machine.

@codecov-io

This comment has been minimized.

Copy link

codecov-io commented Oct 20, 2018

Codecov Report

❗️ No coverage uploaded for pull request base (master@c2a3150). Click here to learn what that means.
The diff coverage is 100%.

Impacted file tree graph

@@            Coverage Diff            @@
##             master    #1655   +/-   ##
=========================================
  Coverage          ?   87.31%           
=========================================
  Files             ?      246           
  Lines             ?    32613           
  Branches          ?     3552           
=========================================
  Hits              ?    28477           
  Misses            ?     3275           
  Partials          ?      861
Impacted Files Coverage Δ
dipy/reconst/forecast.py 92.22% <100%> (ø)

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update c2a3150...f887b77. Read the comment docs.

counter += 1
stop = 2 * l + 1 + counter
diag_lb[counter:stop] = (l * (l + 1)) ** 2
counter = stop

This comment has been minimized.

@arokem

arokem Oct 20, 2018

Member

Nice!

@arokem

This comment has been minimized.

Copy link
Member

arokem commented Oct 20, 2018

I am +1 for merging this. Might be good for others to also take a look.

@Garyfallidis

This comment has been minimized.

Copy link
Member

Garyfallidis commented Oct 22, 2018

Thanks @skoudoro. All PRs should rebase.

@Garyfallidis Garyfallidis merged commit 7e28942 into nipy:master Oct 22, 2018

4 checks passed

codecov/patch No report found to compare against
Details
codecov/project No report found to compare against
Details
continuous-integration/appveyor/pr AppVeyor build succeeded
Details
continuous-integration/travis-ci/pr The Travis CI build passed
Details

@skoudoro skoudoro deleted the skoudoro:fix-1654-forecast branch Oct 22, 2018

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment