Bayesian analysis of random effects panel interval-valued data models

1.   Introduction

Nowadays, big data, artificial intelligence, 5G, and other technologies have gradually entered people's vision, and have a huge impact on our lives. These also promote the total amount of data to show exponential growth. In the face of massive data, how to quickly and accurately obtain, process, mine, and integrate the required information becomes particularly important. According to different actual situations and purpose requirements, the data type is no longer limited to the definite point value data, but also gradually evolves to the uncertain data. As time goes on, many scholars at home and abroad have proposed interval numbers to describe uncertain data. Generally speaking, interval-valued data ^[1] arise due to one of the following two reasons: (ⅰ) imprecise observations of quantities, resulting in the translation of the measured value into an interval of possible values, and (ⅱ) information aggregation. The past three decades have witnessed enormous developments in statistical methods for interval-valued data. In particular, methods of regression focusing on interval-valued data have been extensively developed. For instance, Billard and Diday ^[2] proposed the center method (CM) for linear regression analysis of the interval-valued data, assuming that the interval-valued data had a uniform distribution. Billard and Diday ^[3] proposed the min-max method, taking the lower and upper limits of interval data as special point values, and establishing linear regression models for each of them. Billard and Diday ^[4] proposed the binary center and range method (BCRM), in which the explanatory variable contains information about both the center and range in the regression process. Then, LimaNeto and De Carvalho ^[5] proposed the centre and range method (CRM), which implied that the intervals were converted to centers and ranges. Based on the CRM, Lima Neto and Carvalho ^[6] introduced the constrained center and range method (CCRM), which added non-negative constraints of range coefficients. Giordani ^[7] adapted the CRM and presented a Lasso-based method for the interval-valued regression model. In addition, some scholars tried to establish interval-valued data regression models by other methods instead of the CRM. Maia and De Carvalho ^[8] proposed a multi-layer perceptron neural network based on interval-valued data and the Holt exponential smoothing method. Souza et al. ^[9] considered a parametrized approach that automatically extracted the best reference points from intervals. Kong et al. ^[10] proposed an interval local linear method (ILLM) to fit a regression model with interval-valued explanatory and response variables, which has no restrictions on the form of the regression function. Kong and Gao ^[11] studied the method of moments (MM) estimation method for interval-valued regression models. However, there is almost no literature studying the Bayesian modeling problem of interval panel data models.

As is known to all, a panel data set offers a certain number of advantages over traditional pure cross-section or pure time series data sets and statistical models combining cross-section and time series real-valued data become increasingly popular in economic research. For example, Nuroglu and Kunst ^[12] discussed the impact of exchange rate fluctuations on international trade flows by using panel data analysis and fuzzy data analysis methods. He et al. ^[13] constructed an interval slacks-based measure (SBM) of non-expected output, which analyzed China's environmental technology efficiency based on panel data of various provinces. But so far, the data in panel data models are all real-valued ^{[14,15,16,17,18,19]}, so it is very necessary to build a panel data model for interval-valued data. Recently, Ji et al. ^[20] introduced a panel data regression model for interval-valued data and constructed three kinds of panel interval-valued data regression models, and this is the first attempt to discuss panel interval-valued data models.

