Thin-layer effects in glaciological seismic amplitude-versus-angle (AVA) analysis: implications for characterising a subglacial till unit, Russell Glacier, West Greenland

. Seismic amplitude-versus-angle (AVA) methods are a powerful means of quantifying the physical properties of subglacial material, but serious interpretative errors can arise when AVA is measured over a thinly-layered substrate. A substrate layer with a thickness less than 1/4 of the seismic wavelength, λ , is considered “thin”, and reﬂections from its bounding interfaces superpose and appear in seismic data as a single reﬂection event. AVA interpretation of subglacial till can be vulnerable to such thin-layer effects, since a lodged (non-deforming) till can be overlain by a thin (metre-scale) cap of dilatant (deforming) till. We assess the potential for misinterpretation by simulating seismic data for a stratiﬁed subglacial till unit, with an upper dilatant layer between 0.1– 5.0 m thick ( λ / 120 to > λ / 4, with λ = 12 m). For dilatant layers less than λ / 6 thick, conventional AVA analysis yields acoustic impedance and Poisson’s ratio that indicate contradictory water saturation. A thin-layer interpretation strategy is proposed, that accurately characterises the model properties of the till unit. The method is applied to example seismic AVA data from Russell Glacier, West Greenland, in which characteristics of thin-layer responses are evident. A subglacial till deposit is interpreted, having lodged (acous-tic impedance . . underlying (thickness


Introduction
Seismic reflection methods provide a powerful means of imaging the bed of glaciers and ice masses, and are increasingly used for quantifying substrate material properties (e.g. Smith, 1997;Nolan and Echelmeyer, 1999;Anandakrishan, 2003;Peters et al., 2007. Characterising the physical properties of subglacial material is a key research goal for improving the representation of basal dynamics in predictive ice-flow models (e.g. Pattyn, 1996;Payne, 1999;Truffer et al., 2001;Pimentel et al., 2010;Sergienko and Hulbe, 2011). A hard-bedded glacier is associated with different ice-flow regimes to a soft, sediment-bedded glacier, with the presence of englacial debris and/or subglacial water significantly impacting flow (Iverson et al., 2003;Cohen et al., 2005;Emerson and Rempel, 2007). Changes in the hydrological state of the subglacial system are invoked in the onset of glacier surging Woodward et al., 2003;Burgess et al., 2012). Furthermore, the style and efficiency of subglacial drainage is also influenced by substrate properties, with Darcian and canal flow regimes proposed for sediment substrates (Clarke, 1987;Walder and Fowler, 1994;Ng, 2000), and linked cavities and ice channels for hard-bedded glaciers (Röthlisberger, 1972;Iken and Bindschadler, 1986;Willis et al., 2012). Measurements of the hardness of previously-glaciated surfaces have also been used to inform palaeo-flow reconstructions (e.g. Reinardy et al., 2011).
Amplitude-versus-angle (AVA) analysis of seismic data is particularly useful, delivering at least two mechanical properties of a material, specifically acoustic impedance and Poisson's ratio (Aki and Richards, 2002;Gretener, 2003). AVA concerns the same seismic characteristics and principles to amplitude-versus-offset, AVO, but expresses the variation of reflectivity with actual incidence angle rather than its offset proxy. The method has had notable use in glaciology for identifying subglacial lakes , thick sequences of dilatant (deforming) till (Anandakrishnan, 2003;Peters et al., 2007) and transient changes in subglacial hydrology (Nolan and Echelmeyer, 1999;Kulessa et al., 2008).
Seismic AVA analysis has several key assumptions, including the requirement that a reflective interface separates two otherwise infinite half-spaces. The presence of thinlayering at that interface is a significant complicating factor in AVA interpretation (Swan, 1991;Bakke and Ursin, 1998;Nolan and Echelmeyer, 1999), as this condition can often be violated when surveying over a layered till. A till deposit can be structurally complex, with abrupt variations (both vertical and lateral) in physical properties (Evans et al., 2006) that, critically, are on a smaller spatial scale than the seismic wavelength (∼ 10 m in glaciology, depending on source characteristics and ice thickness; Smith, 2007). We therefore present an investigation of the effect of thin-layering on AVA analysis.
Although AVA practitioners within glaciology may be aware of the potential for thin-layer problems (e.g. Richards, 1988;Nolan and Echelmeyer, 1999), there have been few reported investigations of alternate strategies for interpreting the resulting AVA responses. We therefore review AVA theory and consider the definition of "thin", from a seismic perspective. We then show how thin-layers distort an AVA response, using synthetic seismic data that simulate an acquisition over dilatant (deforming) and lodged (nondeforming) subglacial till layers (e.g. Evans et al., 2006;Peters et al., 2007). We propose an interpretation strategy that honours the structure of the till deposit, and apply this to seismic data acquired on the Russell Glacier outlet of the West Greenland Ice Sheet. By acknowledging the potential for a thinly-layered structure, we are able to recover more information from the AVA response than is possible with conventional interpretation methods.

Reflection coefficients and amplitude-versus-angle responses
A seismic wavelet undergoes partial reflection when it encounters an interface between contrasting acoustic properties. The fraction of wavelet amplitude reflected, termed the "reflection coefficient", R, is influenced both by the magnitude of that contrast (specifically in terms of the velocities of P-(pressure) and S-(shear) waves, density and, to a lesser extent, seismic quality factor, Q), and the incidence angle at which the wavelet arrives. For a wavelet propagating at normal incidence (i.e. θ = 0) from the ith to the j th layer, the P-wave reflection coefficient is a simple function of the acoustic impedance within each layer: where Z denotes acoustic impedance (the product of density and P-wave velocity, v P ). Seismic quality factor, Q, also contributes to the observable reflectivity of an interface, although only significantly where the Q-contrast exceeds an order of magnitude (Bourbié and Nur, 1984;Odebeatu et al., 2006). Two proposed mechanisms for Q to contribute to the apparent reflectivity of an interface are (a) frequency dependency of the reflection coefficient, since different frequency components of a wavelet propagate at different velocities in dispersive materials (Xu et al., 2011), and (b) introduction of a phase lag into the recorded waveform (Lines et al., 2008), since the reaction time of a high-to-low-Q interface is frequency-dependant (Quintal et al., 2009;Morozov, 2011). Q-based reflectivity is undoubtedly an important consideration in quantitative seismic interpetatoin but, at this stage of analysis, we assume small Q-contrasts allowing the associated reflectivity contributions to be neglected. In the general case of non-normal incidence (i.e. θ > 0), a fraction of the incident P-wave is also converted to S-wave energy, hence reflectivity is also influenced by Poisson's ratio, σ , a measure of the stiffness of a material (Mavko et al., 2009): where v S is S-wave velocity. As v S approaches zero, its velocity through a fluid, Poisson's ratio reaches its theoretical maximum of 0.5. The general expression for reflection coefficient with incidence angle, R(θ ), is the "Knott-Zoeppritz equations", a complicated, non-linear system of angle and acoustic properties (Aki and Richards, 2002) and the variation of reflection coefficient with an angle is termed the "AVA response". Figure 1 shows AVA responses modelled at glacier beds underlain by material of varying (a) acoustic impedance and (b) Poisson's ratio (see Table 1). Throughout, the overlying ice has an acoustic impedance of Z ice = 3.50 × 10 6 kg m −2 s −1 ) and a Poisson's ratio of σ ice = 0.33 (Anandakrishnan, 2003). Reference curves for ice-rock and ice-water interfaces are also included as end-member contrasts in glaciological settings. Where the substrate acoustic impedance exceeds that of ice (Fig. 1a, models ii-v), the zero-incidence reflection coefficient is positive and, for these models, switches to negative (i.e. a polarity reversal) between 40-55 • incidence. These polarity reversals can be highly diagnostic of material properties (Anandakrishnan, 2003). In Fig.1b, the acoustic impedance of the substrate is fixed, hence all models express equal R(0), but AVA gradients (i.e. whether R(θ) increases or decreases) are strongly negative in model i, and strongly positive in model v (at least for θ < 55 • ). Contrasts in Poisson's ratio therefore exert the strongest influence on AVA gradient, with a positive AVA gradient associated with an increase in σ across the interface. In the glaciological context, strongly positive AVA gradients imply the presence of water in the substrate, either as a pore-fluid or as a free body (i.e. a subglacial lake).
The characteristics of numerous AVA responses can be summarised on an "AVA cross-plot" (Simm et al., 2000) ( Fig. 1c for previous model curves). Shuey (1985) linearises the Knott-Zoeppritz equations, simplifying the angular variation in reflectivity as where A and B are functions of velocity and density contrasts, estimated by linear regression of the AVA curve where sin 2 (θ) is plotted against R(θ). An AVA cross-plot is therefore a representation of the best-fit Shuey A-and B-terms (e.g. Fig. 1c). Under certain circumstances, A and B can be used to derive quantitative information about subsurface physical properties, but Shuey's linearisation is only strictly valid for cases of negligible refraction across an interface (i.e. for θ up to ∼ 30 • , and for small velocity contrasts). Velocity contrasts across the bed of a glacier, including those applied in our models, typically violate the latter condition, hence A and B cannot be used quantitatively. Nonetheless, the crossplot remains an efficient means of highlighting qualitative characteristics of AVA curves and we use it as such throughout this paper. In Fig. 1c, A and B are derived and cross-plotted for each curve of Fig. 1a and b (cross and circle symbols, respectively). Term A approximates R(0), hence the effect of changing the Poisson's ratio of the substrate is only expressed in B (expected, since all curves in Fig. 1b  common zero-incidence reflectivity). Variations in substrate acoustic impedance modify both A and B, although A shows the greatest change. Uncertainties in either A and B are negligibly small. The cross-plot representation of AVA curves will be revisited in the later sections of this paper.

Thin-layers in glaciological AVA analysis
The Knott-Zoeppritz equations are strictly defined for the interface between two otherwise infinite, homogeneous, isotropic, half-spaces, and the presence of thin-layering between these half-spaces violates this condition (Bakke and Ursin, 1998;Aki and Richards, 2002;Nolan and Echelmeyer, 1999). In seismic terms, a layer is "thin" if it is thinner than one-quarter of the dominant wavelength, λ, of the seismic wavelet (Widess, 1973). At this threshold, termed the "tuning thickness" (Widess, 1973), reflections from the layer's bounding interfaces superpose and a single interface is perceived. For thinner layers, with thickness around λ / 8, the composite response approximates the derivative of the original signal and further thinning also reduces its amplitude (Bakke and Ursin, 1998). For such ultra-thin layers, the composite response may also be influenced by intrabed multiples and mode conversions (i.e. respectively, wavelets that reverberate within the thin layer, and those that convert between P-and S-wave modes). Thin bed effects introduce significant potential for misinterpretation, and Swan (1991) cautions that their impact can be stronger than the underlying lithological AVA effect. How common are "thin-layers" in glaciological seismic settings? With a typical seismic wavelength of ∼ 10 m (Smith, 2007;Anandakrishnan, 2003;Horgan et al., 2008), vertical stratifications spaced more closely than some 2-3 m would be considered seismically thin. Although glaciers themselves are clearly good approximations to half-spaces, subglacial tills frequently contain metre-(or sub-metre-) scale contrasts (Smith, 2007), such as the transitional boundaries between dilatant and lodged till (e.g. Clarke et al., 1984;Echelmeyer and Wang, 1987;Iverson et al., 1988Iverson et al., , 1994Truffer et al., 2001;Boulton et al., 2001;Porter and Murray, 2001;Evans et al., 2006;Iverson, 2011;Reinardy et al., 2011). It would therefore be prudent to consider how thin layering affects a seismic AVA response whenever surveys are conducted over till substrates.
Numerous strategies, developed principally in the hydrocarbon sector but also in glaciology, are available for thinlayer AVA interpretations, although certain simplifications and/or assumptions invalidate them in glaciological settings. Lin and Phair (1993) quantify the properties of a thin layer located in an otherwise homogeneous half-space, assuming that its upper and lower interfaces therefore have equal magnitude reflection coefficients. This is clearly unsuitable where a thin dilatant till layer may be located between ice and lodged till, with markedly different acoustic con-trasts at either interface. Richards (1988) and Nolan and Echelmeyer (1999) derive generalised analytic expressions for thin-layer AVA responses by summing infinite series of intrabed multiples, but this approach neglects wavelet attenuation (effectively imposing an infinite quality factor). More realistically, a finite quality factor would more-rapidly attenuate intrabed multiples; furthermore, the longer a wavelet spends within an attenuating thin-layer, the greater the distortion of its waveform, hence the interference pattern between primaries and intrabed multiples varies with attenuation.
Here, we decompose the composite AVA response to a thin layer geometry using three approaches: a. travel-time analysis of the reflected wavelets that contribute to the composite response, b. assessment of the effective reflectivity of those reflections using the Knott-Zoeppritz equations, and c. forward modelling of seismic data for thin-layer geometries, and derivation of the corresponding AVA response.

Ray-tracing of reflected travel-times
Using Mathworks Matlab ® , we ray-traced travel-times for an isotropic, homogeneous glacier, overlying a horizontal, isotropic layer of dilatant till and a lodged till half-space ( Table 2). An ice thickness of 1 km thick was modelled, although the AVA response is independent of ice thickness provided that the relative range of offsets is maintained. The acoustic impedance of the dilatant till, Z dil , is 3.42 kg m −2 s −1 (Atre and Bentley, 1993; Vaughan et al., 2003) and its Poisson's ratio, σ dil , is 0.494, typical for water-saturated sediment (Gercek, 2007). Lodged till acoustic impedance, Z lod , is 3.90 kg m −2 s −1 (Atre and Bentley, 1993; Vaughan et al., 2003) and its Poisson's ratio, σ lod , is 0.32, representative of average values for granular silty soil (Gercek, 2007). In Sect. 3.3, our synthetic wavelet has a dominant frequency of 150 Hz (consistent with that observed in later analysis of field data) and a nominal wavelength in the dilatant till of 12 m (with v P = 1800 m s −1 ). We therefore vary the thickness of the dilatant till layer, h d , from 0.1 to 5.0 m (λ / 120 Travel-times for nine reflected ray paths ( Fig. 2) were calculated for a common midpoint (CMP) configuration of sources and receivers, offset from 0 to 5000 m (θ varies between 0-68 • ). Critical refractions are not considered at any interface since incidence angles are sub-critical for each wavelet mode. The nomenclature of each ray path describes the mode (i.e. P-or S-wave) of the wavelet as it crosses an interface, while upper-and lower-case letters respectively denote whether the section of the travel-path is a primary or intrabed multiple. We only consider ray paths that arrive at the surface as P-waves, since S-waves travel more slowly through ice (v P ≈ 2v S ) and do not interfere with the primary base-ice reflection, PP ( converts between S-and P-wave modes in the thin-layer, and PSSppP ( Fig. 2i) is a mode-converting intrabed multiple. Example travel-times for h d = λ / 4 = 3.0 m are derived, and plotted ( Fig. 3a) for clarity as a lag behind the arrival time of PP. The strong velocity reduction across the ice-till boundary implies that wavelets are strongly refracted towards the vertical (at PP's maximum incidence of 68 • , PPPP and PSSP are refracted to 27.8 • and 3.0 • , respectively), hence the lag between PP any other event is almost constant with incidence angle. The first 3.4 ms lag behind PP is shaded grey in Fig. 3a, and corresponds to the normal-incidence traveltime of a P-wave through the thin-layer. Given that λ / 4 is the limit of resolution, no events arrive within this lag, hence we would conclude that the basal reflection can be resolved form successive events. However, Widess (1973) states that although resolution criteria are fulfilled at the λ / 4 threshold, this pertains specifically to the first-break of a wavelet and factors including waveform and wavelet duration can still cause interference for successive wavelet half-cycles. In later analysis, our synthetic pulse has several half-cycles and amplitude picks are located within the second of these following the first-break. We therefore expect to observe interference effects beyond h d = λ / 4, until all events lag the second halfcycle of PP by more than 3.4 ms.

Effective reflectivity
Deriving the individual AVA responses of these reflections allows their specific contributions to be identified within the composite response. We use the Knott-Zoeppritz equations to obtain the "effective" reflectivity, R eff , of the nine events in Fig. 2 (Fig. 3b). Effective reflectivity differs from the reflectivity defined in Sect. 2 since it includes all reflection and transmission losses accrued by a wavelet along its propagation path. For the shallowest interface in a model there are no overlying interfaces; hence, a wavelet expresses R eff identical to the reflectivity defined in Eq. (1). However, a primary reflection from the nth interface expresses where T i (= 1 − R i ) is the transmission coefficient across an interface (squared since a wavelet crosses each interface twice as it propagates down from and up to the surface). If no correction for these losses is applied, effective reflectivity is the quantity derived from amplitude analysis of a seismic reflection. For this analysis, it is the effective reflectivity of an interface that governs how much each arrival contributes to the composite wavelet. In Fig. 3b, R eff for each arrival is plotted against the incidence angle of PP at the glacier bed, as would be performed if only a single reflection could be perceived.
The base ice reflection, PP, shows AVA characteristics expected from a transition from low-to-high σ and highto-low Z, specifically negative R(0) (= −0.011), positive AVA gradient (here, for θ < 50 • ) and a polarity reversal (here, at θ ≈ 10 • ). The P-wave reflection from the base of the thin-layer, PPPP, has the opposite character (positive R(0) = +0.066 and negative AVA gradient), although P-wave reflectivity within the dilatant till is low and its intrabed multiples are negligible. The equivalent reflectivity for S-waves is much higher: PSSP is very strong at 60 • incidence (some 4-times stronger than the maximum magnitude of PP) and its multiple wavetrain remains significant even after two intrabed reverberations. However, both P-and S-wave intrabed multiples show successive polarity reversals causing them to interfere destructively in the composite wavelet. This is especially relevant for the S-wave, and their contribution to the composite is less than the large-magnitude reflectivity would suggest.
Together, travel-time and reflectivity models suggest that the strongest contributions to the composite AVA response will be PP and PPPP, with weaker contributions from PSPP www.the-cryosphere.net/6/909/2012/ The Cryosphere, 6, 909-922, 2012 Table 2. Ice and till properties as used in "SKB2" forwardmodelling of seismic data (with reference to Vaughan et al., 2003;. Ice properties are fixed throughout, and S-wave quality factors are equated to reference values for P-waves.  Ray paths that feature intrabed P-wave multiples have negligible reflectivity. and PSSP. Intrabed multiples of the latter may be significant but, given their slower propagation velocity through the thin-layer, they only interfere with the basal reflection in the thinnest h d cases where their lag behind PP is sufficiently small.

Forward modelling of seismic data
The previous modelling approaches provide valuable insight into the character of composite AVA responses, but neglect wavelet attenuation. We therefore forward-modelled seismic data for the models in Table 2, using the software "SKB2" (Laboratoire de Géophysique Interne et Tectonophysique, University of Grenoble; Kennett and Kerry, 1979;Bouchon, 1981;Kennett, 1983;Mallick and Frazer, 1987), which delivers the three-component (i.e. vertical, radial and transverse) seismic response to a one-dimensional reflectivity model (homogenous, isotropic, horizontal layering) for any prescribed range of CMP offsets. Here, we consider only the vertical component of the seismic wavelet since the geophones in our field acquisition (Sect. 5) are also verticalcomponent only. We model a source pulse with 150 Hz frequency, and CMP offsets from 0 to 5000 m.
The P-wave quality factor (Q P ≈ 430) is measured from real seismic data (Sect. 5). The S-wave quality factor, Q S , is equated to Q P (e.g. Abercrombie and Rice, 2005) in the absence of available englacial S-wave attenuation measurements. Since there appear to be no reported examples of quality factors for in situ subglacial till, we allocate a representative value of Q P = Q S = 256 (Sain et al., 2009) obtained for unconsolidated sea-floor mud. This analogue is deemed appropriate since the pore-pressure regime in seafloor sediment is probably more representative of that beneath an ice sheet than a measurement made on a surface till exposure. Although lower quality factors have been recorded for sea-floor sediment (e.g. Best et al., 1994;Ayres and Theilen, 2001;Riedel and Theilen, 2001;Sain and Singh, 2011), we allocate the highest values to minimise the Q-contrast (and hence its reflectivity contribution) across the glacier bed.
Selected forward-modelled seismic responses are shown in Fig. 4, together with a control response corresponding to ice overlying a dilatant till half-space. Trace amplitudes are scaled for geometric spreading and attenuation losses (corrections defined for the dominant frequency of 150 Hz), assuming propagation in ice only. Travel-times correspond to the lag behind the first-break of PP, following the application of non-stretch normal-moveout (NMO) corrections (Perroud and Tygel, 2004) using v P = 3800 m s −1 as measured from field data in Sect. 5. Although we cannot explain the residual moveout that is present in each gather, we consider it negligible since it is less than 0.6 % of the reflection moveout originally displayed (∼ 3 ms vs. ∼ 500 ms). Included in each panel is the corresponding ray-traced travel-time model from Sect. 3.1. Since our amplitude picks (grey triangles; see next section) are located at the absolute maximum of the second wavelet half-cycle, travel-time models are referenced to the amplitude pick in zero-incidence traces.
The characteristics of PP are well-represented in the control response, which has a polarity reversal close to 10 • and increasing reflectivity up to 50 • incidence. In the thinnest case (h d = λ / 120), the most significant contributions to the composite follow within 2 ms of PP and a single reflection, featuring no polarity reversal, is perceived. The λ / 24 case also shows a single reflection, but the interference pattern is different compared to the previous example. At λ / 12, there is sufficient lag between PP and the S-wave modes that the strong PSSP phase can be resolved, yet the basal reflection still appears as a single event. At λ / 4, the theoretical limit of resolution, characteristics of PP can be distinguished, although there is still interference with PPPP. A weak event arrives later than these events: this could be PSPP, implying that the remaining interference is exclusively between P-wave modes. Finally, as λ / 4 is exceeded, PP and PPPP are sufficiently separated that interference is only apparent close to the low-amplitude polarity reversal, and the amplitude variation of PP is similar to the control response.

Composite AVA curves
Composite AVA curves derived from all modelled data are shown in Fig. 5 a and b showing h d less and greater than λ / 8, respectively). For reference, the effective reflectivities of PP and PPPP are also shown. However, since model curves only consider the vertical component of the synthetic wavelet, the equivalent component of the effective reflectivities of PP and PPPP is obtained by applying a cos(θ) obliquity factor (e.g. Yilmaz, 2001). As stated in Sect. 3.1, amplitudes of synthetic wavelets are picked at the absolute maximum of the second wavelet half-cycle. The first-half cycle has very-low amplitude and is barely perceptible even in our noise-free synthetics, hence would most-likely be indistinct given noise conditions in real data. The composite R(0) was obtained by comparing each zero-incidence amplitude to its equivalent in the control model, in which R(0) is explicitly known (i.e. from substituting acoustic impedances for ice, Z ice , and dilatant till, Z dil , into Eq. 1). Since wavelet amplitude scales linearly with reflection coefficient, the ratio of amplitudes is equivalent to the ratio of reflection coefficients.
For h d < λ / 8, the observed AVA responses appear to be hybrids of PP and PPPP curves, with the positive AVA gradient of the former (for θ < 40 • ) but the positive R(0) of the latter. For h d between λ / 8 and λ / 6, AVA curves are more complicated and feature an abrupt switch in gradient near 20 • incidence. In Fig. 3b, this is close to the angle where the effective reflectivity of PP first exceeds that of PPPP, hence we suggest that the character of the composite is influenced by the characteristics of PPPP at small incidence angles and by those of PP at θ > 20 • . Furthermore, Fig. 4 implied that the high reflectivity PSSP phase would not contribute to the composite at h d > λ / 8, hence the apparently complex AVA response actually results from a simpler interference pattern between fewer reflections and mode conversions. As h d approaches λ / 4, the characteristics of PP emerge although, as suggested in Sect. 3.1, some interference is still observed: the zero-incidence reflectivity is shifted to more-negative values with respect to the model curve. Finally, for h d > λ / 4 (= 5.0 m), the characteristics of PP are resolved.
Composite curves are added to the AVA cross-plot (Fig. 5c, which also includes example models from Fig. 1c). As before, best-fit Shuey terms are derived for incidence angles between 0-30 • . The cross-plot reaffirms the observation that AVA curves derived from models of h d ≤ λ / 8 are hybrids of PP and PPPP, since their best-fit terms occupy a region of the cross-plot between the two model curves. Since the two events become distinct for the thickest h d models, A and B tend towards those of PP. The intermediate-thickness h d models (λ / 8 < h d < λ / 4) have B-terms close to zero that are associated with greater uncertainties (uncertainties for all other regressions, and in the A-terms throughout, are negligibly small). Recall that the B-term is associated with the AVA gradient, hence the increased uncertainty here relates to the gradient switch in the associated AVA curves.

Glaciological interpretation of thin-layer AVA responses
Having shown the AVA variability caused by a change in h d , we now consider its significance for a glaciological interpretation. The key interpretative issue is if a thin-layer AVA curve resembles a plausible AVA response to a single interface, it could be misdiagnosed as such. Whilst this is unlikely for the intermediate-thickness models (i.e. λ / 8 < h d ≤ λ / 6), as their abrupt gradient switch is implausible as a singleinterface AVA curve, the other curves are vulnerable to misinterpretation and the incorrect acoustic impedance and/or Poisson's ratio could be allocated to the subsurface.

Conventional interpretation of AVA responses
If it was assumed that the composite AVA curves actually are responses to single-interface geometries, Table 3 shows the apparent acoustic impedance, Z app , and apparent Poisson's ratio, σ app , that they define (also comparing them to model dilatant and lodged till values). To estimate Z app , we rearrange Eq. (1) for Z j and substitute the observed zero-incidence reflectivity, R obs (0), and the acoustic impedance allocated for ice (Z ice = 3.5 kg m −2 s −1 ). For σ app , we define the best-fit AVA curve to the observed response, and record the Poisson's ratio that this defines. Curve fitting was performed by simulating many manifestations of AVA curves, for θ < 30 • , with many trial values of density, v P and v S substituted into the Knott-Zoeppritz; the best-fit curve is defined by the parameter set that minimises the root-mean-square (RMS) error between the trial and observed AVA curve. Tabulated Z app and σ app values suggest that it is only as h d approaches its resolvable that the AVA analysis uniquely identified either till type, specifically the dilatant till. For the thinnest geometries, Z app is indicative of lodged till, whereas σ app may be more suggestive of dilatant till throughout. Although these parameters could be considered to deliver the bulk properties of the substrate, they are incompatible when allocated to a single till type: a subglacial till with a high acoustic impedance would not be expected to also have a high Poisson's ratio, since these quantities point towards opposing water saturation (e.g. Smith, 2007;Vaughan et al., 2003).
The observation of such "paradoxical" quantities could therefore indicate the presence of thin-layers in the substrate. We therefore suggest an interpretative strategy to extract more information about subglacial thin-layering from a composite AVA response.

Thin-layer interpretation for AVA responses
In Fig. 5c, the best-fit terms for all models with h d ≤ λ / 6 plot between those of PP and PPPP. Likewise, despite changes in the AVA gradient, the observed zero-incidence reflectivity in Fig. 5a is always between those of PP and PPPP. The AVA curves for h d of λ / 120, λ / 40 and λ / 24 each express R obs (0) of 0.043, whereas that for h d = λ / 12 expresses R obs (0) = 0.052. These are close to the sum of the effective zero-incidence reflectivity at the upper and lower interfaces of the thin-layer (−0.011 and +0.066, respectively, with a sum of +0.055). With this observation, we use Eqs. 1 and 4 to formulate a method for interpreting thin-layer characteristics from the composite AVA curve.
Since R obs (0) approximates the summed effective reflectivity of the two interfaces, it follows that: and, by substituting Eq. (4), assuming that the upper and lower interfaces of the thin-layer are the first two horizons that are encountered by the seismic wavelet. We next assume that the acoustic impedances of ice, Z ice , and dilatant till, Z dil , can either be measured or assumed, and rearrange Eq. (6) to estimate the acoustic impedance of the lodged till, Z lod . Substituting these into Eq. (1), we fix R 1 (0) and express R 2 (0) in terms of Z lod and Z dil , such that The Cryosphere, 6, 909-922, 2012 www.the-cryosphere.net/6/909/2012/ Rearranging for Z lod therefore gives (note: the rearrangement of these equations is trivial if Z dil is alternatively required).
For our models, R 1 (0) = −0.011, Z ice = 3.5 × 10 6 kg m −2 s −1 and Z dil = 3.42 × 10 6 kg m −2 s −1 . With the zero-incidence reflectivities observed in Fig. 5a, Eq. (8) predicts Z lod of 3.81 × 10 6 kg m −2 s −1 (h d = λ / 120, λ / 40 and λ / 24) and 3.88 × 10 6 kg m −2 s −1 (h d = λ / 12). Models were simulated using Z lod = 3.90 × 10 6 kg m −2 s −1 , hence these estimates are within −2.3 % of the model quantity. The composite AVA response can therefore be interpreted in terms of two acoustic impedances and h d , which must have a maximum value of λ / 8. It is difficult to resolve the Poisson's ratio for the lodged till, although the Poisson's ratio of the dilatant till can be assumed to be at least that expressed by the AVA curve.
For the responses to the thicker layers, where h d = λ / 8 and h d = λ / 6, the observed zero-incidence reflectivity is +0.058 and +0.063, respectively. These imply Z lod is 3.93 and 3.97 × 10 6 kg m −2 s −1 , with respective errors of +0.5 % and +1.5 %. Table 3 suggests that the Poisson's ratio expressed by these AVA responses is neither representative of dilatant nor lodged till, and it may therefore be difficult to constrain a representative value in this case. However, the observation of a gradient switch in an AVA response may indicate that the thickness of the thin-layer is between ∼ λ / 8 and ∼ λ / 6. For the cases of h d approaching λ / 4, Table 3 suggests that both the acoustic impedance and Poisson's ratio of the dilatant till can be accurately estimated despite the remaining interference between PP and PPPP. However, this suggests that the AVA analysis is insensitive to the properties of the underlying lodged till, and it is therefore unlikely that the composite AVA curve would be recognised as the response to a thin-layer geometry. This represents a significant inter-pretative risk, since dilatant till properties may be allocated to the whole substrate instead of the vertical thickness of the thin-layer. A careful interpretation of the actual seismic traces would be recommended, to check for events following immediately behind the basal reflection that could be identified as PPPP.

Application to real data
We now consider thin-layer effects in the interpretation of actual seismic AVA data, acquired on the Russell Glacier outlet of the West Greenland ice sheet (∼ 70 km inland of its western margin) during summer 2010. At the position of the AVA analysis, the ice is 1.08 ± 0.01 km thick (v P = 3800 ± 40 m s −1 , from prestack migration velocity analysis, Sheriff and Geldart, 1999;Bradford et al., 2009) and the bed reflection is planar over ∼ 750 m with an in-line dip less than 3 • . An orthogonal seismic profile suggested there were no potential sources of out-of-plane-reflectivity in this area. The ice surface was non-uniform over the 4 kmlong profile, varying from dry and solid in places to wet and channelised, and friable in others. Consequently, the interference between the primary arrival and its free-surface source ghost are unpredictable from shot-to-shot, hence ghosting represents a source of random error in this analysis. The topographic variation along the line was less than 40 m, and never more than 3 m between successive shot locations, 80 m apart. Sources were 250 g Pentex charges, installed at ∼ 3 m depth, and data were recorded with a Geometrics GEODE system at 48 concrete-mounted, 100 Hz, vertical-component geophones installed at 10 m intervals. Figure. 6a shows a series of traces from a CMP supergather (a combination of successive CMP records). For the measured velocity and ice thickness, the 550 m offset range of this record corresponds to an incidence range of 4-18 • at the glacier bed. For these angles, and elsewhere in the record (up to a maximum angle of ∼ 23 • ), the basal reflection shows no polarity reversal, either with respect to the direct wave or within the offset range of a gather. In conventional AVA analysis, the absence of a polarity reversal would immediately exclude a dilatant till substrate.

Data pre-processing
Minimal data processing was imposed to avoid the potential introduction of amplitude artefacts. We applied a static correction (25 ms) to remove the detonator delay from our traces, a bandpass filter to suppress noise outside of the useful signal bandwidth (centred on 150 Hz), and muting to suppress residual groundroll. Trace amplitudes were corrected for geometric spreading and attenuation losses. To correct for the latter, we obtained Q P by considering the spectral ratio (Dasgupta and Clark, 1998)  involving spectral ratios of direct (Gusmeroli et al., 2010) and/or reflected waves , our approach delivers a depth-averaged quality factor for the whole ice column using wavelets that have sampled virtually the same volume of ice. The multiple arrivals we consider reflect twice from the glacier bed and once from the underside of its surface. This implies that our analysis is more vulnerable to effects of frequency-dependent reflectivity than that of Peters et al. (2012), who consider a primary and a direct arrival, but we consider than the benefit of deriving Q P from wavelets that share a common propagation path to outweigh any such detrimental impacts. The quality factor we measure at this location was Q P = 436 ± 140. Primary and multiple amplitudes were also compared to estimate the zero-offset reflection coefficient of the glacier bed (e.g. Roethlisberger, 1972;King et al., 2003), using primary-multiple pairs within a 0-10 • incidence range. Be-yond this range there is too great a contrast in the AVA response of a primary and its corresponding multiple to make this analysis reliable. The median of the R obs (0) values is +0.110, and observations at successive angles are normalised to this. Figure 6b shows the AVA response we derive at the glacier bed. The reflection points that contribute to this response span a 300 m-wide section of the bed, hence reflectivity is spatially-averaged. However, our wavelet has a Fresnel zone (Lindsey, 1989) of 165 m, hence this distance is comparable to the intrinsic spatial resolution in our data. The signal-tonoise ratio of response is enhanced by averaging bed reflectivity in 2 • incidence angle bins. Uncertainty bars are centred on the median value within each bin, and span the interquartile range  in both the reflectivity and incidence angle directions. A key source of uncertainty is shot-to-shot variability in the coupling of sources and geophones, which was minimised by normalising each trace by the root-mean-square direct wave amplitude in common shot and common receiver gathers. However, given the additional imprecision of v P , Q P and ice thickness, we accept a wide range of potential error in this analysis.

AVA analysis
A Bayesian statistical analysis was applied to derive the best-fit AVA curve (solid line, Fig. 6c) to the observed AVA data, honouring its uncertainties, and to establish the uncertainty in derivative parameters. Given a set of observations and their associated uncertainties, the Bayesian analysis assesses the probability that a given set of model parameters are explained by those observations. In our case, the observed data is R obs , denoting the mean reflectivity in any angle bin, model parameters, M, are the densities and seismic velocities to be substituted into the Knott-Zoeppritz equations, and the probability is denoted as P (M | R obs ). If the uncertainty within a given angle bin, θ i , has a Gaussian distribution, the probability, P (R obs | M), that a set of observations is explained by the model is where ψ i is the standard deviation within each angle bin and g(θ i , M) is the model value output by given θ i and M values, with g here representing the system of Knott-Zoeppritz equations. P (M | R obs ) and P (R obs | M) are related through Bayes' rule that states where P (M) is the "prior" probability for the model parameters. Since we impose no preferred values on the model, we assume that P (M) is constant and P (M | R obs ) = P (R obs | M).
The best-fit M-values are therefore those which maximise P (M | R obs ) or P (R obs | M). However, the distribution of The Cryosphere, 6, 909-922, 2012 www.the-cryosphere.net/6/909/2012/ data in our angle bins has a double-Gaussian form. For the 12-14 • angle bin, the RMS difference between the observed data and best-fit double-Gaussian distribution is 2.5 × 10 −5 , rising to 1.9 × 10 −3 for the best-fit single-Gaussian distribution (a hundred-fold increase in fitting error). To use a double-Gaussian distribution, Eq. (9) is modified to where the subscripts 1 and 2 in the R obs (θ i ) and ψ i terms indicate the mean and standard deviation of the two components in the double-Gaussian distribution. Aside from any earlier analysis so far conducted in this paper, we note at this point that the double-Gaussian distribution may be providing an independent indication of a thinlayer problem, and specifically a subglacial thin-layer that is laterally variable. If a proportion of our reflection points portray a thin-layer of a certain thickness, and that thickness changes beneath nearby reflection points, it is likely that the distribution of reflectivity will be bimodal within a given angle bin. Such statistical analyses may therefore be useful as an interpretative tool for thin-layer analysis, but we will revisit this in future research.
The best-fit AVA curve suggests that the substrate has an acoustic impedance of 4.36 ± 0.18 ×10 6 kg m −2 s −1 and Poisson's ratio close to 0.5 (the precision is difficult to establish since σ is insensitive to v S at low S-wave velocity). The acoustic impedance is strongly indicative of a lodged till deposit with low porosity (e.g. Vaughan et al., 2003), yet Poisson's ratio suggests a highly water-saturated substrate. Accordingly, both A-and B-terms in the Shuey linearisation of this curve are also positive (cross-plotted in Fig. 6c). This contradictory set of observations may diagnose the presence of a thinly-stratified substrate.

Thin-layer analysis
We invoke a thin-layer argument in order to reinterpret the AVA response as a layered till geometry, using Eq. (8). We assume that there is negligible Q-contrast across the interface, and that the maximum thickness of the thin-layer is λ / 6, given the observation of positive A-and B-terms in the Shuey linearisation.
Assuming that maximum v P through the thin-layer of dilatant till is 1800 m s −1 , the maximum wavelength of a 150 Hz wavelet is 12.0 m. Our λ / 6 criterion therefore suggests that the thickness of the thin-layer does not exceed 2 m. For the lower velocity of 1600 m s −1 , as quoted in Peters et al. (2008), the implied wavelength is 10.5 m and the maximum thickness of the thin-layer is 1.8 m.
We therefore conclude that our study site on Russell Glacier is underlain by a stratified till deposit, having dilatant till overlying a lodged till unit. The dilatant till has a maximum thickness of 2 m and a Poisson's ratio that approaches 0.5, indicating high water saturation and non-cohesive sediment (Gercek, 2007). The underlying lodged till has an acoustic impedance of 4.26 ± 0.59 × 10 6 kg m −2 s −1 , suggesting both low porosity and water saturation. However, its acoustic impedance is unlikely to exceed 5×10 6 kg m −2 s −1 , hence the till is probably unlithified (Vaughan et al., 2003). Evidence from statistical analysis gives an initial indication that the thickness of the thin-layer may be laterally variable over a 300 m range.

Discussion and implications
By recognising and invoking a thin-layer AVA argument, we were able to interpret our data with a greater degree of complexity and consistency than with a conventional procedure. Although our approach is based on a qualitative similarity to previous model outputs, a meaningful numerical inversion would require measurement or assumption of additional properties, particularly Q P and Q S for subglacial till (and glacier ice in the latter case). We acknowledge that, throughout, we assume homogenous and isotropic layer properties, having measured both v P and Q P from depth-averaged observations. AVA analysis is complicated in anisotropic cases (e.g. Tsvankin, 2001), and numerous research papers have shown detectable englacial velocity and attenuation contrasts corresponding to changes in ice temperature (Kohnen, 1974;Peters and Anandakrishnan, 2010) and/or crystal orientation (Horgan et al., 2008;Gusmeroli et al., 2012). However, we consider that initial recognition of thin layer issues provides a significant improvement over a conventional interpretation approach. Alternative quantitative approaches to deriving the actual thickness of the thin layer include cepstral relationships, in which the periodicity within an amplitude spectrum can be used to interpret interference patterns (e.g. Hall, 2006;Rubino and Velis, 2009) but trials of these suggested that our data had an insufficient signal-to-noise ratio.
In some cases, thin-layer considerations may be irrelevant. For a glacier underlain by intact bedrock, thin-layer geometries may be neglected since the substrate is unlikely to have strong stratifications (although anisotropic AVA effects may be important if the bedrock is fractured; Xu and Tsvankin, 2001). Likewise, AVA studies of subglacial lakes  are likely to be free from thinlayer effects, provided that the free-water thickness is sufficient to avoid interference. Thin-layer effects should always be considered, however, where AVA analysis is conducted over a subglacial till, given the potential for these units to be thinly stratified. This is particularly relevant where the www.the-cryosphere.net/6/909/2012/ The Cryosphere, 6, 909-922, 2012 interpretation is supplied to predictive models of glacier flow, given the sensitive relationship between the substrate of a glacier and its flow regime (Pattyn, 1996;Truffer et al., 2001;Pimentel et al., 2011;Sergienko and Hulbe, 2011). The risk of neglecting thin-layer considerations is particularly acute where the thickness of a dilatant till approaches its nominally resolvable threshold. For these geometries (i.e. h d ≈ λ / 4), Fig. 5 suggested that the composite AVA response was dominated by the reflection from the upper boundary of the thin-layer, and that from its lower boundary was masked. It is therefore possible that the properties of the underlying unit would be omitted from an interpretation. For example, Peters et al. (2007) use AVA responses to suggest that till deposits beneath Bindschadler Ice Stream (West Antarctica), between 5-20 m thick, are laterally variable and feature cyclic transitions between wet and stiff regimes. Thin-layer considerations allow an alternate model to be suggested, in which the substrate is extensively composed of stiff till, but with a wet till cap whose thickness varies laterally. Where the cap is absent, conventional AVA analysis is valid, and the properties of the stiff till can be recovered. However, where the thickness of the wet till approaches λ / 4 (5 m for Peters et al., 2007), the reflection from the underlying stiff till is masked but the interface is nonetheless present in the substrate.
These interpretations are both plausible substrate geometries, but represent end-members of analysis. Conventional AVA interpretation accommodates all of the seismic variability in a bulk change in till properties and allocates a single set of quantities to the whole thickness of the substrate, whereas a thin-layer interpretation fixes those properties and accommodates any variability in a lateral change in the position of an additional interface. Both interpretations allow the presence of wet and stiff tills, but disagree on the nature of transitions between them. A useful progression from this research would therefore be determining a means of extracting the characteristics of two till types when the thin-layer is close to its resolvable threshold.

Conclusions
Seismic AVA analysis is a powerful method for quantifying the physical properties of subglacial material, although serious misinterpretations can result where thin stratifications are present in the substrate. Here, we have shown how thin layer effects manifest themselves in glaciological AVA responses, and how they can then be interpreted in terms of the thickness, acoustic impedance and Poisson's ratio of the subglacial material. The recognition of thin layer AVA issues is an important step forward in improving the potential to image and characterise the subglacial environment, a key aspect given the importance of subglacial processes in predictive ice-flow models and palaeo-reconstruction. We would therefore recommend that thin-layer AVA considerations become a routine part of interpretation, particularly where stratified subglacial deposits are anticipated.