Spatial Interaction E ﬀ ect of Population Density Patterns in Sub-Districts of Northeastern Thailand

: The north-eastern region in Thailand is the largest in area and population. Its average income per capita is, however, the lowest in Thailand. This phenomenon leads to migration to big cities, which are considered economic centres. We investigated the e ﬀ ect of spatial interaction on the population density pattern in 20 provinces in north-eastern Thailand. Data was obtained from the compilation and preparation of the demographic data of 2676 sub-districts for 2002–2017. A ﬁeld survey was conducted through GPS at educational institutions, hospitals, airports, government o ﬃ ces, and shopping malls. The data was analysed using spatial autocorrelation analysis by a global indicator (global Moran’s I) and a local indicator (local Moran’s I and Getis–Ord G i *). Eight Mueang districts exhibited the high-high (H-H) cluster pattern or hot spot at an increasing yearly rate. In addition, the area with the highest gravity was located near service sources and was found to have the largest population. Moreover, gravity interaction with service sources had a strong positive correlation with migration patterns. Thus, the cluster of areas with the greatest population density is located within the Mueang district in one of the provinces with most service sources, as these places attract people and consequently industrial factories and service trades.


Introduction
The north-eastern region of Isan is the largest in terms of area in Thailand and also hosts the largest population. According to the 12th National Economic and Social Development Plan (2017-2021), the north-eastern region is the economic centre of the Greater Mekong sub-region (GMS) with emphasis on the development of the economic area, tourism, and utilisation of inholdings for the accommodation of the new production bases in the industrial sector and the development of a link from Thailand to Lao and Vietnam. A high speed train project was initiated to connect to the People's Republic of China through Vietnam [1,2]. However, the north-eastern region is still poverty-stricken and features a higher income inequality than other region. This leads to the migration of populations from rural areas to cities, which feature relatively higher economic and industrial growth compared to other regions [3,4], thus, increasing the population density in big cities. The management of urban area expansion and the development of neighbouring towns have, thus, been established to accommodate population growth. Population distribution information is an important variable whose pattern is used at the lowest level of urban development planning of large areas [5,6].
Cluster analysis is one of the methods used to identify patterns in population density. It is also used to drive the creation of economic prosperity, which corresponds to the labour production sector, land, utilities, infrastructure, transportation, and academic institutions [7,8], as support for areas that

Study Area
The study area was the north-eastern part of Thailand (hereafter referred to as the Northeast), a region located on the Korat plateau that looks like a basin, which drains westward from a high altitude into an eastward lowland. The Northeast is bordered by the Mekong River (along the border with Laos) to the north and east, by Cambodia to the southeast, and Sankamphaeng Range south of Nakhon Ratchasima to the south. It has the Phu Phan range as the internal range. As the largest region in Thailand, the Northeast covers 20

Data Acquisition and Preparation
Data from 2676 sub-districts in the Northeast of Thailand for the period ranging from 2002 to 2017 were obtained from the National Statistical Office and Ministry of Interior for use as the database for the analysis of the relationships among population size, public utilities, and public assistance and excluded the illegal population in the analysis. The National Statistic Office obtained the data online

Data Acquisition and Preparation
Data from 2676 sub-districts in the Northeast of Thailand for the period ranging from 2002 to 2017 were obtained from the National Statistical Office and Ministry of Interior for use as the database for the analysis of the relationships among population size, public utilities, and public assistance and excluded the illegal population in the analysis. The National Statistic Office obtained the data online through their network, collected by relevant sub-district administrative organisations. For accuracy, we compared the collected data to the census registration from the sub-district administrative organization. Therefore, the population data were similar when cross-checking the data of the National Statistic Office and the sub-district administrative organisations (SAO) through field surveys and interviews with the authorities. The next step involved the surveyance and data collection of geographic coordinates using global positioning system (GPS) and high-resolution satellite images. The data included five types of public amenities, including schools, hospitals, government offices, airports, and shopping malls. Schools were classified into three categories: large-sized schools with 1500 students or more, medium-sized schools with 500-1499 students, and small-sized schools with less than 500 students. These were used for the analyses of an accessibility and gravity model, which reveals population interactions [29][30][31]. Moreover, data on public utilities, such as main and secondary roads, were also collected to determine factors related to speed, distance, and time according to the standard speed of the Department of Rural Roads. The data were then input into the geographic information system software ArcGIS (ESRI, Canada) for the analysis of road networks and the shortest travel time from sub-districts to the centre of population gravity.

