Several explorations on how to construct an early warning system for local government debt risk in China

This paper aims to explore several ways to construct a scientific and comprehensive early warning system (EWS) for local government debt risk in China. In order to achieve this goal, this paper studies the local government debt risk from multiple perspectives, i.e., individual risk, contagion risk, static risk and dynamic risk. Firstly, taking China’s 30 provinces over the period of 2010~ 2018 as a sample, this paper establishes early warning indicators for individual risk of local government debt, and uses the network model to establish early warning indicators for contagion risk of local government debt. Then, this paper applies the criteria importance though intercrieria correlation (CRITIC) method and coefficient of variation method to obtain the proxy variable Ⅰ, which combines the above two risks. Secondly, based on the proxy variable Ⅰ, both the Markov-switching autoregressive (MS-AR) model and coefficient of variation method are used to obtain the proxy variable Ⅱ, which comprehensively considers the individual risk, contagion risk, static risk and dynamic risk of local government debt. Finally, machine learning algorithms are adopted to generalize the EWS designed in this paper. The results show that: (1) From different perspectives of local government debt risk, the list of provinces that require early warning is different; (2) The support vector machines can well generalize our EWS.


Introduction
Nowadays, the growing public debt and on-going debt crisis in many countries, such as the US, China, Spain, Italy, Brazil and Mexico, are alarming policymakers and scholars. Similarly, ever-growing sovereign debt at sub-national levels like US states and China local governments as well as at the municipality level has caused worldwide concern. Among the highly indebted countries, China is gaining particular attention. This is because since China granted local governments the right to issue bonds in 2009, the balance of local government debt in China has surprisingly increased by 200% in a short decade, i.e., from 10717.491 billion yuan in 2010 to 21330 billion yuan in 2019 (China's National Bureau of Statistics and National Audit Office). Such a rapid growth has spurred an academic debate on how to strengthen the supervision of public debt and how to avoid the outbreak of debt crisis. As a result, the construction of early warning system (EWS) for public debt risk has become one of the research hotspots in academia. An effective EWS can help governments to assess the local government debt risk and take timely measures to prevent the risk. Therefore, scholars have adopted various models to construct the EWS, such as the stochastic model [1][2][3], general equilibrium model [4] and econometric model [5,6]. Although the existing EWSs have great enlightenments in assessing and predicting local government debt risk, further questions need to be explored. Do the existing EWSs consider various types of local government debt as much as possible? Do they take the contagion risk of local government debt into account? Do they explore the regime-switching risk of local government debt? In view of these questions, this paper takes China's 30 provinces over the period of 2010 to 2018 as a sample, and intends to go deep into the EWS for local government debt risk.
The marginal contributions of this paper are as follows: (1) Besides the explicit local government debt, the implicit local government debt will be included in this paper, making our EWS cover more types of local government debt. (2) Based on the spatial spillover effect of China's local government debt risk [7][8][9], the contagion risk of local government debt will not be ignored in our EWS. (3) Rather than using the subjective weighting methods, the objective weighting methods (the criteria importance though intercrieria correlation method and coefficient of variation method) will be used to determine the weight of early warning indicators, thus ensuring the fairness of our indicator weighting. (4) Considering that the data of local government debt risk of China's each province are typically non-linear time series and the Markov-switching autoregressive (MS-AR) model can effectively capture different regimes of time series data [10][11][12], the MS-AR model will be used to study the static risk and dynamic risk of local government debt for each province. (5) The machine learning algorithms will be introduced to debug an optimal classifier to generalize the EWS, thus improving the convenience and practicality of our EWS. However, because of the small data size, it is necessary to follow the changes in China's local government debt risk and update the sample data in future research.
The rest of the paper is organized as follows. Section 2 provides the literature review. Section 3 describes the construction idea of our EWS. Based on the spatial spillover effect of China's local government debt risk, Section 4 calculates the proxy variable Ⅰ which combines the individual risk and contagion risk of local government debt, and obtains the proxy variable Ⅱ which comprehensively considers the individual risk, contagion risk, static risk and dynamic risk of local government debt. Section 5 debugs an optimal classifier based on machine learning algorithms and Section 6 concludes.

Literature review
As western countries granted governments the right to raise debt earlier than China, western scholars made pioneering attempts in constructing the EWS. The study on the EWS for public debt can be tracked back to 1980s [13], and has sprung up in the US, Italy and Colombia. Recently, new models and new techniques are constantly introduced into enhancing the function of EWS, and the research object is extending from sovereign debt risk to sub-national debt risk.
Among the existing EWSs constructed by western scholars, the multinomial logistic regression belonging to the econometric models has been the most widely used. Pioneering studies by Ciarlone and Trebeschi [14] and Fuertes and Kalotychou [15] have applied the multinomial logistic regression to estimate the sovereign debt risk of one country and several countries, respectively. Since then, following in the footsteps of Ciarlone and Trebeschi [14] and Fuertes and Kalotychou [15], western scholars have repeatedly confirmed the predictive value of multinomial logistic regression, such as Fournier and Bétin [16] (see Table 1). Besides the multinomial logistic regression, other models have enriched the construction of EWS. For example, the financial particle theory and objective decomposition method have been comprehensively used by Turay et al. [17] to develop an EWS for local government debt risk. The parametric proportional hazards regression, conventional logistic regression and Bayesian model averaging have been employed by Kamra [18] to explain the uncertainty of EWS for debt crisis in emerging markets. A new regression-tree based approach along with a signal model has been adopted by Savona and Vezzoli [19] to achieve a balance between in-and out-of-sample predictability of sovereign defaults. A more powerful dynamic-recursive early warning model has been established by Dawood et al. [20], who sends more accurate out-of-sample warning signals of sovereign debt crisis.
With the development of software and econometric models, a new trend in the latest work is to introduce the machine learning algorithms into the EWS for municipal debt risk. Representatives are Antulov-Fantulin et al. [21], Alaminos et al. [22] and Zahariev et al. [23]. To be specific, an optimal machine learning model has been selected by Antulov-Fantulin et al. [21] from the gradient boosting machine, random forest, lasso and neural network, finding that the gradient boosting machine performs the best in predicting bankruptcy of local government in Italy. The fuzzy decision trees, AdaBoost, extreme gradient boosting and deep learning neural decision trees all have been proved a good early warning performance for sovereign debt crisis by Alaminos et al. [22]. The support vector machine (SVM) has been confirmed applicable to model the dependence of the debt ratio of Italy and Greece by Zahariev et al. [23].
Faced with less than 15 years of local government debt data, Chinese scholars still made valuable explorations on the construction of EWS. Representatives are the BP neural network by Shi et al. [24] and Hong and Liu [25], the principal component analysis along with a multivariate discriminate analysis by Tao [26], the CRITIC method combined with the grey relational analysis and TOPSIS method by Liu and Lu [27], the TOPSIS method with the Delphi method and SVM algorithm by Li et al. [28], the fuzzy evaluation by Xu et al. [29], the analytic hierarchy process (AHP) with an entropy method by Shen and Jin [30], and the grey prediction model along with the theory of risk energy release by Gao and Zhang [31].
Although previous studies have great enlightenments in constructing EWS for public debt risk, there are still some research gaps at the local government debt level, especially for China's local government debt risk. The research gaps are as follows: (1) In terms of local government debt statistics, previous studies mainly construct the EWS for explicit debt of local government; and limit the research object to local government bonds [1,4,29,32], or debt from local government financial vehicles [26,27], or both the two ones [17,28,30,31,33,34]. (2) In terms of early warning indicators, previous studies only select indicators from single and static perspectives. Each local government is regarded as a single debt risk carrier. Both the contagion effect and the probability of deterioration or turnaround of debt risk are not taken into account. Although few studies (see Wang and Chen [35]) have considered the regime-switching risk of local government debt, the contagion risk of local government debt is unfortunately ignored. (3) In terms of weighting methods, previous studies tend to use subjective methods or a subjective and objective integrated methods [25,[27][28][29][30]. Decision maker needs to compare the criteria subjectively, thus reducing the objectivity of conclusions. (4) In terms of machine learning algorithms, there exists the following problems for China's case: 1) The applicability of machine learning algorithms for China's small sample data. Fundamentally, many machine learning algorithms have requirements for sample size. Taking the BP neural network model as an example, it is likely to fall into the trap of local optimum when applied to small sample.
2) How to optimize the parameters of machine learning algorithms and avoid overfitting in face of China's small sample data. Although the SVM has proved a good performance in generalizing the EWS for government debt risk by western scholars, e.g., Alaminos et al. [22] and Zahariev et al. [23], the SVM has been deliberately simplified when applied to China's case. Some scholars

