Hydrothermal anomalies of the Earth's surface and crustal seismicity related to Ms8.0 Wenchuan EQ

Hydrothermal anomalies related to the Ms8.0 Wenchuan earthquake (EQ) on May 12, 2009, have been widely reported. However, the reported anomalies have not been associated with multi-geosphere analysis, and space–time analysis with crustal seismicity is lacking. In this paper, the space–time variation of hydrothermal parameters, including soil moisture, soil temperature, near-surface relative humidity (RHsig995) and air temperature (TMPsfc), was first extracted and analyzed with the NCEP-FNL reanalysis dataset. The b-value (a seismic parameter from the Gutenberg–Richter law) was calculated and mapped to unravel the crustal stress and rock rupture. Our results reveal a similar time window for hydrothermal anomalies on April 20 and April 30, 2008, and these anomalies are mainly distributed along the southern and middle parts of the Longmenshan fault zone. The surface temperature anomalies lag behind the humidity anomalies, and the accelerating stress accumulation started since June 2007 and lasted for eight to nine months before the mainshock. The b-value mapping shows a segmented difference along strike of the Longmenshan fault, and that regional stress accumulated mainly in the southern parts of the F2 and F3 faults. We propose the occurrence of a complex coupling process led by crustal stress buildup before the Wenchuan EQ. The anomalies are concentrated in the southern part of the surface rupture zone. The prolonged crustal stress accumulation corresponds to the short intermittent hydrothermal response on the Earth’s surface before the Wenchuan EQ. Our findings reveal new hydrothermal anomalies in the Earth’s surface and atmosphere and explore direct link with seismogenic processes in the crust.


Introduction
Monitoring and early warning of tectonic earthquake (EQ) has become a research hot spot in seismology. With the advancement of remote sensing technology, a variety of precursory anomalies related to strong EQ were reported based on integrated analysis of ground and satellite observations. Meanwhile, various hypotheses or models, such as the lithosphere-atmosphere-ionosphere coupling (LAIC) model and/or the lithosphere-coversphere-atmosphere coupling (LCAI) model and positive hole-type charge carriers in the stressed rock (p-hole), have been proposed to reveal these anomalies (Pulinets and Ozounov 2011;Freund 2010Freund , 2011Wu et al. 2012). A catastrophic Ms8.0 EQ had hit Wenchuan country (Sichuan province) in southwestern China, at 06:28 (UTC) on May 12, 2008. The epicenter (31°N, 103.4°E) is located at the middle part of the Longmenshan fold belt (eastern Tibetan Plateau) with the hypocenter depth of about 19 km (Fig. 1a). This earthquake was one of the worst in modern Chinese history, causing over 80,000 fatalities and severe damages to infrastructure (Fu et al. 2011a, b). It is noteworthy that this EQ happened on the active fault boundary around the Bayan Har fault block in the Tibetan plateau (Fig. 1b). With the northward intrusion of India plate into the Eurasian interior, the Bayan Har fault block is generally considered to have  (Tapponnier et al. 2001;Zhang et al. 2010).
Thermal anomalies observed in the tectonic fault intersecting zone or around the epicenter before the Wenchuan mainshock include surface temperature, surface latent heat flux, thermal infrared radiation and outgoing long-wave radiation Jing et al. 2013;Li et al. 2010;Singh et al. 2010;Wu et al. 2012;Xiong et al. 2009;Yang and Mi 2009;Liu et al. 2020). Meanwhile, short-term ionospheric disturbances, such as the total electron content (TEC), critical frequency and maximum electron density of ionospheric F2 layer (foF2) and (NmF2), were widely considered to be possibly related to the Wenchuan EQ (Ding et al. 2010;Hsiao et al. 2010;Jhuang et al. 2010;Liu et al. 2009;Xiong et al. 2008;Xu et al. 2010;Zhao et al. 2010;Zhu et al. 2008;Zhang et al. 2009;Zhou et al. 2008). However, the space-time correlation between these surficial/atmospheric anomalies and seismicity in the crust during the earthquake preparatory phase is still not clear, under the premise of stress increase and energy accumulation release for triggering earthquake. This study aims to identify the hydrothermal anomaly responses to the Wenchuan EQ and analyze the space-time correlation with seismicity to understand the geosphere coupling mechanism of earthquake precursors.