In addition, due to the recent dramatic evolution in advanced computational technologies, Bayesian inference has also received a lot of attention in recent years. In Bayesian regression fields, Park and Casella ^[21] proposed a Bayesian Lasso method for linear models. By using spline approximation, Xu and Zhang ^[22] introduced a Bayesian method for the partially linear model with heteroscedasticity based on the variance modeling technique. Castillo et al. ^[23] studied a high-dimensional linear regression with a sparse prior, which is a mixture of point masses at zero and continuous distributions. Pfarrhofer and Piribauer ^[24] proposed two shrink age priors to make Bayesian variable selection for high-dimensional spatial autoregressive models. Wang and Tang ^[25] considered Bayesian inference on a quantile regression model in the presence of nonignorable missing covariates. Zhang et al. ^[26] considered Bayesian quantile regression analysis for semiparametric mixed-effects double regression models based on the asymmetric Laplace distribution for the errors. Tang et al.^[27] for the first time used Elastic Net, a penalized method, for Bayesian quantile regression of panel data and derived the posterior distributions of all parameters based on the asymmetric Laplace prior distributions, and then constructed Gibbs sampling. Tao et al.^[28] proposed a Bayesian adaptive Lasso quantile regression method based on asymmetric exponential power distribution and applied the method to panel data. However, few works are devoted to constructing the Bayesian framework for interval-valued data. Zhang et al. ^[29] proposed the Bayesian nonparametric regression models by assuming that the upper and lower of the interval were distributed as an asymmetric Laplace distribution. Xu and Qin ^[30] extended the CRM for interval-valued regression models to the Bayesian framework for the first time, and proposed a bivariate Bayesian regression model based on the CRM with known and unknown covariance matrices.

However, to the best of our knowledge, there is little work done for constructing the Bayesian framework for panel interval-valued data. Hence, a Bayesian model for random effects panel interval-valued data is developed on the basis of the center and range method in this paper, and is compared with the Bayesian estimation based on the center method and the Bayesian estimation based on the minimum and maximum method. Finally, our proposed panel interval-valued data models are applied in forecasting of the Air Quality Index.

The outline of the paper is as follows. In Section 2, we introduce a Bayesian model of random effects panel interval-valued data based on the center and range method, and give the likelihood based on this model. In Section 3, the prior distribution and Bayesian posterior inference of the model are given in detail. In Section 4, the results of parameter estimation are obtained through simulation and compared with other methods to illustrate the feasibility of the proposed method. In Section 5, we apply the model to real data, and in Section 6, the paper is concluded with a brief discussion.

2.   Random effects panel interval-valued data models

For panel interval-valued data set $S = \{(X_{it}, y_{it})|i = 1, 2, \ldots, n;t = 1, 2, \ldots, T\}$, $y_{it} = [y_{it}^l, y_{it}^u]$ is assumed as the observed interval-valued dependent variables, where the superscripts $l$ and $u$ represent the lower and upper bounds of the interval, respectively, and let $X_{it} = (x_{it1}, x_{it2}, \ldots, x_{itp})^{T}$ be $p \times 1$ interval-valued independent vectors with $x_{itj} = [x_{itj}^l, x_{itj}^u], i = 1, 2, \ldots, n;t = 1, 2, \ldots, T;j = 1, 2, \ldots, p$.

In the following, we consider Bayesian analysis of random effects panel interval-valued data models:

where $y_{it}^c = \frac{y_{it}^u + y_{it}^l}{2}$, $y_{it}^r = \frac{y_{it}^u - y_{it}^l}{2}$, $X_{it}^c = \frac{X_{it}^u + X_{it}^l}{2}$, $X_{it}^r = \frac{X_{it}^u - X_{it}^l}{2}$. $\beta^c = (\beta_1^c, \ldots, \beta_p^c)^\text{T}$ and $\beta^r = (\beta_1^r, \ldots, \beta_p^r)^\text{T}$ are vectors of p-dimensional unknown parameters. In addition, $\varepsilon_{it}^c$ and $\varepsilon_{it}^r$ are mutually independent and identically distributed normal random variables with zero mean and variances $\sigma_c^2$ and $\sigma_r^2$, respectively. $\alpha_i^c$ and $\alpha_i^r$ denote the random effects associated with individual $i$ and assume that the random effects $\alpha_i^c\sim N(0, \phi_c^2)$, $\alpha_i^r\sim N(0, \phi_r^2)$. Then for time t, model (2.1) can be written as

Then for $s\in \{c, r\}$, the center and range models take on the following form:

