Statistical analyses of reservoir and fracturing parameters for a multifractured shale oil reservoir in Mississippi

The use of multifractured horizontal wells has improved the efficiency of hydrocarbon extraction from shale gas and oil plays. It is highly desirable to estimate the characteristics of the reservoir and well fracturing through production data analysis. Production rate transient data from 16 wells in the Mississippi sections were analyzed to estimate local reservoir permeability and hydraulic fracture parameters. Key factors affecting well productivity have been identified for optimizing future well completion. Based on the theory of distance of investigation during transient linear flow, the reservoir matrix permeability of TMS was estimated to be between 53 nd and 210 nd, averaging at 116 nd. The matrix permeability follows the lognormal distribution and is considered homogeneous according to the coefficient of variation. The fracture half‐length was estimated using the matrix permeability data from the rate transient analysis. The fracture half‐length was found to have a mean value of 234 ft with a standard deviation of 66 ft in a normal distribution. The fracture conductivity was back‐calculated by matching pseudosteady production rate data to Guo et al's (SPE Reserve Eval Eng.12, 2009, 879) productivity model for boundary‐dominated flow. The fracture conductivity was estimated to range from 0.4 md‐ft to 3.2 md‐ft, averaging at 1.4 md‐ft. Based on Guo et al's (SPE Reserve Eval Eng.12, 2009, 879) well productivity model applied to TMS condition, well productivity can be improved significantly by increasing fracture length and conductivity. Future TMS wells should be completed with conductivity values of greater than 10 md‐ft.


| INTRODUCTION
The Tuscaloosa Marine Shale (TMS) across Louisiana and Mississippi is a sedimentary formation that deposited during the Upper Cretaceous. 1,2 It consists of Upper Tuscaloosa, Middle Marine Shale, and Lower Tuscaloosa. 3,4 The central part of Tuscaloosa is composed of dark gray, fissile, and sandy marine shale. 2 Multifractured horizontal wells (MFHW) have enabled commercial production of oil from shale reservoirs. 5,6 More than 80 multifractured horizontal wells were drilled in the TMS play during 2012 and 2015. 4,7 Some wells were reported to have appealing initial production rates of greater than 1000 stb/day. However, drilling activities slowed down due to the high cost of drilling and low oil prices. It is of considerable value to investigate the characteristics of the reservoir and well stimulation through historical production | 617 YANG ANd GUO data analysis, which will provide a useful method for assessment of the production potential of TMS wells.
A broadly used method for assessment of stimulation effectiveness is production data analysis (PDA). PDA begins with the identification of flow regimes, which can be further analyzed for estimating reservoir and/or fracture properties, such as permeability, original oil in place (OOIP), EUR, and fracture half-length. [8][9][10][11][12] For multifractured horizontal wells, two general flow regimes may be observed: transient linear flow and boundary-dominated flow. Transient linear flow is the dominant flow regime in low-permeability reservoirs and may last for a long time, depending on the matrix permeability and fracture spacing, etc The most common method to identify these two flow regimes is the use of a log-log diagnostic plot of production rate versus time when only production data are available. 7,13 For hydraulically fractured wells with infinite fracture conductivity, the log-log plot follows a straight line with a slope of −0.5 during the transient linear flow period. The slope gradually decreases to −1.0 over time. After the pressure transient reaches the extent of the fracture's drainage area, the production data versus time data plotted on the log-log plot should show a slope of −1.
For multifractured horizontal wells with finite conductivity fractures, the shape of a type curve may be different from the shape of infinite-conductivity fractures. In order to account for the influence of flow convergence or skin effect caused by damage to finite conductivity fractures, the square root of time plot was introduced to help identify transient linear flow. [14][15][16] During the linear flow period, the square root of time plot follows a straight line on a Cartesian coordinate system. The slope of this line can be used to calculate the product of well drainage (or fracture half-length) and the square root of permeability. 17,18 Estimates of fracture conductivity have not been well studied in the past. The short-term fracture conductivity is usually estimated by pressure transient analysis of fracture linear flow or bilinear flow. 11,19 However, these two flow periods may not be observed due to wellbore storage. 11 Furthermore, this method is suitable for constant productivity conditions. 19 Another way of estimating fracture conductivity is to match production data. Guo et al 20 proposed an analytical solution for simulating the productivity of multifractured horizontal wells under pseudosteady state conditions. The oil production rate is a function of long-term fracture conductivity. Therefore, an estimate of fracture conductivity can be obtained by solving the pseudosteady state production equation.
Previous studies have investigated the relationship between production rate and underlying factors such as reservoir characteristics, geological parameters, stimulation parameters, and completion parameters. [21][22][23] Esmaili et al 21 presented a data-driven analysis of completion and hydraulic fracturing parameters, including stimulated lateral length, number of clusters, cluster spacing, and minimum distance to the offset well. Temizel et al 22 illustrated the significance of parameters that can affect production rate or cumulative production. They investigated the major parameters affecting the performance of horizontal wells in tight formations. Panja and Deo 23 performed statistical analysis to identify important factors that influence the production potential of low-permeability reservoirs. Li et al 24 stated that the productivity of multifractured wells increases with increased fracture length, increased fracture width, and decreased fracture spacing. There is an urgent need to estimate the reservoir characteristics and stimulation effects through production data analysis, which provides a useful method for assessing the production potential of TMS wells.
The objective of this study was to perform statistical analysis to identify critical parameters that affect TMS well productivity from production data. A method for estimating the fracture conductivity was proposed and verified by field data. Correlation analysis was conducted to investigate the association between cumulative production and other underlying factors, such as permeability, fracture half-length, and effective lateral length. Finally, the production potential of the TMS wells was analyzed.

