Pronounced Changes in Thermal Signals Associated with the Pronounced Changes in Thermal Signals Associated with the Madoi (China) M 7.3 Earthquake from Passive Microwave and Madoi (China) M 7.3 Earthquake from Passive Microwave and Infrared Satellite Data Infrared Satellite Data

: Thermal variations in surface and atmosphere observed from multiple satellites prior to strong earthquakes have been widely reported ever since seismic thermal anomalies were discovered three decades ago. These thermal changes are related to stress accumulation caused by the tectonic activities in the ﬁnal stage of earthquake preparation. In the present paper, we focused on the thermal changes associated with the 2021 Madoi M 7.3 earthquake in China and analyzed the temporal and spatial evolution of the Index of Microwave Radiation Anomaly (IMRA) and the Index of Longwave Radiation Anomaly (ILRA) based on 8-year microwave brightness temperature (MWBT) and 14-year outgoing longwave radiation (OLR) data collected by satellites. We also explored their responses in different tectonic units (seismogenic fault zone and active tectonic block). Our results indicated that the enhanced IMRA was distributed along the seismogenic fault since mid-February and reappeared for a longer time and with stronger intensity in March and April 2021. The pronounced enhancement in the ILRA was observed within one month over Bayan Har tectonic and adjacent blocks. The higher ILRA over the tectonic blocks in the southern Tibet Plateau at the beginning of 2021 could be associated with the regional stress accumulation, as proven by the occurrences of two moderate earthquakes during this period.


Introduction
Thermal changes are observed in tandem with changes in stress in seismically active regions.The anomalous thermal variations triggered by strong earthquakes have been widely observed through satellite thermal infrared images following the work of Gornyi et al. [1], as satellite observation was starting to be considered to be a new tool for the investigation of seismically active regions [2].With the development of satellite observation technology, efforts have been made to use thermal-related satellite data to detect anomalous signals associated with earthquakes, such as infrared/microwave brightness temperature [3-5], surface temperature [6][7][8][9], air temperature [10], outgoing longwave radiation (OLR) [11,12], and surface latent heat flux (SLHF) [13,14], in different parts of the globe.These observations have led to belief in the existence of anomalous thermal signals associated with the earthquake preparation phase.These anomalous thermal signals have been detected months and days prior to strong earthquakes in the epicentral areas and along faults.Thermal anomalies prior to earthquakes are associated with regional geodynamics (e.g., types of seismogenic fault) [15,16].Laboratory experiments [17] have demonstrated the existence of infrared emissions when positive-holes (p-holes) recombine on the surface.Prior to major seismic activity, p-holes are activated in rocks via the rupture of peroxy bonds and recombination on the Earth's surface due to the enhancement of tectonic stress.Such changes could induce a series of physical and chemical processes, e.g., thermal anomalies, air ionization, and increased carbon monoxide [18].The models of Lithosphere-Atmosphere-Ionosphere Coupling (LAIC) [19] and Lithosphere-Coversphere-Atmosphere Coupling (LCAC) [20] have been proposed to understand the physical mechanism.Recently, a method for the comprehensive multi-parameter analysis of the lithosphere, atmosphere and ionosphere that combines thermal-related parameters with other meteorological, atmospheric and ionospheric parameters has been used to study the geosphere coupling during the earthquake preparation stage [21][22][23], and this methodology has proven to be an effective method to weaken the false signals from single parameter models [24,25].In addition, a bottomup pattern of anomalous signals (from the lithosphere to the atmosphere or even in the ionosphere) prior to impending earthquakes has been hypothesized [21,26,27].
Changes in microwave brightness temperature (MWBT) are directly dependent on the surface and subsurface temperature and emissivity [28].In contrast to thermal infrared signals, microwave signals can penetrate clouds and provide information on the surface even under poor weather conditions.In our earlier study, MWBT data were used to detect anomalous signals associated with the major earthquakes that occurred on the boundary of the Bayan Har tectonic block, i.e., the 1997 Manyi M 7. 5, 2001 Kokoxili M 7.8 [4], 2008 Wenchuan M 7.9, 2013 Lushan M 6.6, and 2017 Jiuzhaigou M 6.5 earthquakes [29].The results showed reasonable consistency with the increased MWBT along the seismogenic fault from a few months to a few days prior to the earthquakes, depending on the focal mechanism and seismogenic environment.Outgoing longwave radiation (OLR) provides the radiation information from the top of the atmosphere.An anomalous OLR variation in the epicentral region has been considered to be an important indicator of an impending earthquake [11,12,30], although it is difficult to relate to seismogenic faults due to low spatial resolution.Zhang et al. [31] presented statistical evidence supporting the relation between earthquakes and OLR anomalies.In the present study, two thermal-related parameters (MWBT and OLR) with different spatial resolutions observed from different satellites were combined to analyze the thermal variations associated with the 2021 Madoi M 7.3 earthquake in China.The thermal variations in tectonic units with different spatial scales (seismogenic fault zone and active tectonic block) were investigated according to the different spatial resolutions of the two datasets.Our detailed analysis demonstrated a pronounced change in the thermal signals prior to the earthquake; such changes could be used in the future to help alert people about impending earthquakes in the region.

