9  Seahorse Analysis

print("WebR is working!")
_webr_editor_1 = Object {code: null, options: Object, indicator: Ke}

9.1 Data

The experimental data is derived from a paper we published a couple of years ago. Please check-out the paper where this data was published at:

Janssen et al. (2021) Scientific Reports

Exercise x1

Look up in the paper: - what type of cells have been used? - what compounds were injected during the assays?

For most of the experiments PBMCs have been used. For some experiments RAW264.7 cells were used or monocytes. The dataset we use in this tutorial was an experiment with PBMCs.

We used injections of FCCP, AM/rot, and monensin. Also Hoechst was included, but this was only used for nuclei staining

Let’s load the data first. We use the revive_xfplate() function form the seahtrue package and the excel file that is included with the package.


library(tidyverse)

seahtrue_output_donor_A <-
system.file("extdata",
"20191219_SciRep_PBMCs_donor_A.xlsx",
package = "seahtrue") %>%
revive_xfplate()

_webr_editor_2 = Object {code: null, options: Object, indicator: Ke}
Exercise x2

Explore the data that we just loaded.


seahtrue_output_donor_A %>%
glimpse()

_webr_editor_3 = Object {code: null, options: Object, indicator: Ke}

What is the structure of seahtrue_output_donor_A?


seahtrue_output_donor_A %>%
pluck("rate_data", 1) %>%
glimpse()
_webr_editor_4 = Object {code: null, options: Object, indicator: Ke}

seahtrue_output_donor_A %>%
pluck("raw_data", 1) %>%
head()
_webr_editor_5 = Object {code: null, options: Object, indicator: Ke}

And what is the structure of the rate_data and raw_data tibbles?

The data is a nested tibble. The nested rate_data tibble contains the ocr and ecar data, with for each well the measurements and group annotation and the OCR_wave_bc and ECAR_wave_bc data columns “…._bc” in the columns names means the data is allready background corrected (bc)). The raw_data contains the O2 and pH readings in the experiment, those are used to calculate the OCR and ECAR).

9.2 Plate map

To dive a bit deeper into a single seahtrue experiment, we will first generate an overview of what the experimental set-up was.

Explore the experimental set-up. Remember that a Seahorse experiment is a 96-well plate based experiment, so we will see 96 wells in the plate map.


seahtrue_output_donor_A %>%
pluck("raw_data", 1) %>%
seahtrue::sketch_plate(reorder_legend = TRUE)

_webr_editor_6 = Object {code: null, options: Object, indicator: Ke}
Exercise x3

In the plate map, you see that the corners are colored dark brown The background wells do not contain cells. What is the purpose of these background wells?

There is drift in O2 and pH signals even when there is no sample (cells) in a well. This signal should be substracted from the wells that do contain samples.

Exercise x4

What do the numbers 50.000 to 300.000 in the legend of the plate map?

These are the number of cells that have been plated in each well, so 50.000 cells/well to 300.000 cells/well

Exercise x5

What is the reason why we performed an experiment like this?

In this way, we can check the linearity of the signals and determine the best plating density for our following experiments where we want to compare different conditions or interventions.

9.3 Oxygen consumption rate (OCR)

Now we will plot the OCRs for this experiment.


#before running this, make the seahtrue_output_donor_A
#is assigned in the above code

seahtrue_output_donor_A %>%
purrr::pluck("rate_data", 1) |>
sketch_rate(
param = "OCR",
normalize = FALSE,
take_group_mean = TRUE,
reorder_legend = TRUE
)

_webr_editor_7 = Object {code: null, options: Object, indicator: Ke}
Exercise x6

If we assume that the basal OCR should not be lower than 20 pmol/min (below this value we reach the detection limit), which cell densities are suitable for a Seahorse experiment in this case.

from 150.000 cells/well

Exercise x7

Now check out the parameters that we can set for the sketch_rate() function by using the below code:


?seahtrue::sketch_rate()
_webr_editor_8 = Object {code: null, options: Object, indicator: Ke}

plot the ECAR using this sketch_rate() functions. What do you observe after monensin is injected in the ECAR plot? Look up online how monensin acts on cells why it is doing this with ECAR. What did we observe in OCR? Can you explain that?


seahtrue_output_donor_A %>%
purrr::pluck("rate_data", 1) #.....
#.....
#.....
_webr_editor_9 = Object {code: null, options: Object, indicator: Ke}

seahtrue_output_donor_A %>%
purrr::pluck("rate_data", 1) |>
sketch_rate(
param = "ECAR",
normalize = FALSE,
take_group_mean = FALSE,
reorder_legend = TRUE
)
_webr_editor_10 = Object {code: null, options: Object, indicator: Ke}

