Conditional simulations of max-stable processes C. Dombry , . - - PowerPoint PPT Presentation

Conditional simulations of max-stable processes Mathieu Ribatet – 1 / 24

Conditional simulations of max-stable processes

• C. Dombry†,

F . Éyi-Minko†,

• M. Ribatet‡

† Laboratoire de Mathématiques et Application, Université de Poitiers ‡ Institut de Mathématiques et de Modélisation, Université Montpellier 2

Goal

Conditional simulations of max-stable processes Mathieu Ribatet – 2 / 24

• Max-stable processes are widely used for modelling spatial extremes since

they arise as the only possible (non degenerate) limit of pointwise maxima

• ver independent replicates, i.e.,

max

i=1,...,n

Xi(x)−bn(x) an(x) − → Z(x), n → ∞, x ∈ X ⊂ Rd, for some normalizing functions an > 0 and bn and where Xi are independent copies of a stochastic process X .

Goal

Conditional simulations of max-stable processes Mathieu Ribatet – 2 / 24

• Max-stable processes are widely used for modelling spatial extremes since

they arise as the only possible (non degenerate) limit of pointwise maxima

• ver independent replicates, i.e.,

max

i=1,...,n

Xi(x)−bn(x) an(x) − → Z(x), n → ∞, x ∈ X ⊂ Rd, for some normalizing functions an > 0 and bn and where Xi are independent copies of a stochastic process X .

processes (with continuous spectral measure)?

Setup

Conditional simulations of max-stable processes Mathieu Ribatet – 3 / 24

• Recall that any max-stable process has the spectral characterization

Z(·) = max

i≥1 ζiYi(·),

where

–

Yi(·) are independent copies of a non negative stochastic process such that E[Y (x)] = 1 for all x ∈ X ;

–

{ζi}i≥1 are the points of a Poisson process on (0,∞) with intensity dΛ(ζ) = ζ−2dζ.

Z(·) | {Z(x1) = z1,...,Z(xk) = zk}, for some z1,...,zk > 0 and k conditioning locations x1,...,xk ∈ X .

• 1. Conditional distributions of max-stable

processes

⊲

• 1. Conditional

distributions Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 4 / 24

Decomposition of Φ

• 1. Conditional

distributions Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 5 / 24

• Let Φ a point process on (0,∞)k whose atoms are

ϕi(x) = ζiYi(x), x = (x1,...,xk).

• Consider the two following point processes

Φ− =

• ϕ ∈ Φ: ϕ(xi) < zi, for all i ∈ {1,...,k}
• ,(sub-extremal functions)

Φ+ =

• ϕ ∈ Φ: ϕ(xi) = zi, for some i ∈ {1,...,k}
• .(extremal functions)
• Clearly Φ = Φ− ∪Φ+ and

Φ+ = {ϕ+

1 ,...,ϕ+ k } = {ϕ+ 1 ,...,ϕ+ ℓ },

a.s. (1 ≤ ℓ ≤ k).

Decomposition of Φ

• 1. Conditional

distributions Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 5 / 24

• Let Φ a point process on (0,∞)k whose atoms are

ϕi(x) = ζiYi(x), x = (x1,...,xk).

• Consider the two following point processes

Φ− =

• ϕ ∈ Φ: ϕ(xi) < zi, for all i ∈ {1,...,k}
• ,(sub-extremal functions)

Φ+ =

• ϕ ∈ Φ: ϕ(xi) = zi, for some i ∈ {1,...,k}
• .(extremal functions)
• Clearly Φ = Φ− ∪Φ+ and

Φ+ = {ϕ+

1 ,...,ϕ+ k } = {ϕ+ 1 ,...,ϕ+ ℓ },

a.s. (1 ≤ ℓ ≤ k).

pendent.

Conditional intensity function

• 1. Conditional

distributions Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 6 / 24

Z(x) = max

i≥1 ζiYi(x) = max i≥1 ϕi(x)

• The Poisson point process {ϕi(x)}i≥1 has intensity measure

Λx(A) = ∞ Pr{ζY (x) ∈ A}ζ−2dζ, Borel set A ⊂ Rk.

• We assume that Φ is regular, i.e., Λx(dz) = λx(z)dz, for all

x ∈ X k.

Conditional intensity function

• 1. Conditional

distributions Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 6 / 24

Z(x) = max

i≥1 ζiYi(x) = max i≥1 ϕi(x)

• The Poisson point process {ϕi(x)}i≥1 has intensity measure

Λx(A) = ∞ Pr{ζY (x) ∈ A}ζ−2dζ, Borel set A ⊂ Rk.

• We assume that Φ is regular, i.e., Λx(dz) = λx(z)dz, for all

x ∈ X k.

