Category: Uncategorized

Project Focus:

Data Sci­ence • R

Project Overview

As both a hob­by and to bet­ter my data sci­ence skills, I com­pete in Kag­gle’s Data Sci­ence com­pe­ti­tions. This par­tic­u­lar com­pe­ti­tion stands out as it forced me to crit­i­cal­ly exam­ine the meth­ods I was using, and because I end­ed up peak­ing at num­ber 2 on the leader­board (final­ly fin­ish­ing in the 31st, 89th percentile). 

Skills Used

The crux of the script was a com­bi­na­tion of sim­pler sequence guess­ing mod­els (mode & fre­quen­cy) and a more robust tun­ing mod­el. The dif­fi­cul­ty was in tun­ing the script to select the cor­rect mod­el (mode, fre­quen­cy or lin­ear regres­sion) to solve a sequence. 

---

Project

I recent­ly sub­mit­ted an entry into Kag­gle’s Inte­ger Sequence Learn­ing Com­pe­ti­tion, and as of right now, it’s 2nd on the leader­board. I want­ed to go over my method­ol­o­gy, and how I got to my submission.

The Com­pe­ti­tion

For the com­pe­ti­tion, we’re asked to pre­dict the last term of a sequence. The descrip­tion giv­en by the competition:

7. You read that cor­rect­ly. That’s the start to a real inte­ger sequence, the pow­ers of primes. Want some­thing eas­i­er? How about the next num­ber in 0, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55? If you answered 89, you may enjoy this chal­lenge. Your com­put­er may find it con­sid­er­ably less enjoyable.

The On-Line Ency­clo­pe­dia of Inte­ger Sequences is a 50+ year effort by math­e­mati­cians the world over to cat­a­log sequences of inte­gers. If it has a pat­tern, it’s prob­a­bly in the OEIS, and prob­a­bly described with amaz­ing detail. This com­pe­ti­tion chal­lenges you cre­ate a machine learn­ing algo­rithm capa­ble of guess­ing the next num­ber in an inte­ger sequence. While this sounds like pat­tern recog­ni­tion in its most basic form, a quick look at the data will con­vince you this is any­thing but basic!

The Method­ol­o­gy

The sequences are made up of a lin­ear sequences, log­a­rith­mic sequences, sequences with a mod­u­lus, and many oth­er odd­i­ties. In the below graph, gen­er­at­ed through this script by Kag­gle User Calin Uiore­anu, illus­trates some of these sequences and gives an idea of what we’re work­ing with.

Now that we have a bet­ter idea of what we’re work­ing with, lets get into the methodology.

We’re giv­en 3 files to start:

• train.csv
• test.csv
• sample_submission.csv

All fair­ly nor­mal as far as Kag­gle com­pe­ti­tions go. I won’t be using the train­ing data, and will be jump­ing straight into the test data.

library(plyr)
library(dplyr)
library(readr)
library(stringr)
options(scipen=999) #Prevents Scientific Notation for Large Numbers

Here I’m just load­ing the nec­es­sary libaries we’ll need for the code. Note the “options(scipen=999),” that’s impor­tant to pre­vent R from for­mat­ting the larg­er num­bers into sci­en­tif­ic notation.

Also note that read­ing test.csv pro­duces a data frame of two vari­ables, Id and Sequence. The Id vari­ables are inte­gers, and are exact­ly how we want them. The Sequence vari­able is in strings, so we’ll need to con­vert that to a list a of num­bers. It’s a pret­ty sim­ple use of the str_split() func­tion; I par­tic­u­lar­ly liked how Kag­gle User William Cukier­s­ki wrote his method, so I imple­ment­ed that in my code as well.

parseSequence <- function(str) {
  return(as.numeric(str_split(str, ",")[[1]]))
}

So now into how we’re actu­al­ly going to pre­dict what the last term of each sequence is. That’ll be done using the fol­low­ing 3 meth­ods in this order (i.e. if method 1 fails, move onto method 2):

