Calibration of Gardner coefficient for density–velocity relationships of tertiary sediments in Niger Delta Basin

The relationships between density and velocity are important in many geological analyses that involve rock property parameters. However, the accuracy of the results is often limited when there is a verbatim application of generalized rock property relationships. Many studies in the Niger Delta Basin suffer from this limitation. As a way forward, well logs in the Niger Delta were used to derive lithology-specific coefficients that can be applied in density–velocity transforms that make use of the Gardner equation. Whereas the default coefficient (α) and exponent (β) in the original Gardner equation are 0.31 and 0.25, respectively, fitting the Gardner curve to local data results in a coefficient value of 0.33 and 0.29 for shales and sands, respectively, when the exponent is kept at the default value. Comparing measured density data with estimates obtained from sonic velocities using the original Gardner equation gives a regional mean absolute deviation of 0.13 g/cc while those of the newly derived local coefficients do not exceed 0.05 g/cc giving an improvement of over 60% in the accuracy of estimated rock properties.


Introduction
Many geological interpretations and analyses involve the conventional relationships between density and velocity of rocks. As a rock property, density is required in the identification of lithologies, pore fluids, estimation of porosity, overburden stress estimation and some methods of pore pressure prediction. This makes it a significant requirement in exploration, reservoir characterization and well planning. In seismic analysis, density is needed in the determination of the elastic impedance of an interface and in the improvement of the accuracy of synthetic seismograms. While it is possible to obtain predrill velocity data using seismic technology, density data is hardly obtained along a well path until the well is drilled. As a result of this, density is commonly estimated from seismic velocities using a set of relationships such as Birch (1961), Gardner et al. (1974) and Lindseth (1979). Some of these methodologies are discussed in more detail in Domenico (1984), Castagna et al. (1985) and Han et al. (1986). In particular, Birch (1961) presented a relationship between seismic velocities, densities and atomic mass (an index of rock composition). Following several modifications (e.g. Simmons 1964;Anderson 1967) to the Birch relationships, Gardner et al. (1974) subsequently correlated the variation of bulk density with P-wave velocities for common sedimentary rocks sourced from a wide range of basins, ages and depths from which a density-velocity relationship was developed. The Gardner et al. (1974) relation is: where q = bulk density given in g/cm 3 , V = P-wave velocity, a = 0.31 for V (m/s) and 0.23 for V(ft/s) while b = 0.25.
Although the Gardner calibration considered a variety of sedimentary rocks, the final model is approximately a midline that averaged several lithologies including shale, sandstone, limestone and dolomite (Fig. 1).
Thus, it is common to observe a considerable scattering of real data around the Gardner curve in many plots. This generalization is often a source of errors in some empirical studies (e.g. pore pressure analysis, reservoir characterization, top seal identification and geomechanical studies) that require to work with lithology-specific functions.
In another study, Krasovsky (1981) showed that the Gardner et al. (1974) equation is non-unique for certain datasets from many sedimentary basins across the world as there are several cases where the density-velocity data would not conform to the original Gardner curve but rather deviate consistently from the general rule. The implication of all these is the systematic miscalculation of densities from velocities and vice versa. Therefore, local calibrations are required to improve the accuracy of rock property predictions especially for different rock types and sedimentary basins (such as the Niger Delta) that probably were not captured in the original Gardner calibration.