λs|x,z(u) = λ(s,x)(u,z) λx(z) is the (regular) conditional distribution of Z(x)—if we integrate w.r.t. all possible partitions of x. But not that of Z(·)!!!

Random partitions?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 7 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4

Random partitions?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 7 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4 ϕ1

+

Random partitions?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 7 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4 ϕ1

+

ϕ2

+

Random partitions?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 7 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4 ϕ1

+

ϕ2

+

ϕ3

+

Random partitions?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 7 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4 ϕ1

+

ϕ2

+

ϕ3

+

ϕ4

+

Random partitions?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 7 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4 ϕ1

+

ϕ2

+

ϕ3

+

ϕ4

+

Here the set {x1,...,x5} is partitioned into ({x1,x3},{x2},{x4},{x5})

Random partitions?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 7 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4 ϕ1

+

ϕ2

+

ϕ3

+

ϕ4

+

Here the set {x1,...,x5} is partitioned into ({x1,x3},{x2},{x4},{x5})

• The hitting bounds {zi}i=1,...,k might be reached by several

extremal functions, i.e., Φ+ = {ϕ+

1 ,...,ϕ+ k } = {ϕ+ 1 ,...,ϕ+ ℓ } a.s.,

1 ≤ ℓ ≤ k.

• So we need to take into account all possible ways these hitting

bounds are reached: the hitting scenarios

Why should we bother about Φ−?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 8 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4

Why should we bother about Φ−?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 8 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4 maxΦ+

Why should we bother about Φ−?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 8 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4 maxΦ+ maxΦ−

Why should we bother about Φ−?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 8 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4 maxΦ+ maxΦ−

Why should we bother about Φ−?

• 1. Conditional

distributions

⊲

Random partitions Sampling scheme Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 8 / 24

0.0 0.5 1.0 1.5 2.0 2.5 x Z(x) x2 x5 x1 x3 x4 maxΦ+ maxΦ−

• The atoms of Φ+ are only of interest if we restrict our attention

to the conditioning points x;

• But most often one would like to get realizations at s = x.

maxΦ−(s) > maxΦ+(s)!

A three step procedure

• 1. Conditional

distributions Random partitions

⊲ Sampling scheme

Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 9 / 24

The above key points suggest a three step sampling scheme: Step 1 Draw a random partition τ, i.e., a hitting scenario; Step 2 Given τ of size ℓ, draw the extremal functions ϕ+

1 ,...,ϕ+ ℓ

independently; Step 3 Independently from Steps 1 & 2, draw the sub-extremal functions ϕ−

i , i ≥ 1.

Step 1: The random partitions

• 1. Conditional

distributions Random partitions

⊲ Sampling scheme

Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 10 / 24

• Let Pk the set of all possible partitions of the set {x1,...,xk}.
• Draw a random partition τ ∈ Pk with distribution

πx(z,τ) = 1 C(x,z)

|τ|

• j=1

λxτj (zτj )

• density that some

bounds are reached, i.e., the zτj

• {u<zτc

j }

λxτc

j |xτj ,zτj (u)du

• probability to lie below

the remaining bounds, i.e., below the zτc j

, where the normalization constant C(x,z) is given by C(x,z) =

• θ∈Pk

|θ|

• j=1

λxθj (zθj )

• {u<zθc

j }

λxθc

j |xθj ,zθj (u)du,

and |τ| is the “size” of the partition τ.

Step 2: The extremal functions

• 1. Conditional

distributions Random partitions

⊲ Sampling scheme

Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 11 / 24

• Given τ = (τ1,...,τℓ), draw ℓ independent random vectors

ϕ+

1 (s),...,ϕ+ ℓ (s) from the distribution

Pr

• ϕ+

j (s) ∈ dvj

• = 1

C j

• 1{u<zτc

j }λ(s,xτc j )|xτj ,zτj (vj,u)

• density of ϕ ∈ Φ

such that ϕ(xτj ) = zτj

du

• dvj,

where 1{·} is the indicator function and C j =

• 1{u<zτc

j }λ(s,xτc j )|xτj ,zτj (vj,u)dudvj.

• Define the random vector

Z +(s) = max

j=1,...,ℓϕ+ j (s),

s ∈ X m.

Step 3: The sub-extremal functions

• 1. Conditional

distributions Random partitions

⊲ Sampling scheme

Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 12 / 24

• Independently draw {ζi}i≥1 a Poisson point process on (0,∞)

with intensity ζ−2dζ and {Yi(·)}i≥1 independent copies of Y (·)

• Define the random vector

Z −(s) = max

i≥1 ζiYi(s)1{ζiYi(x)<z},

s ∈ X m.

Step 3: The sub-extremal functions

• 1. Conditional