where ${{\boldsymbol{Y}}_{t}}^s = (y_{1t}^s, y_{2t}^s, \ldots, y_{nt}^s)^\text{T}$, ${\boldsymbol{\alpha}}^s = (\alpha_1^s, \alpha_2^s, \ldots, \alpha_n^s)^\text{T}$, ${{\boldsymbol{X}}_{t}}^s = (X_{1t}^s, X_{2t}^s, \ldots, X_{nt}^s)^\text{T}$, ${\boldsymbol{\beta}}^s = (\beta_1^s, \beta_2^s, \ldots, \beta_p^s)^\text{T}$, ${\boldsymbol{\varepsilon}_t}^s = (\varepsilon_{1t}^s, \varepsilon_{2t}^s, \ldots, \varepsilon_{nt}^s)^\text{T}$.

For convenience, matrix notation and vector notation are used to represent variables and models. Let ${{\boldsymbol{Y}}}^s = (({{\boldsymbol{Y}}_1}^s)^\text{T}, ({{\boldsymbol{Y}}_2}^s)^\text{T}, \ldots, ({{\boldsymbol{Y}}_T}^s)^\text{T})^\text{T}$, ${{\boldsymbol{X}}}^s = (({{\boldsymbol{X}}_1}^s)^\text{T}, ({{\boldsymbol{X}}_2}^s)^\text{T}, \ldots, ({{\boldsymbol{X}}_T}^s)^\text{T})^\text{T}$, ${\boldsymbol{\varepsilon}}^s = (({\boldsymbol{\varepsilon}_1}^s)^\text{T}, ({\boldsymbol{\varepsilon}_2}^s)^\text{T}, \cdots, ({\boldsymbol{\varepsilon}_T}^s)^\text{T})^\text{T}$, ${{\boldsymbol{D}}}^s = {{\mathbf{1}}_T}^s \otimes {{\boldsymbol{I}}_n}^s$, where "$\otimes$" denotes the Kronecker product, ${{\mathbf{1}}_T}^s$ is a $T \times 1$ dimensional vector with all elements 1, and ${{\boldsymbol{I}}_n}^s$ is an $n \times n$ dimensional unit matrix. Then, model (2.3) can also be written in matrix form,

where ${\boldsymbol{\alpha}}^s \sim N(0, \phi_s^2 I_n)$, ${\boldsymbol{\varepsilon}}^s \sim N(0, \sigma_s^2 I_N)$, $N = n\times T$.

Based on the above model, the likelihood function of the model parameters is defined as follows:

3.   Bayesian analysis of the models

To estimate the unknown parameters ${\boldsymbol{\beta}}^s$, $\sigma_s^2$, and $\phi_s^2$, we implement a Bayesian approach ^[31]. Therefore, we should appoint prior distributions for the parameters of the models. First, we assume ${\boldsymbol{\beta}}^s$ has normal prior distributions, i.e., ${\boldsymbol{\beta}}^s \sim N({\boldsymbol{\beta}_0}^s, {\boldsymbol{\Sigma}_0}^s)$. It is assumed that the hyperparameters $\beta_0^s$, $\Sigma_0^s$ are known. In addition, the prior information for other unknown parameters is $\sigma_s^2\sim IG(a_0^s, b_0^s)$, $\phi_s^2\sim IG(c_0^s, d_0^s)$, where $a_0^s$, $b_0^s$, $c_0^s$, $d_0^s$ are hyperparameters to be given, and $"IG"$ means the inverse gamma distribution. Thus, the joint priors of all of the unknown parameter are given by

From the likelihood function (2.5) and prior distributions (3.1), we can obtain the following theorems and give a brief proof of Theorem 1. The proofs of the other theorems are similar and will not be written in detail here.

Theorem 1. Suppose that the parameter ${\boldsymbol{\beta}}^s$ follows a normal prior distribution, i.e., ${\boldsymbol{\beta}}^s\sim N(\beta_0^s, \Sigma_0^s)$, and then the posterior distribution of ${\boldsymbol{\beta}}^s$ follows the normal distribution,

