Complexity
Volume 2020 (2020), Article ID 6746303, 15 pages
https://doi.org/10.1155/2020/6746303

Research Article

Bayesian Estimation of Archimedean Copula-Based SUR Quantile Models

Nachatchapong Kaewsompong, Paravee Maneejuk, and Woraphon Yamaka

Center of Excellence in Econometrics, Faculty of Economics, Chiang Mai University, Chiang Mai 50200, Thailand

Correspondence should be addressed to Woraphon Yamaka; woraphon.econ@gmail.com

Received 3 March 2020; Revised 31 May 2020; Accepted 25 June 2020; Published 16 July 2020

Academic Editor: Zhile Yang

Copyright © 2020 Nachatchapong Kaewsompong et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

We propose a high-dimensional copula to model the dependence structure of the seemingly unrelated quantile regression. As the conventional model faces with the strong assumption of the multivariate normal distribution and the linear dependence structure, thus, we apply the multivariate exchangeable copula function to relax this assumption. As there are many parameters to be estimated, we consider the Bayesian Markov chain Monte Carlo approach to estimate the parameter interests in the model. Four simulation studies are conducted to assess the performance of our proposed model and Bayesian estimation. Satisfactory results from simulation studies are obtained suggesting the good performance and reliability of the Bayesian method used in our proposed model. The real data analysis is also provided, and the empirical comparison indicates our proposed model outperforms the conventional models in all considered quantile levels.

1. Introduction

Typically, the seemingly unrelated regression (SUR) model, developed by Zellner [1], is a system of many structural equations, each having its own dependent variable and probably different sets of independent variables. However, these equations may be correlated through their error terms assuming the existence of the multivariate normal distribution. This model has been applied in many works such as White and Hewings [2]; Adelegan [3]; and Frankel and Poonawala [4].

Recently, the performance of this model has been questioned by many scholars, such as Jun and Pinkse [5] and Waldmann and Kneib [6], as the regression lines are fitted through the means of the independent variables leaving aside the outliers which might carry meaningful information. In general, the data may exhibit a heavy tail. Hence, the normal distribution assumption may fail to provide reliable estimated results. To deal with these issues, Koenker and Bassett [7] introduced the quantile regression. Such an extended linear model has several advantages. For example, the parameter estimates are more robust against outliers, and the extensive analysis of the relationship between the dependent and independent variables is uncovered [6]. Recognizing these advantages of the quantile regression approach, some studies made use of the seemingly unrelated quantile regression (SUQR) model (see [5, 6]). They extended the conventional SUR model by allowing the error term of each quantile equation to be correlated at all quantiles but still with the same normal distribution assumption regarding the marginals. However, the assumption of the multivariate normal distribution may not be appropriate for joining the equation errors which actually have an asymmetric distribution. To relax this assumption, the copula functions are employed to capture the nonlinear dependence structure as well as join the asymmetric error of equations in the multivariate SUQR model.

In the literature, the copula-based models have already been introduced by many studies of which findings demonstrated a higher accuracy in parameter estimation. The pioneering work of Wichitaksorn and Choy [8] suggested using the copula function to join the error terms of the linear regression equation and binary choice equation. They showed that their proposed model is superior to the conventional model in terms of estimation accuracy. Pastpipatkul et al. [9] used the copula to join the errors of the SUR model. Louzada and Ferreira [10] also suggested using the Clayton copula to join the error terms of the bivariate seemingly unrelated regression (SUN) tobit model. They mentioned that their model has the ability to capture the lower tail dependence of the SUN model. Ivanov et al. [11] extended a copula approach to model the dependence of unobserved multivariate factors in the dynamic factor models. They confirmed that the copula-based approach is more general and applicable to several factor models. Recently, Zou et al. [12] suggested using the copula to establish the link between the occurrence of wildlife-vehicle collisions and the underreporting probability. They revealed that the Gaussian copula-based empirical Bayes method is superior to the traditional EB method.

Among various studies using the copula-based models to successfully confirm the role of copula functions in improving the efficiency of the system equation, none of them used the copulas to join the errors of the SUQR model except for Tansuchat et al. [13] who used multivariate elliptical copulas consisting of Gaussian and Student-t to construct the joint equations of the SUQR model. This paper, therefore, attempts to extend the approach of Tansuchat et al. [13] by applying the multivariate Archimedean copulas to join the model errors. Using this copula class has several additional advantages over the elliptical copulas. First, this copula class can model the asymmetric structure of extreme dependence between the errors, which are particularly important in financial modeling. Second, it provides more flexibility to model the dependence structure of the errors as there are many copula families introduced in the class.

Consequently, in this study, we apply various multivariate exchangeable Archimedean copulas to join the equations in the SUQR model. We believe that our model will become more flexible and applicable to investigate the entire conditional distribution of the dependent variable and becomes more robust against outliers. In the estimation aspect, as our proposed model contains large parameter estimates and the full likelihood function is quite complicated, thus, the Bayesian estimation is employed in this study. To confirm the accuracy and reliability of our model and estimation, we conduct the simulation study and real data analysis to evaluate the performance of our model and the Bayesian estimation. To the best of our knowledge, no research introduced the multivariate exchangeable Archimedean copulas to join the error of equations in the SUQR model. To this end, we suggest a Bayesian approach to the multivariate exchangeable Archimedean copula-based SUQR.

The outline of the remaining sections is as follows. In Section 2, we explain the copula-based SUQR model and multivariate copula functions. The posterior of the Bayesian is provided in Section 3. Four simulation studies are provided in Section 4. In Section 5, the real data example is used to show the performance of our model. Finally, Section 6 is the conclusion.

2. Copula-Based Seemingly Unrelated Quantile Regression Model

2.1. Seemingly Unrelated Quantile Regression under Asymmetric Laplace Distribution

Jun and Pinkse [5] and Waldmann and Kneib [6] extended the SUR model of Zellner [1] to quantile in feature to gain more robustness against outliers in the response measurements and allow explaining the entire conditional distribution of the outcome variables in the system equation. Thus, the asymmetric Laplace distribution (ALD) is generally specified as the likelihood of the model. Tu et al. [14] mentioned that the ALD has several mixture representations, for example, a scale mixture of normal with exponential distribution (2015) and a scale mixture of uniform with Gamma distribution according to Wichitaksorn et al. [15]. In this study, we consider the ALD as a mixture of normal as it was proved to be more efficient than as a scale mixture of uniform.

Let and be the dependent and independent variables of equation at time for and . The model is formulated as

where is the vector of parameters of equation i at given quantile is the unobserved error term of equation at time which is assumed to follow the ALD with mean zero, variance , and quantile level or quantile We note that the errors are allowed to be correlated across equations to gain more efficiency. The quantile level has the range . Thus, conditional quantile of given is simply

Typically, are joined through a multivariate distribution, especially the multivariate ALD. Nevertheless, this joint distribution is based on a linear relationship and the same error distributions. To relax this restriction, the copulas are suggested to model the nonnormal and nonlinear dependence structure in the multivariate SUQR model.

2.2. Multivariate Copulas

According to Sklar’s theorem [16], the n continuous marginals can be joined by copula function . Let be an -dimensional joint distribution with marginals . Thus, the n-dimensional joint distribution function can be defined as

where are the realization of random variables (or the standardized residuals). If the marginals , are continuous, then the copula associated to is unique, and the expression in equation (3) can be rewritten as

where is the uniform distribution. To construct the copula density distribution, it can be obtained by

In this study, we consider two copula classes, namely, elliptical copulas and Archimedean copulas. The explicit form of these copula classes is presented in the next section.

2.3. Posterior Distribution of Multivariate Copula Families

Two classes of copulas, namely, elliptical and Archimedean copulas, are presented. Elliptical copulas consist of two families, i.e., Gaussian and Student-t copulas. These copulas are of symmetric dependence structure. Another class is Archimedean copulas consisting of many families, but in this study, we consider only Clayton, Gumbel, Joe, and Frank copulas, which are prominent families and mostly used in many previous studies. To construct the posterior distribution of the copula, we multiply the copula density with the prior distribution.

2.3.1. Elliptical Copulas

Following Smith [17] and Smith et al. [18], the uniform prior is assumed for the copula parameter; thus, the posterior of elliptical copulas can be formulated as in the following:

(1)Posterior Gaussian copula: Let be the standard normal cumulative distribution and be the quantile function of the standard normal distribution; the posterior distribution of the Gaussian copula can be written as

 where is the dependence matrix of the Gaussian copula on quantile . The prior is assumed to have uniform distribution [0, 1]. Thus, the posterior is solely dependent on the copula density. As the Bayesian estimation is used to draw the parameter, we thus draw the updated dependence parameters from the truncated normal distribution [−1, 1] interval. (2) Posterior Student-t copula:

where and are the uniform prior density for and the exponential prior density for , respectively. is the copula density of the cumulative Student-t-distributed random vector. According to Demarta and McNeil [19], this copula density is

where is the dependence matrix of the Student-t copula on quantile and is the degree of freedom of each marginal. is the vector of the inverse Student-t distribution function. is the gamma distribution. In the proposal distribution, we randomly select and from the truncated normal distribution [−1, 1] interval and distribution [0, ] interval, respectively.

2.3.2. Archimedean Copulas

Different from the elliptical copula case, we apply the uninformative prior to these density functions; thus, the posterior density for each of the Archimedean copula is as follows:

The density of this class’ copulas, namely, Frank, Clayton, Gumbel, and Joe, varies as proposed in Hofert et al. [20]. We can write these copulas as follows:

(1)Frank copula:

 where and denotes the polylogarithm of order at .(2)Clayton copula:

(3)Gumbel copula:

 where , with ; and are the starting numbers of the first kind and the second kind, respectively (see Hofert et al. [20]).(4)Joe copula:

where

and , , and denotes the falling factorial.

To draw the updated dependence parameters of these copula functions, we randomly select the dependency parameter () from the proposal truncated normal distribution [0, ] for Clayton, [1, ] for Gumbel and Joe, and [, ] for Frank.

3. Bayesian Inference

We consider the Bayesian approach for estimating all unknown parameters in our proposed model. The Bayesian estimation requires the specification of the likelihood function and the prior distribution for all the estimated parameters. Hence, the posterior density of our proposed model is constructed by multiplying the full likelihood function of the model (ALD densities and copula density) with a prior density of the parameters. Let us consider the first part of the likelihood; according to Yu and Moyeed [21], the density of the asymmetric Laplace distribution is given by

where is the so-called check function defined by , with denoting the usual indicator function. Then, the posterior distribution for and can be written as

where is the prior distribution for and . Note that the copula function is used to join the errors in SUQR equations; hence, the full conditional posterior distribution of this model can be formulated as follows.

For elliptical copulas,

For Archimedean copulas,

where is the parameter set in the model.

In the estimation aspect, all parameters are drawn by an iterative Gibbs sampler with the Metropolis–Hastings algorithm over a partition of parameter blocks: (i) the unknown parameter ; (ii) the variance of the model ; and (iii) the copula dependence parameter . Concerning the prior specification in these three groups, we assume the following priors:

where is the vector of prior means, is the vector of hyperparameters for variance of B, and and are the positive hyperparameters for . We select these three priors since the sign of can be either positive or negative and asymmetric. For the case of copula dependence parameter , we assume it to be uniformly distributed with minimum zero and maximum one as there are various copula families considered in this model, and the prior information of the copula parameter is generally unknown as well as difficult to specify. To simplify our conditional posterior distribution, the uniform prior for the copula parameter is assumed.

As there are large parameter estimates in our model, the block Gibbs sampler with the MH sampling is considered to sample the parameters in the chain. The adaptive sampler algorithm can be explained as follows [22]:

(1)Starting at an initial parameter value : in this step, the initial value of is computed from the traditional quantile regression estimation while and , for elliptical copulas and Archimedean copulas, respectively.(2)Updating the candidate parameter based on the proposal function; the proposal function is simulated as follows:

(i)Simulate from the normal distribution with

  and , where , , and is the number of coefficients in each equation.

(ii)Simulate from the inverse Gaussian distribution, , with

  and , where is the random generation for the gamma distribution.

(iii)Simulate from the inverse-gamma distribution, , with

(iv)Simulate from the truncated uniform distribution, , where and are the lower bound and upper bound of the copula parameter, respectively.(3)At the j-th iteration, the acceptance function is employed to strike a balance between the following two constraints: (a) the sampler should tend to approach higher probability areas under the full posterior distribution and (b) the sampler should explore the space and avoid getting stuck at one site. Each candidate in each iteration is considered to be a proposal with acceptance probability equal to one, with the “proposal function” selected appropriately. In this computation aspect, the acceptance function can be defined as the ratio between posterior based on candidate parameters and posterior based on previously updated parameters . Thus, we can calculate the acceptance probability:

 Then, set with probability . Otherwise, set .(4)Repeat steps 2 and 3 for in order to obtain samples .

