Julia GLM.jl and MixedModels.jl based implementation of the cluster-based permutation test for time series data, powered by JuliaConnectoR.
As of March 2025, Google Colab supports Julia. This means
{jlmerclusterperm} just works out of the box. Try it out in a demo
notebook
that runs some of the code from the Ito et al. 2018 case study
vignette.
Install the released version of jlmerclusterperm from CRAN:
install.packages("jlmerclusterperm")
Or install the development version from GitHub with:
# install.packages("remotes")
# remotes::install_github("yjunechoe/jlmerclusterperm")
install.packages("jlmerclusterperm", repos = c('https://yjunechoe.r-universe.dev', 'https://cloud.r-project.org'))
Using jlmerclusterperm requires a prior installation of the Julia
programming language, which can be downloaded from either the official
website or using the command line utility
juliaup. Julia version >=1.8 is
required and
1.9
or higher is preferred for the substantial speed improvements.
Before using functions from jlmerclusterperm, an initial setup is
required via calling jlmerclusterperm_setup(). The very first call on
a system will install necessary dependencies (this only happens once and
takes around 10-15 minutes).
Subsequent calls to jlmerclusterperm_setup() incur a small overhead of
around 30 seconds, plus slight delays for first-time function calls. You
pay up front for start-up and warm-up costs and get blazingly-fast
functions from the package.
# Both lines must be run at the start of each new session
library(jlmerclusterperm)
jlmerclusterperm_setup()
See the Get Started page on the package website for background and tutorials.
clusterpermute()
A time series data: vwp_sim, a simulated eyetracking experiment of
looks to a target (elog, the empirical logit of Fixations) by Age
(adult vs. child) and Condition (a related vs. unrelated competitor).
matplot(
unique(vwp_sim$Time), tapply(vwp_sim$elog, list(vwp_sim$Time, interaction(vwp_sim$Age, vwp_sim$Condition)), mean),
type = "l", col = 1:2, lty = rep(1:2, each = 2), lwd = 3, ylab = "Looks to target", xlab = "Time"
)
legend("bottomright", c("Adult Unrelated", "Adult Related", "Child Unrelated", "Child Related"), col = c(1, 1, 2, 2), lty = c(2, 1, 2, 1), lwd = 3)
Preparing a specification object with make_jlmer_spec():
vwp_spec <- make_jlmer_spec(
formula = elog ~ 1 + Condition * Age,
data = vwp_sim,
subject = "Subject", trial = "Item", time = "Time"
)
vwp_spec
Cluster-based permutation test with clusterpermute():
set_rng_state(123L)
CPA <- clusterpermute(
vwp_spec,
threshold = 2,
nsim = 100
)
CPA
Collecting results as data frames with tidy(), e.g., to plot clusters:
clusters <- tidy(CPA$empirical_clusters)
clusters
#> # A tibble: 3 × 7
#> predictor id start end length sum_statistic pvalue
#> <chr> <fct> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 Condition1 1 450 1350 19 -114. 0.00990
#> 2 Age1 1 100 2000 39 595. 0.00990
#> 3 Condition1__Age1 1 650 1250 13 55.9 0.00990
par(mar = c(4, 10, 1, 1))
plot(NA, xlim = range(vwp_sim$Time), ylim = c(0.5, 3.5), yaxt = "n", xlab = "Time", ylab = "")
segments(clusters$start, as.integer(factor(clusters$predictor)), clusters$end, lwd = 15, lend = "butt",
col = ifelse(clusters$pvalue < 0.05, "steelblue", "grey70"))
axis(2, at = 1:3, labels = levels(factor(clusters$predictor)), las = 1)
Including random effects:
vwp_re_spec <- make_jlmer_spec(
formula = elog ~ 1 + Condition * Age +
(1 + Condition | Subject) + (1 + Condition | Item),
data = vwp_sim,
subject = "Subject", trial = "Item", time = "Time"
)
set_rng_state(123L)
clusterpermute(
vwp_re_spec,
threshold = 2,
nsim = 100
)$empirical_clusters
Computing time-wise statistics of the observed data:
empirical_statistics <- compute_timewise_statistics(vwp_spec)
matplot(unique(vwp_sim$Time), t(empirical_statistics), type = "l", lty = 1, lwd = 3, ylab = "t-statistic", xlab = "Time")
abline(h = c(-2, 2), lty = 3)
legend("topright", rownames(empirical_statistics), col = 1:3, lwd = 3)
Identifying empirical clusters:
empirical_clusters <- extract_empirical_clusters(empirical_statistics, threshold = 2)
empirical_clusters
Simulating the null distribution:
set_rng_state(123L)
null_statistics <- permute_timewise_statistics(vwp_spec, nsim = 100)
null_cluster_dists <- extract_null_cluster_dists(null_statistics, threshold = 2)
null_cluster_dists
Significance testing the cluster-mass statistic:
calculate_clusters_pvalues(empirical_clusters, null_cluster_dists, add1 = TRUE)
Iterating over a range of threshold values:
walk_threshold_steps(empirical_statistics, null_statistics, steps = c(1.5, 2, 2.5))
The paper Maris & Oostenveld (2007) which originally proposed the cluster-based permutation analysis.
The JuliaConnectoR package for powering the R interface to Julia.
The Julia packages GLM.jl and MixedModels.jl for fast implementations of (mixed effects) regression models.
Existing implementations of CPA in R (permuco, permutes, etc.) whose designs inspired the CPA interface in jlmerclusterperm.
If you use jlmerclusterperm for cluster-based permutation test with mixed-effects models in your research, please cite one (or more) of the following as you see fit.
To cite jlmerclusterperm:
To cite the cluster-based permutation test:
To cite the Julia programming language:
To cite the GLM.jl and MixedModels.jl Julia libraries, consult their Zenodo pages: