Interpretation of Falling-Head Tests in Presence of Random Measurement Error

Field data are tainted by random and several types of systematic errors. The paper presents a review of interpretation methods for falling-head tests. The statistical robustness of each method is then evaluated through the use of synthetic data tainted by random error. Six synthetic datasets are used for this evaluation. Each dataset has an average relative error for water elevation Z, respectively, of 0.04%, 0.11%, 0.22%, 0.34%, 0.45%, and 0.90% (absolute errors on elevation are, respectively, 0.10, 0.25, 0.50, 1.0, and 2.0 mm for a range of water elevation change of 150 mm during test). Each synthetic dataset is composed of 40 synthetic tests (each test consisting of 18 data couples of synthetic falling-head measurements). Results show that the Z-t method is the most accurate and precise, followed by the Hvorslev method when a correction is applied and the velocity method when appropriately interpreted. Advice on how to interpret falling-head tests tainted by random error concludes the study.


Introduction
A number of methods exist to measure saturated hydraulic conductivity of soils. The constant head test [1,2] is one of these methods. In a laboratory setting, water inlet and outlet need to be maintained at constant elevations. Measurements of water flow rate in function of hydraulic gradient (head loss per flow path distance) allow the computation of the hydraulic conductivity of the soil using the following equation: where Q is the flow rate, L is the distance separating the two piezometer measuring tips, h j are measured heads at measuring tips 1 and 2, and A is cross-section of the permeameter. ASTM [1,2] is a good source for permeameter design and testing methodologies. When performing hydraulic conductivity tests in the laboratory, Chapuis et al. [3] stress the importance of good saturation. They describe a testing protocol to reach and measure saturation levels of tested samples. Constant head tests can also be performed in the field. CAN/BNQ [4,5] describe methods where water is injected in a cased borehole at a constant head in an aquifer or aquitard under steady state conditions [6]. Hydraulic conductivity can then be computed by using where C is a shape factor that is, for example, 2.75 times the inner diameter D for an end of casing test, h 1 is the total head of water in the cased borehole, and h 2 is the total head of the free surface (or piezometric level) of the tested soil layer. This method is based on the hypothesis that the injected volume of water will have negligible influence on the piezometric level (PL) of the soil surrounding the injection zone. When testing low hydraulic conductivity soils, the main drawbacks of constant head tests are that they are extremely lengthy in time, and, in the case of field tests, the PL of the soil layer may be unknown or difficult to measure. Fallinghead tests are effective answers to both these drawbacks. They can be completed in a short-time span and Chapuis et al. [7] demonstrated that their interpretation does not necessitate the PL of the surrounding soil layer.
Field (and laboratory) data are generally tainted with measurement errors. These are of two types, systematic and random errors. A systematic error is, for example, created to the total head by incorrect calibration of a piezometer. This type of error can be constant in time or can follow some sort of time drift. Random errors are characterized by a zero mean and some standard deviation. These errors are the result of a number of phenomena. Human reading of a standpipe scale and white noise due to electrical interference in piezoelectric transducers are examples. Piezoelectric transducers also yield step-wise data due to conversion of analog values to digital code. Measurements are then higher or lower than true values. This paper will review the theoretical background for the interpretation of falling-head tests when in presence of random errors. Five (5) interpretation methods are presented: log (or Hvorslev), traditional interpretation of velocity (or Chapuis), alternate interpretation of velocity, Z-t and optimised log (or corrected Hvorslev). Each method will be applied on synthetic datasets. Data within each set includes a random error component of predefined variance. The object will be to evaluate how each interpretation method is sensitive to random measurement error. A comparison between methods and a discussion on the interpretation of falling-head tests will follow.

