Spatiotemporal Analysis of Environmental Factors on the Birdstrike Risk in High Plateau Airport with Multi-Scale Research

: The aircrafts’ engine performance deteriorates sharply during the take-o ﬀ and landing at high plateau airport. This situation increases the take-o ﬀ or landing distance, aggravating the hidden danger of birdstrikes at high plateau airport. This paper ﬁrst used GIS to classify and rasterize the bird data and calculated the monthly Birdstrike Risk Index (BRI) within 6, 13, and 25 km radii of Lhasa Airport, based on the bird observation data of Tibet and the birdstrike data of Lhasa Airport from 2015 to 2019. The spatiotemporal relationships between the BRI and the environmental factors around Lhasa Airport were compared by the Geographically or Temporally Weighted Regression (GWR or TWR) model and Geographically and Temporally Weighted Regression (GTWR) model. The results showed that the temporal nonstationary e ﬀ ect of environmental factors was more signiﬁcant than that of spatial nonstationary at Lhasa Airport. Besides, the composition of land types had positive impacts on birdstrike risk within the 6 km radius, and this scope was broader than that of the plain airport. Within the 13 km and 25 km ranges, the water distribution and the altitude during dry season also positively impacted birdstrike risk. Moreover, the key factor to birdstrike risk was the water distribution in December.


Introduction
High plateau airport is an airport with an altitude of more than 2438 m (8000 feet). Due to the complex and harsh terrain and meteorological conditions, the aircraft's engine performance is reduced, such as difficulty to start, thrust reduction, and other problems. Therefore, the runway length, thickness, and width of the high plateau airport are higher than those of the plain airport [1,2]. As a result, aircraft occupy runways and airways for a long time during take-off and landing at high plateau airports, resulting in a higher birdstrike risk than plain airports. Thus, the hidden danger of birdstrike in the high plateau airport is one of the important factors affecting the airport's safe operation.
At present, the assessment methods of birdstrike risk at the airport are mainly based on the risk matrix with the combination of occurrence probability and historical severity. For example, the airport specific data was used to evaluate birdstrike risk by creating a simple probability-time-severity matrix [3]. Furthermore, Ning et al. [4] used the real-time bird radar data of Beijing Capital International Airport in 2014 to evaluate the birdstrike hazard degree by combining the Delphi Method and Analytic Hierarchy Process. The probability estimation model was used to determine the birdstrike probability based on the relative position of the flock and the runway, the flock and the height of the aircraft, the aircraft route, and the bird avoidance ability to evaluate birdstrike risk in the airfield area. However, goodness-of-fit than the global Ordinary Least Squares and GWR model and can visually analyze the results in time and space [16]. The analysis of the factors that affect the birdstrike risk is also a typical spatiotemporal geography issue. Therefore, this paper will use the GTWR model to deeply explore the spatiotemporal influence mechanism of environmental factors on birdstrike risk in high altitude airports.
In this paper, Lhasa Kongga International Airport is used as the research airport. We use the bird situation analysis report of Lhasa Airport from 2015 to 2019 and the bird observation data in Tibet to evaluate the birdstrike risk due to the lack of bird detection radar and the inability to achieve artificial bird observation. The bird observation data includes observation location, time, bird species, and number. First of all, bird observation data are statistically classified and rasterized by GIS technology. The BRI existing research is improved to generate the mBRI suitable for Lhasa Airport. Secondly, taking the calculated mBRI as the dependent variable and the selected five environmental factors as independent variables, the traditional GWR model, TWR model, and GTWR model are used to compare and analyze the relationship. Finally, the spatiotemporal effect of environmental factors affecting birdstrike risk is visualized and analyzed based on the GTWR model results. The definition of acronyms frequently in this paper are shown in Appendix A (Table A1).

Study Area and Data Collection
Lhasa Kongga International Airport is located in Gongga County, Shannan City, southwest of Lhasa City, Tibet Autonomous Region, China. It is about 85 km away from the center of Lhasa City. The airport grade is 4E, the runway is 4000 m long and 45 m wide, and operates on a single runway. With an altitude of 3600 m, the airport is one of the world's highest civil airports and the largest aviation hub in the Tibet Autonomous Region. Lhasa Airport is located on the south side of the Brahmaputra River valley, with complex topography, high mountains on both sides of the north and south, and poor clearance conditions. Within the range of the go-around point (about 13 km from the end of the runway), there are continuous mountains of more than 4000 m, affecting take-off and landing significantly. The annual average temperature of Lhasa Airport is 9 • C, the annual average relative humidity is 46.2%, the easterly wind prevails throughout the year, the yearly average wind speed is 2.5 m/s, and there are no distinct four seasons.
According to the Lhasa Airport's flight procedures, which are different from those for an ordinary plain airport, the study area in this paper is a circular area within a 25 km radius from the airport center, and the whole study area is divided into 1*1 km grids. Moreover, the birdstrike risk of Lhasa Airport is estimated in three different circular buffers (6 km, 13 km, and 25 km). On this basis, the temporal and spatial distribution of the influencing factors is studied, and the schematic diagram of the study area is shown in Figure 1. Figure 1 was generated using ArcGIS software to extract the remote sensing monitoring data of land use status in China.
Since Lhasa Airport is located between river valleys and surrounded by mountains, the geographical environment is complex. So it is difficult to observe birds. We used the bird situation analysis report of the Lhasa Airport airfield area from 2015 to 2019 and the bird observations of the entire Tibet area (eBird Observation Dataset [17]). Moreover, bird observation data covers the observation location, time, bird species, and number in Tibet from 2015 to 2019 (Table S1). It is rasterized using GIS technology to facilitate the following calculation. To simplify the data management, based on the bird population composition proposed by Soldatini et al. [6], we reclassified the bird data combined with the high altitude birds' unique ecological features (habitat, foraging, weight, and clustering or not). The specific grouping is shown in Table 1. At the same time, we compared the bird situation analysis report of Lhasa Airport from 2015 to 2019. In the past five years, the most expelled birds were Turtle Dove, Ruddy Shelduck, Passer montanus, Bar-headed Goose, and Small Skylark, accounting for 21.59%, 16.99%, 14.16%, 12.91%, and 8.81% of all expelled birds, respectively. The number of Sustainability 2020, 12, 9357 4 of 18 manual monitoring and evictions peaked in August at the airfield area, showing a fluctuating trend as a whole. Considering the breeding and migration features of birds, according to the high plateau's unique meteorological characteristics and our investigation, the research time is divided into two stages: the dry season (November-April) and the rainy season (May-October).
Sustainability 2020, 12, x FOR PEER REVIEW  4 of 18 birds, according to the high plateau's unique meteorological characteristics and our investigation, the research time is divided into two stages: the dry season (November-April) and the rainy season (May-October).  Notes: The bold indicates that the monthly average bird species risk of these groups greater than 5 (See Figure 2).  Notes: The bold indicates that the monthly average bird species risk of these groups greater than 5 (See Figure 2).

