Group Projects

Your final is a group project. You will also have two “mini” projects. They comprise a large part of your grade.

You need to start planning soon.

To aid in your planning, here are the required elements of your final project.

1. You must find existing data to analyze. Aggregating and merging data from multiple sources is encouraged.

2. You must visualize 3 intersting features of that data.

3. You must come up with some analysis—using tools from this course—which relates your data to either a prediction or a policy conclusion.

4. You must think critically about your analysis and be able to identify potential issues.

5. You must present your analysis as if presenting to a C-suite executive.

Your mini-projects along the way will be more structured, but will serve to guide you towards the final project.

Teams

Please form teams of 2-3 people. Once all agree to be on a team, have one person email our TA, Nick, armst427@msu.edu and cc all of the members of the team so that nobody is surprised to be included on a team. Title the email [SSC442] - Group Formation. Tell us your team name (be creative), and list in the email the names of all of the team members and their email address (in addition to cc-ing those team members on the email).

If you opt to not form a team, you will be automatically added to the “willing to be randomly assigned” pool and will be paired with another 1-2 persons.

Send this email by September 21st and we will assign un-teamed folks at the beginning of the following week. Project 1 is due in mid-October. See schedule for all the important project dates.

Guiding Question

For future lectures, the guiding questions will be more pointed and at a higher level to help steer your thinking. Here, we want to ensure you remember some basics and accordingly the questions are straightforward.

• Why do we want tidy data?
• What are the challenges associated with shaping things into a tidy format?

Learning From Data

The following are the baisc requirements for statistical learning

1. A pattern exists

2. This pattern is not easily expressed in a closed mathematical form

3. You have data

Formalization

We think of our outcome-of-interest as a reponse or target that we wish to predict or wish to learn something about.

We generically refer to the response as \(Y\)

Other aspects of the data are known as features, inputs, predictors, or regressors. We call one of these \(X_i\).

• The subscript \(i\) indicates that we have an \(X\) realized for every individual observation in our data

We can refer to the input vector collectively as:

\[X = \begin{bmatrix} x_{11} & x_{12} \\ x_{21} & x_{22} \\ \vdots & \vdots \\ x_{N1} & x_{N2} \end{bmatrix}\]

This is sometimes called a “rectangular array” of data. Sounds complicated, but you’ve definitely encountered it before. Excel holds “rectangular arrays” of data.

We are seeking some unknown function that maps \(X\) to \(Y\)

Put another way, we are seeking to explain \(Y\) as follows:

\[Y = f(X) + e\]

The Target Function

We call the function \(f: \mathcal{X} \rightarrow \mathcal{Y}\) the target function

How do we find the function? We don’t! We get as close as we can, though:

• Observe data \((\mathbf{x}_1, y_1), \cdots, (\mathbf{x}_N, y_N)\)

• Use some algorithm to approximate \(f\)

• Produce final hypothesis function \(g \approx f\)

• Evaluate how well \(g\) approximates \(f\) and iterate as needed.

Why Estimate an Unknown Function?

With a good estimate of \(f\) we can make predictions of \(Y\) at new points \(X = x\)

We can also understand which components of \(X = (X_1, X_2, \cdots, X_m)\) are important in explaining \(Y\), and which are (potentially) irrelevant

• e.g. GDP and yearsindustrialized have a big impact on emissions but hydroutilization typically does not.

Depending on the complexity of \(f\), we may be able to meaningfully understand how each component of \(X\) affects \(Y\).

(But we should be careful about assigning causal interpretations, more on this later)

The Parable of the Marbles

(Courtesy of Prof. Bushong)

Imagine a bag of marbles with two types of marbles: ♣️ and ♦️.

We are going to pick a sample of \(n\) marbles (with replacement).

We want to learn something about \(\mu\), the objective probability to pick a ♣️.

In addition to defining the objective probability of picking a ♣️, we have an observed fraction \(\eta\), which will define as the fraction of ♣️ in the sample.

Question: Can we say anything exact and for-sure about \(\mu\) (outside the data) after observing \(\eta\) (the data)?