where $\mu_{\beta^s} = \Sigma_{\beta^s} (\frac{({{\boldsymbol{X}}}^s)^\mathit{\text{T}}({{\boldsymbol{Y}}}^s-{{\boldsymbol{D}}}^s{\boldsymbol{\alpha}}^s)}{\sigma_s^2} + ({\boldsymbol{\Sigma}_0}^s)^{-1}{\boldsymbol{\beta}_0}^s)$, $\Sigma_{\beta^s} = (\frac{1}{\sigma_s^2}({{\boldsymbol{X}}}^s)^\mathit{\text{T}} {{\boldsymbol{X}}}^s+({\boldsymbol{\Sigma}_0}^s)^{-1})^{-1}$.

The posterior distribution of ${\boldsymbol{\beta}}^s$ is as follows, and we make a transformation of its form:

Obviously, the conditional distribution of ${\boldsymbol{\beta}}^s$ agrees with the form of the probability density of the multivariate normal distribution, so ${\boldsymbol{\beta}}^s$ obeys the normal distribution, thus this paper obtains Theorem 1.

Theorem 2. Based on model assumptions, we can obtain the posterior distribution of ${\boldsymbol{\alpha}}^s$ as follows, that is

where $\mu_{\alpha^s} = \Sigma_{\alpha^s}\frac{1}{\sigma_s^2}({{\boldsymbol{D}}}^s)^\mathit{\text{T}}({{\boldsymbol{Y}}}^s-{{\boldsymbol{X}}}^s{\boldsymbol{\beta}}^s)$, $\Sigma_{\alpha^s} = (\frac{({{\boldsymbol{D}}}^s)^\mathit{\text{T}}{{\boldsymbol{D}}}^s}{\sigma_s^2}+\frac{{{\boldsymbol{I}}_n}^s}{\phi_s^2})^{-1}$.

Theorem 3. Suppose that the prior distribution of $\sigma_s^2$ is known, i.e., $\sigma_s^2\sim IG(a_0^s, b_0^s)$, and then the posterior distribution of $\sigma_s^2$ follows the inverse gamma distribution

where $a_*^s = \frac{N}{2}+a_0^s$, $b_*^s = \frac{1}{2}({{\boldsymbol{Y}}}^s-{{\boldsymbol{D}}}^s{\boldsymbol{\alpha}}^s-{{\boldsymbol{X}}}^s{\boldsymbol{\beta}}^s)^\mathit{\text{T}}({{\boldsymbol{Y}}}^s-{{\boldsymbol{D}}}^s{\boldsymbol{\alpha}}^s-{{\boldsymbol{X}}}^s{\boldsymbol{\beta}}^s)+b_0^s$.

Theorem 4. Suppose that the prior distribution of $\phi_s^2$ is known, i.e., $\phi_s^2\sim IG(c_0^s, d_0^s)$, and then the posterior distribution of $\phi_s^2$ follows the inverse gamma distribution

where $c_*^s = \frac{n}{2}+c_0^s$, $d_*^s = \frac{({\boldsymbol{\alpha}}^s)^\mathit{\text{T}}{\boldsymbol{\alpha}}^s}{2} +d_0^s$.

From Theorems 1 to 4, the posteriors of ${\boldsymbol{\beta}}^s$, $\sigma_s^2$, $\phi_s^2$ are all familiar distributions and can be sampled directly. The specific algorithm is shown in Table 1. Thus based on the above MCMC algorithm, a converged posteriori sample can be collected and this sample is noted as ${\boldsymbol{\theta}}^{(sim)} = ({{\boldsymbol{\beta}}^s}^{(sim)}, {\sigma_s^2}^{(sim)}, {\phi_s^2}^{(sim)})$, $sim = 1, 2, \ldots, M, M < Sim$. As such, the posteriori estimation of the parameters $(\hat{{\boldsymbol{\beta}}^s}, \hat{\sigma_s^2}, \hat{\phi_s^2})$ can be respectively estimated as follows:

During the implementation process, the algorithm is repeated 100 times, and the final parameter estimation is the result based on the 100 times average.

4.   Simulation study

In this section, we investigate the performance of the proposed model and Bayesian estimation method via Monte Carlo simulation. We compare Bayesian analysis of random effects panel interval-valued data models based on the center and range method (BCRM) with Bayesian analysis of random effects panel interval-valued data models based on the center method (BCM) and Bayesian analysis of random effects panel interval-valued data based on the minimum and maximum method (BMMM). In order to demonstrate the quality of Bayesian estimation and prediction, Section 4.1 introduces the measurement methods. Section 4.2 shows the data generation process and parameter details. Then we present all simulation results and conclusions in Section 4.3.

4.1.   Measurements

There are three measurements for evaluating the performances of different models:

1) The upper and lower bounds root mean-square errors, i.e., $RMSE_U$ and $RMSE_L$ ^[5] are

2) The root mean-square error $RMSE_H$ ^[32] is

3) The rate of different intervals(RI) defined by Hu and He ^[33] is

where $\omega(\cdot)$ represents the width of the intervals.

4.2.   Data generating process

In this subsection, based on the center and range method (BCRM), four configurations C1, C2, C3, and C4 are generated, as shown in Table 2.

The data sets are generated from the following models:

Referring to Xu and Qin ^[30], first we fixed the regression coefficients ${\boldsymbol{\beta}}^c = (8, 7, 6)^T$ and ${\boldsymbol{\beta}}^r = (4, 6, 5)^T$. Second, ${{\boldsymbol{X}}}^s$ is generated from the uniform distribution $U(0, 2)$ and ${\boldsymbol{\alpha}}^s$ is generated from $N(0, 0.5I_n)$ and fixed for each replication. We generate ${{\boldsymbol{Y}}}^s$ from the above model, where ${\boldsymbol{\varepsilon}}^s \sim N(0, \sigma_s^2I_N)$. In addition, we consider the noninformative prior information type of hyperparameter values for unknown parameters ${\boldsymbol{\beta}}^s$, $\sigma_s^2$, and $\phi_s^2$ in the simulation: ${\boldsymbol{\beta}_0}^s = {{\mathbf{0}}_3}$, ${\boldsymbol{\Sigma}_0}^s = 10\times {{\boldsymbol{I}}_3}$, $a_0^s = b_0^s = c_0^s = d_0^s = 0.01$, where ${{\mathbf{0}}_3}$ is a 3-dimensional vector with all elements being 0. Further, we choose n to be 50,100 and T to be 12, and therefore the sample sizes are $N = n\times T$ = 600 and 1200.

For configurations C1 and C2, we assume that there is a linear relationship hidden in ${\boldsymbol{\varepsilon}}^c$ and ${\boldsymbol{\varepsilon}}^r$ that ${\boldsymbol{\varepsilon}}^r = {\boldsymbol{\varepsilon}}^c+ {{\boldsymbol{e}}}$, where ${{\boldsymbol{e}}}$ is the random error generated from $N(0, \sigma^2 I_N)$ and $\sigma^2$ = 0.5. In configurations C3 and C4 it is assumed that there is no linear relationship between ${\boldsymbol{\varepsilon}}^c$ and ${\boldsymbol{\varepsilon}}^r$.

Based on the above settings and the generated data sets, the preceding proposed algorithm is used to evaluate the Bayesian estimation of unknown parameters based on 100 replications. In order to obtain better and more accurate results, we collect the observation results of the following J = 2000 for statistical inference by discarding the first 3000 burn-in iterations.

4.3.   Results