Representatives
Model/ Technique Contribution Limitation

Ciarlone and Trebeschi [14]
Multinomial logistic regression It discovers which macroeconomic factors determine the debt crisis, and obtains crisis signals from the default parameters.
Due to the limited sample size, the two aspects of EWS for debt crisis are not well integrated.

Fuertes and Kalotychou [15];
It controls for country, regional and time heterogeneity to improve the forecast power of sovereign default models.
The EWS is based on binary dependent variable model, with an inherent endogeneity problem.

Fournier and Bétin [16]
It complements the literature on debt limits in advanced economies.
The sample does not include low-income countries.
The Delphi method (a subjective weighting method) affects the objectivity of weights of indicators.

Kamra [18]
Parametric proportional hazards regression+ Conventional logistic regression+ Bayesian model averaging It expands the data coverage to emerging markets, and uses alternate empirical techniques to select early warning indicators.
The EWS is too complex. The warning speed and convenience need to be improved.

Savona and Vezzoli [19]
Regression-tree based approach + Signal model It reports different risk indicators used for four main purposes, and helps local governments to deal with exogenous financial crises.
Any claim of generalizability beyond the reviewed material needs to be further verified.

Dawood et al. [20]
Dynamic-recursive early warning model It proposes a different specification of crisis variable, and develops a more powerful dynamic-recursive forecasting technique.
The sample data of developing countries are very limited.

Antulov-Fantulin et al. [21]
Machine learning algorithms It enriches the literature on EWS by using the recent machine learning algorithms.
More explanations about the role which social, human and cultural capital play in the municipal defaulting process need to be elaborated.

Alaminos et al. [22]
It shows the superiority of computational techniques over statistics in terms of the precision in the EWS for sovereign debt risk.
The country strength model and financial strength model need to be modified to increase the generalization of EWS.

Zahariev et al. [23]
It constructs a good visible EWS by using Python software.
The selection of early warning indicators needs to be elaborated.

Shi et al. [24]
BP neural network It fills the research gap in the EWS for China's government financial risk.
The BP neural network has slow convergence speed and low training efficiency, and tends to fall into the trap of local optimum due to the small sample size.

Hong and Liu [25]
It adopts a subjective and objective integrated method in the indicator weighting, and constructs a nonlinear simulation EWS for China's local government debt risk.
The subjectivity of AHP a method.

Tao [26]
Principal component analysis + Multivariate discriminate analysis It assesses the default risk and debt-bearing boundaries of local governments in China.
The sample is limited to implicit debt of local government.

Liu and Lu [27]
CRITIC method + Grey relational analysis+ TOPSIS method It uses the CRITIC b method to weight indicators. The weight derived from TOPSIS c method leads to the problem of pseudo weight. The sample is limited to implicit debt of local government.

Li et al. [28]
TOPSIS method+ Delphi method + SVM algorithm It introduces machine learning algorithms into the construction of EWS for China's local government debt risk.
The subjectivity of TOPSIS and Delphi method.
The penalty parameter C is ignored in the parameter optimization of SVM.

Xu et al. [29]
Fuzzy evaluation It uses the AHP and fuzzy comprehensive evaluation to construct the EWS for China's H province.
The subjectivity of AHP and fuzzy evaluation method. The implicit debt is ignored in the sample.