Interpretation of Falling-Head Tests
Data from falling head tests, be they obtained in the laboratory or in the field, may be interpreted by a number of methods. Chapuis [8] indicates that when the deformations of the soil can be neglected, falling-head tests are governed by the Laplace equation. Its solutions, the harmonic functions, have several properties. One of them relates flux into the soil (Q soil ) to flow into the pipe (Q inj ) through a mass-balance equation: where c is a shape factor that depends on the geometry of the injection zone and on the hydraulic boundaries of the problem, H is the applied hydraulic head difference, and k is the hydraulic conductivity. This equation is the starting point of the Hvorslev, velocity (Chapuis), and Z-t methods.
Another equation is the starting point of another method for cases where the soil deformation is assumed to be elastic and not negligible [9]. However, this method contains physical and mathematical confusions according to the mathematical, physical, and numerical proofs by Chapuis [8] and the experimental proofs of Chapuis and Chenaf [10]. According to the equations of Chapuis [8], the effect of soil deformation can be neglected when the soil is an aquifer or an overconsolidated aquitard. It is no longer negligible for compressible aquitards when they are tested using either a falling-head test with a very small injection pipe or a pulse test between packers. Chapuis and Cazaux [11] gave suggestions on how to correctly handle the instantaneous (elastic) and delayed deformations in such cases. In a falling-head test, Q inj is the flow through the inflow pipe (often a standpipe connected to the borehole casing) of internal cross-section S inj : where t is time. Equations (3) and (4) yield Rearranging gives Integrating leads to the solution proposed by Hvorslev [12]: where H j is the instantaneous head loss at time t j obtained from the difference between the total heads in the inlet standpipe and PL at the boundary of the surrounding soil, and C = c/S inj is a shape factor that depends on the inlet/outlet geometry. Hvorslev [12] compiled shape factors for a number of inlet/outlet geometries. For field conditions illustrated in Figure 1, shape factor C is where D and d are, respectively, the internal diameters of the casing and of the standpipe. For laboratory falling-head tests, C is given by where A is the area of the permeameter cross-section, "a" is the cross-section of the inlet standpipe, and L is the flow distance through the soil sample. In the laboratory, the outlet is a constant head basin controlled by an overflow weir. In the field, the outlet head is set by the geometry and the boundary conditions of the surrounding soil, that is, the PL. In the laboratory, it is relatively easy to accurately and precisely measure the outlet piezometric head. In the field, an accurate measure of the PL of the surrounding soil is not necessarily trivial. An error in the PL value will introduce a bias or systematic error of constant value in head difference H j of (7). Chapuis [13] showed that this systematic error will produce a curved ln(H) versus time plot. Since the hydraulic conductivity k is the slope of this plot, a curved plot wrongly suggests that k changes with time. A concave downward curve suggests that k increases during the duration of the test. If the test is repeated assuming the same systematic error in the PL of the surrounding soil, the same curved plot will be observed. Chiasson [14] mathematically demonstrated that when the systematic error in the PL of the surrounding soil is not zero (i.e., H 0 / = 0), the relationship between hydraulic conductivity and the slope of the plot as expressed by (7) is no longer valid. Thus the direct use of (7) as proposed by Hvorslev without questioning the PL value is faulty practice [13]. Chiasson [15], furthermore, concludes that Hvorslev's equation is incomplete and proposes a correction to this method and an interpretation through log optimization. This is described later in the paper.

Traditional Interpretation of the Velocity Plot.
The velocity method proposed by Chapuis et al. [7] is one where the unknown PL of the soil is not needed for hydraulic conductivity determination. Therefore, a systematic error in the assumed PL has no consequence in computations ( [16] as quoted by [13]). For this method, the following definition is first introduced: where H 0 is the vertical distance between Z(t), the elevation above ground of the inlet falling water level within the standpipe, and the PL of the surrounding soil ( Figure 1). If Z(t) is erroneously assumed as the total head loss between inlet and the PL of the soil, H 0 can be seen as the systematic error, or bias. Rearranging (5) using (10) then gives By plotting y = Z as a function of x = v = dZ/dt, a straight line should be obtained with slope m v = −1/Ck and intercept H 0 . Thus, the slope of this plot is related to k by Choosing Y = Z and X = v defines the traditional velocity plot. Modern spreadsheets permit to easily compute the slope m v using the least squares method, a best linear unbiased estimator (BLUE).
The velocity plot has many advantages. A part from estimating hydraulic conductivity, it can be used to graphically obtain the systematic error (if any) of the assumed PL, identify hydraulic fracturing during the test ( [7,17,18] in piezometers) or consolidation at the beginning of the test of the soil surrounding the tip of the borehole [8,19], and even identify leakage in buried pipes [20].

