Instrumentation and field methods

We use a high-power impulse transmitter (Reference Narod and ClarkeNarod and Clarke, 1994) in conjunction with a resistively loaded transmitting antenna of 4 m half-length, producing a centre frequency of 12 MHz. Geometry of the equipment during data collection is shown in Figure 2. The direct air wave is used to trigger the receiver, and data are transmitted via an optical cable to a battery-powered Fluke Model 97 Scopemeter. A serial interface connects the oscilloscope to an HP200LX palmtop computer that controls data collection, processing and storage.

For case of travel, antennas are oriented parallel to the propagation direction of the air wave and to the given transect. In accessible areas of the glacier, we attempt to collect continuous lines of data generally transverse to the direction of ice flow to ensure good lateral resolution. Where possible, spacing between soundings is either 12.5 or 25.0 m, with surface-elevation survey data collected less frequently (usually every 50.0 m). We implicitly assume that surface elevation gradients are smaller than basal elevation gradients so that surface coordinate data can be linearly interpolated between survey points.

Fig. 2. Geometry and components of the data-collection equipment. The transmitting end comprises the miniature high-power impulse transmitter and power supply between two antennas of 4 m half-length. Separated from the transmitter by 25 m, the receiving antenna is connected to a digital scope-meter and HP palmtop computer. The sounding location is assumed to be vertically beneath the midpoint of the two antennas, hence the position of the survey target.

Picking direct and reflected arrivals

Figure 3 shows a series of typical radar-sounding traces observed on the HP palmtop computer where the direct and reflected waves are indicated. We attempt to pick the time of the first detectable energy in the arrivals, and therefore the incidence of the wave as opposed to the peak. This strategy avoids the potential error introduced with data of variable amplitude. The results of this processing are travel times for the direct and reflected waves for each radar sounding. Ice thickness is readily computed from the travel-time difference, assuming a homogeneous ice velocity and simple parallel-planar geometry of the ice surface and bed in the vicinity of the measurement. Data migration is not necessary in this case because much of Trapridge Glacier flows unconfined over its bed; steeply dipping valley walls, usually the hallmark of glaciated terrain, are not an obvious feature of this basin. Migration techniques have been developed to handle this problem (Reference Welch, Pfeffer, Harper and HumphreyWelch and others, 1998) common to many other alpine glaciers.

Fig. 3. Example radar traces from three locations showing the direct air wave (d) and bed reflections (r). Traces A–C are presented in order of decreasing data quality and increasing uncertainty in arrival picks, demonstrating the need to quantify pick confidence in the analysis. Englacially reflected wave arrivals (e) are not uncommon in areas of complex geometry, in which case we use the signal-to-noise ratio or travel-time comparisons with neighbouring points to determine which arrival is associated with the bed.

For each radar trace, data quality is assessed and recorded by quantifying two attributes; the first pertains to the amplitude of the direct arrival, and the second to the clarity of the time-domain location of the reflected wave. The reflection index is quantitatively dependent on the number of oscilloscope pixels of uncertainly in the arrival location. Depending on oscilloscope parameters, each pixel represents the travel time through 0.68–1.36 m of ice, so an uncertainly of one pixel can introduce errors in the final ice-thickness estimates of 0.035–0.07 m. This index is thus a confidence factor for the pick that we exploit to preferentially value high-quality data in the interpolation.

Data preconditioning

Because the data are relatively dense in some areas and sparse in others, we combine measurements that occur within the same 5 m gridcell. For the ice-thickness data, we compute a weighted average where each datum is favoured according to its reflection quality index. In some cases, very poor data are excluded if they are spatially redundant. For the ice-surface elevation dataset, all values are weighted equally since there is no processing step that introduces additional subjective errors.

Both datasets under consideration, ice-surface elevation and ice thickness, are statistically corrupted by spatial trends that are apparent when either variable is plotted as a function of x, the easting coordinate (approximately down-flow). We anticipate spatial trends in this direction because of the ice-surface slope and the tendency for the glacier to thin toward the terminus, but in practice these trends must be empirically identified. Fortunately, in both cases the variation with x is approximately linear. To precondition the raw surface-elevation data z _s , we search for a best-fit line which, when subtracted from the data, results in a mean approximately equal to zero. Let us denote the transformed surface-elevation data _zs . The result of removing this trend from the data is shown in Figure 4a and b where comparative histograms of z _s and _zs qualitatively demonstrate that the transformed variable is normally distributed about zero.

