dilution_series.qmd

# Regressions for many lipids

```{r}
#| label: setup
#| echo: false
#| message: false
source("_common.R")
library(knitr)

```

## Libraries

```{r}
#| label: load-lib
#| echo: true
#| message: false
library(here)
library(tidyverse)
library(broom)
library(ggpmisc)

```

## Overview lm and 'broom' package

```{r}
#| label: intro
#| echo: true
#| message: false

d_test <- tibble(InjVol = c(0,0.2, 0.4, 0.6, 0.8,1),
                  Response  = c(12, 23, 34,44, 89, 101)) 
                 
d_test

# Linear model
model <- lm(formula = Response ~ InjVol, data = d_test)

# Get result summary
summary(model)

# Get r^2 only
summary(model)$r.squared

# Using broom functions to summarize model results into a table
broom::glance(model) 
broom::tidy(model) 
 
```

## Import Datasets

```{r}
#| label: read-data
#| echo: true
#| message: false

d_orig <- read_csv(here("data/Testdata_Lipidomics_flat_wide_annotated_V1.csv"))

```

## Prepare Data

```{r select-rename-columns}
#| label: prepare-data
#| echo: true
#| message: false

# Convert to long format

# Convert to long format
d_long <- d_orig |> 
  pivot_longer(cols = -DataFileName:-InjVol, 
               names_to = "Lipid" , 
               values_to = "Area")

# Get a table with RQCs only and sort by Lipid
d_rqc <- d_long |> 
  filter(QCtype == "RQC") |> 
  arrange(Lipid)

```

## Run regression for each lipid

In this example a logistic regression is used. The output of `glm()` is converted to a tidy table using the `broom::tidy()` function.

```{r filter-tables}
#| label: reg-analysis
#| echo: true
#| message: false
#| tbl-cap: "Model results"
#| tbl-colwidths: [60,40]



model <- as.formula("Area ~ InjVol")

d_res <- d_rqc %>%
  group_by(Lipid) %>%
  nest() %>%
  mutate(
    models = map(data, function(x) lm(model, data = x)), 
    #mandel = map(data, \(x) DCVtestkit::calculate_mandel(x, "InjVol", "Area")),
    #ppa = map(data, \(x) DCVtestkit::calculate_pra_linear(x, "InjVol", "Area")),
    tidy = map(models, function(x) broom::glance(x))) |> 
  unnest(c(tidy)) |> 
  dplyr::select(-data, -models)

# Fix DCVtestkit::calculate_pra_linear currently returning a list instead vector
# d_res$ppa <- unlist(d_res$ppa)

```

The results contain the combined estimates, errors, and *P* values for each term for each lipid species.

```{r}
#| label: tbl_results
#| echo: false
kable(head(d_res, n = 10),digits = 6)
```