Learning about Lengths of Baseball Games

One complaint about baseball is that it is a relatively long game, compared to American football or basketball.

To understand more about the lengths of baseball games, I collected the lengths (in minutes) of all games played during the 2007 season.  Here are my questions:

1.  Are the times of baseball games normally distributed?

2.  If not, how do the times differ from a normal curve?

I'll illustrate the R work. You'll need a couple of packages: LearnEDA and vcd (for the rootogram function).

1.  First I'll read in the baseball game times, construct bins and the rootogram (this is essentially a histogram where you plot the roots of the counts instead of the counts).

library(LearnEDA)

data=read.table("game.time.07.txt",header=T)

attach(data)

# set up bins and construct a rootogram

bins=seq(100,300,by=20)

bin.mids=(bins[-1]+bins[-length(bins)])/2

h=hist(time,breaks=bins)

h$counts=sqrt(h$counts)

plot(h,xlab="TIME",ylab="ROOT FREQUENCY",main="")

2.  Now I want to fit a normal comparison curve.  I figure out the mean and standard deviation of the normal curve by using letter values.

# find mean and standard deviation of matching Gaussian curve

S=lval(time)

f=as.vector(S[2,2:3])

m=as.integer((f[1]+f[2])/2)

sd=as.integer((f[2]-f[1])/1.349)

I get that the times are approximately N(175, 24).

3.  Is this normal approximation ok?  The function fit.gaussian computes the expected counts for each bins and computes the residuals sqrt(observed) - sqrt(expected).

# function fit.gaussian.R fits Gaussian curve to the counts

s=fit.gaussian(time,bins,m,sd)

# output observed and expected counts

output=round(data.frame(s$counts, sqrt(s$counts), s$probs, s$expected, 

     sqrt(s$expected), sqrt(s$counts)-sqrt(s$expected)),2)

names(output)=c("count","root","prob","fit","root fit","residual")

output

   count root prob   fit root fit residual

1      1 1.00 0.01  1.00     1.00     0.00

2      8 2.83 0.06  6.08     2.47     0.36

3     22 4.69 0.19 19.17     4.38     0.31

4     27 5.20 0.32 31.34     5.60    -0.40

5     20 4.47 0.27 26.60     5.16    -0.69

6     12 3.46 0.12 11.72     3.42     0.04

7      6 2.45 0.03  2.67     1.64     0.81

8      1 1.00 0.00  0.32     0.56     0.44

9      1 1.00 0.00  0.02     0.14     0.86

10     1 1.00 0.00  0.00     0.02     0.98

# place root expected counts on top of rootogram

lines(bin.mids,sqrt(s$expected),lwd=4,col="red")

3.  The rootogram function plots the residuals.  Both plots below show the same info, but in different ways.

# load library vcd that contains rootogram function

library(vcd)

# again save histogram object

h=hist(time,breaks=bins,plot=F)

# illustrate "hanging" style of rootogram

r=rootogram(h$counts,s$expected,type="hanging")  # this is the hanging style 

# illustrate "deviation" style of rootogram

r=rootogram(h$counts,s$expected,type="deviation") 

What have we learned?  Although the lengths of baseball games are roughly normal in shape, we see that there are some large negative residuals on the right tail.

What does this mean?  Baseball game lengths have a long right-tail which means one tends to see many long games.  If you know anything about baseball, you know that baseball can go into extra-innings when there the game is tied after 9 innings, and these extra-inning games cause the long right tail.

Median Polishing Student Grades

Our department recently was interested in exploring the relationship between student class attendance and performance in a number of 100-level math classes.   We all know that missing math classes will have an adverse effect on grades, but we wanted to learn more about this relationship.

For many students in four 100-level math classes, we collected the student's attendance (percentage of classes attended) and their performance (a percentage where a 90% corresponds to A work, 80-90% a B, etc).   We decided to categorize attendance using the categories 0-50%, 50-70%, 70-90%, over 90%) and then we found the mean performance for students in each categorization of attendance and course.  We got the following two-way table.

          Attendance Pct

