Detecting Tourism Typologies of Regional Destinations Based on Their Spatio-Temporal and Socioeconomic Performance: A Correlation-Based Complex Network Approach for the Case of Greece

Within the context that tourism-seasonality is a composite phenomenon described by temporal, geographical, and socio-economic aspects, this article develops a multilevel method for studying time patterns of tourism-seasonality in conjunction with its spatial dimension and socio-economic dimension. The study aims to classify the temporal patterns of seasonality into regional groups and to configure distinguishable seasonal profiles facilitating tourism policy and development. The study applies a multilevel pattern recognition approach incorporating time-series assessment, correlation, and complex network analysis based on community detection with the use of the modularity optimization algorithm, on data of overnight-stays recorded for the time-period 1998–2018. The analysis reveals four groups of seasonality, which are described by distinct seasonal, geographical, and socio-economic profiles. Overall, the analysis supports multidisciplinary and synthetic research in the modeling of tourism research and promotes complex network analysis in the study of socio-economic systems, by providing insights into the physical conceptualization that the community detection based on the modularity optimization algorithm can enjoy to the real-world applications.


Introduction
A variety of theories and models have been developed over time aiming to study the uneven developmental dynamics that emerge between regions and to contribute to the academic dialogue towards the direction of seeking a solution to the regional problem that concerns the reduction of spatial inequalities [1][2][3][4][5][6][7][8][9]. In general, the uneven spatial development affects the economic growth [2,3], the quality of the environment, the culture, and the social conceptualization [4][5][6], and, therefore, it can be considered as a major driver for the sustainable economic development [7,8]. This aspect highlights the importance of studying spatial inequalities in the research of regional economics [1,9], where tourism is a major field due to its socio-economic configuration that affects the regional, national, and worldwide economic interaction and thus the economic development [8][9][10]. Provided that tourism is a time-variant socio-economic phenomenon that is by default equipped with a spatial reference, the study of tourism within a common temporal and spatial context [6][7][8] can provide insights towards the demand for integration, which can contribute to the multidisciplinary conceptualization of regional economics [9,11]. Amongst the many others, a major property of tourism [12][13][14][15][16], which diversely affects almost every destination [9,17], Despite that seasonality affects almost every tourism-destination (column C1, Table A1, Appendix A), in the current literature, the symbiotic relationship between seasonality and space does not seem to be studied in a comprehensive quantitative context. This observation is verified by the bi-dimensionality (column B1, Table A1, Appendix A) of the relevant research that was previously emerged [9], one pillar of which deals with the study of multivariate time (temporal) aspects of seasonality [6,14,17,[30][31][32][33]45] that is mainly applied within the context of a single destination, whereas the second one, which studies geographical aspects of seasonality, overlooks the effect that time induces to the spatial patterns of tourism-demand [42,43]. Some bright exceptions to this dual conceptualization of seasonality can be found in the works of [4][5][6]9,17,44,49], which attempt to integrate the spatial and temporal dimension of this phenomenon through different methodological aspects, such as cluster analysis, econometric modeling, tourist area life cycle assessment, classification and dimension reduction modeling. However, such approaches rather highlight the demand for integration than adequately configure a comprehensive context for integrating the time (temporal) and geographical (spatial) dimensions of tourism-seasonality [9].
Aiming at serving the demand for integration, this paper proposes a new method for studying temporal patterns of tourism-seasonality in conjunction with their spatial information. The proposed method builds on a joint consideration of correlation and complex network analysis to convert into a complex network the temporal information of the seasonality patterns of the 51 (NUTS III) Greek regions, which are configured for the period 1998-2018, and then to classify these available patterns into communities, according to their level of connectivity in the network model. Complex network (or network) analysis is considered by many scholars as a separate discipline, called network science [51][52][53][54], uses the network paradigm to represent complex communication (or other interacting) systems as graphs, namely as pair-sets of nodes and edges that are converted into connectivity matrices, where many topological measures can be computed. In tourism research, network structures prevail [55][56][57][58] because tourism is by definition an act of interaction (or mobility) and thus the conceptualization of connectivity is by default applicable in the modeling of tourism systems. However, the network paradigm has not yet been sufficiently incorporated in tourism research and therefore there are many ongoing avenues of further research towards this direction [55][56][57]. For instance, a major benefit of the network conceptualization in the description of interconnected systems regards its ability to incorporate, within a single model, structural and functional (flow) information, at any scale (e.g., local, neighborhood, regional, or global), and thus to provide a dual microeconomic and macroeconomic description of real-world systems [58]. Within this context, this paper uses the network paradigm to model correlations between seasonality patterns as a complex network, which (according to the dual conceptualization of networks) allows incorporating topological information related to the connectivity of correlated data with the scale and variability information captured by correlation analysis. The study examines the case of Greece, which enjoys a considerable tourism specialization, where almost 20% of its GDP is related to direct tourism activities [46][47][48][49][50][51], it is a top-ranked in the EE27 in terms of the number of visitors [46], amongst which the 50-60% of foreign arrivals is implemented in the high-season period of July to September [47]. Greece has more than 16,000 km of coastline and more than 1350 islands, islets, and rocky-islands [47], and its rich geomorphological configuration has led to a composite tourism product [50], which is unevenly distributed throughout the various tourism-destinations of the country [49][50][51].
Within this context, the proposed methodological framework allows classifying the Greek prefectures according to their seasonality patterns, aiming to be insightful about the spatial configuration of tourism-seasonality in Greece and how it can be managed in conjunction with the dynamics of the regional space. The main research question addressed in this paper is configured within an exploratory context and it concerns whether temporal patterns of tourism seasonality can be related to the spatial distribution and the socioeconomic attributes of the regional destinations in Greece. The overall approach can drive to more efficient and sustainable management, planning, and hospitality strategies in tourism, which can contribute the Greek tourism-destinations to the goal of balanced regional development.
The rest of this paper is structured as follows: Section 2 describes the methodology of the study, the available data, and the variables that are included in the analysis. Section 3 shows the results of the analysis and provides a discussion within the context of tourism management and regional economics. Finally, in Section 4 the conclusions are summarized.

