SpiecEasi icon indicating copy to clipboard operation
SpiecEasi copied to clipboard

How to change the output of SpiecEasi into a correlation matrix similar to the downloaded from MENAP ?

Open YuZhang-learner opened this issue 2 years ago • 5 comments

I hope that change the output of SpiecEasi into a correlation matrix similar to the downloaded from MENAP will eventually allow me to calculate robustness with the code in this web web page. How do I get this matrix?

YuZhang-learner avatar Sep 08 '22 19:09 YuZhang-learner

Using glasso, you could convert the sparse inverse covariance matrix to a correlation matrix. An example of this in the README using American Gut data

secor  <- cov2cor(getOptCov(se.gl.amgut))

However, it looks like that MENAP script is just thresholding a correlation matrix anyway, so why not feed in the SpiecEasi network as the 'network.raw' object instead. We've done essentially this robustness analysis on SpiecEasi networks in a few papers now.

zdk123 avatar Sep 08 '22 20:09 zdk123

Thank you, actually I wanted to use SpiecEasi Network, but I got an error when running the following code.

rand.remov.once<-function(netRaw, rm.percent, sp.ra, abundance.weighted=T){ id.rm<-sample(1:nrow(netRaw), round(nrow(netRaw)rm.percent)) net.Raw=netRaw #don't want change netRaw net.Raw[id.rm,]=0; net.Raw[,id.rm]=0; ##remove all the links to these species if (abundance.weighted){ net.stength= net.Rawsp.ra } else { net.stength= net.Raw }

sp.meanInteration<-colMeans(net.stength)

id.rm2<- which(sp.meanInteration<=0) ##remove species have negative interaction or no interaction with others remain.percent<-(nrow(netRaw)-length(id.rm2))/nrow(netRaw)

remain.percent }

rm.p.list=seq(0.05,0.2,by=0.05) rmsimu<-function(netRaw, rm.p.list, sp.ra, abundance.weighted=T,nperm=100){ t(sapply(rm.p.list,function(x){ remains=sapply(1:nperm,function(i){ rand.remov.once(netRaw=netRaw, rm.percent=x, sp.ra=sp.ra, abundance.weighted=abundance.weighted) }) remain.mean=mean(remains) remain.sd=sd(remains) remain.se=sd(remains)/(nperm^0.5) result<-c(remain.mean,remain.sd,remain.se) names(result)<-c("remain.mean","remain.sd","remain.se") result })) }

library(SpiecEasi) library(phyloseq)

data('amgut2.filt.phy') se.gl.amgut2 <- spiec.easi(amgut2.filt.phy, method = 'glasso', lambda.min.ratio = 0.01, nlambda = 20, pulsar.params = list(rep.num = 50))

Weighted.simu<-rmsimu(netRaw=se.gl.amgut2 , rm.p.list=seq(0.05,1,by=0.05), sp.ra=otu_table(amgut2.filt.phy), abundance.weighted=T,nperm=100)`

Error in h(simpleError(msg, call)) : error in evaluating the argument 'x' in selecting a method for function 't': argument of length 0

How should I do?

YuZhang-learner avatar Sep 09 '22 08:09 YuZhang-learner

you're passing in the spiec.easi object, not the re-fit adjacency matrix. Perhaps try getRefit(se.gl.amgut2)

zdk123 avatar Sep 09 '22 16:09 zdk123

There is still an error: image

zdk123

YuZhang-learner avatar Sep 11 '22 08:09 YuZhang-learner

maybe it doesn't like the sparse matrix output? Try as.matrix?

zdk123 avatar Jan 17 '23 23:01 zdk123