--- title: "Transition Networks from Time Series: Estimation and Inference with Nestimate" author: - name: "Mohammed Saqr" affiliation: "University of Eastern Finland" - name: "Sonsoles López-Pernas" affiliation: "University of Eastern Finland" - name: "Manuel J. Gómez" affiliation: "University of Murcia" output: rmarkdown::html_vignette: toc: true number_sections: true vignette: > %\VignetteIndexEntry{Transition Networks from Time Series: Estimation and Inference with Nestimate} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} --- ```{r setup, include=FALSE} knitr::opts_chunk$set( collapse = TRUE, comment = "#>", fig.align = "center", fig.width = 7, fig.height = 4.5, dpi = 110, out.width = "100%" ) library(tsn) has_nestimate <- requireNamespace("Nestimate", quietly = TRUE) has_cograph <- requireNamespace("cograph", quietly = TRUE) && requireNamespace("igraph", quietly = TRUE) ``` # Introduction A student rates their momentary anxiety on a 0-100 scale, again and again over a term. That record is a numeric time series: one quantity, measured repeatedly in time. A transition network asks a different question of the same record. Once the continuous scale is reduced to a few qualitative *states* such as `low`, `mid`, and `high`, the network's nodes are those states and its directed edges the probability of passing from one to the next. The series says how high a measurement was; the network says what a measurement is likely to be followed by. Recovering the network from the series takes two steps: reducing the continuous series to a sequence of states, and estimating the state-to-state transitions. There is a third step that is easy to skip, and it is the one this vignette keeps in view: deciding what the fitted model may legitimately claim. The **tsn** package performs the first two steps in a single call and hands the result to **Nestimate**, which supplies the third. `ts_tna()` discretizes a series, estimates the transition network, and returns a native Nestimate model, so Nestimate's descriptive and inferential verbs apply to it with no conversion step in between. This vignette follows one analysis through the whole path: building the model, reading its structure, and, where the data support it, testing that structure with bootstrap confidence intervals, centrality stability, a Markov-order test, and a permutation comparison of two periods. The distinction between description and inference organizes everything below. Descriptive summaries such as centralities, state distributions, entropy, and pruning ask what a fitted network looks like, and a network fitted to a single series answers them. Inference asks how far the estimate would move under a different sample, and a single series is a single sample: there is nothing to resample. The dividing line is the unit of analysis, and knowing where it falls is what keeps a transition model from being read for more than it can support. Two bundled data sets carry the analysis. `esm_srl` is an experience-sampling study in which 41 students rated momentary self-regulation, motivation, and anxiety on 0-100 scales; `srl` is a balanced daily panel in which 36 students reported nine self-regulated-learning indicators for 156 occasions each. ```{r libs, eval=has_nestimate} library(tsn) library(Nestimate) data(esm_srl) data(srl) ``` # From measurements to a model `ts_tna()` takes a data frame and the name of one numeric column, cuts that column into states, and estimates the transition network in one call. Here it runs on the 79 momentary `anxiety` reports of one student, Jamal, from the `esm_srl` study. ```{r build, eval=has_nestimate} jamal <- subset(esm_srl, name == "Jamal") model <- ts_tna( jamal, series = "anxiety", labels = c("low", "mid", "high") ) model ``` The printed object is the complete model: a row-stochastic weight matrix whose entry *(i, j)* is the probability of moving to state *j* given that the system is currently in state *i*, together with the initial-state distribution. Each row sums to one because from any state the process must go somewhere. One entry already rewards reading: after a high-anxiety report, the most likely next state is `low` (0.423), so for this student high anxiety resolves sharply rather than escalating. ## The bridge: this object *is* a Nestimate model Everything that follows rests on one property of the returned object. ```{r class, eval=has_nestimate} class(model) ``` `ts_tna()` does not return a private structure that must be converted before Nestimate can read it. It returns an object that inherits from **`netobject`**, Nestimate's core network class, and Nestimate's verbs dispatch on that class directly. There is no export step and no adapter call between the model and the tests applied to it, which is what lets a single analysis move from construction to inference without leaving the object behind. ## Four constructors, one interface The same discretized sequence supports four definitions of an edge, and the four constructors take identical arguments and return the same class. `ts_tna()` gives relative frequencies, so an edge is the conditional probability of the next state. `ts_ftna()` gives the raw transition counts behind those probabilities. `ts_atna()` applies attention weighting, and `ts_cna()` counts undirected co-occurrence rather than directed movement. Changing the edge definition means changing only the constructor. ```{r builders-run, eval=has_nestimate} ts_ftna(jamal, series = "anxiety", labels = c("low", "mid", "high")) ``` The remainder of the vignette uses `ts_tna()`, because the inferential tests are written for probabilities. # One series, and what it can support Asking a single-series model for a bootstrap returns a warning rather than a result: > A network with one long sequence is not recommended and can't be validated > using bootstrap and other confirmatory testings. The warning is correct. A bootstrap resamples *sequences* to estimate how far a quantity would move under a different sample of them; one series is one sequence, so there is nothing to resample and no sampling distribution to form. The permutation and stability tests rest on the same logic. This is a statement about the unit of analysis, not a dead end. A long series contains many transitions; what it lacks is a way to group them into units that can be resampled. The `segment` argument supplies that unit by cutting the series into shorter blocks, each treated as its own sequence. This is the block bootstrap, and tsn provides its two standard forms. The series below is one student's 156 daily `effort` reports from the `srl` panel. Two features of the call matter before the segmenting itself. States are cut at 40 and 70 on the 0-100 scale with `discretization = "threshold"` rather than at the default tertiles, so the states mean the same thing outside this sample; the next section returns to why that matters. And segmentation is applied *after* discretization, so every block inherits the same cut points and no block is ever discretized against itself. ```{r segment-setup, eval=has_nestimate} students <- subset(srl, !is.na(effort)) erik <- subset(students, name == "Erik") whole <- ts_tna( erik, value = "effort", time = "day", discretization = "threshold", breaks = c(40, 70), labels = c("low", "moderate", "high") ) ``` Partitioning cuts the 156 observations into consecutive, non-overlapping blocks. At `segment = 10` there are sixteen of them. ```{r segment-blocks, eval=has_nestimate} blocks <- ts_tna( erik, value = "effort", time = "day", segment = 10, discretization = "threshold", breaks = c(40, 70), labels = c("low", "moderate", "high") ) bootstrap_network(blocks, iter = 500, seed = 2026) ``` Sixteen sequences, and the bootstrap now runs. Partitioning is not free: every cut destroys the transition that straddles it, so sixteen blocks cost fifteen of the 155 transitions, about ten percent, and the estimate shifts slightly. Sliding the window instead of partitioning keeps all of them. At `segment = 2, overlap = TRUE` the blocks are the consecutive lag-1 pairs, so every transition appears exactly once and the network is *identical* to the unsegmented one. ```{r segment-lag1, eval=has_nestimate} lag_one <- ts_tna( erik, value = "effort", time = "day", segment = 2, overlap = TRUE, discretization = "threshold", breaks = c(40, 70), labels = c("low", "moderate", "high") ) all.equal(as.matrix(lag_one), as.matrix(whole)) ``` ```{r segment-lag1-boot, eval=has_nestimate} bootstrap_network(lag_one, iter = 500, seed = 2026) ``` The two schemes trade the point estimate against independence. Sliding blocks preserve the estimate exactly and yield the most units; partitioned blocks respect whatever dependence exists inside a block, at the price of the boundary transitions lost to the cuts. On this series the two run neck and neck: mean interval widths of 0.253 for the partitioned blocks against 0.254 for the lag-1 pairs, with a single edge crossing the stability criterion under partitioning and none under sliding. That is the real lesson of the exercise: 156 volatile observations of one person yield intervals a quarter of a probability unit wide however they are resampled. Where the two schemes do differ, prefer partitioning, because resampling lag-1 pairs assumes exactly the first-order Markov property the model itself assumes, and the Markov-order test below rejects that property for these data. Blocks of 10 to 30 retain roughly 90% to 97% of the transitions, which is usually the right trade; blocks of 2 retain only half and are better avoided. None of this manufactures information. Segmenting a single series asks how stable one person's own dynamics are, and says nothing about anyone else. A claim about people in general needs people in general, so the rest of the vignette uses all 36 of them. ## Choosing states that mean something The default `discretization = "quantile"` cuts at the tertiles, forcing each of the three states to hold a third of the observations, so the state sizes carry no information, by construction. That is why every model built from the `srl` panel uses fixed thresholds of 40 and 70 on the 0-100 effort scale instead. Fixed breaks make the states interpretable outside the sample, since below 40 is low effort wherever it occurs, and they give every model one shared alphabet, which is what makes the period comparison below legitimate. (The grouped `esm_srl` model at the end keeps the quantile default; grouping pools the series before cutting it, so its groups share an alphabet without breaks being named.) ```{r many, eval=has_nestimate} effort_model <- ts_tna( students, value = "effort", id = "name", time = "day", discretization = "threshold", breaks = c(40, 70), labels = c("low", "moderate", "high") ) effort_model ``` The states now mean something outside the data, the same cut points apply to every comparison made later, and the state distribution is worth reading. ```{r distribution, eval=has_nestimate} state_distribution(effort_model) ``` High effort is the panel's most common state at 41.7% of the days, and it is also the most persistent: the `high -> high` probability of 0.575 is the largest entry in the matrix. Effort, once high, tends to continue, and the asymmetry runs one way: a low-effort day returns to `low` at only 0.415. ## One network per person The pooled model estimates a single network across everyone. `series_networks()` splits it into one complete model per student, each rebuilt on the same alphabet so node sets stay aligned and comparable. ```{r per-person, eval=has_nestimate} per_person <- series_networks(effort_model) head(summary(per_person)) ``` Each row indexes a full `ts_tna` model that prints and plots on its own: the models are the list elements, and this is the tidy index over them. Because the cut points are fixed scale positions rather than sample quantiles, a student's personal network describes their own effort, not their rank within the panel. The index already shows individual structure: one of the first six students never makes two of the nine transitions at all, so their network carries seven edges rather than nine. # Describing the model ## Centrality ```{r centrality, eval=has_nestimate} net_centrality(effort_model) ``` The message notes that centralities exclude self-transitions by default, which matters for a transition network: the diagonal is usually the largest entry, so a state can be highly persistent while receiving little traffic from elsewhere. Excluding the diagonal measures how central a state is to *movement* between states rather than to staying put. Here `high` and `moderate` receive similar between-state traffic (in-strengths of 0.60 and 0.57) while `low` receives markedly less (0.38); the betweenness column is uniformly zero because in a fully connected three-state network no state is an obligatory bridge. Passing `loops = TRUE` includes the diagonal. ## How much of the sequence is structure, and how much is noise ```{r entropy, eval=has_nestimate} transition_entropy(effort_model) ``` The entropy rate is normalized against `log2(n_states)`, so 0 is a perfectly deterministic process and 1 means every next state is equally likely from wherever the process currently sits. A value of 1 requires more than an uninformative current state: a lopsided but state-independent process still scores well below it. Cut at the tertiles, these same effort reports score 0.96; cut at 40 and 70 they score 0.93. The fixed thresholds recover slightly more structure, and either way the message is sobering: day-to-day effort is mostly noise around a weak first-order signal, which is worth knowing before any edge is over-interpreted. # Pruning Not every edge deserves a place in the picture. `net_prune()` returns a model of the same class with weak edges removed, so its result flows straight back into any verb above. ```{r prune, eval=has_nestimate} net_prune(effort_model, method = "threshold", threshold = 0.3) ``` Four edges fall below the threshold, and they are not a random four: every *descending* route (`moderate -> low`, `high -> low`, `high -> moderate`) disappears, together with the direct `low -> high` leap. What survives is the ascending ladder (`low -> moderate` at 0.314, `moderate -> high` at 0.329) and the three persistence loops. Effort climbs stepwise and decays only weakly. Prune after testing rather than before, however. A threshold ranks edges by size, and size is not the same as interest: a small but reliably estimated edge can matter more than a large but uncertain one. The bootstrap in the next section shows exactly this happening. # Bootstrap: which edges survive resampling? `bootstrap_network()` resamples the sequences with replacement, re-estimates the network on each resample, and reports a confidence interval for every edge. ```{r bootstrap, eval=has_nestimate} boot <- bootstrap_network(effort_model, iter = 200, seed = 2026) boot ``` `summary()` returns the tidy, one-row-per-edge table. ```{r bootstrap-summary, eval=has_nestimate} summary(boot) ``` Five of the nine edges clear the stability criterion, and four of them are edges pruning also kept: the `moderate` and `high` persistence loops and the two ascending steps, `low -> moderate` and `moderate -> high`. The fifth is where the two verdicts part. The bootstrap keeps `high -> moderate` (0.260), a single step down from the top of the scale that pruning had discarded as too small, and it drops `low -> low` (0.415), the large loop pruning kept. `low -> low` is the outgoing edge of the rarest state (`low` fills 24.4% of the days), so it rests on the fewest transitions and just misses the criterion (p = 0.060). This is the pruning caution made concrete: ranking edges by size and ranking them by stability need not agree, and here they trade one edge in each direction. The four edges that fail resampling are the long-range moves, the `low -> high` leap and the two descents into `low` (`moderate -> low`, `high -> low`), together with the `low -> low` near miss. This is what a moderate panel buys: the ascending backbone and the strong loops are firm, while the rarest-state and long-range edges are not. # Stability: would the same states rank first in a smaller study? `centrality_stability()` repeatedly drops a proportion of the sequences, re-computes the centralities, and reports the largest proportion that can be dropped while the ordering still tracks the original, a correlation-stability coefficient. ```{r stability, eval=has_nestimate} set.seed(2026) centrality_stability(effort_model, iter = 100) ``` The in-strength ordering holds to a 0.60 drop, acceptable against the conventional bars of 0.5 for acceptable and 0.7 for excellent, but the out-strength ordering fails the bar at 0.40. With 36 sequences, claims about which state *sends* the most between-state traffic should not be made from this panel; claims about which state receives it are on firmer ground. This is the stability check doing its job, and the honest report includes the number that failed. # Markov order: is a first-order model enough? Every transition network here assumes the Markov property, that the next state depends on the current one and nothing earlier. That assumption is testable, by comparing the fit of models of increasing order. ```{r markov, eval=has_nestimate} set.seed(2026) markov_order_test(effort_model) ``` The test disagrees with the first-order model. AIC and the permutation likelihood-ratio test prefer order 3 and BIC order 2, meaning yesterday carries information about tomorrow beyond what today already conveys. This does not invalidate the first-order network, but it bounds the claim: the network is a summary of one-step behaviour, not a full description of the process. The natural follow-up is a higher-order network. ```{r hon, eval=has_nestimate} build_hon(effort_model, max_order = 2) ``` The higher-order construction keeps every node at order 1: it finds no individual context that shifts the next-step distribution enough to justify splitting the node. The two criteria disagree because they ask different questions. The likelihood test detects aggregate improvement pooled across all sequences, while the higher-order rule applies a stricter, per-node test. Both results are real, and both are worth reporting: there is higher-order structure in aggregate, but not localized in any single state. # Comparing two periods The tests so far describe one network. `permutation()` compares two, shuffling the group labels to build a null distribution for every edge. The `srl` panel covers 156 occasions per student; splitting at occasion 78 contrasts the first half of the course with the second, with the same 36 students on both sides. Both models are built with the same fixed `breaks`, which matters more than it appears: had each period been discretized on its own quantiles, the two networks would rest on different cut points, and any difference between them could reflect the thresholds rather than behaviour. ```{r two-models, eval=has_nestimate} first_half <- ts_tna( subset(students, day <= 78), value = "effort", id = "name", time = "day", discretization = "threshold", breaks = c(40, 70), labels = c("low", "moderate", "high") ) second_half <- ts_tna( subset(students, day > 78), value = "effort", id = "name", time = "day", discretization = "threshold", breaks = c(40, 70), labels = c("low", "moderate", "high") ) comparison <- permutation(first_half, second_half, iter = 2000, seed = 2026) comparison ``` ```{r perm-table, eval=has_nestimate} summary(comparison) ``` Nothing separates the halves: no edge comparison comes near significance (the smallest p-value is 0.49) and the largest estimated difference on any edge is 0.033. The compositions agree just as closely: ```{r composition, eval=has_nestimate} state_distribution(first_half) state_distribution(second_half) ``` High-effort days make up 42.0% of the first half and 41.5% of the second; every proportion moves by less than a percentage point. Level and dynamics are both stationary across the course: whatever effort habits these students have, they are set before the midpoint and do not drift. A null this uniform is a finding: it licenses pooling the whole course into the single model above, which a real first-to-second-half shift would have forbidden. Two cautions bound how far that carries. Failing to reject nine tests is not the same as proving the dynamics identical, and the balanced design helps the comparison: the same 36 students stand on both sides of the cut, so the null is not an artefact of different people in different periods. # More than two groups: one model per condition Building two models by hand is manageable when there are two. When the variable that separates the groups is already a column in the data, the `group` argument builds one network per level in a single call, and it scales past two, which the pairwise form does not. The `esm_srl` study carries a `day_type` column marking weekday and weekend reports, so momentary effort can be split by day type. ```{r grouped, eval=has_nestimate} by_day <- ts_tna( subset(esm_srl, !is.na(effort)), value = "effort", id = "name", time = "occasion", group = "day_type", labels = c("low", "mid", "high") ) by_day ``` Two things happen here that would be tedious to arrange by hand. The states are cut from the *pooled* series before the split, so both networks rest on one set of thresholds and share a node set, the same discipline the fixed `breaks` enforced earlier, applied automatically. And the result is a Nestimate `netobject_group`, so the compatible grouped verbs take it directly: `net_prune()`, `state_distribution()`, `net_centrality()`, `compare_model()`, `permutation()`, `mosaic_plot()`, and their neighbours, while `plot()` with a `group` argument draws one selected network. Not every verb applies: those needing a precision matrix or a clustering attribute reject a transition model. `as.data.frame()` gives the tidy view of how much data stands behind each network: ```{r grouped-index, eval=has_nestimate} as.data.frame(by_day, what = "groups") ``` The sequence counts repay a second look: the weekday network holds 2,033 observations in 238 sequences, the weekend network 783 in 233. That is not padding. A run is cut at every change of student *and* at every change of day type, and the two cuts are not the same. A change of student separates measurements that were never consecutive; a change of day type falls between two genuinely adjacent measurements, and the transition across it is excluded from both networks because it belongs to neither day type. A working week of roughly two reports a day therefore forms one sequence of eight to ten observations, a weekend one of three to four, and the Friday-to-Saturday transition is counted in neither. ```{r grouped-distribution, eval=has_nestimate} state_distribution(by_day) ``` The compositions are nearly identical: each state holds close to a third of the reports on both day types, which pooled quantile cutting guarantees for the whole panel, though not for each group taken alone. Whatever distinguishes weekends, it is not how much high-effort time they contain. The difference, if any, must live in the dynamics, and the shared alphabet is what makes that question answerable: ```{r grouped-permutation, eval=has_nestimate} day_comparison <- permutation(by_day, iter = 1000, seed = 2026) summary(day_comparison) ``` One difference holds up under the null: `high -> mid` runs at 0.254 on weekdays and 0.342 on weekends (p = 0.023). Its mirror image, weekend `high -> high` retention falling to 0.534 from 0.631, sits exactly on the boundary (p = 0.046) and should be read as suggestive rather than established. Read together, the pattern is coherent: on weekends a high-effort spell is more likely to soften into the mid range and less likely to hold. This is the reverse of the two-period result above: there, composition was free to move and nothing moved; here, composition is pinned by construction, so the dynamics carry the difference. A mosaic plot shows both transition structures at a glance, one panel per day type, with the weekend `high -> mid` cell visibly the larger of the pair: ```{r grouped-mosaic, eval=has_nestimate, fig.width=9, fig.height=5, fig.alt="Mosaic panels of the weekday and weekend transition networks, in which tile area encodes transition probability, showing the weekend high-to-mid cell visibly larger than its weekday counterpart."} mosaic_plot(by_day) ``` Because both networks share an alphabet, the pair can also be set side by side and compared as matrices: ```{r grouped-compare, eval=has_nestimate} compare_model(by_day, i = "weekday", j = "weekend") ``` A Pearson correlation of 0.966 with a maximum absolute difference of 0.097, the `high -> high` retention gap, confirms the same reading: the two networks are close, and where they part is in how long high effort lasts. Pruning applies group-wise and returns a grouped model, so the collection can be carried through a workflow without being taken apart: ```{r grouped-prune, eval=has_nestimate} net_prune(by_day, method = "threshold", threshold = 0.3) ``` One caution generalizes beyond this example. Every level of the grouping column gets a network, however little data stands behind it, so read the `sequences` and `observations` columns of the index before reading any edges. A network estimated from a handful of transitions is a level of a column, not a result. # Plotting a model Any single network has a `plot()` method, and a grouped model draws one selected group at a time via its `group` argument. Here is the pooled effort model: three states as nodes, transition probabilities as directed edges, and self-loops for persistence. ```{r plot, eval=has_nestimate && has_cograph, fig.width=9, fig.height=6, fig.alt="Transition network of daily effort states, with low, moderate, and high nodes joined by directed edges whose thickness reflects transition probability and self-loops that mark each state's persistence."} plot(effort_model) ``` The heavy self-loops confirm what the state distribution reported numerically: each state mostly returns to itself, and movement between states is the smaller part of the process. The plotting vignette covers the plot types and their arguments in full. # Where the line falls Everything descriptive in this vignette, from centralities, state distributions, entropy, and pruning to per-student splits and the grouped models, runs on a model built from a single series. Everything inferential, meaning the bootstrap, stability, the Markov-order test, and the permutation test, requires many sequences, because resampling and permutation both treat the sequence as the unit of analysis. Three practices follow. Build with `id` from the start if you intend to test anything, so the sequences exist when a test asks for them. Fix your `breaks` before you compare anything, or pass `group` and let the pooled discretization fix them for you, so a difference between networks is a difference in behaviour and not in thresholds. And when a comparison rests on few sequences, say so: the `sequences` column exists to be read, and so does a stability coefficient that comes back at 0.40.