Olink QC Report

Analysis Details

This is a QC report based on the following data:

• Olink Data: /home/ec2-user/pipelines/proteomics/54gene-olink-qc/results/olink_data_postqc/example_54gene.tsv.gz
• Metadata/Technical Variables: /home/ec2-user/pipelines/proteomics/54gene-olink-qc/results/dataset_metadata/example_54gene.metadata.tsv.gz

Overall structure across the dataset

Here we visualize the overall PCs based on just the data matrix of probes by samples and the NPX values contained. We have colored each point arbitrarily by the first metadata variable that is available (to limit having hard requirements on the metadata). For a more detailed breakdown of technical metadata that correlate with proteomic variation, please see the last section of this report for a statistical test of significant covariates.

fig.pca <- plot_ly(merged.meta.df,
  x = ~PC1, y = ~PC2, text = ~SampleID,
  type ="scatter",
  mode = "markers",
  #name = as.formula(paste("~", meta.variables[1])),
  color = as.formula(paste("~", meta.variables[1])),
  marker = list(size = 10)
  ) %>% layout(
    title = "PCA Overview of Olink Dataset",
    legend = list(title = list(text = paste("<b>", meta.variables[1], "</b>"))),
    xaxis = list(title = "PC1"),
    yaxis = list(title = "PC2")
)

fig.pca

Missingness across samples

Missingness across samples are measured by the proportion of assays that fail Olink’s internal QC measures. Full details on how to interpret their QC metrics and warnings can be found here: https://www.olink.com/faq/how-is-quality-control-of-the-data-performed/.

QC Warnings from Olink are provided based on internal controls and detection controls with respect to a given 96-well plate. The rule for determining if a QC Warning is emitted is if any of these control values deviate more than 0.3 NPX from the plate median. As a general rule, if more than 1/6th of the samples fail the QC then the run is deemed to be unreliable and this should be noted and sent as immediate feedback to the Olink team.

olink.merged.df <- fread(snakemake@input[["merged_olink_data"]], sep = "\t")
qc.warning.samples.df <- olink.merged.df %>%
  group_by(SampleID) %>%
  summarise(QC_Warning = sum(QC_Warning == "WARN"), N = n()) %>%
  mutate(Fraction_QC_Warning = QC_Warning / N)

reactable(qc.warning.samples.df)

Missingness across Olink Panels

The Olink Explore platform allows one to evaluate proteomic variation across different panels (e.g. Neurology or Oncology). Since each 96-well plate has a limitation of 384 probes, each panel is effectively a different experiment. As these are the same samples, but different panels, each panel may exhibit different properties across samples, including missingness. Here we quantify the extent to which each sample has probes not passing QC on each of the panels to get a more granular view at how this varies across multiple panels.

qc.warning.samples.panel.df <- olink.merged.df %>%
  group_by(SampleID, Panel) %>%
  summarise(QC_Warning = sum(QC_Warning == "WARN"), N = n()) %>%
  mutate(Fraction_QC_Warning = QC_Warning / N)

reactable(qc.warning.samples.panel.df)

We also display this missingness below across all panels to visualize individuals whose missingness may be driven largely by a single panel or if missingness is well-distributed across panels:

warning.samples.panel.df <- qc.warning.samples.panel.df %>%
  pivot_wider(
    id_cols = "SampleID",
    names_from = Panel,
    values_from = QC_Warning,
    values_fn = mean
  )

fig <- reshape2::melt(warning.samples.panel.df, id.vars = "SampleID") %>%
  plot_ly(
    x = ~SampleID, y = ~value, type = "bar",
    name = ~variable, color = ~variable
  ) %>%
  layout(
    yaxis = list(title = "Count"),
    xaxis = list(categoryorder = "total descending"),
    legend = list(title = list(text = "<b> Panel </b>")),
    barmode = "stack"
  )
fig

Per-Panel Sample Outliers

Outliers (either in terms of missingness or measurement) can appear on a per-panel basis and we may want to filter specific individuals per-panel. In this setting we are recreating the olink_qc_plot function from the OlinkAnalyze package but for a more interactive view per-panel. Here we set the threshold for outlier detection as being outside:

• outside 3 standard deviations on either the sample median NPX values
• outside 3 standard deviations outside the inter-quartile range

Dynamic-Range Outliers