Methodology and Data
This study configures a methodological framework for analyzing temporal patterns of tourism-seasonality amongst the 51 Greek regions and for classifying them into regional groups, according to their relevance in connectivity that is defined after applying complex network analysis to the seasonal data. The purpose of the analysis is to detect linear and non-linear commonalities in the temporal patterns of tourism-seasonality of the Greek regions, to categorize them into groups, and to detect distinctive seasonal profiles. Towards this direction, the proposed methodology builds on complex network analysis [51][52][53][54] that is applied to the seasonal data after converting them to graphs through a correlation-based criterion. In particular, the methodological framework of the study consists of the five steps shown in Figure 1. At the first step, the seasonality-variables expressing the seasonality patterns of the Greek prefectures are computed, on data expressing the monthly number of total overnightstays per prefecture (NUTS III), which are recorded for the years from 1998 up to 2018. The available data were granted upon request by the Hellenic Statistical Authority [59] to be used under an exclusive license for this study. In particular, 51 seasonality-variables are configured, each one representing the time-series of the monthly overnight-stays recorded per prefecture (codes and variable-names are shown in Table A2, in the Appendix A).
At the second step, the seasonality-variables are first included in a correlation analysis, where the Pearson's bivariate coefficients of correlations [60] are computed for all possible pairs of the available variables. The square table of correlation coefficients is then filtered according to the level of significance, where only significant coefficients are preserved. Therefore, for significances 1%, 5%, and 10%, this allows producing three correlation tables corresponding to a certain significance level. At next, each table is converted to a network (graph) model, each one expressing a weighted graph of 51 nodes, where the links have weights equal to the significant correlation coefficients, for this case.
Further, a community detection analysis [61] that builds on the modularity optimization algorithm of [62] is applied to each network model (corresponding to 1%, 5%, and 10% levels of significance). This algorithm generates communities by maximizing the intracommunity connectivity (i.e., the connectivity of the nodes included within a community), which simultaneously minimizes the inter-community (i.e., between the communities) connectivity, as is described by the mathematical expression [61]: where Q is the modularity-function [61] and m is the number of network edges (links). The community detection [61,62] is implemented through a greedy algorithm that is applied in two steps [54,61]. At the first step, all nodes in the graph are assigned to different communities. Afterward, the nodes are step-by-step being swept and placed into collectivecommunities if the assignment of a node into an adjacent community increases the weighted modularity function of the original model. At the second step, the collective-communities are replaced by new (super-) nodes, and the process is repeated up to the point that the modularity function cannot increase anymore [61]. Generally, greedy algorithms are sensitive to regional conditions, where the algorithm can be "trapped" in the area of a local-optimum and thus cannot move to a better (global) solution [54,61,62]. Therefore, the generated communities do not a priori describe the (globally) optimum partition of the graph-model into communities but a quite satisfactory division, which ensures that nodes within the communities are well-enough densely connected and well-enough sparsely connected between the communities. According to the comparison provided by the author of [61], amongst 12 in number modularity optimization algorithms, the algorithm of Blondel et al. [62], which is chosen in this paper to run the community detection analysis, belongs to the pair of cases with the smallest order of computational complexity O(m). This property gives this algorithm a privilege to be considered as one of those with optimum performance and, also, its availability to be built-in in the open-source software of [55] makes it a popular choice for community detection analysis. However, the community classification outcome is inevitably dependent on the algorithm that is chosen to perform the analysis, and therefore the results of the community detection cannot be considered unique. For a detailed review about community detection in networks and the optimization methods, the interested reader is referred to the reference work of [61]. Within this context, to manage the restriction due to the greedy architecture of the modularity optimization algorithm, which brings uncertainty both at the step of the algorithm's choice (amongst several others) and at the step of applying the specific algorithm to the network dataset (where the initialization of the algorithm is by default submitted to randomness), we apply some repairing techniques, as follows: First, the community detection analysis is designed to be applied to three correlation-based complex networks (i.e., those converted at 0.01, 0.05, and 0.10 level of significance tables) instead of one, so that the final communities to result by the synthesis of these three aspects incorporating different levels of uncertainty and variability in network connectivity. Secondly, we run the procedure of community detection for a single significance-network several times and we then configure the final community outcome by the mode (i.e., most frequent solution) of all runs. Finally, the results of the analysis in this paper are compared with those obtained in previous work of [9], which applied a grouping to the available seasonal patterns of the Greek prefectures based on the principal components analysis method. Interpreting the results of the comparison in terms of regional dynamics and tourism development will provide an evaluation of the proposed in this paper methodology. All these three repairing approaches are expected to reduce the sensitivity of the modularity optimization algorithm due to its greedy architecture. At the third step of the methodology, the different patterns referring to the modularity communities resulted for 1%, 5%, and 10% levels of significance are jointed to a single classification of seasonal modularity-patterns. This joint-classification takes into account the unique modularity-community patterns of the form (Q 1% (i), Q 5% (i), Q 10% (i)), where each node (i) belongs to.
At the fourth step, the patterns resulted from the classification of the seasonality modularity-patterns are evaluated by their geographical and socio-economic features that are approximated by a set of 51 variables, as shown in Table A3 in the Appendix A. This analysis is based on the formulation of error-bars expressing 95% confidence-intervals (CIs), which are computed for the mean values [54,60] of the available socio-economic and geographical (SEG) variables amongst the groups referring to the available seasonal modularity-patterns. When the error-bars do not overlay, the mean values referring to modularity-groups are considered as statistically different, under a 95% certainty [9,54,60]. Thus, comparisons amongst the modularity-groups applied for every SEG variable (Table A2, Appendix A) can reveal which group outperforms (i.e., has a maximum performance) or underperforms (i.e., has a minimum performance) in terms of a single SEG attribute. Within this context, the comparison-results are further tabulated to formulate the total profiles (including seasonality, socio-economic, and geographical information) of the available modularity-groups [9]. This part of the analysis develops linkages between the temporal (as expressed by the time-series patterns of overnight-stays on which the analysis is built), socio-economic, and geographical dimension in the analysis of tourism-seasonality and thus to contribute to the conceptualization of integrating modeling of the complex phenomenon of tourism-seasonality.
The overall approach provides a tool of quantitative analysis that can be proven useful for the scholars and practitioners activating on regional planning, tourism policy, and hospitality management. The findings and the results of the analysis are discussed at the fifth step of the methodology, where an overall assessment and conclusions are formulated.