• No. It is possible for the sample to be all ♣️, ♣️, ♣️, ♣️, ♣️ even when the bag is is 50/50 ♦️

  • No matter what we draw, we can’t (based on that draw alone) eliminate the possibility of drawing a ♦️.

  • And unless we assume that the only two values in the world are ♦️ and ♣️, we can’t rule out 💩¹

Question: If we really can’t say anything exact and for-sure, then why do we do things like polling (e.g. to predict the outcome of a presidential election)?

• The bad case, that we draw something that has is completely misleading, is possible but not probable.

Outside the Data

Put another way, since \(f\) is unknown, it can take on any value outside the data we have, no matter how large the data.

• This is called No Free Lunch

You cannot know anything for sure about \(f\) outside the data without making assumptions.

Is there any hope to know anything about \(f\) outside the data set without making assumptions about \(f\)?

Yes, if we are willing to give up the “for sure”.

Instead of

\[ Y = f(X) \]

we can work with

\[ Y = f(X) + e \]

and if we guess at an \(f\) function that gets pretty close, \(e\) will be pretty small. If \(e\) is pretty small, then we’re in good shape for prediction!

Hoeffding’s Inequality

Hoeffding’s Inequality states, loosely, that \(\eta\) (the fraction of our sample that is ♣️) cannot be too far from \(\mu\).

\[\mathbb{P}\left[|\eta - \mu| > \epsilon \right] \leq 2e^{-2\epsilon^2n}\] \(\eta \approx \mu\) is called probably approximately correct (PAC) learning.

An example of Hoeffding’s Inequality

Example: n = 1,000. Draw a sample and observe \(\eta\)

• \(\sim\) 99% of the time, \(\mu - .05 \leq \eta \leq \mu+.05\)

  • This is implied by setting \(\epsilon = 0.05\) and using \(n=1,000\)

• 99.9999996% of the time \(\mu - .10 \leq \eta \leq \mu + .10\)

What does this mean?

If I repeatedly pick a sample of size 1,000, observe \(\eta\) and claim that \(\mu \in \left[\eta - .05, \eta + .05\right]\) (or that the error bar is \(\pm 0.05\)), I will be right 99% of the time.

On any particular sample you may be wrong, but not often.

Tidy data

We say that a data table is in tidy format if each row represents one observation and columns represent the different variables available for each of these observations. The murders dataset is an example of a tidy data frame.

library(dslabs)
data(murders)
head(murders)

##        state abb region population total
## 1    Alabama  AL  South    4779736   135
## 2     Alaska  AK   West     710231    19
## 3    Arizona  AZ   West    6392017   232
## 4   Arkansas  AR  South    2915918    93
## 5 California  CA   West   37253956  1257
## 6   Colorado  CO   West    5029196    65

Each row represent a state with each of the five columns providing a different variable related to these states: name, abbreviation, region, population, and total murders.

To see how the same information can be provided in different formats, consider the following example:

library(dslabs)
data("gapminder") # gapminder will now be a data.frame in your "environment" (memory)
tidy_data <- gapminder %>%
  filter(country %in% c("South Korea", "Germany") & !is.na(fertility)) %>%
  select(country, year, fertility)
head(tidy_data, 6)

##       country year fertility
## 1     Germany 1960      2.41
## 2 South Korea 1960      6.16
## 3     Germany 1961      2.44
## 4 South Korea 1961      5.99
## 5     Germany 1962      2.47
## 6 South Korea 1962      5.79

This tidy dataset provides fertility rates for two countries across the years. This is a tidy dataset because each row presents one observation with the three variables being country, year, and fertility rate. However, this dataset originally came in another format and was reshaped for the dslabs package. Originally, the data was in the following format:

##       country 1960 1961 1962
## 1     Germany 2.41 2.44 2.47
## 2 South Korea 6.16 5.99 5.79

The same information is provided, but there are two important differences in the format: 1) each row includes several observations and 2) one of the variables’ values, year, is stored in the header. For the tidyverse packages to be optimally used, data need to be reshaped into tidy format, which you will learn to do throughout this course. For starters, though, we will use example datasets that are already in tidy format.

Although not immediately obvious, as you go through the book you will start to appreciate the advantages of working in a framework in which functions use tidy formats for both inputs and outputs. You will see how this permits the data analyst to focus on more important aspects of the analysis rather than the format of the data.