Reported anomalies of Wenchuan EQ
After the Wenchuan earthquake, many reports of pre-earthquake heat radiation, gas, temperature and electromagnetic anomalies were generated from satellite remote sensing, assimilation dataset and ground monitoring stations, with these many parameters belonging to different Earth layers. For the study of the Earth's surficial and atmospheric parameters, Singh et al. (2010) used the data of sensor brightness temperature aboard the polar orbit satellite, temperature, skin temperature and relative humidity from the Chongqing meteorological station and found that the epicenter's air temperature in each horizontal pressure gradient has gradually risen on May 8. The skin temperature rose with the air temperature on May 6 and peaked on May 1. The relative humidity of the epicenter dropped significantly on May 8. Meanwhile, the vertical relative humidity profile also shows that the lower atmosphere (700 hpa) dropped significantly and the upper atmosphere (400 hPa) rose significantly on May 5. Jing et al. (2013) used OLR accumulative daily and monthly data on temperature, air pressure and relative humidity from the polar-orbiting satellite for anomaly identification. The results show that OLR presented an anomaly in the wider areas in the Sichuan Basin from April 29 to April 30, and it turned into an anomalous stripshaped area concentrating on the Longmenshan fault zone on May 5. On the day before the earthquake (May 11), SLHF presented a slightly anomalous area in the northern part of the Longmenshan fault zone. On the day of the earthquake, it turned into a cluster of strong anomalies south of the epicenter. Yang and Mi (2009) used the surface long-wave radiant flux, temperature data of underground 10-20 cm and the AT2m parameter to identify, respectively, anomalies of the surface long-wave radiant flux on March 10-15 and May 6-8, and the suspected anomaly of the underground 10-20 cm temperature on the March 14 and AT2m on March 15, April 28 and May 6. Wu et al. (2012), based on assimilation and remote sensing data, identified five parameters of the multiperiodic quasi-synchronous anomaly from April 26 to May 11. Latest studies also indicated a chain of precursor changes in four parameters (aerosol optical depth (AOD); skin temperature (SKT); surface latent heat flux (SLHF); and total column water vapor (TCWV)) within two months before the Wenchuan EQ (Liu et al. 2020). The reported ionospheric anomalies before the Wenchuan EQ were unprecedented. The ionospheric parameters show an anomaly from April 29 to the day of the earthquake and the maximum anomaly on May 9 (3 days before the earthquake), and the main spatial features of the various ionospheric anomalies are located in the vicinity and south of the epicenter (Ding et al. 2010;Hsiao et al. 2010;Jhuang et al. 2010;Liu et al. 2009;Xiong et al. 2008;Xu et al. 2010;Zhao et al. 2010;Zhu et al. 2008;Zhang et al. 2009;Zhou et al. 2008).
However, published studies on the Wenchuan EQ lack hydrothermal anomalies of the Earth's surface (i.e., lithospheric-atmospheric transition layer), which pose difficulties in revealing the coupling process of the geosphere anomalies. This is supplemented by our new multi-parametric study on hydrothermal anomalies.

