Estimation of Heat Flow Using a Bottom Simulating Reflection Based on 3D Seismic Data, West Africa

- A Bottom Simulating Reflection (BSR) has been identified using 3D seismic data from the Cameroon continental slope margin. The BSR covers an area of c. 350 km 2 in water depths ranging between 940 m to 1,750 m across an area characterized by high- and low-gradient slopes, gullies, scours, and fans. The thickness of the Gas Hydrate Stability Zone (GHSZ) is ~100 - 250 m, assuming an average velocity of 1,800 m/s. Pockmarks are inten-sively developed across the slope and most of them are observed in the BSR area. Geothermal estimation is based on hydrate stability conditions for pure methane - seawater system, hydrostatic pressure model, and a range of P-wave velocity models for the GHSZ, ranging from 1,600 -1,800 to 2,000 - 2,200 m/s. Geothermal gradient is calculated showing the range and distribution of thermal gradients in the BSR area from 0.046 o C/m to 0.094 o C/m with assumed GHSZ velocity of 1,800 m/s. Thermal gradient anomalies have been observed in association with gullies, vertically stacked channels and in some individual pockmarks. The highest anomalies of 0.08 o C/m – 0.094 o C/m are found in the depression areas of pockmark trains, within seafloor gullies. These positive anomalies are most likely controlled by active or recently active fluid advection and expulsion through the Cameroon slope.

BSR-derived geothermal gradients can be used to detect thermal gradient anomalies. BSR generally corresponds to the base of the gas hydrate stability zone (GHSZ). It therefore will follow the surface of constant temperature (isotherm) at a constant water depth (Dillon et al., 1994). Hence, the variation of the BSR depths from a predicted background model will reflect local heat flow variations cause by geothermal anomalies, sediment conductivity variations, fluid flow, and nonsteady-state sedimentation/ erosion or faulting (Shankar et al., 2010). Positive anomaly heat flows are generally interpreted as a consequence of fluid migration along unsealed fault planes, fluid venting sites, and large fluid circulation that could have radically perturbed P-T regimes (Minshull and White, 1989;Foucher et al., 1990).
In petroleum systems, heat flow and thermal history have important implications for maturation of basins with respect to oil and gas generation. Measured geothermal gradients along the West African continental margin; in Congo, Gabon, Angola, and Namibia vary between 45 o C/ km -57 o C/km (from ODP Leg 175 data, Figure  1). Deep-water setting offshore Cameroon ( Figure  2) is a new exploration area where, due to a lack of well data, the geothermal regime is poorly constrained. However, newly acquired highresolution 3D seismic data reveals the existence of a well developed BSR in deep-water setting offshore Cameroon (Figure 3) (Le et al., 2015), suggesting the widespread occurrence of natural gas hydrates. Indeed seismically, well-defined BSRs can be used to estimate thermal gradient in the shallow subsurface (<500 mbsf), providing useful information to constrain the thermal regime along the Cameroon margin. This study aims to estimate the geothermal gradient/heat flow based on the BSRs depth, and to investigate the linkages between fluid migration and thermal evolution across the continental slope system offshore southern Cameroon.