Madoi (China) M 7.3 Earthquake
On 21 May 2021 (18:04:13 UTC), a strong earthquake (Mw 7.3, USGS) hit the northern Tibetan Plateau and the Bayan Har tectonic block (hereafter referred to as the Madoi M 7.3 earthquake).The epicenter was located at 98.34 • E, 34.59 • N, with a focal depth of 17 km according to China Earthquake Networks Center (CENC).The focal mechanism of the main earthquake event presented an NW-striking, left-lateral strike-slip fault with a normal faulting component (Kunlunshankou-Jiangcuo Fault) [32,33].The main event ruptured ~156 km on the main fault according to InSAR and GPS measurements [34].The Madoi M 7.3 earthquake was the largest earthquake within the Bayan Har block since 1947 (according to USGS earthquake catalogue), although a few major earthquakes (the 1997 Manyi M 7.5, 2001 Kokoxili M 7.8, 2008 Wenchuan M 7.9, 2010 Yushu M 6.9, 2013 Lushan M 6.6, and 2017 Jiuzhaigou M 6.5 earthquakes) occurred in the boundary of the Bayan Har tectonic block over the past 30 years (Figure 1).A number of sand liquefactions in different locations were observed in the epicentral region after the Madoi M 7.3 earthquake [35] associated with changes in the stress regime.Such changes are the cause of increasing turbidity in the epicentral region water bodies that have been detected by Sentinel SAR data [36].The thermal-related signals associated with this earthquake have been reported using the MODIS infrared brightness temperature (IRBT) and AMSR2 microwave brightness temperature (MWBT) data [37,38].The infrared variation has shown that an increased IRBT has covered most areas of Tibetan Plateau since October 2020 [37], and an enhancement in MWBT within the Bayan Har block, which is more pronounced compared to the outside of this block, has been observed since January 2021 [38]; this may be associated with the localized activities in this block.

Satellite Data
In the present work, MWBT data collected by the Advanced Microwave Scanning Radiometer 2 (AMSR2) instrument, onboard the Japanese Space Exploration Agency Global Change Observation Mission 1st-Water (GCOM-W1) satellite, were used.AMSR2 is the successor of the Advanced Microwave Scanning Radiometer (AMSR) and the Advanced Microwave Scanning Radiometer for EOS Aqua (AMSR-E), which provides global observation data at multiple frequencies (6.9, 7.3, 10.7, 18.7, 23.8, 36.5 and 89.0 GHz) and multiple polarizations (vertical and horizontal).Descending brightness temperature data with a spatial resolution of 10 km from AMSR2 Level 3 products were considered to avoid solar radiation during the daytime.
The 1.0-degree OLR dataset computed from infrared channel of the NOAA polar orbit satellite was used to study thermal variations associated with strong earthquake activities in view of its relationship with the thermal dynamic process prior to strong earthquakes [11,12].Additionally, in this case, only night observation was considered to reduce the influence of solar radiation.

Methods
The methods of anomalous signals detection for MWBT and OLR were derived from the Robust Satellite Technique (RST) proposed by Tramutoli et al. [40].The RST approach was applied to different satellite data (thermal and microwave; polar-orbiting and geostationary) to identify the possible anomalous signals related to natural disaster events, such as earthquakes [41,42], volcanoes [43,44], floods [45,46] and dust storms [47][48][49].
The RST approach is focused on the changes between the value at the target period (e.g., day week, and month) and the multi-year background values.The background values, including the temporal means µ t (x, y) and the standard deviation σ t (x, y) for a target period t (e.g., day, a few days, a week, or a month), were computed.µ t (x, y) and σ t (x, y) are defined by the following Equations ( 1) and (2): where V t,n (x, y) is the parameter value observed using a satellite at the target period t in the year n and at the location with latitude x and longitude y.N is the number of the years used as the background.
The anomalous variation S t (x, y) was computed with Equation (3): The IMRA (Index of Microwave Radiation Anomaly) and ILRA (Index of Longwave Radiation Anomaly) were defined for MWBT and OLR data based on Equations ( 1)-(3).According to the data availability, the years from 2007 to 2020 and from 2013 to 2020 were considered background years for OLR and MWBT data, respectively.
As microwave signals are strongly dependent on soil moisture and surface roughness, we retrieved AMSR2 10 km soil moisture and MODIS 0.05-degree normalized difference vegetation index (NDVI) products in Southwest China.We computed time-series data on the monthly area-averaged soil moisture and NDVI in the epicentral region (5 • × 5 • ) during 2021 (Figure 2a).The soil moisture and NDVI prior to the Madoi earthquake were found to be lower than those after the event, although soil moisture slightly increased in May.Furthermore, we computed the 5-month-averaged (January to May in 2021) soil moisture and NDVI around the epicenter (Figure 2b,c).According to the spatial distributions for the soil moisture and NDVI, we used the region with relatively lower soil moisture and less vegetation for the spatial analysis of MWBT variations (yellow box in Figure 2b,c).An anomalous MWBT signal prior to the earthquake was observed in low-soil-moisture and NDVI environments using DMSP SSM/I 19 GHz with horizontal polarization [4].Here, we used AMSR2 MWBT data at 18.7 GHz horizontal polarization to analyze the possible IMRA variation associated with the Madoi earthquake.

Microwave Brightness Temperature Variation
Three-day-averaged IMRA values were computed to avoid missing data due to th satellite orbital gaps.Earlier, we found high IMRA values mainly distributed along th

