### Install ggforestplot from GitHub

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/ggforestplot.html

Instructions for installing the ggforestplot package directly from GitHub. It also mentions how to install with vignette building enabled.

```R
# install.packages("devtools")
devtools::install_github("NightingaleHealth/ggforestplot")

# To install with vignettes:
devtools::install_github("NightingaleHealth/ggforestplot", build_vignettes = TRUE)
```

--------------------------------

### Install ggforestplot from GitHub

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/README.md

Installs the ggforestplot R package directly from its GitHub repository. Requires the 'devtools' package to be installed first.

```R
# install.packages("devtools")
devtools::install_github("NightingaleHealth/ggforestplot")
```

--------------------------------

### Install and Load Tidyverse Packages

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

Installs the complete tidyverse collection of R packages, which provides a suite of tools for data manipulation, exploration, and visualization. It then loads these packages into the current R session.

```R
# Install the tidyverse packages
install.packages("tidyverse")

# Attach the tidyverse packages into the current session
library(tidyverse)
```

--------------------------------

### Install ggforestplot from GitHub (R)

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/index.html

Instructions for installing the ggforestplot R package directly from GitHub using `devtools`. It mentions the prerequisite of having `devtools` installed and the option to build vignettes.

```R
# install.packages("devtools")
devtools::install_github("NightingaleHealth/ggforestplot")
```

--------------------------------

### R Example: Metabolic Data Analysis

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/discovery_regression.html

An R code example demonstrating the use of `discovery_regression` with simulated metabolic data. It shows data preparation steps using `dplyr` and `tidyr`, including selecting, transforming, and gathering data into a long format before calling `discovery_regression`.

```R
library(magrittr)

# Linear Regression Example

# We will use the simulated demo data that come with the package,
# ggforestplot::df_demo_metabolic_data

# Extract the names of the NMR biomarkers for discovery analysis
nmr_biomarkers <- dplyr::intersect(
  ggforestplot::df_NG_biomarker_metadata$machine_readable_name,
  colnames(df_demo_metabolic_data)
)

# Select only variables to be used for the model and collapse to a long
# format
df_long <-
  df_demo_metabolic_data %>%
  # Select only model variables
  dplyr::select(tidyselect::all_of(nmr_biomarkers), gender, BMI) %>%
  # Log-transform and scale biomarkers
  dplyr::mutate_at(
    .vars = dplyr::vars(tidyselect::all_of(nmr_biomarkers)),
    .funs = ~ .x %>% log1p() %>% scale() %>% as.numeric()
  ) %>%
  # Collapse to a long format
  tidyr::gather(
    key = machine_readable_name,
    value = biomarkervalue,
    tidyselect::all_of(nmr_biomarkers)
  )

df_assoc_per_biomarker <-
  discovery_regression(
    df_long = df_long,
```

--------------------------------

### GGForestPlot: R Code Examples

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/plot_all_NG_biomarkers.html

Illustrative R code examples demonstrating the usage of the `plot_all_NG_biomarkers` function for visualizing biomarker associations, including data preparation and custom layout creation.

```R
if (FALSE) {
# Join the built-in association demo dataset with a variable that contains
# the machine readable names of Nightingale biomarkers. (Note: if you
# have built your association data frame using the Nightingale CSV result file,
# then your data frame should already contain machine readable names.)
df <-
  df_linear_associations %>%
  left_join(
    select(df_NG_biomarker_metadata,
      name,
      machine_readable_name
    ),
    by = "name"
  )

# Print effect sizes for Nightingale biomarkers in a 2-page pdf
plot_all_NG_biomarkers(
  df = df,
  machine_readable_name = machine_readable_name,
  # Notice that when name is not defined explicitly, names from
  # df_NG_biomarker_metadata are used
  estimate = beta,
  se = se,
  pvalue = pvalue,
  colour = trait,
  filename = "biomarker_linear_associations.pdf",
  xlab = "1-SD increment in BMI
per 1-SD increment in biomarker concentration",
  layout = "2016"
)

# Custom layout can also be provided
layout <- df_NG_biomarker_metadata %>%
  dplyr::filter(
    .data$group == "Fatty acids",
    .data$machine_readable_name %in% df$machine_readable_name
  ) %>%
  dplyr::mutate(
    group_custom = .data$subgroup,
    column = dplyr::case_when(
      .data$group_custom == "Fatty acids" ~ 1,
      .data$group_custom == "Fatty acid ratios" ~ 2
    ),
    page = 1
  ) %>%
  dplyr::select(
    .data$machine_readable_name,
    .data$group_custom,
    .data$column,
    .data$page
  )

plot_all_NG_biomarkers(
  df = df,
  machine_readable_name = machine_readable_name,
  # Notice that when name is not defined explicitly, names from
  # df_NG_biomarker_metadata are used
  estimate = beta,
  se = se,
  pvalue = pvalue,
  colour = trait,
  xlab = "1-SD increment in BMI
per 1-SD increment in biomarker concentration",
  layout = layout
)

# log odds for type 2 diabetes
df <-
  df_logodds_associations %>%
  left_join(
    select(df_NG_biomarker_metadata,
      name,
      machine_readable_name
    ),
    by = "name"
  ) %>%
  # Set the study variable to a factor to preserve order of appearance
  # Set class to factor to set order of display.
  dplyr::mutate(
    study = factor(
      study,
      levels = c("Meta-analysis", "NFBC-1997", "DILGOM", "FINRISK-1997", "YFS")
    )
  )

# Print effect sizes for Nightingale biomarkers in a 2-page pdf
plot_all_NG_biomarkers(
  df = df,
  machine_readable_name = machine_readable_name,
  # Notice that when name is not defined explicitly, names from
  # df_NG_biomarker_metadata are used
  estimate = beta,
  se = se,
  pvalue = pvalue,
  colour = study,
  logodds = TRUE,
  filename = "biomarker_t2d_associations.pdf",
  xlab = "Odds ratio for incident type 2 diabetes (95% CI)
per 1−SD increment in metabolite concentration",
  layout = "2016",
  # Restrict limits as some studies are very weak and they take over the
  # overall range.
  xlims = c(0.5, 3.2)
)
}
```