Dependent Variable-Monthly Birdstrike Risk (mBRI) at Lhasa Airport
Based on the Birdstrike Risk Index (BRI) and the Attraction Risk Index of the airport by Cecilia [6] and Francesca et al. [7], the mBRI is proposed to evaluate the birdstrike risk within a 25 km radius from Lhasa Airport. The mBRI consists of the following components: (1) bird species factors, including species weight and population aggregation; (2) aircraft factors, including the effect of takeoffs and landings of aircrafts' number on the birdstrike (Table S2); (3) historical birdstrike factors, including the impact degree of historical birdstrike on flights and the historical number of birdstrikes caused by different bird species (Table S3). For details, please refer to the classification criteria of Effect On Flight (EOF) proposed by Cecilia et al. (4) Height factors, including the flying height of birds and aircraft, may interfere with aircraft only when a group of species can reach the height of the aircraft passing over the study area. Therefore, we also introduce a height indicative variable h according to the bird flight altitude [18,19] and flight program of high plateau airport, which determines whether or not the species can reach the height of the aircraft to interfere with the aircraft. Multiplying h with its corresponding other evaluation indicators can exclude those bird species that are unlikely to cause interference to aircraft [7]. Based on the above factors, determine the bird species risk GF i of the ith group, as shown in Equation (1).
where W i is the average weight of the ith bird species group, g. Ag i is the aggregation factor of the ith bird species group. BS i is the monthly average number of birdstrike of the ith bird species group.
TLS is the monthly average takeoffs and landings of the airport. EOF (95) i is the 95th percentile that the effect of a birdstrike on the flight in ith bird species group, and the specific bird information is shown in Table 1. Thus, the actual monthly bird species risk GSR i (n) of group i in the nth grid can be calculated, as shown in Equation (2).
where GF i (n) is the bird species risk of group i in the nth grid, whether there is bird species group i in the nth grid is determined by nuclear density analysis and rasterization. If the density value of bird species group i in the nth grid is greater than 1, it is deemed to be present, otherwise considered to be absent. MB i (n) is the number of monthly occurrences of group i in the nth grid, this parameter of bird species group i in the nth grid is obtained by Kriging and rasterization of the observed numbers over the past five years in this paper. Finally, combined with birds' flying altitude, the altitude indicative variable h is introduced to obtain the mBRI in the nth grid, as shown in Equation (3).
where h i (n) determines whether or not the ith group of bird species could reach the height of flight of aircraft over the nth grid area, it is 1 if it can be achieved, otherwise is 0. Figure 2 shows the average activity height of different bird species and aircraft flight altitude of the aircraft in the three circular buffers around the airport. The yellow-marked group is the monthly average bird species risk greater than 5 (See Table S4). That is, these bird species have a higher birdstrike risk.
Sustainability 2020, 12, x FOR PEER REVIEW 6 of 18 Figure 2. Schematic diagram of bird species activity height and aircraft height. The yellow-marked group indicates these bird species have a higher birdstrike risk than other groups. These groups were bold in Table 1.

Independent Variable-Environmental Factors
(1) Vegetation coverage: Since vegetation provides food and habitat for herbivorous birds, the further analysis of vegetation coverage mainly considers the impact of vegetation coverage changes on the birdstrike risk. The Normalized Difference Vegetation Index (NDVI) is recognized as an effective parameter to characterize vegetation changes, and it contains useful information about vegetation coverage [20]. (2) Probability of birds appearing in different land types (PB): This paper makes a statistical analysis of bird observation data from 2015 to 2019 in Tibet. We found that different birds will forage in distinct habitats at different seasons, indicating that different land types Figure 2. Schematic diagram of bird species activity height and aircraft height. The yellow-marked group indicates these bird species have a higher birdstrike risk than other groups. These groups were bold in Table 1.