A two-step transformation is applied to precondition the ice-thickness data. To cluster the data more effectively about a line, we first take the natural logarithm of ice thickness h _i; this also ensures that the interpolation will produce positive estimates for ice thickness, because the recovered data will be the exponential of the interpolated estimate. The second step is to remove the linear trend from the logarithmic data to obtain the transformed variable _hi. Figure 4c and d show histograms for the raw and final transformed data, demonstrating that the two-step transformation not only eliminates the bimodality of the original data but again yields an approximately normal distribution.

Fig. 4. Histograms of raw and transformed data. (a) Raw ice-surface elevation z_s. (b) Transformed ice-surface elevation _zs. (c) Raw ice thickness h_i. (d) Transformed ice thickness _hi. For both datasets the transformations have resulted in a data mean near zero and hare produced approximately normal distributions from bimodal raw data. The bimodal nature of distributions in (a) and (c) reflects distinct environmental and geometric factors associated with the upper and lower basins.

Geostatistical analysis

To quantify the spatial characteristics of our datasets, we use a well-known tool of geostatistical analysis, the experimental semivariogram. A variogram 2 γ(h) for a particular variable Y, considered in one dimension for simplicity, characterizes the variability between values of Y at different spatial locations (Y(x), Y(x + h)), or the autocovariance of this variable as a function of lag h (separation) (Reference Journel and HuijbregtsJournel and Huijbregts, 1978; Reference CarrCarr, 1995). Following Carr, the semivariogram γ(h) = 2 γ(h)/2 is defined mathematically as

(1)

where N is the number of data points under consideration and Y(x _i +h) is the notation for a variable that is a distance h away from Y(x _i ), so h = |x _i+h – x _i |. The semivariogram, presented here for a single variable type (Y), is a measure of similarity within the dataset and relies on the fact that similarity is inversely proportional to distance. Formally, this assumption is known as the intrinsic hypothesis and it states that variable similarity is exclusively a function of lag, not of absolute spatial location.

The utility of this analysis is in identifying the appropriate semivariogram model and quantifying its statistical parameters. Among the most common variogram models are those termed Gaussian and spherical. These are shown in Figure 5 along with visual explanations of their corresponding characteristic parameters. The equation of an idealized Gaussian semivariogram is (Reference Journel and HuijbregtsJournel and Huijbregts, 1978)

(2)

with γ(0) = 0 where CO, C and a are the standard representations of the statistical parameters. The “nugget” CO is the value of the semivariogram when the lag is extrapolated back to zero. Ideally CO = 0 because each data point should correlate exactly with itself, so a non-zero nugget value is indicative of random noise or spatial autocorrelation on a scale smaller than the minimum lag. Parameter C = sill – CO where the sill is defined as the constant value that γ(h) approaches at large lags and is related to the statistical variance of the data. For the Gaussian model, the “range’’ a, or transitional correlation lag, is determined by , where γ(a′) = 0.95(sill). An idealized spherical model is written as (Reference DavidDavid, 1977)

(3)

with γ(0) = 0 for CO and C as previously described, and a defined as the lag at which the semivariogram attains the sill value. For spherical models the range is determined graphically rather than computed.

Fig. 5. Two common idealized variogram models: (a) Gaussian and (b) spherical. Increasing values of γ(h) for small values of h indicate spatial correlation of the data at small lags, diminishing as separation between data pairs increases. Parameters germane to the interpolation are indicated for both models. The sill relates to the statistical variance of the data, the nugget is an indication of data errors (if it is non-zero as shown here) and the range delimits the maximum lag at which correlation is significant. The range corresponds to different values of γ in (a) and (b) because it must be uniquely calculated for each model type in order to represent the same geostatistical property.

Directional semivariograms

Further insight into the nature of spatial autocorrelation of the data can be gained by computing directional semivariograms, which use a subset of the data chosen based on data-pair orientation. Results can reveal the imperfection of simple data transformations in removing small-scale trends, and anisotropy of the transformed data can then he quantified as a further constraint on the interpolation. We compute directional semivariograms for eight angles for each of the two datasets. As expected from the omnidirectional results, nearly all of the directional semivariograms for _zs are fitted well with Gaussian models, and models for _hi are unanimously spherical. The parameters that resulted from fitting models to these semivariograms are listed in Table 2 for both datasets.