Alternate Interpretation of the Velocity Plot.
Interpretation of the velocity method using the method of least squares should yield the best unbiased linear fit. Chiasson [15] underlined that statistical estimation by least squares is theoretically based on one dependant variable being a function of another that is independent. The independent variable is a controlled variable; that is, it is the user that decides at which value a measurement of the dependant variable will be made. Thus, by definition, the independent variable has no measurement error.
In the traditional interpretation of the velocity method, velocity v during a time increment Δt is considered as the independent variable and the average elevation Z m = (Z j + Z j+1 )/2 during the same time increment is the dependent variable. When data has low scatter in Z (or t), this has little effect on the result. When data has some scatter in Z (or t), Chiasson [14,15] shows that interpretation problems arise! The act of choosing v as the controlled variable when it has high statistical scatter clearly departs from least square estimation theory.
Between variables v and Z m , Z m displays the least measurement error. One could thus choose to consider Z as the control variable and time t could be measured at a certain value of Z. With Z m as the control variable, velocity v = ΔZ/Δt is the dependant variable. This makes more sense since by definition velocity v is a function of elevation Z and time t, that is, it is dependent on Z and t. Rearranging (11) to isolate v on the left hand side gives From (13), the hydraulic conductivity k and the bias H 0 are:

Z-t Method.
A third method proposed by Chiasson [14] plots raw elevation data Z as a function of time t. The solution for (6) using (10) is as follows (see Chiasson for demonstration): where H i is the head difference at initial time t j = t 0 = 0 and "a" is an exponent parameter where Let Z j be the measurement of Z at time t j , for j = {0, . . . n}, and let Z * (t j ) be the estimated water level in the standpipe at time t j , using estimated parameters H i * , a * , and H * 0 . The best unbiased estimator will then be obtained by numerically minimising the following equation: In doing this, Hvorslev made two implicit suppositions that H 0 is known a priori and that the initial reading for H i at time t 0 = 0 has no measurement error. Clearly, this can never be perfectly the case. These two assumptions will only be approximately true, the extent of the approximation being a function of error. Chapuis [13] and Chiasson [14] show how making these implicit suppositions will adversely affect the interpretation of the test and the value of the hydraulic conductivity k. The interpretation process could use the velocity method to determine H 0 and interpret the fallinghead test by rearranging (7) using (10) to obtain and to plot y = ln(Z j + H 0 ) in function of x = t. Unfortunately, this approach is incomplete since it implicitly makes the supposition that the initial reading at time t 0 = 0 has no measurement error. Since the error on the first reading (16) is approximately equal to Z 0 + H 0 . (where Z 0 is measured and H 0 is estimated by the velocity method) Thus, the hydraulic conductivity computed this way will generally not be adversely affected. This interpretation is found to yield hydraulic conductivity values that are always close to being equal to the value obtained by the velocity method.
Since two methods give the same result, one could conclude that this confirms the validity of the hydraulic conductivity value. Chiasson [14] shows that if H 0 and H i are obtained by another interpretation method, that is, the Z-t method, interpretation of a falling-head test using (20) will give values equal to those obtained by the Z-t method! This brings Chiasson to conclude that (20) (Hvorslev's method) with H 0 estimated from another method (velocity or Z-t) cannot be used to confirm the validity of the k value obtained by the same other method (velocity or Z-t), since the value of k that is obtained is dependant of the method used to estimate H 0 . A remedy to this is to correct Hvorslev's method by estimating H 0 , H i , and k by a least squares optimisation technique similar to the one used in the Z-t method. Rewriting (20) with Let then y * = ln(Z * (t j ) + H 0 ) be the estimated natural logarithm of the total head in the standpipe at time t j , with estimated parameters b * ln = ln(H * i ), m * ln = −k * C, and H * o . The best unbiased estimator will then be obtained by numerically minimising the following equation: while being subjected to the unbiased condition n j=0 y j − y * t j = 0.
This variant of Hvorslev's log method independently estimates all unknown parameters. Theoretically, it can be used to separately evaluate k and confirm values estimated by both velocity and Z-t methods.

Performance of Interpretation Methods or Their Statistical Robustness
The synthetic tests presented here correspond to conditions encountered in overconsolidated aquitards such as compacted clay liners (CCLs) used in solid waste cells. These have hydraulic conductivities in the order of 1 × 10 −6 mm/s with a typical thickness of one metre. CCL are also very rigid in nature due to the applied high energy compaction. Both high rigidity and limited thickness minimize effects on the local PL of immediate and delayed aquitard deformation [21]. Furthermore, the velocity plot will yield a straight line when pore water volume change is negligible [8] and straight velocity plots are observed in tests performed on CCLs (see, e.g., [14]). The solution for differential equation (6) expressed by (16) with H i = 701 mm, H 0 = 400 mm, and k = 5×10 −7 mm/s and a test setup giving a shape factor C = 235.3 mm −1 was used to generate synthetic data as follows: where ε is a normal law distributed random fluctuation with zero mean and standard deviation σ = ΔZ/1.96. The random fluctuation component corresponds to a synthetic random measurement error. By definition, ΔZ is the absolute error of the synthetic dataset composed of measurements Z. Six synthetic datasets of increasing absolute error on Z were generated this way. Absolute errors for each synthetic dataset will be in increasing order ΔZ =: ±0. 10 [14]). Each method reviewed earlier is applied to these six synthetic datasets to investigate their sensitivity to data tainted by random measurement error.