Independent Variable-Environmental Factors
(1) Vegetation coverage: Since vegetation provides food and habitat for herbivorous birds, the further analysis of vegetation coverage mainly considers the impact of vegetation coverage changes on the birdstrike risk. The Normalized Difference Vegetation Index (NDVI) is recognized as an effective parameter to characterize vegetation changes, and it contains useful information about vegetation coverage [20]. (2) Probability of birds appearing in different land types (PB): This paper makes a statistical analysis of bird observation data from 2015 to 2019 in Tibet. We found that different birds will forage in distinct habitats at different seasons, indicating that different land types have different attraction levels to birds. We have counted the proportions of 11 groups of birds in different land types in different months over the past five years, as the attraction of different land types to birds. Therefore, the sum of the attractiveness in different land types to 11 groups of birds is an explanatory variable to analyze the effect of the land composition around Lhasa Airport on the birdstrike risk. (3) Human population density (POP): This paper uses the population per square kilometer during the studied period as an explanatory variable to analyze the impact of human activities around the airport on the birdstrike risk. (4) Water distribution (WD): According to the well-known Tobler's first law of geography [21], the assumption is that nearby water has a more significant effect on birds than water farther apart. Therefore, the reciprocal of the nearest Euclidean distance between different grids and the water body is used to describe the water distribution variable, and the influence of the water's spatial distance on the birdstrike risk in different areas is considered. (5) Altitude (DEM): The Digital Elevation Model (DEM) information is directly used as the studied area's altitude. The descriptive statistics of the independent variable data in 2018 used in the paper are shown in Table 2.

Geographical and Temporal Weighted Regression (GTWR) Model
Firstly, the definition of the Geographical Weighted Regression (GWR) model will be more helpful for understanding the Geographical and Temporal Weighted Regression (GTWR) model. The GWR model is an extension based on the traditional global regression model, which embeds the data's spatial characteristics into the model. It is a common spatial statistical method used to analyze the local changes of spatial phenomena. The fundamental equation of the GWR model [22] is described as Equation (4).
is the kth regression parameter of the BRI of the ith grid, X ik is the kth environment variable, ε i is the random error term. The GWR model's salient feature allows the estimated values of parameters to change with different spaces to measure its spatial nonstationarity. However, the GWR model cannot solve the temporal nonstationarity, which is also an essential factor affecting the birdstrike risk [13]. The treatment of spatial variation in the GTWR model is similar to that in the GWR model. The GTWR model uses space-time distance to construct a space-time weight Sustainability 2020, 12, 9357 7 of 18 matrix to measure the relationship between birdstrike risk and different regions and time. The GTWR model is described as Equation (5).
where t i is the temporal attribute of the ith grid. Based on the Weighted Least Square regression, the unknown parameters of grid i are estimated, and the objective function is shown in Equation (6). (7).
It is still assumed that the data points "close" to the point i in the space-time coordinate system have a greater impact on the parameter estimation of the point i than the data points slightly farther away. In this sense, "closeness" here includes the proximity of space and time. Due to the different measurement units of space distance and time distance (in our research, location in meters and time in months), it is more practical to use the elliptical coordinate system to measure the "closeness" between the regression point i and its adjacent point j.
hence, the spatiotemporal distance d ST between points i and j as Equation (8) [14].
where µ is a scale factor used to balance the different influences caused by the different measurement units of the distance of space and time. In fact, the essential effect of µ is to enlarge or reduce the temporal distance effect to match with spatial distance. Therefore, if appropriate adjustments are made to the parameter µ, the spatiotemporal distance d ST can be used to measure the degree of "closeness" in spacetime. As shown in Equation (9), the Gaussian kernel function based on distance attenuation is used to construct the spatiotemporal weight matrix W ij .
where h ST is the spatiotemporal bandwidth [15]. h S is the spatial bandwidth, h T is the temporal bandwidth. Figure 3 shows the primary process of parameter estimation by the GTWR model [22]. First, the traditional GWR model is specified based on the input data of Equation (4). The goodness index is then fitted based on the cross-validation to define the optimal spatial bandwidth for each period. Secondly, the optimal time bandwidth is determined based on cross-validation again. Once the optimal spatial and temporal bandwidths are obtained, they can construct the spatiotemporal weight matrix W i (x i , y i , t i ). Then the Equation (7) can be used to estimate the local parameters β(x i , y i , t i ). Figure 3 shows the primary process of parameter estimation by the GTWR model [22]. First, the traditional GWR model is specified based on the input data of Equation (4). The goodness index is then fitted based on the cross-validation to define the optimal spatial bandwidth for each period. Secondly, the optimal time bandwidth is determined based on cross-validation again. Once the optimal spatial and temporal bandwidths are obtained, they can construct the spatiotemporal weight matrix ( , , ) i i i i W x y t . Then the Equation (7) can be used to estimate the local parameters ( ) , , the sum of the squares of these errors can be expressed as Equation (10). Equation (10) is minimized to determine the optimal bandwidth. The selection of bandwidth h S is also an important step in revising the GTWR model. We use the cross-validation (CV) method to determine the bandwidth h S . When the bandwidth is h S , assuming thatŶ −i (h S ) is the fitted value for Y i with point i excluded from the GTWR model process, the sum of the squares of these errors can be expressed as Equation (10). Equation (10) is minimized to determine the optimal bandwidth.