Data Visualization
The available 51 seasonality-variables of overnight-stays participating in the analysis are first presented to the diagrams shown in Figure 2. This is done for visually detecting possible similarities amongst the seasonality patterns of the 51 Greek prefectures. Each seasonality-variable refers to a Greek region (prefecture) and has 252 monthly records of tourism overnight-stays (i.e., is of length 252), which are recorded for the period from January 1998 to December 2018 [59]. As is evident, almost all available variables have a distinguishable periodicity, which is evident by the oscillating curves of the line-plots. In these patterns, differences in scale (i.e., the oscillation height) and in trend (e.g., some seasonality patterns are described by an increasing trend) can be observed, which are more distinguishable in the metric rather than in the semi-log scale. In particular, in the metric representation, three scale-levels (zones) can be observed in the overlaid line plots of Figure 2, where approximately the first zone ranges from 0-500,000, the second one from 500,000-2 million, and the third one from 2-4 million (measured in overnight-stays). These levels correspond to different oscillation heights that express different levels of magnitude in tourism overnight-stays. However, in the semi-log representation, these zones observed in the metric scale do not seem to propose a distinguishable grouping, similarly to the metric case. Instead, in the semi-log scale, one concentration of "noisy" (i.e., not following clear periodical) patterns can be observed at the range of 10-20,000, another cluster of distinguishable periodical patterns can be observed at the range of 100-1 million, and a final one at the range of 100,000-3 million (all expressing overnight-stays). The previous visual examination of the available data's time-series curves indicates that pattern recognition on data that are jointly submitted to the temporal and geographical configuration is a task of high complexity, where supervision based on intuition can be easily lost, especially on big datasets. Within this context, this paper applies a multilevel approach based on correlation and complex network analysis to jointly manage the temporal and spatial information describing the available dataset. This multi-level approach deals with complexity by integrating linear and non-linear methods and aims to contribute to pattern recognition in seasonality and generally in tourism research.

