Jeans Anisotropic Modeling (JAM): Overview

• Jeans Anisotropic Modeling (JAM) is a dynamical modeling approach that uses the Jeans equations and Multi-Gaussian Expansion to derive mass profiles and anisotropy from kinematic maps.
• It employs symmetry assumptions, deprojection techniques, and anisotropy prescriptions to predict projected second velocity moments for direct comparison with observations.
• JAM is applied to diverse systems from galaxies to clusters, though its accuracy relies on equilibrium and axisymmetry, which can be challenged by triaxiality or dynamic feedback.

Jeans Anisotropic Modeling (JAM) is a family of dynamical modeling methods that predicts projected second velocity moments from parametrized tracer and mass distributions by solving the steady-state Jeans equations under symmetry and anisotropy assumptions. In its standard axisymmetric form, JAM combines a deprojected light model—usually expressed as a Multi-Gaussian Expansion (MGE)—with a total gravitational potential and an anisotropy prescription, then fits the model prediction for $V_{\rm rms}\equiv\sqrt{V^2+\sigma^2}$ to observed kinematic maps in order to infer the mass distribution, stellar mass-to-light ratio, dark halo structure, inclination, and anisotropy (Cappellari, 2012, Li et al., 2017). Subsequent developments extended the formalism to all six projected second moments, to spherical and spherically aligned axisymmetric variants, and to applications involving multiple stellar populations, dwarf galaxies, galaxy clusters, strong-lensing galaxies, and the Milky Way halo (Cappellari, 2015, Cappellari, 2019).

1. Formal definition and dynamical assumptions

The canonical JAM framework assumes a stationary, axisymmetric, collisionless system in cylindrical coordinates $(R,\phi,z)$, with the velocity ellipsoid aligned with the cylindrical axes and the mixed term $\sigma_{Rz}$ neglected. Under these assumptions, the Jeans equations reduce to

$$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$

with a meridional anisotropy parameter

$\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$

The observable used in most implementations is the line-of-sight second moment,

$V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$

which is compared to spatially resolved spectroscopy after projection and, when appropriate, PSF convolution (Li et al., 2017).

This cylindrically aligned formulation is the one most often associated with JAM in integral-field studies. It presumes that the tracer kinematics reflect an equilibrium solution in a smooth, time-independent potential, and that the geometry is sufficiently close to axisymmetry for an oblate deprojection to be meaningful. These assumptions make the method computationally efficient and analytically tractable, but they also define the main regimes in which systematic error enters.

A distinct extension, JAM$_{\rm sph}$, assumes that the velocity ellipsoid is aligned with spherical polar coordinates rather than cylindrical ones. In that case the vanishing cross term is $\langle v_r v_\theta\rangle=0$, the anisotropy is written

$$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$

and the axisymmetric Jeans equations are recast in $(r,\theta)$. This variant was developed to model systems for which spherical alignment is empirically or theoretically better motivated, and it complements the original cylindrically aligned formalism rather than replacing it (Cappellari, 2019).

An additional formal extension provides all six projected second moments, not only the line-of-sight term. In observed coordinates $(R,\phi,z)$0, JAM can therefore predict radial-velocity, proper-motion, and cross-moment observables, enabling applications to resolved-star systems and nearby galaxies for which proper motions are available (Cappellari, 2012, Cappellari, 2015).

2. Multi-Gaussian Expansion, deprojection, and observable prediction

A defining feature of JAM is its use of the Multi-Gaussian Expansion to represent projected light, tracer density, and often the total mass. In a standard photometric parameterization,

$(R,\phi,z)$1

where each Gaussian has peak surface brightness $(R,\phi,z)$2, dispersion $(R,\phi,z)$3, and projected axial ratio $(R,\phi,z)$4. For a given inclination, the MGE can be deprojected into a three-dimensional luminosity density, which then serves as the tracer density or, after scaling by a stellar mass-to-light ratio, as the stellar mass density (Li et al., 2017).

The practical workflow is correspondingly standardized. One supplies photometry for the MGE light model, kinematic maps for $(R,\phi,z)$5 and $(R,\phi,z)$6, and parametric forms for the total mass profile and anisotropy. JAM then computes the gravitational potential from the stellar component plus a dark halo or other additional mass components, integrates the intrinsic second moments along the line of sight, and predicts the two-dimensional $(R,\phi,z)$7 field for direct comparison with the data (El-Badry et al., 2016).

