Joint Projected Normal & Skew-Normal (JPSN)

• JPSN is a flexible multivariate distribution that unifies circular data via the projected normal and linear data via the skew-normal frameworks.
• It employs latent variable augmentation to facilitate Bayesian inference and efficiently capture complex interdependencies between angular and linear components.
• Its effectiveness is demonstrated in animal movement studies, where it improves state detection and predictive accuracy compared to simpler independent models.

The Joint Projected Normal and Skew-Normal (JPSN) distribution is a class of flexible multivariate models for joint circular-linear or poly-cylindrical data. It combines the projected normal (PN) distribution for circular components and the skew-normal (SN) distribution for linear components, extending multivariate modeling in applications such as animal movement studies where angular (directional) and linear (e.g., speed or length) variables co-occur. The construction leverages latent variable augmentations to facilitate Bayesian inference and Gibbs sampling, addresses non-identifiability in the projected normal component via post-processing, and naturally supports quantification of arbitrary multivariate dependence structures (Mastrantonio, 2015, Mastrantonio, 2017).

1. Constructive Definition

Let $p$ denote the number of circular variables, $$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$, and $q$ the number of linear variables, $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$. The latent variable construction introduces $W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$ for each angular component, with $R_i = \|W_i\| > 0$ and $\Theta_i = \text{atan2}(W_{i2}, W_{i1})$. The linear variables incorporate a skew-normal structure via a latent random effect $d \sim \mathcal{HN}_q(0, I_q)$, so that conditional on $d$, $Y \mid d \sim N_q(\mu_y + \Lambda d, \Sigma_y)$, where $$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$0.

Stacking $$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$1, the augmented vector $$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$2 is modeled jointly as: $$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$3 where

$$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$4

The full (augmented) likelihood is: $$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$5 where $$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$6 is determined by $$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$7. The marginal joint density $$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$8, which integrates out the latent $$\Theta = (\Theta_1, \ldots, \Theta_p)^\top \in [0,2\pi)^p$$9, is intractable; hence inference proceeds using the augmented representation (Mastrantonio, 2015, Mastrantonio, 2017).

2. Relation to Projected Normal and Skew-Normal Distributions

The JPSN construction unifies two well-established frameworks:

• Projected Normal (PN): For $q$0, the circular component $q$1 arises by projection: $q$2 for $q$3. Marginalization over the radius recovers the PN density for angles (Mastrantonio, 2017).
• Sahu-Sahu Skew-Normal (SSN): If $q$4 and $q$5, the linear component $q$6 has a multivariate skew-normal distribution, with augmented density $q$7.

The JPSN replaces the product of separate circular and linear densities for $q$8 with a single $q$9-variate normal, coupling circular and linear components via $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$0, thus encoding full multivariate dependence (Mastrantonio, 2017).

3. Parameters, Identifiability, and Posterior Representation

The parameters consist of:

• $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$1 and $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$2 for circular latent normals;
• $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$3 and $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$4 for the linear component;
• $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$5 capturing block covariance;
• $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$6, the skew-normal parameter vector for linear part.

A structural non-identifiability arises from the PN: scaling $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$7 by any $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$8 leaves $Y = (Y_1, \ldots, Y_q)^\top \in \mathbb{R}^q$9 invariant. Identifiability is often enforced post hoc by setting $W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$0 for each $W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$1 (white noise constraint), so after sampling unconstrained $W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$2, draws are mapped to an identifiable space: $W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$3

$W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$4

All posterior summaries are then reported for $W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$5 (Mastrantonio, 2017, Mastrantonio, 2015).

4. Bayesian Inference and Augmented Gibbs Sampling

The Bayesian model employs conjugate priors: $W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$6 follows a normal-inverse-Wishart prior, and $W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$7 is Gaussian a priori. Inference leverages the tractability of the augmented framework, enabling Gibbs updates as if for standard multivariate normal models:

• Given the latent radii $W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$8 and skew-normals $W_i = (W_{i1}, W_{i2})^\top \in \mathbb{R}^2$9, updating $R_i = \|W_i\| > 0$0 reduces to standard normal-inverse-Wishart conjugate forms.
• Full conditional for $R_i = \|W_i\| > 0$1 is Gaussian, due to its linear regression idiom in $R_i = \|W_i\| > 0$2 versus $R_i = \|W_i\| > 0$3.
• Each $R_i = \|W_i\| > 0$4 is sampled from a truncated normal with mean and covariance determined by the linear block of $R_i = \|W_i\| > 0$5 and $R_i = \|W_i\| > 0$6.
• $R_i = \|W_i\| > 0$7 is sampled by slice sampling or Metropolis-Hastings, exploiting the $R_i = \|W_i\| > 0$8-dimensional structure.