| Reservoir permeability
The reservoir permeability can be estimated by using the concept of the distance of investigation, [25][26][27] where k m is the permeability in md, ϕ is the porosity, μ o is the viscosity in cp, c t is the total compressibility in psi -1 , t elf is the time at the end of linear flow in days, y e is the distance from fracture to the boundary (ie, half of fracture spacing, 0.5S f ) in ft, and S f is the fracture spacing in ft.

| Fracture half-length
Under the constant flowing pressure condition, the square root of time plot can be used to describe the transient flow behavior of hydraulically fractured wells, 11,14,17 where q o is the oil production rate in stb/day, t is the time in days, m is the slope of the square root of time plot, and b' is the corresponding intercept representing the apparent skin.
The product of well drainage area and the square root of permeability can be calculated from the slope of the square root of time plot, 8,10,17,18 (1) where A cm is the well drainage area in ft, 2 B o the oil formation factor in rb/stb, p i is the initial formation pressure in psi, and p w is the bottom-hole flowing pressure in psi.

The drainage area for multifractured wells is
where x f is the fracture half-length in ft, n f is the number of fractures, and h is the fracture height in ft. If we assume that all the fractures are identical in fracture half-length (ie, A cm = 4n f x f h), the average fracture half-length can be computed,

| Fracture conductivity
The oil production rate in the pseudosteady state flow period of multifractured horizontal wells can be estimated using the following equation, 20,28 and where p a is the average formation pressure in psi, and C f is the fracture conductivity in md-ft. The average formation pressure can be calculated by means of the following equation if the formation pressure is known, 20,28 where p e is the reservoir pressure in psi.
The average formation pressure can also be estimated using OOIP, 7,29,30 where N p is the cumulative production in stb, and N i is the OOIP in stb.
OOIP can be calculated using the following equation, where S w is water saturation. It is evident that the oil production rate is a function of fracture conductivity. Solving Equation (5), one can get the average fracture conductivity of multifractured horizontal wells.