Tables 3 and 4 show the Bayesian estimates when n = 50,100 and T is fixed at 12 based on 100 repetitions in all configurations. The accuracy of the estimation is expressed by BIAS and standard deviation SD, and the following conclusions are obtained:

1) Comparing configurations C1, C2, and C4, it can be found that when n is the same, SD of C1 are basically the smallest among the three configurations, while SD of C4 are basically the largest, which indicates that the smaller the variance of the error term, the smaller the estimated error.

2) As the sample size increases, the accuracy of Bayesian estimation for the BCRM gradually gets better. For example, in configuration C1, when n = 50 and 100, the SD of $\beta_1^c$ is 0.0476 and 0.0336, respectively.

In general, the BIAS and SD of all parameters in the four configurations are very small, indicating that the Bayesian estimation effect is good.

To investigate the convergence of the proposed algorithm, we visualize the values of 5000 iteration updates in a loop for each configuration when n = 50, and the results are shown in Figure 1. It can be seen that all the parameters are roughly in a straight line, indicating that the algorithm has a considerable convergence speed.

The next step is to do a predictive study of the models. We compare the Bayesian method based on the center and range method with the other two methods mentioned above, and all results are listed in Tables 5–8, where the standard deviation is in parentheses. 75% of the data is selected as the training set and 25% of the data as the test set. Each case is repeated 100 times and the results of the prediction error are averaged. The following conclusions can be drawn:

1) When n is the same, the BCRM is better than the BCM and the BMMM. For example, in C1 configuration, when n = 50, the $RMSE_L$, $RMSE_U$, and $RMSE_H$ of the BCRM are 0.723, 1.585, and 1.631, respectively, and RI is 0.931. The $RMSE_L$, $RMSE_U$, $RMSE_H$, and RI of the BCM and BMMM methods are 6.531, 6.664, 7.094, 0.714 and 6.296, 2.360, 6.320, 0.768, respectively. The first three data are larger than in the BCRM, while RI is smaller than in the BCRM. This shows that the BCRM has better prediction effect than the BCM and the BMMM.

2) Also in the case of n = 50, the $RMSE_L$ of the BCRM method in C1 configuration is 0.723 and the standard deviation is 0.041, and that of the BCRM method in C2 configuration is 0.738, 0.044, indicating that the smaller the variance of the random error, the higher the accuracy of the prediction.

3) As the sample size increases, the effectiveness of Bayesian prediction gradually improves. For example, in C2 configuration, the standard deviation of the $RMSE_L$ for the BCRM decreases with increasing the sample size, changing from 0.044 to 0.029.

In general, the prediction effect of the BCRM based on the center and range method is satisfactory in different configurations and different sample sizes.

5.   Empirical analysis

This section applies the proposed model to the estimation and prediction of the AQI and compares it with other methods. The concentration of all kinds of pollutants changes with space and time. The panel interval-valued data can be used to describe this variation, and this paper aims to construct panel interval-valued data models for the AQI. Based on the AQI-related data in Beijing, Tianjin, Shijiazhuang, and Chongqing in China, this study selected AQI-related data from 4 representative cities for 40 consecutive days (2023.7.20–2023.8.28). Among them, the data of the first 30 days are used to train the models, and the remaining data are used to test the models.

For the panel interval-valued data set $S = (X_{it}, y_{it})$, $i = 1, 2, \ldots, N$; $t = 1, 2, \ldots, T$, $y_{it} = [y_{it}^L, y_{it}^U]$ is considered to be the observed interval-valued dependent variable, which is the AQI. $y_{it}^L$ represents the minimum value of the AQI in the ith city on date t and $y_{it}^U$ represents the maximum value of the AQI in the ith city on date t. $X_{it} = (X_{it1}, X_{it2}, \ldots, X_{it6})^{T}$ are interval-valued independent vectors, which represent the values of $CO$, $NO_2$, $O_3$, $PM_{10}$, $PM_{2.5}$, and $SO_2$, respectively. $X_{itj} = [a_{itj}, b_{itj}]$, $i = 1, 2, 3, 4$, $t = 1, 2, \ldots, 40$, $j = 1, 2, \ldots, 6$, $a_{itj}$ indicates the minimum value of the jth pollutant in the ith city on date t, $b_{itj}$ indicates the maximum value of the kth pollutant in the ith city on date t.