Correlation Analysis
In the second step of the analysis, the seasonality-variables participate in a correlation analysis to detect linear relations between the possible variable pairs. The results of the analysis are shown in Figure 3, where the correlation tables for each level of significance (i.e., 1%, 5%, and 10%) are illustrated in heat-maps. In all cases of significance, the analysis first illustrates a distinguishable sub-   (50), on the island of Crete. As it can be observed in Figure 4, the geographical arrangement of these regions is based on the state of the adjacency (where prefectures are mainly located in central and south Greece), have geographical relevance, and their spatial distribution configures a "U"-shaped pattern.  Table A2  Another distinguishable area that can be observed in Figure 3 is the cross-junction area A 2 , which includes the seasonality patterns of the prefectures of Xanthi (5), Thessaloniki (6), Hmathia (7), Kilkis (8), and Pella (9). All these regions have geographical relevance, as they are located in north-east Greece (Figure 4), as well as tourism-seasonality (functional) relevance, to the extent that they are in their majority positively correlated to all the other seasonality patterns of the Greek regions. A final distinguishable area that can be observed in Figure 3 is the cross-junction area A 3 , which includes the seasonality patterns of the neighbor prefectures of Kozani (13), Grevena (14), Kastoria (15), and Florina (16). These regions have also geographical relevance due to their neighboring location in north Greece ( Figure 4) and their functional relevance regards their either negative or insignificant correlations to all other seasonality patterns of the Greek regions.