--------------------------------

### Inspect Example Dataset Structure

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/ggforestplot.html

Shows how to view the structure and first few rows of the `df_linear_associations` dataset, which is used as input for the forestplot function.

```R
# Linear associations of blood biomarkers to Body Mass Index (BMI), insulin
# resistance (log(HOMA-IR)) and fasting glucose
ggforestplot::df_linear_associations %>%
  print()

```

--------------------------------

### Basic Forest Plot Usage

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/ggforestplot.html

Demonstrates how to load the ggforestplot and tidyverse packages, filter the example dataset for BMI associations, and generate a basic forest plot using the forestplot() function.

```R
# Load and attach the package
library(ggforestplot)

# Load and attach other useful packages
# install.packages("tidyverse")
library(tidyverse)

# Filter only associations to BMI for the first 30 biomarkers of the example
# dataset
df <-
  ggforestplot::df_linear_associations %>%
  filter(
    trait == "BMI",
    dplyr::row_number() <= 30
  )

# Draw a forestplot of cross-sectional, linear associations
ggforestplot::forestplot(
  df = df,
  name = name,
  estimate = beta,
  se = se
)
```

--------------------------------

### Discovery Regression Function Usage (GLM Example)

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/discovery_regression.html

Shows how to use the discovery_regression function for logistic regression (GLM) to model incident diabetes risk. It involves data transformation, model formula specification, and parameter mapping.

```R
# Select only variables to be used for the model and
# collapse to a long format
df_long <-
  df_demo_metabolic_data %>%
  dplyr::select(
    tidyselect::all_of(nmr_biomarkers),
    gender,
    incident_diabetes,
    BMI,
    baseline_age
  ) %>%
  dplyr::mutate_at(
    .vars = dplyr::vars(tidyselect::all_of(nmr_biomarkers)),
    .funs = ~ .x %>% log1p() %>% scale() %>% as.numeric()
  ) %>%
  tidyr::gather(
    key = machine_readable_name,
    value = biomarkervalue,
    tidyselect::all_of(nmr_biomarkers)
  )

df_assoc_per_biomarker_gender <-
  discovery_regression(
    df_long = df_long,
    model = "glm",
    formula =
      formula(
        incident_diabetes ~ biomarkervalue + factor(gender) + BMI + baseline_age
      ),
    key = machine_readable_name,
    predictor = biomarkervalue
  )
```

--------------------------------

### ggforestplot geom_effect Example

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/geom_effect.html

Demonstrates how to use geom_effect with ggplot2 to plot study effects and their confidence intervals, including mapping aesthetics like color, shape, and significance.

```R
library(ggplot2)
library(magrittr)
df <-
  # Use built-in demo dataset
  df_linear_associations %>%
    # Arrange by name in order to filter the first few biomarkers for more
    # than one studies
    dplyr::arrange(name) %>%
    # Estimate confidence intervals
    dplyr::mutate(
      xmin = beta - qnorm(1 - (1 - 0.95) / 2) * se,
      xmax = beta + qnorm(1 - (1 - 0.95) / 2) * se
    ) %>%
    # Select only first 30 rows (10 biomarkers)
    dplyr::filter(dplyr::row_number() <= 30) %>%
    # Add a logical variable for statistical significance
    dplyr::mutate(filled = pvalue < 0.001)

g <-
  ggplot(data = df, aes(x = beta, y = name)) +
  # And point+errorbars
  geom_effect(
    ggplot2::aes(
      xmin = xmin,
      xmax = xmax,
      colour = trait,
      shape = trait,
      filled = filled
    ),
    position = ggstance::position_dodgev(height = 0.5)
  )
print(g)
```

