R Linear Regression: Linear Models (Advanced Tibbles Data Structure)

Advanced dplyr: Tibbles

• Tibbles can be regarded as a modern version of data frames and are the default data structure in the tidyverse.
• Some functions that do not work properly with data frames do work with tibbles.

Code: Tibbles vs DF

# stratify by HR
dat <- Teams %>% filter(yearID %in% 1961:2001) %>%
  mutate(HR = round(HR/G, 1), 
         BB = BB/G,
         R = R/G) %>%
  select(HR, BB, R) %>%
  filter(HR >= 0.4 & HR<=1.2)
  
# calculate slope of regression lines to predict runs by BB in different HR strata
dat %>%  
  group_by(HR) %>%
  summarize(slope = cor(BB,R)*sd(R)/sd(BB))
  
# use lm to get estimated slopes - lm does not work with grouped tibbles
dat %>%  
  group_by(HR) %>%
  lm(R ~ BB, data = .) %>%
  .$coef
  
# inspect a grouped tibble
dat %>% group_by(HR) %>% head()
dat %>% group_by(HR) %>% class()

Tibbles vs DataFrames: The Differences

• Tibbles are more readable than data frames.
• If you subset a data frame, you may not get a data frame. If you subset a tibble, you always get a tibble.
• Tibbles can hold more complex objects such as lists or functions.
• Tibbles can be grouped.

Code: Tibbles

# inspect data frame and tibble
Teams
as.tibble(Teams)

# subsetting a data frame sometimes generates vectors
class(Teams[,20])

# subsetting a tibble always generates tibbles
class(as.tibble(Teams[,20]))

# pulling a vector out of a tibble
class(as.tibble(Teams)$HR)

# access a non-existing column in a data frame or a tibble
Teams$hr
as.tibble(Teams)$hr

# create a tibble with complex objects
tibble(id = c(1, 2, 3), func = c(mean, median, sd))

do Function

• The do() function serves as a bridge between R functions, such as lm(), and the tidyverse.
• We have to specify a column when using the do() function, otherwise we will get an error.
• If the data frame being returned has more than one row, the rows will be concatenated appropriately.

Code: do Function

# use do to fit a regression line to each HR stratum
dat %>%  
    group_by(HR) %>%
    do(fit = lm(R ~ BB, data = .))
    
# using do without a column name gives an error
dat %>%
    group_by(HR) %>%
    do(lm(R ~ BB, data = .))

# define a function to extract slope from lm
get_slope <- function(data){
  fit <- lm(R ~ BB, data = data)
  data.frame(slope = fit$coefficients[2], 
             se = summary(fit)$coefficient[2,2])
}

# return the desired data frame
dat %>%  
  group_by(HR) %>%
  do(get_slope(.))
  
# not the desired output: a column containing data frames
dat %>%  
  group_by(HR) %>%
  do(slope = get_slope(.))
  
# data frames with multiple rows will be concatenated appropriately
get_lse <- function(data){
  fit <- lm(R ~ BB, data = data)
  data.frame(term = names(fit$coefficients),
    slope = fit$coefficients, 
    se = summary(fit)$coefficient[,2])
}

dat %>%  
  group_by(HR) %>%
  do(get_lse(.))

broom Function

• The broom package has three main functions, all of which extract information from the object returned by lm and return it in a tidyverse friendly data frame.
• The tidy() function returns estimates and related information as a data frame.
• The functions glance() and augment() relate to model specific and observation specific outcomes respectively.

Code: broom Function

# use tidy to return lm estimates and related information as a data frame
library(broom)
fit <- lm(R ~ BB, data = dat)
tidy(fit)

# add confidence intervals with tidy
tidy(fit, conf.int = TRUE)

# pipeline with lm, do, tidy
dat %>%  
  group_by(HR) %>%
  do(tidy(lm(R ~ BB, data = .), conf.int = TRUE)) %>%
  filter(term == "BB") %>%
  select(HR, estimate, conf.low, conf.high)
  
# make ggplots
dat %>%  
  group_by(HR) %>%
  do(tidy(lm(R ~ BB, data = .), conf.int = TRUE)) %>%
  filter(term == "BB") %>%
  select(HR, estimate, conf.low, conf.high) %>%
  ggplot(aes(HR, y = estimate, ymin = conf.low, ymax = conf.high)) +
  geom_errorbar() +
  geom_point()
  
# inspect with glance
glance(fit)