Shen and Jin [30]
AHP method+ Entropy method It selects indicators based on the inherent logic of local government debt risks, and constructs a multilevel EWS.
The subjectivity of AHP method.

Gao and Zhang [31]
Grey prediction model+ Theory of risk energy release It constructs a dynamic model of local government debt risk assessment by Vensim PLE software, and obtains a grey prediction model.
The program design is too complicated, with slow operation speed.
(Continued ) only optimize the parameter gamma in Gaussian kernel function of SVM, without considering the penalty parameter C (see Li et al. [28]). 3) Application reason is not convincing. Different from the real cases of local government bankruptcy in western countries, China has not witnessed a case of local government bankruptcy due to some policy factors. Therefore, Chinese scholars cannot define the local governments as default or non-default by machine learning classifiers like western scholars. From the few attempts made by Hong and Liu [25] and Li et al. [28], the machine learning algorithms is adopted to extract experts' subjective experience-because the subjective weighting methods such as AHP and Delphi are widely used in their studies. However, whether the experts' experience is valuable and professional enough has not been demonstrated, making their conclusions less interpretable and objective.
In view of the above research gaps, this paper intends to construct an EWS from multiple perspectives of local government debt risk, and debugs a machine learning classifier to generalize the EWS. We believe that our EWS is of some reference significance to China and other countries with similar risk problems. That is, China and other countries can learn from the ideas of this paper to construct a comprehensive EWS to quickly grasp the local government debt risk signals.

Conceptual framework
There are four types of local government debt risks in this paper: The individual risk of local government debt: due to the fact that local government as a subnational administrative agency cannot easily go bankruptcy, the individual risk of local government debt can be defined as "too big to fail" of financial risk. Although the Ministry of Finance and the People's Bank of China are exploring the local government bankruptcy, financial institutions and the public still adhere to the firm beliefs in the government. Then, the moral hazard of "too big to fail" encourages an aggressive expansion of local government debt, thus accumulating a great amount of local government debt risk. Nowadays, scholars have attached great importance to the individual risk of local government debt, and established various EWSs to predict the risk. Representatives are Wijayanti and Rachmanira [6], Turay et al. [17], Hong and Liu [25], Li et al. [28], Gao and Zhang [31], Jin et al. [36] and Zhang [37].
The contagion risk of local government debt: it can be defined as "too tightly correlated to fail" of financial risk. Referring to the empirical conclusions of Li et al. [9], there exists spatial correlations and spillover effect of local government debt risks in China's 30 provinces. Such spatial characteristics are likely to cause contagion of local government debt risks between provinces.
The static risk of local government debt: it refers to how risky the local government debt is when it maintains its original risk state. As the local government debt risk is likely to transit from low-risk state to high-risk state in some China's provinces (see Wang and Chen [35]), it is necessary to pay attention to the state changes of local government debt risk. Obviously, if the local government debt of a province is in high-risk state for a long time, the local government debt risk of that province is greater; conversely, if the local government debt of a province is in lowrisk state for a long time, the local government debt risk of that province is smaller. The dynamic risk of local government debt: it describes the regime-switching risk of local government debt when it switches between low-risk state and high-risk state. If the local government debt of a province switches from a low-risk state to a high-risk state, the province is experiencing a deterioration of local government debt risk; conversely, if the local government debt of a province switches from a high-risk state to a low-risk state, the province is experiencing an easing of local government debt risk.

Hypothesis
Risks are not single and static, but contagious and dynamic. The contagion effect and dynamics can be described by empirical models. Regarding the contagion risk, Bianchi et al. [38] take the network structure perspective and use the standard eigenvector centrality to model contagion in financial market; Anagnostou et al. [39] incorporate contagion in portfolio credit risk model by using network theory; Berloco et al. [40] use the network model to capture firms' fragility to shocks. Regarding the dynamic risk, Jutasompakorn et al. [41] identify the banking crisis dates via the MS-AR model; Xaba et al. [42] explore the performance of MS-AR model to forecast the quarterly exchange rate of South Africa; Makatjane and Kagiso [43] realize a dynamic early warning of the Johannesburg stock exchange all-share index through a two regime MS-AR model. Referring to the above practices, this paper will apply the network model and MS-AR model to construct a comprehensive EWS local government debt risk in China. Based on the above, this paper proposes the following hypothesis: Hypothesis 1: A proxy variable that comprehensively considers the individual risk, contagion risk, static risk and dynamic risk of local government debt of each province can be obtained by network model and MS-AR model.
In risk research, there is a heated discussion on foundational issues about concepts and perspectives. The development of well-founded risk perspectives is a crucial issue to intensify the practice of risk analysis. Fundamentally, different risk perspectives may lead to different early warning results [44]. Taking China's data as sample, Shen et al. [30] construct an EWS from the individual risk perspective, finding that Qinghai, Hunan, Guizhou, Heilongjiang and Jilin have greater local government debt risk. Wang and Chen [35] study from the static and dynamic perspectives, and conclude that Jiangsu, Hebei, Shandong, Chongqing, Xinjiang and Henan provinces have greater local government debt risk. Different from previous studies, this paper has four risk perspectives, i.e., individual risk, contagion risk, static risk and dynamic risk. Such multiple risk perspectives may lead to different early warning results. Thus, this paper proposes the following hypothesis: Hypothesis 2: Different risk perspectives may lead to different list of provinces that require early warnings.
Moreover, due to the numerous advances in software, recent works increasingly use the machine learning algorithms to finalize prediction. Scholars have showed how machine learning algorithms outperform the traditional econometric models, and recognized the superiority of machine learning algorithms in constructing EWS (see Antulov-Fantulin et al. [21], Alaminos et al. [22] and Zahariev et al. [23]). Based on this, this paper will use machine learning algorithms to improve the generalization performance of our EWS. Therefore, the following hypothesis is proposed: Hypothesis 3: An optimal machine learning classifier can be debugged to generalize our EWS.

