Nothing
context("RACS simulation")
test_that("rbdr produces simulations with the correct coverage and covariance", {
grain <- owin(xrange = c(-5, 5), yrange = c(-5, 5))
win <- owin(xrange = c(0, 20), c(0, 20))
if (identical(Sys.getenv("NOT_CRAN"), "true")){
win <- owin(xrange = c(0, 200), c(0, 200))
}
lambda <- 4.2064E-3
#calculate theoretical values of the model
truecoveragefrac <- bdrcoverageprob(lambda, grain)
xy <- as.mask(dilationAny(win, win), eps = c(0.2, 0.2))
#eps of 1 (a 10 grain width) is too large! eps of 0.1
thcovariance <- bdrcovar(lambda, grain, xy)
## use covariance to form approximate confidence interval
setcovB <- setcov(win, eps = c(thcovariance$xstep, thcovariance$ystep))
harmims <- harmonise.im(thcovar = thcovariance, setcovwin = setcovB)
integrand <- eval.im((thcovariance - truecoveragefrac ^ 2) * setcovB, envir = list(thcovariance = harmims$thcovar, setcovB = harmims$setcovwin))
exactvariance <- (1 / (area.owin(win))^2) * sum(integrand) * integrand$xstep * integrand$ystep
##assume that window large enough that distribution of estimator is Gaussian
##So half-width of 0.05 conf interval is
confint_halfwidth <- -qnorm(0.025, mean = 0, sd = sqrt(exactvariance))
xi <- rbdr(lambda, grain, win)
phat <- coverageprob(xi, win)
# plot(xi, col = "black")
# plot(win, add = TRUE)
expect_equal(phat, truecoveragefrac, tolerance = confint_halfwidth)
#now on to checking covariance
cvchat <- racscovariance(as.mask(xi, eps = c(0.2, 0.2)), obswin = win, estimators = list("pickaint"), drop = TRUE)
truecvc.iso <- rotmean(thcovariance[disc(radius = 50), drop = FALSE], padzero = FALSE)
cvchat.iso <- rotmean(cvchat[disc(radius = 50), drop = FALSE], padzero = FALSE)
cvciso <- collapse.fv(harmonise.fv(True = truecvc.iso, Est = cvchat.iso), different = "f")
expect_lt(with.fv(cvciso, max(True - Est)), 10 * confint_halfwidth) #the 10 is just a guess based on 4th order properties being much harder to estimate
})
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.