The sensitivity of GNSS measurements in Fennoscandia to distinct three-dimensional upper-mantle structures

The sensitivity of global navigation satellite system (GNSS) measurements in Fennoscandia to nearby viscosity variations in the upper mantle is investigated using a 3D finite element model of glacial isostatic adjustment (GIA). Based on the lateral viscosity structure inferred from seismic tomography and the location of the ice margin at the last glacial maximum (LGM), the GIA earth model is subdivided into four layers, where each of them contains an amalgamation of about 20 blocks of different shapes in the central area. The sensitivity kernels of the three velocity components at 10 selected GNSS stations are then computed for all the blocks. We find that GNSS stations within the formerly glaciated area are most sensitive to mantle viscosities below and in its near proximity, i.e., within about 250 km in general. However, this can be as large as 1000 km if the stations lie near the center of uplift. The sensitivity of all stations to regions outside the ice margin during the LGM is generally negligible. In addition, it is shown that prominent structures in the second (250–450 km depth) and third layers (450–550 km depth) of the upper mantle may be readily detected by GNSS measurements, while the viscosity in the first mantle layer below the lithosphere (70–250 km depth) along the Norwegian coast, which is related to lateral lithospheric thickness variation there, can also be detected but with limited sensitivity. For future investigations on the lateral viscosity structure, preference should be on GNSS stations within the LGM ice margin. But these stations can be grouped into clusters to improve the inference of viscosity in a specific area. However, the GNSS measurements used in such inversion should be weighted according to their sensitivity. Such weighting should also be applied when they are used in combination with other GIA data (e.g., relative sea-level and gravity data) for the inference of mantle viscosity.


