# NOT RUN {
#
# Trying out on a simulated data set
# Generate observation locations (50 is really too small but this is
# just to make this run quickly)
n <- 50
set.seed(124)
x = matrix(runif(2*n), nrow=n)
#generate observations at the locations
trueTheta = .1
trueSigma = .01
covMat = exp( -rdist(x,x) /trueTheta )
# As rcode: y = t(covMat)%*% (rnorm(n)) + trueSigma * rnorm( n)
y <- t(covMat) %*% (rnorm(n)) + trueSigma * rnorm( n)
# Use exponential covariance estimate constant function for mean
out <- MLESpatialProcess(x, y,
smoothness = .5,
mKrig.args = list( m = 1),
)
out$summary
# }
# NOT RUN {
#Now a (small) grid search over the smoothness
# add replicated data fields ( i.e. independent copies drawn from same covariance model)
# to make this stable
set.seed(124)
# Y = t(covMat)<!-- %*% matrix(rnorm(n*200), n,200) + trueSigma * matrix((rnorm(n*200)), n,200) -->
Y = t(covMat)%*% matrix(rnorm(n*200), n,200) +
trueSigma * matrix((rnorm(n*200)), n,200)
#This may take a few seconds
testSmoothness = c(.5, 1.5, 2.0)
for( nu in testSmoothness){
out = MLESpatialProcess(x, Y, cov.args = list(Covariance="Matern"),
smoothness = nu, cov.params.start= list( lambda=.5))
print( out$MLEJoint$summary)
}
# }
# NOT RUN {
# example with a covariate
# }
# NOT RUN {
data(COmonthlyMet)
ind<- !is.na( CO.tmean.MAM.climate)
x<- CO.loc[ind,]
y<- CO.tmean.MAM.climate[ind]
elev<- CO.elev[ind]
obj2<- MLESpatialProcess( x,y)
obj3<- MLESpatialProcess( x,y, Z=elev)
# elevation makes a difference
obj2$MLEJoint$summary
obj3$MLEJoint$summary
# }
# NOT RUN {
# }
# NOT RUN {
# fits for first 10 days from ozone data
data( ozone2)
NDays<- 5
O3MLE<- NULL
for( day in 1: NDays){
cat( day, " ")
ind<- !is.na(ozone2$y[day,] )
x<- ozone2$lon.lat[ind,]
y<- ozone2$y[day,ind]
print( length( y))
out<- MLESpatialProcess( x,y,
Distance="rdist.earth")$MLEJoint$summary
O3MLE<- rbind( O3MLE, out)
}
# NOTE: names from summary:
#[1] "lnProfileLike.FULL" "lambda" "theta"
#[4] "sigmaMLE" "rhoMLE" "funEval"
#[7] "gradEval" "eff.df" "GCV"
plot( log(O3MLE[,"lambda"]), log(O3MLE[,"theta"]))
# }
Run the code above in your browser using DataCamp Workspace