Manipulating data frames

The dplyr package from the tidyverse introduces functions that perform some of the most common operations when working with data frames and uses names for these functions that are relatively easy to remember. For instance, to change the data table by adding a new column, we use mutate. To filter the data table to a subset of rows, we use filter. Finally, to subset the data by selecting specific columns, we use select.

library(dslabs)
data("murders")
murders <- mutate(murders, rate = total / population * 100000)

head(murders)

##        state abb region population total     rate
## 1    Alabama  AL  South    4779736   135 2.824424
## 2     Alaska  AK   West     710231    19 2.675186
## 3    Arizona  AZ   West    6392017   232 3.629527
## 4   Arkansas  AR  South    2915918    93 3.189390
## 5 California  CA   West   37253956  1257 3.374138
## 6   Colorado  CO   West    5029196    65 1.292453

filter(murders, rate <= 0.71)

##           state abb        region population total      rate
## 1        Hawaii  HI          West    1360301     7 0.5145920
## 2          Iowa  IA North Central    3046355    21 0.6893484
## 3 New Hampshire  NH     Northeast    1316470     5 0.3798036
## 4  North Dakota  ND North Central     672591     4 0.5947151
## 5       Vermont  VT     Northeast     625741     2 0.3196211

new_table <- select(murders, state, region, rate)
filter(new_table, rate <= 0.71)

##           state        region      rate
## 1        Hawaii          West 0.5145920
## 2          Iowa North Central 0.6893484
## 3 New Hampshire     Northeast 0.3798036
## 4  North Dakota North Central 0.5947151
## 5       Vermont     Northeast 0.3196211

library(dplyr)
library(dslabs)
data(murders)

murders <- mutate(murders, population_in_millions = population / 10^6)

select(murders, state, population) %>% head()

filter(murders, state == "New York")

no_florida <- filter(murders, state != "Florida")

filter(murders, state %in% c("New York", "Texas"))

filter(murders, population < 5000000 & region == "Northeast")

The pipe: %>%

With dplyr we can perform a series of operations, for example select and then filter, by sending the results of one function to another using what is called the pipe operator: %>%. Some details are included below.

We wrote code above to show three variables (state, region, rate) for states that have murder rates below 0.71. To do this, we defined the intermediate object new_table. In dplyr we can write code that looks more like a description of what we want to do without intermediate objects:

\[ \mbox{original data } \rightarrow \mbox{ select } \rightarrow \mbox{ filter } \]

For such an operation, we can use the pipe %>%. The code looks like this:

murders %>% select(state, region, rate) %>% filter(rate <= 0.71)

##           state        region      rate
## 1        Hawaii          West 0.5145920
## 2          Iowa North Central 0.6893484
## 3 New Hampshire     Northeast 0.3798036
## 4  North Dakota North Central 0.5947151
## 5       Vermont     Northeast 0.3196211

This line of code is equivalent to the two lines of code above. What is going on here?

In general, the pipe sends the result of the left side of the pipe to be the first argument of the function on the right side of the pipe. Here is a very simple example:

16 %>% sqrt()

## [1] 4

We can continue to pipe values along:

16 %>% sqrt() %>% log2()

## [1] 2

The above statement is equivalent to log2(sqrt(16)).

Remember that the pipe sends values to the first argument, so we can define other arguments as if the first argument is already defined:

16 %>% sqrt() %>% log(base = 2)

## [1] 2

Therefore, when using the pipe with data frames and dplyr, we no longer need to specify the required first argument since the dplyr functions we have described all take the data as the first argument. In the code we wrote:

murders %>% select(state, region, rate) %>% filter(rate <= 0.71)

murders is the first argument of the select function, and the new data frame (formerly new_table) is the first argument of the filter function.

Note that the pipe works well with functions where the first argument is the input data. Functions in tidyverse packages like dplyr have this format and can be used easily with the pipe. It’s worth noting that as of R 4.1, there is a base-R version of the pipe |>, though this has its disadvantages. We’ll stick with %>% for now.