Introduction
It is well known that observations of the glacial isostatic adjustment (GIA) process allow us to determine the earth's viscosity structure, especially that beneath formerly glaciated areas such as Fennoscandia and North America.So far, the most frequently employed GIA data are relative sea-levels, global navigation satellite system (GNSS) measurements and the gravity-rate-of-change data from the Gravity Recovery and Climate Experiment (GRACE) satellite mission (e.g., Peltier, 2004;Wu et al., 2013).These data are commonly used to infer lithospheric thickness and radial variations of mantle viscosity.Owing to recent improvement in modeling techniques and advances in computational power, lateral variations of both lithospheric thickness and mantle viscosity can also be inferred.In view of that, it is important to understand the capability of the many GIA observations for the determination of lateral lithospheric and mantle variations.This study will analyze how sensitive class "A" GNSS stations of the EUREF Permanent Network (EPN, Bruyninx et al., 2013) are to distinct areas in the upper mantle beneath Fennoscandia.(Note that class "A" stations are the best and well-maintained stations of the EPN -they have position accuracy of 1 cm at all epochs of the time span of the used observations (Bruyninx et al., 2013).) Sensitivity (or Fréchet) kernels of GIA observations at a specific GNSS station show how sensitive they are to viscosity variations in a specific region of the mantle in comparison Published by Copernicus Publications on behalf of the European Geosciences Union.
to another region.The methodology in computing sensitivity kernels for a laterally homogeneous earth can be found in Mitrovica and Peltier (1991) and Peltier (1998), while that for a laterally heterogeneous earth is described in Wu (2006).The fundamentals of Fréchet kernel sensitivity are also developed in seismology (see Dahlen and Tromp, 1998, for a review).
Sensitivity kernels of GIA observations to radial changes in viscosity have been calculated in Mitrovica and Peltier (1991, 1993, 1995), Peltier and Jiang (1996a, b), and Peltier (1998).These studies showed that sensitivity is generally higher in the upper mantle than in the lower mantle.This is especially true for data in Fennoscandia, where the resolving power of GIA observations is too low to provide accurate inference of lower-mantle viscosity (Mitrovica and Peltier, 1993;Steffen and Wu, 2011).However, the data in North America can see deeper because the load is wider and the sensitivity is higher down to the shallow part of the lower mantle (Mitrovica and Peltier, 1995).
The sensitivity kernels for selected stations of the BIFROST (Baseline Inferences for Fennoscandian Rebound Observations, Sea Level, and Tectonics) project to radial viscosity variations have been studied by Milne et al. (2004).Interestingly, they found sufficiently high sensitivities for the lower mantle.This was not supported by Steffen et al. (2006), who showed with a 3-D model, that lateral variations in lower-mantle viscosity do not affect the GNSS velocity field in Fennoscandia.As pointed out in Wu (2006), the sensitivity of the Fennoscandian data to the lower mantle may actually be due to Laurentia.
When lateral viscosity is included in the earth model, the normal mode formulation of Mitrovica and Peltier (1991) and Peltier (1998) no longer applies or becomes impractical to apply due to mode coupling (Wu, 2002).To overcome this, Wu (2006) showed that the sensitivity kernel can equivalently be computed from the difference in response between a model with a small but fixed perturbed viscosity in a single mantle block (or layer) at the location of interest and the response of the reference model without the perturbation.In Wu (2006), an axisymmetric (2-D) earth model and simplistic ice load, with size comparable to the Laurentian ice load, are used to study the sensitivity of global GIA data.Later, Steffen et al. (2007) employed a 3-D earth model with a realistic (4-D) ice model to study the sensitivity of GNSS stations in Fennoscandia.An advantage of the latter study is that all three horizontal components can be investigated.In addition, the spatial resolution is much higher, with element block sizes of 600 km × 600 km or 1000 km × 1000 km.In any case, both studies showed that sensitivity is highest within the formerly glaciated area.The radial or vertical variation of sensitivity kernels for uplift rate peaks around 300-450 km depth but becomes small below the upper mantle.Also, both the load distribution and the deglaciation history strongly affect the magnitude of sensitivity (Steffen et al., 2007).Furthermore, if one is interested in the lateral viscosity outside the former ice margin, then the horizontal velocities should be analyzed, although that also depends on the size of the perturbed mantle region.
It should be noted that in the models of Wu (2006) and Steffen et al. (2007), the perturbed viscosity has a fixed magnitude and lies in a single rectangular block in an otherwise laterally homogeneous mantle.When lateral heterogeneity in the mantle is taken into account, the shape of the blocks must be modified to reflect the shape of the lateral heterogeneity.Also, the magnitude of the viscosity in the block must reflect the true viscosity value there.This study will include such changes for a 3-D viscosity model and a realistic ice load history will also be used.
The GNSS stations where sensitivity kernels are computed in this study (Fig. 1) belong to the class "A" stations, except for station Vaasa, which is of class "B" (positions with an accuracy of 1 cm at the epoch of minimal variance of the station).The selected stations are also used in BIFROST investigations to GIA (see Lidberg et al., 2010).
The aims of this study are (i) to investigate the sensitivity of velocity fields at selected GNSS stations to certain regions of the mantle with similar viscosity and location relative to the former ice margin; and (ii) outline where (future) GNSS stations in Fennoscandia would be helpful to identify lateral viscosity changes.The next section describes the model in more detail.This is followed by the presentation and discussion of the results.Finally, the conclusion is presented in Sect. 5.