Research design
The research design is shown in Fig 1. Firstly, early warning indicators of individual risk and contagion risk of local government debt will be selected, respectively. Secondly, both the CRITIC method and coefficient of variation method will be used to calculate the proxy variable Ⅰ of each province. The proxy variable Ⅰ reflects the individual risk and contagion risk of local government debt. Thirdly, based on the proxy variable Ⅰ, the MS-AR model will be applied to investigate the regime-switching risk of local government debt for each province. Then, the proxy variable Ⅱ which comprehensively considers the individual risk, contagion risk, static risk and dynamic risk of local government debt can be obtained by the MS-AR estimation. Finally, machine learning algorithms will be introduced to debug an optimal classifier to generalize our EWS.

Methodology
Due to the availability of data, this paper chooses a provincial-level unit as the research area [45][46][47][48]. Additionally, in order to cover as many types of local government debt as possible, this paper classifies the local government debt into explicit debt and implicit debt according to Polackova's definitions [49]. Then, referring to the statistical caliber of Mao and Huang [50] and Wang [51], the explicit debt is composed of local government bonds and debt re-loans, and the implicit debt is composed of urban investment bonds. Table 2. The data come from the Wind

PLOS ONE
Several explorations on how to construct an early warning system for local government debt risk in China For negative values of X3: if a province's current new debt > 0 with current fiscal revenue increment <0, the province will be classified in the severe warning area; if a province's current new debt < 0 with current fiscal revenue increment > 0, the province will be classified in the mild warning area. b For negative values of X4: if a province's current GDP growth > 0 with current debt balance growth < 0, the province will be classified in the mild warning area; if a province's current GDP growth < 0 with current debt balance growth > 0, the province will be classified in the severe warning area.
c For negative values of X11: if a province's current fiscal revenue growth > 0 with current fiscal expenditure growth < 0, the province will be classified as in the mild warning area; if a province's current fiscal revenue growth < 0 with current fiscal expenditure growth > 0, the province will be classified in the severe warning area. https://doi.org/10.1371/journal.pone.0263391.t002 database and the empirical results of Li et al. [9]. However, Li et al. [9] do not provide the statistics of principal repayment and interest payable for local government bonds, urban investment bonds and debt re-loans. This paper supplements the relevant calculation as follows: For the local government bonds: (1) The principal repayment of local government bonds of one province at the end of year t equals the sum of matured local government bonds of that province in the t year. (2) Since each bond has its own ways of paying interest, terms in the t year, interest rates and amounts; the interest of local government bonds should be calculated one bond by one bond for each province. Then, these interests will be aggregated as the total interest of that province at the end of year t. The calculation are conducted by Stata 14.0.
For the urban investment bonds: the debt principal repayment and interest are calculated in the same way as those of local government bonds.
For the debt re-loans: referring to the Notice (No. 479 [1999]) of "the Ministry of Finance on Certain Issues Concerning the Repayment of Principals and Interests of Debt Re-loans": (1) The principal repayment of debt re-loans of one coastal developed province at the end of year t is measured by the sum of debt re-loans obtained in the t−6 year. The total interest of that province in the t year equals the sum of debt re-loan balance and principal repayment multiplied by an annual interest rate of 5.5%. (2) The principal repayment of debt re-loans of one central or western province at the end of year t is measured by the sum of debt re-loans obtained in the t−10 year. The interest payable of that province in the t year equals the sum of debt re-loan balance and principal repayment multiplied by an annual interest rate of 5%.

Early warning indicators of contagion risk.
According to social network theory, centrality can be used to express the importance of edges and nodes. The degree of influence of the node on other nodes can be estimated by centrality analysis [53].
Based on this, the centrality indicators from Li et al. [9] will be used to measure the contagion risk of local government debt, as shown in Table 3. However, the mild warning, moderate warning and severe warning are not distinguished in Table 3. In order to determine the warning degree of each province, these indicators will be divided into three equal shares in terms of the year and value. The division rules are: (1) For positive indicators, those provinces ranking in the top 1/3 will be classified in the severe warning area, those in the bottom 1/3 will be in the mild warning area, and those in the middle 1/3 will be in the moderate warning area; (2) For negative indicators, the opposite is true.

Individual risk index and contagion risk index.
There exists various weighting methods, i.e., the subjective weighting represented by AHP, the objective weighting represented by CRITIC method and coefficient of variation method, the subjective-objective integrated weighting, and the newly-emerging weighting represented by BWM (Best-Worst Method) [56], SWARA (Simultaneous Evaluation of Criteria and Alternatives) [57], SECA (Simultaneous Evaluation of Criteria and Alternatives) [58] and MEREC (Method Based on the Removal Effects of Criteria) [59]. Considering that subjective weighting incorporates subjective judgments of decision maker while some objective weighting needs fussy work with much calculation, this paper adopts the CRITIC method and coefficient of variation method to give objective weights to indicators.
The CRITIC method extracts information from a decision matrix to determine the objective weights of indicators [60]. The quantification between indicator j and other contradicting indicators is represented by X n i¼1 ð1 À r ij Þ. Where r ij refers to the correlation coefficient between indicator i and indicator j, σ j serves as the standard deviation of indicator j, and n is defined as the number of indicators. Theoretically, the objective weight of each indicator results from the comparison between the importance and contradiction of indicators. Therefore, C j , as the amount of information contained in the indicator j, can be expressed as Eq (1): In Eq (1), the larger the C j is, the more information the indicator j contains, and the greater weight will be given to indicator j. Then, W j , as the objective weight of indicator j, can be where i is a given node. If there is a tie between node i and node j, then a(i,j) =1.
Degree centrality equals the sum of ties of a given node. The larger the degree centrality of a province is, the more ties the province shares with other provinces in the spatial network of local government debt risk; vice versa.

Y1
Closeness Centrality C C ðiÞ ¼ ½ P N j¼1 dði; jÞ� À 1 . The farness of a node is the sum of the lengths of the geodesics to every other node. The reciprocal of farness is exactly the closeness centrality. In addition, the larger the farness of a node is, the more marginal the node will be.
Closeness centrality is a measure for a given node that is not controlled by other nodes. According to Freeman [54], the node with the farthest geodesic distance from the central node controls the least information resources, power, prestige and influence. In this paper, the larger the closeness centrality of a province is, the less important the province will be in the spatial network of local government debt risk.

