EM-SMC-RKHS Procedure

• The EM-SMC-RKHS procedure is a nonparametric inference framework that models continuous-time dynamics by encoding unknown functions in an RKHS and applying EM for latent variable optimization.
• It adapts to different settings by employing kernel regression with spike-and-slab selection for CTMCs and particle-based SMC filtering for sparse, noisy SDE observations.
• The framework ensures convergence and practical performance through finite-dimensional approximations, demonstrated by metrics like MSE and Kolmogorov distance in applications from biomedical models to complex dynamical systems.

The EM-SMC-RKHS procedure denotes a family of inference schemes that couple reproducing kernel Hilbert space parameterizations of unknown continuous-time dynamics with expectation-maximization updates and, when latent trajectories make the likelihood intractable, Sequential Monte Carlo approximations of the E-step. In the literature summarized here, the construction appears in two closely related forms. For covariate-based continuous-time Markov chains, the primary development is an RKHS model with Frequentist penalization and a Bayesian EMVS posterior-mode algorithm, while an EM-SMC-RKHS variant is presented as a conceptual extension for latent-path settings (Han et al., 6 May 2025). For sparse and noisy stochastic differential equations, the EM-SMC-RKHS procedure is developed explicitly as a nonparametric drift-learning algorithm whose M-step is reduced to finite-dimensional kernel regression by a generalized representer theorem (Ganguly et al., 15 Aug 2025). Its Bayesian interpretation is underwritten by the RKHS–Gaussian random field correspondence, under which regularized RKHS estimators coincide with posterior means in the quadratic-loss case and with finite-dimensional MAP estimators under general losses (Aravkin et al., 2013).

1. Conceptual architecture

At its core, the procedure has three layers. The RKHS layer represents unknown transition mechanisms or drift fields as functions in a reproducing kernel Hilbert space; the EM layer treats either latent trajectories or latent inclusion variables as missing data; and the SMC layer supplies particle-based approximations when the relevant filtering distribution is analytically intractable. In the continuous-time Markov chain formulation, the unknown object is the covariate-dependent generator $Q(z)$, with off-diagonal intensities represented through unconstrained log-intensity functions in an RKHS. In the stochastic differential equation formulation, the unknown object is the entire drift function $b:\mathbb{R}^d\to\mathbb{R}^d$, again learned in a vector-valued RKHS (Han et al., 6 May 2025).

Component	Role	Explicit status
RKHS	Finite kernel expansion via generalized representer theorem	Explicit in CTMC and SDE
EM	Iterative optimization with latent variables	Explicit in EMVS and SDE EM
SMC	Approximate filtering distribution and E-step expectations	Explicit in SDE; conceptual in CTMC

A central point of terminology is that “EM-SMC-RKHS” is not used in a single uniform sense across these works. In the CTMC setting, Sequential Monte Carlo is not part of the estimation procedure actually implemented; the paper states that it “does not explicitly employ Sequential Monte Carlo (SMC), but its CTMC–RKHS–EMVS framework naturally suggests an extension for latent-path or more complex observation models” (Han et al., 6 May 2025). By contrast, in the SDE setting the procedure is an implemented algorithm: sparse observations render the marginal likelihood intractable, EM decomposes the optimization, SMC approximates the filtering law, and the M-step becomes penalized empirical risk minimization in RKHS (Ganguly et al., 15 Aug 2025).

This distinction also clarifies a common misconception. The phrase does not identify a single canonical algorithmic template with fixed update equations; rather, it identifies a design pattern whose precise latent variables, filtering targets, and regularization mechanisms depend on the underlying continuous-time model.

2. RKHS representation and Bayesian meaning

The RKHS component provides the nonparametric function space in which the unknown dynamics are estimated. In the CTMC model, off-diagonal rates are parameterized by

$q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$

with diagonal entries determined by

$q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$