Microwave Brightness Temperature Variation
Three-day-averaged IMRA values were computed to avoid missing data due to the satellite orbital gaps.Earlier, we found high IMRA values mainly distributed along the seismogenic fault in other major earthquakes in the Bayan Har tectonic block [4,29].Here, we firstly focused on the IMRA variation within the seismogenic fault (the region surrounded by green lines in Figure 2b).To obtain the extent of the enhanced IMRA in the seismogenic fault, the percentages of the IMRA pixels with various levels (IMRA > 1, IMRA > 2, IMRA > 3 and IMRA > 4) were computed within the fault zone from January to July 2021 (Figure 3a).By comparing the different IMRA limits, we found IMRA > 2 criteria to be an ideal choice for analyzing anomalous signals associated with the Madoi earthquake.At a lower threshold value (i.e., IMRA > 1), high percentages occurred frequently but were difficult to relate to earthquake activities.At higher threshold values (i.e., IMRA > 3 and IMRA > 4), possibly useful information was masked, which can be seen in the limited data indicated by blue and black bars in Figure 3a.It should be noted that threshold values could vary from place to place due to different seismogenic environments [50].
Remote Sens. 2022, 14, x FOR PEER REVIEW 6 of 19 one.The most notable increase in the IMRA was observed in the third phase (phase III), i.e., 21-26 April 2021.Different from phases I to III, the fourth phase (phase IV) in early-May lasted for a short time and immediately disappeared.After the occurrence of the main earthquake event, four slight enhancements of the IMAR were found.Considering that the MWBT increases with decreases in soil moisture, area-averaged Global Land Data Assimilation System (GLDAS) soil moisture data with a 0.25-degree spatial resolution in this region were used to understand the increased IMRA (Figure 3c).It was clearly demonstrated that the increased IMRA after the earthquake and the increase in phase IV coincided with the decrease in soil moisture (yellow bars in Figure 3), although the decrease in soil moisture during phase IV was not huge.The increase in the IMRA may have lagged behind the decreased soil moisture by 1-2 days.However, phases I to III with an enhanced IMRA prior to the Madoi earthquake occurred with the stable soil moisture environment (pink bars in Figure 3), demonstrating that the IMRA was enhanced independently of soil moisture during these periods.Finally, we considered the IMRA variations with the IMRA > 2 values.The enhanced intensity was defined by using the percentage multiplied by the average value for the pixels with IMRA > 2 to remove the random enhancement of individual pixels.For example, instead of only analyzing the average value, if only one pixel has a high IMRA value (that could be random noise), it can be removed by using this method.Figure 3b shows the time series of IMRA variations within the seismogenic fault zone before and after the Madoi M 7.3 earthquake (highlighted by the red star) using IMRA > 2 values.Four phases with the enhanced IMRA were observed prior to the main event on 21 May 2021 (Figure 3b).The first and second phases (phase I and phase II) occurred in mid-February and mid-March 2021, respectively.The second enhancement lasted for longer time than the first one.The most notable increase in the IMRA was observed in the third phase (phase III), i.e., 21-26 April 2021.Different from phases I to III, the fourth phase (phase IV) in early-May lasted for a short time and immediately disappeared.After the occurrence of the main earthquake event, four slight enhancements of the IMAR were found.Considering that the MWBT increases with decreases in soil moisture, area-averaged Global Land Data Assimilation System (GLDAS) soil moisture data with a 0.25-degree spatial resolution in this region were used to understand the increased IMRA (Figure 3c).It was clearly demonstrated that the increased IMRA after the earthquake and the increase in phase IV coincided with the decrease in soil moisture (yellow bars in Figure 3), although the decrease in soil moisture during phase IV was not huge.The increase in the IMRA may have lagged behind the decreased soil moisture by 1-2 days.However, phases I to III with an enhanced IMRA prior to the Madoi earthquake occurred with the stable soil moisture environment (pink bars in Figure 3), demonstrating that the IMRA was enhanced independently of soil moisture during these periods.
To further confirm the relation between the increased IMRA and the Madoi earthquake in 2021, the IMRA variations in the seismogenic fault zone during the same period in 2019 and 2020, during which no earthquakes with M > 4 occurred in the study region, were analyzed.No similar high IMRA values were observed in the two years (Figure 4a).IMRA variations caused by soil moisture during the months from May to July in 2021 were not observed due to more stable changes in the soil moisture in 2019 and 2020 compared to 2021 in the same period (Figure 4b).To further confirm the relation between the increased IMRA and the Madoi earthquake in 2021, the IMRA variations in the seismogenic fault zone during the same period in 2019 and 2020, during which no earthquakes with M > 4 occurred in the study region, were analyzed.No similar high IMRA values were observed in the two years (Figure 4a).IMRA variations caused by soil moisture during the months from May to July in 2021 were not observed due to more stable changes in the soil moisture in 2019 and 2020 compared to 2021 in the same period (Figure 4b).The spatial distributions of the IMRA around the epicentral region during phases I-III were also obtained for spatial analysis (Figure 5).During phase I (in mid-February 2021), high IMRA values were distributed along not only the seismogenic fault zone but also the boundary of the Bayan Har tectonic block (Figure 5a,b).However, during phases The spatial distributions of the IMRA around the epicentral region during phases I-III were also obtained for spatial analysis (Figure 5).During phase I (in mid-February 2021), high IMRA values were distributed along not only the seismogenic fault zone but also the boundary of the Bayan Har tectonic block (Figure 5a,b).However, during phases II and III (in mid-April and late-May), high IMRA values were distributed along the seismogenic fault zone (Figure 5c-h).The enhancement in the IMRA was observed along the faults in the adjacent block of Bayan Har (Figure 5a,d-h).It is worth noting that the distribution of the high IMRA values was closer to the seismogenic fault during phase III (Figure 5g,h).

Outgoing Longwave Radiation Variation
Compared to the 10 km (0.1 degree) AMSR2 MWBT data, NOAA daily 1.0-degree OLR data have a lower spatial resolution, which makes them more suitable to study thermal variations in large-scale tectonic units.In the present study, the long-term OLR variations in five tectonic blocks (see their locations in Figure 1) were analyzed.The area-

