Shoubiao Zhu
Current Organizations / Institutes
The organizations/institutes from the scholar's latest article
Research Subjects
Shoubiao Zhu has a total of 40 co-authors. Here are the top co-authors by number of shared publications:
Miao Miao
Department of Earth and Space Sciences Southern University of Science and Technology Shenzhen China
4 shared publicationsJie Yuan
National Institute of Natural Hazards Ministry of Emergency Management of China Beijing China
2 shared publications- Jie Yuan
School of Resources and Environment, Xingtai University, Xingtai, Hebei, China
2 shared publications Yao-Lin Shi
Key Lab of Calculation Geodynamics; Graduate School of Chinese Academy of Sciences; Beijing; 100039; China
2 shared publicationsYaoLin Shi
Chinese Academy of Sciences
2 shared publicationsYaolin Shi
2 shared publications
Latest Publications
- Open Access
- 18 May 2026
- Published by MDPI AG in Remote Sensing
- Vol. 18 (10) , 1618
Precise Contemporary Crustal Strain and Rotation Rates Derived from GNSS Measurements in the Pamir–Tian Shan Region
The Pamir–Tian Shan domain constitutes one of the most actively deforming intracontinental orogenic systems associated with continued India–Eurasia convergence. Characterizing present-day deformation in this region is fundamental to deciphering its geodynamic evolution and assessing seismic risk. Existing strain rate models based on GNSS measurements display noticeable discrepancies, largely attributable to variations in analytical strategies and uneven station distribution. In this study, we determine the present crustal strain and rotation fields across the Pamir–Tian Shan area using the most updated GNSS velocity solution referenced to stable Eurasia. To address the issues of inconsistent strain rate field results and lack of reliability verification in previous studies based on GNSS data, this paper computes the crustal strain rate field (principal strain rate, maximum shear strain rate, dilatation strain rate, and rotational strain rate) with a grid spacing of 0.75° × 0.75° in the study area, followed by numerical validation of the results’ reliability. The derived strain field is characterized by dominant NNW–SSE shortening throughout much of the orogenic system, with peak compressional strain rates (~1.0 × 10−7 yr−1) concentrated along the Pamir Frontal Thrust. By contrast, the interior of the Pamir Plateau exhibits clear EW extension, consistent with areas affected by normal-faulting earthquakes. High values of shear strain rates are primarily localized along major active fault systems, whereas negative dilatational components indicate overall contraction within the Tian Shan. The rotation-rate distribution reveals clockwise rotation of the Tarim Basin (approximately 0.6°/Myr) together with counterclockwise rotation affecting the Pamir and Tian Shan blocks, accommodated by prominent strike–slip fault networks. The close spatial agreement between the modeled strain patterns, active tectonic structures, and focal mechanism solutions supports the reliability of the inferred deformation field. The research results of this paper are of great scientific significance for in-depth study of the tectonic evolution and earthquake disaster assessment in the Pamir–Tian Shan region.
- 18 May 2026
- Published by MDPI AG in Remote Sensing
- Vol. 18 (10) , 1618
Open AccessPrecise Contemporary Crustal Strain and Rotation Rates Derived from GNSS Measurements in the Pamir–Tian Shan Region
The Pamir–Tian Shan domain constitutes one of the most actively deforming intracontinental orogenic systems associated with continued India–Eurasia convergence. Characterizing present-day deformation in this region is fundamental to deciphering its geodynamic evolution and assessing seismic risk. Existing strain rate models based on GNSS measurements display noticeable discrepancies, largely attributable to variations in analytical strategies and uneven station distribution. In this study, we determine the present crustal strain and rotation fields across the Pamir–Tian Shan area using the most updated GNSS velocity solution referenced to stable Eurasia. To address the issues of inconsistent strain rate field results and lack of reliability verification in previous studies based on GNSS data, this paper computes the crustal strain rate field (principal strain rate, maximum shear strain rate, dilatation strain rate, and rotational strain rate) with a grid spacing of 0.75° × 0.75° in the study area, followed by numerical validation of the results’ reliability. The derived strain field is characterized by dominant NNW–SSE shortening throughout much of the orogenic system, with peak compressional strain rates (~1.0 × 10−7 yr−1) concentrated along the Pamir Frontal Thrust. By contrast, the interior of the Pamir Plateau exhibits clear EW extension, consistent with areas affected by normal-faulting earthquakes. High values of shear strain rates are primarily localized along major active fault systems, whereas negative dilatational components indicate overall contraction within the Tian Shan. The rotation-rate distribution reveals clockwise rotation of the Tarim Basin (approximately 0.6°/Myr) together with counterclockwise rotation affecting the Pamir and Tian Shan blocks, accommodated by prominent strike–slip fault networks. The close spatial agreement between the modeled strain patterns, active tectonic structures, and focal mechanism solutions supports the reliability of the inferred deformation field. The research results of this paper are of great scientific significance for in-depth study of the tectonic evolution and earthquake disaster assessment in the Pamir–Tian Shan region.
- 21 October 2024
- Published by Elsevier in Journal of Asian Earth Sciences
- Vol. 277
Contemporary crustal strain rates derived from GPS measurements in southwestern China
- 21 October 2024
- Published by Elsevier in Journal of Asian Earth Sciences
- Vol. 277
Contemporary crustal strain rates derived from GPS measurements in southwestern China
- Open Access
- 7 February 2022
- Published by MDPI AG in Remote Sensing
- Vol. 14 (3) , 779
Strain-Rates from GPS Measurements in the Ordos Block, China: Implications for Geodynamics and Seismic Hazards
A number of devastating earthquakes have occurred around the Ordos Block in recent history. For the purpose of studying where the next major event will occur surrounding the Ordos Block, much work has been done, particularly in the investigation of the Earth’s surface strain rates based on GPS measurements. However, there exist striking differences between the results from different authors although they used almost the same GPS data. Therefore, we validated the method for the calculation of GPS strain rates developed by Zhu et al. (2005, 2006) and found that the method is feasible and has high precision. With this approach and the updated GPS data, we calculated the strain rates in the region around the Ordos Block. The computed results show that the total strain rates in the interior of the Block are very small, and the high values are mainly concentrated on the peripheral zones of the Ordos Block and along the large-scale active faults, such as the Haiyuan fault, which are closely aligned to the results by geological and geophysical observations. Additionally, the strain rate results demonstrated that all rifted grabens on the margin of the Ordos Block exhibit extensional deformation. Finally, based on the strain rate, seismicity, and tectonic structures, we present some areas of high earthquake risk surrounding the Ordos Block in the future, which are located on the westernmost of the Weihe Graben, both the east and westernmost of the Hetao Graben, and in the middle of the Shanxi Graben. Hence, this work is significant in contributing to a better understanding of the geodynamics and seismic hazard assessment.
- 7 February 2022
- Published by MDPI AG in Remote Sensing
- Vol. 14 (3) , 779
Open AccessStrain-Rates from GPS Measurements in the Ordos Block, China: Implications for Geodynamics and Seismic Hazards
A number of devastating earthquakes have occurred around the Ordos Block in recent history. For the purpose of studying where the next major event will occur surrounding the Ordos Block, much work has been done, particularly in the investigation of the Earth’s surface strain rates based on GPS measurements. However, there exist striking differences between the results from different authors although they used almost the same GPS data. Therefore, we validated the method for the calculation of GPS strain rates developed by Zhu et al. (2005, 2006) and found that the method is feasible and has high precision. With this approach and the updated GPS data, we calculated the strain rates in the region around the Ordos Block. The computed results show that the total strain rates in the interior of the Block are very small, and the high values are mainly concentrated on the peripheral zones of the Ordos Block and along the large-scale active faults, such as the Haiyuan fault, which are closely aligned to the results by geological and geophysical observations. Additionally, the strain rate results demonstrated that all rifted grabens on the margin of the Ordos Block exhibit extensional deformation. Finally, based on the strain rate, seismicity, and tectonic structures, we present some areas of high earthquake risk surrounding the Ordos Block in the future, which are located on the westernmost of the Weihe Graben, both the east and westernmost of the Hetao Graben, and in the middle of the Shanxi Graben. Hence, this work is significant in contributing to a better understanding of the geodynamics and seismic hazard assessment.
- Open Access
- 29 November 2021
- Published by American Geophysical Union (AGU) in Journal of Geophysical Research: Solid Earth
- Vol. 126 (12)
Dynamic Rupture Simulations of the 2008 7.9 Wenchuan Earthquake: Implication for Heterogeneous Initial Stress and Complex Multifault Geometry
The 2008 Mw7.9 Wenchuan earthquake ruptured northwest‐dipping imbricate oblique reverse faults along the Longmenshan thrust belt at the eastern margin of the Tibetan Plateau, and developed one of the most complex coseismic rupture patterns for reverse faulting events ever reported in intraplate settings. To evaluate to what extent complex fault geometry and heterogeneities in initial stress can explain the characteristics of coseismic observations, we simulate spontaneous dynamic rupture propagation governed by slip‐weakening friction law on the geometrically complex multifault system loaded by heterogeneous stress. Our model reproduces many observed features including the multifault ruptures, locations of the maximum slip, Interferometric Synthetic Aperture Radar observations, source time functions and peak ground velocities. Our results suggest that the maximum horizontal principal stress in the northern part of Beichuan Fault (BFC) is rotated approximately counterclockwise compared to that in the southern part of BFC, yet this rotation is not sufficient to produce a supershear rupture. We infer that the fault core zone in the southwestern Beichuan Fault (BCF) may have been more severely damaged than the northeastern counterpart, which provides an explanation that the rupture in the northern BCF propagates faster than in the southern portion. Besides, we find that Coulomb failure stress changes on Wenchuan‐Maoxian Fault (WMF) may be counteracted by the aftershock slips of the Lixian Fault, which suggests Wenchuan earthquake could not significantly push the WMF closer to failure.
- 29 November 2021
- Published by American Geophysical Union (AGU) in Journal of Geophysical Research: Solid Earth
- Vol. 126 (12)
Open AccessDynamic Rupture Simulations of the 2008 7.9 Wenchuan Earthquake: Implication for Heterogeneous Initial Stress and Complex Multifault Geometry
The 2008 Mw7.9 Wenchuan earthquake ruptured northwest‐dipping imbricate oblique reverse faults along the Longmenshan thrust belt at the eastern margin of the Tibetan Plateau, and developed one of the most complex coseismic rupture patterns for reverse faulting events ever reported in intraplate settings. To evaluate to what extent complex fault geometry and heterogeneities in initial stress can explain the characteristics of coseismic observations, we simulate spontaneous dynamic rupture propagation governed by slip‐weakening friction law on the geometrically complex multifault system loaded by heterogeneous stress. Our model reproduces many observed features including the multifault ruptures, locations of the maximum slip, Interferometric Synthetic Aperture Radar observations, source time functions and peak ground velocities. Our results suggest that the maximum horizontal principal stress in the northern part of Beichuan Fault (BFC) is rotated approximately counterclockwise compared to that in the southern part of BFC, yet this rotation is not sufficient to produce a supershear rupture. We infer that the fault core zone in the southwestern Beichuan Fault (BCF) may have been more severely damaged than the northeastern counterpart, which provides an explanation that the rupture in the northern BCF propagates faster than in the southern portion. Besides, we find that Coulomb failure stress changes on Wenchuan‐Maoxian Fault (WMF) may be counteracted by the aftershock slips of the Lixian Fault, which suggests Wenchuan earthquake could not significantly push the WMF closer to failure.
- Open Access
- 22 November 2021
- Published by American Geophysical Union (AGU) in Earth and Space Science
- Vol. 8 (11)
Spatiotemporal Evolution of Pore Pressure Changes and Coulomb Failure Stress in a Poroelastic Medium for Different Faulting Regimes
Spatial distribution of aftershocks has been explained by the co‐seismic static Coulomb Failure Stress changes (ΔCFS) to some extent. In practice, the earth's crust is porous medium saturated with fluid filling pores, cavities, cracks, and faults. The pore fluid flow plays an important role in the evolution of the stress field in porous medium of the earth. The ΔCFS can be influenced by the co‐seismic pore pressure changes and the post‐seismic pore fluid flow. In order to understand such impacts, we built 3D fault‐slip finite element models to study the evolution of the pore pressure changes and the ΔCFS based on the fully coupled poroelastic theory. The strike‐slip, reverse, and normal faulting models were investigated separately. The patterns of co‐seismic ΔCFS were similar to the elastic stress changes on a homogeneous half‐space, but there were considerable differences in details especially in the near‐fault area when considering co‐seismic pore pressure change. Due to post‐seismic pore fluid flow, the near‐field stress shadow narrowed gradually in the strike‐slip model. In the reverse faulting model, both the near‐field stress shadow and enhanced areas expanded. In the normal faulting model, the stress shadow near the epicenter reduced. The pore pressure changes decayed sharply at the beginning and then gradually approached to zero because the system evolves to drained conditions. These results can help us to understand the effect of pore fluid flow on Coulomb Failure Stress evolution and might improve the understanding of earthquake triggering and aftershock forecasting.
- 22 November 2021
- Published by American Geophysical Union (AGU) in Earth and Space Science
- Vol. 8 (11)
Open AccessSpatiotemporal Evolution of Pore Pressure Changes and Coulomb Failure Stress in a Poroelastic Medium for Different Faulting Regimes
Spatial distribution of aftershocks has been explained by the co‐seismic static Coulomb Failure Stress changes (ΔCFS) to some extent. In practice, the earth's crust is porous medium saturated with fluid filling pores, cavities, cracks, and faults. The pore fluid flow plays an important role in the evolution of the stress field in porous medium of the earth. The ΔCFS can be influenced by the co‐seismic pore pressure changes and the post‐seismic pore fluid flow. In order to understand such impacts, we built 3D fault‐slip finite element models to study the evolution of the pore pressure changes and the ΔCFS based on the fully coupled poroelastic theory. The strike‐slip, reverse, and normal faulting models were investigated separately. The patterns of co‐seismic ΔCFS were similar to the elastic stress changes on a homogeneous half‐space, but there were considerable differences in details especially in the near‐fault area when considering co‐seismic pore pressure change. Due to post‐seismic pore fluid flow, the near‐field stress shadow narrowed gradually in the strike‐slip model. In the reverse faulting model, both the near‐field stress shadow and enhanced areas expanded. In the normal faulting model, the stress shadow near the epicenter reduced. The pore pressure changes decayed sharply at the beginning and then gradually approached to zero because the system evolves to drained conditions. These results can help us to understand the effect of pore fluid flow on Coulomb Failure Stress evolution and might improve the understanding of earthquake triggering and aftershock forecasting.
- Open Access
- 1 January 2021
- Published by Taylor & Francis in Geomatics, Natural Hazards and Risk
- Vol. 12 (1) , 2048-2069
Estimation of seismic hazard around the Ordos Block of China based on spatial and temporal variations of b-values
Ruinous earthquakes have frequently occurred on the border of the Ordos Block, historically. Thus, where is a future disastrous earthquake around the Ordos? This problem is becoming more and more urgent both in the scientific community and the public, particularly after the 2008 Wenchuan earthquake (Ms = 8.0). To this end, we at first compute b-values in and around the Ordos Block based on the updated catalogue with earthquakes (M ≥ 0.0) from 1 January 2010 to 31 December 2020, and presented the spatial distribution of b-values. By examining the distribution feature of b-values, we chose 14 locations on the border of the Ordos Block where the b-value is relatively low and neo-tectonic faults are active. After that, we calculate the temporal variations of b-values in these 14 places, respectively. Comparing the 14 curves of b-variations with seismic precursors on b-values, we found that four of fourteen b-variations have characteristics of earthquake precursors. Therefore, these four places, located in the western end of the Weihe Graben, both the eastern and western ends of Hetao Graben, and the northeastern front of the Tibet, are considered to be seismically dangerous in the future. Naturally, the seismic damage in these areas is more serious than other places. Hence, the study is helpful in earthquake prediction and also can provide important clues for earthquake hazard assessment around the Ordos Block.
- 1 January 2021
- Published by Taylor & Francis in Geomatics, Natural Hazards and Risk
- Vol. 12 (1) , 2048-2069
Open AccessEstimation of seismic hazard around the Ordos Block of China based on spatial and temporal variations of b-values
Ruinous earthquakes have frequently occurred on the border of the Ordos Block, historically. Thus, where is a future disastrous earthquake around the Ordos? This problem is becoming more and more urgent both in the scientific community and the public, particularly after the 2008 Wenchuan earthquake (Ms = 8.0). To this end, we at first compute b-values in and around the Ordos Block based on the updated catalogue with earthquakes (M ≥ 0.0) from 1 January 2010 to 31 December 2020, and presented the spatial distribution of b-values. By examining the distribution feature of b-values, we chose 14 locations on the border of the Ordos Block where the b-value is relatively low and neo-tectonic faults are active. After that, we calculate the temporal variations of b-values in these 14 places, respectively. Comparing the 14 curves of b-variations with seismic precursors on b-values, we found that four of fourteen b-variations have characteristics of earthquake precursors. Therefore, these four places, located in the western end of the Weihe Graben, both the eastern and western ends of Hetao Graben, and the northeastern front of the Tibet, are considered to be seismically dangerous in the future. Naturally, the seismic damage in these areas is more serious than other places. Hence, the study is helpful in earthquake prediction and also can provide important clues for earthquake hazard assessment around the Ordos Block.
- Open Access
- 2 March 2020
- Published by MDPI AG in Applied Sciences
- Vol. 10 (5) , 1687
Effects of Barriers on Fault Rupture Process and Strong Ground Motion Based on Various Friction Laws
A barrier may induce a supershear rupture on a fault. This paper focuses on two questions: One is whether the existence of a barrier accelerates the propagation speed of a whole fault rupture, and the other is what are the effects of friction laws and strength of a barrier on the rupture propagation process. For these purposes, classical slip-weakening, rate-state, and modified slip-weakening friction laws are employed to simulate the effect of a barrier on the fault rupture process. The simulation results showed that the rupture speed of the fault obviously decreases when the rupture front propagates to the barriers, and the rupture speed obviously increases when the rupture front leaves barriers. It was also found that a barrier on a fault may induce a supershear rupture via the rate-state friction law. The simulation results also showed that with the increase of barrier strength, the rupture speed near barriers fluctuates more and more; when the barrier strength exceeds a certain level, a supershear rupture area appears on the fault; with the increase of barrier strength, the propagation distance of the rupture at supershear wave velocity correspondingly increases. In addition, with the increase of barrier strength, the overall rupture duration of the fault slightly increases. This indicates that a barrier cannot shorten the total duration of a fault rupture. Though a barrier will lead to a supershear rupture, it just regulates the distribution of the rupture speed on the fault surface. Moreover, with the increase of barrier strength, the peak ground acceleration caused by rupture through the barrier also increases, indicating that the existence of a barrier may lead to the intensification of seismic hazards.
- 2 March 2020
- Published by MDPI AG in Applied Sciences
- Vol. 10 (5) , 1687
Open AccessEffects of Barriers on Fault Rupture Process and Strong Ground Motion Based on Various Friction Laws
A barrier may induce a supershear rupture on a fault. This paper focuses on two questions: One is whether the existence of a barrier accelerates the propagation speed of a whole fault rupture, and the other is what are the effects of friction laws and strength of a barrier on the rupture propagation process. For these purposes, classical slip-weakening, rate-state, and modified slip-weakening friction laws are employed to simulate the effect of a barrier on the fault rupture process. The simulation results showed that the rupture speed of the fault obviously decreases when the rupture front propagates to the barriers, and the rupture speed obviously increases when the rupture front leaves barriers. It was also found that a barrier on a fault may induce a supershear rupture via the rate-state friction law. The simulation results also showed that with the increase of barrier strength, the rupture speed near barriers fluctuates more and more; when the barrier strength exceeds a certain level, a supershear rupture area appears on the fault; with the increase of barrier strength, the propagation distance of the rupture at supershear wave velocity correspondingly increases. In addition, with the increase of barrier strength, the overall rupture duration of the fault slightly increases. This indicates that a barrier cannot shorten the total duration of a fault rupture. Though a barrier will lead to a supershear rupture, it just regulates the distribution of the rupture speed on the fault surface. Moreover, with the increase of barrier strength, the peak ground acceleration caused by rupture through the barrier also increases, indicating that the existence of a barrier may lead to the intensification of seismic hazards.
- Open Access
- 1 January 2020
- Published by Taylor & Francis in Geomatics, Natural Hazards and Risk
- Vol. 11 (1) , 2123-2146
Physical mechanism for severe seismic hazard in the 2010 Yushu, China, earthquake (Mw= 6.9): insights from FEM simulations
The moderate-size 2010 Yushu earthquake (Mw 6.9) caused severe seismic damage around the source region which is situated on the sparsely populated hinterland of the Tibetan Plateau. But until now the mechanisms have not been well understood. To this end, we constructed the model to simulate the fault spontaneous rupture propagation, in which the realistic fault with a curved bend is imbedded. Our modeling results show that the special fault geometry controlled the rupture behavior, which encouraged rupture propagation speed transition from subshear to supershear with the speed of 5.17 km/s larger than the local shear wave velocity. Moreover, calculation results demonstrated that strong ground motion acceleration was greatly intensified by supershear ruptures, leading to widespread destruction. This may be the main reason why serious earthquake damage happened in the Yushu earthquake. In particular, we can see from numerous numerical experiments that rupture styles will be different if geometries of fault are varied. It is confirmed that it is the special geometry of the seismogenic fault of the Yushu mainshock produced grave seismic hazard. Thus, deeply investigating fault geometry will be helpful to better understand seismic source process and seismic hazard assessment.
- 1 January 2020
- Published by Taylor & Francis in Geomatics, Natural Hazards and Risk
- Vol. 11 (1) , 2123-2146
Open AccessPhysical mechanism for severe seismic hazard in the 2010 Yushu, China, earthquake (Mw= 6.9): insights from FEM simulations
The moderate-size 2010 Yushu earthquake (Mw 6.9) caused severe seismic damage around the source region which is situated on the sparsely populated hinterland of the Tibetan Plateau. But until now the mechanisms have not been well understood. To this end, we constructed the model to simulate the fault spontaneous rupture propagation, in which the realistic fault with a curved bend is imbedded. Our modeling results show that the special fault geometry controlled the rupture behavior, which encouraged rupture propagation speed transition from subshear to supershear with the speed of 5.17 km/s larger than the local shear wave velocity. Moreover, calculation results demonstrated that strong ground motion acceleration was greatly intensified by supershear ruptures, leading to widespread destruction. This may be the main reason why serious earthquake damage happened in the Yushu earthquake. In particular, we can see from numerous numerical experiments that rupture styles will be different if geometries of fault are varied. It is confirmed that it is the special geometry of the seismogenic fault of the Yushu mainshock produced grave seismic hazard. Thus, deeply investigating fault geometry will be helpful to better understand seismic source process and seismic hazard assessment.
- Open Access
- 1 January 2020
- Published by Taylor & Francis in Geomatics, Natural Hazards and Risk
- Vol. 11 (1) , 1196-1211
Spatial and temporal b-value precursors preceding the 2008 Wenchuan, China, earthquake (Mw = 7.9): implications for earthquake prediction
The sudden and unexpected Wenchuan earthquake (Mw = 7.9) occurred on the Longmen Shan Fault without any definite precursors. Therefore, a new computing strategy of b-values was applied in the study area around the epicenter of the Wenchuan mainshock, which is divided into nine circular sub-regions with a radius of 50 km uniformly distributed in space, ignoring geological structures. Then the temporal mapping of b-values is performed at each sub-region. The simulated results show that, at one sub-region where the Wenchuan shock is located, the b-value increased to ∼1.25 in 2002 and was almost unchanged for over five years until 2007. But, the b-value dropped to ∼0.5 in Nov. 2007 and remained the same for about three months. Then, b-value increased again since March 2008, and on May 12, 2008, the Wenchuan mainshock occurred. Indeed, according to previous results, this pattern of b-variations is a typical precursor for major earthquakes, which has not been discovered by many retrospect studies for the Wenchuan earthquake. Therefore, this study, which may capture subtle precursors on b-variations for major earthquakes, is significant for earthquake prediction and seismic hazard assessment.
- 1 January 2020
- Published by Taylor & Francis in Geomatics, Natural Hazards and Risk
- Vol. 11 (1) , 1196-1211
Open AccessSpatial and temporal b-value precursors preceding the 2008 Wenchuan, China, earthquake (Mw = 7.9): implications for earthquake prediction
The sudden and unexpected Wenchuan earthquake (Mw = 7.9) occurred on the Longmen Shan Fault without any definite precursors. Therefore, a new computing strategy of b-values was applied in the study area around the epicenter of the Wenchuan mainshock, which is divided into nine circular sub-regions with a radius of 50 km uniformly distributed in space, ignoring geological structures. Then the temporal mapping of b-values is performed at each sub-region. The simulated results show that, at one sub-region where the Wenchuan shock is located, the b-value increased to ∼1.25 in 2002 and was almost unchanged for over five years until 2007. But, the b-value dropped to ∼0.5 in Nov. 2007 and remained the same for about three months. Then, b-value increased again since March 2008, and on May 12, 2008, the Wenchuan mainshock occurred. Indeed, according to previous results, this pattern of b-variations is a typical precursor for major earthquakes, which has not been discovered by many retrospect studies for the Wenchuan earthquake. Therefore, this study, which may capture subtle precursors on b-variations for major earthquakes, is significant for earthquake prediction and seismic hazard assessment.