Cemoody/bhtsne Barnes-Hut t-SNE

## Introduction

This PR presents the Barnes-Hut implementation of t-SNE. t-SNE is used to visualize high-dimensional data in a low dimensional space that attempts preserve the pairwise high-dimensional similarities in a low-dimensional embedding. The Barnes-Hut algorithm, which is used by astrophysicists to perform N-body simulations, allows the calculation of the t-SNE embedding in O(NlogN) time instead of O(N^2). 

This effectively allows us to learn embeddings of data sets with millions of elements instead of tens of thousands.

## Example

An iPython Notebook is available where we've put together simple install instructions & a demo of both methods:
http://nbviewer.ipython.org/urls/gist.githubusercontent.com/cemoody/01135ef2f26837548360/raw/76ce7f0bac916a516501ead719b513b22430cad0/Barnes-Hut%20t-SNE%20Demo.ipynb

![bhtsne](https://cloud.githubusercontent.com/assets/3419930/5566392/24d8b7dc-8edb-11e4-9e67-4f319fc6dc59.png)

## Performance
This compares the timings for the Barnes-Hut approximation against the current method (note the log y-axis).

![timing](https://cloud.githubusercontent.com/assets/3419930/5566592/9f8ba956-8ee2-11e4-88e9-2504485f06c2.png)

The following diagram shows the relative speedup between both versions:
![speedup](https://cloud.githubusercontent.com/assets/3419930/5566360/f4802f94-8ed9-11e4-891c-da999a0426e8.png)

## TODO 

- [ ] Create fit & transform methods to update using new data
- [ ] Change sparse inputs to use method='standard'
- [ ] Change dimension > 4 to use method='standard'
- [ ] Warnings for 64bit being transformed to 32bit
- [ ] Do we have memory leaks?
- [ ] Segfaults on data larger than 50k elements
- [x] PEP8 the code
- [x] Include usage documentation
- [x] Remove extra imports, print statements, pdb statements
- [x] Ensure python 3 works
- [x] Ensure output dimensionality works for 2D or 3D
- [x] Am I using memviews or returning full arrays appropriately?
- [x] Incorporate into SKLearn
- [x] Remove GIL as much as possible in Cython code
- [x] Add answer tests
- [x] Ensure old t-SNE tests works
- [x] Changed perplexity calculation to nearest neighbors
- [x] Changed positive gradient calc to nearest neighbors

## Learn

- You can read more about the technique in general here:

http://lvdmaaten.github.io/tsne/

- The Barnes-Hut approximation is here:

http://lvdmaaten.github.io/publications/papers/JMLR_2014.pdf

Initially I wrote it from scratch, following the paper. I used to do
astrocomputing, where I'd written oct-tree codes before (which are similar
to this Barnes-Hut approximation). I then compared results with the LvdM's
code and used those answers to debug my own. Where his code uses
Vantage-Point trees we use BallTree, and there are other major structural
differences. But I also changed the some of names of my methods and
variables to mirror his (for example `compute_non_edge_forces` or `sum_Q`)
for easier comparison & debugging. This PR is not a 'direct translation': I
did not take a chunk of C++ code and try to rewrite it in Cython. But at
the same time, I certainly read & understood LvdM's code and it helped
speed up progress on this one.

So I'm not sure if that means it's derivative or not. To be safe, I'll
contact LvdM and ask for authorization.

chris

On Mon, Dec 29, 2014 at 12:51 AM, Olivier Grisel <notifications@github.com>
wrote:

> In sklearn/manifold/bhtsne.pyx
> <https://github.com/scikit-learn/scikit-learn/pull/4025#discussion-diff-22305672>
> :
>
> > +# distutils: extra_compile_args = -fopenmp
> > +# distutils: extra_link_args = -fopenmp
> > +# cython: boundscheck=False
> > +# cython: wraparound=False
> > +# cython: cdivision=True
> > +import numpy as np
> > +from cython.parallel cimport prange, parallel
> > +from libc.stdio cimport printf
> > +cimport numpy as np
> > +cimport cython
> > +cimport openmp
> > +
> > +# Implementation by Chris Moody & Nick Travers
> > +# original code by Laurens van der Maaten
> > +# for reference implemenations and papers describing the technique:
> > +# http://homepage.tudelft.nl/19j49/t-SNE.html
>
> When you say: "original code" you mean that this implementation is a
> direct translation in Cython of the original C++ code?
>
> If so this is a violation of the of the legal statement of the original
> implementation that says in the README:
>
> You are free to use, modify, or redistribute this software in any way you want, but only for non-commercial purposes. The use of the software is at your own risk; the authors are not responsible for any damage as a result from errors in the software.
>
> scikit-learn is a BSD licensed project, it should only include compatibly
> licensed code. In particular the non-commercial exclusion is not compatible.
>
> If the implementation in this PR is not a derivative work there is no
> problem. If it was based on the original implementation instead of
> reimplementing from the paper they we need to ask the original for
> authorization before being able to consider the inclusion of this
> contribution in scikit-learn.
>
> —
> Reply to this email directly or view it on GitHub
> <https://github.com/scikit-learn/scikit-learn/pull/4025/files#r22305672>.
>