Faster NNLS algorithm

[rdeits]

Hi! Thanks for putting this package together.

I've been working on building a fast Julia implementation of the NNLS algorithm from Lawson and Hansen. It's really just a brain-dead port of the original Fortran code (extracted from the copy in scipy.optimize), but having it in Julia means I'll be able to experiment with numeric types other than double.

Some quick benchmarking shows that the ported code is 100 to 1000 times faster than NonNegLeastSquares.nnls(), and its memory allocation is reduced by a similar factor. Would you be interested in bringing that implementation into NonNegLeastSquares.jl? I'm afraid the code is longer and messier than what's currently here (the original source is a 400 line Fortran subroutine full of gotos), but I'm working on cleaning it up. And I think the performance benefit could be worth it.

I've pushed my code to a new repo here: https://github.com/rdeits/NNLS.jl if you'd like to take a look.

[ahwillia]

Nice! I'd be happy to include it. Algorithmically speaking, is it the same as NNLS? Just a faster implementation?

[rdeits]

As far as I understand, yes. It's precisely this algorithm: https://github.com/scipy/scipy/blob/v0.18.1/scipy/optimize/nnls/nnls.f which references the same book as https://github.com/ahwillia/NonNegLeastSquares.jl/blob/master/src/nnls.jl#L12

I think that a lot of the performance comes from the fact that the Fortran code has been optimized to do almost everything in-place, and I've kept that benefit in the Julia implementation. In fact, you can even construct an NNLSWorkspace object for the solver to use, which lets the entire NNLS algorithm run with nearly zero memory allocation (32 bytes currently. Hopefully I can get that to 0).

[ahwillia]

Very cool. Would love to see if this could benefit the other algorithms that have been implemented​ here. Go ahead and open a PR!

On Mar 25, 2017 9:07 PM, "Robin Deits" notifications@github.com wrote:

As far as I understand, yes. It's precisely this algorithm: https://github.com/scipy/scipy/blob/v0.18.1/scipy/optimize/nnls/nnls.f which references the same book as https://github.com/ahwillia/ NonNegLeastSquares.jl/blob/master/src/nnls.jl#L12

I think that a lot of the performance comes from the fact that the Fortran code has been optimized to do almost everything in-place, and I've kept that benefit in the Julia implementation. In fact, you can even construct an NNLSWorkspace object for the solver to use, which lets the entire NNLS algorithm run with nearly zero memory allocation (32 bytes currently. Hopefully I can get that to 0).

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/ahwillia/NonNegLeastSquares.jl/issues/2#issuecomment-289256959, or mute the thread https://github.com/notifications/unsubscribe-auth/AAm20ToklhgP9O40rary5nDxAtAEw4g0ks5rpeSXgaJpZM4MpSAQ .

[rdeits]

Ok, will do. One quick question: do you want to continue supporting the tol argument to nnls()? I don't think that tolerance actually appears in Lawson's book, and there's no numeric tolerance in the Fortran implementation. So far, I haven't encountered any cases where I've wanted that tolerance, but perhaps you have?

[ahwillia]

It was only there for numerical reasons, but if it that is unnecessary I'm happy to take it out.

[gsaxena888]

@rdeits I came accross this entry when searching for a fast nnls solver in java. By any chance, would you know if this Julia code could be written in Java with more-or-less comparable performance or is there something "special" about either Julia or Fortran that Java would have a hard (or impossible) time emulating?

On a completely different but related note, does this NNLS julia implementation you wrote use some of the techniques described in the Fast NNLS paper, https://pdfs.semanticscholar.org/c524/cab2ec917c2adbbb453088479c7f800d9f76.pdf? (If I understood the basic premise of that paper, the fast NNLS version would further speed up NNLS by a factor of up to 20x.)

Finally, if Java could work (and as a bonus if the fast NNLS would work too), you wouldn't know anyone (yourself/grad student/colleague etc) interested in prototyping a Java solution for either the original NNLS or the Fast NNLS?

Thanks in advance, @rdeits.

[rdeits]

I'm not a Java expert, but I believe that it should be possible to implement this algorithm efficiently in Java. Julia does have an advantage in that it can create new numeric types with no overhead (just like C++), whereas Java Objects have a small but significant overhead that makes them hard to use to create custom numeric types. I suspect that would mean that Julia and Java would behave similarly for primitive types like float or double, but Julia would be better for custom numeric types.

However, given that this code already exists in the form of a well-established Fortran library, as shown here: https://github.com/scipy/scipy/blob/v0.18.1/scipy/optimize/nnls/nnls.f I suspect you'd be better off simply calling that Fortran library from Java rather than re-implementing it.

I'm not familiar with the Fast NNLS paper, so I do not believe its results are implemented here.