Skew selection for factor stochastic volatility models

TL;DR

This paper introduces factor stochastic volatility models with skew error distributions, employing Bayesian sparsity to efficiently model skewness in high-dimensional financial data, improving prediction and portfolio performance.

Contribution

It develops a novel Bayesian sparse skew structure for factor stochastic volatility models using the generalized hyperbolic skew t-distribution.

Findings

Skewness is significant in common-factor processes.

Sparse skew structure enhances predictive accuracy.

Improves portfolio performance.

Abstract

This paper proposes factor stochastic volatility models with skew error distributions. The generalized hyperbolic skew t-distribution is employed for common-factor processes and idiosyncratic shocks. Using a Bayesian sparsity modeling strategy for the skewness parameter provides a parsimonious skew structure for possibly high-dimensional stochastic volatility models. Analyses of daily stock returns are provided. Empirical results show that the skewness is important for common-factor processes but less for idiosyncratic shocks. The sparse skew structure improves prediction and portfolio performance.

Tables

Table 1. Table 1: S&P500 Sector Index.

1	INDU	Industrials
2	MATR	Materials
3	ENRS	Energy
4	CONS	Consumer Staples
5	INFT	Information Technology

Table 2. Table 2: Cumulative log predictive density ratios ( LPDR LPDR \mathrm{LPDR} ), relative to Model S0 , for the US-sector indices.

	Horizon ($h$ days)
Model	1	2	3	4	5	Total
SY	4.1	28.6	90.9	26.3	54.4	204.4
SF	172.7	113.0	220.9	146.2	201.2	854.1
SYF	71.6	50.2	75.6	46.8	79.2	323.3
SSYF	195.5	204.2	231.4	146.9	216.8	994.9

Table 3. Table 3: VaR results: the number of violations for the US-sector indices. The value with an asterisk “ * ” shows that the p 𝑝 p -value of the likelihood ratio test is less than 10%. The null hypothesis of the test is that the expected ratio of violations is equal to α 𝛼 \alpha .

		Target return ($\overline{m}$)			Target
Model	$\alpha$	0.005%	0.01%	0.02%	-free
S0	0.5%	5	4	6	6*
	1%	9	8	7	8
	5%	31	29	29	30
SY	0.5%	7*	7*	8*	5
	1%	11*	10*	8	8
	5%	26	27	31	29
SF	0.5%	5	5	7*	3
	1%	6	7	10*	8
	5%	30	29	30	27
SYF	0.5%	4	5	6*	5
	1%	5	6	7	8
	5%	26	26	25	26
SSYF	0.5%	4	3	2	3
	1%	9	7	6	7
	5%	27	28	29	26

Table 4. Table 4: List of US individual stocks. The sector categories: Financials (FI), Industrials (ID), Information Technology (IT), Consumer Discretionary (CD), Consumer Staples (CS), Telecommunication Services (TL), Real Estate (RE), Energy (EN), and Utilities (UT).

No.	Ticker	Name	Sector	No.	Ticker	Name	Sector
1	MS	Morgan Stanley	FI	11	TGT	Target	CD
2	BA	Boeing	ID	12	KO	Coca-Cola	CD
3	STI	Sun Trust Banks	FI	13	GIS	General Mills	CS
4	WFC	Wells Fargo	FI	14	K	Kellogg	CS
5	BLK	BlackRock	FI	15	T	AT&T	TL
6	MET	Met Life	FI	16	AMT	American Tower	RE
7	AXP	American Express	FI	17	COG	Cabot Oil & Gas	EN
8	L	Loews	FI	18	EIX	Edison International	UT
9	IBM	IBM	IT	19	FE	FirstEnergy	UT
10	DIS	Walt Disney	CS	20	AEE	Ameren	UT

Table 5. Table 5: Cumulative returns over 250 business days for portfolios obtained from competing models, with the number of factors p = 1 𝑝 1 p=1 and 2 2 2 .

		Target return ($\overline{m}$)
Model	$p$	0.005%	0.01%	0.02%
S0	1	0.063	0.058	0.047
SF	1	0.092	0.101	0.118
SYF	1	0.120	0.122	0.125
SSYF	1	0.147	0.153	0.165
S0	2	0.109	0.103	0.089
SF	2	0.102	0.089	0.062
SYF	2	0.114	0.116	0.120
SSYF	2	0.121	0.130	0.147

Equations

$$w_{it}=m_i+{\beta }_i{z}_{it}+\sqrt{z_{it}}{\epsilon }_{it},{\epsilon }_{it}\sim N(0,1),$$

$$\displaystyle\mbox{\boldmath$y$}_{t}$$

$$\displaystyle\mbox{\boldmath$f$}_{t}$$

$$\displaystyle\mbox{\boldmath$y$}_{t}$$

$$w_{it}$$

$$\displaystyle\mbox{\boldmath$f$}_{t}$$

$$v_{jt}$$

$$\displaystyle\mbox{\boldmath$h$}_{t+1}\,=\,\mbox{\boldmath$\mu$}+\mbox{\boldmath$\Phi$}(\mbox{\boldmath$h$}_{t}-\mbox{\boldmath$\mu$})+\mbox{\boldmath$\eta$}_{t},$$

$$\displaystyle\left(\begin{array}[]{c}\mbox{\boldmath$\varepsilon$}_{t}\\ \mbox{\boldmath$\eta$}_{t}\end{array}\right)\,\sim\,N\left(\mbox{\boldmath$0$},\left[\begin{array}[]{cc}\mbox{\boldmath$I$}&\mbox{\boldmath$R$}\\ \mbox{\boldmath$R$}&\mbox{\boldmath$S$}\end{array}\right]\right),$$

