Gravity Variation in Siberia : GRACE Observation and Possible Causes

We report the finding, from the GRACE observation, of an increasing trend in the gravity anomaly in Siberia at the rate of up to 0.5 μgal yr-1 during 2003/1 2009/12, in the backdrop of a negative anomaly of magnitude on the order of ~-10 mgal. In consideration of the non-uniqueness of the gravitational inverse problem, we examine in some detail the various possible geophysical causes to explain the increasing gravity signal. We find two geophysical mechanisms being the most plausible, namely the melting of permafrost and the GIA post-glacial rebound. We conclude that these two mechanisms cannot be ruled out as causes for the regional gravity increase in Siberia, based on gravity data and in want of ancillary geophysical data in the region. More definitive identification of the contributions of the various causes awaits further studies.


IntROduCtIOn
The twin-satellite Gravity Recovery and Climate Experiment (GRACE) mission has been yielding high-precision measurements for the Earth's gravity field since its launch in 2002.A positive "secular" increase of gravity in the Siberia proper was detected in the GRACE time-variable gravity data (wherein, here the term "secular" means the linear trending of the observed signal during the period of observation).In this paper we report and analyze this gravity signal, compare it with other relevant numerical models, and raise and discuss the possible geophysical causes.In consideration of the non-uniqueness of the gravitational inverse problem, we raise several candidate geophysical causes for the secular gravity increase.We conclude with two of these causes being the most plausible.
One main candidate, as raised recently by Muskett and Romanovsky (2009) and Ogawa (2010), pertains to the melting of Siberia's permafrost layer as a result of recent global warming.Together with an increase in precipitation in the studied region during the studied period, it causes the water holding capacity of the soil to increase, and hence an increase in gravity.However, the dynamic processes or scenarios through which this happens have yet to be thoroughly investigated.
Another plausible cause we consider is the glacial isostatic adjustment (GIA).Beginning ~120 K years before present, the last (late Pleistocene) ice-age glacial period underwent rapid melt-away of continental ice sheets around 10 K years ago.Based on present-day crustal uplift and other geophysical observations, accepted models of GIA, such as ICE-5G (VM2) (Peltier 2004), put the largest post-glacial rebound in Laurentia in North America with secondary contributions from Fennoscandia in the northern Europe and (to a lesser certainty) the Antarctica, but little direct GIA activity in Siberia which lies in the same latitude band.On the other hand, there have been conflicting hypotheses about the existence and actual extent of possible Siberia glaciation, ranging from the "restricted" scenario of moderate glaciation in Arctic and eastern Siberia only (Brigham-Grette et al. 2001) to the maximum reconstruction of extensive ice sheet coverage across the Arctic Eurasia (Grosswald and Hughes 2002).The uncertainty stems from the lack of ground surveys, geophysical observations and systematic chronology in the area (Stauch and Gualtieri 2008).Terr. Atmos. Ocean. Sci., Vol. 22, No. 2, 149-155, April 2011

GRAVIty FORmulAtIOn And dAtA
The Earth's disturbing potential of gravity (relative to the ellipsoidal Earth) is customarily expressed in terms of its spherical harmonic decomposition (Kaula 1966;Heiskanen and Moritz 1967): where M is the Earth's mass, and the scale-length a is often chosen conveniently to be the mean Equatorial radius.The dimensionless coefficients C nm and S nm are the (normalized) Stokes coefficients of degree n and order m, and P nm is the 4π-normalized associated Legendre function.An alternative way of expressing the gravity field is by the free-air gravity anomaly (relative to the ellipsoidal Earth): Note the emphasis on the higher degree or shorter wavelength components by the factor n -1 in Eq. ( 2).Physically the Stokes coefficients are simply the normalized multipoles of the density distribution of the Earth r t^h (e.g., Chao and Gross 1987;Chao 2005): where dV = r 2 sinθ dθ dλ dr is the element of the volume integral over V 0 .Equation ( 3) is how the gravity observation relates to the mass distribution of the Earth, the (3-D) inversion of which, however, is notoriously non-unique (e.g., Chao 2005).Launched in 2002, the twin-satellite GRACE mission yields accurate measurements for the Earth's gravity field via high-precision determination of the satellite orbits and the satellite-to-satellite perturbations thereof.The standard GRACE data are released in the form of spherical harmonic coefficients typically up to degree n = 120 at monthly intervals (Tapley et al. 2004; see http://podaac.jpl.nasa.gov/grace/documentation.html).Corrections for solid and ocean tides (including pole tides) are applied during the solution; the atmospheric and oceanic mass transport effects are also corrected for beforehand using general circulation model products (ECMWF for the atmosphere and the mass-conserving baroclinic ocean model which is a further develop-ment of Hamburg Ocean Primitive Equation, "HOPE", see GRACE AOD1B product description document for detail) in order to eliminate the associated high-frequency aliasing errors.The C 20 solutions are replaced, as is the common practice for better quality, with those obtained independently from the satellite laser ranging observations (Cheng and Ries 2007).Note also that the degree-1 terms are excluded from the above equations because they vanish by definition as long as the Earth's center of mass is chosen to be the origin of the reference system, as is the case here.
The utility of the GRACE data lies primarily in their (tiny) temporal variations, or the time-variable gravity (TVG) signals, as the observed TVG in any given region signifies mass redistribution in the region over time (subject to the non-uniqueness of the inversion).Thus, removing the temporal average (or the "static") gravity field from each monthly solution reveals TVG in the form of monthly "maps."