Spatial Autocorrelation Analysis
This study used SAA at two levels: global and local indicators. Global Moran's I was used for global indicator spatial autocorrelation analysis, whereas local Moran's I and Getis-Ord G i * were used for local indicator spatial autocorrelation. The relationship between population density in each of the sub-districts and in 8-16 nearby sub-districts was examined to analyse the continuity weight using Queen contiguity. Using Euclidean distance, the next step involved a selection of a distance metric via the GeoDa software for calculation. The confidence level for the result was 95% [32][33][34].
To analyse the population density using one variable, global Moran's I was used to examine the global indicator spatial autocorrelation. For the years 2002, 2007, 2012, and 2017, the data concerning population density were used to identify the overall patterns of the Northeast in one unit at the regional level by comparing the resulting point value in the position under investigation to other positions using a spatial weight matrix [18,35], as shown in Equation (1).
where I is the spatial correlation based on global Moran's I, x i is the population per area i, x j is population per area j which is adjacent to area i, x is the average population density, n is the total number of sub-districts in the Northeast, and W ij is the spatial weight matrix of distance from the sub-districts' centre point. Local Moran's I was used to analyse the local indicator of spatial autocorrelation. Such statistics were used to measure the population density correlation of areas at the sub-district level and nearby sub-districts in the pattern of agglomeration from a total of 2676 sub-districts, as shown in Equation (2).
where I i is the population density correlation based on local Moran's I, Z i is the deviation from the population per area i which is close to area j population per area the deviation from the is Z j , i, and W ij is the spatial weight matrix of distance from the sub-districts' centre point.
The coefficient of correlation based on Moran's I is between −1 and +1. A positive spatial autocorrelation indicates a significant spatially aggregated distribution of a high or low population density (cluster). In contrast, a negative value indicates a significant spatial difference between the area of each sub-district in terms of population density (disperse), whereas zero indicates a non-significant spatial autocorrelation (random) [36,37]. Local Moran's I illustrates the patterns of population density clusters in four classes: high-high (H-H), low-low (L-L), high-low (H-L), and low-high (L-H) [15,17,18].
Getis-Ord G i * was used to analyse local indicators of spatial autocorrelation for the identification of clusters of data, called hot spot analysis. The areas with the highest and lowest population densities yielded Z-score values that were higher than +1.96, called hot spots, and lower than −1.96, called cold spots, respectively. A positive Z-score indicates a high-value cluster, while a negative one indicates a low-value cluster, determined at a confidence of 95%, Getis-Ord G i * was obtained from Equation (3).
is the distribution of the population density (Z-scores), d is the distance from the central point of one sub-district to that of all neighbouring sub-districts, x i is the population per area i, x j is the population per area j which is close to area i is the spatial weight matrix of the W ij distance from the central point of sub-districts, and n is the total number of sub-districts in the Northeast.

Gravity Centre Model
The gravity model was used to analyse the interaction between the areas that affect the activities for area development. This model is a well-known method of analysing economic phenomena that affect the population. In this work, the clustered districts from the LISA analysis were specifically selected and created mean centre points of sub-districts for applying the model. This study analysed access to amenities (e.g., schools, hospitals, government offices, airports, and shopping malls), which affects population clusters, and the effectiveness of access to the road network. Applying the gravity centre model, the equation of Huff [38,39] and Lakshmanan and Hansen was set and applied [40]. Two factors affecting the analysis included the population in the area and the potential accessibility to the centre sources [14,31,41].
To analyse the gravity index of living in the area surrounding the centre of services with the spatial interaction, the calculation was performed based on the human resources index and the accessibility index using the gravity model, as expressed below.
G ij is the index of access to the service area. PH i is the human resources index and TH ij is the accessibility index. The next step involved the classification of the gravity of access to the service area into different levels using natural breaks (Jenks). The groups with G ij at the same level were grouped together; the variance within the same group was lowest but the variance between groups was highest. This study classified G ij into four groups, namely, higher, high, medium, and low. The human resources index (PH i ) is an attraction factor and used to analyse the population per area, as expressed in Equation (5) [8].
where PH i is the human resources index, P i is the population in the sub-district i, P max and P min are the maximum and minimum of 2676 sub-districts in the north-eastern part of Thailand. The accessibility index (TH ij ) is a factor of travelling to the central sources via the shortest path obtained from the analysis of the spatial network structure. A network analysis tool (closed facility function) was used with time and distance to process the shortest path route for providing the travel time. It was used to obtain data concerning the built road line network to determine limited speed by road type, the centre of the sub-district, and the centre of business quarters, with which the shortest path based on road type for travelling to the centre was determined, as stated in Equation (6).
where TH ij is the accessibility index (probability of area j being visited by the inhabitant of area i). T ij is the drive time measured between the gravity point and sub-district central point of the clustering pattern. λ is the exponent of the drive time (equals 1), and F j is the attractiveness of area j (equals 1).
After the survey and the collection of data on the geographic coordinates of public utilities, the centre of attraction for services and the destination point of travel to reach the central source was determined using mean centre analysis to determine the weight value of all positions of public utilities, which equal 1, to indicate the central source of the centre of city clusters.