Outgoing Longwave Radiation Variation
Compared to the 10 km (0.1 degree) AMSR2 MWBT data, NOAA daily 1.0-degree OLR data have a lower spatial resolution, which makes them more suitable to study thermal variations in large-scale tectonic units.In the present study, the long-term OLR variations in five tectonic blocks (see their locations in Figure 1) were analyzed.The area-averaged daily ILRA values in five tectonic blocks (Bayan Har, Qiangtang, Qaidam, Lhasa, and Sichuan-Yunnan) during the months from May to June in 2021 were computed using the method discussed in the Section 3.2 (Figure 6).To eliminate regional differences, a modified ILRA was obtained for comparison by subtracting the time-averaged ILRA in each block.We used the mean ILRA plus two standard deviations as a threshold to identify the unusual ILRA variations in each block (dashed lines in Figure 6).The maximum ILRA enhancement (red dashed box in Figure 6) was observed on 4 May 2021 (17 days prior to the main earthquake event) in the Bayan Har (BH) block, which is the tectonic block where the epicenter is located.Around this time, enhanced ILRA values were also observed in the two adjacent blocks of the Bayan Har block, i.e., the Qiangtang (QT) and Qaidam (QD) blocks (blue and green dashed boxes, respectively, in Figure 6), while the intensity in the adjacent blocks was less than that in the Bayan Har block (as can be seen from the length of dashed box).We also observed increased ILRA values in the Lhasa (LS) and Qiangtang (QT) blocks in early March, and these may be associated with two moderate earthquakes that occurred in the Qiangtang block in the same month (M 5. each block.We used the mean ILRA plus two standard deviations as a threshold to identify the unusual ILRA variations in each block (dashed lines in Figure 6).The maximum ILRA enhancement (red dashed box in Figure 6) was observed on 4 May 2021 (17 days prior to the main earthquake event) in the Bayan Har (BH) block, which is the tectonic block where the epicenter is located.Around this time, enhanced ILRA values were also observed in the two adjacent blocks of the Bayan Har block, i.e., the Qiangtang (QT) and Qaidam (QD) blocks (blue and green dashed boxes, respectively, in Figure 6), while the intensity in the adjacent blocks was less than that in the Bayan Har block (as can be seen from the length of dashed box).We also observed increased ILRA values in the Lhasa (LS) and Qiangtang (QT) blocks in early March, and these may be associated with two moderate earthquakes that occurred in the Qiangtang block in the same month (M 5.  Furthermore, we obtained ILRA spatial variations during the period, and these showed distinct enhancements prior to the Madoi M 7.3 earthquake in ILRA time-series variations (red dash boxes in Figure 6).We clearly observed an increase in the ILRA on 29 April 2021 (22 days prior to the main event) and distributions around the epicentral region and its southwest region.Three days later, the distribution of the high ILRA extended to west and gradually covered most area of the Bayan Har block on 4 May 2021, which was the maximum enhancement in the Bayan Har block, as shown in Figure 6.Afterwards, the high ILRA values moved close to the epicenter and diffused to the north on 7-8 May 2021; then, they disappeared (Figure 7).It should be noted that the ILRA increase during this period occurred just after phase III of the increased IMRA (Figure 3b).Furthermore, we obtained ILRA spatial variations during the period, and these showed distinct enhancements prior to the Madoi M 7.3 earthquake in ILRA time-series variations (red dash boxes in Figure 6).We clearly observed an increase in the ILRA on 29 April 2021 (22 days prior to the main event) and distributions around the epicentral region and its southwest region.Three days later, the distribution of the high ILRA extended to west and gradually covered most area of the Bayan Har block on 4 May 2021, which was the maximum enhancement in the Bayan Har block, as shown in Figure 6.Afterwards, the high ILRA values moved close to the epicenter and diffused to the north on 7-8 May 2021; then, they disappeared (Figure 7).It should be noted that the ILRA increase during this period occurred just after phase III of the increased IMRA (Figure 3b).

Discussion
In the present study, we analyzed MWBT and OLR data before and after the 2021 Madoi M 7.3 earthquake to find the possible relationship between the thermal anomalous variations and active tectonic units with different scales using the different spatial resolutions of the observation data from different satellites.Our results clearly demonstrate MWBT variation in the seismogenic areas (Figure 3) and OLR variation in the Bayan Har tectonic block (Figure 6) in good response to the main earthquake event on 21 May 2021.
Qi et al. [38] found that increases in MWBT in the Bayan Har tectonic block were more pronounced than its outside based on the monthly mean background, which was different than our method.Here, we further computed the area-averaged IMRA > 2 with time in the Bayan Har block and its eastern part (the part of the Bayan Har block in the yellow box of Figure 2b) and compared these values with the IMRA variations within the seismogenic fault zone (Figure 8).The area in the Bayan Har block but outside the yellow box in Figure 2 was not considered for computing due to a large fluctuation of soil moisture.We observed that slight increases in the area-averaged IMRA appeared intermittently in the whole block from late January to early May in 2021, but it showed lower values after the main event (Figure 8a).For the eastern Bayan Har block, the IMRA increased in the months of January and late-March, and the increase in late-March showed the maximum variation compared to other days (green bar in Figure 8b).The high IMRA trends in the Bayan Har block and its eastern part seemed to be synchronized, but the latter was more distinct.However, for the seismogenic fault zone, the most obvious IMRA variation occurred in late-April (red bar in Figure 8c), and no high IMRA values were observed in the month of January 2021 (yellow bar in Figure 8c) compared to the changes in the Bayan Har block and its eastern part (yellow bar in Figure 8a,b).The slight increase in the IMRA during January 2021 observed in the Bayan Har block and its eastern part could be associated with the stress accumulation in this tectonic block during this period.In late March, we observed the maximum IMRA enhancement within the eastern Bayan Har tectonic block, which may reflect the enhancement of regional tectonic activity.The two moderate earthquakes (M 5.5 and M 5.7) in March 2021 occurred in its southern block, i.e., Qiangtang block, a few days following this maximum IMRA enhancement.The clockwise rotation of eastern Tibetan Plateau was caused by the collision between the Indian plate and the Eurasian plate, as proven with GPS observations [51], which could make the eastern Bayan Har block more vulnerable to the seismic activities in its southern block due

Discussion
In the present study, we analyzed MWBT and OLR data before and after the 2021 Madoi M 7.3 earthquake to find the possible relationship between the thermal anomalous variations and active tectonic units with different scales using the different spatial resolutions of the observation data from different satellites.Our results clearly demonstrate MWBT variation in the seismogenic areas (Figure 3) and OLR variation in the Bayan Har tectonic block (Figure 6) in good response to the main earthquake event on 21 May 2021.
Qi et al. [38] found that increases in MWBT in the Bayan Har tectonic block were more pronounced than its outside based on the monthly mean background, which was different than our method.Here, we further computed the area-averaged IMRA > 2 with time in the Bayan Har block and its eastern part (the part of the Bayan Har block in the yellow box of Figure 2b) and compared these values with the IMRA variations within the seismogenic fault zone (Figure 8).The area in the Bayan Har block but outside the yellow box in Figure 2 was not considered for computing due to a large fluctuation of soil moisture.We observed that slight increases in the area-averaged IMRA appeared intermittently in the whole block from late January to early May in 2021, but it showed lower values after the main event (Figure 8a).For the eastern Bayan Har block, the IMRA increased in the months of January and late-March, and the increase in late-March showed the maximum variation compared to other days (green bar in Figure 8b).The high IMRA trends in the Bayan Har block and its eastern part seemed to be synchronized, but the latter was more distinct.However, for the seismogenic fault zone, the most obvious IMRA variation occurred in late-April (red bar in Figure 8c), and no high IMRA values were observed in the month of January 2021 (yellow bar in Figure 8c) compared to the changes in the Bayan Har block and its eastern part (yellow bar in Figure 8a,b).The slight increase in the IMRA during January 2021 observed in the Bayan Har block and its eastern part could be associated with the stress accumulation in this tectonic block during this period.In late March, we observed the maximum IMRA enhancement within the eastern Bayan Har tectonic block, which may reflect the enhancement of regional tectonic activity.The two moderate earthquakes (M 5.5 and M 5.7) in March 2021 occurred in its southern block, i.e., Qiangtang block, a few days following this maximum IMRA enhancement.The clockwise rotation of eastern Tibetan Plateau was caused by the collision between the Indian plate and the Eurasian plate, as proven with GPS observations [51], which could make the eastern Bayan Har block more vulnerable to the seismic activities in its southern block due to interactions between adjacent tectonic blocks.Comparing the IMRA changes within the Bayan Har block and its eastern part shows that the maximum IMRA observed in the seismogenic fault occurred in late April 2021, which was closer to the date of the main earthquake (Figure 8c).The occurrence of maximum enhanced intensity (IMRA > 2) gradually approached the day of the main event from large-scale (tectonic block) to small-scale (seismogenic fault) tectonic units (blue arrows in Figure 8).to interactions between adjacent tectonic blocks.Comparing the IMRA changes within the Bayan Har block and its eastern part shows that the maximum IMRA observed in the seismogenic fault occurred in late April 2021, which was closer to the date of the main earthquake (Figure 8c).The occurrence of maximum enhanced intensity (IMRA > 2) gradually approached the day of the main event from large-scale (tectonic block) to small-scale (seismogenic fault) tectonic units (blue arrows in Figure 8).Tectonic stress changes can enhance the emission of carbonaceous gas (e.g., CO or CO2) from shallow or deep crust.Positive spatial correlations between Earth degassing and tectonic regimes and seismically active zones were proven by, respectively, Tamburello et al. [52] and Chiodini et al. [53].Anomalous variations in CO prior to earthquakes have been observed from satellite sensors [21,54,55].Many researches have reported that in increasing stress environments, water fluctuates can cause emissions of gases near the Earth's crust that enhance CO and other gases in epicentral areas [26,54,56].In the present work, we analyzed the near-surface CO variations from the AIRS (Atmospheric Infrared Tectonic stress changes can enhance the emission of carbonaceous gas (e.g., CO or CO 2 ) from shallow or deep crust.Positive spatial correlations between Earth degassing and tectonic regimes and seismically active zones were proven by, respectively, Tamburello et al. [52] and Chiodini et al. [53].Anomalous variations in CO prior to earthquakes have been observed from satellite sensors [21,54,55].Many researches have reported that in increasing stress environments, water fluctuates can cause emissions of gases near the Earth's crust that enhance CO and other gases in epicentral areas [26,54,56].In the present work, we analyzed the near-surface CO variations from the AIRS (Atmospheric Infrared Sounder) onboard on AQUA satellite using the method proposed by Jing and Singh [25].The Index of CO Anomaly (ICOA) at the single pressure level was computed using Equations ( 1)-(3).The differences in the ICOA between the high pressure level corresponding to the low altitude and the low pressure corresponding to high altitude were obtained and showed the anomalous variations from lithospheric contribution.Here, we considered 500 and 200 hPa pressure levels since the epicentral region is in a plateau area.It was interesting to observe two distinct enhancements in the ICOA (exceeding the threshold line of ICOA = 3, as marked with red circle in Figure 9) on 23 February and 23 March 2021 over the epicentral region (Figure 9) following two phases with an increased IMRA within the seismogenic fault zone (phases I and II in Figure 3b, respectively).The enhancement in the ICOA may indicate high emissions of gas from the Earth's crust due to the tectonic activities along the fault during this period.We did not observe an increase in CO during phase III in Figure 3b, which could have been because the residence time of gases in the atmosphere depend on the weather condition at that time (e.g., wind speed, air temperature, and air humidity).Geochemical process changes during the different phases prior to the earthquake commencement could be another reason.
Sounder) onboard on AQUA satellite using the method proposed by Jing and Singh [25].The Index of CO Anomaly (ICOA) at the single pressure level was computed using Equations ( 1)-(3).The differences in the ICOA between the high pressure level corresponding to the low altitude and the low pressure corresponding to high altitude were obtained and showed the anomalous variations from lithospheric contribution.Here, we considered 500 and 200 hPa pressure levels since the epicentral region is in a plateau area.It was interesting to observe two distinct enhancements in the ICOA (exceeding the threshold line of ICOA = 3, as marked with red circle in Figure 9) on 23 February and 23 March 2021 over the epicentral region (Figure 9) following two phases with an increased IMRA within the seismogenic fault zone (phases I and II in Figure 3b, respectively).The enhancement in the ICOA may indicate high emissions of gas from the Earth's crust due to the tectonic activities along the fault during this period.We did not observe an increase in CO during phase III in Figure 3b, which could have been because the residence time of gases in the atmosphere depend on the weather condition at that time (e.g., wind speed, air temperature, and air humidity).Geochemical process changes during the different phases prior to the earthquake commencement could be another reason.Additionally, the daily area-averaged ILRA in 2021 for the five tectonic blocks in the Tibetan Plateau was computed and compared to analyze the relation between the OLR variations and earthquake activities in different tectonic blocks.The ILRA values in 2019 were also computed for comparison since no earthquake with a magnitude of greater than 5.5 occurred in the Tibetan Plateau.The ILRA was modified by subtracting the average annual ILRA value to remove possible regional variations.An increased ILRA was observed in the Bayan Har (BH) and its adjacent blocks, i.e., the Qiangtang (QT) and Qaidam (QD) tectonic blocks, about 2 weeks prior to the Madoi M 7.3 earthquake (red arrows in Figure 10a-c), but no similar behaviors were observed in the Lhasa (LS) and Sichuan-Yunnan (SY) blocks (Figure 10d,e).For the other two moderate earthquakes (M 5.7 and M 5.5) in Qiangtang (QT) block in March 2021, increased ILRA values were observed in late February within Qiangtang (QT) block and its adjacent blocks-the Lasa (LS) and Sichuan-Yunnan (SY) blocks (blue arrows in Figure 10b,d,e).At the beginning of 2021 Additionally, the daily area-averaged ILRA in 2021 for the five tectonic blocks in the Tibetan Plateau was computed and compared to analyze the relation between the OLR variations and earthquake activities in different tectonic blocks.The ILRA values in 2019 were also computed for comparison since no earthquake with a magnitude of greater than 5.5 occurred in the Tibetan Plateau.The ILRA was modified by subtracting the average annual ILRA value to remove possible regional variations.An increased ILRA was observed in the Bayan Har (BH) and its adjacent blocks, i.e., the Qiangtang (QT) and Qaidam (QD) tectonic blocks, about 2 weeks prior to the Madoi M 7.3 earthquake (red arrows in Figure 10a-c), but no similar behaviors were observed in the Lhasa (LS) and Sichuan-Yunnan (SY) blocks (Figure 10d,e).For the other two moderate earthquakes (M 5.7 and M 5.5) in Qiangtang (QT) block in March 2021, increased ILRA values were observed in late February within Qiangtang (QT) block and its adjacent blocks-the Lasa (LS) and Sichuan-Yunnan (SY) blocks (blue arrows in Figure 10b,d,e).At the beginning of 2021 (green bar in Figure 10), the most distinct ILRA enhancements occurred in the Lasa (LS) and Qiangtang (QT) blocks, although all five tectonic blocks showed high ILRA values in varying degrees, which could have been related to the stress accumulation in the Tibet Plateau during this period.It is worth noting that a slight increase in the IMRA in the month of January was also observed in the Bayan Har block (Figure 8a).In comparison, lower ILRA values were observed within the five tectonic blocks in 2019 when earthquakes of M > 5.5 did not occur (blue lines in Figure 10).Using MODIS monthly brightness temperature data, Yang et al. [37] observed two increased thermal anomaly processes in the Tibetan Plateau prior to the Madoi earthquake, which started from June 2020 and January 2021, and they found that the anomalies moved from west to east with time.Our results based on MWBT and OLR data with higher temporal resolutions shed further light on the thermal evolution within the different scale tectonic units (seismogenic fault and active tectonic block) during the second thermal anomaly process since January 2021.
(green bar in Figure 10), the most distinct ILRA enhancements occurred in the Lasa (LS) and Qiangtang (QT) blocks, although all five tectonic blocks showed high ILRA values in varying degrees, which could have been related to the stress accumulation in the Tibet Plateau during this period.It is worth noting that a slight increase in the IMRA in the month of January was also observed in the Bayan Har block (Figure 8a).In comparison, lower ILRA values were observed within the five tectonic blocks in 2019 when earthquakes of M > 5.5 did not occur (blue lines in Figure 10).Using MODIS monthly brightness temperature data, Yang et al. [37] observed two increased thermal anomaly processes in the Tibetan Plateau prior to the Madoi earthquake, which started from June 2020 and January 2021, and they found that the anomalies moved from west to east with time.Our results based on MWBT and OLR data with higher temporal resolutions shed further light on the thermal evolution within the different scale tectonic units (seismogenic fault and active tectonic block) during the second thermal anomaly process since January 2021.In view of the good response of the ILRA to the two moderate earthquakes that occurred in the Qiangtang block in March 2021 (Figure 10), we further produced and analyzed daily ILRA spatial distributions from January to March 2021.The data clearly showed the increased ILRA values frequently appeared in the epicentral region and its surrounding areas more than two months prior to the two moderate earthquakes in March 2021 (Figure 11).Within one month prior to the two events, e.g., 26 February, 5-7 March, and 15 March 2021, high ILRA values were repeatedly distributed along the India-Eurasia collision zone (Figure 11).Moreover, we computed the accumulated ILRA prior to the two moderate earthquakes and the M 7.3 event (Figure 12).For the M 5.5 and M 5.7 earthquakes in March 2021, the accumulative value of daily ILRA in February was computed.For the Madoi M 7.3 earthquake, the ILRA values during the period from 28 April to 8 May 2021 were used for the accumulation calculation.Prior to the two moderate earthquakes in March 2021, the high values of the accumulated ILRA were concentrated in between the epicenters of two earthquakes and along the India-Eurasia collision zone (Figure 12a).However, prior to the Madoi M 7.3 earthquake, high accumulative ILRA values were distributed in the epicentral area and its west region, which were in the Bayan Har block and its southern adjacent block (Figure 12b).Our results show a good correspondence between the high accumulative ILRA values and the seismogenic block and/or the epicenter of the impending strong earthquake.

Conclusions
The thermal variations associated with the 2021 Madoi M 7.3 earthquake in China were investigated using microwave and infrared satellite data.Microwave brightness temperature (MWBT) and outgoing longwave radiation (OLR) data with differences in spatial resolution were used to demonstrate the response to the different active tectonic units (seismogenic fault and active tectonic block) from a few months to a few days prior to the main earthquake event.An increased IMRA along the seismogenic fault since mid-February 2021 and the distinct enhancement in the ILRA within one month prior to the main event in the Bayan Har tectonic block and its adjacent blocks were observed.The occurrence of maximum IMRA enhancement gradually approached the day of the main event from large-scale (active tectonic block) to small-scale (seismogenic fault zone) tectonic units.The increased ILRA located in the southern Tibet Plateau during the beginning of 2021 may have been related to the stress accumulation in the Tibet Plateau during this period, as proven by two moderate earthquakes that occurred in March 2021 in this region.Our results show that pronounced changes in thermal signals could be an indicator of tectonic activity associated with the emissions of carbonaceous [55] and radon [57] gases from the deep crust and/or the p-holes effect [58] caused by changes in regional stress during the earthquake preparation stage.

Figure 1 .
Figure 1.Distribution map of elevation, tectonic blocks in Tibetan Plateau, and the Madoi M 7.3 earthquake on 21 May 2021.The boundaries of active tectonic blocks originate from the work of Zhang et al. [39].The 2021 Madoi earthquake was the largest one within the Bayan Har tectonic block, although a few major earthquakes have occurred along the boundary of this block over the past three decades.BH-Bayan Har block; QT-Qiangtang block; QD-Qaidam block; LS-Lhasa block; SY-Sichuan-Yunnan block.

1 Figure 2 .
Figure 2. Spatio-temporal distribution of soil moisture and NDVI.(a) Monthly area-averaged (5° 5°) variations in the soil moisture and NDVI around the epicentral region in 2021.The red triangl indicates the month in which Madoi earthquake occurred.(b) Five-month-averaged soil moistur in epicentral region from January to May 2021.(c) Five-month-averaged NDVI in epicentral regio from January to May 2021.

Figure 2 .
Figure 2. Spatio-temporal distribution of soil moisture and NDVI.(a) Monthly area-averaged (5 • × 5 • ) variations in the soil moisture and NDVI around the epicentral region in 2021.The red triangle indicates the month in which Madoi earthquake occurred.(b) Five-month-averaged soil moisture in epicentral region from January to May 2021.(c) Five-month-averaged NDVI in epicentral region from January to May 2021.

Figure 3 .
Figure 3.Time series of IMRA and soil moisture during the period of January to July 2021 in the seismogenic fault zone.(a) Percentages of IMRA with different scales; (b) temporal variation in IMRA > 2 around the occurrence of Madoi earthquake; (c) temporal variation in soil moisture.Roman numerals indicate four phases with increased IMRA prior to the main event.Pink bars indicate the variations independent of soil moisture, and yellow bars indicate the variations caused by the decrease in soil moisture.Cyan arrows indicate the behavior of the decrease in soil moisture.

Figure 3 .
Figure 3.Time series of IMRA and soil moisture during the period of January to July 2021 in the seismogenic fault zone.(a) Percentages of IMRA with different scales; (b) temporal variation in IMRA > 2 around the occurrence of Madoi earthquake; (c) temporal variation in soil moisture.Roman numerals indicate four phases with increased IMRA prior to the main event.Pink bars indicate the variations independent of soil moisture, and yellow bars indicate the variations caused by the decrease in soil moisture.Cyan arrows indicate the behavior of the decrease in soil moisture.

Figure 4 .
Figure 4. Time series of IMRA > 2 in the seismogenic fault zone (a) and soil moisture (b) during the period of January-July in 2019 and 2020.

Figure 4 .
Figure 4. Time series of IMRA > 2 in the seismogenic fault zone (a) and soil moisture (b) during the period of January-July in 2019 and 2020.

19 Figure 5 .
Figure 5. Spatial distribution of the IMRA during the three enhanced IMRA phases: (a,b) Phase I, (c-f) phase II, and (g,h) phase III.Black lines indicate the boundaries of tectonic blocks, and red indicate the main active faults.

Figure 5 .
Figure 5. Spatial distribution of the IMRA during the three enhanced IMRA phases: (a,b) Phase I, (c-f) phase II, and (g,h) phase III.Black lines indicate the boundaries of tectonic blocks, and red lines indicate the main active faults.

Figure 6 .
Figure 6.Time series of ILRA variations in different tectonic blocks from March to June in 2021.BH-Bayan Har block; QT-Qiangtang block; QD-Qaidam block; LS-Lhasa block; SY-Sichuan-Yunnan block.Colored dashed lines indicate the threshold (mean ILRA plus two standard deviations) in tectonic blocks with corresponding color.Colored dashed boxes indicate the parts exceeding the corresponding threshold.The red star indicates the day that the Madoi M 7.3 earthquake occurred, and the blue stars indicate two moderate earthquakes in Qiangtang block in March 2021.

Figure 6 .
Figure 6.Time series of ILRA variations in different tectonic blocks from March to June in 2021.BH-Bayan Har block; QT-Qiangtang block; QD-Qaidam block; LS-Lhasa block; SY-Sichuan-Yunnan block.Colored dashed lines indicate the threshold (mean ILRA plus two standard deviations) in tectonic blocks with corresponding color.Colored dashed boxes indicate the parts exceeding the corresponding threshold.The red star indicates the day that the Madoi M 7.3 earthquake occurred, and the blue stars indicate two moderate earthquakes in Qiangtang block in March 2021.

Figure 7 .
Figure 7. Spatial distribution of ILRA variation prior to the Madoi M 7.3 earthquake.Black lines indicate the boundaries of tectonic blocks, and red stars indicate the epicenter of the Madoi M 7.3 earthquake.

Figure 7 .
Figure 7. Spatial distribution of ILRA variation prior to the Madoi M 7.3 earthquake.Black lines indicate the boundaries of tectonic blocks, and red stars indicate the epicenter of the Madoi M 7.3 earthquake.

Figure 8 .
Figure 8.Time series of IMRA > 2 in the different scale tectonic units associated with the Madoi M7.3 earthquake: (a) Bayan Har tectonic block, (b) eastern Bayan Har tectonic block, and (c) seismogenic fault zone.Blue arrows indicate the maximum enhanced IMRA intensity in different tectonic units.

Figure 8 .
Figure 8.Time series of IMRA > 2 in the different scale tectonic units associated with the Madoi M 7.3 earthquake: (a) Bayan Har tectonic block, (b) eastern Bayan Har tectonic block, and (c) seismogenic fault zone.Blue arrows indicate the maximum enhanced IMRA intensity in different tectonic units.

Figure 9 .
Figure 9. Spatial-temporal variation in the ICOA associated with the Madoi M 7.3 earthquake.The red dashed line indicates the threshold line with ICOA = 3.The red circles indicate the values exceeding the threshold line.

Figure 9 .
Figure 9. Spatial-temporal variation in the ICOA associated with the Madoi M 7.3 earthquake.The red dashed line indicates the threshold line with ICOA = 3.The red circles indicate the values exceeding the threshold line.

Figure 10 .
Figure 10.Time series of ILRA variations in different tectonic blocks in 2019 and 2021: (a) Bayan Har block-BH; (b) Qiangtang block-QT; (c) Qaidam block-QD; (d) Lhasa block-LS; and (e) Sichuan-Yunnan block-SY.The red star represents the Madoi M 7.3 earthquake, and black stars indicate two moderate earthquakes that occurred in the Qiangtang tectonic block.

Figure 10 .
Figure 10.Time series of ILRA variations in different tectonic blocks in 2019 and 2021: (a) Bayan Har block-BH; (b) Qiangtang block-QT; (c) Qaidam block-QD; (d) Lhasa block-LS; and (e) Sichuan-Yunnan block-SY.The red star represents the Madoi M 7.3 earthquake, and black stars indicate two moderate earthquakes that occurred in the Qiangtang tectonic block.

19 Figure 11 .
Figure 11.Spatial distribution of ILRA variation from January to March 2021.Black lines indicate the boundaries of tectonic blocks.The red stars indicate the epicenter of the Madoi M 7.3 earthquake, and magenta dots indicate two moderate earthquakes in March 2021.

Figure 11 .
Figure 11.Spatial distribution of ILRA variation from January to March 2021.Black lines indicate the boundaries of tectonic blocks.The red stars indicate the epicenter of the Madoi M 7.3 earthquake, and magenta dots indicate two moderate earthquakes in March 2021.

Figure 12 .
Figure 12.Spatial distribution of accumulative ILRA variation prior to (a) the two moderate earthquakes in the Qiangtang tectonic block and (b) the Madoi M 7.3 earthquake.