Skip to contents

Performs linear trend tests on selected Principal Coordinates (PCs) of beta diversity distance matrices over time, for different groups. Allows using PCoA, NMDS, t-SNE, UMAP for dimension reduction, with PCoA as default.

Usage

generate_beta_pc_trend_test_long(
  data.obj = NULL,
  dist.obj = NULL,
  pc.obj = NULL,
  pc.ind = c(1, 2),
  subject.var,
  time.var,
  group.var = NULL,
  adj.vars = NULL,
  dist.name = c("BC"),
  ...
)

Arguments

data.obj

A list object in a format specific to MicrobiomeStat, which can include components such as feature.tab (matrix), feature.ann (matrix), meta.dat (data.frame), tree, and feature.agg.list (list). The data.obj can be converted from other formats using several functions from the MicrobiomeStat package, including: 'mStat_convert_DGEList_to_data_obj', 'mStat_convert_DESeqDataSet_to_data_obj', 'mStat_convert_phyloseq_to_data_obj', 'mStat_convert_SummarizedExperiment_to_data_obj', 'mStat_import_qiime2_as_data_obj', 'mStat_import_mothur_as_data_obj', 'mStat_import_dada2_as_data_obj', and 'mStat_import_biom_as_data_obj'. Alternatively, users can construct their own data.obj. Note that not all components of data.obj may be required for all functions in the MicrobiomeStat package.

dist.obj

Distance matrix between samples, usually calculated using mStat_calculate_beta_diversity function. If NULL, beta diversity will be automatically computed from data.obj using mStat_calculate_beta_diversity.

pc.obj

A list containing the results of dimension reduction/Principal Component Analysis. This should be the output from functions like mStat_calculate_PC, containing the PC coordinates and other metadata. If NULL (default), dimension reduction will be automatically performed using metric multidimensional scaling (MDS) via mStat_calculate_PC. The pc.obj list structure should contain:

points

A matrix with samples as rows and PCs as columns containing the coordinates.

eig

Eigenvalues for each PC dimension.

vectors

Loadings vectors for features onto each PC.

Other metadata

like method, dist.name, etc.

See mStat_calculate_PC function for details on output format.

pc.ind

Numeric vector specifying which principal coordinate (PC) axes to include in the trend test, e.g. c(1,2) for PC1 and PC2. Should not exceed the number of PCs calculated in pc.obj.

subject.var

Character string specifying the column name in metadata containing unique subject IDs. This should uniquely identify each subject in the study. Required to fit mixed effects models over time.

time.var

Character string specifying the column in metadata containing the numeric time variable. Should contain ordered time points for trend test. Required.

group.var

Character string specifying the column in metadata containing a grouping variable. Separate models will be fitted for each group. Optional, can be left NULL.

adj.vars

Character vector specifying columns in metadata containing covariates to adjust distance matrices for prior to ordination. Optional, can be left NULL.

dist.name

A character vector specifying which beta diversity indices to calculate. Supported indices are "BC" (Bray-Curtis), "Jaccard", "UniFrac" (unweighted UniFrac), "GUniFrac" (generalized UniFrac), "WUniFrac" (weighted UniFrac), and "JS" (Jensen-Shannon divergence). If a name is provided but the corresponding object does not exist within dist.obj, it will be computed internally. If the specific index is not supported, an error message will be returned. Default is c('BC', 'Jaccard').

...

Additional named arguments to pass to lmer():

  • control: Control parameters for lmer.

  • weights: Prior weights for the residuals.

Value

A nested list by distance > PC. Each element contains:

  • coef: Data frame of coefficients from mixed effects model.

  • model: Fitted lmer model object.

Details

This function allows performing linear trend tests on PCs of beta diversity distances over time, across groups, with adjustments.

It checks for pre-calculated distances and PCs, generating them from data if needed. Sufficient PCs should be calculated to cover indices specified.

For each distance, PCs are extracted for the selected indices. The metadata is joined and formatted for mixed effects modeling with time as numeric.

The model formula is created using the response, time, group, and subject variables. Mixed effects models are fitted with lmer() and coefficients extracted.

See also

mStat_calculate_beta_diversity to generate distance matrices.

mStat_calculate_PC to generate PCs from distances.

Examples

if (FALSE) { # \dontrun{
library(vegan)
data("ecam.obj")
generate_beta_pc_trend_test_long(
  data.obj = ecam.obj,
  dist.obj = NULL,
  pc.obj = NULL,
  pc.ind = c(1, 2),
  subject.var = "studyid",
  time.var = "month",
  group.var = "diet",
  adj.vars = "delivery",
  dist.name = c('BC')
)

data("subset_T2D.obj")
generate_beta_pc_trend_test_long(
  data.obj = subset_T2D.obj,
  dist.obj = NULL,
  pc.obj = NULL,
  pc.ind = c(1, 2),
  subject.var = "subject_id",
  time.var = "visit_number_num",
  group.var = "subject_race",
  adj.vars = "subject_gender",
  dist.name = c('BC')
)
} # }