Estimate bias using the compositional least-squares approach described in McLaren, Willis, and Callahan (2019).
estimate_bias(observed, actual, ...) # S3 method for matrix estimate_bias(observed, actual, margin, boot = FALSE, times = 1000) # S3 method for otu_table estimate_bias(observed, actual, ...) # S3 method for phyloseq estimate_bias(observed, actual, ...)
| observed | Abundance matrix of observed compositions. |
|---|---|
| actual | Abundance matrix of actual or reference compositions for the
same samples and taxa in |
| ... | Arguments passed to the matrix method. |
| margin | Matrix margin that corresponds to observations (samples);
|
| boot | Whether to perform bootstrapping. |
| times | Number of bootstrap replicates. |
A mc_bias_fit object with coef(), fitted(), residuals(), and
summary() methods.
Bias is estimated by applying center() to the compositional error matrix
defined by observed/actual, which requires that observed and actual
are non-zero for the same sample-taxa pairs. For convenience, this
function will automatically set values in observed to 0 whose
corresponding entries are 0 in actual, but it is up to you to replace 0
values in observed with a non-zero value (such as a pseudocount).
Requirements for observed and actual: The row and column names (for
matrices) or taxa and sample names (for phyloseq objects) must match, but
can be in different orders. Any taxa and samples in observed but not in
actual will be dropped prior to estimation.
# Load data from the cellular mock communities of Brooks et al 2015 dr <- system.file("extdata", package = "metacal") list.files(dr)#> [1] "brooks2015-actual.csv" "brooks2015-observed.csv" #> [3] "brooks2015-sample-data.csv"actual <- file.path(dr, "brooks2015-actual.csv") |> read.csv(row.names = "Sample") |> as("matrix") observed <- file.path(dr, "brooks2015-observed.csv") |> read.csv(row.names = "Sample") |> subset(select = - Other) |> as("matrix") sam <- file.path(dr, "brooks2015-sample-data.csv") |> read.csv() # Estimate bias with bootstrapping for error estimation mc_fit <- estimate_bias(observed, actual, margin = 1, boot = TRUE)#>#> Summary of a metacal bias fit. #> #> Estimated relative efficiencies: #> # A tibble: 7 × 4 #> taxon estimate gm_mean gm_se #> <chr> <dbl> <dbl> <dbl> #> 1 Atopobium_vaginae 0.285 0.285 1.03 #> 2 Gardnerella_vaginalis 0.160 0.159 1.05 #> 3 Lactobacillus_crispatus 2.29 2.29 1.03 #> 4 Lactobacillus_iners 4.68 4.69 1.02 #> 5 Prevotella_bivia 1.79 1.79 1.04 #> 6 Sneathia_amnii 4.59 4.59 1.03 #> 7 Streptococcus_agalactiae 0.250 0.250 1.03 #> #> Geometric standard error estimated from 1000 bootstrap replicates.