Network neuroscience models the brain as interacting elements. However, a large number of elements imply a vast number of interactions, making it difficult to assess which connections are relevant and which are spurious. Zalesky et al. (2010) proposed the Network-Based Statistics (NBS), which identifies clusters of connections and tests their likelihood via permutation tests. This framework shows a better trade-off of Type I and II errors compared to conventional multiple comparison corrections. NBS uses General Linear Hypothesis Testing (GLHT), which may underestimate the within-subject variance structure when dealing with longitudinal samples with a varying number of observations (unbalanced samples). We implemented NBR, an R-package that extends the NBS framework adding (non)linear mixed-effects (LME) models. LME models the within-subject variance in more detail, and deals with missing values more flexibly. To illustrate its advantages, we used a public dataset of 333 human participants (188/145 females/males; age range: 17.0-28.4 y.o.) with two (n=212) or three (n=121) sessions each. Sessions include a resting-state fMRI scan and psychometric data. State anxiety scores and connectivity matrices between brain lobes were extracted. We tested their relationship using GLHT and LME models for balanced and unbalanced datasets, respectively. Only the LME approach found a significant association between state anxiety and a subnetwork that includes the cingulum, frontal, parietal, occipital, and cerebellum. Given that missing data is very common in longitudinal studies, we expect that NBR will be very useful to explore unbalanced samples.Significant StatementLongitudinal studies are increasing in neuroscience, providing new insights into the brain under treatment, development, or aging. Nevertheless, missing data is highly frequent in those studies, and conventional designs may discard incomplete observations or underestimate the within-subject variance. We developed a publicly available software (R package: NBR) that implements mixed-effect models into every possible connection in a sample of networks, and it can find significant subsets of connections using non-parametric permutation tests. We demonstrate that using NBR on larger unbalanced samples has higher statistical power than when exploring the balanced subsamples. Although this method is applicable in general network analysis, we anticipate this method being potentially useful in systems neuroscience considering the increase of longitudinal samples in the field.