Geological setting
Prolific oil and gas activities have led to detailed studies of the geology of the Niger Delta Basin and these abound in the literature (e.g. Short and Stäuble 1967;Knox and Omatsola 1988;Doust and Omatsola 1990;Nwajide 2013). Specifically, Nwajide (2013) discussed the geology of the delta as an integral part of the Benue Trough and other genetically related sedimentary basins in Nigeria. According to this author, the Tertiary Niger Delta that developed as a prograding extensional complex covers a surface area of over 100,000 km 2 and is composed of an overall regressive clastic sequence that is as much as 12 kilometres in thickness. From the northern area of the delta towards the south, the sediments successively range in age from Eocene to Quaternary (Fig. 2). The overall configuration of the basin suggests that the delta prograded southwestward, successively forming depobelts that represent the most active portion of the delta at each stage of its development (Doust and Omatsola 1990). These depobelts include the Northern Delta, Greater Ughelli, Central and Coastal swamps as well as the Offshore and Distal Offshore.
The stratigraphic succession (Fig. 3) comprises three major lithofacies that are commonly distinguished based on their sand-shale ratios. At the upper part is the Benin Formation which is a massive continental deposit of alluvial and coastal plain sands, with thicknesses often in excess of 2000 m (Avbovbo 1978). This is underlain by the Agbada Formation which is an alternation of paralic sands, shales and clays including the Opuama Channel Complex and Qua Iboe Collapse System (Nwajide 2013). In many  Doust and Omatsola 1990). Inset is the delta outbuilding from the Eocene Northern Delta towards the Ocean parts of the delta, the thickness of the sand-dominated Agbada Formation can be as much as 3700 m (Avbovbo 1978). The lowermost sediments of Agbada Formation grade downwards into the more deeply buried thick marine shale sequences, clays, silts and often turbidite sands of the Akata Formation that were deposited in shallow marine shelf to deepwater environments (Nwajide 2013). The thickness of the shale-dominated Akata Formation is estimated to be over 7000 m (Doust and Omatsola 1990). The formation is predominantly overpressured and is the major hydrocarbon source rock in the basin while the paralic shore face sands of the Agbada Formation are the major reservoirs. The sand and shale strata are so distinct that they can be significantly correlated across many wells and seismic sections. Considering this distinctiveness and dominance, the unique properties of the two major rock types need to be understood for proper evaluation of lithospecific parameters.

Materials and methods
A suite of well logs that included gamma ray (GR), sonic velocity (Vp), bulk density (Rho), neutron porosity (NPHI), density porosity (DPHI) and caliper logs were sourced from 13 wells. The wells were chosen because they have good quality log data and wide coverage of the basin, making it possible for each of the depobelts to be represented in the study. The logs were quality checked with caliper logs for consistency of wellbore diameter. This makes it possible to identify and edit out sections of the logs that were compromised by borehole rugosity. Log headers were checked for depth datum (i.e. whether the quoted depths were relative to Kelly Bushing, subsea or along hole, etc.) and further verified from composite logs and well files. Where necessary, conversions were made from any of the differing depth references to the true vertical depth subsea (TVDSS) format that was uniformly used in the study. Because of the need to work with consistent sand and shale properties, the top and base of main sand and shale units were also delineated from the composite logs and well reports. Considering that the average geothermal gradient of the basin is 25°C/km, data were taken from well sections that were not expected to be greatly affected by diagenetic mineralogical alterations and late overpressure processes which are mostly active when temperatures tend to exceed 70°C. Moreover, velocity and density logs of shales were visually examined for consistency of mechanical compaction behaviour (i.e. a general trend of similar density and velocity behaviours with depth).

Lithology discrimination
Log data had to be isolated for shales and sands that are of particular interest to the study. Routine identification of lithologies often relies on values of gamma ray logs. Typically, values of gamma ray across mudstones are higher than in other sedimentary rocks because of the high potassium content (and presence of adsorbed uranium and thorium) of clay minerals. Hence, shales can often be distinguished from other lithologies using a simple gamma ray cut-off. In the Niger Delta, this cut-off is usually 65 API such that gamma ray values in excess of this threshold are assigned to shales while sands are assigned to lower values. However, high gamma ray sands also exist in the Niger Delta and could be wrongly picked as mudstones. Picking the cleanest sand and shale units therefore requires a combination of tools.