Fig. 6. Experimental semivariograms (points) plotted with best-fit models (lines). (a) _zs. A Gaussian model is chosen to fit this semivariogram for ice-surface elevation. (b) _hi. The experimental calculation and model are plotted only for lags less than 600 m due to an insufficient number of ice-thickness data pairs at large lags. A spherical model produces the best fit to the data. Both experimental semivariograms are generated with a minimum lag of 40 m, although a choice of 20 or 60 m yields nearly identical results.

Table 1. Geostatistical parameters extracted from experimental omnidirectional semivariogram analysis of _zs and _hi

Of direct importance is the comparison of the ranges determined from the directional calculations, which gives the orientation and magnitude of the spatial autocorrelation anisotropy. If vectors with orientations parallel to the individual semivariogram angles, and magnitudes proportional to the corresponding ranges, are plotted with common midpoints, this anisotropy is revealed as an ellipse. From the associated angles and ranges for _zs, an ellipse striking approximately northwest–southeast is definitive. It has an aspect ratio of 1.8, meaning that the range is nearly twice as great along its major axis as perpendicular to it. The directional ranges for _hi vary from 150 to 375 m, and the best-fit ellipse yields an anisotropy angle of 45^° west of north with a ratio of 1.4. This orientation is very close to that for _zs, so residual spatial trends are similarly present in both preconditioned datasets. Having quantified the orientation and magnitude of spatial anisotropy, the required geostatistical analysis is complete, and these constraints can now be applied in the interpolation.

Table 2. Geostatistical parameters extracted from directional semivariogram analysis of _zs and _hi. (For a linear model, C is defined as the slope of the semivariogram)

Kriging and post-processing

A published routine (Reference CarrCarr, 1995) is used to krig both datasets with the input information discussed above, and the true recovered data obtained from the raw kriging results by reversing the transformations performed initially. All ice-thickness values are positive as a result of the logarithmic data transformation. Our kriging results are relatively insensitive to the exact values of the geostatistical input parameters but are extremely unreasonable if an unrepresentative variogram model or random parameter values are chosen. Because much of the Trapridge Glacier terminus is an ice cliff, not a gently thinning tongue, we choose to suppress information about the ice extern in the interpolation. Furthermore, introducing values of zero ice thickness into t he dataset disrupts the distribution and requires more complex data transformations. Following this approach we must artificially impose the glacier margin on the interpolated dataset, so we combine survey data and an aerial photograph to generate 25 points that adequately define its position. We then set ice-thickness values to zero that occur outside the margin and specify corresponding ice-surface elevation values equal to the local bed elevation.

Piezometric surface

Water moves in response to potential gradients, so a useful hydrological predictor is a map of the piezometric or potential surface. Such a map cannot be constructed without knowledge of the presence and pressure of water everywhere, but one can predict water-flow paths by assuming a spatially uniform flotation fraction (P _w/P _i) of water pressure. In these calculations it is often assumed that the water pressure is everywhere equal to the ice-overburden pressure, perhaps an acceptable approximation for winter when an efficient subglacial drainage system is not developed. It is well established from borehole water-pressure measurements that summer meltwater can induce excursions in excess of flotation, but this condition cannot be maintained across the entire glacier bed. Therefore we consider two situations: P _w = 0.5 P _i and P _w = P _i, characteristic limits of diurnal summer water-pressure variations.

Total hydraulic potential φ is the sum of pressure potential and elevation potential (Reference ShreveShreve, 1972),

(4)

where P _w is water pressure, ρ _w = 1000 kg m^-3 is the density of water, g = 9.81 m s^-2 is the acceleration due to gravity, and z _b is bed elevation. The water pressure in terms of the ice-overburden pressure is P _w = f ρ _i gh _i, where f = P _w/P _i is the flotation fraction. ρ _i = 917 kg m^-3 is the density of ice, and h _i is ice thickness. Thus we have

(5)

Figure 9 shows the results of this calculation for f = 0.5 and f = 1, visualized as contours of equal hydraulic potential. Contours that appear beyond the glacier margin are expressions of topographic potential only.