In the sampling method, we specify the number iteration to be 50,000, whereas the first 20,000 iterations are discarded as burn-in. Then, we can obtain the estimated parameter by averaging the remaining 30,000 simulated sets of parameters . As we consider many copula families, we employ the deviance information criterion (DIC) to compare the performance of different copula families.

4. Simulation Study

To learn about the performance of the Bayesian estimation for fitting our proposed model, we conduct four simulation studies:

(1)First, we examine and evaluate the accuracy of the Bayesian estimation on our model under Clayton, Gumbel, and Frank copulas(2)In the second part, we examine the performance of our model when misspecified copula is assumed(3)In simulation study 3, we evaluate the finite sample performance of the Bayesian estimation under various sample sizes(4)Finally, we investigate the performance of the Bayesian method in the high-dimension setting

4.1. Simulation Study 1: Accuracy in Parameter Estimation

In the simulation study, Archimedean copula families, namely, Clayton, Gumbel, and Frank, are considered to model the dependence structure of the SUQR. In this study, the simulation is the realization of SUQR with three equations. Thus, our simulated equation can be written as

Note that the three-dimension error terms, , are assumed to have asymmetric Laplace distribution with skewness or quantile . In this simulation, we first simulate the uniform margins from the three-dimension copula model. We set the true value of the correlation coefficient of Gumbel, Clayton, and Frank copula parameters to be 2.5, 0.5, and 2, respectively. Then, the obtained uniform data are transformed to be errors , and using the quantile function of the asymmetric Laplace distribution with variance and skewness or quantile parameter . The independent variables and are randomly simulated from a standard normal distribution. The true parameters in this simulation model are provided in Tables 1–3. We generate 1,000 data sets each with  = 500. To assess the accuracy of the Bayesian estimation, we consider the average of the estimated parameters and the average of the estimated standard deviations. The simulations are conducted for three quantile levels, say  = 0.25, 0.5, and 0.75, and the results are provided, respectively, in Tables 1–3, showing the average of the estimated parameters and their standard errors. We can see that our Bayesian estimation produces the unbiased parameter estimates for the SUQR model. It is observed that the estimated parameters are close to the true values in all cases, and the average standard deviation from the parameter is reliable. This simulation result indicates the adequacy and reliability of our Bayesian estimation of all unknown parameters in the copula-based SUQR model. In this simulation study, the hyperparameters a, b, c, and d are specified as 0.1 to reflect weak prior information.

Table 1: Mean and standard deviation of the Bayesian estimate of the Gumbel copula at three quantile levels.

Parameter
True
 = 0.25
True
 = 0.50
True
 = 0.75

1
0.918 (0.114)
1
0.828 (0.058)
1
1.258 (0.052)

1
1.187 (0.105)
5
5.024 (0.036)
2
1.697 (0.045)

1
1.035 (0.100)
1
0.949 (0.053)
1
1.225 (0.060)

4
3.580 (0.092)
2
2.384 (0.395)
4
3.570 (0.210)

−3
−3.009 (0.101)
−2
−1.819 (0.152)
−0.2
−0.198 (0.070)

2
2.091 (0.204)
2
1.984 (0.088)
2
3.242 (0.051)

2
2.308 (0.101)
5
4.883 (0.088)
1
0.904 (0.089)

3
3.736 (0.101)
2
1.938 (0.196)
3
2.659 (0.162)

3
3.125 (0.310)
3
3.020 (0.302)
3
3.147 (0.313)

2.5
2.458 (0.076)
2.5
2.535 (0.072)
2.5
2.401 (0.085)

Note: ( ) is the average standard deviation of the parameters.

Table 2: Mean and standard deviation of the Bayesian estimate of the Clayton copula at three quantile levels.

Parameter
True
 = 0.25
True
 = 0.50
True
 = 0.75

