Modeling spatial uncertainty
Modeling Uncertainty in the Earth Sciences Jef Caers Stanford University
A reminder: models of uncertainty
Models of uncertainty P(A) P(A|B b) f X ( X Y y) Samples: x1 , x2 , x3 ,..., xn
A set of samples drawn by Monte Carlo simulation are a valid model of uncertainty
Motivation
uncertain uncertain certain or uncertain uncertain
Spatial Input parameters
Spatial Stochastic model
Physical model
response
Forecast and decision model
uncertain
Datasets
Physical input parameters
uncertain
Raw observations
uncertain/error
Spatial stochastic simulation
Modeling spatial uncertainty
Spatial stochastic simulation
Seed 1 Seed 2 Earth models (input uncertainty)
input
Spatial Stochastic simulator
(spatial uncertainty)
Data
(samples, well, geophysics)
Conditional simulation = constrained to data Unconditional simulation = not constrained to any data
Hard data
Definition: Hard data = direct, exact information at
the scale of modeling
All the rest is soft data
If you dig out a volume of this size at that location and measure the variable of interest than you have hard data
Example
Size of the core = 5 miles Grid cells 6 inch x 6 inch x 12 inch 1 grid cell = 200ft x 200 ft x 3 ft
If we assume the core to be hard data, i.e. assign it to a grid cell then we basically assume that there is NO spatial variation in that grid cell
Example spatial stochastic simulation
Input uncertainty: Range horizontal= [15,25] Range vertical = 5 Isotropy or horizontal anisotropy max range = [20,25] min range = [10,15] Azm=45
Example spatial stochastic simulation
If I take a fixed input Range horizontal= 25, Range vertical = 5 Horizontal Isotropy
Result
horizontal
vertical
The Earth model reflects our knowledge about The sampled values at their location The interpreted spatial continuity model
Spatial uncertainty
and infinitely more if you want
Input and spatial uncertainty
and infinitely more if you want
Object-based stochastic simulation
Modeling spatial uncertainty
Unconditional simulation
Principle mostly used : rejection sampler
Place an object drawn from the pdfs of the various
parameters defining the object Place the object, either randomly or according to trend If object violates interaction or other rules
Then reject Else accept
This will make sure that the 3D Earth model created
reflects exactly the Boolean model specified
Conditional simulation
Use rejection sampler: slow Use Metropolis sampler
Create an initial 3D Earth model: it probably violates data
constraints Propose a perturbation
Move an object Remove an object Add an object
Accept the change with a probability a, where a is dependent
on how much improvement was made
These methods are iterative methods, hence slow and
may not converge
Example: conditional simulation
After 198 iterations Some well data still not honored After 200 iterations
Final result: constrained to all wells
Channel position changes After 202 iterations After 204 iterations
Channel added
Same channel removed
Courtesy: Norwegian Computing Center
Conditional Boolean simulation
Too slow for practical Earth modeling involving
uncertainty
Any conflict between model and data makes it even
slower
Can only be constrained to
Limited amount of well data Not to geophysical data or more complex data
Training-image-based simulation
Modeling spatial uncertainty
Conditional Earth model
Idea
Generate a single unconditional Boolean Earth model
Anchor patterns to data
Actual Data Extract 3D Patterns
Principle of sequential simulation
A training image A reservoir with a 2x2 grid
Another realization
Step 1. Pick a cell
Step 2. Assign probability 50%
Step 3. Assign color
Step 4. Pick a cell
100%
Step 5. Assign color
Step 6. Final result
In general
= data
data event
P(A|B)?
Algorithm
Generic sequential simulation algorithm (1) Assign any hard data to grid cells if required (2) Define a random path (3) Loop over all grid cells (1) Determine P(A|B) B=any data and previously simulated values (2) Draw from P(A|B) a value (3) Add that value to the data set
How to use the training image ?
Simulation grid with some data points
u2 u4
u?
u1
u3
Using a training image
Training image
P(A|B)=1/4 A = blue
Scanning is CPU-demanding
Simulation grid with some data points
u2
u?
u4 u1
u3
Template
Create a data-base
Training image
Search tree
14 11
Construction requires scanning training image one single time
5 3
Minimizes memory demand
3 1 2 5 3 0 1 1
Allows retrieving all training probabilities for the template adopted
0 2
Spatial continuity at large scale
Training image Freeze coarse grid nodes and use them as conditioning data to simulate finer grid nodes Finer simulation grid
Coarse simulation grid
fine template
Coarse template
Example
Facies type Tidal bars Conceptual description Elongated ellipses w/ upper sigmoidal cross-section Sheets (rectangles) Stratigraphy Length (m) 2000 to 4000 Width (m) 500 Thickness (ft) 3 to 7 Anywhere
Tidal sand flats Estuarine sands Transgr. Lags
Anywhere, eroded by sand bars Top of reservoir Top of estuarines
2000
1000
Sheets (rectangles) Sheets (rectangles)
4000 3000
2000 1000
8 4
Training image Background shales Tidal sand flats Transgressive lags Tidal bars Estuarine sands
Example
N
Aerial proportion maps
Estuarine sands
Plan view of stratigraphic grid with location of the 140 wells
Sand bars
Background shale Vertical proportion curves Background shale Estuarine sands 1
Facies model
Sand bars
Variogram-based simulation
Modeling spatial uncertainty
Introduction to spatial estimation
Data point
uj ?
Spatial estimation = What is the best guess for the value at the location where no data was taken ?
There is only one single guess that is the best Depends on what you determine as best
What is best?
Best = as close as possible to the unknown truth Consider a situation where you want to estimate the
total amount of pollution of Pb at a specific location. You have two methodologies do so. Consider that you apply these methodologies to 10 sites
Estimation: principles
site 1 2 3 4 5 6 7 8 9 10 estimation estimation unknown error method 1 (1) 1 m (1) m 2 (1) m 3 (1) m 4 (1) m 5 (1) m 6 (1) m 7 (1) m 8 (1) m 9 (1) 10 m method 2 (2) 1 m (2) m 2 (2) m 3 (2) m 4 (2) m 5 (2) m 6 (2) m 7 (2) m 8 (2) m 9 (2) 10 m real# Pb m1 m2 m3 m4 m5 m6 m7 m8 m9 m10 error method 1 method 2 (1) (2) e1 e1
(1) e2 (1) e3 (1) e4 (1) e5 (1) e6 (1) e7 (1) e8 (1) e9 (1) e10 (2) e2 (2) e3 (2) e4 (2) e5 (2) e6 (2) e7 (2) e8 (2) e9 (2) e10
Estimation: principles
Unbiased: the average error is zero (it is a property
measured over many trials Best ?
Average square error is zero ? Absolute value of error is zero ?
loss
-error
+error
Introduction to Kriging
What is kriging ? Is an estimation method Finds the best (Least Square) linear estimate of the unknown Accounts for the variogram
* Spatial correlation between unknown and data * Redundancy between data But is mostly used in sequential simulation
Linear estimation
Problem Inverse distance solution
1 1 e.g. 3
2 i 2
z1=0.5 d
1
z izi
i1 3
z2=0.9 d
2
d 1
2 i
1d 1d
1 2
2 d3
1d 3
d
3
z3=1.5
h13 h12 1
Layered system
Inverse distance mapping
d i i 1 di 1
1 2 1 2
3 Situation 1 2=1/3 1=1/3 1=1/4
3 Situation 2 2=1/4
3=1/3 Inverse distance
3=1/2
kriging
Linear estimation
To be estimated
Z* SK (uj )
,n
Data Z(ua ), a 1,
uj ?
Direction of major continuity
Z (uj ) a Z(ua ), what is a ?
* SK a1
Principle 1:
a datum close in geological distance to the unknown should get a large weight
Principle 2:
Data close together are redundant and should share their weight
How does kriging do it?
z(u1)=0.5 h12
z(u) izi
i1 3
z(u2)=0.9
h13
h23
z(u3)=1.5
Record all the distance between Data locations Data location vs location of unknown Calculate the covariance function for those distances
How does kriging do it?
Var ( z ) C (h12 ) C (h13 ) 1 C (h 01 ) C (h ) Var ( z ) C (h ) C ( h ) 12 23 02 2 C ( h13 ) C ( h 23 ) Var ( z ) 3 C (h 03 )
Solve this linear system of equations
u j ) a z(ua ) is the kriging estimate z(
a1 n
Note: mean = assumed zero
What is the average error we make?
Var ( z ) C (h12 ) C (h13 ) 1 C (h 01 ) C (h ) Var ( z ) C (h ) C ( h ) 12 23 02 2 C ( h13 ) C ( h 23 ) Var ( z ) 3 C (h 03 )
Solve this linear system of equations
2 (u j ) 2 a C (h0a ) is the kriging variance
a1
Note: mean = assumed zero
Example
e d ? d
Var (z) Var (z) C (2d) 1 C (d) Var (z) Var (z) C (2d) C (d) 2 C (2d) C (2d) Var (z) 3 C (d)
1 1 / 4 1 / 4 2 3 1 / 2
1 1 1 2 (u) var(z) C (d) C (d) C (d) var(z) C (d) 4 4 2
Note: mean = assumed zero
Various flavors of kriging
Simple kriging
You assume there is no trend You assume you know the mean of the variable over the
domain
Ordinary kriging
You dont want to assume anything explicitly
Kriging with locally varying mean
You assume there is a trend and you know that trend exactly
Kriging with trend
You assume there is trend but you only know the type of
trend (dont know it exactly its magnitude)
Example
Back to sequential simulation
(1) Assign any hard data to grid cells if required (2) Define a random path (3) Loop over all grid cells (1) Determine P(A|B) B=any data and previously simulated values (2) Draw from P(A|B) a value (3) Add that value to the data set
Kriging is not simulation !
Kriging: what is the best guess?
uj ?
uj ?
Simulation: what is the uncertainty as expressed through a probability or probability distribution or a set of Possible outcomes
Gaussian simulation
Variance: 2
uj ?
Mean: m
mean: Z (u j ) a Z (ua ) is the simple kriging mean
* SK a1
variance: (u j ) 1 a Cov(u j ua ) is the kriging variance
2 SK a1
Transformation of the data
Gaussian simulation assumes the data is Gaussian
Problem ?
Prior to simulation, perform a normal score transform of the hard data
After sequential visit of all cells, perform a back-transform which is the exact reverse of the normal score transform
Normal score transform
Unit free
Complete SGS algorithm
(1) Transform the data into normal score domain
(2) Assign any hard data to grid cells if required
(3) Define a random path
(4) Loop over all grid cells (1) Determine, using kriging, the distribution P(A|B) B=any data and previously simulated values (2) Draw from P(A|B) a value (3) Add that value to the data set
(5) Back-transform all simulated values (requires extra/interpolation)
Example