2.3.2.1 Modify the number of samples and their library sizes

To change the library sizes, set the lib.size argument equal to a vector with the target library sizes. The length of this vector becomes the new number of samples. The advantage of the parametric mode is that it handles coordinated changes automatically. As a result, when changing library sizes in parametric mode, there’s no need for the user to provide the marginal totals of non-zero cells explicitly.

For example, to generate data having library sizes that are uniformly distributed between 1,000 and 10,000 while changing the number of samples to 700, we execute the code

new.lib.size <- sample(1000:10000, size=700, replace=TRUE)

count.ibd.modified <- MIDASim.modify(
  count.ibd.setup, 
  lib.size = new.lib.size   # other arguments left NULL
  )

2.3.2.2 Modify taxa features by mean relative abundances or location parameters

Specifying new values for mean relative abundances (mean.rel.abund) or location parameters (gengamma.mu) can modify taxa features.

You can alter taxon‑level characteristics in either of two equivalent ways:

• mean.rel.abund provide supply a new vector of target mean relative abundances. MIDASim will calculate the corresponding location parameters of the generalized‑gamma model, then updates all downstream probabilities and sparsity counts.
• gengamma.mu - provide the mean specification and shift \(\mu_j\) directly. This is mathematically identical to changing the means on the log scale.

Because the two arguments modify the same quantity, specify one or the other, not both.

For example, if we want to increase \(log(p_j)\) by \(\beta=0.1\) for the first 10 taxa (e.g. modify these taxa according to a compositional model), the following is the code for generating a set of target relative abundances as described,

beta <- 0.1
new.mean.rel.abund <- count.ibd.setup$mean.rel.abund
new.mean.rel.abund[1:10] <- exp(beta) * new.mean.rel.abund[1:10]
new.mean.rel.abund <- new.mean.rel.abund / sum(new.mean.rel.abund)  # renormalize

then the target taxa relative abundances can be used in MIDASim.modify() as the following,

count.ibd.modified <- MIDASim.modify(
  count.ibd.setup, 
  mean.rel.abund = new.mean.rel.abund
  )

Alternatively, the users can also directly change the location parameters \(\mu_j\) by

\[\mu_j \rightarrow \mu_j + \beta \]

count.ibd.modified <- MIDASim.modify(
  count.ibd.setup, 
  gengamma.mu = count.ibd.setup$mu.est + beta
  )

2.3.2.3 Modify taxa features by individual‑specific relative abundances

MIDASim v2.0 allows you to pass an individual-by-taxon matrix of expected relative abundances (individual.rel.abund). This is the most straightforward way to encode subject‑level covariate effects (binary or continuous, one or many) in a single step.

Illustration with two covariates:

• \(Y1\) is binary: 50 controls (\(Y=0\)) and 50 cases (\(Y=1\)). Cases get a uniform 0.1 log-fold increase in the first 10 taxa.
• \(Y2\) is continuous: 100 draws from a standard normal. Each unit increase in \(Y2\) affects the second 10 taxa (11-20) by taxon-specific log-fold changes sampled from uniform distribution, Unif(−1, 1).

Y1 <- rep(c(0, 1), each = 50)  # binary case-control
Y2 <- rnorm(100)   # continuous
Y.all <- cbind(Y1, Y2)

n.taxa <- length(count.ibd.setup$mean.rel.abund)
beta1 <- beta2 <- rep(0, n.taxa) # effect on each taxon
beta1[1:10] <- rep(0.1, 10)
beta2[11:20] <- runif(10, -1, 1)
beta.all <- cbind(beta1, beta2)

new.individual.rel.abund <- matrix(count.ibd.setup$mean.rel.abund, 
                                   nrow = 100, ncol = n.taxa, byrow = T)
new.individual.rel.abund <- exp(Y.all %*% t(beta.all))  * new.individual.rel.abund  # introduce signals
new.individual.rel.abund <- new.individual.rel.abund/rowSums(new.individual.rel.abund)  # normalize

Feed this matrix to MIDASim.modify() as the following,

count.ibd.modified <- MIDASim.modify(
  count.ibd.setup, 
  individual.rel.abund = new.individual.rel.abund
  )

The resulting setup object now carries subject‑specific compositions that reflect the desired subject-level signals and is ready for data generation with MIDASim().

2.3.3.1 Modify the number of samples and their library sizes

To change the library sizes, set the lib.size argument equal to a vector with the target library sizes. The number of samples can be modified implicitly by specifying library sizes with a vector length equal to the desired sample size.

Note that, if we wish to change the library sizes, the proportion of non-zero taxa would also change. Typically, larger library sizes are accompanied by fewer zero cells. However, since the model is nonparametric, there is no simple parametric relationship between the changes in zero cells and relative abundances. MIDASim.modify() requires the user to specify the marginal totals of non-zeros together when a change in lib.size is made. This can be achieved by providing the marginal proportions using the sample.1.prop and taxa.1.prop arguments within the function.”