The method is flexible with respect to the mass model. Stellar mass is commonly taken to follow the deprojected MGE light up to a factor $(R,\phi,z)$8, while dark matter may be represented either by additional MGE components or by an analytic halo such as a generalized NFW profile. Supermassive black holes can be included as point masses or very compact Gaussian components. Because Gaussian components admit analytic or semi-analytic expressions for the potential and force, the MGE formalism is the main computational reason JAM remains fast across large samples (Cappellari, 2012, Li et al., 2015).

Projection is central to the method. The intrinsic second moments are transformed to the observed frame, integrated along the line of sight, and matched to the sampling of the actual data. For IFU data this ordinarily means PSF convolution and integration over Voronoi bins. For resolved-star applications, the same projected second moments enter discrete likelihoods for line-of-sight velocities and proper motions (Cappellari, 2012, 2206.12121).

3. Parameters inferred, multiple components, and structural degeneracies

In routine use, JAM returns a best-fit total mass model, a stellar mass-to-light ratio, structural halo parameters, anisotropy parameters, and often inclination. In galaxy applications these outputs are frequently interpreted as enclosed mass profiles, dark matter fractions, central density slopes, or IMF-related dynamical mass normalizations; in cluster applications they are used to infer $(R,\phi,z)$9, $\sigma_{Rz}$0, and the full density profile shape (El-Badry et al., 2016, Shi et al., 2024).

A persistent feature of the method is the mass–anisotropy degeneracy. Changes in anisotropy can partially mimic changes in mass normalization or halo structure in the projected second moments, especially when the radial coverage is limited or the inclination is poorly constrained. Large simulation-based assessments show that total mass is generally more robust than the decomposition into stellar and dark components. In the Illustris-based benchmark, the enclosed total mass within $\sigma_{Rz}$1 was recovered with typical accuracy $\sigma_{Rz}$2–$\sigma_{Rz}$3 and negligible median bias, while the recovered $\sigma_{Rz}$4 had a $\sigma_{Rz}$5 scatter of $\sigma_{Rz}$6–$\sigma_{Rz}$7, and the enclosed stellar and dark masses had $\sigma_{Rz}$8 scatters of $\sigma_{Rz}$9–$$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$0 (Li et al., 2015).

Geometry introduces a second major degeneracy. Axisymmetric oblate deprojection is accurate for many regular fast rotators, but prolate or triaxial systems violate the adopted geometry. In the same Illustris analysis, there was no significant bias for oblate galaxies, whereas for prolate galaxies the JAM-recovered stellar mass was on average $$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$1 higher than the input value and the dark matter mass $$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$2 lower, even though the total mass remained relatively well constrained (Li et al., 2015).

The framework has also been generalized to systems with multiple stellar populations. In that setting each stellar component can have its own density profile, flattening, mass, scale length, anisotropy, rotation, age, metallicity, IMF, and mass-to-light ratio. A generalized rotational decomposition was introduced for cases in which the standard Satoh decomposition cannot be applied, specifically when the azimuthal second-moment offset $$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$3 becomes negative. The generalized form,

$$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$4

extends Jeans-based modeling to counter-rotating or multi-component configurations that are awkward to encode in the original single-component JAM setup (Caravita et al., 2021).

These degeneracies imply a methodological principle that recurs throughout the literature: second moments constrain the total potential more directly than they constrain the detailed partition between stars and dark matter, and more directly than they constrain the anisotropy profile.

4. Numerical validation and benchmark performance

Quantitative validation of JAM has been carried out with cosmological simulations, mock clusters, dwarf satellites, and strong-lensing samples, revealing a sharp contrast between cases that respect the equilibrium-axisymmetry assumptions and cases that violate them through triaxiality, non-equilibrium, or feedback-driven potential fluctuations (Li et al., 2015, El-Badry et al., 2016, 2206.12121, Shi et al., 2024, Huang et al., 28 Feb 2025).

