Enhancing Insurer Portfolio Resilience and Capital Efficiency with Green Bonds: A Framework Combining Dynamic R-vine Copulas and Tail Risk Modeling

1. Introduction

Robust risk modeling is crucial for modern insurance portfolio management, especially for life and non-life insurers navigating volatile markets. Traditional approaches often fail to capture tail risks and nonlinear dynamics, particularly during stress periods. As regulatory demands tighten under solvency frameworks, insurers increasingly require advanced modeling techniques to accurately assess volatility and extreme losses (Bollerslev, 1986; Daly, 2008).

Insurance companies are increasingly integrating Environmental, Social, and Governance (ESG) considerations into their investment strategies to strengthen portfolio resilience and support long-term sustainability goals. Among ESG instruments, green bonds have emerged as a compelling asset class, channeling funds into projects such as renewable energy and low-carbon infrastructure (Climate Bonds Initiative, 2022). Their dual potential to generate financial returns and promote sustainability has attracted insurers seeking diversification and alignment with ESG mandates (Ferrer et al., 2021; Han et al., 2024; Karim et al., 2024; Park et al., 2020; Papavassiliou et al., 2025). Green bonds have also been associated with improved Sharpe ratios and reduced capital charges under solvency frameworks (Taghizadeh-Hesary et al., 2021; Zhang et al., 2021). In Thailand, where insurers manage sizable and duration-sensitive portfolios, green bond inclusion supports both financial stability and national sustainability objectives (Fitrah and Soemitra, 2022; Ramadhan, 2020). Regulatory support is central to enabling this transition. Tools such as green taxonomies, favorable capital treatment, and targeted incentives can accelerate ESG adoption across the sector (Flammer, 2020; Huang and Lin, 2023; Okeke et al., 2024). Recognizing the strong risk-adjusted performance of green bonds, the Office of Insurance Commission (OIC) could facilitate their broader adoption by incorporating sustainability-related investment guidance into supervisory frameworks. This would help insurers align investment strategies with solvency goals while advancing Thailand’s sustainable finance agenda.

Advanced financial econometrics provides crucial tools for modeling the complex risks faced by insurers. The Autoregressive Moving Average Generalized Autoregressive Conditional Heteroskedasticity (ARMA–GARCH) model, especially with the Glosten–Jagannathan–Runkle (GJR) extension, is valuable for capturing volatility clustering and asymmetric shock responses (Adegboyo and Sarwar, 2025; Glosten et al., 1993; Hidayana et al., 2021; Ma et al., 2024; Wei et al., 2025). The GJR–GARCH model accounts for the leverage effect, where negative shocks induce greater volatility than positive ones (Liu and Hung, 2010). Financial returns often exhibit skewness and fat tails, which, if unaccounted for, can lead to underestimation of risk. To address this, skewed Student-t innovations are widely adopted to improve Value-at-Risk (VaR) accuracy and better capture these distributional properties (Akanbi et al., 2025; Al-Khasawneh et al., 2024; Hansen, 1994; Harvey and Siddique, 1999; Lambert and Laurent, 2001; Patra and Gupta, 2025). However, even these enhanced models may underestimate extreme losses. Extreme Value Theory (EVT), particularly through the Generalized Pareto Distribution (GPD), enhances tail risk modeling by focusing on exceedances beyond a high threshold (Braione and Scholtes, 2016; Majumder, 2018). When applied to GARCH-filtered residuals, EVT improves VaR estimation under stress conditions, which is crucial for insurers operating under solvency requirements (Chen et al., 2024; Okou and Amar, 2023).

Insurers require methods that capture interdependence among assets, moving beyond univariate modeling. Correlation matrices often fail to reflect nonlinear and time-varying relationships, particularly during periods of financial distress. Copula models address this limitation by constructing flexible joint distributions. R-vine copulas are well suited for modeling complex dependency structures using sequences of bivariate copulas (Brechmann and Czado, 2013). Their dynamic extensions allow dependencies to evolve over time, which is essential for capturing contagion effects and shifts in market conditions (Raza et al., 2025; Zhou and Ji, 2021). When integrated with GARCH-based marginal models, this approach offers significant advantages for multivariate risk modeling. It enhances portfolio-level risk assessment and supports the design of capital-efficient, ESG-aligned investment strategies.

This study integrates ARMA–GJR–GARCH models with skewed Student-t innovations, EVT, and dynamic R-vine copulas into a unified framework tailored for insurer portfolios. Capturing volatility dynamics, tail risks, and evolving interdependencies, the framework improves VaR estimation and Sharpe ratio assessment under both normal and stressed conditions (Ahmadi et al., 2023; Han and Li, 2022; Jeleskovic et al., 2024). Using daily data on Thai financial assets and green bonds from 2014 to 2024, the study evaluates VaR forecasts and rolling Sharpe ratios across optimized portfolios. Results reveal that green bond inclusion enhances insurer portfolio performance, strengthens capital adequacy, improves risk-adjusted returns through diversification, and supports alignment with Thailand’s sustainability goals. This research offers practical insights for risk practitioners and regulators in emerging markets. The proposed framework enables more accurate capital requirement calculations and supports ESG integration within supervisory investment guidelines, crucial as Thai regulatory guidance, market maturity, and sustainability commitments rapidly evolve.

The remainder of this paper is structured as follows: Section 2 describes the dataset and methodology, including model specifications and validation procedures; Section 3 presents empirical findings and discussions on VaR forecasting, backtesting, and portfolio performance analysis; and Section 4 concludes with key findings, policy implications, and suggestions for future research directions.

2. Domain of Experiment and Methodology

This study analyzes a comprehensive dataset of 2,869 daily observations spanning January 1, 2014, to December 31, 2024, using logarithmic returns. Logarithmic returns were calculated using the formula $r_t=\mathrm{l}\mathrm{n}(P_t/P_{t-1})$, where $P_t$ is the price at time t and $P_{t-1}$ is the price at time t$-$ 1. The dataset encompasses a variety of significant Thai financial assets: the SET index, the Dubai crude oil, the Thai bullion gold, the 3–7-year government bond index, the 7–10-year government bond index, the JPY/THB exchange rate, the property sector index, and the Bloomberg Barclays MSCI US green bond index. Data obtained from Datastream International and Bloomberg exhibit statistical properties typical of financial time series. With the exception of the exchange rate, the series display negative skewness, indicating a higher likelihood of large negative returns. High kurtosis values confirm leptokurtosis in most series, indicating a greater concentration around the mean and fatter tails than expected under normality, suggestive of increased extreme events. The Jarque–Bera (JB) test indicates that gold returns are normally distributed, unlike the other assets. The Augmented Dickey–Fuller (ADF) test results confirm that all series are stationary. Table 1 presents the summary statistics of daily returns for all eight assets examined in this study.

2.1. Autoregressive Moving Average Glosten–Jagannathan–Runkle Generalized Autoregressive Conditional Heteroskedasticity with Skewed Student-t Innovations

This study addresses key complexities in financial time series, such as volatility clustering, fat tails, and asymmetric shock responses, across eight asset classes, including green bonds. To model individual asset dynamics, it adopts the Autoregressive Moving Average–Glosten–Jagannathan–Runkle Generalized Autoregressive Conditional Heteroskedasticity (ARMA–GJR–GARCH) model with skewed Student-t innovations, jointly modeling the conditional mean and variance of returns. Portfolios of five and six assets (including green bonds) are then constructed by maximizing the Sharpe ratio using Dynamic R-vine copulas. This integrated approach offers a flexible and robust structure for volatility modeling and risk forecasting. To capture excess kurtosis and skewness common in financial return distributions, the innovation term follows the skewed Student-t distribution introduced by Hansen (1994). This distribution accounts for both fat tails and asymmetry, enhancing the precision of Value-at-Risk (VaR) estimates under both normal and extreme market conditions.

The GJR component is central to capturing asymmetric volatility in response to shocks, an effect often observed in financial markets. In the general ARMA(m,n)–GJR–GARCH(p,q) model, p represents the number of ARCH terms (lagged squared residuals) and asymmetry terms, while q denotes the number of GARCH terms (lagged conditional variances). Estimating mean and variance with skewed Student-t innovations, the model is given by the following:

$$\left\{\begin{array}{c}{r}_t = \mu +{\sum }_{i=1}^m{\varnothing }_i r_{t-i}+{\sum }_{j=1}^n{\theta }_{j }{\epsilon }_{t-j}+ {\epsilon }_t\\ {\epsilon }_t = {\sigma }_t{z}_t, z_t~Skewed-t\left(\mathrm{0,1},\nu ,\eta \right) \end{array}\right.$$

(1)

In this model, $r_t$ represents the return of the time series at time t, and $\mu$ is the constant mean. The autoregressive coefficients ${\varnothing }_i$ capture the influence of past returns (AR part of order m), while the moving average coefficients ${\theta }_{j }$account for the impact of past shocks (MA part of order n). The innovation term ${\epsilon }_t$ represents unpredictable shocks, modeled as the product of the conditional standard deviation ${\sigma }_t$ and a shock $z_t$, which follows a skewed Student-t distribution with zero mean and unit variance.

The skewed Student-t distribution is parameterized by degrees of freedom $\nu$> 2 to ensure finite variance and a skewness parameter $\eta \in (-1, 1)$, where $\eta =0$ corresponds to the symmetric Student-t distribution. Its density function $f(z_t;\nu ,\eta )$ is defined as follows:

$$f\left(z_t;\nu ,\eta \right)=\left\{\begin{array}{c}bc{\left[1+{\displaystyle \frac{1}{\nu -2}}{\left({\displaystyle \frac{b{z}_t+a}{1-\eta }}\right)}^2\right]}^{-{\displaystyle \frac{\nu +1}{2}}},z_t<-{\displaystyle \frac{a}{b}}\\ bc{\left[1+{{\displaystyle \frac{1}{\nu -2}}\left({\displaystyle \frac{{bz}_t+a}{1+\eta }}\right)}^2\right]}^{-{\displaystyle \frac{\nu +1}{2}}}, z_t\ge -{\displaystyle \frac{a}{b}}\end{array}\right.$$

(2)

where constants $a, b$, and $c$ are given by the following:

$$\left\{\begin{array}{c}a=4\eta c{\displaystyle \frac{\nu -2}{\nu -1}} \\ b=\sqrt{1+3{\eta }^2-a^2} \\ c={\displaystyle \frac{\mathsf{\Gamma }\left(\frac{\nu +1}{2}\right)}{\sqrt{\pi (\nu -2)}\cdot \mathsf{\Gamma }\left(\frac{\nu }{2}\right)}} \end{array}\right.$$

(3)

The conditional variance ${\sigma }_t^2$ equation under the GJR–GARCH(p,q) specification is as follows:

$${ \sigma }_t^2= \omega + {\sum }_{i=1}^p{\alpha }_i{\epsilon }_{t-i}^2+{\sum }_{i=1}^p{\gamma }_i{\epsilon }_{t-i}^2{I}_{\left\{{\epsilon }_{t-i<0}\right\}}+ {\sum }_{j=1}^q{\beta }_j{\sigma }_{t-j}^2$$

(4)

The constant $\omega$ helps ensure the variance remains positive. The ARCH terms ${\alpha }_i$ represent the effects of lagged squared shocks, while the GARCH terms ${\beta }_j$ capture the persistence of volatility through lagged conditional variances. The asymmetry terms ${\gamma }_i$, together with the indicator function $I_{\left\{{\epsilon }_{t-i<0}\right\}}$, allow the model to distinguish between the effects of positive and negative shocks. This indicator equals 1 when the lagged shock is negative, and 0 otherwise.

Under this specification, a positive shock contributes ${\alpha }_i{\epsilon }_{t-i}^2$ to the conditional variance, whereas a negative shock contributes ${(\alpha }_i+{\gamma }_i){\epsilon }_{t-i}^2$. A positive and statistically significant ${\gamma }_i$ indicates a leverage effect, meaning that negative shocks have a greater impact on volatility than positive shocks of the same magnitude. The GJR–GARCH model is particularly relevant for analyzing assets like green bonds, which may react asymmetrically to environmental developments or policy shifts. While the ARMA component captures serial dependence in returns, GJR–GARCH effectively models both volatility persistence and its asymmetric responses.

Model selection involved a systematic search over ARMA–GJR–GARCH specifications by varying the parameters m, n, p, and q. Each model was estimated using pseudo-maximum likelihood and validated through residual diagnostics to ensure zero mean, homoscedasticity, and no autocorrelation. The model with the lowest Akaike Information Criterion (AIC) was selected, reflecting the best trade-off between fit and simplicity.

2.2. Extreme Value Theory

In financial risk management, accurately capturing extreme tail behavior is essential for stress testing and capital adequacy. This study first filters serial dependence and conditional heteroskedasticity using ARMA(m,n)-GJR-GARCH(p,q) models with skewed Student-t innovations. The resulting standardized residuals $z_t={\epsilon }_t/{\sigma }_t$ are then used for tail modeling. Therefore, this study adopts Extreme Value Theory (EVT), specifically the Peak Over Threshold (POT) method, to model the tail distribution of these standardized residuals.

EVT offers a more refined framework for analyzing extreme outcomes and estimating tail-related risk measures such as VaR, especially under stressed conditions (Degen and Embrechts, 2008; Reiss and Thomas, 2007). Accurately modeling such extreme movements is crucial for reliable risk assessment and capital determination. The effectiveness of EVT in financial applications is well-documented, with notable contributions by Ayusuk and Sriboonchitta (2015), Melina et al. (2024), Muela et al. (2023), Roy (2022), Singh et al. (2013), and Uluceviz (2025).

Building on this foundation, the present study incorporates the POT method alongside the ARMA–GJR–GARCH modeling framework to analyze the tail behavior of standardized residuals using the Generalized Pareto Distribution (GPD). This integrated approach enhances the understanding of extreme market events in light of underlying volatility dynamics and supports more accurate estimation of capital needs.

2.2.1. Peak Over Threshold

To estimate one-day VaR under stressed conditions, this study applies the POT method to eight financial assets, including green bonds. Following Nortey et al. (2015), asset returns are first filtered using ARMA–GJR–GARCH models to capture time-varying volatility. Standardized residuals are then used as input for tail risk modeling. This approach is consistent with established EVT applications in finance (Rosso, 2015) and enhances the robustness of stress testing and capital requirement assessments.

Let $X$ denote a sequence of financial returns with cumulative distribution function $F\left(x\right)$, and let $u$ be a predefined high threshold. POT focuses on the distribution of exceedances above this threshold, specifically the conditional distribution of $X$ given that $X>u$. An exceedance occurs when $X>u$. Defining the excess as $y=X-u$, the conditional distribution function $F_u\left(y\right)$ is given by the following:

$$F_u\left(y\right)=\mathrm{Pr}(\left(X-u\right)\le \left.y\right|X>u)={\displaystyle \frac{F\left(u+y\right)-F\left(u\right)}{1-F\left(u\right)}}, 0\le y\le x_F-u$$

(5)

where $x_F$ is the right endpoint of $F$.

The POT method assumes a sequence of independent and identically distributed (i.i.d.) losses $X_i$. It models the conditional excess distribution $F_u\left(y\right)$ for values exceeding the threshold $u$ (Magnou, 2017). As noted by Gharib et al. (2017), such exceedances are often modeled using the Generalized Pareto Distribution (GPD), which is well-suited for capturing extreme tail risk and is critical in stress testing and VaR estimation (Basu, 2011). The GPD is widely used in finance to model low probability events (Sharpe and Juarez, 2019).

Following the ARMA–GJR–GARCH modeling, the next step is to characterize the distribution of extreme values for each asset. The goal is to find a suitable parametric distribution to model $F_u\left(y\right)$, with threshold $u$ chosen to ensure that excesses are well-approximated by the GPD. Given the sensitivity of GPD quantiles to threshold choice (Muela et al., 2023), a sample mean excess function can be used to determine the threshold (Nurhadi, 2016). Accurate tail modeling is particularly critical in this study, as it involves financial assets such as green bonds. Underestimating risk for these assets could lead to significant misjudgments.

The foundational work of Balkema and De Haan (1974) and Pickands (1975) demonstrates that given a sufficiently high threshold $u$, the distribution of exceedances $y=X-u$, conditional on $X>u$, converges to the GPD. This convergence is mathematically expressed as follows:

$$G_{\xi ,\beta }\left(y\right)=\left\{\begin{array}{c}1-{\left(1+ \xi {\displaystyle \frac{y}{\beta }}\right)}^{-{\displaystyle \frac{1}{\xi }}} if \xi \ne 0 \\ 1- e^{-{\displaystyle \frac{y}{\beta }}} if \xi =0 \end{array}\right.$$

(6)

where $\xi$ and $\beta$ represent the shape and scale parameters, respectively. The flexibility of the GPD lies in its ability to model various tail behaviors: when $\xi >0$, it models heavy-tailed distributions (e.g., Pareto distributions); when $\xi \to 0,$ it simplifies to the exponential distribution; and when $\xi <0$, it reflects a bounded, short-tailed distribution, sometimes referred to as the Pareto Type II.

2.2.2. Estimation of Value-at-Risk

Traditional risk measures often fail to adequately capture the dynamics of financial time series during extreme market conditions (Omari et al., 2020). EVT provides a robust framework for modeling the tails of distributions, making it especially valuable for estimating VaR (Zhang and Zhang, 2016).

Given that the exceedances over a sufficiently high threshold $u$ converge in distribution to the GPD, the cumulative distribution function $F\left(x\right)$ for values above $u$ can be expressed as follows:

$$F\left(x\right)=\left[1-F\left(u\right)\right]F_u\left(y\right)+F\left(u\right), X>u$$

(7)

where $F_u\left(y\right)$ is the distribution of exceedances $y=X-u$ given $X>u$.

Substituting the GPD approximation for $F_u\left(y\right)$, the tail of $F\left(x\right)$ becomes:

$$F\left(x\right)=\left[1-F\left(u\right)\right]G_{\xi ,\beta }\left(x-u\right)+F\left(u\right), X>u$$

(8)

The empirical estimator of the tail distribution is then:

$$\widehat{F}\left(x\right)=1-{\displaystyle \frac{N_u}{n}}{\left(1+\widehat{\xi }{\displaystyle \frac{\left(x-u\right)}{\widehat{\beta }}}\right)}^{-{\displaystyle \frac{1}{\widehat{\xi }}}}, X>u$$

(9)

where $n$ is the total number of observations, $N_u$ is the number of exceedances above the threshold $u$, and $\widehat{\xi }$ and $\widehat{\beta }$ are the maximum likelihood estimates of the shape and scale parameters, respectively.

To compute the unconditional VaR at a confidence level $q>F\left(u\right)$, Equation (9) can be inverted to obtain the following:

$${\mathrm{V}\mathrm{a}\mathrm{R}}_q=u+{\displaystyle \frac{\widehat{\beta }}{\widehat{\xi }}}\left[{\left({\displaystyle \frac{n}{N_u}}\left(1-q\right)\right)}^{-\widehat{\xi }}-1\right]$$

(10)

This formula provides the estimated VaR at quantile $q$, incorporating the tail behavior beyond the threshold $u$. The POT approach under EVT is particularly suited for modeling such tail risks in financial stress testing.

2.3. Dynamic R-Vine Copulas: Capturing Time-Varying Dependencies

This section introduces the dynamic Regular vine (R-vine) copula model to capture evolving dependence structures among financial assets. Unlike static copulas, the dynamic R-vine framework reflects time-varying interdependencies, which are essential for modeling volatility clustering and asymmetric co-movements in financial markets (Zhou and Ji, 2021). R-vine copulas, introduced by Bedford and Cooke (2001, 2002), decompose high-dimensional dependence into structured pairs of bivariate copulas, enabling flexibility and modularity. They outperform traditional copulas in representing conditional and asymmetric dependencies.

The Student-t copula is used due to its ability to capture symmetric tail dependence, critical for modeling joint extreme events like market crashes. Unlike the Gaussian copula, it incorporates both correlation ($\rho$) and degrees of freedom ($\nu$) parameters, enhancing its relevance in capital adequacy and risk management. In the dynamic copula framework, copula parameters evolve over time using a rolling window maximum likelihood estimation (MLE) method. The model recalibrates copula parameters over a moving window of observations, instead of specifying a time-series process. This approach allows each bivariate pair to exhibit unique, time-varying dependence patterns (Yang et al., 2021), enabling a granular and adaptive representation of joint dynamics.

In this study, the modeling process begins by estimating marginal distributions using ARMA–GJR–GARCH models with skewed Student-t innovations (see Section 2.1). Standardized residuals from each model are tested to ensure they are i.i.d., a prerequisite for copula modeling (Brechmann and Czado, 2013). These residuals are then transformed via the probability integral transform to produce uniform margins for the R-vine copulas. To capture time-varying dependence, copula parameters are subsequently estimated using a rolling window approach.

2.3.1. Sklar’s Theorem and the Dynamic R-Vine Framework

Sklar’s theorem (Sklar, 1959) underpins copula theory by decoupling marginal distributions from their joint dependence structure. For asset returns at time $t$, the joint distribution $H_t\left(r_{1,t}, r_{2,t},\dots ,r_{d,t} \right)$ is expressed as follows:

$$H_t\left(r_{1,t}, r_{2,t},\dots ,r_{d,t} \right)=C_t\left(F_{1,t}\left(r_{1,t}\right), \dots ,F_{d,t}\left(r_{d,t}\right);{\theta }_t\right)$$

(11)

where $F_{i,t}\left(r_{i,t}\right)$ represents the time-varying marginal distribution function of asset $i$ at time $t$, and ${\theta }_t$ denotes the copula parameter vector at time $t$, capturing the evolving dependence structure.

The returns $r_{i,t}$ are modeled using ARMA–GJR–GARCH processes with skewed Student-t innovations. The standardized residuals $z_{i,t}$ are transformed into uniform variables $u_{i,t}$ via $u_{i,t}=F_{i,t}\left(z_{i,t}\right)$.

These transformed values are used as inputs to the R-vine copulas, which models their conditional dependence through a series of bivariate copulas structured in hierarchical trees. Parameters ${\theta }_t$ are re-estimated over time using a 600-day rolling window, allowing the model to adjust to market dynamics.

2.3.2. R-Vine Copula Structure and Hierarchical Decomposition

The R-vine copula decomposes the joint copula $C_t$ into a product of conditional bivariate copulas, enabling flexible modeling of high-dimensional dependencies. The joint copula density $c_t\left(u_{1,t}, ..., u_{d,t}\right)$ is decomposed into conditional bivariate copulas as follows:

$$c_t\left(u_{1,t}, ..., u_{d,t}\right)=\prod _{k=1}^{d-1}\prod _{e\in E_k}{c}_{j\left(e\right),\left.l\left(e\right)\right|D\left(e\right)}(F\left(\left.u_{j\left(e\right),t}\right|u_{D\left(e\right),t}\right),F\left(\left.u_{l\left(e\right),t}\right|u_{D\left(e\right),t}\right);{\theta }_{j\left(e\right),\left.l\left(e\right)\right|D\left(e\right),t})$$

(12)

Each term $c_{j\left(e\right),\left.l\left(e\right)\right|D\left(e\right)}(\cdot )$ denotes a bivariate copula density for variables $j\left(e\right)$ and $l\left(e\right)$, conditional on the set $D$(e), where ${\theta }_{j\left(e\right),\left.l\left(e\right)\right|D\left(e\right),t}$ is the time-varying copula parameter at time $t$, and $E_k$ denotes the set of edges in tree $k$ of the vine structure, with each level representing conditional dependencies. The conditional distributions $F\left(\left.u_{j\left(e\right),t}\right|u_{D\left(e\right),t}\right)$ and $F\left(\left.u_{l\left(e\right),t}\right|u_{D\left(e\right),t}\right)$ are computed recursively.

The hierarchical structure is encoded in a sequence of trees $V_t=\left(T_{1,t},T_{2,t}, \dots , T_{d-1,t}\right)$, where each tree $T_{k,t}$ is constructed from edges in $T_{k-1,t}$, and all trees obey the proximity condition to ensure model validity. The proximity condition requires that, for each tree $T_{k,t}$, an edge can be formed only between nodes that share a common node in the previous tree $T_{k-1,t}$, thereby ensuring consistent conditional dependence modeling. This hierarchical structure enables tailored modeling of both pairwise and conditional dependencies, enhancing the flexibility and accuracy of multivariate dependence representation (Bedford and Cooke, 2001, 2002).

2.3.3. Copula Family

All bivariate copulas $c_{j\left(e\right),\left.l\left(e\right)\right|D\left(e\right)}(\cdot )$ used in the R-vine structure are specified as bivariate Student-t copulas due to their ability to capture symmetric tail dependence, which is crucial in modeling joint extreme market movements. The Student-t copula with correlation parameter $\rho$ and degrees of freedom $\nu$ is defined as follows:

$$C_{\rho ,\nu }^t\left(u_1,u_2\right)=t_{\rho ,\nu }\left(t_{\nu }^{-1}\left(u_1\right),t_{\nu }^{-1}\left(u_2\right)\right)$$

(13)

where $t_{\nu }^{-1}\left(\cdot \right)$ is the inverse cumulative distribution function (quantile function) of the univariate Student-t distribution with $\nu$ degrees of freedom, and $t_{\rho ,\nu }\left(\cdot ,\cdot \right)$ is the cumulative distribution function of the bivariate Student-t distribution with correlation $\rho$ and degrees of freedom $\nu$.

The corresponding copula density function is as follows:

$$c_{\rho ,\nu }^t\left(u_1,u_2\right)= {\displaystyle \frac{f_{\rho ,\nu }\left(t_{\nu }^{-1}\left(u_1\right),t_{\nu }^{-1}\left(u_2\right)\right)}{f_{\nu }\left(t_{\nu }^{-1}\left(u_1\right)\right)\cdot f_{\nu }\left(t_{\nu }^{-1}\left(u_2\right)\right)}}$$

(14)

where $f_{\rho ,\nu }\left(\cdot ,\cdot \right)$ is the bivariate Student-t density, and $f_{\nu }\left(\cdot \right)$ is the univariate Student-t density. These properties enable the Student-t copula to jointly model extreme outcomes in both tails of the distribution, making it particularly suitable for financial risk modeling, including VaR forecasts and capital adequacy analyses, as discussed in the empirical section (Demarta and McNeil, 2005).

2.3.4. Kendall’s Tau and Tail Dependence

To initialize and evaluate dependencies within the R-vine framework, Kendall’s tau is used for its robustness and its direct relationship with copula parameters. For a bivariate copula $C\left(u_1,u_2\right)$, Kendall’s tau $\left(\tau \right)$ defined as follows:

$$\tau =4{\int }_0^1{\int }_0^{1}C\left(u_1,u_2\right)·c\left(u_1,u_2\right)\mathrm{d}{u}_1{du}_2-1$$

(15)

where $c\left(u_1,u_2\right)={\displaystyle \frac{{\partial }^{2}C\left(u_1,u_2\right)}{\partial u_1\partial u_2}}$ is the copula density, assuming it exists.

Tail dependence coefficients further quantify the likelihood of simultaneous extreme outcomes. The upper and lower tail dependence coefficients, ${\lambda }_U$ and ${\lambda }_L$, are defined as follows:

$$\left\{\begin{array}{c}{\lambda }_U=\underset{t\to 1^{-}}{\mathit{lim}}{\displaystyle \frac{1-2t+C\left(t,t\right)}{1-t}}\\ {\lambda }_L= \underset{t\to 0^{+}}{\mathit{lim}}{\displaystyle \frac{C\left(t,t\right)}{t}} \end{array}\right.$$

(16)

The Student-t copula employed in this study exhibits symmetric, non-zero tail dependence, offering robustness in modeling extreme co-movements (Joe, 2014; Nelsen, 2006).

2.3.5. Estimation of Time-Varying Copula Parameters

To capture the evolving dependence structure among financial assets, time-varying copula parameters ${\theta }_t$ are estimated using a rolling window maximum likelihood estimation (MLE) procedure. At each time point $t$, the copula parameters are estimated over a moving window of size $W$ (e.g., 2,269 trading days), which advances forward by one day throughout the sample period. The local log-likelihood function at time $t$ is defined as follows:

$$L\left({\theta }_t\right)=\sum _{k=1}^{d-1}\sum _{e\in E_k}\sum _{s=t-W+1}^{t}log\left[c_{j\left(e\right),\left.l\left(e\right)\right|D\left(e\right)}(F\left(\left.u_{j\left(e\right),s}\right|u_{D\left(e\right),s}\right),F\left(\left.u_{l\left(e\right),s}\right|u_{D\left(e\right),s}\right);{\theta }_{j\left(e\right),\left.l\left(e\right)\right|D\left(e\right),s})\right]$$

(17)

where $c_{j\left(e\right),\left.l\left(e\right)\right|D\left(e\right)}\left(\cdot ,·;\theta \right)$ is the bivariate copula density for edge $e$ in tree $k$, parameterized by $\theta$; ${\theta }_{j\left(e\right),\left.k\left(e\right)\right|D\left(e\right),s}$ denotes the copula parameter vector at time $s$ for the conditional pair $\left(j\left(e\right),l\left(e\right)|D\left(e\right)\right)$; $E_k$ denotes the set of edges in tree $k$ of the vine structure; and $F\left(\left.u_{j\left(e\right),s}\right|u_{D\left(e\right),s}\right)$ and $F\left(\left.u_{l\left(e\right),s}\right|u_{D\left(e\right),s}\right)$ are the recursively computed conditional marginal distributions of the pseudo-observations (see Section 2.3.2).

This framework accommodates a wide range of bivariate copula families, including Gaussian, Clayton, Gumbel, Frank, and Student-t copulas. The copula family and associated parameters for each pair are selected based on model fit criteria such as the Akaike Information Criterion (AIC) or Bayesian Information Criterion (BIC).

If the Student-t copula is selected for a given bivariate pair, the log-likelihood function becomes the following:

$$L\left({\theta }_t\right)=\sum _{k=1}^{d-1}\sum _{e\in E_k}\sum _{s=t-W+1}^{t}log\left[c_t\left(F\left(\left.u_{j\left(e\right),s}\right|u_{D\left(e\right),s}\right),F\left(\left.u_{l\left(e\right),s}\right|u_{D\left(e\right),s}\right);{\rho }_{j\left(e\right),\left.l\left(e\right)\right|D\left(e\right),s},{\nu }_{\left.j\left(e\right),l\left(e\right)\right|D\left(e\right),s}\right)\right]$$

(18)

where ${\rho }_{j\left(e\right),\left.l\left(e\right)\right|D\left(e\right),s}$ is the time-varying correlation parameter; ${\nu }_{\left.j\left(e\right),l\left(e\right)\right|D\left(e\right),s}$ is the degrees of freedom; $c_t\left(\cdot ,\cdot ,\rho ,\nu \right)$ is the bivariate Student-t copula density function, defined in Section 2.3.3; and the conditional marginals $F\left(\left.u_{j\left(e\right),s}\right|u_{D\left(e\right),s}\right)$ and $F\left(\left.u_{l\left(e\right),s}\right|u_{D\left(e\right),s}\right)$ are computed recursively as before. This methodological framework and estimation strategy closely follow Aas et al. (2009) and Brechmann and Czado (2013).

2.3.6. Model Selection and Vine Structure

The vine structure is selected using the Dißmann et al. (2013) algorithm, which uses empirical Kendall’s tau and Akaike Information Criterion (AIC) to guide pair-copula selection and tree construction. While the flexible R-vine is the default, a D-vine may be chosen when a natural economic ordering exists. This data-driven approach ensures both statistical fit and economic interpretability.

2.4. Forecasting Method

This study adopts a rolling window approach with a fixed size of 2,269 daily observations to estimate ARMA–GJR–GARCH, ARMA–GJR–GARCH–EVT, and dynamic R-vine copula models. At each step, the window advances by one day to re-estimate parameters and produce one-step-ahead forecasts. The marginal models generate forecasts of conditional means and variances for all eight assets. Based on these forecasts, two types of portfolios comprising $d$ assets, where $d$ = 5 and 6 (with the six-asset portfolio including a green bond), are constructed from the full set of eight. Concurrently, the dynamic R-vine copula is estimated using standardized residuals, and S pseudo-random samples are drawn via inverse Rosenblatt sampling to simulate joint returns. These simulations yield time-varying forecasts of portfolio means and covariances, dynamically updating the risk and return structure over 600 out-of-sample observations.

2.5. Value-at-Risk Measures

This section introduces the Value-at-Risk (VaR) framework used to quantify potential portfolio losses under normal market conditions. VaR estimates are based on ARMA–GJR–GARCH and ARMA–GJR–GARCH–EVT forecasts, together with simulations from a dynamic R-vine copula, forming the basis for forward-looking risk and capital adequacy assessment.

2.5.1. ARMA-GJR-GARCH VaR (with skewed Student-t)

To quantify market risk, one-day-ahead VaR estimates are obtained from the conditional mean and variance forecasts produced by the ARMA-GJR-GARCH model. Let ${\mu }_{t+1}$ and ${\sigma }_{t+1}^2$ denote the conditional mean and variance of asset returns at time $t+1$. Assuming the standardized residuals follow a skewed Student-t distribution (Hansen, 1994), the VaR at the $(1-\alpha$) confidence level is computed as follows:

$${\mathrm{V}\mathrm{a}\mathrm{R}}_{t+1}^{\alpha }={\mu }_{t+1}+{\sigma }_{t+1}{q}_{\alpha , skewed-t}$$

(19)

where ${\mathrm{V}\mathrm{a}\mathrm{R}}_{t+1}^{\alpha }$ is the VaR for the next period $t+1$ at a confidence level of $\alpha$, ${\mu }_{t+1}$ is the conditional mean return forecast from the ARMA model for the next period $t+1$, ${\sigma }_{t+1}$ is the conditional standard deviation forecast from the GJR-GARCH model for the next period $t+1$, and $q_{\alpha , skewed-t}$ is the $\alpha$-quantile of the standardized skewed Student’s t-distribution with estimated skewness $\left(\eta \right)$ and degrees of freedom $\left(\nu \right)$ parameters. A 95% confidence level is applied for one-day VaR, in line with Thai regulations. For 10-day VaR, consistent with Basel’s 99% confidence level, a rolling window generates 10 one-step-ahead forecasts from the ARMA–GJR–GARCH model.

2.5.2. ARMA-GJR-GARCH-EVT VaR

To estimate one-day-ahead VaR under stressed conditions, this study combines ARMA–GJR–GARCH filtering with EVT. Returns are first modeled using ARMA–GJR–GARCH with skewed Student-t innovations, and the resulting standardized residuals are fitted using the POT method with a GPD to capture tail extremes. To assess capital requirements during market stress, the one-day VaR at time $t+1$ is estimated at a 97.5% confidence level. This tail-focused approach enhances risk quantification by better modeling the distribution of extreme losses. The one-day VaR at time $t+1$ derived from the conditional EVT model is given by the following:

$${\mathrm{V}\mathrm{a}\mathrm{R}}_{t+1}^{\alpha }={\mu }_{t+1}+ {\sigma }_{t+1}\left[u+{\displaystyle \frac{\widehat{\beta }}{\widehat{\xi }}}\left[{\left({\displaystyle \frac{n}{N_u}}\left(1-q\right)\right)}^{-\widehat{\xi }}-1\right]\right]$$

(20)

where ${\mu }_{t+1}$ and ${\sigma }_{t+1}$are the one-step-ahead forecasts of the conditional mean and standard deviation, respectively; $n$ is the total number of observations; $N_u$ is the number of observations exceeding the threshold $u$; and $\widehat{\xi }$ and $\widehat{\beta }$ are the estimated shape and scale parameters. This two-stage filtering and tail-modeling framework follows the extreme value approach of McNeil and Frey (2000).

2.5.3. Dynamic R-Vine Copula VaR

To estimate VaR under evolving market conditions, this study applies a dynamic R-vine copula to model time-varying dependencies among portfolios of five and six assets, each modeled using ARMA–GJR–GARCH with skewed Student-t innovations. The one-day portfolio VaR at a confidence level $\alpha$ is computed using a Monte Carlo simulation. At each time $t$, a set of $S$ pseudo-random return vectors, ${\left\{r_{t+1}^{\left(s\right)}\right\}}_{s=1}^S$, is generated from the dynamic R-vine copula, conditioned on the current marginal forecasts. Each $r_{t+1}^{\left(s\right)}$ represents the $s^{th}$ simulated return vector drawn from the copula model. The simulated portfolio returns are then computed as follows:

$$R_{t+1}^{\left(s\right)}={w^{T}r}_{t+1}^{\left(s\right)}$$

(21)

where $w$ is the vector of portfolio weights. The empirical $\alpha$-quantile of the simulated portfolio returns ${\left\{R_{t+1}^{\left(s\right)}\right\}}_{s=1}^S$ provides the one-day portfolio VaR estimate at time $t+1$:

$${\mathrm{V}\mathrm{a}\mathrm{R}}_{t+1}^{\alpha }=-{Quantile}_{\alpha }\left({\left\{R_{t+1}^{\left(s\right)}\right\}}_{s=1}^S\right)$$

(22)

where ${Quantile}_{\alpha }$ denotes the empirical quantile at level $\alpha$, capturing left-tail portfolio risk. This simulation-based dynamic copula framework builds on the methodology developed by Aas et al. (2009) and Brechmann and Czado (2013).

This study constructs four portfolios—two for life insurers and two for non-life insurers—each comprising five conventional assets or six assets including green bonds. All portfolios are designed with two key objectives: duration matching, to align asset durations with insurer liabilities, and Sharpe ratio maximization, to optimize risk-adjusted returns. The optimization respects Thai OIC investment limits and standard constraints $\sum _{i=1}^n{w}_i=1$ and $w_i\ge 0$, where $w_i$ is the weight of the ith asset.

Using the optimal weights, the one-day 95% VaR is estimated via a dynamic R-vine copula, which captures time-varying dependencies and non-Gaussian features. These estimates serve as proxies for capital requirements. Comparing five- and six-asset portfolios for both life and non-life insurers allows this study to assess whether green bond inclusion reduces capital needs and improves Sharpe ratios under evolving market conditions. This analysis is consistent with Markowitz portfolio theory (1952), adapted to the regulatory and liability constraints specific to the insurance sector.

2.6. Backtesting

Backtesting evaluates the accuracy of VaR models by comparing predicted losses with actual outcomes. This study applies two widely used methods: Kupiec’s unconditional coverage test and Christoffersen’s conditional coverage test, which assess model reliability from different perspectives (Ziggel et al., 2014).

2.6.1. Kupiec’s Unconditional Coverage Test

Kupiec’s test assesses whether the observed number of VaR exceedances aligns with the expected frequency implied by the model’s confidence level (Kupiec, 1995). Let $N$ be the number of exceedances over $T$ trading days, and let $\alpha$ be the probability of an exceedance, dictated by the confidence level of the VaR model. Under the null hypothesis $H_0$, $N$ follows a binomial distribution with parameters $\left(T,\alpha \right)$. The likelihood ratio test statistic is calculated as follows:

$${LR}_{UC}=-2\ln \left[{\displaystyle \frac{{\left(1-\alpha \right)}^{T-N}{\left(\alpha \right)}^N}{{\left(1-\frac{N}{T}\right)}^{T-N} {\left(\frac{N}{T}\right)}^N}} \right] \sim {\chi }_{1 }^2$$

(23)

If the calculated ${LR}_{UC}$ exceeds the critical value from the chi-squared distribution with 1 degree of freedom, the null hypothesis is rejected. This indicates that the model’s predicted exceedance frequency is inaccurate and may either underestimate or overestimate risk.

2.6.6. Christoffersen’s Conditional Coverage Test

While Kupiec’s test evaluates whether the overall frequency of VaR exceedances conforms to the model’s confidence level, it does not consider the timing or sequence of those exceedances. Christoffersen’s conditional coverage test addresses this limitation by assessing whether exceedances are independently distributed over time (Christoffersen, 1998). The null hypothesis states that the likelihood of an exceedance on day $t+1$ is unaffected by whether one occurred on day $t$. To test this, the number of transitions between exceedance states across consecutive days is recorded. Specifically, $n_{ij}$ represents the number of transitions from state $i$ on day $t$ to state $j$ on day $t+1$, where $i$ and $j$ can be either 0 (no exceedance) or 1 (exceedance). This structure enables the estimation of conditional probabilities to determine whether exceedances occur independently over time. Based on these transitions, the following conditional probabilities are defined:

$${\pi }_{01}=P\left(I_{t+1}=1|I_t=0\right)=\mathrm{P}\mathrm{r}\left(Exceedance tomorrow|No exceedance today\right)$$

$${\pi }_{11}=P\left(I_{t+1}=1|I_t=1\right)=\mathrm{P}\mathrm{r}\left(Exceedance tomorrow|Exceedance today\right)$$

Under the assumption of independence, ${\pi }_{01}$ and ${\pi }_{11}$ should be equal and approximately equal to $\alpha$, the VaR model’s stated exceedance probability.

Christoffersen’s test uses a likelihood ratio statistic to compare the likelihood of the observed data under the independence assumption with the likelihood without this constraint:

$${LR}_{CC}=2ln\left[{\left(1-{\pi }_{01}\right)}^{n_{00}}{\pi }_{01}^{n_{01}}{\left(1-{\pi }_{11}\right)}^{n_{10}}{\pi }_{11}^{n_{11}}\right]-2ln\left[{\left(1-\alpha \right)}^{T-N}{\alpha }^N\right] ~ {\chi }_2^2$$

(24)

where $N$ is the total number of exceedances in a sample of $T$ observations. If the ${LR}_{CC}$ statistic exceeds the critical value from the chi-squared distribution with two degrees of freedom, the null of independence is rejected. This indicates clustered exceedances and potential shortcomings in the VaR model.

3. Empirical Results and Discussions

This section details the empirical findings of the study, beginning with the estimation results of the ARMA-GJR-GARCH models incorporating a skewed Student-t distribution for each of the eight assets: the SET index, the Dubai crude oil, the Thai bullion gold, the 3–7-year government bond index, the 7–10-year government bond index, the JPY/THB exchange rate, the property sector index, and the Bloomberg Barclays MSCI US green bond index. Subsequently, the analysis transitions to the ARMA-GJR-GARCH-EVT model, which leverages extreme value theory to model tail risk. Finally, the dynamic R-vine copula model is employed to evaluate dynamic correlations among the assets and their subsequent impact on portfolio VaR. The results for each model are discussed in detail, highlighting key findings and their implications for risk management within the Thai insurance industry.

3.1. ARMA-GJR-GARCH Estimation (with Skewed Student-t Distribution)

Table 2 presents the in-sample parameter estimates derived from the selected ARMA-GJR-GARCH models, employing a skewed Student-t distribution, for each asset. Model selection, guided by the AIC, resulted in distinct ARMA orders for different assets. The results show that ARMA(3,3)-GJR-GARCH(1,1) for the SET, ARMA(3,2)-GJR-GARCH(1,1) for the crude oil, ARMA(1,0)-GJR-GARCH(1,1) for the gold, ARMA(3,2)-GJR-GARCH(1,1) for the 3–7-year government bond, ARMA(1,0)-GJR-GARCH(1,1) for the 7–10-year government bond, ARMA(1,1)-GJR-GARCH(1,1) for the exchange rate, ARMA(5,3)-GJR-GARCH(1,1) for the property, and ARMA(1,0)-GJR-GARCH(1,1) for the green bond. The diverse ARMA orders suggest varying degrees of short-term dependencies in the return series of the assets.

Note: *** significant at 0.01, ** significant at 0.05, * significant at 0.1.

A consistent observation across all assets is the presence of significant volatility clustering, evidenced by the near-unity sum of the ARCH and GARCH coefficients in the selected models. This reinforces the necessity of employing time-varying volatility models for effective risk management. Furthermore, the statistically significant and positive coefficients of lagged squared returns confirm strong ARCH effects, and the significant coefficients on the lagged conditional variance confirm strong GARCH effects, indicating that both past shock and past volatility information are critical for forecasting future volatility. These findings align with previous research that has documented persistent volatility in both developed and emerging markets (Floros, 2007; Lin et al., 2020). This study builds upon these findings by demonstrating the applicability of GARCH models, specifically the GJR-GARCH, for capturing volatility dynamics relevant to investment risk management for insurers within the Thai insurance industry.

The GJR-GARCH specification captures asymmetries in volatility responses to shocks via the gamma coefficient. All assets displayed statistically significant gamma estimates, though the sign and magnitude varied. A positive gamma, as seen in the SET, the crude oil, and the property, implies an inverse leverage effect, where positive shocks increase volatility more than negative ones. Conversely, a negative gamma, observed in the gold, the 3–7-year government bond, the 7–10 government bond, the exchange rate, and the green bond, aligns with the traditional leverage effect, where negative shocks have a stronger impact on volatility. Moreover, the statistically significant skewness and shape parameters support the use of the skewed Student-t distribution, which accommodates both asymmetry and excess kurtosis in the standardized residuals. This specification improves model fit and yields more realistic volatility and risk estimates (Nugroho et al., 2021).

Diagnostic tests, including the Ljung–Box test on standardized squared residuals and the ARCH–LM test, generally support the adequacy of the selected ARMA–GJR–GARCH models for capturing volatility dynamics, as most assets exhibit no significant autocorrelation or residual ARCH effects (p-values > 0.05). An exception is the property, where the Ljung–Box test indicates mild autocorrelation at higher lags, and the ARCH–LM test reveals significant residual ARCH effects. While such autocorrelation may be tolerable, the persistence of this heteroskedasticity suggests that additional modeling is needed to fully capture extreme volatility behavior.

Although the ARMA–GJR–GARCH model effectively captures time-varying volatility and asymmetric responses to shocks, it may fall short in modeling extreme tail events in this sector. To address this limitation, the next section introduces an EVT framework, which complements the GARCH model by explicitly focusing on the distribution of extreme returns. This two-stage approach enhances the ability to assess tail-related risks and improves the reliability of Value-at-Risk (VaR) estimates for risk management in the Thai insurance sector. By incorporating EVT, the model is better equipped to capture recent and extreme fluctuations in returns, thereby producing more accurate VaR forecasts across a range of confidence levels.

3.2. ARMA-GJR-GARCH-EVT Estimation for Tail Risk Assessment

Accurate tail risk modeling is crucial for capital adequacy, especially during market downturns. To enhance tail risk estimation, this study integrates the ARMA–GJR–GARCH model with EVT, focusing on the left tail due to observed negative skewness and potential downside risk. The framework incorporates the skewed Student-t distribution to account for key stylized facts in financial time series such as asymmetry and leptokurtosis, inadequately captured by models assuming normality (Huang et al., 2014). By modeling asymmetric volatility via the GJR-GARCH structure and capturing skewness and excess kurtosis through skewed Student-t innovations, the approach offers a more nuanced and realistic representation of volatility and tail behavior.

After estimating the conditional mean and variance using the ARMA–GJR–GARCH specification, the POT approach from EVT is applied to model extreme losses in the standardized residuals. Selecting an appropriate threshold is a critical step: if set too low, non-extreme observations may be included, distorting tail estimates; if too high, the number of exceedances may be too small, increasing estimation variance and reducing reliability (Coles, 2001). In line with the recent literature (Eita and Djemo, 2022; Huang et al., 2017; Li, 2017; McNeil and Frey, 2000), this study adopts the 93rd-percentile as the threshold level. This higher threshold is particularly suitable for stress-testing in regulatory and solvency contexts that demand a focus on rare but impactful losses.

Across all asset classes, the selected threshold yields an average of approximately 159 exceedances, providing a stable sample size for tail modeling. Despite differences in absolute threshold values, the consistent exceedance count supports the robustness of the chosen threshold. As shown in Table 3, the estimated shape parameter ($\xi$) is close to zero for most asset classes, suggesting that their tails are well-approximated by an exponential distribution. However, the property exhibits distinctly heavier tails, with a positive and statistically significant shape parameter; this implies a Pareto-type distribution and greater exposure to extreme downside risk. The Kolmogorov–Smirnov goodness-of-fit test supports these findings, with p-values above 0.05 for the exponential fit in most cases, while a heavier-tailed model provides a better fit for the property. These results confirm the effectiveness of the ARMA–GJR–GARCH–EVT framework in capturing heterogeneous tail risks across asset classes, which is vital for robust investment risk assessment.

This is particularly important in light of Section 3.1, where the property exhibited residual ARCH effects even after modeling with ARMA–GJR–GARCH and skewed Student-t distribution. The incorporation of EVT directly addresses these residual tail risks by explicitly modeling extreme losses that conventional GARCH structures may not fully capture. The ARMA–GJR–GARCH–EVT framework is thus employed specifically to estimate capital requirements under stress scenarios, focusing on the extreme downside risks faced by insurers.

Notably, the findings underscore the inadequacy of assuming normally distributed returns, as such assumptions can substantially underestimate capital requirements for Thai insurers during market stress, especially given the fat tails and skewness typical of financial return distributions. The GJR-GARCH component captures the asymmetric volatility response to negative shocks, a well-documented phenomenon in financial markets, while the ARMA component accounts for autocorrelation in return series. This two-stage modeling framework improves the accuracy of VaR estimates and facilitates more prudent and risk-sensitive capital allocation.

In the next section, the analysis is extended by incorporating dynamic dependence structures through R-vine copulas, enabling realistic modeling of time-varying and nonlinear relationships across asset classes in a diversified investment environment.

3.3. ARMA-GJR-GARCH-EVT Estimation for Tail Risk Assessment

This section presents empirical results from R-vine copula models for life and non-life insurer portfolios. Building on the ARMA–GJR–GARCH marginal models discussed in Section 3.1, the analysis captures nonlinear and asymmetric dependencies among financial assets and evaluates the impact of green bond inclusion on diversification, capital adequacy, and risk-adjusted returns. R-vine copulas with Student-t pairings effectively model tail dependencies, critical for portfolio risk management across diverse financial assets in extreme markets.

3.3.1. Portfolio Analysis for Life Insurers

This subsection presents R-vine copula results for life insurer portfolios, comparing a five-asset baseline with a six-asset version that includes green bonds. Both portfolios are optimized for the Sharpe ratio, meet capital requirements, and align with life insurers’ long-term liabilities. The analysis evaluates whether green bond inclusion enhances diversification, reduces tail risk, and improves overall portfolio performance.

The results of the copula estimations are summarized in Table 4. Both configurations use D-vine Student-t copulas, crucial for managing extreme market risks by capturing nonlinear, asymmetric, and tail-dependent relationships. The upper panel of Table 4 presents the five-asset portfolio, comprising the SET index (1), the Thai bullion gold (2), the property sector index (3), the 7–10-year government bond index (4), and the JPY/THB exchange rate (5). The strongest unconditional dependence is observed between the SET and the property ($\tau$ = 0.87), reflecting equity market concentration. A strong positive dependence between the gold and the exchange rate ($\tau$ = 0.80, ${\lambda }_1$ = ${\lambda }_2$ = 0.75) reflects safe-haven behavior, while a notable negative dependence between the SET and the exchange rate ($\tau$ = –0.72) suggests diversification potential. As the R-vine structure progresses, tail dependence and conditional relationships weaken, indicating that extreme co-movements are concentrated among a few key asset pairs. The estimated degrees of freedom ($\nu$) range from 3.38 to 30, reflecting varying tail heaviness across the dependence structure. Overall, the five-asset model exhibits strong statistical performance, with the log-likelihood at 7,733.54, the AIC at –15,427.08, and the BIC at –15,312.54, forming a robust baseline for evaluating the added value of green bond inclusion in the six-asset configuration.

The lower panel of Table 4 presents the six-asset portfolio, which includes the Bloomberg Barclays MSCI US green bBond index (6). The green bond significantly alters the dependence structure. In Tree 1, the green bond exhibits strong unconditional positive dependence with the gold ($\tau$ = 0.75) and the exchange rate ($\tau$ = 0.74), along with substantial symmetric tail dependence (${\lambda }_1$ = ${\lambda }_2$ = 0.62 and 0.53, respectively), positioning them alongside traditional safe-haven assets and reinforcing their stabilizing role during market stress. In higher-order trees, the green bond displays moderate to weak conditional dependencies, such as a modest positive link with the SET ($\tau$ = 0.34), a weak negative association with the 7–10-year government bond ($\tau$ = –0.21), and a slightly stronger negative relationship with the property ($\tau$ = –0.25). These relationships suggest that the green bond provides differentiated exposure to macroeconomic risks and may help mediate cross-asset risk transmission. By Tree 5, dependencies are negligible, highlighting their role in diffusing residual co-movement. Statistically, the six-asset model outperforms the five-asset configuration, with the log-likelihood of 9,545.86, the AIC of –19,031.73, and the BIC of –18,859.92. These improvements indicate that green bonds enrich the dependence structure, enhance diversification, and strengthen tail-risk modeling. This clearly supports capital-efficient and resilient portfolio design aligned with sustainability-oriented investment objectives.

Table 5 presents the optimal portfolio weights for life insurers under Sharpe ratio maximization, comparing a traditional five-asset configuration with a six-asset portfolio that includes the green bond. In the five-asset portfolio, the 7–10-year government bond receives the highest allocation (50%) due to its low volatility and strong duration-matching benefits. The SET accounts for 30%, while the gold and the property are each allocated 5%. The exchange rate receives a 10% weight, offering diversification and safe-haven characteristics. The introduction of the green bond in the six-asset portfolio reshapes the asset mix while preserving portfolio stability. The 7–10-year government bond remains at 50%, while the SET is reduced to 20%. The green bond and the property are each allocated 10%, with the gold and the exchange rate maintaining their 5% weights. This reallocation underscores green bonds’ contribution to diversification and resilience, while supporting ESG-aligned objectives.

These optimized in-sample weights provide a baseline for the subsequent dynamic analysis. The next step evaluates time-varying portfolio performance using 600-day rolling-window forecasts of VaR and the Sharpe ratio for life insurers.

3.3.2. Portfolio Analysis for Non-Life Insurers

This subsection presents R-vine copula results for non-life insurer portfolios, built using the same Sharpe ratio optimization and duration-matching approach as applied to life insurers. Due to their shorter-duration liabilities, non-life insurers require greater sensitivity to liquidity and volatility. While return optimization remains relevant, capital preservation and regulatory compliance take precedence.

Two portfolios are analyzed: a five-asset baseline and a six-asset version including green bonds. The five-asset portfolio uses a D-vine structure, while the six-asset R-vine captures increased interdependencies with green bond inclusion. Both utilize Student-t copulas to model nonlinear, asymmetric, tail-dependent relationships, crucial under market stress. The upper panel of Table 6 presents the five-asset portfolio, comprising the SET index (1), the Dubai crude oil (2), the property sector index (3), the 3–7-year government bond Index (4), and the JPY/THB exchange rate (5). The D-vine structure reveals the strongest unconditional dependence between the SET and the property ($\tau$ = 0.87, ${\lambda }_1$ = ${\lambda }_2$ = 0.64), reflecting high co-movement within the equity market. Notable negative dependencies are observed between the exchange rate and both the SET ($\tau$ = –0.72) and the crude oil ($\tau$ = –0.83), highlighting the diversification potential of currency exposure. A moderate negative relationship between the crude oil and the 3–7-year government bond ($\tau$ = –0.61) suggests opposing responses to macroeconomic shocks such as inflation or interest rate changes. As the vine progresses, dependencies weaken. For instance, Tree 3 reports a weak positive link between the crude oil and the property ($\tau$ = 0.21) and a moderately negative relationship between the SET and the 3–7-year government bond ($\tau$ = –0.60), indicating differentiated behavior under stress. Tail dependence is concentrated in the first tree. The model demonstrates strong statistical performance, with the log-likelihood of 8,037.66, the AIC of –16,035.31, and the BIC of –15,920.77.

The lower panel of Table 6 presents the six-asset portfolio, incorporating the Bloomberg Barclays MSCI US green bond index (6). With increased interdependencies, a general R-vine structure is adopted. In Tree 1, the green bond exhibits strong positive dependence with the exchange rate ($\tau$ = 0.74, ${\lambda }_1$ = ${\lambda }_2$ = 0.54) and a strong negative dependence with the 3–7-year government bond ($\tau$ = –0.74), indicating their potential as a hedge against interest-rate-sensitive assets. In higher-order trees, the green bond shows moderate conditional dependence with the SET ($\tau$ = 0.42) and weak negative relationships with the property ($\tau$ = –0.17) and the crude oil ($\tau$ = –0.18), suggesting differentiated macroeconomic exposure. The six-asset model improves the statistical fit, with the log-likelihood of 9,797.16, the AIC of –19,534.31, and the BIC of –19,362.50. These results confirm that green bonds enhance the dependence structure, improve diversification, and contribute to tail-risk mitigation. This supports the construction of capital-efficient and resilient portfolios for non-life insurers while aligning with sustainable investment strategies (Abakah et al., 2022).

Table 7 presents the optimal portfolio weights for non-life insurers under Sharpe ratio maximization, comparing a traditional five-asset configuration with a six-asset portfolio that includes green bonds. In the five-asset portfolio, the 3–7-year government bond receives the highest allocation (55%), reflecting their role as low-risk, liquid instruments aligned with non-life insurers’ short-duration liabilities. The exchange rate accounts for 30%, offering diversification and safe-haven benefits under volatile conditions. The remaining 15% is equally distributed across the SET, the crude oil, and the property, reflecting measured exposure to higher-volatility assets. In the six-asset portfolio, green bond inclusion reshapes the allocation while preserving the conservative structure. The 3–7-year government bond is reduced to 40%, while the green bond is allocated 10%, highlighting their role in enhancing stability and supporting ESG mandates. The exchange rate remains at 30%, the SET increases to 10%, and the crude oil and the property retain 5% each. This reallocation reflects a strategic response to the risk-return dynamics introduced by green bonds and supports capital efficiency under regulation-sensitive investment policies.

Consistent with the life insurer analysis, these in-sample weights are used as inputs for a 600-day rolling-window estimation. This dynamic approach captures the evolving risk–return profile of non-life insurer portfolios through forecasts of VaR and the Sharpe ratio.

3.4. Value-at-Risk Backtesting

This section evaluates the predictive performance of the VaR forecasts using two established backtesting procedures: the Kupiec unconditional coverage test and the Christoffersen conditional coverage test. Backtesting serves as a critical validation step, assessing whether the VaR estimates adequately capture potential losses and comply with solvency capital requirements (Smolović et al., 2017). The analysis employs a rolling-window estimation approach, as described in Section 2.4, to generate 600 daily out-of-sample VaR forecasts for each model under review.

Three modeling frameworks are evaluated across different confidence levels and horizons: (1) ARMA–GJR–GARCH with skewed Student-t innovations: Captures volatility clustering, asymmetric responses, and fat-tailed return distributions. It provides one-day-ahead VaR forecasts at the 95% confidence level, consistent with the requirement of the Thai OIC, and ten-day-ahead forecasts at the 99% level, aligned with Basel standards. These estimates reflect the capital buffers that life and non-life insurers are required to hold at specified confidence levels under financial solvency frameworks. (2) ARMA–GJR–GARCH–EVT hybrid model: Combines conditional volatility modeling with EVT to improve tail risk estimation, especially under extreme market conditions. It generates one-day-ahead VaR forecasts at the 97.5% confidence level, a threshold commonly used in insurance sector stress testing (Paraschiv et al., 2020), offering more conservative and robust estimates for capital adequacy. (3) Dynamic R-vine copula model: Captures time-varying, nonlinear dependencies among asset returns. Dynamically estimated copula parameters are used to forecast one-day-ahead VaR at the 95% confidence level and to compute Sharpe ratios based on optimized weights for life and non-life insurer portfolios.

3.4.1. ARMA–GJR–GARCH Model (with Skewed Student-t Innovations) and ARMA–GJR–GARCH–EVT Performance

This subsection presents backtesting results for two VaR forecasting models: the ARMA–GJR–GARCH with skewed Student-t innovations and its tail-risk-enhanced extension incorporating EVT. Both models are estimated using a 600-day rolling window and evaluated using the Kupiec unconditional coverage and Christoffersen conditional coverage tests to ensure forecast reliability for regulatory compliance and capital adequacy.

The ARMA–GJR–GARCH model generates one-day-ahead VaR forecasts at the 95% confidence level and ten-day-ahead forecasts at the 99% level. It effectively captures asymmetric volatility and fat tails, producing accurate forecasts for most assets. As shown in Table 8 and Figure 1 and Figure 2, most assets pass both backtests. Minor exceptions include the 7–10-year government bond, with p-values below 0.1, suggesting mild underestimation of risk. For the ten-day-ahead horizon, Kupiec test violations are observed for the crude oil and the 3–7-year government bond, while the Christoffersen test flags the crude oil, the exchange rate, and the green bond, indicating areas for improvement.

To enhance tail-risk estimation, the ARMA–GJR–GARCH–EVT model forecasts one-day-ahead VaR at the 97.5% confidence level. Table 9 and Figure 3 present results that demonstrate strong model performance, with exceedance rates generally aligning with expectations. Notable exceptions include the exchange rate (Kupiec p-value of 0.087) and the 3–7-year government bond (Christoffersen p-value of 0.037), indicating some inconsistency in the timing of exceedances.

The green bond exhibits one-day-ahead VaR values of −0.65% under the ARMA–GJR–GARCH model and −0.84% under the EVT-enhanced specification. These values exceed those of the 3–7-year (−0.14%, −0.21%) and the 7–10-year (−0.31%, −0.43%) Thai sovereign bonds but remain lower than the SET (−1.16%, −1.56%). Compared to other asset classes, the green bond also demonstrates lower risk than the crude oil (−3.38%, −4.46%), the gold (−1.34%, −1.64%), the exchange rate (−0.94%, −1.09%), and the property (−1.49%, −1.91%). Over a ten-day-ahead horizon at the 99% confidence level, the green bond shows a VaR of −2.86%, again falling between the 3–7-year (−0.58%) and the 7–10-year (−1.34%) government bonds, and riskier assets such as the SET (−5.19%), the gold (−5.92%), and the crude oil (−14.93%).

Despite some discrepancies, both models provide valuable insights. The ARMA–GJR–GARCH model is well-suited for routine risk monitoring, while its EVT-enhanced extension improves reliability under stressed market conditions. One-day-ahead VaR estimates from the GARCH model range from 0.13% to 3.37%, and ten-day-ahead estimates range from 0.57% to 14.93%. The SET’s average ten-day-ahead VaR is 5.19%, exceeding the typical range in developed markets (Degiannakis et al., 2014), supporting the case for localized capital adequacy standards. The EVT model yields higher one-day-ahead VaR estimates—from 0.21% to 4.46%—and shows the SET at 1.56%, below the 2.60%–3.99% range observed in developed economies (Echaust and Just, 2020), indicating differences in volatility or capital regulation.

The relatively higher VaR of the green bond in this study, compared to conventional Thai government bonds, indicates that it may not appear to be a lower-risk asset when considered in isolation. This contrasts with findings from developed markets, such as Liaw (2020), which describe a “greenium” effect, characterized by lower yields and reduced perceived risk driven by sustained ESG demand. However, such results are typically reported within the same market, credit tier, or issuer type. In contrast, this study compares an international green bond with Thai sovereign bonds, representing a cross-jurisdictional and cross-credit-category assessment. The elevated VaR observed for green bonds may therefore reflect differences in market structure, issuer composition, duration exposure, and currency denomination, rather than a true absence of the greenium effect.

While the green bond may not appear lower risk than Thai government bonds on a standalone basis, it exhibits lower capital requirements than other asset classes discussed earlier. Its predominantly negative dependence with other assets, as revealed by the R vine copula analysis, suggests potential diversification benefits. The next subsection examines this potential by evaluating portfolio-level risk and return through rolling window Sharpe ratios and VaR forecasts.

3.4.2. Performance Evaluation of the Dynamic R-Vine Copula Model: Rolling-Window Sharpe Ratios and VaR Forecasts

This section assesses the dynamic R-vine copula model’s performance in forecasting one-day-ahead VaR at the 95% confidence level and calculating Sharpe ratios using a 600-day rolling window. The model combines ARMA–GJR–GARCH marginal models with skewed Student-t innovations and a dynamic copula to capture time-varying, asymmetric dependencies among assets. This approach is valuable for insurers managing portfolio risk in volatile markets while aligning investments with solvency and sustainability goals. Unlike univariate models, R-vine copulas capture evolving nonlinear interdependencies, making them suitable for life and non-life insurers with distinct asset-liability management and risk-adjusted return objectives.

Table 10 presents rolling performance results for life and non-life portfolios with and without green bonds. For life insurers, green bond inclusion significantly enhances portfolio efficiency. The average Sharpe ratio rises from –0.0432 to 0.0063, while the average one-day-ahead VaR decreases from –0.1576% to –0.1428%, indicating improved risk-adjusted returns and reduced capital requirements. The six-asset portfolio achieves a higher Sharpe ratio in 80.50% of rolling windows and delivers a lower VaR in 60.33%, demonstrating meaningful diversification benefits. These findings align with Casal et al. (2025) and Gupta et al. (2025), supporting the role of green bonds as dual contributors to portfolio performance and solvency for life insurers. For non-life insurers, the six-asset portfolio delivers a moderate improvement in the average Sharpe ratio, increasing from −0.0928 to −0.0441. Additionally, 60.67% of rolling windows exhibit higher Sharpe ratios, underscoring the presence of diversification benefits. This improvement is accompanied by a slight rise in downside risk, as the average one-day-ahead VaR increases from −0.1004% to −0.1264%. For insurers that prioritize short-term liquidity and volatility control, the trade-off may still be justifiable, particularly for those with a higher risk appetite who are seeking enhanced returns. These findings align with Jareño et al. (2024), reinforcing that ESG-focused assets can enhance portfolio dynamics even when risk reduction is not uniform.

In summary, the dynamic R-vine copula model demonstrates that green bonds fundamentally alter inter-asset relationships, boosting portfolio efficiency, especially for life insurers, despite their higher standalone VaR. This restructuring reduces tail risk, improves Sharpe ratios for both life and non-life insurers, and lowers portfolio-level VaR, notably for life insurers. These findings align with Mensi et al. (2022), who highlight the effectiveness of copula-based models in capturing dynamic co-movements under market volatility. The results also support Tsoukala and Tsiotas’s (2021) emphasis on integrating ESG assets into complex financial systems. These benefits are more pronounced for life insurers, whose long-term investment horizons and solvency-driven mandates make them especially responsive to diversification gains. While Bouri et al. (2023) and Pham and Nguyen (2022) report green bond vulnerabilities to oil price volatility, policy uncertainty, and limited hedging capacity, the current findings highlight the importance of evaluating ESG assets through a dependence-aware, multivariate framework. Although the extent of benefit varies by insurer type, green bond integration supports both portfolio stability and ESG-aligned performance when assessed using advanced models such as dynamic copulas, which are able to capture nonlinear dependencies. These results, which clarify the nuanced role of green bonds in insurance portfolio construction, are further synthesized in the concluding section.

4. Conclusions

This study presents a robust framework for assessing capital requirements and optimizing insurer portfolios, particularly those that include green bonds, through advanced risk modeling. By integrating ARMA–GJR–GARCH models with skewed Student-t innovations, extreme value theory, and dynamic R-vine copulas, the methodology captures essential risk features such as volatility clustering, tail risk, and time-varying dependencies. This approach is well-suited for regulatory stress testing and solvency assessment, offering a comprehensive perspective under both normal and adverse market conditions. These advanced models can support regulators in establishing appropriate capital requirements by providing more accurate risk evaluations. Backtesting results further validate the models’ effectiveness and support their application in proactive risk oversight and regulatory supervision.

The analysis reveals that although green bonds may not exhibit a “greenium” effect in isolation compared to Thai sovereign bonds, the risk level associated with green bonds is lower than that of more volatile asset classes such as equities, oil, property, gold, and exchange rates. Incorporating green bonds into diversified insurer portfolios enhances risk-adjusted returns by increasing Sharpe ratios and reducing portfolio-level value-at-risk, in most rolling periods. This can potentially lead to reduced capital charges. These findings underscore the importance of evaluating capital requirements at the portfolio level, rather than focusing solely on individual assets. These benefits are especially evident among life insurers, whose long-term horizons and solvency mandates align with the characteristics of green bonds. While non-life insurers also experience diversification, short-term risk requires careful consideration. As ESG-aligned instruments, green bonds offer insurers a pathway to strengthen financial resilience while supporting national sustainability priorities.

Therefore, financial regulators in Thailand should proactively integrate green bonds into supervisory investment frameworks for insurer portfolios, actively encouraging their adoption to enhance both financial stability and sustainability. Given the observed improvements in risk-adjusted returns and potential reductions in capital requirements associated with green bond allocations, regulators should explicitly recognize the benefits of ESG assets within solvency-oriented investment policies. This should include conducting capital assessments at the portfolio level, utilizing advanced models like dynamic copula-based approaches that account for nonlinear and time-varying dependencies to accurately reflect the risk-reducing properties of green bonds. By strengthening regulatory guidance in this way, regulators can not only enhance insurers’ financial resilience and potentially reduce capital charges under the Office of Insurance Commission of Thailand regulatory frameworks, but also actively steer the insurance sector towards supporting Thailand’s sustainable finance agenda and commitments under the United Nations Sustainable Development Goals (SDGs).

Future research could further enhance forecasting accuracy by exploring alternative GARCH family models, such as EGARCH, TGARCH, or GJR GARCH M, and regime switching frameworks to better capture market asymmetries and structural shifts. Integrating EVT copulas models may improve tail dependence modeling under extreme conditions, while advanced distributions such as generalized hyperbolic or finite mixtures could more precisely represent skewness and fat tails. These extensions would reinforce insurer risk management and promote adaptive, sustainability-aligned investment strategies.