The Way Toward Growth: A Time-series Factor Decomposition of Socioeconomic Impulses and Urbanization Trends in a Pre-crisis European Region

The present study investigates long-term urbanization and suburbanization trends - and the consequent impact on economic expansion and social change - in a divided region of Mediterranean Europe (Attica, Greece) by performing a time series (1965–2008) dynamic factor analysis of 14 socioeconomic indicators that re�ect different aspects of metropolitan growth. Attica was partitioned in two spatial domains, the ‘Greater Athens’ area (hereafter the ‘core’ district) and the rest of the region (hereafter the ‘ring’ district) with the aim at quantifying the (possible bi-directional) spatio-temporal propagation of socioeconomic impulses to metropolitan growth. The exploratory scheme, integrating Multi-way Factor Analysis (MFA) with Continuous Wavelength Transform (CWT) and rapidity-of-change metrics grounded on complex thinking, delineates latent mechanisms of urban expansion, indicating substantial divergences in the development path of the two districts. While the ‘core’ district experienced population increase and settlement densi�cation, the ‘ring’ district underwent a suburbanization process resulting in a moderate – and slower – concentration of economic functions. Re�ecting – at least in part – growth impulses’ propagation from urban to rural areas, the economic interplay between ‘core’ and ‘ring’ areas delineates a complex development path accelerating spatial polarization in central and peripheral locations. Our results de�nitely highlight the importance of ‘system thinking’ in regional studies and applied economics.