distributions Random partitions

⊲ Sampling scheme

Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 12 / 24

• Independently draw {ζi}i≥1 a Poisson point process on (0,∞)

with intensity ζ−2dζ and {Yi(·)}i≥1 independent copies of Y (·)

• Define the random vector

Z −(s) = max

i≥1 ζiYi(s)1{ζiYi(x)<z},

s ∈ X m.

˜ Z(s) = max

• Z +(s),Z −(s)
• follows the conditional distribution of Z(s) given Z(x) = z.

As an aside

• 1. Conditional

distributions Random partitions

⊲ Sampling scheme

Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 13 / 24

• The conditional cumulative distribution function is

Pr{Z(s) ≤ a | Z(x) = z} =

• τ∈Pk

πx(z,τ)

|τ|

• j=1

Fτ,j (a)

• Steps 1 & 2

Pr[Z(s) ≤ a,Z(x) ≤ z] Pr[Z(x) ≤ z]

• Step 3

,

where Fτ,j(a) =

• {y<zτc

j ,u<a} λ(s,xτc j )|xτj ,zτj (u,y)dydu

• {y<zτc

j } λtτc j |xτj ,zτj (y)dy

.

As an aside

• 1. Conditional

distributions Random partitions

⊲ Sampling scheme

Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 13 / 24

• The conditional cumulative distribution function is

Pr{Z(s) ≤ a | Z(x) = z} =

• τ∈Pk

πx(z,τ)

|τ|

• j=1

Fτ,j (a)

• Steps 1 & 2

Pr[Z(s) ≤ a,Z(x) ≤ z] Pr[Z(x) ≤ z]

• Step 3

,

where Fτ,j(a) =

• {y<zτc

j ,u<a} λ(s,xτc j )|xτj ,zτj (u,y)dydu

• {y<zτc

j } λtτc j |xτj ,zτj (y)dy

.

• Remark. It is “clear” that Z(·) | {Z(x) = z} is not max-stable.

Brown–Resnick processes

• 1. Conditional

distributions Random partitions Sampling scheme

⊲ Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 14 / 24

Example 1. If Z is a Brown–Resnick process, i.e., Z(x) = max

i≥1 ζi exp{εi(x)−γ(x)},

x ∈ X , then the intensity function is λx(z) = Cx exp

• −1

2 logzTQx logz+Lx logz k

• i=1

z−1

i ,

z ∈ (0,∞)k, and the conditional intensity function is

λs|x,z(u) = (2π)−m/2|Σs|x|−1/2 exp

• − 1

2 (logu−µs|x,z)T Σ−1

s|x(logu−µs|x,z)

m

• i=1

u−1

i

,

i.e., the extremal functions are log-Normal processes.

Schlather processes

• 1. Conditional

distributions Random partitions Sampling scheme

⊲ Examples

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 15 / 24

Example 2. If Z is a Schlather process, i.e., Z(x) =

• 2πmax

i≥1 ζi max{0,εi(x)},

x ∈ X , then the intensity function is λx(z) = π−(k−1)/2|Σx|−1/2ax(z)−(k+1)/2Γ k +1 2

• ,

z ∈ Rk, where ax(z) = zT Σ−1

x z, and the conditional intensity function is

λs|x,z(u) = π−m/2(k +1)−m/2|˜ Σ|−1/2

• 1+ (u−µ)T ˜

Σ−1(u−µ) k +1 −(m+k+1)/2 Γ

• m+k+1

2

• Γ
• k+1

2

• ,

i.e., the extremal functions are Student processes.

• 2. Markov chain Monte–Carlo sampler

(for Step 1)

• 1. Conditional

distributions

⊲ 2. MCMC sampler

Computational burden Full conditional distributions The full conditional distributions are nice !

• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 16 / 24

Do you recognize these numbers?

• 1. Conditional

distributions

• 2. MCMC sampler

Computational burden Full conditional distributions The full conditional distributions are nice !

• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 17 / 24

1 1 2 5 15 52 203 877 4140 21147 115975 678570 4213597 27644437 190899322 1382958545 10480142147 82864869804 682076806159 5832742205057 ...

Do you recognize these numbers?

• 1. Conditional

distributions

• 2. MCMC sampler

Computational burden Full conditional distributions The full conditional distributions are nice !

• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 17 / 24

1 1 2 5 15 52 203 877 4140 21147 115975 678570 4213597 27644437 190899322 1382958545 10480142147 82864869804 682076806159 5832742205057 ...

• Remark. Recall that Bell(k) is the number of partitions of a set

with k elements. Hence with our terminology we have # hitting scenarios = Card(Pk) = Bell(k).

Computational burden

• 1. Conditional

distributions

