Showing posts with label clustering. Show all posts
Showing posts with label clustering. Show all posts

Saturday, May 09, 2020

Fun and Learn with Manning LiveProjects


The pandemic has forced most people indoors. With it, there has been a corresponding rise in online education companies offering courses to help you update your skills. Most of them follow the so-called "freemium" model, where you can watch the course videos and do the exercises, but if you want certification or support, you have to pay. In the past, I have aggressively taken advantage of these free offers, and have learned a lot in the process, so I am very grateful for the freemium model and hope it continues to exist, but nowadays I find myself being a bit more selective than I used to be. However, recently I came across a very interesting product being offered by Manning.com -- a "liveProject" on Discovering Disease Outbreaks from News Headlines, that promises to provide hands-on exposure to the consumer about Pandas, Scikit-Learn, text extraction, KMeans and DBScan clustering, as they do the project.

Although, in all fairness, while the idea is somewhat uncommon, it is not completely novel. Kaggle was there first, with their Beginner Datasets and Machine Learning Projects. However, there is one important difference -- a Manning liveProject is broken into steps, each of which has high level instructions on the prescribed approach to solve that step, but supplemented by educational material excerpted from one of Manning's books. I thought that it was a really cool idea to repurposing existing content and opening it up to a potentially different demographic. In that sense, it reminds me of the Google Places API, created by combining maps that powered Google maps and the location feedback from consumers using it.

In any case, the project setup is to discover one or more disease outbreaks from newspaper headlines collected over some time frame, and plot them on a map to discover clusters. If the cluster is over multiple geographical areas, it can be classified as a pandemic. I signed up primarily because (a) most of the clustering I have done so far involve topics and terms in text, so geographic clustering seemed new and cool to me, and (b) my son is an aspiring data scientist, and I figured that maybe we could do a bit of pair programming and learn together. However, the project turned out to be quite interesting and I got sucked in, and I ended up optimizing for (a) more than for (b). Oh well :-).

I forked the project template provided by one of the instructors, and implemented the steps of the project as Jupyter notebooks, and finally wrote up my project report (mandatory deliverable for the liveProject) as the README.md file for my fork. Steps are listed under the Methods section. At a high level, the transition from a list of newspaper headlines to disease clusters on a map (World and US) involved the following steps:

The project provides around 650 news paper headlines captured from various news agencies over an unspecified time interval, so it reflects the state of the world for some snapshot. We tag the country and city in the headlines using regular expression. Specifically, we build regexes out of the list of countries and cities in the GeoNamesCache library, and run them against the headlines, capturing the city and country names found in each headline. Of the 650 headlines, 634 could be fully resolved with both country and city names, 1 with only country name, and 15 for which neither country nor city could be found. The resolved city and country names are used to look up the latitude and longitude coordinates for each of the 634 cities, again using the GeoNamesCache. The other headlines are dropped from further analysis.

The coordinates of the cities are then plotted on a world map (Figure 1), and it looks like there are disease outbreaks all over the place during that time frame. Note that the project also additionally asks to look specifically at the United States, but in order to keep the blog post short, we don't talk about it here. But you can find those visualizations in the notebooks.


Clustering them using the K-Means algorithm helps somewhat, but basically clusters the points by longitude -- the first cluster is the Americas, the second is Europe, Africa and West Asia, and the third is South Asia and Australia.


Clustering the headlines the density based method DBSCAN produces more fine grained clusters.


The distance measure used in the clustering above was standard Euclidean distance, which is more suitable for a flat earth. For a spherical earth, a better distance measure would be the Great Circle Distance. Using that distance measure, and standard hyperparameters for DBSCAN, we get a cluster which is even more fine grained.


At this point, it looks like most of the United States and Western Europe is afflicted by one major disease or another. Given that the visualizations clearly indicated disease clusters, we wanted to find if these were all about the same disease or different diseases. We then extracted and manually looked at the "most representative" newspaper headlines (i.e., headlines that had coordinates closest to the centroids of each cluster), looking for readily identifiable diseases, then looking at the surrounding words, then using these words to look for more diseases. Using this strategy, we were able to get a count of headlines for each disease. It turned out that even though different diseases were being talked about, the dominant one was the Zika virus.

So, we filtered out the newspaper headlines for the Zika virus (around 200 of them), and reclustered them using their latitude and longitude using DBSCAN and the Great Circle Distance, and we got this.


Based on this visualization, we see that the biggest outbreak seems to be in the central part of the Americas, with two big clusters in Southern United States, Mexico, and Ecuador in South America. There is also a significant cluster in South East Asia, and one in North India, and smaller outbreaks in Western Asia. Since our pretend client is the World Health Organization (WHO), we are supposed to make a recommendation, and the recommendation is that this is a pandemic since the outbreak is across multiple countries.

This was a fun exercise, and I learned about map based clustering and visualization, which was relatively new to me, since I had never used it before. I think the liveProject idea is very powerful and has lots of potential. If you are curious about the code, the notebooks are here.

Sunday, October 06, 2019

Efficient Reverse Image Search on Lucene using Approximate Nearest Neighbors


I wrote last week about our Search Summit, an internal conferences where engineers and product managers from Elsevier, LexisNexis, and various other organizations that make up the RELX Group, get together to share ideas and best practices on search. The Search Summit is in its third year, and given that its a single track conference and runs for two days, there are limited speaking slots available, and competition for these slots is quite fierce. Our acceptance rate this year was around the 40-50% mark, which I guess puts us on par with some big external conferences in terms of exclusivity. In any case, I guess that was my long-winded way of saying that this year I didn't get accepted for a speaking slot, and was offered the opportunity to present a poster instead, which I accepted.

This was the first time I made a poster for work, and quite honestly, I didn't know where to begin. My only experience with posters was back when I was at school, where I remember that we used crayons a lot, and then helping my kids build posters for their school, which involved a lot of cutting out pictures and glueing them on to a posterboard. Fortunately, as I learned from some googling, technology has progressed since then, and now one can build the poster as a single Microsoft Powerpoint (or Apple Keynote or Google slides) slide (I used Powerpoint). Quite a few sites provide free templates for various standard poster sizes as well, I ended up choosing one from Genigraphics. One you are done designing the poster, you save it as a PDF, and upload the PDF to sites such as Staples Print Services for printing. And that's pretty much all there is to it. Here is the poster I ended up presenting at the search summit.


Basically, the poster describes a Proof of Concept (POC) application that attempts to do reverse image search, i.e., given an image, it returns a ranked list of images, similar to the given image, in the corpus. My dataset for my experiments was a set of approximately 20,000 images from the ImageCLEF 2016 competition. The notion of similarity that I use is the one learned by a RESNET network pre-trained on ImageNet. The idea is that this sort of network, trained on a large generic image corpus such as ImageNet, learns how to recognize colors, edges, shapes, and other more complex features in images, and these features can be used to provide a notion of similarity across images, much like tokens in a document.

RESNET (and other pre-trained networks like it) generate dense and (comparatively) low-dimensional vectors of features. For example, the last layer of a trained RESNET model will generate feature vectors of size 2048. On the other hand, text is typically vectorized for Information Retrieval as a vector over all the tokens (and sometimes token combinations as well) in its vocabulary, so feature vector sizes of 50,000 or 100,000 are not uncommon. In addition, such vectors tend to be quite sparse as well, since a given token is relatively unlikely to occur in many documents. Lucene (and probably other search engines) leverage the high-dimensionality and sparseness of text feature vectors using its Inverted Index data structure (aka Postings List), where the key is the token and the value is a priority queue of document ids that the token occurs in.

