test_differential_abundance {tidybulk}R Documentation

Add differential transcription information to a tbl using edgeR.

Description

test_differential_abundance() takes as imput a 'tbl' formatted as | <SAMPLE> | <TRANSCRIPT> | <COUNT> | <...> | and returns a 'tbl' with additional columns for the statistics from the hypothesis test.

Usage

test_differential_abundance(
  .data,
  .formula,
  .sample = NULL,
  .transcript = NULL,
  .abundance = NULL,
  .contrasts = NULL,
  method = "edgeR_quasi_likelihood",
  significance_threshold = 0.05,
  minimum_counts = 10,
  minimum_proportion = 0.7,
  fill_missing_values = FALSE,
  scaling_method = "TMM",
  omit_contrast_in_colnames = FALSE,
  action = "add"
)

## S4 method for signature 'spec_tbl_df'
test_differential_abundance(
  .data,
  .formula,
  .sample = NULL,
  .transcript = NULL,
  .abundance = NULL,
  .contrasts = NULL,
  method = "edgeR_quasi_likelihood",
  significance_threshold = 0.05,
  minimum_counts = 10,
  minimum_proportion = 0.7,
  fill_missing_values = FALSE,
  scaling_method = "TMM",
  omit_contrast_in_colnames = FALSE,
  action = "add"
)

## S4 method for signature 'tbl_df'
test_differential_abundance(
  .data,
  .formula,
  .sample = NULL,
  .transcript = NULL,
  .abundance = NULL,
  .contrasts = NULL,
  method = "edgeR_quasi_likelihood",
  significance_threshold = 0.05,
  minimum_counts = 10,
  minimum_proportion = 0.7,
  fill_missing_values = FALSE,
  scaling_method = "TMM",
  omit_contrast_in_colnames = FALSE,
  action = "add"
)

## S4 method for signature 'tidybulk'
test_differential_abundance(
  .data,
  .formula,
  .sample = NULL,
  .transcript = NULL,
  .abundance = NULL,
  .contrasts = NULL,
  method = "edgeR_quasi_likelihood",
  significance_threshold = 0.05,
  minimum_counts = 10,
  minimum_proportion = 0.7,
  fill_missing_values = FALSE,
  scaling_method = "TMM",
  omit_contrast_in_colnames = FALSE,
  action = "add"
)

## S4 method for signature 'SummarizedExperiment'
test_differential_abundance(
  .data,
  .formula,
  .sample = NULL,
  .transcript = NULL,
  .abundance = NULL,
  .contrasts = NULL,
  method = "edgeR_quasi_likelihood",
  significance_threshold = 0.05,
  minimum_counts = 10,
  minimum_proportion = 0.7,
  fill_missing_values = FALSE,
  scaling_method = "TMM",
  omit_contrast_in_colnames = FALSE,
  action = "add"
)

## S4 method for signature 'RangedSummarizedExperiment'
test_differential_abundance(
  .data,
  .formula,
  .sample = NULL,
  .transcript = NULL,
  .abundance = NULL,
  .contrasts = NULL,
  method = "edgeR_quasi_likelihood",
  significance_threshold = 0.05,
  minimum_counts = 10,
  minimum_proportion = 0.7,
  fill_missing_values = FALSE,
  scaling_method = "TMM",
  omit_contrast_in_colnames = FALSE,
  action = "add"
)

Arguments

.data

A 'tbl' formatted as | <SAMPLE> | <TRANSCRIPT> | <COUNT> | <...> |

.formula

A formula with no response variable, representing the desired linear model

.sample

The name of the sample column

.transcript

The name of the transcript/gene column

.abundance

The name of the transcript/gene abundance column

.contrasts

A character vector. See edgeR makeContrasts specification for the parameter 'contrasts'. If contrasts are not present the first covariate is the one the model is tested against (e.g., ~ factor_of_interest)

method

A string character. Either "edgeR_quasi_likelihood" (i.e., QLF), "edgeR_likelihood_ratio" (i.e., LRT)

significance_threshold

A real between 0 and 1 (usually 0.05).

minimum_counts

A real positive number. It is the threshold of count per million that is used to filter transcripts/genes out from the scaling procedure.

minimum_proportion

A real positive number between 0 and 1. It is the threshold of proportion of samples for each transcripts/genes that have to be characterised by a cmp bigger than the threshold to be included for scaling procedure.

fill_missing_values

A boolean. Whether to fill missing sample/transcript values with the median of the transcript. This is rarely needed.

scaling_method

A character string. The scaling method passed to the backend function (i.e., edgeR::calcNormFactors; "TMM","TMMwsp","RLE","upperquartile")

omit_contrast_in_colnames

If just one contrast is specified you can choose to omit the contrast label in the colnames.

action

A character string. Whether to join the new information to the input tbl (add), or just get the non-redundant tbl with the new information (get).

Details

Maturing lifecycle

At the moment this function uses edgeR only, but other inference algorithms will be added in the near future.

Value

A 'tbl' with additional columns for the statistics from the hypothesis test (e.g., log fold change, p-value and false discovery rate).

A 'tbl' with additional columns for the statistics from the hypothesis test (e.g., log fold change, p-value and false discovery rate).

A 'tbl' with additional columns for the statistics from the hypothesis test (e.g., log fold change, p-value and false discovery rate).

A 'tbl' with additional columns for the statistics from the hypothesis test (e.g., log fold change, p-value and false discovery rate).

A 'SummarizedExperiment' object

A 'SummarizedExperiment' object

Examples



	test_differential_abundance(
	 tidybulk::counts_mini,
	    ~ condition,
	    sample,
	    transcript,
	    `count`
	)

	# The functon `test_differential_abundance` operated with contrasts too

 test_differential_abundance(
	    tidybulk::counts_mini,
	    ~ 0 + condition,
	    sample,
	    transcript,
	    `count`,
	    .contrasts = c( "conditionTRUE - conditionFALSE")
 )



[Package tidybulk version 1.0.2 Index]