Modeling
The finite-element method is used to model the GIA process in Fennoscandia.The earth model used is flat with isotropic, compressible, Maxwell-viscoelastic layers, but the finite-element model allows both vertical and lateral variations to be taken into account.This model is described in Steffen et al. (2006) and is based on the approach of Wu (2004) which has been used successfully in many GIA investigations in North America (e.g., Wu, 2005), Fennoscandia (e.g., Steffen et al., 2006), the Barents Sea (e.g., Kaufmann and Wu., 1998), Antarctica (e.g., Kaufmann et al., 2005) and Iceland (e.g., Schmidt et al., 2012).The model consists of a central area of 3000 km × 3000 km size, where each element has a horizontal length of 100 km.The ice-load history model FBKS8 (Lambeck et al., 1998) is applied to the surface in the central area.Surrounding the central area is a 60 000 km wide peripheral area, that is connected to infinity horizontally with semi-infinite elements.
In this preliminary study, where we are not interested in the sensitivity of small-scale features in the mantle, we continue to use the laterally heterogeneous model U3L1_V1 in Steffen et al. (2006) to define the viscosity and shape of the blocks.This model has a uniform 70 km thick lithosphere on top of a laterally heterogeneous mantle which consists of  four layers in the upper mantle and another four in the lower mantle.The lateral viscosity variations in each layer of the upper mantle are converted from the SH shear wave tomography model S20A (Ekström and Dziewonski, 1998) using the scaling relationship derived from Ivins and Sammis (1995), but modified by Steffen et al. (2006).The viscosity structure within the four upper-mantle layers is shown in Fig. 2 with solid black contour lines.Model U1L1_V1 of Steffen et al. (2006) is used as the reference model, on which the blocks with lateral viscosity are superposed.It has an uppermantle viscosity value of 4 × 10 20 Pa s, which is a good average value of upper-mantle viscosity beneath Fennoscandia (Steffen and Wu, 2011).Lower-mantle viscosity is set to 2 × 10 22 Pa s (Steffen and Kaufmann, 2005).Elastic parameters (density ρ, shear modulus µ and bulk modulus κ) of the model are obtained by volume-averaging the values in the Preliminary Reference Earth Model (PREM; Dziewonski and Anderson, 1981) and are shown in Table 1 of Steffen et al. (2006).
To investigate the sensitivity of lateral heterogeneity in the mantle, finite elements with similar viscosities (i.e., within one order of magnitude in the first mantle layer and within half an order of magnitude in the three other layers) are grouped together to form blocks, so that the blocks reflect the lateral viscosity structure within each layer (see red lines in Fig. 2).These blocks of similar viscosity are further subdivided into groups of blocks that lie inside the former ice margin and those lying outside in order to study the effect of location relative to the former ice margin.In addition, we design three groups of blocks in the center of uplift that have the same shape in all four layers.This will help us investigate the sensitivity as a function of depth without the shapeeffects of the blocks.The number of blocks in the first mantle layer is 22, while that in the second and bottom layers is 19 and the third layer has 18 blocks.That means 78 models are generated -where each model has a different block of lateral heterogeneity included in Model U1L1_V1.One should note that different ice-load models (e.g., ICE-5G; Peltier, 2004) or different seismic tomography models (e.g., Grand et al., 1997) will give different shapes of the subdivisions.But, that should not significantly change the major conclusions of this investigation which is quite general in nature.