murders <- mutate(murders, rate =  total / population * 100000,
                  rank = rank(-rate))

my_states <- filter(murders, region %in% c("Northeast", "West") &
                      rate < 1)

select(my_states, state, rate, rank)

mutate(murders, rate =  total / population * 100000,
       rank = rank(-rate)) %>%
  select(state, rate, rank)

my_states <- murders %>%
  mutate SOMETHING %>%
  filter SOMETHING %>%
  select SOMETHING

Summarizing data

An important part of exploratory data analysis is summarizing data. The average and standard deviation are two examples of widely used summary statistics. More informative summaries can often be achieved by first splitting data into groups. In this section, we cover two new dplyr verbs that make these computations easier: summarize and group_by. We learn to access resulting values using the pull function.

library(dplyr)
library(dslabs)
data(heights)
head(heights)

##      sex height
## 1   Male     75
## 2   Male     70
## 3   Male     68
## 4   Male     74
## 5   Male     61
## 6 Female     65

s <- heights %>%
  filter(sex == "Female") %>%
  summarize(average = mean(height), standard_deviation = sd(height))
s

##    average standard_deviation
## 1 64.93942           3.760656

s$average

## [1] 64.93942

s$standard_deviation

## [1] 3.760656

heights %>%
  filter(sex == "Female") %>%
  summarize(median = median(height), minimum = min(height),
            maximum = max(height))

##     median minimum maximum
## 1 64.98031      51      79

heights %>%
  filter(sex == "Female") %>%
  summarize(range = quantile(height, c(0, 0.5, 1)))

murders <- murders %>% mutate(rate = total/population*100000)

summarize(murders, mean(rate))

##   mean(rate)
## 1   2.779125

us_murder_rate <- murders %>%
  summarize(rate = sum(total) / sum(population) * 100000)
us_murder_rate

##       rate
## 1 3.034555

class(us_murder_rate)

## [1] "data.frame"

us_murder_rate %>% pull(rate)

## [1] 3.034555

us_murder_rate <- murders %>%
  summarize(rate = sum(total) / sum(population) * 100000) %>%
  pull(rate)

us_murder_rate

## [1] 3.034555

class(us_murder_rate)

## [1] "numeric"

heights %>% group_by(sex)

## # A tibble: 1,050 × 2
## # Groups:   sex [2]
##    sex    height
##    <fct>   <dbl>
##  1 Male       75
##  2 Male       70
##  3 Male       68
##  4 Male       74
##  5 Male       61
##  6 Female     65
##  7 Female     66
##  8 Female     62
##  9 Female     66
## 10 Male       67
## # ℹ 1,040 more rows

heights %>%
  group_by(sex) %>%
  summarize(average = mean(height), standard_deviation = sd(height))

## # A tibble: 2 × 3
##   sex    average standard_deviation
##   <fct>    <dbl>              <dbl>
## 1 Female    64.9               3.76
## 2 Male      69.3               3.61

murders %>%
  group_by(region) %>%
  summarize(median_rate = median(rate))

## # A tibble: 4 × 2
##   region        median_rate
##   <fct>               <dbl>
## 1 Northeast            1.80
## 2 South                3.40
## 3 North Central        1.97
## 4 West                 1.29

Sorting data frames

When examining a dataset, it is often convenient to sort the table by the different columns. We know about the order and sort function, but for ordering entire tables, the dplyr function arrange is useful. For example, here we order the states by population size:

murders %>%
  arrange(population) %>%
  head()

##                  state abb        region population total       rate
## 1              Wyoming  WY          West     563626     5  0.8871131
## 2 District of Columbia  DC         South     601723    99 16.4527532
## 3              Vermont  VT     Northeast     625741     2  0.3196211
## 4         North Dakota  ND North Central     672591     4  0.5947151
## 5               Alaska  AK          West     710231    19  2.6751860
## 6         South Dakota  SD North Central     814180     8  0.9825837

With arrange we get to decide which column to sort by. To see the states by murder rate, from lowest to highest, we arrange by rate instead:

murders %>%
  arrange(rate) %>%
  head()