--------------------------------

### Nightingale Biomarker Grouping Data Frame Structure

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/df_grouping_all_NG_biomarkers.html

Describes the structure of the df_grouping_all_NG_biomarkers data frame, which contains custom groupings for Nightingale Health biomarkers. This data frame is used in conjunction with the plot_all_NG_biomarkers function.

```APIDOC
df_grouping_all_NG_biomarkers:
  A data frame (tibble) with example custom groupings for Nightingale Health biomarkers.

  Structure:
  - version: Character. Biomarker platform version (e.g., '2016', '2020').
  - layout: A data frame (tibble) defining the layout for plotting.
    - machine_readable_name: Character. Biomarker machine readable name.
    - group_custom: Character. Custom group titles.
    - column: Integer. Column number for layout.
    - page: Integer. Page number for layout.

  Related Functions:
  - plot_all_NG_biomarkers: Function that utilizes this data frame for plotting.
  - df_NG_biomarker_metadata: Related data frame with biomarker metadata.
```

--------------------------------

### Inspect Variable Names

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

Prints the names of all variables (biomarkers and sample identifiers) present in the loaded NMR dataset. This is a common step to understand the structure and available data columns.

```R
names(df_nmr_results)
```

--------------------------------

### Get Nightingale Colour Hex Code

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/ng_colour.html

Retrieves the hex code for a specific Nightingale brand colour by its name. This function is useful for programmatically applying brand colours in visualizations or other design elements.

```R
ng_colour("dark pesto")

```

--------------------------------

### Display Log Odds Associations Data

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/ggforestplot.html

Shows how to print the `df_logodds_associations` dataset, which contains log odds ratios for blood biomarkers and associated metadata like standard errors and p-values.

```R
# Odds Ratios of Blood Biomarkers with Incident Type 2 Diabetes
ggforestplot::df_logodds_associations %>%
  print()
#> # A tibble: 1,135 x 6
#>    name          study            beta     se       pvalue     n
#>    <chr>         <chr>           <dbl>  <dbl>        <dbl> <dbl>
#>  1 Isoleucine    Meta-analysis  0.285  0.0500 0.0000000126 11777
#>  2 Leucine       Meta-analysis  0.286  0.0506 0.0000000165 11783
#>  3 Valine        Meta-analysis  0.182  0.0558 0.00110      11777
#>  4 Phenylalanine Meta-analysis  0.271  0.0538 0.000000453  11782
#>  5 Tyrosine      Meta-analysis  0.168  0.0541 0.00189      11779
#>  6 Alanine       Meta-analysis  0.122  0.0540 0.0237       11784
#>  7 Glutamine     Meta-analysis -0.172  0.0575 0.00284      11744
#>  8 Glycine       Meta-analysis -0.0773 0.0617 0.211        11496
#>  9 Histidine     Meta-analysis  0.0473 0.0553 0.393        11777
#> 10 Lactate       Meta-analysis  0.119  0.0529 0.0245       11781
#> # … with 1,125 more rows
```

--------------------------------

### Perform PCA and Calculate Variance Explained in R

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

This snippet demonstrates how to perform Principal Component Analysis (PCA) on metabolic data, scale and center it, and then calculate the cumulative variance explained by each principal component. It identifies the number of components explaining 99% of the variance, which is used for Bonferroni correction.

```R
# Perform principal component analysis on metabolic data
df_pca <- 
  df_full_data %>%
  select(tidyselect::all_of(nmr_biomarkers)) %>%
  nest(data = dplyr::everything()) %>%
  mutate(
    # Perform PCA on the data above, after scaling and centering to 0 
    pca = map(data, ~ stats::prcomp(.x, center = TRUE, scale = TRUE)),
    # Augment the original data with columns containing each observation's 
    # projection into the PCA space. Each PCA-space dimension is signified 
    # with the .fitted prefix in its name 
    pca_aug = map2(pca, data, ~broom::augment(.x, data = .y))
  )

# Estimate amount of variance explained by each principal component
df_pca_variance <- 
  df_pca %>%
  unnest(pca_aug) %>%
  # Estimate variance for the PCA projected variables
  summarize_at(.vars = vars(starts_with(".fittedPC")), .funs = ~ var(.x)) %>%
  # Gather data in a long format 
  gather(key = PC, value = variance) %>%
  # Estimate cumulative normalized variance
  mutate(cumvar = cumsum(variance / sum(variance)),
         PC = str_replace(PC, ".fitted", ""))

# Find number of principal components that explain 99% 
pc_99 <- 
  df_pca_variance %>%
  filter(cumvar <= 0.99) %>%
  nrow()

print(pc_99)
#> [1] 41

# Corrected significance threshold
psignif <- signif(0.05 / pc_99, 1)
```