Traditional Velocity Plot.
As underlined earlier, when measurement errors on Z are small, interpretation by least squares of the traditional velocity plot will yield good results. An illustration of this statement is given in Figure 2. In this  When relative measurement error on Z is significant (in the order of ΔZ/[Z max − Z min ] = 0.7%), the velocity plot displays considerable scatter (Figure 3). Such a plot will suggest that data from the test are of questionable quality. With the same shape factor and (12), one finds by least squares m v = −2643.6 sec −1 giving a hydraulic conductivity k = 1.61 × 10 −6 mm/s and H 0 = −28.8 mm. This is a 222% relative error on k and −107% relative error on H 0 . Hence, relatively small measurement errors in elevation Z yield appreciable estimation errors for k and H o .
Chiasson [15] observes that interpretation by least squares of the traditional velocity plot systematically gives higher k values than the methods earlier presented. The six synthetic datasets, generated by (24) with ΔZ: ±0.10 mm, ±0.25 mm, ±0.50 mm, ±0.75 mm, ±1.0 mm, and ±2.0 mm, confirm this systematic bias ( Figure 4). Also, as absolute measurement error increases on Z, so does the average and statistical scatter of the hydraulic conductivity obtained by the velocity method. Thus, this method of interpretation, for data subjected to measurement error, produces a systematic bias in the interpreted k value.
There is a clear trend between measurement error on Z and the hydraulic conductivity obtained from least squares of the traditional velocity plot (average trend in Figure 4). Furthermore, a good correlation is observed between the relative error (Δk/k) on the value of k obtained by the velocity plot and the coefficient of determination of the same  plot ( Figure 5). This further demonstrates that interpretation by least squares of the traditional velocity plot is not statistically robust. Based on this work and on Chiasson [14], least squares of the traditional velocity plot must not be used in a straightforward manner on raw data when scatter is observed. That is when R 2 is less than 0.92. Coefficients of determination higher than this threshold will yield relative errors on k below 10% ( Figure 5). The above analysis shows that increasing statistical scatter in the velocity plot decreases the slope of the least squares line (and increases the interpreted k value). The value of intercept H 0 will as a consequence decrease in a concurrent fashion. This result for the estimated slope value is in accordance with a relationship between the true slope and the slope estimated by least squares as presented by Snedecor and ([22] unfortunately, the coefficient of the relationship is often difficult to determine). Thus, least squares of the traditional velocity plot must not be used to extract the correct (or true) piezometric level at the injection point when the coefficient of determination R 2 is below 0.92.
Velocity plot scatter can be reduced in cases where the number of data couples is sufficiently high. In such cases, Chapuis [21] proposes a five-step interpretation procedure to reach this goal. When datasets are too small, least squares of the traditional velocity plot must not be used to extract k and H 0 . Nevertheless, It remains warranted to prepare a velocity plot to confirm, as earlier underlined, that it is linear and that no other phenomena intervenes during the test.