Results
The normalized sensitivity kernels K lj (r i ) of block j in layer i at location l is computed based on the approach of Wu (2006), which is modified from Peltier (1998): where δp l is the difference between the prediction p 3−D l of the perturbed 3-D model and the prediction p 1−D l of the reference model U1L1_V1 at GNSS location l. (Here, the prediction p l is one of the three velocity components.)δm j (r i ) is the (dimensionless) viscosity perturbation of block j in layer i (i.e., the difference between the log of the viscosity of the block in model U3L1_V1 and that in model U1L1_V1 which is log of 4 × 10 20 Pa s).Also, V j (r i ) is the fractional volume of this particular block.The latter is given by where V j (r i ) is the block volume, and V model the volume of the entire central area in the model, which includes the upper and lower mantle.(For example, V 1 (r 2 ) refers to the volume of block 1 in the second layer shown in Fig. 2.) Normalization by this term is useful as we are only interested in the relative amplitude of the sensitivity kernels, i.e., which viscosity block is comparatively more sensitive to the particular measurement at a GNSS station.Unlike the approach of Wu (2006), we introduce the extra dimensionless normalization factor {V max (r i )}, which is the value of the maximum fractional volume of the four layers investigated.This normalization is introduced here mainly for plotting purpose.The kernels for the three velocity components are calculated for the location of each of the 10 selected EPN stations.Thus, we are able to analyze the relative sensitivity of the station to every block.Figure 3 presents two typical sensitivity kernels for the three velocity components (EW, NS, Z) to all the different viscosity blocks in the model at the two stations Kiruna and Brussels.Kiruna (Fig. 3a) is located above block 1 and also not too far away from the center of rebound.Figure 3 clearly shows that sensitivity of any velocity component to a block in the first mantle layer (70-250 km depth) is small.However, sensitivity is highest for the blocks right underneath the station or close to it, provided that the blocks lie within the former ice margin.The largest sensitivity for the vertical velocity is in block 1 of the second mantle layer (250-450 km depth).This block also has the highest sensitivity for any velocity component at any one of the 10 selected stations.Sensitivity is smaller in the two neighboring blocks in the north (block 4) and south (block 2), and become almost negligible for all other blocks in layer 2. In blocks 1, 2 and 5 in the third mantle layer (450-550 km depth) as well as blocks 1 and 4 in the fourth mantle layer (550-670 km depth), sensitivity of the vertical component is also larger than in other blocks of that layer.Sensitivities for horizontal velocities are generally smaller, but they also peak in the second layer and their amplitudes decrease with depth.The exception is the EW component in block 5of the third layer, which shows the largest sensitivity for a horizontal component of all stations.
For the station in Brussels (Fig. 3b), which is located outside the former ice margin, sensitivities of all velocities to any block in any layer become almost negligible.This is true even for the viscosity blocks directly underneath the station.Any sensitivity that shows up marginally is related to horizontal velocities and they are located within the former ice margin.
The presentation of kernels as in Fig. 3 does not allow us to see visually where the blocks with significant sensitivity lie in the map.To overcome this, we set an arbitrary threshold value for the normalized sensitivity and plot only the blocks with sensitivity above the threshold value in Figs.4-13.(Likewise, in the rest of the paper, when we say the data is sensitive to a block, we mean that the sensitivity of the block is above the threshold.)Among the 10 stations, Kiruna shows the largest kernel values, and Brussels has the lowest.After testing, we find a threshold of 3 mm yr −1 for the vertical velocity.A threshold of 1 mm yr −1 is found for  For example, GNSS velocities measured at Kiruna (Fig. 4) have sensitivity to several blocks in each layer.However, not all components of the velocity field are sensitive to the same viscosity block.In the first mantle layer, vertical velocity is insensitive to any viscosity block but horizontal velocities show sensitivity to blocks west of the station at the edge of the former ice sheet.In the second mantle layer, the vertical and NS components have sensitivity to the block below the station as well as to the block south of it.In addition, the NS component is also sensitive to the block north of it while blocks in the west and east are detectable by the EW component.This shows that horizontal velocities may provide information for adjacent blocks.In the third layer, there are more viscosity blocks with sufficiently large sensitivity than in the second layer.The vertical component is sensitive to the underlying block and blocks surrounding it, the EW component to blocks west and east of the station, and the NS component to the underlying block and all three blocks south of the station within the former ice margin.Both NS and EW component are also sensitive to a block southeast of  the station.At the bottom of the upper-mantle layer, the underlying block and the one north of it influence the vertical component.Several blocks have sensitivity in the horizontal components, which have similar characteristics as found for the third mantle layer.
The reader can continue to look at how sensitivity varies with depth at each station in Figs.5-13.However, it is more profitable to compare the block sensitivities for all stations at a certain layer of the upper mantle.This is what we will do next.

