Nothing
## ---- include = FALSE---------------------------------------------------------
knitr::opts_chunk$set(
collapse = TRUE,
comment = "#>"
)
## ----setup--------------------------------------------------------------------
library(NBR)
data("voles")
brain_labs <- NBR:::voles_roi
dim(voles)
head(voles)[1:8]
## ----upper triangle-----------------------------------------------------------
nnodes <- length(brain_labs)
tri_pos <- which(upper.tri(matrix(nrow = nnodes, ncol = nnodes)), arr.ind = T)
head(tri_pos)
## ----input network, fig.align = "center"--------------------------------------
library(lattice)
avg_mx <- matrix(0, nrow = nnodes, ncol = nnodes)
avg_mx[upper.tri(avg_mx)] <- apply(voles[-(1:3)], 2, function(x) mean(x, na.rm=TRUE))
avg_mx <- avg_mx + t(avg_mx)
# Set max-absolute value in order to set a color range centered in zero.
flim <- max(abs(avg_mx))
levelplot(avg_mx, main = "Average", ylab = "ROI", xlab = "ROI",
at = seq(-flim, flim, length.out = 100))
## ----NBR-LME, eval = FALSE----------------------------------------------------
# set.seed(18900217)
# before <- Sys.time()
# library(nlme)
# nbr_result <- nbr_lme_aov(net = voles[,-(1:3)],
# nnodes = 16,
# idata = voles[,1:3],
# nperm = 5,
# mod = "~ Session*Sex",
# rdm = "~ 1+Session|id",
# na.action = na.exclude)
# after <- Sys.time()
# show(after-before)
## ----multicore NBR-LME, eval = FALSE------------------------------------------
# set.seed(18900217)
# before <- Sys.time()
# library(nlme)
# library(parallel)
# nbr_result <- nbr_lme_aov(
# net = voles[,-(1:3)],
# nnodes = 16,
# idata = voles[,1:3],
# nperm = 1000,
# nudist = T,
# mod = "~ Session*Sex",
# rdm = "~ 1+Session|id",
# cores = detectCores(),
# na.action = na.exclude
# )
# after <- Sys.time()
# show(after-before)
## ----NBR-LME results----------------------------------------------------------
nbr_result <- NBR:::voles_nbr
show(nbr_result$fwe)
## ----component display, fig.align = "center"----------------------------------
# Plot significant edges
edge_mat <- array(0, dim(avg_mx))
edge_mat[nbr_result$components$Session[,2:3]] <- 1
levelplot(edge_mat, col.regions = rev(heat.colors(100)),
main = "Component", ylab = "ROI", xlab = "ROI")
## ----component cum-pval, fig.height = 3, fig.width = 5, fig.align = "center"----
null_ses_str <- nbr_result$nudist[,2] # Null distribution for Session strength
obs_ses_str <- nbr_result$fwe$Session[,4] # Observed Session strength
nperm <- length(null_ses_str)
cumpval <- cumsum(null_ses_str >= obs_ses_str)/(1:nperm)
# Plot p-value stability
plot(cumpval, type="l", ylim = c(0,0.06), las = 1,
xlab = "Permutation index", ylab = "p-value",
main = "Cumulative p-value for Session strength")
abline(h=0.05, col="red", lty=2)
# Add binomial marginal error
mepval <- 2*sqrt(cumpval*(1-cumpval)/1:nperm)
lines(cumpval+mepval, col = "chartreuse4")
lines(cumpval-mepval, col = "chartreuse4")
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.