Geological setting
The Kribi-Campo Subbasin is one of the several divergent basins formed along the West African continental margin along with the Rio-Muni, North Gabon, and Sergipe Alagoas Basins known as the rift branch of the proto-Antlantic (Teisserenc and Villemin, 1989;De Matos, 1992) (South America). The Kribi-Campo Subbasin is located between 2º20' N to 3º00' N latitude, 9º00'E to 9º50' E longitude, extends over 6,150 km 2 offshore and   Figure. 1. Geothermal gradients estimated from Ocean Drilling Program (ODP Leg 175) data, using Adara tool to measure in-situ formation temperature in three points basic. Penetration depths are less than 600 mbsf, most of depths are less 200 m below seafloor (mbsf). Geothermal gradient ranges from 45 o C/ km to 57 o C/ km. The costal line map is extracted from National Geophysical Data Centre (NGDC) (http://rimmer.ngdc.noaa.gov/mgg/coast/getcoast.html). 45 km 2 in triangular onshore area (Ntamak-Nida et al., 2008). The basin is bounded by Pre-Cambrian basement outcropping close to the present day shoreline to the east, the southern margin is delimited by the Kribi Fracture Zone (KFZ) ( Figure  1), separating the Kribi-Campo Subbasin from the Equatorial Guinea Rio-Muni Basin. The northern boundary of the basin is unclearly determined and a previous study suggested that it is possibly the northern extension end of Aptian Salt (Davison, 2007). However, based on our data analysis, salt seems to be absent in the deep water part of the Kribi-Campo Subbasin. The northwest part of the basin is defined by the Cameroon Volcanic Line (CVL), trending NE -SW and extending onshore.
The Kribi-Campo Subbasin contains three successions: a lower syn-rift sequence (Barremian-Aptian), a rift to drift transition sequence (Late Aptian), and a post-rift sequence (Albian to Present) (Ntamak-Nida et al., 2010). The thickness of the sedimentary section in the Kribi-Campo Subbasin varies from 8 to 10 km (Pauken, 1992). During its evolution, the basin experienced several regional tectonic events resulting in uplift, deformation, and erosion, including the Albian-Aptian break up unconformity (115 Ma), Santonian uplift (85 Ma), and resumed Eocene-Oligocene uplift (c. 45 Ma) (Ntamak-Nida et al., 2010), marked by major unconformities (Lawrence et al., 2002). Turner et al. (2008) linked these uplift events with the cooling episodes identified from the apatite fission track analysis, vitrinite reflectance data, and stratigraphic observations. According to these authors, heat flow is not constant over time; the basin has experienced higher temperatures in the past. The evidence of higher temperature is observed in the case of paleo-geothermal gradient in the well Rio-Muni-1 decaying between 58 o C/ km during the mid-Cretaceous and 21.5 o C/km at the present day. The thermal gradient of the basin should take into account the effect of the Cameroon Volcanic Line to the northwest. The growth of the Cameroon Volcanic Line began in the Oligocene with peak activity in mid-Miocene, however still active with Mount Cameroon. (Meyers et al., 1998). Volcanic activities may cause anomalous temperature enhancing thermal maturity of shallow source rocks.

Dataset
The studied area is located offshore Cameroon in water depths ranging from 600 m to 1,900 m. The dataset is a 3D seismic volume which was acquired in a northeast−southwest orientation, covering an area of 1,500 km2, with a bin spacing of 25 m and a total record length of 6.6 sTWT. The interval focused on this study is down to c. 1,000 ms below the seafloor corresponding to Pliocene to Holocene sequence ( Figure 3). The seismic data are displayed as a zero phase, SEG normal polarity, i.e. orange/brown peak indicating an increase in acoustic impedance. In this shallow section, dominant frequency in the 44 -48 Hz range give a vertical resolution (λ/4) of 8 -9 N te m R. Figure 2. Bathymetry map of the continental slope on the Cameroon margin shows location of the studied area, number of wells has been drilled within the basin (Pauken, 1992) and the fracture zone system (Meyers et al., 1996). The costal line map is extracted from National Geophysical Data Centre (NGDC) (http://rimmer.ngdc.noaa.gov/mgg/ coast/getcoast.html).
m, assuming seismic velocity in the shallow sediment of c. 1,600 m/s. Additional information was obtained from the Ocean Drilling Program (ODP) and published data from surrounding areas. The additional data includes: (i) Velocity-depth relationship to convert BSR reflection from time (Ferguson et al., 1993;Yuan et al., 1996;Pecher et al., 2001;Rao et al., 2001;Ojha and Sain, 2007;Minshull and Keddie, 2010); (ii) The seafloor temperature from two sources, Sterling Energy Company (2002) and the National Oceanographic Data Centre (NODC) ; (iii) Salinity from NODC ; (iv) Thermal conductivity of the sediment interval between BSR and seafloor from ODP drill cores site 1075 -Congo Basin, site 1078 -Angola Basin and site 1081 -Namibia Basin with the average number (Wefer et al., 1998).

Methodology
Geothermal gradient/surface heat flow in the deep-water setting offshore Cameroon is calculated based on the existence of BSR which is known as a surface of constant temperature at a certain depth. Estimation method has been applied to calculate the geothermal gradient. The method is a two-point approximation of the geothermal gradient introduced by Yamano et al. (1982). Regarding this method, Grevemeyer and Villinger (2001) noted that it is generally simplistic, but it can be applied to place additional constraints on the thermal gradient along continental margins. According to these authors, the two-point approximation is only appropriate for homogeneous sediments with a constant thermal conductivity. This method depends on both estimated temperature at BSR depth and seafloor temperature as used to estimate temperature change between seafloor and BSR depth. Geothermal gradient can be displayed in a map view, showing the variation and anomalies in the shallow geothermal gradient.

Seafloor And BSR Mapping
In order to calculate the geothermal gradient, a surface of BSR needs to be mapped. Seismic interpretation has been carried out using the Petrel software version (2009-2) with a fine grid basis for quality control; auto-tracking was used for the seafloor resulting in a good quality continuous surface. BSR is characterized by a negative impedance contrast reflection, opposite polarity to the seafloor reflection. In some areas, the BSRs have been mapped with a strong confidence, i.e. all the criteria mentioned by Berndt et al. (2004) are satisfied.
To calculate the geothermal gradients and heat flow systematically, BSR and seafloor surfaces were converted into random 14,000 point grid for each surface and exported into Excel format. The number points are enough to represent the data and it is below the point limitation in Excel software.

-Temperature at the BSR depth
In this study, the formula used to compute the geothermal gradient is extracted from Dickens and Quinby Hunt (1994) for the "pure methaneseawater" system at the given pressure between 2.5 -10 MPa. However, it can be applied for a higher pressure (Dickens and Quinby Hunt, 1994).  1) where T is temperature (K) and P is pressure (MPa) An assumption of pure methane has been applied for the present study. A study conducted by Charlou et al. (2004) indicated that the gas hydrates recovered from Congo -Angola at 3,000 m water depth are 100% methane hydrate. Samples from deep to ultra-deep-water Nigeria also indicated for the hydrate that are primary composed of biogenic methane, more than 99% of hydrocarbon gases (Brooks et al., 2000). Further analysis of the carbon and hydrogen isotope ratio in methane in Congo-Angola and Nigeria of both studies proved the biogenic origin of methane, resulting from bacterial transformation of CO 2 which is a common process in shallow marine sediments (Claypool and Kaplan, 1974;Whiticar, 1999). However, small amounts of thermogenic gas can be observed in natural gas hydrate samples (Brooks et al., 2000).

-Geothermal gradient from BSR depth
Geothermal gradient is calculated by using formula (2) below based on the difference of the temperatures at the BSR and the seafloor, and the depth of the BSR beneath the seafloor (Shankar et al., 2010). Heat flow estimated for the GHSZ The heat flow per an area unit will be achieved by using formula (3) below as the product of thermal conductivity (K) and temperature gradient (Lee, 1963) , assuming a linear temperature gradient and using an average thermal conductivity from the thermal conductivity range which taken from one well in the adjacent area.
where Q = heat flow (W/m 2 ) K = thermal conductivity (W/mK -Watts per meter per Kelvin)

∆G = Geothermal gradient (K/m)
Geothermal gradient and heat flow calculation requires a number of parameters. Those parameters availability and low uncertainty will bring the estimation close to the real value.

Results
Seafloor and BSR mapping accompaning with temperature and heat flow estimation has been broken into seven steps applied in software 'Pe- Step 1. Seismic Mapping: Sea floor and BSR mapping Seafloor has been mapped as a positive reflection using auto-tracking. Bottom simulating reflection is mapped as a negative polarity which satisfied all of the criteria defined by Berndt et al. (2004). BSR mapping used fine interpretation grid 5 × 5 inline and crossline spacing. Subsequently, surfaces are generated for seafloor and BSR area shown in the Figures 3a and b and Table 1.
Step 2. Depth Conversion: Velocity-depth relationship to convert reflection times to depths. P-wave velocity (Vp) in the seawater is chosen a standard value of 1500 m/s. The P-wave velocity in the interval between seafloor and the BSR (in the GHSZ) varies in calculations conducted by different authors, ranges from 1,600 m/s to 2,300 m/s (Ferguson et al., 1993;Yuan et al., 1996;Pecher et al., 2001;Rao et al., 2001;Ojha and Sain, 2007;Minshull and Keddie, 2010). Some authors concluded the effect of gas hydrate on velocity of the gas hydrate bearing sediment is significant Ojha and Sain, 2007), but others ignored it as minimal (Shankar et al., 2010). This is important because the increase in velocity will Step Seismic velocity ( result in increasing BSR depth and the thickness of GHSZ. It leads to a higher temperature at the BSR depth, and produce higher geothermal gradient. Gangully et al. (2002) indicated that an increase of 180 m/s velocity for the GHSZ would change the calculated heat flow by only 2 -3 %.
In this study, as there is no direct estimate of the velocity, it is assumed a constant velocity in the GHSZ and run calculations for discrete velocity values between 1,600 m/s and 2,200 m/s in order to examine the variation of the geothermal gradients and heat flows with the variation of velocity in the GHSZ.
Conversion from two-way-travel time (s) to depth (m) basically uses the equation below: The result is shown in Table 1 for the water depth (D seafloor in m), the thickness of GHSZ (Z GHSZ in m), and the BSR depth (D BSR in m).

..(6)
Step 3. Pressure at the BSR Depth Hydrostatic pressure was applied for estimation pressure at the BSR depth. There is no clear statements for which is more accurate pressure model, hydrostatic, or lithostatic (Ganguly et al., 2000). The difference in pressure for the two models is considerably large, leading to a significant change in thermal gradient results. Geothermal gradient derived from the hydrostatic model is c. 8−12 % smaller than that of lithostatic model (Ganguly et al., 2000). Lithostatic model is only applied if overpressure occurs (Ferguson et al., 1993). Even though overpressure are unlikely at the shallow depth (Hyndman et al., 1993), in the area where active focused fluid flows are observed, overpressure may locally occur at very shallow sediments, such as in the VØring Plateau (Brown et al., 2006) or Congo Basin (Gay et al., 2007).
where P BSR_Hydro = Hydostatic pressure at the BSR depth (MPa) ρ seawater = Density of the seawater= 1034 kg/ m 3 in which Salt concentration (S)= 34 ‰ , achieved from National Oceanographic Data Centre (NODC) . g = Gravity at sea level = 9.8 m/s 2 D BSR = Total depth to the BSR from the sea level (m) The results are shown in Table 1 Step 4. Seafloor Temperature Equation Seafloor temperature is derived from a number of sources including well and temperature profile through the water column acquired during seismic acquisition campaign of Sterling Energy Company (2002) (Figure 4). Seafloor temperature from Sterling Energy Company has been used for the calculation since it has been locally measured, thus suggested for more precisely. Seafloor temperatures over the BSR-mapped area, range from 4.11 o C to 6.67 o C, generally defined by the following equation: Step 5. Geothermal Gradient Geothermal gradients were computed by using the Equation (2) where the thickness of the GHSZ and temperatures at the BSR depth and seafloor temperature were achieved in the Step 2, Step 4, and 5 respectively. The results are shown in Table  1 and Figure 5 for Vp of 1,800 m/s. Geothermal gradients calculated in two-point basic, provide the geothermal gradient distribution for the whole area that BSR occurred. The result shows geother- Step 6. Thermal Conductivity: Thermal conductivity from probe measurements to obtain heat flow Thermal conductivity difference between sediments and hydrate -bearing sediments is usually small (Grevemeyer and Villinger, 2001), and likely not detectable in borehole measurements . Here thermal conductivity ranging between 0.7 -0.8 W/mK are used to compute the heat flow. The value was taken from needle probe measurements on ODP drill cores site 1075 -Congo Basin (0.6 -0.8 W/mK), site 1078 -Angola Basin (0.8 -1.1 W/mK), and site 1081 -Namibia Basin (0.7 -1.1 W/mK). All sites have shallow penetration depth c. 200m (Wefer et al., 1998).

Step 7. Heat Flow
Heat flow was estimated by using Equation (3) in which thermal conductivity and temperature gradient (∆G) was gained in Step 6 and Step 7. The results are shown in Table 1 and Figure 5 (P-wave velocity of 1,800 m/s).

Geothermal Gradient Variation with The Thickness of GHSZ
Calculation has been carried out based on some values of P-wave velocities, 1,600 m/s, 1,800 m/s, 2,000 m/s, and 2,200 m/s. The results are illustrated in Table 1. It is apparent from this table that the geothermal gradient/heat flow is the highest when considering low P-wave velocities. There is a clear relationship between geothermal gradient/heat flow and the thickness of the GHSZ; geothermal gradient decreases with the thickening of the GHSZ (Table  1). The increase of velocity linked with the greater of thickness of the GHSZ leads to reducing geothermal gradient. In general, an increase in P-wave velocity by 200 m/s at the interval between seafloor and BSR produces a decrease in the calculated geothermal gradient /heat flow by c. 9 -10%.
The relationship between velocity and heat flow was also mentioned by Gangully et al. (2000). According to these authors, heat flow is

Geothermal Gradient/Heat Flow Estimations
In order to calculate geothermal gradient/heat flow, one important factor is the velocity within the GHSZ to convert the thickness of the GHSZ into depth and finally the BSR depth. However, the velocity of the GHSZ is unclear as there are no borehole data in the studied area. Therefore, the velocity of the GHSZ will be assumed based on the published data and the shallow sediment velocity from adjacent ODP wells.
Different authors have measured the thickness of the GHSZ by using different P-wave velocities. Winters et al. (2007) investigated the effect of gas hydrate and ice on acoustic velocity in different sediment types. In their study, the P wave velocity in fine-grained sediments was 1970 m/s. However, this value can only be applied for the interval of gas hydrate bearing sediments. From seismic data, it is difficult to determine clearly the gas hydrate distribution and concentration within the sedimentary column. It therefore used an average value of P-wave velocity between velocity in shallow marine sediment (>1,600 m/s) (ODP Leg 175) and gas hydrate bearing sediment of 1,800 m/s for the velocity of GHSZ. This value will be used for discussion and illustration for the results.
One of the factors caused thermal gradient high possibly is thermal conductivity. Low thermal conductivity in shallow sediments is interpreted to occur in the studied area due to Pleistocene mud-rich sediments, which is interpreted, based on seismic facies, widespreading over the whole area. The occurrence of gas hydrate tends to decrease the thermal conductivity of sediments (Stoll and Bryan, 1979), but not significant (Grevemeyer and Villinger, 2001;Henninges et al., 2005). The occurrence of mud-rich sediments in shallow depth caps the fan lobe-channel sands below could have been abruptly prevented the heat transmitting through it, leading to raising the thermal gradient in the shallow sediment.

Geothermal Gradient Variation along The Slope and Its Anomalies
Geothermal gradient varies significantly across the studied area between 0.046 o C/m and 0.094 o C/m ( Figure 5). In general, the geothermal gradients increase with thickenings of the GHSZ;

I J O G
Indonesian Journal on Geoscience, Vol. 8 No. 3 December 2021: 297-311 with low thermal gradient upslope and high thermal gradient downslope.
The interesting aspect of this dataset is the geothermal gradient anomalies which coincide with some single pockmarks, pockmark trains ( Figure  6), gully 3 and the vertically stacked channels complexes (Figure 7) (see more descriptions in Le et al. (2015). In this study, the term geothermal gradient anomaly is used to describe the higher geothermal gradient compared with the surrounding areas.
Thermal gradients at some single pockmarks which are associated with seismic chimneys, pull up BSR also displayed as anomalies with circular shapes, the value is c. 0.075 o C/m which is about 0.01 o C, higher than those in surrounding areas. An example of this is shown in Figure 6 in  which the geotherm is slightly higher and gives a closure circular contour corresponding to a single pockmark on the seafloor. Furthermore, thermal gradients are slightly higher along the vertically stacked channels by c.0.005 o C/m (Figure 7). Pockmark has been known as a direct and violent migration of fluid from the overpressured zone below. Once it occurred, it may bring the high heat from deep to the shallower section and caused the anomalies high geothermal gradient.
All the anomalies observed show a strong relationship with underlying sediment wave boundaries (Figure 8), fan/channels complexes, and seabed fluid flow features, characterizing a shallow plumbing system enhancing warm fluid advection from depth.

Conclusions
Thermal gradient/heat flow along the Cameroon margin has been calculated using BSR depths from 3D seismic data. A wide range of thermal gradient and heat flow has been observed with consistent trend as well as some local anomalies. The results could be highlighted as follows: • The geothermal gradient in the shallow sediments in the studied area ranges from 0.046 o C/m to 0.094 o C/m for the case Vp = 1,800 m/s. The addition of c. 8% -10 % is possibly necessary if lithostatic pressure environment is present. Estimated thermal gradients are higher than the measured value from deep wells suggested for the cause of fluid migrations and low thermal conductivity in shallow gas hydrate bearing -sediments. • Positive heat flow anomalies are detected in gullies, vertically stacked channels, and some single pockmarks which are associated with pull-up BSRs and seismic chimneys.
Anomalies are significantly high in the gully, and reach a peak in the depression areas of pockmark trains of 0.08 o C/m -0.094 o C/m indicating recent activity of the fluid migration network in the studied area. • Heat flow/geothermal gradient are very sensitive to the thickness of the GHSZ or the changes of P-wave velocity. The values will decrease c. 9 -10% for the greater of 200 m/s P-wave velocity.
The studied area has a high thermal gradient in the shallow sediments which can be lower at depth. Many assumptions make it difficult to constrain uncertainties; hence future investigations might need to re-assess the estimated value with the availability of well data.