The main optimization at that stage seems to be preallocating the weights. I don't know pandas but such a thing would have been possible without dropping any of

The main optimization at that stage seems to be preallocating the weights. I don't know pandas but such a thing would have been possible without dropping any of the linalg libraries I do know how to use. I doubt the author's C++ implementations beat BLAS/LAPACK, but since they're not shown I can only guess. I've done stuff like this before but the tooling is really no fun, somewhere between 2 and 3 I'd just write it all in C++. Changing the interface just to get parallelism out seems not great - give it to the user for free if the array is long enough - but maybe it was more reasonable for the non-trivial real problem.