Spatial and Temporal Feature of BRI
The monthly average Birdstrike Risk Index (BRI) of Lhasa Airport during the dry and rainy season is shown in Figure 4. The areas of high BRI value in the dry season are larger than that of the rainy season. Moreover, the dry season's high-value regions clustered in the 13 km radius of the airport and approached the runway. However, in the rainy season, the high BRI areas are only distributed in a small part of the whole studied area and far away from the airport. When the aircraft lands at a horizontal distance of about 13 km from the runway, its flight height is about 730 m. As the aircraft continues to descend, it gradually enters the area of high birdstrike risk, which is most apparent in the dry season. The opposite is true during the take-off. Moreover, the fields with low BRI value are all distributed outside the 13 km radius of the airport, showing that the birdstrike risk is closely related to aircraft height.
Also, the minimum (Min.), the first quartile (1st Qu.), the median (Median), the mean (Mean), the third quartile (3rd Qu), and the maximum (Max.) six statistics are selected to describe the BRI of Airport based on the calculation results in Section 2.2. Furthermore, all the calculation results of BRI are shown in Table S4. Due to the lack of bird observation data in January and February (the reason may be that the coldest months are January and February on the Tibetan Plateau, with almost no bird activity), January and February will be ignored in the subsequent calculation. It can be seen from Figure 5 that, in general, the highest birdstrike risk occurred in December at Lhasa airport, which is different from Sustainability 2020, 12, 9357 9 of 18 the period with the plain airport [7,23]. Next, we will continue to explore why the difference in the pattern of manifestation over birdstrike risk between the high plateau airport and the plain airport.
The monthly average Birdstrike Risk Index (BRI) of Lhasa Airport during the dry and rainy season is shown in Figure 4. The areas of high BRI value in the dry season are larger than that of the rainy season. Moreover, the dry season's high-value regions clustered in the 13 km radius of the airport and approached the runway. However, in the rainy season, the high BRI areas are only distributed in a small part of the whole studied area and far away from the airport. When the aircraft lands at a horizontal distance of about 13 km from the runway, its flight height is about 730 m. As the aircraft continues to descend, it gradually enters the area of high birdstrike risk, which is most apparent in the dry season. The opposite is true during the take-off. Moreover, the fields with low BRI value are all distributed outside the 13 km radius of the airport, showing that the birdstrike risk is closely related to aircraft height. Also, the minimum (Min.), the first quartile (1st Qu.), the median (Median), the mean (Mean), the third quartile (3rd Qu), and the maximum (Max.) six statistics are selected to describe the BRI of Airport based on the calculation results in Section 2.2. Furthermore, all the calculation results of BRI are shown in Table S4. Due to the lack of bird observation data in January and February (the reason may be that the coldest months are January and February on the Tibetan Plateau, with almost no bird activity), January and February will be ignored in the subsequent calculation. It can be seen from Figure 5 that, in general, the highest birdstrike risk occurred in December at Lhasa airport, which is different from the period with the plain airport [7,23]. Next, we will continue to explore why the difference in the pattern of manifestation over birdstrike risk between the high plateau airport and the plain airport.    Also, the minimum (Min.), the first quartile (1st Qu.), the median (Median), the mean (Mean), the third quartile (3rd Qu), and the maximum (Max.) six statistics are selected to describe the BRI of Airport based on the calculation results in Section 2.2. Furthermore, all the calculation results of BRI are shown in Table S4. Due to the lack of bird observation data in January and February (the reason may be that the coldest months are January and February on the Tibetan Plateau, with almost no bird activity), January and February will be ignored in the subsequent calculation. It can be seen from Figure 5 that, in general, the highest birdstrike risk occurred in December at Lhasa airport, which is different from the period with the plain airport [7,23]. Next, we will continue to explore why the difference in the pattern of manifestation over birdstrike risk between the high plateau airport and the plain airport.

Model Comparison
Before establishing the model, the multicollinearity between explanatory variables should be tested. It will lead to model estimation distortion or inconsistent if there is serious multicollinearity in the model. The variance inflation factor of each variable in two periods is calculated respectively in this paper. The results showed that these five variables' variance inflation factors were much less than 10 and floated at about 2 [24]. It indicated that these variables had no multicollinearity, and the spatiotemporal analysis could be continued.
The same data were used for GWR, TWR, and GTWR models analysis. As the effects of different environmental factors on birdstrike risk vary in space and time, the output data of the estimated coefficients of this series of local regression models are voluminous (all calculation results can be seen in Tables S5-S7). Moreover, the evaluation indexes of different models are also given in Figure 6, such as goodness-of-fit, residual sum of squares, and Akaike information criterion at the same time.
It is worth noting that in the dry season, the goodness-of-fit of GTWR was 63.9% and 46.6% better than GWR and TWR models. It also increased by 44.9% and 25.1%, respectively, in the rainy season.
coefficients of this series of local regression models are voluminous (all calculation results can be seen in Tables S5-7). Moreover, the evaluation indexes of different models are also given in Figure 6, such as goodness-of-fit, residual sum of squares, and Akaike information criterion at the same time. It is worth noting that in the dry season, the goodness-of-fit of GTWR was 63.9% and 46.6% better than GWR and TWR models. It also increased by 44.9% and 25.1%, respectively, in the rainy season. Furthermore, by comparing the residual sum of squares and Akaike information criterion's diminution, it is further shown that the GTWR model provided better fitting than the GWR and TWR models. We believe that this result is because the GTWR model can handle both spatial and temporal heterogeneity. On the other hand, it can be seen from Figure 6 that the TWR model showed better goodness-of-fit than the GWR model. It showed that the spatial nonstationary effect was not as significant as the temporal nonstationarity.
Besides, it is necessary to detect further the nonstationarity of each group of parameters in the GTWR model, that is, whether or not the parameters have significant differences in the study area. An effective method is to compare the twice standard error of global Ordinary Least Squares with the Inter Quartile Range of the GTWR model and the latter's larger values indicating significant nonstationary [15,22]. As shown in Table 3, there were no extra local variations in the coefficients of DEM, POP, and WD variables in the TWR model. There were extra local variations in all GWR and GTWR models parameters, indicating that the estimation results of the GTWR model had significant spatiotemporal nonstationarity. Therefore, it is suggested to use the GTWR model to reveal the spatiotemporal variation of environmental factors on birdstrike risk at Lhasa Airport.  Furthermore, by comparing the residual sum of squares and Akaike information criterion's diminution, it is further shown that the GTWR model provided better fitting than the GWR and TWR models. We believe that this result is because the GTWR model can handle both spatial and temporal heterogeneity. On the other hand, it can be seen from Figure 6 that the TWR model showed better goodness-of-fit than the GWR model. It showed that the spatial nonstationary effect was not as significant as the temporal nonstationarity.
Besides, it is necessary to detect further the nonstationarity of each group of parameters in the GTWR model, that is, whether or not the parameters have significant differences in the study area. An effective method is to compare the twice standard error of global Ordinary Least Squares with the Inter Quartile Range of the GTWR model and the latter's larger values indicating significant nonstationary [15,22]. As shown in Table 3, there were no extra local variations in the coefficients of DEM, POP, and WD variables in the TWR model. There were extra local variations in all GWR and GTWR models parameters, indicating that the estimation results of the GTWR model had significant spatiotemporal nonstationarity. Therefore, it is suggested to use the GTWR model to reveal the spatiotemporal variation of environmental factors on birdstrike risk at Lhasa Airport.