In order to present a ranked list of similar images, a naive (and exhaustive) solution would be to compute some measure of vector similarity (such as Cosine similarity or Euclidean similarity) between the input image and every other image in the corpus, a O(N) operation. A better solution might be to partition the similarity space using something like a KD-tree such that comparison is done against a subset of most similar images, and that is a O(log N) operation. However, the Lucene inverted index makes the search an almost O(1) operation -- first looking up by token, and then navigating down a priority queue of a subset of document IDs containing that token. My reasoning went as follows -- if I could somehow project the dense, low-dimensional vectors from RESNET back into a sparse, high-dimensional space similar to that used for text vectors without losing too much information, I would be able to do image search as efficiently as text search.

I experimented with two such projection techniques -- Random K-Means and Random Projections.

  • Random K-Means involves running a bunch of K-Means clusterings over the image set, typically 50-200 of them, each with some random value for k. Each image (point in similarity space) is represented by a sequence of its cluster memberships across all the clusterings. The intuition here is that similar images, represented by two points that are close in the similarity space, will share cluster memberships across more clusterings than dissimilar images, represented by points that are further apart. A 2D schematic on the bottom center of the poster illustrates the idea.
  • Random Projections is similar, but instead of partitioning the space into a random number of circular clusters, we partition the space itself using a (D-1) dimensional hyperplane, where D is the size of the dense feature vector (2048). The number of random projections used is typically much higher than the number of K-Means, typically 1,000-10,000. This effectively results in binary clusters where points behind the hyperplane get assigned a 0 and points in front get assigned a 1. As before, each image is represented by a sequence of its membership across all the random projections. The intuition here is similar to random K-Means. The 2D schematic on the bottom right of the poster illustrates how random projections work.

My baseline for experiments was caption based text search, i.e., similar images were found by doing text search on their captions. I took a random set of 20 images, and manually evaluated the relevance for top 10 similar images for each on a 4-point scale. Its only one person, so its not very objective, but I needed a place to start. The guideline I followed was to give it the highest rating (4) if the image is identical, or almost identical, to the source image -- possibly either the same image, or a composite containing the same image, or a minor modification. The second highest rating (3) was for images that were basically the same thing -- for example, if the source image was an X-ray of the left lung, I would rate an X-ray of the right as very similar. The third rating (2) were for images that I couldn't rate 4 or 3, but which I would expect to see in the set of similar images. Finally, the lowest rating (1) is for images that I don't expect to see in the results, i.e., irrelevant images.

I used the ratings to compute Discounted Cumulative Gain (DCG) for k=1, 3, 5, and 10, where k is the number of top results taken from the ranked list of similar images for each image, and then averaged them across the 20 images. I then computed the Ideal DCG (IDCG) for k=10 by sorting the ranked list based on my ratings, and then divided the computed DCG by the IDCG to get the Normalized Discounted Cumulative Gain (NDCG). These numbers are shown in the table under the Results section in the poster. The bolded numbers indicate the best performers. Results for k=1 serve as a sanity check, since all the searches in my evaluation returned the source image as the top result in the ranked list of similar images. Overall, the numbers seem quite good, and Random K-Means with 250 clusters produced the best results.

Random K-Means provided the best results in qualitative terms as well. The block of images at the center of the poster show a brain MRI as the source image, and the top 3 results from the baseline, the best random K-Means (250 clusters), and the best Random Projections (10000 projections). A green box indicates relevant results, and a red box indicates irrelevant results. As you can see, the caption text search based baseline produces more diverse results, but doesn't do very well otherwise. Random K-Means and Random Projections produce similar images that are more visually similar to the source image, and between the two, Random K-Means tends to produce better results.

In the future, I hope to extend this approach to text vectors as well. Specifically, I want to capitalize on the increase in recall that I noticed with caption based text search, by building embeddings (word2vec, ELMo, BERT, etc.) from the caption and treat these embeddings as additional feature vectors for the image. On the other hand, I do notice that the trend has been to use special purpose backends for vector search, such as the Non-Metric Space Library (NMSLib) based on Numpy, or the GPU-based FAISS Library for efficient similarity search, or hybrid search engines such as Vespa. So while it might be intellectually satisfying to try to make a single platform serve all your search needs, it might more worthwhile to see if one of these other toolkits might provide better results.

I also wanted to acknowledge the work of Simon Hughes. The approaches used in this poster are based heavily on ideas originally proposed by him on his DiceTechJobs/VectorsInSearch project, as well as ideas and suggestions proposed by him and others on the Relevancy and Matching Technology #vectors-in-search slack channel.



Saturday, March 03, 2018

An implementation of the Silhouette Score metric on Spark


While clustering some data on Spark recently, I needed a quantitative metric to evaluate the quality of the clustering. Couldn't find anything built-in, so (predictably) went looking on Google, where I found this Stack Overflow page discussing this very thing. However, as you can see from the accepted answer, the Silhouette score by definition is not very scalable, since it requires measuring the average distance from each point to every other point in its cluster, as well as the average distance from each point to every other point in all the other clusters. This is clearly an O(N2) operation and not likely to scale for large N. What caught my eye, however, was Geoffrey Anderson's answer, who suggested using the Simplified Silhouette score, which requires only the distance from every point to its own centroid and is thus an O(N) operation and therefore quite scalable.

However, while this is indeed the algorithm suggested in the cited paper (Multiple K Means++ Clustering of Satellite Image Using Hadoop MapReduce and Spark by Sharma, Shokeen and Mathur), this simplified metric only reports on how tight each cluster is. The intent of the original Silhouette score is to report both on tightness of individual clusters as well as separation between clusters. Turns out that such a simplified Silhouette metric does exist, and is defined in detail in this paper titled An Analysis of the Application of Simplified Silhouette to the Evaluation of k-means Clustering Validity (PDF) by Wang, et al. Interestingly, this is also the definition used in the implementation of Silhouette score in Scikit-Learn.

The latter formulation of the Simplified Silhouette Index (SSI) is shown below. For each point i, call the distance to its own cluster centroid ai, and call the distance to the nearest neighboring centroid bi. The Silhouette score for the i-th point is given by SSIi as shown below. Here ai is the indicator of cluster tightness and bi is the indicator of cluster separation. One thing to notice is that in order to compute bi you will need to compute the distance from point i to all the other centroids except its own, so the complexity of the algorithm is O(Nk) where k is the number of clusters — so still linear, but can be high for large k.


Values for SSI can vary in the range [-1, 1]. Values near 0 indicate overlapping clusters, and negative values generally indicate that the point has been wrongly clustered, since it is closer to a different cluster than it's own. Best values of SSI are close to 1. The average SSI across all the points in the corpus gives us an indication of how good the cluster is.

In addition, the distribution of SSI values in each cluster can be histogrammed as shown in this KMeans clustering and Silhouette analysis example on Scikit-Learn. A Clustering where the points are distributd approximately equally across clusters tend to be better, given similar values of average SSI.

A flow diagram for my implementation is shown below. There are 3 inputs needed, an RDD of point vectors keyed by an sequential record ID (rid), an RDD of predictions consisting of a record ID and cluster ID pair, and an RDD of centroids consisting of the clusterID and centroid vector. The RDD of points is just the input vectors with an additional sequential record ID, which you can easily provide with a zipWithIndex call on the original input RDD of vectors to cluster. The prediction RDD is the output of model.predict on the clustering model. The centroids RDD is the output of model.clusterCenters, with an additional zipWithIndex to get the cluster IDs. The point and prediction RDD are joined on the record ID, and the centroids RDD is converted to a lookup dictionary of cluster ID to cluster vector and broadcasted to the workers. The joined RDD and the broadcasted lookup table are used to compute SSI for each point. We retain the cluster ID in case we want to compute the histograms by cluster. I didn't do this because my data is too large for this information to be visually meaningful, I ended up plotting the distribution of number of points in each cluster instead.


