shot_noise.qmd

---
title: "Shot Noise Analysis"
author: "Kantnerova et al. (2023)"
date: "`r format(Sys.Date(), '%d %b %Y')`"
number-sections: false
toc: true
toc-depth: 2
fig-width: 8
fig-height: 6
df-print: tibble
embed-resources: true
format:
  html:
    code-tools: true
    code-fold: false
    code-summary: "Show the code"
    toc-float: true
crossref:
  fig-prefix: Fig.
  tbl-prefix: Table
  ref-hyperlink: true
editor: source
editor_options:
  chunk_output_type: console
---

# Setup

Using `r R.version.string` , tidyverse version `r packageVersion("tidyverse")`, and isoorbi version `r packageVersion("isoorbi")`.

```{r}
#| label: setup
#| message: false

# load packages
library(isoorbi) # for orbitrap function
library(dplyr) # for mutating data frames
library(forcats) # for recoding factors
library(ggplot2) # for data visualization
library(cowplot) # arrange multipanel plots
```

```{r, include=FALSE}
# global knitting options for code rendering
knitr::opts_chunk$set(
  collapse = TRUE, comment = "#>",
  dev = c("png", "pdf"), fig.keep = "all",
  dev.args = list(pdf = list(encoding = "WinAnsi", useDingbats = FALSE)),
  fig.path = sprintf("plots/%s_", gsub("\\.rmarkdown", "", knitr::current_input()))
)
```

# Trifluoroacetate (TFA)

## Data

```{r}
# load, process, and export tfa data
data_tfa <-
  "data/shot_noise/tfa" |>
  orbi_find_isox() |>
  orbi_read_isox() |>
  orbi_simplify_isox() |>
  orbi_flag_satellite_peaks() |>
  orbi_flag_weak_isotopocules(min_percent = 10) |>
  orbi_flag_outliers(agc_fold_cutoff = 2) |>
  orbi_define_basepeak("M0")

shot_noise_tfa <- 
  data_tfa |>
  orbi_analyze_shot_noise() |>
  orbi_export_data_to_excel("output/shot_noise_tfa.xlsx")
```

## Bonus Figure: TFA satellite peaks

```{r}
#| label: fig-bonus-tfa-satellite-peaks
#| fig-cap: Satellite peaks
#| fig-width: 15
#| fig-height: 5
#| warning: false
data_tfa |> orbi_plot_satellite_peaks()
```

## Bonus Figure: TFA isotopocule coverage

```{r}
#| label: fig-bonus-tfa-isotopocule-coverage
#| fig-cap: Isotopocule coverage
#| fig-width: 15
#| fig-height: 3
#| warning: false
data_tfa |> orbi_plot_isotopocule_coverage()
```

## Bonus Figure: TFA ratios

```{r}
#| label: fig-bonus-tfa-ratios
#| fig-cap: Isotopocule ratios vs. M0
#| fig-width: 15
#| fig-height: 5
#| warning: false
data_tfa |> orbi_plot_raw_data(y = ratio)
```

## Figure 10: TFA shotnoise vs counts/time

```{r}
#| label: fig-10-shotnoise-tfa
#| fig-cap: Shot noise vs counts and vs time.
#| fig-width: 10
#| fig-height: 5
#| warning: false

# individual plots
tfa_vs_ions <-
  shot_noise_tfa |>
  orbi_filter_isox(filename = "TFA_M0_1uscan_15kRes") |>
  orbi_plot_shot_noise(x = "n_effective_ions") +
  labs(title = "vs counts") +
  theme(plot.title = element_text(hjust = 0.5))

tfa_vs_time <-
  shot_noise_tfa |>
  orbi_filter_isox(filename = "TFA_M0_1uscan_15kRes") |>
  orbi_plot_shot_noise(permil_target = 1) +
  labs(title = "vs time") +
  theme(plot.title = element_text(hjust = 0.5))

# combine
plot_grid(
  tfa_vs_ions + theme(legend.position = c(0.85, 0.74)), 
  tfa_vs_time + theme(legend.position = "none"),
  align = "h", nrow = 1, axis = "tb"
)
```

## Extended Data Figure 6: TFA shotnoise at different IT

