Exploring Brexit with dynamic spatial panel models: some possible outcomes for employment across the EU regions

Starting with a reduced form derived from standard urban economics theory, this paper estimates the possible job-shortfall across UK and EU regions using a time-space dynamic panel data model with a spatial moving average random effects structure of the disturbances. The paper provides a logical rational for the presence of spatial and temporal dependencies involving the endogenous variable, leading to estimates based on a state-of-the-art dynamic spatial generalized moments estimator proposed by Baltagi et al. (Reg Sci Urban Econ, 2018. https://doi.org/10.1016/j.regsciurbeco.2018.04.013). Given reliable interregional trade estimates, the simulations are based on a linear predictor which utilizes different regional interdependency matrices according to assumptions about interregional trade post-Brexit. The results indicate that heightened barriers to trade will evidently cause job-shortfalls both in the UK and across the EU, but it is also shown that there is a considerable amount of asymmetry in the outcomes across regions and sectors.


Introduction
The paper aims to simulate potential outcomes for employment across 255 EU regions as a result of the impending and, at the time of writing, probable departure of the UK from the EU, commonly known as Brexit. Most analysts agree that Brexit will have momentous consequences for the UK and (remaining) EU economies, but there is very minimal analysis at the regional scale, and analysis typically fails to account for interconnectivity at the regional level. Some regional impact studies have been carried out by Dhingra et al. (2017a, b), Los et al. (2017) and McCann (2018), and the current paper complements or contrasts with this research by applying a state-of-the-art dynamic spatial panel data model, in which a pan-European 1 3 approach is adopted involving the majority of EU regions 1 and all UK regions. This modelling approach is ideally suited to capturing the impact of spatial interconnectivity of the European regions and projecting the long-run consequences of Brexit across EU and UK regions, thus enabling comparison of the impact on both sides of the Channel and the Irish Sea. To this end, we use the dynamic spatial panel data model and prediction equation recently introduced into the literature by Baltagi et al. (2018) and applied in different contexts by  and Fingleton and Szumilo (2019). The model assumes that employment in a given region depends on the levels of production and investment within that region, as shown in the basic economic model which underpins the estimating equation, and it also depends on demand coming from all the regions of the EU and UK, as determined by interregional trade flows. Additionally, the model proposes that employment levels in any regional are closely linked to employment levels in the region in the previous period and on employment levels in trade-connected EU and UK regions in the previous period. Following this literature, a rational basis for the presence of spatial and temporal lags is introduced which more typically is ad hoc in the spatial econometrics literature. In addition, the model takes account of unobserved factors which also affect the level of employment. These are captured by region-specific random effects which are also spatially interdependent. An additional feature of the approach adopted is the way in which endogeneity is handled, with the application of internal instruments in the spirit of Arellano and Bond (1991), thus eliminating the often difficult search for valid external instruments.
The focus of analysis is the so-called job-shortfall which could arise due to Brexit. In other words, the intention is not to forecast what happens to the actual levels of employment in each region, which would be the predicted change in the number of jobs, but to simulate what the impact of Brexit would be assuming no consequential responses such as jobs created by new trade links formed post-Brexit, changes to the UK's competitivity and consequences for demand and employment due to changes in exchange rates and prices, changes to migration flows in and out of the UK, changes in the competitivity of firms if trade barriers are increased and regulations relaxed, and possible changes in levels of inward and outward investment and capital stock if capital relocates. Contemplating these and other possible consequences enhance uncertainty regarding what might be the actual change in the levels of employment in the UK regions, so in this paper the focus is on attempting to simulate the job-shortfall due to Brexit per se.
Stated more explicitly, the empirical analysis bases the spatial interdependence of levels of employment across different regions on how closely they are connected in terms of trade. We assume that employment levels partly reflect demand for a given region's good and services coming from the UK and EU regions. Naturally, since about 50% of the UK's trade in 2019 is with countries outside the EU, demand coming from these non-EU countries will also affect the levels of employment. For both UK and EU regions generally, we assume that the non-EU component of demand

The model
The reduced form used as a basis to simulate the Brexit effect assumes that employment partly depends on level of output, as measured by GVA (Gross Value Added), denoted by t , and (a proxy for) the level of capital within the region, based on GFCF (Gross Fixed Capital Formation), which is denoted by t . To show this, we start with the theoretical model given as Eq. (1), which is based on Eq. (30) given in "Appendix". The N by 1 vector d t is the density of employment per unit area, and t is the level of efficiency of labour at time t, so that the product d t t is the number of labour efficiency units. This is related to ̃ t , which is a measure of output in the competitive final goods and services sector in each region at time t, via the constant parameters and ̃ , thus In order to obtain total output t , it is assumed that ̃ t = t , in which is an N by 1 vector giving the share of total output in each region that is competitive final goods and services output. For simplicity of estimation, it is assumed that is constant over time. Also the employment levels are t = d t in which is the area of land in each region. Taking logs gives Rearranging (2) gives To obtain (3), I assume that labour efficiency t = t t , with more efficient labour having a higher level of output per unit of capital ̃ t . As shown below in Eq. (14), an approximation to the log level of capital is lñ t = − lnã +b ln t , hence ln t = ln t + lnã −b ln t , and from this Collecting together constants as c with an N by 1 vector of ones, and reorganizing gives in which the error term t captures the time-invariant regional heterogeneity in land and in shares , which are unobserved, as given in Eq. (12).
In the dynamic context, it is reasonable to assume that disparities in employment levels across locations will persist as an equilibrium outcome to unchanging and fundamental causes. We therefore proceed, following in parallel the exposition in Baltagi et al. (2018), to assume that (log) employment levels across regions, denoted by the N by 1 vector ln t at time t will persist at dynamically stable levels so that ln t = ln t−1 unless there are changes in the levels of t or t , or changes in common factors, interregional trade, or unobserved effects. If such a disturbance occurs at time t and is ephemeral, then ln t ≠ ln t−1 but over a subsequent period of quiescence t → T then once again we expect employment levels to converge on a new equilibrium, at which ln T = ln T−1 . Assume data are observed where ln t ≠ ln t−1 but tending to converge, so that ln t = f (ln t−1 ) , and an autoregressive process is assumed, hence in which is an N by 1 vector and is a scalar parameter. In the long-run with abs( ) < 1 , and with no subsequent disturbances, the process converges to ln T = (1− ) .
Consider next connectivity between regions in the form of a matrix * N , which is a time-invariant N by N matrix, where N is the number of regions. Spatial interdependence between regions is a feature of many different situations, and can be modelled either via an autoregressive process involving the dependent variable, or via spatial interdependence of the errors, or by both as in this paper. The problem of how to model dependence between N regions is typically resolved by application of (2) ln + ln t = ln +̃ ln t +̃ ln t −̃ ln Exploring Brexit with dynamic spatial panel models: some… an N by N matrix of constant quantitative values or weights assigned to the cells of * N which indicate the existence and importance of a link between each pair of regions. In many spatial econometric applications, connectivity of the N regions will be some function of the distance between them, be it geographical distance or some measure of economic distance. In the paper, we proxy economic distance by the level of trade between each pair of regions, more trade equals shorter economic distance. Usefully, * N provides a parsimonious parametrization for interdependence between, in this case, employment levels in different regions. As explained by LeSage and Pace (2009), once we allow for dependence relations between a set of N by N entities on a single variable, for example as represented by the N by 1 vector ln t , there are potentially N 2 −N parameters that define individual interdependence, such as the relation between ln it and ln jt , having excluded dependence of an observation on itself. This leads to an over-parametrization problem, which can be solved by imposing an a priori structure, or weights matrix * N , on the interdependence relations, thus reducing the number of parameters to be estimated from N 2 −N to one, denoted here by 1 . For purposes of interpreting parameter estimates, we normalize by dividing * N by the maximum eigenvalue of * N to give 2 N . Using this normalization, the maximum eigenvalue of N is 1, and the continuous range for which N − 1 N is non-singular is 1 min(eig) < 1 < 1 max(eig) = 1,in which 1 is a scalar spatial autoregressive parameter.
Multiplying (5) by 1 N gives Subtracting (6) from (5) leads to another logically consistent expression, in which the spatial dependence implied by (6) can be seen in (7) as an explicit cause of variation in ln t . Thus, and N is an identity matrix of order N. In order to solve Eq. (7), given appropriate parameter restrictions, Eq. (7) converges to ln T = N − N −1 N . Introducing additional covariates by writing N = (c + ) , in which c is a constant N by 1 vector, is an N by k matrix and is a k by 1 vector, gives The matrix N retains (scaled) absolute levels rather than shares as the basis of interregional connectivity, and we make the standard assumptions for a weights matrix, that it comprises fixed (non-stochastic) non-negative values with zeros on the leading diagonal and its row and column sums are uniformly bounded in absolute value, and maintain the same assumption for −1 N = N − 1 N −1 (Elhorst 2014, p. 99).
In order to maintain dynamically stable simulations, following Elhorst (2001Elhorst ( , 2014, Parent and LeSage (2011, p. 478;2012, p. 731) and Debarsy et al. (2012, p. 162), the largest characteristic root e max of −1 N N should be less than 1. This restriction ensures that employment converges to equilibrium levels Additional realism is introduced in three ways. First, the restriction that = − 1 is removed since this greatly simplifies estimation. However, we anticipate that ̂ ≈ −̂ 1̂ . Second, taking account of the variables in Eq. (4), the time-invariant matrix is replaced by time-varying matrix 3 t . Third, spatially dependent unobservables are represented by the error term t . Although the system may, depending on −1 N N , still tend towards equilibrium, equilibrium will be continuously disturbed and new equilibrium levels established as t varies. For simplicity of estimation, interregional connectivity is assumed to remain constant over the estimation period. These considerations lead to the model of employment levels 4 given in Eqs. (9,10,11,12), which is a time-space dynamic panel data model, thus Given 1t = ln t , 2t = ln t , 3t = ln t , t = 1t 2t 3t and = 1 2 3 T , Eq. (9) can be stated more explicitly as The presence of the district-invariant mean of the dependent variable t attempts to allow for the presence of observed or unobserved common factors affecting all districts at each point in time. This approach is motivated by Pesaran (2015) who provides a major treatise on the different approaches to modelling dynamic spatial panel data with common factors, and by Bailey et al. (2016) who ask, 'to what extent are the observed dependencies between different spatial units due to common factorsfor example, aggregate shocks-that affect different units rather than being the result of local interactions that generate spatial spillover effects?'. They propose the use of cross-unit averages to extract common factors, an approach that has also been Exploring Brexit with dynamic spatial panel models: some… applied by  and Fingleton and Szumilo (2019). The introduction of common factors to spatial econometric models has also been considered by among others Vega and Elhorst (2016) and Ertur and Musolesi (2017).
The disturbances t capture the effects of the spatially dependent unobserved variables, with a compound structure (12) comprising time-invariant unobserved unitspecific interregional heterogeneity represented by i with i = 1, … , N and unobserved idiosyncratic shocks represented by it ;i = 1, … , N, t = 1, … , T . These are assumed to be independent of each other and are collectively represented by u it . It is important to recognize that the i s represent the net effect of unobserved variables which in the short run can be treated as time invariant.
Most usually, the assumption is that spatial dependence is an autoregressive (SAR-RE) process, such that t = 2 N t + t . However, in this paper the assumption for the error process is a spatial moving average process (SMA-RE) as in Eq. (11), thus t = N t , where N = N − 2 N . This means that the error process is such that a shock in a region affects only neighbouring regions as defined by a row-standardized interregional contiguity matrix 5 N . In contrast, an SAR-RE process would entail shocks affecting all regions. There are two reasons for this. First, assuming SMA-RE rather than SAR-RE errors improves the predictive performance of the estimator, as described in Sect. 6. Second, SMA-RE errors might proxy for omitted spillovers, which otherwise might be captured by the spatial lags N t . This is pertinent since the presence of N t on the right-hand side of (9) could adversely affect estimation. As explained by Fingleton et al. (2017) and Baltagi et al. (2018), an SMA-RE error specification 'mitigates against the problem for instrumental variable estimation identified by Pace et al. (2012)'. In two-stage least squares (2SLS) estimation, the instrument set should comprise the 'exogenous' variables ( t ) and their spatial lags ( N t ), and kept to a low order to avoid linear dependence and retain full column rank for the matrix of instruments (1998, Kelejian andPrucha 1999). The performance of the estimation procedure could be suboptimal, as explained by Pace et al. (2012), by including N t among the set of explanatory variables. This is because with spatial lags ( N t ) among the set of regressors, then spatial lags of the spatial lags ( 2 N t , 3 N t , . . .) feature among the instruments, and this could lead to a weak instrument problem. To avoid this, SMA-RE errors are adopted as an alternative way to capture local spillovers.

Data
In estimating Eq. (10), data for employment ( t ), output as measured by Gross Value Added (GVA, t ) and capital as proxied by a function of Gross Fixed Capital Formation (GFCF, t ), both denominated in €2005m, are taken from the Cambridge Econometrics European Regional Economic database, with observations over the 5 The matrix N has ẽ −1 max = 1 , where ẽ is the vector of purely real characteristic roots of N . Often N = N , otherwise we assume that N has the same properties as N , and with the restriction that � � e 10-year period from 2001 to 2010 used to estimate the model. Data are also available for 2011−2012 , but are held back to allow out-of-sample prediction tests of the model and some rivals. t is used to reflect capital stock ̃ t , for which data are unavailable, on the basis of a simple relationship which is assumed to exist between the two variables. t measures gross net investment (acquisitions minus disposals of produced fixed assets) in fixed capital assets and so provides an indicator of changes to the stock of capital. The assumption is that t is a nonlinear function of a constant fraction ã of ̃ t so that hence As a test of the viability of this approximation, assume a standard model for the evolution of capital stock which is depreciating at a constant rate d so that in which T is a large number. One problem with (15) is that it requires the initial capital stock at time t = 1 , i.e. ̃ 1 . However, given arbitrary values for ̃ 1 and d , values for ã and b can be found, whereby (14) provides a reasonable approximation to the outcome of iterations (15). A more realistic test is provided by the existence of both (albeit experimental estimates of) capital stock 6 (Derbyshire et al. 2010) and of wellfounded GFCF data. Using the latest available data for both t and ̃ t , which is for the year, t = 2008 , and taking logs of (14), leads to a loglinear regression of lñ t on ln t which gives OLS estimates of the constant lnã −1 = 2.4546 (t ratio = 13.5628 ) and slope b = 1.0195 (t ratio = 50.8118) , with R 2 = 0.8888. The plot of ln t against lñ t shows a significant linear relationship and no evidence of outliers or of heteroscedasticity It thus appears that the model given as Eq. (13) provides a good approximation. The estimated ã = 0.0859 suggests the approximate proportion of the capital stock that is invested, and, by comparison, ∑ t ∕ ∑̃ t = 0.0686. The matrix N is based on estimated bilateral trade flows between EU NUTS2 regions. The data come from the PBL (the Netherlands Environmental Assessment Agency) 7 who developed a new methodology which is close to that of Simini et al. (2012). Details of the methodology are given in Thissen et al. (2013Thissen et al. ( , 2013a, see also Gianelle (2014). The method follows a top-down approach and therefore is consistent with the national accounts of the different countries. Given the total international exports and imports on the country level, interregional trade flows are derived using data on business travel (services) and on freight transport (goods).
6 I am grateful to Cambridge Econometrics for providing these data. 7 We are grateful to Mark Thiessen, who kindly provided the data. The data can be visualized at http:// thema sites .pbl.nl/eu-trade /index 2.html?vis=net-score s.

3
Exploring Brexit with dynamic spatial panel models: some… Additionally, exports that went to EU destination countries' final demand 8 were also included. Trade flows involving regions of non-EU countries such as Switzerland and Norway were obtained on the basis of interregional trade flows estimated by the best linear disaggregation method of Chow and Lin (1971), which was initially used to break down annual time series into quarterly series (see Abeysinghe and Lee 1998;Doran and Fingleton 2014). In this, commencing with aggregate trade values 9 between 21 EU counties, these were allocated to the NUTS2 regions. A parallel approach has been used by Polasek et al. (2010), Vidoli and Mazziotta (2010), and Fingleton et al. (2015). More details of the method are provided in "Appendix". Finally, OLS regression of the log PBL trade flows on log Chow-Lin trade flows produced parameters used to predict the missing PBL regional trade flows for Switzerland and Norway using the values for these regions obtained via the Chow-Lin approach. For estimation, the start-of-period trade flows for the year 2000 are used. This year is chosen because it is the earliest available, so it is treated as exogenous to t , t and t , for t = 2001 to 2010. Prediction is based on the 2010 trade flows supplemented in the same way by Chow-Lin data. Estimates are also given in "Appendix" Table A3 based on a N matrix constructed entirely from the Chow-Lin trade flows. These simply use great circle distances and year 2000 GVA levels, and so are also assumed to be exogenous. The comparative predictive performance of each set of estimates is discussed in Sect. 6.

Estimator for the time-space dynamic panel data model
Comprehensive overviews of spatial panel econometrics are given by Pesaran (2015, Chapters 29 and 30) and Baltagi (2013, Chapter 13) which highlight its growing importance for the applied econometrician. The estimator used in this paper, introduced by Baltagi et al. (2018), adds to the available methodology by allowing a wider range of spatial interaction effects which include the spatial lag of the temporal lag of the dependent variable W N ln t−1 , thus avoiding bias due to constraints necessary for dynamic stability and stationarity, and also by allowing spatial moving average compound error dependence rather than the usual autoregressive compound error process found in the majority of spatial econometric models. The estimator, which is applied to Eq. (9), is based on the earlier paper by Baltagi et al. (2014), which extends the approach of Arellano and Bond (1991) by the introduction of extra moments in line with the presence and availability of spatial lags (see also Bouayad-Agha and Védrine 2010). Since the estimator is described elsewhere, a simple outline sketch is provided here focussing on the treatment of regressors 8 1. Final consumption expenditure by households and non-profit organizations 2. Final consumption expenditure by government 3. Net capital formation 4. Inventory adjustment 9 They are downloadable from http://cid.econ.ucdav is.edu/data/undat a/undat a.html, see also Feenstra et al. (2005). 1 3 as predetermined rather than exogenous. 10 Hence in Eq. (10), ln t and ln t are considered to be predetermined alongside endogenous right-hand side variables ln t−1 , N ln t−1 and N ln t and ln t .
Focussing on the endogenous dependent variable ln t , the instruments include ln t lagged by two periods, and its spatial lag N ln t also lagged by two periods, so that the moments equations (16 and 17) hold assuming it is serially uncorrelated and E(Δ it , Δ it−2 ) = 0 . Thus, following Baltagi et al. (2018), with, we have in which E denotes the expectation. Also, if we were to assume exogenous rather than predetermined regressors 1 , 2 this leads to (18) for t = 3, … , T . Given that in (18) the regressors 1 , 2 are exogenous, the moments equations are satisfied including the entire set 11 , … , 1T , 21 , … , 2T , N 11 , … , N 1T and N 21 , … , N 2T regardless of time t. As explained in Baltagi et al. (2018), additional instruments can be generated via the matrix 2 N , but for simplicity these are omitted from the estimators used in the current paper.
Strict exogeneity rules out any feedback from past shocks to current values of the variable, and the need to accommodate feedback leads to the preferred estimator based on predetermined regressors (see Bond 2002;Pesaran 2015). Predetermined regressors are contemporaneously uncorrelated, so that corr( t , t ) = 0, but do depend on earlier shocks so that, for example, corr( t , t−1 ) ≠ 0 . This means that an adjustment to ln ,which embodies , at time t does not have an instantaneous effect on output and capital investment time t but takes effect at t + 1 and later. This allows an extension to the set of instruments (compared with assuming endogeneity, where all endogenous variables are lagged by two periods), by the inclusion of 1t−1 , 2t−1 , N 1t−1 , and N 2t−1 so that Given the set of instruments as in Eq. (19), these are used to obtain initial estimates of , 1 , , 1 , 2 and 3 , having first differenced the data to eliminate the time-invariant Exploring Brexit with dynamic spatial panel models: some… individual effects which are correlated with the time and space-lagged dependent variables. The resulting estimates are then used to give estimated errors which lead to estimates of the parameters of the spatial moving average error process, namely 2 , 2 and 2 using moments equations given in Fingleton (2008). Given these, preliminary one-stage consistent spatial GM estimates are obtained, followed by the two-stage Spatial GM estimates of , 2 , and based on a robust version of the variance-covariance matrix. Table 1 shows that the estimate for the spatial lag of the temporal lag ( N ln t−1 ) is not dissimilar to −̂ ̂ 1 , in line with expectation stemming from an equilibrium process. Also, Table 1 estimates are stationary and dynamically stable, as shown by the largest characteristic root of −1 N N which is equal to 0.6874, and the stationary bounds for 2 are � e −1 min = −1.1239 < 2 < � e −1 max = 1. Observe that the negative values of ̂ 2 imply positive spatial dependence among the errors. Among the instrument set, we have several endogenous variables, one is the dependent variable lagged by two periods, ln t−2 , its spatial lag lagged by two periods N ln t−2 , ln t−2 and N ln t−2 . To satisfy the orthogonality conditions and moments equations for these instruments, we require a lack of serial correlation in the it , in other words we need to satisfy the assumption that E Δ it , Δ it−2 = 0 . Arellano and Bond (1991) give a test m 2 = cov Δ it , Δ it−2 ∕s.e . which is asymptotically N(0, 1) under the null of no serial correlation. In our case m 2 = − 0.9389 with two-tailed p value equal to 0.3478. Thus, we assume that there is an absence of serial correlation as required. Note also that m 1 = cov Δ it , Δ it−1 ∕s.e = − 6.22 indicating significant first-order serial correlation as one would expect, since if the it are serially uncorrelated, Δ it has first order moving average serial correlation. A second complementary approach to testing the validity of the instrument set is via the application of the Sargan-Hansen test of over-identifying restrictions, which is equal to 253.3. This is insignificant when referred to the 2 314 distribution, and while this evidently supports the moments conditions implied by our dynamic spatial panel model, one should be cautious because it may have low power, given the presence of many moments conditions (Bowsher 2002;Pesaran 2015).

Estimates
"Appendix" Table 4 gives the estimates of some rival estimators, including one with SMA-RE errors but assuming exogenous regressors (Table A1), and with  (Table A2). As noted in Sect. 6, the predictive ability of these rivals is not as good as obtained via the preferred estimates summarized in Table 1.

Prediction
In order to support the preferred model summarized by Table 1, a cross-validation strategy is employed to assess the performance of competing estimators 'by comparing their predictive ability on data which have not been used in model estimation' (Anselin 1988). Out-of-sample predictions of the level of employment across regions are obtained for the years 2011 and 2012 using 2011 and 2012 data combined with the parameter estimates obtained for data over the estimation period from 2001 to 2010.
Following Chamberlain (1984), Sevestre and Trognon (1996), and Baltagi et al. (2014Baltagi et al. ( , 2018, the linear predictor is in which E[.] denotes the expectation, so this can be seen to be identical to Eq. (9) but with expectations. With regard to the estimate of the time-invariant component of the error term , assuming a spatial moving average error process gives Eq. (9) rewritten thus (20) In order to obtain estimates ̂ (t) estimates ̂ N = N −̂ 2 N ,̂ N = N −̂ 1 N , N = ̂ +̂ N and ĉ and ̂ are used along with random draws from t ∼ N(0,̂ 2 ). We then take the mean over time of the ̂ (t) s for t = 2, … , T , subsequently scaling so that it has variance equal to ̂ 2 , thus giving the estimate ̂ of the time-invariant error component. The outcome is the prediction Eq. (23) for T + 1 = 2011, in which 1T+1 = ln T+1 , 2T+1 = ln T+1 , and 3T+1 = ln T+1 , t = 1, … , T.

3
For two-step ahead, 11 1T+2 = ln T+2 , 2T+2 = ln T+2 and 3T+2 = ln T+2 . Figure 1 shows a close correlation between predicted log employment ln̂ T+1 and observed log employment, suggesting that the preferred estimator giving Table 1 estimates would be a good basis for simulating the impact on employment following Brexit.
The preference for Table 1 estimates is based on the mean of the Exploring Brexit with dynamic spatial panel models: some… on SMA-RE errors and predetermined (and exogenous) regressors, but with the additional variables N 1 , the spatial lag of ln , and N 2 , the spatial lag of ln k , with N given by the PBL trade data. This is thus a form of spatial Durbin specification, but with regressors t = 1t , 2t , N 1t , N 2t the additional covariates evidently cause a problem of weak instruments, giving dynamically unstable nonstationary estimates, as reflected by the largest characteristic root of −1 N N equal to 1.0663 (1.9041) and, with in Eqs. (22) and (23) RMSE = 7.4403 (3.0918). The same spatial Durbin specification again assuming predetermined regressors but with 2 restricted to zero gives a largest characteristic root equal to 1.1127 and RMSE = 3.3017 . The same spatial Durbin specification assuming exogenous regressors and with a spatial autoregressive (SAR) error process gives a largest characteristic root equal to 2.489 and RMSE = 20.9333.These results point to the viability of Table 1 estimates for prediction purposes.

Simulating the Brexit effect
The approach adopted is to use the parameters estimates in Table 1 to predict the impact on employment of presumably reduced trade between the UK and the remaining EU regions in the year 2020 and beyond. Attention is focussed on 2020 and later, given that the UK's formal exit from the EU is scheduled for the first half of 2019, so 2020 will be the first full year outside the EU. Given a lack of appropriate and accessible data, for instance with the same geography as up to 2011, beyond 2011 employment could be predicted on the basis of assumptions about the level of , and in 2020. From = 2020 onwards, there are two scenarios, one based on the trade flows assuming no-Brexit effect, and the other assuming a Brexit effect on trade flows, and the difference between them is taken as the Brexit effect. Regarding the no-Brexit effect scenario, this applies matrix N ,which is based on the latest available trade flows pertaining to the year 2010. The prediction is then given by the solution to Also is an (N by 3) matrix containing the forward projections ln , ln and ln , thus The second scenario is to assume that bilateral trade between the UK regions and the (remaining) EU regions is, for example, 2% lower than it would otherwise be. Thus, the % job-shortfall at time is lñ − ln̂ . Using the equilibrium solution of (8), but also taking into account at time = T and ̂ N , employment converges to

Fig. 9 % employment shortfall UK and Ireland due to industry
Similarly Thus, the % job-shortfall with long-run convergence at time T is ln̂ T − lñ T , hence One assumption might be that , and in 2020 are at the same level as observed in each region in 2011. An alternative assumption could be that from 2011 onwards they grow at their historical rates, taken over the period from 1991 to 2011 in each region. On this basis on average the level of and in 2025 is approximately 25% more than the 2011 levels. However, Table 2 gives simulation outcomes for the examples of Inner London and Paris which illustrate the relative insensitivity of Fig. 10 Frequency distribution from Fig. 9 1 3 Exploring Brexit with dynamic spatial panel models: some… job-shortfall to assumed regressor levels. The table shows that doubling the level of and gives equilibrium job-shortfalls that are only about 6% higher, as a result of the same increases applying to both Brexit and no-Brexit outcomes. Table 2 also shows that doubling the trade reduction in each region in effect doubles the jobshortfall. It shows that doubling the % trade reduction, say from 2 to 4% , has the effect of doubling the job-shortfall in each region. Increasing the % reduction by a factor of 8, going from 2 to 16% , increases each region's job-shortfall by a factor of 8. This means that ratio of Inner London to Paris remains stable (which for this pair of regions is approximately 1.98) regardless of what is assumed for % trade reduction. This stability of the outcome ratio exists for any pair of regions so that maps of job-shortfalls would be in a sense identical-identifying the same regions with large or small levels giving constant outcome ratios-irrespective of the assumed % trade reduction. This geographical stability is a result of the assumptions made within the simulation exercise, with trade for all UK to EU trade flows reduced by the same %. This means that the subsequent map patterns are immune to the assumed reduction in trade, although the scales would differ, where we focus on a trade reduction other than 2% . So in this way, we see we have an element of robustness in our simulations. Of course in an ideal world, one might wish to make changes to trade on an individual region by region and sector by sector basis rather than assume that trade reduces by the same amount across all regions and all sectors. However, this is very much the unknown, although some sectorally specific estimates are given elsewhere. Sectorally specific Brexit impacts are obtained by assuming that trade in specific sectors alone is restricted. While this is unrealistic, it is likely that there will be sectorally differentiated impacts but it is difficult to know by how much trade in manufactures, for example, will be reduced compared with trade in services. 12 A simple approach is therefore to assume that a specific sector is impacted by Brexit, but that there is zero impact on other sectors. This highlights the geography of the sectorspecific trade impacts, because the sectoral trade patterns have different geographies and therefore the impacts have different geographical distributions to the outcomes assuming a global reduction across all sectors (Table 3). Exploring Brexit with dynamic spatial panel models: some…

Results
The initial outcomes relate to a reduction in EU-UK trade of 2% across all sectors. The predicted % changes in employment across the EU and UK regions assuming the 2011 levels for , and are shown by Fig. 2. This shows the dynamic paths for each region to 2050, with convergence to steady state occurring after 2030. From this, it is evident that the maximum equilibrium job-shortfall is − 2.56%, in the case of Inner London, with most other regions falling below 1% . Figure 3 shows the geographical pattern of the Brexit impact equal to lñ − ln̂ for = 2025 , indicating a maximum shortfall by 2025 of − 2.34% (Inner London). The picture which emerges from the simulation is that the negative Brexit impact is diverse across regions and bilateral, with both UK regions and EU regions likely to see a job-shortfall. Figure 3 shows larger negative impacts in regions with strong trading links to the UK,  Figure 4 gives the frequency distribution of Fig. 3 data, highlighting the fact that despite some large impacts, for about 160 of the 255 regions, Brexit is likely to have close to zero effect on employment. Figure 5 shows that within the UK, Inner and Outer London ( −1.32% ) are expected to have the biggest % shortfall by 2025, with impacts generally higher along the Thames valley in Berkshire, Bucks and Oxfordshire ( −1.09% ) towards Gloucestershire, Wiltshire and North Somerset ( − 1.19% ). Generally, %s are higher around the Greater South East and in some of the large conurbations (Birmingham − 0.78% , Manchester − 0.76% , West Yorkshire − 0.67% ) than in more rural and peripheral regions. Figure 6 gives the frequency distribution of Fig. 5 data, emphasizing the Inner London outlier, with many regions having a job-shortfall of less than − 0.5% . As noted above, if one were to assume different reductions in trade other than 2% , the outcomes for employment would be different, but proportional to the 2% impact, so that the ratio of impacts in different regions and the geographical pattern would be identical.
Next, consider the separate impacts on employment of restricted trade in the manufacturing sector, defined as the production of food, beverages and tobacco, textiles and leather, coke, refined petroleum, nuclear fuel and chemicals, electrical and optical equipment and transport equipment and other manufacturing. Simulating on the same basis, region paths converge to equilibrium levels as with Fig. 2, although the equilibrium levels differ from those of Fig. 2. We take a snapshot across the dynamic paths in Figs. 7, 8, 9 and 10, showing the % shortfall in employment in manufactures by region for the year 2025. Figure 7 shows that the geography of the impact due to 2% less industry trade is very similar to the overall pattern shown in Fig. 3. However, a comparison of Figs. 4 and 8 emphasizes the differences in the levels of impact, with the maximum level in Inner London. Figures 9 and 10 show the % job-shortfall in the UK and Ireland. Again the impact in the South and East of Britain and, especially, London is clear, and the effect on Ireland remains pronounced.
We also estimate the impact of reduced trade within the manufacturing sector. Of particular interest is the group of industries defined for trade purposes as comprising 'electrical and optical equipment and transport equipment', which includes the all-important production of vehicles. With technological development, we have seen the geographical fragmentation of production processes involved in vehicle manufacture, with the development of spatially dispersed value chains as different elements of the production processes optimally located in different regions or countries. Also just in time processes means that quick and easy access to parts and components used in manufacturing vehicles is important, so interregional connectivity is important, and any disruption of it due to increased barriers to trade will have a significant impact. The impacts of a 2% trade reduction are summarized in Figs. 11 and 12, which picks out some of the production hot spots. In Figs. 13 and 14, the integrated nature of production and dispersed knock-on effects of reduced trade are evident from the relatively even spread of impacts across regions.
Compared with industry, the impact of reduced trade in services 13 is much less symmetrical, with the bulk of the job-shortfall occurring in Britain and Ireland. This is clear from Figs. 15 and 16. Inner London clearly stands out with the most significant projected job-shortfall compared with all other EU regions, and apart from the South and Eastern region of Ireland, almost all nonzero job-shortfalls occur in Great Britain. Figures 17 and 18 emphasize the polarized effect of service trade reduction, with Inner London standing out as an outlier. Outer London and Southern and Eastern Ireland see comparable effects, ahead of all the other UK regions. Focussing  Fig. 15 13 Distribution, Hotels and restaurants, Transport storage and communication, Financial intermediation and Real estate renting and business activities on the highly important Financial Intermediation sector, Figs. 19 and 20 emphasize even more strongly the asymmetric impact on Brexit, with more than 200 non-UK-EU regions having almost zero job-shortfall. The most affected non-UK-EU region is the South and East of Ireland, followed by Luxembourg. Figures 21 and 22 illustrate the role of Inner London in particular as a centre for financial intermediation, but overall for the UK increased trade barriers for this element of services seems to have a less profound impact on the job-shortfall than the more geographically widespread and deeper impact of reduced trade in transport equipment.

Conclusion
The paper shows negative Brexit-induced impacts on employment which affect not only the UK regions but also employment levels in EU regions, especially those which are close trading partners. This pan-European interregional interdependency is captured in the state-of-the-art model by spatial and temporal interactions based on the best available trade flow estimates which determine the strength of interdependence. This means that employment within a region not only depends on the levels of output and capital within the region, but also on demand coming from other regions which are trading partners. The impacts will be multi-way, what happens to employment in London depends partly on what happens in Paris, which depends on what happens in Munich, which depends on what happens in London, etc. The approach adopted has been to assume a reduction in trade between EU and UK regions which gives a corresponding reduction in demand for jobs. The predicted job-shortfalls depend on the assumed global % reduction in trade between UK and EU regions, but the modelling assumptions ensured stability in the geography of Brexit impact.
In the paper, the impact of Brexit is measured in terms of the 2025 job-shortfall, which is the reduction in the number of jobs in each region due to Brexit assuming no alternative sources of employment are put in place. This of course might be a false assumption, as the pro-Brexit lobby has consistently emphasized the potential stimulus of new trade deals with other non-EU countries. Therefore, the Brexit impact as reflected in the maps of job-shortfall indicates those regions which could be in the greatest need of alternative compensating sources of employment. Thus, the paper is not predicting a job-loss per se, simply a potential job-loss without successful alternative trade arrangements post-Brexit. Additional employment due to trade diversion effects due to higher UK-EU trade barriers (Ortiz Valverde and Latorre 2018, Dhingra et al. 2017a, Krueger 1999) could possibly be captured within the current modelling set-up via changes to the levels of output and capital in each region, but these would be difficult to estimate and there is some empirical evidence that they might be quite small (Krueger 1999;Magee 2008).
Key outcomes are as follows. First, both UK and EU regions are negatively affected, this is a lose-lose scenario. Second, the deepest most concentrated impact is for the UK, many EU regions are barely affected, but the South and East regions of Ireland are the worst affected EU region, with impacts on a par with the worst affected UK regions. In addition the Ile de France, Oberbayern, Stuggart and Dusseldorf stand out as regions likely to see significant job-shortfalls. Third, in the UK, Southern regions, especially in and around London, and big cities, are expected to see the largest job-shortfall. Fourth, the effect of manufacturing trade reduction will be evidently larger than for services, but the service impact is more asymmetric, with the bulk of the job-shortfall focussed on UK regions, especially London. This is even more the case when comparing the effect of trade reductions in vehicles and financial intermediation.
Overall, the simulations suggest that the biggest Brexit impact on UK regions will occur in the richer South East and urban areas, which is in line with work from LSE based on GVA, which shows that 'areas in the South of England, and urban areas, are harder hit by Brexit...the areas that were most likely to vote remain are those that are predicted to be most negatively impacted by Brexit' (Dhingra et al. 2017b). This interpretation is in direct contrast to other work which maintains that 'the regions which voted Leave also tended to be more dependent on Europe for their prosperity than the regions which voted Remain' (Los et al. 2017).
Clearly, Brexit is a complex phenomenon leading to diverse interpretations of outcomes, as evident in the special issue of Papers in Regional Science (McCann 2018). The outcomes presented in this paper are based on model assumptions, but it is argued that the main driver of the results is the data, not imposed assumptions. Nevertheless, great caution is needed in interpreting the validity and value of any 'prediction' effort. It is worth recalling the words of Box and Draper (1986), 'Essentially, all models are wrong but some are useful'. David Spiegelhalter, Professor of the Public Understanding of Risk at the University of Cambridge, refers to Donald Rumsfeld as the patron saint of Risk Analysis, who will be remembered for famously saying that 'but there are also unknown unknowns. There are things we do not know we don't know'. We should therefore put forward predictions with all due humility, but clearly and without fear, because we don't want to come across as 'dithering scientists'. In defence of the approach adopted, there is support from the words of Pesaran (1990), who points out that 'Econometric models are important tools for forecasting and policy analysis, and it is unlikely that they will be discarded in the future. The challenge is to recognize their limitations and to work towards turning them into more reliable and effective tools. There seem to be no viable alternatives'.

Appendices
See Table 4. Exploring Brexit with dynamic spatial panel models: some…

Theoretical background
The theoretical background to the model specification commencing with Eq. (1) is derived from standard urban economics as given by Abdel-Rahman and Fujita (1990), Ciccone and Hall (1996) and Fujita and Thisse (2002), among others. The theory commences with a constant returns to scale Cobb-Douglas production function for the output q of the competitive final goods and services sector in which m denotes sector-specific labour efficiency units and the level of composite services is given by i, and h is area of land. Assume h = 1 theñ Assume that the equilibrium output of each service firm is i e and there are g firms, depending on the total services effective labour. We obtain g by dividing the total services effective labour by the effective labour per firm, thus In (26), (1 −̃ ) equals the services labour share of total effective labour ẽ under competitive equilibrium in the labour market. For each services firm, we have internal returns to scale, with s denoting the fixed labour requirement and a the marginal labour requirement, so that ai e + s is the effective labour per firm at equilibrium. From the CES production function, we obtain where ̃ is a measure of monopoly power in the monopolistically competitive services sector, and ̃ ̃ −1 is the constant elasticity of substitution. Substituting i into equation (25) gives Substituting for g in Eq. (28) gives which simplifies to This shows that with = 1 there are increasing returns ( � > 1) if service firms are relevant to output in the competitive sector ( � < 1) and also possess monopoly power � > 1 . However, < 1 indicates that land is also a relevant factor, and depending on the value of a tendency to increasing returns could be offset by the congestion effect caused by the restriction of production to a unit of land, leading to � < 1 hence diminishing returns.

Other estimates The Chow-Lin approach
This commences with aggregate trade values between Ñ = 21 EU counties (denoted by the square (Ñ ×Ñ) matrix Ñ , in which the subscript Ñ means national). There are Ñ = 21 unobserved intra-country trade flows, thus giving p = 420 observations for the year 2000. The Ñ 2 = 441 'observations' are the dependent variable in a Weighted Least Squares regression model. All observed trade flows are given the weight 1, and those unobserved weighted zero. In terms of parameter estimates, an entirely equivalent procedure is to estimate the 420 observed trade flows by OLS. In the regression, the explanatory regressors are great circle distances between country pairs ( Ñ ), and the product of each pair of country's national GVA level ( ̃ Ñ ) in the year 2000, so that given the (Ñ × 1) vector ̃ Ñ , � Subsequently, Ñ , Ñ and Ñ are reshaped as (Ñ 2 × 1) vectors Ñ , Ñ and Ñ , and together with Ñ which is an (Ñ 2 × 1) vector of ones, these variables provide the input for the regression denoted by Eq. (31), thus giving the estimates ̂ Ñ . Interregional trade estimates are based on the (national-level) regression parameter estimates ̂ Ñ and on the estimated regression residuals ̂ Ñ . Thus, trade flows between regions, denoted by R , are obtained by applying the national-level estimates ̂ Ñ to the regressors measured at the regional level, denoted by R and R . Also, an equal share of the national-level residuals ̂ Ñ is added to regions corresponding to country pairs. This process is summarized by the two equations: with log bilateral interregional trade flows ln R then being obtained using To obtain the (R 2 × 1) vector ̂ R , we calculate the (R ×R) matrix � � R = � � N ⊤ where = ⊤ ( ⊤ ) −1 and is an (Ñ ×R) indicator matrix with ones in each country's row indicating which regions are within that country, and zeros indicating those which are not. The (Ñ ×Ñ) matrix of residuals ̂ Ñ is formed by reshaping the (Ñ 2 × 1) vector ̂ Ñ so that ̂ Ñ (1 ∶Ñ, 1) =̂ Ñ (1 ∶Ñ),̂ Ñ (1 ∶Ñ, 2) =̂ Ñ (Ñ + 1 ∶ 2Ñ), … , Ñ (1 ∶Ñ,Ñ) =̂ Ñ (Ñ 2 −Ñ + 1 ∶Ñ 2 ). The outcome ̂ R is an (R ×R) matrix containing the inter-country residuals ̂ Ñ allocated equally to all regions corresponding to country pairs. Also, there are block diagonal zeros as a result of the unobserved intra-country trade flows, and zeros across four rows representing four regions equal to the countries Estonia, Lithuania, Latvia and Luxembourg which were excluded from the initial regression. Finally, ̂ R is reshaped as the (R 2 × 1) vector ̂ R of Eq. (32) to give ln R . The resulting (R ×R) matrix of interregional trade flows is * N = R .