The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.
A complete workflow: from sequences to a
group comparison
Mohammed Saqr
Sonsoles López-Pernas
This vignette shows a complete applied workflow: define the sequence
data, fit an LSA model, read the fitted transitions, visualise the main
structure, add uncertainty evidence, start from a raw event log, and
compare groups.
It assumes the method overview in
vignette("intro", package = "lagdynamics"). It uses only a
small plotting subset; the full gallery is in
vignette("plotting", package = "lagdynamics"). The full
confirmatory battery is in
vignette("confirmatory", package = "lagdynamics").
Define the sequence data
engagement contains 138 students observed over 15 weekly
time points. Each row is one student sequence. Each cell is a weekly
engagement state: Active, Average, or
Disengaged.
head(engagement)
#> 1 2 3 4 5 6 7
#> 1 Active Disengaged Disengaged Disengaged Active Active Active
#> 2 Average Average Average Average Average Average Average
#> 3 Average Active Active Active Active Active Active
#> 4 Active Active Active Active Active Average Average
#> 5 Active Active Active Average Active Average Active
#> 6 Average Average Disengaged Average Disengaged Disengaged Average
#> 8 9 10 11 12 13 14 15
#> 1 Active Average Active Average Active <NA> <NA> <NA>
#> 2 Average Average Active Average Average Average Average Average
#> 3 Active Active Average Active Active Active Active Active
#> 4 Active Active Active Average Active Average Active Average
#> 5 Active Active Disengaged Active Active Active Average Average
#> 6 Average Average Disengaged Average Disengaged Average Average Average
This first step fixes the scope of the claim: the empirical units are
students, the events are weekly states, and the ordering is column
order.
Fit the model
lsa() estimates the transition table, the counts
expected under independence, and the adjusted residual for each
transition. The adjusted residual is the local test statistic: positive
values mark over-represented transitions and negative values mark
avoided transitions.
fit <- lsa(engagement)
fit
#> Lag Sequential Analysis - classical (lag 1, directed)
#> 3 states | 1734 transitions | 1870 events | 136 sequences
#> states: Active, Average, Disengaged
#> independence: G² = 618.3, df = 4, p <2e-16
#>
#> Significant transitions (p < 0.05): 7 of 9
#> strongest over-represented (of 3):
#> Active -> Active z = +21.7 ***
#> Disengaged -> Disengaged z = +15.4 ***
#> Average -> Average z = +12.5 ***
#>
#> Initial states:
#> Active 0.382 ████████████████████████
#> Average 0.368 ███████████████████████
#> Disengaged 0.250 ████████████████
Read the fitted transitions
The fitted object is read through verbs. transitions()
returns one row per from -> to edge.
transitions(fit)
#> from to lag count expected prob prob_col adj_res p
#> 1 Active Active 1 459 247 0.698 0.7051 21.663 4.60e-104
#> 2 Average Active 1 153 282 0.204 0.2350 -12.906 4.16e-38
#> 3 Disengaged Active 1 39 122 0.120 0.0599 -10.549 5.11e-26
#> 4 Active Average 1 176 290 0.267 0.2307 -11.319 1.06e-29
#> 5 Average Average 1 458 330 0.610 0.6003 12.453 1.35e-35
#> 6 Disengaged Average 1 129 143 0.397 0.1691 -1.736 8.25e-02
#> 7 Active Disengaged 1 23 121 0.035 0.0719 -12.557 3.64e-36
#> 8 Average Disengaged 1 140 139 0.186 0.4375 0.176 8.60e-01
#> 9 Disengaged Disengaged 1 157 60 0.483 0.4906 15.390 1.90e-53
#> yules_q kappa kappa_z kappa_p lift sign significant
#> 1 0.828 0.4438 19.19 4.68e-82 1.858 over TRUE
#> 2 -0.601 -0.4994 -14.70 6.81e-49 0.543 under TRUE
#> 3 -0.698 -0.7069 -11.46 2.03e-30 0.320 under TRUE
#> 4 -0.534 -0.4242 -12.48 9.39e-36 0.608 under TRUE
#> 5 0.553 0.2275 9.83 7.97e-23 1.386 over TRUE
#> 6 -0.108 -0.1618 -2.96 3.03e-03 0.902 under FALSE
#> 7 -0.826 -0.8272 -13.41 5.14e-41 0.189 under TRUE
#> 8 0.011 -0.0903 -1.66 9.79e-02 1.010 over FALSE
#> 9 0.754 0.3136 13.56 7.34e-42 2.618 over TRUE
For interpretation, it is usually clearer to begin with transitions
that depart from chance.
transitions(fit, significant = TRUE)
#> from to lag count expected prob prob_col adj_res p
#> 1 Active Active 1 459 247 0.698 0.7051 21.7 4.60e-104
#> 2 Average Active 1 153 282 0.204 0.2350 -12.9 4.16e-38
#> 3 Disengaged Active 1 39 122 0.120 0.0599 -10.5 5.11e-26
#> 4 Active Average 1 176 290 0.267 0.2307 -11.3 1.06e-29
#> 5 Average Average 1 458 330 0.610 0.6003 12.5 1.35e-35
#> 6 Active Disengaged 1 23 121 0.035 0.0719 -12.6 3.64e-36
#> 7 Disengaged Disengaged 1 157 60 0.483 0.4906 15.4 1.90e-53
#> yules_q kappa kappa_z kappa_p lift sign significant
#> 1 0.828 0.444 19.19 4.68e-82 1.858 over TRUE
#> 2 -0.601 -0.499 -14.70 6.81e-49 0.543 under TRUE
#> 3 -0.698 -0.707 -11.46 2.03e-30 0.320 under TRUE
#> 4 -0.534 -0.424 -12.48 9.39e-36 0.608 under TRUE
#> 5 0.553 0.227 9.83 7.97e-23 1.386 over TRUE
#> 6 -0.826 -0.827 -13.41 5.14e-41 0.189 under TRUE
#> 7 0.754 0.314 13.56 7.34e-42 2.618 over TRUE
The companion tables keep the transition interpretation grounded in
node volume, starting states, and tablewise tests.
nodes(fit)
#> state outgoing incoming
#> 1 Active 658 651
#> 2 Average 751 763
#> 3 Disengaged 325 320
tests(fit)
#> test statistic df p
#> 1 lrx2 618 4 1.69e-132
#> 2 x2 629 4 7.41e-135
initial(fit)
#> state init_prob
#> 1 Active 0.382
#> 2 Average 0.368
#> 3 Disengaged 0.250
Visualise the fitted model
The default plot is the residual heatmap. It displays the LSA model
directly: rows are source states, columns are target states, and colour
is the adjusted residual.
The residual network shows the same departure-from-chance model as
directed edges.
plot(fit, type = "network")
A TNA-style transition network is a different view of the same
process: edge weights are conditional transition probabilities rather
than residuals.
plot(fit, type = "network", weights = "tna")
Use the plotting vignette for chord diagrams, sunbursts, comparison
heatmaps, and uncertainty plots.
Add uncertainty evidence
An edge claim needs more than a fitted residual.
certainty_lsa() gives analytic credible intervals for
transition probabilities.
cert <- certainty_lsa(fit)
cert
#> <lsa_certainty> (analytic Dirichlet-Multinomial)
#> engine: classical
#> prior: Dirichlet(0.50)
#> CI level: 95% | inference: stability
#> certain edges: 7 of 9
transitions(cert) |> head(4)
#> from to observed prob_observed prob_mean prob_se prob_ci_low
#> 1 Active Active 459 0.698 0.697 0.0179 0.6611
#> 2 Average Active 153 0.204 0.204 0.0147 0.1760
#> 3 Disengaged Active 39 0.120 0.121 0.0180 0.0879
#> 4 Active Average 176 0.267 0.268 0.0172 0.2345
#> prob_ci_high p_value stable adj_res_observed adj_res_stable
#> 1 0.731 5.56e-20 TRUE 21.7 TRUE
#> 2 0.233 6.06e-04 TRUE -12.9 TRUE
#> 3 0.158 9.45e-02 FALSE -10.5 FALSE
#> 4 0.302 1.21e-04 TRUE -11.3 TRUE
bootstrap_lsa() resamples sequences, refits the model,
and reports edge-level variability. It is useful when the population may
be a mixture of different sequence types.
boot <- bootstrap_lsa(fit, R = 100)
boot
#> <lsa_bootstrap>
#> engine: classical
#> level: sequence
#> replicates: 100
#> CI level: 95%
#> stable edges: 7 of 9
transitions(boot) |> head(4)
#> from to observed count_mean count_se count_ci_low count_ci_high
#> 1 Active Active 459 465.2 51.55 378 560
#> 2 Average Active 153 153.9 16.30 126 187
#> 3 Disengaged Active 39 39.1 6.37 27 52
#> 4 Active Average 176 178.4 17.68 148 217
#> adj_res_observed adj_res_mean adj_res_se adj_res_ci_low adj_res_ci_high
#> 1 21.7 21.5 1.46 18.2 23.92
#> 2 -12.9 -13.0 1.41 -15.6 -10.60
#> 3 -10.5 -10.4 1.26 -12.5 -8.16
#> 4 -11.3 -11.3 1.53 -14.1 -8.63
#> adj_res_p_boot adj_res_stable prob_observed prob_mean prob_ci_low
#> 1 0 TRUE 0.698 0.697 0.633
#> 2 0 TRUE 0.204 0.206 0.162
#> 3 0 TRUE 0.120 0.124 0.089
#> 4 0 TRUE 0.267 0.269 0.223
#> prob_ci_high yules_q_observed yules_q_mean yules_q_ci_low yules_q_ci_high
#> 1 0.746 0.828 0.824 0.743 0.871
#> 2 0.256 -0.601 -0.603 -0.690 -0.508
#> 3 0.166 -0.698 -0.690 -0.798 -0.575
#> 4 0.324 -0.534 -0.530 -0.633 -0.418
reliability_lsa() moves from edge evidence to
whole-network evidence by repeated split-half refitting.
reliability_lsa(fit, R = 30)
#> <lsa_reliability>
#> engine: classical
#> replicates: 30
#> weights: prob
#> method: pearson
#> n sequences: 136
#> split-half r: 0.969 (sd = 0.021)
#> 95% CI: [0.925, 0.993]
Start from a raw event log
Real process data often arrives as a long event log. The bundled
group_regulation_long data has one timestamped regulation
action per row, plus an achievement-group label.
head(group_regulation_long)
#> Actor Achiever Group Course Time Action
#> 1 1 High 1 A 2025-01-01 10:27:07 cohesion
#> 2 1 High 1 A 2025-01-01 10:35:20 consensus
#> 3 1 High 1 A 2025-01-01 10:42:18 discuss
#> 4 1 High 1 A 2025-01-01 10:50:00 synthesis
#> 5 1 High 1 A 2025-01-01 10:52:25 adapt
#> 6 1 High 1 A 2025-01-01 10:57:31 consensus
lsa() can recover sequences and fit the model in one
call. actor identifies the sequence, action
identifies the event state, and time orders events. If
there is no actor column, action alone is valid and the log
is treated as one long sequence.
fit_log <- lsa(group_regulation_long, actor = "Actor",
action = "Action", time = "Time")
fit_log
#> Lag Sequential Analysis - classical (lag 1, directed)
#> 9 states | 25533 transitions | 27533 events | 2000 sequences
#> states: adapt, cohesion, consensus, coregulate, discuss, emotion, monitor, plan, synthesis
#> independence: G² = 13203.8, df = 64, p <2e-16
#>
#> Significant transitions (p < 0.05): 72 of 81
#> strongest over-represented (of 23):
#> emotion -> cohesion z = +58.2 ***
#> discuss -> synthesis z = +48.0 ***
#> synthesis -> adapt z = +38.8 ***
#> consensus -> coregulate z = +35.3 ***
#> consensus -> plan z = +32.6 ***
#> ... and 18 more
#>
#> Initial states:
#> consensus 0.214 ████████████████████████
#> plan 0.204 ███████████████████████
#> discuss 0.175 ████████████████████
#> emotion 0.151 █████████████████
#> monitor 0.144 ████████████████
#> cohesion 0.060 ███████
#> synthesis 0.019 ██
#> coregulate 0.019 ██
#> adapt 0.011 █
The fitted probabilities are exposed for downstream
transition-network workflows.
transition_probabilities(fit_log)
#> adapt cohesion consensus coregulate discuss emotion monitor
#> adapt 0.000000 0.2731 0.477 0.0216 0.0589 0.1198 0.0334
#> cohesion 0.002950 0.0271 0.498 0.1192 0.0596 0.1156 0.0330
#> consensus 0.004740 0.0149 0.082 0.1877 0.1880 0.0727 0.0466
#> coregulate 0.016244 0.0360 0.135 0.0234 0.2736 0.1721 0.0863
#> discuss 0.071374 0.0476 0.321 0.0843 0.1949 0.1058 0.0223
#> emotion 0.002467 0.3253 0.320 0.0342 0.1019 0.0768 0.0363
#> monitor 0.011165 0.0558 0.159 0.0579 0.3754 0.0907 0.0181
#> plan 0.000975 0.0252 0.290 0.0172 0.0679 0.1468 0.0755
#> synthesis 0.234663 0.0337 0.466 0.0445 0.0629 0.0706 0.0123
#> plan synthesis
#> adapt 0.0157 0.00000
#> cohesion 0.1410 0.00354
#> consensus 0.3958 0.00758
#> coregulate 0.2391 0.01878
#> discuss 0.0116 0.14098
#> emotion 0.0998 0.00282
#> monitor 0.2156 0.01605
#> plan 0.3742 0.00179
#> synthesis 0.0752 0.00000
initial(fit_log)
#> state init_prob
#> 1 adapt 0.0115
#> 2 cohesion 0.0605
#> 3 consensus 0.2140
#> 4 coregulate 0.0190
#> 5 discuss 0.1755
#> 6 emotion 0.1515
#> 7 monitor 0.1440
#> 8 plan 0.2045
#> 9 synthesis 0.0195
Fit groups
For group-specific models, name the grouping column.
lsa() estimates one model per group under the same state
labels.
gfit <- lsa(group_regulation_long, actor = "Actor", action = "Action",
time = "Time", group = "Achiever")
gfit
#> <lsa_group>
#> engine: classical
#> states: 9 (adapt, cohesion, consensus, coregulate, discuss, emotion, monitor, plan, synthesis)
#> groups: 2
#> - High: 1000 sequences
#> - Low: 1000 sequences
The reading verbs now return a leading group column.
transitions(gfit, significant = TRUE) |> head(6)
#> group from to lag count expected prob prob_col adj_res p
#> 1 High consensus adapt 1 14 38.13 0.00413 0.0979 -4.59 4.45e-06
#> 2 High coregulate adapt 1 20 10.05 0.02237 0.1399 3.27 1.06e-03
#> 3 High discuss adapt 1 48 22.52 0.02396 0.3357 5.88 4.00e-09
#> 4 High emotion adapt 1 5 17.42 0.00323 0.0350 -3.19 1.40e-03
#> 5 High plan adapt 1 4 32.51 0.00138 0.0280 -5.72 1.06e-08
#> 6 High synthesis adapt 1 40 3.13 0.14388 0.2797 21.21 7.61e-100
#> yules_q kappa kappa_z kappa_p lift sign significant
#> 1 -0.544 -0.6606 -4.98 6.36e-07 0.367 under TRUE
#> 2 0.371 0.0636 2.90 3.68e-03 1.990 over TRUE
#> 3 0.466 0.1803 5.21 1.86e-07 2.132 over TRUE
#> 4 -0.589 -0.7375 -3.46 5.48e-04 0.287 under TRUE
#> 5 -0.824 -0.8859 -5.99 2.04e-09 0.123 under TRUE
#> 6 0.905 0.2406 19.62 1.08e-85 12.800 over TRUE
nodes(gfit) |> head(6)
#> group state outgoing incoming
#> 1 High adapt 141 143
#> 2 High cohesion 938 936
#> 3 High consensus 3392 3439
#> 4 High coregulate 894 954
#> 5 High discuss 2003 1986
#> 6 High emotion 1550 1517
Compare groups
compare_lsa() tests whether groups differ by permuting
group labels at the sequence level. The result includes edge-level
differences and an omnibus comparison.
cmp <- compare_lsa(gfit, R = 100, adjust = "BH")
cmp
#> <lsa_comparison>
#> groups: High vs Low
#> measure: log_or difference (High - Low)
#> R: 100 label permutations
#> edges: 35 significant of 78 tested (adjust = BH)
#> omnibus: statistic = 79.48, p = 0.009901
transitions(cmp) |> head(6)
#> from to log_or_a log_or_b diff p_perm p_adj significant
#> 1 adapt adapt -1.183 -3.17 1.984 NA NA FALSE
#> 2 cohesion adapt -0.794 -3.92 3.127 0.0099 0.0249 TRUE
#> 3 consensus adapt -1.219 -1.97 0.748 0.0594 0.1030 FALSE
#> 4 coregulate adapt 0.778 -1.08 1.855 0.0099 0.0249 TRUE
#> 5 discuss adapt 1.010 2.25 -1.241 0.0099 0.0249 TRUE
#> 6 emotion adapt -1.352 -3.10 1.751 0.0891 0.1390 FALSE
The default comparison plot is a back-to-back barrel display. It
shows direction, magnitude, and the permutation p-value for each
edge.
In short
The complete analysis is a claim-to-evidence chain:
fit <- lsa(engagement)
transitions(fit, significant = TRUE)
plot(fit)
plot(fit, type = "network", weights = "tna")
certainty_lsa(fit)
bootstrap_lsa(fit)
reliability_lsa(fit)
log_fit <- lsa(log_data, actor = "Actor", action = "Action", time = "Time")
group_fit <- lsa(log_data, actor = "Actor", action = "Action",
time = "Time", group = "Group")
compare_lsa(group_fit)
Descriptive claims use the fitted transition tables and plots. Edge
claims add uncertainty. Whole-network claims add reliability or
stability. Group-difference claims use permutation comparison.
These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.