Skip to contents

Differential analysis of proximity scores

DifferentialProximityAnalysis.Rd

Runs differential analysis (Running Wilcoxon rank-sum test) on proximity scores calculated from Proximity Network Assay (PNA) data generated with the Pixelator data processing pipeline.

Usage

DifferentialProximityAnalysis(object, ...)

# S3 method for class 'data.frame'
DifferentialProximityAnalysis(
  object,
  contrast_column,
  reference,
  targets = NULL,
  group_vars = NULL,
  proximity_metric = "log2_ratio",
  metric_type = c("all", "self", "co"),
  backend = c("dplyr", "data.table"),
  p_adjust_method = c("bonferroni", "holm", "hochberg", "hommel", "BH", "BY", "fdr"),
  min_cells_per_group = 10,
  verbose = TRUE,
  ...
)

# S3 method for class 'Matrix'
DifferentialProximityAnalysis(
  object,
  group_data,
  contrast_column,
  reference,
  targets = NULL,
  group_vars = NULL,
  proximity_metric = "log2_ratio",
  p_adjust_method = c("bonferroni", "holm", "hochberg", "hommel", "BH", "BY", "fdr"),
  diff_threshold = 0.01,
  min_pct = 0,
  min_diff_pct = -Inf,
  min_cells_per_group = 10,
  verbose = TRUE,
  ...
)

# S3 method for class 'Seurat'
DifferentialProximityAnalysis(
  object,
  contrast_column,
  reference,
  targets = NULL,
  assay = NULL,
  group_vars = NULL,
  lazy = FALSE,
  min_exp_join_count = 0,
  min_cells_per_group = 10,
  diff_threshold = 0.01,
  min_pct = 0,
  min_diff_pct = -Inf,
  proximity_metric = "log2_ratio",
  metric_type = c("all", "self", "co"),
  backend = c("dplyr", "data.table"),
  method = c("seurat", "legacy"),
  p_adjust_method = c("bonferroni", "holm", "hochberg", "hommel", "BH", "BY", "fdr"),
  verbose = TRUE,
  ...
)

Arguments

object

An object containing PNA proximity scores

...

Not yet implemented

contrast_column

The name of the column where the group labels are stored. This column must include target and reference.

reference

The name of the reference group

targets

The names of the target groups. These groups will be compared to the reference group. If the value is set to NULL (default), all groups available in contrast_column except reference will be compared to the reference group.

group_vars

An optional character vector with column names to group the tests by.

proximity_metric

The proximity metric to use. Any numeric data column in the proximity score table can be selected. The default is "pearson_z".

metric_type

One of "all", "self" or "cross". If "all", all pairwise comparisons are considered. If "self", only protein pairs of the same type are considered. If "cross", protein pairs of different type are considered.

backend

One of "dplyr" or "data.table". The latter requires the dtplyr package to be installed.

p_adjust_method

One of "bonferroni", "holm", "hochberg", "hommel", "BH", "BY" or "fdr". (see ?p.adjust for details)

min_cells_per_group

Minimum number of cells required per group. Target groups with less observations than min_cells_per_group will be skipped.

verbose

Print messages

group_data

A tibble with a column for the contrast and optional group variables. The rownames of this tibble should correspond to the columns names of the matrix object.

diff_threshold

Minimum difference in proximity metric to consider a pair of groups for testing. Default is 0.1. This parameter is only used when the method argument is set to "seurat".

min_pct

Minimum percentage of cells in either group that must express a marker pair for it to be considered for testing. Default is 0. This parameter is only used when the method argument is set to "seurat".

min_diff_pct

Minimum difference in percentage of cells expressing a marker pair between the two groups for it to be considered for testing. Default is -Inf. This parameter is only used when the method argument is set to "seurat".

assay

Name of assay to use

lazy

If TRUE, the proximity scores will be loaded lazily and filtered using the duckdb backend.

min_exp_join_count

Minimum number of join counts required for a marker pair to be included in the analysis. Dropped protein pairs (those with fewer than min_exp_join_count counts) will be treated as missing entries. With method = "seurat", these are treated as having a proximity score of 0.

method

One of "seurat" or "legacy". The former uses the Seurat framework for differential testing, while the latter uses a custom implementation. The main difference between the two methods is that missing observations are handled differently. With the "seurat" method, all missing values are set to 0, while the "legacy" method ignores missing values. For the former, this means that the number of observations per group is always equal to the number of cells in that group, while for the latter ("legacy"), the number of observations per group can be less than the number of cells in that group. Ignoring missing values can lead to confusing estimates of group statistics and misses comparisons where one of the two test groups have no observations. The "legacy" method is provided for backward compatibility and may be removed in future versions.