Interaction between Gravity Index and Related Factors
The spatial correlation analysis between the gravity index, population density, migration, local Moran's I, and Getis-Ord G i * was conducted as expressed in Equation (7).
where ρ is Pearson's correlation coefficient. x, y are variables, and µ x , µ y are the mean values of variables x and y.
The correlation value is between −1 and +1: a correlation of 0 is an indication of the lack of correlation between variables x and y. A correlation value that is non-zero indicates a correlation between variables x and y. A positive correlation indicates that the two variables change in the same direction, while a negative one means that the variables change in opposite directions. The approach presented here was examined for appropriateness and approved by the Mahasarakham University Ethics Committee in Human Research No. 131/2019.

Spatial Pattern Analysis Using Spatial Autocorrelation
An analysis of the distribution of population using spatial correlation based on global Moran's I for the period ranging from 2002 to 2017 yielded a positive global Moran's I value (p < 0.05) every year. This indicated that the area with a similar density (high or low) was clustered and the value of global Moran' s I increased every year (Table 1). The spatial autocorrelation analysis of the population density based on local Moran's I showed that in the year 2002, most areas had no correlation (2053 sub-districts from a total of 2676 sub-districts). Among the areas that were correlated (p < 0.05), most featured the L-L pattern (485 districts) followed by H-H (97 sub-districts), L-H (17 sub-districts) and H-L (24 sub-districts) with a similar amount, which yielded similar characteristics in the results obtained for years 2007, 2012, and 2017 (Table 2). From 2002 to 2017, the area with the H-H pattern increased every year, whereas the area with the L-L pattern decreased. This indicated an increase in clustered sub-districts with a high population density, while there was a decrease in clustered sub-districts with a low population density. The area with the H-H pattern was mostly observed in Mueng districts and an area roughly resembling Khorat, Mun River, and Chee River basins, whereas the L-L pattern was mostly observed in the mountain range and forests in the west of the Northeast. This area covers the Phetchabun, Dong Phaya Yen, and Dong Sankampang mountain ranges, which drain into the east of the Phu Phan mountain range. The south was the Khao Phanom Dong Rak mountain range, which runs along the entire extent from east to west, and this area is adjacent to the border with Cambodia ( Figure 3).
The data of each province for the period ranging from 2002 to 2007 show that only eight provinces were clustered (from a total of 20 provinces), with the correlation being unequal (p < 0.05) as illustrated in Figure 2. Five provinces were initially found to have the H-H pattern (Buriram, Roi Et, Ubon Ratchathani, Nong Khai, and Surin). Three more provinces were found to have the H-H pattern (Nakhon Ratchasima, Khon Kaen, and Udon Thani), as illustrated in Figures 3 and 4. For the last 15 years, many areas have expanded from the centre of this region, such as Nakhon Ratchasima, being the capital of the lower north-eastern region, the gateway linking the central region and the eastern region. This city is the hub of economic, trade, industrial and agricultural advancements, Khmer cultural tourism, sports, medicine and education, geographically linking the economic advancement with other provinces in the region, including Buriram, Surin, and Ubon Ratchathani. Khonkaen is the capital of the central north-eastern region, which clustered and expanded to nearby districts and created new cities in the area. Udon Thani is the capital of the higher north-eastern region, being the hub of services, agricultural industries, and sugarcane mills and links to Nong Khai, which is the economic hub around the Mekong River, allowing trade across the border. Moreover, the population settled along the main transportation network connecting Nakhon Ratchasima, Khon Kaen, Udon Thani, and Nong Khai.   The analysis through Getis-Ord Gi* indicated that in 2002, clusters of high population density areas or hot spots were in 122 sub-districts. Contrarily, clusters of low population density areas or cold spots were observed in 506 sub-districts. There were 2048 sub-districts that were not clustered (p < 0.05) remaining. These results from 2002 to 2017 also suggested an increase in the number of subdistricts with a hot spot every year, whereas those with a cold spot decreased as illustrated in Table  2. This table shows the result of the analysis based on Getis-Ord Gi* (Figure 4), which corresponded with the result of the analysis based on local Moran's I.