Posterior draws are mapped to the identifiable space. Convergence and mixing are evaluated via trace plots, effective sample size (ESS), and Geweke’s diagnostic (Mastrantonio, 2015, Mastrantonio, 2017).

5. Dependence, Closure, and Interpretability

The JPSN enjoys closure under marginalization: any subset of $R_i = \|W_i\| > 0$9 is again JPSN-distributed for corresponding parameter sub-blocks. This property is inherited from the stability of both the projected normal and skew-normal distributions under subsetting.

Multivariate dependence can be quantified as follows:

• Circular–circular correlation (Jammalamadaka–Sarma):

$\Theta_i = \text{atan2}(W_{i2}, W_{i1})$0

where $\Theta_i = \text{atan2}(W_{i2}, W_{i1})$1 is an i.i.d. copy of $\Theta_i = \text{atan2}(W_{i2}, W_{i1})$2.

• Circular–linear dependence (Mardia’s measure):

$\Theta_i = \text{atan2}(W_{i2}, W_{i1})$3

• Linear–linear dependence is characterized by Pearson correlation.

This allows comprehensive interpretability of joint circular-linear structures (Mastrantonio, 2017).

6. Extensions to Time Series: Hidden Markov Models

Time-dependent extensions model sequences $\Theta_i = \text{atan2}(W_{i2}, W_{i1})$4 as emissions of a $\Theta_i = \text{atan2}(W_{i2}, W_{i1})$5-state hidden Markov model (HMM), where each state $\Theta_i = \text{atan2}(W_{i2}, W_{i1})$6 admits a state-specific JPSN$\Theta_i = \text{atan2}(W_{i2}, W_{i1})$7 emission law. The latent discrete state sequence $\Theta_i = \text{atan2}(W_{i2}, W_{i1})$8 evolves by a Markov chain with transition matrix $\Theta_i = \text{atan2}(W_{i2}, W_{i1})$9 (finite $d \sim \mathcal{HN}_q(0, I_q)$0) or as a stick-breaking process (infinite-state sticky HDP-HMM; sHDP-HMM).

Inference employs beam sampling (slice sampling for the state sequence), combining the JPSN Gibbs steps for within-state parameters. This enables state-dependent modeling of heterogeneous and dependent circular and linear behavior in, e.g., animal movement data, supporting automatic (Bayesian nonparametric) inference for the effective number of states (Mastrantonio, 2015).

7. Applications and Empirical Performance

Applications include modeling movement of free-ranging Maremma Sheepdogs (six animals, $d \sim \mathcal{HN}_q(0, I_q)$1 angles and $d \sim \mathcal{HN}_q(0, I_q)$2 log-step lengths, $d \sim \mathcal{HN}_q(0, I_q)$3 observations) with the sHDP–JPSN–HMM. The posterior mode for the number of behavioral states was three, corresponding to interpretable behaviors: (1) low-speed, tortuous motion (rest or livestock attending), (2) intermediate speed with moderate alignment, and (3) high-speed, nearly straight boundary patrol. The JPSN–HMM reveals intra- and inter-individual dependence in both turning and step-length, whereas simpler HMM models with independent marginals tend to over-segment the behavior and obscure social structure.

Out-of-sample predictive scores (Continuous Ranked Probability Score, CRPS) demonstrated that the JPSN–HMM outperformed alternative models using independent von Mises or wrapped Cauchy (for angles) and Gamma or Weibull (for lengths): mean circular CRPS $d \sim \mathcal{HN}_q(0, I_q)$4 versus $d \sim \mathcal{HN}_q(0, I_q)$5–$d \sim \mathcal{HN}_q(0, I_q)$6, mean linear CRPS $d \sim \mathcal{HN}_q(0, I_q)$7 versus $d \sim \mathcal{HN}_q(0, I_q)$8–$d \sim \mathcal{HN}_q(0, I_q)$9 for the alternatives (Mastrantonio, 2015).

A comparable analysis on zebra movement data ($d$0, $d$1, $d$2) confirmed improved recovery of empirical dependence and predictive accuracy, with the JPSN attaining the lowest CRPS across held-out circular and linear segments (Mastrantonio, 2017).

8. Computational and Diagnostic Considerations

The principal computational bottleneck is the per-iteration cost: $d$3 for normal-inverse-Wishart updates in HMM extensions, driven by the covariance structure. Metropolis steps for $d$4 require tuning. The model’s structure allows for block-wise sparse or diagonal approximations in high dimension. Mixing and convergence are assessed via trace plots and ESS; identifiability is reliably enforced via post-processing (Mastrantonio, 2015, Mastrantonio, 2017).

---

Primary sources:

• (Mastrantonio, 2015) The Joint Projected and Skew Normal
• (Mastrantonio, 2017) The joint projected normal and skew-normal: a distribution for poly-cylindrical data