Context	Setup	Main result
Illustris early-type galaxies	1,413 simulated galaxies with $$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$5	$$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$6 recovered to $$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$7–$$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$8; $$\frac{\partial (\rho\,\sigma_R^2)}{\partial R} + \rho\,\frac{\sigma_R^2 - \sigma_\phi^2 + \overline{v}_\phi^2}{R} = -\rho\,\frac{\partial \Phi}{\partial R}, \qquad \frac{\partial (\rho\,\sigma_z^2)}{\partial z} = -\rho\,\frac{\partial \Phi}{\partial z},$$9 scatter $\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$0–$\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$1; prolate stellar mass $\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$2 and dark matter $\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$3
FIRE low-mass galaxies	Bursty dwarfs with feedback-driven potential fluctuations	Typical Jeans mass errors $\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$4; excursions $\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$5; quenched gas-stripped rerun recovers masses within $\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$6
Auriga dwarf satellites	28 dwarfs with realistic discrete tracers	$\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$7 ensemble-unbiased with $\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$8 dex scatter; LOS-only 2,000-star case increases scatter to $\beta_z \equiv 1 - \frac{\sigma_z^2}{\sigma_R^2}.$9 dex
TNG300 galaxy clusters	28 massive clusters traced by satellites	Bound-satellite case: $V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$0 bias $V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$1 dex, scatter $V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$2 dex; $V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$3 bias $V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$4 dex, scatter $V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$5 dex
Time-delay lenses	Axisymmetric and triaxial ETG mocks	Axisymmetric JAM recovery $V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$6 in $V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$7 for axisymmetric galaxies and $V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$8 for triaxial galaxies modeled with axisymmetric JAM

Two conclusions follow directly from these tests. First, the total mass near a characteristic tracer radius is usually the most stable quantity. This is explicit in the cluster analysis, where $V_{\rm rms} \equiv \sqrt{V^2 + \sigma^2},$9 is more tightly constrained than $_{\rm sph}$0, and in the dwarf analyses, where mass near the tracer half-mass radius is markedly better recovered than the full inner density profile (Shi et al., 2024, 2206.12121). Second, a formally successful $_{\rm sph}$1 fit does not guarantee a correct physical decomposition. In the FIRE dwarfs, models that fit $_{\rm sph}$2 most tightly did not necessarily yield the correct mass or anisotropy, and in feedback-dominated systems non-equilibrium effects exceeded the classic mass–anisotropy degeneracy when anisotropy was not known a priori (El-Badry et al., 2016).

The strong-lensing benchmark adds an important counterexample. For the slow-rotator ETG regime characteristic of strong-lens samples, spherically aligned axisymmetric JAM produced essentially unbiased aperture-dispersion recovery, and the residual uncertainty associated with triaxiality was subdominant to current measurement errors. In that restricted regime, the alignment and projection approximations are unusually benign (Huang et al., 28 Feb 2025).

5. Scientific applications

JAM has become a general-purpose inference tool for stellar and satellite kinematics. In the MaNGA DR13 analysis of 816 galaxies, including 403 early types and 413 spirals, axisymmetric JAM was combined with stellar population synthesis to measure dynamical mass-to-light ratios and an IMF mismatch parameter. The result was a clear increase of dynamical-to-SPS $_{\rm sph}$3 with $_{\rm sph}$4: low-dispersion early types were consistent with a Chabrier-like IMF, while high-dispersion systems approached a Salpeter-like IMF, and spirals showed a similar trend with larger scatter. Incorporating stellar $_{\rm sph}$5 gradients strengthened the correlation (Li et al., 2017).

In Milky Way work, the spherically aligned axisymmetric variant JAM$_{\rm sph}$6 has been used to fit Gaia-based three-dimensional halo kinematics for K giants and blue horizontal branch stars. The near-spherical alignment of the observed velocity ellipsoid motivated the use of JAM$_{\rm sph}$7, and the preferred model allowed the halo flattening $_{\rm sph}$8 to vary with Galactocentric radius. In the joint KG+BHB fit, the halo became more oblate around $_{\rm sph}$9–$\langle v_r v_\theta\rangle=0$0 kpc and rounder both inside and outside that range, with $\langle v_r v_\theta\rangle=0$1, $\langle v_r v_\theta\rangle=0$2 kpc, and $\langle v_r v_\theta\rangle=0$3 (Zhang et al., 14 May 2025).

At the scale of galaxy clusters, JAM has been adapted from stellar tracers to satellite galaxies. Applied to 28 TNG300-1 clusters with $\langle v_r v_\theta\rangle=0$4, the method recovered $\langle v_r v_\theta\rangle=0$5 and $\langle v_r v_\theta\rangle=0$6 with small ensemble biases even in redshift-space selections that included $\langle v_r v_\theta\rangle=0$7 contamination, and DESI-like fiber incompleteness barely modified the velocity-dispersion profiles relevant for the fit (Shi et al., 2024).