COURSE   0-50 50-70 70-90 90-100  

---------------------------------       

MATH 112  61.9  68.0  68.8  75.5

MATH 115  54.3  71.8  76.6  83.3

MATH 122  56.1  64.8  73.1  78.1

MATH 126  50.6  66.3  73.9  78.1

I'll demonstrate median polish with this data.

1.  First, I have to get this data into a matrix form into R.  I'll first use the matrix command to form the matrix (by default, one enters data column by column) and then I'll add row and column labels.

grade=matrix(c(61.9,54.3,56.1,50.6,

               68.0,71.8,64.8,66.3,

               68.8,76.6,73.1,73.9,

               75.5,83.3,78.1,78.1),c(4,4))

dimnames(grade)=list(c("MATH 112","MATH 115","MATH 122","MATH 126"),

                     c("0-50","50-70","70-90","90-100"))

To check that I've read in the data correctly, I'll display "grade":

> grade

         0-50 50-70 70-90 90-100

MATH 112 61.9  68.0  68.8   75.5

MATH 115 54.3  71.8  76.6   83.3

MATH 122 56.1  64.8  73.1   78.1

MATH 126 50.6  66.3  73.9   78.1

2.  Now I can implement median polish to get an additive fit.  I'll store the output in the variable "my.fit" and then I'll display the different components.

> my.fit=medpolish(grade)

Here's the common value.

> my.fit$overall

[1] 69.7

Here are the row effects.

> my.fit$row

MATH 112 MATH 115 MATH 122 MATH 126 

 -0.7375   4.5000   0.2875  -0.2875 

Here are the column effects.

> my.fit$col

     0-50     50-70     70-90    90-100 

-16.35000  -2.75625   2.75625   8.40000 

I'll interpret this additive fit.

• The average performance of these students is 69.7%.
  
  
• Looking at the row effects, we see that MATH 115 students get grades that are 4.5 - (0.2875) approx 4.2 points higher than MATH 122 students.  MATH 122 students tend to be 0.2875 - (-0.2875) approx 0.57 points higher than MATH 126 students, and MATH 112 students are about a half percentage point lower than MATH 126 students.
  
  
• Looking at the column effects, we see a clear relationship between attendance and performance.  The best (90-100%) attenders do 8.4 - 2.75 = 5.65 points better (on average) than the 70-90 attendance group, do 8.4 - (-2.75) = 11.15 points better than the 50-70 attendance group, and a whopping 8.4 - (-16.35) = 24.75 points better than the "no shows" (the under 50% attending group).

3.  It might be helpful to plot this fit.  I wrote a function plot2way in the LearnEDA package:

> library(LearnEDA)

> plot2way(my.fit$overall+my.fit$row,my.fit$col,dimnames(grade)[[1]],

  dimnames(grade)[[2]])

Note that the plot2way function has four arguments:  the row part, the column part, the vector of names of the rows, and the vector of names of the columns.  Here's the figure.