Gravity Centre Index
Two factors were used to analyse the gravity-centre index: gravity index and the travel time it took to reach the central source. An analysis of eight cities with cluster patterns indicates that the sub-   The analysis through Getis-Ord Gi* indicated that in 2002, clusters of high population density areas or hot spots were in 122 sub-districts. Contrarily, clusters of low population density areas or cold spots were observed in 506 sub-districts. There were 2048 sub-districts that were not clustered (p < 0.05) remaining. These results from 2002 to 2017 also suggested an increase in the number of subdistricts with a hot spot every year, whereas those with a cold spot decreased as illustrated in Table  2. This table shows the result of the analysis based on Getis-Ord Gi* (Figure 4), which corresponded with the result of the analysis based on local Moran's I.

Gravity Centre Index
Two factors were used to analyse the gravity-centre index: gravity index and the travel time it took to reach the central source. An analysis of eight cities with cluster patterns indicates that the sub- The analysis through Getis-Ord G i * indicated that in 2002, clusters of high population density areas or hot spots were in 122 sub-districts. Contrarily, clusters of low population density areas or cold spots were observed in 506 sub-districts. There were 2048 sub-districts that were not clustered (p < 0.05) remaining. These results from 2002 to 2017 also suggested an increase in the number of sub-districts with a hot spot every year, whereas those with a cold spot decreased as illustrated in Table 2. This table shows the result of the analysis based on Getis-Ord G i * (Figure 4), which corresponded with the result of the analysis based on local Moran's I.

Gravity Centre Index
Two factors were used to analyse the gravity-centre index: gravity index and the travel time it took to reach the central source. An analysis of eight cities with cluster patterns indicates that the sub-districts were classified at the highest level (dark red). They were located near the service centre. However, some sub-districts with a gravity index at the highest level were found to be located far from the central source of the city. This yielded different results when compared to those sub-districts with a gravity index at a high level (light red), as these areas had a large population size. Ubon Ratchathani in Rai Noi sub-district, Meung District, is such an example; although it had a population size of 25,190 people and a distance of 11,957 m from the centre, which is quite far, it was still classified at the highest level (dark red). In contrast, while the sub-districts in the cities had a population size of 11,291, and the distance from the centre was around 1649 m, they were classified at a high level (light red). The area within a 30 km radius around the centre of the city was designated as the point of maximum accessibility, where business quarters and services are found, eventually leading to a high population density. The population density, however, decreases when such areas are far from the centre. The centres of regional capitals are service hubs, making them high gravity index centres, including Nakhon Ratchasima, Udon Thani, and Khon Kaen, while Buriram, Roi Et, Surin, and Nong Khai are in the centre of the city, and Ubon Ratchathani is in the educational centre area as illustrated in Figure 5. The gravity model can be used to plan economics and locations, such as regional trade agreements near provinces of geographic indicators. Meanwhile, the effects to economic growth are internal, and it is external economies that create jobs, thus labour movement grows.

