makeContrastsDream {variancePartition}R Documentation

Construct Matrix of Custom Contrasts

Description

Construct the contrast matrix corresponding to specified contrasts of a set of parameters.

Usage

makeContrastsDream(
  formula,
  data,
  ...,
  contrasts = NULL,
  suppressWarnings = FALSE
)

Arguments

formula

specifies variables for the linear (mixed) model. Must only specify covariates, since the rows of exprObj are automatically used a a response. e.g.: ~ a + b + (1|c) Formulas with only fixed effects also work

data

data.frame with columns corresponding to formula

...

expressions, or character strings which can be parsed to expressions, specifying contrasts

contrasts

character vector specifying contrasts

suppressWarnings

(default FALSE). suppress warnings for univariate contrasts

Details

This function expresses contrasts between a set of parameters as a numeric matrix. The parameters are usually the coefficients from a linear (mixed) model fit, so the matrix specifies which comparisons between the coefficients are to be extracted from the fit. The output from this function is usually used as input to dream().

This function is inspired by limma::makeContrasts() but is designed to be compatible with linear mixed models for dream()

Names in ... and contrasts will be used as column names in the returned value.

Value

matrix of linear contrasts between regression coefficients

Examples

# load library
# library(variancePartition)

# Intialize parallel backend with 4 cores
library(BiocParallel)
register(SnowParam(4))

# load simulated data:
# geneExpr: matrix of gene expression values
# info: information/metadata about each sample
data(varPartData)

form <- ~ 0 + Batch + (1|Individual) + (1|Tissue) 

# Define contrasts
L = makeContrastsDream( form, info, contrasts = c(Batch1_vs_2 = "Batch1 - Batch2", Batch3_vs_4 = "Batch3 - Batch4", Batch1_vs_34 = "Batch1 - (Batch3 + Batch4)/2"))

# show contrasts matrix
L

# Plot to visualize contrasts matrix
plotContrasts(L)

# Fit linear mixed model for each gene
# run on just 10 genes for time
fit = dream( geneExpr[1:10,], form, info, L=L)

# examine contrasts after fitting
head(coef(fit))

# show results from first contrast
topTable(fit, coef="Batch1_vs_2")

# show results from second contrast
topTable(fit, coef="Batch3_vs_4")

# show results from third contrast
topTable(fit, coef="Batch1_vs_34")


[Package variancePartition version 1.23.3 Index]