Data and method
Four parameters characterizing the hydrothermal conditions, i.e., volumetric soil moisture level 1 (SML1) at 0-7 cm b.g.l. (below ground level), soil temperature level 1 (STL1) at 0-7 cm b.g.l. from ERA-Interim reanalysis dataset, near-surface relative humidity (RHsig995) and air temperature (TMPsfc) from NCEP-FNL reanalysis dataset, were selected operationally for every 6 h/day in long-term intervals and within a month before the mainshock for anomaly identification. ERA-Interim reanalysis dataset is a new global atmosphere reanalysis product of the European Centre for Medium-Range Weather Forecasts (ECMWF), which covers the data-rich period from 1979 to now. The gridded data were transformed into a regular N128 Gaussian grid with 0.71° × 0.71° spatial resolution (https ://apps.ecmwf .int/datas ets/data). NCEP-FNL reanalysis dataset is produced with the same National Center for Environmental Prediction (NCEP) model as the one used for the Global Forecast System (https ://rda.ucar.edu/datas ets) and is also represented by a Gaussian grid with 1° × 1° spatial resolution.
In the time-series analysis of the epicenter pixel (black rectangles in Fig. 3), the historical data (year 2000-2008) since March 1 before the earthquake were extracted to calculate the mean (μ) and standard deviation (σ) as background year of 2000-2007. Here, we defined a deviation with over-quantity > 1.5 σ as a threshold and thus as a potential anomaly. Consider the following weather system conditions: (1) Compared with noon time, the sunrise time and sunset time in April and May were greatly different due to the rugged terrain in Wenchuan. When the sun is at low altitude, it is greatly affected by the terrain and the net radiation flux is unstable; (2) the southern trough hardly affects the Wenchuan area, which is blocked by the Hengduan Mountains and the Tibetan Plateau. Thus, the region experienced less from weather system, and local net radiation and moisture conditions directly influence the surface meteorological elements in most of the day. Hence, all the ERA-Interim and NCEP-FNL data at 06:00 UTC (14:00 LT) were selected for information extraction and anomaly recognition. The daily average and the maximum/minimum values based on the ground-based station data were then analyzed.
Considering the possible uncertainties arisen from different seasons, terrains and latitudes, we computed the differential images of the changed parameters (ΔI) by subtracting the 2008 daily value from the multiple annual mean (2000)(2001)(2002)(2003)(2004)(2005)(2006)(2007) to determine the spatial anomalous distributions. The other differential images between the 2007 (background year) and the multi-annual history (2000)(2001)(2002)(2003)(2004)(2005)(2006) were also compared with the normal background. The image calculation is expressed as follows: where I t = daily value of a parameter in 2008 and μ t = corresponding daily mean estimated over the years 2000-2007. ΔI t images on the same day in 2007 were applied for comparison, and I i was adjusted to the means of 2000-2006. In addition, a collection of the daily average temperature, relative humidity and precipitation data from Dujiangyan meteorological monitoring station in Sichuan province was adopted to supplement the time-series analysis (yellow dot in Fig. 1). The station is adjacent to the Wenchuan epicenter with elevation of 698.5 m.

Multi-parameter anomaly
The anomalous temperature and humidity of the soil do not show quasi-synchronization in time series of the four parameters of the epicenter pixels ( Fig. 2a, b). The anomalous time of SML1 was mainly on April 17-22, and the maximum appeared on April 20 (0.42 m 3 / m 3 ). While the STL1 anomaly obviously lagged behind SML1, the abnormal intensity was weak. It is only slightly greater than the μ + 1.5σ background critical threshold of 0.06 K on April 30 (Fig. 2b). It shows that the soil properties of the epicenter pixels had changed anomalously before the earthquake, and the humidity change was more severe and earlier than that of temperature.
Besides, temperature and humidity curves of the Earth's surface fluctuated significantly, and the anomalous drop and rise of RHsig995 and TMPsfc had occurred many times (Fig. 2c,d). Among these, the strongest negative anomalies of RHsig995 occur on March 26 and April 30, which are below the critical threshold of μ − 2σ background by over 4%, while the values of most other time were basically above 50%. It was in contrary to the humid climatic conditions of perennial April and May in Sichuan. Except for April 23, the TMPsfc value of the month before the 2008 earthquake exceeded the historical average (black line in Fig. 2d). Therefore, the anomalous surface temperature rise appeared most frequently and (mostly) exceeded the upper limit of the μ + 2σ background threshold and reached a peak on May 6. However, the maximum abnormality occurred on May 5, exceeding the μ + 2σ background threshold of 3.5 K. In fact, the overall value of the epicenter pixel started its anomalous change from March 16 to May 7. It shows that the anomalous surface temperature rise in the epicenter area was significant, and it occurred frequently and lasted longer before the earthquake. It is thus indicated the perturbation of the surface temperature and humidity after the anomalous changes of the soil layer. Very recently, Liu et al. (2020) used the CAPRI algorithm to detect anomalies in the time series of SKT, SLHF and TCWV from the ERA-Interim datasets before the Wenchuan EQ. Their results show that the earlier precursor's time (mainly from April 7 to April 14) is also broadly coeval with SML1, STL1 and TMPsfc.
According to the spatial proximity criterion, the different images with typical anomalous spatial distribution features were selected for analysis after preliminary visual screening. The area ΔSML1 with an anomalous rise was mainly inside the Sichuan Basin (Fig. 3a). The anomalous area on April 20 was concentrated in the southern and central parts of the Longmenshan fault zone and was consistent with the NNE trend of the Longmenshan fault zone.
(1) 1 3 The main difference of the anomaly field exceeds 0.11 m 3 /m 3 . The spatial distribution in the same period in 2007 shows consistent low values (Fig. 3b), which supports the significant spatial correlation between the anomalous characteristics of the day and the seismogenic structure. The ΔSML1 spatial characteristics indicate that the topmost soil moisture on April 20 suddenly increased before the earthquake in the Wenchuan seismogenic fault zone, but did not spread to the northern part of the fault. Compared with ΔSML1, ΔSTL1 did not exhibit strong anomalous spatial characteristics for many days. The anomalous area on April 30 was mainly located in the epicenter and to the south of the area (Fig. 3c). The ΔSTL1 of the entire Longmenshan fault zone had spread to a slight strip (2.5-3.5 K). The same period in 2007 also shows a "peaceful" image, indicating that the soil temperature along the Longmenshan fault zone had risen anomalously before the mainshock (Fig. 3d).
Anomalous changes in soil properties would continue to affect the (near-) surface atmosphere. ΔRHsig995 of the Longmenshan fault zone suddenly dropped on April 28, and its anomalous spatial morphology was completely consistent with the F1, F2 and F3 locations (also located in the future surface rupture zone) (Fig. 3e). Compared with the surrounding areas, the anomaly of ΔRHsig995 on April 28 was spatially independent, while the same period in 2007 shows substantial increase in the eastern margin of the Qinghai-Tibetan Plateau. The negative anomalous field of ΔRHsig995 was not seen either in the same period of 2007 (Fig. 3f). It is shown that after the soil moisture change (April 20), the relative humidity of the upper fault-zone surface had also changed. The ΔTMPsfc spatial anomaly was mainly on April 30 and was concentrated to the west of the epicenter. This shows a spatial pattern consistent with ΔRHsig995 for the same period (Fig. 3g). This spatial anomalous warming feature is absent in the same period of 2007 (Fig. 3h).
Dujiangyan city is located on the F3 fault in the southeastern part of the Wenchuan EQ epicenter (Fig. 1). Surface rupture occurred after the mainshock. Here, we focused on the changes of Dujiangyan ground-based meteorological data in the period of the above hydrothermal parameter anomalies as the background of historic data from 1998 to 2007. It shows that the overall surface temperature change at the Earth's surface was similar to the TMPsfc parameters and peaked on April 30 and May 11 above threshold of μ + σ background (Fig. 4a). In particular, the surface temperature was up to 28.5 °C on April 30, while there is weak anomalous change of air temperature but not remarkable (< μ + σ). Meanwhile, the relative humidity drops on April 28-30 were quasi-synchronous with RHsig995 of NCEP-FNL data (Fig. 4b). The minimum value on April 28 reached 49%, which was higher than the RHsig995 epicentral pixel value (30%), but was largely the same as the abnormal change time. In addition, the Dujiangyan meteorological monitoring station has recorded almost no atmospheric precipitation on April 20 and April 30 (Fig. 4c), further supporting that the soil temperature and humidity changes on both days were unrelated to meteorological conditions. Compared with the anomalous changes in soil humidity, the soil temperature had responded to the Wenchuan EQ later but at a higher intensity. In short, data from the Dujiangyan meteorological station verify the reliability of the assimilation parameter results and confirmed that the main hydrothermal anomaly changes on ground surface probably occurred on April 20 and April 30 before the Wenchuan EQ. 1 3 Fig. 4 Time series of air temperature (a), relative humidity (b) and precipitation (c) during the month before the mainshock and historic records at the Dujiangyan meteorological stations

Seismicity of b-value
The b-value is a key coefficient to express the fundamental relationship between the magnitude and frequency of earthquakes in the Gutenberg-Richter law (Gutenberg and Richter 1944). The b-value is widely used to study the size distributions of small-/medium-magnitude seismic events and associated with several geophysical properties, including regional stress, material homogeneity and temperature gradient (Gulia and Wiemer 2010;Wiemer 2004, 2005;Tormann et al. 2015;Urbancic et al. 1992;Wiemer and Wyss 2002;Wyss and Wiemer 2000). Therefore, the b-value can serve as a proxy of crustal stress conditions and seismicity. Wang et al. (2012) analyzed 19 strong earthquakes (M6.5 or above) since 1981 in the Sichuan-Yunnan region and analyzed the b-value spatial distributions between every two strong earthquakes (i.e., the preearthquake phase of the strongest earthquake) (Wang et al. 2012). Although the relationship between the b-value spatial distribution and strong earthquake activity is unobvious (many strong earthquakes do not occur in low b-value areas), there were low-value zones in the southern part of the Longmenshan fault zone. Besides, all b-value changes with time before strong earthquakes and shows a congruous feature of systematic degradation. Tian (2013) selected the earthquake history from 1970 to 2009 before the Wenchuan EQ and analyzed the b-value spatial distribution of the Longmenshan fault zone. It is found that the Wenchuan-Maoxian and the Beichuan-Pingwu sections of the Longmenshan fault zone show an anomalous low b-value before the mainshock, which had returned to normal after the earthquake. Due to the different choices of the space-time range of the microseismic catalogues, existing studies often yielded different results.
In this paper, the seismic subsets used for b-value calculation were strictly selected according to the statistical characteristics of the microseismic catalogue. The calculation method based on the least squares was used to analyze the temporal variation and map spatial pattern of the b-value before and after the Wenchuan EQ in the Longmenshan fault zone. We further unraveled the anomalous variation features of the environment (e.g., lithospheric tectonic stress and rock rupture) to better understand the lithosphere-Earth's surface, even up to atmosphere, anomaly coupling process.
To ensure the seismic completeness of each subset, the minimum completeness magnitude (Mc) is firstly calculated by using the best fitting method with a coverage over the whole region for each subset. Due to the strong aftershocks after the Wenchuan EQ, the Mc at P2 was as high as 3.0, and the average b-value also dropped from normal 0.98 ± 0.03 at P1-2 (before the EQ) to 0.77 ± 0.15 (Fig. 5b, c). Accordingly, any later records of Mc are discarded from the two subsets for the subsequent b-value calculation and mapping. In the b-value space mapping, the grid spatial resolution of 0.1° × 0.1° and the scanning radius of 40 km are set, and the minimum number of events for calculating the b-value is 20. Meanwhile, Mc is determined based on the magnitude-cumulative frequency relationship of each subset. Once the Mc is determined, all the subsets would be selected, and then, the MATLAB software package ZMAP (Wiemer 2001)

b-value anomaly
Magnitudes changes with time of P1-2 and P2 phases in the study area reveal high frequency seismicity before and after the Wenchuan EQ (Fig. 6a, b). The strong aftershocks lasted longer and did not subside until July 2008 (Fig. 6b). In spatial distribution during the P1-2 phase, the epicenter was widely distributed in the high-altitude area of the Tibetan Plateau on the western side of the Longmenshan fault zone, and a few scattered in the southeastern corner of the Sichuan Basin (Fig. 6c). The seismic intensity was dominated by Ml < 3.0 events, and only one Ml > 5.0 event occurred in Gansu province on June 20, 2006 (Ml = 5.3). In/around the Longmenshan fault zone, seismic events were concentrated between F2 and F3 faults and were mostly Ml < 4.0 (Fig. 6c). The P2 phase was dominated by strong (Ml > 5.0) NE-striking aftershocks were concentrated in the middle and northern parts of the Longmenshan fault zone (Fig. 6d). In the P1-2  the P1-2 phase are slightly above the average of the study area (Fig. 6e). In the P2 phase, the b-values are distinctly different between the southern, middle and northern parts of the Longmenshan fault zone. Due to the surface rupture at the epicenter and along the F2 and F3 faults, the stress was fully released and the b-values are significantly higher than those of the northern part. Due to the surface rupture development and the stress release, the b-values around the mainshock have returned to the normal level (~ 1.0), and the stress accumulation with higher b-values had transferred to the unruptured area north of the Longmenshan fault zone (Fig. 6f). Such spatial characteristics of b-values indicate certain stress accumulation enhancement in the southern Longmenshan fault zone before the mainshock, while b-values of the foreshocks and the Wenchuan main rupture zone are relatively high, and the stress accumulation is unobvious. Based on the known locations of the surface rupture zone after the earthquake (Fu et al. 2011a, b), a subset of the buffer zone (30 km radius, from P1-2 to P2 phase) of the post-earthquake F2 and F3 fault surface rupture zones was extracted for the time-series b-value analysis, and the number of recorded events (n) is 5412 and 2548, respectively. A subset of the southern F2 fault (where no surface rupture occurred) at a 30 km radius in the same phase was also extracted (n = 1290) for comparison. The time-cumulative frequency relationship of the three buffer subsets shows a substantially consistent preand post-earthquake segmentation (Fig. 7). According to the magnitude-accumulation frequency relationship, the corresponding complete subsets are filtered out by Mc = 1.3, 1.4 and 1.0. The time-series calculation of the b-value also adopts the optimal fitting method and sets 200 and 300 event-sliding windows and ten time-superposition windows.
The b-value dropped systematically in three buffer subsets and experienced three stages of "stationary-decline-recovery" (Fig. 7). However, the drop peak of b-value in each subset does not coincide with the time when the minimum was reached. For example, the main Wenchuan EQ rupture zone along the F2 fault has encountered rapid drop in b-value since June 2007 and reached the minimum under calculations of 200 and 300 event-sliding windows during February 20-March 20, 2008, i.e., the stress accumulation was up to nine months. The b-values then rebounded rapidly, and the Ms8.0 Wenchuan EQ occurred in the middle of the rebound (Fig. 7a). The similar scenario is present in the other rupture zone along the F3 fault, but the takeoff time was later (since July 2007) and the stress acceleration was shorter when b-value reaches the minimum (over eight months) (Fig. 7c). Mean of the b-value time-series distribution of the two rupture zones is 1 to 1.1 prior to the Wenchuan EQ, which may indicate stronger seismicity of future rupture zones than where there is no surface rupture in the southern F2 fault (mean < 1.02) (Fig. 7b). Therefore, we infer that stress accumulation in the crust lasted for at least 8-9 months and the release accelerated probably in February-March 2008 before the mainshock.
To further compare the above seismicity changes with the b-value, the other parameters reported by Huang (2008) are cited to show the spatiotemporal characteristics in Fig. 8, including RTL (region-time-length, which describes the decrease/increase in the seismicity measured from the epicentral distance), time and rupture length compared to the background level and Q-value (average RTL of same time window). Figure 8a shows a seismic gap of RTL parameter during 2006-2007. The Ms8.0 Wenchuan EQ occurred after RTL returned to the background level (~ 0). Meanwhile, the Q-map during January-June 2007 shows a clear seismic quiescence around the mainshock epicenter (Fig. 8b). Obviously, the "decrease-rebound" time trend of RTL from 2006 to 2008 is largely consistent with our b-value time series (Fig. 7), and the negative Q-value region is also distributed along the Longmenshan fault zone, similar to the low b-value area in Fig. 6a.

Anomaly response of the Earth's surface
Based on the above spatiotemporal evolution analysis of hydrothermal parameters at the Earth's surface, there is a similar time window for multi-parameter anomalies of before/ after April 20 and April 30 (22-12 days before the Wenchuan EQ). Temperature anomalies of the upper soil and the ground surface lag behind the humidity anomalies. Moreover, anomalies of each parameter have also multi-phase variation on April 20 and/or April 30. Spatial distribution of the surface hydrothermal anomaly is mainly along the southern and middle parts of the Longmenshan fault zone and is NE striking (ΔSML1, ΔRHsig995 and ΔTMPsfc). Our results indicate that the hydrothermal transfer process and changes from the lithosphere to the Earth's surface (even to the atmosphere) are continuous and stable. Due to the more open-space environment on the Earth's surface, surface temperature anomalies occurred in successive periods until the mainshock. This agrees with the published results of thermal parameter anomalies (e.g., SLHF and AT2m), supporting that the atmospheric thermal anomalies were strongly fluctuated before the Wenchuan EQ. Moreover, broadly coeval anomalies of more parameters (including SML1, STL1 and TMPsfc) we analyzed during earlier April can verify the physical changes related to the Wenchuan EQ on the Earth's surface. In addition, it indicates various anomalies from multi-source observation occurred a month before the Ms8.0 mainshock.

Anomaly of crustal rupture
During the P1-2 before the mainshock, the overall b-values in the study area are relatively high, but show segmented differences along strike in the Longmenshan fault zone. F2 fault corresponds to the high b-value (> 1.3) area in the middle-northern parts, while F3 fault  ( modified from Huang 2008) corresponds to the relatively low b-value (0.85-0.92) area in the southern parts of the two faults. The b-values in the middle-southern Longmenshan fault zone are relatively high and returned to normal level in the P2 phase, due to the surface rupture development and substantial stress-energy release. It is noteworthy that major b-value anomalies are present in the northern section, which indicates that the post-earthquake stress accumulated and transferred to the other area (without surface rupture) of the Longmenshan fault zone. The b-value time series of three buffer subsets of surface rupture zone show the "b-value drop" phenomenon. However, compared with the southern F2 fault (without surface rupture), the timing of abrupt b-value drop in the two surface rupture zones was close to the mainshock, with shorter duration (8-9 months) and more drastic fluctuation. This indicates that the crustal stress environment below the rupture surface was more variable in the P1-2 phase.
The accelerating stress accumulation started since June 2007, which is clearly earlier and longer than the two surface rupture zones. This indicates significant segmented difference between the stress accumulation and lithospheric rupture development in the Longmenshan fault zone. The sudden b-value change shortly before the earthquake is more sensitive to the surface rupture zone, i.e., crustal stress changed more drastically in the P1-2 phase below the surface rupture zone. In addition, the time-series results of RTL, as a synthetic seismic activity parameter (Huang 2008), also show that the "placid" (RTL < 0) anomaly phase of seismic activity is also present in 2006-2007, which agrees with our calculation results. This further supports the anomalous significance of seismic activity in the Longmenshan area in 2006-2007.

Coupling effects of the geosphere
The discussion above demonstrates that the anomalies before the Wenchuan EQ experienced a geosphere's coupling process. Long before the earthquake, the material (magma) from deep mantle ascended and passed through the equilibrium rebound area under the Longmenshan fault zone and entered the crust. The crustal fracture underwent stress accumulation and fault activity driven by plate compression and deep geodynamics. In particular, in some areas in the southern F2 and F3 faults, the crustal stress accumulation was gradually intensified. As a result, the underground gas was extruded out from rocks with the thermal energy and was transported up to the Earth's surface along faults and fractures. Meanwhile, the rocks were strongly deformed due to compression, resulting in the crystal lattice deformation. The activated p-hole propagated widely in the rocks and eventually reached the surface.
On the Earth's surface, the moisture content and soil layer temperature were all changed by the gas from the lithosphere. In addition to the water-level change caused by the stress itself, the surrounding groundwater would also participate in the changes of soil properties in the form of fluid. After the deep gas is thoroughly mixed with the soil gas, other tracer gases (notably Rn) were transported to the surface. Then, a series of coupling effects progressed to result in the multi-parametric anomalies. The Rn decay and the ionization, ion hydration and local greenhouse gas effects (caused by p-hole) may still be fundamental effects. In addition, the discovery of OLR anomalies provides direct support for the coupling of thermal radiation energy changes at upper atmosphere. The local additional electric field would eventually cause ionospheric disturbances with the global electromagnetic cycle. The b-value of the lithosphere/hydrothermal anomaly analysis of the Earth's surface can provide a solid basis for the geosphere anomalous coupling mechanism.
In conclusion, this study reveals new hydrothermal anomalies in the Earth's surface and atmosphere and explores direct link with crustal seismogenic processes. Our results have explained the quasi-synchronism of the appearance of hydrothermal parametric and connected to seismic anomalies.