Lecture 14 Covariance Functions 3/08/2018 1 More on Covariance - PowerPoint PPT Presentation

Lecture 14 Covariance Functions 3/08/2018 1

More on Covariance Functions 2

Nugget Covariance 3 𝐷𝑝𝑤(𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = 𝜏 2 1 {ℎ=0} where ℎ = |𝑢 𝑗 − 𝑢 𝑘 | 2 1.00 1 0.75 draw 0 Draw 1 C 0.50 y Draw 2 −1 0.25 −2 0.00 0 5 10 15 20 0 5 10 15 20 h x

(- / Power / Square) Exponential Covariance 4 𝐷𝑝𝑤(𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = 𝜏 2 exp (−(ℎ 𝑚) 𝑞 ) where ℎ = |𝑢 𝑗 − 𝑢 𝑘 | Covariance − l=12, sigma2=1 Exponential 1.00 Cov 2 0.75 Exp 0 C 0.50 y Pow Exp (p=1.5) 0.25 Sq Exp −2 0.00 −4 0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00 h x Powered Exponential (p=1.5) Square Exponential 1 2 0 0 y y −1 −2 −2 −3 0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00 x x

Matern Covariance √ 2𝜉 ℎ ⋅ 𝑚) where ℎ = |𝑢 𝑗 −𝑢 𝑘 | 5 2𝜉 ℎ ⋅ 𝑚) √ 𝜉 𝐿 𝜉 ( 𝐷𝑝𝑤(𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = 𝜏 2 2 1−𝜉 Γ(𝜉) ( Covariance − l=2, sigma2=1 Matern − v=1/2 1.00 2 0.75 1 v=1/2 C 0.50 y 0 v=3/2 −1 v=5/2 0.25 −2 0.00 0 2 4 6 0 2 4 6 h x Matern − v=3/2 Matern − v=5/2 1 2 0 1 y y −1 0 −2 −1 −3 0 2 4 6 0 2 4 6 x x

Matern Covariance • A Gaussian process with Matérn covariance has sample functions that are ⌈𝜉 − 1⌉ times differentiable. (product of an exponential and a polynomial of order 𝑞 ). • When 𝜉 = 1/2 the Matern is equivalent to the exponential covariance. • As 𝜉 → ∞ the Matern converges to the square exponential covariance. • A Gaussian process with Matérn covariance has paths that are ⌈𝜉⌉ − 1 times differentiable. 6 • 𝐿 𝜉 is the modified Bessel function of the second kind. • When 𝜉 = 1/2 + 𝑞 for 𝑞 ∈ N + then the Matern has a simplified form

Rational Quadratic Covariance −𝛽 7 ) 𝛽 𝐷𝑝𝑤(𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = 𝜏 2 (1 + ℎ 2 𝑚 2 where ℎ = |𝑢 𝑗 − 𝑢 𝑘 | Covariance − l=12, sigma2=1 Rational Quadratic − alpha=1 1.00 1 0.75 0 alpha=1 0.50 y alpha=3 −1 alpha=10 0.25 −2 0.00 0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00 h x Rational Quadratic − alpha=10 Rational Quadratic − alpha=100 2 1 1 y y 0 0 −1 −1 0.00 0.25 0.50 0.75 1.00 0.00 0.25 0.50 0.75 1.00 x x

Rational Quadratic Covariance • is a scaled mixture of squared exponential covariance functions with different characteristic length-scales ( 𝑚 ). • As 𝛽 → ∞ the rational quadratic converges to the square exponential covariance. • Has sample functions that are infinitely differentiable for any value of 𝛽 8

Spherical Covariance 0 where ℎ = |𝑢 𝑗 −𝑢 𝑘 | otherwise 9 2 (ℎ ⋅ 𝑚) 3 )) 𝐷𝑝𝑤(𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = {𝜏 2 (1 − 3 2 ℎ ⋅ 𝑚 + 1 if 0 < ℎ < 1/𝑚 Covariance − sigma2=1 Spherical − l=1 1.00 1 0.75 0 l=1 0.50 y l=3 −1 0.25 l=10 −2 −3 0.00 0.0 0.3 0.6 0.9 0.0 0.3 0.6 0.9 h x Spherical − l=3 Spherical − l=10 2 2 1 y y 0 0 −1 −2 −2 0.0 0.3 0.6 0.9 0.0 0.3 0.6 0.9 x x

Periodic Covariance 10 𝐷𝑝𝑤(𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = 𝜏 2 exp (−2 𝑚 2 sin 2 (𝜌ℎ 𝑞 )) where ℎ = |𝑢 𝑗 − 𝑢 𝑘 | Covariance − l=2, sigma2=1 Periodic − p=1 1.00 1 forcats::as_factor(Cov) 0.75 p=1 0 0.50 y p=2 0.25 p=3 −1 0.00 0 1 2 3 4 0 2 4 6 h x Periodic − p=2 Periodic − p=3 1 2 1 0 y y 0 −1 −1 −2 −2 0 2 4 6 0 2 4 6 x x

Linear Covariance 𝐷𝑝𝑤(𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = 𝜏 2 11 𝑐 + 𝜏 2 𝑤 (𝑢 𝑗 − 𝑑)(𝑢 𝑘 − 𝑑) 1.0 0.5 0.0 y −0.5 −1.0 0.00 0.25 0.50 0.75 1.00 x

Combining Covariances If we definite two valid covariance functions, 𝐷𝑝𝑤 𝑏 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) and 𝐷𝑝𝑤 𝑐 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) then the following are also valid covariance functions, 𝐷𝑝𝑤 𝑏 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) + 𝐷𝑝𝑤 𝑐 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) 𝐷𝑝𝑤 𝑏 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) × 𝐷𝑝𝑤 𝑐 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) 12

Linear × Linear → Quadratic 13 𝐷𝑝𝑤 𝑏 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = 1 + 2 (𝑢 𝑗 × 𝑢 𝑘 ) 𝐷𝑝𝑤 𝑐 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = 2 + 1 (𝑢 𝑗 × 𝑢 𝑘 ) Cov_a * Cov_b 5 0 y −5 −10 −2 −1 0 1 2 x

14 Linear × Periodic 𝐷𝑝𝑤 𝑏 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = 1 + 1 (𝑢 𝑗 × 𝑢 𝑘 ) 𝐷𝑝𝑤 𝑐 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = exp (−2 sin 2 (2𝜌 ℎ)) Cov_a * Cov_b 2 0 y −2 −4 0 1 2 3 x

Linear + Periodic 15 𝐷𝑝𝑤 𝑏 (𝑧 𝑢 𝑗 , 𝑧 𝑢 𝑘 ) = 1 + 1 (𝑢 𝑗 × 𝑢 𝑘 ) 𝐷𝑝𝑤 𝑐 (ℎ = |𝑢 𝑗 − 𝑢 𝑘 |) = exp (−2 sin 2 (2𝜌 ℎ)) Cov_a + Cov_b 0 −1 draw −2 Draw 1 y Draw 2 −3 −4 −5 0 1 2 3 x

Sq Exp × Periodic → Locally Periodic 16 𝐷𝑝𝑤 𝑏 (ℎ = |𝑢 𝑗 − 𝑢 𝑘 |) = exp (−(1/3)ℎ 2 ) 𝐷𝑝𝑤 𝑐 (ℎ = |𝑢 𝑗 − 𝑢 𝑘 |) = exp (−2 sin 2 (𝜌 ℎ)) Cov_a * Cov_b 2 1 0 y −1 −2 0 2 4 6 x

Sq Exp (short) + Sq Exp (long) √ 3/2)ℎ 2 ) 17 3ℎ 2 ) √ 𝐷𝑝𝑤 𝑏 (ℎ = |𝑢 𝑗 − 𝑢 𝑘 |) = (1/4) exp (−4 𝐷𝑝𝑤 𝑐 (ℎ = |𝑢 𝑗 − 𝑢 𝑘 |) = exp (−( Cov_a + Cov_b 1 0 y −1 −2 0.0 2.5 5.0 7.5 10.0 x

Sq Exp (short) + Sq Exp (long) (Seen another way) 18 Cov_A (short) Cov_B (long) Cov_A + Cov_B 2 1 0 y −1 −2 −3 0.0 2.5 5.0 7.5 10.0 0.0 2.5 5.0 7.5 10.0 0.0 2.5 5.0 7.5 10.0 x

BDA3 example 19

BDA3 http://research.cs.aalto.fi/pml/software/gpstuff/demo_births.shtml 20

Births (one year) 1. Smooth long term trend ( sq exp cov ) 2. Seven day periodic trend with decay ( periodic × sq exp cov ) 3. Constant mean 21

Component Contributions We can view our GP in the following ways, but with appropriate conditioning we can also think of 𝐳 as being the sum of multipe independent GPs 𝐳 = 𝝂 + 𝑥 1 (𝐮) + 𝑥 2 (𝐮) + 𝑥 3 (𝐮) where 𝑥 1 (𝐮) ∼ 𝒪(0, 𝚻 1 ) 𝑥 2 (𝐮) ∼ 𝒪(0, 𝚻 2 ) 𝑥 3 (𝐮) ∼ 𝒪(0, 𝜏 2 𝐉 ) 22 𝐳 ∼ 𝒪(𝝂, 𝚻 1 + 𝚻 2 + 𝜏 2 𝐉 )

Decomposition of Covariance Components 0 ⎣ Σ 1 Σ 2 Σ 1 Σ 1 0 Σ 2 Σ 2 ⎦ ⎤ ⎥ ⎦ ⎞ ⎟ ⎠ therefore 𝑢 ⎡ ⎢ ⎥ ⎜ ⎢ ⎣ 𝑧 𝑥 1 𝑥 2 ⎤ ⎥ ⎦ 23 ⎝ ⎡ ⎢ ⎣ 𝝂 0 0 ⎤ Σ 1 + Σ 2 + 𝜏 2 𝐉 ∼ 𝒪 ⎛ , ⎡ 𝑥 1 | 𝐳, 𝝂, 𝜾 ∼ 𝒪(𝝂 𝑑𝑝𝑜𝑒 , 𝚻 𝑑𝑝𝑜𝑒 ) 𝝂 𝑑𝑝𝑜𝑒 = 0 + Σ 1 (Σ 1 + Σ 2 + 𝜏 2 𝐽) −1 (𝐳 − 𝝂) 𝚻 𝑑𝑝𝑜𝑒 = Σ 1 − Σ 1 (Σ 1 + Σ 2 + 𝜏 2 𝐉) −1 Σ 1

Births (multiple years) 1. slowly changing trend ( sq exp cov ) 2. small time scale correlating noise ( sq exp cov ) 3. 7 day periodical component capturing day of week effect ( periodic × sq exp cov ) 4. 365.25 day periodical component capturing day of year effect ( periodic × sq exp cov ) 5. component to take into account the special days and interaction with weekends ( linear cov ) 6. independent Gaussian noise ( nugget cov ) 7. constant mean 24

Mauna Loa Exampel 25

26 Atmospheric CO 2 390 Source NOAA y 360 Scripps (co2 in R) 330 1960 1980 2000 x

GP Model Based on Rasmussen 5.4.3 (we are using slightly different data and −𝛽 ) 𝛽 27 𝑧 ̄ parameterization) 𝐳 ∼ 𝒪(𝝂, 𝚻 1 + 𝚻 2 + 𝚻 3 + 𝚻 4 + 𝜏 2 I ) {𝝂} 𝑗 = {𝚻 1 } 𝑗𝑘 = 𝜏 2 1 exp (−(𝑚 1 ⋅ 𝑒 𝑗𝑘 ) 2 ) 2 exp (−(𝑚 2 ⋅ 𝑒 𝑗𝑘 ) 2 ) exp (−2 (𝑚 3 ) 2 sin 2 (𝜌 𝑒 𝑗𝑘 /𝑞)) {𝚻 2 } 𝑗𝑘 = 𝜏 2 3 (1 + (𝑚 4 ⋅ 𝑒 𝑗𝑘 ) 2 {𝚻 3 } 𝑗𝑘 = 𝜏 2 {𝚻 4 } 𝑗𝑘 = 𝜏 2 4 exp (−(𝑚 5 ⋅ 𝑒 𝑗𝑘 ) 2 )

JAGS Model } }” alpha ~ dt(0, 2.5, 1) T(0,) } l[i] ~ dt(0, 2.5, 1) T(0,) sigma2[i] ~ dt(0, 2.5, 1) T(0,) for(i in 1:5){ } Sigma[i,i] <- sigma2[1] + sigma2[2] + sigma2[3] + sigma2[4] + sigma2[5] for (i in 1:length(y)) { } ml_model = ”model{ Sigma[j,i] <- Sigma[i,j] Sigma[i,j] <- k1[i,j] + k2[i,j] + k3[i,j] + k4[i,j] k4[i,j] <- sigma2[4] * exp(- pow(l[5] * d[i,j],2)) k3[i,j] <- sigma2[3] * pow(1+pow(l[4] * d[i,j],2)/alpha, -alpha) k2[i,j] <- sigma2[2] * exp(- pow(l[2] * d[i,j],2) - 2 * pow(l[3] * sin(pi*d[i,j] / per), 2)) k1[i,j] <- sigma2[1] * exp(- pow(l[1] * d[i,j],2)) for (j in (i+1):length(y)) { for (i in 1:(length(y)-1)) { y ~ dmnorm(mu, inverse(Sigma)) 28

Diagnostics 29 sigma2[1] sigma2[2] sigma2[3] sigma2[4] sigma2[5] 0.8 40 0.04 2.0 6000 0.6 30 1.5 4000 0.03 0.4 20 1.0 2000 0.2 10 0.5 0.02 0 0.0 0.0 0 l[1] l[2] l[3] l[4] l[5] 0.020 1.2 6 0.06 0.9 estimate 0.015 1.0 0.6 4 0.04 0.010 0.8 0.3 2 0.005 0.02 0.6 0.0 0 250500750 1000 0 250500750 1000 0 250500750 1000 0 250500750 1000 alpha 8 6 4 2 0 0 250500750 1000 .iteration