Actually, this plot would look better if I had figured out how to rotate the figure so that the FIT lines are horizontal.  (I'll give extra credit points to anyone who can fix my function to do that.)

From this graph we see that the best performances are the MATH 115 students who attend over 90% of the classes; the worst performances are the MATH 112 students who have under 50% attendance.

4.  Are we done?  Not quite.  We have looked at the fit, but have not looked at the residuals -- the differences between the performance and the fit.  They are stored in my.fit$residual -- I will round the values so they are easier to view.

> round(my.fit$residual)

         0-50 50-70 70-90 90-100

MATH 112    9     2    -3     -2

MATH 115   -4     0     0      1

MATH 122    2    -2     0      0

MATH 126   -2     0     2      0

I see one large residual that I have highlighted in red.  MATH 112 students who don't come to class (under 50% attendance) seem to do 9 points better than one would expect based on the additive fit.   This might deserve further study.

My Smoothing Tribute to the Phillies

Last night was an exciting night.  My team, the Philadelphia Phillies, won the World Series for only the second time in their history that started in 1883.

So it seems appropriate to give a statistical tribute to the Phillies.

I collected the Phillies winning percentage (percentage of games won) for each season in their history.  Below I plot the winning percentage against the season year.

It is hard to see the general pattern in this graph, so it makes sense to smooth.  I could use Tukey's resistant smooth described in the EDA notes, but I'll illustrate the use of an alternative smoothing method called lowess.  Here is the R code to graph the scatterplot and overlay the lowess smooth.

plot(Year,Win.Pct,pch=19,col="red")

lines(lowess(Year,Win.Pct,f=.2),lwd=3,col="red")

title("PHILLIES WINNING PERCENTAGES",col="red")

abline(h=50,lwd=2)

I added the horizontal line at WIN.PCT = 50 that corresponds to an average season.

Looking at the graph, you should see why the Philadelphia fans are so excited about the Phillies winning this season.

1.  The Phillies teams generally have been crummy, especially between 1920 and 1950.  The team hasn't experienced much success.

2.  But there have been two recent periods where the Phillies have been successful.  One period was in the 1970's and the climax of this success was the Phillies first World Series win in 1980.  The second period is in recent history and of course the Phillies won their second World Series in 2008.

Comments on Plotting Homework

Here are some general and specific comments on your efforts on Homework 4 on plotting and straightening.

First, most of you did well on this homework.  You were good in describing the fit and the pattern in the residuals.  Also you made reasonable choices at reexpressing the x and/or the y variable so that the graph looked pretty straight.

But there were some things that caused you to lose points.

1.  What are you looking for in a residual plot?  In this assignment, the focus was looking for nonlinear patterns.    For example, is it appropriate to fit a line to the (year, population) data for the England and Wales dataset?  We can answer this question by looking for a nonlinear pattern in the plot of residuals against year.  There is significant curvature (that is, a quadratic pattern) in the residual plot which tells you that the population growth is not linear.

By the way, some of you plotted log population against log year -- why did you take the log of year?  It doesn't make any sense to me.  

2.  Always talk about the fit and the residuals in the context of the data.  Someone plotted x against y without telling me the variables.  The fun part of statistics is that you can always talk about the application.

3.  Remember that funny problem where the scatterplot shows two group of points with a clear separation?  This is one type of nonlinear pattern that you won't be able to straighten with a single choice of power transformation.  But since the graph clearly divides into two parts, it makes sense to treat this as two independent problems and try to straighten each part.

4.  Should one fit a line by least-squares or by a resistant line?  In many situations, it won't make a difference -- either fit will work.  But least-squares can give you relatively poor fits when there are outliers.  

How can you tell if least squares isn't the best fit?  Look at the residual plot.  If you still see some increasing or descreasing pattern, then this tells you that least-squares hasn't explained all of the "tilt" pattern in the graph.

Straightening

This week, the main topic is reexpressing either the x or y variable to straighten a nonlinear pattern that we see in a scatterplot. Although the manipulation may be straightforward, it is possible to miss the main message in this material. Here are some questions that may help in your explanations in the homework.

Question 1: What is a simple description of the pattern in a scatterplot?

Question 2: Why do we prefer to fit lines instead of more complicated curves like quadratic or cubic?

Question 3: Is it possible to straighten all nonlinear patterns by power transformations?

Question 4: Can you think of a situation or an example where it is not possible to straighten by power transformations?

Don't forget to look at your data first -- you may see rightaway that it is not possible to straighten the graph.

Example of Resistant Fitting

In baseball, the objective is to win games and a team wins a game by scoring more runs than its opponent. An interesting question is "how important is a single run" towards the goal of winning a game? Suppose one collects the runs scored, the runs allowed, the number of wins and the number of losses for a group of teams. Bill James (a famous guy who works on baseball data) discovered the empirical relationship

log(wins/losses) = 2 log(runs scored/runs allowed)

He called this the Pythagorean Relationship.

Let's try to demonstrate this relationship by use of a resistant fit.

1. First, I collected data for all baseball teams in the 2008 season. The dataset teams2008. txt contains for each of the 30 teams ...

Team -- the name of the team
Wins -- the number of wins
Losses -- the number of losses
Runs.Scored -- the total number of runs scored
Runs.Allowed -- the total number of runs allowed

2. I read this dataset into R and compute the variables log.RR and log.WL.

data=read.table("http://bayes.bgsu.edu/eda/data/teams2008.txt",header=T)
attach(data)

log.RR=log(Runs.Scored/Runs.Allowed)
log.WL=log(Wins/Losses)

3. I graph log.RR against log.WL and add team labels to the graph. As we hoped, the relationship looks pretty linear.

plot(log.RR,log.WL,pch=19)
text(log.RR,log.WL,Team,pos=2)



















4. I next fit a resistant line using the rline function in the LearnEDA package. I add the
fitted line to the graph.

the.fit=rline(log.RR,log.WL,iter=4)
curve(the.fit$a+the.fit$b*(x-the.fit$xC),add=TRUE)

5. If Bill James' relationship holds, the slope of the resistant line should be close to 2.

the.fit
$a
[1] 0.01006079
$b
[1] 1.801718

It approximately holds since the slope of 1.8 is close to 2.

6. To see if this is a reasonable fit, we compute the fit and the residuals.

FIT=the.fit$a+the.fit$b*(log.RR-the.fit$xC)
RESIDUAL=log.WL-FIT

and then plot the residuals, adding the team labels.

plot(log.RR,RESIDUAL,pch=19)
abline(h=0)
text(log.RR,RESIDUAL,Team,pos=2)



















7. What are we looking for in the residual plot? First, we look for general patterns that we didn't see earlier in the first plot. I don't see any trend, so it appears that we removed the tilt by fitting a line.

Also we are looking for unusually small or large residual. Here a "lucky team" corresponds to a team who seemed to win more games than one would expect based on their wins and losses.

Which team was unusually lucky in the 2008 season? A hint: they were a "heavenly" team from the west coast.

Reexpressing for Symmetry

Your instructor is currently in "Phillies Heaven".  His baseball team rarely has a chance to win the World Series (they have won one World Series in over 120 seasons of competing) and they currently in the National League Championship against the Dodgers!

Oh right -- I'm supposed to talk about the EDA class.

I finished the grading on your "symmetry" homework and Fathom assignment.    You generally did fine, but I'll explain some issues that may have caused you to lose points.

What was I looking for?

In the notes, we talked about several methods for assessing symmetry of a batch, including looking a the sequence of midsummaries, using a symmetry plot, and Hinkley's quick method.   For each dataset you consider, here's a outline of what you should do:

1.  Demonstrate (using some method) that your raw data looks nonsymmetric.

2.  Experiment with power transformations with different choices of p to try to make the reexpressed data approximately symmetric.  Use one of our methods to see if the "p-reexpression" works.

3.  Convince me that you have found a reasonable transformation by graphing the reexpressed data (say by a histogram or a stemplot).

It is best if you explain (with words) your process of finding the best reexpression.  I'm much more interested in your thought process than your computer output.

Here are a few other pitfalls.

1.  Some of you were confused when you considered reexpressions with negative values of p.

p = 1  -- data is really right-skewed (big positive value of Hinkley's d)

p = 0 -- data is slightly right-skewed (smaller positive value of d)

p = -1 -- data is right skewed (positive value of d)

What is going on?  The problem is that you were defining a reexpression like p = -1/2 as

data^(-1/2)

when you should have used

-data^(-1/2)

By adding the negative sign, all of your transformations are increasing functions, and then

you can make better sense of the reexpressed graphs and the methods.

2.  Reexpression only works when there is sufficient spread in your data.  Suppose you have data that ranges from 50 to 60 -- here the range is only 60/50 = 1.2.  Reexpression using any value of p won't work -- that is, it won't change the shape of the data.

3.  One of you used a normal probability plot to determine if the graph was symmetric.  What's wrong with this?  Well, it is not one of the methods we discussed in the notes.  Second, checking for normality is different (but related) than checking for symmetry.  Data can be symmetric but not normally distributed.