Alternate Interpretation of the Velocity
Plot. An appropriate interpretation of the velocity method should be one that uses as the X variable the measurement with minimal error, reserving the role of Y to the other. The same data that was plotted in Figure 3 is used to illustrate this alternative interpretation. The permuted plot, where velocity v is in function of average elevation of falling-head Z m , does not yield a better correlation coefficient, since it measures scatter and this is unchanged ( Figure 6). The hydraulic conductivity computed from this plot is on the other hand significantly different (it is lower). Using (14), the same shape factor C and the least square slope m Z (Figure 6) yields  With the permuted plot, relative error on k has decreased to −11.6% and to 21.0% for relative error on H 0 . This is a considerable improvement in relation to the 222% relative error for k and −107% relative error for H 0 that was obtained earlier with least squares of the traditional velocity plot. Computed hydraulic conductivity values from this "alternate interpretation" of the velocity method show no correlation with test data scatter as characterised by the coefficient of determination (Figure 7). On average, the alternate interpretation of the velocity method will yield good hydraulic conductivity values. On the other hand, scatter in measurements, although less problematic than with the traditional interpretation of the velocity method, will still yield rather high relative errors (Δk/k) for k (Figure 8). The alternate interpretation can thus be qualified as being unbiased, that is, on average accurate but not precise. Results from this study indicate that the alternate interpretation should be used with caution when data scatter in the velocity plot yields R 2 less than 0.87. Otherwise, relative error on the k value may be greater than ±10% (Figure 7). The same observations apply for intercept H 0 (Figure 9). When R 2 is less than 0.87, it is recommended to follow the five steps interpretation procedure proposed by Chapuis [21] but using least squares on X-Y variables of the alternate velocity plot.

Z-t Method.
Chiasson [14] introduced this method after observing that scatter in the velocity plot (likewise with the permuted velocity plot) is inherent to the computation of the velocity ( [21] gives the equation of the relative error). Chiasson thus proposes to use raw data, that is, [t j , Z j ] data couples, and directly plot them on a Z-t graph. Difference  in scatter amplitude is evident when comparing scatter in a velocity plot (or permuted velocity) with scatter in a Z-t plot (compare Figures 3 and 6 with Figure 10).
Applying the Z-t method to the same dataset earlier used with traditional and alternate interpretations of the velocity method gives k = 4.68 × 10 −7 mm/s, H 0 = 439.5 mm and H i = 739.9 mm. This corresponds to a relative error for k of −6.4% and of 9.9% for relative error on H 0 . This is an improvement comparative to earlier presented methods, be it the traditional or alternate interpretation of the velocity plot.  The Z-t method is unbiased; that is, there is no significant correlation with scatter intensity (Figure 11). It yielded high coefficients of determination for the complete suite of studied absolute errors, meaning that (16) well explains the relationship between Z and t.
The Z-t method yields good hydraulic conductivity values and it is less sensitive to data scatter (Figures 11 and 12).  It is thus an accurate interpretation method. It is also a precise method since relative error on k is of only 3.6% when synthetic data has ΔZ = 0.25 mm and increases by approximately the same increment for each 0.25 mm increment to ΔZ. Results from this study indicate that the Zt method can be used even when absolute errors on Z are of the order of 2.0 mm (ΔZ/Z = 1.3%) for which relative error on k will be of 28.4%.
The same behavior is also observed for intercept values H 0 , that is, no systematic error (accurate value). Relative errors ΔH 0 /H 0 are smaller than in both earlier discussed methods although they are 1.6 times higher than hydraulic conductivity relative errors Δk/k. Thus, the Z-t method gives an accurate value for intercept H 0 and gives the highest precision.

Corrected Hvorslev's Method or Optimised Log [Z + H 0 ]
Method. Again, for the same dataset used for the other three methods, the corrected Hvorslev method gives k = 4.40 × 10 −7 mm/s, H 0 = 482.6 mm, and H i = 783.1 mm. This result for k is approximately equivalent in accuracy to the Z-t value obtained earlier. Investigating the complete dataset of synthetic data, results for k show some sensitivity to data scatter, more than in the case of the Z-t method ( Figure 13). The method thus appears acceptable when ΔZ < 1.0 mm. Higher absolute errors in Z appear to destabilise the method. This is due to numerical optimisation which is less stable for this method.