We see that after monensin injection, ECAR increases. This is because glycolysis is increased because cell membrane transporters are activated that increase ATP demand. OCR does not change much, because the mitochondrial activity is allready blocked by AM/rot

9.4 Plotting basal and maximal respiration

Pluck the injection_info table from the seahtrue_output_donor_A dataset to see what how the injections were defined in the experimental set-up before running the seahorse.


seahtrue_output_donor_A %>%
pluck("injection_info",1)

_webr_editor_11 = Object {code: null, options: Object, indicator: Ke}

This is not a typical mito-stress test experiment, where we inject oligomycin, FCCP and antimycinA/rotenone sequentially. Instead we inject only FCCP and antimycinA/rotenone.

To get the maximal and basal respiration out of the ocr rate table, we need to do some calculations. We first make some assumptions and defintions:

We call each interval between two injections or between start and an injection or between an injection and the end a phase

Each phase has a unique name that is named after the injection that was last. The first phase after the start is called init_ocr and we also typically have the phases om_ocr, fccp_ocr and amrot_ocr. Phases are marked with either _ocr or _ecar, because these are distinct parameters.

To calculate respiration parameters, like basal respiration (= basal_ocr), we define the following:

  • basal_ocr = init_ocr - amrot_ocr.
  • max_ocr = fccp_ocr - amrot_ocr
  • spare_ocr = fccp_ocr - init_ocr
  • proton_leak = om_ocr = amrot_ocr
  • atp_linked = init_ocr - om_ocr

We also use indices to have relative parameters:

  • spare_ocr_index = (spare_ocr / basal_ocr)*100
  • basal_ocr_index = (basal_ocr / max_ocr)*100
  • leak_index = (proton_leak / basal_ocr)*100
  • coupling_index = (atp_linked / basal_ocr)*100) %>%

Another important assumption is that we are not using average values to represent each phase, but intead we use a specific measurement. The reason for this is that we assume that for all phases, except FCCP, three measurements are needed in time to get to steady-state. For FCCP injection, we assume that it reaches steady-state fast, or at least its maximal ocr, so we take the first measurement after injection as the measurement representing the FCCP phase.

Let’s now put that into R code. We call the type of experiment we did in this dataset a maximal capacity (maxcap) test.

We also injected monensin, which can maximize ECAR, but we don’t need it for OCR calculations.


# first define which timepoints are what
param_set_maxcap_ocr <- c(init_ocr = "m3",
fccp_ocr = "m4",
amrot_ocr = "m9",
mon_ocr = "m12"
)

#next do some pivoting, renaming and selecting
seahtrue_output_donor_A %>%
pluck("rate_data",1) %>%
select(well,
measurement,
group,
ocr = OCR_wave_bc) %>%
pivot_wider(names_from = measurement,
names_prefix = "m",
values_from = ocr) %>%
rename(all_of(param_set_maxcap_ocr)) %>%
select(contains(c("wel"," group", "ocr")))
_webr_editor_12 = Object {code: null, options: Object, indicator: Ke}

Now for each well we have the parameters related to the phases that we defined in the parameter set.

Next we want to calculate the respiration parameters and indices.


seahtrue_output_donor_A %>%
pluck("rate_data",1) %>%
select(well,
measurement,
group,
ocr = OCR_wave_bc) %>%
pivot_wider(names_from = measurement,
names_prefix = "m",
values_from = ocr) %>%
rename(all_of(param_set_maxcap_ocr)) %>%
select(contains(c("well","group", "ocr"))) %>%
#piped here to the mutate
mutate(non_mito_ocr = amrot_ocr,
basal_ocr = init_ocr - non_mito_ocr,
max_ocr = fccp_ocr - amrot_ocr,
spare_ocr = max_ocr - basal_ocr,
spare_ocr_index = (spare_ocr / max_ocr)*100,
basal_ocr_index = (basal_ocr/max_ocr)*100)
_webr_editor_13 = Object {code: null, options: Object, indicator: Ke}

With this data, we can plot our typical basal and maximal bar/scatter plots that we see in papers, presentations and theses.


webr::install("ggdist")

param_set_maxcap_ocr <- c(init_ocr = "m3",
fccp_ocr = "m4",
amrot_ocr = "m9",
mon_ocr = "m12"
)

