Nothing
m2LL.stream <- function(theta, m2LLdata, X, dist.hydro = NULL,
a.mat = NULL, b.mat = NULL, net.zero = NULL,
weight = NULL, x.dat = NULL, y.dat = NULL,
Del.i = NULL, A.5 = NULL,
CorModels, use.nugget, use.anisotropy, useTailDownWeight,
EstMeth = "REML",loglik.environment=loglik.environment,
REs=NULL, scale, maxrang = NULL)
{
if((max(theta) > 20) | (min(theta) < -20)) return(1e+32)
theta1 <- untrans.theta(theta = theta, scale = scale)
if(!is.null(maxrang))
if (any(theta1[!is.na(maxrang)] > maxrang[!is.na(maxrang)])) return(1e+32)
z <- m2LLdata
n <- length(z)
p <- length(X[1,])
if(length(CorModels) == 0) V <- diag(rep(theta1, times = n))
else V <- makeCovMat(theta1, dist.hydro = dist.hydro, w.matrix = weight,
a.mat = a.mat, b.mat = b.mat, net.zero = net.zero,
x.row = x.dat, y.row = y.dat,
x.col = x.dat, y.col = y.dat,
CorModels = CorModels, useTailDownWeight = useTailDownWeight,
use.nugget = use.nugget, use.anisotropy = use.anisotropy, REs)
V.eigenvalues <- eigen(V, symmetric=TRUE, only.values=TRUE)
if(any(V.eigenvalues$values <=0)) stop("covariance matrix is not positive definite")
if(!is.null(Del.i)) V <- Del.i*A.5*t((Del.i*A.5) * V)
qrV <- qr(V)
ViX <- solve(qrV,X)
covbi <- crossprod(X,ViX) ## Computationally more efficient than covbi <- t(X) %*% ViX
covb <- solve(covbi)
b.hat <- covb %*% crossprod(ViX,z) ##b.hat <- covb %*% t(ViX) %*% z
r <- z - X %*% b.hat
f1 <- t(r) %*% solve(qrV,r) + sum(log(abs(diag(qr.R(qrV)))))
if(EstMeth == "REML") f1 <- f1 + sum(log(svd(covbi)$d))
nmult <- (n - p)*log(2*pi)
if(EstMeth == "ML") nmult <- n*log(2*pi)
result <- f1 + nmult
## Add this value of theta and the -loglik to an environment for
## examination later
RESULT <- get("RESULT",loglik.environment)
RESULT <- rbind(RESULT,c(theta1,result))
assign("RESULT",RESULT,loglik.environment)
result
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.