Lecture 12 Gaussian Process Models 10/16/2018 1 Multivariate - - PowerPoint PPT Presentation

Lecture 12

Gaussian Process Models

10/16/2018

1

Multivariate Normal

Multivariate Normal Distribution

For an 𝑜-dimension multivate normal distribution with covariance 𝚻 (positive semidefinite) can be written as

𝐙

𝑜×1 ∼ 𝑂( 𝝂 𝑜×1

, 𝚻

𝑜×𝑜) where {𝚻}𝑗𝑘 = 𝜏2 𝑗𝑘 = 𝜍𝑗𝑘 𝜏𝑗 𝜏𝑘

⎛ ⎜ ⎝ 𝑍1 ⋮ 𝑍𝑜 ⎞ ⎟ ⎠ ∼ 𝑂 ⎛ ⎜ ⎝ ⎛ ⎜ ⎝ 𝜈1 ⋮ 𝜈𝑜 ⎞ ⎟ ⎠ , ⎛ ⎜ ⎝ 𝜍11𝜏1𝜏1 ⋯ 𝜍1𝑜𝜏1𝜏𝑜 ⋮ ⋱ ⋮ 𝜍𝑜1𝜏𝑜𝜏1 ⋯ 𝜍𝑜𝑜𝜏𝑜𝜏𝑜 ⎞ ⎟ ⎠ ⎞ ⎟ ⎠

2

Density

For the 𝑜 dimensional multivate normal given on the last slide, its density is given by

(2𝜌)−𝑜/2 det(𝚻)−1/2 exp (−1 2(𝐙 − 𝝂)′

1×𝑜

𝚻−1

𝑜×𝑜(𝐙 − 𝝂) 𝑜×1

)

and its log density is given by

−𝑜 2 log 2𝜌 − 1 2 log det(𝚻) − −1 2(𝐙 − 𝝂)′

1×𝑜

𝚻−1

𝑜×𝑜(𝐙 − 𝝂) 𝑜×1

3

Sampling

To generate draws from an 𝑜-dimensional multivate normal with mean 𝝂 and covariance matrix 𝚻,

• Find a matrix 𝐁 such that 𝚻 = 𝐁 𝐁𝑢, most often we use

𝐁 = Chol(𝚻) where 𝐁 is a lower triangular matrix.

• Draw 𝑜 iid unit normals (𝒪(0, 1)) as 𝐴
• Obtain multivariate normal draws using

𝐙 = 𝝂 + 𝐁 𝐴

4

Sampling

To generate draws from an 𝑜-dimensional multivate normal with mean 𝝂 and covariance matrix 𝚻,

• Find a matrix 𝐁 such that 𝚻 = 𝐁 𝐁𝑢, most often we use

𝐁 = Chol(𝚻) where 𝐁 is a lower triangular matrix.

• Draw 𝑜 iid unit normals (𝒪(0, 1)) as 𝐴
• Obtain multivariate normal draws using

𝐙 = 𝝂 + 𝐁 𝐴

4

Sampling

To generate draws from an 𝑜-dimensional multivate normal with mean 𝝂 and covariance matrix 𝚻,

• Find a matrix 𝐁 such that 𝚻 = 𝐁 𝐁𝑢, most often we use

𝐁 = Chol(𝚻) where 𝐁 is a lower triangular matrix.

• Draw 𝑜 iid unit normals (𝒪(0, 1)) as 𝐴
• Obtain multivariate normal draws using

𝐙 = 𝝂 + 𝐁 𝐴

4

Sampling

To generate draws from an 𝑜-dimensional multivate normal with mean 𝝂 and covariance matrix 𝚻,

• Find a matrix 𝐁 such that 𝚻 = 𝐁 𝐁𝑢, most often we use

𝐁 = Chol(𝚻) where 𝐁 is a lower triangular matrix.

• Draw 𝑜 iid unit normals (𝒪(0, 1)) as 𝐴
• Obtain multivariate normal draws using

𝐙 = 𝝂 + 𝐁 𝐴

4

Bivariate Example

𝝂 = (0 0) 𝚻 = (1 𝜍 𝜍 1)

rho=−0.9 rho=−0.7 rho=−0.5 rho=−0.1 rho=0.9 rho=0.7 rho=0.5 rho=0.1 −2.5 0.0 2.5 −2.5 0.0 2.5 −2.5 0.0 2.5 −2.5 0.0 2.5 −4 −2 2 −4 −2 2

x y 5

Marginal distributions