seahtrue_output_donor_A %>%
pluck("rate_data",1) %>%
select(well, measurement, group, ocr = OCR_wave_bc) %>%
pivot_wider(names_from = measurement,
names_prefix = "m",
values_from = ocr) %>%
rename(all_of(param_set_maxcap_ocr)) %>%
select(contains(c("well","group", "ocr"))) %>%
#piped here to the mutate
mutate(non_mito_ocr = amrot_ocr,
basal_ocr = init_ocr - non_mito_ocr,
max_ocr = fccp_ocr - amrot_ocr,
spare_ocr = max_ocr - basal_ocr,
spare_ocr_index = (spare_ocr / max_ocr)*100,
basal_ocr_index = (basal_ocr/max_ocr)*100) %>%
filter(group %in% c("150.000", "250.000")) %>%
ggplot(aes(x = group, y = max_ocr, color = group))+
geom_bar(data = . %>%
summarize(median_max_ocr = median(max_ocr),
.by = group),
mapping = aes(
x = forcats::fct_reorder(
group,
parse_number(group)),
y = median_max_ocr,
fill = group),
stat="identity",
_webr_editor_14 = Object {code: null, options: Object, indicator: Ke}
Exercise 8

The plot above shows maximal ocr. Now make your own plot with 1) basal_ocr and 2) all groups except background. Make sure to order the group legend tidyly and have reable x-axis labels.


#use the code block above and change that code

_webr_editor_15 = Object {code: null, options: Object, indicator: Ke}

#solution 1 and 2 together
#lines with changes marked with #VB

webr::install("ggdist")

param_set_maxcap_ocr <- c(init_ocr = "m3",
fccp_ocr = "m4",
amrot_ocr = "m9",
mon_ocr = "m12"
)

group_order <- seahtrue_output_donor_A %>% #VB
pluck("rate_data",1) %>%
pull(group) %>% unique()

