While calculating team ratings for a machine learning-based March Madness prediction, I ran into a couple of situations where my code got slow as I expanded it to include all the teams over all the seasons. By slow I mean several hours, and that was longer than I was willing to wait. I needed to be able to recalculate ratings on-demand, in a few minutes at most.
For my offensive-defensive rating, I started with an implementation similar to the one in Offensive and defensive team ratings for the Premier League 2014-2015:
It's looping over all the rows and columns in the matrix many times. With college basketball having over three-hundred teams, this wasn't going to work for me. I figured out that it could be vectorized and that numpy could handle it more efficiently than me:
Not only is it faster, but I would argue the more concise code is easier to understand as well.
For my Markov Chain ranking, I started with an implementation similar to the one in A Markov Chain ranking of Premier League teams (14/15 season):
I don't mean to disparage these two posts in any way. They are awesome and really helped me. I just had a different situation and had to worry about performance, and here I figured out that I could get rid of both loops:
Eigenvectors to the rescue. The eigenvector for the largest eigenvalue is the stationary distribution we are trying to find and we don't need to do the 100,000 step random walk.
Showing posts with label linear algebra. Show all posts
Showing posts with label linear algebra. Show all posts
Thursday, February 23, 2017
Sunday, February 8, 2015
Eigenvectors and PageRank
Eigenvectors are another topic that I used in a math class or two but never really understood what they were for. As I mentioned in my last post Bipartite Graphs and Google Adwords, I like how the professors of the Coursera class Mining Massive Datasets gave great examples of how to apply some academic concepts to real world problems. Their PageRank material is no different as it explains how eigenvectors are used in the algorithm. Here is my summary of the class' link analysis chapter which I'm sure is greatly simplified from what Google really does, but I think it's interesting nonetheless.
Early search engines used an inverted index which I covered in my tf-idf post. That approach is vulnerable to term spam. People would put invisible text on their pages (same color as the background) to trick web crawlers into thinking that the page was about a topic it wasn't. The PageRank algorithm solves that problem by giving more importance to pages that have many in-links, and especially important in-links. Even if someone creates millions of fake pages that link to their real page, it's importance will not be that high because the fake pages will have a low importance as no one is linking to them.
The web can be thought of as a directed graph where web pages are the vertexes and links are the arcs. The basic idea of PageRank is that you simulate web surfers starting at random spots in the web graph and then following random out-links. The most important pages, the pages with the most in-links, will end up with the most surfers.
To implement this simulation, consider a transition matrix M which has n rows and columns where n is the number of pages crawled. Mij has the value 1/k if page j has k arcs out and one of them is to page i, otherwise it has the value 0. Start with a vector v0 which has all elements set to 1/n. The first step is to multiply v0 by M, the second step is to multiple by M2, and so on (this is an example of a Markov process which means you can make predictions on its future state based only on the present state and do not need the history).
If the web graph is one strongly connected component and doesn't have any dead ends, then v = Mv is the limiting distribution. The limit would be reached when multiplying the distribution by M another time doesn't change the distribution. Vector v is an eigenvector of M because v = λMv. Vector v is also the principal eigenvector because M is a stochastic matrix and the eigenvalue associated with the principal eigenvector is 1. In practice it takes 50-75 iterations to reach this limit.
M is too big for us to use Gaussian elimination. M is also very sparse, so it makes sense to store only the non-zero elements. An entire column can be represented by the out-degree of a page and the row number of the non-zero elements (since those element's values will be 1 divided by the out-degree). That means we only need to store a little over 4 bytes for each non-zero element. Then we can solve for v using MapReduce and v tells us the "page rank" of each page.
In reality the web is not strongly connected. There are dead ends and there are spider traps, or cycles, in the graph. This problem is solved by the concept of taxation. Each surfer has a small probability of teleporting to a random page instead of following an out-link. There are also other algorithm variations for dealing with link spam. For example, besides page content, you also consider the link text or words near the link, so you are getting other people's take on what the page is about instead of relying solely on the page owner.
Early search engines used an inverted index which I covered in my tf-idf post. That approach is vulnerable to term spam. People would put invisible text on their pages (same color as the background) to trick web crawlers into thinking that the page was about a topic it wasn't. The PageRank algorithm solves that problem by giving more importance to pages that have many in-links, and especially important in-links. Even if someone creates millions of fake pages that link to their real page, it's importance will not be that high because the fake pages will have a low importance as no one is linking to them.
The web can be thought of as a directed graph where web pages are the vertexes and links are the arcs. The basic idea of PageRank is that you simulate web surfers starting at random spots in the web graph and then following random out-links. The most important pages, the pages with the most in-links, will end up with the most surfers.
To implement this simulation, consider a transition matrix M which has n rows and columns where n is the number of pages crawled. Mij has the value 1/k if page j has k arcs out and one of them is to page i, otherwise it has the value 0. Start with a vector v0 which has all elements set to 1/n. The first step is to multiply v0 by M, the second step is to multiple by M2, and so on (this is an example of a Markov process which means you can make predictions on its future state based only on the present state and do not need the history).
If the web graph is one strongly connected component and doesn't have any dead ends, then v = Mv is the limiting distribution. The limit would be reached when multiplying the distribution by M another time doesn't change the distribution. Vector v is an eigenvector of M because v = λMv. Vector v is also the principal eigenvector because M is a stochastic matrix and the eigenvalue associated with the principal eigenvector is 1. In practice it takes 50-75 iterations to reach this limit.
M is too big for us to use Gaussian elimination. M is also very sparse, so it makes sense to store only the non-zero elements. An entire column can be represented by the out-degree of a page and the row number of the non-zero elements (since those element's values will be 1 divided by the out-degree). That means we only need to store a little over 4 bytes for each non-zero element. Then we can solve for v using MapReduce and v tells us the "page rank" of each page.
In reality the web is not strongly connected. There are dead ends and there are spider traps, or cycles, in the graph. This problem is solved by the concept of taxation. Each surfer has a small probability of teleporting to a random page instead of following an out-link. There are also other algorithm variations for dealing with link spam. For example, besides page content, you also consider the link text or words near the link, so you are getting other people's take on what the page is about instead of relying solely on the page owner.
UPDATE:
This same idea can be applied to sports ratings where links are replaced by something like goals scored. My Machine Learning for NCAA Basketball Prediction - Performance Edition post has more details and some code.
Labels:
big data,
coursera,
linear algebra
Saturday, August 9, 2014
Vectorization: A Matrix Multiplication Example
It was in the Coursera Machine Learning class that I first learned about vectorization. When dealing with vectors and matrices, instead of writing loops and computing operations on scalars, you can instead just perform the operations directly on the higher dimensional data structures. It is a type of parallelism where one processor performs operations on multiple data simultaneously, but does not involve any concurrency.
Vectorization is made possible by array programming languages and libraries like Octave, R, and NumPy. It can be a compiler optimization or some kind of interface can be made available to the programmer so they can indicate the operations to vectorize. At a lower level, this is implemented with SIMD processor instructions or using, for example, 32-bit instructions to simulate vector computations of 16 or 8-bit types.
A simple example to demonstrate this, I think, is matrix multiplication. In Java (no vectorization possible without using the JNI) a simple (but partially optimized) implementation looks like this:
Multiplying 1000 x 1000 matrices containing random numbers takes about 6 seconds on my laptop. You can see how this will not work for machine learning calculations on matrices with millions of elements.
In Python, I can do the same multiplication using NumPy like this:
It takes about .02 seconds. Obviously this isn't a proper performance comparison, but the two order of magnitude difference is still illustrative of the power of vectorization.
Vectorization is made possible by array programming languages and libraries like Octave, R, and NumPy. It can be a compiler optimization or some kind of interface can be made available to the programmer so they can indicate the operations to vectorize. At a lower level, this is implemented with SIMD processor instructions or using, for example, 32-bit instructions to simulate vector computations of 16 or 8-bit types.
A simple example to demonstrate this, I think, is matrix multiplication. In Java (no vectorization possible without using the JNI) a simple (but partially optimized) implementation looks like this:
Multiplying 1000 x 1000 matrices containing random numbers takes about 6 seconds on my laptop. You can see how this will not work for machine learning calculations on matrices with millions of elements.
In Python, I can do the same multiplication using NumPy like this:
It takes about .02 seconds. Obviously this isn't a proper performance comparison, but the two order of magnitude difference is still illustrative of the power of vectorization.
Labels:
java,
linear algebra,
parallelism,
python