Discussion
Let us begin with the first mantle layer (70-250 km depth).Horizontal components of all stations except Brussels are sensitive to at least one block of the area along the Norwegian Atlantic coast (blocks 6-8 in UM1 in Fig. 2), where there is a strong gradient in the viscosity of the first mantle layer.Since lateral thickening of the lithosphere from west to east under Fennoscandia (Steffen and Wu, 2011) can appear as a strong viscosity gradient in the first mantle layer, sensitivity to viscosity blocks in the first mantle layer along the Norwegian Atlantic coast also implies sensitivity to lithospheric thickness variations there.The station of Vaasa is the only one where the vertical component is sensitive to a block (the underlying one).Therefore, the vertical component does most likely not contain sufficient information of the first mantle layer.
For the second mantle layer (250-450 km depth), the GNSS stations in the center of rebound or in the southern     west of the station reflect sensitivity in the EW component.

Solid
The two stations near the center of rebound, Skellefteå and Vaasa, are sensitive to the largest number of blocks in the second layer.The number of blocks in the second layer is reduced as the station moves away from the center of rebound.Also, these blocks are generally located in the vicinity of the GNSS station, i.e., within a radius of about 500 km.Svetloe   and Brussels show no sensitivity to any block in this layer, while Riga has sensitivity (for all components) to block number 7 in the second layer (Fig. 2).Block number 7, which lies below the southern part of the Scandinavian Peninsula, is of special interest.Many stations show sensitivity of at least one velocity component to this particular block.In other words, the viscosity of this block   affects the measurements of these stations within the former ice margin, and strongly affects the value of viscosity inverted from such GNSS measurements.Also, by comparing the viscosity inverted from the vertical component of GNSS stations within the former ice margin and the viscosity inverted without the stations in the center, one should get a feeling of the accuracy of the viscosity inverted for this area.
For the third mantle layer (450-550 km depth), the number of blocks with enough sensitivity (for any velocity component) is larger than that for the second layer.In general, the same characteristics as for the second layer apply, but there are also some additional findings.Vertical velocity of the stations near the center of rebound and in the southern part of the Scandinavian Peninsula has sensitivity to surrounding blocks and partly (Kiruna, Skellefteå and Vaasa) to blocks about 1000 km away from the station.Horizontal components at the stations Svetloe and Brussels are now sensitive to certain blocks which are located more than 800 km away from the station.Riga and Svetloe, which both lie at intermediate distance from the center, have sensitivity for one of the horizontal components to one block along the Norwegian coast, which is on the opposite side of the glaciated area.However, note that the sensitivity values for most of these additional findings are close to the chosen threshold value and so these results should be interpreted with care.
The sensitivity of the fourth mantle layer (550-670 km depth) has similar characteristics as that of the third mantle layer, but the vertical component with enough sensitivity is only found in the underlying block and/or blocks within 400 km distance from the station.The stations in Riga and Svetloe still have sensitivity to the Norwegian coast area.Brussels has sensitivity to an area closer to the station than in the third mantle layer.Horizontal velocities of most stations are sensitive to blocks with distances up to 1000 km away, and thus provide viscosity information of the lowest part of the upper mantle underneath the Scandinavian Peninsula.
Our results strongly support the usage of stations near the center of rebound (e.g., Skellefteå or Vaasa) to investigate the viscosity structure in the upper mantle below Fennoscandia.The vertical component gives information of the viscosity structure in an area of 500 to 1000 km around the station from about 250 to 670 km depth.Horizontal velocities may enlarge this area to more than 1000 km, especially in the third and fourth layer.The farther one goes away from the center, the less information can be obtained.
An interesting result is that the horizontal components at many GNSS stations -even those on the other side of the former ice margin (e.g., in Riga and Svetloe) -have enough sensitivity in almost all the layers in the upper mantle below the Norwegian Atlantic coast.Thus, thorough analysis of the horizontal velocities in Fennoscandia can probably result in better estimation of viscosity or lithospheric thickness variations there.
Stations outside the former glaciated area are of limited value, e.g., sensitivities found for Brussels are very close to the threshold value and thus are quite small compared to values found for stations in the center of rebound.
Oceanic areas far off the coast, i.e., the ones that were never affected by ice load on top, do not show any significant sensitivity at any GNSS station.Similarly, blocks in the southwest do not have enough sensitivity.However, if ice load on the British Isles is investigated, then the situation may be different.Future investigations with British GNSS stations should analyze their potential sensitivity for the area.