Spatial Feature of Coefficients
One of the essential features of this series of weighted regression models is that the local estimated parameters of spatial relationships are mappable and can be visualized. We grouped the regression coefficients into several intervals. We used different colors to visualize the spatial variations of these environmental factors' effects on the birdstrike risk at Lhasa Airport. In this paper, the mean values of regression coefficients in different months at the same location were used to show the influencing effect's spatial variation pattern. The GTWR model considers both temporal and spatial heterogeneity of each variable's coefficients to vary with time. According to the GTWR model results, the human population density had almost no effect on the spatial distribution of birdstrike risk at Lhasa Airport, so the remaining four factors' spatial variation was only depicted in the paper. Figure 7 shows the spatial effects of vegetation coverage (NDVI). Given the spatial impact of vegetation coverage, the spatial influence patterns are very different in the two research periods. The NDVI had almost no positive effects on birdstrike risk in the entire study area in the dry season. However, in the rainy season, the positive effects' areas were mainly clustered in the southern area between 6 km and 25 km away from the runway, which indicated that vegetation coverage in the rainy season had a more significant effect on the birdstrike risk. In the high plateau, the vegetation coverage in the rainy season is higher than that in the dry season. Correspondingly, Lucas et al. [25] also found that NDVI is the most crucial predictor of bird species richness. It shows that the high vegetation coverage in the rainy season attracts more birds in the high plateau airport and aggravates the birdstrike risk.

