WEBVTT Kind: captions Language: en-US 00:00:03.590 --> 00:00:08.060 Well, good morning. And thank you for joining us for today’s ESC seminar. 00:00:08.060 --> 00:00:13.820 Next week, we have a talk by Buka Nweke from UCLA. 00:00:13.820 --> 00:00:16.900 Just a quick reminder about today’s seminar. 00:00:16.900 --> 00:00:22.140 Please keep your cameras turned off and keep your microphones on mute 00:00:22.140 --> 00:00:27.000 until the end of the presentation and when it’s time to take questions. 00:00:27.000 --> 00:00:31.760 And, with that, I pass it on to Elizabeth Cochran to introduce today’s speaker. 00:00:31.779 --> 00:00:37.720 - All right. Thanks, Noha. So today’s speaker is Iason Grigoratos, 00:00:37.720 --> 00:00:42.920 who studied civil engineering at the University of Patras in Greece. 00:00:42.920 --> 00:00:46.880 He later obtained his master’s degree in earthquake engineering and 00:00:46.880 --> 00:00:50.570 engineering seismology in the University of Grenoble Alpes 00:00:50.570 --> 00:00:57.760 in France and IUSS in Pavia in Italy. He stayed in Italy to start his Ph.D. 00:00:57.760 --> 00:01:00.530 on induced seismicity but then moved to the U.S. to spend 00:01:00.530 --> 00:01:04.680 the last two years of his Ph.D. as a visiting research fellow at 00:01:04.680 --> 00:01:09.090 the University of Texas at Austin, where he is now a postdoctoral fellow 00:01:09.090 --> 00:01:14.759 assessing the elevated seismic risk and hazard related to human activities. 00:01:14.759 --> 00:01:21.409 So I’ll pass it over to Iason, and his title today is Oil and Gas 00:01:21.409 --> 00:01:24.509 Activities in Oklahoma Have Led to Unprecedented 00:01:24.509 --> 00:01:27.860 Seismic Hazard and Risk Level. Thanks. 00:01:28.460 --> 00:01:31.460 - Thank you very much for the introduction. 00:01:31.469 --> 00:01:36.219 So a little story before I start. When I decided to work on the topic 00:01:36.219 --> 00:01:40.419 of induced seismicity, I didn’t know where to start, what to read exactly, 00:01:40.419 --> 00:01:45.469 so my adviser, Paolo, told me to check what the USGS was doing. 00:01:45.469 --> 00:01:48.640 So I quickly found this exact seminar series, and I saw a presentation by 00:01:48.640 --> 00:01:52.759 Mark McClure about Oklahoma. And fast forward to four years later. 00:01:52.760 --> 00:01:58.460 I’m very happy to be here today, even remotely, presenting my results to you. 00:01:59.490 --> 00:02:03.600 So this study was conducted in close collaboration with Ellen Rathje from 00:02:03.620 --> 00:02:08.300 the University of Texas at Austin, Alexan Savvaidis from the Bureau of 00:02:08.300 --> 00:02:12.380 Economic Geology and the TexNet Seismological Network at UT-Austin, 00:02:12.380 --> 00:02:15.739 and Paolo Bazzurro from the University School of 00:02:15.739 --> 00:02:19.060 Advanced Studies, IUSS Pavia. 00:02:19.060 --> 00:02:22.689 So what’s special about Oklahoma? Well, as many of you know, 00:02:22.689 --> 00:02:25.980 over the last decade, Oklahoma has experienced unprecedented levels 00:02:25.980 --> 00:02:32.060 of seismicity about 100 times higher than the historical average. 00:02:32.060 --> 00:02:37.500 In 2015 alone, it experienced more than 700 earthquakes above magnitude 3 – 00:02:37.500 --> 00:02:42.200 more per unit area than California, to put things into perspective. 00:02:42.200 --> 00:02:45.790 And Oklahoma also had three events above magnitude 5 00:02:45.790 --> 00:02:50.379 in 2016 alone, causing some non-structural damages. 00:02:50.379 --> 00:02:55.099 The lack of seismic loads in the public – the lack of seismic loads 00:02:55.099 --> 00:02:59.370 in the building code has led to – has led to unrest in the public. 00:02:59.370 --> 00:03:04.709 And also housing prices have also dropped in high-seismicity areas. 00:03:04.709 --> 00:03:08.680 Here you can see the time history of the monthly seismicity rates, 00:03:08.680 --> 00:03:13.480 declustered and non-declustered, with almost no events before 2009 00:03:13.480 --> 00:03:17.080 and a rapid acceleration after 2014. 00:03:17.080 --> 00:03:20.810 The scientific community linked this unprecedented seismicity to the 00:03:20.810 --> 00:03:24.579 increasing wastewater disposal rates. I show them in blue here. 00:03:24.579 --> 00:03:27.749 And the regulators mandated their reduction in 2015, 00:03:27.749 --> 00:03:31.140 with the seismicity rates decreasing since then. 00:03:32.820 --> 00:03:36.379 Now, you can see that there’s a big time lag between the onset 00:03:36.379 --> 00:03:39.800 of the injection rates in the year 2000s and the onset 00:03:39.800 --> 00:03:43.049 of seismicity many, many years later. 00:03:43.049 --> 00:03:46.680 We will come back to that, but it’s important to point that out. 00:03:46.680 --> 00:03:49.520 Apart from wastewater disposal, hydraulic fracturing has also 00:03:49.530 --> 00:03:52.650 been linked with pockets of seismicity in Oklahoma. 00:03:52.650 --> 00:03:57.879 But the main driver is assumed to be wastewater disposal. 00:03:57.879 --> 00:04:01.079 So the outline of my talk is, introduction into the underlying 00:04:01.079 --> 00:04:04.300 physical processes, proposed earthquake recurrence model for 00:04:04.300 --> 00:04:09.150 wastewater disposal, hindcasting results of the earthquake recurrence model, 00:04:09.150 --> 00:04:13.029 hypothesis testing for the causal factors – in this case, wastewater 00:04:13.029 --> 00:04:16.600 disposal, forecasting results of the earthquake recurrence model. 00:04:16.600 --> 00:04:20.380 I’ll talk about some issues with the b values. 00:04:20.380 --> 00:04:23.520 Then I’ll show the results of a time-dependent seismic risk model – 00:04:23.520 --> 00:04:29.540 so, seismic hazard maps and other things. 00:04:29.540 --> 00:04:34.600 Sorry. I’m getting some control requests. 00:04:34.600 --> 00:04:39.520 And then I’ll talk about the seismic risk results, some loss maps, 00:04:39.530 --> 00:04:44.310 and some direct loss estimates. And finally conclusions, limitations, 00:04:44.310 --> 00:04:48.260 future research, and acknowledgements. 00:04:49.600 --> 00:04:53.340 Given that this is probably not the only Microsoft Teams event in your calendar 00:04:53.340 --> 00:04:57.710 today, I will not talk about quality control and injection data, 00:04:57.710 --> 00:05:00.870 catalog homogenization – we merged a lot of different catalogs 00:05:00.870 --> 00:05:04.139 together, and we homogenized the moment magnitudes. 00:05:04.139 --> 00:05:06.620 And I will not talk about declustering issues, even though 00:05:06.620 --> 00:05:10.370 they are very important, variations in magnitude of 00:05:10.370 --> 00:05:14.420 completeness, issues with M-max, and disaggregation of seismic hazard. 00:05:14.420 --> 00:05:16.480 Anyone who’s interested in these topics, 00:05:16.480 --> 00:05:20.060 they are covered in the publications. 00:05:20.060 --> 00:05:23.180 So the two main oil and gas activities related to induced seismicity 00:05:23.180 --> 00:05:27.280 in Oklahoma are hydraulic fracturing and wastewater disposal. 00:05:27.280 --> 00:05:29.980 Hydraulic fracturing is a drilling technique used to create 00:05:29.980 --> 00:05:34.720 new pathways for extracting oil and gas from deep underground. 00:05:37.460 --> 00:05:42.100 A new recent method in hydraulic fracturing applied to previously 00:05:42.110 --> 00:05:47.229 unproductive formations lead to massive quantities of coproduced 00:05:47.229 --> 00:05:50.960 flowback fluids, which we call wastewater, being pumped up 00:05:50.960 --> 00:05:55.389 together with the crude oil. These volumes are then disposed 00:05:55.389 --> 00:05:59.740 nearby through injection into deep underground reservoirs. 00:05:59.740 --> 00:06:05.310 The ratio of pumped wastewater to produced oil can be as high as 10 to 1. 00:06:06.400 --> 00:06:13.520 So, in the past 20 years, about 2.2 cubic kilometers of wastewater, or about a 00:06:13.530 --> 00:06:19.930 third of the volume of the Salton Sea, have been disposed under Oklahoma. 00:06:19.930 --> 00:06:25.490 With original stresses being near the strength limit of the crust, 00:06:25.490 --> 00:06:29.050 small stress perturbations like this can affect fault stability 00:06:29.050 --> 00:06:32.980 and trigger earthquakes, as you can see in this sketch. 00:06:32.980 --> 00:06:38.280 Increased fluid pore pressure due to wastewater injection reduces the 00:06:38.280 --> 00:06:42.830 normal stresses on the fault plane, bringing it closer to failure. 00:06:42.830 --> 00:06:47.419 So humans cause small pressure perturbations that destabilize faults, 00:06:47.420 --> 00:06:50.680 releasing energy that has accumulated over decades. 00:06:50.680 --> 00:06:53.780 Now, I should note that, once a rupture is initiated, 00:06:53.780 --> 00:06:59.640 all its characteristics are very similar to typical tectonic events of similar depth. 00:07:00.800 --> 00:07:04.520 Now, let me present the earthquake recurrence model that we developed 00:07:04.520 --> 00:07:07.620 for wastewater disposal. We modified the Gutenberg-Richter 00:07:07.620 --> 00:07:11.200 relation, allowing the a value to be time-dependent, 00:07:11.200 --> 00:07:15.879 as you can see here, or, rather, injected volume-dependent. 00:07:15.879 --> 00:07:19.700 To do that, we employed Shapiro’s physics-based seismogenic index 00:07:19.700 --> 00:07:23.460 model, the capital sigma, to quantify the number of earthquakes induced 00:07:23.460 --> 00:07:30.660 by injection for a unit volume of fluid. So, in blue is the injection-dependent 00:07:30.660 --> 00:07:35.380 term of the a value. The a value, in a way, is the whole parentheses here. 00:07:35.380 --> 00:07:39.620 Natural tectonic seismicity has, of course, not stopped, and that is 00:07:39.620 --> 00:07:44.260 why we added this a-tectonic term in Shapiro’s original equation. 00:07:45.520 --> 00:07:51.479 So the additional seismicity rates scale linearly with the injected volume, 00:07:51.479 --> 00:07:55.129 which is considered a proxy for pore pressure perturbation. 00:07:55.129 --> 00:07:58.900 And any references are always at the bottom of the slide. 00:07:59.980 --> 00:08:03.100 So the monthly injection and seismicity data between 00:08:03.100 --> 00:08:08.000 2006 and 2018 are aggregated into blocks. 00:08:08.000 --> 00:08:11.960 Lambda is the number of earthquakes above magnitude 3 00:08:11.960 --> 00:08:17.660 in that block and for that specific month, j. 00:08:17.660 --> 00:08:22.919 And v is the injected volume within that block for that specific month. 00:08:22.919 --> 00:08:25.300 And we do this process – we aggregate everything for 00:08:25.300 --> 00:08:28.939 every month and for every block. 00:08:28.939 --> 00:08:33.450 Acknowledging that injection in one block could affect the stress field in 00:08:33.450 --> 00:08:36.850 a neighboring block, as you can see from the sketch bubble here, 00:08:36.850 --> 00:08:43.840 we distribute the injected volumes of each well in space and time. 00:08:43.840 --> 00:08:46.540 We did that following Theis’ equation for pore pressure diffusion 00:08:46.540 --> 00:08:50.840 in non-leaky confined aquifers. So, for example, we can see 00:08:50.840 --> 00:08:56.510 how much of the volume from this well falls within the gray block with time. 00:08:56.510 --> 00:08:58.440 This is this function. 00:08:58.440 --> 00:09:02.340 So, as the pore pressure front reaches the block, 00:09:02.340 --> 00:09:05.520 moving away from the well, the percentage increases. 00:09:05.520 --> 00:09:09.230 And, when the pore pressure front has exceeded the middle of the block 00:09:09.230 --> 00:09:16.060 and is moving to the right, in a way, this function is decreasing again. 00:09:16.060 --> 00:09:21.390 The diffusivity value was chosen to be 0.3 square meters per second, 00:09:21.390 --> 00:09:25.010 following other studies. This value represents large-scale 00:09:25.010 --> 00:09:30.980 diffusivity that is often derived from seismicity migration analysis. 00:09:30.980 --> 00:09:35.370 It is a middle ground in sorts of ways between the diffusivity 00:09:35.370 --> 00:09:38.780 of the basement and the diffusivity of the aquifer. 00:09:40.230 --> 00:09:44.900 So the maps here show the effect of this process into the cumulative volumes 00:09:44.900 --> 00:09:51.100 with a 5-kilometer grid before and after we have done the distribution process. 00:09:51.100 --> 00:09:54.380 You can see that the distributed volumes at the bottom here are 00:09:54.380 --> 00:09:58.030 better correlated with the observed seismicity that is these 00:09:58.030 --> 00:10:02.700 black circles overlapped here. These are the earthquakes. 00:10:02.700 --> 00:10:06.970 So this approach works quite well producing reasonable estimates, 00:10:06.970 --> 00:10:10.210 as you can see from this map, even though it only needs 00:10:10.210 --> 00:10:14.290 one model parameter – this large-scale diffusivity. 00:10:14.290 --> 00:10:18.000 Alternative options include complex hydrologic models to 00:10:18.000 --> 00:10:22.130 do this pore pressure diffusion approach, which need more 00:10:22.130 --> 00:10:26.300 model parameters that are harder to constrain sometimes. 00:10:27.580 --> 00:10:31.080 As you can see from this plot, I have used all the disposal wells, 00:10:31.080 --> 00:10:35.980 targeting any formation, and not only the ones injecting into the Arbuckle. 00:10:37.180 --> 00:10:40.840 So, to address the [inaudible] time lag between injection seismicity that 00:10:40.840 --> 00:10:46.180 I mentioned earlier, we looked at rate-and-state friction, which dictates 00:10:46.190 --> 00:10:50.470 that changes in seismicity rate occur over characteristics time – 00:10:50.470 --> 00:10:56.160 over characteristic time scales that decrease with increasing stressing rate. 00:10:56.160 --> 00:11:01.440 For simplification, we assumed that the time lag here is inversely proportional 00:11:01.440 --> 00:11:06.330 to the monthly distributed injection rate at any given time – 00:11:06.330 --> 00:11:09.420 so this is time-dependent as well – and is controlled by this 00:11:09.420 --> 00:11:12.540 free model parameter, theta. 00:11:12.540 --> 00:11:14.940 V is the injection rate in each block, 00:11:14.940 --> 00:11:19.380 distributed spatiotemporally according to the previous slide. 00:11:19.380 --> 00:11:24.400 So, for this simplified linear case – this inversely proportional linear case, 00:11:24.410 --> 00:11:28.100 which is an assumption, a simple reservoir model by Norbeck 00:11:28.100 --> 00:11:32.590 and Rubinstein has theta being the product of the direct effect – 00:11:32.590 --> 00:11:35.870 the effective normal stress at seismogenic depth, the rock porosity, 00:11:35.870 --> 00:11:39.270 the thickness of the reservoir, the a of the investigated zone – 00:11:39.270 --> 00:11:43.370 so our block in this case, and the total compressibility. 00:11:43.370 --> 00:11:47.780 So there is some physical basis for this parameter, theta, according to 00:11:47.780 --> 00:11:53.160 some assumptions. But, even so, we will calibrate it empirically. 00:11:53.160 --> 00:11:56.960 So we’ll get this theta value from the data directly. 00:11:58.580 --> 00:12:01.120 So here’s an example for implementing the timeline. 00:12:01.120 --> 00:12:07.220 If the injection for the first volume for a given value of theta leads to 00:12:07.220 --> 00:12:11.260 a time lag of two months, then the third month, which has 00:12:11.260 --> 00:12:16.170 double the volume, will have half the time lag, so one month. 00:12:16.170 --> 00:12:22.070 So applying the time lag to the first row gives us this lag volume time history. 00:12:22.070 --> 00:12:25.900 So this volume moves two months to the future, two months, 00:12:25.900 --> 00:12:29.120 one month, one month, and we get this new time history. 00:12:29.120 --> 00:12:32.370 And now, this is what we actually use in the main equation. 00:12:32.370 --> 00:12:36.040 V is replaced by this lagged distributed injection rate. 00:12:36.040 --> 00:12:38.330 So it’s the same thing, but now we have done some 00:12:38.330 --> 00:12:42.000 extra computations to get this volume for this block. 00:12:43.000 --> 00:12:46.920 Now, seismicity in crystalline rocks is absent when the stress field 00:12:46.920 --> 00:12:50.640 is smaller than the maximum stress previously applied. 00:12:50.640 --> 00:12:53.540 As a result, when simulating seismicity, we ignore the 00:12:53.540 --> 00:12:57.580 lagged injection rate that followed this global maximum. 00:12:57.580 --> 00:13:02.000 Temporary injection reversals like these ones here are 00:13:02.010 --> 00:13:04.460 not considered long enough and strong enough for the 00:13:04.460 --> 00:13:08.660 pore pressure to reduce significantly across the entire block. 00:13:10.300 --> 00:13:13.640 So this constraint that I just said can help us explain why, 00:13:13.640 --> 00:13:18.240 after 2015, the seismicity dropped by 65% even though 00:13:18.240 --> 00:13:22.000 the injection rates dropped only by 40%. 00:13:23.640 --> 00:13:27.420 Now, we have declustered the earthquake catalog using 00:13:27.420 --> 00:13:30.920 the Reasenberg method. This is important to note. 00:13:30.920 --> 00:13:35.120 But, even though the aftershocks are, in principle, removed, we are still 00:13:35.120 --> 00:13:39.560 dealing with a non-homogeneous Poisson process because each month 00:13:39.560 --> 00:13:43.460 the earthquake rate depends on the injection rate. 00:13:43.460 --> 00:13:47.220 Since we have this equation linking injection with seismicity, we can 00:13:47.230 --> 00:13:52.200 model the whole process as a homogeneous Poisson process with 00:13:52.200 --> 00:13:56.610 varying mean rate, lambda of t. So lambda of t is changing every 00:13:56.610 --> 00:14:01.500 month, but within each month, it is a Poisson process incrementally. 00:14:02.540 --> 00:14:05.840 Having defined the statistical distribution for lambda, 00:14:05.840 --> 00:14:09.600 we can use maximum likelihood to calibrate our model given the 00:14:09.600 --> 00:14:12.560 observed injection and seismicity data. 00:14:12.560 --> 00:14:17.490 We essentially have two – three parameters – theta and sigma that we 00:14:17.490 --> 00:14:24.790 need to calibrate using the seismicity and injection data after 2006. 00:14:24.790 --> 00:14:27.890 The background tectonic rate, a-tectonic, is estimated from 00:14:27.890 --> 00:14:31.970 the original catalog before 2009, and it’s very low. 00:14:31.970 --> 00:14:34.530 And the b value is not important for now since we are using 00:14:34.530 --> 00:14:39.260 all the events above magnitude 3. We’ll come back to the b value later. 00:14:40.100 --> 00:14:44.460 So some early results. First, we defined two randomly 00:14:44.470 --> 00:14:48.860 placed grids of 20-kilometer blocks – this is a sketch – and looked at 00:14:48.860 --> 00:14:52.860 the results of our simulations statewide. 00:14:53.390 --> 00:14:57.780 The blue lines represent the wastewater injection rates, distributed and not. 00:14:57.790 --> 00:15:01.780 The red line is our simulated catalog. And the black line is the observed 00:15:01.780 --> 00:15:06.160 declustered catalog. As you can see, the results, the red lines, 00:15:06.160 --> 00:15:10.860 are quite sensitive to the grid selection. So changing the grid, shifting it 00:15:10.860 --> 00:15:16.380 10 kilometers to the north, to the south, or to the east, affects our fit. 00:15:16.380 --> 00:15:21.140 To ensure stable results, we resorted to a novel two-scale fitting procedure 00:15:21.140 --> 00:15:25.040 because larger grids provide richer data sets for model calibration. 00:15:25.040 --> 00:15:27.800 So we like this 20-kilometer scale. 00:15:27.800 --> 00:15:31.170 But also the assumptions behind the seismogenic index model are 00:15:31.170 --> 00:15:37.070 more applicable to smaller scales. Because it’s for a single injection well. 00:15:37.070 --> 00:15:39.810 That’s what the original development. 00:15:39.810 --> 00:15:45.790 So, to address these instabilities, we defined a dense 5-kilometer grid, 00:15:45.790 --> 00:15:51.830 and for each 5-kilometer block, we defined 16 different 20-kilometer 00:15:51.830 --> 00:15:55.200 blocks that enclosed it. Here is some examples of the 00:15:55.200 --> 00:15:58.040 20-kilometer blocks surrounding the 5-kilometer block. 00:15:58.040 --> 00:16:01.880 There are seven more examples, but I think you get the idea of how 00:16:01.880 --> 00:16:08.880 we created these 16 20-kilometer blocks around the 5-kilometer block. 00:16:08.880 --> 00:16:11.560 So now, bear with me for this slide. 00:16:11.560 --> 00:16:17.320 The regression of the model parameters is always done at the 20-kilometer scale, 00:16:17.320 --> 00:16:22.980 maximizing the total likelihood, L-i, of these 20-kilometer blocks. 00:16:22.990 --> 00:16:25.750 So maximizing this, we get sigma and theta values 00:16:25.750 --> 00:16:28.120 for each 20-kilometer block. 00:16:28.120 --> 00:16:32.660 We then apply these 16 pairs of sigma and theta to the single 00:16:32.670 --> 00:16:37.300 5-kilometer block in question. And we obtain 16 likelihood values 00:16:37.300 --> 00:16:41.820 for this 5-kilometer block, which we call L-i of 5. 00:16:42.620 --> 00:16:46.380 Now, we have 16 pairs. We have to choose one of them. 00:16:46.380 --> 00:16:51.100 The chosen pair is the one maximizing this function. 00:16:51.100 --> 00:16:56.660 I told you what L-i of 5 is, L-i of 20. And this g function is a correction 00:16:56.670 --> 00:17:00.070 factor because some of the 20-kilometer blocks have more 00:17:00.070 --> 00:17:04.790 earthquakes than the others, biasing the Poisson likelihoods to lower values. 00:17:04.790 --> 00:17:07.720 So we have to correct for the sample size. 00:17:10.020 --> 00:17:13.780 So the regression is always done at the 20-kilometer scale, 00:17:13.790 --> 00:17:17.079 but the model calibration takes into account both scales. 00:17:17.079 --> 00:17:20.470 And you end up with a final spatial resolution of 5 kilometers. 00:17:20.470 --> 00:17:23.169 Because, in the end, we apply the final sigma 00:17:23.169 --> 00:17:26.740 and theta into every 5-kilometer block. 00:17:28.620 --> 00:17:33.140 So let’s look at the goodness of fit of our simulated time histories 00:17:33.140 --> 00:17:38.220 with a two-scale – with a multi-scale calibration 00:17:38.230 --> 00:17:44.230 using all the data between 2006 and the end of the 2017 to do this calibration. 00:17:44.230 --> 00:17:47.170 So the blue lines are the injection rate. Same idea as before. 00:17:47.170 --> 00:17:50.450 Red is the simulation. Black is the declustered catalog. 00:17:50.450 --> 00:17:55.200 These dotted lines here are the actual monthly values – the black and red. 00:17:55.200 --> 00:18:00.250 But, from now on, I will only show this smooth time history here just for 00:18:00.250 --> 00:18:06.060 illustration purposes. It’s easier to follow the solid red and black lines. 00:18:07.169 --> 00:18:10.820 So the simulation, as you can see, the red line is in very good agreement 00:18:10.820 --> 00:18:14.180 with observed earthquake rates, especially after 2014. 00:18:14.180 --> 00:18:16.140 We get a very good agreement here. 00:18:16.140 --> 00:18:20.570 Where we need to capture both the rise and the decay of the seismicity 00:18:20.570 --> 00:18:23.730 very well, even though the original seismogenic index model 00:18:23.730 --> 00:18:27.400 is applicable only to increasing injection rates. 00:18:28.440 --> 00:18:31.620 The two free parameters, sigma and theta, span several orders 00:18:31.629 --> 00:18:37.450 of magnitude, you can see here. And indicating a high degree of spatial 00:18:37.450 --> 00:18:40.120 heterogeneity in agreement with other studies. 00:18:40.120 --> 00:18:43.980 So it’s important to have this 5-kilometer resolution – 00:18:43.980 --> 00:18:47.780 this fine resolution because things are changing spatially. 00:18:47.780 --> 00:18:52.100 And, if you remember, I said that theta – theta – I’m Greek. 00:18:52.100 --> 00:18:56.450 I call it theta – can be approximated as the product 00:18:56.450 --> 00:18:59.680 of six hydromechanical parameters. 00:18:59.680 --> 00:19:04.440 What’s interesting is that the median theta value across all the blocks is 00:19:04.440 --> 00:19:07.920 consistent with the commonly adopted regional values of these 00:19:07.920 --> 00:19:11.679 six hydromechanical parameters that comprise it. 00:19:11.679 --> 00:19:15.330 So literature values agree with what we get here for these 00:19:15.330 --> 00:19:18.429 hydromechanical parameters. As a product, it agrees. 00:19:18.429 --> 00:19:21.080 So this is some sort of encouragement that our 00:19:21.080 --> 00:19:26.220 time lag approach has some merit. 00:19:26.220 --> 00:19:28.860 Now let’s look at the result of a more regional scale. 00:19:28.860 --> 00:19:31.600 You can see three regions here – northwest or north-central 00:19:31.600 --> 00:19:34.090 and central Oklahoma. 00:19:34.090 --> 00:19:37.600 The model performs very well in the two northern regions, 00:19:37.600 --> 00:19:41.340 while it underestimates a bit the seismicity between 00:19:41.340 --> 00:19:45.100 2014 and 2016 here in central Oklahoma. 00:19:46.250 --> 00:19:51.000 We will come back to that later with some ideas, but one explanation, 00:19:51.000 --> 00:19:55.460 of course, is that this area is affected by complex earthquake interactions 00:19:55.460 --> 00:19:57.480 that our declustering cannot capture. 00:19:57.480 --> 00:20:02.260 There are some papers indicating that this area is affecting by remote 00:20:02.260 --> 00:20:07.919 triggering from very large events thousands of kilometers away. 00:20:07.920 --> 00:20:11.380 So this area needs your attention. 00:20:12.420 --> 00:20:16.400 So this slide is exactly the same as the one before, but I have now 00:20:16.400 --> 00:20:22.420 added a simulation that I did with my injection data using the model 00:20:22.429 --> 00:20:27.169 of Norbeck and Rubinstein, which also uses disposal rates. 00:20:27.169 --> 00:20:33.549 So the green is this model, and it’s not able to capture the 00:20:33.549 --> 00:20:39.080 rapid rise and decay of the seismicity that well, at least for 00:20:39.080 --> 00:20:43.980 this specific case and their literature values. 00:20:45.220 --> 00:20:49.200 Now let’s look at three smaller scales to demonstrate 00:20:49.200 --> 00:20:53.100 the consistent performance of the model across spatial scales. 00:20:53.100 --> 00:20:57.880 So you see the three polygons here – the three regions that we show. 00:20:57.880 --> 00:21:03.100 So Oklahoma City, which is the capital, I think – it’s the biggest city, for sure – 00:21:03.100 --> 00:21:07.480 you can see these – barely see this light blue line here, which is the actual 00:21:07.480 --> 00:21:12.220 injected volume within Polygon 1, and it’s, like, almost negligible. 00:21:12.220 --> 00:21:21.080 But wells further south of Polygon 1 contribute to these distributed injection 00:21:21.080 --> 00:21:27.169 rates in the solid blue, and they help us do the successful simulation. 00:21:27.169 --> 00:21:30.440 So it’s important to do this volume distribution because sometimes 00:21:30.440 --> 00:21:34.080 neighboring wells are affecting what’s going on. 00:21:35.980 --> 00:21:38.580 Now, in Cushing and Pawnee, where earthquakes above 00:21:38.580 --> 00:21:44.720 magnitude 5 occurred, our simulation has, again, a rather good fit. 00:21:44.720 --> 00:21:48.080 So, at small scales, regional scales, statewide, 00:21:48.080 --> 00:21:51.480 the performance is consistent. 00:21:51.480 --> 00:21:55.130 Now, we want to test the statistical significance of the prevailing 00:21:55.130 --> 00:21:57.610 hypothesis that wastewater injection is, indeed, 00:21:57.610 --> 00:22:01.250 the driver of the seismicity rates. 00:22:01.250 --> 00:22:05.559 So far, we have assumed that this is true everywhere in the state. 00:22:05.559 --> 00:22:09.480 Two hypotheses are constructed to do that – a null hypothesis that 00:22:09.480 --> 00:22:12.669 assumes that there is no relationship between injection and seismicity – 00:22:12.669 --> 00:22:14.649 everything is tectonic. 00:22:14.649 --> 00:22:18.549 And we have the alternative hypothesis, with total likelihood, L-1, 00:22:18.549 --> 00:22:21.909 which does assume a relationship, is the equation that I showed before 00:22:21.909 --> 00:22:24.340 with the injection rates here. 00:22:25.320 --> 00:22:29.580 So, as I said, we assumed an incrementally homogeneous 00:22:29.580 --> 00:22:35.820 Poisson process for each month. So, for each block, for each month, 00:22:35.820 --> 00:22:41.379 we get a likelihood by comparing the observed seismicity rates 00:22:41.380 --> 00:22:45.160 from the catalog with the model seismicity rate, 00:22:45.160 --> 00:22:50.480 either for the null hypothesis or for the alternative hypothesis. 00:22:50.480 --> 00:22:54.539 So we can compare the observed with the model for every month. 00:22:54.539 --> 00:22:58.779 We get likelihoods for every month. And, when we get the product of 00:22:58.779 --> 00:23:04.029 all these monthly likelihoods for each block, we get the total likelihood, 00:23:04.029 --> 00:23:09.159 L-zero, if we compare this pair, and L-1 if we compare this pair. 00:23:09.160 --> 00:23:16.640 We also define R as L-1 over L-zero. This is the likelihood ratio. 00:23:16.640 --> 00:23:21.309 So, although our values larger than 1 indicate that the alternative hypothesis 00:23:21.309 --> 00:23:26.740 is, indeed, more likely, this criterion is not sufficient to statistically reject 00:23:26.740 --> 00:23:31.020 the null hypothesis. We have to do more work to achieve that. 00:23:31.020 --> 00:23:36.340 So we need a reference statistical distribution for the ratio, R, 00:23:36.350 --> 00:23:41.409 for which the null hypothesis is true. The reference distribution enables 00:23:41.409 --> 00:23:44.320 a direct comparison of the two hypotheses, even though the 00:23:44.320 --> 00:23:48.280 alternative hypothesis has higher model complexity. 00:23:49.780 --> 00:23:53.499 So to generate this reference distribution, we followed the 00:23:53.500 --> 00:24:00.820 framework of McClure et al. For each block, we fed our proposed – 00:24:00.820 --> 00:24:06.080 sorry, again, I have some requests for control. Okay. 00:24:06.080 --> 00:24:11.039 So, for each block, we fed our proposed model with real seismicity 00:24:11.040 --> 00:24:18.140 data and fake, but realistic, injection data, producing 100 synthetic 00:24:18.140 --> 00:24:20.520 R values that fulfill the null hypothesis. 00:24:20.520 --> 00:24:23.900 They fulfill the null hypothesis because the injection data are fake. 00:24:23.909 --> 00:24:29.230 So they couldn’t have affected these seismicity rates in real life. 00:24:29.230 --> 00:24:33.890 So we call these R values, the ones from the synthetic data, 00:24:33.890 --> 00:24:37.710 the R-null distribution. This distribution reflects how likely 00:24:37.710 --> 00:24:43.540 it is for our earthquake recurrence model to find purely coincidental 00:24:43.540 --> 00:24:49.020 correlation between the observed seismicity and random injection data. 00:24:49.020 --> 00:24:53.200 And, by comparing this R-null distribution for each block with 00:24:53.210 --> 00:24:59.970 the unique real R value that we get with the real injection data 00:24:59.970 --> 00:25:03.980 for that block, we can compute the p values. 00:25:03.980 --> 00:25:11.539 So, for example, if four of the 100 R-null values are larger than 00:25:11.540 --> 00:25:18.889 the one real value, R-real, then the p value for that block is 0.04. 00:25:20.260 --> 00:25:23.680 Now, we do that for each 20-kilometer block. 00:25:23.680 --> 00:25:26.400 Because we start with the 20-kilometer blocks. 00:25:26.400 --> 00:25:31.049 And then post-processing the results to get a finer resolution of 5 kilometers, 00:25:31.049 --> 00:25:32.850 and we end up with this map. 00:25:32.850 --> 00:25:38.100 So let me explain that a bit more. The p value of each seismogenic 00:25:38.100 --> 00:25:44.360 5-kilometer block here is the median of the p values that we get from the 16 00:25:44.360 --> 00:25:48.820 20-kilometer blocks that enclose each 5-kilometer block. 00:25:50.620 --> 00:25:57.820 So this is the legend for the colors. 76% of the blocks have a p value less 00:25:57.820 --> 00:26:02.600 than 0.05, indicating that they are – we are very confident that the 00:26:02.610 --> 00:26:07.759 wastewater disposal explains very well the observed seismicity rates and much, 00:26:07.759 --> 00:26:11.269 much better than the tectonic assumption. 00:26:11.269 --> 00:26:16.149 So these red blocks hosted 84% of events above magnitude 3 00:26:16.149 --> 00:26:20.980 and 91% of events above magnitude 4, including the four largest earthquakes. 00:26:20.980 --> 00:26:24.240 This observation strongly supports the consensus of the scientific 00:26:24.240 --> 00:26:28.740 community that the main cause of the rise in seismicity is, indeed, 00:26:28.740 --> 00:26:33.180 wastewater disposal. But now we finally know, 00:26:33.180 --> 00:26:39.420 and we have a quantitative measure to support this claim. 00:26:39.420 --> 00:26:43.860 Our elevated b values in the south and southeastern Oklahoma – you 00:26:43.860 --> 00:26:46.880 can see here in the SCOOP/STACK and the Anadarko Basin. 00:26:46.880 --> 00:26:50.419 I’m in agreement with Skoumal et al., who link the seismicity there – 00:26:50.419 --> 00:26:54.860 these are their boxes – directly with hydraulic fracturing activities. 00:26:54.860 --> 00:27:00.100 So we shouldn’t see red colors here because there the seismicity is caused 00:27:00.100 --> 00:27:03.639 by hydraulic fracturing and not wastewater disposal. 00:27:03.639 --> 00:27:07.549 This provides an external validation for our map from 00:27:07.549 --> 00:27:11.280 a study with different input data sets and methods. 00:27:11.900 --> 00:27:16.340 Now, there are some high p values in the center of the state here. 00:27:16.350 --> 00:27:19.720 And these p values remain similar if I limit the analysis to 00:27:19.720 --> 00:27:22.980 only the Arbuckle wells. So this is not the issue. 00:27:22.980 --> 00:27:26.059 So what’s going on there? Well, triggering via cascading 00:27:26.059 --> 00:27:28.600 [inaudible] interactions, as I said before, might be 00:27:28.600 --> 00:27:31.600 an alternative cause of the seismicity there. 00:27:31.600 --> 00:27:37.980 Alternatively, this area might require a detailed hydrogeologic model 00:27:37.980 --> 00:27:41.880 like the one from Langenbruch, Zoback, and Weingarten. 00:27:41.960 --> 00:27:45.760 to explain some [inaudible] pore pressure changes. 00:27:45.760 --> 00:27:49.040 However, another scenario worth investigating in the future is that the 00:27:49.049 --> 00:27:54.990 seismicity in that area is at least partly associated with fluid extraction. 00:27:54.990 --> 00:28:00.239 So the abundance of wastewater in that area has led to significant portions 00:28:00.239 --> 00:28:04.769 of it being disposed in neighboring areas dozens of kilometers away. 00:28:04.769 --> 00:28:07.850 So they take wastewater from here. They cannot dispose it there. 00:28:07.850 --> 00:28:12.240 The formations are not favorable. And they move it dozens of kilometers 00:28:12.240 --> 00:28:16.620 away. So they take the fluids from there and put it elsewhere. 00:28:16.630 --> 00:28:19.350 So, over time, this is a fluid extraction process 00:28:19.350 --> 00:28:21.879 which can itself cause earthquake activity. 00:28:21.880 --> 00:28:26.640 So I’m not saying this is the answer, but we should at least consider it. 00:28:27.480 --> 00:28:31.419 We consider this white polygon enclosed in the red blocks as the 00:28:31.420 --> 00:28:34.400 area of interest for wastewater disposal. 00:28:34.400 --> 00:28:39.659 Wastewater disposal is not the driver of seismicity for the rest of the state. 00:28:39.659 --> 00:28:46.070 And this is important to understand because applying induced models 00:28:46.070 --> 00:28:49.529 for wastewater injection where the seismicity is driven by other 00:28:49.529 --> 00:28:54.489 processes leads to biases. So before, I showed you this plot 00:28:54.489 --> 00:28:57.179 for the simulation according to Norbeck and Rubinstein. 00:28:57.179 --> 00:29:02.940 This plot is a bit unfair because I simulated the whole state seismicity 00:29:02.940 --> 00:29:07.799 using all the volumes, even here, where hydraulic fracturing is the answer. 00:29:07.799 --> 00:29:20.149 So that’s why we get this fit. But when we actually limit the area of 00:29:20.149 --> 00:29:24.789 investigation within the AOI, we see a much better fit because this model is 00:29:24.789 --> 00:29:30.080 more applicable within the AOI where wastewater injection is the driver. 00:29:30.080 --> 00:29:34.080 So it’s important to quantify the causal factor before we apply 00:29:34.080 --> 00:29:38.480 a model to do forecasting or hindcasting. 00:29:40.100 --> 00:29:43.620 Now let’s see how the model performs when we calibrate it using the data 00:29:43.629 --> 00:29:47.759 until the end of 2014. We can see the calibration periods here. 00:29:47.759 --> 00:29:51.669 And then simulate the earthquake rates for the next three years, 00:29:51.669 --> 00:29:56.830 given that we know that the injection rates that followed are these ones. 00:29:56.830 --> 00:30:02.549 We see that the simulation – the yellow line here, here, here, and here – is still 00:30:02.549 --> 00:30:06.309 in good agreement with the observations both statewide and regionally. 00:30:06.309 --> 00:30:08.759 We can see the regions here. 00:30:08.760 --> 00:30:15.620 I also show in magenta, here, here, and here, the annually updated 00:30:15.620 --> 00:30:23.540 forecasts of the USGS for induced seismicity. They simply lag the mean 00:30:23.549 --> 00:30:27.119 seismicity rate of the previous year to the next. 00:30:27.119 --> 00:30:31.529 But I applied them to my declustered catalog. 00:30:31.529 --> 00:30:35.539 So I applied the concept of the USGS to my catalog. 00:30:35.539 --> 00:30:38.090 Also for the Norbeck approach, for example, I applied it to 00:30:38.090 --> 00:30:41.460 my injection data. So this USGS approach 00:30:41.460 --> 00:30:44.060 remains agnostic to variations in the disposal rates. 00:30:44.060 --> 00:30:47.540 Because it only looks at the seismicity of the previous year. 00:30:50.200 --> 00:30:55.620 Now here you can see how sigma and theta – this is sigma, this is theta – 00:30:55.620 --> 00:30:59.440 vary when we calibrate them either through the end of 2014 00:30:59.450 --> 00:31:01.539 or through the end of 2017. 00:31:01.539 --> 00:31:04.539 So you can see that they are fairly stable when you see red colors here, 00:31:04.539 --> 00:31:08.389 you see them here, blue here – you see the blue here. 00:31:08.389 --> 00:31:10.840 Of course, this has more blocks because there are more aspects 00:31:10.840 --> 00:31:14.850 around the state by the end of 2017. 00:31:14.850 --> 00:31:18.820 But overall, we see stable parameters, as we should expect. 00:31:21.500 --> 00:31:24.879 So now, could the model have predicted the rapid change 00:31:24.880 --> 00:31:30.020 in seismicity early in advance before the issue became widespread, 00:31:30.020 --> 00:31:32.660 say at the end of 2011? 00:31:32.660 --> 00:31:36.960 So, if you look at this map, there are much fewer data to work with back 00:31:36.970 --> 00:31:41.559 in 2011 – these blue earthquakes are what we had back in 2011. 00:31:41.559 --> 00:31:45.330 And several areas in the north had no seismicity back then 00:31:45.330 --> 00:31:49.710 in 2011 but later hosted a lot of earthquakes in gray. 00:31:49.710 --> 00:31:54.789 So we have to fill gaps here that we have no sigma and theta to work with. 00:31:54.789 --> 00:31:57.580 And we do spatial extrapolation of sigma and theta to be able to 00:31:57.580 --> 00:31:59.929 capture this northern region. So we take the [inaudible] 00:31:59.929 --> 00:32:05.960 from here and we apply it also in the north. 00:32:05.960 --> 00:32:09.309 So despite the scarcity of data – this is the calibration period – 00:32:09.309 --> 00:32:14.269 our forecast captures both the timing and the amplitude of the seismicity rate 00:32:14.269 --> 00:32:19.450 peak remarkably well. Because we did the spatial extrapolation, which is 00:32:19.450 --> 00:32:23.860 a crude assumption, we have some overestimated rates before 2014. 00:32:25.520 --> 00:32:30.480 Now, we also calibrated the model through the end of 2017 00:32:30.490 --> 00:32:34.619 to do some forecasting. So we forecasted the seismicity rates, 00:32:34.619 --> 00:32:38.169 assuming that the injection rates will remain stable and equal to 00:32:38.169 --> 00:32:42.359 their December 2017 value. In reality, the injection rates are 00:32:42.360 --> 00:32:46.200 still in decline, so therefore our forecast is conservative. 00:32:46.900 --> 00:32:51.600 So, contrary to the model of Norbeck and Rubinstein, our model predicts 00:32:51.600 --> 00:32:56.320 decreasing seismicity in red under stable injection rates, 00:32:56.320 --> 00:33:02.289 with the seismicity reaching pre-2009 levels after 2025. 00:33:02.289 --> 00:33:04.239 So this is the forecast. 00:33:05.220 --> 00:33:09.559 Now, all the results so far were based on the seismicity rate 00:33:09.559 --> 00:33:11.710 of earthquakes above magnitude 3. 00:33:11.710 --> 00:33:16.489 Usually earthquakes below magnitude 4.5 do not cause damages, 00:33:16.489 --> 00:33:18.549 and therefore we need to use the Gutenberg-Richter 00:33:18.549 --> 00:33:23.320 b value to extrapolate to larger magnitudes. 00:33:23.320 --> 00:33:26.040 If we fit the Gutenberg-Richter curve, as you can see here, 00:33:26.040 --> 00:33:33.120 to the entire declustered catalog, the fit of the larger magnitudes is very poor. 00:33:34.500 --> 00:33:38.520 However, not all declustered events are created equal. 00:33:38.520 --> 00:33:41.880 What do I mean? If we look at this Prague sequence, 00:33:41.899 --> 00:33:46.360 we can see that there are several isolated events in these orange 00:33:46.360 --> 00:33:50.860 polygons. And there is one blue main shock – the biggest event. 00:33:50.860 --> 00:33:55.370 And these are the sequence – this is the sequence of the main shock in black. 00:33:55.370 --> 00:34:01.350 So all the black events within the blue cloud are removed, 00:34:01.350 --> 00:34:04.409 and only the main shock remains during the declustering. 00:34:04.409 --> 00:34:06.639 So, when we decluster, all the black events are removed. 00:34:06.640 --> 00:34:11.780 We stay with the blue main shock and all the isolated events. 00:34:13.320 --> 00:34:18.100 However, if we plot the isolated events and all of the main shocks – 00:34:18.110 --> 00:34:20.470 the clustered main shocks – so the clustered main shock 00:34:20.470 --> 00:34:24.580 is the largest event of a full cluster. 00:34:26.190 --> 00:34:31.000 So you can see that the b value of the clustered main shocks is close to 1, 00:34:31.010 --> 00:34:33.300 in agreement with the historical data, 00:34:33.300 --> 00:34:36.660 and also the larger magnitudes are well-covered. 00:34:37.640 --> 00:34:42.420 So discarding the isolated events, and just focusing on the main shocks, 00:34:42.420 --> 00:34:45.760 it is justified because we – because the isolated events 00:34:45.760 --> 00:34:49.400 do not grow to magnitudes above 4.2. 00:34:49.400 --> 00:34:51.780 And we are interested in potentially damaging events, 00:34:51.780 --> 00:34:55.360 which tend to be larger than 4.2. 00:34:56.360 --> 00:35:00.140 So these are the current probabilities of potentially damaging events 00:35:00.140 --> 00:35:04.080 for different years. And you can see that, for 2020, 00:35:04.080 --> 00:35:09.550 they are about half of 2015, while in 2025, we expect them 00:35:09.550 --> 00:35:13.180 to be an order of magnitude above the tectonic rates. 00:35:13.180 --> 00:35:15.260 But, again, these are conservative estimates 00:35:15.260 --> 00:35:18.140 because the injection rates, now they are decreasing. 00:35:18.140 --> 00:35:21.681 So, before I move into the seismic hazard risk results, I would like 00:35:21.681 --> 00:35:25.590 to show this table summarizing the six available models in the 00:35:25.590 --> 00:35:31.180 literature for wastewater injection seismicity in Oklahoma. 00:35:31.180 --> 00:35:34.760 This is our model, and these are the other three. 00:35:34.760 --> 00:35:39.300 All these models use the injection data. We believe that our model is 00:35:39.300 --> 00:35:42.070 an attractive alternative because of its relative simplicity. 00:35:42.070 --> 00:35:48.510 It has fewer free parameters, its fine resolution down to 5 kilometers. 00:35:48.510 --> 00:35:55.190 and its ability to model the lack of seismicity before 2009, before arbitrary 00:35:55.190 --> 00:35:59.510 thresholds for earthquake triggering. So this time lag issue that I talked about 00:35:59.510 --> 00:36:03.480 is important, and there is no consensus on how to address it. 00:36:03.480 --> 00:36:07.600 Some models assume that there is a volume or pressure threshold – 00:36:07.600 --> 00:36:11.660 so, until you reach the pressure threshold, you see no earthquakes. 00:36:11.660 --> 00:36:17.420 Others use nonlinear functions of the seismogenic index equation. 00:36:17.420 --> 00:36:22.280 Others use time lag. Others have no time lag. 00:36:22.280 --> 00:36:28.500 So it’s important to figure out why this time lag exists. 00:36:28.500 --> 00:36:33.130 And these thresholds are very sensitive to the magnitude of completeness 00:36:33.130 --> 00:36:36.660 because, okay, you might see no earthquakes up to 2009 00:36:36.660 --> 00:36:41.200 for magnitude 3 and above, but if, in some magic way, you could see 00:36:41.210 --> 00:36:45.130 magnitude zero, you would of course see them way before 2009. 00:36:45.130 --> 00:36:50.810 So having these thresholds arbitrarily set to a value is quite sensitive 00:36:50.810 --> 00:36:53.310 to the magnitude of completeness. And when we extrapolate to 00:36:53.310 --> 00:36:56.570 larger magnitudes, this becomes a very big issue. 00:36:56.570 --> 00:36:59.160 And also the b values change. The declustering approaches 00:36:59.160 --> 00:37:04.260 are very different. You can see a lot of differences in this table. 00:37:04.260 --> 00:37:07.930 I’m not saying that one is exactly for sure correct, 00:37:07.930 --> 00:37:12.240 but it’s important to highlight these differences. 00:37:12.240 --> 00:37:17.190 Now, now that we know the current rates of potentially damaging events, 00:37:17.190 --> 00:37:20.350 we can assess the economic losses in a probabilistic way. 00:37:20.350 --> 00:37:24.810 So, for the simple case of a single house next to a fault, the risk 00:37:24.810 --> 00:37:28.360 calculation would be as follows. We forecast the location and 00:37:28.360 --> 00:37:32.240 the size of an earthquake. We convert that to ground motion 00:37:32.240 --> 00:37:35.920 at the house’s location using a suitable ground motion model 00:37:35.930 --> 00:37:40.040 and a Vs30 value of the site. And this is the seismic hazard [inaudible]. 00:37:40.040 --> 00:37:44.740 We get the ground motion at the site. This is their seismic hazard component. 00:37:44.740 --> 00:37:48.600 Then, for this given ground motion we just computed, we know, 00:37:48.600 --> 00:37:51.810 from another empirical model, what damages we should expect 00:37:51.810 --> 00:37:54.900 at this house and how much they will cost to repair. 00:37:54.900 --> 00:37:58.620 And when we do that for all the ruptures and all the buildings, 00:37:58.620 --> 00:38:02.390 we get the total losses for our region. 00:38:02.390 --> 00:38:04.570 The proposed time-dependent seismic hazard model has 00:38:04.570 --> 00:38:08.550 four main components. The spatially distributed seismicity 00:38:08.550 --> 00:38:11.000 rates coming out of the earthquake recurrence model. 00:38:11.000 --> 00:38:16.200 A statewide b value of 1.05. A Vs30 model of Oklahoma. 00:38:16.200 --> 00:38:20.180 And then region-specific ground motion model calibrated on 00:38:20.180 --> 00:38:24.350 regional induced earthquakes. You can see the references here. 00:38:24.350 --> 00:38:28.570 So M-min was set to 4.5 given the lack of seismic revisions in the code. 00:38:28.570 --> 00:38:33.120 M-max was set to 7 given the extent of the regional fold network. 00:38:33.120 --> 00:38:38.940 And, since we’re aiming at one-year estimate, the mean annual rate 00:38:38.940 --> 00:38:44.710 in each block was simply taken as the sum of the [inaudible] seismicity rates 00:38:44.710 --> 00:38:50.800 in that block over that specific year. This is the function. This is a 00:38:50.800 --> 00:38:55.110 commonly used approximation. And, of course, the mean annual rate 00:38:55.110 --> 00:39:00.410 is changing every year because the injection is changing every month. 00:39:00.410 --> 00:39:04.270 Here we show one-year hazard maps at 10% annual probability of exceedance 00:39:04.270 --> 00:39:08.840 for the most relevant intensity measure, SA 0.3 seconds. 00:39:09.320 --> 00:39:15.420 Compared to 2015, the forecasted values for 2020 are about 3 times 00:39:15.420 --> 00:39:19.120 smaller, even though the injection rates are assumed to be only 00:39:19.120 --> 00:39:21.700 50% smaller than in 2015. 00:39:21.700 --> 00:39:27.200 So this is an important observation. The color bar is the same here. 00:39:28.390 --> 00:39:32.200 To be able to compare our results with the long-term hazard maps 00:39:32.210 --> 00:39:37.610 of the USGS used in the building code, I computed maps for SA 0.2 seconds 00:39:37.610 --> 00:39:43.690 and 475 years return period. So this is the map that I computed. 00:39:43.690 --> 00:39:49.190 In 2015, the seismic hazard estimates at several sites in Oklahoma were 00:39:49.190 --> 00:39:58.260 higher than in most parts of the San Andreas Fault, which is quite striking. 00:39:59.320 --> 00:40:05.860 And, for 2012, the hazard in most of the state is estimated to be double 00:40:05.870 --> 00:40:10.110 the regional long-term average according to the latest national maps. 00:40:10.110 --> 00:40:14.340 Please note that this color bar is different here, but these are the same. 00:40:16.960 --> 00:40:19.740 To the risk model now – the exposure data. 00:40:19.740 --> 00:40:23.320 The total replacement cost of the building inventory in Oklahoma 00:40:23.320 --> 00:40:28.380 is 640 billion USD, which is about 3 times the annual GDP, 00:40:28.380 --> 00:40:32.780 with only 87 billion attributed to structural elements. 00:40:32.780 --> 00:40:36.730 The rest is non-structural elements and contents. 00:40:37.460 --> 00:40:42.380 61% of the buildings are made from wood, 26 are unreinforced masonry, 00:40:42.380 --> 00:40:47.570 and 10% are mobile homes. Almost all buildings are residential. 00:40:47.570 --> 00:40:51.690 And all buildings were designed with little to no seismic provisions. 00:40:51.690 --> 00:40:54.040 The exposure data were given to us by GEM – 00:40:54.040 --> 00:40:57.140 the Global Earthquake Model in Pavia. 00:40:57.140 --> 00:41:00.610 Each building type has three vulnerability functions. 00:41:00.610 --> 00:41:03.590 One for the structural elements – you can see here some examples. 00:41:03.590 --> 00:41:08.280 One for the non-structural elements, and one for the contents. 00:41:08.280 --> 00:41:11.090 These functions provide the loss percentage for a given level 00:41:11.090 --> 00:41:15.670 of spectral acceleration. For example, if any wood family 00:41:15.670 --> 00:41:19.920 home in Oklahoma experiences SA of 0.3 seconds equal to 1 g, 00:41:19.920 --> 00:41:25.680 it will lose 10% of its value in the structural elements. 00:41:25.690 --> 00:41:28.930 Unreinforced masonry are the most valuable in terms of structural elements. 00:41:28.930 --> 00:41:32.200 And mobile homes in terms of non-structural elements. 00:41:32.200 --> 00:41:36.030 These [inaudible] building codes were given to us, again, by GEM. 00:41:36.030 --> 00:41:38.980 So let’s look at some loss exceedance curves now. 00:41:38.980 --> 00:41:44.050 The ratio of the 2015 loss exceedance curve over those coming from the 00:41:44.050 --> 00:41:48.050 next year to reach values larger than one order of magnitude – you can see 00:41:48.050 --> 00:41:52.510 here the difference, about 3 times the difference observed in the 00:41:52.510 --> 00:41:55.770 hazard values. So in the hazard, we’d see 3 times difference. 00:41:55.770 --> 00:41:59.740 Here we see 10 times difference. These demonstrate the importance 00:41:59.740 --> 00:42:03.001 of seismic risk modeling. As a comparison, we show the model 00:42:03.001 --> 00:42:07.480 by Gupta and Baker for the year 2015, which leads to much higher losses, 00:42:07.480 --> 00:42:10.170 which even they consider an over-estimation. 00:42:10.170 --> 00:42:13.990 They used the PGA – PGA as their intensity measure and different 00:42:13.990 --> 00:42:16.580 [inaudible] building codes. We think that our model is able to 00:42:16.580 --> 00:42:21.330 yield more realistic results, as you will see in the following slides. 00:42:21.330 --> 00:42:25.470 So, [inaudible] here, 66% of the losses come from non-structural elements, 00:42:25.470 --> 00:42:30.610 20% come from contents, and only 14% come from structural elements. 00:42:30.610 --> 00:42:34.160 This is the evolution of the average annual losses with time. 00:42:34.160 --> 00:42:39.440 And it follows roughly the time history of the seismicity rates. 00:42:39.440 --> 00:42:41.460 So it’s pretty close to what you would expect 00:42:41.460 --> 00:42:44.280 just by looking at the seismicity rates. 00:42:48.920 --> 00:42:52.180 Now, in this map, we show the spatial distribution of the 00:42:52.180 --> 00:42:58.600 forecasted average annual losses for the year 2023. 00:42:58.600 --> 00:43:01.940 The largest risk is identified in Oklahoma City, where most of the 00:43:01.940 --> 00:43:05.550 exposed assets are co-located with significant levels of 00:43:05.550 --> 00:43:08.810 forecasted ground shaking. you see the bull’s-eye here. 00:43:08.810 --> 00:43:12.910 However, the elevated loss estimates are not limited within the area of 00:43:12.910 --> 00:43:16.820 interest for wastewater disposal. The second-largest city, Tulsa, 00:43:16.820 --> 00:43:21.960 is affected by elevated ground motions 50 kilometers away. 00:43:21.960 --> 00:43:24.860 This indicates the interplay between hazard and exposure 00:43:24.860 --> 00:43:29.040 but ultimately controls the seismic risk in Oklahoma. 00:43:29.040 --> 00:43:35.530 We believe that the risk-based metrics, such as this map, could replace 00:43:35.530 --> 00:43:38.690 over-simplified decision variables such as the maximum magnitude, 00:43:38.690 --> 00:43:44.720 or a certain IM level threshold, that are commonly used in traffic light systems. 00:43:44.720 --> 00:43:48.540 Some sensitivity results for the year 2015. 00:43:48.540 --> 00:43:51.140 So the M-max is not very important, as expected. 00:43:51.140 --> 00:43:53.980 The M-min is significant for low return periods. 00:43:53.980 --> 00:43:57.870 But the ground motion model is what really affects our results. 00:43:57.870 --> 00:44:01.740 You can see here, we used Zalachoris and Rathje here. 00:44:01.740 --> 00:44:06.700 But you can use other models that are, in theory, applicable for eastern North 00:44:06.700 --> 00:44:10.960 America and induced seismicity, and you will get very different answers. 00:44:12.500 --> 00:44:16.880 And this is because the actual spectral accelerations from these 00:44:16.890 --> 00:44:20.160 models are very different. This is 5.8 magnitude for 10 kilometers. 00:44:20.160 --> 00:44:23.310 You can see very big differences. And this is kind of surprising because, 00:44:23.310 --> 00:44:28.580 for example, Novakovic and Zalachoris have very similar data sets. 00:44:28.580 --> 00:44:33.820 So – and this model is for eastern North America tectonic earthquakes. 00:44:33.820 --> 00:44:36.820 So this is a big difference that we should highlight. 00:44:36.820 --> 00:44:39.900 And for the b value, it’s not very important unless – 00:44:39.900 --> 00:44:43.680 if you do seismic moment balance. So if you don’t do seismic moment 00:44:43.680 --> 00:44:47.820 balance when you change the b value, then you get very different answers. 00:44:47.820 --> 00:44:50.590 So now let’s compare our modeled losses with observed ones that we 00:44:50.590 --> 00:44:54.920 infer from the insurance data for the period between 2010 and 2016. 00:44:54.920 --> 00:44:58.140 Looking at the insurance data, we estimate the total observed losses 00:44:58.140 --> 00:45:02.720 at 85 million, while the risk model for M-max equal to the observed one of 00:45:02.720 --> 00:45:08.760 5.8 yields 100 million, which is a very good agreement given the uncertainties. 00:45:08.770 --> 00:45:14.100 So we are happy with this result. It validates, in a way, our risk model. 00:45:14.100 --> 00:45:18.600 To put things into perspective, these observed direct losses are only 3% 00:45:18.600 --> 00:45:22.720 of the oil and gas production tax revenue over that time period. 00:45:22.720 --> 00:45:27.140 This is why Oklahoma tolerates these high seismicity rates 00:45:27.140 --> 00:45:30.020 for a whole decade, perhaps. 00:45:30.020 --> 00:45:32.240 So, in conclusion, we developed a semi-empirical 00:45:32.240 --> 00:45:35.320 model linking wastewater injection to seismicity rates. 00:45:35.320 --> 00:45:38.660 The model can be calibrated during the early rise in seismicity 00:45:38.660 --> 00:45:42.080 and forecast earthquake rates given an injection scenario. 00:45:42.080 --> 00:45:44.500 The results are in good agreement with the observed seismicity 00:45:44.500 --> 00:45:48.200 at different scales, outperforming some models in the literature. 00:45:48.200 --> 00:45:51.840 Scientific consensus over fundamental modeling decisions is, 00:45:51.840 --> 00:45:53.560 unfortunately, lacking. 00:45:53.560 --> 00:45:57.690 84% of the felt seismicity in Oklahoma, including the four largest earthquakes, 00:45:57.690 --> 00:46:02.810 can be associated with wastewater disposal at a 95% confidence level. 00:46:02.810 --> 00:46:05.320 While pockets of seismicity in central Oklahoma cannot be 00:46:05.320 --> 00:46:08.890 associated with either wastewater disposal or hydraulic fracturing. 00:46:08.890 --> 00:46:11.500 And we should investigate that. 00:46:11.500 --> 00:46:16.200 In terms of hazard and risk, in 2015, the seismic hazard in several places 00:46:16.210 --> 00:46:20.100 in Oklahoma was higher than the long-term average across 00:46:20.100 --> 00:46:23.130 the San Andreas Fault. Again, striking. 00:46:23.130 --> 00:46:25.640 And the seismic hazard maps illustrate the incompatibility 00:46:25.640 --> 00:46:28.480 of the regional seismic provisions with the current seismicity. 00:46:28.480 --> 00:46:33.020 Now, we are not saying that we should update the building code every year 00:46:33.030 --> 00:46:37.490 because the seismicity changes, but we just highlight this discrepancy. 00:46:37.490 --> 00:46:40.030 Our lost estimates are in good agreement with the paid insurance 00:46:40.030 --> 00:46:44.460 claims. Only 14% of the modeled direct economic losses come from 00:46:44.460 --> 00:46:47.960 structural elements. And our proposed methodology is generic enough 00:46:47.960 --> 00:46:51.020 to be applicable to other regions and can also be used as 00:46:51.020 --> 00:46:54.240 a monitoring tool as I just showed with the last one. 00:46:54.240 --> 00:46:59.280 Now, there is no – there is currently no scientific consensus regarding 00:46:59.280 --> 00:47:02.060 the source of the apparent time lag between the start of injection and 00:47:02.060 --> 00:47:04.950 the onset of seismicity, and that requires future research. 00:47:04.950 --> 00:47:08.900 There is a clear need for a new declustering method. 00:47:08.900 --> 00:47:13.050 More research is needed to explain the variations in the observed b values. 00:47:13.050 --> 00:47:17.500 A hydrogeologic model for the pore pressure would lead to better results, 00:47:17.500 --> 00:47:21.140 perhaps, for our model. And our seismic hazard and risk 00:47:21.140 --> 00:47:24.450 models do not account for epistemic uncertainty. 00:47:24.450 --> 00:47:27.480 And now we could also apply the model to hydraulic fracturing and 00:47:27.480 --> 00:47:29.810 see how it behaves. We can actually – we have 00:47:29.810 --> 00:47:33.340 actually done that, and we have some preliminary results. 00:47:33.340 --> 00:47:38.570 We have modeled all the seismicity in Oklahoma after 2017 that is outside 00:47:38.570 --> 00:47:42.800 the area of interest for wastewater disposal, and we see a very good fit. 00:47:42.800 --> 00:47:46.220 We used a variable time lag – all the things that I showed, 00:47:46.220 --> 00:47:49.030 and we get a good fit. So we can expand to calculating 00:47:49.030 --> 00:47:53.220 p values by hydraulic fracturing. We can do that in Texas, in Oklahoma, 00:47:53.220 --> 00:47:57.840 in Canada, maybe. There are a lot of potential pathways that we can follow. 00:47:57.840 --> 00:47:59.240 These are the main references. 00:47:59.240 --> 00:48:03.200 These three papers contain what I showed you in this presentation. 00:48:03.200 --> 00:48:08.190 And you can read them if you want to see more information about that. 00:48:08.190 --> 00:48:13.060 And these are my acknowledgements. Funding from IUSS, TexNet CISR 00:48:13.060 --> 00:48:17.090 project, the University of Texas at Austin, and also I used resources 00:48:17.090 --> 00:48:22.040 from the DesignSafe, the TACC, and GEM gave me a lot of data. 00:48:22.040 --> 00:48:26.420 And this is my email if you have any questions. Thank you. 00:48:29.860 --> 00:48:34.260 - Okay. Thank you so much for an interesting talk. 00:48:34.270 --> 00:48:38.560 If anyone has questions, feel free to unmute yourself 00:48:38.560 --> 00:48:43.680 and ask a question, or ask the question in the chat. 00:48:46.840 --> 00:48:49.780 - Hey, Jason. Really nice talk. I really enjoyed it. 00:48:49.780 --> 00:48:52.640 You’ve done a ton of work here. 00:48:52.640 --> 00:48:54.980 I had a question. I wanted to go back to the slide 00:48:54.980 --> 00:48:59.050 where you were looking at the p values. I was curious – I don’t know exactly 00:48:59.050 --> 00:49:02.500 where that is. What are – what are the local towns around there? 00:49:02.500 --> 00:49:06.500 I’m just wondering what the – what earthquakes were there. 00:49:08.100 --> 00:49:11.760 - About the center of the state. - Right. Yeah. What are – 00:49:11.770 --> 00:49:13.800 what are the – what are the cities around there? 00:49:13.800 --> 00:49:16.680 And were there, you know, decent-size earthquakes? 00:49:16.680 --> 00:49:23.090 - So this is Oklahoma City. I think this is a magnitude 5, 00:49:23.090 --> 00:49:27.080 and this is the Prague one. So it’s between Prague, Cushing, 00:49:27.080 --> 00:49:31.980 and Oklahoma City, roughly. - Okay, yeah. Yeah. So it’s close to – 00:49:31.980 --> 00:49:33.930 it’s sort of close to Cushing in the northeast. 00:49:33.930 --> 00:49:37.310 I’m just wondering, was there – was there even a lot of seismicity 00:49:37.310 --> 00:49:41.480 in the area? So would it – would you expect it to be well-constrained? 00:49:42.520 --> 00:49:47.740 - Well, I think there is enough to – for the model to work. 00:49:47.740 --> 00:49:52.600 It worked in, for example, here, where you have similar levels of seismicity. 00:49:53.460 --> 00:49:56.060 But, yeah, I think we can look into that. 00:49:56.070 --> 00:50:00.589 I looked at some InSAR data, and there is uplift in that area. 00:50:00.589 --> 00:50:04.800 So this kind of is in favor of the extraction. 00:50:06.080 --> 00:50:09.480 Or, either it’s downlift, sorry. 00:50:09.490 --> 00:50:13.700 But I don’t know. Of course, the depths are different for the extraction 00:50:13.700 --> 00:50:16.710 and where the earthquakes are. So I don’t know what’s the answer, 00:50:16.710 --> 00:50:24.480 but I think – I’m curious why the model works everywhere else and not there. 00:50:25.420 --> 00:50:29.020 I don’t know if you have any insights, but … 00:50:29.020 --> 00:50:31.960 - No, no. I’m just – I’m just curious. Maybe you and I can talk offline, 00:50:31.960 --> 00:50:36.200 but thanks. - Yeah. Yeah. Thanks a lot. 00:50:37.600 --> 00:50:39.780 Thank you. 00:50:39.780 --> 00:50:43.630 - Okay. And we have a question from Andy Barbour, who asks, 00:50:43.630 --> 00:50:46.920 I’m wondering if you can comment on any differences between your 00:50:46.920 --> 00:50:51.540 method for assessing loss and what the PAGER uses. 00:50:52.440 --> 00:50:54.900 - Hm, that’s a good question. 00:50:56.260 --> 00:51:03.700 So, first of all, the input data for the hazard – well, the PAGER, I guess, 00:51:03.700 --> 00:51:06.980 works for every new earthquake that you can run the simulation. 00:51:06.980 --> 00:51:13.260 I think the vulnerability curves are different – I think PAGER 00:51:13.260 --> 00:51:18.260 has all different vulnerability curves. The one we got from GEM are the 00:51:18.260 --> 00:51:25.500 latest one that they used for the Global Risk Map. 00:51:25.500 --> 00:51:29.400 And the ground motion models are different. 00:51:29.410 --> 00:51:31.980 We used the Zalachoris and Rathje one. 00:51:31.980 --> 00:51:41.130 And the exposure – I also think PAGER has older replacement costs. 00:51:41.130 --> 00:51:49.000 We have 2018 replacement costs, so these affect the exposure data as well. 00:51:49.000 --> 00:51:52.550 So I think there are differences. I’m not sure – maybe PAGER 00:51:52.550 --> 00:51:56.310 is updated to newer models that I’m not aware of. 00:51:56.310 --> 00:52:03.140 Maybe someone knows more, but I think key models are different. 00:52:06.700 --> 00:52:10.480 - Thanks. Elizabeth, do you want to ask a question? 00:52:10.480 --> 00:52:18.180 - Yeah. So I was wondering, for the – when you showed the full state 00:52:18.180 --> 00:52:26.600 forecast, prior to the peak in 2014-2015, the number of – 00:52:26.600 --> 00:52:30.080 or, the forecast number of events exceeds the number observed. 00:52:30.090 --> 00:52:33.360 And that’s true for Norbeck and Rubinstein, and it’s true 00:52:33.360 --> 00:52:38.730 for your forecast as well. Do you have any sense of 00:52:38.730 --> 00:52:45.660 why that seismicity rate is overestimated early on? 00:52:46.540 --> 00:52:49.280 Yeah. There. - Yes. Sorry, I was looking 00:52:49.280 --> 00:52:51.680 through the slides. - Yeah. 00:52:52.360 --> 00:52:55.940 - So we discussed that a bit in the paper. 00:52:58.680 --> 00:53:03.740 So one assumption could be that the magnitude of completeness 00:53:03.740 --> 00:53:10.460 is optimistic. So maybe we were not capturing all the magnitude 3s before 00:53:10.460 --> 00:53:16.640 2010 because the network of Oklahoma really got updated in 2010. 00:53:17.460 --> 00:53:23.440 So one assumption is that the model – the red line is correct, and the 00:53:23.440 --> 00:53:26.490 black line should have more earthquakes early on because 00:53:26.490 --> 00:53:29.340 we have missed earthquakes due to the [inaudible]. 00:53:29.340 --> 00:53:33.680 So this is the optimistic answer. 00:53:34.780 --> 00:53:42.120 Now, the other potential reason is that the model really tries to fit the 00:53:42.120 --> 00:53:47.740 high seismicity rates after 2014. And, because it doesn’t have poroelastic 00:53:47.740 --> 00:53:52.890 effects, that might explain rapid changes in seismicity for gradual changes in the 00:53:52.890 --> 00:53:59.070 injection rates, it kind of forces itself to overestimate earlier seismicity. 00:53:59.070 --> 00:54:03.370 So there are some people that talk about poroelastic effects, 00:54:03.370 --> 00:54:08.200 and perhaps they could help. But I think they need so many 00:54:08.200 --> 00:54:13.830 parameters to model them, that we prefer to stay with this 00:54:13.830 --> 00:54:16.400 first order pore pressure effects. 00:54:18.020 --> 00:54:19.520 - Thanks. 00:54:21.360 --> 00:54:23.120 - Annemarie? 00:54:23.780 --> 00:54:28.960 - Hi. Thanks for a great talk. This is maybe outside of, like, 00:54:28.970 --> 00:54:32.020 what you focused on, but I’m wondering if you have an opinion on, 00:54:32.020 --> 00:54:35.820 like, the state of the ground motion models for this region. 00:54:35.820 --> 00:54:40.620 Like, at least what you showed, there’s a lot of variability between the models. 00:54:40.620 --> 00:54:44.670 Do you have an opinion – not necessarily if one of those is better 00:54:44.670 --> 00:54:48.130 or not, but do we need – do we still need more work on these models? 00:54:48.130 --> 00:54:51.840 Or is this just capturing the range, sort of, of uncertainty in the state 00:54:51.840 --> 00:54:55.630 of knowledge for the ground motion modeling? 00:54:56.160 --> 00:54:59.800 - Yeah. That’s a good question. 00:55:01.760 --> 00:55:06.220 I think – I think we should focus more on this issue 00:55:06.220 --> 00:55:09.060 because it really affects all the results. 00:55:09.060 --> 00:55:15.310 And I would like to point out, as well, that, if you really want to do what I did 00:55:15.310 --> 00:55:21.030 using the Vs30 model for the eastern part of the United States, the only 00:55:21.030 --> 00:55:26.360 ground motion model right now that has proper Vs30 scaling for eastern 00:55:26.360 --> 00:55:31.050 North America that it’s it not derived from west North America data 00:55:31.050 --> 00:55:34.700 is the one in Zalachoris and Rathje, and that’s why we used it. 00:55:34.700 --> 00:55:40.290 So it’s important, first of all, to get more Vs30 scaling relations 00:55:40.290 --> 00:55:42.760 for all the other ground motion models. 00:55:42.760 --> 00:55:47.670 Apart from that, to me, it was very striking that Novakovic and Zalachoris 00:55:47.670 --> 00:55:52.730 have so different spectral acceleration values predicted. 00:55:52.730 --> 00:55:56.110 Because their data sets are induced earthquakes in Oklahoma. 00:55:56.110 --> 00:55:58.380 So they are very similar. 00:55:59.320 --> 00:56:05.160 And also the tectonic one for eastern North America is even higher. 00:56:05.160 --> 00:56:12.300 So I don’t think the stress drops of induced seismicity are that different 00:56:12.300 --> 00:56:16.900 from tectonic ones to really change our approach to the ground motion models. 00:56:16.900 --> 00:56:21.800 But I think we should investigate why they’re so different. 00:56:21.800 --> 00:56:25.750 And also, the shape is different as well. - Mm-hmm. 00:56:25.750 --> 00:56:32.720 - And so, yeah, I think – I think there is more work to be done here. 00:56:32.720 --> 00:56:36.380 I don’t have the perfect answer, though, on what exactly we should do. 00:56:36.380 --> 00:56:42.670 - No, just curious if you think there’s still room for things to investigate or – 00:56:42.670 --> 00:56:49.200 yeah. Great. - [inaudible] scaling. Yeah. 00:56:49.200 --> 00:56:51.820 - Thanks. - Of course, all these have limited 00:56:51.820 --> 00:56:56.440 data sets because we don’t have enough earthquakes between 4 and 6. 00:56:56.440 --> 00:57:03.400 - Okay. Are there – are there more questions? I have one, which is, 00:57:03.410 --> 00:57:11.220 do you think there is – there is evidence that injection rates may play a role 00:57:11.220 --> 00:57:16.200 in the induced seismicity on shorter time scales too? 00:57:18.100 --> 00:57:23.440 - What kind of time scales? - Well, what’s the – what are the 00:57:23.440 --> 00:57:27.880 time scales of the – I forget the time scale of the data that you were using … 00:57:27.880 --> 00:57:30.860 - [inaudible] - Monthly? 00:57:30.860 --> 00:57:35.540 - I was using monthly aggregated data. 00:57:35.540 --> 00:57:43.340 I would assume that, yeah, there is a – there is an influence in shorter time 00:57:43.340 --> 00:57:52.560 scales. But we don’t have the data daily. We have after 2016, for the area of 00:57:52.560 --> 00:57:59.700 interest, for the Arbuckle-only wells. So you could do daily data 00:57:59.700 --> 00:58:02.660 aggregation – with the daily data aggregation for the hydraulic 00:58:02.660 --> 00:58:10.270 fracturing, for example, so the model works, even for daily inputs. 00:58:11.320 --> 00:58:15.360 But, if we want to capture this whole time period for the whole Oklahoma 00:58:15.360 --> 00:58:18.610 for all the target formations, we had to use monthly values. 00:58:18.610 --> 00:58:23.200 But I agree. We do under-sample, in a way, the data. 00:58:23.200 --> 00:58:27.040 Because we do the monthly aggregation. 00:58:27.040 --> 00:58:29.420 It is very important for the hydraulic fracturing, this is – 00:58:29.420 --> 00:58:32.600 you cannot do monthly for hydraulic fracturing. 00:58:33.000 --> 00:58:34.400 - Okay. 00:58:36.260 --> 00:58:39.920 Are there any other questions? Okay, well, if not, we know – 00:58:39.920 --> 00:58:46.050 we have Iason’s email address. We can always reach out to him 00:58:46.050 --> 00:58:50.360 if we have other questions too. And, with that, let’s thank our speaker, 00:58:50.360 --> 00:58:53.460 and see you guys next week. - Thank you very much. 00:58:53.460 --> 00:58:54.860 - Thank you. 00:58:54.860 --> 00:58:56.780 - Thanks, Iason. - Great talk. 00:58:56.780 --> 00:58:58.640 It’s a great body of work. Thanks. 00:58:58.640 --> 00:59:02.080 - That was great. Thank you. - Good job, Iason. 00:59:02.080 --> 00:59:04.460 - Thank you. Thank you, Paul. 00:59:04.460 --> 00:59:08.960 - Good job, Iason. - Thanks. Thanks, Alexander.