Conclusions
Unlike previous studies, this paper includes realistic structures of lateral viscosity variation under Fennoscandia to investigate the sensitivity of GNSS measurements in 10 selected stations.These GNSS stations are backbones of the EPN and the BIFROST project and thus represent excellent and well-maintained stations of high accuracy.We employed a 3-D finite element model that has been commonly used in the last two decades.A realistic ice load of the ice model FBKS8 (Lambeck et al., 1998) was also applied.
Our results confirm previous findings (see e.g., Milne et al., 2004;Steffen et al., 2007) that GNSS stations are most sensitive to viscosity changes underneath a station, but mainly at a depth between 250 and 550 km.Both horizontal and vertical velocities show significant sensitivities provided that the GNSS station and the block are located within the former ice margin.The depth of sensitivity depends on the ice thickness -thinner ice gives less information on the fourth layer, which confirms the resolving power of GIA data in general.
For stations closer to the center of rebound or mid-distance between the center and the ice margin, the sensitivity is largest for the viscosity blocks right underneath (thus to a lateral extent of about 250 km), but also a few other blocks nearby have sensitivity if these blocks or parts of such a block are located within a lateral distance of about 500 to 1000 km.The latter is in contrast to the findings by Steffen et al. (2007), who showed that the sensitivity of neighboring blocks is mainly negligible.This difference is related to the regular block structure used in Steffen et al. (2007).To test or confirm this conclusion, future studies should use a different block structure, which is based, for example, on a different seismic tomography model.
Stations outside the former glaciated area do not have sufficient sensitivity to viscosity directly underneath.This is different to the findings by Steffen et al. (2007), who found that horizontal velocities of such stations might be helpful.This is probably due to the approach of averaging the kernels of a block they used.This may have increased the kernel value for blocks that covered glaciated and non-glaciated areas.It should be noted that Steffen et al. (2007) already suggested using a more realistic viscosity block structure for a sophisticated analysis to find out if their result is correct or not.In turn, it is indicated that horizontal velocities at stations outside the former glaciated area have sensitivity to certain regions within the former glaciated area at 450-670 km depth, which should be further investigated in the future as well.
Regarding the planning of future GNSS stations for GIA research, Wu et al. (2010)  States and NW Russia.Our results clearly support this argument as both regions are located within the former glaciated area.Furthermore, the dense network installed in the countries of Norway, Sweden, Finland and Denmark will be further densified in the next few years.In Sweden, for example, the network will consist of 400 stations by 2020 (Lantmäteriet, 2011).Recently, 20 stations of the existing network have been proposed as new EPN stations (Engberg et al., 2013), which means that the quality of the observed data will increase.Together with the existing network of GNSS stations, they should allow a thorough investigation of lateral viscosity structure under Fennoscandia.
The results from this study are helpful in future investigations on lateral variations of mantle viscosity and lithospheric thickness.We recommend a careful grouping of GNSS velocity measurements from selected areas, e.g., from the uplift center or the Baltic States, to investigate the vertical viscosity profile underneath the center of rebound or the viscosity structure of the Norwegian Atlantic coast, respectively.In this or in combined analyses with other GIA observations, such as relative sea-level data or gravity on ground and in space (Steffen et al., 2012(Steffen et al., , 2014)), the observations should be properly weighted according to their sensitivity to a specific region.

Figure 1 .
Figure 1.Location of EPN/BIFROST station used in our investigation.The blue line marks the location of the ice margin at the Last Glacial Maximum.
investigated optimal locations in Fennoscandia and suggested more stations in the Baltic www.solid-earth.net/5