SECulAR tVG In SIbERIA
Figure 1 gives the GRACE map of the static free-air gravity anomaly g in mgal [Eq.( 2)] evaluated at r = a up to spherical harmonic degree 60.As stated, this is the reference field to be removed from the GRACE data to reveal the TVG signals.We shall come back to examine Fig. 1 later.
Our studied time period for TVG is 7 years long from 2003/1 to 2009/12.The GRACE TVG data obtained after removal of the static field are plagued with short-wavelength noises on regional to local scales including the "striping" noise resulting from the satellite orbit artifacts, so spatial filtering is necessary (Wahr et al. 1998).In this paper we apply the anisotropic fan filter of Zhang et al. (2009), combined with the correlated-error filter [with degree-4 polynomial fit for coefficients of spherical harmonic order 6 and higher, see Swenson and Wahr (2006)], with equivalent filter length (hence spatial resolution or half-wavelength) of 350 km.
We then analyze the temporal signal at each grid point on the map by fitting the (monthly) time series with seasonal (annual + semi-annual) sinusoids plus a linear term.The coefficient of the linear term is a measure of the secular rate of the increasing or decreasing trend of TVG at the given location.Figure 2 shows the distribution of the resultant secular rate of g (μgal yr -1 ) around the Arctic region.As GRACE TVG data have already had the atmospheric and oceanic gravity signals corrected out, many remaining TVG signals are evident, especially on land, which have been linked to secular mass transports on or in the Earth.For example, recent studies have identified the secular decreasing of mass due to the melting away of coastal Greenland icesheet (e.g., Luthcke et al. 2006) and southern Alaskan glaciers (e.g., Chen et al. 2006), as well as the increasing TVG associated with the GIA post-glacial rebound of Laurentide (e.g., Tamisiea et al. 2007) and Fennoscandia (e.g., Steffen et al. 2008).Also evident, but less conspicuous, is the TVG increase in the Siberia proper.Below we shall search for the possible geophysical causes for this gravity signal.