All the input RDDs are available to the algorithm as structured text files from the clustering process. Here is the Scala code I used to implement this Simplified Silhouette index. We use Databricks Notebooks for most of our Spark analytics work, so the code below doesn't contain the boilerplate that you need to start Spark jobs directly on EMR or a standard Spark cluster. But it should be relatively simple to add that stuff in if needed.

1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
import breeze.linalg._

val OUTPUT_POINTS = "/path/to/points"
val OUTPUT_PREDICTIONS = "/path/to/predictions"
val OUTPUT_CENTROIDS = "/path/to/centroids"

// read centroids file and convert to lookup table, then broadcast to workers
val centroidsRDD = sc.textFile(OUTPUT_CENTROIDS)
  .map(line => {
    val Array(cid, cstr) = line.split('\t')
    val cvec = DenseVector(cstr.split(',').map(x => x.toDouble))
    (cid.toLong, cvec)
  })

val cid2centroid = centroidsRDD.collect.toMap
val b_cid2centroid = sc.broadcast(cid2centroid)

// read points file
val pointsRDD = sc.textFile(OUTPUT_POINTS)
  .map(line => {
    val Array(rid, pstr) = line.split('\t')
    val pvec = DenseVector(pstr.split(',').map(x => x.toDouble))
    (rid.toLong, pvec)
  })
pointsRDD.persist()

// read predictions file
val predictionsRDD = sc.textFile(OUTPUT_PREDICTIONS)
  .map(line => {
    val Array(rid, cid) = line.split('\t')
    (rid.toLong, cid.toLong)
  })
predictionsRDD.persist()

// join pointsRDD and predictionsRDD, look up centroid vectors from broadcast
// and compute a, b and SSI for all points, group by cluster ID
def euclideanDist(v1: DenseVector[Double], v2: DenseVector[Double]): Double = norm(v1 - v2, 2)

val predictedPointsRDD = pointsRDD.join(predictionsRDD)    // (rid, (pvec, cid))
  .map(rec => {
    val rid = rec._1
    val pvec = rec._2._1
    val cid = rec._2._2
    val cvec = b_cid2centroid.value(cid)
    val aDist = euclideanDist(pvec, cvec)
    val bDists = b_cid2centroid.value.toList
      .filter(cc => cc._1 != cid)
      .map(cc => {
        val otherCid = cc._1
        val otherCvec = cc._2
        val otherDist = euclideanDist(pvec, otherCvec)
        (otherCid, otherDist)
      })
    val bDist = bDists.sortWith((a, b) => a._2 < b._2).head._2
    val ssi = if (aDist == 0.0 && bDist == 0.0) 0.0D 
              else (bDist - aDist) / max(List(aDist, bDist))
    (cid, ssi)
  })
predictedPointsRDD.persist()

// compute mean SSI
val acc = sc.doubleAccumulator("acc_ssi")
val sumSSI = predictedPointsRDD.foreach(rec => acc.add(rec._2))
val meanSSI = acc.value.toDouble / predictedPointsRDD.count
print("mean SSI: %.5f\n".format(meanSSI))

predictedPointsRDD.unpersist()
pointsRDD.unpersist()
predictionsRDD.unpersist()

I have used this code to evaluate two clustering algorithms (KMeans and Bisecting KMeans) on my data, using two different approaches to vectorizing the data, and with various values of K (number of clusters). The mean SSI metric provided me the ability to reduce an entire operation to a single number that I could compare across runs. I thought this was very helpful and helped me decide which outputs to keep and which to discard without having to physically scan each output. I hope this code is useful to others who might need to evaluate their clustering algorithms on Spark.


Friday, July 15, 2016

Trip Report: Data Science Summit 2016 @ San Francisco


Earlier this week, I was at the Data Science Summit 2016 at San Francisco. This post is my trip report. The event was organized by Turi, better known as the people behind GraphLab Create. This OReilly article provides a quick backstory about the evolution of the company from a group of students and professors at Carnegie-Mellon University (CMU).

While the talks spanned a wide variety of subjects, there were a few unifying themes across the conference. They were, in no particular order, Distributed Systems, General Machine Learning (ML), Deep Learning (DL), Natural Language Processing (NLP), Recommendations, Unsupervised Learning, Online Learning, Visualization, Explainability and Graph Theory. I initially thought of classifying the talks along these lines, but then realized that a talk can span multiple categories. So I am going to cover them chronologically, and tag them with the themes that I think they belong to.

Day 1


Keynote 1 - by Pradeep Dubey, Intel Labs

[machine learning], [deep learning]

Pradeep Dubey explains how computing is moving from Inside-Out problems to Outside-In problems. Inside-Out problems are those which we understand analytically, such as a car moving up an incline. Outside-In problems are those for which we can observe the behavior but don't know how it works. An example given was trying to predict a person's social network from their purchasing behavior. Outside-In problems require lots of computing power, and Pradeep goes on to describe all the work being done at Intel to support such large ML/DL workloads on Intel CPUs using the MKL 2017 Beta toolkit (available now). His blog post provides some (quite impressive IMO) benchmark information.

Keynote 2 - by Carlos Guestrin, CEO Turi and Prof of ML, University of Washington (UoW)

[distributed systems], [machine learning], [deep learning], [unsupervised], [online learning], [visualization], [explainability]

Prof Guestrin runs through the themes of the conference using demos using GraphLab Create to demonstrate each theme, and also to show off the capabilities of the tool. His list (different from mine) included Data Distribution, Online Learning, Explanations, Automated Feature Engineering and Visualization. He is one of the authors of XGBoost: A Scalable Tree Boosting System, which he hinted could be used for automatic feature engineering. He also showed some very interesting demos around how GraphLab Create explains decisions from DL based visual recognition systems.

Deep Dive: A Dark Data System - by Chris Re

[machine learning], [unsupervised]

Chris Re describes Lattice.IO, a commercial system that automatically extracts entities from text. Initial training is by a process called Data Programming, where the trainer specifies functions that allow the system to detect patterns. The system is then fed a large body of text, and it uses the patterns and the text to teach itself about other related entities and extracts them, giving each of the extractions a probability of being correct. Humans can then tell it whether its right or wrong, and the system updates itself. Lattice.IO is based on the DeepDive Project from Stanford University. The HazyResearch/snorkel project has a number of IPython noteboosk with examples of data programming.

Petuum+ for Big ML: what next after the Parameter Server - by Prof Eric Xing, CMU

[distributed systems], [machine learning]

Prof Xing describes an alternative to the centralized Parameter Server that is usually found in Distributed ML systems. One alternative is to make it more of a peer-to-peer (P2P) system, but the communication overhead between peers becomes quite high. So the idea is to factorize out and broadcast the Sufficient Factors of the model to all workers and reconstruct the update matrix at each worker. Because Sufficient Factors are much smaller than the actual matrix, communication costs go down and such a P2P system becomes feasible. More information in this paper.

CoCoA: A communication-efficient primal-dual framework for distributed optimization, by Prof Mike Jordan, University of California, Berkeley (UCB)

[distributed systems], [machine learning]