##           state abb        region population total      rate
## 1       Vermont  VT     Northeast     625741     2 0.3196211
## 2 New Hampshire  NH     Northeast    1316470     5 0.3798036
## 3        Hawaii  HI          West    1360301     7 0.5145920
## 4  North Dakota  ND North Central     672591     4 0.5947151
## 5          Iowa  IA North Central    3046355    21 0.6893484
## 6         Idaho  ID          West    1567582    12 0.7655102

Note that the default behavior is to order in ascending order. In dplyr, the function desc transforms a vector so that it is in descending order. To sort the table in descending order, we can type:

murders %>%
  arrange(desc(rate))

murders %>%
  arrange(region, rate) %>%
  head()

##           state abb    region population total      rate
## 1       Vermont  VT Northeast     625741     2 0.3196211
## 2 New Hampshire  NH Northeast    1316470     5 0.3798036
## 3         Maine  ME Northeast    1328361    11 0.8280881
## 4  Rhode Island  RI Northeast    1052567    16 1.5200933
## 5 Massachusetts  MA Northeast    6547629   118 1.8021791
## 6      New York  NY Northeast   19378102   517 2.6679599

murders %>% top_n(5, rate)

##                  state abb        region population total      rate
## 1 District of Columbia  DC         South     601723    99 16.452753
## 2            Louisiana  LA         South    4533372   351  7.742581
## 3             Maryland  MD         South    5773552   293  5.074866
## 4             Missouri  MO North Central    5988927   321  5.359892
## 5       South Carolina  SC         South    4625364   207  4.475323

library(NHANES)

## Warning: package 'NHANES' was built under R version 4.2.3

data(NHANES)

library(dslabs)
data(na_example)
mean(na_example)

## [1] NA

sd(na_example)

## [1] NA

mean(na_example, na.rm = TRUE)

## [1] 2.301754

sd(na_example, na.rm = TRUE)

## [1] 1.22338

Tibbles

Tidy data must be stored in data frames. We introduced the data frame in our section on data.frames and have been using the murders data frame throughout the unit. In an earlier section we introduced the group_by function, which permits stratifying data before computing summary statistics. But where is the group information stored in the data frame?

murders %>% group_by(region)

## # A tibble: 51 × 6
## # Groups:   region [4]
##    state                abb   region    population total  rate
##    <chr>                <chr> <fct>          <dbl> <dbl> <dbl>
##  1 Alabama              AL    South        4779736   135  2.82
##  2 Alaska               AK    West          710231    19  2.68
##  3 Arizona              AZ    West         6392017   232  3.63
##  4 Arkansas             AR    South        2915918    93  3.19
##  5 California           CA    West        37253956  1257  3.37
##  6 Colorado             CO    West         5029196    65  1.29
##  7 Connecticut          CT    Northeast    3574097    97  2.71
##  8 Delaware             DE    South         897934    38  4.23
##  9 District of Columbia DC    South         601723    99 16.5 
## 10 Florida              FL    South       19687653   669  3.40
## # ℹ 41 more rows

Notice that there are no columns with this information. But, if you look closely at the output above, you see the line A tibble followed by dimensions. We can learn the class of the returned object using:

murders %>% group_by(region) %>% class()

## [1] "grouped_df" "tbl_df"     "tbl"        "data.frame"

The tbl, pronounced tibble, is a special kind of data frame. The functions group_by and summarize always return this type of data frame. The group_by function returns a special kind of tbl, the grouped_df. We will say more about these later. For consistency, the dplyr manipulation verbs (select, filter, mutate, and arrange) preserve the class of the input: if they receive a regular data frame they return a regular data frame, while if they receive a tibble they return a tibble. But tibbles are the preferred format in the tidyverse and as a result tidyverse functions that produce a data frame from scratch return a tibble.

Tibbles are very similar to data frames. In fact, you can think of them as a modern version of data frames. Nonetheless there are three important differences which we describe next.

class(murders[,4])

## [1] "numeric"

class(as_tibble(murders)[,4])

## [1] "tbl_df"     "tbl"        "data.frame"

class(as_tibble(murders)$population)

## [1] "numeric"

murders$Population

## NULL

as_tibble(murders)$Population

## Warning: Unknown or uninitialised column: `Population`.

## NULL