1.  Lin­ear Fit
2.  Fre­quen­cy Table
3.  Mode

Let’s start with the Mode method­ol­o­gy since it’s the sim­plest. For this, you sim­ply find the mode in a giv­en sequence, and that’ll be your guess for the last term in the sequence. The Mode Bench­mark seen on the Leader­board has an accu­ra­cy of 5.7%, so it’s a decent fall-back if our first two options fail.

Mode <- function(x) {
  ux <- unique(x)
  ux[which.max(tabulate(match(x, ux)))]
}

A lot of these sequences are sub-sequences of oth­er sequences with­in our dataset. So we’ll be look­ing at the last term in each sequence (i.e. say a sequence ends 185), find oth­er sequences that con­tain that num­ber (i.e. say 10 sequences con­tain 185), and then find what fol­lows that num­ber (i.e. if 287 fol­lows 185 in 7/10 sequences, that’ll be our guess). This idea was inspired by this script on Kaggle

buildFrequencyTable <- function(sequences)
{
  # Collate all pairs of consecutive integers across all sequences
  pairs <- ldply(sequences,
                 function(sequence)
                 {                 
                   if(length(sequence) <= 2)
                   {
                     return(NULL)
                   }
                   data.frame(last=head(sequence,-1),
                              following=tail(sequence,-1))
                 })
  
  # For each unique integer across all sequences, find the most common
  # integer that comes next
  frequencyTable <- ddply(pairs,
                          "last",
                          function(x)
                          {
                            data.frame(last=x$last[1],
                                       following=Mode(x$following))
                          })
  
  return(frequencyTable)
}

That’s the method itself. It’ll gen­er­ate a Fre­quen­cy Table that lists the inte­ger most com­mon­ly fol­low­ing every unique inte­ger. To inte­grate into code simply:

test<-read.csv("test.csv")
test$Sequence<- sapply(test$Sequence, FUN = parseSequence)  
frequencyTable <- buildFrequencyTable(test$Sequence)
test$last <- sapply(test.sample$Sequence, 
                           function(sequence) tail(unlist(sequence),1))
test <- join(test,frequencyTable,type="left")

This gives a a data frame, test, with 4 vari­ables — Id (ID vari­able), Sequence (the sequence), last (the last term in the sequence) & fol­low­ing (the most com­mon inte­ger fol­low­ing the last term in the sequence).  We’ll be using this in our Lin­ear Fit Prediction.

The Lin­ear Fit pre­dic­tion code is below

predictor<-function(seq.data,freq.soln){
  seq<-unlist(seq.data)
  if(length(seq)<2){
    return(tail(seq,1))
  }

  for(numberOfPoints in 1:(length(seq)-1)){
    df <- data.frame(y=tail(seq,-numberOfPoints))
    
    formulaString <- "y~"
    
    for(i in 1:(numberOfPoints))
    {
      df[[paste0("x",i)]] <- seq[i:(length(seq)-numberOfPoints+i-1)]
      formulaString <- paste0(formulaString,"+x",i)
    }
    
    fit <- lm(formula(formulaString),df)
    mae <- max(abs(fit$residuals))
    
    df <- list()
    for(i in 1:numberOfPoints)
    {
      df[[paste0("x",i)]] <- seq[length(seq)-numberOfPoints+i]
    }
    prediction<-predict(fit,df) if(mae>0 && mae<1){
      return(round(prediction))
    }
  }
  if(!is.na(freq.soln)){
    return(freq.soln)
  }


    return(Mode(seq))
  
}

I’ll break down the indi­vid­ual com­po­nents of that code block.

Here all we’re doing is unlist­ing the sequence that came in so we can with it. The first check is to see if the sequence is just one inte­ger — if so, we obvi­ous­ly can’t do any­thing with it and return the first ele­ment. I sup­pose this could be inte­grat­ed into the Mode fall­back, but if in the future I want to track how many short sequences there are, it helps.

  seq<-unlist(seq.data)
  iterations<<-iterations+1
  if(length(seq)<2){
    print("Returned Short Seq")
    return(tail(seq,1))
  }
  }