Yet another approach to reducing communication overhead between nodes in distributed systems. Prof Mike Jordan advocates CoCoA, that uses local computation to reduce the communication overhead. As a result, experiments using CoCoA converged to the same solution 25x faster than comparable experiments without CoCoA. More information in this paper.

Evolution of the SFrame: Scalable Data Structure for ML - by Sethu Raman, Turi

[distributed systems], [machine learning], [graph theory]

This was something I have been curious about ever since I started using GraphLab Create as a student at the first course of the Coursera ML Specialization (I got a 1-year student license). SFrames allow you to treat your dataset on disk as if it was in memory, thus allowing you to work with data that is the size of your disk rather than the size of your RAM. It is open source and comes with its Python API. Unfortunately, about the only ML toolkit that uses it is GraphLab Create, so if you want to use it with Scikit-Learn, you would have to figure out how to do the vector arithmetic yourself. A graph abstraction of SFrame is SGraph, which is just a pair of SFrames (one for nodes and one for edges).

Tools for explorers, explaining and evaluating your recommender system - by Dr Chris DuBois, Turi

[recommendations], [explainability]

This was a demo of various features built into GraphLab Create to explain the behavior of a recommender system.

Advancing the Python Data Stack with Apache Arrow - by Wes McKinney, Cloudera

[distributed systems], [machine learning]

This talk was more about interoperability between different systems, using Apache Arrow as the data middleware. Apache Arrow supports an intermediate format and converters to write into different formats. Wes McKinney (creator of Pandas) has collaborated with Hadley Wickham (creator of many R packages) to seamlessly transfer data back and forth between Pandas dataframes and R dataframes. For users of Spark and Python (non-Spark), Apache Arrow also promises to one day allow reading Parquet files from standalone Python programs.

Using Graphs for improving recommendations, Amit Bhattacharya, Teachers pay teachers.

[machine learning], [recommendations], [unsupervised], [graph theory]

This is a very interesting application of graph theory to build a recommendation system. The user population are teachers who purchase books from the site, and the recommender's job is to suggest new books to purchase. A graph was built based on the user and purchase data currently available - teachers who purchase the same books are linked by an edge with weight proportional to the number of books they have in common. A few central users in this graph are labelled manually (Elementary School, High School Math, etc), and Label Propagation (functionality built into GraphLab Create) used to label the other teachers into the most probable cluster they belong to. Finally, a user is recommended books that other teachers in that cluster are purchasing.

Lessons learned from 2MM machine learning models - Dr Anthony Goldbloom, Kaggle

[machine learning]

Nice overview of general strategies used by Kaggle contestants. Popular algorithms used include XGBoost, Random Forests and Deep Learning (CNN) for Image competitions, RNN/LSTM used to a lesser extent. Also a quick glimpse into Kaggle's future plans of providing a more rounded metric of a contestant's data science skills as a whole.

Design for X - Amanda Cesari, Concur Labs

[distributed systems], [machine learning]

Amanda Cesari provides an overview of Concur Labs Data Science Stack, which includes Apache Spark and GraphLab Create. The general approach is to use Spark to analyze large volumes of data and reduce it to a medium size, then process it with GraphLab Create. This is quite pragmatic given that GraphLab Create can handle medium sized data thanks to SFrames, and has more choices in terms of algorithms compared to MLLib. She also covers a case study using Anomaly Detection using the above stack.

DSSTNE - A new deep learning framework for large sparse datasets - by Scott Le Grand, Teza Technologies

[deep learning], [recommendations]

Scott Le Grand describes DSSTNE (pronounced Destiny), Amazon's DL framework for handling super-sparse matrices. Amazon's catalogs are very large, and off-the-shelf DL packages could not handle the degree of sparseness they required. The network described looks conceptually like an Autoencoder that takes a 1-hot encoding of an item as input and generates embeddings (similarity scores, recommendations) as output. DSSTNE currently supports only the fully connected model, but has plans to support CNN and RNN in the future. Scott suggests using Nervana's Neon for building CNN and RNNs. Scott has reported that DSSTNE is 15% faster than Tensorflow, more details on his blog post.

Developing customer insights at Microsoft Visual Studio - Sai Tulasi Neppali, Microsoft

[recommendations], [unsupervised]

Sai Tulasi Neppali describes how her group used telemetry data collected from Visual Studio users to categorize them into three classes based on their usage data. The categorization helps to drive email campaigns designed to retain these users in different ways.

Deep Learning and Machine Learning: A view from the trenches - Supratim Banerjee, India Equity Partners

[machine learning], [deep learning], [recommendations], [unsupervised]

The use case here is to make sure trucks in the company's fleet are optimally loaded to maximize the revenue per kilometer. Photos of various sizes of loads were taken and image vectors extracted from them (most likely by using GraphLab Create's built in functionality using a CIFAR-10 CNN model). A k-nearest neighbors job was run on the vectors with k=10, and a graph was created with each node connected to its 10 neighbors. PageRank was run on the graph to find the top N important nodes, and these nodes were manually classified as full or empty. For a new picture, it is converted to an image vector and cosine similarity computed against these N vectors - the label of the closest vector is assigned to the new image.

The Data Science behind Bot blocking - William Cox, Distil Networks

[machine learning], [recommendations]

A very entertaining and riveting talk about the constant one-upmanship between a web bot and its human defenders. I liked the idea of fingerprinting users based on their activity, so it is easy to detect anomalous behavior against baselines with similar fingerprints.

DAY 2


Natural Language Understanding Pipelines: from keywords and grammar to inference and prediction - by Dr David Talby, Atigeo

[machine learning], [deep learning], [natural language processing]

Dr Talby describes Atigeo's Natural Language Understanding (NLU) pipeline. Their input corpus is the MIMIC dataset and their technology stack contains Apache Spark, Elasticsearch and UIMA. They started with simple dictionary based attributes and off the shelf NER, but have since created drug-disease knowledge graphs using word2vec embeddings on the critical care notes (reduced to a bag of UMLS concepts). He provides notebooks at Atigeo/nlp_demo that show some aspects of what they are doing and how.

The power of geospatial graph visualization - by Corry Lanum, Cambridge Intelligence

[visualization]

Nice presentation show how geospatial graphs (ie, overlaying data on top of maps) can increase the understanding of the data.

Exploratory Data Analysis 2.0 - Jock MacKinlay, Tableau

[visualization]

Very nice and detailed demo of Tableau features with 2 example datasets. Demonstrates the power and flexibility of the Tableau tool to develop insights from data.

Staying shallow and lean in a deep learning world - Dr Xavier Amatriain, Quora

[machine learning]

Dr Xavier Amatriain talks about the pitfalls of using DL indiscriminately. He mentions several other algorithms which are as deserving of Data Scientist's attention but which are not being considered because of DL's popularity, such as Factor Methods, Non-parameteric Bayesian Models, Online Learning, Reinforcement Lerning and Learning to Rank. He also talks about how DL models are hard to explain and mentions the Why should I trust you? paper, which lays out a technique for explanation that should be adopted by all ML models.

Matrix Factorization at scale: a comparison of scientific data analytics on Spark and MPI using three case studies - Prof Michael Mahoney, UCB

[distributed systems]

Prof Mahoney describes 3 matrix factorization techniques (NMF, PCA and CX) on Spark and Cray, and shows how using MPI locally can result in speedups. More details in his paper.

Deep Personalization - by Prof Alex Smola, CMU

[machine learning], [deep learning], [recommendations]

Prof Alex Smola talks about how to capture implicit recommendations that vary with time. Most recommendation systems do not consider how user preferences change over time. He uses survival analysis to model this change. User and Time embeddings are fed into an LSTM to produce time varying recommendations.