```{r}
#| label: fig-ext-data-6-shotnoise-tfa-ITs
#| fig-cap: Shot noise at different ITs.
#| fig-width: 10
#| fig-height: 5
#| warning: false
shot_noise_tfa |>
  orbi_filter_isox(filename = c("TFA_1E4_AGC", "TFA_M0_1uscan_15kRes")) |>
  # change color legend
  mutate(
    IT_info = filename |>
      fct_recode(
        "IT = 0.03 ms" = "TFA_1E4_AGC",
        "IT = 0.75 ms" = "TFA_M0_1uscan_15kRes"
      )
  ) |>
  orbi_plot_shot_noise(color = "IT_info") +
  # wrap by the ratio label
  facet_wrap(~ratio_label) +
  theme(legend.position = c(0.41, 0.75))
```

## Extended Data Figure 7: TFA shotnoise at different resolutions

```{r}
#| label: fig-ext-data-7-shotnoise-tfa-resolution
#| fig-cap: Shot noise at different resolutions.
#| fig-width: 10
#| fig-height: 5
#| warning: false
shot_noise_tfa |>
  orbi_filter_isox(filename = c("TFA_M0_1uscan_120kRes", "TFA_M0_1uscan_15kRes")) |>
  # change color legend
  mutate(
    res_info = filename |>
      fct_recode(
        "res. = 120k" = "TFA_M0_1uscan_120kRes",
        "res. = 15k" = "TFA_M0_1uscan_15kRes"
      )
  ) |>
  orbi_plot_shot_noise(color = "res_info") +
  # wrap by the ratio label
  facet_wrap(~ratio_label) +
  theme(legend.position = c(0.42, 0.75))
```

# Model Peptide (MRFA)

## Data

```{r}
# load, process, and export mrfa data
data_aas <- 
  "data/shot_noise/mrfa" |>
  orbi_find_isox() |>
  orbi_read_isox() |>
  orbi_simplify_isox() |>
  orbi_flag_satellite_peaks() |>
  orbi_flag_weak_isotopocules(min_percent = 10) |> 
  orbi_flag_outliers(agc_fold_cutoff = 2) |> 
  orbi_define_basepeak("M0") |>
  mutate(
    compound = compound |>
      fct_recode(
        "Alanine" = "Ala_imm",
        "Arginine" = "Arg_imm",
        "Methionine" = "Met_imm",
        "Phenylalanine" = "Phe_imm"
      )
  )

shot_noise_aas <-
  data_aas |>
  orbi_analyze_shot_noise() |> 
  orbi_export_data_to_excel("output/shot_noise_MRFA.xlsx") 
```

## Table

```{r}
# example of first few rows
shot_noise_aas |>
  arrange(compound, isotopocule, scan.no) |>
  select(compound, scan.no, time.min, isotopocule,
         ratio, ratio_rel_se.permil, shot_noise.permil) |>
  head(10) |>
  knitr::kable()
```

## Bonus Figure: Amino Acids satellite peaks

```{r}
#| label: fig-bonus-mrfa-satellite-peaks
#| fig-cap: Satellite peaks
#| fig-width: 10
#| fig-height: 9
#| warning: false
data_aas |> orbi_plot_satellite_peaks()
```

## Bonus Figure: Amino Acids isotopocule coverage

```{r}
#| label: fig-bonus-mrfa-isotopocule-coverage
#| fig-cap: Isotopocule coverage
#| fig-width: 10
#| fig-height: 9
#| warning: false
data_aas |> orbi_plot_isotopocule_coverage()
```

## Bonus Figure: Amino Acids ratios

```{r}
#| label: fig-bonus-mrfa-ratios
#| fig-cap: Isotopocule ratios vs. M0
#| fig-width: 10
#| fig-height: 9
#| warning: false
data_aas |> orbi_plot_raw_data(y = ratio)
```

## Extended Data Figure 8: Amino Acids Shot Noise

```{r}
#| label: fig-ext-data-8-shotnoise-mrfa
#| fig-cap: amino acids shotnoise
#| warning: false
#| fig-width: 8
#| fig-height: 10
shot_noise_aas |>
  filter(compound == "Methionine" | isotopocule != "2H") |>
  orbi_plot_shot_noise()
```