Fun with Kaggle Leaderboards

Our team recently fought hard to earn our way into the top 10% of Kagglers in the AXA telematics challenge. I even wrote a post about it here where I talk about some of the feature generation I contributed. It was our first time in a real data science competition, and I’ve got to say, it was quite a rush. But, after it was all over, it was hard to describe the heat of the race to the top. Ranks on the public leaderboard and ROC curve scores did little to provide a realistic feel of the competitors around us.

I wanted to change that.

I was reading a blog post about plateaus in the 2013 Titanic Kaggle challenge. Trevor, the author, plotted the score vs. rank for all competitors on the public leaderboard. In this way, he was able to see plateaus where there were larger or smaller changes in rank for a given change in score.

Static Plot: Any given day

I wanted to test some web scraping skills using SelectorGadget and the XML and rvest packages, so the first approach to the problem was just get a quick look at the current leaderboard.


# AXA Standings
axaROC<-html("") %>%
  html_nodes("td:nth-child(4) , .abbr")
axaRank<-html("") %>%

In the above snippet, I used SelectorGadget to isolate the html nodes that stored the scores and ranks for all competitors. I can’t stress how valuable the extension is for scraping well-organized webpages.

The output of html_nodes() is unusable as-is as numeric data. With xmlValue, it was a relatively simple matter of trimming the fat around them to get at the scores and ranks themselves.

score<-as.numeric(unlist(lapply(axaROC,function(x) xmlValue(x,trim=T))))
rank<-as.numeric(unlist(lapply(axaRank,function(x) xmlValue(x,trim=T))))

Final AXA Leaderboard

Finally, ggplot2 created this view of the public leaderboard on the final day. I also manually entered the benchmarks that the Kaggle Admins had put in place as well as 0.66 which was a popular post that provided R code to achieve a leaderboard score of 0.66 using simple logistic regression.

AXA Static

Score vs. Rank for the AXA Telematics Challenge. Our team, Vivi’s Angels, is highlighted in red.

The ggplot2 code was straightforward, but I include it here in case you are interested.


  annotate("text",x=350,y=0.92,label="Vivi's Angels")+
  geom_hline(yintercept=0.5,color="dark blue")+
  annotate("text",x=500,y=0.485,label="All 1's Benchmark")+
  geom_hline(yintercept=0.53154, color="blue")+
  annotate("text",x=500,y=0.545,label="Trip Length Benchmark")+
  geom_hline(yintercept=0.66,color="light blue")+
  annotate("text",x=500,y=0.675,label="Speed Quantiles")+
  xlab("Private Leaderboard Rank")+
  ylab("Private Leaderboard Score")+
  ggtitle("AXA Telematics: Final Standings")+

This gave some sense of climbing a mountain. We were able to quickly see not only that we had traversed a steep climb from some of the benchmarks, but also that ahead of us in the top 5 and 10% was a very steep climb.

Animated Plot: What about during the competition?

We had a sense of where we ended up, but I wanted to see the changing landscape and our path to the final standings. Well, this was easy enough to plot in R after downloading the AXA public leaderboard history as a .csv, which Kaggle conveniently provides at the bottom of each competition’s public leaderboard page. This file contained all timestamped and identified updates made to the leaderboard since the beginning of the competition. Juicy data!

The animation package let me send it a function that produced a series of plots and then stitch those plots together as an animated gif. This was a technique a classmate of mine covered at the NYCDSA Bootcamp, so I was eager to try it myself.

The first step was to point to the downloaded .csv file in the working directory. It would have been nice to do this automatically, however, this required cookies which would be different for all users. I thought it wasn’t really a huge trade-off in usability to have the user download and unzip the csv themselves.

# Info about the competition
### download the zip file and place it in the working directory
### go to the public leaderboard page and download the file at the end

# What is your user or team name?
teamName<-"Vivi's Angels"

# Start of the competition

# End of the competition