| Verification of the method
The performance of the proposed method was verified using data from the Eagle Ford Shale (EFS) well #6 estimated by Alotaibi et al 18 The initial formation pressure is 8,000 psi; the flow pressure is 2000 psi; the oil viscosity is 0.4 cp; the porosity is 0.09; the total compressibility is 1.14 × 10 -5 /psi; the oil formation factor is 1.4 rb/stb; the number of fractures is 12; the fracture spacing is 402 ft; the slope of the square root of time is 0.000136; and the time at the end of the half-slope is about 40 days. The results are listed in Table 1.
The time at the end of the linear flow can be estimated by using the time at the end of half-slope. Previous studies have shown that pressure propagation reaches the boundary before the end of half-slope. 17,22 The time at the end of the linear flow of a multifracture horizontal well is approximately 0.32 times at the end of the half-slope (ie, t elf = 0.32t ehs ) according to Tabatabaie et al 26 Therefore, the time at which the liner flow ends is about 12.8 days. The matrix permeability is computed to be 0.051 md, according to Equation (1). The average fracture half-length is then calculated to be 176 ft, which is 2.27% higher than the estimated 172 ft fracture half-length by Alotaibi et al 18 At the beginning of the boundary-dominated flow, EFS well #6 had an oil production rate of approximately 1,236 stb. Solving Equations (5) and (7) simultaneously gives a fracture conductivity of 23 md-ft. It should be mentioned that the average formation pressure is computed using Equation 7,29,30 because we do not know the cumulative oil production. It is assuming that the formation pressure is equal to the initial value. The fracture conductivity of well #6 estimated by Li et al 24 using the effect of fracture-closure stress on resin-coated-sand permeability is about 30 md-ft. Obviously, by using Equation (5), the underprediction of fracture conductivity is estimated to be −23.33%. If the oil production rate remains the same, the higher the