Classification of Seasonality Patterns Based on Community Detection
As being evident by the previous correlation analysis, the temporal dimension of tourism-seasonality can enjoy a spatial conceptualization, which is expressed by the state of neighborhood or adjacency amongst regions with correlated seasonality patterns. However, pattern recognition based exclusively on correlation analysis cannot be that insightful because the geographical patterns that emerged by correlations are quite primary and selective, to the extent that they do not include the total cases of the available population. Within this context, the correlation tables of Figure 3 are converted to graph-models, where a link expresses the value of a significant correlation coefficient between two nodes (prefectures). These graph-models are afterward divided into communities, by applying the modularity optimization algorithm [54,61,62], which generates communities with dense connectivity within and sparse connectivity between them.
The results of the community detection are shown in Figure 5, where regions are colored according to the community they belong to. As it can be observed, for the 1% and 5% level of significance the graph-models are divided into 3 communities, whereas for the 10% level of significance the graph-models are divided into 2 communities. Communities for the graph-models of 1% and 5% levels of significance include the same number of members, the first one consisting of 13 prefectures, the second one of 6, and the third (and bigger) one of 32 elements. On the other hand, the first community for the graph-model corresponding to the 10% level of significance includes 18 members (prefectures), whereas the second one includes 33 elements. As it can be observed, the community grouping corresponding to the case of 10% significance level results by merging the 13-element and 6-element communities of the 1% and 5% significance levels and by shifting one element, the prefecture of Kozani (13), to the 32-element community of the 1% and 5% significance levels. This observation implies first that the 13-element and 6-element communities of the 1% and 5% significance levels become relevant when uncertainty increases (and thus they are linked in a weak-ties context) and secondly that the prefecture of Kozani (13), shows an interesting outlier performance due to its community shift, suggesting a special case to study.
To evaluate the information provided by the community detection of the available (interregional seasonal correlation) graph-models and to manage the results of the community detection algorithm within a context of repairing their heuristic configuration, a joint assessment is applied to the available modularity classes (communities). In particular, based on the modularity class that each prefecture belongs to, we construct ordered triplets of the form p = (p 1 ,p 2 ,p 3 ), where the first coordinate (p 1 , , ) expresses the modularity class for the 1% level of significance, the second one ( ,p 2 , ) for 5% sig., and the third one ( , ,p 3 ) for 10% significance. Therefore, based on the available ordered triplets p(i), which correspond to the i = 1, . . . , 51 modularity (seasonal) classes of the Greek prefectures, four (4) in number unique patterns (clusters) can be distinguished, which are expressed by the triplet-groups (0,0,1), (1,1,0), (1,1,1), and (2,2,0), as is shown in Figure 5. Amongst these patterns (clusters) of modularity classes, the first one configures a sigmoid ("S"-shaped) spatial pattern (Figure 6b) that consists of the prefectures of Drama (2), Xanthi (5), Thessaloniki (6), Imathia (7), Kilkis (8), Pella (9), Serres (11), Kastoria (15), Ioannina (17), at northern Greece, the prefectures of Karditsa (22), Trikala (24), and Viotia (26), at central Greece, and Chania (51), at the island of Crete. The second one includes a single prefecture, Kozani (13), which is located in north Greece (Figures 5 and 6) and thus it configures a point or dot ("•"-shaped) spatial pattern (Figure 6c) at the central Macedonia region. As is evident, this community consists of the single element with the most unstable connectivity, as denoted by the shifting captured by the intersection consideration of communities between the 1%, 5%, and 10% significance levels. The third cluster configures a linear spatial pattern (Figure 6d) consisting of the prefectures of Grevena (14), Florina (16), at north Greece, Arta (18) and Evrytania (28), at central Greece, and the prefecture of Arkadia (37), at the region of Peloponnesus.  According to the community detection analysis, the modularity communities can relate the temporal with the spatial dimension of tourism-seasonality within a more perplexing context than this shown in Figure 4, where the modularity groups (p 1 , p 2 , p 3 ) that were generated configure diverse but also distinguishable spatial patterns in the geographical map of Figure 6. Within this framework of complexity, to further decompose and examine the multiplex information included in these modularity groups, we illustrate the lineplots of the overnight stay time-series per modularity group, as is shown in Figure 7. As it can be observed, modularity groups (0,0,1) and cluster (2,2,0) include seasonality patterns with more discrete periodical patterns, whereas groups (1,1,0) and (1,1,1) are described as more noisy patterns. In particular, the first modularity group (0,0,1), which corresponds to the "S"-shaped spatial pattern of Figure 6, includes cases that are described by a periodic pattern with an increasing trend. The case included in the second modularity group (1,1,0), which refers to the "dot"-shaped spatial pattern of Figure 6, is described by a noisy seasonal pattern with a decreasing trend, where linearity is more obvious in the time-series curve than periodicity. Also, the third modularity group (1,1,1), which corresponds to the "I"-shaped spatial pattern of Figure 6, includes cases that are described by a noisy periodical pattern that configures (on average) a bell-shaped curve with its convex part defined at the period 2003-2012. Finally, the fourth (and larger) modularity group (2,2,0), which corresponds to the "O"-shaped spatial pattern of Figure 6, includes cases that are described by an approximately cyclical pattern (i.e., a periodical with the same oscillation amplitude) that configures (on average) a U-shaped curve with its concave part also defined at the period 2003-2012. Overall, this analysis allows claiming that the four available modularity groups have, along with their spatial patterns, distinguishable seasonality patterns, on average.  Figure 6) to the triplet-groups: (a) (0,0,1) of the "S"-shaped spatial pattern, (b) (1,1,0) of the "dot"-shaped spatial pattern, (c) (1,1,1) of the "I"-shaped spatial pattern, and (d) (2,2,0) of the "O"-shaped spatial pattern.
To detect differences at the scale of seasonality amongst the available modularity groups, we apply further analysis based on error-bars of confidence intervals (CIs) [54,60]. In particular, we construct 95% CIs for the mean values of the Gini coefficient [9,37] and the RSI (Relative Seasonal Index) [14,36], which are composite indicators capturing aspects of seasonality on a set of data. The results are shown in Figure 8, where it can be observed that the modularity group (2,2,0), the one with the "O"-shaped spatial pattern, has on average the maximum score of seasonality, which is statistically different than these of all the other groups. On the other hand, the modularity group (1,1,0), described by the "dot"-shaped spatial pattern, has on average the minimum score of seasonality, which is statistically different than this of (1,1,0) and (2,2,0) groups. However, no statistical difference can be detected between the (0,0,1) and (1,1,1) groups, implying that their seasonality patterns lie on the same scale level. Overall, this analysis highlights that community detection resulted from the modularity optimization algorithm sufficed to configure communities separating the minimum and maximum performance in tourism-seasonality, as it can be captured by the Gini coefficient and the Relative Seasonal Index.

Socio-Economic Determination of the Modularity Seasonal Groups
In the fourth step of the methodology, the available modularity groups are interpreted in terms of their geographical and socio-economic information. To do so, the analysis builds on comparisons amongst error-bars of 95% confidence intervals, which are computed for a set of 51 (in number) available variables that are divided into thematic categories (geographical, seasonality, transport, demographic, productivity, tourism, environmental, and cultural), as is shown in the Appendix A (Table A3). Based on the error-bars comparisons, the statistically significant minimum and maximum performances per variable are tabulated and shown in the Appendix A (Table A4). According to this analysis, the prefectures belonging to the modularity group (0,0,1) have high urbanization, winter tourism activity specialization, the capacity of "Natura" (i.e., biodiversity protected) areas, and the capacity of woodland parks. On the other hand, the prefectures included in this group have not rich geomorphological configuration, relatively low seasonality performance, and a small number of ports, are not specialized in the secondary sector, not equipped with industrial areas, and have a small number of camping areas and beaches. The prefecture of Kozani (13), which configures the modularity group (1,1,0), is located in northwest Greece, is a mainland prefecture with rich geomorphological configuration, high income, tourism specialization, and specialization in winter tourism activities.
On the other hand, this single-member group has low seasonality, tourism capacity in terms of accommodation, and level of environmental and cultural resources. The prefectures belonging to the modularity group (1,1,1), which are mainly located on the mainland, are cases with high urbanization, but with low geomorphological configuration, seasonality, sea geomorphology, and tourism profile, as well as demographic and income background. Finally, the prefectures belonging to the modularity group (2,2,0) are mainly located in southeast Greece, where their sea geomorphologic configuration is dominant. Prefectures of this group have high seasonality, urbanization, primary and tertiary sector specialization, tourism profile, and environmental and cultural resources. On the contrary, cases of this modularity group have low mountainous geomorphology. All previous interpretation is tabulated in Table 1, where the socio-economic and geographical profiles of the modularity groups are illustrated. Table 1. The semiology of the modularity groups (see Figure 5) resulted from the analysis.
Coastal or island area; seasonality; poor road density and roads; low primary and tertiary sector specialization; low tourism profile; low environmental wealth, low cultural resources profile.
Coastal or island area; poor geomorphological configuration; seasonality; poor ports and airports configuration; population and human capital; low income; low secondary sector specialization; low tourism profile; low beach environment; low cultural resources profile.