--------------------------------

### Prepare and Order Biomarker Data for Forest Plot

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

This R code snippet defines a specific order for biomarker groups and then filters and arranges the Nightingale biomarker metadata accordingly. This prepares the data for creating a forest plot with a custom biomarker order.

```R
group_order <- c(
  "Branched-chain amino acids",
  "Aromatic amino acids",
  "Amino acids",
  "Fluid balance",
  "Inflammation",
  "Fatty acids",
  "Triglycerides",
  "Other lipids",
  "Cholesterol"
)

df_with_groups <- 
  ggforestplot::df_NG_biomarker_metadata %>%
  select(name = name,
         group) %>%
  filter(group %in% group_order) %>%
  arrange(factor(group, levels = group_order))
```

--------------------------------

### df_demo_metabolic_data Dataset Description

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/df_demo_metabolic_data.html

This entry describes the structure and content of the df_demo_metabolic_data dataset. It details the number of rows and columns, the types of variables included (e.g., id, gender, baseline_age, BMI, incident_diabetes, age_at_diabetes), and the nature of the biomarker data.

```R
## df_demo_metabolic_data

A data frame (tibble) with 1887 rows and 256 columns.
`id` is a character variable with the ID number of fictional individual;
`gender` is a character variable with the gender information on each individual;
`baseline_age` is a numeric variable with the age at the time of blood draw;
`BMI` is a numeric variable with the BMI of each individual at the time of blood draw;
`incident_diabetes` is a numeric variable with values 1 or 0 for whether the diabetes event occured or not during the observed time, respectively;
`age_at_diabetes` is a numeric variable with the age at the end of the study and the rest of the variables are numeric containing the machine readable names of Nightingale NMR blood biomarkers.
```

--------------------------------

### Read and Harmonize Clinical Data

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

Reads the clinical data CSV file, renames the 'identifier' column to 'Sample_id' to match the NMR data, and harmonizes the ID format by adding the 'ResearchInstitutionProjectLabId_' prefix. This prepares the clinical data for merging with the NMR results.

```R
# Read the clinical_data.csv data file
df_clinical_data <- readr::read_csv(
  # Enter the correct location for your file below
  file = "/path/to/file/clinical_data.csv") %>%
  # Rename the identifier column to Sample_id as in the NMR result file
  rename(Sample_id = identifier) %>%
  # Harmonize the id entries in clinical data with the ids in the NMR data 
  mutate(Sample_id = paste0(
    "ResearchInstitutionProjectLabId_", 
    as.numeric(Sample_id)
  ))

# Inspect the first 10 entries
print(df_clinical_data)
#> # A tibble: 1,887 x 6
#>    Sample_id          gender baseline_age   BMI incident_diabet… age_at_diabetes
#>    <chr>              <chr>         <dbl> <dbl>            <dbl>           <dbl>
#>  1 ResearchInstituti… male           57.6  25.2                1            59.9
#>  2 ResearchInstituti… male           36.6  25.5                0            52.2
#>  3 ResearchInstituti… female         37.9  31.4                0            56.6
#>  4 ResearchInstituti… male           57.9  30.5                0            70.9
#>  5 ResearchInstituti… male           70.9  28.6                0            86.5
#>  6 ResearchInstituti… female         41.1  26.9                0            63.8
#>  7 ResearchInstituti… male           47.1  24.2                0            59.7
#>  8 ResearchInstituti… female         32.4  23.7                0            47.2
#>  9 ResearchInstituti… female         58.1  19.3                1            72.4
#> 10 ResearchInstituti… male           52.4  31.6                0            67.4
#> # … with 1,877 more rows
```

--------------------------------

### Join NMR and Clinical Data

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

Performs a left join operation to combine the NMR results data (assumed to be in `df_nmr_results`) with the processed clinical data (`df_clinical_data`). The join is performed based on the harmonized `Sample_id` column, integrating patient demographics and outcomes with biomarker information.

```R
# Join NMR result file with the clinical data using column "Sample_id"
df_full_data <- dplyr::left_join(
  x = df_nmr_results,
  y = df_clinical_data,
  by = "Sample_id"
)
```

--------------------------------

### Plot All Nightingale Biomarkers

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/ggforestplot.html

