update.homTP {hetGP}

Fast homTP-update

Description

Update existing homTP model with new observations

Usage

## S3 method for class 'homTP'
update(
  object,
  Xnew,
  Znew = NULL,
  lower = NULL,
  upper = NULL,
  noiseControl = NULL,
  known = NULL,
  maxit = 100,
  ...
)

Arguments

object

initial model of class homTP

Xnew

matrix of new design locations; ncol(Xnew) must match the input dimension encoded in object

Znew

vector new observations at those new design locations, of length nrow(X). NAs can be passed, see Details

lower, upper, noiseControl, known

optional bounds for MLE optimization, see mleHomTP. If not provided, they are extracted from the existing model

maxit

maximum number of iterations for the internal L-BFGS-B optimization method; see optim for more details

...

no other argument for this method.

Details

In case hyperparameters need not be updated, maxit can be set to 0. In this case it is possible to pass NAs in Znew, then the model can still be used to provide updated variance predictions.

Examples

## Not run: 
##------------------------------------------------------------
## Example : Sequential Homoskedastic TP moding 
##------------------------------------------------------------
set.seed(42)

## Spatially varying noise function
noisefun <- function(x, coef = 1){
  return(coef * (0.05 + sqrt(abs(x)*20/(2*pi))/10))
}

df_noise <- 3
nvar <- 1
n <- 10
X <- matrix(seq(0, 2 * pi, length=n), ncol = 1)
mult <- sample(1:50, n, replace = TRUE)
X <- rep(X, mult)
Z <- sin(X) + noisefun(X) * rt(length(X), df = df_noise)

testpts <- matrix(seq(0, 2*pi, length = 10*n), ncol = 1)
mod <- mod_init <- mleHomTP(X = X, Z = Z, covtype = "Matern5_2",
                                lower = rep(0.1, nvar), upper = rep(50, nvar))
preds <- predict(x = testpts, object = mod_init) 
plot(X, Z)
lines(testpts, preds$mean, col = "red")


nsteps <- 10
for(i in 1:nsteps){
  newIds <- sort(sample(1:(10*n), 5))
  
  newX <- testpts[rep(newIds, times = sample(1:50, length(newIds), replace = TRUE)), drop = FALSE] 
  newZ <- sin(newX) + noisefun(newX) * rt(length(newX), df = df_noise)
  points(newX, newZ, col = "blue", pch = 20)
  mod <- update(object = mod, newX, newZ)
  X <- c(X, newX)
  Z <- c(Z, newZ)
  plot(X, Z)
  print(mod$nit_opt)
}
p_fin <- predict(x = testpts, object = mod) 
lines(testpts, p_fin$mean, col = "blue")
lines(testpts, p_fin$mean + sqrt(p_fin$sd2) * qt(0.05, df = mod$nu + length(Z)),
      col = "blue", lty = 2)
lines(testpts, p_fin$mean + sqrt(p_fin$sd2) * qt(0.95, df = mod$nu + length(Z)),
      col = "blue", lty = 2)
lines(testpts, p_fin$mean + sqrt(p_fin$sd2 + p_fin$nugs) * qt(0.05, df = mod$nu + length(Z)),
      col = "blue", lty = 3)
lines(testpts, p_fin$mean + sqrt(p_fin$sd2 + p_fin$nugs) * qt(0.95, df = mod$nu + length(Z)),
      col = "blue", lty = 3)

mod_dir <-  mleHomTP(X = X, Z = Z, covtype = "Matern5_2",
                          lower = rep(0.1, nvar), upper = rep(50, nvar))
p_dir <- predict(x = testpts, object = mod_dir)
print(mod_dir$nit_opt)
lines(testpts, p_dir$mean, col = "green")
lines(testpts, p_dir$mean + sqrt(p_dir$sd2) * qt(0.05, df = mod_dir$nu + length(Z)),
      col = "green", lty = 2)
lines(testpts, p_dir$mean + sqrt(p_dir$sd2) * qt(0.95, df = mod_dir$nu + length(Z)),
      col = "green", lty = 2)
lines(testpts, p_dir$mean + sqrt(p_dir$sd2 + p_dir$nugs) * qt(0.05, df = mod_dir$nu + length(Z)),
      col = "green", lty = 3)
lines(testpts, p_dir$mean + sqrt(p_dir$sd2 + p_dir$nugs) * qt(0.95, df = mod_dir$nu + length(Z)),
      col = "green", lty = 3)

lines(testpts, sin(testpts), col = "red", lty = 2)

## Compare outputs
summary(mod_init)
summary(mod)
summary(mod_dir)



## End(Not run)