Pith. sign in

REVIEW 3 major objections 4 minor 2 cited by

HOMC: A MATLAB Package for Higher Order Markov Chains

T0 review · 3 major / 4 minor · reviewed 2026-08-04 · deepseek-v4-flash

Challenge this review Re-run · record.json Download PDF Read on arXiv ↗

Pith's one-line read The paper introduces HOMC, the first dedicated MATLAB package for higher-order Markov chains, and shows how to compute their key tensors and classifications.

desk verdict A legitimate software paper with a real artifact, but the erp stopping rule needs an error bound before the state-classification output can be trusted. read the letter →

arxiv 2510.02664 v1 pith:ENSXWRCP submitted 2025-10-03 stat.CO math.PR

Jianhong Xu This is my paper
classification stat.COmath.PR MSC 60J1015A69
keywords higherorderMarkovchainsMATLABtensorstransitionprobabilitiesever-reachinglimitingprobabilitydistributionsmeanfirstpassagetimesclassificationofstates
verification ladder T0 review T1 audit T2 compute T3 formal

The pith

A machine-rendered reading of the paper's core claim, the machinery that carries it, and where it could break.

The reading

The paper introduces HOMC, a MATLAB package for higher-order Markov chains—chains whose next state depends on several past states, not just the present one. Its claim is that this is the first dedicated software of its kind, and that it can compute the main numerical quantities from the recent theory of such chains: k-step transition tensors, limiting probability distributions, ever-reaching probability tensors, and mean first passage time tensors, plus regularity and ergodicity checks and state classification. If that is right, an entire toolkit of tensor-based Markov chain results becomes directly usable for numerical experiments and applications.

What carries the argument

The box product ⊠ of two m-th order n-dimensional tensors A and B is the tensor C with c_{i1...im} = Σ_j a_{i1 j i2 ... i_{m-1}} b_{j i2 ... im}. It reduces to ordinary matrix multiplication when m=2, but is non-associative for m≥3, which is why the package provides dedicated powering (bpow) rather than relying on matrix-style identities. This product is the algebraic engine: tensor powers give k-step transition probabilities, and the recurrences for first-passage and mean-first-passage tensors are box-product equations. Memory efficiency comes from generating index tuples on the fly with MATLAB's ind2sub rather than storing all index combinations.

What would settle it

Run erp on the two-state first-order chain with transition matrix P = [[1−α, 0.5], [α, 0.5]] for α = 10^−10 and default tolerance 10^−6. The routine stops when max(|F[k]|) < 10^−6, which for this chain happens around k = 21; the returned F_{12} is about 2×10^−9, while the true ever-reaching probability from state 2 to state 1 is 1 (the sum of the long-tailed geometric series). This divergence would show the truncation criterion is not a valid error bound.

Watch

Extended reading notes

Core claim

The central claim is that a single MATLAB package, built around a non-associative tensor product called the box product (⊠), can replace bespoke code for almost all numerical work with higher-order Markov chains. For an m-th order, n-dimensional transition tensor P, the box product defines tensor powers whose entries are exactly k-step transition probabilities; the same operation yields the ever-reaching probability tensor as a series of box products and the mean first passage time tensor as the solution of a tensor linear equation. The package also constructs the reduced first-order chain's transition matrix via a Khatri-Rao product, so limiting distributions can be obtained from a dominant

Load-bearing premise

The ever-reaching probability function stops summing the infinite series as soon as the largest term falls below a tolerance, without a proof that the omitted tail is that small; for chains with slowly decaying first-passage probabilities this can silently misreport the ever-reaching probability by far more than the tolerance.

Editorial extensions

If this is right

  • Users can compute the limiting distribution of a higher-order chain by forming the reduced first-order transition matrix and taking a dominant eigenvector, without manual bookkeeping.
  • The ever-reaching probability tensor gives a direct state classification: a state is recurrent iff its diagonal entries in F are all 1, transient otherwise.
  • The package makes the recent theoretical results on higher-order chains available as off-the-shelf routines, lowering the barrier to numerical experimentation.
  • Because it also handles first-order chains (m=2) as a special case, it can serve as a single tool for a broad range of Markov chain models.
  • The O(n^m) memory footprint is the theoretical minimum for dense tensors, so the package scales as far as the underlying data structure can reasonably go.