Introduces the `plot_all_NG_biomarkers()` function from the ggforestplot package, designed to plot and print all Nightingale blood biomarkers into a single PDF file, supporting different layout versions.

```R
# Join the built-in association demo dataset with a variable that contains
# the machine readable names of Nightingale biomarkers. (Note: if you
```

--------------------------------

### CoxPH Regression with discovery_regression

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

Prepares data by log-transforming, scaling, and reshaping biomarker data into a long format. It then fits a Cox proportional hazards model using `ggforestplot::discovery_regression` to estimate gender-, age-, and BMI-adjusted hazard ratios for incident type 2 diabetes. The function requires a formula specifying the survival outcome and predictors, along with data keys.

```R
# Select only variables to be used for the model and collapse to a long data 
# format
df_long <-
  df_full_data %>%
  # Select only model variables
  dplyr::select(tidyselect::all_of(nmr_biomarkers), gender, baseline_age, BMI, age_at_diabetes, incident_diabetes) %>%
  # Log-transform and scale biomarkers
  dplyr::mutate_at(
    .vars = dplyr::vars(tidyselect::all_of(nmr_biomarkers)),
    .funs = ~ .x %>% log1p() %>% scale %>% as.numeric()
  ) %>%
  # Collapse to a long format
  tidyr::gather(
    key = biomarkerid,
    value = biomarkervalue,
    tidyselect::all_of(nmr_biomarkers)
  )

# Estimate sex-adjusted associations of metabolite to T2D
# Notice that parameter model below is set to 'coxph'
df_assoc_per_biomarker_diab <-
  discovery_regression(
    df_long = df_long, 
    model = "coxph",
    formula =
      formula(
        survival::Surv(
          time = baseline_age,
          time2 = age_at_diabetes,
          event = incident_diabetes
        ) ~ biomarkervalue + factor(gender) + BMI
      ),
    key = biomarkerid,
    predictor = biomarkervalue
  ) %>%
  # Join this dataset with the grouping data in order to choose a different 
  # biomarker naming option
  left_join(
    select(
      df_NG_biomarker_metadata, 
      name,
      biomarkerid = machine_readable_name
    ), 
    by = "biomarkerid")

head(df_assoc_per_biomarker_diab)
#> # A tibble: 6 x 5
#>   biomarkerid    estimate     se  pvalue name          
#>   <chr>             <dbl>  <dbl>   <dbl> <chr>         
#> 1 Total_C          0.0161 0.0690 0.815   Total-C       
#> 2 non_HDL_C        0.0524 0.0715 0.464   non-HDL-C     
#> 3 Remnant_C        0.0876 0.0708 0.216   Remnant-C     
#> 4 VLDL_C           0.238  0.0725 0.00101 VLDL-C        
#> 5 Clinical_LDL_C  -0.0229 0.0700 0.744   Clinical LDL-C
#> 6 LDL_C            0.0207 0.0707 0.770   LDL-C
```

--------------------------------

### Plot GlycA Distribution by Obesity Status

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

Generates boxplots to visualize the distribution of GlycA levels for obese and non-obese subjects. It involves adding an obesity status column based on BMI and then plotting using ggplot2, removing the x-axis title for clarity.

```R
df_full_data %>
  # Add a column to mark obese/non-obese subjects
  mutate(obesity = ifelse(BMI >= 30, yes = "Obese", no = "Not obese")) %>
  # Filter out subjects with missing values (if any)
  filter(!is.na(obesity)) %>
  # Box plots for each group
  ggplot(aes(y = GlycA, x = obesity, color = obesity)) +
  geom_boxplot() +
  # Remove x axis title 
  theme(axis.title.x = element_blank())
```

--------------------------------

### Display Biomarker Groups with ggforestplot

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

Shows how to extract and display unique blood biomarker groups from the df_NG_biomarker_metadata data frame using R and the ggforestplot package. This helps in understanding the available categories for plotting.

```R
ggforestplot::df_NG_biomarker_metadata %>%
  pull(group) %>%
  unique()
```

--------------------------------

### Document discovery_regression Function

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/discovery_regression.html

Documents the `discovery_regression` function, which fits multiple regression models (lm, glm, coxph) to data. It details the function's arguments, including data format, model types, formula specification, key and predictor variables, and verbosity. It also describes the expected return data frame containing estimates, standard errors, and p-values.