The vector of off-diagonal log-intensity functions is assumed to belong to a vector-valued RKHS $\mathcal{H}_{\mathbf K}$, and a generalized Representer Theorem implies that the minimizer lies in the finite span of kernel sections at the observed covariates: $$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$ With a block-diagonal kernel, each transition arm has its own scalar kernel and coefficient sequence, so that

$$\tilde q(i,j;z)=\sum_{l=1}^L K_{(i,j)}(z,z^{(l)})\,\tilde\beta_l^{(i,j)}.$$

This turns the infinite-dimensional CTMC likelihood optimization into a finite-dimensional problem over kernel coefficients (Han et al., 6 May 2025).

In the SDE formulation, the unknown drift $b$ belongs to a vector-valued RKHS $\mathcal H_{\mathsf K}$, and the estimation target is the penalized negative log-likelihood functional

$\loss(b)=-\ell(b\mid \by_{1:M_0})+\lambda\|b\|_{\mathcal H_{\mathsf K}}^2.$

After the E-step has been approximated, the generalized representer theorem shows that the M-step minimizer admits a finite expansion supported at particle locations: $b:\mathbb{R}^d\to\mathbb{R}^d$0 The support of the kernel expansion is therefore data-adaptive at the level of latent-state reconstruction rather than fixed solely by observed covariates (Ganguly et al., 15 Aug 2025).

The Bayesian interpretation of this regularization is given by the Gaussian random field analysis in “The connection between Bayesian estimation of a Gaussian random field and RKHS” (Aravkin et al., 2013). For the estimator

$b:\mathbb{R}^d\to\mathbb{R}^d$1

if $b:\mathbb{R}^d\to\mathbb{R}^d$2 is assigned a zero-mean Gaussian random field prior with covariance $b:\mathbb{R}^d\to\mathbb{R}^d$3 and the likelihood is

$b:\mathbb{R}^d\to\mathbb{R}^d$4

then $b:\mathbb{R}^d\to\mathbb{R}^d$5. Under quadratic loss, the RKHS estimator is the posterior mean. Under general loss, the posterior over the infinite-dimensional function is not a proper density, but for any finite set of sampling locations the MAP estimate of the finite vector of samples coincides with the RKHS estimator evaluated at those locations. This corrects another common misunderstanding: under non-Gaussian losses the RKHS solution is not the MAP of an infinite-dimensional field, but it is the MAP for every finite collection of function values (Aravkin et al., 2013).

3. EM formulations: latent paths and latent inclusions

The EM layer appears in two distinct forms.

In the SDE problem, the hidden variable is the full latent trajectory on a fine discretization grid,

$b:\mathbb{R}^d\to\mathbb{R}^d$6

The joint model couples the latent path density with a sparse noisy observation process, and the marginal likelihood requires integrating over all unobserved path segments. The exact EM map is defined through

$b:\mathbb{R}^d\to\mathbb{R}^d$7

with $b:\mathbb{R}^d\to\mathbb{R}^d$8 taken as the filtering distribution of the latent path under the current drift. The exact EM sequence is

$b:\mathbb{R}^d\to\mathbb{R}^d$9

Here the E-step computes the filtering distribution $q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$0, and the M-step minimizes the expected complete-data negative log-likelihood plus the RKHS penalty (Ganguly et al., 15 Aug 2025).

In the CTMC Bayesian formulation, EM is used differently. The data consist of subject-specific covariates together with jump times and states, and the model introduces spike-and-slab latent indicators $q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$1 for the kernel coefficients $q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$2. The joint posterior is non-conjugate because the CTMC likelihood is non-Gaussian in the coefficients, direct MCMC is difficult in high dimension, and the paper therefore adopts an EMVS-style posterior-mode search. The latent variables in the E-step are not trajectories but inclusion indicators. Given $q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$3, the E-step computes posterior inclusion probabilities $q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$4; the M-step then maximizes an expected log-posterior in which the spike-and-slab hierarchy induces a weighted quadratic penalty on the coefficients (Han et al., 6 May 2025).