Reading between the lines

Editorial extensions of the paper, not claims the author makes directly.

  • The truncation rule in the ever-reaching probability routine is a practical heuristic, not a certified numerical method; for chains with heavy-tailed first-passage times, the returned tensor can be wrong by far more than the tolerance. A tail bound or an alternative linear-system solve would make this routine reliable.
  • The non-associativity of the box product is a genuine conceptual difference from matrix Markov chains: standard matrix-power shortcuts do not transfer, so users must think in terms of the package's bpow rather than ordinary algebra.
  • The dense O(n^m) storage limits practical state spaces; extending the package to exploit sparse or structured transition tensors would open up higher-order chain models with much larger state spaces.
  • The routines assume the theory in the related papers (regularity/ergodicity for uniqueness); a useful diagnostic addition would be warnings when a chain fails the ergodicity check before computing mean first passage times.
Share X Bluesky LinkedIn Reddit HN

Editorial analysis

A structured set of objections, weighed in public.

Desk editor's note, referee report, and a circularity audit.

Referee Report

3 major / 4 minor

Summary. The manuscript presents HOMC, a MATLAB package for higher-order Markov chains. Its functions implement a tensor 'box' product and powers (bprod, bpow), linear indexing (lind), diagonal tensors (diagt), matricization/tensorization (t2mat, mat2t), construction of the reduced first-order chain (rcmat, krprod), ever-reaching probabilities (erp), and mean first passage times (mfptd, mfpti). The paper claims these tools can compute the k-step transition tensor, limiting distributions, ever-reaching probability tensors, and MFPT tensors, and can check ergodicity/regularity and classify states as recurrent or transient. The manuscript includes worked examples, including a residual check of the MFPT equation and first-order reductions to standard matrix results.

Significance. If the package performs as claimed, it would fill a genuine gap: no other dedicated MATLAB toolkit for higher-order Markov chains is known to the community, and the tensor-based formulation is natural for multidimensional transition data. The paper's internal checks are a real strength: Example 3.2 verifies the MFPT equation to zero residual, and Examples 2.4 and 3.3 reproduce standard first-order matrix results. The package builds on the authors' prior theoretical work [8,9,11], which is appropriate for a software presentation. The main weakness is the erp truncation rule, which has no tail error bound and directly affects the state-classification feature; additionally, the 'first of its kind' claim is not benchmarked against existing software, and the manuscript does not include a test suite or code listing.

major comments (3)
  1. [§3.1, Eq. (5) and erp stopping rule] The stopping test max(abs(F[k])) < tol does not control the omitted tail sum_{l>k} F[l]. Since all entries are nonnegative, a small per-step increment does not imply a small remaining tail; for a nearly decomposable chain the tail can be as large as tol/(1−ρ), where ρ is the spectral radius of the substochastic matrix for transitions avoiding the target state. With ρ close to 1, the tail can exceed tol by orders of magnitude. This is load-bearing because the same section uses the exact condition f_{ii3...im}=1 to classify recurrent states. Please provide a tail bound, or compute F by solving the linear tensor equation F = P + (F−F_d)⊠P (analogous to Eq. (7)), or at least clearly document the limitation and require the user to supply a much stricter tol after an explicit convergence check.
  2. [§3.1, state classification] The recurrent/transient classifier hinges on exact equality f_{ii3...im}=1. With the approximate erp output, a recurrent state may be returned as 0.999997 and a transient state with true probability 0.99999 may be rounded to 1. The default tol=1e-6 gives no assurance in either direction. This is not a purely theoretical concern: the classification routine is a headline feature of the package. The manuscript should either return the approximate probability together with a confidence statement, or use an exact method for the binary recurrent/transient decision.
  3. [General, §2–§3, reproducibility] The paper's verification consists of four examples, but no automated test suite, code listing, or systematic function-reference mapping is included. For a software paper, this is insufficient to establish that the package behaves correctly beyond the displayed inputs. I recommend adding a small regression test suite and a table listing each function, its signature, and the examples that exercise it. This would also make the claim that HOMC computes 'all important quantities' in the cited works independently checkable.