```R
discovery_regression(
  df_long,
  model = c("lm", "glm", "coxph"),
  formula,
  key = key,
  predictor = predictor,
  verbose = FALSE
)

Arguments
---------
df_long
a data frame in a long format that contains a `key` column with the names of the variables for which to estimate the measures of effect, e.g. Nightingale NMR biomarker names, a `value` column with the numeric values of the respective keys and other columns with the response and other variables to adjust for. Note: You may use `tidyr::gather` to collapse your data frame to key-value pairs; it is recommended to `dplyr::select` only the variables that you want to go into your model, in order to avoid memory overheads in case of large datasets.

model
a character, setting the type of model to fit on the input data frame. Must be one of `'lm'`, `'glm'` or `'coxph'`, for linear, logistic and cox proportional hazards regression, respectively.

formula
a formula object. For the case of models `'lm'` and `'glm'` it must be an object of class `stats::formula` (or one that can be coerced to that class): a symbolic description of the model to be fitted. The details of model specification are given under ‘Details’ in the documentation of `stats::formula`. For the case of model `'coxph'`, the formula object must have the response on the left of a ~ operator, and the rest of the terms on the right. The response must be a survival object, as returned by the `Surv` function.

key
the name of the `key` variable.

predictor
the name of the variable for which (adjusted) univariate associations are estimated. This is stated here in order to individuate the predictor of interest over many, possible cofactors that are also present in `formula`.

verbose
logical (default FALSE). If TRUE it prints a message with the names of the predictor and outcome. This may come in handy when, for example, fitting multiple outcomes.

Value
-------
A data frame with the following columns: a character variable with the same name as the `key` parameter and numeric variables `estimate`, `se` and `pvalue` with the values of the respective variables of the linear model. If `predictor` is factor, additional column `term` is returned indicating factor level of `predictor`
```

--------------------------------

### Estimate Biomarker Associations with BMI via Linear Regression

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

Estimates sex- and age-adjusted linear associations between numerous biomarkers and BMI using the `discovery_regression` function. This involves preparing the data by log-transforming and scaling biomarkers, collapsing it into a long format, and then running the regression model.

```R
# Extract names of NMR biomarkers
nmr_biomarkers <- dplyr::intersect(
  ggforestplot::df_NG_biomarker_metadata$machine_readable_name,
  colnames(df_nmr_results)
)

# NMR biomarkers here should be 250
stopifnot(length(nmr_biomarkers) == 250)

# Select only variables to be used for the model and collapse to a long data 
# format
df_long <-
  df_full_data %>
  # Select only model variables
  dplyr::select(tidyselect::all_of(nmr_biomarkers), gender, baseline_age, BMI) %>
  # Log-transform and scale biomarkers
  dplyr::mutate_at(
    .vars = dplyr::vars(tidyselect::all_of(nmr_biomarkers)),
    .funs = ~ .x %>% log1p() %>% scale %>% as.numeric()
  ) %>
  # Collapse to a long format
  tidyr::gather(
    key = biomarkerid,
    value = biomarkervalue,
    tidyselect::all_of(nmr_biomarkers)
  )

# Estimate sex- and age-adjusted associations of metabolite to BMI
df_assoc_per_biomarker_bmi <-
  ggforestplot::discovery_regression(
    df_long = df_long,
    model = "lm",
    formula =
      formula(
       biomarkervalue ~ BMI + factor(gender) + baseline_age
      ),
    key = biomarkerid,
    predictor = BMI
  ) %>
  # Join this dataset with the grouping data in order to choose a different 
  # biomarker naming option
  left_join(
    select(
      df_NG_biomarker_metadata, 
      name,
      biomarkerid = machine_readable_name
    ), 
    by = "biomarkerid")

head(df_assoc_per_biomarker_bmi)
```

--------------------------------

### Read NMR Metabolomics Data

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

Reads a CSV file containing Nightingale NMR metabolomics data using the `readr::read_csv` function. It specifies custom NA values ('NA', 'TAG') and column types, ensuring that 'Sample_id' is treated as a character and all other columns as doubles.

```R
# Read the biomarker concentration file
df_nmr_results <- readr::read_csv(
  # Enter the correct location for your file below
  file = "/path/to/file/12345-Results.csv",
  # Set not only NA but TAG string as <NA>
  na = c("NA", "TAG"), 
  col_types = cols(.default = col_double(), 
                   Sample_id = col_character())
)
```

--------------------------------

### Inspect Biomarker Order

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/ggforestplot.html

Prints the first 10 biomarker names from the grouping data frame (`df_grouping`) to inspect the order of biomarkers as determined by the metadata.

```r
# Print the first 10 rows of df_grouping to inspect the order of biomarkers
df_grouping %>%
  pull(name) %>%
  print()

```

--------------------------------

### Draw Forest Plot with Grouping and Significance in R

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

This snippet prepares data for plotting by joining association data with grouping information and then draws a forest plot. It uses the `forestplot` function, specifying the data frame, p-value column, corrected significance threshold (`psignif`), and x-axis label. Faceting by group is applied using `ggforce::facet_col`.