Next, I had to create the function that would generate each plot, for any given day. Here, genPlot takes several inputs:
daysIn describing the days from the beginning of the competition (the day for which a plot will be generated)
yMin and yMax which define the range of scores to show on the generated plot
higherBetter which determines which direction to plot the axes. For most competitions, higher is better
PLHist is the public leaderboard dataframe loaded above
startDate which defines “t=0” for the plots of interest.

  thisPL<-leaderboard %>%
    filter(SubmissionDate<thisTime) %>%
    group_by(TeamId) %>%

  } else{

  # Add the day's ranks
  thisPL<-thisPL %>%

If a team of interest is specified outside of the function (teamID), this part of genPlot will get the coordinates (score and rank) for the team, if it has made a submission. It will also set up a timeStamp to float in the upper left corner of the plot and a team stamp to identify the team name, if it is present (this uses the annotate function in the grid package).

  teamScore<-ifelse(teamID %in% thisPL$TeamId,as.numeric(thisPL[which(thisPL$TeamId==teamID),2]),NA)
  teamRank<-ifelse(teamID %in% thisPL$TeamId,as.numeric(thisPL[which(thisPL$TeamId==teamID),3]),NA)

  timeStamp <- grobTree(textGrob(as.character(thisTime), x=0.1,  y=0.95, hjust=0,
                                 gp=gpar(col="red", fontsize=18, fontface="bold")))
  teamStamp <- grobTree(textGrob(teamName, x=0.1,  y=0.9, hjust=0,
                                 gp=gpar(col="blue", fontsize=18, fontface="bold")))

The last part of the function sets up the output plot. If a team of interest has been supplied, it will add blue crosshairs to signify the team’s location and it will annotate the plot with the teamStamp.

    #   xlim(c(totalTeams,0))+

  if( print(p)

Cool! Now we have a code that will plot the team’s position on any given day, specified by days from the start date. Now we want to animate this to show progress. Luckily, the animate package makes this really simple. I defined the pl.animate function which, in addition to passing backdoor arguments to genPlot, allows the user to specify the animation interval (animInt) and the start day of the animation. Finally, the saveGIF() function calls the pl.animate function, specifying an interval (in seconds) for each frame to appear. I thought 4 frames per second was just the right speed for a ~90 day competition.

pl.animate <- function(animInt=1,dayzero=start,...) {
  lapply(seq(0,as.numeric(latest-dayzero),animInt), function(i) {

saveGIF(pl.animate(yMin=0.4,yMax=1), interval = .25,

The Climb

Animated Public Leaderboard for the AXA Telematics Competition. Go, team, go! Click on the image to (re)load the animated GIF as WordPress tends to bug out.

And so there it is, if you click on the image above, you can observe the changing layout of scores and ranks throughout 90 days of competition and follow the crosshair of Vivi’s angels as we climb up the hill.


I’ll definitely be using this before entering and during competitions to see how things are changing. Even in the static image, plateaus other than the benchmark may provide some clue as to a useful forum post or idea someone has discovered. The benefit of the animated plot shows just how quickly teams make it into the top scoring tiers.

Finally, I invite you to fork this code from my github page and use it for your own team or to see any leaderboards you may be interested in.

Kaggle AXA Telematics: Trip Matching

The following is adapted from my RPub. The code can be found on my Github repo.

I’m proud to say I was part of the NYCDSA Bootcamp Team, Vivi’s Angels, for the AXA Telematics Kaggle competition. Now that the competition is over and the scores have been tallied, we are all learning so much from those who have started to share their approaches to solving the problem of identifying the primary owner of a car merely from the x-y data of the trip he or she took.

While no one here at the bootcamp is walking away with the $30,000 prize money, we do want to share our approach and code to such an intriguing and fun challenge. I was part of the feature engineering team and spent most of my time developing the trip-matching feature for use in our GBM modeling. This post will discuss how (and why) we wanted to get this feature up and running.

Early in the competition, we came across posts like this one that spoke of the success of trip matching in improving model results. Some domain experts even posted in the forums, saying that in the real world, GPS data was the best indicator of whether a driver was the owner of the car or not. The problem was that AXA rotated and flipped some trips so that they couldn’t be matched in this way! The first order of business was to re-align the trips, then compare them to one another to test for matches.

Once the team decided that trip matching should be incorporated into the features of our model, there were three challenges:

  • Translating the approach from Python to R
  • Optimizing computational time (the original codes took days to run on multiple cores)
  • Finding the best way to implement the data as a feature

Coding for Trip Matching

This code was written in R to do two things:

  • Rotate/flip trips so that most of the trip was in the first quadrant.
  • Compare trips of roughly similar shape using euclidean distance.

With 200 trips to compare for each of over 2700 drivers, parallelization was a must. library(doParallel) enabled the use of foreach loops that greatly sped up computation, even on only three cores on my Macbook Air. The flag, .combine=rbind takes the output from each loop and rbinds it to the dataframe in a manner.

cl <- makeCluster(1)

dataDir <- "data/"
drivers <- list.files(dataDir)

similarTrips<-foreach(driver=drivers,.combine=rbind) %dopar%{
 # Each driver's 200 trips were loaded as a binary file before being compared with one another


Rotating/flipping trips

Before the comparison for loops, transformations were applied via mutate() in the dplyr package. These transformations rotate each trip so that the last point is on the positive x-axis by changing the coordinate system. The trips were then flipped if more than half of the x or y coordinates were negative. This was useful later when trips were compared roughly before calculating their euclidean distances. When the euclidean distances were calculated, the rotated y values were used.

test<-drives %>%
 group_by(tripID) %>%
 mutate(x = x, y = y,
 rot.x.flip=ifelse(sum(rot.x<0)>floor(rows/2), -rot.x, rot.x),
 rot.y.flip=ifelse(sum(rot.y<0)>floor(rows/2), -rot.y, rot.y)

To give a visual idea of this rotation, Figure 1 and Figure 2 show the raw and rotated trips, respectively.

Raw x-y data for driver 1, all 200 trips.

Fig. 1: Raw x-y data for driver 1, all 200 trips.

Fig. 2 Trips rotated such that the last point is on the positive x-axis, then flipped such that most points are in the first quadrant.

Fig. 2: Trips rotated such that the last point is on the positive x-axis, then flipped such that most points are in the first quadrant.

Even at this point, the naked eye can see some similarity among the rotated, spaghetti-esque mess. The next steps will automate the trip comparison to ultimately reveal similar trips.

Comparing individual trips

Only unique pairs of trips were considered, greatly reducing computational time.

 # begin a nested loop to check all UNIQUE combinations of trips
 for(i in 1:199){
 focus<-select(test[test$tripID==i,], foc.x=rot.x, foc.y=rot.y)
 for(k in (i+1):200) {
 compare<-select(test[test$tripID==k,], cmp.x=rot.x, cmp.y=rot.y)

Once two trips were selected for potential comparison, their “footprint” was compared. Trips that were more than 20% different in x- or y- range were not compared, since that large a difference would probably indicate a poor match. This selection by if statements reduced computational time by over 30% over checking every pair, generally resulting in approximately 1000 comparisons.

 if(!(diff(range(focus$foc.x))<0.8*diff(range(compare$cmp.x))) &
 !(diff(range(compare$cmp.x))<0.8*diff(range(focus$foc.x))) &
 !(diff(range(focus$foc.y))<0.8*diff(range(compare$cmp.y))) &

In order to compare trips with euclidean distance, their vectors must be the same length. The forum post that inspired this post imputed extra values by repeating values from the shorter of the two trips. Initially, that was my approach; however, it was faster (by about 10%) and more conservative (no imputed values) to truncate the longer of the two vectors.

Finally, the euclidean distance was calculated between the focusTrim and compareTrim y-values. I had to apply a normalization in order for the euclidean distance to be an acurate metric for similarity of trips of all sizes. Without such a normalization, trips with larger y- values would have greater euclidean distance, even if they were just as good a match. I chose to normalize by the mean of the largest y-values from each of the two trips.

 # Trim Trips
 trimLength<- min(nrow(compare),nrow(focus))

 #Calculate the Euclidean Distance
 data.frame(driver=driverID, tripA=i,tripB=k,
 } # end if loop
 } # end comparison loop
 } # end focus loop

Visualizing the result

Arranging the output (similarTrips) by eucDist, we get a “rank” of trip similarity, with the most similar trips towards the top. Keep in mind that eucDist doesn’t hold much value, but serves to indicate which trips are relatively similar.

> head(arrange(similarTrips,eucDist),n=10)
 driver tripA tripB eucDist
1     1    76   182 0.6139725
2     1    46    86 0.7146952
3     1    33    68 0.8082965
4     1   139   160 0.8460641
5     1    60    68 0.9580528
6     1    31   160 1.0072537
7     1    17   170 1.0351923
8     1    68   154 1.0593948
9     1    33    91 1.0752706
10    1    43   129 1.1583642

After looking at the spaghetti plot (Fig. 1) for so long, plotting the unique trips from the tripA and tripB column produces a gratifyingly organized verification of some success in matching trips:

Unique trips from the top 10 matched pairs for Driver 1.  Note this is the plot of rotated trips.

Fig. 3: Unique trips from the top 10 matched pairs for Driver 1. Note this is the plot of rotated trips.

Our team now had a trip matching algorithm of our own (written in R instead of Python), and it was relatively fast at about 50 seconds per driver per cpu (about 38 cpu hours for the whole set).

Which eucDist should we choose?

This was not an easy question to answer. Choosing several drivers at random and plotting the above facets along with the ranked list of eucDist values provided grounds to produce a clustering of trips by the same techniques used in social network analysis. For instance, in the above output for Driver 1, trips 33, 60, 68, 91, and 154 belong to the same undirected network cluster. This script on my github covers how to produce plots like the ones below, which shows clusters for all linked trips for Driver 1 and Driver 10. Here, green edges indicate particularly close euclidean distances between trips and I have filtered out vertices (trips) with fewer than 3 connections.


If all of this manual verification sounds like it was mindnumbing, don’t worry: the monotony of verifying a few of these top matches for several drivers was outweighed by the sheer joy of finding order amidst the chaos.

Trip matches as features

Now that we could get a rough idea of which trips were similar, we needed to somehow code these trips for our gradient boosting machines model. We kicked around a lot of ideas:

Maybe any trip that matches another trip is the primary owner.

This turned out to produce bad results. Intuitively, we should have realized that even non-owners of the vehicles might take the same routes from time to time.

Maybe trips that have a match should “round” our initial prediction to 0 or 1.

Again, this reduced our performance on the public leaderboard. Although it would have been so nice to apply such a manual post-processing boost to our predictions!

Maybe trips should receive a score based on their similarity at different eucDist thresholds.

This was what ultimately produced an improvement in our public leaderboard score. I chose 5 thresholds, ranging from very conservative to very lax and assigned scores to them such that trips that matched at conservative thresholds (and more lax ones as a result) received a higher score than those that only matched at the “liberal thresholds.”

A Happy Ending

We had a blast in our first-ever Kaggle competition, finally able to edge our way into the top 10%. More importantly, we learned a TON about team organization and workflow and this was my first project working in the Agile/Spring Development cycle. Admittedly, the trip matching algorithm was successful weeks before its results could be turned into useful predictors. This is where having a great team helped. Conversations over coffee or lunch often produced new insights. I was lucky to be part of such a dedicated and hard-working team.

Further Reading

If you liked this blog post, please check out work by my teammates (links upcoming):

Julian Asano – Team Organization
Tim Schmeier – Visualization and Dynamic Time Warping
Sylvie Lardeux – Trip Attributes
Jason Liu – GBM Mechanic

Kaggle Forum Machine Learning Cloud


I had been wanting to explore some text scraping and manipulation in R, so I decided to try my hand at a word cloud. I know they’re not especially quantitative, but people like them; they’re cute.

With Kaggle on the brain, I had been reading a lot of the forums. After each competition, the top three teams are (usually) interviewed by the Kaggle folks and are asked things like:

  • “What was your background prior to entering the Challenge?
  • “What was your most important insight into the dataset?”
  • “Which tools did you use?”
  • “What preprocessing and supervised learning methods did you use?”

Those last two were what I was after. My goal was to get a word cloud of the tools and methods being employed by the top performers on Kaggle.

So I got to work gently scraping the Kaggle blogs.

Finding Kaggle blog links (and titles)

I constructed a blog page URL that held several posts per page. With no visible limit to the number of pages (that might otherwise signal a place for this loop to end), I checked to see the number of post links on the next page. If it was 1 or more, the loop proceeded.

On each page, I used the CSS markers to grab the post titles and corresponding links, saving the result into a data frame.


# For constructing the blog links

## Cycle through blog pages and get all post links
while(nextLinkCount>=1){ # keep going until there are 0 links on the page

  titleXML<-html(paste0(baseURL,page,pNum)) %>%
    html_nodes("#homepage h1 a:nth-child(1)")
  # Get the post titles and links on the current page
  postTitles<-unlist(lapply(titleXML,function(x) xmlValue(x)))
  postLinks<-unlist(lapply(titleXML,function(x) xmlGetAttr(x,"href")))

  # Check the number of post links on the next page
  nextLinkCount<-length(html(paste0(baseURL,page,pNum+1)) %>%
                            html_nodes("#homepage h1 a:nth-child(1)"))


Just the Winners, please.

Kaggle is a great resource for seeing how folks are using data science techniques. However, they also blog about quite a few other things. The three main categories are:
Data Science News and Editorials (33 posts as of 3/28)
Kaggle News (118 posts as of 3/28)
Tutorials and Winners’ Interviews (116 posts as of 3/28)

Obviously, I only wanted the last of these, so for each of the blog links, I checked the category, later filtering so that all I was left with was the Tutorials & Winners’ Interviews.

for(i in 1:length(blogDF$link)){
  catXML<-html(blogDF$link[i]) %>%
    html_nodes(".categories a")


# get just the posts about winning teams
winPosts<-filter(blogDF,category=="Tutorials and Winners\' Interviews")

Idea: Check the tags!

Now I had the relevant posts. My first thought was to check the tags on each page. Some posts were tagged in a helpful way, i.e. “GBM, Python, scikitlearn”. However, I quickly saw that this was not helpful: Relatively few posts were tagged and even then, some posts were tagged oddly. One post was even tagged “Celine Dion.” A dead end here, but My Heart Will Go On.

Grabbing the (important) Text

There was a sea of text, and as much as Kagglers like their SVMs and GBMs, those mentions would be easily swamped out by competition-specific words in the word cloud. So, I needed to isolate the answers to just the important questions. I had seen that a few posts had interview questions marked with helpful node ids tags like “question” or even “h1” but that others were simply strong or bold. I resorted to something fairly brute-force, just grabbing all the paragraph text.

# get the text from the winning pages
for(i in 1:nrow(winPosts)){
  pXML<-html(winPosts$link[i]) %>%
    html_nodes("div p, h3")
  postText<-unlist(lapply(pXML,function(x) xmlValue(x)))

What was the question?

With all the text from the interviews scraped, I had to pick out the questions (and specifically, the questions that asked about methods and tools). I needed all of the questions because I would later be extracting all the text in the post between the important question and the *next* question. Spot checking the question list, I saw that most that asked about approaches involved the words “tools, “methods” or “learning.”

# Tag an entry as a question

# Tag an entry as involving tools/methods/learning i.e. "Which machine learning methods did you use?"
entrytxt$techniques<-grepl(x=entrytxt$text,pattern="tools|methods|learning") &
  (sapply(gregexpr("\\W+", entrytxt$text), length) + 1)<11 # Filter out lengthy erroneous labels

Some answers

I set an index of questions and an index of machine learning questions, looping through the text data frame to extract only the text that was in response to the methods/tools questions of interest.

# Index for all questions
Qindex <- which(entrytxt$isQuestion)

# Index for ML questions
MLindex <- which(entrytxt$techniques)

# Find the paragraphs between a ML question and the next question
for(i in 1:length(Qindex)){
  if(Qindex[i] %in% MLindex){
    startMLAns <- Qindex[i]+1
    endMLAns <- Qindex[i+1]-1
    MLpara <- rbind(MLpara,data.frame(start=startMLAns,end=endMLAns))

# Take the actual text from the index list
for(i in 1:nrow(MLpara)){
  MLAns<-rbind(MLAns,paste(as.character(entrytxt$text[MLpara$start[i]:MLpara$end[i]]),collapse=" "))

answerText<-paste(MLAns,collapse=" ")

# Write the text file to a temporary directory for the wordcloud

Finally, it was time for a word cloud. Using tm and wordcloud, I followed this tutorial to get up and running. I had to remove certain “meaningless” words before I was able to see the random forest through the trees.


meaninglessWords<-c("data","used", "using","features","code","different","models","table","problem","like","category", "first","second","third","one","two","three","approach","many","number","blog","model","feature", "also","learning","top","score","competition","create","user","much")

winSpeak<-Corpus(DirSource("tmp/")) %>%
  tm_map(removeWords, stopwords("english")) %>%
  tm_map(stemDocument) %>%
  tm_map(removeWords,meaninglessWords) %>%
wordcloud(winSpeak, scale=c(3,0.5),
          colors=brewer.pal(8, "Dark2")

The right tool for the right job… and a giant toolbox

Sort of unsurprisingly, no particular machine learning method or coding language pops out of the background. “Random,” is gigantic, but this is not necessarily because of the mention of random trees or random forest. The cloud shows that there is no one tool we should all be using and no clear-cut method to win Kaggle.

On a cool note, I recently stumbled across Kaggle’s beta wiki and posts like this one on Random Forests that link to relevant information and even competition winners who used that random forests primarily. This will be an awesome tool for even the casual Kaggler (or an aspiring data scientist seeking to expand his toolset).