Spatial Feature of Coefficients
One of the essential features of this series of weighted regression models is that the local estimated parameters of spatial relationships are mappable and can be visualized. We grouped the regression coefficients into several intervals. We used different colors to visualize the spatial variations of these environmental factors' effects on the birdstrike risk at Lhasa Airport. In this paper, the mean values of regression coefficients in different months at the same location were used to show the influencing effect's spatial variation pattern. The GTWR model considers both temporal and spatial heterogeneity of each variable's coefficients to vary with time. According to the GTWR model results, the human population density had almost no effect on the spatial distribution of birdstrike risk at Lhasa Airport, so the remaining four factors' spatial variation was only depicted in the paper. Figure 7 shows the spatial effects of vegetation coverage (NDVI). Given the spatial impact of vegetation coverage, the spatial influence patterns are very different in the two research periods. The NDVI had almost no positive effects on birdstrike risk in the entire study area in the dry season. However, in the rainy season, the positive effects' areas were mainly clustered in the southern area between 6 km and 25 km away from the runway, which indicated that vegetation coverage in the rainy season had a more significant effect on the birdstrike risk. In the high plateau, the vegetation coverage in the rainy season is higher than that in the dry season. Correspondingly, Lucas et al. [25] also found that NDVI is the most crucial predictor of bird species richness. It shows that the high vegetation coverage in the rainy season attracts more birds in the high plateau airport and aggravates the birdstrike risk.   Figure 8 shows the spatial effects of birds' probability of appearing on different land types (PB). There were regions with significant positive effects, whether in the dry or the rainy season. In the dry season, the positive impacts' regions were located in the southwest, extending to 25 km from the airport. In contrast, in the rainy season, they were mainly distributed in the airport's 13 km radius. However, Morgan et al. [10] found that the increase in the diversity of land types in the 3 km radius of the plain airport will increase the birdstrike rate. Some experts have proposed to vigorously eliminate the bird attractants in the airport's 5 km radius to reduce birdstrike occurrence in the plain airport [13]. However, this region expanded to 13 km in the rainy season at Lhasa Airport, even broadened to 25 km in the dry season. The reason for this phenomenon may consist of the following two parts: at first, the special geographical environment and bird species composition at the high plateau airport are very different from those in the plain airport. Secondly, the aircraft's performance sharply declines when entering and leaving the Lhasa Airport due to the special geographical environment and harsh meteorological conditions of the high plateau. Therefore, when the aircraft lands at Lhasa Airport, it must enter the final approach fix ahead of time than those plain airports, and it also will take longer to complete the whole take-off process. Hence, the aircraft's time to occupy the runway and low-altitude airway is longer, so the corresponding areas with high birdstrike risk are larger than that of the plain airport.
sharply declines when entering and leaving the Lhasa Airport due to the special geographical environment and harsh meteorological conditions of the high plateau. Therefore, when the aircraft lands at Lhasa Airport, it must enter the final approach fix ahead of time than those plain airports, and it also will take longer to complete the whole take-off process. Hence, the aircraft's time to occupy the runway and low-altitude airway is longer, so the corresponding areas with high birdstrike risk are larger than that of the plain airport. Besides, in the dry season of Lhasa Airport, the grassland area accounts for 82.6% of the total regions where PB had positive effects on birdstrike risk, followed by 6.9% of the forest land. In the rainy season, the largest areas still were grassland, accounting for 57.1% of the total areas where PB had positive effects on birdstrike risk, but there was a sharp decline compared with the dry season. It was cultivated land, which accounted for 20.2%, and the water area also accounted for 16.6%. This change shows that in Lhasa Airport, the land type that needs to reduce birds' attractiveness is grasslands. Other studies have shown that farmland and wetlands around plain airports are the primary reasons for the increase in the number of birds and birdstrike [7,10]. Concisely, due to the high plateau airport's particular geographical environment, the birdstrike defense and measures to reduce birds' attraction are different from those in the plain airport. The scope of prevention is also larger than that in the plain airport. Figure 9 shows the spatial effects of water distribution (WD) affecting the birdstrike risk in different seasons. There was almost no effect on the birdstrike risk within the 6 km radius of Lhasa Airport for the water distribution. However, the areas that the WD had positive impacts were distributed on the west outside the 6 km scope, whether in the dry or rainy season. In the areas within the 13 km to 25 km radius, the positive influences' regions gradually extended to the study area's edge. The most positive significant areas were clustered in the southwest of the airport in the dry season; inversely, they were distributed in the northwest of the airport in the rainy season. This phenomenon may be related to the activity area of different bird species around Lhasa Airport. For Besides, in the dry season of Lhasa Airport, the grassland area accounts for 82.6% of the total regions where PB had positive effects on birdstrike risk, followed by 6.9% of the forest land. In the rainy season, the largest areas still were grassland, accounting for 57.1% of the total areas where PB had positive effects on birdstrike risk, but there was a sharp decline compared with the dry season. It was cultivated land, which accounted for 20.2%, and the water area also accounted for 16.6%. This change shows that in Lhasa Airport, the land type that needs to reduce birds' attractiveness is grasslands. Other studies have shown that farmland and wetlands around plain airports are the primary reasons for the increase in the number of birds and birdstrike [7,10]. Concisely, due to the high plateau airport's particular geographical environment, the birdstrike defense and measures to reduce birds' attraction are different from those in the plain airport. The scope of prevention is also larger than that in the plain airport. Figure 9 shows the spatial effects of water distribution (WD) affecting the birdstrike risk in different seasons. There was almost no effect on the birdstrike risk within the 6 km radius of Lhasa Airport for the water distribution. However, the areas that the WD had positive impacts were distributed on the west outside the 6 km scope, whether in the dry or rainy season. In the areas within the 13 km to 25 km radius, the positive influences' regions gradually extended to the study area's edge. The most positive significant areas were clustered in the southwest of the airport in the dry season; inversely, they were distributed in the northwest of the airport in the rainy season. This phenomenon may be related to the activity area of different bird species around Lhasa Airport. For example, among the high-risk bird species, the number of Laridae and Anatidae in the rainy season was significantly higher than that in the dry season at Lhasa Airport. Their active areas are in the vicinity of waters and wetlands. The Lhasa River, a tributary located in the northwest of the airport, is surrounded by wetlands and is more attractive to these waterfowls than the Brahmaputra River's mainstream. Also, the most significant effect of WD in the dry season was higher than that in the rainy season. The reason probably may be that the YamdrokTso (one of the three sacred lakes in Tibet, which is the largest habitat for waterfowl in southern Tibet), is located about 25 km southwest of Lhasa Airport, is more attractive to waterfowls during the dry season.
In the previous study, some experts have proposed that the incidence of adverse birdstrike was closely affected by large areas of wetlands, waters, and crops within the 8 km and 13 km radius of the airport [10]. Given the WD, the areas with positive impacts on the birdstrike risk are still larger and more significant than that in the plain airport. Moreover, there are also studies that the possibility of birds using a pond other than its 8 km is zero [26]. On the other hand, this shows that no matter if they are in the plateau or the plain, the waters near the airport do have a close relationship with bird activities, which increases the probability of birdstrike. Therefore, the scope of water influence near Lhasa airport was much broader than that of the plain airport. vicinity of waters and wetlands. The Lhasa River, a tributary located in the northwest of the airport, is surrounded by wetlands and is more attractive to these waterfowls than the Brahmaputra River's mainstream. Also, the most significant effect of WD in the dry season was higher than that in the rainy season. The reason probably may be that the YamdrokTso (one of the three sacred lakes in Tibet, which is the largest habitat for waterfowl in southern Tibet), is located about 25 km southwest of Lhasa Airport, is more attractive to waterfowls during the dry season. In the previous study, some experts have proposed that the incidence of adverse birdstrike was closely affected by large areas of wetlands, waters, and crops within the 8 km and 13 km radius of the airport [10]. Given the WD, the areas with positive impacts on the birdstrike risk are still larger and more significant than that in the plain airport. Moreover, there are also studies that the possibility of birds using a pond other than its 8 km is zero [26]. On the other hand, this shows that no matter if they are in the plateau or the plain, the waters near the airport do have a close relationship with bird activities, which increases the probability of birdstrike. Therefore, the scope of water influence near Lhasa airport was much broader than that of the plain airport. Figure 10 shows the spatial effects of altitude (DEM) affecting the birdstrike risk in different seasons. For the DEM, in the dry season, the most significant positive impact areas were distributed in the north and northwest of the airport, almost entirely within the airport's 13 km radius. In the range from 13 km to 25 km radius of the airport, the positive influence degree decreased, but the area increased. However, in the rainy season, the DEM had almost no positive effect on birdstrike risk in the whole study area. The spatial variation between different seasons is also related that the activity regions and heights of birds are disparate in different seasons. Most of the birds that are active in areas with significant positive impacts on birdstrike risk in the dry season are small raptors (such as Falco tinnunculus) and flocking passerine (such as Montifringilla ruficollis and Eremophila alpestris). Nevertheless, during the rainy season, the numbers of passerine have declined. However, birds of the Falconiformes have a wide range of activities and are not restricted by terrain heterogeneity. The range of Passeriformes is relatively narrow, so the terrain plays a vital role in distributing passerine birds [12]. Another part of the reason may be that due to uneven solar illumination, there is a significant difference in the distribution of woodlands on both sides of the airport, which indirectly affects the distribution of passerine birds.
On the other hand, the National Wildlife Strike Database for civil aviation in the USA have been statistically analyzed. For birdstrike that causes severe damage to the aircraft at the plain airport, showed that 66% occurred at less than or equal to 500 ft (153 m), 29% of them between 501 ft and 3500 ft (153-1067 m), the height here is the height above the ground level [27]. However, there are peaks of more than 4000 m on both sides of Lhasa airport. Most of the aircraft is less than 500 ft away from  Figure 10 shows the spatial effects of altitude (DEM) affecting the birdstrike risk in different seasons. For the DEM, in the dry season, the most significant positive impact areas were distributed in the north and northwest of the airport, almost entirely within the airport's 13 km radius. In the range from 13 km to 25 km radius of the airport, the positive influence degree decreased, but the area increased. However, in the rainy season, the DEM had almost no positive effect on birdstrike risk in the whole study area. The spatial variation between different seasons is also related that the activity regions and heights of birds are disparate in different seasons. Most of the birds that are active in areas with significant positive impacts on birdstrike risk in the dry season are small raptors (such as Falco tinnunculus) and flocking passerine (such as Montifringilla ruficollis and Eremophila alpestris). Nevertheless, during the rainy season, the numbers of passerine have declined. However, birds of the Falconiformes have a wide range of activities and are not restricted by terrain heterogeneity. The range of Passeriformes is relatively narrow, so the terrain plays a vital role in distributing passerine birds [12]. Another part of the reason may be that due to uneven solar illumination, there is a significant difference in the distribution of woodlands on both sides of the airport, which indirectly affects the distribution of passerine birds.

