R for BiologistsWeek 5 – Data manipulation with dplyr
Learning objectives
- Chain operations together into a workflow using pipes
- Learn the key dplyr functions for manipulating your data
- Learn about faceting in ggplot2 to split your data into separate categories and create a series of sub-plots arranged in a grid
Data manipulation with dplyr
dplyr is one of the packages that gets loaded as part of the tidyverse.
library(tidyverse)## ── Attaching core tidyverse packages ────────────────────────────────────────────────────────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.5
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.5.1 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.1
## ✔ purrr 1.0.2
## ── Conflicts ──────────────────────────────────────────────────────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorsdplyr is the Swiss army knife in the tidyverse, providing many useful functions for manipulating tabular data in data frames or tibbles. We’re going to look at the key functions for filtering our data, modifying the contents and computing summary statistics.
We’ll also introduce the pipe operator, %>%, for chaining operations together into mini workflows in a way that makes for more readable and maintainable code.
Finally, we’ll return to plotting and look at a powerful feature of ggplot2, faceting, that allows you to divide your plots into subplots by splitting the observations based on one or more categorical variables.
We’ll again use the METABRIC data set to illustrate how these operations work.
metabric <- read_csv("data/metabric_clinical_and_expression_data.csv")
metabric## # A tibble: 1,904 × 32
## Patient_ID Cohort Age_at_diagnosis Survival_time Survival_status Vital_status
## <chr> <dbl> <dbl> <dbl> <chr> <chr>
## 1 MB-0000 1 75.6 140. LIVING Living
## 2 MB-0002 1 43.2 84.6 LIVING Living
## 3 MB-0005 1 48.9 164. DECEASED Died of Dis…
## 4 MB-0006 1 47.7 165. LIVING Living
## 5 MB-0008 1 77.0 41.4 DECEASED Died of Dis…
## 6 MB-0010 1 78.8 7.8 DECEASED Died of Dis…
## 7 MB-0014 1 56.4 164. LIVING Living
## 8 MB-0022 1 89.1 99.5 DECEASED Died of Oth…
## 9 MB-0028 1 86.4 36.6 DECEASED Died of Oth…
## 10 MB-0035 1 84.2 36.3 DECEASED Died of Dis…
## # ℹ 1,894 more rows
## # ℹ 26 more variables: Chemotherapy <chr>, Radiotherapy <chr>,
## # Tumour_size <dbl>, Tumour_stage <dbl>, Neoplasm_histologic_grade <dbl>,
## # Lymph_nodes_examined_positive <dbl>, Lymph_node_status <dbl>,
## # Cancer_type <chr>, ER_status <chr>, PR_status <chr>, HER2_status <chr>,
## # HER2_status_measured_by_SNP6 <chr>, PAM50 <chr>, `3-gene_classifier` <chr>,
## # Nottingham_prognostic_index <dbl>, Cellularity <chr>, …dplyr verbs
The dplyr operations are commonly referred to as “verbs” in a data manipulation grammar. These verbs have a common syntax and work together in a consistent and uniform manner. They all have the following shared behaviours:
The first argument in each function is a data frame (or tibble)
Any additional arguments describe what operation to perform on the data frame
Variable names, i.e. column names, are referred to without using quotes
The result of an operation is a new data frame
Chaining operations using %>%
Let’s consider again an earlier example in which we filtered the METABRIC data set to retain just the patients who were still alive at the time of the study and had survived for more than 10 years (120 months). We use filter() to select the rows corresponding to the patients meeting these criteria and can then use select() to only display the variables (columns) we’re most interested in.
patients <- filter(metabric,
Survival_status == "LIVING" & Survival_time > 120)
patient_details <- select(patients,
Patient_ID,
Survival_time,
Tumour_stage,
Nottingham_prognostic_index)
patient_details## # A tibble: 545 × 4
## Patient_ID Survival_time Tumour_stage Nottingham_prognostic_index
## <chr> <dbl> <dbl> <dbl>
## 1 MB-0000 140. 2 6.04
## 2 MB-0006 165. 2 4.05
## 3 MB-0014 164. 2 4.02
## 4 MB-0039 164. 1 2.04
## 5 MB-0045 165. 2 5.04
## 6 MB-0053 161. 2 3.05
## 7 MB-0054 160. 2 4.07
## 8 MB-0060 141. 2 4.05
## 9 MB-0062 154. 1 4.03
## 10 MB-0066 157. 2 4.03
## # ℹ 535 more rowsHere we’ve used an intermediate variable, patients, which we only needed in order to get to the final result. We could just have used the same name to avoid cluttering our environment and overwritten the results from the filter() operation with those of the select() operation.
patients <- select(patients,
Patient_ID,
Survival_time,
Tumour_stage,
Nottingham_prognostic_index)In order to avoid creating intermediate objects, another less readable way of writing this code is to nest the filter() function call inside the select(). Nested function calls are very common in a lot of R code you may come across. For example, imagine we have a vector of numbers and we want to find the square root of the mean. We could do it two steps:
x <- c(24, 53, 352, 643, 46)
myMean <- mean(x)
sqrt(myMean)## [1] 14.95326Or we can nest the commands inside each other:
sqrt(mean(x))## [1] 14.95326This is convenient for simple commands, but is unwieldy and not easy to read with more complex commands:
patients <- select(filter(metabric, Survival_status == "LIVING", Survival_time > 120), Patient_ID, Survival_time, Tumour_stage, Nottingham_prognostic_index)
nrow(patients)## [1] 545However, there is another way chaining together a series of operations into a mini workflow that is elegant, intuitive and makes for very readable R code. For that we need to introduce a new operator, the pipe operator, %>%.
The pipe operator %>%
The pipe operator takes the output from one operation, i.e. whatever is on the left-hand side of %>% and passes it in as the first argument to the second operation, or function, on the right-hand side.
x %>% f(y) is equivalent to f(x, y)
For example:
select(starwars, name, height, mass)
can be rewritten as
starwars %>% select(name, height, mass)
This allows for chaining of operations into workflows, e.g.
starwars %>%
filter(species == “Droid”) %>%
select(name, height, mass)
The %>% operator comes from the magrittr package and is available when we load the tidyverse using library(tidyverse).
Piping in R was motivated by the Unix pipe, |, in which the output from one process is redirected to be the input for the next. This is so named because the flow from one process or operation to the next resembles a pipeline.
NOTE: The %>% operator is from the packages magrittr. R now has a native pipe operator, |>, but we’ll stick with the %>% operator as it is more widely used and has a few advantages over the native pipe.
We can rewrite the code for our filtering and column selection operations as follows.
patients <- metabric %>%
filter(Survival_status == "LIVING" & Survival_time > 120) %>%
select(Patient_ID, Survival_time, Tumour_stage, Nottingham_prognostic_index)Note how each operation takes the output from the previous operation as its first argument. This way of coding is embraced wholeheartedly in the tidyverse hence almost every tidyverse function that works on data frames has the data frame as its first argument. It is also the reason why tidyverse functions return a data frame regardless of whether the output could be recast as a vector or a single value.
“Piping”, the act of chaining operations together, becomes really useful when there are several steps involved in filtering and transforming a data set.
The usual way of developing a workflow is to build it up one step at a time, testing the output produced at each stage. Let’s do that for this case.
We start by considering just the patients who are living.
patients <- metabric %>%
filter(Survival_status == "LIVING")We then add another filter for the survival time.
patients <- metabric %>%
filter(Survival_status == "LIVING") %>%
filter(Survival_time > 120)At each stage we look at the resulting patients data frame to check we’re happy with the result.
Finally we only want certain columns, so we add a select() operation.
patients <- metabric %>%
filter(Survival_status == "LIVING") %>%
filter(Survival_time > 120) %>%
select(Patient_ID, Survival_time, Tumour_stage, Nottingham_prognostic_index)
# print out the result
patients## # A tibble: 545 × 4
## Patient_ID Survival_time Tumour_stage Nottingham_prognostic_index
## <chr> <dbl> <dbl> <dbl>
## 1 MB-0000 140. 2 6.04
## 2 MB-0006 165. 2 4.05
## 3 MB-0014 164. 2 4.02
## 4 MB-0039 164. 1 2.04
## 5 MB-0045 165. 2 5.04
## 6 MB-0053 161. 2 3.05
## 7 MB-0054 160. 2 4.07
## 8 MB-0060 141. 2 4.05
## 9 MB-0062 154. 1 4.03
## 10 MB-0066 157. 2 4.03
## # ℹ 535 more rowsWhen continuing our workflow across multiple lines, we need to be careful to ensure the %>% is at the end of the line. If we try to place this at the start of the next line, R will think we’ve finished the workflow prematurely and will report an error at what it considers the next statement, i.e. the line that begins with %>%.
# R considers the following to be 2 separate commands, the first of which ends
# with the first filter operation and runs successfully.
# The second statement is the third and last line, is not a valid commmand and
# so you'll get an error message
patients <- metabric %>%
filter(Survival_status == "LIVING")
%>% filter(Survival_time > 120)This is very similar to what we saw with adding layers and other components to a ggplot using the + operator, where we needed the + to be at the end of a line.
We’ll be using the pipe %>% operator throughout the rest of the course so you’d better get used to it. But actually we think you’ll come to love it as much as we do.
Modifying data using mutate()
In one of the examples of filtering observations using filter() above, we created a new logical variable called Deceased.
metabric$Deceased <- metabric$Survival_status == "DECEASED"This is an example of a rather common type of data manipulation in which we create a new column based on the values contained in one or more other columns. The dplyr package provides the mutate() function for this purpose.
metabric <- mutate(metabric, Deceased = Survival_status == "DECEASED")Both of these methods adds the new column to the end.
Note that variables names in the mutate() function call do not need to be prefixed by metabric$ just as they don’t in any of the dplyr functions.
We can overwrite a column and this is quite commonly done when we want to change the units. For example, our tumour size values are in millimetres but the formula for the Nottingham prognostic index needs the tumour size to be in centimetres. We can convert to Tumour_size to centimetres by dividing by 10.
metabric <- mutate(metabric, Tumour_size = Tumour_size / 10)Nottingham Prognostic Index
The Nottingham prognostic index (NPI) is used to determine prognosis following surgery for breast cancer. Its value is calculated using three pathological criteria: the size of the tumour, the number of lymph nodes involved, and the grade of the tumour.
The formula for the Nottingham prognostic index is:
NPI = (0.2 * S) + N + G
where
- S is the size of the tumour in centimetres
- N is the node status (0 nodes = 1, 1-3 nodes = 2, >3 nodes = 3)
- G is the grade of tumour (Grade I = 1, Grade II = 2, Grade III = 3)
We’ll recalculate the Nottingham prognostic index using the values Tumour_size, Neoplasm_histologic_grade and Lymph_node_status columns in our METABRIC data frame.
metabric <- metabric %>%
mutate(NPI = (0.2 * Tumour_size) +
Lymph_node_status +
Neoplasm_histologic_grade)
select(metabric,
Tumour_size,
Lymph_node_status,
Neoplasm_histologic_grade,
NPI,
Nottingham_prognostic_index)## # A tibble: 1,904 × 5
## Tumour_size Lymph_node_status Neoplasm_histologic_grade NPI
## <dbl> <dbl> <dbl> <dbl>
## 1 2.2 3 3 6.44
## 2 1 1 3 4.2
## 3 1.5 2 2 4.3
## 4 2.5 2 2 4.5
## 5 4 3 3 6.8
## 6 3.1 1 3 4.62
## 7 1 2 2 4.2
## 8 2.9 2 2 4.58
## 9 1.6 2 3 5.32
## 10 2.8 1 2 3.56
## # ℹ 1,894 more rows
## # ℹ 1 more variable: Nottingham_prognostic_index <dbl>There is a discrepency. Can you see what the problem is?
It appears that the tumour size wasn’t correctly converted into centimetres in the original NPI calculation or that the wrong scaling factor for tumour size was used.
ggplot(data = metabric) +
geom_point(mapping = aes(x = Age_at_diagnosis, y = NPI), na.rm = TRUE)
The round() function is really useful for rounding numerical values to a specified number of decimal places. We’ll read in the METABRIC data again and create a small workflow that carries out the tumour size conversion, computes the NPI, rounds the tumour size and the resulting NPI value to 1 decimal place and displays the results in decreasing order of NPI.
read_csv("data/metabric_clinical_and_expression_data.csv") %>%
mutate(Tumour_size = Tumour_size / 10) %>%
mutate(NPI = (0.2 * Tumour_size) +
Lymph_node_status +
Neoplasm_histologic_grade) %>%
mutate(Tumour_size = round(Tumour_size, digits = 1)) %>%
mutate(NPI = round(NPI, digits = 1)) %>%
arrange(desc(NPI)) %>%
select(Tumour_size, Lymph_node_status, Neoplasm_histologic_grade, NPI)## # A tibble: 1,904 × 4
## Tumour_size Lymph_node_status Neoplasm_histologic_grade NPI
## <dbl> <dbl> <dbl> <dbl>
## 1 18 3 3 9.6
## 2 16 3 3 9.2
## 3 15 3 3 9
## 4 13 3 3 8.6
## 5 18.2 2 3 8.6
## 6 12 3 3 8.4
## 7 10 3 3 8
## 8 9.9 3 3 8
## 9 9.9 3 3 8
## 10 8.4 3 3 7.7
## # ℹ 1,894 more rowsMutating multiple columns
In that last workflow we included the same rounding operation applied to two different variables. It would be nice to be able to carry out just the one mutate() but apply it to both Tumour_size and NPI columns and we can using across().
across is used inside mutate and takes two arguments, first a selection of columns to operate across and second a function to apply to those columns.
So if we just wanted to round to the the nearest integer:
metabric %>%
mutate(across(c(Tumour_size, NPI), round)) %>%
select(Patient_ID, Tumour_size, NPI)## # A tibble: 1,904 × 3
## Patient_ID Tumour_size NPI
## <chr> <dbl> <dbl>
## 1 MB-0000 2 6
## 2 MB-0002 1 4
## 3 MB-0005 2 4
## 4 MB-0006 2 4
## 5 MB-0008 4 7
## 6 MB-0010 3 5
## 7 MB-0014 1 4
## 8 MB-0022 3 5
## 9 MB-0028 2 5
## 10 MB-0035 3 4
## # ℹ 1,894 more rowsThis applies the function round across all columns specified in the first argument for across.
This is nice and simple but things can get a little more complex if we want more than just a straightforward function.
We may come across situations where we’d like to apply the same operation to multiple columns but where there is no available function in R to do what we want or situations we need to pass extra arguments to our function.
In these cases the syntax get a little trickier. We need to write what is known as an “anonymous function”.
Let’s say we want to convert the petal and sepal measurements in the iris data set from centimetres to millimetres. We’d either need to create a new function to do this conversion or we could use what is known as an anonymous function.
There is no ‘multiply by 10’ function and it seems a bit pointless to create one just for this conversion so we’ll use an anonymous function instead – anonymous because it has no name, it’s an in situ function only used in our mutate() function call.
iris %>%
as_tibble() %>%
mutate(across(Sepal.Length:Petal.Width, ~ .x * 10))## # A tibble: 150 × 5
## Sepal.Length Sepal.Width Petal.Length Petal.Width Species
## <dbl> <dbl> <dbl> <dbl> <fct>
## 1 51 35 14 2 setosa
## 2 49 30 14 2 setosa
## 3 47 32 13 2 setosa
## 4 46 31 15 2 setosa
## 5 50 36 14 2 setosa
## 6 54 39 17 4 setosa
## 7 46 34 14 3 setosa
## 8 50 34 15 2 setosa
## 9 44 29 14 2 setosa
## 10 49 31 15 1 setosa
## # ℹ 140 more rowsThe ~ denotes that we’re using an anonymous function (it is the symbol for formulae in R) and the .x is a placeholder for the column being operated on. In this case, we’re multiplying each of the columns between Sepal.Length and Petal.Width inclusive by 10.
If we wish to use an existing function but provide additional arguments, then the syntax is:
~func(.x, ...)
Where .x is going to be replaced with the column names and ... is the extra arguments we need to add.
So if we want to round to 1 decimal place:
metabric %>%
mutate(across(c(Tumour_size, NPI), ~round(.x, digits = 1))) %>%
select(Patient_ID, Tumour_size, NPI)## # A tibble: 1,904 × 3
## Patient_ID Tumour_size NPI
## <chr> <dbl> <dbl>
## 1 MB-0000 2.2 6.4
## 2 MB-0002 1 4.2
## 3 MB-0005 1.5 4.3
## 4 MB-0006 2.5 4.5
## 5 MB-0008 4 6.8
## 6 MB-0010 3.1 4.6
## 7 MB-0014 1 4.2
## 8 MB-0022 2.9 4.6
## 9 MB-0028 1.6 5.3
## 10 MB-0035 2.8 3.6
## # ℹ 1,894 more rowsWe can use the range operator and the same helper functions as we did for selecting columns using select() inside across().
For example, we might decide that our expression values are given to a much higher degree of precision than is strictly necessary.
metabric %>%
mutate(across(ESR1:MLPH, ~round(.x, digits = 2))) %>%
select(Patient_ID, ESR1:MLPH)## # A tibble: 1,904 × 9
## Patient_ID ESR1 ERBB2 PGR TP53 PIK3CA GATA3 FOXA1 MLPH
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 MB-0000 8.93 9.33 5.68 6.34 5.7 6.93 7.95 9.73
## 2 MB-0002 10.0 9.73 7.51 6.19 5.76 11.2 11.8 12.5
## 3 MB-0005 10.0 9.73 7.38 6.4 6.75 9.29 11.7 10.3
## 4 MB-0006 10.4 10.3 6.82 6.87 7.22 8.67 11.9 10.5
## 5 MB-0008 11.3 9.96 7.33 6.34 5.82 9.72 11.6 12.2
## 6 MB-0010 11.2 9.74 5.95 5.42 6.12 9.79 12.1 11.4
## 7 MB-0014 10.8 9.28 7.72 5.99 7.48 8.37 11.5 10.8
## 8 MB-0022 10.4 8.61 5.59 6.17 7.59 7.87 10.7 9.95
## 9 MB-0028 12.5 10.7 5.33 6.22 6.25 10.3 12.2 10.9
## 10 MB-0035 7.54 11.5 5.59 6.41 5.99 10.2 12.8 13.5
## # ℹ 1,894 more rowsOr we could decide that all the columns whose names end with “_status” are in fact categorical variables and should be converted to factors.
metabric %>%
mutate(across(ends_with("_status"), as.factor)) %>%
select(Patient_ID, ends_with("_status"))## # A tibble: 1,904 × 7
## Patient_ID Survival_status Vital_status Lymph_node_status ER_status PR_status
## <chr> <fct> <fct> <fct> <fct> <fct>
## 1 MB-0000 LIVING Living 3 Positive Negative
## 2 MB-0002 LIVING Living 1 Positive Positive
## 3 MB-0005 DECEASED Died of Dis… 2 Positive Positive
## 4 MB-0006 LIVING Living 2 Positive Positive
## 5 MB-0008 DECEASED Died of Dis… 3 Positive Positive
## 6 MB-0010 DECEASED Died of Dis… 1 Positive Positive
## 7 MB-0014 LIVING Living 2 Positive Positive
## 8 MB-0022 DECEASED Died of Oth… 2 Positive Negative
## 9 MB-0028 DECEASED Died of Oth… 2 Positive Negative
## 10 MB-0035 DECEASED Died of Dis… 1 Positive Negative
## # ℹ 1,894 more rows
## # ℹ 1 more variable: HER2_status <fct>Computing summary values using summarise()
We’ll cover the fifth of the main dplyr ‘verb’ functions, summarise(), only briefly here. This function computes summary values for one or more variables (columns) in a table. Here we will summarise values for the entire table but this function is much more useful in combination with group_by() in working on groups of observations within the data set. We will look at summarizing groups of observations next week.
Any function that calculates a single scalar value from a vector can be used with summarise(). For example, the mean() function calculates the arithmetic mean of a numeric vector. Let’s calculate the average ESR1 expression for tumour samples in the METABRIC data set.
mean(metabric$ESR1)## [1] 9.607824The equivalent operation using summarise() is:
summarise(metabric, mean(ESR1))## # A tibble: 1 × 1
## `mean(ESR1)`
## <dbl>
## 1 9.61If you prefer Oxford spelling, in which -ize is preferred to -ise, you’re in luck as dplyr accommodates the alternative spelling.
Both of the above statements gave the same average expression value but these were output in differing formats. The mean() function collapses a vector to single scalar value, which as we know is in fact a vector of length 1. The summarise() function, as with most tidyverse functions, returns another data frame, albeit one in which there is a single row and a single column.
Returning a data frame might be quite useful, particularly if we’re summarising multiple columns or using more than one function, for example computing the average and standard deviation.
summarise(metabric, ESR1_mean = mean(ESR1), ESR1_sd = sd(ESR1))## # A tibble: 1 × 2
## ESR1_mean ESR1_sd
## <dbl> <dbl>
## 1 9.61 2.13Notice how we also named the output columns in this last example.
summarise()
summarise() collapses a data frame into a single row by calculating summary values of one or more of the columns.
It can take any function that takes a vector of values and returns a single value. Some of the more useful functions include:
-
Centre:
mean(),median() -
Spread:
sd(),mad() -
Range:
min(),max(),quantile() -
Position:
first(),last() -
Count:
n()
Note the first(), last() and n() are only really useful when working on groups of observations using group_by().
n() is a special function that returns the number of observations; it doesn’t take a vector argument, i.e. a column.
It is also possible to summarise using a function that takes more than one value, i.e. from multiple columns. For example, we could compute the correlation between the expression of FOXA1 and MLPH.
summarise(metabric, correlation = cor(FOXA1, MLPH))## # A tibble: 1 × 1
## correlation
## <dbl>
## 1 0.898Summarizing multiple columns
As with mutate() we can use across() to summarise across multiple columns using the same function.
summarise(metabric, across(c(FOXA1, MLPH), mean))## # A tibble: 1 × 2
## FOXA1 MLPH
## <dbl> <dbl>
## 1 10.8 11.4We can use the everything() selector to summarise across all the columns in the table.
metabric %>%
select(ESR1:MLPH) %>%
summarise(across(everything(), mean))## # A tibble: 1 × 8
## ESR1 ERBB2 PGR TP53 PIK3CA GATA3 FOXA1 MLPH
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 9.61 10.8 6.24 6.20 5.97 9.50 10.8 11.4But, you have to be careful that all columns can be summarised with the given summary function. For example, what happens if we try to compute an average of a set of character values?
metabric %>%
summarise(across(everything(), ~mean(.x, na.rm = TRUE)))## Warning: There were 14 warnings in `summarise()`.
## The first warning was:
## ℹ In argument: `across(everything(), ~mean(.x, na.rm = TRUE))`.
## Caused by warning in `mean.default()`:
## ! argument is not numeric or logical: returning NA
## ℹ Run `dplyr::last_dplyr_warnings()` to see the 13 remaining warnings.## # A tibble: 1 × 34
## Patient_ID Cohort Age_at_diagnosis Survival_time Survival_status Vital_status
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 NA 2.64 61.1 125. NA NA
## # ℹ 28 more variables: Chemotherapy <dbl>, Radiotherapy <dbl>,
## # Tumour_size <dbl>, Tumour_stage <dbl>, Neoplasm_histologic_grade <dbl>,
## # Lymph_nodes_examined_positive <dbl>, Lymph_node_status <dbl>,
## # Cancer_type <dbl>, ER_status <dbl>, PR_status <dbl>, HER2_status <dbl>,
## # HER2_status_measured_by_SNP6 <dbl>, PAM50 <dbl>, `3-gene_classifier` <dbl>,
## # Nottingham_prognostic_index <dbl>, Cellularity <dbl>,
## # Integrative_cluster <dbl>, Mutation_count <dbl>, ESR1 <dbl>, ERBB2 <dbl>, …We get a lot of warning messages and NA values for those columns for which computing an average does not make sense.
We can combine across() with where() and a logical test (e.g. is.numeric) to select those columns for which a summarization function is appropriate.
metabric %>%
summarise(across(where(is.numeric), ~mean(.x, na.rm = TRUE)))## # A tibble: 1 × 19
## Cohort Age_at_diagnosis Survival_time Tumour_size Tumour_stage
## <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 2.64 61.1 125. 2.62 1.75
## # ℹ 14 more variables: Neoplasm_histologic_grade <dbl>,
## # Lymph_nodes_examined_positive <dbl>, Lymph_node_status <dbl>,
## # Nottingham_prognostic_index <dbl>, Mutation_count <dbl>, ESR1 <dbl>,
## # ERBB2 <dbl>, PGR <dbl>, TP53 <dbl>, PIK3CA <dbl>, GATA3 <dbl>, FOXA1 <dbl>,
## # MLPH <dbl>, NPI <dbl>It is possible to summarise using more than one function in which case a list of functions needs to be provided.
metabric %>%
summarise(across(c(ESR1, ERBB2, PGR), list(mean, sd)))## # A tibble: 1 × 6
## ESR1_1 ESR1_2 ERBB2_1 ERBB2_2 PGR_1 PGR_2
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 9.61 2.13 10.8 1.36 6.24 1.02Pretty neat but I’m not sure about those column headings in the output – fortunately we have some control over these.
metabric %>%
summarise(across(c(ESR1, ERBB2, PGR), list(Average = mean, stdev = sd)))## # A tibble: 1 × 6
## ESR1_Average ESR1_stdev ERBB2_Average ERBB2_stdev PGR_Average PGR_stdev
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 9.61 2.13 10.8 1.36 6.24 1.02Anonymous functions can also be used with summarise.
For example, suppose we want to compute the correlation of the expression for FOXA1 against all other genes to see which was most strongly correlated. Here is how we could do this in a single summarise() statement using an anonymous function.
metabric %>%
summarise(across(c(ESR1:MLPH, -FOXA1), ~cor(.x, FOXA1)))## # A tibble: 1 × 7
## ESR1 ERBB2 PGR TP53 PIK3CA GATA3 MLPH
## <dbl> <dbl> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 0.724 0.280 0.390 -0.0700 -0.149 0.781 0.898Notice how we selected all genes between ESR1, the first gene column in our data frame, and MLPH, the last gene column, but then excluded FOXA1 as we’re not all that interested in the correlation of FOXA1 with itself (we know the answer is 1).
Faceting with ggplot2
Finally, let’s change track completely and take a look at a very useful feature of ggplot2 – faceting.
Faceting allows you to split your plot into subplots, or facets, based on one or more categorical variables. Each of the subplots displays a subset of the data.
There are two faceting functions, facet_wrap() and facet_grid().
Let’s create a scatter plot of GATA3 and ESR1 expression values where we’re displaying the PR positive and PR negative patients using different colours. This is a very similar to a plot we created last week.
ggplot(data = metabric,
mapping = aes(x = GATA3, y = ESR1, colour = PR_status)) +
geom_point(size = 0.5, alpha = 0.5)
An alternative is to use faceting with facet_wrap().
ggplot(data = metabric, mapping = aes(x = GATA3, y = ESR1)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_wrap(vars(PR_status))
This produces two plots, side-by-side, one for each of the categories in the PR_status variable, with a banner across the top of each for the category.
The variable(s) used for faceting are specified using vars() in a similar way to the selection of variables for mutate_at() and summarise_at().
We can still use separate colours if we prefer things to be, well, colourful.
ggplot(data = metabric,
mapping = aes(x = GATA3, y = ESR1, colour = PR_status)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_wrap(vars(PR_status))
Faceting is usually better than displaying groups using different colours when there are more than two or three groups when it can be difficult to really tell which points belong to each group. A case in point is for the 3-gene classification in the GATA3 vs ESR1 scatter plot we created last week. Let’s create a faceted version of that plot.
ggplot(data = metabric,
mapping = aes(x = GATA3, y = ESR1, colour = `3-gene_classifier`)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_wrap(vars(`3-gene_classifier`))
This helps explain why the function is called facet_wrap(). When it has too many subplots to fit across the page, it wraps around to another row. We can control how many rows or columns to use with the nrow and ncol arguments.
ggplot(data = metabric,
mapping = aes(x = GATA3, y = ESR1, colour = `3-gene_classifier`)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_wrap(vars(`3-gene_classifier`), nrow = 1)
ggplot(data = metabric,
mapping = aes(x = GATA3, y = ESR1, colour = `3-gene_classifier`)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_wrap(vars(`3-gene_classifier`), ncol = 2)
We can combine faceting on one variable with a colour aesthetic for another variable. For example, let’s show the tumour stage status using faceting and the HER2 status using colours.
ggplot(data = metabric,
mapping = aes(x = GATA3, y = ESR1, colour = HER2_status)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_wrap(vars(Neoplasm_histologic_grade))
Instead of this we could facet on more than variable.
ggplot(data = metabric, mapping = aes(x = GATA3, y = ESR1)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_wrap(vars(Neoplasm_histologic_grade, HER2_status))
Faceting on two variables is usually better done using the other faceting function, facet_grid(). Note the change in how the formula is written.
ggplot(data = metabric, mapping = aes(x = GATA3, y = ESR1)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_grid(vars(Neoplasm_histologic_grade), vars(HER2_status))
Again we can use colour aesthetics alongside faceting to add further information to our visualization.
ggplot(data = metabric, mapping = aes(x = GATA3, y = ESR1, colour = PAM50)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_grid(rows = vars(Neoplasm_histologic_grade),
cols = vars(HER2_status))
Finally, we can use a labeller to change the labels for each of the categorical values so that these are more meaningful in the context of this plot.
grade_labels <- c("1" = "Grade I", "2" = "Grade II", "3" = "Grade III")
her2_status_labels <- c("Positive" = "HER2 positive",
"Negative" = "HER2 negative")
ggplot(data = metabric, mapping = aes(x = GATA3, y = ESR1, colour = PAM50)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_grid(rows = vars(Neoplasm_histologic_grade),
cols = vars(HER2_status),
labeller = labeller(Neoplasm_histologic_grade = grade_labels,
HER2_status = her2_status_labels))
This would certainly be necessary if we were to use ER and HER2 status on one side of the grid.
er_status_labels <- c("Positive" = "ER positive", "Negative" = "ER negative")
#
ggplot(data = metabric, mapping = aes(x = GATA3, y = ESR1, colour = PAM50)) +
geom_point(size = 0.5, alpha = 0.5) +
facet_grid(rows = vars(Neoplasm_histologic_grade),
cols = vars(ER_status, HER2_status),
labeller = labeller(Neoplasm_histologic_grade = grade_labels,
ER_status = er_status_labels,
HER2_status = her2_status_labels))
Summary
In this session we have covered the following concepts:
- Building up workflows by chaining operations together using the pipe operator
- Sorting rows based on values in one or more columns
- Modifying a data frame by either adding new columns or modifying existing ones
- Summarizing the values in one or more columns
- Faceting of ggplot2 visualizations
Exercises
These exercises will test your ability to manipulate the METABRIC data set by changing the values of existing columns or adding new columns by computing new variables from existing ones.
We are expecting you to use the 5 main dplyr ‘verb’ functions: select(), filter(), arrange(), mutate() and summarize(). Please use the pipe operator, %>%, in cases where more than one operation is required to achieve the outcome requested.
library(tidyverse)
metabric <- read_csv("data/metabric_clinical_and_expression_data.csv")1. Investigate the subset of long-surviving breast cancer patients that didn’t receive chemo or radiotherapy
First obtain the subset of patients that received neither chemotherapy or radiotherapy and survived for more than 20 years.
Answer
patients <- metabric %>%
mutate(Survival_time_years = Survival_time / 12) %>%
filter(Chemotherapy == "NO" &
Radiotherapy == "NO" &
Survival_time_years > 20)
patients## # A tibble: 70 × 33
## Patient_ID Cohort Age_at_diagnosis Survival_time Survival_status Vital_status
## <chr> <dbl> <dbl> <dbl> <chr> <chr>
## 1 MB-0270 1 30.0 337. LIVING Living
## 2 MB-2517 2 60.3 268. DECEASED Died of Oth…
## 3 MB-2610 2 59.0 281. LIVING Living
## 4 MB-2763 2 70.3 276. LIVING Living
## 5 MB-2770 2 60.3 274. LIVING Living
## 6 MB-2793 2 67.6 267. LIVING Living
## 7 MB-2797 2 50.1 253. LIVING Living
## 8 MB-2819 2 62.0 271. LIVING Living
## 9 MB-2838 2 67.1 270. LIVING Living
## 10 MB-2840 2 60.5 269 LIVING Living
## # ℹ 60 more rows
## # ℹ 27 more variables: Chemotherapy <chr>, Radiotherapy <chr>,
## # Tumour_size <dbl>, Tumour_stage <dbl>, Neoplasm_histologic_grade <dbl>,
## # Lymph_nodes_examined_positive <dbl>, Lymph_node_status <dbl>,
## # Cancer_type <chr>, ER_status <chr>, PR_status <chr>, HER2_status <chr>,
## # HER2_status_measured_by_SNP6 <chr>, PAM50 <chr>, `3-gene_classifier` <chr>,
## # Nottingham_prognostic_index <dbl>, Cellularity <chr>, …Now look at the breakdown of these patients in terms of ER status. Count the numbers of ER positive and ER negative patients in this subset. Calculate the proportion that are ER positive.
Answer
table(patients$ER_status)##
## Negative Positive
## 8 62mean(patients$ER_status == "Positive")## [1] 0.8857143What does this tell us? Calculate the proportion of ER positive patients in the whole cohort by way of a comparison.
Answer
table(metabric$ER_status)##
## Negative Positive
## 445 1459mean(metabric$ER_status == "Positive")## [1] 0.76628152. Which patients have spent the largest proportion of their lives dealing with breast cancer? Note that the current Survival_time column is months. Present the final results as a table showing Patient_ID, Age, Age_at_diagnosis, Survival_time_in_years and Proportion_of_life_since_diagnosis. Round all numbers in the final table to 2 decimal places and sort it so that the patient with the highest proportion is at the top.
Hint
- Create a column of Survival_time_years
- Add Survival_time_years it to Age_at_diagnosis to get final Age
- Calculate the proportion: Survival_time_years/Age
Answer
metabric %>%
mutate(Survival_time_years = Survival_time / 12) %>%
mutate(Age = Age_at_diagnosis + Survival_time_years) %>%
mutate(Proportion_of_life_since_diagnosis = Survival_time_years / Age) %>%
select(Patient_ID,
Age,
Age_at_diagnosis,
Survival_time_years,
Proportion_of_life_since_diagnosis) %>%
arrange(desc(Proportion_of_life_since_diagnosis)) %>%
mutate(across(where(is.numeric), ~round(.x, digits = 2)))## # A tibble: 1,904 × 5
## Patient_ID Age Age_at_diagnosis Survival_time_years Proportion_of_life_si…¹
## <chr> <dbl> <dbl> <dbl> <dbl>
## 1 MB-0270 58.1 30.0 28.1 0.48
## 2 MB-4332 60.0 34.4 25.7 0.43
## 3 MB-6062 55.1 31.6 23.5 0.43
## 4 MB-2753 53.8 31.0 22.9 0.42
## 5 MB-2957 52.9 31.0 21.8 0.41
## 6 MB-2823 58.5 36.0 22.5 0.38
## 7 MB-2779 59.7 36.8 22.9 0.38
## 8 MB-2853 58.7 36.2 22.5 0.38
## 9 MB-4881 59.9 37.0 22.8 0.38
## 10 MB-4375 64.3 40.0 24.3 0.38
## # ℹ 1,894 more rows
## # ℹ abbreviated name: ¹Proportion_of_life_since_diagnosis