The lin­ear mod­el works back­wards from a sequence. I got the inspi­ra­tion (and most of the code) for this sec­tion from this script. The last X terms (defined by num­berOf­Points) from the end of the sequence are used to mod­el the sequence. The first for loop with­in the over­ar­ch­ing for loop is used to cre­ate the for­mu­la string (i.e. y ~ x1 + x2 + x3…) that’ll be plugged into the lin­ear mod­el func­tion, lm(). We also cal­cu­lat­ed what the max­i­mum resid­ual is.

Now that we have our mod­el, we cre­ate a data frame of the val­ues we need to use the pre­dict() func­tion and gen­er­ate our pre­dic­tion. To save some com­pu­ta­tion time, this code block cna be placed in the fol­low­ing if statement.

Recall the max­i­mum resid­ual we cal­cu­lat­ed ear­li­er. If that val­ue is 0, we like­ly have 0 degrees of free­dom (more pre­dic­tor vari­ables than sam­ples). This Stack Over­flow post goes into it nice­ly. I used a max­i­mum resid­ual of less than 1 as a cut-off point to deter­mine if a mod­el was good enough. That can (& should be tuned) for max­i­mum accu­ra­cy. If the lin­ear mod­el we’ve made for a giv­en sequence meets that cri­te­ria we return the pre­dic­tion from the lin­ear fit.

  for(numberOfPoints in 1:(length(seq)-1)){
    df <- data.frame(y=tail(seq,-numberOfPoints))
    
    formulaString <- "y~"
    
    for(i in 1:(numberOfPoints))
    {
      df[[paste0("x",i)]] <- seq[i:(length(seq)-numberOfPoints+i-1)]
      formulaString <- paste0(formulaString,"+x",i)
    }
    
    fit <- lm(formula(formulaString),df)
    mae <- max(abs(fit$residuals))
    
    df <- list()
    for(i in 1:numberOfPoints)
    {
      df[[paste0("x",i)]] <- seq[length(seq)-numberOfPoints+i]
    }
    prediction<-predict(fit,df) if(mae>0 && mae<1){
      return(round(prediction))
    }
  }

If our lin­ear mod­el does not meet the cri­te­ria, we move on to our Fre­quent Solu­tion method, and if that fails, Mode is the fail-safe.

 if(!is.na(freq.soln)){
    return(freq.soln)
  }

    return(Mode(seq))

To get the actu­al solu­tion file out it’s simply:

results<-mapply(test.sample$Sequence, test.sample$following, FUN=predictor)
solution.df<-data.frame(test.sample$Id,results)
colnames(solution.df)<-c('Id','Last')
write.csv(solution.df,file="Soln.csv",row.names = F)

This yields a solu­tion that’s 18.63% accu­rate, good for 2nd on the leader­board as of June 20th. I’d like to acknowl­edge the fol­low­ing scripts:

• Explor­ing Sequences
• Lin­ear Models
• Fre­quen­cies
• Mode

Those scripts did the brunt of the work, and taught me a lot. I high­ly rec­om­mend them.

The Method­ol­o­gy

Obvi­ous­ly there’s a lot that can still be added to this script. Cur­rent­ly ~27% of pre­dic­tions use the lin­ear mod­el, ~38% use fre­quen­cy pre­dic­tion and ~34% use the mode pre­dic­tion. Look­ing through some of the sequences by hand, this mod­el miss­es on a lot of geo­met­ric sequences and sequences with a mod­u­lus term. I’m not entire­ly sure to attack sequences with a mod­u­lar term, but I’d imag­ine sequences with a geo­met­ric term can be mod­eled. The imme­di­ate next step is to add an addi­tion­al mod­el between the lin­ear mod­el and fre­quen­cy mod­el step. hope­ful­ly low­er­ing the amount of sequences that require the Mode fallback.