```R
# Join the association data frame with group data
df_to_plot <-
  df_assoc_per_biomarker_bmi %>%
  # use right_join, with df_grouping on the right, to preserve the order of 
  # biomarkers it specifies. 
  dplyr::right_join(., df_with_groups, by = "name") %>%
  dplyr::mutate(
    group = factor(.data$group, levels = group_order)
  ) %>%
  tidyr::drop_na(.data$estimate)

# Draw a forestplot of cross-sectional, linear associations.
forestplot(
  df = df_to_plot,
  pvalue = pvalue,
  psignif = psignif,
  xlab = "1-SD increment in biomarker concentration\nper unit increment in BMI"
) + 
  ggforce::facet_col(
    facets = ~group,
    scales = "free_y",
    space = "free"
  )
```

--------------------------------

### Plot all Nightingale Biomarkers to PDF

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/nmr-data-analysis-tutorial.html

The `plot_all_NG_biomarkers()` function from the `ggforestplot` package plots all Nightingale blood biomarkers into a single PDF file. It requires a data frame and specifies columns for biomarker ID, name, estimate, standard error, and p-value. The function supports different layouts and allows specifying an output filename.

```R
plot_all_NG_biomarkers(
  df = df_assoc_per_biomarker_bmi,
  machine_readable_name = biomarkerid,
  name = name,
  estimate = estimate,
  se = se,
  pvalue = pvalue,
  xlab = "1-SD increment in biomarker concentration\nper unit increment in BMI",
  filename = "all_ng_associations_bmi.pdf"
)
```

--------------------------------

### Nightingale Color Scale Constructors (R)

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/ng_scale.html

Defines discrete and continuous color scale functions for ggforestplot, allowing users to apply Nightingale Health color palettes to ggplot2 plots. Includes options for palette selection, reversal, and aesthetic mapping.

```APIDOC
scale_colour_ng_d(..., palette = "all", reverse = FALSE, aesthetics = "colour")
scale_fill_ng_d(..., palette = "all", reverse = FALSE, aesthetics = "fill")
scale_colour_ng_c(
  ...,
  palette = "magma",
  reverse = FALSE,
  values = NULL,
  space = "Lab",
  na.value = "grey50",
  guide = "colourbar",
  aesthetics = "colour"
)
scale_fill_ng_c(
  ...,
  palette = "magma",
  reverse = FALSE,
  values = NULL,
  space = "Lab",
  na.value = "grey50",
  guide = "colourbar",
  aesthetics = "fill"
)

Arguments:
  ...: Additional arguments passed to discrete_scale() or continuous_scale() to control name, limits, breaks, labels, and so forth.
  palette: Character name of the Nightingale (or viridis) colour palette. Defaults to "all" for discrete, "magma" for continuous.
  reverse: Boolean indicating whether the palette should be reversed. Defaults to FALSE.
  aesthetics: Character string or vector of character strings listing the name(s) of the aesthetic(s) that this scale works with (e.g., "colour", "fill").
  values: If colours should not be evenly positioned along the gradient, this vector gives the position (between 0 and 1) for each colour in the colours vector.
  space: Colour space in which to calculate gradient. Must be "Lab".
  na.value: Missing values will be replaced with this value. Defaults to "grey50".
  guide: A function used to create a guide or its name. Defaults to "colourbar".

Example:
  # Example taken from ggplot2::scale_colour_discrete()
  dsamp <- ggplot2::diamonds[sample(nrow(ggplot2::diamonds), 1000), ]
  d <- ggplot2::ggplot(dsamp, ggplot2::aes(carat, price)) +
    ggplot2::geom_point(ggplot2::aes(colour = clarity)) +
    ggforestplot::scale_colour_ng_d()
  print(d)
```

--------------------------------

### Prepare Data and Plot Nightingale Biomarkers

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/articles/ggforestplot.html

This R code snippet demonstrates how to join biomarker metadata with association data and then generate a custom forest plot for all Nightingale biomarkers. It uses `left_join` for data merging and `plot_all_NG_biomarkers` for visualization.

```R
df <-
  df_linear_associations %>%
  left_join(
    select(
      df_NG_biomarker_metadata,
      name,
      machine_readable_name
    ),
    by = "name"
  )

p <- plot_all_NG_biomarkers(
  df = df,
  machine_readable_name = machine_readable_name,
  estimate = beta,
  se = se,
  pvalue = pvalue,
  colour = trait,
  layout = "2016",
  xlab = "1-SD increment in cardiometabolic trait
      per 1-SD increment in biomarker concentration"
)
```

--------------------------------

### Basic Forest Plot with Log Odds

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/forestplot.html

