"CpermSigKnow.MBparams.MBnetwork.MBdataDE.MBgibbsPars" <-
function(params,graph,data,gibbsPars,quiet=FALSE,pcut=.05,nperm=100,contrasts=NULL,...)
{
# Compatibility checks
    isCompat1 <- compatObjects(params,graph)
    isCompat2 <- compatObjects(graph,data)
    isValid <- validObject(gibbsPars)

# Run Gibbs sampler on all data
    if(!quiet) cat("\nRunning gibbs sampler on all data...\n\n")
    graph <- calcCPDs(x=graph,data=data,contrasts)
    graph <- CgibbsSampler(x=graph,params=params,gibbsPars=gibbsPars,quiet=quiet,...)
    lods.orig <- calcLods(graph)

# Save old Sim
    Sim.old <- params@Sim
    
# Run permutation gibbs sampler
    lods.rand <- matrix(nrow=length(lods.orig),ncol=nperm)
    for(i in 1:nperm) {
      if(!quiet) cat(paste("\nRunning gibbs sampler on permutation ",i,"...\n\n",sep=""))
      permute <- sample(1:params@N)
      for(j in 1:params@m)
	params@Sim[[j]] <- params@Sim[[j]][permute,permute]
      graph <- CgibbsSampler(x=graph,params=params,gibbsPars=gibbsPars,quiet=quiet,...)
      lods.rand[,i] <- calcLods(graph)
      params@Sim <- Sim.old
    }

# Calculate permutation significance
    pval.perm <- apply(cbind(lods.orig,lods.rand),1,function(x) sum(x[2:(nperm+1)][!is.na(x[2:(nperm+1)])] >= x[1])/sum(!is.na(x[2:(nperm+1)])))
    nsig <- sum(pval.perm <= pcut,na.rm=TRUE)

    return(nsig)
}