$$\displaystyle\tilde{y}_{it}\,=\,\left\{\begin{array}[]{ll}y_{it}-\mbox{\boldmath$B$}_{(i)}\mbox{\boldmath$f$}_{t}&(i\leq k),\\ f_{jt}&(i>k;\,j=i-k),\end{array}\right.i=1,\ldots,q,$$

$${\tilde{y}}_{it}$$

$$h_{i,t+1}$$

$$\displaystyle\left(\begin{array}[]{c}\varepsilon_{it}\\ \eta_{it}\\ \end{array}\right)$$

$${\beta }_i\sim \kappa N({\beta }_i|0,{\tau }_0^2)+(1-\kappa ){\delta }_0({\beta }_i),i=1,\dots ,q,$$

$$y_{it}$$

$$f_{jt}$$

$${\alpha }_{\ell t}$$

$${\hat{u}}_{\ell t}$$

$$\displaystyle\hat{y}_{it}\,=\,d_{it}+\mbox{\boldmath$g$}_{it}^{\prime}\mbox{\boldmath$b$}_{i}+\hat{u}_{it},\quad\hat{u}_{it}\sim N(0,\hat{\sigma}_{it}^{2}),$$

$$\displaystyle d_{it}\,=\,\left\{\begin{array}[]{ll}f_{it}&(\mathrm{for}\ i=2,\ldots,p),\\ 0&(\mathrm{for}\ i=p+1,\ldots,k).\end{array}\right.$$

$$\displaystyle\hat{\mbox{\boldmath$W$}}_{i}\,=\,\left(\mbox{\boldmath$W$}_{i0}^{-1}+\sum_{t=1}^{T}\frac{\mbox{\boldmath$g$}_{it}\mbox{\boldmath$g$}_{it}^{\prime}}{\hat{\sigma}_{it}^{2}}\right)^{-1},\quad\hat{\mbox{\boldmath$b$}}_{i}\,=\,\hat{\mbox{\boldmath$W$}}_{i}\left(\mbox{\boldmath$W$}_{i0}^{-1}\mbox{\boldmath$b$}_{i0}+\sum_{t=1}^{T}\frac{(\hat{y}_{it}-d_{it})\mbox{\boldmath$g$}_{it}}{\hat{\sigma}_{it}^{2}}\right).$$

$$\displaystyle\hat{\mbox{\boldmath$y$}}_{t}$$

$$\displaystyle\mbox{\boldmath$f$}_{t}$$

$${\beta }_i|\cdot \sim {\hat{\kappa }}_{i}N({\beta }_i|{\hat{\beta }}_i,{\hat{\tau }}_i^2)+(1-{\hat{\kappa }}_i){\delta }_0({\beta }_i),$$

$$\displaystyle\mathrm{LPDR}_{h}(M_{i})\,=\,\sum_{j=1}^{J}\log\left\{\frac{p_{M_{i}}(\mbox{\boldmath$y$}_{T_{j}+h}|\mbox{\boldmath$y$}_{1:T_{j}})}{p_{M_{0}}(\mbox{\boldmath$y$}_{T_{j}+h}|\mbox{\boldmath$y$}_{1:T_{j}})}\right\},$$

$$2\log \left\{{\left(\frac{n}{N}\right)}^n{\left(1-\frac{n}{N}\right)}^{N-n}\right\}-2\log \left\{{\alpha }^n(1-\alpha {)}^{N-n}\right\}.$$

Full text

Key words: Factor stochastic volatility; Generalized hyperbolic skew t-distribution; Portfolio allocation; Skew selection; Stock returns.

1 Introduction

Multivariate analyses of stochastic volatility (SV) have been widely discussed in both the theoretical and empirical contexts of financial econometrics literature (see e.g., Chib et al., 2009). Among them, the class of factor stochastic volatility (FSV) models is a theoretically and practically effective approach to understanding dynamic variances and covariances in multivariate series of financial variables (Geweke and Zhou, 1996; Aguilar et al., 1999; Pitt and Shephard, 1999; Aguilar and West, 2000; Han, 2005; Chib et al., 2006). Summarizing common dynamics among the series by latent factors yields an efficient modeling strategy for multivariate volatility, and for high-dimensional applications in particular. Previous studies have proposed generalizations on the basis of a standard FSV model, for example, with time-varying factor loadings (Lopes and Carvalho, 2007; Carvalho et al., 2011; Nakajima and West, 2013b; Zhou et al., 2014), and heavy-tailed error distributions (e.g., Chib et al., 2006; Ishihara and Omori, 2017; Li and Scharth, 2018).

The skewness for a return distribution of the financial variables has been recognized as an important aspect of model fit and predictive performance (e.g. Hansen, 1994; Eberlein et al., 1998; Aas and Haff, 2006). The skew error distributions have been incorporated into univariate SV models (Nakajima and Omori, 2012; Abanto-Valle et al., 2015; Langrock et al., 2015; Kobayashi, 2016; Takahashi et al., 2016) and a Cholesky-type multivariate SV (Nakajima, 2017). Little is known, however, about the skew error distribution of the FSV models.

To fill this gap, this paper introduces an FSV model with skew error distribution. The generalized hyperbolic (GH) skew $t$-distribution (see e.g., Barndorff-Nielsen and Shephard, 2001; Aas and Haff, 2006) is incorporated into common-factor processes and the idiosyncratic shocks in the FSV model. It is simple and straightforward to assume that each of the factor processes and idiosyncratic shocks is associated with a different skewness parameter. Given this full model, a key question is whether all of the skewness parameters are necessary to describe the skewness in the multivariate series of financial variables. Because the factor processes summarize the common dynamics in multivariate SV series, some of the skewness parameters in the idiosyncratic shocks may not be relevant, and furthermore, not all of the latent processes may need the skew distribution. If so, then, assuming all possible skewness parameters can lead to more uncertainty in model-fit estimation and prediction, especially for high-dimensional applications.

To address this point, the paper proposes a skew selection strategy for the proposed FSV model that includes the skew error distribution. It employs Bayesian sparsity modeling, which has been widely exploited in the literature (e.g., George and McCulloch, 1993, 1997; West, 2003; Clyde and George, 2004). A so-called spike-and-slab prior selects a zero or non-zero skewness parameter for each factor process and each idiosyncratic shock. This approach allows the proposed FSV model to select-out redundant skewness from the full model and to provide a parsimonious multivariate skew structure. The proposed approach can easily scale up, which has practical relevance for high-dimensional problems. The paper is related to the literature on Bayesian sparsity modeling in the context of multivariate SV for reducing the number of parameters. Loddo et al. (2011) propose a Bayesian stochastic search approach; Kastner (2018) develops a global-local shrinkage prior; and Nakajima and West (2013b) and Zhou et al. (2014) incorporate time-varying sparsity with a latent thresholding process.

This paper provides two empirical analyses using multivariate daily stock returns to illustrate the proposed modeling framework. The first analysis illustrates the key characteristics of the proposed model with a concise dataset of five US-sector indices. A forecasting exercise and a Value-at-Risk (VaR) analysis for portfolio are examined. In the second analysis, a larger-scale dataset of 20 US individual stock returns is analyzed, which highlights the advantages of the skew FSV model and the sparsity skewness modeling in a high-dimensional, practical setting. Portfolio performance in an out-of-sample forecasting exercise is provided. Empirical results from these analyses show that the skewness is important for common-factor processes but less for idiosyncratic shocks. The sparse skew structure contributes to improving stock returns and VaR forecasts as well as portfolio performance.

Li and Scharth (2018) develop the same type of the skew FSV model, incorporating the sparsity prior for the skewness parameter, as well as for the leverage-effect parameter. Their empirical analysis using daily stock returns shows that the FSV model with skewness performs better than a no-skew model for VaR forecasts. The focus of the current paper is different; it concentrates on how a sparse skew structure plays a role in improving the forecasting ability of the skew FSV model. The empirical analyses use several diagnostics of forecasting ability in addition to the VaR measure and reveal that the skew FSV model with the skew selection mechanism performs better than the one without it.

Section 2 formulates the FSV models that includes the skew distribution. Section 3 explains a Bayesian estimation strategy for model fitting. Section 4 provides the first empirical analysis using five US-sector stock indices. Section 5 shows the second analysis using 20 individual stock returns. Section 6 concludes.

2 The model

This section describes the GH skew $t$-distribution in a univariate case (Section 2.1), and then defines the FSV models with skew distribution (Section 2.2). The identification issue for the factor loadings is addressed (Section 2.3), and the sparsity structure for the skewness is introduced (Section 2.4).

2.1 GH skew $t$-distribution

The random variable that follows the GH skew $t$-distribution, defined by $w_{it}$ ($t=1,2,\ldots$) for a series $i$, can be written in the form of normal variance-mean mixture:

[TABLE]

where $z_{it}$ follows the generalized inverse gaussian (GIG) distribution (see e.g., Barndorff-Nielsen and Shephard, 2001; McNeil et al., 2005; Aas and Haff, 2006). This original form of the GH skew $t$-distribution is quite flexible for determining a shape of distribution. Previous studies (Prause, 1999; Aas and Haff, 2006; Nakajima and Omori, 2012) discuss that the parameters in the GH skew-$t$ distribution are typically difficult to jointly estimate. Following Nakajima and Omori (2012), to derive a parsimonious structure of the skewed distribution for the SV models, the current paper assumes that (i) $z_{it}\sim IG(\nu_{i}/2,\nu_{i}/2)$, where $IG$ denotes the inverse gamma distribution; (ii) $m_{i}=-\beta_{i}c_{i}$, where $c_{i}\equiv\mathrm{E}(z_{it})=\nu_{i}/(\nu_{i}-2)$, for $\mathrm{E}(w_{it})=0$; and (iii) $\nu_{i}>4$ for the finite variance of $w_{it}$. The resulting density function describes the shape of distribution with two parameters $\beta_{i}$ and $\nu_{i}$. Although the skewness and kurtosis of the distribution are functions of both these two parameters (see Aas and Haff, 2006), the $\beta_{i}$ and $\nu_{i}$ essentially represent the degree of skewness and kurtosis, respectively.

The literature discusses several alternatives of skew $t$-distributions (Hansen, 1994; Fernández and Steel, 1998; Prause, 1999; Jones and Faddy, 2003; Azzalini and Capitanio, 2003). Specifications differ among the distributions, which yield different shapes of the density, tail behaviors in particular. This paper uses the GH skew $t$-distribution defined above to exploit the form of the normal variance-mean mixture, as it makes Bayesian computation for model-fitting easy and efficient as described later in the paper. While it is of interest to compare different specifications of the skew $t$-distributions in terms of model fit and forecasting, this paper sticks to the GH skew $t$-distribution because the main purpose of the analysis is to investigate the skew structure of financial time series in the context of FSV models. Running a horse race among the different skew $t$-distributions is left for future work.

2.2 Factor stochastic volatility with skew distribution

Define a $k\times 1$ vector of response, \mbox{\boldmathy}_{t}=(y_{1t},\ldots,y_{kt})^{\prime}, for $t=1,\ldots,T$. A standard FSV model is formulated as

[TABLE]

where \mbox{\boldmathf}_{t}=(f_{1t},\ldots,f_{pt})^{\prime} is a $p\times 1$ vector of the latent, common factors; $B$ is a $k\times p$ matrix of factor loadings; \mbox{\boldmath\Lambda}_{t}=\mathrm{diag}(\exp(h_{1t}),\ldots,\exp(h_{kt})), and \mbox{\boldmath\Psi}_{t}=\mathrm{diag}(\exp(h_{k+1,t}),\ldots, $\exp(h_{qt}))$, with $q=k+p$. This form of the model arises from the idea of dimension reduction in the unconditional variance of \mbox{\boldmathy}_{t}, denoted by \mbox{\boldmath\Sigma}_{t}, as \mbox{\boldmathy}_{t}\sim N(\mbox{\boldmath0},\mbox{\boldmath\Sigma}_{t}). With a large $k$, the matrices in eqns. (2) and (3) yield an effective reduction of dimensionality in \mbox{\boldmath\Sigma}_{t}, with \mbox{\boldmath\Sigma}_{t}=\mbox{\boldmathB}\mbox{\boldmath\Psi}_{t}\mbox{\boldmathB}^{\prime}+\mbox{\boldmath\Lambda}_{t}. Because the matrices \mbox{\boldmath\Psi}_{t} and \mbox{\boldmath\Lambda}_{t} are diagonal and a typical choice of $p$ is small, the number of parameters can be reduced effectively.

We build the proposed skew $t$-distribution on the simple form of the FSV model defined in eqns. (2) and (3). A variety of generalizations can be easily introduced, however. A local trend can be added to the observation equation (2); the factor loading $B$ can be time-varying; and the process of \mbox{\boldmathf}_{t} can be dynamic, for example, as it follows a vector autoregressive (VAR) model (see e.g., Chib et al., 2009). The following discussion on incorporating the skewness and its estimation method can be straightforwardly applied to the more elaborated FSV models.

The new class of the FSV models with the GH skew $t$-distribution is defined by

[TABLE]

where each element of \mbox{\boldmathw}_{t} and \mbox{\boldmathv}_{t} follows the GH skew $t$-distribution specified by eqn. (1). Define \mbox{\boldmathh}_{t}=(h_{1t},\ldots,h_{qt})^{\prime}, and \mbox{\boldmath\varepsilon}_{t}=(\varepsilon_{1t},\ldots,\varepsilon_{qt})^{\prime}. The latent SV process is defined as the standard autoregressive form:

[TABLE]

where \mbox{\boldmath\mu}=(\mu_{1},\ldots,\mu_{q})^{\prime}, \mbox{\boldmath\Phi}=\mathrm{diag}(\phi_{1},\ldots,\phi_{q}), \mbox{\boldmath\eta}_{t}=(\eta_{1t},\ldots,\eta_{qt}), and

[TABLE]

where \mbox{\boldmathS}=\mathrm{diag}(\sigma_{1}^{2},\ldots,\sigma_{q}^{2}), and \mbox{\boldmathR}=\mathrm{diag}(\rho_{1}\sigma_{1},\ldots,\rho_{q}\sigma_{q}). This formulation makes each of the processes $\{h_{it}\}_{t=1}^{T}$ follow the univariate autoregressive process. The parameter $\rho$ measures the so-called leverage effect, which is widely addressed in financial econometrics (e.g., Black, 1976; Nelson, 1991; Yu, 2005; Omori et al., 2007). In a stock market, for example, a large price drop tends to lead to a large increase in the volatility the next day. This behavior can be described using a negative correlation ($\rho<0$) between the two shocks $(\varepsilon_{t},\eta_{t})$.

A key aspect of the proposed model structure is that, conditional on (\mbox{\boldmathB},\mbox{\boldmathf}_{t}), the model reduces to $q$-fold state-space representation of the univariate SV model. To see this, define

[TABLE]

where \mbox{\boldmathB}_{(i)} denotes a $1\times p$ vector of the $i$-th row of $B$. Because of the diagonal structure of the matrices (\mbox{\boldmath\Lambda}_{t},\mbox{\boldmath\Psi}_{t},\mbox{\boldmathS},\mbox{\boldmathR}), we obtain the univariate SV model with leverage effect, that is,

[TABLE]

for $i=1,\ldots,q$. In a Bayesian estimation framework, this structure enables us to easily scale up the dimension of the responses for the computation of parameter estimates. With the conditionally independent priors described below, we can obtain an efficient and fast parallel computation.

2.3 Identification

A general modeling issue on the latent factors is restriction for identification. The factor models require a proper parameter restriction on the loading matrix or another valid approach. The most common approach, the one taken here, is the upper triangular form of the factor-loading matrix to order factors according to the selected first $p$ time-series variables (e.g., Aguilar and West, 2000; Lopes and West, 2004). Specifically, we assume that $b_{ii}=1$, for $i=1,\ldots,p$; and $b_{ij}=0$, for $i<j$, and $i\leq p$. Although this approach constrains the model because it depends on the selection and ordering of the first $p$ variables, it is practical to track it with a robustness check that examines different sets of the first $p$ series. Alternative strategies can be found in the existing studies (e.g., Früwirth-Schnatter and Lopes, 2018).

2.4 Sparsity prior for the skewness

The Bayesian sparsity modeling strategy is incorporated into the skewness parameters to develop the skew selection in the FSV model. The sparsity prior for $\beta_{i}$ is given by

[TABLE]

where $\delta_{0}$ denotes the Dirac delta function at zero, and $0<\kappa<1$. This prior assigns the normal distribution for $\beta_{i}$ with probability $\kappa$ and the point mass at $\beta_{i}=0$ with probability $1-\kappa$. This is the so-called spike-and-slab prior that is often used in the literature (e.g., George and McCulloch, 1993, 1997; West, 2003; Clyde and George, 2004). This approach yields an effective skew selection for the FSV model, selecting out redundant skewness parameters from the full set of skewness and providing parsimonious multivariate skew structure.

It is not easy to theoretically derive the skewness of $y_{it}$ in the resulting FSV model. The following brief simulation study addresses how the skewness of $y_{it}$ depends on the skewness parameters for the idiosyncratic shocks and the common-factor processes. A sample of size $T=1{,}000$, $k=3$, $p=2$ is simulated from the skew FSV model. The parameter values are set as $\phi_{i}=0.995$, $\sigma_{i}=0.05$, $\rho_{i}=-0.5$, and $\nu_{i}=8$ for $i=1,\ldots,q$. These values are selected following previous studies (e.g. Nakajima and Omori, 2012; Nakajima, 2017). For the level of volatility, we set $\mu_{i}=-11$, for $i=1,2,3$, while $\mu_{i}=-10$ for $i=4,5$; to assume that the factor volatilities are higher than volatilities of idiosyncratic shocks, which are typically observed in practice. The free parameters in the factor loading matrix $B$ are randomly generated from the uniform distribution $U[0.5,1.5]$.

Figure 1 exhibits summaries of the skewness of $y_{it}$ from 1,000 sets of simulated data for different values of the skewness parameter: \mbox{\boldmath\beta}\equiv(\beta_{1},\ldots,\beta_{q})= (i) $(-1,-1,-1,0,0)$, (ii) $(0,0,0,-1,0)$, (iii) $(0,0,0,-1,-1)$, and (iv) -\mbox{\boldmath1}_{1\times 5}. Case (i) assumes that only the idiosyncratic shocks follow the skew distribution; the resulting skewness of \mbox{\boldmathy}_{t} is not so large as Case (iv) which assumes both the idiosyncratic shocks and the factor processes follow the skew distribution. This reflects the structure of the factor model: $y_{it}$ is the sum of the factor and the idiosyncratic shock. The level of the factor volatility is assumed to be higher than that of idiosyncratic shock in this simulation exercise, the symmetry of the factor processes weakens the degree of skewness caused from the idiosyncratic shock.

Case (ii) assumes that only the first factor process follows the skew distribution. The degree of skewness of the first series (Y1) is larger than that in Case (i), because the skewness in the first factor process dominates the idiosyncratic shock. In contrast, the degree of skewness in the second (Y2) and the third (Y3) series is slightly smaller than the first series, because the second factor, which follows the symmetric distribution but does not appear in the first series, weakens the degree of skewness. The level of volatility for the second factor, however, is typically smaller than that for the first factor in practice, which leads to a smaller difference between the skewness for the first series and the others. The key finding here is that the resulting skewness of \mbox{\boldmathy}_{it} also depends on $\mu$.

Case (iii) assumes that both two factor processes follow the skew distribution, which exhibits a notably similar degree of skewness as Case (iv). From this result, the assumption of $\beta_{i}\neq 0$ for all $i=1,\ldots,q$, as in Case (iv), can be redundant to describe the skewness of \mbox{\boldmathy}_{t}. The skew selection structure with the sparsity prior for $\beta$ is expected to effectively address this point for the skew FSV model, exploring the best parsimonious structure to select out the redundant parameters and improve the prediction performance.

3 Bayesian estimation

The Bayesian analysis and computation are developed for the proposed FSV model based on the Markov chain Monte Carlo (MCMC) method. The sampling estimation scheme can be built on popular efficient conditional samplers for univariate SV models with a leverage effect (Omori et al., 2007; Omori and Watanabe, 2008; Nakajima and Omori, 2012), for the state space dynamic models (e.g., Prado and West, 2010), and for the multivariate FSV models (e.g., Lopes and Carvalho, 2007; Chib et al., 2009; Ishihara and Omori, 2017). Li and Scharth (2018) utilize particle Gibbs with an ancestor sampling algorithm to sample the volatility process in the skew FSV model. In contrast, the current paper shows that the full posterior analysis can be implemented by using the simpler algorithm of multi-move sampler for the volatility without any particle sampling.

Based on observations \mbox{\boldmathy}=\{\mbox{\boldmathy}_{1},\ldots,\mbox{\boldmathy}_{T}\}, the conditional samplers are developed for each of the following latent variables and model parameters:

The $q$-fold univariate SV processes $\{h_{i1},\ldots,h_{iT}\}$, their parameters \mbox{\boldmath\theta}_{i}\equiv\{\mu_{i},\phi_{i},\sigma_{i},\rho_{i}, $\nu_{i}\}$, and mixing latent processes $\{z_{i1},\ldots,z_{iT}\}$ for $i=1,\ldots,q$. 2. 2.

The factor process \mbox{\boldmathf}_{t} and the loading matrix $B$. 3. 3.

The skewness parameter $\beta$ and the sparsity parameter $\kappa$.

A key computational strategy here is decoupling the univariate SV processes in the conditional samplers. Generally, a sampler for the SV process requires a relatively heavier computational burden and makes the mixing of the MCMC algorithm slower. Its simple multivariate extension and, further, with the skewness distribution can be more challenging. However, as mentioned before, the representation of the normal variance-mean mixture leads to the efficient computation with those conditional samplers for the SV processes and their parameters. In addition, the conditional independence of the univariate SV leads to a conditional sampler for $\beta_{i}$ under the sparsity prior.

Each of the conditional sampler is outlined as follows.

Generation of SV and mixing latent process

Based on eqns. (6)–(12), each univariate SV process is generated using the conditional sampler for the SV with leverage effect. Omori and Watanabe (2008) develop the efficient multi-move sampler for this type of non-linear state space model and Nakajima and Omori (2012) utilize it on the SV model with the GH skew $t$-distribution, which is taken in the current analysis. For the generation of the mixing latent process $z_{it}$, the conditional posterior distribution is the product of the kernel of the inverse gamma distribution and another part. A standard Metropolis-Hasting algorithm is employed with a candidate generated from the inverse gamma distribution.

To develop the conditional sampler for the parameters defining each of the SV process, we simply assume traditional priors that are independent across $i$: normal priors for $\mu_{i}$, shifted beta priors for $\phi_{i}$ and $\rho_{i}$, inverse gamma priors for $\sigma_{i}^{2}$, and truncated gamma priors for $\nu_{i}$ with $\nu_{i}>4$. Then, the conditional samplers are implemented with standard posterior distributions for $\mu_{i}$ and the MH algorithm for ($\sigma_{i},\rho_{i}$) and $\nu_{i}$ (see Nakajima and Omori, 2012).

Note that because these latent processes $\{h_{i1},\ldots,h_{iT},z_{i1},\ldots,z_{iT}\}$ and the parameter \mbox{\boldmath\theta}_{i} are conditionally independent across $i$ in the conditional posterior inference, we can implement parallel computing for the generation of $q$-fold SV processes.

Generation of factor process and loadings

Define $\hat{\eta}_{it}=(h_{i,t+1}-\mu_{i})-\phi_{i}(h_{it}-\mu_{i})$, and $\hat{\sigma}_{it}^{2}=z_{it}(1-\rho^{2}_{i})e^{h_{it}}$, for $i=1,\ldots,q$. Conditional on the mixing latent process, SV process and its parameters, we obtain

[TABLE]

where

[TABLE]

for $\ell=1,\ldots,q$.

Define \mbox{\boldmathb}_{i}=(B_{i1},\ldots,B_{ir})^{\prime} as the $r\times 1$ vector of free parameters in the $i$-th row of $B$, with $r=\mathrm{min}(i-1,p)$, for $i=2,\ldots,k$. We generate the sample of \mbox{\boldmathb}_{i}, for $i=2,\ldots,k$, based on the following regression model:

[TABLE]

where $\hat{y}_{it}=y_{it}-\alpha_{it}$, \mbox{\boldmathg}_{it}=(f_{1t},\ldots,f_{rt})^{\prime}, and

[TABLE]

We set a prior as \mbox{\boldmathb}_{i}\sim N(\mbox{\boldmathb}_{i0},\mbox{\boldmathW}_{i0}). Then, the conditional posterior distribution of \mbox{\boldmathb}_{i} is N(\hat{\mbox{\boldmathb}}_{i},\hat{\mbox{\boldmathW}}_{i}), where

[TABLE]

For the generation of \mbox{\boldmathf}_{t}, define \hat{\mbox{\boldmathy}}_{t}=(\hat{y}_{1t},\ldots,\hat{y}_{kt})^{\prime}, and \hat{\mbox{\boldmathu}}_{t}=(\hat{u}_{1t},\ldots,\hat{u}_{kt})^{\prime}. We obtain the following form:

[TABLE]

where \mbox{\boldmath\alpha}_{t}=(\alpha_{k+1,t},\ldots,\alpha_{qt})^{\prime}, and \hat{\mbox{\boldmath\Psi}}_{t}=\mathrm{diag}(\hat{\sigma}_{k+1,t}^{2},\ldots,\hat{\sigma}_{qt}^{2}). We generate \mbox{\boldmathf}_{t} from the conditional posterior distribution, N(\hat{\mbox{\boldmatha}}_{t},\hat{\mbox{\boldmathV}}_{t}), where \hat{\mbox{\boldmathV}}_{t}=(\hat{\mbox{\boldmath\Psi}}_{t}^{-1}+\mbox{\boldmathB}^{\prime}\hat{\mbox{\boldmath\Sigma}}^{-1}\mbox{\boldmathB})^{-1}, and \hat{\mbox{\boldmatha}}_{t}=\hat{\mbox{\boldmathV}}_{t}(\hat{\mbox{\boldmath\Psi}}_{t}^{-1}\mbox{\boldmath\alpha}_{t}+\mbox{\boldmathB}^{\prime}\hat{\mbox{\boldmath\Sigma}}_{t}^{-1}\hat{\mbox{\boldmathy}}_{t}), with \hat{\mbox{\boldmath\Sigma}}_{t}=\mathrm{diag}(\hat{\sigma}_{1t}^{2},\ldots,\hat{\sigma}_{kt}^{2}), for $t=1,\ldots,T$.

Generation of skewness parameters

Before we derive the posterior distribution of the $\beta_{i}$ under the sparsity prior, consider a simple normal prior $N(0,\tau_{0}^{2})$ for $\beta_{i}$. Then, the conditional posterior distribution of $\beta_{i}$ appears to be the normal distribution; we denote it as $N(\hat{\beta}_{i},\hat{\tau}_{i}^{2})$. Turning into the sparsity prior, the conditional posterior distribution is derived as conjugate as

[TABLE]

where $\hat{\kappa}_{i}=\kappa\gamma_{i}/(\kappa\gamma_{i}+1-\kappa)$, and $\gamma_{i}=\exp(\hat{\beta}_{i}^{2}/2\hat{\tau}_{i}^{2})\hat{\tau}_{i}/\tau_{0}$. In this conditional posterior, $\beta_{i}$ follows the normal distribution with probability $\hat{\kappa}$, and shrunk at zero otherwise.

For the generation of $\kappa$, we set a beta prior and derive the conjugate posterior distribution conditional on the number of $\beta_{i}$’s such that $\beta_{i}\neq 0$.

4 A study of US-sector indices

The first real-data analysis applies the proposed FSV model to a series of $k=5$ daily stock returns of US-sector indices. This analysis demonstrates how the skew error distribution and sparsity modeling for the skewness parameter play the role of improving the model fit and forecasting ability of the FSV model.

4.1 Data and setup

The data consist of five selected S&P500 Sector indices, listed in Table 1. The sample period spans from January 4, 2006 to December 29, 2017, resulting in $T=3{,}019$ business days. The returns are computed as the log difference of the daily closing price, and demeaned prior to estimation.

The number of factors is determined as $p=1$ from the conventional factor analysis on the $k=5$ series. The first conventional factor explains 83% of variations, and the second one only an additional 8%. A screen-plot shows a significant drop from the first to the second factor. Further, an analysis of predictive performance, which is explained below, confirms that the model with $p=1$ factors performs better than the one with $p=2$.

In theory, there is no ideal method to determine a selection of the first $p$ variables in \mbox{\boldmathy}_{t}, which technically induce the factor processes under the triangular identification structure embedded in the factor-loading matrix $B$. Ideally, those variables are strongly correlated to other series to efficiently extract the common factors. A simple practical approach taken here is to compute a correlation matrix of \mbox{\boldmathy}_{t} and to sum up each row of the matrix, which represents a degree of correlation with other series. In the dataset, the Industrials sector is the series with the highest degree of correlation, and is selected as the first variable in \mbox{\boldmathy}_{t}, as shown in Table 1.

For the computation, the following priors are assumed: $\mu_{i}\sim N(-11,1)$, $(\phi_{i}+1)/2\sim B(20,1.5)$, $\sigma_{i}^{-2}\sim G(20,0.01)$, $(\rho_{i}+1)/2\sim B(1,1)$, $\nu_{i}\sim G(24,0.8)I[\nu_{i}>4]$, $\beta_{i}\sim\kappa N(\beta_{i}|0,10)+(1-\kappa)\delta_{0}(\beta_{i})$, and $\kappa\sim B(2,2)$, where $B$ and $G$ denotes the beta and gamma distributions, respectively. These priors mainly reflect those assumed in the existing studies (e.g., Nakajima and Omori, 2012; Nakajima, 2017). The size of MCMC iteration is 50,000 samples after the burn-in period of 5,000 samples. To check a convergence of the MCMC draws, the convergence diagnostic of Geweke (1992) is computed. In the following analyses, the null that the Markov chain converges is not rejected at 5 percent significance level.

To measure the efficiency of the MCMC algorithm, the relative numerical efficiency measure, or so-called effective sample size is computed (e.g., Chib, 2001). The effective sample size for selected parameters appears to be mostly less than 100, which is comparable to the standard univariate SV models in the existing studies. This figure shows that the developed MCMC algorithm mixes well and yields a practically efficient estimation scheme.

The analysis focuses on the following five FSV specifications:

• •

Model S0: Symmetric $t$-distribution for \mbox{\boldmathy}_{t} and \mbox{\boldmathf}_{t} ($\beta_{i}\equiv 0$).

• •

Model SY: Skew $t$-distribution for \mbox{\boldmathy}_{t}, symmetric $t$-distribution for \mbox{\boldmathf}_{t} ($\beta_{i}=0$, for $i>p$), no sparsity on $\beta_{i}$ for all $i$ (i.e., $\kappa\equiv 1$).

• •

Model SF: Skew $t$-distribution for \mbox{\boldmathf}_{t}, symmetric $t$-distribution for \mbox{\boldmathy}_{t} ($\beta_{i}=0$, for $i\leq p$), no sparsity on $\beta_{i}$.

• •

Model SYF: Skew $t$-distribution for both \mbox{\boldmathy}_{t} and \mbox{\boldmathf}_{t}, no sparsity on $\beta_{i}$.

• •

Model SSYF: Skew $t$-distribution for both \mbox{\boldmathy}_{t} and \mbox{\boldmathf}_{t}, with the sparsity prior on $\beta_{i}$, for all $i$.

A key focus here is on which skewness structure– for only the idiosyncratic shocks, latent processes, or both of them– is desirable and how the sparsity skew selection works in forecasting the multivariate returns.

4.2 Posterior estimates

Prior to the forecasting exercise, posterior estimates based on a fit of the SSYF model are reported. Figure 2 plots posterior means of the SV process for the factor process and idiosyncratic shocks. The SV of the factor process (F1) exhibits the higher level of volatility and the largest variation among the SV processes, which indicates the variation explained by the latent factor dominates the one by each idiosyncratic shock. The volatilities of five stock indices are highly correlated and their common variation is summarized well by the latent factor. The factor SV increases largely at the time of market turmoil in the global financial crisis around 2009 and the European debt crisis around 2011. The SV of idiosyncratic shock to the Energy sector (Y3) exhibits a different trajectory from other sectors. In particular, it increases relatively higher than the other series in late 2014 and afterwards. This kind of decoupling in the SV process could reflect an industry-specific business environment. The oil price and other commodity price dropped significantly from late 2014 to 2016, which may have driven up the volatility of stock prices in the energy sector.

Figure 3 shows a result of posterior estimates for selected parameters obtained from the SSYF model, with posterior median exhibited in dots and solid line, the 50% and 90% credible intervals in dark and light filled area, respectively. The posterior median of $\phi_{i}$ is quite high, around $0.985$ to $0.998$, which suggests a strong persistence in the SV process. The parameter estimates for the factor is notably different from those for the idiosyncratic shocks. The factor process is characterized by a lower value of $\phi_{i}$ and a higher value of $\sigma_{i}$ and $\mu_{i}$ than the idiosyncratic shocks, which implies a larger variation and a higher level of the factor SV process. In addition, the posterior median of $\rho_{i}$ for the factor is about $-0.4$ with its 90% credible intervals excluding zero. In contrast, those of all the idiosyncratic shocks include zero. This result suggests that the leverage effect is observed only in the factor process, but not in the idiosyncratic shocks. The degree $\nu_{i}$ of freedom in the GH skew-$t$ distribution is lower for the Consumer Staples (Y4) and Information Technology (Y5) sectors than the others, which indicates heavier-tails of idiosyncratic shocks in those sectors.

Figure 4 plots the posterior estimates of the skewness parameter $\beta_{i}$ from Models SYF and SSYF. For both models, the posterior median of $\beta_{i}$ for the factor (F1) is about $-1.0$ and the 90% credible intervals do not include zero. In contrast, $\beta_{i}$’s for idiosyncratic shocks are mostly less relevant. In Model SYF, the 90% credible intervals for Y1 to Y3 include zero. The posterior distributions for those sectors are remarkably shrunk towards zero in Model SSYF. While the 90% credible intervals for Y5 excludes zero in Model SYF, the posterior distribution is partly shrunk in Model SSYF with its posterior median at zero. The posterior probability of $\beta_{i}=0$ is about 80–95% for Y1 to Y3 and 70% for Y5. The Consumer Staples sector (Y4) appears to have relevant skewness in its idiosyncratic shocks with the posterior probability of $\beta_{i}=0$ only 1%.

4.3 Forecasting performance

A recursive forecasting exercise is examined to study the model comparison among the competing specifications. The analysis is based on a sequence of updating data and forecasting returns one to five days ahead, in order to illustrate a practical use of the model in the financial market. Specifically, first, the model parameters are estimated for a period, over $t=1\ldots,T_{1}$, using the posterior computation with the MCMC algorithm. Through the MCMC algorithm, samples of returns $(y_{T_{1}+1},\ldots,y_{T_{1}+5})$ are generated by simulating the latent process from one to five days ahead to obtain the posterior predictive density $\pi(y_{T_{1}+h}|y_{1},\ldots,y_{T_{1}})$ for $h=1,\ldots,5$. We then update the returns for five business days and estimate the model for the returns ($y_{1},\ldots,y_{T_{1}+5}$) using the full MCMC analysis with the computation of posterior predictive density for the next five days, $(y_{T_{2}+1},\ldots,y_{T_{2}+5})$, where $T_{2}=T_{1}+5$. We repeat this prediction exercise recursively for $T_{j}=T_{j-1}+5$, $j=2,\ldots,J$, until we obtain forecasts for 500 business days, that is, $J=100$.

We set $T_{1}=2{,}519$, which corresponds to the first pre-sample period, to January 6, 2016; the forecasting exercise period goes from January 7, 2016 to December 29, 2017. This experiment reflects practical forecasting activity in the financial market, given the historical data over about 10 years at the starting point, recursively forecasting returns in the following week with the updated data. The exercise period is quite long and includes periods of both bull and bear market sentiments.

The forecasting performance of each model is compared based on two diagnostics. The first one is out-of-sample predictive density. Specifically, the log predictive density ratio (LPDR) for Model $i$’s $h$-day ahead forecast is computed as

[TABLE]

where $p_{M}$ denotes the predictive density under model $M$, and $T_{1}$ denotes the set of time points for the recursive forecasting. This quantity measures how high the likelihood is, when evaluated at the realized return (\mbox{\boldmathy}_{T_{j}+h}), that is, how accurate the predictive density is, relative to the base model $M_{0}$.

Table 2 reports the result of the LPDR for five models relative to Model S0, showing the LPDR at each horizon $h=1,\ldots,5$ and the sum of all the horizons. Because all the LPDRs are positive, the forecasting performance of the newly introduced skew-$t$ FSV models is superior to the standard FSV model with the symmetric $t$-distribution. Model SF performs better than Model SY, which suggests that the skewness in the factor process is more important for prediction than the ones in the idiosyncratic shocks. Interestingly, Model SYF performs worse than Model SF. This result suggests that the statistically irrelevant skewness parameters for the idiosyncratic shocks, as observed in the left panel of Figure 4, worsens the predictive ability of the skew FSV model, unless the shrinkage method is applied. The LPDRs of Model SSYF, which is the highest among the competing models for all the horizons, confirms that the shrinkage method for the skewness parameters effectively selects statistically relevant ones and consequently improves the forecasting performance.

The second forecasting exercise is based on Value-at-Risk (VaR). The VaR is a standard tool to measure a risk of possible shortfalls in portfolio. We compute the VaR for portfolios that consist of the five stock indices using the forecasts from each competing model, and evaluate the accuracy of the VaR forecast by comparing it with the realized loss of the portfolio return. This exercise focuses on the prediction ability of the skew-$t$ FSV models for a tail risk of the stock returns, which has practical importance in portfolio management.

We consider a standard portfolio allocation strategy based on a mean-variance optimization following Markowitz (1959). The analysis uses standard Bayesian mean-variance optimization. Let \mbox{\boldmathm}_{t} and \mbox{\boldmathD}_{t} denote the forecast mean vector and variance matrix of \mbox{\boldmathy}_{t}, respectively. Note that these variables are computed from the posterior predictive distribution obtained in the recursive forecasting, as illustrated in the previous exercise. Define \mbox{\boldmath\omega}_{t} as a vector of portfolio weights, and $\bar{m}$ as a target level of return. We solve the optimization problem with respect to the portfolio weights \mbox{\boldmath\omega}_{t} as it minimizes the forecast variance of the portfolio return among the restricted portfolios whose expectation is equal to $\bar{m}$. We also assume that the amount of money invested into the portfolio is fixed during the sample period, and that the resource can be freely reallocated to arbitrary long or short positions without any transaction cost.

The mathematical form of the optimization problem is minimizing \mbox{\boldmath\omega}_{t}^{\prime}\mbox{\boldmathD}_{t}\mbox{\boldmath\omega}_{t}, subject to \mbox{\boldmath\omega}_{t}^{\prime}\mbox{\boldmathm}_{t}=\bar{m}, and \mbox{\boldmath\omega}_{t}^{\prime}\mbox{\boldmath1}=1. The solution is \mbox{\boldmath\omega}_{t}^{(\bar{m})}=\mbox{\boldmathK}_{t}(\mbox{\boldmath1}^{\prime}\mbox{\boldmathK}_{t}\mbox{\boldmathq}_{t}\mbox{\boldmathm}_{t}-\mbox{\boldmathm}_{t}^{\prime}\mbox{\boldmathK}_{t}\mbox{\boldmathq}_{t}\mbox{\boldmath1}), where \mbox{\boldmathq}_{t}=(\mbox{\boldmath1}\bar{m}-\mbox{\boldmathm}_{t})/d_{t}, and d_{t}=(\mbox{\boldmath1}^{\prime}\mbox{\boldmathK}_{t}\mbox{\boldmath1})(\mbox{\boldmathm}_{t}^{\prime}\mbox{\boldmathK}_{t}\mbox{\boldmathm}_{t})-(\mbox{\boldmath1}^{\prime}\mbox{\boldmathK}_{t}\mbox{\boldmathm}_{t})^{2}, with \mbox{\boldmathK}_{t}=\mbox{\boldmathD}_{t}^{-1}. Note that the portfolio is reallocated every business day. This experiment assumes a practical situation in which investors allocate their resource every day based forecasts updated weekly. The analysis uses a range of daily target returns of $\bar{m}=0.005\%$, $0.01\%$, and $0.02\%$, implying a yearly return of approximately $1.25\%$, $2.5\%$, and $5.0\%$, respectively. In addition, we consider a target-free portfolio that is implemented by minimizing the optimization problem only with \mbox{\boldmath\omega}_{t}^{\prime}\mbox{\boldmath1}=1.

For a diagnostic of model comparison in this VaR forecast exercise, we employ the test statistic of Kupiec (1995). Let $n$ denote the number of VaR violations, and $N$ the total number of forecasting days. Note that the expected number of violations for $\alpha$ quantile is $\alpha N$. Under the null hypothesis that the expected ratio of violations is equal to $\alpha$, the likelihood ratio statistic is given by

[TABLE]

As this statistic is asymptotically distributed as $\chi^{2}(1)$, we use it to test whether the VaR forecast is accurate (see Kupiec, 1995). The analysis examines the likelihood ratio test for $\alpha=0.5\%$, $1\%$, and $5\%$ levels.

Table 3 reports the number of VaR violations and results of the likelihood ratio test. The value with the asterisk (“*****”) shows that the null hypothesis is rejected at a 10% significance level. For Model S0, the target-free portfolio does not pass the test. For Models SY, SF, and SYF, some of the target-return portfolios exhibit a rejection of the test. In particular, the number of violations in Model SY indicates its optimistic VaR forecast, possibly due to the lack of skewness parameter in the factor process. Finally, Model SSYF provides no portfolio that is rejected in the selected range of portfolio allocation rules and VaR levels. This confirms that the proposed skew selection strategy works remarkably well for forecasting VaR in portfolio risk management.

5 A higher-dimensional study: US individual stocks

The second empirical analysis uses daily returns of 20 individual stocks to illustrate a high-dimensional application of the skew FSV models. The stocks are randomly selected from the compositions of S&P500 index, and listed in Table 4. The selected set of securities covers major industrial sectors. The sample period is $T=751$ business days from January 3, 2014 to December 30, 2016. As in the previous analysis, the returns are computed as the log difference of the daily closing price.

We examine the number of factors $p=1$ and $2$ from the conventional factor analysis. The first and second conventional factors explain more than 50% of variations in total. In addition, an additional factor does not appear to improve the return performance in the following forecasting exercise. The Morgan Stanley and Boeing are selected for the first two series of \mbox{\boldmathy}_{t} for the $p=2$ factors. These two countries from the Financials and Industrials sectors are considered to parsimoniously summarize common movements among the individual stock series. The model with $p=1$ loads Morgan Stanley for the first factor. Note that a robustness check with randomly selected stocks for the first two series revealed that this selection of the first two series do not significantly affect the main result and its implication in this analysis.

A model comparison is examined based on out-of-sample forecasting performance in the recursive forecasting exercise in the same manner as the previous analysis. The posterior predictive distributions are computed for one- to five-day horizons. The first set of data spans to January 5, 2016 and the forecasting period is from January 5, 2016 to December 30, 2016, to obtain the forecasts over 250 business days. Cumulative returns from the optimized portfolios allocations with the target returns ($\bar{m}$) are compared among the competing models. Model SY is excluded from the model comparison because it performs relatively less well in the previous analysis. For the posterior computation, all the prior settings are the same as the previous analysis.

Table 5 reports the cumulative returns over the forecasting period for portfolios with three target returns. Model SSYF provides the uniformly highest cumulative returns both for $p=1$ and $p=2$ among the competing models for each target return. The return from the $p=1$ SSYF model performs slightly better than the $p=2$ model. This result confirms that the skew FSV model with parsimonious skew selection structure improves the forecasting performance. Interestingly, Model SF is outperformed by Model S0 for $p=2$, while Model SYF outperforms Model S0, which indicates that the skewness for idiosyncratic shocks is partly relevant to describe dynamics of the high-dimensional multivariate returns.

Figure 5 shows posterior estimates of selected parameters from Model SSYF with $p=1$. The posterior estimates of $\phi_{i}$ indicates quite high persistence for the SV processes of both the factors and idiosyncratic shocks, while some of the idiosyncratic shocks exhibit relatively wider credible intervals. For the parameter $\rho_{i}$, the posterior median of the factor is about $-0.7$, suggesting that the leverage effect exists as a common behavior among the individual stocks. Additionally, several idiosyncratic shocks (Y1, Y7, Y17) exhibit relevant leverage effect. The posterior estimates of the degree $\nu_{i}$ of the freedom are quite low for several stocks, with some of the posterior medians being below 10. This estimate indicates that the price movement of individual stocks is relatively more volatile and the return distribution is more heavily-tailed than the sector indices examined in the previous analysis.

The estimates of $\beta_{i}$ show considerable shrinkage for most of the idiosyncratic shocks. While some of them exhibit that the posterior medians or credible intervals depart from zero, 16 of 20 idiosyncratic shocks exhibit a posterior median at zero. For Y18 to Y20, the skewness is relatively more relevant than the other series with the 50% credible intervals excluding zero. The result in the forecasting performance partly reflects this relevant skewness for idiosyncratic shocks. The $\beta_{i}$ for the factor is around $-1.0$, which indicates the importance of the skewness in the factor process.

6 Concluding remarks

This paper develops the FSV model with GH skew-$t$ distribution. The Bayesian sparsity prior is embedded to the skewness parameters, which provides a parsimonious skew structure for the FSV model. Two empirical analyses using daily stock returns of US-sector indices and individual stocks show that the advantage of the skew selection structure that effectively shrinks irrelevant skewness parameters and consequently improves forecasting performance.

Possible extensions can be explored in several directions. The factor loadings can be time-varying and the factor process can be more complicated dynamics, such as VAR and VARMA processes. From another perspective, it is of interest to model a time-varying skewness that take non-zero values when relevant and zero otherwise using time-varying sparsity techniques (e.g. Nakajima and West, 2013a; Huber et al., 2019).