Group (2,2,0) (32 prefectures)
Southern and eastern location; coastal or island area, high seasonality; rich road density and ports configuration; high urbanization; high primary and tertiary sector specialization; high tourism profile; high environmental quality; high capacity of cultural resources.
Poor geomorphological configuration; poor rail configuration; low secondary sector specialization; low mountainous activities.
(a) as defined by the max values of As it can be observed from Table 1, the first modularity group (0,0,1), which is a group with a sigmoid geographical distribution and a seasonal pattern described by periodicity and an increasing trend, is a group where its relatively low seasonality pattern seems to be driven by forces related to winter tourism activities and MICE tourism (i.e., meetings, incentives, conferencing, exhibitions) attractions due to urbanization configuration. The second modularity group (1,1,0), concerning a single-prefecture group described by a noisy low seasonality pattern with a decreasing trend, has its seasonality related to the mainland geomorphology, high rail, and airport integration (but poor road density), and a mountainous tourism specialization. The third modularity group (1,1,1), which configures a linear geographical pattern and has a noisy periodical pattern with a bell-shaped curve running a stage of decline, is a group owing to its seasonal pattern to high environmental wealth and forces related to tertiary specialization and low transport integration and tourism profile. Finally, the fourth modularity group (2,2,0), which is the largest group described by a circular geographical configuration and an approximately cyclical pattern of seasonality, has its seasonality related to its rich sea geomorphology, environmental, and cultural wealth, as long as transport integration. This group is mainly related to the 3S (i.e., sea, sun, and sand) tourism activities that are the dominant model describing the nature of the tourism product in Greece.

Further Analysis and Overall Assessment
To further evaluate the effectiveness of the overall correlation-based complex network analysis applied in this paper, we compare the results of the modularity grouping with those presented in a previous work [9] of the authors, which were produced by applying a Principal Components Analysis (PCA) method to the same overnight-stays database. This comparative approach, on the one hand, is expected to provide insights into the uncertainty due to the modularity optimization algorithm and, on the other hand, to evaluate the effectiveness of the method proposed in this paper. In particular, the authors of [9] applied a PCA to the same set of seasonality-variables (see Table A2 Appendix A) by choosing the covariance matrix extraction and the Varimax rotation methods. The PCA resulted in seven (7) principal components explaining a high percentage (84.86%) of the total variance.
Provided that the PCA's component matrix includes component loadings (which are defined by correlations) between the principal components and all the available variables, further min-max filtering was applied in this approach to singularly assign one variable in a unique PCA component. The PCA yielded two groupings of the available 51 Greek prefectures, where the first one (which was generated by considering only the maximum -M-component loadings per variable) consisted of five (5) groups including 38 (PCM#1), 7 (PCM#2), 4 (PCM#3), 1 (PCM#4), and 1 (PCM#5) cases (prefectures) respectively, while the second one (which was generated by considering only the minimum -m-component loadings per variable) consisted of seven (7) groups including 1 (PC#1), 7 (PC#2), 13 (PC#3), 14 (PC#4), 8 (PC#5), 5 (PC#6), and 3 (PC#7) cases, respectively. The labeling of both PCA groupings was according to the ranking of the principal components resulted from their eigenvalues shown in the scree plot.
Within this context, Table 2 is constructed for the comparison and measurement of the relevance between the modularity grouping (MOD) applied to this paper with the PCA grouping applied in the previous paper [9] of the same authors. This comparison table has a multilayer block-structure consisting of a window (block) named "Intersection Frequencies" that displays the set-theoretic intersections PC#i ∩ MOD#j between the i-th PCA and j-th modularity groups, two other blocks named "Relevance to PCA/MOD Groups" that display the percentage (%) of coverage that the PCA/MOD groups have for each intersection, and a final window displaying relevance differences that imply which (PCA/MOD) group is more relevant (dominant) to each intersection.
Next, according to the relevance-assessment shown in Table 2, the 6 cases (prefectures) composing the complementary set C#1 = PCM#1−PCM#1 ∩ (2, 2, 0) participate in an amount of 38.5% in the configuration of the modularity group (0,0,1), 20% in the configuration of the modularity group (1,1,1), and 13.2 + 2.6 = 15.8% in the configuration of the PCM#1 group. These shares imply that set C#1 is more relevant in the configuration of the modularity than the PCA groups and therefore they can provide another indication of better performance of the MOD than the PCA grouping.
Finally, as is evident by the block-structure of Table 2, the PCA grouping is by default described by a higher level of complexity than the MOD grouping because of the dual (min-filtering and max-filtering) structure describing the PCA groups. Therefore, in terms of resolution of the grouping outcome, the modularity grouping can be considered more reliable than the PCA's because it results in a singular instead of dual grouping outcome. Moreover, amongst the two PCM and PCm groups, the modularity grouping appears more relevant to the max-filtering (PCM) than to the min-filtering (PCm) outcome. This relevance is obvious in terms of groups' number (i.e., the number of MOD groups is closer to the number of PCM than PCm groups), but it can be further supported by the rank-size distribution of the groups, which follows a power-law pattern in the case of PCM and MOD groupings whereas a Gaussian-like for the PCm grouping.  Frequencies) between i-th PCA and j-th modularity groups, with j = (2,2,0), (0,0,1), (1,1,1), or (1,1,0), their percentage (%) relevance to the PCA/MOD groups (window: Relevance to PCA/MOD Groups), and a metric test (window: Differences in Relevance) implying which (PCA/MOD) group is dominant for each intersection.