Interaction of Gravity Index and Factors Related to Population
To examine the interaction gravity index and related factors (population density, migration based on local Moran's I, and Getis-Ord Gi*), this study analysed the area data of the sub-districts

Interaction of Gravity Index and Factors Related to Population
To examine the interaction gravity index and related factors (population density, migration based on local Moran's I, and Getis-Ord G i *), this study analysed the area data of the sub-districts with a statistically significant gravity index only. These sub-districts were located in eight provinces, namely Buriram, Khon Kaen, Nakhon Ratchasima, Nong Khai, Roi Et, Surin, Ubon Ratchithani, and Udon Thani.
The results of the analysis of the gravity index and population density interaction suggested that the correlation (r) was positive, between 0.23 and 0.99 (Figure 6), which indicated a gravity index and population density at different levels, except for provinces such as Ubon Ratchathani. An r = 0.23 indicates that the area with the highest gravity was slightly associated with population density. namely Buriram, Khon Kaen, Nakhon Ratchasima, Nong Khai, Roi Et, Surin, Ubon Ratchithani, and Udon Thani. The results of the analysis of the gravity index and population density interaction suggested that the correlation (r) was positive, between 0.23 and 0.99 (Figure 6), which indicated a gravity index and population density at different levels, except for provinces such as Ubon Ratchathani. An r = 0.23 indicates that the area with the highest gravity was slightly associated with population density.
The gravity index was positively associated with the migration of the population with a value between 0.61 and 0.98. Upon considering all provinces, the r value was higher between the gravity index and other factors (population density, local Moran's I and Getis-Ord Gi*) The gravity index was positively associated with local Moran's I and Getis-Ord Gi*, with a value between 0.54 and 0.98 (local Moran's I) and 0.46 and 0.93 (Getis-Ord Gi*). Local Moran's I had a higher value due to the analysis of the standard division value of all the sub-districts in this region, combined with the nearby areas, while Getis-Ord Gi* measured the nearby areas only. Therefore, the gravity index was highly correlated with migration, population density, local Moran's I, and Getis-Ord Gi*.

Discussion
An analysis of spatial data distribution using the spatial autocorrelation method was practically applied to study medium and small-scale areas [13,23,42,43]. Meanwhile, an analysis of large scale areas in different time periods was infrequently performed [44,45]. This study analysed regional distributions at the sub-district level in the years 2002, 2007, 2012, and 2017 using spatial autocorrelation analysis. There are two levels of the spatial autocorrelation analysis: 1) global indicator and 2) local indicator. The global indicator tool is Moran's I, while the local indicator tools are local Moran's I and Getis-Ord Gi*. In addition, the analysis of the gravity model of the service centre and the interaction between the gravity index and related factors were performed.
The analysis of the spatial population distribution in the Northeast using global Moran's I yielded a positive index in every year of study. This indicated the presence of a clustering of high population density sub-districts. Furthermore, the index of analysis increased throughout the 15 years, indicating positive clustering and the presence of relatively stronger positive spatial autocorrelations among provinces over time. The results of development policies in the Northeast over many years, along with the status of the region as an important connecting area for neighbouring countries, has caused the invigoration of the economy for some provinces along with labour related migrations for education also contributed to this trend. The gravity index was positively associated with the migration of the population with a value between 0.61 and 0.98. Upon considering all provinces, the r value was higher between the gravity index and other factors (population density, local Moran's I and Getis-Ord G i *).
The gravity index was positively associated with local Moran's I and Getis-Ord G i *, with a value between 0.54 and 0.98 (local Moran's I) and 0.46 and 0.93 (Getis-Ord G i *). Local Moran's I had a higher value due to the analysis of the standard division value of all the sub-districts in this region, combined with the nearby areas, while Getis-Ord G i * measured the nearby areas only. Therefore, the gravity index was highly correlated with migration, population density, local Moran's I, and Getis-Ord G i *.