| Statistical analysis of TMS wells
In this section, the characteristics of hydraulic fractures in the TMS trend are statistically analyzed. Production data as of March 2019 came from 55 TMS wells in Mississippi were collected. These data were downloaded from the Mississippi State Oil & Gas Board (www.ogb.state.ms.us). This Web site reports information on TMS wells, such as well location, monthly oil production, water production rate, natural gas production, and well production days. Drilling activities were centered in Amite and Wilkinson Counties. However, drilling activities in the region slowed down in the summer of 2015 due to high drilling costs and low commodity prices. The statistical results of 16 TMS wells were analyzed because the pseudosteady state flow state was identified by flow regime diagnosis analysis. In fact, we collected the stage information for 23 wells. However, the transition from the transient flow to the boundary-dominated flow of seven wells cannot be determined due to abrupt changes in the production rate over time. The location of these wells is presented in Figure 1. Table 2 lists the statistical results of 16 wells. Two TMS wells were analyzed as examples to illustrate the process in detail. TMS well #1 has a true vertical depth (TVD) of 12 328 ft and an effective lateral (L) of 5204 ft. This well was completed with 19 stages. Figure 2(A) presents the square root of time plot for well #1. The estimated time at the end of half-slope was about 522 days, and the corresponding production rate and cumulative oil production were 176 stb/day and 189 Mstb, respectively. Therefore, the estimated time at the end of the linear flow was approximately 167 days. This means that the pressure had reached the virtual boundary at t = 167 days. Porosity, oil viscosity, and total compressibility are 0.08, 0.5 cp, and 1.0 × 10 -5 /psi, respectively. 4 Substituting these data into Equation (1) gives a matrix permeability of 111 nd. The initial formation pressure gradient and flowing pressure were 0.52 psi/ft and 0.3 psi/ft, respectively. Then, the fracture half-length was computed to be 360 ft according to Equation (4). Then, the OOIP was calculated as 17.08 MMstb. Solving Equation (5) Figure 2B.
The time determined at the end of half-slope was 576 days, and the corresponding oil production rate and cumulative production were 144 stb/day and 168 Mstb, respectively. The proposed method gives matrix permeability, fracture halflength, and fracture conductivity of 86 nd, 290 ft, and 1.3 md-ft, respectively. The probability plot is an excellent way of identifying the possible distribution of data. In this study, the distribution of data shown in Table 2 was determined using the quantilequantile (Q-Q) plot. Figure 3 shows a typical plot of the time quantiles and theoretical values at the end of linear flow. The time at the end of the linear flow of 16 wells ranged from 48 to 252 days. The Anderson-Darling (AD) test was performed to evaluate the null hypothesis that the time was from a normal distribution (CI: 95%). A lower AD value indicates a better fit, while a lower P-value (eg, <0.05) indicates that the data do not follow the assumed distribution. The P-value was found to be 0.67, and the data points approximated a straight line, indicating that the data were normally distributed. The average value of t elf was about 160 days, suggesting that these wells experienced long-term transient linear flow. According to Equation (4), we know that the time at the end of linear flow is related to matrix permeability and fracture half-length. Lower matrix permeability and larger fracture half-length tend to increase the time at the end of the linear flow. with an average of 116 nd. Lohr and Hackley 31 tested the permeability of four core samples from two TMS wells in Mississippi. The permeability they tested ranged from 59 to 133 nd, averaging at 99 nd. Obviously, the experimental values were of a similar order of magnitude as the values predicted using production data. The P-value was found to be .74, which was higher than the significance level (α = 0.05). Therefore, we fail to reject the null hypothesis that the permeability data follow a lognormal distribution. From the lognormal probability plot, we can see that the permeability data form an approximately straight line along the fitted distribution line, indicating that the lognormal distribution appears to be a good fit for the permeability data of 16 TMS wells. The coefficient of variation can be used to measure the heterogeneity of permeability. For a normal distribution, the heterogeneity of permeability is defined as the ratio between the sample variation and its mean value, 32 where C v is the coefficient of variation, σ is the standard deviation, and μ is the mean value.
In many cases, the logarithm of permeability is found to be normally distributed. 28 In this case, the coefficient of variation is expressed as, 33 where s is the standard deviation of ln(k m ), and k m is in md.
If C v < 0.5, the permeability is considered homogeneous, and if C v > 1, the permeability is assumed to be very heterogeneous. 32 The standard deviation of ln(k m ) was found to be 0.41, which yielded a coefficient of variation of 0.42. Because C v was <0.5, the permeability of 16 wells was considered to be uniform. Figure 5 shows the normal Q-Q plot comparing estimated fracture half-length on the vertical axis to the standard normal data on the horizontal axis. The fracture half-length of 16 wells ranged from 97 to 360 ft, averaging 234 ft. The AD value and P-value were found to be 0.19 and 0.88, respectively. Therefore, the fracture half-length may follow a normal distribution. Besides, the linear relationship of the data points indicated that the fracture half-length data were  Figure 6 presents the normal and lognormal Q-Q plots of fracture conductivity. The fracture conductivity ranges from 0.4 to 3.2 md-ft, with an average of 1.4 md, according to Table 2. The average fracture conductivity of EFS wells was less than 5 md-ft. 34 Evans et al 34 indicated that the effective fracture half-length was about 200 ft because the fracture conductivity was reduced from just 200 ft from the wellbore to <0.1 md-ft. Due to a similar fracturing design, the fracture conductivity of TMS wells drilled between 2012 and 2014 was estimated to be close to EFS wells. As can be seen from Figure 6, P-values for both cases are higher than 0.05. Therefore, we cannot reject the null hypothesis that the matrix permeability comes from a normal distribution or a logarithmic distribution. As mentioned earlier, lower AD values are associated with higher P-values. Given the higher P-value, the normal distribution was chosen as the best fit for the fracture conductivity data. Figure 7 shows the probability plot of dimensionless fracture conductivity. According to Table 2, the dimensionless fracture conductivity is between 13 and 149, with an average of 59. The logarithm of dimensionless fracture conductivity follows a straight line, indicating that the dimensionless fracture conductivity is a lognormal distribution. Figure 8 shows the Q-Q plot comparing the standard normal distribution to the estimated OOIP. The P-value was higher than the confidence value of 0.05. The linearity of the data suggests that OOIP follows a normal distribution with mean and standard values of 13.6 and 5.0 MMstb, respectively.
Correlation analysis was conducted to investigate the association between cumulative production as of March 2019 (N p ) and other potential factors, such as permeability, fracture half-length, and effective lateral length (L). If the absolute value of the correlation coefficient is less than 0.65, we cannot determine the strong relationship between different variables. 35 Furthermore, if the P-value is less than the significance level (0.05 in this study), the correlation is significant. As can be seen from Table 3 and Figure  9, higher cumulative production levels are associated with higher fracture half-length, effective lateral length, and OOIP values. Pearson's correlation (r) between cumulative production and fracture half-length was found to be 0.73 (>0.65), indicating a strong correlation between these two variables. If the effective lateral length increases, the cumulative production increases. Therefore, the cumulative oil production also increases with the increase in OOIP, as shown in Figure 9C. The P-value between cumulative oil production and matrix permeability, fracture spacing, and fracture conductivity is much >0.05, which provides inconclusive evidence for the correlation between these variables. It should be mentioned that other variables not  present in our study may affect the correlation analysis. In addition, statistical analysis has only 16 observations, and more data are needed to understand the correlation between variables better.