tectonic Signal?
The shield of the Siberia Craton sits far away from any known plate boundaries and away from active seismicity or volcanism.The scheme of tectonics of Siberia has many versions but presently no agreeable consensus (e.g., Glebovitsky et al. 2008).Few time-variable geodetic measurements were made in a systematic manner in Siberia barring GRACE.Available ground data exist only for the ~2000 km region around Lake Baikal; they are sporadic, rather erratic and noisy, and lying too far south to the area in question (San'kov et al. 2003).There the conventional leveling before the 1980s and more recent GPS campaign data seem to suggest general uplift (particularly on the eastern part).So do the four regional GPS/GIS continuous stations around Baikal (station KSTU, IRKT, ULAZ, and ULAB in Mongolia) relative to station URUM in Xingjiang, China, thereby a general slow vertical uplift on the order of 1 mm yr -1 can be deduced (cf.San'kov et al. 2003).It is not clear whether this crustal displacement is related to the regional tectonics (Tatavian and San'kov 2008, personal communication), perhaps to the existence of the Baikal rift valley.At this point we cannot completely rule out the tectonics as one source of the observed TVG in Siberia.On the other hand however, the general, albeit slow, uplift is consistent with a weak GIA post-glacial rebound to be discussed later.

land Hydrological Signal?
Next we turn attention to climate-related mass variations.If one makes the assumption that the observed TVG comes only from mass redistributions occurring on the surface of the Earth, for example the redistribution of oceanic or continental hydrological masses, then the 3-D Eq. ( 3) reduces to a form where the 2-D inversion for surface mass change becomes unique (Chao 2005): Here are the harmonic coefficients of the surface density variation / , , cos t tP exp imm ^h, kn l is the load Love number of degree n, and the factor 1 kn + l ^h is to "undo" the mass loading effect.Upon division by the density of pure water 1 g cm -3 , Eq. ( 4) gives the convenient surface "equivalent water thickness" ETH.Note that ETH sees an even higher emphasis than the TVG anomaly in Eq. ( 2) on the higher degree or shorter wavelength components by the factor 2n + 1.Actually, their map  representations are very similar to each other (apart from their units).In this section we shall only show TVG and not ETH for simplicity and consistency with the rest of the paper.
We first examined the atmospheric contribution (which has already been removed from the GRACE data) according to the ECMWF output data, and those of the NCEP model for comparison.They show no secular variability in Siberia during the studied period 2003/1 -2009/12.We can thus readily rule out the possibility that the observed secular TVG comes from errors in the atmospheric pressure variation.
On the other hand, there may be some secular landhydrological anomaly in Siberia presumably related to climate change during the studied period.We take the GLDAS (Global Land Data Assimilation System) model (Rodell et al. 2004) and compute the global gravity anomaly of the monthly soil moisture mass distribution according to Eq. ( 2) and subject it under the same filtering processing as with GRACE data.Figure 3a gives the secular hydrological TVG (in μgal yr -1 ) in Siberia showing a large region of negative trend lying toward the west, and a weaker positive trend toward the east, of the region of the GRACE-observed secular TVG increase.
We then take the difference of GRACE -GLDAS (i.e., Fig. 2 minus Fig. 3a), the resulting map is shown as Fig. 3b, signifying what Siberia "should" see without the GLDAS land hydrology -a more extensive and more uniform posi-tive TVG in central Siberia.We have also taken another land hydrology model, namely the NCEP-CPC (Climate Prediction Center), and found according to that model a TVG trend pattern similar to, but weaker than, GLDAS in general, and hence predicting even less influence.Overall, we can conclude that, to the extent of uncertainty of the land hydrological models, the land hydrological contribution in TVG trend in Siberia is insignificant during the studied period.
On the other hand, conventional hydrological models do not treat or include permafrost which actually plays a major role in the water budget in Siberia's changing climate, as to be discussed in the next section.

Permafrost Signal?
Global warming induces melting of permafrost and thickening of the active layer (the layer near permafrost surface that melts in summer and freezes in winter) over a great expanse of area in the Arctic region.From east to west, Lena, Yenisei and Ob are the three main Eurasian Arctic rivers running through the Siberian permafrost, with thickness reaching a few hundred of meters.At first glance, the melting of permafrost would lead to a higher run-off, where an increasing discharge to the Arctic Ocean (~2.0 gigaton yr -1 ) was indeed observed (Peterson et al. 2002); a resultant reduction, rather than increase, of the water storage would be expected.However, an inter-annual increase of precipita-  tion has also been observed -the CMAP (Climate Prediction Center Merged Analysis of Precipitation) data show a clear increasing trend in the Lena river basin over the last few years in these dry regions with normal annual precipitation of only ~200 mm (Sugimoto et al. 2003).In addition, the peculiar dynamical behavior of the continuous permafrost region as a consequence of the impermeability of the frozen soil, which prohibits melt water from migrating into depth, should be considered.Moreover, excess water from summer precipitation does not always run off, but can be stored as ice lens during winter.Such water could be used in the next dry summer, stabilizing the Taiga forest.Repeated freezing and melting of the permafrost due to warming would create and expand large depressions, or thermokarst "alas" lakes, on the soil surface, as suggested by visual observations (Iijima et al. 2010), enabling larger storage of water.These behaviors are not modeled in conventional land hydrology models.Thus, suppose the GRACE-observed secular TVG increase is all caused by the increase of net water storage of water-holding capacity in the permafrost region.The GRACE TVG data translate [by Eq. ( 4)] into ETH of 10.7 ± 1.2 mm yr -1 for the Lena basin, 5.1 ± 1.1 mm yr -1 for the Yenisey basin, and 6.8 ± 2 mm yr -1 for the Ob basin, equivalent to mass increase of 30.7, 12.9, and 9.6 gigaton yr -1 , respectively.The detail of the above is reported in Ogawa et al. (2010) (see also Muskett and Romanovsky 2009), which also studied the relations with the regional temperature changes.
We can compare the above GRACE-observed values of the Lena and Ob basins (where reliable run-off data are available) with those estimated from the CMAP data, namely 11.8 ± 1.1 and 14.8 ± 1.7 mm yr -1 .The latter are obtained by integrating the budget P-E-R, where precipitation (P), evapo-transpiration (E) and river run-off (R) were derived by observation and model calculations.Such estimates give the net budget regardless of any mechanism or role played by the permafrost.Thus, we see that the GRACE-observed TVG increase agrees well with the hydrological modelderived water mass increase for the Lena basin, but only accounts for just over a half of the latter for the Ob basin.The latter result suggests a negative mass contribution has also been in action.

GIA Post-Glacial Rebound Signal?
Under GIA, Laurentia and Fennoscandia represent mass void following the relieve of the last ice-age icesheet load around 10 K years ago but since only partially refilled by a post-glacial rebound of the mantle.Such regions are therefore characterized by a negative static gravity signature, as indeed shown in Fig. 1.At the same time, these regions should also exhibit a positive (increasing) trend in the secular rate of TVG, as shown in Fig. 2.
The Laurentide peak static g anomaly is about -34 mgal, of which about 40% is attributed to GIA, or -14 mgal.The corresponding peak secular rate of g from GRACE is ~+1.5 μgal yr -1 (Tamisiea et al. 2007).Fennoscandia peak static g anomaly is about -10 mgal; the GRACE peak secular rate is 1 -1.2 μgal yr -1 (Steffen et al. 2008).In both cases, the timescale of the GIA recovery is ~10000 years.
The corresponding peak values for the Siberia anomaly are ~-10 -15 mgal and ~0.5 μgal yr -1 .It should be noted that the GIA behavior applies well to the Siberia gravity anomaly in a qualitative sense but with a somewhat longer timescale, even though the GIA activity is considered to be absent in Siberia [e.g., the currently accepted past ice model ICE-5G (VM2) of Peltier 2004].In fact, Siberia along with Laurentia and Fennoscandia are the three regions in the high-latitude northern Hemisphere that manifest the same gravitational pattern and behavior said above.However, the latter is of course only a necessary, rather than sufficient, condition to equate a region with post-glacial rebound.
We shall point out that the TVG signal of the postglacial rebound if observed on the ground would be about twice as big as the GRACE-observed but of the opposite sign, that is, of a larger negative (decreasing) trend, because the gravity effect of the geometrical uplift of the ground surface would offset the gravity effect of increasing regional mass.We await such gravity measurements in the Siberia region to be reported.

COnCluSIOnS
We report the GRACE TVG observation in Siberia during 2003/1 -2009/12 -an increasing trend in the gravity anomaly g at the rate of up to 0.5 μgal yr -1 , in the backdrop of a negative anomaly of magnitude on the order of ~-10 mgal.We then examine in some detail the various possible geophysical causes to explain the increasing gravity signal.We find two geophysical mechanisms being the most plausible, namely the melting of permafrost and the GIA post-glacial rebound; more definitive identification of their contributions awaits further studies.In the least we conclude that these two mechanisms cannot be ruled out as causes for the regional gravity increase in Siberia, based on gravity data and in want of ancillary geophysical data in the region.
The need to raise various possible causes is traced back to the fundamental non-uniqueness of the gravitational inverse problem.Besides the fact that the observed data are the sum total of all contributing sources, just knowing the Stokes' harmonics alone (as in the case of the GRACE data) is inadequate to determine the depths at which the observed TVG originates (Chao 2005).Thus, to further this investigation, one obviously needs independent geophysical and geodetic data sources, such as ground-based long-term gravity, vertical crustal motion from GPS or leveling or possibly satellite altimetry (Hwang 2009, personal communication), extensive geological and geomorphologic survey, as well as ancillary meteorological and land hydrological observations, particularly those of the changing arctic permafrost environment which has yet to be better monitored and understood.

Fig. 1 .
Fig. 1.Static field of free-air gravity anomaly g (mgal) according to the GRACE data up to maximum harmonic degree of 60 over 2003/1 -2009/12, in a polar projection down to 50 degree N. Note the negative patterns in Siberia.