Betweenness
Centrality g jk ðiÞ g jk , where g jk equals the sum of geodesics linking node j and node k, and g jk (i) represents the sum of geodesics which pass through node i.
Betweenness centrality is a measure of the number of times a node occurs on a geodesic, reflecting the degree to which a node controls the flow of information between other nodes. The larger the betweenness centrality of a province is, the stronger ability the province has to promote the spillover effect of local government debt risk.

Eigenvector
Centrality where N(g) refers to a given network; g stands for the sum of nodes. A = a i,j is the adjacency matrix corresponding to the network N(g); if node i is adjacent to node j, then a i,j = 1, otherwise a i,j = 0. M(i) represents the set of adjacent nodes of node i. λ illustrates the maximum eigenvalue of A = a i,j by the Perron-Frobenius theorem.
According to the eigenvector centrality, the importance of a node is not only determined by the sum of its adjacent nodes, but also by the importance of the adjacent nodes. Thus, the larger the eigenvector centrality of a province is, the closer the province is to the province of source of contagion.

Bonacich's Power
Given an adjacency matrix A, the centrality of node i (described by ci) is given by ci = ∑A ij (α+βcj), where α and β are parameters. The centrality of each node is therefore determined by the centrality of the node to which it is connected. In addition, the value of α is used to normalize the measure; and the value of β is defined as an attenuation factor, which gives the amount of dependence of each node's centrality on the centralities of the nodes it is adjacent to.
According to Bonacich [55], the greater the Bonacich power of a node, the stronger the node's ability to possess resources through connections with other nodes. In this paper, the greater Bonacich power of a province has, the stronger the contagion ability the province has in the spatial network of local government debt risk.