| Production potential of TMS wells
According to Equation (5), the oil production rate is mainly related to matrix permeability, fracture conductivity, fracture half-length, and fracture spacing. This section provides sensitivity analysis to identify the production potential of TMS wells. We attempt to determine the extent to which these factors affect productivity and rank them in descending order of oil productivity. The average parameters of 16 TMS wells were used unless otherwise stated. The pay-zone thickness is 137 ft; the matrix permeability is 116 nd; the fracture spacing is 64 ft; the average formation pressure is 4980 psi; the flowing pressure is 3600 psi; the fracture half-length is 234 ft; the fracture conductivity is 1.4 md-ft; the oil viscosity is 0.5 cp; the oil formation factor is 1.3 rb/stb. Substituting these data into Equation (5) produces productivity at the beginning of the boundary-dominated flow of 127 stb/day. Figure 10 shows the effect of matrix permeability on the oil production rate. It is evident that oil productivity increases as the matrix permeability increases. For example, if the matrix permeability is 200 nd, the estimated oil production rate is about 194 stb/day, provided that other parameters are constant. Therefore, areas with higher permeability are more likely to be the sweet spot for future infill development. Figure 11 presents the influence of fracture conductivity on oil production rate. The oil production rate increases with the increase in fracture conductivity in the low-conductivity region, but it is not proportional. If the fracture conductivity can be increased by 10 times from 1.4 to 14 md-ft, the oil production rate will increase to 171 stb/day, an increase of 35%. Previous studies have shown that channel fracturing leads to a 46% increase in production in EFS trend because it increases the fracture conductivity by 1.5-2.5 orders of magnitude. 36 However, it is evident that the increase in oil production slows down when the fracture conductivity exceeds a specific value. It seems that reasonable fracture conductivity for this case is about 10 md-ft. Figure 12 shows the influence of fracture half-length on oil production rate. If the fracture half-length is increased, the oil production rate increases but less than proportionally. For example, if the fracture half-length is increased from 234 ft to 300 ft and 400 ft, the oil production will increase by 15% and 24%, respectively. It should be mentioned that the optimization of completion method is not only focused on maximizing productivity but also on additional costs. The pursuit of long fractures or small fracture spacing means that more fracturing fluid and proppant are needed, which in turn increases the completion cost. Therefore, an economic analysis should be performed to consider the additional cost and determine which design is the most cost-effective. It can be seen from Figures  10-12 that the matrix permeability is most sensitive to productivity, and the fracture half-length is the least sensitive. Therefore, it is recommended to identify core areas with a permeability bubble map in TMS trend.