Fig. 9. Contours of equal hydraulic potential. (a) P_w/P_i = 0.5. (b) P_w/P_i = 1. Both cases are realistic for behaviours observed during the melt season. Water moves perpendicular to equipotential lines and is channeled into areas where contours are diverted to the west. Obvious outlets for water are the terminal reaches of the central and southern ice lobes. The box indicates our study area.

Both assumptions give similar results, with, f = 0.5 producing slightly more structure, a reflection of the greater influence of bed topography which varies over shorter spatial scales. Best seen for f = 1, the strongest suggestions of water channelization at the bed (indicated as contours diverted to the west) are from the upper basin into the transitional region, and near the terminus. For both f = 0.5 and f = 1, the northwest–southeast-trending potential low emerging from the upper basin appears to relieve most of the upper catchment. This potential trough continues with northern concavity into the lower basin, though it heroines much less distinct. Water flow to the glacier terminus is likely divided beneath the centre and southern ice lobes. We expect a more distributed drainage network for f = 0.5 than for f = 1 due to the slightly greater structure in its potential distribution.

Terrain analysis

Quantitative description of terrain attributes has long been used in geomorphology to characterize drainage in river basins, and is explained along with its history by Reference Zevenbergen and ThorneZevenbergen and Thorne (1987). Glaciological application of terrain analysis is not unprecedented: Reference SharpSharp and others (1993) computed upstream area distributions for Haut Glacier d’Arolla, Switzerland, based on DEMs to substantiate inferences drawn from experimental dye-tracing results; Reference Marshall, Clarke, Dyke and FisherMarshall and others (1996) have used aspects of terrain analysis to determine controls on fast-flowing ice; and Reference Bahr and PeckhamBahr and Peckham (1996) investigated the viability of describing glacier networks using statistical topology models.

For our purposes, the most relevant terrain characteristic described by Zevenbergen and Thorne is upstream area which, for a particular gridcell, is the sum of all gridcell areas that are upstream and connected. In a landscape model with known water thicknesses in each cell, computing upstream area would allow total runoff volume passing through any cell to be ascertained directly. The term “upstream” usually refers to higher topographic elevations, but to adapt this for our purposes, we use hydraulic potential differences rather than elevation gradients alone.

Much labour has gone into creating rules to calculate upstream area accurately (e.g. Reference Costa-Cabral and BurgesCosta-Cabral and Burges, 1994; Reference TarbotonTarboton, 1997). Algorithms of varying sophistication have been developed for these calculations, each attended by certain drawbacks related to grid dependence, artificial dispersion, computational memory requirements and numerical complexity. The simplest method, referred to as D8 or steepest descent (Reference TarbotonTarboton, 1997), uses a transfer rule where one gridcell gives its own area plus its upstream area to a single neighbour with the largest gradient between it and the donor. This is a numerically simple calculation but has an obvious grid bias. In a glaciological context, D8 might be useful for representing a network of discrete conduits, although it still suffers from limited flow-orientation possibilities. Reference Quinn, Beven, Chevallier and PlanchonQuinn and others (1991) improved on this by introducing a multiple-direction method that incorporates weighted area transfer where all downstream neighbours among the eight nearest receive area in proportion to relative gradients. This method largely alleviates the grid bias encountered using D8, even though only eight discrete flow directions are possible, but introduces numerical dispersion and performs poorly at boundaries. Despite this, Reference Costa-Cabral and BurgesCosta-Cabral and Burges (1994) advocate this method over the steepest-descent algorithm. A clever way of approximating continuous-flow-direction possibilities called D∞ was developed by Reference TarbotonTarboton (1997); in this approach each gridcell is divided into triangular facets, and upstream area is partitioned between two receiving cells according to the angle of steepest descent. The method performs similarly to the multiple-direction method with reduced dispersion. Tarboton’s results demonstrate this point: the drainage morphology derived from both is similar, but images resulting from D∞ calculations are sharper. In one comparative test, though, D∞ failed to reproduce proper symmetry that the multiple-direction method captured (Reference TarbotonTarboton, 1997). For our purposes, we choose the multiple-direction method as a satisfactory compromise in that it clearly outperforms D8 and is substantially easier to program than D∞. Because we are interested in qualitative drainage morphology rather than quantitative discharge calculations, this method is adequate. We have performed the same calculations using the D8 method and find that our conclusions remain unaltered, despite the fact that we believe it is an inferior representation of a hydraulically diffusive glacier bed.