Value

A tbl_df object with test results

Details

If you are working with a Seurat object created with pixelatorR, the proximity scores can be accessed with ProximityScores.

The input object should contain a contrast_column (character vector or factor) that includes information about the groups to compare. A typical example is a column with sample labels, for instance: "control", "stimulated1", "stimulated2". If the input object is a Seurat object, the contrast_column should be available in the meta.data slot. For those familiar with FindMarkers from Seurat, contrast_column is comparable to the group.by parameter.

The targets parameter specifies a character vector with the names of the groups to compare reference. targets can be a single group name or a vector of group names while reference can only refer to a single group. Both targets and reference should be present in the contrast_column. These parameters are similar to the ident.1 and ident.2 parameters in FindMarkers.

Additional groups

The tests are always computed between targets and reference, but it is possible to add additional grouping variables with group_vars. If group_vars is used, each comparison is split into groups defined by the group_vars. For instance, if we have annotated cells into cell type populations and saved these annotations in a meta.data column called "cell_type", we can pass "cell_type" to group_vars="cell_type" to split tests across each cell type.

Types of comparisons

Consider a scenario where we have a Seurat object (seurat_object) with Proximity Network Assay (PNA) data. seurat_object contains a meta.data column called "sampleID" that holds information about what samples the components originated from. This column could have three sample IDs: "control", "stimulated1" and "stimulated2". In addition, we have a column called "cell_type" that holds information about the cell type identity of each component.

  1. If we want to compare the "stimulated1" group to the "control" group:

    
dp_markers <- DifferentialProximityAnalysis(
   object = seurat_object,
   contrast_column = "sampleID",
   reference = "control",
   targets = "stimulated1"
)

  2. If we want to compare the "stimulated1" and "stimulated2" groups to the "control" group:

    
dp_markers <- DifferentialProximityAnalysis(
  object = seurat_object,
  contrast_column = "sampleID",
  reference = "control",
  targets = c("stimulated1", "stimulated2")
)

  3. If we want to compare the "stimulated1" and "stimulated2" groups to the "control" group within each cell type:

    
dp_markers <- DifferentialProximityAnalysis(
   object = seurat_object,
   contrast_column = "sampleID",
   reference = "control",
   targets = c("stimulated1", "stimulated2"),
   group_vars = "cell_type"
)

Method details

By default, the function uses FindMarkers under the hood to compute the tests (method = "Seurat"). The proximity scores are stored in long format in a table (either as a table in the Seurat object or in the associated PXL file). The function reshapes the table to wide format and uses FindMarkers to compute the tests.

Note that missing observations in the proximity score table are replaced with 0's when reshaping to wide format. This means that if a marker pair is not present in a component, it will be treated as having a proximity score of 0 for that component. With this strategy, we cannot distinguish between a marker pair that is truly absent in a component and a marker pair that is present but has a proximity score of 0. Marker pairs rarely have a proximity score of 0, except for extremely low abundant pairs which are usually filtered out in the pre-processing steps.

With the "legacy" method, the missing observations are ignored and the tests are computed only on the marker pairs that are present in both groups. This strategy has two important limitations: 1) it tends to focus the results on subsets of cells which are not necessarily representative of the entire population of interest, and 2) if one of the two conditions (target or reference) has no observations for a pair, that comparison will be skipped.

See also

Other DA-methods: RunDAA(), RunDCA(), RunDPA()

Examples

library(dplyr)
example_data <- tidyr::expand_grid(
  marker_1 = c("HLA-ABC", "B2M", "CD4", "CD8", "CD20", "CD19", "CD45", "CD43") %>%
    rep(each = 50),
  marker_2 = c("HLA-ABC", "B2M", "CD4", "CD8", "CD20", "CD19", "CD45", "CD43") %>%
    rep(each = 50)
) %>%
  mutate(
    join_count_z = rnorm(n(), sd = 10)
  )

example_data <- example_data %>%
  mutate(sampleID = "ctrl") %>%
  bind_rows(
    example_data %>% mutate(join_count_z = join_count_z + 1) %>%
      mutate(sampleID = "treatment")
  )

# Compute statistics
dp_results <- DifferentialProximityAnalysis(
  example_data,
  contrast_column = "sampleID",
  reference = "ctrl",
  proximity_metric = "join_count_z",
  metric_type = "self"
)
#> ℹ Running Wilcoxon rank-sum test across the following comparisons:
#> 
#>   • treatment vs ctrl