Calibration Error
In cases where the elevation of falling head tests is monitored with a piezoelectric transducer, a bias will be introduced if it is faultily calibrated (error in calibration is a type of systematic error). Piezoelectric transducers are linear voltage devices that convert voltage into a corresponding water pressure. In such case, the relationship between faultily calibrated readings and true elevations is where Z and Z r are, respectively, the true and faulty calibrated (biased) elevation of the water column, Z B and m B This defines the relationship between true velocity (at left) and faultily calibrated velocity (at right): According to (11), using (25), the relationship between faulty velocity and faulty elevation of the water column is: This later equation shows that the plot of faulty elevation in function of faulty velocity has the same slope as one obtained from true elevation and velocity, that is, (11). Then, according to (11), the faulty plot will give the correct hydraulic conductivity. This is not the case for the intercept value as can be viewed by comparing (28) and (11). Since the intercept defines the piezometric level (Figure 1), the plot of faulty calibrated readings will not yield a correct PL.
Since the same equations are used to prepare the alternate velocity plot, the same conclusions will also apply. Developing (15) with (25) will find the same results for the Z-t method and (21) with (25) will yield the same for the corrected Hvorslev method.
When the PL is of importance, and it is certainly often the case, a verification of calibration needs to be performed. If when performing the falling head test, the elevation of the measuring tip of the piezoelectric transducer is noted, it is then possible to back calculate the PL value after calibration of the transducer in the laboratory. For the case of piezoelectric transducers in open piezometers, a field comparison between transducer readings and direct measurements of water elevation will permit to obtain correct calibration.

Conclusion
Interpretation methods for falling-head tests were evaluated for their sensitivity to measurement error in elevations measurements of the falling water column. Least squares of the traditional velocity plot are found to be the less appropriate approach for interpreting falling-head tests to evaluate hydraulic conductivity, even when measurement errors are relatively small. It tends to systematically overestimate the true hydraulic conductivity (i.e., the soil is interpreted as being considerably more permeable than actual). The introduction of an error on the inlet stand pipe elevation Z systematically has a greater impact on falling head velocity v than on Z m = (Z j + Z j+1 )/2 [14,15]. This thereby always increases the scatter range of v values more than the scatter range of Z m values. As a consequence, the slope of the velocity graph will flatten, lowering the slope and thus yielding the interpretation of a higher k value.
A simple corrective measure is to assign the measurement with the less error as the independent X variable, and the other as the dependant variable Y . This yields a more appropriate interpretation method (which is just a permutation between independent and dependent variables used in the traditional velocity plot), where Z m is assigned to the X axis and v to the Y axis. This interpretation is found to be accurate (i.e., on average, it yields the correct k). The alternate interpretation of the velocity method also displays considerably less scatter in computed k values. It is thus recommended to discontinue the use of the traditional velocity plot and replace it by the more statistically robust alternate interpretation of the velocity method and its permuted plot. In cases where the coefficient of determination R 2 is less than 0.87, it is recommended to follow the five steps procedure proposed by Chapuis [21]. This procedure reduces statistical scatter of the velocity plot (traditional or alternate) and thus increases the coefficient of determination. Least squares of the alternate velocity plot can then be applied to obtain more accurate hydraulic conductivity and piezometric level of the surrounding soil.
The best method is found to be Z-t, with the corrected Hvorslev trailing not too far behind. Both these methods are accurate and the Z-t is particularly precise when compared to the other studied methods. It is important to underline here that the Z-t method has some weaknesses. It cannot identify a number of phenomena that may develop during a test. These phenomena can be identified with a velocity plot [7,8,17,19,20]. It is also stressed that the Z-t method is only valid when the velocity plot follows a linear trend. Thus, it should always be used in conjunction with a velocity plot (preferably the alternate one).
By using more than one method, it is possible to better evaluate the accurateness and precision of interpreted hydraulic conductivity values. As a rule of thumb, if the difference between the alternate interpretation of the velocity plot and the Z-t method is small (and if the velocity plot follows a linear trend), the Z-t value can be considered accurate and precise. If absolute errors on Z are less than ±1.0 mm (error relative to the range of elevation change ΔZ/[Z max − Z min ] = 0.7%), this study shows that it may be concluded that computed hydraulic conductivity values are accurate and precise.