These two EM constructions are formally related but statistically different. One integrates over latent paths; the other integrates over latent sparsity indicators. This suggests that “EM” in the EM-SMC-RKHS label should be read structurally rather than narrowly: it denotes iterative surrogate optimization driven by conditional expectations, not a single fixed missing-data formulation.

4. Sequential Monte Carlo as the E-step engine

Sequential Monte Carlo becomes necessary when the E-step requires filtering distributions that are analytically intractable. In the sparse-observation SDE setting, the filtering law

$q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$5

is not Gaussian, and the marginal likelihood has no analytic expression. The SMC approximation therefore represents the filtering distribution by a weighted empirical measure

$q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$6

Particle proposals are factorized over observation intervals, and the importance weights are updated by the ratio of the observation density and true latent transition densities to the proposal density. To avoid degeneracy, effective sample size

$q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$7

is monitored and resampling is performed when $q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$8 (Ganguly et al., 15 Aug 2025).

A distinctive feature of the SDE construction is the proposal design. Rather than using a drift-free modified diffusion bridge, the proposal uses a first-order linear SDE approximation in which the drift is expanded through its Jacobian. This yields Gaussian bridge proposals with mean and covariance determined by an ODE system for $q(i,j;z)=\exp(\tilde q(i,j;z)),\qquad i\neq j,$9 and $q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$0, and finally a Gaussian proposal

$q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$1

The proposal therefore depends on the current drift estimate and its Jacobian, which is important when the drift is highly nonlinear (Ganguly et al., 15 Aug 2025).

The CTMC paper does not implement SMC, but it specifies how SMC would enter in a latent-path extension. Given current $q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$2, an SMC algorithm would approximate $q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$3 by weighted particles, and the E-step would compute expected sufficient statistics such as the expected number of transitions $q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$4 and the expected total intensity exposure $q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$5. These quantities would replace the complete-data terms in the CTMC likelihood before the RKHS/EMVS M-step (Han et al., 6 May 2025).

Because this CTMC use of SMC is explicitly presented as conceptual, it should not be conflated with an implemented algorithm. A plausible implication is that the EM-SMC-RKHS label is most literal in the SDE setting and more schematic in the CTMC setting.

5. M-step structure, sparsity, and shrinkage

Once the E-step produces either expected latent-path statistics or posterior inclusion probabilities, the M-step becomes an optimization over RKHS coefficients.

For the CTMC model, the Frequentist estimator minimizes the negative log-likelihood plus a normed square penalty in RKHS. The paper states that, unlike the SDE case, no closed-form solution exists because of the CTMC likelihood structure, and quasi-Newton gradient descent algorithms such as Nelder–Mead or BFGS via R’s optim() are used. In the Bayesian spike-and-slab version, the E-step computes

$q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$6

and the M-step maximizes a CTMC log-likelihood plus a weighted $q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$7 penalty whose weights are

$q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$8

Hyperparameter selection for $q(i,i;z)=-\sum_{j\neq i}q(i,j;z).$9 is based not on grid search but on clustering Frequentist estimates $\mathcal{H}_{\mathbf K}$0 with mixtools in R; in simulations the paper fixes $\mathcal{H}_{\mathbf K}$1 (Han et al., 6 May 2025).

For the SDE model, the approximate M-step after SMC has a weighted least-squares form. Writing particle increments as

$\mathcal{H}_{\mathbf K}$2

the objective becomes a quadratic form in $\mathcal{H}_{\mathbf K}$3 plus the RKHS penalty. The representer theorem reduces this to a finite linear system for the stacked coefficient vector $\mathcal{H}_{\mathbf K}$4: $\mathcal{H}_{\mathbf K}$5 with

$\mathcal{H}_{\mathbf K}$6

The M-step is therefore explicit once the particle system has been fixed (Ganguly et al., 15 Aug 2025).

The SDE paper also proposes a hybrid Bayesian EM variant in which kernel coefficients are given Gaussian priors with inverse-gamma hyperpriors on local variances,