Determination of Katahara-type sand/shale cut-off
Katahara, (2008) method uses differential porosities (i.e. the difference between neutron porosity, NPHI and porosity derived from density, DPHI) to discriminate shales from other lithologies. This is because bulk density log is a function of the total density of both rock matrix and pore fluids and can be used to calculate a rock porosity with assumed water and grain densities. On the other hand, the neutron logs respond to hydrogen ion concentration in the rocks, which derives from pore-held fluids (PHF) and claybound water (CBW). Porosities can be calculated from neutron logs by assuming that all the hydrogen ions detected are contained in the pore water. When plotted, porosities from both tools overlap across clean fluid-filled sands but become separated from each other across shales (Fig. 4). If the neutron tool is not understood, the additional CBW it detects could be erroneously attributed to pore volume. Katahara (2008) used this approach to discriminate between sands and shales by plotting the porosity difference from these two logs against density. The method uses both the change in slope and density difference to distinguish shales from other lithologies. Following this procedure, wireline log data were used to determine an appropriate shale cut-off on a typical Katahara plot of Niger Delta data. Figure 5 is a cross-plot of bulk density (q) data and the difference between neutron porosity (NPHI) and density porosity (DPHI). DPHI was calculated from density logs thus: where q ma = assumed matrix density = 2.7 g/cc; q bulk = value in density log; q f = average water density = 1.03 g/cc.
Results obtained by this method were checked against gamma ray logs and lithology descriptions in the composite logs and well files. The identified shale intervals correspond to high gamma ray. However, in the sand section of the well, only the Katahara (2008) method could discriminate the lithology as the gamma ray indicated high values that could be falsely taken as shales. Based on the slopechange threshold of 0.2 established in Fig. 5, density and velocity logs were filtered for the two lithologies of sands and shales (using NPHI-DPHI C0.2 for shales and NPHI-DPHI B0.2 for sands). An example of the discriminated data is shown in a Gardner et al. (1974) type of densityvelocity cross-plot (Fig. 6). The sands are coloured yellow, and grey mass is the shale. Fig. 4 Plot of gamma ray, DPHI and NPHI logs in a succession of sands (orange) and shales (green) in the Niger Delta. Gamma ray values along some of the sand packages are same as those in the shales Fig. 5 Determination of sand/shale cut-off in the Niger Delta using a cross-plot of density and NPHI-DPHI. Data are coloured according to their gamma ray values. The slope change between sand-rich and clay-rich is on 0.2. Gamma ray, however, widely varied in the sandrich part of the curve Fig. 6 Example plot of distinct sand (yellow) and shale (grey) on a Gardner type of density-velocity log plot As expected of normally compacted sediments, the plot in Fig. 6 shows that velocity increases with density. However, the sands generally trend lower densities than shales of corresponding velocity. This is because the mineral composition of the two rock types is not the same and shales are more responsive to compaction than sands. The discrete responses of the two rock types on cross-plots further highlights the differences in their elastic properties and will therefore negate many generalized rock property relationships such as Gardner et al. (1974). This suggests that rock property correlations will be improved when specific rock types are accounted for.

Determination of lithology-specific coefficients
To examine the deviation effect of the lithological differences, the original Gardner curve will be compared with data. From this comparison, the default Gardner coefficient will be perturbed to obtain unique values for the two rock types. Modifying the coefficient in this way helps to preserve the original form of the Gardner equation for a better correlation of subsequent results. The coefficients that bestfit velocity-density relationships for the shales and sands were determined by two methods: 1. Adjusting the coefficients for visual best fitting curves through the dataset. 2. Adjusting the coefficients to minimize mean errors from estimated densities.

Visual curve fitting through sand and shale
This involved plotting two sets of curves on a cross-plot of density and velocity of sands and shales using the default Gardner coefficient value of 0.31. Then adjusting this coefficient individually for each of the curves until one of the curves visually makes a best fit with the shales and the other the sand. For the shales, the coefficient needed to be increased to between 0.31 and 0.35 while the sands required the value to be decreased from 0.31 to between 0.30 and 0.27. In Fig. 7, the red curves are fits through sand, grey dashed line is the default Gardner curve while the green curve is the fit through shales. The plots show that the original Gardner line is at best, the lower bound for the sands and the upper bound for the shales.