Modularity Group
Group Size 32 13 5 1 32 13 5 1 32 13 5 1 32 13 5 1 Group Label (2,2,0) (0,0,1) (1,1,1) (1,1,0) (2,2,0) (0,0,1) (1,1,1) (1,1,0) (2,2,0) (0,0,1) (1,1,1) (1,1,0)  In spatial terms, Figures 9 and 10 show the geographical distribution of the maxfiltering (PCM) and min-filtering (PCm) groupings resulted from the Principal Components Analysis being available in the previous work of [9]. As it can be observed in Figure 9, group PCM#1 configures an "O"-shaped pattern, PCM#2 a "U"-shaped pattern, PCM#3 an "I"-shaped, and PCM#4 and PCM#5 "dot"-shaped patterns. On the other hand, group PCm#1 configures a "dot"-shaped pattern, PCM#2 an "O"-shaped pattern, PCM#3 and PCM#4 "U"-shaped patterns, PCM#5 an "S"-shaped pattern, and PCM#6 and PCM#7 "I"-shaped patterns, as is shown in Figure 10. Between the PCM and PCm groupings, the geographical distribution of the former appears more relevant to this of the MOD grouping, whereas PCm appears less relevant and messier. On the one hand, this relevance (that is observed between the spatial patterns of the PCM and MOD groupings) is consistent with the quantification shown in Table 2 and generally implies that different clustering techniques can lead to converging outcomes when the system being examined is described by considerable underlying causality. On the other hand, measurable differences that are captured between different grouping outcomes can quantify the reliability of each clustering technique. Within this context, the modularity grouping proposed in this paper appears overall more reliable than the PCA, in terms of grouping singularity, size, and configuration. This added value can be due to the complex-network approach applied in this paper, which incorporates topological information related to connectivity of interconnected data instead of just considering their scale and variability. that are computed in the paper of [9], where (a) is the aggregate group, including all 5 available PCM groups, (b) is the PCM#1 group of an "O"-shaped spatial pattern, (c) is the PCM#2 group of a "U"-shaped spatial pattern, (d) is the PCM#3 group of an "I"-shaped spatial pattern, (e) is the PCM#4 group of a "dot"-shaped spatial pattern, and (f) is the PCM#5 group of a "dot"-shaped spatial pattern. that are computed in the paper of [9], where (a) is the aggregate group, including all 7 available PCm groups, (b) is the PCm#1 groups of a "dot"-shaped spatial pattern, (c) is the PCA#2 group of an "O"-shaped spatial pattern, (d) and (e) are the PCA#3 and PCA#4 groups of a "U"-shaped spatial patterns, (f) is the PCA#5 group of an "S"-shaped spatial pattern, and (g) and (h) are the PCA#6 and PCA#7 groups of an "I"-shaped spatial patterns.
However, despite these indications of better performance describing the proposed complex-network-based method in comparison to the PCA grouping, a study that builds on clustering, numerical analysis, and heuristic optimization is inevitably submitted to uncertainty due to the modeling choices and computational constraints. In this paper, the number and the configuration of the generated groups were a result of optimization based on a heuristic algorithm, where iteration suggests a prime repairing technique to consider and manage the variability of results. Towards this direction, several avenues of further research can emerge. One may concern examining the convergence of the results due to iterations implemented within a single grouping method, such as the modularity grouping proposed in this paper. In this avenue, the characteristic parameters (e.g., degrees of freedom, variability, scale, etc.) describing the dataset and those describing the grouping outcome (e.g., number, scale, and variability of groups) can be correlated. Another avenue is examining the convergence between different grouping methods, similarly to the comparisons between the MOD and PCA groupings applied in this paper. In this approach, the variability observed in the grouping results (e.g., in the number and size of groups) should be evaluated in the light of the computational mechanism describing each method, which can be based either on variance, or distance, or scale, or positioning, or similar relevance-like conceptualization. Despite such computational-based approaches, avenues for further research can also emerge at the level of implementation, and particularly they can concern the management of the resulting groups in terms of tourism policy and regional development. For instance, an effective clustering on temporal patterns of tourism demand can be insightful for the tourism policy-makers to better distribute the national budget of tourism development across the regions based on their relevant profiles or to develop more efficient tourism developmental strategies that would be customized to the tourism profile of each group. Overall, the limitations that are applicable in this study can excellently suggest avenues of further research that can contribute to the integration of the spatial, temporal, and socioeconomic aspects in tourism research.