Introduction
Theoretical approaches referring to the City Life Cycle (CLC) notion de ne urbanization as a process when settlements expand at the cost of the surrounding areas (Cross, 1990;Lever, 1993;Felsenstein, 2002), implying an intense movement of population and activities from the countryside towards central locations (Mordridge and Parr, 1997;Van Criekingen, 2010;Nam et al., 2012).The intrinsic concentration of industry and services in core settlements (Klaassen et al., 1981;van den Berg et al., 1982;Economou, 1997) -together with other factors such as the promotion of new life-style patterns (Haase et al., 2010;Kabisch and Haase, 2011;Makarem, 2016), the increased nancial and institutional power of metropolis (Hall, 1997a(Hall, , 1997b;;Gordon and Cox, 2012), and infrastructure development (Hale and Moberg, 2003) -determined favorable conditions for economic growth in central locations (Lucy and Phillips, 1997;Timàr and Varadi, 2001; Lee and Mason, 2010).With the expansion of population and activities (Dahan and Tsiddon, 1998), negative externalities such as tra c, pollution, crimes, and poverty arose in central cities, pushing medium-and high-income households to move to suburbs with better environmental quality and housing conditions (Vaiou, 1997;Maloutas, 2007;Cervellati and Sunde, 2011;Morelli et al., 2014).As a consequence, ring's population increased compared to city cores' population (e.g.Salvati and Sabbi, 2014), resulting in a considerable diversi cation of land-use, economic activities, and social classes in peri-urban locations (Zhou and Ma, 2000;Arapoglou and Sayas, 2009;Heider and Siedentop, 2020).These premises, and other critical issues intrinsically associated with the CLC (e.g.Champion, 1990;Cheshire, 1995;Bourne, 1996;Lestegàs, 2019), justify the development of a comprehensive (multivariate) framework (sensu Coccossis et al., 2005;Couch et al., 2007;Arapoglou and Sayas, 2009) identifying and characterizing the different stages of a metropolitan cycle -considering together the synergic impact of both economic and non-economic (e.g.social, demographic, environmental, territorial) factors of change (Kroll and Kabisch, 2012; Zambon et al., 2019; Torrado et al., 2021).Unfortunately, comprehensive multi-dimensional approaches assessing socioeconomic transformations are relatively scarce (e.g.Morya and Ram, 2020;De Vidovich and Scolari, 2022;Salvati, 2022).Moreover, experimental designs speci cally aimed at distinguishing sign, direction, and timing of the impulses to development paths intrinsically associated with urbanization and suburbanization, are partial and provide sometimes mixed evidence (apart from few exceptions, e.g.Cuberes, 2011;Sànchez-Vidal et al., 2014;Sheng et al., 2014).On the contrary, metropolitan complexity re ecting less simplistic development paths justi es an in-depth analysis of multi-domain indicators (Redding, 2022).This analysis is increasingly required to investigate together socioeconomic and territorial dimensions of regional growth and change (Oppio et al., 2018;Walker, 2018;Wang and Dong, 2022), going beyond the evidence grounded on a quantitative analysis of individual indicators (Zanella et al., 2015;Delmelle et al., 2021;Fernandez and Hartt, 2021).Evidence from such analyses may inform spatial planning reconnecting local speci cities with broader trends in regional development (Giuliano et al., 2019; Shen and Wu, 2020; Bailey, 2021).
Assuming local development as an eminently multivariate concept constituted of several (e.g.social, economic, territorial, political) interacting dimensions (e.g.Rodriguez-Pose and Fratesi, 2007; Di Feliciantonio and Salvati, 2015; Cuberes and Gonzales-Val, 2017), we investigated the latent mechanism at the base of economic impulses' propagation during urbanization and suburbanization in two (neighbouring) spatial domains, the 'Greater Athens' area (hereafter the 'core') and the rest of the Attica region (hereafter the 'ring').A time series  dynamic factor analysis of 14 socioeconomic indicators that re ect different aspects of metropolitan development (e.g.Pili et al., 2017) was integrated with a Continuous Wavelength Transform (CWT) approach (Salvati, 2022) with the aim at deriving 'rapidity-of-change' metrics grounded on complex thinking (e.g.Salvati and Serra, 2016).This operational framework contributes to identify (apparent and latent) mechanisms of metropolitan growth (e.g.Fratesi and Rodriguez-Pose, 2014), quantifying the interrelationships between local development paths (De Rosa and Salvati, 2016).

Methodology Study area
Attica is the most populated administrative region of Greece (Nuts-2 level of the European Territorial Statistic Nomenclature) covering nearly 3,800 km 2 of mainland and insular districts and overlapping, for most of the area, the boundaries of the 'Urban Atlas' region of Athens -the capital city of Greece (Di Feliciantonio et al., 2018).Attica includes more than 100 mainland municipalities -both urban and rural -and various islands in the Saronikos gulf, Aegean Sea (e.g.Salamina, Aegina, Poros, Hydra, Spetses).The region mostly consists of mountains bordering the at area hosting the Greater Athens' agglomeration (administered by nearly 50 municipalities).Three coastal plains (Messoghia, Marathona, Thriasio) are located outside the Greater Athens' area and host the largest proportion of rural population in Attica (Morelli et al., 2014).
To investigate suburbanization patterns and economic (both synchronic and delayed) impulses of urbanization, an inner district ('core') and a 'ring' area with distinctive demographic, economic, and territorial characteristics (respectively corresponding with Greater Athens and the rest of Attica) were identi ed and selected as the elementary spatial units of this study (Economou, 1997).In the mid-1960s, Greater Athens was an expanding urban area (4050 inhabitants/km 2 ) while the rest of Attica played the role of a strictly rural region (60 inhabitants/km 2 ) surrounding the central city.Suburbanization in the 'ring' district was documented in earlier studies as a process depending on socioeconomic dynamics in downtown Athens (Pili et al., 2017).

Data sources and socioeconomic indicators
A total of 14 contextual indicators covering the time period between 1965 and 2008 were made available on a year base from o cial statistics aggregated for two districts: Greater Athens ('core') and the rest of Attica ('ring').Indicators were developed from elementary data derived from the database collected by ELSTAT at prefectural (NUTS-3) level in Greece.Socioeconomic time-series data were also collected from the www.nomoi.grweb-site disseminating a set of information on the economic performance of Greek regions.Indicators were subdivided into four development dimensions covering a number of elementary (socioeconomic) issues (Table 1): (i) population and job market (3 indicators), (ii) land-use and settlements (4), (iii) economic performance (3), and metropolitan functions (4).
The selected indicators were taken as representative of dynamic urban contexts found in Athens between the mid-1960s and the mid-2000s, in accordance with earlier studies (Salvati and Serra, 2016).The most recent period (since the late 2000s) was excluded from the analysis because it is representative of the most intense economic recession hitting Greece in the last half century.Being characterized by austerity urbanism and inherent perturbations in both housing and job markets (Morelli et al., 2014), this stage determined important, and less predictable, impacts on local development (Di Feliciantonio et al., 2018).These impacts depend on exogenous shocks, being unrepresentative of the effects of long-term (or short-term) economic dynamics and structural change (Salvati and Serra, 2016), whose in-depth analysis is the main scope of the present work.

Statistical analysis
In this study, socioeconomic development was assumed as a response of complex adaptive systems to metropolitan transformations (Salvati and Serra, 2016).Being regarded as an intrinsic property of the system, rapidity of change in a comprehensive set of local development descriptors was investigated using a dynamic multidimensional analysis run on 2 data matrices (distinguishing 'core' from 'ring' dynamics) constituted of indicators reported in Table 1 by year (n = 43 observations from 1965 to 2008).A statistical framework was adopted in this study integrating four analysis' steps that include: (i) exploratory time series analysis, (ii) parametric and non-parametric correlations, (iii) a Multiway Factor Analysis (MFA) identifying 'fast' and 'slow' variables and estimating rapidity of change, as well as an (iv) advanced time series analysis in the frequency eld (based on a Continuous Wavelength Transform algorithm).Steps (i) and (ii) prepares data to a re ned, analysis in steps (iii) and (iv).Taken together, this procedure delineates the intimate, multivariate relationship between developmental processes in 'core' and 'ring' districts, with the nal aim at assessing metropolitan transitions, and the importance of selected economic (and non-economic) drivers of urbanization and suburbanization in Athens.More speci cally, MFA was applied to our case study with the aim to select independent, latent variables (factors) managing data redundancy and serial autocorrelation in a multivariate distribution of observations (Coppi, 1994).This methodology highlights complex structures in higher-order datasets, where data have three (or more) dimensions (Coppi et al., 2010).Linking different variables with comparable spatio-temporal patterns on a few signi cant factors (Lavit et al., 1994), this strategy provides an indirect measure of the extent to which a system's characteristics (i.e.development indicators) have substitutes to ensure functioning in the event of a transition or a shock (Esco er and Pagès, 2008).

Exploratory time-series analysis
Table 2 reports the results of a classical (univariate) time-series analysis based on coe cients of temporal autocorrelation (total and partial) by lag.Considering the intrinsic limit of the available time series (n = 43 observations over time), analysis was restricted to ve lags, from 1 to 5. The empirical results of this analysis documents the evident time structure of some variables, especially for 'core' district (Den, Car, Dru).However, results of partial autocorrelations indicate signi cant coe cients only at lag 1.The temporal structure of both data matrices ('Core' and 'Ring') was thus removed through computation of rst-differenced time series, whose temporal autocorrelation resulted to be systematically insigni cant for all variables and districts analyzed.X(I,J,T) = {xijt}, i = 1,…,I; j = 1,…,J; t = 1,…,T (Equation 1) where indexes i, j, and t correspond with variables, units, and times, respectively.From a descriptive (nonprobabilistic) perspective, MFA combines a cross-section Principal Component Analysis (PCA) through (i) spectral decomposition of a correlation matrix computed on the input dataset with (ii) investigation of the time series dimension balancing data matrices' contribution to the overall variability (Coppi et al., 2010).
Equation 1 was extended to a broader speci cation in line with the objective of this study analyzing the following X(I,T,S) array: where indexes i, t, and s correspond with variables, times (i.e.years), and spatial units (i.e.'core' and 'ring' district), respectively.Input variables were standardized prior to analysis.Based on the seminal work by Coppi and Bolasco (1989), MFA decomposes the variance-covariance matrix R relative to X(I,T,S) where the two-way behavior of each local system (i.e.development path) is described with a 'variables-times' representation (I•T).Concerning the I•S observations over the S units, R was decomposed into the sum of three variancecovariance matrices: R = *S T + *S S + *S TS (Equation 3) where *S T is the matrix of the static structure of the spatial units, i.e. the variance-covariance matrix calculated on the variables' average (by year) for each district; *S S is the matrix delineating the average system's dynamics, i.e. the variance-covariance matrix run on the average of spatial units' behavior and re ecting the spatial variability (i.e.across districts) of the temporal units, independently of the temporal dynamic characteristic of the individual districts.
*S TS is the matrix of the differential dynamics between districts, i.e. the variance-covariance matrix of the interactions between time units (years) and space (metropolitan districts), re ecting the variability due to the difference between the dynamic of the overall average of the units (the two districts jointly), and the speci c dynamic of each individual districts.
Considering the X(.,.,S) dimension, if all the sets of variables are introduced as active elements (like in a general factor analysis), i.e. without balancing their in uence, a single set can contribute quite by itself to the construction of the rst axes (Esco er and Pagès, 1994).Based on these considerations, MFA introduces a generalized, weighted factor analysis, with constant weights for the variables of the same group, and varying weights assigned to variables belonging to different groups (in our case, S dimension, see above).Going beyond repeated (partial) factor analyses as the number of variables' groups (Pagés, 1996), a global analysis in a MFA perspective is de ned when two (or more) sets of variables are simultaneously introduced as active ones, and require balancing the in uences of these sets or groups (i.e. the S dimension).
The in uence of one set derives from its structure, in the sense of its inertia distribution in the different space directions (Esco er and Pagès, 2008).For instance, if a set presents a high inertia in one direction, this direction will in uence the rst axis of the global analysis (Coppi, 1994).This suggests to normalize the highest axial inertia of each set.Balancing maximum axis' inertia rather than the total inertia gives the MFA important properties that can be used in several applications (Bove and Di Ciaccio, 1994).Weights are constructed in a way that the maximum axial inertia of a variable's group is equal to unity, i.e. assigning a weight equal to the inverse of the rst eigenvalue of the group' s factor analysis to each within-group variable (Lavit et al., 1994).Technically speaking, it is done by weighting each variable of the set j by with λj 1 denoting the rst eigenvalue of the (partial) factor analysis applied (restrictively) to j variables' set.
Having normalized each cloud of variables-observations making its highest axial inertia equal to 1 (Pagès, 1996), it should be noted that this weighting system does not balance total inertia of the different sets, thus implying that a set having a high dimensionality will have a high global in uence in the sense that this set will contribute to several (extracted) axes (Esco er and Pagès, 1994).
As an exploratory analysis not grounded on hypothesis testing, the selection of signi cant factors in MFA was based on a priori eigenvalue threshold (Esco er and Pagès, 2008), similarly to PCA.At the same time, MFA provides classical outputs of general factor analysis, that is to say, for each axis: (i) coordinates, contributions and squared cosines of individuals, as well as (ii) correlation coe cients between factors and continuous variables (Kroonenberg, 2008).We speci cally considered loadings and scores respectively from dimension (i) and (ii) when de ning independent dimensions of local development in the study area (Morelli et al., 2014).
MFA results thus allow an explicit evaluation of changes over time in the position of each unit and case since they are projected into the same factorial plane (Cecchini et al., 2019).
Investigating 'fast' and 'slow' variables of local development using MFA MFA assures a homogeneous representation of the trajectory over time (i.e. the intrinsic change: sensu Kroonenberg, 2008) for each development indicator moving from a given area to another one, e.g. from 'core' to 'ring' districts.In this perspective, co-ordinates of the center of gravity of variables (i.e.developmental indicators) by observation (i.e.years) and spatial unit (i.e.districts) can be used (González et al., 2009) to derive appropriate metrics of rapidity-of-change.
In other words, we assumed local development as a multivariate, latent dimension of growth evolving over time with separate characteristics among districts in the same metropolitan area (Morelli et al., 2014).MFA was adopted to identify such individual paths ('core' and 'ring' separately) and the resulting equilibrium dynamic at the metropolitan scale (the results of the global analysis).This rationale allows identi cation of 'fast' and 'slow' drivers of change (Ferrara et al., 2016).At the same time, considering all the indicators together, estimation of whole-system rapidity of change over time (i.e.year by year), was carried out with the same logic, i.e. estimating the differential rapidity of developmental dynamics between 'core' and 'ring' districts (Morelli et al., 2014).Following Salvati and Serra (2016), a multidimensional metric of rapidity-ofchange (C' i ) for each developmental indicator was derived as the Euclidean, n-dimensional distance between factor loadings at 'core' and 'ring' districts according to the following equation: C' i = √((x a,1 − x a,0 ) 2 + (x b,1 -x b,0 ) 2 + (x ...,1 − x ...,0 ) 2 + (x n,1 − x n,0 ) 2 ) (Equation 5) where x a,s is the coordinate loading on factor a at a given location s ('ring' = 1 or 'core = '0') and n is the number of extracted factors, i.e. overpassing the eigenvalue threshold mentioned above (Ferrara et al., 2016).
Speci c dynamics in developmental indicators were regarded as 'fast' (or 'slow'), if the related C' i was above (or below) the median value of the overall rapidity of change computed for the respective time interval (Cecchini et al., 2019).
Considering all the indicators together (i.e.re ecting the local development path), estimation of whole-system rapidity-of-change over time (i.e.year by year), was carried out with the same logic, i.e. estimating the differential rapidity of developmental dynamics between 'core' and 'ring' districts (Morelli et al., 2014).Following Salvati and Serra (2016), a multidimensional metric of rapidity of change (M') was calculated as the Euclidean, n-dimensional distance between coordinate scores observed at 'core' and 'ring' separately for each temporal unit (year) according to the following equation: M' (t) = √((y a,1 − y a,0 ) 2 + (y b,1 -y b,0 ) 2 + (y ...,1 − y ...,0 ) 2 + (y n,1 − y n,0 ) 2 ) (Equation 6) where y a,s is the coordinate score on factor a at a given location s ('ring' = 1 or 'core = '0') and n is the number of extracted factors (i.e.overpassing the eigenvalue threshold mentioned above).Rapidity-of-change estimated for each year (t) was illustrated graphically (Salvati and Zitti, 2009); this plot evaluates similarities between all pairs of points along with the time series, representing the degree of similarity between values at the rst and last point of the time series (Gonzales et al., 2009).
Investigating the temporal structure of factors with a cross-correlation analysis Cross-correlation coe cients by lag (from -10 to 10, including the synchronic time at lag 0) were computed on the M' metric time-series considering separately dynamics at 'core' and 'ring' districts and the absolute difference between them.Cross-correlation is a statistical technique implemented in time series analysis with the aim at identifying signi cant lags in the degree of correlation between the same variable observed at two districts (i.e.'core' vs 'ring') over a given time span (Christensen, 1991).In line with the characteristic of the input matrix, we assumed the dominance of synchronic cross-correlations against cross-corrections at any lag.Comparing the empirical results of cross-correlation analysis (both synchronic and lagged) derived from different coe cients (e.g.parametric vs non-parametric) provides indirect proofs of robustness and conceptual representativeness of the exploratory multivariate approach (e.g.Salvati, 2014).
In our case, the results by lag derived from a traditional cross-correlation coe cient (based on productmoment Bravais-Pearson coe cient) were checked for robustness (i.e.linear vs non-linear correlation) using a series of non-parametric, Spearman cross-correlation coe cients calculated at the same lag series.Both coe cients range between 1 (the highest positive cross-correlation over time at a given lag) and -1, representing the highest negative cross-correlation (Reinsel, 2003).Being appropriate for analysis of time series with different statistical characteristics, comparison between different correlation coe cients allows a re ned exploration of data with linear and non-linear temporal structure, under both normal and non-normal statistical conditions (Salvati, 2022).A comparative analysis of cross-correlation structures nally provides an enhanced overview of a given phenomenon, irrespective of the length over time (i.e. the number of observations), complexity in spatial interactions, and the intrinsic statistical features of the input data.

Investigating local development paths with Continuous Wavelength Transform
MFA was complemented with results of a Continuous Wavelet Transform (CWT) approach adopting a Morlet function and running on the complete time series of rapidity of change (M' (t) ) metric by year and district.The CWT is a procedure exploring latent structures of data sets at small, intermediate and large scales simultaneously with the nal objective at identifying (i) periodicities at multiple wavelengths, (ii) self-similarity, and (iii) other (intrinsic) features characteristic of the studied time series.CWT algorithm is based on fast convolution of the signal with the wavelet at different time scales, using the fast Fourier transform (Wei, 2018).Based on CWT, a graphical outcome called biplot was calculated, with the vertical axis (logarithmic size scale) representing the signal recorded at a scale of only two consecutive data points at the bottom, and a scale of one-fourth of the whole sequence at the top.Assuming that one unit on this axis corresponds to a doubling of the size scale, the bottom of the CWT biplot provides a ne-grained representation of short-term time trends, and the top depicts a smoothed overview of trends at broader time scales.

Multiway Factor Analysis: summary results
Considering together the two datasets that illustrate the socioeconomic dynamics in both 'core' and 'ring' districts of the Athens' metropolitan area, the dynamic multi-factor analysis (MFA) has extracted ve axes that explain 47% of the overall variance (Table 4).Generally speaking, the extracted axes document a similar distribution of variables' loadings in both metropolitan districts.However, Axis 1 (12.9% of the total variance) highlights similarities and some differences between 'core' and 'ring' districts as regards speci c (socioeconomic) functions.Natural population balance (Nat) received a positive loading being higher in the 'ring' than in the 'core' district.This nding may re ect a more intense population growth in suburban locations because of the younger demographic structure.The proportion of agricultural areas in the landscape (Cro) assumed signi cant loadings with different signs between 'core' and 'ring' districts, con rming the land-use divide in urban and rural areas typical of post-war compact cities.The density of doctors (Doc) showed a similar loading's distribution, having totalized negative and positive values respectively for 'core' and 'ring' districts.Finally, the density of commercial businesses (Dru) received negative loadings -higher in the 'core' than in the 'ring' -highlighting a preference for central location typical of trade and tourism activities.As a proxy of urban concentration, a positive loading was nally assigned to density of private cars (Car) exclusively for 'core' district; conversely, a negative loading was assigned to household electricity consumption (Ele) in the 'ring' district.This indicator refers to the percentage of residential consumption in total consumption and may represent a proxy of residential suburbanization -con rming the preferential location of industries and other economic activities in central areas.
Axis 2 (11.3% of total variance) summarizes the multivariate dimension of settlement and land-use as opposed to spatial divides in local nance performances and negative externalities of urban concentration.In fact, the number of new dwellings per 100 inhabitants (New) received a positive loading for both districts (higher in the 'core' than in the 'ring' -indicating the latent effects of urban densi cation), while the average number of rooms per apartment (Roo) achieved a positive loading only in the 'core' district.In the 'ring' district, positive loadings were observed for the ratio between expenses and revenues of local nance (Exp) and for the density of medical doctors (Doc), while the density of road accidents (Acc) has received a negative loading.Axis 3 (8.7% of total variance) depicts the opposition between negative externalities of agglomeration (Acc, receiving positive loadings in the 'core' district) and density of primary services (Hos, receiving negative loadings in the 'core' district) which showed, in turn, a clear spatial polarization (Hos received a positive loading in the 'ring' district).
Axis 4 (7.8% of total variance) delineated a gradient of settlement size and territorial divides in urban and suburban areas, attributing positive loadings to the number of rooms per dwelling (Roo), being higher in the 'core', and the percentage of cropland in total landscape (Cro) -being higher in the 'ring'.This axis, speci cally for the 'ring', also highlighted the imbalance between the demand for private transport (Car, negative loadings) and the supply of some primary public services (Hos, positive loadings).Con rming the ndings of Axis 4, Axis 5 (6.9% of the total variance) highlights two dimensions of economic polarization distinguishing 'core' from 'ring' districts.In the 'core' district, there was an imbalance between the supply of primary services (Doc, positive loadings) and the demand for private transport (Car, negative loadings), while in the 'ring' district there was an imbalance between settlement endowment (Roo, positive loadings) and performance local nance (Exp, negative loadings).Results of a cross-correlation analysis by lag run separately on the rst ve axes (A1-A5) extracted from MFA (Table 5) document the existence of signi cant synchronic (i.e.lag 0) and positive cross-correlations for all the tested axes (from 1 to 5), with decreasing correlation coe cients (moving from 0.70 for Axis 1 to 0.47 for Axis 5).Cross-correlations at non-0 lags were insigni cant (p > 0.05) in all cases, excepting for two weak positive cross-correlation signals detected at lags +2 and +4 exclusively for Axis 1.However, correlation intensity in both cases (i.e. the value of the associated coe cient) was systematically below the intensity of the corresponding synchronic (i.e.lag 0) correlation coe cient.Comparable results were derived from nonparametric cross-correlations (data not shown).Based on MFA loadings, the rapidity-of-change characteristic of local development paths in Greater Athens and the rest of Attica was calculated by time interval using a standardized metric that estimates the intrinsic dynamics (i.e.distance) between 'core' and 'ring' districts (Figure 1).Assuming an average value of the distance metric close to 0.7, 'fast' variables included Cro, Irr, Acc, Exp, Hos, and Doc.These variables showed a marked dynamism over time, with increasing differentiation between 'ring' and 'core' districts.This speci c evolution was basically associated with development dimensions such as (i) agriculture (Cro, Irr) -more dynamic at suburban locations than urban locations, (ii) economic performances (Acc) -increasing more evidently at urban locations than suburban locations, and (iii) metropolitan functions (Hos, Doc, Dru)expanding more rapidly in suburbs than in core settlements.Population and settlement features assumed values around the average (Den, Roo) or below average (Nat, Roo), indicating a more static evolution over time.A similar dynamic between 'core' and 'ring' district was observed for job participation (Par), electricity consumption (Ele), and car density (Car), suggesting how these variables were less able to discriminate developmental processes at central and peripheral locations.
Quantifying 'rapidity of change' by district and year Rapidity-of-change characteristic of each sub-system (i.e.district) and year was quanti ed using a standardized metric derived from MFA results (Figure 2 'core' and a bit delayed and prolonged over time (i.e. the 1990s) for 'ring'.Signi cant parts of the spectrums were delineated with black borders.The difference in the standardized metric between 'ring' and 'core' districts was signi cant for both the 1990s and the early 2000s, delineating an evident acceleration of developmental dynamics in the 'ring' in respect to 'core' that may re ect intense suburbanization in the area.

Discussion
A dynamic multi-factor analysis integrated with exploratory time series approaches contributed to reveal the spatio-temporal direction of urbanization and suburbanization, delineating the individual contribution of speci c processes of change to metropolitan dynamics and local development (Salvati and Serra, 2016).The empirical results of our study provides, at the same time, a place-based interpretation of urbanization and suburbanization processes, and insights for a more general analysis of long-term metropolitan growth and local development in advanced -but spatially polarized -economies (Kabisch and Haase, 2011;Di Feliciantonio and Salvati, 2015;Makarem, 2016).Earlier studies have documented the variable impact of a vast ensemble of socioeconomic forces underlying urbanization and suburbanization, distinguishing economic from non-economic drivers of change at different spatial scales (Arapoglou and Sayas, 2009;van Criekingen, 2010;Gordon and Cox, 2012).The extension to which each factor contributes to urban growth depends largely on place-speci c conditions (Felsenstein, 2002;Maloutas, 2007;Nam et al., 2012).
The results of dynamic multi-factor analysis have also indicated how changes over time in both economic and non-economic variables may exert a relevant contribution to suburbanization (Haase et al., 2010) The methodological framework presented in our work appears to be exible enough and oriented towards the exploratory analysis of multivariate datasets characterized by a limited number of observations, which frequently happens when working with time series data concerning disaggregated spatial partitions fully comparable over time.Originally used for the investigation of temporal dynamics (Coppi and Bolasco, 1989), the extension of the multi-way analysis to the spatial dimension allows the veri cation of the latent relationships between relevant variables in a perspective of complex system thinking (Salvati and Serra, 2016).Based on a continuous wavelength transform algorithm, the exploratory analysis of time series in the frequency domain re nes and completes the results of the multi-way analysis (e.g.

Conclusions
Our analytical framework highlights the differential development path between metropolitan districts and the possible (logical and temporal) interconnections between the individual dimensions of the development process.Overcoming the limits of urban economic theory and econometric approaches based on the speci cation of causal relationships between variables, an exploratory framework oriented toward complex thinking seems to provide an original and updated interpretation of latent development dynamics at both urban and metropolitan scale.This analysis provides innovative knowledge informing spatial planning and sustainable development/territorial cohesion policies.

Figures
Figure 1 , panel a) referring separately to 'core' and 'ring' areas.A standardized difference in the metric estimating rapidity-of-change between the two districts was further calculated (Figure2, panel b).Both metrics outlined the existence of two time intervals with a distinctive evolution in the multivariate set of developmental indicators: (i) a stage dominated by compact urbanization with huge population increase (mid-1960s to mid-1980s) in central locations and less intense development of peripheral districts, and (ii) a subsequent stage re ecting spatially dispersed (low-density) settlement expansion in more peripheral locations with stable (or slightly declining) population at inner districts (mid-1980s to mid-2000s).In the second stage, the estimated rapidity-of-change assumed values systematically higher in the 'ring' than in the 'core' district (Figure2, panel a).Additionally, rapidity-of-change was clearly more homogeneous in the rst stage (compact urban growth) than in the second stage (dispersed urban growth), suggesting the existence of mixed and fragmented socioeconomic processes at the base of suburbanization dynamics.These ndings may support the existence of a 'dual growth' mechanism characteristic of metropolitan development in the study area.The resulting difference metric ('ring' minus 'core') and the associated 5-years moving average (Figure2, panel b) indicate a substantial acceleration in the socioeconomic dynamism of suburbs at least since the mid-1980s, after a period of stability or moderate increase.Interestingly, rapidity-of-change slowed down at the end of the study period -being sensitive to the incipient recession.Con rming the empirical results presented earlier in this paper, the temporal structure of the standardized metric of 'rapidity-of-change' comparing 'core' and 'ring' dynamics was insigni cant (Figure2, panel c).The only signi cant coe cient (p < 0.05) was detected at lag 0, documenting a positive and intense synchronic correlation using both Pearson r (0.70) and Spearman r s (0.68) coe cients.Investigating local development paths with a Continuous Wavelength Transform algorithmEmpirical results of a Wavelet analysis are illustrated in Figure3using biplots that report spectrums of longterm and short-term temporal patterns of the standardized measure of rapidity-of-change, re ective of development dynamics respectively at 'core' (left panel) and 'ring' (middle panel) districts.The highest part of the Continuous Wavelength Transform (CWT) biplot indicates a generalized acceleration of developmental dynamics along the investigated time interval in both districts, more evident in central years (i.e. the 1980s) for The differential time pattern between central city and the ring observed for service, transport, and economic performance indicators, suggests how economic development was spatially heterogeneous and mostly divided in urban and rural areas (e.g.Couch et al., 2007;Shen and Wu, 2020;Bailey, 2021).Infrastructural improvements were delayed by nearly one decade in the ring district(Morelli et al., 2014).An expanding tertiary sector, and adjustments in the cost of energy and transportation, made the socioeconomic structure of 'core' and 'ring' districts more similar and balanced since the 1980s (e.g.Vaiou, 1997;Maloutas, 2007;Di Feliciantonio and Salvati, 2015), suggesting how urbanization has indirectly fuelled exurban development during the whole study period(Di Feliciantonio et al., 2018;Zambon et al., 2019;Salvati, 2022).
Hale and Moberg, 2003; Morya and Ram, 2020; Delmelle et al., 2021), providing a graphical outcome as a function of the development path of the two neighboring areas investigated.

Table 1 .
List of variables used in the present study.

Table 2 .
(Bove and Di Ciaccio, 1994nd partial autocorrelation coe cients by lag and district (only signi cant coe cients at p < 0.01 were reported).Pearson) and non-parametric (Spearman) coe cients (i.e.considering together sign, intensity, and signi cance) was run on rst-differentiated time series with the aim at testing linear (or nonlinear) relationships among the individual dimensions of local development.Pearsonr and Spearman r s coe cients range between -1 (negative association) and +1 (positive association) with values close to 0 indicating an uncorrelated pattern.The subsequent analysis assumed a linear process of local development, with linear interactions among individual components of metropolitan growth.Statistical signi cance (testing MFA) deals with data matrices in which a set of individuals is described by several sets of variables with similar characteristics, e.g.being completely quantitative(Lavit et al., 1994).More speci cally, MFA is a multivariate exploratory framework investigating the latent relationship characteristic of a data matrix with 3 (or more) ways(Bove and Di Ciaccio, 1994), e.g. in a simple three-way case, 'units • variables • times' arrays such as: against the null hypothesis of no correlation between indicators) was set up at p < 0.05.Concordance in both sign and intensity of Pearson r and Spearman r s coe cients indicates a linear, pair-wise correlation between indicators.Table3reports the main results of this analysis documenting how every signi cant pair-wise correlation detected within the indicator's set is linear.This allows application of a traditional multi-way analysis of developmental paths in both 'core' and 'ring' districts based on the classical decomposition of correlation matrices via Pearson metric.Table3.Correlation matrix (using both parametric (Pearson) and non-parametric (Spearman) pair-wise coe cients) by district; signi cance tested at p < 0.05 after Bonferroni's correction for multiple comparisons).DruInvestigating 'fast' and 'slow' variables of local development using Multi-way Factor Analysis Multi-way Factor Analysis (

Table 4 .
Loading coordinates to signi cant axes of a Multi-way Factor Analysis (MFA) by metropolitan district in the study area.