In strong-lensing cosmography, JAM supplies the stellar-kinematic prior that helps break the mass-sheet degeneracy. The central result of the recent projection study is that unresolved spherical JAM can bias the inferred mass-sheet parameter by up to $\langle v_r v_\theta\rangle=0$8–$\langle v_r v_\theta\rangle=0$9, whereas spherically aligned axisymmetric JAM largely removes that bias and leaves a residual random uncertainty of $$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$0–$$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$1 in the stellar velocity dispersion, corresponding to $$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$2–$$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$3 in $$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$4 (Huang et al., 28 Feb 2025).

Edge-on disc galaxies define a different use case. In the GECKOS-MUSE analysis of seven highly inclined discs, JAM fits to $$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$5 were intentionally kept simple in order to isolate the impact of dust and non-axisymmetric structure. When dust was masked appropriately, disc regions reached $$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$6, and enclosed dynamical masses remained constant within $$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$7 across all dust-masking levels. Residual velocity fields then revealed structures that JAM itself could not reproduce, including bar-related residuals in NGC 3957 and IC 1711 and an X-shaped morphology consistent with a boxy/peanut bulge in NGC 0522 (Rutherford et al., 10 Sep 2025).

6. Limitations, misconceptions, and practical use

The central limitation of JAM is not numerical but physical: it assumes equilibrium and symmetry. Whenever the potential is time dependent, the galaxy is strongly triaxial, or the tracer population is not phase mixed, the fitted second moments can cease to represent an equilibrium configuration. The FIRE dwarf study makes this point especially sharply. Feedback-driven gas outflows and inflows change the potential on $$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$8 Myr timescales, the stellar radial dispersion can vary by nearly a factor of two, and the sign of the mass bias depends on the gas-flow phase: post-starburst outflows overestimate mass, while net inflows underestimate it (El-Badry et al., 2016).

This directly corrects a common misconception: a good $$\beta \equiv 1 - \frac{\overline{v_\theta^2}}{\overline{v_r^2}},$$9 map fit does not imply that anisotropy has been correctly recovered, nor even that the mass profile is unbiased. In the same low-mass systems, Jeans modeling alone could not reliably constrain $(r,\theta)$0, and the minimum in a $(r,\theta)$1-based $(r,\theta)$2 surface often favored anisotropy profiles that were not the true ones (El-Badry et al., 2016). The Auriga dwarf analysis reaches a related conclusion in the context of the core–cusp problem: inner masses are ensemble-unbiased, but quiescent Sagittarius-like systems and star-forming systems with strong outflows show opposite systematic trends because contraction and expansion drive them out of equilibrium (2206.12121).

Axisymmetry is the second recurrent caveat. The Illustris assessment shows that triaxial or prolate structure biases the stellar-versus-dark decomposition even when the total mass remains relatively stable, while the GECKOS analysis shows that bars, boxy/peanut bulges, and nuclear discs survive as coherent residuals in otherwise acceptable JAM fits (Li et al., 2015, Rutherford et al., 10 Sep 2025). In such regimes, JAM is best understood as a controlled axisymmetric baseline rather than a complete physical description.

The practical recommendations that emerge across the literature are consistent. For star-forming dwarfs, recent star-formation activity, H$(r,\theta)$3/UV indicators, disturbed morphologies, and evidence for gas outflows or inflows should be treated as warnings that the steady-state assumption is likely violated; for those systems, adopting $(r,\theta)$4 systematic uncertainties on enclosed mass and larger uncertainties on inner slope or core size is recommended, whereas quenched gas-poor dwarfs are substantially more reliable targets (El-Badry et al., 2016). For massive galaxies, high-resolution imaging, extended radial coverage, and independent constraints on $(r,\theta)$5 or the halo reduce the stellar–dark degeneracy, and comparisons between cylindrically aligned JAM and spherically aligned JAM$(r,\theta)$6 provide a direct robustness check on inferred density slopes (Li et al., 2015, Cappellari, 2019). For clusters, external priors from weak lensing or X-ray data improve constraints on outer halo structure and inclination (Shi et al., 2024).

Within those limits, JAM remains one of the most efficient moment-based tools in galactic dynamics. Its durability derives from a specific balance: the assumptions are strong enough to make large-sample inference and semi-analytic projection feasible, yet flexible enough to incorporate anisotropy, multiple components, dark matter, black holes, proper motions, and alternative alignment prescriptions. The method is therefore most powerful when its outputs are interpreted as controlled dynamical constraints tied to clearly stated symmetry and equilibrium assumptions, rather than as assumption-free reconstructions of the underlying distribution function.