$\mathcal{H}_{\mathbf K}$7

Marginally this induces a multivariate Student-$\mathcal{H}_{\mathbf K}$8 prior, replacing the global RKHS penalty with coefficient-specific shrinkage. The resulting updates retain the same linear-system structure but with a diagonal local-precision term in place of $\mathcal{H}_{\mathbf K}$9 (Ganguly et al., 15 Aug 2025).

The comparison between these variants is instructive. The CTMC method emphasizes spike-and-slab variable selection and sparse posterior modes; the SDE method emphasizes continuous shrinkage through inverse-gamma mixing. This suggests a methodological spectrum from discrete inclusion control to heavy-tailed local shrinkage, both implemented within the broader RKHS-EM framework.

6. Convergence, empirical assessment, and methodological boundaries

The SDE paper provides explicit convergence results. Under Assumption 2.1, the M-step functional is strictly convex, lower semicontinuous, and has a unique minimizer. For the exact EM sequence, the penalized loss is non-increasing and convergent, the sequence is strongly precompact in $$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$0, and any limit point belongs to the stationary set. For approximate EM, if the filtering approximations converge in KL divergence to the true filtering distribution, one can choose the particle numbers $$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$1 large enough so that the approximate EM sequence retains loss convergence and convergence to the stationary set (Ganguly et al., 15 Aug 2025).

The paper also states an important qualification: a raw SMC approximation does not give KL convergence of the density, although kernel density estimates could in principle satisfy the required condition. This is a methodological boundary rather than a contradiction. The practical algorithm is SMC-based and behaves well empirically, but the strongest abstract convergence statement is formulated for filtering approximations satisfying a KL condition (Ganguly et al., 15 Aug 2025).

Empirical evaluation is reported in both model classes. In the CTMC study, effectiveness is assessed through “the normalized difference between estimated and true nonlinear transition functions” and through “the difference in the probability of getting absorbed in one the final states,” intended to capture long-term behavior. The paper also presents mean squared error of log-intensity functions and reports absorption-probability differences that are “small (on the order of $$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$2 or less) for sufficiently large datasets.” In a follicular cell lymphoma application, the asymptotic distribution of ending states under the RKHS-based method is compared with a Cox PH-based multistate model and the actual observed state; the RKHS method matches the observed trajectory more closely in some cases (Han et al., 6 May 2025).

In the SDE study, reported performance metrics are mean squared error between $$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$3 and $$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$4, and the Kolmogorov distance between stationary distributions of the true and estimated SDEs. The numerical examples include a double-well potential SDE, a multiplicative-noise variant, a Gamma SDE, Michaelis–Menten kinetics, and a near-deterministic SIR model. Reported results include MSE around $$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$5 for Model 3 with $$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$6 of the data, Kolmogorov distances around $$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$7–$$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$8 in one-dimensional examples, MSE around $$\boldsymbol{\tilde q}(z)=\sum_{l=1}^L \mathbf K(z,z^{(l)})\tilde\beta_l.$$9 for Michaelis–Menten kinetics, and MSE $$\tilde q(i,j;z)=\sum_{l=1}^L K_{(i,j)}(z,z^{(l)})\,\tilde\beta_l^{(i,j)}.$$0 for the SIR system (Ganguly et al., 15 Aug 2025).

Taken together, these results establish the procedure as a nonparametric continuous-time modeling framework with two distinct operational modes. In one mode, fully observed or nearly fully observed CTMC trajectories permit RKHS likelihood optimization with EMVS sparsity, and SMC enters only as a proposed extension for partially observed settings. In the other, sparse and noisy SDE observations require a genuine EM-SMC-RKHS algorithm in which particle filtering is indispensable. The broader Bayesian-RKHS theory indicates why these constructions remain finite-dimensional at the optimization stage: regularized function estimation in RKHS is equivalent, on finite sets of locations, to MAP or posterior-mean estimation under Gaussian random field priors (Aravkin et al., 2013).