Surface Urban Heat Island and Thermal Profiles Using Digital Image Analysis of Cities in the El Bajío Industrial Corridor, Mexico, in 2020

The Surface Urban Heat Island (SUHI) effect refers to the difference in Land Surface Temperature (LST) between an urban area and its surrounding non-urban area. LST can provide detailed information on the variations in different types of land cover. This study, therefore, analyzes the behavior of LST and SUHIs in fourteen cities in the El Bajío Industrial Corridor, Mexico, using Landsat satellite images from 2020, with QGIS software. It utilizes thermal profiles to identify the land uses that intensify LST, which are essentially those that are anthropologically altered. The results show that the increases in LST and SUHI are more pronounced in cities with greater urban conglomeration, as well as those where there are few green areas and a sizeable industrial or mixed area, with few or no bodies of water. In addition, the increase in temperature in the SUHI is due to certain crops such as vegetables, red fruits, and basic grains such as corn, wheat, and sorghum that use fallow as part of agricultural practices, located around urban areas, which minimizes natural areas with arboreal vegetation.


Introduction
The transformation of natural space through urbanization, which includes residential, commercial, and industrial developments, is largely responsible for its environmental conditions, including urban climatology. Urbanization processes replace land and natural areas with built surfaces, whose materials (such as concrete, asphalt, and tiles) are largely characterized by being almost impermeable. They are composed of rigid, rough elements with sharp edges, low reflectivity, and high albedo, with decreased water absorption capacity and thermal behavior conducive to heat storage and emissions [1][2][3]. These elements, together with humidity and pollutants, contribute to increasing the air temperature of cities in relation to their less urbanized environment, such as the surrounding suburban and rural areas, through a phenomenon known as the Urban Heat Island effect [3][4][5][6][7][8][9][10][11][12][13][14].
The urban heat phenomenon can also be characterized by surface temperatures. While surface temperatures can be higher and more variable than concurrent air temperatures due to the complexity of surface types in urban settings and variations in urban topography [15,16], they tend to be associated more with surface conditions, as described by [17]. The morphology and intensity of heat islands are conditioned, among other factors, by the distribution and composition of large green areas or water mirrors, which appear as relatively low-temperature areas, in comparison with the built-up surfaces of their surroundings [18][19][20]. Since surfaces heat and cool faster than air, the highest surface temperatures are observed at midday [21,22].