How to analyze 500,000h/day of human to human conversation with bleeding edge Deep Learning models - by Yishay Carmiel, Spoken Labs

[distributed systems], [machine learning], [deep learning]

Yishoy Carmiel describes his learnings when faced with processing large amounts of conversation data. He describes techniques to reduce processing times in DNNs for audio processing, including frame subsampling, using WFST beam search, using Deep Autoencoders to reduce the number of features, binarizing the weights and inputs. His changes resulted in a 35x boost in performance and he was ultimately able to process the volume within the time and expense budgeted.

The exploit-explore dilemna of music recommendation - by Dr Oscar Celma, Pandora

[recommendations], [graph theory]

Dr Oscar Celma talks about how to balance exploit (play songs that user is known to like) vs explore (play songs that the user might like based on past preferences). The decision to switch a given user from an exploit song to an explore song is made by using Markov chains over a graph of songs and user preferences. However, for any major changes in the algorithm, AB testing is done on a control small control group and rolled out to the general user base only if retention and activity metrics indicate that the change was received well.

Understanding cortical principles and building intelligent machines - by Subutai Ahmad, Numenta

[machine learning], [deep learning], [unsupervised], [online learning]

Subutai Ahmad describes the Hierarchical Temporal Memory (HTM) which is a general model of the neocortex. He describes his system as Neuroscience applied to streaming analytics, which is exactly how the human brain learns. Details of HTM are on the numenta/nupic project. He also describes NAB, a streaming anomaly detection benchmark that detects anomalies in real-time with a small amount of initial training and learns adaptively thereafter. The NAB project is available at numenta/NAB.

Product Reviews and NLP analysis and Elasticsearch - Dr Lynn Cherny, Ghostweather R&D

[machine learning], [natural language processing]

Dr Lynn Cherney delivers a nice tutorial about using NLP on Yelp! Product Review dataset. Much of it is about data analysis and loading into Pandas dataframes so it can be loaded into ElasticSearch and queried from it. She has supporting notebooks at arnicas/nlp_elasticsearch_reviews.

Scalable Learning and Recognition - Prof Ali Farhadi, University of Washington (UoW)/Allen AI

[deep learning], [online learning]

Prof Farhadi describes a number of systems that he and his students have built, that attempt to learn visually. The first is Learn EVerything About ANything (LEVAN), which learns by crawling the web for images and looking at image metadata. The second system is Visual Knowledge Extraction (VisKe), that learns interactions between entities and is able to answer questions about them. It uses factor graphs to model relationships between concepts and find the most probable explanation. The third is You Only Look Once (YOLO), which focuses on very fast recognition of images in photographs, and is described more fully in this paper. For YOLO, he uses XNOR-Net, which are CNNs with binary weights, which resulted in 30x boost in performance without loss in accuracy. They are described in this paper.

Towards Transparent AI systems: Do Humans and deep networks look at the same regions while answering visual questions? - by Prof Dhruv Batra, Virginia Tech

[deep learning], [explainability]

Prof Batra discusses how to verify that a deep learning vision system is doing what it appears to be doing. This research is an offshoot of the Visual QA (VQA) project. It generates attention maps of VQA models against human attention using visualizations (visual occlution and partial decomposition) and rank order correlation methods. The work is described in greater detail in this paper.

And here are the talks I wanted to go to but could not because I was attending one in a parallel session that I thought was more interesting. I hope to watch these once Turi publishes all the videos. If the videos are made public (which I hope they will), I will post the link to the videos once I have it.

  • Making data accessible with SQL on everything - Tomer Shiran, Apache Drill
  • Machine Learning for Analyzing complex time series - Prof Emily Fox, UoW
  • MOOCS/Turn 4: what have we learned? - Prof Daphne Koller, Coursera
  • Why did you recommend that? - Delip Rao, Joostware
  • AUC at what cost? - Alex Korbonits, Remitly
  • Next generation image processing - Dr Lukasz Kidzrinski, Deepart
  • Large scale Deep Learning with Tensorflow - Jeff Dean, Google
  • Engineering Open Machine Learning Software - Andreas Mueller, NYU
  • Machine Learning in Production - Dr Yucheng Low, Turi
  • Churn Prediction, Aggregate Features and Visualizations - Dr Srikrishna Sridhar, Turi
  • Active Learning and Human in the loop - Lukasz Biewald, Crowdflower
  • Personalizing image search with feature vectors - Rodrigo Nunes, The Real Self

Overall, I thought it was quite a nice conference. Because it is organized by a for-profit company, there were quite a few talks from employees and interesting user stories from satisfied customers and partners. However, because of the company's roots and connections in academia, there were quite a few talks from highly acclaimed researchers as well. I thought the mix between academic and business focused talks was as perfect as it could be. Looking forward to a few months of digging through all the github repositories and paper links I collected in this conference.


Sunday, May 15, 2016

Feature Reduction in Images