# Defining outliers based on medians and assay IQR
olink.outliers.df <- olink.merged.df %>%
  inner_join(merged.meta.df, by = "SampleID") %>%
  group_by(Panel, SampleID, Index) %>%
  mutate(
    iqr = IQR(NPX, na.rm = TRUE),
    sample_median = median(NPX, na.rm = TRUE),
    pc1_mean = mean(PC1, na.rm = TRUE),
    pc2_mean = mean(PC2, na.rm = TRUE)
  ) %>%
  ungroup() %>%
  select(SampleID, Index, Panel, iqr, sample_median, pc1_mean, pc2_mean, PC1, PC2) %>%
  distinct() %>%
  group_by(Panel) %>%
  mutate(
    pc1_outlier_low = mean(pc1_mean, na.rm = TRUE) - pc.sd * sd(pc1_mean, na.rm = TRUE),
    pc1_outlier_high = mean(pc1_mean, na.rm = TRUE) + pc.sd * sd(pc1_mean, na.rm = TRUE),
    pc2_outlier_low = mean(pc2_mean, na.rm = TRUE) - pc.sd * sd(pc2_mean, na.rm = TRUE),
    pc2_outlier_high = mean(pc2_mean, na.rm = TRUE) + pc.sd * sd(pc2_mean, na.rm = TRUE),
    median_low = mean(sample_median, na.rm = TRUE) - median.sd * sd(sample_median, na.rm = TRUE),
    median_high = mean(sample_median, na.rm = TRUE) + median.sd * sd(sample_median, na.rm = TRUE),
    iqr_low = mean(iqr, na.rm = TRUE) - iqr.sd * sd(iqr, na.rm = TRUE),
    iqr_high = mean(iqr, na.rm = TRUE) + iqr.sd * sd(iqr, na.rm = TRUE),
  ) %>%
  ungroup() %>%
  mutate(valid_outliers = remove.outliers) %>%
  mutate(outlier_range = ifelse(sample_median > median_low & sample_median < median_high & iqr > iqr_low & iqr < iqr_high & valid_outliers, FALSE, TRUE)) %>%
  mutate(pc_outlier = ifelse(pc1_mean > pc1_outlier_low & pc1_mean < pc1_outlier_high & pc2_mean > pc2_outlier_low & pc2_mean < pc2_outlier_high & valid_outliers, FALSE, TRUE)) %>%
  mutate(outlier = (outlier_range | pc_outlier) & valid_outliers)

uniq.panels <- unique(olink.outliers.df$Panel)

fig1 <- plot_ly(olink.outliers.df %>% filter(Panel == "Neurology"),
  x = ~sample_median,
  y = ~iqr,
  type = "scatter",
  mode = "markers",
  text = ~SampleID,
  color = ~outlier,
  marker = list(size = 10)
) %>%
  layout(
    yaxis = list(title = "IQR"),
    xaxis = list(title = "Sample Median"),
    legend = list(title = list(text = "<b> Outlier </b>")),
    title = "Neurology"
  )

fig2 <- plot_ly(olink.outliers.df %>% filter(Panel == "Cardiometabolic"),
  x = ~sample_median,
  y = ~iqr,
  type = "scatter",
  mode = "markers",
  text = ~SampleID,
  color = ~outlier,
  marker = list(size = 10)
) %>%
  layout(
    yaxis = list(title = "IQR"),
    xaxis = list(title = "Sample Median"),
    legend = list(title = list(text = "<b> Outlier </b>")),
    title = "Cardiometabolic"
  )

fig3 <- plot_ly(olink.outliers.df %>% filter(Panel == "Oncology"),
  x = ~sample_median,
  y = ~iqr,
  type = "scatter",
  mode = "markers",
  text = ~SampleID,
  color = ~outlier,
  marker = list(size = 10)
) %>%
  layout(
    yaxis = list(title = "IQR"),
    xaxis = list(title = "Sample Median"),
    legend = list(title = list(text = "<b> Outlier </b>")),
    title = "Oncology"
  )

fig4 <- plot_ly(olink.outliers.df %>% filter(Panel == "Inflammation"),
  x = ~sample_median,
  y = ~iqr,
  type = "scatter",
  mode = "markers",
  text = ~SampleID,
  color = ~outlier,
  marker = list(size = 10)
) %>%
  layout(
    yaxis = list(title = "IQR"),
    xaxis = list(title = "Sample Median"),
    legend = list(title = list(text = "<b> Outlier </b>")),
    title = "Inflammation"
  )

fig1

fig2

fig3

fig4

PCA-Based Outliers

The PCA-based outliers are also visualized in the following plot. The rule to determine outliers is based on a sample being outside 3 outside standard deviations from the means of PC1 and PC2

fig <- plot_ly(olink.outliers.df,
  x = ~PC1,
  y = ~PC2,
  type = "scatter",
  mode = "markers",
  text = ~SampleID,
  color = ~pc_outlier,
  marker = list(size = 10)
) %>%
  layout(
    yaxis = list(title = "PC2"),
    xaxis = list(title = "PC1"),
    legend = list(title = list(text = "<b> Outlier </b>")),
    title = "PCA-based Outliers"
  )

fig

Total Recommended Outlier Filtering

Based on the data presented above, we recommend removing the following Sample / Panel combinations:

olink.outliers.filt.df <- olink.outliers.df %>%
  filter(outlier == TRUE) %>%
  select(SampleID, Panel, outlier_range, pc_outlier)

reactable(olink.outliers.filt.df)

Missingness across Olink Probes

Assay Warnings across Probes

Olink also provides a statement of whether each Assay passes their QC checks across all individuals, based on the Incubation Control and the Amplification Control setups per-plate. If an assay has an “Assay Warning”, then we will want to remove the assay from downstream analysis like differential expression or pQTL analyses. The Olink team recommends that we remove any probes that show an assay warning from downstream analyses.

if ("Assay_Warning" %in% colnames(olink.merged.df)) {
  missingness.olink.id.assays <- olink.merged.df %>%
    group_by(OlinkID) %>%
    summarise(Assay_Warning = sum(Assay_Warning == "WARN"), N = n()) %>%
    mutate(Fraction_Assay_Warning = Assay_Warning / N)

  reactable(missingness.olink.id.assays)
} else {
  cat(paste(getwd(), snakemake@input[["merged_olink_data"]], " does not contain Assay Warnings!"))
}