Proposition - For an 𝑜-dimensional multivate normal with mean 𝝂 and covariance matrix 𝚻, any marginal or conditional distribution of the 𝑧’s will also be (multivariate) normal. For a univariate marginal distribution,

𝑧𝑗 = 𝒪(𝝂𝑗, 𝚻𝑗𝑗)

For a bivariate marginal distribution,

𝐳𝑗𝑘 = 𝒪 ((𝝂𝑗 𝝂𝑘 ) , (𝚻𝑗𝑗 𝚻𝑗𝑘 𝚻𝑘𝑗 𝚻𝑘𝑘 ))

For a 𝑙-dimensional marginal distribution,

𝐳𝑗,⋯,𝑙 = 𝒪 ⎛ ⎜ ⎝ ⎛ ⎜ ⎝ 𝝂𝑗 ⋮ 𝝂𝑙 ⎞ ⎟ ⎠ , ⎛ ⎜ ⎝ 𝚻𝑗𝑗 ⋯ 𝚻𝑗𝑙 ⋮ ⋱ ⋮ 𝚻𝑙𝑗 ⋯ 𝚻𝑙𝑙 ⎞ ⎟ ⎠ ⎞ ⎟ ⎠

6

Marginal distributions

Proposition - For an 𝑜-dimensional multivate normal with mean 𝝂 and covariance matrix 𝚻, any marginal or conditional distribution of the 𝑧’s will also be (multivariate) normal. For a univariate marginal distribution,

𝑧𝑗 = 𝒪(𝝂𝑗, 𝚻𝑗𝑗)

For a bivariate marginal distribution,

𝐳𝑗𝑘 = 𝒪 ((𝝂𝑗 𝝂𝑘 ) , (𝚻𝑗𝑗 𝚻𝑗𝑘 𝚻𝑘𝑗 𝚻𝑘𝑘 ))

For a 𝑙-dimensional marginal distribution,

𝐳𝑗,⋯,𝑙 = 𝒪 ⎛ ⎜ ⎝ ⎛ ⎜ ⎝ 𝝂𝑗 ⋮ 𝝂𝑙 ⎞ ⎟ ⎠ , ⎛ ⎜ ⎝ 𝚻𝑗𝑗 ⋯ 𝚻𝑗𝑙 ⋮ ⋱ ⋮ 𝚻𝑙𝑗 ⋯ 𝚻𝑙𝑙 ⎞ ⎟ ⎠ ⎞ ⎟ ⎠

6

Marginal distributions

Proposition - For an 𝑜-dimensional multivate normal with mean 𝝂 and covariance matrix 𝚻, any marginal or conditional distribution of the 𝑧’s will also be (multivariate) normal. For a univariate marginal distribution,

𝑧𝑗 = 𝒪(𝝂𝑗, 𝚻𝑗𝑗)

For a bivariate marginal distribution,

𝐳𝑗𝑘 = 𝒪 ((𝝂𝑗 𝝂𝑘 ) , (𝚻𝑗𝑗 𝚻𝑗𝑘 𝚻𝑘𝑗 𝚻𝑘𝑘 ))

For a 𝑙-dimensional marginal distribution,

𝐳𝑗,⋯,𝑙 = 𝒪 ⎛ ⎜ ⎝ ⎛ ⎜ ⎝ 𝝂𝑗 ⋮ 𝝂𝑙 ⎞ ⎟ ⎠ , ⎛ ⎜ ⎝ 𝚻𝑗𝑗 ⋯ 𝚻𝑗𝑙 ⋮ ⋱ ⋮ 𝚻𝑙𝑗 ⋯ 𝚻𝑙𝑙 ⎞ ⎟ ⎠ ⎞ ⎟ ⎠

6

Marginal distributions

Proposition - For an 𝑜-dimensional multivate normal with mean 𝝂 and covariance matrix 𝚻, any marginal or conditional distribution of the 𝑧’s will also be (multivariate) normal. For a univariate marginal distribution,

𝑧𝑗 = 𝒪(𝝂𝑗, 𝚻𝑗𝑗)

For a bivariate marginal distribution,

𝐳𝑗𝑘 = 𝒪 ((𝝂𝑗 𝝂𝑘 ) , (𝚻𝑗𝑗 𝚻𝑗𝑘 𝚻𝑘𝑗 𝚻𝑘𝑘 ))

For a 𝑙-dimensional marginal distribution,