1
1.004 (0.036
1
1.141 (0.057)
1
0.929 (0.132)

1
1.114 (0.033)
5
5.071 (0.105)
2
1.944 (0.100)

1
1.027 (0.096)
1
0.8104 (0.053)
1
1.205 (0.073)

4
3.524 (0.033)
2
2.443 (0.035)
4
4.742 (0.055)

−3
−2.967 (0.035)
−2
−2.179 (0.115)
−0.2
−0.1284 (0.0034)

2
2.107 (0.197)
2
2.226 (0.083)
2
2.062 (1.413)

2
2.281 (0.100
5
5.030 (0.098)
1
0.833 (0.078)

3
3.752 (0.083)
2
1.627 (0.153)
3
2.655 (0.157)

3
3.158 (0.300)
3
3.096 (0.292)
3
2.709 (0.256)

0.5
0.3792 (0.095)
0.5
0.370 (0.097)
0.5
0.471 (0.106)

Note: ( ) is the average standard deviation of the parameters.

Table 3: Mean and standard deviation of the Bayesian estimate of the Frank copula at three quantile levels.

Parameter
True
 = 0.25
True
 = 0.50
True
 = 0.75

1
1.016 (0.127)
1
1.012 (0.067)
1
1.509 (0.046)

1
1.196 (0.153)
5
4.795 (0.403)
2
1.738 (0.043)

1
1.599 (0.124)
1
1.179 (0.149)
1
0.739 (0.060)

4
3.038 (0.809)
2
2.389 (1.001)
4
3.640 (0.349)

−3
−2.841 (0.061)
−2
−2.371 (0.578)
−0.2
−0.182 (0.060)

2
1.791 (0.045)
2
2.496 (0.551)
2
2.514 (0.242)

2
2.070 (0.104)
5
5.113 (0.110)
1
1.099 (0.109)

3
3.135 (0.209)
2
1.896 (0.189)
3
3.193 (0.224)

3
3.291 (0.332)
3
3.019 (0.366)
3
3.044 (0.301)

2
2.979 (0.481)
2
3.015 (0.458)
2
1.977 (0.419)

Note: ( ) is the average standard deviation of the parameters.

4.2. Simulation Study 2: Robustness and Kullback–Leibler Divergence

In this section, another simulation study is proposed to measure the performance of the model using the relative entropy, also known as Kullback–Leibler divergence (KLD) [23]. This relative entropy is the measurement of the difference between two probability distributions. Consider the continuous probability distribution, and let and denote the density of probabilities and ; thus, the KLD is given as

In this study, we define as an alternative or approximated posterior and define as a true posterior function when all parameters are known. We simulate the data sets in a similar way as in the previous simulation study (Section 4.1). If the copula family is known, then the proposed estimator and model could provide an accurate result. However, this estimation will be valid only in the simulation study. In practice, we need to select the true copula function to join the error terms in the SUQR model. Noh et al. [24] discussed the copula misspecification and suggested that the selection of the wrong copula function will bring about bias in the estimation of the model.

The purpose of this simulation study is to investigate the distance between the true SUQR function and its approximation when the copula is correctly specified and when the copula is misspecified. In this simulation study, we set the true posterior function to be the Student-t SUQR function for quantile levels 0.25, 0.50, and 0.75. We compare the true SUQR function and its approximation (in terms of posterior function) among the SUQR function family (i.e., Gaussian, Student-t, Joe, Clayton, Gumbel, and Frank) as well as the conventional SUQR function (noncopula-based).

Figure 1 illustrates the three-panel result for quantile levels 0.25, 0.50, and 0.75. As expected, the Student-t copula-based SUQR model achieves its minimum and is close to the true SUQR function line (dashed line) at every quartile level. We also made a performance comparison between our proposed model and the conventional SUQR function (M0) and found that our proposed model performs better than the conventional model since . In the case of misspecified copula function, we observed that the misspecified copula functions bring a larger deviation of the approximated SUQR function from the true function. According to these results, we can say that our proposed model is the robust model, and the incorrectly specified copula function will lead to the low accuracy of the model.

Figure 1: The performance of copula-based SUQR vs. conventional SUQR at different quantile levels, .

4.3. Simulation Study 3: The Finite-Sample Properties of the Bayesian Estimation

In the third simulation study, the finite-sample properties of the Bayesian estimation in our proposed model are investigated upon the calculated absolute Bias and mean squared error (MSE) of the estimator. Again, we simulate the data sets the same way as in the first simulation study. The absolute Bias and MSE can be calculated by

where  = 1,000 is the number of Monte Carlo replications and and are the estimated values and the true values, respectively. The sample sizes are fixed at 100, 500, and 1,000 d for each replication.

Tables 4–6 contain the results of the Bayesian estimation over the 1,000 simulated data sets with three different dimensions for the sample size. The most important finding from Tables 4–6 is that the Bayesian estimation provides reliable parameter estimates as the absolute Biases and MSEs are close to zero. In addition, the Biases and MSEs seem to be lower as the sample size increases. This result indicates that the Bayesian estimator is asymptotically unbiased and consistent for the estimated parameter in our model. This same pattern of convergence to zero is repeatedly obtained, considering different levels of quantile .

Table 4: Absolute Bias and MSE of parameter estimates from the Bayesian estimation of the Gumbel copula with different sample sizes.

Parameter
Absolute Bias
MSE

 = 0.25
 = 100
 = 500
 = 1000
 = 100
 = 500
 = 1000

0.0200
0.0115
0.0028
0.0747
0.0147
0.049

0.0224
0.0073
0.0048
0.0699
0.0133
0.0053

0.0356
0.0043
0.0058
0.0625
0.0120
0.0051

0.0019
0.0016
0.0010
0.1138
0.0136
0.0051

0.0465
0.0087
0.0012
0.0637
0.0108
0.0062

0.0103
0.0101
0.0032
0.0766
0.0105
0.0067

0.0431
0.0281
0.0083
0.0149
0.0094
0.0079

 = 0.50
 = 100
 = 500
 = 1000
 = 100
 = 500
 = 1000

0.0149
0.0073
0.0009
0.0409
0.0083
0.0039

0.0066
0.0074
0.0037
0.0437
0.0100
0.0042

0.0076
0.0082
0.0005
0.0492
0.0075
0.0036

0.0007
0.0005
0.0004
0.0689
0.0085
0.0034

0.0500
0.0117
0.018
0.0477
0.0084
0.0049

0.0018
0.0011
0.0030
0.0467
0.0089
0.0046

0.9005
0.8391
0.4918
0.8443
0.7433
0.1588

 = 0.75
 = 100
 = 500
 = 1000
 = 100
 = 500
 = 1000

0.0079
0.0013
0.0001
0.0499
0.0111
0.0015

0.0085
0.0062
0.0030
0.0665
0.0125
0.0125

0.0060
0.0037
0.0012
0.0464
0.0110
0.0014

0.0032
0.0018
0.0008
0.0795
0.0121
0.0030

0.0142
0.0123
0.0120
0.0544
0.0085
0.0025

0.0099
0.0080
0.0040
0.0855
0.0147
0.0058

0.2271
0.1321
0.1045
0.5613
0.4351
0.0105

Table 5: Absolute Bias and MSE of parameter estimates from the Bayesian estimation of the Clayton copula with different sample sizes.

Parameter
Absolute Bias
MSE

 = 0.25
 = 100
 = 500
 = 1000
 = 100
 = 500
 = 1000

0.0819
0.0172
0.0100
0.0912
0.0123
0.0060

0.0598
0.0033
0.0010
0.0793
0.0141
0.0052

0.0488
0.0130
0.0093
0.0641
0.0132
0.0058

0.0171
0.0105
0.0026
0.0915
0.0153
0.0047

0.0571
0.0114
0.0098
0.0682
0.0145
0.0060

0.0114
0.0131
0.0016
0.0747
0.0162
0.0045

0.7584
0.2990
0.2099
0.7940
0.8015
0.6381

 = 0.50
 = 100
 = 500
 = 1000
 = 100
 = 500
 = 1000

0.0414
0.0148
0.0098
0.0473
0.0101
0.0041

0.0109
0.0010
0.0007
0.0650
0.0104
0.0044

0.0405
0.0087
0.0015
0.0379
0.0091
0.0040

0.0291
0.0087
0.0045
0.0750
0.0091
0.0037

0.0328
0.0121
0.0095
0.0497
0.0104
0.0035

0.0094
0.0016
0.0013
0.0535
0.0089
0.0040

0.8278
0.4235
0.2075
0.7639
0.2002
0.1655

 = 0.75
 = 100
 = 500
 = 1000
 = 100
 = 500
 = 1000

0.0183
0.0026
0.0009
0.0646
0.0095
0.0049

0.0059
0.0012
0.0013
0.0822
0.0158
0.0073

0.0178
0.0048
0.0033
0.0533
0.0103
0.0056

0.0252
0.0037
0.0013
0.0816
0.0104
0.0049

0.0290
0.0190
0.0057
0.0709
0.0107
0.0054

0.0061
0.0010
0.0010
0.0715
0.0164
0.0045

0.5555
0.5387
0.4306
0.3322
0.2926
0.2831

Table 6: Absolute Bias and MSE of parameter estimates from the Bayesian estimation of the Frank copula with different sample sizes.

Parameter
Absolute Bias
MSE

 = 0.25
 = 100
 = 500
 = 1000
 = 100
 = 500
 = 1000

0.0155
0.0151
0.0061
0.0803
0.0100
0.0064

0.0003
0.0071
0.0042
0.0874
0.0127
0.0059

0.0367
0.0058
0.0016
0.0774
0.0123
0.0076

0.0114
0.0083
0.0053
0.0694
0.0134
0.0070

0.0078
0.0038
0.0030
0.0899
0.0135
0.0042

0.0078
0.0012
0.0008
0.0778
0.0108
0.0063

0.7279
0.7099
0.5023
0.6394
0.5238
0.4236

 = 0.50
 = 100
 = 500
 = 1000
 = 100
 = 500
 = 1000

0.0160
0.0070
0.0070
0.0411
0.0082
0.0002

0.0150
0.0057
0.0050
0.0627
0.0102
0.0004

0.0157
0.0068
0.0010
0.0497
0.0086
0.0006

0.0053
0.0015
0.0005
0.4675
0.0096
0.0004

0.0082
0.0001
0.0000
0.0441
0.0106
0.0050

0.0240
0.0002
0.0001
0.0549
0.0075
0.0023

0.0058
0.0095
0.0003
0.1395
0.0284
0.0015

 = 0.75
 = 100
 = 500
 = 1000
 = 100
 = 500
 = 1000

0.0013
0.0001
0.0009
0.0741
0.0108
0.0028

0.0216
0.0173
0.0042
0.0676
0.0126
0.0058

0.0107
0.0057
0.0034
0.0575
0.0101
0.0052

0.0179
0.0019
0.0018
0.0575
0.0159
0.0052

0.0217
0.0060
0.0011
0.0753
0.0131
0.0051

0.0245
0.0122
0.0021
0.0697
0.0119
0.0061

1.3613
1.0375
0.3658
0.9613
0.9184
0.8112

4.4. Simulation Study 4: Evaluation of the Performance of the Bayesian Estimation in the High-Dimension Copula-Based SUQR

Our model is supposed to propose a quite general model with the possibly arbitrary n-dimension copula-based model. So, the last simulation study is conducted to examine the performance of our model when the number of equations is large. We consider the sample size  = 1,000, the dimension n∈{3, 5, 7}, the dependency parameter as 3 for all dimensions, and the one-parameter (Archimedean) copula families of Gumbel, Clayton, and Frank. For each of these combinations, we generate a random sample of the corresponding size and compute . We repeat this procedure R = 1,000 times and compute the absolute Bias and MSE for each of the three quantile levels, . To simplify this simulation study, we set the marginal parameters for each equation as 1, say . Since the Bayesian estimation provides the similar marginal parameter estimates (parameters in each equation), we focus here on the copula parameter estimates .

Once all the data are simulated, we fit the copula-based SUQR model with , and the absolute Bias and MSE are recorded. The results are shown in Figure 2. We can see a pattern of diverging from zero of the absolute Bias and MSE when n increases. As a general statement, we can say that the Bias and MSE tend to deviate from zero when the number of copula dimensions increases, indicating that the estimates based on the Bayesian estimation may not provide a good result. Our result is in line with the study of Embrechts and Hofert [25], which reveals that the estimation of the parameter in the exchangeable copula models becomes increasingly severe in higher dimensions. Another reason is that the estimation cannot be improved with the higher computation cost. This is important in that we may not gain accurate parameter estimates in higher dimensions. However, our Bayesian estimation is still performing promisingly acceptable in high dimensions as the absolute Bias and MSE are not quite high.

Figure 2: Absolute Bias (first column) and MSE (second column) of the estimates of at different quantile levels, , where  = 1,000 repetitions.

5. Real Data Example

5.1. Estimation Results

In this section, we illustrate the applicability of our proposed model and the Bayesian estimation developed in this study, using the same data set as in Tansuchat et al. [13]. The data set consists of several variables measured in 177 months. Here, we focus on three US stock returns in the NASDAQ market, consisting of ADOBE, APPLE, and MICROSOFT (MICRO), and two additional factors consisting of small minus big (SMB) (a proxy for company size) and high minus low (HML) (a proxy for book-to-market values). The summary of the data statistics is presented in Table 7. Note that all data are presented in log-return form.

Table 7: Descriptive statistics.

APPLE
MICRO
ADOBE
NASDAQ
SMB
HML

Mean
0.011
0.001
0.003
0.001
0.299
0.246

Median
0.013
0.002
0.010
0.005
0.240
0.040

Maximum
0.137
0.092
0.138
0.061
6.860
13.910

Minimum
−0.185
−0.096
−0.263
−0.110
−6.540
−9.670

Std. dev.
0.048
0.032
0.050
0.027
2.549
2.681

Skewness
−0.718
−0.018
−1.580
−0.842
0.124
0.395

Kurtosis
5.157
3.937
9.969
4.758
2.747
7.248

Jarque−Bera
49.521
6.485
431.860
43.681
0.924
137.671

ADF test
4.388
6.245
9.545
4.589
6.414
5.879

Note: denotes 1% significant level.

According to the data description, the mean of SMB has the largest mean value. All stock returns show a negative skewness, while SMB and HML show a positive skewness, indicating that they are more likely to have a negative return on stock markets. The kurtosis of all series is greater than 3, except for SMB. This means that the returns are very much volatile and might exhibit nonnormal distribution. The Jarque–Bera test is then conducted to investigate the normal distribution of these data. The result shows that data reject the null hypothesis of normality at 1% significant level. In addition, we also conduct the augmented Dicky–Fuller (ADF) test to examine the stationary our variables and found that all data are stationary at 1% significant level.

Our empirical SUQR model is constructed under the Fama–French approach [26, 27]. The empirical model can be shown as

where is the return of asset i at time t, is the risk-free rate which is measured by the US. Treasury bill, , is the return of the market, is small minus big (a proxy for company size), and is high minus low (a proxy for book-to-market values). The parameters are referred to as beta risk, level of exposure to size risk, and value risk, respectively.

The choice of hyperparameters is very delicate in our Bayesian estimation; it is important to determine an approriate prior informations. In this empirical study, we decide to suggest three priors, which are as follows:

(1)Weak informative prior: the hyperparameters a, b, c, and d are specified as 0.1.(2)Diffusion prior a = 0, b = 100I, c = 1.5, and d = 0.5, where I is the identity matrix. Santos and Bolfarine [28] suggested that this prior has not presented any problem, and the posterior distribution seems to be insensitive to minor changes.(3)Informative prior: a = 2, b = 1, c = 1, and d = 0.1.

Prior to showing the estimated results, we compare the performance of various copula-based models as well as the conventional model of Jun and Pinkse [5] and Waldmann and Kneib [6] (M0) and the multivariate elliptical copula-based SUQR of Tansuchat et al. [13]. As a sensitivity analysis, we have also investigated various hyperparameters. By using the DIC, Table 8 reports the comparison results on 0.25, 0.50, and 0.75 quantiles. Among the trial runs of several alternative copula functions as well as hyperparameter prior, we learn that Clayton presents the lowest DIC for quantile level at , while Frank presents the lowest DIC for This indicates that different copula families can model different quantile models. In addition, we compare the best fit model with the conventional model of Jun and Pinkse [5] and Waldmann and Kneib [6] and with the two elliptical copula-(Gaussian- and Student-t-) based models of Tansuchat et al. [13]. The results reveal that our model is more adequate for this data set at all quantiles than other symmetric copula families. We also conclude that if the same copula is used to construct the joint between errors of the SUQR model, the unreliable result may be obtained. Moreover, we also learn that the model selection result seems to be insensitive to the hyperparameter prior as similar results are obtained. However, when we compare the DIC of these three priors, we find that the informative prior seems to be more valid than the others.

Table 8: Model selection and prior sensitivity check.

DIC
Gaussian
Student-t
Clayton
Gumbel
Joe
Frank
M0

Weak prior

−4028.48
−3602.62
−4441.98
−2301.10
−2980.83
−3316.27
−3018.11

−2501.28
−3313.28
−3293.29
−2238.56
−2218.23
−3517.42
−2458.22

−3074.94
−2982.58
−4075.10
−2430.20
−2405.26
−3352.41
−3001.32

Diffusion prior

−4149.18
−3515.92
−4193.47
−2596.13
−2596.13
−3462.22
−3084.18

−2473.64
−3298.69
−3203.15
−2251.99
2244.39
3500.01
−2444.98

−3000.63
−2814.36
−3978.64
−2221.54
−2211.18
−3354.69
−2987.14

Informative prior

−4052.51
−3687.39
−4512.68
−2347.66
−2001.35
−3398.41
−3021.01

−2641.31
−3343.14
−3300.88
−2288.98
−2210.35
−3601.01
−2500.54

−2942.64
−2848.54
−41871.1
−2200.87
−2209.91
−3358.69
−3005.05

The bold number indicates the lowest DIC for each quantile level.

Table 9 shows the estimated posterior mean and standard deviation parameters for quantiles 0.25, 0.50, and 0.75. We can observe that all estimated parameters of all stocks change when the quantile changes. This indicates that our model is robust against the outliers. In addition, we also find that the values of these copula dependences are different in different quantiles. It is noticed that the estimated copula parameters are quite high for all quantiles indicating a strong positive relationship among the errors of our model. To illustrate the correlation among the 19 marginals (cumulative of the standardized residual of each equation), three-dimension scatter plots among marginals for 0.25, 0.50, and 0.75 are illustrated in Figure 3. Finally, for this application study, we can see that the model estimates are reliable, and acceptable results are obtained. Moreover, as a sensitivity analysis, we have also provided the estimation results of the best fit model based on various hyperparameter priors. We find that the estimated posterior mean for each predictor is quite similar under various priors, and they do agree on the importance of the variables.

Table 9: Estimates of the parameters and their standard errors.

Parameter
Weak prior
Diffusion prior
Informative prior

−0.0057 (0.0042)
−0.0005 (0.0001)
0.0848 (0.0003)
−0.0514 (0.0132)
0.0071 (0.0012)
0.0750 (0.0445)
−0.0420 (0.0102)
0.0091 (0.0022)
−0.0662 (0.0004)

1.3604 (0.0141)
0.3232 (0.0247)
0.5866 (0.0033)
1.2079 (0.1158)
0.9491 (0.1442)
0.4971 (0.1254)
1.7481 (0.0212)
0.7030 (0.0254)
1.6397 (0.0047)

−0.0112 (0.0016)
−0.0011 (0.0001)
−0.0023 (0.0002)
−0.0013 (0.0002)
−0.0021 (0.0010)
−0.0004 (0.0002)
−0.0031 (0.0011)
−0.0053 (0.0001)
−0.006 (0.0002)

−0.0006 (0.0015)
−0.0001 (0.0001)
−0.0007 (0.0001)
−0.0015 (0.0021)
−0.0002 (0.0002)
−0.0001 (0.0010)
−0.0010 (0.0015)
−0.0008 (0.0001)
−0.0021 (0.0002)

−0.0124 (0.0028)
−0.0015 (0.0001)
0.0231 (0.0012)
−0.0538 (0.0144)
−0.0562 (0.0254)
0.0597 (0.0223)
−0.0448 (0.0024)
−0.0020 (0.0002)
0.0577 (0.0012)

0.8528 (0.0011)
0.7995 (0.0022)
0.4134 (0.0049)
1.6114 (0.1482)
1.0801 (0.2443)
0.8571 (0.1840)
1.8110 (0.0022)
1.0664 (0.0065)
0.6558 (0.0053)

−0.0044 (0.0015)
−0.0028 (0.0001)
−0.0124 (0.0001)
−0.0019 (0.0023)
−0.0050 (0.0012)
−0.0013 (0.0001)
−0.0034 (0.0014)
−0.0015 (0.0001)
−0.0043 (0.0001)

−0.0446 (0.0024)
−0.0001 (0.0001)
−0.0003 (0.0001)
−0.0004 (0.0001)
−0.0003 (0.0002)
−0.0009 (0.0003)
−0.0314 (0.0022)
−0.0007 (0.0001)
−0.0003 (0.0001)

−0.0111 (0.0039)
0.0025 (0.0011)
0.0428 (0.0021)
−0.0541 (0.0113)
0.0042 (0.0022)
0.0547 (0.0123)
−0.0431 (0.0042)
0.0030 (0.0012)
0.0546 (0.0020

1.3206 (1.0111)
1.1543 (0.0033)
0.8983 (0.3241)
0.9690 (1.0031)
0.8371 (0.0412)
0.7745 (0.4112)
1.7503 (1.0032)
1.0361 (0.0022)
0.6351 (0.3214)

−0.0011 (0.0020)
−0.0022 (0.0001)
−0.0485 (0.0023)
−0.0045 (0.0022)
−0.0013 (0.0005)
−0.0122 (0.0034)
−0.0078 (0.0032)
−0.0058 (0.0010)
−0.0068 (0.0015)

−0.0058 (0.0023)
−0.00028 (0.0001)
−0.0006 (0.0003)
−0.0021 (0.0015)
0.0001 (0.0001)
−0.0002 (0.0003)
−0.0017 (0.0008)
−0.0003 (0.0002)
−0.0004 (0.0002)

Copula
Clayton
Frank
Clayton
Clayton
Frank
Clayton
Clayton
Frank
Clayton

8.1421 (0.0121)
0.0001 (0.0012)
13.6251 (0.8812)
7.6548 (0.0135)
0.0001 (0.0001)
10.6257 (0.5547)
5.9844 (1.3554)
0.0001 (0.0001)
11.2898 (1.0548)

Note: the brackets ( ) denote the standard deviation.

Figure 3: Scatter plots among standardized residuals at (left panel), (middle panel), and (right panel).

To formally check the convergence of MCMC chains, the marginal posterior distributions can be visualized by plotting the histograms of the simulated parameter draws. Figure 4 displays the posterior distributions for quantiles 0.25, 0.50, and 0.75. We only present the case of copula parameters as similar results are obtained in other parameters. The plots show that the MCMC sampler converges to a normal distribution and mixes very well. We may thus conclude that the Bayesian method estimates the copula parameter quite well.

Figure 4: Posterior histograms for copula parameters for the copula-based SUQR model at (left panel), (middle panel), and (right panel). The vertical red line represents the mean parameter estimate.

5.2. Goodness-of-Fit Tests

To carry out the goodness-of-fit test for our proposed models, Cramer–von Mises (CvM) method is conducted in this section. Genest et al. [29] suggested that the CvM test is the most powerful test to check the goodness of fit of copula models. In this test, significant statistics indicate that the copula models based on the data are rejected. The result is reported in Table 10. For  = 0.25, the result shows that the Clayton copula yields the highest value, indicating Clayton copula-based SUQR offers a better fit for the data than other models at quantile 0.25. In the case of  = 0.50 and  = 0.75, Frank and Clayton copulas, respectively, provide the better fit than other copula-based models.

Table 10: Goodness-of-fit test for the copula-based SUQR model.

DIC
Gaussian
Student-t
Clayton
Gumbel
Joe
Frank

Weak prior

0.1247
0.0134
0.7849
0.0098
0.0148
0.0049

0.0525
0.2733
0.1149
0.0148
0.0849
0.5247

0.0049
0.0049
0.4549
0.0049
0.0443
0.1041

Diffusion prior

0.1244
0.0134
0.7848
0.0095
0.0141
0.0048

0.0521
0.2731
0.1147
0.0147
0.0844
0.5247

0.0042
0.0050
0.4549
0.0048
0.0443
0.1040

Informative prior

0.1240
0.0135
0.7845
0.0097
0.0150
0.0050

0.0521
0.2734
0.1145
0.0147
0.0851
0.5248

0.0047
0.0054
0.4543
0.0048
0.1445
0.1042

Note: this table presents the values of the goodness-of-fit tests for our copula-based models. Bold numbers indicate the highest value, which indicates that the copula model provides the best fit to the data.

6. Conclusion

In this paper, we introduced the multivariate exchangeable Archimedean copula to join the errors of the seemingly unrelated quantile regression (SUQR). The model becomes more accurate and robust against the outlier relationship between the dependent and independent variables. We also introduced the Bayesian Markov chain Monte Carlo approach to estimate the parameter sets of our proposed model. As the posterior distribution of the copula parameter does not appear to be in any form, therefore, we employ a Bayesian estimation together with a Gibbs sampler with the Metropolis–Hastings algorithm to infer the full posterior distribution. To examine the accuracy of our Bayesian estimation and performance of our proposed model, we present the simulation study and real data analysis.

Four simulation studies are conducted. The result of the first simulation study shows the accuracy of the Bayesian estimation. The results confirm that our proposed model is well estimated as reliable estimation results are obtained for every quartile level, and the parameter estimates on average are close to their true values. In the second simulation, we adopt the Kullback–Leibler divergence (KLD) to measure the distance between the true posterior probability distribution and the approximated posterior probability when the copula function is unknown. The result confirms the robustness of our model and presents the closest distance between the correctly specified and the true probability function. The third simulation is proposed for examining the finite-sample properties of the Bayesian estimation in our proposed model. The result shows that the Bayesian estimation provides reliable parameter estimates as the absolute Biases and MSEs converge to zero when the sample size increases. Finally, the performance of the Bayesian estimation in high-dimension copula-based SUQR is investigated, and the result reveals that we may not gain the accurate parameter estimates inhigher dimensions. However, the Bayesian estimation is still performing promisingly acceptable in high dimensions as the absolute Bias and MSE are not quite high.

In the real data application, we apply our proposed model to the data set provided by Tansuchat et al. [13]. We quantify and measure the risk of the stock price through the Fama–French model with three-factor analysis. We compare the performance of our model with the Archimedean copula-based SUQR models of Tansuchat et al. [13] as well as the conventional models of Jun and Pinkse [5] and Waldmann and Kneib [6]. The result shows that our proposed model is flexible and has the potential to capture the extreme market condition.

Although our proposed model presents a good performance in both simulation and real data studies, Gibbs algorithm with MH may not be the most suitable for our proposed model, given the large number of parameters to estimate. In addition, it can be difficult to choose the proposal functions in the Gibbs sampler with the MH sampling algorithm. For future study, the adaptive Gibbs sampler with the MH algorithm, which uses the history of the process to tune the appropriated proposal distribution, can be applied to our proposed model. Moreover, future research should consider applying our proposed model subject to pseudo-cyclical-structural changes since many financial time series exhibit behavioral change over time.

Data Availability

In this study, we use the simulated data to show the performance of our model, and the simulation processes are already explained in the paper. For the real data analysis section, we use the same data of Tansuchat et al. [13]. These data can be freely collected from http://www.investing.com (http://www.investing.com) or from Thomson Reuter Datastream. By the way, the data are available from the corresponding author upon request (woraphon.econ@gmail.com).

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

The authors are grateful for the financial support from the Centre of Excellence in Econometrics, Faculty of Economics, Chiang Mai University. They thank Dr. Laxmi Worachai for valuable comments to improve this paper.

References

A. Zellner, “An efficient method of estimating seemingly unrelated regressions and tests for aggregation bias,” Journal of the American Statistical Association, vol. 57, no. 298, pp. 348–368, 1962.

E. N. White and G. J. Hewings, “SPACE‐Time employment modeling: some results using seemingly unrelated regression estimators,” Journal of Regional Science, vol. 22, no. 3, pp. 283–302, 1982.

J. O. Adelegan, “Foreign direct investment and economic growth in Nigeria: a seemingly unrelated model,” African Review of Money Finance and Banking, vol. 1, pp. 5–25, 2000.

J. Frankel and J. Poonawala, “The forward market in emerging currencies: less biased than in major currencies,” Journal of International Money and Finance, vol. 29, no. 3, pp. 585–598, 2010.

S. J. Jun and J. Pinkse, “Efficient semiparametric seemingly unrelated quantile regression estimation,” Econometric Theory, vol. 25, no. 05, pp. 1392–1414, 2009.

E. Waldmann and T. Kneib, “Bayesian bivariate quantile regression,” Statistical Modelling, vol. 15, no. 4, pp. 326–344, 2015.

R. Koenker and G. Bassett Jr, “Robust tests for heteroscedasticity based on regression quantiles,” Econometrica: Journal of the Econometric Society, vol. 50, pp. 43–61, 1982.

N. Wichitaksorn and S. T. B. Choy, “Modeling dependence of seemingly unrelated Tobit model through copula: a Bayesian analysis,” Thailand Econometrics Society, vol. 3, pp. 6–19, 2011.

P. Pastpipatkul, P. Maneejuk, A. Wiboonpongse, and S. Sriboonchitta, “Seemingly unrelated regression based copula: an application on Thai rice market,” in Causal Inference in Econometrics, pp. 437–450, Springer, Berlin, Germany, 2016.

F. Louzada and P. H. Ferreira, “Modified inference function for margins for the bivariate clayton copula-based SUN Tobit Model,” Journal of Applied Statistics, vol. 43, no. 16, pp. 2956–2976, 2016.

E. Ivanov, A. Min, and F. Ramsauer, “Copula-based factor models for multivariate asset returns,” Econometrics, vol. 5, no. 2, p. 20, 2017.

Y. Zou, X. Zhong, J. Tang et al., “A copula-based approach for accommodating the underreporting effect in wildlife‒vehicle crash analysis,” Sustainability, vol. 11, no. 2, p. 418, 2019.

R. Tansuchat, P. Maneejuk, W. Yamaka, and S. Sriboonchitta, “Copulas based seemingly unrelated quantile regression,” Journal of Physics: Conference Series, vol. 1053, no. 1, p. 012102, 2018.

S. Tu, M. Wang, and X. Sun, “Bayesian variable selection and estimation in maximum entropy quantile regression,” Journal of Applied Statistics, vol. 44, no. 2, pp. 253–269, 2017.

N. Wichitaksorn, J. J. Wang, S. B. Choy, and R. Gerlach, “Analyzing return asymmetry and quantiles through stochastic volatility models using asymmetric Laplace error via uniform scale mixtures,” Applied Stochastic Models in Business and Industry, vol. 31, no. 5, pp. 584–608, 2015.

A. Sklar, Functions de repartition and dimensions et leurs marges, l’Institut de Statistique de L’Université de Paris, Paris, France, 1959.

M. S. Smith, Bayesian Approaches to Copula Modelling, Oxford Scholarship, Oxford, UK, 2011.

M. S. Smith, Q. Gan, and R. J. Kohn, “Modelling dependence using skew t copulas: Bayesian inference and applications,” Journal of Applied Econometrics, vol. 27, no. 3, pp. 500–522, 2012.

S. Demarta and A. J. McNeil, “The t copula and related copulas,” International Statistical Review, vol. 73, no. 1, pp. 111–129, 2005.

M. Hofert, M. Mächler, and A. J. Mcneil, “Likelihood inference for archimedean copulas in high dimensions under known margins,” Journal of Multivariate Analysis, vol. 110, pp. 133–150, 2012.

K. Yu and R. A. Moyeed, “Bayesian quantile regression,” Statistics & Probability Letters, vol. 54, no. 4, pp. 437–447, 2001.

J. Merhi Bleik, “Fully Bayesian estimation of simultaneous regression quantiles under asymmetric laplace distribution specification,” Journal of Probability and Statistics, vol. 2019, 2019.

S. Kullback and R. A. Leibler, “On information and sufficiency,” The Annals of Mathematical Statistics, vol. 22, no. 1, pp. 79–86, 1951.

H. Noh, A. E. Ghouch, and T. Bouezmarni, “Copula-based regression estimation and inference,” Journal of the American Statistical Association, vol. 108, no. 502, pp. 676–688, 2013.

P. Embrechts and M. Hofert, “Statistical inference for copulas in high dimensions: a simulation study,” ASTIN Bulletin: The Journal of the IAA, vol. 43, no. 2, pp. 81–95, 2013.

E. F. Fama and K. R. French, “The cross-section of expected stock returns,” Journal of Finance, vol. 47, no. 2, pp. 427–465, 1992.

E. F. Fama and K. R. French, “Multifactor explanations of asset pricing anomalies,” Journal of Finance, vol. 51, pp. 55–84, 1996.

B. Santos and H. Bolfarine, “Bayesian quantile regression analysis for continuous data with a discrete component at zero,” Statistical Modelling, vol. 18, no. 1, pp. 73–93, 2018.

C. Genest, B. Rémillard, and D. Beaudoin, “Goodness-of-fit tests for copulas: a review and a power study,” Insurance: Mathematics and Economics, vol. 44, no. 2, pp. 199–213, 2009.