An example way determine the appropriate number of non-zero cells for the new library size is to fit a Shape Constrained Additive model (SCAM) model \[ log10(Z_{i.}) = f\left(log_{10}(N_{i})\right)+\epsilon_i\] where \(f\) is a monotone smoothing spline, \(Z_{i\cdot}\) is the number of non-zero cells in each sample, \(\epsilon_i\) is the error term. We could adjust the number of non-zero cells for each taxon in a way that satisfies the constraint that the product of sample.1.prop and the number of taxa should be equal to the product of taxa.1.prop and the number of samples,

\[Z_{i.}\times J = Z_{.j}\times n\]

For example, if we want to generate data with library sizes uniformly distributed between 1,000 and 10,000 while changing the number of samples to 700 and allowing the proportion of non-zero cells to adjust according to the SCAM model, we can start by generating the marginal totals of non-zero cells and transforming the totals to proportions using the following:

new.lib.size <- sample(1000:10000, size=700, replace=TRUE)

obs.sample.1.ct <- rowSums(count.ibd > 0)
xvar <- log10(rowSums(count.ibd))
scamfit.non0 <- scam::scam( log10(obs.sample.1.ct) ~  s( xvar, bs = "mpi" ))
sample.1.ct <- 10^(predict(scamfit.non0,
                           newdata = data.frame(xvar = log10(new.lib.size) )) )

n.taxa <- ncol(count.ibd)
new.sample.1.prop <- sample.1.ct/n.taxa
new.taxa.1.prop <- fitted$taxa.1.prop * (sum(new.sample.1.prop) * n.taxa / 700) / sum(fitted$taxa.1.prop)

Then we can use MIDASim.modify() by

count.ibd.modified <- MIDASim.modify(
  count.ibd.setup, 
  lib.size = new.lib.size,
  sample.1.prop = new.sample.1.prop,
  taxa.1.prop = new.taxa.1.prop
  )

2.3.3.2 Modify features of taxa by mean relative abundances or proportion of non-zero cells

If the mean relative abundances are modified, adjustments to the values of non-zero relative abundances may be necessary. Let \(\hat{p}_j\) represent the mean relative abundances of taxa (mean.rel.abund), \(\hat{p}^{(1)}_j\) denote the mean relative abundances of taxa among non-zero samples (mean.rel.abund.1), and \(\hat{\delta}_j\) represent the proportion of non-zero cells for taxa (taxa.1.prop).

These three quantities are related by \[\hat{p}_j = \hat{\delta}_j\hat{p}^{(1)}_j,\]

For example, if \(\hat{p}_j\) is altered while keeping \(\hat{\delta}j\) constant, MIDASim calculates the mean relative abundance of non-zero cells as \(p^{(1)}_j=p_j/\delta_j\). Subsequently, it determines the value \(\alpha_j\) for each taxon in a way that ensures that \({ \pi^\alpha{i,j} | \pi{ij}>0 }\) has a mean equal to \(p^{(1)}_j\) for each taxon.

Suppose we want to increase \(log(p_j)\) by \(\beta=0.1\) for the first 10 taxa (e.g. modify these taxa according to a compositional model). The following is the code for generating a set of target relative abundances as described,

beta <- 0.1
new.mean.rel.abund <- count.ibd.setup$mean.rel.abund
new.mean.rel.abund[1:10] <- exp(beta) * new.mean.rel.abund[1:10]
new.mean.rel.abund <- new.mean.rel.abund / sum(new.mean.rel.abund)

then the target taxa relative abundances can be used in MIDASim.modify() as the following,

count.ibd.modified <- MIDASim.modify(count.ibd.setup,
                                     mean.rel.abund = new.mean.rel.abund)

Users also have the option to modify the proportion of non-zero cells for all taxa. However, because this modification affects the overall total count of non-zero cells throughout the entire dataset, users are required to specify the marginal proportions of non-zero cells for all samples.

Suppose we wish to increase the proportion of non-zero cells for first 10 taxa by 0.05. We can achieve this by

new.taxa.1.prop <- count.ibd.setup$taxa.1.prop 
new.taxa.1.prop[1:10] <- new.taxa.1.prop[1:10] + 0.05

Then to ensure that appropriate proportions of non-zero cells are provided for all samples, here we adjust them accordingly by

r <- sum(new.taxa.1.prop) / sum(count.ibd.setup$taxa.1.prop)
new.sample.1.prop <- count.ibd.setup$sample.1.prop * r

then the target proportion of non-zero cells for taxa can be used in MIDASim.modify() as the following,

count.ibd.modified <- MIDASim.modify(
  count.ibd.setup,
  taxa.1.prop = new.taxa.1.prop,
  sample.1.prop = new.sample.1.prop
  )

Users can also simultaneously adjust both the mean relative abundances and the proportion of non-zero cells for taxa using the following code by

count.ibd.modified <- MIDASim.modify(
  count.ibd.setup,
  mean.rel.abund = new.mean.rel.abund,
  taxa.1.prop = new.taxa.1.prop,
  sample.1.prop = new.sample.1.prop
  )