𝐳𝑗,⋯,𝑙 = 𝒪 ⎛ ⎜ ⎝ ⎛ ⎜ ⎝ 𝝂𝑗 ⋮ 𝝂𝑙 ⎞ ⎟ ⎠ , ⎛ ⎜ ⎝ 𝚻𝑗𝑗 ⋯ 𝚻𝑗𝑙 ⋮ ⋱ ⋮ 𝚻𝑙𝑗 ⋯ 𝚻𝑙𝑙 ⎞ ⎟ ⎠ ⎞ ⎟ ⎠

6

Conditional Distributions

If we partition the 𝑜-dimensions into two pieces such that

𝐙 = (𝐙1, 𝐙2)𝑢 then

𝐙

𝑜×1 ∼ 𝒪 ⎛

⎜ ⎜ ⎝ (𝝂1 𝝂2 )

𝑜×1

, (𝚻11 𝚻12 𝚻21 𝚻22 )

𝑜×𝑜

⎞ ⎟ ⎟ ⎠ 𝐙1

𝑙×1

∼ 𝒪( 𝝂1

𝑙×1

, 𝚻11

𝑙×𝑙

) 𝐙2

𝑜−𝑙×1

∼ 𝒪( 𝝂2

𝑜−𝑙×1

, 𝚻22

𝑜−𝑙×𝑜−𝑙

) then the conditional distributions are given by 𝐙𝟐 | 𝐙2 = 𝐛 ∼ 𝒪(𝝂𝟐 + 𝚻𝟐𝟑 𝚻−1

𝟑𝟑 (𝐛 − 𝝂𝟑), 𝚻𝟐𝟐 − 𝚻𝟐𝟑 𝚻−1 𝟑𝟑 𝚻𝟑𝟐)

𝐙𝟑 | 𝐙1 = 𝐜 ∼ 𝒪(𝝂𝟑 + 𝚻𝟑𝟐 𝚻−1

𝟐𝟐 (𝐜 − 𝝂𝟐), 𝚻𝟑𝟑 − 𝚻𝟑𝟐 𝚻−1 𝟐𝟐 𝚻𝟑𝟐)

7

Conditional Distributions

If we partition the 𝑜-dimensions into two pieces such that

𝐙 = (𝐙1, 𝐙2)𝑢 then

𝐙

𝑜×1 ∼ 𝒪 ⎛

⎜ ⎜ ⎝ (𝝂1 𝝂2 )

𝑜×1

, (𝚻11 𝚻12 𝚻21 𝚻22 )

𝑜×𝑜

⎞ ⎟ ⎟ ⎠ 𝐙1

𝑙×1

∼ 𝒪( 𝝂1

𝑙×1

, 𝚻11

𝑙×𝑙

) 𝐙2

𝑜−𝑙×1

∼ 𝒪( 𝝂2

𝑜−𝑙×1

, 𝚻22

𝑜−𝑙×𝑜−𝑙

) then the conditional distributions are given by 𝐙𝟐 | 𝐙2 = 𝐛 ∼ 𝒪(𝝂𝟐 + 𝚻𝟐𝟑 𝚻−1

𝟑𝟑 (𝐛 − 𝝂𝟑), 𝚻𝟐𝟐 − 𝚻𝟐𝟑 𝚻−1 𝟑𝟑 𝚻𝟑𝟐)

𝐙𝟑 | 𝐙1 = 𝐜 ∼ 𝒪(𝝂𝟑 + 𝚻𝟑𝟐 𝚻−1

𝟐𝟐 (𝐜 − 𝝂𝟐), 𝚻𝟑𝟑 − 𝚻𝟑𝟐 𝚻−1 𝟐𝟐 𝚻𝟑𝟐)

7

Gaussian Processes

From Shumway, A process, 𝐙 = {𝑍 (𝑢) ∶ 𝑢 ∈ 𝑈}, is said to be a Gaussian process if all possible finite dimensional vectors 𝐳 = (𝑧𝑢1, 𝑧𝑢2, ..., 𝑧𝑢𝑜)𝑢, for every collection of time points 𝑢1, 𝑢2, … , 𝑢𝑜, and every positive integer 𝑜, have a multivariate normal distribution. So far we have only looked at examples of time series where 𝑈 is discete (and evenly spaces & contiguous), it turns out things get a lot more interesting when we explore the case where 𝑈 is defined on a continuous space (e.g.

• r some subset of

).

8

Gaussian Processes