Discussion
An analysis of spatial data distribution using the spatial autocorrelation method was practically applied to study medium and small-scale areas [13,23,42,43]. Meanwhile, an analysis of large scale areas in different time periods was infrequently performed [44,45]. This study analysed regional distributions at the sub-district level in the years 2002, 2007, 2012, and 2017 using spatial autocorrelation analysis. There are two levels of the spatial autocorrelation analysis: (1) global indicator and (2) local indicator. The global indicator tool is Moran's I, while the local indicator tools are local Moran's I and Getis-Ord G i *. In addition, the analysis of the gravity model of the service centre and the interaction between the gravity index and related factors were performed.
The analysis of the spatial population distribution in the Northeast using global Moran's I yielded a positive index in every year of study. This indicated the presence of a clustering of high population density sub-districts. Furthermore, the index of analysis increased throughout the 15 years, indicating positive clustering and the presence of relatively stronger positive spatial autocorrelations among provinces over time. The results of development policies in the Northeast over many years, along with the status of the region as an important connecting area for neighbouring countries, has caused the invigoration of the economy for some provinces along with labour related migrations for education also contributed to this trend.
The results of the local Moran's I and Getis-Ord G i * analysis were consistent. A dense population was identified in clusters in urban areas near the central areas of eight provinces. This indicated that the neighbouring sub-district population density changed in the same direction. For example, if one sub-district has a high population density, the neighbouring sub-district will also have the same densely populated level. In addition, if one district has a low population density, the neighbouring sub-districts will also have a low population density. Throughout the period ranging from 2002 to 2017, hot spots increased, while cold spots decreased.
One of the reasons for the observed patters is that the population in the Northeast increased rapidly from 21,551,434 people in 2002 to 21,989,477 people in 2017. Public utilities and, agricultural processing plants increased according to the Department of Industrial Works (2014) [46]. There were 41,220 registered factories in 2005; meanwhile, 44,143 factories existed in 2014. The transportation or movement between provinces has become more convenient due to the extension of the main roads, including the fact that in 2015 the airport in Roi Et province was upgraded to accommodate long haul flights. In addition, regional airports were renovated and upgraded to international levels in Khon Kaen, Ubon Ratchathani, and Udon Thani provinces [47]. The economic growth resulted in an average per capita income increase for the population during 2005-2014 from 32,308 to 71,286 Baht per capita [48]. Similarly, Luo [27] found that economic fragmentation was an important factor affecting the patterns of the population distribution in the area of Sichuan-Yunnan-Guizhou province, China.
Most of the areas with high-low and low-high clustering were located around the city areas due to the gap between economic development cities and that in urban areas. Therefore, the rates of economic development and urban population density have also increased. Some reasons for this are migration, whereby people move from the suburbs into the city areas. Similarly, in Nanijung, Jiangsu province, China, a high-low clustering has been reported due to the large population in Nanijung downtown. However, the neighbouring cities has small populations [28].
Eight provinces with H-H clustering were different from the 12 provinces in the Northeast without clustering, especially in the border cities, such as Surin, Nong Khai, Ubon Ratchathani. Meanwhile, Buriram province is considered to be the distinctive point of tourism, sporting destinations for football, and world-class racing tracks. Nakhon Ratchasima, Khon Kaen, Udon Thani provinces are defined as the centres of trade and education among neighbouring countries. A cluster has also been found in the Mueang district area, Roi Et province. Although the province is not considered to be a distinct city like the others, the clustering results revealed a dense population in the centre of the city and areas around the city. The reason for this pattern is likely due to the province being the centre of the Northeast, where it represents a hub for the connection of services and product distribution to nearby provinces.
The results of the gravity analysis in the service centres of cities with high indexes such as Nakhon Ratchasima, Udon Thani, and Khon Kaen, were consistent with the reality of the city areas. As these are a large cities in the region and economic centres of trade, health services, and education, which are utilities for supporting population growth, the cities have expanded rapidly from the past to the present [49]. Similar results from Bunea [50] indicated that migration was related to the attraction of the area, especially the population size and distance of migration. If the population was large, there was also a chance of migration. On the other hand, migration opportunities were reduced in instances of far proximity. Meanwhile, Mafi [51] suggested an analysis of accesses to service areas within the shortest time based on speed and types of roads. The results of the analysis were consistent with the shortest time used to access the service centre, speed limited by road types, and the quality of the road affecting the gravity of the local population travels and transportations.
An analysis of spatial data distribution using spatial autocorrelation and the comparison of gravity indexes based on population-related factors showed that the coherent provinces of the sub-district had a high population density and a large population. The main reasons for the increase in clustering were economic and service sources, for instance, hospitals, multifarious, famous educational institutions, being promoted as a world-class sports city, a border town of international trades, et cetera, which can attract populations of migrants in search of job opportunities [50]. Therefore, the development of public utilities and transportations should play a key role in attracting people, which will also affect industrial properties, trade, and services.
For further suggestions in the study of the phenomenon of populations incoherent with accurate patterns to identify economic conditions, in areas with high growth and urbanisation, the sizes and concentrations will be higher than the urbanised communities. Although this study is limited to the collection of the gross domestic product (GDP) at the district levels of Thailand, population density data per area can show the patterns of population distribution, and correlated to gravity and population-related factors, which can effectively solve the problems of local data storages in large areas, regions, and countries. In addition, the approach can further be developed to quickly identify the concentration and distribution patterns of different economic zones. The results of the study can be utilised to help develop an area in line with the population density of each province in the Northeast with efficiency and accuracy in a wide area to a local level.

Conclusions
Funding: This research was supported by the Faculty of Humanities and Social Sciences, Mahasarakham University, Fiscal Year 2019 under unit number 033/2019. The funding agency made no contribution to the study design, data collection, analysis and interpretation of data, writing of the report, and the decision to submit the article for publication.