The hardware and bandwidth for this mirror is donated by METANET, the Webhosting and Full Service-Cloud Provider.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]metanet.ch.

Simulation with non-canonical matrices

Matthew Stephens

2023-10-18

Goal

To try out some simulations that don’t match the canonical covariance matrices and illustrate how the data driven matrices help.

Simple simulation

Here the function simple_sims_2 simulates data in five conditions with just two types of effect:

  1. shared effects only in the first two conditions; and

  2. shared effects only in the last three conditions.

library(ashr)
library(mashr)
set.seed(1)
simdata = simple_sims2(1000,1)
true.U1 = cbind(c(1,1,0,0,0),c(1,1,0,0,0),rep(0,5),rep(0,5),rep(0,5))
true.U2 = cbind(rep(0,5),rep(0,5),c(0,0,1,1,1),c(0,0,1,1,1),c(0,0,1,1,1))
U.true  = list(true.U1 = true.U1, true.U2 = true.U2)

Simple simulation

Run 1-by-1 to add the strong signals and ED covariances.

data   = mash_set_data(simdata$Bhat, simdata$Shat)
m.1by1 = mash_1by1(data)
strong = get_significant_results(m.1by1)
U.c    = cov_canonical(data)
U.pca  = cov_pca(data,5,strong)
U.ed   = cov_ed(data,U.pca,strong)

# Computes covariance matrices based on extreme deconvolution,
# initialized from PCA.
m.c    = mash(data, U.c)
m.ed   = mash(data, U.ed)
m.c.ed = mash(data, c(U.c,U.ed))
m.true = mash(data, U.true)
  
print(get_loglik(m.c),digits = 10)
print(get_loglik(m.ed),digits = 10)
print(get_loglik(m.c.ed),digits = 10)
print(get_loglik(m.true),digits = 10)

The log-likelihood is much better from data-driven than canonical covariances. This is good! Indeed, here the data-driven fit is very slightly better fit than the true matrices, but only very slightly.

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.