Daddy U: Sabermetrics

Okay, let's see if we can use Julia for some sabermetrics.

First, head on over to juliabox.org and register. I have a secret code.

Load some packages

We need to add on some external packages. Sadly these can be a bit slow to load, so let's do it know:

using Gadfly, DataFrames

(These are pre-installed on juliabox.)

Download some data

f = "http://seanlahman.com/files/database/lahman-csv_2015-01-24.zip"
d = "/tmp/lahman"
if !isdir(d)
  mkdir(d)
  fout = "$d/lahman-csv.zip"
  download(f, fout)
  run(`unzip -n -d $d $fout`)
end

Read in some data

Julia had a DataFrame type to hold tabular data. (Akin to Python's Pandas or R's data.frame types.) We can use readtable to read a csv file into a data frame

d = readtable("$d/Batting.csv")

Data frames store variables in columns. To access a single column we use the notation d[:name]. The : makes the name a symbol. This works fine unless there are spaces in the names, in which case use :(name with spaces) as a pattern.

d[:RBI]

99846-element DataArray{Int64,1}:
  0
 13
 19
 27
 16
  5
  2
 34
  1
 11
 18
  0
  1
  ⋮
 13
  0
 28
 10
 30
 17
  0
  0
  1
 38
 52
 60

Subsetting

There are values with no RBI info. These are stored as NA (not 0). To see them we need to

• subset the data to show where RBIs are NA

Indexing of d can be done by column (as above) or by row and column. The notation is d[rowindex,columnindex].

The row index can be specified by number (1 or by range 1:3 or by vector [1,3,5], use : to refer to all of them). Similarly the column index.

More commonly we can index by a logical vector with a length that matches the number of rows. These are generated by logical operators or by special "is" functions, as here with isna:

d[isna(d[:RBI]), :yearID]

5573-element DataArray{Int64,1}:
 1882
 1882
 1882
 1882
 1882
 1882
 1882
 1882
 1882
 1882
 1882
 1882
 1882
    ⋮
 1999
 1999
 1999
 1999
 1999
 1999
 1999
 1999
 1999
 1999
 1999
 1999

We see the early years, but also some in 1999. What happened there?

d[(isna(d[:RBI])) & (d[:yearID] .== 1999), :]

I'm guessing American League pitchers with no atbats.

That syntax is pretty clunky. You need to put the expressions in (), and use .== instead of just == to get comparison as we are comparing each element of d[:yearID] to 1999. The DataFramesMeta package can make it easier. We have to install this one thoughL

Pkg.add("DataFramesMeta")
using DataFramesMeta

Then the above becomes something like:

@linq d |>
where(isna(:RBI)) |>
where(:yearID .== 1999)

Splitting data

Okay, let's get at some anomalies. Which players have the most home runs?

The data is stored by player so we need to split up the data by each player and add up the :HR values. When we split the data, we will get a smaller data frame. This function will add up the HRs ( paying attention to NA values):

function add_hrs(d)
  hrs = d[:HR]
  hrs = hrs[!isna(hrs)]
  sum(hrs)
end

add_hrs (generic function with 1 method)

Here we split the data and apply the function

a = by(d, :playerID, add_hrs)

We want to sort by the numbers. We can sort as follows:

sort(a[:x1])

18405-element DataArray{Int64,1}:
   0
   0
   0
   0
   0
   0
   0
   0
   0
   0
   0
   0
   0
   ⋮
 569
 573
 583
 586
 609
 612
 630
 654
 660
 714
 755
 762

But to sort the data frame we want to shuffle the items around. The sortperm function helps here. It returns the indices of the permutation that will sort the data. Here we pass in rev=true to get reverse order:

a[sortperm(a[:x1], rev=true),:]

Go Bonds!

On base percentage

What function will compute the on base percentage?

$$ OBP = \frac{H + BB + HBP}{AB + BB + HBP + SF} $$

Do we have these in the data?

names(d)

22-element Array{Symbol,1}:
 :playerID
 :yearID  
 :stint   
 :teamID  
 :lgID    
 :G       
 :AB      
 :R       
 :H       
 :x2B     
 :x3B     
 :HR      
 :RBI     
 :SB      
 :CS      
 :BB      
 :SO      
 :IBB     
 :HBP     
 :SH      
 :SF      
 :GIDP    

We are good to go. Here we define a function to compute the OBP. The only issue is the NA values must be addressed and those who barely played are weeded out (and given a 0.0 value)

function compute_obp(d)
u = d[:, [:H, :BB, :HBP, :AB, :SF]]
## make NA -> 0
for nm in names(u)
    u[isna(u[nm]), nm] = 0
end
sum(u[:AB]) > 100 || return(0.0)	
top =  @with(u, :H + :BB + :HBP)
bottom = @with(u, :AB + :BB + :HBP + :SF)
sum(top) / sum(bottom)
end

compute_obp (generic function with 1 method)

The pattern is as before:

OBP = by(d, :playerID, compute_obp)
OBP = OBP[sortperm(OBP[:x1], rev=true), :]