From Shumway, A process, 𝐙 = {𝑍 (𝑢) ∶ 𝑢 ∈ 𝑈}, is said to be a Gaussian process if all possible finite dimensional vectors 𝐳 = (𝑧𝑢1, 𝑧𝑢2, ..., 𝑧𝑢𝑜)𝑢, for every collection of time points 𝑢1, 𝑢2, … , 𝑢𝑜, and every positive integer 𝑜, have a multivariate normal distribution. So far we have only looked at examples of time series where 𝑈 is discete (and evenly spaces & contiguous), it turns out things get a lot more interesting when we explore the case where 𝑈 is defined on a continuous space (e.g. R or some subset of R).

8

Gaussian Process Regression

Parameterizing a Gaussian Process

Imagine we have a Gaussian process defined such that

𝐙 = {𝑍 (𝑢) ∶ 𝑢 ∈ [0, 1]},

• We now have an uncountably infinite set of possible 𝑢’s and 𝑍 (𝑢)s.
• We will only have a (small) finite number of observations

𝑍 (𝑢1), … , 𝑍 (𝑢𝑜) with which to say something useful about this

infinite dimensional process.

• The unconstrained covariance matrix for the observed data can have

up to 𝑜(𝑜 + 1)/2 unique values∗

• Necessary to make some simplifying assumptions:
• Stationarity
• Simple parameterization of Σ

9

Parameterizing a Gaussian Process

Imagine we have a Gaussian process defined such that

𝐙 = {𝑍 (𝑢) ∶ 𝑢 ∈ [0, 1]},

• We now have an uncountably infinite set of possible 𝑢’s and 𝑍 (𝑢)s.
• We will only have a (small) finite number of observations

𝑍 (𝑢1), … , 𝑍 (𝑢𝑜) with which to say something useful about this

infinite dimensional process.

• The unconstrained covariance matrix for the observed data can have

up to 𝑜(𝑜 + 1)/2 unique values∗

• Necessary to make some simplifying assumptions:
• Stationarity
• Simple parameterization of Σ

9

Parameterizing a Gaussian Process

Imagine we have a Gaussian process defined such that

𝐙 = {𝑍 (𝑢) ∶ 𝑢 ∈ [0, 1]},

• We now have an uncountably infinite set of possible 𝑢’s and 𝑍 (𝑢)s.
• We will only have a (small) finite number of observations

𝑍 (𝑢1), … , 𝑍 (𝑢𝑜) with which to say something useful about this

infinite dimensional process.

• The unconstrained covariance matrix for the observed data can have

up to 𝑜(𝑜 + 1)/2 unique values∗

• Necessary to make some simplifying assumptions:
• Stationarity
• Simple parameterization of Σ

9

Parameterizing a Gaussian Process

Imagine we have a Gaussian process defined such that

𝐙 = {𝑍 (𝑢) ∶ 𝑢 ∈ [0, 1]},

• We now have an uncountably infinite set of possible 𝑢’s and 𝑍 (𝑢)s.
• We will only have a (small) finite number of observations

𝑍 (𝑢1), … , 𝑍 (𝑢𝑜) with which to say something useful about this

infinite dimensional process.

• The unconstrained covariance matrix for the observed data can have

up to 𝑜(𝑜 + 1)/2 unique values∗

• Necessary to make some simplifying assumptions:
• Stationarity
• Simple parameterization of Σ

9

Parameterizing a Gaussian Process

Imagine we have a Gaussian process defined such that

𝐙 = {𝑍 (𝑢) ∶ 𝑢 ∈ [0, 1]},

• We now have an uncountably infinite set of possible 𝑢’s and 𝑍 (𝑢)s.
• We will only have a (small) finite number of observations

𝑍 (𝑢1), … , 𝑍 (𝑢𝑜) with which to say something useful about this

infinite dimensional process.

• The unconstrained covariance matrix for the observed data can have

up to 𝑜(𝑜 + 1)/2 unique values∗

• Necessary to make some simplifying assumptions:
• Stationarity
• Simple parameterization of Σ

9

Covariance Functions

More on these next week, but for now some simple / common examples Exponential Covariance:

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − |𝑢 − 𝑢′| 𝑚 )

Squared Exponential Covariance:

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − (|𝑢 − 𝑢′| 𝑚 )2)

Powered Exponential Covariance (𝑞 ∈ (0, 2]):

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − (|𝑢 − 𝑢′| 𝑚 )𝑞)

10

Covariance Functions

More on these next week, but for now some simple / common examples Exponential Covariance:

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − |𝑢 − 𝑢′| 𝑚 )

Squared Exponential Covariance:

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − (|𝑢 − 𝑢′| 𝑚 )2)