Conclusions
This paper developed a multilevel methodology for analyzing time (temporal) patterns of tourism-seasonality and for grouping them into regional classes and therefore to configure distinguishable seasonality profiles facilitating tourism policy and development. The proposed method was built on a multilevel pattern recognition incorporating time-series assessment, correlation (linear assessment), and complex network analysis (non-linear assessment) to deal jointly with the temporal and spatial dimension of tourism-seasonality and to interpret the results of the analysis in socio-economic terms. The available data were registrations of overnight-stays for the years from 1998 up to 2018, which were provided by the Hellenic Statistical Authority [59] under an exclusive license. The analysis revealed four groups of seasonality, which were configured by a connectivity conceptualization based on the modularity optimization algorithm, where each group showed to be described by distinguishable spatial and seasonality patterns and to enjoy distinct socio-economic semiology. In particular, the first modularity group has a relatively low seasonality, with a pattern described on average by periodicity and an increasing trend, and a sigmoid geographical distribution. The seasonal pattern of this group is driven by forces related to winter tourism activities and MICE tourism attractions. The second, single-prefecture, modularity group has the lowest seasonality and is described by a noisy pattern with a decreasing trend, where its seasonality is related to rich geomorphological mountainous configuration and a mountainous tourism specialization. The third modularity group has a relatively low seasonality with a pattern described by a linear geographical distribution and on average by a noisy periodical configuration with bell-shaped time evolution. The seasonal pattern of this group is driven by forces related to high environmental wealth and tertiary specialization but low transport integration and tourism profile. Finally, the fourth and largest in size modularity group has the highest seasonality, with a pattern described on average by an approximately cyclical pattern and a circular geographical configuration. The seasonal pattern of this group is driven by forces related to rich sea geomorphology, environmental and cultural wealth, and transport integration, all of which describe the dominant 3S model of the tourism product in Greece. Overall, the analysis supports multidisciplinary and synthetic research in the modeling of tourism-seasonality, which is a composite socio-economic phenomenon of high complexity, and also promotes complex network analysis in the study of socio-economic systems by providing insights into the physical conceptualization of community detection in real-world applications.  9  Table A3. The socio-economic and geographical (SEG) variables * used to determine the socio-economic profiles of the modularity groups of Figure 5. SEG.20 C SEC The specialization (% of the GDP) in the tertiary (C) sector.
[66, 69] SEG.21 TOURISM GDP The specialization (% of the GDP) in the tourism sector. [66] SEG.22 TILLING LAND The proportion of the tilling-land areas to the total regional area. [67] SEG.23 FORESTS The proportion of forest-areas to the total regional area. [67] SEG.24 INLAND WATERS The proportion of the inland-water-areas to the total regional area. [67] SEG. 25 INDUSTRIAL AREA The proportion of the industrial-areas to the total regional area.
[67, 69] SEG.26 LAND AREA The proportion of (non-mountainous) land-areas to the total regional area. [67] SEG.27 SEMI MOUNTAIN AREA The proportion of the semi-mountain areas to the total regional area. [67] SEG.28 MOUNTAIN AREA The proportion of the mountain-areas to the total regional area. [67] SEG. 29 MOUNT ACTIVITIES The number of mountain-activities (e.g., walking paths, mount sports, climb fields). [67] SEG.30 CLIMB FIELDS The number of climb-fields. [67] SEG.31 MOUNT ROUTES The number of mountain-routes.