minor comments (4)
  1. [Abstract and §1] 'First of its kind' is a strong claim. Either report a brief search of existing MATLAB toolboxes for Markov chains/HMMs/tensor methods, or soften the wording to 'to our knowledge, no dedicated package...'.
  2. [§2.2, identity tensor] When defining I, specify explicitly that it is a left identity for the box product but not a right identity. This prevents confusion with standard tensor identity tensors and explains the non-associativity warning.
  3. [Example 3.1] The text states that erp used the first 67 terms, but the exact MATLAB call and formatting are not shown. Include the full command sequence and the MATLAB version/rounding behavior so that the displayed outputs can be reproduced bit-for-bit.
  4. [§3.2, mfpti] The iterative stopping criterion is standard, but no error estimate is given for the fixed-point iteration. Please state whether the convergence proof in [9] provides a contraction rate that would justify interpreting tol as an a posteriori bound on the error.

Circularity Check

0 steps flagged · score 2.0 of 10

No structural circularity; package implements displayed recurrences and matrix-level checks, with only a minor self-citation for the limiting-distribution theorem.

full rationale

This is a software paper rather than a derivation chain. The package implements formulas that are either stated explicitly in the paper (the box-product definition, the P^k recurrence, the ever-reaching sum (5), and the mean-first-passage linear system (6)) or taken from the author's own prior mathematics. No parameter is fitted to data and then renamed as a prediction, and no output quantity is defined in terms of the algorithm's stopping rule. The limiting distribution is produced by constructing the reduced-chain matrix Q and calling MATLAB's eig; the ever-reaching tensor is an actual partial sum of the displayed series; the mean first passage time routines solve or iterate the displayed linear system. The one point touching the self-citation boundary is Section 2.5, where the uniqueness of the limiting distribution is imported from [11], an under-review preprint by the same author; however, this is an externally checkable mathematical theorem about regular higher-order chains, not an ansatz or a fitted value, and the rest of the package is independently checkable against first-order matrix special cases. The erp truncation criterion in Section 3.1 (stop when max(abs(F^{[k]})) < tol) lacks a stated tail bound, which is an omitted numerical proof and a robustness concern, but it is not circular reasoning. The minor self-citation to [11] does not make the central claim circular, hence the low score.

Assumptions & free parameters 2 free parameters · 4 assumptions · 0 invented entities

The package adds no free physical or mathematical constants; the numerics rest on prior theoretical results (mostly same-author) about box-product powers, reduced first-order chains, and mean first passage systems. The tolerance parameters govern approximation quality and are not fitted to data. No new entities are introduced.

free parameters (2)
  • erp tolerance = 1e-6 (default)
    Stopping criterion for truncating the ever-reaching probability series (5); chosen by hand, no tail-error bound is derived in the paper.
  • mfpti tolerance = 1e-6 (default)
    Stopping criterion for the iterative mean first passage time solver in Section 3.2; affects accuracy of the reported approximation.
assumptions (4)
  • domain assumption The box product I⊠A=A and the power recursion A^{k+1}=A^k⊠A provide the correct k-step transition probabilities for higher-order Markov chains.
    Invoked in Section 2.2 to define P^k and to infer regularity from P^10>0; these results are cited from the author's earlier papers [9,11], not proved in this preprint.
  • domain assumption For an ergodic higher-order chain, the mean first passage time tensor is the unique solution of μ=E+(μ−μ_d)⊠P and the fixed-point iteration converges.
    Used in Section 3.2 for both mfptd and mfpti; cited to [9] and not rederived.
  • domain assumption The reduced first-order chain with transition matrix Q=G*P (columnwise Khatri-Rao product) correctly governs the higher-order chain's limiting distribution.
    Used in Section 2.5 and Example 2.6 to obtain π; cited to [11], an arXiv preprint under review.
  • domain assumption If some power P^k is entrywise positive, the chain is regular and has a unique limiting probability distribution.
    Used in Example 2.2; cited to [11].