tibble(id = c(1, 2, 3), func = c(mean, median, sd))

## # A tibble: 3 × 2
##      id func  
##   <dbl> <list>
## 1     1 <fn>  
## 2     2 <fn>  
## 3     3 <fn>

grades <- tibble(names = c("John", "Juan", "Jean", "Yao"),
                     exam_1 = c(95, 80, 90, 85),
                     exam_2 = c(90, 85, 85, 90))

grades <- data.frame(names = c("John", "Juan", "Jean", "Yao"),
                     exam_1 = c(95, 80, 90, 85),
                     exam_2 = c(90, 85, 85, 90))
class(grades$names)

## [1] "character"

grades <- data.frame(names = c("John", "Juan", "Jean", "Yao"),
                     exam_1 = c(95, 80, 90, 85),
                     exam_2 = c(90, 85, 85, 90),
                     stringsAsFactors = FALSE)
class(grades$names)

## [1] "character"

as_tibble(grades) %>% class()

## [1] "tbl_df"     "tbl"        "data.frame"

The dot operator

One of the advantages of using the pipe %>% is that we do not have to keep naming new objects as we manipulate the data frame. As a quick reminder, if we want to compute the median murder rate for states in the southern states, instead of typing:

tab_1 <- filter(murders, region == "South")
tab_2 <- mutate(tab_1, rate = total / population * 10^5)
rates <- tab_2$rate
median(rates)

## [1] 3.398069

We can avoid defining any new intermediate objects by instead typing:

filter(murders, region == "South") %>%
  mutate(rate = total / population * 10^5) %>%
  summarize(median = median(rate)) %>%
  pull(median)

## [1] 3.398069

We can do this because each of these functions takes a data frame as the first argument. But what if we want to access a component of the data frame. For example, what if the pull function was not available and we wanted to access tab_2$rate? What data frame name would we use? The answer is the dot operator.

For example to access the rate vector without the pull function we could use

rates <-   filter(murders, region == "South") %>%
  mutate(rate = total / population * 10^5) %>%
  .$rate
median(rates)

## [1] 3.398069

In the next section, we will see other instances in which using the . is useful.

do

The tidyverse functions know how to interpret grouped tibbles. Furthermore, to facilitate stringing commands through the pipe %>%, tidyverse functions consistently take data frames and return data frames, since this assures that the output of a function is accepted as the input of another. But most R functions do not recognize grouped tibbles nor do they return data frames. The quantile function is an example we described earlier. The do function serves as a bridge between R functions such as quantile and the tidyverse. The do function understands grouped tibbles and always returns a data frame.

In the summarize section (above), we noted that if we attempt to use quantile to obtain the min, median and max in one call, we will receive something unexpected. Prior to R 4.1, we would receive an error. After R 4.1, we actually get:

data(heights)
heights %>%
  filter(sex == "Female") %>%
  summarize(range = quantile(height, c(0, 0.5, 1)))

We probably wanted three columns: min, median, and max. We can use the do function to fix this.

First we have to write a function that fits into the tidyverse approach: that is, it receives a data frame and returns a data frame. Note that it returns a single-row data frame.

my_summary <- function(dat){
  x <- quantile(dat$height, c(0, 0.5, 1))
  tibble(min = x[1], median = x[2], max = x[3])
}

We can now apply the function to the heights dataset to obtain the summaries:

heights %>%
  group_by(sex) %>%
  my_summary

## # A tibble: 1 × 3
##     min median   max
##   <dbl>  <dbl> <dbl>
## 1    50   68.5  82.7

But this is not what we want. We want a summary for each sex and the code returned just one summary. This is because my_summary is not part of the tidyverse and does not know how to handled grouped tibbles. do makes this connection:

heights %>%
  group_by(sex) %>%
  do(my_summary(.))

## # A tibble: 2 × 4
## # Groups:   sex [2]
##   sex      min median   max
##   <fct>  <dbl>  <dbl> <dbl>
## 1 Female    51   65.0  79  
## 2 Male      50   69    82.7

Note that here we need to use the dot operator. The tibble created by group_by is piped to do. Within the call to do, the name of this tibble is . and we want to send it to my_summary. If you do not use the dot, then my_summary has no argument and returns an error telling us that argument "dat" is missing. You can see the error by typing:

heights %>%
  group_by(sex) %>%
  do(my_summary())

If you do not use the parenthesis, then the function is not executed and instead do tries to return the function. This gives an error because do must always return a data frame. You can see the error by typing:

heights %>%
  group_by(sex) %>%
  do(my_summary)

So do serves as a bridge between non-tidyverse functions and the tidyverse.

The purrr package

In previous sections (and labs) we learned about the sapply function, which permitted us to apply the same function to each element of a vector. We constructed a function and used sapply to compute the sum of the first n integers for several values of n like this:

compute_s_n <- function(n){
  x <- 1:n
  sum(x)
}
n <- 1:25
s_n <- sapply(n, compute_s_n)
s_n

##  [1]   1   3   6  10  15  21  28  36  45  55  66  78  91 105 120 136 153 171 190
## [20] 210 231 253 276 300 325

This type of operation, applying the same function or procedure to elements of an object, is quite common in data analysis. The purrr package includes functions similar to sapply but that better interact with other tidyverse functions. The main advantage is that we can better control the output type of functions. In contrast, sapply can return several different object types; for example, we might expect a numeric result from a line of code, but sapply might convert our result to character under some circumstances. purrr functions will never do this: they will return objects of a specified type or return an error if this is not possible.

The first purrr function we will learn is map, which works very similar to sapply but always, without exception, returns a list:

library(purrr) # or library(tidyverse)
n <- 1:25
s_n <- map(n, compute_s_n)
class(s_n)

## [1] "list"

If we want a numeric vector, we can instead use map_dbl which always returns a vector of numeric values.

s_n <- map_dbl(n, compute_s_n)
class(s_n)

## [1] "numeric"

This produces the same results as the sapply call shown above.

A particularly useful purrr function for interacting with the rest of the tidyverse is map_df, which always returns a tibble data frame. However, the function being called needs to return a vector or a list with names. For this reason, the following code would result in a Argument 1 must have names error:

s_n <- map_df(n, compute_s_n)

We need to change the function to make this work:

compute_s_n <- function(n){
  x <- 1:n
  tibble(sum = sum(x))
}
s_n <- map_df(n, compute_s_n)
head(s_n)

## # A tibble: 6 × 1
##     sum
##   <int>
## 1     1
## 2     3
## 3     6
## 4    10
## 5    15
## 6    21

Because map_df returns a tibble, we can have more columns defined in our function and returned.

compute_s_n2 <- function(n){
  x <- 1:n
  tibble(sum = sum(x), sumSquared = sum(x^2))
}
s_n <- map_df(n, compute_s_n2)
head(s_n)

## # A tibble: 6 × 2
##     sum sumSquared
##   <int>      <dbl>
## 1     1          1
## 2     3          5
## 3     6         14
## 4    10         30
## 5    15         55
## 6    21         91

The purrr package provides much more functionality not covered here. For more details you can consult this online resource.

Tidyverse conditionals

A typical data analysis will often involve one or more conditional operations. In the section on Conditionals, we described the ifelse function, which we will use extensively in this book. In this section we present two dplyr functions that provide further functionality for performing conditional operations.

x <- c(-2, -1, 0, 1, 2)
case_when(x < 0 ~ "Negative",
          x > 0 ~ "Positive",
          x == 0  ~ "Zero")

## [1] "Negative" "Negative" "Zero"     "Positive" "Positive"

murders %>%
  mutate(group = case_when(
    abb %in% c("ME", "NH", "VT", "MA", "RI", "CT") ~ "New England",
    abb %in% c("WA", "OR", "CA") ~ "West Coast",
    region == "South" ~ "South",
    TRUE ~ "Other")) %>%
  group_by(group) %>%
  summarize(rate = sum(total) / sum(population) * 10^5)

## # A tibble: 4 × 2
##   group        rate
##   <chr>       <dbl>
## 1 New England  1.72
## 2 Other        2.71
## 3 South        3.63
## 4 West Coast   2.90

x >= a & x <= b

between(x, a, b)

exp(mean(log(murders$population)))