Upstream area is computed for a range of flotation fractions by repeated sweeps in four directions, allowing areas to cascade down from gridcells of high potential to those of lower potential. Results for P _w/P _i = 0.5 and P _w/P _i = 1 are shown in Figure 10. DEMs with square 20 m gridcells are used in both cases. Additional tests were done with DEMs of 5 and 10 m resolution to elucidate any grid sensitivity in our method; spatial patterns obtained were visually indistinguishable in all three cases, demonstrating the robustness of this method under our working grid geometry and resolution. In describing these results, we refer to water transport and drainage patterns as interpreted from the realizations of upstream area. This terminology is justified insofar as these results represent water-flow probability distributions. We do not use terms such as discharge, which imply quantitative knowledge, since we cannot, in general, compute flow volumes; this is possible only if a constitutive equation is invoked to relate water-equivalent thickness to subglacial water pressure, a relationship which necessarily depends on the subglacial medium but is not in the spirit of the present geometrical analysis.

Fig. 10. Upstream area distributions governed by hydraulic potential. (a) P_w/P_i = 0.5. (b) P_w/P_i = 1. Plots are presented as logarithms of upstream area to remove the overwhelming downslope trend. Upstream area is possibly an indication of basin drainage morphology, or, for a particular location, a proxy for hydraulic connection probability. Bright areas represent preferential water-flow paths or areas of high connection probability. The box indicating our study area is detailed in Figure 11.

In Figure 10, brightness indicates high values of upstream area, and thus hydraulically favourable pathways in our interpretation. Results are presented as logarithms of the original calculations to remove the overwhelming downslope trend which obscures details. Shaded areas are hydraulically resistive. In the case where P _w/P _i = 0.5 (Fig. 10a), both ice and bed topography exert strong controls over the predicted drainage morphology. The result is a system dendritic in certain areas, diverting flow around relatively small-scale obstacles and spreading out into a uniform distribution in others. This arises because the bed is predisposed to unconfined distributed drainage, rather than channelization. One notable exception is the confluence suggested at the eastern margin where a major outlet develops, the result of a depression between the two bedrock ridges. Important to note is the discontinuity of drainage between the upper and lower basins (west and cast). The potential for ponded water is indicated by light areas surrounded by dark (hydraulic barriers). Such an area is easily visible in the upper basin where water is collected from several tributaries but flow is interrupted by prohibitive topography. These areas have potential glaciological significance as water-storage locations. Without an efficient connection between the upper and lower basins in this case, there is no obvious outlet for water entering these depressions. A further interesting feature of this upstream area distribution is the fact that no significant channelization occurs in the vicinity of the study area, suggesting distributed, inefficient flow if water pressure is low.

For the case of full flotation (P _w/P _i = 1), the upstream area distribution changes significantly. Water is channelized beneath the ice along a down-flow axis with a southeast excursion as it exits the upper basin avoiding a large bedrock obstacle. This axis then bends northeast back into the central area of the lower basin through the study area. Most importantly, it is the ice-surface slope that is sufficient to drive water from the overdeepened upper basin to the lower reaches of the glacier and establish the hydraulic connection that was missing at low water pressures. The qualitative differences between these two cases (high and low water pressures) is in accord with our speculations based on the hydraulic potential distributions (Fig. 9); drainage tends to be distributed and less confined at low water pressures, while high water pressures carry the influence of the ice surface and lead to more channelized morphologies. As the water enters the lower basin where bed and surface slopes are mild, the intense channelization is interrupted by two areas of distributed water. The first interruption diverts some water flow to the southern terminal lobe. This water anastomoses before converging again into a single channel that exits the glacier. The second hydraulic interruption occurs further down-flow in the centre of the basin very near our study area. Here again the flow paths eventually reconverge with a few tributaries between the bedrock ridges and emerge at the extreme end of the centre lobe. These results could be substantiated or refuted if there were subglacial outlets where water discharge could be measured, but because sediments at the margin of Trapridge Glacier are frozen, there is no observational evidence of subglacial drainage issuing from the ice-bed contact. Thus we infer that discharge generally occurs through the subsurface where it is difficult to detect.