how reviews work

0 comments
Cite this review

Pith. "Pith review of HOMC: A MATLAB Package for Higher Order Markov Chains." pith.science (2026). https://pith.science/paper/ENSXWRCP

@misc{pith2026251002664,
  author       = {Pith},
  title        = {Pith review of: HOMC: A MATLAB Package for Higher Order Markov Chains},
  year         = {2026},
  howpublished = {\url{https://pith.science/paper/ENSXWRCP}},
  note         = {Machine review of arXiv:2510.02664}
}
abstract

We present a MATLAB package, which is the first of its kind, for Higher Order Markov Chains (HOMC). It can be used to easily compute all important quantities in our recent works relevant to higher order Markov chains, such as the $k$-step transition tensor, limiting probability distribution, ever-reaching probability tensor, and mean first passage time tensor. It can also be used to check whether a higher order chain is ergodic or regular, to construct the transition matrix of the associated reduced first order chain, and to determine whether a state is recurrent or transient. A key function in the package is an implementation of the tensor ``box'' product which has a probabilistic interpretation and is different from other tensor products in the literature. This HOMC package is useful to researchers and practitioners alike for tasks such as numerical experimentation and algorithm prototyping involving higher order Markov chains.

Discussion (0). Sign in to comment.

Forward citations

Cited by 2 Pith papers

Reviewed papers in the Pith corpus that reference this work. Sorted by Pith novelty score. Full citation record

  1. An Alternative Framework for Irreducibility and Primitivity of Nonnegative Tensors

    math.RA 2026-06 unverdicted novelty 6.0 of 10

    Introduces s-irreducibility and s-primitivity for nonnegative tensors as a framework that includes matrix results as special cases and expands prior tensor characterizations, with implications for higher-order Markov chains.

  2. Can a Higher Order Markov Chain Be Treated as a First Order Markov Chain?

    math.PR 2025-12 conditional novelty 4.0 of 10

    A higher-order Markov chain cannot in general be treated through its reduced first-order chain: ergodicity, ever-reaching probabilities, mean first passage times, and state classification can all fail to transfer, tho...

Reference graph

Works this paper leans on

28 extracted references · 1 linked inside Pith · cited by 2 Pith papers

  1. [1]

    Baena-Mirabete, P

    S. Baena-Mirabete, P. Puig, Parsimonious higher order Markov models for rating transitions, Journal of the Royal Statistical Society Series A: Statistics in Society 181 (2018): 107–131

    2018

  2. [2]

    Chang, T

    K. Chang, T. Zhang, On the uniqueness and non-uniqueness of the pos- itiveZ-eigenvector for transition probability tensors, Journal of Mathe- matical Analysis and Applications 408 (2013): 525–540

    2013

  3. [3]

    Deneshkumar, M

    V. Deneshkumar, M. Manoprabha, K. Annan, Markov chain modeling for immunological status of HIV patients, Communications in Statistics: Case Studies, Data Analysis and Applications 5 (2019): 85–91

    2019

  4. [4]

    Doob, Stochastic Processes, Wiley Classic Library Edition, Wiley, New York, 1990

    J. Doob, Stochastic Processes, Wiley Classic Library Edition, Wiley, New York, 1990

    1990

  5. [5]

    Flett, N

    G. Flett, N. Kelly, An occupant-differentiated, higher-order Markov chain method for prediction of domestic occupancy, Energy and Buildings 125 (2016): 219–230

    2016

  6. [6]

    Geiger, A sufficient condition for a unique invariant distribution of a higher-order Markov chain, Statistics and Probability Letters 130 (2017): 49–56

    B. Geiger, A sufficient condition for a unique invariant distribution of a higher-order Markov chain, Statistics and Probability Letters 130 (2017): 49–56

    2017

  7. [7]

    Gleich, L

    D. Gleich, L. Lim, Y. Yu, Multilinear pagerank, SIAM Journal on Matrix Analysis and Applications 36 (2015): 1507–1541

    2015

  8. [8]

    L. Han, K. Wang, J. Xu, Higher order ergodic Markov chains and first passage times, Linear and Multilinear Algebra 70 (2022): 6772–6779

    2022

Show all 28 references
  1. [9]

    L. Han, J. Xu, Ever-reaching probabilities and mean first passage times of higher order ergodic Markov chains, Linear and Multilinear Algebra 72 (2024): 59–75. 18

    2024

  2. [10]

    L. Han, J. Xu, On classification of states in higher order Markov chains, Linear Algebra and Its Applications 685 (2024): 24–45

    2024

  3. [11]

    L. Han, J. Xu, On limiting probability distributions of higher order Markov chains, under review,arXiv.2506.08874v2

    arXiv

  4. [12]

    S. Hu, L. Qi, Convergence of a second order Markov chain, Applied Mathematics and Computation 241 (2014): 183–192

    2014

  5. [13]

    Islam, R

    M. Islam, R. Chowdhury, A higher order Markov model for analyzing covariate dependence, Applied Mathematical Modelling 30 (2006): 477– 488

    2006

  6. [14]

    Kemeny, J

    J. Kemeny, J. Snell, Finite Markov Chains, Springer-Verlag, New York, 1960

    1960

  7. [15]

    Kolda, B

    T. Kolda, B. Bader, Tensor decompositions and applications, SIAM Review 51 (2009): 455–500

    2009

  8. [16]

    D. Ky, L. Tuyen, A higher order Markov model for time series fore- casting, International Journal of Applied Mathematics and Statistics 57 (2018): 1–18

    2018

  9. [17]

    J. Lan, X. Li, V. Jilkov, C. Mu, Second-order Markov chain based multiple-model algorithm for maneuvering target tracking, IEEE Trans- actions on Aerospace and Electronic Systems 49 (2013): 3–19

    2013

  10. [18]

    C. Li, S. Zhang, Stationary probability vectors of higher-order Markov chains, Linear Algebra and Its Applications 473 (2016): 114–125

    2016

  11. [19]

    W. Li, M. Ng, On the limiting probability distribution of a transition probability tensor, Linear and Multilinear Algebra 62 (2014): 362–385

    2014

  12. [20]

    Z. Liu, Y. Luo, Y. Zhu, State estimation for linear dynamic system with multiple-step random delays using high-order Markov chain, IEEE Access 8 (2020): 76218–76227

    2020

  13. [21]

    Martin, R

    C. Martin, R. Shafer, B. Larue, An order-ptensor factorization with applications in imaging, SIAM Journal on Scientific Computing 35 (2013): 474–490. 19

    2013

  14. [22]

    Masseran, Markov chain model for the stochastic behaviors of wind- direction data, Energy Conversion and Management 92 (2015): 266–274

    N. Masseran, Markov chain model for the stochastic behaviors of wind- direction data, Energy Conversion and Management 92 (2015): 266–274

    2015

  15. [23]

    L. Qi, Z. Luo, Tensor Analysis: Spectral Theory and Special Tensors, SIAM, Philadelphia, 2017

    2017

  16. [24]

    Sanjari, H

    M. Sanjari, H. Gooi, Probabilistic forecast of PV power generation based on higher order Markov chain, IEEE Transactions on Power Systems 32 (2017): 2942–2952

    2017

  17. [25]

    Smilde, R

    A. Smilde, R. Bro, P. Geladi, Multi-Way Analysis: Applications in the Chemical Sciences, Wiley, West Sussex, UK, 2004

    2004

  18. [26]

    S. Wu, M. Chu, Markov chains with memory, tensor formulation, and the dynamics of power iteration, Applied Mathematics and Computation 303 (2017): 226–239

    2017

  19. [27]

    Xiong, R

    H. Xiong, R. Mamon, A higher-order Markov chain-modulated model for electricity spot-price dynamics, Applied Energy 233 (2019): 495–515

    2019

  20. [28]

    Y. Yang, H. Jang, B. Kim, A hybrid recommender system for sequential recommendation: combining similarity models with Markov chains, IEEE Access 8 (2020): 190136–190146. 20

    2020

Pith tools

Reviewed August 4, 2026 · model on record in the stance chip above.