Fit the model

formula <- y ~ 0 + b0 + covtemp + covprec + f(s, model = spde)
res <- inla(
  formula,
  family = "gaussian",
  data = inla.stack.data(stk.full),
  control.predictor = list(
    compute = TRUE,
    A = inla.stack.A(stk.full)
  ),
  control.compute = list(return.marginals.predictor = TRUE)
)