Materials
The limit of the cities in the study area was defined by the National Urban System (NUS) version 2018 [27], which identified the municipalities comprising the ten Metropol itan Areas: Aguascalientes, Celaya, Guanajuato, La Piedad-Pénjamo, León, Moroleón-Uri angato, Querétaro, Rioverde, San Francisco del Rincón, and San Luis Potosí, and the fou conurbations: Salamanca, Irapuato, San Juan del Río, and San Miguel de Allende, all o which are located in El Bajío, Mexico, and subsequently delimited the urban areas through the Geostatistical Framework [28]. (The Geostatistical Framework is a unique and nationa system designed by the INEGI, which presents the division of the national territory into different levels of disaggregation to refer geographically to the statistical information o the institutional censuses and surveys and of the state units, each element of these layers has attributes of name and geostatistical key [28]).
Images from the Operational Land Imager (OLI) and Thermal Infrared Sensor (TIRS aboard the Landsat 8 were selected for this study since, according to [29], nearly 53% o researchers have used one or more Landsat images in their SUHI studies, which under lines the importance of the Landsat constellation for SUHI research. The main advantages

Materials
The limit of the cities in the study area was defined by the National Urban System (NUS) version 2018 [27], which identified the municipalities comprising the ten Metropolitan Areas: Aguascalientes, Celaya, Guanajuato, La Piedad-Pénjamo, León, Moroleón-Uriangato, Querétaro, Rioverde, San Francisco del Rincón, and San Luis Potosí, and the four conurbations: Salamanca, Irapuato, San Juan del Río, and San Miguel de Allende, all of which are located in El Bajío, Mexico, and subsequently delimited the urban areas through the Geostatistical Framework [28]. (The Geostatistical Framework is a unique and national system designed by the INEGI, which presents the division of the national territory into different levels of disaggregation to refer geographically to the statistical information of the institutional censuses and surveys and of the state units, each element of these layers has attributes of name and geostatistical key [28]).
Images from the Operational Land Imager (OLI) and Thermal Infrared Sensor (TIRS) aboard the Landsat 8 were selected for this study since, according to [29], nearly 53% of researchers have used one or more Landsat images in their SUHI studies, which underlines the importance of the Landsat constellation for SUHI research. The main advantages of using Landsat images for conducting studies related to obtaining LST and SUHIs, as noted by [29], are that they make it possible to achieve consistent, reliable, voluminous, and freely downloadable data files, in addition to the fact that they have a comparatively high resolution of 60-120 m.
Since the swath coverage of the images is 185 km x 185 km, five images were downloaded to cover the entire area of El Bajío, Mexico. Although the repeated cycle of these images is sixteen days, just one cutoff point was used for 2020. To select the download dates, the study considered the average Maximum Temperatures by State and Nationwide (2020), obtained from the National Water Commission (Comisión Nacional del Agua, Spanish acronym CONAGUA) platform, which shows that the months with Maximum Temperature for cities in El Bajío are April and May. This made it possible to visualize the images of both months and those with the lowest percentage of cloudiness were selected.
The Landsat images were downloaded from the EarthExplorer platform of the United States Geological Survey (USGS) (EarthExplorer data are free and in the public domain; data available from U.S. Geological Survey, National Geospatial Program), while the methodology for obtaining the LST and SUHIs of the Landsat images was developed using QGIS 3.18.3 a free and open-source software.
To visualize the changes in LST and SUHI by type of cover and current land use, the GlobeLand30 Earth land-cover map (the GlobeLand30 2020 is an open access product available at http://www.globallandcover.com/ (accessed on 28 December 2022)) was used as an input, making use of its most recent version created in 2020 in tiles N13_20 and N14_20. This consists mainly of multispectral 30 m images, including TM5, ETM+, and OLI multispectral images from Landsat (USA), as well as China's Environment and Disaster Reduction Satellite HJ-1 imagery and GF-1 16-m resolution multispectral imagery [30], and it is classified into ten land cover classes: crops, trees and forests, grasslands, bushes/shrubs, wetlands/accumulated vegetation, water bodies, tundra, built-up areas, bare soil, and permanent snow and ice. Additionally, the coverage of built-up areas of GlobeLand 30 was recategorized to integrate more categories within the cities, complementing it with current land uses.
Information was extracted on the limits of the blocks comprising each of the cities in the study area from the shapefiles of the massive download page of the National Housing Inventory (   Obtaining inputs: As depicted in Table 1, the Landsat products were downloaded from the EarthExplorer platform. The search criteria used to download the dataset for this study correspond to the weather information provided by the Comisión Nacional del Agua (Comisión Nacional del Agua, Spanish acronym, CONAGUA, 2020) database [31]. Where it is observed that the hottest months corresponding to the Bajío area are April and May, thus selecting the images closest in time and best visualized, i.e., with little or no cloud cover. The acquired data belong to level 1, which are orthorectified products. They consist of quantified data and scaled, calibrated Digital Numbers (DN) representing multispectral image data obtained through both the Operational Land Imager (OLI) and the Thermal Infrared Sensor (TIRS). The products are presented in a 16-bit positive integer format (scaled to 55,000 gray levels) and can be rescaled to Top of Atmosphere (TOA) reflectance and/or radiance using radiometric change of scale coefficients provided in the product metadata file (MTL file). The MTL file also contains the thermal constants required to convert thermal infrared sensor (TIRS) data to the brightness temperature in the satellite [32]. Obtaining inputs: As depicted in Table 1, the Landsat products were downloaded from the EarthExplorer platform. The search criteria used to download the dataset for this study correspond to the weather information provided by the Comisión Nacional del Agua (Comisión Nacional del Agua, Spanish acronym, CONAGUA, 2020) database [31]. Where it is observed that the hottest months corresponding to the Bajío area are April and May, thus selecting the images closest in time and best visualized, i.e., with little or no cloud cover. The acquired data belong to level 1, which are orthorectified products. They consist of quantified data and scaled, calibrated Digital Numbers (DN) representing multispectral image data obtained through both the Operational Land Imager (OLI) and the Thermal Infrared Sensor (TIRS). The products are presented in a 16-bit positive integer format (scaled to 55,000 gray levels) and can be rescaled to Top of Atmosphere (TOA) reflectance and/or radiance using radiometric change of scale coefficients provided in the product metadata file (MTL file). The MTL file also contains the thermal constants required to convert thermal infrared sensor (TIRS) data to the brightness temperature in the satellite [32]. Spectral Radiance in OLI and TIRS Sensor: Images are processed in absolute radiance units using 32-bit floating-point calculations. These values are converted to 16-bit integer values in the finished Level 1 product. They can then be converted to spectral radiance using the reflectance scale factors provided in the metadata file [33]: where L λ = The spectral radiance value of (TOA) measured in values of (Watts/(m 2 * srad * µm)). M L = The specific multiplicative scaling factor (RADIANCE_MULT_BAND_x, where x is the number of the band). A L = The specific additive scaling factor obtained from the metadatum (RADIANCE_ADD_BAND_x, where x is the number of the band). Q cal = Standard product quantified and calibrated by pixel values (DN). This value refers to each of the bands in the image.
Top of Atmosphere Reflectance (TOA) in OLI sensor: The entire 16-bit values in the Level 1 product can also be transformed into TOA reflectance. The following equation is used to turn Level 1 DN values into reflectance [33]: where ρλ = TOA planetary reflectance, without correcting for the solar angle. M ρ = Multiplicative scaling factor for metadata scaling (REFLECTANCE_MULT_BAND_x, where x is the band number). A ρ = Additive scaling factor based on metadata (REFLECTANCE_ADD_BAND_x, where x is the band number). Q cal = Quantified and calibrated product pixel values (DN).
One should consider that ρλ does not have a correction for the solar angle. Therefore, to obtain TOA reflectance in the sun angle, the following conversion should be used: where ρλ = Planetary TOA reflectance. θ SE = Local solar elevation angle; the sun elevation angle of the center of the scene in degree is provided by the metadata as (SUN_ELEVATION). θ SZ = Local solar zenith angle; Obtaining the Brightness Temperature: To calculate the Land Surface Temperature, it is necessary to obtain the brightness temperature of the thermal bands based on the radiance Earth 2023, 4 99 result. The brightness temperature is the effective temperature seen by the satellite under the assumption of unit emissivity [33]. The conversion equation is as follows: where T = Brightness temperature of Top of Atmosphere. L λ = TOA spectral radiance (Watts/(m 2 * srad * µm)). k 1 = Specific thermal conversion content of metadata band(K1_CONSTANT_BAND_x, where x is the thermal band number). k 2 = Specific thermal conversion content of metadata band (K2_CONSTANT_BAND_x, where x is the thermal band number).
Threshold Method based on NDVI to calculate Emissivity: Surface emissivity is a measure of the inherent efficiency of the surface, since it converts thermal energy into radiant energy over the surface. It depends on the composition, roughness, and moisture content of the surface and the observation conditions, in other words, the wavelength, pixel resolution, and viewing angle [34]. Thermal infrared emissivity is therefore an essential parameter for both surface characterization and atmospheric correction methods. The mapping of emissivity from satellite data is therefore a key issue to resolve. The main problem is temperature coupling and the effect of emissivity on thermal radiances [35].
To calculate emissivity, [36] found a high correlation between the emissivity measurement and the Normalized Difference Vegetation Index (NDVI) derived from reflectance in the red and near-infrared bands. This relationship is of potential use in sensors with infrared thermal bands because the NDVI can easily be derived from radiometric spectral measurements found on board many of the operational meteorological satellites. Accordingly, various approaches have been used to predict the emissivity of the Earth's surface from NDVI values [35][36][37]. Therefore, in this study, the Normalized Difference Vegetation Index (NDVI) Threshold method used to calculate the emissivity of the Earth's surface was chosen. It is worth recalling that NDVI is defined as equation 5 [38].
where ρ 2 : the reflectance measured in near-infrared, and ρ 1 : the reflectance corresponding to red wave bands.
For Landsat OLI 8, it is: Once the NDVI has been obtained, the threshold method uses certain NDVI values (Thresholds) to distinguish between soil pixels (NDVI<NDVIs) and complete vegetation pixels (NDVI>NDVIv). Therefore, for pixels composed of soil and vegetation (mixed pixels, NDVIs ≤ NDVI ≤ NDVIv), the method uses the following equation [37].
where ε υ and ε s are, respectively, soil emissions and vegetation; P υ is the proportion of vegetation (also called fractional vegetation cover, FVC); and C is a term that considers the Earth 2023, 4 100 cavity effect due to surface roughness (C = 0 for flat surfaces). The proportion of vegetation is obtained through the equation [39]: where NDV I S = 0.2 y NDV I V = 0.5 for global conditions [40]. However, these values can be obtained from the NDVI histogram [41]. Using the histogram methodology designed by [41], the thresholds for each urban area were obtained in the El Bajío area, as shown in Table 2. Finally, using the geometric model proposed by authors [42], the cavity term (C λ ) for a mixed area and a close view of the nadir is given by the value: where F is a geometric factor that varies between zero and one, depending on the geometric distribution of the surface. Given that F cannot be estimated from remote VNIR/TIR perception data, a mean value is usually chosen. Land Surface Temperature: The Land Surface Temperature was calculated using the equation proposed by [43]: where λ = 10.895 is equal to the wavelength of the radiance emitted; ε λ is the emissivity calculated and ρ = h c σ = 1.438 × 10 −2 mk, whereby σ is Boltzmann's constant (1.3810 23 J/K), h is Planck (6.62610 34 Js); and c is the speed of light (2.99810 8 m/s) [44]. The result is obtained in Kelvin degrees, as a result of which the value of −273.15 is subtracted to convert it to Celsius degrees.
Surface Urban Heat Islands: From the LST calculated in equation [10] and converted to degrees Celsius in the Landsat Diurnal images, the SUHI was obtained using the following equation: LSTI ( • C)= Corresponds to the value obtained from the LST of each pixel within the metropolitan or conurbation area.
LSTO ( • C)= The average LST outside the metropolitan or conurbation area in a defined area of 10 km around each city.
The SUHIs were obtained by the difference in each pixel of the land surface temperature within the metropolitan or conurbation area minus the average LST in a 10 km zone of influence outside the urban area, obtaining a positive SUHI where the Surface Temperature has increased inside the cities. Subsequently, the results were plotted, and each city's minimum, maximum, and average values were tabulated.

Transects
The use of transects in the various types of cover and soil is proposed since, as mentioned earlier [44], the profiles will help to explain the factors that shape the thermal landscape of the city. Figure 3 shows the steps that will be taken throughout this subsection.
Earth 2023, 4, FOR PEER REVIEW 9 (°) = Corresponds to the value obtained from the LST of each pixel within the metropolitan or conurbation area.
(°) = The average LST outside the metropolitan or conurbation area in a defined area of 10 km around each city.
The SUHIs were obtained by the difference in each pixel of the land surface temperature within the metropolitan or conurbation area minus the average LST in a 10 km zone of influence outside the urban area, obtaining a positive SUHI where the Surface Temperature has increased inside the cities. Subsequently, the results were plotted, and each city's minimum, maximum, and average values were tabulated.

Transects
The use of transects in the various types of cover and soil is proposed since, as mentioned earlier [44], the profiles will help to explain the factors that shape the thermal landscape of the city. Figure 3 shows the steps that will be taken throughout this subsection.

Description of Categorization of El Bajío Cities for Transect Plotting
For the transect plotting, straight lines were drawn over each city, considering the 10 km buffer around the urban localities within the metropolitan areas and conurbations. Straight lines were drawn, attempting to go through extreme temperature points and diverse types of land cover and current land use. Distances were subsequently established by city size to place points along the lines and obtain the resulting LST and SUHI values as well as the type of cover coinciding with those points.
Cities were grouped according to the length of the lines considered for the transects and depending on the separation distance between the points placed on the straight lines, taking three groups into account: three kilometers for the metropolitan cities of Querétaro, León, Celaya, and La Piedad-Pénjamo, two kilometers for the metropolitan cities of San Luis Potosí, Aguascalientes, and Guanajuato, and one and a half kilometers for the metropolitan cities of Rioverde, Moroleón-Uriangato, San Francisco del Rincón, and the conurbated cities of San Juan del Río, Irapuato, Salamanca, and San Miguel de Allende. These distances between points were selected to obtain between 15 and 30 points per line to make the transect graphs, to ensure that there was sufficient information without saturating the graphs or hampering their interpretation.

Description of Categorization of El Bajío Cities for Transect Plotting
For the transect plotting, straight lines were drawn over each city, considering the 10 km buffer around the urban localities within the metropolitan areas and conurbations. Straight lines were drawn, attempting to go through extreme temperature points and diverse types of land cover and current land use. Distances were subsequently established by city size to place points along the lines and obtain the resulting LST and SUHI values as well as the type of cover coinciding with those points.
Cities were grouped according to the length of the lines considered for the transects and depending on the separation distance between the points placed on the straight lines, taking three groups into account: three kilometers for the metropolitan cities of Querétaro, León, Celaya, and La Piedad-Pénjamo, two kilometers for the metropolitan cities of San Luis Potosí, Aguascalientes, and Guanajuato, and one and a half kilometers for the metropolitan cities of Rioverde, Moroleón-Uriangato, San Francisco del Rincón, and the conurbated cities of San Juan del Río, Irapuato, Salamanca, and San Miguel de Allende. These distances between points were selected to obtain between 15 and 30 points per line to make the Earth 2023, 4 102 transect graphs, to ensure that there was sufficient information without saturating the graphs or hampering their interpretation.

Land Cover and Current Land Use
To obtain a sample of the transects, it was necessary to identify the land cover and current land use in the study area to determine the behavior of the daytime LST in Landsat in the various categories, and thereby identify how daytime LST interacts with current land cover and use in El Bajío, Mexico. To this end, each state's official census information from the National Housing Inventory (Inventario Nacional de Viviendas, Spanish acronym INV, 2010) was extracted and used. Information on the blocks was subsequently combined with the number of dwellings, restricting it to an area containing both metropolitan and suburban and buffer areas. It was applied to regions in the El Bajío Industrial Corridor, highlighting the VIV0 column, referring to the number of houses for residential use in its respective table. A zero in the box indicates that there is no dwelling assigned as a residence. GlobeLand30 types of land cover were combined with the current land uses to identify each of the changes in LST of the landscapes present inside and outside the metropolitan and suburban areas. The following steps were taken to incorporate the inputs:

•
The DENUE points were incorporated by an analysis of proximity to the blocks. They were therefore converted from points to polygons to identify the blocks in which economic activities were located.

•
The blocks were subsequently classified according to the number of residential dwellings as well as the economic activity or activities that take place within them, Table 3. This category includes tertiary activities involving the distribution of goods such as trade, transportation, services, and information management A category of roads was also constructed by adapting the GlobeLand30 land cover to the size of urban locations, and block areas were subtracted to obtain the roads in polygons. In addition, green areas available in the 2017 Geostatistical Framework were included, comprising objects such as airstrips and medians. Finally, the GlobeLand30 types of cover were integrated with current land use classes and road and green area categories. Once the data had been incorporated, the resulting layer was rasterized.

Transect Design
Once the data resulting from the Landsat LST had been obtained, as well as the layer resulting from the integration of land cover and current land uses, the layers were superimposed over the study area, including the buffer area, and four transects were drawn per city. Every attempt was made to go through the diverse types of cover and certain points of interest within each urban area, such as airports, heat-insulating roofs, bodies of water, and green areas.
Subsequently, obtaining points on the transects separated by the estimated distance of 3, 2, and 1.5 km depending on the city was automated within the processing of the QGIS tool modeler to achieve the sample of raster values of LST, types of cover, and current land use. Figure 4 shows a diagram with the procedure and the tools used for generating points in the transects and extracting values from the rasters.
Earth 2023, 4, FOR PEER REVIEW 11 6 Trade and services This category includes tertiary activities involving the distribution of goods such as trade, transportation, services, and information management A category of roads was also constructed by adapting the GlobeLand30 land cover to the size of urban locations, and block areas were subtracted to obtain the roads in polygons. In addition, green areas available in the 2017 Geostatistical Framework were included, comprising objects such as airstrips and medians. Finally, the GlobeLand30 types of cover were integrated with current land use classes and road and green area categories. Once the data had been incorporated, the resulting layer was rasterized.

Transect Design
Once the data resulting from the Landsat LST had been obtained, as well as the layer resulting from the integration of land cover and current land uses, the layers were superimposed over the study area, including the buffer area, and four transects were drawn per city. Every attempt was made to go through the diverse types of cover and certain points of interest within each urban area, such as airports, heat-insulating roofs, bodies of water, and green areas.
Subsequently, obtaining points on the transects separated by the estimated distance of 3, 2, and 1.5 km depending on the city was automated within the processing of the QGIS tool modeler to achieve the sample of raster values of LST, types of cover, and current land use. Figure 4 shows a diagram with the procedure and the tools used for generating points in the transects and extracting values from the rasters. Once the tools had been automated, this made it possible to process the cities to which the same distances were going to be applied and reduce processing times. Finally, the samples of values were exported in tables to graph the results of the data obtained.

Landsat Daytime Land Surface Temperature 2020
Land Surface Temperatures obtained with Landsat images for each metropolitan area and conurbation are detailed in Appendix A.1; this shows the behavior of the daytime LST obtained from the Landsat images of the metropolitan cities and conurbations of El Bajío, Mexico. The Landsat images have higher resolution and show the areas with high Once the tools had been automated, this made it possible to process the cities to which the same distances were going to be applied and reduce processing times. Finally, the samples of values were exported in tables to graph the results of the data obtained.

Landsat Daytime Land Surface Temperature 2020
Land Surface Temperatures obtained with Landsat images for each metropolitan area and conurbation are detailed in Appendix A.1; this shows the behavior of the daytime LST obtained from the Landsat images of the metropolitan cities and conurbations of El Bajío, Mexico. The Landsat images have higher resolution and show the areas with high or low LST in greater detail. Within the urban area of the cities, those with the maximum temperatures are the conurbation of Salamanca (59.7 • C) and the conurbations of León (53.1 • C) and Irapuato (52.7 • C), while Querétaro (21.6 • C) and León (24.3 • C) are the cities with the minimum LST. It is worth noting that, although the metropolitan area of Querétaro is one of the largest, most industrial cities in El Bajío, it has the lowest temperatures, due to the bodies of water inside the city. The opposite is true of the conurbation of Salamanca which, although the Lerma River runs through it, has a high LST identified in crop areas. In regard to the outside of the urban area of each city, the LST was obtained within a 10 km buffer.
The exterior of the urban area records lower daytime LST in all cities, since, it does not exceed the minimum LST of 30 • C, unlike the interior. It also has lower maximum LSTs than the interior area. As one would expect, areas close to bodies of water have minimum temperatures, as in the city of Querétaro, both inside and outside the urban area; although, the presence of thermal construction materials or specific buildings also causes minimum LSTs as in the metropolitan area of Aguascalientes, both inside and outside the urban area.

Daytime Landsat Surface Urban Heat Island 2020
An analysis of the daytime SUHIs in Appendix A.2 shows that positive values indicate that LST is higher inside the urban area than outside in the buffer, while negative values indicate that the LST is higher in the buffer than in the urban area. As well as, in most of the El Bajío cities, the highest LST occurs within the city, creating a positive SUHI. The Landsat images show a minimum to maximum range of SUHIs from −18.1 to 18.9 • C, whereas for the fourteen cities, the average minimum is −3.4 • C and the average maximum of 2.5 • C, with an average difference of −1 • C. In other words, during the daytime in El Bajío in May 2020, the LST was 1 • C higher outside the urban area than inside. The behavior of SUHIs with Landsat in the El Bajío region can be seen in Figure 5 in their minimum, maximum, and average values by city. It is also possible to see the cities with negative SUHIs and the way some had extremely low minimums, as in the case of León with −18.1 • C or Querétaro with −16.8 • C, whereas Salamanca had a maximum of 18.9 • C. or low LST in greater detail. Within the urban area of the cities, those with the maximum temperatures are the conurbation of Salamanca (59.7 °C) and the conurbations of León (53.1 °C) and Irapuato (52.7 °C), while Querétaro (21.6 °C) and León (24.3 °C) are the cities with the minimum LST. It is worth noting that, although the metropolitan area of Querétaro is one of the largest, most industrial cities in El Bajío, it has the lowest temperatures, due to the bodies of water inside the city. The opposite is true of the conurbation of Salamanca which, although the Lerma River runs through it, has a high LST identified in crop areas. In regard to the outside of the urban area of each city, the LST was obtained within a 10 km buffer.
The exterior of the urban area records lower daytime LST in all cities, since, it does not exceed the minimum LST of 30 °C, unlike the interior. It also has lower maximum LSTs than the interior area. As one would expect, areas close to bodies of water have minimum temperatures, as in the city of Querétaro, both inside and outside the urban area; although, the presence of thermal construction materials or specific buildings also causes minimum LSTs as in the metropolitan area of Aguascalientes, both inside and outside the urban area.

Daytime Landsat Surface Urban Heat Island 2020
An analysis of the daytime SUHIs in Appendix A.2 shows that positive values indicate that LST is higher inside the urban area than outside in the buffer, while negative values indicate that the LST is higher in the buffer than in the urban area. As well as, in most of the El Bajío cities, the highest LST occurs within the city, creating a positive SUHI. The Landsat images show a minimum to maximum range of SUHIs from −18.1 to 18.9 °C, whereas for the fourteen cities, the average minimum is −3.4 °C and the average maximum of 2.5 °C, with an average difference of −1 °C. In other words, during the daytime in El Bajío in May 2020, the LST was 1 °C higher outside the urban area than inside. The behavior of SUHIs with Landsat in the El Bajío region can be seen in Figure 5 in their minimum, maximum, and average values by city. It is also possible to see the cities with negative SUHIs and the way some had extremely low minimums, as in the case of León with −18.1 °C or Querétaro with −16.8 °C, whereas Salamanca had a maximum of 18.9 °C.  The analysis of the Land Surface Temperature and the Surface Urban Heat Islands with daytime Landsat images shows that the LST not only varies according to the urban limit, but that due to its spatial resolution, it gives rise to exploration depending on the types of cover and soils inside and outside the urban limit. In this way, transect plotting is proposed to identify the LST according to the corresponding land cover or use and analyze its behavior.

Transects
We used transects of different distances to measure the SUHI effect in different classes of land cover and current land use. Four transects were built for each city in the study. Points were subsequently calculated every 3, 2, and 1.5 km according to the size of the city on the transects, and samples of the LST, land cover, and current land use were obtained; the results are shown in Appendix B.
In the largest cities with separation points of three kilometers, the transects show a temperature increase during the day in urban areas. The most temperate zones correspond to forests, followed by scrub and crops. In some cases, the transects tend to show SUHI and high LST points in urban areas. However, temperature anomalies are also found in certain crop areas, possibly due to the lack of bodies of water and the corridors of construction areas formed outside the urban limits.
Regarding Querétaro, the daytime images show the high surface temperatures associated with crop areas. There is a need to analyze the agricultural practices and the species planted in these areas. Moreover, it is essential to identify the materials used to build urban residential buildings due to the low temperatures they reflect. However, the behavior is as expected since the existing bodies of water do not significantly reduce the temperature of the surroundings, or at least this is not noticeable in the transects. Moreover, the SUHIs are concentrated in the western part of the city due to the concentration of manufacturing industries in that area ( Figure 6).
Regarding the cities with separation points of two kilometers in their transects, as in the city of Guanajuato, the types of cover with low LST are trees and forests, bushes and scrub, and bodies of water. On the other hand, grasslands, crops, and residential land use increase LST. Due to the type of vegetation surrounding the city of Guanajuato, LST and SUHIs are more noticeable in the south since this area is predominantly under cultivation (Figure 7).
Finally, in the case of Rioverde, high surface temperatures are associated with bushes and scrub in cities with transects that have points every 1.5 km. In contrast, low temperatures are found in wooded areas. The urban area of Rioverde is a transition zone between the two, with no temperature anomalies. High LST values can be seen north of the city of Rioverde, whereas low temperatures are associated with the southwest.
The city of San Miguel de Allende has a significant percentage of water cover outside the urban area, demonstrating how bodies of water keep the surrounding areas of the urban area cool. Another factor in keeping the surrounding areas cool is the advantage of retaining a high percentage of the natural ecosystem around the urban area.
As shown in Figure 8, Transect 2 highlights the importance of water bodies in decreasing temperature and their function as heat sinks for neighboring areas. The SUHIs are therefore concentrated in the northwestern part of San Miguel de Allende; this region is the furthest from an extensive body of water and wooded areas.
The LST values summarized in the form of box plots in Figure 9 are especially representative of the land-cover categories that take up a high percentage of the cities in the El Bajío Industrial Corridor, Mexico. In this case, the crop areas and other land use related to exposed soils devoid of vegetation stand out as they present a higher LST than the green areas found within all the cities included in this study.
We associate the month of May with seasonal crops, which are of greater scope and profitability for many farmers in the country. In that case, there is a lower production of these foods since it also goes through a dry season, generating the presence of bare soils or with little vegetation, which consequently generates an increase in the LST. On the other hand, the dynamics within the cities are different. Many vegetation areas correspond to gardens, parks, and park areas, which maintain constant hydration by the neighborhoods or governmental maintenance programs, keeping the areas with continuous vegetation.
We assume that the crop zones correspond to the Panorama Agrolimentario 2020, as in the rest of the country, where grain corn, beans, sorghum, wheat, vegetables, and greens are the main crops [45]. On the other hand, the variety of techniques and conditions of cultivation must be considered, such as fallow, irrigation, rainfed, or/and greenhouse.    or with little vegetation, which consequently generates an increase in the LST. On the other hand, the dynamics within the cities are different. Many vegetation areas correspond to gardens, parks, and park areas, which maintain constant hydration by the neighborhoods or governmental maintenance programs, keeping the areas with continuous vegetation.
We assume that the crop zones correspond to the Panorama Agrolimentario 2020, as in the rest of the country, where grain corn, beans, sorghum, wheat, vegetables, and greens are the main crops [45]. On the other hand, the variety of techniques and conditions of cultivation must be considered, such as fallow, irrigation, rainfed, or/and greenhouse.

Discussion
Changes in land use are related to the context of anthropogenic climate change and global warming, meaning that their understanding and active management are increasingly important [46]. As noted in studies such as those by [47,48], the occurrence of Urban Surface Heat Islands depends on a range of factors, including large-scale and mesoscale climate conditions, microclimate, geographic location, the biophysical characteristics of the land surface, urban structure, population size, wind corridors, human lifestyle, and energy consumption. Apparently, zoning, as an urban planning tool, profoundly affects the physical characteristics of urban landscapes by imposing restrictions such as maximum building height and density, the extent of impervious surface and open space, and types of land use and activities. These variables, in turn, control surface energy exchange,

Discussion
Changes in land use are related to the context of anthropogenic climate change and global warming, meaning that their understanding and active management are increasingly important [46]. As noted in studies such as those by [47,48], the occurrence of Urban Surface Heat Islands depends on a range of factors, including large-scale and mesoscale climate conditions, microclimate, geographic location, the biophysical characteristics of the land surface, urban structure, population size, wind corridors, human lifestyle, and energy consumption. Apparently, zoning, as an urban planning tool, profoundly affects the physical characteristics of urban landscapes by imposing restrictions such as maximum building height and density, the extent of impervious surface and open space, and types of land use and activities. These variables, in turn, control surface energy exchange, surface and subsurface hydrology, microto mesoscale weather and climate systems, and other environmental processes [49].
The heat island effect is mainly due to a combination of factors, including the concentration of heat-absorbing surfaces, the lack of vegetation, and the presence of heat-generating sources. The view of profiles or transects in the different cities, therefore, enables one to observe the impact of these zonings on well-defined transects that show changes in land use and their association with changes in surface temperature, such as works in [50,51].
Within Mexican cities, [50] present a study of a city in northeastern Mexico. Characterized by large surfaces exposed to the sun, the study of this city shows a predominance of the city's high temperatures due to pavements and buildings with scarce vegetation within the urban area. All measurements were obtained during the warm period in transects in the main roads of the urban layout, to which there is a significant presence of heavy industry, light industry, and large commercial corridors. In a similar study in Toluca city in central Mexico, [25] concluded that the high temperatures correspond to the area with the most extensive urban infrastructure located towards the city's center. The lack of vegetation, the scarcity of bodies of water, multiple buildings, and construction with heat-absorbing materials, and a large concentration of vehicles characterize this area. Outside of cities, the heat island effect is generally less pronounced. The lower density of buildings and pavement and the presence of vegetation and agricultural activity can help dissipate heat and cool the air. However, fallow agricultural areas can also have the opposite effect and increase the urban heat island effect [52,53].
A study by Gough in 2020 [54] found that if fallow land is without proper management and the vegetation is allowed to dry out or die, the urban heat island effect can increase. In this case, the land can become a significant heat source rather than a heat sink and contribute to the urban heat island effect. Therefore, proper management and planning are crucial to mitigate potential adverse effects on the urban heat island [53][54][55].
Urban Heat Island studies are generally carried out mainly in two ways: first, by measuring UHI in air temperature using car transects and weather station networks; and second, by measuring UHI in surface temperature using remote sensing. The vast majority of climatological studies of UHI in Mexico have been carried out using in situ data, which have the advantage of high temporal resolution and an extensive data record; however, they have a poor spatial resolution and involve high costs due to measurement instruments acquisition and regular field trips. Whereas satellite thermal data, as in this case, have more excellent spatial coverage, allowing identification of UHI in a more versatile way, with the advantage of having more spatial and even temporal information at different coverage scales, such as local and regional.
Utilizing satellite data has become increasingly important for researchers and decision makers to understand the climatic effects of urbanization and, more importantly, how to contribute to more sustainable urban development. In such a way, the results of UHI studies between cities with different climatic and socioeconomic environments, such as ours, use data that allow for the identification of areas vulnerable to UHI, which, in turn, can be related to other fields of study, such as health, by identifying diseases that have developed in specific areas of the cities or sectors of the population derived from UHI. This research work is a precedent for other studies in the methodology employed and how the results are presented in terms of comparing the same region in different temporalities or regions with similar landcover characteristics. In addition, it can serve as theoretical support material to forecast future behavior and propose measures that will allow, if not diminish, the effect of the UHI to be controlled or redistributed through implementing public policy strategies.

Conclusions
Through the incorporation of thermal profiles, the proposed methodology makes it possible to identify the land uses that intensify the LST, providing detailed information on the variation with distinct types of land cover. It also details the importance, for example, of green and blue areas, as well as crop areas within urban environments. In the present study, we used box plots to show the distribution of LST transects of different distances in different land use types. We analyzed the urban heat island effect in cities of the El Bajío Industrial Corridor in Mexico.
With respect to the case study, the peri-urban expansion that increases the use of urban land in a gradient ranging from natural ecosystems to crop areas and built areas (such as buildings, industries, housing, and roads), significantly affects the temperature change in cities.
The characteristics of the cities in El Bajío that most intensify the effect of the SUHIs have greater urban conglomeration; in other words, the built areas are not scattered. Other characteristics include having limited vegetation or green areas and a sizable industrial or mixed area, with few or no bodies of water, which have been seen to operate as temperature regulators. Furthermore, located around urban areas, crops such as vegetables, red fruits, corn, wheat, and sorghum that use fallow as part of agricultural practices minimize natural areas and directly contribute to the increase in temperature in the SUHI. This warrants a more specific study to identify the species causing this to find a way to mitigate this effect.
Finally, processing the images with the tools described proved helpful in automating the process. Moreover, since the tools used in this study are affordable for everyone, and given that the software is free and open source, the proposed procedure can be widely used in other cities. Earth 2023, 4, FOR PEER REVIEW 29 San Miguel de Allende Irapuato Salamanca San Juan del Río  Guanajuato Figure A3. Transects of the Ten Metropolitan Cities and Four Conurbations, with the SUN code. Figure A3. Transects of the Ten Metropolitan Cities and Four Conurbations, with the SUN code.