Given the purely geometric nature of this analysis, we acknowledge that it is not generally applicable to all subglacial environments. We begin with a simple assumption (uniform water-pressure distribution) that we believe adequately represents a given hydrological state at the bed, but do not address the mechanisms that may have brought about its development. Many processes in an open system are likely to invalidate this assumption. For example, it has been shown that englacial routing can deliver large volumes of surface melt to the bed, and in some cases transport the water over significant horizontal distances before releasing it (Reference Fountain and WalderFountain and Walder, 1998); for a system where englacial transport of surface water is sufficient to maintain spatial pressure gradients at the bed, this type of analysis may be inappropriate. Furthermore, while results in Figure 10b suggest a uniaxial drainage morphology and perhaps a potential location for the development of a conduit, this analysis takes no account of detailed subglacial physics such as is appropriate for Röthlisberger or Nye channels or for canals as proposed by Reference Walder and FowlerWalder and Fowler (1994). Indeed the interaction of a subglacial water sheet and conduits would certainly preclude the possibility of sustaining uniform water pressures at the bed over short time-scales. The channel suggested in Figure 10b is 60–80 m wide and interpreted as a preferential flow path through permeable sediments rather than a channel incised in rock or ice, because the drainage of Trapridge Glacier is better characterized as a diffusive Darcian system than as a network of conduits (Reference StoneStone, 1993). In essence, this analysis is a way of isolating the hydrological effects of ice-surface and bed geometries that will undoubtedly not reflect the true complexity of the drainage network but will nevertheless modulate its development.

Despite these caveats, we attempt to predict changes in connection probability in our study area with water pressure, if upstream area can be used as an adequate proxy for this. Such behaviour has been documented at Trapridge Glacier by Reference Murray and ClarkeMurray and Clarke (1995) who identified areas of the bed that were hydraulically connected during periods of high pressure and became isolated at low pressures. We now examine our borehole-drilling records in an effort to characterize the spatial connection probability of the glacier bed.

Borehole connection record

Hot-water drilling has proceeded at Trapridge Glacier since 1980 and a log has been maintained to record drilling observations, including the connection status of each hole (whether hydraulic communication was established at the base of the glacier, indicated by a rapid drop in borehole water level). We use this information from boreholes drilled between 1989 and 1997 to construct an average connection-probability map of the study area. Figure 11a maps the borehole locations within the study area and indicates the connection status of each. From the intermingling of connected and unconnected holes, the extreme spatial and temporal heterogeneity of bed hydrology is evident. From these data we generate an interpolated connection-probability map (Fig. 11b) using the statistical kriging method we have described in this paper. For comparison with terrain analysis results, we choose 20 m gridcells; such low resolution obscures the sharp heterogeneity of Figure 11a, but indicates a region of high connection probability (bright) in the centre of the study area. This could be due to factors other than ice and bed topography, such as laterally extensive cracks near the base of the ice or continuous presence of permeable sediments.

Fig. 11. Comparison of borehole connection probability and upstream area distributions for the study area. (a) Drilling observations, 1989–97, used to derive connection probabilities, (b) Connection probability distribution for 20 m gridcells obtained by kriging. The box indicates the area of highest connection probability predicted from borehole data, and appears in (c) and (d) for spatial comparison, (c) Upstream area distribution for P_w/P_i = 0.5. (d) Upstream area distribution for P_w/P_i = 1. Note the difference in shading scales in (c) and (d).

Not surprisingly, neither upstream area calculation predicts the detailed connection probability very well. Bearing in mind that the study area is the only location where we have successfully drilled connected holes on Trapridge Glacier, the fact that most of the water flow is predicted to be routed through this area (especially in the case P4w/P _i = 1) is remarkable in itself. Furthermore, we note that the connection “hot spot” in the centre of the region is slightly better predicted by the upstream area distribution for P _w/P _i = 1 than for any other distribution resulting from P _w/P _i = 0.1–0.9, including P _w/P _i = 0.5 shown here. The difference in brightness scales in Figure 11 between panels c and d, chosen to emphasize contrasts within each panel, suppresses the evidence that the high-water-pressure regime focuses a greater fraction of the total flow through the study-area than does the low-pressure regime. Thus absolute values of upstream area in the study region and in the hot spot are much higher when P _w/P _i = 1. Because our drilling records are based on observations made in the midst of the melt season, it is sensible that P _w/P _i = 1 comes closer to capturing this feature than P _w/P _i = 0.5.