Linear Regression 3 Ways in Julia | Hacker News Reader

Linear Regression 3 Ways in Julia

boss-level.com

81 pointsluu20 comments

every time i see a Julia post i want to convert my Matlab code and start anew - then i remember seeing a blog some time ago where the vectorizations i performed to speed up my Matlab code would be faster as for-loops in Julia. Is this still the case? If not, i'll probably consider the task much more seriously

That is still the case.

The following blog post goes into the reasons, and the remedies. In particular, there are ways to automatically devectorise your code.

http://julialang.org/blog/2013/09/fast-numeric/

I'm a little surprised they haven't done loop fusion optimizations like Haskell's "vector". Do they plan to?

I'm not sure what you mean. The Haskell vector library doesn't use any compiler side intelligence.

If you mean loop merging, the Repa 4 library does something pretty interesting called series fusion. But it also needs a compiler plugin so that it can do that optimization internally.

I'm actually working on my own vectorization / fusion friendly numerical tools in Haskell. One philosophical point that's pretty darn important Is making it obvious and clear how to get good performance in a robust way. Many auto vectorization codes (eg try auto vectorization in clang) require a bit of work to tickle correctly. So in many ways, I think the best approach is to encourage idioms that don't require compiler smarts.

I think you mean stream fusion. While cool, it doesn't actually work in real cases yet. Indeed the ticket to implement it has been closed as invalid because it needs more research.

http://ghc.haskell.org/trac/ghc/ticket/915

It's still one of the most impressive moves to a sufficiently smart compiler.

Ummmm. False. No compiler support is needed for Stream fusion. Stream fusion is used in the follow widely used Haskell libraries: Vector, Text, and ByteString. In fact, for any fixed k, I can manufacture natural look Haskell code in these libraries that would be roughly k times slower sans stream fusion.

The different fusion strategies have different trade offs, so it could also be that I'm not reading your with the correct nuance.

Don't get me wrong, I'm not criticising what it _can_ do. It's just that the really exciting idea is the SSC idea: where you can just write clean code and have the compiler work out the best execution strategy.

At the moment, SF seems to support natural looking code, but not always truly natural code. It's tantalisingly close, but it doesn't seem to be getting any closer.

could you unpack your remarks in a bit more detail?

devectorizing is a terrible step back in expressiveness. It shouldn't be hard to make a compiler that can compute the example without a bunch of temporary arrays. I would be shocked if heavy duty packages like Matlab make all those temporary arrays.

On the other hand, I've seen (and written) plenty of code that went through contortions to use vector operations in IDL (http://en.wikipedia.org/wiki/IDL_(programming_language)). It would be great to have a language where for loops and vector expressions are fast.

Matlab does in fact make all those temporary arrays. So does R. See my comment above [1] – the idea is not that you should devectorize everything, but that you can if you need to. Automatic devectorization of general purpose code is not something anyone knows how to do – even in Haskell, which is the best candidate language for such things – see moomin's comment below [2].

[1] https://news.ycombinator.com/item?id=6341406

[2] https://news.ycombinator.com/item?id=6341187

In many cases, our vectorized code is as fast as Matlab or R. In some cases, we need to optimize the way our memory allocation and GC works to reduce the memory pressure created by temporaries. With that, julia should be as good as anything out there on vectorized codes.

There are certain things that are still slower when you write vectorized code than they would be if you manually devectorized them. The only reason that's not true in Matlab or Python is that loops are so slow.

The real point of Dahua's "Fast Numeric" blog post [1] is that you can devectorize in the high-level language if you need to, not that you always should. Most code is not a performance bottleneck and can use the most convenient form. In cases where you identify a real hotspot [2], in Matlab you'd have to resort to writing C code via MEX. In Julia, you do one or more of the following: (a) change some copying function to a mutating version by sticking a "!" on the end of it; (b) move an array allocation out of a loop and pass the pre-allocated array into a mutating function inside the loop; (c) manually devectorize some code if absolutely necessary – or use the Devectorize package to do it automatically. At no point do you ever have to resort to writing code in any other language.

Also keep in mind that slowness is not the only problem one encounters with numerical code – I've personally have much more frequent trouble with memory usage. I've written far too much insanely convoluted Matlab code to try to avoid blowing out the memory on a 100GB machine (or 1TB these days). It's generally just barely doable using find, sparse, sub2ind and ind2sub. If I had to resort to writing C code every time I encountered one of these memory problems, I would just write the damned thing in C and scrap Matlab altogether. On the other hand in Julia for such situations, I can just write a few for loops and not only do I avoid blowing out my memory, but the code gets faster too, and I'm never required to leave my high-level environment.

[1] https://news.ycombinator.com/item?id=6329237

[2] http://docs.julialang.org/en/latest/stdlib/profile/

Does Julia do any of the copy-on-write tricks that other numerical packages do? I.e. when I make a copy of an array does it immediately need space for a full copy of the original?

No, there is no copy-on-write. We try to keep the behavior and performance as transparent as possible and cow is distinctly non-transparent. If you want a copy, you have to explicitly make a copy.

The julia-users mailing list discussion of this blog post:

https://groups.google.com/forum/#!topic/julia-users/T5YYKQzl...

I'm impressed that these characters can be used in code:

Σ₁,Σ₀, σ²

Many moder programming languages allow arbitrary unicode symbols as variable names. I am now thinking about doing this in Python.

Fair enough. I use Sublime and write in python on Ubuntu. How might I go about incorporating Unicode?

Hmm, here's a thought: Instead of naming a function "sum", call it Σ, and similarly use √ instead of "sqrt". This would enable you to write things like

rms = √(1/len(x)*Σ(x))

Which is really quite cool from a mathematical view point.

Might be nice to read, but I wouldn't want to edit code written like that, since the symbols could well be harder and slower to write on a normal keyboard / editor.

Reply on news.ycombinator.com