seahtrue_output_donor_A %>%
pluck("rate_data",1) %>%
select(well, measurement, group, ocr = OCR_wave_bc) %>%
pivot_wider(names_from = measurement,
names_prefix = "m",
values_from = ocr) %>%
rename(all_of(param_set_maxcap_ocr)) %>%
select(contains(c("well","group", "ocr"))) %>%
#piped here to the mutate
mutate(non_mito_ocr = amrot_ocr,
basal_ocr = init_ocr - non_mito_ocr,
max_ocr = fccp_ocr - amrot_ocr,
spare_ocr = max_ocr - basal_ocr,
spare_ocr_index = (spare_ocr / max_ocr)*100,
basal_ocr_index = (basal_ocr/max_ocr)*100) %>%
filter(group!= c("Background")) %>% #VB
ggplot(aes(x = group, y = basal_ocr, color = group))+ #VB
geom_bar(data = . %>%
summarize(
median_basal_ocr = median(basal_ocr), #VB
_webr_editor_16 = Object {code: null, options: Object, indicator: Ke}

9.5 Functional omics with Seahorse analysis - ATP interpretations

Next, we will continue with the transformation of OCR and ECAR to ATP units. The calculations are based on the Mookerjee et al paper in J Biol Chem:

Mookerjee et al. (2017) J Biol Chem

The calculations are performed under the hood in the calculate_space() function of the seahtrue package.

For these excercises we will work with a similar dataset as the seahtrue_output_donor_A, but now from a different individual. This is seahtrue_output_donor_C. Since there is a difference in the Seahorse tracjectory this is a nice example to learn how to interpret seahorse experiments.

First, load the new data.


library(tidyverse)

# raw path on github
url <- paste0(
"https://raw.githubusercontent.com/vcjdeboer/seahtrue/develop-gerwin/inst/extdata/",
"20200110_SciRep_PBMCs_donor_C.xlsx"
)

# Save it locally
download.file(
url,
destfile = "20200110_SciRep_PBMCs_donor_C.xlsx",
mode = "wb")

# Load it into seahtrue
seahtrue_output_donor_C <-
"20200110_SciRep_PBMCs_donor_C.xlsx" %>%
seahtrue::revive_xfplate()

_webr_editor_17 = Object {code: null, options: Object, indicator: Ke}

seahtrue_output_donor_C
_webr_editor_18 = Object {code: null, options: Object, indicator: Ke}

Now lets calculate the atp space for both donors and combine the data for plotting:


df_space_C <-
#first pluck the rate_data
seahtrue_output_donor_C %>%
pluck("rate_data", 1) %>%
#only take the 200.000 cells/well group
filter(group %in% c("200.000")) %>%

#next call the calculate_space function
seahtrue::calculate_space(
# provide the injection info
param_set_ocr = c(init_ocr = "m3",
fccp_ocr = "m4",
amrot_ocr = "m9",
mon_ocr = "m12"),
param_set_ecar = c(init_ecar = "m3",
fccp_ecar = "m4",
amrot_ecar = "m9",
mon_ecar = "m12")
) %>%
mutate(group = "donor_C")

df_space_A <-
#first pluck the rate_data
seahtrue_output_donor_A %>%
pluck("rate_data", 1) %>%
#only take the 200.000 cells/well group
filter(group %in% c("200.000")) %>%

#next call the calculate_space function
seahtrue::calculate_space(
# provide the injection info
_webr_editor_19 = Object {code: null, options: Object, indicator: Ke}
Exercise x9

Verify that the correct measurements were put as arguments for the calculate_space() function

Since we do not take the average of the three measurements for each phase, because we think that the cells reach a steady state after injection after a certain time interval, we take one measurement for each phase. In this experiment, we injected fccp, AM/rot and monension consecutively. The fccp injection typically give the highest value allready after one measurement so we take measurement 4 (m4). for the others we take the third measurement in the phase (m3, m9 and m12)

param_set_ocr = c(init_ocr = “m3”, fccp_ocr = “m4”, amrot_ocr = “m9”, mon_ocr = “m12”),

Now plot the trajectory and compare the two trajectories.


my_palette<- c("#1f78b4","#33a02c")

df_space_total %>%
seahtrue::plot_bioenergetic_trajectory(palette = my_palette)


_webr_editor_20 = Object {code: null, options: Object, indicator: Ke}
Exercise x10

Now that we generated the trajectory through ATP space, explain what you observe. How does the trajectory change for each individual? What is the biggest difference between donor A and C

please discuss with one of your lecturers

Next, we will plot the ATP space. For this we will use the plot_bioenergetic_space() function. The space plot shows the total allowed ATP generation at the condition that the cells are in. The dot is the position in the space in the initial (basal) phase.


my_palette<- c("#1f78b4","#33a02c")

df_space_total %>%
seahtrue::plot_bioenergetic_space(palette = my_palette)


_webr_editor_21 = Object {code: null, options: Object, indicator: Ke}
Exercise x11

Explain the difference between the two spaces between donor A and B. Also, what does the dashed line represent?

please discuss with one of your lecturers

The following exercise is a bit more advanced, it requires some R skills, but gives insights in the relevance of the bioenergetic space.

Exercise ADVANCED x12

We have another dataset where we isolated T cells from an individual and stimulated the cells with CD3/CD28 (= IC). Also, we injected a number of chemical compounds. Please focus for now on the “IC + 1.5 uM BAM15” and “1.5 uM BAM15” group. Instead of the uncoupler FCCP we now used BAM15 instead which is more gentle on the cell membrane.

Keep the fccp_ocr and fccp_ecar nomenclature for your param_sets. The function does not recognize bam15_ocr or bam15_ecar.

Find out which measurements you should use to calculate the space.

Filter out well C10. This is a technical outlier.

Use “agilent” as the atp conversion model.

Plot the ATP space for the activated and non-activated T cells.

  • What do you observe?
  • What does it mean when the basal J_ATP (the dot in the space) is more closer to the dashed line?
#load the data
# raw path on github
url <- paste0(
"https://raw.githubusercontent.com/vcjdeboer/seahtrue/develop-gerwin/inst/extdata/",
"20251102_Tcell_Activation_VP.xlsx"
)

# Save it locally
download.file(
url,
destfile = "20251102_Tcell_Activation_VP.xlsx",
mode = "wb")

# Load it into seahtrue
df_t_cell_activation <-
"20251102_Tcell_Activation_VP.xlsx" %>%
seahtrue::revive_xfplate()
_webr_editor_22 = Object {code: null, options: Object, indicator: Ke}

We want to see the difference between activate T cells and not activated, so we should take the measurement 9 where the cells were activated or not (at measurement three both groups are still the same!).


#use this param sets
my_param_set_ocr = c(init_ocr = "m9", om_ocr = "m12",
fccp_ocr = "m15", amrot_ocr = "m18")
my_param_set_ecar = c(init_ecar = "m9", om_ecar = "m12",
fccp_ecar = "m15", amrot_ecar = "m18")

my_rate <- df_t_cell_activation %>%
purrr::pluck("rate_data", 1) %>%
dplyr::filter(well != "C10") %>%
dplyr::filter(stringr::str_detect(group, "1.5 uM BAM15"))

df_space_tcell <- calculate_space(rate = my_rate,
param_set_ocr = my_param_set_ocr,
param_set_ecar = my_param_set_ecar,
conversion_model = "agilent"
)


my_palette <- c("#e41a1c", "#984ea3")

df_space_tcell %>%
seahtrue::plot_bioenergetic_space(palette = my_palette)
_webr_editor_23 = Object {code: null, options: Object, indicator: Ke}
Downloading package: seahtrue