Minimizing mean errors in estimated density
This involves predicting density from velocity separately for shale and sand using the default Gardner equation. Then, the coefficient in the prediction equation is continuously adjusted until the lowest possible error is achieved. The errors are calculated as the average of the differences between measured density data and the predicted. The example in Fig. 8 shows that the original Gardner (red plot) systematically underestimates the density hence the magnitudes of the absolute errors are greater compared to the local fits (green plot). Figure 9 shows that in the presented case, this difference can on the average be as much as 0.15 g/cc with the original Gardner relationship.
Using the described workflow, the mean absolute errors for both default Gardner coefficient and the local best fit for sands and shales in each of the wells are summarized in Table 1. From these individual well results, regional mean values were then calculated.
Previously in Fig. 7, it was shown that the derived density-velocity curves for shales and sands were quite different from the original Gardner line. The fitting coefficient decreases towards the sand and increases towards the shales. Consequently, the errors in estimated densities were minimized with the data-specific constraining of the coefficients (Table 1). For the shales, the best-fitting coefficients range from 0.31 to upwards of 0.34 while in the sands these range from 0.28 to 0.3. This results in a calculated regional mean coefficient values of 0.33 and 0.29 for shale and sand, respectively. These values would suggest a new set of distinct Gardner type of relationships. In the final model shown in Fig. 10, the sand curve which is a sinistral shift from the original Gardner line has a coefficient value of 0.29 while the shale curve, the dextral shift from the original Gardner line has a coefficient value of 0.33. Based on the deductions from the previous discussions, the new Gardner type of equations is presented in the same form as Eq. (1).
For the sands, this is: The relationship in the shales is:

Discussion
The original Gardner equation and locally derived Eqs. (3) and (4) were used to estimate the average densities for shales and sands in the wells. The estimated average densities were then compared to the measured average densities to determine which method is more accurate.  Considering the geology of the study area, the comparison is made for sediments that are expected to be in the domain of mechanical sediment compaction with minimal diagenetic alterations. Figure 11 shows there is substantial variation in the densities estimated for each of shale and sand in all the wells using the original Gardner transform. As an example, the variation can be as much as 0.19 g/cc as shown by shales in well AS and sands in well SE. Contrary to this, variations between estimated average densities in a single well and actual data using the local calibrations is not up to 0.08 g/cc for the shales (see well GB) and sands (see well AS). The range of errors associated with the local derivatives was used to determine the percentage improvement of the estimates over the original Gardner technique (see Fig. 12).
In well KC for example, the mean absolute error in estimated density of shales using the original Gardner equation is 0.15 g/cc but this error is reduced to 0.03 g/cc with the local calibration (Fig. 11). This implies that an improved accuracy of almost 80% is achieved with the modified relationships (Fig. 12). Also in the sands, the original Gardner gives a mean absolute error value of 0.18 g/cc in well SE but the local transform reduces this error to 0.02 g/cc resulting in an improved accuracy of almost 90%. Therefore, the re-calibrated Gardner equations are preferred for density-velocity transforms because of their high ability to improve the accuracy of estimates.

Summary and conclusion
The empirical relationship of Gardner et al. (1974) does not appear to give good estimates of density from velocity data in the Niger Delta. For this reason, the equation needs to be modified in order to satisfy local conditions and improve the accuracy of ensuing analyses. More realistic models of the density and better fits to observed data are possible by adjusting the coefficient in the equation so that minimal error is achieved in the analysis of named rock types. For example, as far as shales and sands are concerned, which are of greater interest to this study, the use of a default coefficient value of 0.31 in the Gardner equation is unsatisfactory. The uncertainty in the estimates is minimized when the coefficient is adjusted to provide individual reasonable fits to the shales and sands. Consequently, better estimates of density from velocity data of sands is obtained with Eq. (3) while Eq. (4) is best suited for the shales. The difference between the local modifications and the original Gardner relationship can be attributed to the differences in the composition of the originating rock types and probably because the Niger Delta is not genetically related to the sedimentary basins that Gardner et al. (1974) used in their empirical calibrations. In the Niger Delta, the shales and sands are of massive thickness and can be worked predrill because the density contrasts between both rock types makes it easier to distinguish them on most seismic data. Accepted that this separation could be difficult to realize in thin cyclic interbedding of shales and sands, in such diffuse cases, it is suggested that the default Gardner equation which is empirically an average for the two rock types may be incorporated and used with the local calibrations discussed in this study. Fig. 11 Errors associated with the estimation of densities from velocities using the derived best-fit coefficients in shales (green) and sands (yellow) are less than 0.08 g/cc. Errors from default Gardner in shales and sands are shown as grey and pink bars, respectively, and can be as much as 0.19 g/cc Fig. 12 Range of improvement in accuracy of estimates based on new coefficients over the default Gardner value. Accuracy of estimates can be improved by almost 90% using the locally calibrated coefficients