Powered Exponential Covariance (𝑞 ∈ (0, 2]):

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − (|𝑢 − 𝑢′| 𝑚 )𝑞)

10

Covariance Functions

More on these next week, but for now some simple / common examples Exponential Covariance:

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − |𝑢 − 𝑢′| 𝑚 )

Squared Exponential Covariance:

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − (|𝑢 − 𝑢′| 𝑚 )2)

Powered Exponential Covariance (𝑞 ∈ (0, 2]):

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − (|𝑢 − 𝑢′| 𝑚 )𝑞)

10

Covariance Functions

More on these next week, but for now some simple / common examples Exponential Covariance:

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − |𝑢 − 𝑢′| 𝑚 )

Squared Exponential Covariance:

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − (|𝑢 − 𝑢′| 𝑚 )2)

Powered Exponential Covariance (𝑞 ∈ (0, 2]):

Σ(𝑧𝑢, 𝑧𝑢′) = 𝜏2 exp ( − (|𝑢 − 𝑢′| 𝑚 )𝑞)

10

Covariance Function - Correlation Decay

exp cov sq exp cov 0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00

d corr l

1 2 3 4 5 6 7 8 9 10

11

Correlation Decay - AR(1)

Recall that for a stationary AR(1) process:

𝛿(ℎ) = 𝜏2

𝑥𝜚|ℎ| and 𝜍(ℎ) = 𝜚|ℎ|

therefore we can draw a somewhat similar picture about the decay of 𝜍 as a function of distance.

0.25 0.50 0.75 1.00

rho phi

0.1 0.3 0.5 0.7 0.9

12

Example

−2 −1 1 0.00 0.25 0.50 0.75 1.00

t y 13

Prediction

Our example has 15 observations which we would like to use as the basis for predicting 𝑍 (𝑢) at other values of 𝑢 (say a sequence of values from 0 to 1). For now lets use a square exponential covariance with 𝜏2 = 10 and 𝑚 = 5 We therefore want to sample from 𝐙𝑞𝑠𝑓𝑒|𝐙𝑝𝑐𝑡.

𝐙𝑞𝑠𝑓𝑒 | 𝐙𝑝𝑐𝑡 = 𝐳 ∼ 𝒪(𝚻𝑞𝑝 𝚻−1

𝑝𝑐𝑡 𝐳, 𝚻𝐪𝐬𝐟𝐞 − 𝚻𝑞𝑝 𝚻−1 𝑞𝑠𝑓𝑒 𝚻𝑝𝑞)

14

Prediction

Our example has 15 observations which we would like to use as the basis for predicting 𝑍 (𝑢) at other values of 𝑢 (say a sequence of values from 0 to 1). For now lets use a square exponential covariance with 𝜏2 = 10 and 𝑚 = 5 We therefore want to sample from 𝐙𝑞𝑠𝑓𝑒|𝐙𝑝𝑐𝑡.

𝐙𝑞𝑠𝑓𝑒 | 𝐙𝑝𝑐𝑡 = 𝐳 ∼ 𝒪(𝚻𝑞𝑝 𝚻−1

𝑝𝑐𝑡 𝐳, 𝚻𝐪𝐬𝐟𝐞 − 𝚻𝑞𝑝 𝚻−1 𝑞𝑠𝑓𝑒 𝚻𝑝𝑞)

14

Prediction

Our example has 15 observations which we would like to use as the basis for predicting 𝑍 (𝑢) at other values of 𝑢 (say a sequence of values from 0 to 1). For now lets use a square exponential covariance with 𝜏2 = 10 and 𝑚 = 5 We therefore want to sample from 𝐙𝑞𝑠𝑓𝑒|𝐙𝑝𝑐𝑡.

𝐙𝑞𝑠𝑓𝑒 | 𝐙𝑝𝑐𝑡 = 𝐳 ∼ 𝒪(𝚻𝑞𝑝 𝚻−1

𝑝𝑐𝑡 𝐳, 𝚻𝐪𝐬𝐟𝐞 − 𝚻𝑞𝑝 𝚻−1 𝑞𝑠𝑓𝑒 𝚻𝑝𝑞)

14

Draw 1

−4 −2 0.00 0.25 0.50 0.75 1.00

t y 15

Draw 2

−4 −2 0.00 0.25 0.50 0.75 1.00

t y 16

Draw 3

−4 −2 2 0.00 0.25 0.50 0.75 1.00

t y 17

Draw 4

−4 −2 2 0.00 0.25 0.50 0.75 1.00

t y 18