• 2. MCMC sampler

⊲

Computational burden Full conditional distributions The full conditional distributions are nice !

• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 18 / 24

• In Step 1, we need to sample from a discrete distribution

whose state space is Pk. Combinatorial explosion Hence we cannot compute C(x,z) in πx(z,τ) = 1 C(x,z)

|τ|

• j=1

λxτj (zτj )

• {u<zτc

j }

λxτc

j |xτj ,zτj (u)du.

Computational burden

• 1. Conditional

distributions

• 2. MCMC sampler

⊲

Computational burden Full conditional distributions The full conditional distributions are nice !

• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 18 / 24

• In Step 1, we need to sample from a discrete distribution

whose state space is Pk. Combinatorial explosion Hence we cannot compute C(x,z) in πx(z,τ) = 1 C(x,z)

|τ|

• j=1

λxτj (zτj )

• {u<zτc

j }

λxτc

j |xτj ,zτj (u)du.

• Remark. We will use a Gibbs sampler since the full conditional

distributions are especially convenient.

Full conditional distributions

• 1. Conditional

distributions

• 2. MCMC sampler

Computational burden

⊲

Full conditional distributions The full conditional distributions are nice !

• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 19 / 24

• Want to sample from Pr[θ ∈ · | θ−j = τ−j], θ ∼ πx(z,·) where τ−j

is the restriction of τ ∈ Pk to the set {x1,...,xk}\{xj}.

• The number of possible states for θ is

b+ =

• ℓ

if {xj} is a partitioning set of τ, ℓ+1

• therwise.

Example 3. For τ = ({x1,x2},{x3}) we have ℓ = 2 and Restriction τ∗

−2 = τ−2

τ∗

−3 = τ−3

Possible states ({x1,x2},{x3}) ({x1,x2},{x3}) ({x1},{x2,x3}) ({x1,x2,x3}) ({x1},{x2},{x3}) ——

The full conditional distributions are nice !

• 1. Conditional

distributions

• 2. MCMC sampler

Computational burden Full conditional distributions

⊲

The full conditional distributions are nice !

• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 20 / 24

• For all τ∗ ∈ Pk such that τ∗

−j = τ−j,

Pr[θ = τ∗ | θ−j = τ−j ] = πx(z,τ∗)

• ˜

τ∈Pk

πx(z, ˜ τ)1{˜

τ−j =τ−j }

∝ |τ∗|

j=1 wτ∗,j

|τ|

j=1 wτ,j

,

where wτ,j = λxτj (zτj )

• {u<zτc

j } λxτc j |xτj ,zτj (u)du.

the Gibbs sampler is especially convenient!

• 3. Application
• 1. Conditional

distributions

• 2. MCMC sampler

⊲ 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 21 / 24

Precipitation around Zurich in 2000

• 1. Conditional

distributions

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 22 / 24

20 40 60 80 30 40 50 60 70 80 90 40 60 80 100 1.00 1.02 1.04 1.06 1.08 1.10

Figure 1: From left to right, maps on a 50 × 50 grid of the pointwise 0·025, 0·5 and 0·975 sample quantiles for rainfall (mm) obtained from 10000 conditional simulations of Brown– Resnick processes having semi variogram γ(h) = (h/38)0·69. The rightmost panel plots the ratio of the width of the pointwise confidence intervals with and without taking estimation uncertainties into account. The squares show the conditional locations.

Temperature anomalies for the 2003 European heatwave

• 1. Conditional

distributions

• 2. MCMC sampler
• 3. Application

Conditional simulations of max-stable processes Mathieu Ribatet – 23 / 24

1000 2000 3000 4000

Arosa (1840) Bad Ragaz (496) Basel (316) Bern (565) Chateau d’Oex (985) Davos (1590) Engelberg (1035) Gd−St−Bernard (2472) Locarno−Monti (366) Lugano (273) Montana (1508) Montreux (405) Neuchatel (485) Oeschberg (483) Santis (2490) Zurich (556)

25 50 75 100

(km) (m)

2.5 3 3.5 4 4.5

(°C)

Figure 2: Left: Topographical map of Switzerland showing the sites and altitudes in metres above sea level of 16 weather stations for which annual maxima temperature data are avail-

• able. Right: Map of temperature anomalies (◦C), i.e., the difference between the pointwise

medians obtained from 10000 conditional simulations and unconditional medians esti- mated from the fitted Schlather process.

• As expected the largest deviations occur in the plateau region
• f Switzerland
• The differences range between 2·5◦C and 4·75◦C

THANK YOU! ANY QUESTIONS?

Dombry, C. Éyi-Minko, F . and Ribatet, M. (2012) Conditional simulation of max-stable processes To appear in Biometrika.