Temporal Feature of Coefficients
By introducing the time dimension into the traditional GWR model to use the GTWR model, we can get the time series in which each factor's coefficients vary with the month from the GTWR model results. In this context, the temporal variation of coefficients about different factors can be calculated and compared.
The estimated value of the birdstrike risk index from the GTWR model and the study area's actual value is compared in Figure 11, which shows that the GTWR model has a better fitting advantage. Moreover, the months with relatively high birdstrike risk are April, August, and On the other hand, the National Wildlife Strike Database for civil aviation in the USA have been statistically analyzed. For birdstrike that causes severe damage to the aircraft at the plain airport, showed that 66% occurred at less than or equal to 500 ft (153 m), 29% of them between 501 ft and 3500 ft (153-1067 m), the height here is the height above the ground level [27]. However, there are peaks of more than 4000 m on both sides of Lhasa airport. Most of the aircraft is less than 500 ft away from these peaks when entering and leaving the airport, which increases the risk of some alpine birds (such as small raptors and Leucosticte brandti) colliding with aircraft, especially in the dry season.

Temporal Feature of Coefficients
By introducing the time dimension into the traditional GWR model to use the GTWR model, we can get the time series in which each factor's coefficients vary with the month from the GTWR model results. In this context, the temporal variation of coefficients about different factors can be calculated and compared.
The estimated value of the birdstrike risk index from the GTWR model and the study area's actual value is compared in Figure 11, which shows that the GTWR model has a better fitting advantage. Moreover, the months with relatively high birdstrike risk are April, August, and December at Lhasa Airport, which is also different from birdstrike at the plain airport. However, at ordinary plain airports, birds' risk of colliding with aircraft is the lowest in winter [23]. More birdstrikes are caused by increased bird activity result from bird migration in spring and autumn [28,29]. Therefore, it is necessary to analyze the birdstrike risk at high plateau airport separately to enhance birdstrike defense capability at high plateau airport.