Draw 5

−4 −2 2 0.00 0.25 0.50 0.75 1.00

t y 19

Many draws later

−5.0 −2.5 0.0 2.5 5.0 0.00 0.25 0.50 0.75 1.00

t y 20

Exponential Covariance

−7.5 −5.0 −2.5 0.0 2.5 5.0 0.00 0.25 0.50 0.75 1.00

t y 21

Exponential Covariance - Draw 2

−7.5 −5.0 −2.5 0.0 2.5 5.0 0.00 0.25 0.50 0.75 1.00

t y 22

Exponential Covariance - Draw 3

−4 4 0.00 0.25 0.50 0.75 1.00

t y 23

Exponential Covariance - Posterior

−5.0 −2.5 0.0 2.5 5.0 0.00 0.25 0.50 0.75 1.00

t y 24

Powered Exponential Covariance (𝑞 = 1.5)

−6 −3 3 6 0.00 0.25 0.50 0.75 1.00

t y 25

Back to the square exponential

−5.0 −2.5 0.0 2.5 5.0 0.00 0.25 0.50 0.75 1.00

t y 26

Changing the range (𝑚)

Sq Exp Cov − sigma2=10, l=15 Sq Exp Cov − sigma2=10, l=20 Sq Exp Cov − sigma2=10, l=5 Sq Exp Cov − sigma2=10, l=10 0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00 −5 5 −5 5

t y 27

Effective Range

For the square exponential covariance

𝐷𝑝𝑤(𝑒) = 𝜏2 exp (−(𝑚 ⋅ 𝑒)2) 𝐷𝑝𝑠𝑠(𝑒) = exp (−(𝑚 ⋅ 𝑒)2)

we would like to know, for a given value of 𝑚, beyond what distance apart must observations be to have a correlation less than 0.05?

28

Changing the scale (𝜏2)

Sq Exp Cov − sigma2=5, l=20 Sq Exp Cov − sigma2=15, l=20 Sq Exp Cov − sigma2=5, l=15 Sq Exp Cov − sigma2=15, l=15 0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00 −5 5 −5 5

t y 29

Fitting

gp_sq_exp_model = ”model{ y ~ dmnorm(mu, inverse(Sigma)) for (i in 1:N) { mu[i] <- 0 } for (i in 1:(N-1)) { for (j in (i+1):N) { Sigma[i,j] <- sigma2 * exp(- pow(l*d[i,j],2)) Sigma[j,i] <- Sigma[i,j] } } for (k in 1:N) { Sigma[k,k] <- sigma2 + 0.00001 } sigma2 ~ dlnorm(0, 1.5) l ~ dt(0, 2.5, 1) T(0,) # Half-cauchy(0,2.5) }”

30

Trace plots

l sigma2 250 500 750 1000 20 30 40 50 60 1 2 3 4

.iteration estimate term

l sigma2

param post_mean post_med post_lower post_upper l 30.20 28.70 20.63 51.51 sigma2 1.44 1.33 0.72 2.78

31

Fitted models

−2 2 0.00 0.25 0.50 0.75 1.00

t y

Post Mean Model

32

Forcasting

−3 −2 −1 1 2 0.0 0.5 1.0 1.5

t y

Post Mean Model

33

Improving the model

gp_sq_exp_model2 = ”model{ y ~ dmnorm(mu, inverse(Sigma)) for (i in 1:N) { mu[i] <- 0 } for (i in 1:(N-1)) { for (j in (i+1):N) { Sigma[i,j] <- sigma2 * exp(- pow(l*d[i,j],2)) Sigma[j,i] <- Sigma[i,j] } } for (k in 1:N) { Sigma[k,k] <- sigma2 + nugget } sigma2 ~ dlnorm(0, 1.5) l ~ dt(0, 2.5, 1) T(0,) # Half-cauchy(0,2.5) nugget ~ dlnorm(0, 1) }”

34

Trace plots

l nugget sigma2 250 500 750 1000 5 10 15 20 0.0 0.5 1.0 2 4 6 8

.iteration estimate term

l nugget sigma2

param post_mean post_med post_lower post_upper l 7.01 6.75 2.17 11.79 nugget 0.13 0.09 0.03 0.57 sigma2 1.73 1.53 0.64 4.04

35

Fitted models

Post Mean Model 0.00 0.25 0.50 0.75 1.00 −3 −2 −1 1

t y 36

Forcasting

Post Mean Model 0.0 0.5 1.0 1.5 −3 −2 −1 1 2 3

t y 37