Generates a basic forest plot using the `forestplot` function from the ggforestplot package. It maps study to color and sets the x-axis label for log odds ratios, after preparing the data with dplyr.

```R
df_logodds <-
  df_logodds_associations %>%
  dplyr::arrange(name) %>%
  dplyr::filter(dplyr::row_number() <= 30) %>%
  dplyr::mutate(
    study = factor(
      study,
      levels = c("Meta-analysis", "NFBC-1997", "DILGOM", "FINRISK-1997", "YFS")
    )
  )

forestplot(
  df = df_logodds,
  estimate = beta,
  logodds = TRUE,
  colour = study,
  xlab = "Odds ratio for incident type 2 diabetes (95% CI)
  per 1-SD increment in biomarker concentration"
)
```

--------------------------------

### Plot Odds Ratios with ggforestplot

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/README.md

Illustrates plotting odds ratios for type 2 diabetes risk using the `forestplot` function. It processes data to include multiple cohorts and their meta-analysis, customizing appearance with shapes.

```R
# Get subset of example, log odds ratios, data frame
df_logodds <-
  df_logodds_associations %>%
  dplyr::arrange(name) %>%
  dplyr::left_join(ggforestplot::df_NG_biomarker_metadata, by = "name") %>%
  dplyr::filter(group == "Amino acids") %>%
  # Set the study variable to a factor to preserve order of appearance
  # Set class to factor to set order of display.
  dplyr::mutate(
    study = factor(
      study,
      levels = c("Meta-analysis", "NFBC-1997", "DILGOM", "FINRISK-1997", "YFS")
    )
  )

# Forestplot
forestplot(
  df = df_logodds,
  estimate = beta,
  logodds = TRUE,
  colour = study,
  shape = study,
  title = "Associations to type 2 diabetes",
  xlab = "Odds ratio for incident type 2 diabetes (95% CI)
  per 1−SD increment in metabolite concentration"
) +
  # You may also want to add a manual shape scale to mark meta-analysis with a
  # diamond shape
  ggplot2::scale_shape_manual(
    values = c(23L, 21L, 21L, 21L, 21L),
    labels = c("Meta-analysis", "NFBC-1997", "DILGOM", "FINRISK-1997", "YFS")
  )
```

--------------------------------

### plot_all_NG_biomarkers Function Documentation

Source: https://github.com/nightingalehealth/ggforestplot/blob/master/docs/reference/plot_all_NG_biomarkers.html

Documentation for the `plot_all_NG_biomarkers` function, detailing its parameters, return value, and usage for creating forest plots of biomarker associations.

```APIDOC
plot_all_NG_biomarkers(df, machine_readable_name, estimate, se, pvalue, colour, filename = NULL, logodds = FALSE, paperwidth = 7, paperheight = 7, xlims = NULL, layout = "2016", ...)

  Generates forest plots for biomarker associations from Nightingale Health data.

  Arguments:
    df: A data frame containing the association data. Must include columns specified by machine_readable_name, estimate, se, and pvalue.
    machine_readable_name: The name of the column in `df` that contains the machine-readable names of Nightingale biomarkers.
    estimate: The name of the column in `df` containing the effect size estimates (e.g., beta coefficients).
    se: The name of the column in `df` containing the standard errors for the estimates.
    pvalue: The name of the column in `df` containing the p-values for the estimates.
    colour: The name of the column in `df` to be used for coloring the points (e.g., trait, study).
    filename: Optional. If provided, the plot will be saved to this file path. Defaults to NULL (plot is not saved).
    logodds: Logical. If TRUE, the plot will be formatted for log-odds ratios (e.g., for type 2 diabetes associations). Defaults to FALSE.
    paperwidth: The width of the output plot in inches. Defaults to 7.
    paperheight: The height of the output plot in inches. Defaults to 7.
    xlims: NULL or a numeric vector of length 2 specifying the common x-axis limits across all biomarker subgroups. Defaults to NULL.
    layout: Specifies the plot layout. Can be one of the predefined layouts (e.g., "2016") or a custom layout tibble.
    ...: Additional `ggplot2` graphical parameters such as `title`, `ylab`, `xlab`, `xtickbreaks`, etc., to be passed along.

  Value:
    If filename is NULL, a list of plot objects (one for each page in the layout) is returned.

  Details:
    The function uses a custom grouping specified by `df_grouping_all_NG_biomarkers`. The input `df` and `df_grouping_all_NG_biomarkers` are joined by `machine_readable_name`. Another `df` variable may be used for y-axis labels, defined in the `name` input parameter (though `name` is not explicitly listed as a parameter here, it's mentioned in the details section of the source text).

  Related Functions:
    `df_grouping_all_NG_biomarkers`: Used for defining custom layouts.
```