Temporal Feature of Coefficients
By introducing the time dimension into the traditional GWR model to use the GTWR model, we can get the time series in which each factor's coefficients vary with the month from the GTWR model results. In this context, the temporal variation of coefficients about different factors can be calculated and compared.
The estimated value of the birdstrike risk index from the GTWR model and the study area's actual value is compared in Figure 11, which shows that the GTWR model has a better fitting advantage. Moreover, the months with relatively high birdstrike risk are April, August, and December at Lhasa Airport, which is also different from birdstrike at the plain airport. However, at ordinary plain airports, birds' risk of colliding with aircraft is the lowest in winter [23]. More birdstrikes are caused by increased bird activity result from bird migration in spring and autumn [28,29]. Therefore, it is necessary to analyze the birdstrike risk at high plateau airport separately to enhance birdstrike defense capability at high plateau airport. It can be seen from Figure 12 that the human population density had almost no effect on the birdstrike risk at Lhasa Airport as a whole but only showed a positive impact in May, and the coefficient value was only 0.1. The reason may be that the Lhasa airport is sparsely populated and far away from the urban area. Compared with other environmental factors, the POP has little influence It can be seen from Figure 12 that the human population density had almost no effect on the birdstrike risk at Lhasa Airport as a whole but only showed a positive impact in May, and the coefficient value was only 0.1. The reason may be that the Lhasa airport is sparsely populated and far away from the urban area. Compared with other environmental factors, the POP has little influence on the birdstrike risk at Lhasa airport in the time dimension. Furthermore, the other four factors' coefficients fluctuate greatly over time. Most of the monthly average coefficients swing in a negative direction. Vegetation coverage, the probability of birds appearing in different land types, and water distribution factors all showed positive effects in September and October. The reason may be that with the rise in temperature, the number of waterbirds and herbivorous birds has increased significantly. The risk of birdstrike has also increased compared with other months. on the birdstrike risk at Lhasa airport in the time dimension. Furthermore, the other four factors' coefficients fluctuate greatly over time. Most of the monthly average coefficients swing in a negative direction. Vegetation coverage, the probability of birds appearing in different land types, and water distribution factors all showed positive effects in September and October. The reason may be that with the rise in temperature, the number of waterbirds and herbivorous birds has increased significantly. The risk of birdstrike has also increased compared with other months. Furthermore, the water distribution shows the most significant positive effect on birdstrike risk in December, which is related to the high plateau's special environment. The water area in the dry season has a high aggregation of birds. Besides, altitude showed the most significant positive effect on birdstrike in August. However, the degree of influence was not significant compared with other factors. This change may be related to the number and distribution of small raptors and flocking passerine and perhaps also have something to do with the inconsistent solar illumination on both sides of Lhasa Airport.
In short, the temporal effects of different environmental factors have a specific influence on the birdstrike risk at Lhasa Airport, such as the factors' greatest positive impact at Lhasa Airport is that the human population density (POP) in May, the vegetation coverage (NDVI) in September, the probability of birds looking in different land types (PB) in October, the water distribution (WD) in December, and the altitude (DEM) in August. Still, their significant degree is different, and the biggest is the effect of water distribution. According to the results, the Lhasa airport managers can make targeted arrangements and deployments for birdstrike prevention measures in advance in different months.

Conclusions
In this paper, the GTWR model was first used to analyze the relationship between the birdstrike risk and its environmental factors within the 25 km radius of Lhasa Airport, and further revealed the influence of environmental factors around the airport on the birdstrike risk in space and time. The bird observation data of Tibet and the birdstrike statistical data of Lhasa Airport from 2015 to 2019 were used to identify the spatiotemporal influence of environmental factors on the birdstrike risk. The main conclusions are as follows: Firstly, the GTWR model has a better fitting and analytical ability than the GWR and TWR models, indicating that the effect of different environmental factors on birdstrike risk has significant spatial and temporal characteristics. Secondly, the temporal nonstationary effect was more significant than the spatial. Then, for the feature of spatial influence, the factors of altitude in the dry season and land type in the rainy season had significant positive impacts on birdstrike risk within the 6 km radius Furthermore, the water distribution shows the most significant positive effect on birdstrike risk in December, which is related to the high plateau's special environment. The water area in the dry season has a high aggregation of birds. Besides, altitude showed the most significant positive effect on birdstrike in August. However, the degree of influence was not significant compared with other factors. This change may be related to the number and distribution of small raptors and flocking passerine and perhaps also have something to do with the inconsistent solar illumination on both sides of Lhasa Airport.
In short, the temporal effects of different environmental factors have a specific influence on the birdstrike risk at Lhasa Airport, such as the factors' greatest positive impact at Lhasa Airport is that the human population density (POP) in May, the vegetation coverage (NDVI) in September, the probability of birds looking in different land types (PB) in October, the water distribution (WD) in December, and the altitude (DEM) in August. Still, their significant degree is different, and the biggest is the effect of water distribution. According to the results, the Lhasa airport managers can make targeted arrangements and deployments for birdstrike prevention measures in advance in different months.

Conclusions
In this paper, the GTWR model was first used to analyze the relationship between the birdstrike risk and its environmental factors within the 25 km radius of Lhasa Airport, and further revealed the influence of environmental factors around the airport on the birdstrike risk in space and time. The bird observation data of Tibet and the birdstrike statistical data of Lhasa Airport from 2015 to 2019 were used to identify the spatiotemporal influence of environmental factors on the birdstrike risk. The main conclusions are as follows: Firstly, the GTWR model has a better fitting and analytical ability than the GWR and TWR models, indicating that the effect of different environmental factors on birdstrike risk has significant spatial and temporal characteristics. Secondly, the temporal nonstationary effect was more significant than the spatial. Then, for the feature of spatial influence, the factors of altitude in the dry season and land type in the rainy season had significant positive impacts on birdstrike risk within the 6 km radius around Lhasa Airport. The influencing regions of this land type were larger than that of the plain airport. Moreover, the most significant positive effects were the water distribution within the 6 km to 25 km radius of Lhasa Airport. Lastly, the factors that positively affected birdstrike risk at Lhasa Airport vary significantly in different periods. The key factor was the water distribution in December, followed by the human population density in May.
According to the researched results, the Lhasa Airport planners and managers can apply targeted and suitable bird-repelling measures and preventive measures to create a safe airport environment with a low birdstrike risk. Moreover, other high plateau airports can also appropriately adjust their airfield's precautionary measures to reduce the birdstrike risk based on our study. In the current GTWR model in this paper, we used the linear weighting method to combine spatial and temporal distances, and the advanced nonlinear weighting algorithm needs to be further studied.

Conflicts of Interest:
The authors declare no conflict of interest.