First of all, the center data and range data of all dependent variables are made into a Q-Q plot and the results are shown in Figure 2. It is considered that the data is approximately normally distributed.

$\alpha_i$ is a random effect and $\beta^s = (\beta_1^s, \beta_2^s, \ldots, \beta_6^s)$, where the prior information $\beta_0^s$ is generated from the linear least squares method. Other prior information is consistent with the numerical simulation. We take the mean value of the last 2000 observations as the estimated value of the parameters based on 5000 iterations, and the results are shown in Table 9.

It can be seen from the parameters of the central data that $NO_2$, $O_3$, $PM_{10}$, and $PM_{2.5}$ have a positive impact on the AQI, that is, the greater the concentration of the pollutants, the greater the AQI and the more severe the air pollution, while $CO$ and $SO_2$ have a negative impact on the AQI. Similarly, in the range data, $CO$, $O_3$, $PM_{10}$, and $PM_{2.5}$ have a positive impact on the AQI, while $NO_2$ and $SO_2$ have a negative impact on the AQI. The first 75% data (the first 30 days) are selected as the training set, and the last 25% data (the last 10 days) as the test set, and the AQI is predicted by the three methods mentioned above. The prediction results are shown in Table 10.

From the table, we can see that the BCRM and the BMMM are relatively similar, but consistently better than the BCM. Specifically, the results obtained for the BCRM method are as follows. First, the BCRM method has relatively small values in the indicators of $RMSE_L$, $RMSE_U$, and $RMSE_H$, which proves the accuracy of BCRM method estimation. It reached 0.5820 on $RI$, which is relatively good. In addition, the BMMM method also predicts well and even better than the BCRM method in some indicators, which may be due to the fact that the interval-valued data in the empirical study of this paper are presented in the form of minimum and maximum values. Then, comparison charts of different methods based on real data and forecast data are listed in Figure 3, where a solid line represents the predicted data, a dotted line represents the real data, the horizontal coordinate represents the predicted frequency, and the vertical coordinate represents the predicted value of the AQI. As can be seen from Figure 3, the trend of the forecast data based on the BCRM is generally consistent with the real data and the distance between the real data and forecast data is small, which also explains the accuracy of the estimation of the BCRM method to a certain extent. The distance between real data and forecast data is larger for the BCM method, while the distance between real data and forecast data is also smaller for the BMMM method.

6.   Conclusions and discussion

In this paper, a Bayesian model based on the center and range method (BCRM) of random effects panel interval-valued data is proposed and compared with the Bayesian model based on the center method (BCM) of random effects panel interval-valued data and the Bayesian model based on the maximum and minimum method (BMMM) of random effects panel interval-valued data. The results show that: on the one hand, the deviation and standard deviation of parameters estimated by the BCRM model are very small; and on the other hand, in the case of different configurations and different sample sizes, the prediction errors $RMSE_L$, $RMSE_U$, and $RMSE_H$ of the BCRM are basically the smallest, while the interval coverage $RI$ is basically the largest, that is, the BCRM has the best prediction effect. The model is applied to the empirical study of AQI prediction, and the results show the effectiveness of the proposed Bayesian model.

In addition, the model deserves further study. Specific considerations are as follows: Quantile regression has greater flexibility in the distribution of random errors and is therefore robust to non-normal errors and outliers. Therefore, how to combine quantile regression with a Bayesian model of random effect panel interval-valued data is worthy of further study.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

This research was supported by the National Social Science Fund of China (Grant No. 23BTJ069).

Conflict of interest

The authors declare there are no conflicts of interest.