Average Reciprocal Distance
. The average reciprocal distance is the reciprocal of the average geodesic distance of d(i,j). Where d(i,j) refers to the shortest geodesic distance between node i and all other reachable nodes, and N represents the total amount of reachable nodes of node i in the network.
The average reciprocal distance emphasizes the value of nodes in the network. In this paper, the larger the ARD(i) of a province is, the larger probability the province has to accept or transfer the local government debt risks. expressed as Eq (2): Before the application of CRITIC method, the early warning indicators in Tables 2 and 3 require data processing in the SPSSAU software. Considering that the indicators in Tables 2 and  3 are both positive and negative, the negative ones should take the "reciprocal" to be consistent with the positive ones. In addition, since the dimensions of indicators are different and not comparable, all of the indicators should be "standardized". After the data processing, the CRITIC method can be used to obtain the weights, as shown in the last column of Tables 2 and 3.
After data processing, the values of 0, 1 and 2 will be respectively assigned to the indicators in the mild warning area, moderate warning area and severe warning area. Then, the values will be multiplied by corresponding CRITIC method weights, to obtain the individual risk index and contagion risk index for each province. Due to the limited space, the two indexes are presented in Appendices A and B in S1 Appendix (see all appendices at https://doi.org/10.7910/DVN/1EOZXG).
4.1.4. The proxy variable I: Individual and contagion risk. "Too big to fail" and "too tightly correlated to fail" are equally important in risk prevention. Regulatory authorities should be vigilant against the two risks at the same time. Therefore, it is necessary to calculate the proxy variable Ⅰ to combine the above individual risk index and contagion risk index.
The coefficient of variation method is suitable for researchers to conduct a comprehensive evaluation of the research object. The coefficient of variation of the i-th index (i = 1,2,. . .,n) is denoted as CV i , shown in Eq (3): In Eq (3), σ i represents the standard deviation of the i-th index, and � x i serves as the mean value of the i-th index. The weight W i of the i-th index can be obtained by Eq (4): Before the application of coefficient of variation method, the individual risk index and contagion risk index need to be "standardized", in order to solve the dimensional problem. In this case, the "interval" method is adopted by SPSSAU.
After data processing, the weights of individual risk index and contagion risk index are calculated to be 41.27% and 58.73%, respectively. Based on the two weights, the proxy variable Ⅰ can be obtained, as shown in Table 4. It can be seen from Table 4 that provinces with the larger proxy variable Ⅰ are Hunan, Jiangsu, Hubei, Shandong and Tianjin provinces (These provinces have repeatedly ranked among the top 5 provinces in terms of the proxy variable Ⅰ during the sample period. Specifically, Hunan appeared 9 times, Jiangsu 7 times, Hubei 6 times, Shandong and Tianjin each 5 times).
Consequently, on the basis of Table 4, Appendices A and B in S1 Appendix, China's 30 provinces can be classified into four types without regard for static risk and dynamic risk of local government debt: (1) The first type are the provinces with high-risk, i.e., the provinces ranking the top 10 in both the individual risk index and contagion risk index. It is found from Fig 2 that the number distribution of the four types is relatively balanced. Regulatory authorities can implement different policies based on the Fig 2. (1) For the first type of provinces, regulatory authorities should try to reduce their local government debt risk as much as possible. (2) For the fourth type of provinces, regulatory authorities can temporarily not list them as the risky provinces; however, they should pay more attention to the dynamic trend of debt risk. (3) For the second and third types of provinces, which are the provinces inconsistent with the severity of individual risk and contagion risk, regulatory authorities do not need to balance these two risks at the same time, thus reducing the difficulty of local government debt management. However, regulatory authorities can still take some risk management measures. For the second type of provinces, regulatory authorities should conduct one-to-one debt governance; for the third type of provinces, regulatory authorities should pay more attention to the competition among the provinces, so as to curb the spatial spillover effect of local government debt risk.

The proxy variable II
It can be seen from Table 4 that the proxy variable Ⅰ of each province is not static, but has fluctuations and structural changes during the sample period. Therefore, it is imperative to conduct a dynamic monitoring of China's local government debt risk. Considering that the MS-AR model can effectively capture different regimes of time series data [10][11][12], this paper applies the MS-AR model to explore the transition probabilities and expected durations of the proxy variable Ⅰ for each province.

The MS-AR model. The MS-AR model is initially introduced by Hamilton [10]
to study the regime-switching of market-related time series. The model is shown in Eq (5): In Eq (5), {x t } represents a time series, which is characterized by the p-order auto-regressive of M; S t refers to a state variable, which is a first-order Markov process with a state space of J = (1,2,. . .,M); y ¼ ðC S t ; φ S t k ; s S t ; P ij jS t 2 J; k ¼ 1; 2; :::; p; i; j 2 JÞ denotes a set of undetermined parameters, where the C S t , φ S t k and s S t respectively stand for the intercept, coefficient and residual mean square error of the autoregressive model under the state of S t ; {ε t ,t = 1,2,. . .,T} is defined as a white noise sequence with zero mean following a normal distribution; and P ij describes the one-step transition probability of regime i switching to regime j.
Identifying the P ij requires selecting the number of states, i.e., the number of regimes. After many attempts, this paper finds that the two regimes of high-risk state and low-risk state is the best. Therefore, P ij has a total of 2×2 elements in this paper, with the transition probabilities matrix p shown in Eq (6): According to the conclusion of Hamilton [10], the θ in Eq (5) can be obtained by using maximum likelihood estimation. Consequently, the conditional logarithmic likelihood function of θ can be expressed as Eq (7): In Eq (7), O t−1 describes the set of all series samples ( where the probability variable ξ i,t−1 and the state-mixture density function η j,t are defined as Eq (8) and Eq (9), respectively: Obviously, Eq (7) can be derived from Eq (8), Eq (9) and Bayes formula. Thus, the θ can be obtained by iterative estimation.
Before the application of MS-AR model, the following definitions should be given: The regime-switching probability P ij (i6 ¼j): the probability of being in regime j at t depends on the regime i at t−1, and the probability P ij (i6 ¼j) of regime i switching to regime j is defined as Eq (10). Then, the regime-switching matrix P can be derived from P ij (i6 ¼j).
The steady-state probability P ij (i6 ¼j): different from the regime switching probability P ij (i6 ¼j), the steady-state probability P ij (i = j) describes the probability that the regime remains unchanged. Then, the steady-state probabilities matrix Q can be derived from P ij (i = j). Since QP = Q, i.e., P T Q T = Q T , Q can be obtained by finding the eigenvector of P T with eigenvalue 1.
The duration D i : D i is denoted as the duration that the regime i can last, and satisfies Eq (11):

The MS-AR estimation of the proxy variable I.
Due to the limited space, the MS-AR estimation is detailed in Appendix C in S1 Appendix, and Chongqing is taken as an example to interpret Appendix C in S1 Appendix, as shown in Table 5.
It can be seen from Table 5 that: (1) Both the coefficient of high-risk state (1.223327) and the coefficient of low-risk state (0.997505) have passed the 1% significance test, showing that the MS-AR model of Chongqing province is quite reasonable. (2) Although the coefficient of high-risk state is relatively large, the steady-state probability of high-risk state (0.60986) is less than the steady-state probability of low-risk state (0.846999), indicating that Chongqing's local government debt has a relatively strong stability in the low-risk state. (3) The probability of transiting from high-risk state to low-risk state is 0.39014, while the probability from low-risk state to high-risk state is 0.153001, indicating that Chongqing's local government debt changes from a high-risk state to a low-risk state more frequently. In other words, Chongqing's local government debt risk is more likely to ease than worsen in the future. (4) The duration of high-risk state (2.563184) is much less than that of low-risk state (6.535899), reflecting that Chongqing's local government debt lasts longer in the low-risk state.

The proxy variable II: Individual, contagion, static and dynamic risk.
Based on the Table 4, Appendices A-C in S1 Appendix, the proxy variable Ⅱ which comprehensively considers the individual risk, contagion risk, static risk and dynamic risk of local government debt can be obtained by using the coefficient of variation method.
Referring to the practice of Wang and Chen [35], the proxy variable Ⅱ includes the following seven indicators: the average of the proxy variable Ⅰ in Table 4; the average of the individual risk index in Appendix A in S1 Appendix; the average of the contagion risk index in Appendix B in S1 Appendix; the steady-state probability of low-risk state, steady-state probability of high-risk state, duration of low-risk state and duration of high-risk state in Appendix C in S1 Appendix.
Before the application of coefficient of variation method, some data processing is needed. Concretely, the above negative indicators should take the "reciprocal" to be consistent with the positive ones. In addition, the seven indicators need to be standardized by the "interval" method, so as to solve the dimension problem.
After data processing, the weights of the seven indicators can be calculated. The results are shown in Table 6. Then, the indicators are multiplied by the corresponding weights to obtain the proxy variable Ⅱ of each province. The results are shown in Table 7. In Table 7, the larger the proxy variable Ⅱ of a province, the higher the local government debt risk in the province; and the province with the largest proxy variable Ⅱ ranks the 1 st , indicating that the province is the most risky province of local government debt; et cetera.
In order to show Table 7 more intuitively, China's 30 provinces have been divided into five equal shares according to their risk ranking, as shown in Fig 3. In Fig 3, the five equal shares correspond to highest-risk, higher-risk, medium-risk, lower-risk and lowest-risk provinces, respectively. It should be emphasized that Fig 3 is the final classification result of this paper, which comprehensively considers the individual risk, contagion risk, static risk and dynamic risk of local government debt. Consequently, Hypothesis 1 is verified. A comparison between Figs 2 and 3 finds that different risk perspectives lead to different list of provinces that require early warnings: (1) From the perspective of individual risk and contagion risk in Fig 2, regulatory authorities should send different early warning signals to different provinces: high-risk signals to Hunan, Anhui, Tianjin, Shananxi, Hubei, Jiangsu and Guizhou provinces; vulnerability signals to Qinghai, Guangxi, Gansu, Xinjiang, Heilongjiang, Ningxia, Hainan and Yunnan provinces; and high-correlation risk signals to Shandong, Henan, Hebei, Zhejiang, Shanxi, Liaoning and Jilin provinces. (2) From the perspective of individual risk, contagion risk, static risk and dynamic risk in Fig 3, regulatory authorities should send highest-risk signals to Hunan, Jiangsu, Hubei, Henan, Shandong and Shananxi provinces; and higher-risk signals to Anhui, Guangdong, Guizhou, Guangxi, Sichuan and Beijing provinces. Therefore, the above different list of risky provinces confirms the Hypothesis 2. The larger the average of the proxy variable I, the higher the individual risk and contagion risk of local government debt.

Average of Individual Risk Index
The larger the average of individual risk index, the higher the individual risk of local government debt.

Average of Contagion Risk Index
The larger the average of contagion risk index, the higher the contagion risk of local government debt.

Steady-state Probability of Low-risk State
The larger the steady-state probability of low-risk state, the larger the probability that the local government debt risk will remain in a low-risk state, thus the lower the local government debt risk is.

Steady-state Probability of High-risk State
The larger the steady-state probability of high-risk state, the larger the probability that the local government debt risk will remain in a high-risk state, thus the higher the local government debt risk is.

Duration of Lowrisk State
The larger the duration of low-risk state, the longer the local government debt risk will remain in a low-risk state, thus the lower the local government debt risk is.

Z6
Duration of Highrisk State The larger the duration of high-risk state, the longer the local government debt risk will remain in a high-risk state, thus the higher the local government debt risk is.

Sensitivity analysis.
A sensitivity analysis is conducted based on changing attribute weights to show the stability of results. To be specific, individual risk and contagion risk index are respectively obtained by the coefficient of variation method. Then, the proxy variable Ⅰ is obtained by the independent weight method, and the proxy variable Ⅱ is obtained by the CRITIC method (see S1 Dataset for more details). The results are shown in Appendix D in S1 Appendix.
Comparing the Appendix D in S1 Appendix and Table 7, it is found that the list of provinces receiving early warning signals in Appendix D in S1 Appendix is consistent with that of Table 7 except for Inner Mongolia and Hebei. The comparison demonstrates that changing attribute weights does not have a significant impact on the robustness of the results of this paper.

Generalization based on machine learning algorithms
Although the EWS designed in this paper considers the local government debt risk from multiple perspectives, the calculation is complex and cumbersome. It is necessary to improve the generalization ability of our EWS. Therefore, the machine learning algorithms are adopted to debug a classifier to realize a rapid early warning of our EWS. Machine learning has three algorithms, i.e., supervised learning, unsupervised learning and reinforcement learning. Considering that only provinces with the highest and higher risks will receive early warning signals in this paper, the provinces with the highest and higher risks can thereupon be labeled as 1, and the remaining provinces are labeled as 0. Then, labels 0 and label 1 are used as output variables, and the X 1~X11 and Y 1~Y6 are used as input variables. Obviously, it is a binary classification problem, which is most suitable for using the supervised learning.
Among many algorithms of supervised learning, k-neighbors, linear function, neural network, random forest, gradient boosting machine and support vector machine have been widely used by scholars. Considering China's small sample data, it is necessary to exclude the inappropriate supervised learning algorithms first. According to Andreas and Sarah [61], kneighbors has slow prediction and cannot deal with a dataset of many feature vectors; linear classifier is limited in low dimensions in face of small sample data; neural network takes a long training time and is easily trapped in the local optimum in the small sample data. However, the random forest, gradient boosting machine and support vector machine have advantages of statistical, computational and representative, without high requirements for sample size [21][22][23]. Therefore, this paper selects random forest, gradient boosting machine and support vector machine, as detailed in Table 8. Compared to previous studies, (1) This paper fully considers the characteristics of small sample data of local government debt in China. (2) This paper does not deliberately omit any parameter in the tuning process. (3) This paper uses the machine learning algorithms to generalize the EWS rather than extract experience from experts for the subjective weights. That is, there is no "black box" in this paper, thus maintaining the objectivity of conclusions.

Random forest
Random forest (RF) uses the bootstrap resampling to repeatedly select n samples (usually 2/3) from the original training set T to generate a new training set, each of which is used to train a tree independently. Then, the n decision trees form a forest, in which each tree has the same distribution. Theoretically, the classification error depends on the classification ability of each tree and the correlation between them. In addition, the unselected sample is called out of bag data (OOB), whose error is an unbiased estimate that can be used to validate the performance of the model to prevent overfitting [62].
The generalization error P � of RF is defined as Eq (12): where ρ refers to the average of the correlation of decision trees; S represents the average strength of decision trees. From Eq (12), it is necessary to reduce the correlation between decision trees or increase the strength of decision trees, so as to improve the generalization of RF. To achieve this goal, random disturbance term of feature variables is introduced, resulting in different split nodes of each tree. Finally, the training set and each node variable are randomly generated in the forest.
In the process of growing each classification tree, the splitting of each node is determined by the "purity" of split sample. The "purity" has two criteria, i.e., Gini index and entropy, as shown in Eq (13) and Eq (14), respectively: where p i is denoted as the probability of classification i. The smaller the Gini index or entropy, the higher the purity of the sample, thus the better the performance of tree splits.  Tables 2 and 3.

Test Set
The data of year 2018 As shown in Table 7.

Training Set
The data of year 2016 and 2017 Taking 2010~2016 and 2010~2017 as the sample period, the proxy variable II and risk ranking of local government debt risk in 2016 and 2017 are calculated respectively. The results are shown in Appendices E and F in S1 Appendix.

Gradient boosting machine
Different from RF, gradient boosting machine (GBM) uses a continuous method to grow decision trees. The algorithm of GBM is as follows [63]: Suppose that A = (x 1 ,x 2 ,. . .,x n ) serves as the independent variable, B = {y i } refers to the dependent variable and i = 1,2,. . .,n. For a given dataset, the variables in A need to be mapped to B through a mapping function f � (x). In addition, the difference between the mapping function and the real function is represented by the loss function L(y,f(x)). Obviously, the prediction model should minimize the loss function L(y,f(x)) and initialize the mapping function f � (x) as Eq (15): According to the GBM, the direction of L(y,f(x)) declines along with the gradient direction, due to the fact that L(y,f(x)) changes most significantly in that direction. Thereupon, the negative gradient value of L(y,f(x)) is approximately to the residual and can be defined as Eq (16): The pseudo-residual derived from Eq (16) should be matched with a base classifier g m (x), which includes various parameters and is trained by the training set fðx i ; r im Þg n i¼1 . Then, the δ m can be obtained by the following optimization: The residual coefficient y m is derived from Eq (17): After iteration of m = 1,2,. . .,M, the optimized prediction function can be obtained, as shown in Eq (18): To sum up, the GBM mainly includes the following parameters: learning rate, tree depth and number of iterations.

Support vector machine
Support vector machine (SVM) attempts to construct an optimal separation hyperplane to classify the positive data and negative data [64]. It starts with mapping the sample data from original space to high-dimensional space through the nonlinear mapping φ(�), and optimizes the classification decision function shown in Eq (19). For data points distributed between the classification hyperplane and the support vector, a relaxation variable should be given, with some penalty imposed for the wrong classification. In addition, the optimal hyperplane should satisfy the constraints in Eq (20). Then, the decision function obtained can be expressed as Eq (21), along with the Gaussian radial basis function (RBF) as shown in Eq (22): Where x i refers to the input sample; y i represents the expected output vector; ω serves as the weight vector; C is denoted as the penalty parameter and C2(0,+1); ξ i is defined as the slack variable and ξ i �0,(i = 1,2,. . .); b describes the bias vector; K(x i ,x) is exactly the kernel function, i.e., Gaussian radial basis function; and a i means the Lagrange multiplier. In brief, ω and b determine the position of optimal separation hyperplane, and γ and C are the parameters that should be optimized.

Performance of RF/GBM/SVM classifiers
This paper uses the 5-fold cross validation and learning curve to debug the RF/GBM/SVM classifiers. The results are shown in Table 9. It can be seen from Table 9 that the roc_auc scores of RF/GBM/SVM classifiers all approach 90%, which are high level in small sample data and exceed their previous roc_auc scores. In addition, all of the three classifiers run very fast. Therefore, regulatory authorities only need to input the original data of X 1~X11 and Y 1~Y6 in 2019 and later years into the classifiers, and then the output labeled 1 or 0 will be given by our EWS. In order to present the performance of RF/GBM/SVM classifiers more intuitively, this paper provides the ROC curves of the three classifiers, as shown in Fig 4. Obviously, the inflection points of the three ROC curves are all located in the upper left corner, verifying the good performance of the three classifiers again.
Furthermore, the prediction results for the 2018 test set are listed in Table 10. Comparing the true labels with the predicted labels, it is found that both the GBM classifier and SVM classifier correctly warn 26 provinces; while the RF classifier only correctly warns 25 provinces. Therefore, the GBM classifier and SVM classifier are the best-performing classifiers in terms of their high prediction accuracy. However, the GBM classifier is inferior to SVM classifier in terms of the roc-auc score, as shown in Table 9. Taken together, the SVM classifier has a better fitting effect on our test set. Thus, the Hypothesis 3 is verified.

Conclusions and policy recommendations
According to the spatial correlations and spillover effect of China's local government debt risk, this paper constructs a scientific and comprehensive early warning system (EWS) for local government debt risk in China. The data of China's 30 provinces over the period of 2010 to 2018 are taken as a sample data. Through the selection of early warning indicators for individual risk and contagion risk of local government debt, the proxy variable I which combines the two risks are obtained by the CRITIC method and coefficient of variation method for each province. Then, based on the proxy variable Ⅰ, the proxy variable Ⅱ which comprehensively considers the individual, contagion, static and dynamic risk of local government debt are estimated by the MS-AR model and coefficient of variation method. Finally, machine learning algorithms are applied to generalize our EWS. The results show that: (1) From different risk perspectives, the list of provinces that require early warning is different. Specifically, from the individual and contagion risk perspectives, Hunan, Anhui, Tianjin, Shananxi, Hubei, Jiangsu and Guizhou provinces receive high-risk signals; Qinghai, Guangxi, Gansu, Xinjiang, Heilongjiang, Ningxia, Hainan and Yunnan provinces receive vulnerability signals; while Shandong, Henan, Hebei, Zhejiang, Shanxi, Liaoning and Jilin provinces receive high-correlation risk signals. From the individual, contagion, static and dynamic risk perspectives, Hunan, Jiangsu, Hubei, Henan, Shandong and Shananxi provinces receive highest-risk signals; while Anhui, Guangdong, Guizhou, Guangxi, Sichuan and Beijing provinces receive higher-risk signals. (2) The SVM classifier can accurately send early warning signals to highest-risk and higher-risk provinces. Once the original data of X 1~X11 and Y 1~Y6 are input to the SVM classifier, high-accuracy signals can be output, thus realizing a rapid and comprehensive early warning of China's local government debt risk.
Based on the above conclusions, it is recommended that China's regulatory authorities broaden and innovate regulatory ideas, and strengthen the comprehensive supervision of individual risk, contagion risk, static risk and dynamic risk of local government debt.
However, due to data limitations, only ten years of data is available for the analysis, and only the MS-AR model is suitable to investigate the trend of local government debt risk. Further research is expected to focus on the following points: (1) We will keep track of China's data to expand the sample size. (2) Considering that local government debt risk can be studied from multiple perspectives with multiple risk indicators, once the sample size is expanded, we  will build a Markov switching mixed frequency dynamic factor model to characterize the trend of local government debt risk. Then, we can analyze the dynamic transformation path between different risk regimes, thus forming a more in-depth understanding of local government debt risk.