# True parameter values:
truebeta <- c(1,2,.2) # beta = (intercept, linear in x, linear in y)
truetheta <- c(.5,2,.02) # theta = (range, sill, nugget)
# Auxiliary covariance function, parameterized by
# theta = (range, sill, nugget)
spherical <- function(distmat, theta, eps = 1e-06) {
Sigma <- distmat
d <- Sigma@entries/theta[1]
Sigma@entries <- ifelse(d < eps,
theta[3]+ theta[2],
ifelse(d < 1, theta[2]*(1 - 1.5*d + 0.5*d^3), 0))
return( Sigma)
}
# We now define a grid, distance matrix, and a sample:
x <- seq(0,1,l=5)
locs <- expand.grid( x, x)
X <- as.matrix( cbind(1,locs)) # design matrix
Covariance <- 'spherical' # covariance function
distmat <- nearest.dist( locs, upper=NULL) # distance matrix
Sigma <- spherical( distmat, truetheta) # true covariance matrix
set.seed(15)
y <- c(rmvnorm.spam(1,X %*% truebeta,Sigma)) # construct sample
# Here is the negative 2 log likelihood:
neg2loglikelihood.spam( y, X, distmat, Covariance,
truebeta, truetheta)
# We pass now to the mle:
res <- mle.spam(y, X, distmat, Covariance,
truebeta, truetheta,thetalower=c(0,0,0),thetaupper=c(1,Inf,Inf))
# Similar parameter estimates here, of course:
mle.nomean.spam(y-X%*%res$par[1:3], distmat, Covariance,
truetheta, thetalower=c(0,0,0), thetaupper=c(1,Inf,Inf))Run the code above in your browser using DataLab