At ElasticON 2016 earlier this year, Doug Turnbull spoke about Image Search in his talk The Ghost in The Search Machine. The idea is simple - treating an image as a bag of pixel intensity "words". Thus a 640 x 480 color image of size is represented by 307,200 words, each word representing a single pixel. Since the image is in color, each pixel has 3 numbers, one for red, green and blue respectively. In our representation, an example word might be 251_253_248 representing (251, 253, 248) for RGB respectively, and each image would be represented by 307,200 such words. Now that the image is just a bunch of words, one can use Lucene (Doug's example was ElasticSearch) to look up similar images given an image.

A similar approach was demonstrated by Adrian Rosenbrock on his blog post Hobbits and Histograms - A How-To Guide to Building your First Image Search Engine in Python, where he uses histograms to detect similar images. Instead of a Lucene index, he uses Chi-squared distance between the histogram for the query image and the other images in the index.

I have been thinking of a hybrid approach where I use a Lucene index to store images as a sequence of pixel words, but where I increase the recall somewhat by doing a bit of feature reduction to have fewer words. So, looking at my original example, since my pixel word is composed of three values in the range 0-255, it can have 16,777,216 possible values. However if I bin them into 25 possible values per channel, my vocabulary size reduces to just 625. This post explores a couple of ideas around binning the pixels.

For my test, I used this image of a butterfly I got from the free stock photos at PhotoRack. The charts on the right represent the histograms for the pixel counts for different intensities for each of the RGB channels.






The code to generate the two images are as follows:

1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
import matplotlib.pyplot as plt

img = plt.imread("/path/to/1132868439-12104.jpg")

# original image
plt.imshow(img)
plt.show()

# channel histograms
channels = ["red", "green", "blue"]
for ix, ch in enumerate(channels):
    plt.subplot(3, 1, ix+1)
    plt.hist(img[:, :, ix].reshape(img.shape[0] * img.shape[1],), color=ch[0])
    plt.grid(True)
plt.tight_layout()
plt.show()

One way to reduce the number of features is to bin the pixel intensities for each channel in equally spaced blocks. Since I want only 25 values in each channel, I create 25 bins. I then create a new image where each pixel is replaced with the mid-point of the range for the bin it belongs to. The code to do this is shown below.

1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
import numpy as np

# binning each matrix into 25 bins and replacing with binned value
bins = np.linspace(0, 255, 25)
binned_img = ((10 * np.digitize(img, bins)) + 5).astype("uint8")
plt.imshow(binned_img)
plt.show()

# channel histograms (binned image)
for ix, ch in enumerate(channels):
    plt.subplot(3, 1, ix+1)
    plt.hist(binned_img[:, :, ix].reshape(img.shape[0] * img.shape[1],), 
             color=ch[0])
    plt.grid(True)
plt.tight_layout()
plt.show()

And here is the resulting image and the histogram. Notice that parts of the red and green channels are slightly taller and parts of the blue channel slightly shorter than the original histogram.






Finally, instead of equal sized bins, I use K-Means to cluster the pixel intensities in each channel into 25 clusters, and replace the pixel intensities with the value of the centroid of the cluster it belongs to. Here is the code to do this:

1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
from sklearn.cluster import KMeans

temp_img = np.zeros((img.shape[0] * img.shape[1], 3), dtype="uint8")
for ix, ch in enumerate(channels):
    kmeans = KMeans(n_clusters=25)
    kmeans.fit(img[:, :, ix].reshape(img.shape[0] * img.shape[1], 1))
    centroids = kmeans.cluster_centers_
    labels = kmeans.labels_
    for lx, label in enumerate(labels):
        temp_img[lx, ix] = int(centroids[label])
clustered_img = temp_img.reshape((img.shape[0], img.shape[1], img.shape[2])) 
plt.imshow(clustered_img)
plt.show()

# channel histogram (clustered image)
for ix, ch in enumerate(channels):
    plt.subplot(3, 1, ix+1)
    plt.hist(clustered_img[:, :, ix].reshape(img.shape[0] * img.shape[1],), 
             color=ch[0])
    plt.grid(True)
plt.tight_layout()
plt.show()

And here are the corresponding image and channel histograms. The histograms seem to be closer to the original than the equal sized bins above.






The nice thing about the two approaches is that the change in image quality seems to be quite minor to the human eye (at least my human eye), but these images can be represented by much smaller vectors. Hopefully this will translate to faster performance and/or smaller Lucene term indexes.

Sunday, August 30, 2015

Categorizing Medical Transcripts using DMOZ and Word2Vec


Sometime back, a colleague mentioned during a conversation that his PhD dissertation involved using DMOZ for clustering query terms. The technique seemed interesting in relation to a problem we were trying to solve at the time, but also got me thinking that perhaps this idea could be useful for categorizing documents as well. Since DMOZ provides a comprehensive hierarchy of web categories, and links to representative documents in each of these categories, the category (or categories) of an unseen document can be determined by computing the similarity of this document against the representative documents. Word2Vec vectors can be used to reduce both reference and test documents to the same latent space for comparison. To test this idea out, I used a subset of DMOZ categories under medical specialties to categorize a small collection of Medical Transcription documents I had crawled some time back. This post describes the effort.

DMOZ data is available for download as a pair of large RDF files - structure.rdf and content.rdf. The structure.rdf file contains the categories defined as path like strings. Categories are defined using the Topic tag and nested using narrow and symbolic tags. The content.rdf file also contains these path-like categories defined using the Topic tag and provides links to representative web pages for these categories using the link tag. Since we are only interested in Medical Specialties, we used the following code to parse the two RDFs in a streaming manner (SAX), and extract the links from them. Each link is then sent to Boilerpipe's Article Extractor, which strips out the HTML markup and removes irrelevant text from the page (using a combination of rules and machine learning).

1
  2
  3
  4
  5
  6
  7
  8
  9
 10
 11
 12
 13
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
// Source: src/main/scala/dmozcat/preprocessing/ReferencePageExtractor.scala
package dmozcat.preprocessing

import java.io.File
import java.io.FileWriter
import java.io.PrintWriter
import java.net.URL
import java.util.concurrent.FutureTask
import java.util.concurrent.Callable
import java.util.concurrent.TimeUnit
import java.util.concurrent.TimeoutException

import javax.xml.parsers.SAXParserFactory

import org.xml.sax.Attributes
import org.xml.sax.helpers.DefaultHandler

import de.l3s.boilerpipe.extractors.ArticleExtractor

import scala.collection.mutable.ArrayBuffer

object ReferencePageExtractor extends App {
    val dataDir = "data"
    val structXml = new File(dataDir, "structure.rdf.u8")
    val contentXml = new File(dataDir, "content.rdf.u8")
    val outputCsv = new File(dataDir, "ref-pages.csv")
    val extractor = new ReferencePageExtractor(structXml, contentXml, outputCsv)
    extractor.extract()
}

class ReferencePageExtractor(structXml: File, contentXml: File, outputCsv: File) {
    
    def extract(): Unit = {
        val factory = SAXParserFactory.newInstance()
        // parse structure RDF to get list of topics
        val structParser = factory.newSAXParser()
        val structHandler = new DmozStructHandler()
        structParser.parse(structXml, structHandler)
        val topicSet = structHandler.topicSet()
        // parse content RDF to get list of URLs for each topic
        val contentParser = factory.newSAXParser()
        val contentHandler = new DmozContentHandler(topicSet)
        contentParser.parse(contentXml, contentHandler)
        val contentUrls = contentHandler.contentUrls()
        // download pages and write to file
        val writer = new PrintWriter(new FileWriter(outputCsv), true)
        contentUrls.foreach(topicUrl => {
            val topic = topicUrl._1
            val url = topicUrl._2
            val downloadTask = new FutureTask(new Callable[String]() {
                def call(): String = {
                    try {
                        val text = ArticleExtractor.INSTANCE
                            .getText(new URL(url))
                            .replaceAll("\n", " ")
                        Console.println("Downloading %s for topic: %s"
                            .format(url, topic))
                        text
                    } catch {
                        case e: Exception => "__ERROR__"
                    }
                }
            })
            new Thread(downloadTask).start()
            try {
                val text = downloadTask.get(60, TimeUnit.SECONDS)
                if (!text.equals("__ERROR__")) {
                    writer.println("%s\t%s\t%s".format(topic, url, text))
                } else {
                    Console.println("Download Error, skipping")
                }
            } catch {
                case e: TimeoutException => Console.println("Timed out, skipping")
            }
        })
        writer.flush()
        writer.close()
    }
}

class DmozStructHandler extends DefaultHandler {
    
    val contentTopics = Set("narrow", "symbolic")
    
    var isRelevant = false
    val topics = ArrayBuffer[String]()
    
    def topicSet(): Set[String] = topics.toSet
    
    override def startElement(uri: String, localName: String, 
            qName: String, attrs: Attributes): Unit = {
        if (!isRelevant && qName.equals("Topic")) {
            val numAttrs = attrs.getLength()
            val topicName = (0 until numAttrs)
                .filter(i => attrs.getQName(i).equals("r:id"))
                .map(i => attrs.getValue(i))
                .head
            if (topicName.equals("Top/Health/Medicine/Medical_Specialties"))
                isRelevant = true
        }
        if (isRelevant && contentTopics.contains(qName)) {
            val numAttrs = attrs.getLength()
            val contentTopicName = (0 until numAttrs)
                .filter(i => attrs.getQName(i).equals("r:resource"))
                .map(i => attrs.getValue(i))
                .map(v => if (v.indexOf(':') > -1)
                    v.substring(v.indexOf(':') + 1) else v)
                .head
            topics += contentTopicName
        }
    }
    
    override def endElement(uri: String, localName: String, 
            qName: String): Unit = {
        if (isRelevant && qName.equals("Topic")) isRelevant = false
    }
}

class DmozContentHandler(topics: Set[String]) extends DefaultHandler {
    
    var isRelevant = false
    var currentTopicName: String = null
    val contents = ArrayBuffer[(String, String)]()
    
    def contentUrls(): List[(String, String)] = contents.toList
    
    override def startElement(uri: String, localName: String, 
            qName: String, attrs: Attributes): Unit = {
        if (!isRelevant && qName.equals("Topic")) {
            val numAttrs = attrs.getLength()
            val topicName = (0 until numAttrs)
                .filter(i => attrs.getQName(i).equals("r:id"))
                .map(i => attrs.getValue(i))
                .head
            if (topics.contains(topicName)) {
                isRelevant = true
                currentTopicName = topicName
            }
        }
        if (isRelevant && qName.equals("link")) {
            val numAttrs = attrs.getLength()
            val link = (0 until numAttrs)
                .filter(i => attrs.getQName(i).equals("r:resource"))
                .map(i => attrs.getValue(i))
                .head
            contents += ((currentTopicName, link))
        }
    }
    
    override def endElement(uri: String, localName: String, 
            qName: String): Unit = {
        if (isRelevant && qName.equals("Topic")) {
            isRelevant = false
            currentTopicName = null
        }
    }
}

Output of this code is a tab separated file containing the category name, file URL, and the text. There are 464 records in all and 44 unique categories. The data looks something like this.

1
2
3
Top/Health/Medicine/Medical_Specialties/Aerospace_Medicine  http://www.civilavmed.com/  Home How are you as an AME adapting to the new ...
Top/Health/Medicine/Medical_Specialties/Aerospace_Medicine  http://www.asma.org/        Sunday, September 20, 2015 8:00 AM ICASM 2015 The Congress ...
...

The next step is to vectorize the reference text (the third column in the file shown above). We use OpenNLP's models to segment the text into sentences and extract noun phrases, using the excellent Scala interface provided by the NLP Tools Project from the University of Washington.

One of the artifacts published by the Word2Vec project is the set of vectors for around 3 million words and phrases trained on the Google News dataset. This information is released as a binary file, and we use Gensim's Word2Vec module to convert to TSV format that our code can use. Using the noun phrase retrieved using OpenNLP, we construct unigrams, bigrams and trigrams and look them up against this word vector dictionary. The vectors for the ngrams so found are added and normalized to build a single vector for each category.

1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
// Source: src/main/scala/dmozcat/vectorize/TextVectorizer.scala
package dmozcat.vectorize

import scala.Array.canBuildFrom

import org.apache.log4j.Level
import org.apache.log4j.Logger
import org.apache.spark.SparkConf
import org.apache.spark.SparkContext
import org.apache.spark.rdd.RDD.rddToPairRDDFunctions

import breeze.linalg.DenseVector
import breeze.linalg.InjectNumericOps
import breeze.linalg.norm

object TextVectorizer {

    def main(args: Array[String]): Unit = {
        Logger.getLogger("org.apache.spark").setLevel(Level.ERROR)
        Logger.getLogger("org.apache.spark.storage.BlockManager").setLevel(Level.ERROR)
        
        // arguments
        val awsAccessKey = args(0)
        val awsSecretKey = args(1)
        val inputFile = args(2)
        val word2vecFile = args(3)
        val outputDir = args(4)
        
        val conf = new SparkConf()
        conf.setAppName("TextVectorizer")
        
        val sc = new SparkContext(conf)
        
        sc.hadoopConfiguration.set("fs.s3n.awsAccessKeyId", awsAccessKey)
        sc.hadoopConfiguration.set("fs.s3n.awsSecretAccessKey", awsSecretKey)
        
        val input = sc.textFile(inputFile)
            .map(line => parseLine(line))
            .mapPartitions(p => extractNGrams(p)) // (key, ngram)
            .map(kv => (kv._2, kv._1))            // (ngram, key)
            
        val wordVectors = sc.textFile(word2vecFile)
            .map(line => {
                val Array(word, vecstr) = line.split("\t")
                val vector = new DenseVector(vecstr.split(",").map(_.toDouble))
                (word.toLowerCase, vector)        // (ngram, vector)
            })
        
        // join input to wordVectors by word
        val inputVectors = input.join(wordVectors)  // (ngram, (key, vector))
            .map(nkv => (nkv._2._1, nkv._2._2))     // (key, vector)
            .aggregateByKey((0, DenseVector.zeros[Double](300)))(
                (acc, value) => (acc._1 + 1, acc._2 + value),
                (acc1, acc2) => (acc1._1 + acc2._1, acc1._2 + acc2._2))
            .mapValues(countVec => 
                (1.0D / countVec._1) * countVec._2) // (key, mean(vector))
            
        // save document (id, vector) pair as flat file
        inputVectors.map(kvec => {
            val key = kvec._1
            val value = (kvec._2 / norm(kvec._2, 2))
                .toArray
                .map("%.5f".format(_))
                .mkString(",")
            "%s\t%s".format(key, value)
        }).saveAsTextFile(outputDir)
    }
        
    def parseLine(line: String): (String, String) = {
        val cols = line.split("\t")
        val key = cols.head
        val text = cols.last
        (key, text)
    }
    
    def extractNGrams(p: Iterator[(String, String)]): 
            Iterator[(String, String)] = {
        val t2n = new NGramExtractor()
        p.flatMap(keyText => t2n.ngrams(keyText))
    }
}

The actual vectorization logic is encapsulated inside the NGramExtractor, whose code is shown below. If you have used the OpenNLP API directly, you will appreciate the convenience of using the NLPTools API. Plus, since the models are captured as dependencies, you no longer have to explicitly deal with the OpenNLP model files. The only downside of the NLPTools API is the number of dependencies you must specify in your build file - this is because NLPTools unify access to a bunch of underlying tools, and each of them have different licenses, so having different JAR files for each is how NLPTools deals with this.

1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
// Source: src/main/scala/dmozcat/vectorize/NGramExtractor.scala
package dmozcat.vectorize

import edu.knowitall.tool.sentence.OpenNlpSentencer
import edu.knowitall.tool.postag.OpenNlpPostagger
import edu.knowitall.tool.tokenize.OpenNlpTokenizer
import edu.knowitall.tool.chunk.OpenNlpChunker

import scala.collection.mutable.ArrayBuffer

class NGramExtractor {

    val sentencer = new OpenNlpSentencer
    val postagger = new OpenNlpPostagger
    val tokenizer = new OpenNlpTokenizer
    val chunker = new OpenNlpChunker
    
    def ngrams(keyText: (String, String)): List[(String, String)] = {
        val key = keyText._1
        val text = keyText._2
        // segment text into sentences
        val sentences = sentencer.segment(text)
        // extract noun phrases from sentences
        val nounPhrases = sentences.flatMap(segment => {
            val    sentence = segment.text
            val chunks = chunker.chunk(sentence)
            chunks.filter(chunk => chunk.chunk.endsWith("-NP"))
                .map(chunk => (chunk.string, chunk.chunk))
                .foldLeft(List.empty[String])((acc, x) => x match {
                    case (s, "B-NP") => s :: acc
                    case (s, "I-NP") => acc.head + " " + s :: acc.tail
                }).reverse
        })
        // extract ngrams (n=1,2,3) from noun phrases
        val ngrams = nounPhrases.flatMap(nounPhrase => {
            val words = nounPhrase.toLowerCase.split(" ")
            words.size match {
                case 0 => List()
                case 1 => words
                case 2 => words ++ words.sliding(2).map(_.mkString("_"))
                case _ => words ++ 
                    words.sliding(2).map(_.mkString("_")) ++
                    words.sliding(3).map(_.mkString("_"))
            }
        })
        ngrams.map(ngram => (key, ngram)).toList
    }
}

The text from the documents to be categorized are similarly vectorized by extracting the noun phrases and looking up vectors for its unigrams, bigrams and trigrams against the Word2Vec dictionary. In order to use the same code to do this, we preprocess our documents from a set of files in a directory to a single file containing tab-separated filename and corresponding text.

1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
// Source: src/main/scala/dmozcat/preprocessing/TestDataReformatter.scala
package dmozcat.preprocessing

import java.io.File
import scala.io.Source
import java.io.PrintWriter
import java.io.FileWriter

object TestDataReformatter extends App {
    
    val inputDir = new File("/Users/palsujit/Projects/med_data/mtcrawler/texts")
    val outputFile = new PrintWriter(new FileWriter(new File("/tmp/testdata.csv")), true)
    inputDir.listFiles().foreach(inputFile => {
        val filename = inputFile.getName()
        val text = Source.fromFile(inputFile)
            .getLines
            .map(line => if (line.endsWith(".")) line else line + ".")
            .mkString(" ")
            .replaceAll("\\s+", " ")
        outputFile.println("%s\t%s".format(filename, text))
    })
    outputFile.flush()
    outputFile.close()
}

We vectorize the reference and test datasets in two separate Amazon EMR jobs. Amazon officially announced support for Spark on EMR couple of months ago (at the Spark Summit in San Francisco in June this year), and I've been meaning to try it out. Their instructions for running Spark on EMR is very detailed and worked perfectly for me, the only thing you need to know is to switch to the "advanced option" when defining your cluster.

Essentially Spark on EMR works in client mode. The first step invokes an Amazon custom job to execute a Hadoop dfs command to copy your JAR from S3 to the driver box. The next step actually runs your JAR against the Spark cluster. The JAR file needs to have (in this case) the OpenNLP models baked in, which can be done by invoking "sbt assembly".

The output of the job run against the file of reference text is a file of (category name, category vector) pairs, and a file of (filename, vector) pairs when run against the test dataset.

Next, we compute the cosine similarity of each document vector (in word2vec space) against each of the category vectors (also in word2vec space). Since there are only 44 category vectors, we convert this to a 44x300 matrix (the word2vec vectors provided are 300 dimensional) and broadcast it to the workers. Figuring out the best category is simply a matter of multiplying this matrix by the document vector. Both the matrix entries and the vector are normalized by their L2 norms, so the result is a vector of cosine similarities. We then find the top N (3) similarities, and their associated category names.

1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
// Source: src/main/scala/dmozcat/vectorize/KNearestCategories.scala
package dmozcat.vectorize

import scala.Array.canBuildFrom
import org.apache.log4j.Level
import org.apache.log4j.Logger
import org.apache.spark.SparkConf
import org.apache.spark.SparkContext
import breeze.linalg.DenseVector
import breeze.linalg.norm
import breeze.linalg.DenseMatrix

object KNearestCategories {
    
    def main(args: Array[String]): Unit = {
        
        Logger.getLogger("org.apache.spark").setLevel(Level.ERROR)
        Logger.getLogger("org.apache.spark.storage.BlockManager").setLevel(Level.ERROR)
        
        // arguments
        val awsAccessKey = args(0)
        val awsSecretKey = args(1)
        val refVectorsDir = args(2)
        val testVectorsDir = args(3)
        val bestcatsFile = args(4)
        
        val conf = new SparkConf()
        conf.setAppName("KNearestCategories")
        
        val sc = new SparkContext(conf)
        
        sc.hadoopConfiguration.set("fs.s3n.awsAccessKeyId", awsAccessKey)
        sc.hadoopConfiguration.set("fs.s3n.awsSecretAccessKey", awsSecretKey)
    
        // read reference vectors and broadcast them to workers for
        // replicated join
        val refVectors = sc.textFile(refVectorsDir)
            .map(line => {
                val Array(catStr, vecStr) = line.split("\t")
                val cat = catStr.split("/").last
                val vec = new DenseVector(vecStr.split(",").map(_.toDouble))
                val l2 = norm(vec, 2.0)
                (cat, vec / l2)
            }).collect
        val categories = refVectors.map(_._1)
        // we want to take the Array[DenseVector[Double]] and convert
        // it to DenseMatrix[Double] so we can do matrix-vector multiplication
        // for computing similarities later
        val nrows = categories.size
        val ncols = refVectors.map(_._2).head.length
        val catVectors = refVectors.map(_._2)
            .reduce((a, b) => DenseVector.vertcat(a, b))
            .toArray
        val catMatrix = new DenseMatrix[Double](ncols, nrows, catVectors).t
        // broadcast it
        val bCategories = sc.broadcast(categories)
        val bCatMatrix = sc.broadcast(catMatrix) 
        
        // read test vectors representing each test document
        val testVectors = sc.textFile(testVectorsDir)
            .map(line => {
                val Array(filename, vecStr) = line.split("\t")
                val vec = DenseVector(vecStr.split(",").map(_.toDouble))
                val l2 = norm(vec, 2.0)
                (filename, vec / l2)
            })
            .map(fileVec => 
                (fileVec._1, topCategories(fileVec._2, 3,
                    bCatMatrix.value, bCategories.value)))
            
        testVectors.map(fileCats => "%s\t%s".format(
                fileCats._1, fileCats._2.mkString(", ")))
            .coalesce(1)
            .saveAsTextFile(bestcatsFile)
    }

    def topCategories(vec: DenseVector[Double], n: Int,
            catMatrix: DenseMatrix[Double],
            categories: Array[String]): List[String] = {
        val cosims = catMatrix * vec
        cosims.toArray.zipWithIndex
            .sortWith((a, b) => a._1 > b._1)
            .take(n)                           // argmax(n)
            .map(simIdx => (categories(simIdx._2), simIdx._1)) // (cat,sim)
            .map(catSim => "%s (%.3f)".format(catSim._1, catSim._2))
            .toList
    }
}

This is also a Spark job executed on EMR. One thing I learned was that you cannot use broadcast variables with Scala objects that extend App because of scoping issues, as discussed in this StackOverflow page and this Spark JIRA. The final output looks like this:

1
 2
 3
 4
 5
 6
 7
 8
 9
10
2022.txt  Pain_Management (0.875), Podiatry (0.831), Sports_Medicine (0.826)
1928.txt  Pain_Management (0.858), Behavioral_Medicine (0.843), Podiatry (0.840)
1296.txt  Surgery (0.864), Podiatry (0.849), Radiotherapy (0.832)
0996.txt  Cardiology (0.864), Pain_Management (0.847), Radiotherapy (0.846)
2000.txt  Pain_Management (0.751), Cardiology (0.736), Radiotherapy (0.735)
0853.txt  Radiotherapy (0.773), Cardiology (0.767), Podiatry (0.734)
2228.txt  Pain_Management (0.918), Podiatry (0.904), Surgery (0.889)
0361.txt  Surgery (0.859), Podiatry (0.848), Pain_Management (0.843)
0916.txt  Cardiology (0.823), Radiotherapy (0.820), Pain_Management (0.801)
3078.txt  Palliative_Care (0.903), Surgery (0.887), Critical_Care (0.885)

Thats all I have for today. It was a lot of fun playing with Word2Vec and the KnowItAll APIs for OpenNLP, as well as using Spark on EMR. Nowadays people consider DMOZ slightly outdated compared to Wikipedia, but it provides a nice hierarchical topic taxonomy that can be used effectively for this kind of work. I hope you enjoyed reading it, hopefully it gave you some ideas to categorize your own data. All the code for this post is available on GitHub here.