#' #' Header for all (concatenated) test files #' #' Require spatstat.model #' Obtain environment variable controlling tests. #' #' $Revision: 1.5 $ $Date: 2020/04/30 05:31:37 $ require(spatstat.model) FULLTEST <- (nchar(Sys.getenv("SPATSTAT_TEST", unset="")) > 0) ALWAYS <- TRUE cat(paste("--------- Executing", if(FULLTEST) "** ALL **" else "**RESTRICTED** subset of", "test code -----------\n")) ## ## tests/marcelino.R ## ## $Revision: 1.4 $ $Date: 2020/04/30 02:18:23 $ ## local({ if(FULLTEST) { Y <- split(urkiola) B <- Y$birch O <- Y$oak B.lam <- predict (ppm(B ~polynom(x,y,2)), type="trend") O.lam <- predict (ppm(O ~polynom(x,y,2)), type="trend") Kinhom(B, lambda=B.lam, correction="iso") Kinhom(B, lambda=B.lam, correction="border") Kcross.inhom(urkiola, i="birch", j="oak", B.lam, O.lam) Kcross.inhom(urkiola, i="birch", j="oak", B.lam, O.lam, correction = "iso") Kcross.inhom(urkiola, i="birch", j="oak", B.lam, O.lam, correction = "border") } }) #' tests/mctests.R #' Monte Carlo tests #' (mad.test, dclf.test, envelopeTest, hasenvelope) #' $Revision: 1.5 $ $Date: 2022/05/23 04:09:49 $ local({ if(FULLTEST) { fitx <- ppm(redwood~x) ax <- envelopeTest(fitx, exponent=2, nsim=9, savefuns=TRUE) print(ax) } }) # # tests/mppm.R # # Basic tests of mppm # # $Revision: 1.23 $ $Date: 2022/07/08 04:55:13 $ # if(!FULLTEST) spatstat.options(npixel=32, ndummy.min=16) local({ ## test interaction formulae and subfits fit1 <- mppm(Points ~ group, simba, hyperframe(po=Poisson(), str=Strauss(0.1)), iformula=~ifelse(group=="control", po, str)) fit2 <- mppm(Points ~ group, simba, hyperframe(po=Poisson(), str=Strauss(0.1)), iformula=~str/id) fit2w <- mppm(Points ~ group, simba, hyperframe(po=Poisson(), str=Strauss(0.1)), iformula=~str/id, weights=runif(nrow(simba))) # currently invalid # fit3 <- mppm(Points ~ group, simba, # hyperframe(po=Poisson(), pie=PairPiece(c(0.05,0.1))), # iformula=~I((group=="control") * po) + I((group=="treatment") * pie)) fit1 fit2 fit2w # fit3 if(FULLTEST) { ## run summary.mppm which currently sits in spatstat-internal.Rd summary(fit1) summary(fit2) summary(fit2w) # summary(fit3) } ## test vcov algorithm vcov(fit1) vcov(fit2) # vcov(fit3) if(FULLTEST) { fit4 <- mppm(Points ~ group, simba, hyperframe(str=Strauss(0.1)), iformula=~str/group) fit4 summary(fit4) vcov(fit4) fit0 <- mppm(Points ~ group, simba) anova(fit0, fit4, test="Chi") ## [bug from Fernando Milesi] fit5 <- mppm(Wat ~ id, data=hyperframe(Wat=waterstriders), interaction=StraussHard(4.5, 1.5), iformula=~Interaction:id) fit5 summary(fit5) vcov(fit5) } ## test subfits algorithm if(FULLTEST) { s1 <- subfits(fit1) s2 <- subfits(fit2) # s3 <- subfits(fit3) s4 <- subfits(fit4) s5 <- subfits(fit5) ## validity of results of subfits() p1 <- solapply(s1, predict) p2 <- solapply(s2, predict) # p3 <- solapply(s3, predict) p4 <- solapply(s4, predict) p5 <- solapply(s5, predict) } }) local({ if(FULLTEST) { ## cases of predict.mppm W <- solapply(waterstriders, Window) Fakes <- solapply(W, runifpoint, n=30) FakeDist <- solapply(Fakes, distfun) H <- hyperframe(Bugs=waterstriders, D=FakeDist) fit <- mppm(Bugs ~ D, data=H) p1 <- predict(fit) p2 <- predict(fit, locations=Fakes) p3 <- predict(fit, locations=solapply(W, erosion, r=4)) locn <- as.data.frame(do.call(cbind, lapply(Fakes, coords))) df <- data.frame(id=sample(1:3, nrow(locn), replace=TRUE), D=runif(nrow(locn))) p4 <- predict(fit, locations=locn, newdata=df) fitG <- mppm(Bugs ~ D, data=H, use.gam=TRUE) p1G <- predict(fitG) p2G <- predict(fitG, locations=Fakes) p3G <- predict(fitG, locations=solapply(W, erosion, r=4)) p4G <- predict(fitG, locations=locn, newdata=df) } }) local({ ## [thanks to Sven Wagner] ## factor covariate, with some levels unused in some rows if(FULLTEST) { set.seed(14921788) H <- hyperframe(X=replicate(3, runifpoint(20), simplify=FALSE), Z=solist(as.im(function(x,y){x}, owin()), as.im(function(x,y){y}, owin()), as.im(function(x,y){x+y}, owin()))) H$Z <- solapply(H$Z, cut, breaks=(0:4)/2) fit6 <- mppm(X ~ Z, H) v6 <- vcov(fit6) s6 <- subfits(fit6) p6 <- solapply(s6, predict) ## random effects fit7 <- mppm(X ~ Z, H, random=~1|id) v7 <- vcov(fit7) s7 <- subfits(fit7) p7 <- solapply(s7, predict) fit7a <- mppm(X ~ Z, H, random=~x|id) v7a <- vcov(fit7a) s7a <- subfits(fit7a) p7a <- solapply(s7a, predict) ## multitype: collisions in vcov.ppm, predict.ppm H$X <- lapply(H$X, rlabel, labels=factor(c("a","b")), permute=FALSE) M <- MultiStrauss(matrix(0.1, 2, 2), c("a","b")) fit8 <- mppm(X ~ Z, H, M) v8 <- vcov(fit8, fine=TRUE) s8 <- subfits(fit8) p8 <- lapply(s8, predict) c8 <- lapply(s8, predict, type="cif") fit9 <- mppm(X ~ Z, H, M, iformula=~Interaction * id) v9 <- vcov(fit9, fine=TRUE) s9 <- subfits(fit9) p9 <- lapply(s9, predict) c9 <- lapply(s9, predict, type="cif") ## and a simple error in recognising 'marks' fit10 <- mppm(X ~ marks, H) } }) local({ if(FULLTEST) { ## test handling of offsets and zero cif values in mppm H <- hyperframe(Y = waterstriders) (fit1 <- mppm(Y ~ 1, data=H, Hardcore(1.5))) (fit2 <- mppm(Y ~ 1, data=H, StraussHard(7, 1.5))) (fit3 <- mppm(Y ~ 1, data=H, Hybrid(S=Strauss(7), H=Hardcore(1.5)))) s1 <- subfits(fit1) s2 <- subfits(fit2) s3 <- subfits(fit3) ## prediction, in training/testing context ## (example from Markus Herrmann and Ege Rubak) X <- waterstriders dist <- solapply(waterstriders, function(z) distfun(runifpoint(1, Window(z)))) i <- 3 train <- hyperframe(pattern = X[-i], dist = dist[-i]) test <- hyperframe(pattern = X[i], dist = dist[i]) fit <- mppm(pattern ~ dist, data = train) pred <- predict(fit, type="cif", newdata=test, verbose=TRUE) ## examples from Robert Aue GH <- Hybrid(G=Geyer(r=0.1, sat=3), H=Hardcore(0.01)) res <- mppm(Points ~ 1, interaction = GH, data=demohyper) print(summary(res)) sub <- subfits(res, verbose=TRUE) print(sub) } }) local({ if(FULLTEST) { ## test handling of interaction coefficients in multitype case set.seed(42) XX <- as.solist(replicate(3, rthin(amacrine, 0.8), simplify=FALSE)) H <- hyperframe(X=XX) M <- MultiStrauss(matrix(0.1, 2, 2), levels(marks(amacrine))) fit <- mppm(X ~ 1, H, M) co <- coef(fit) subco <- sapply(subfits(fit), coef) if(max(abs(subco - co)) > 0.001) stop("Wrong coefficient values in subfits, for multitype interaction") } }) local({ if(FULLTEST) { ## test lurking.mppm ## example from 'mppm' n <- 7 H <- hyperframe(V=1:n, U=runif(n, min=-1, max=1)) H$Z <- setcov(square(1)) H$U <- with(H, as.im(U, as.rectangle(Z))) H$Y <- with(H, rpoispp(eval.im(exp(2+3*Z)))) fit <- mppm(Y ~ Z + U + V, data=H) lurking(fit, expression(Z), type="P") lurking(fit, expression(V), type="raw") # design covariate lurking(fit, expression(U), type="raw") # image, constant in each row lurking(fit, H$Z, type="P") # list of images } }) local({ if(FULLTEST) { ## test anova.mppm code blocks and scoping problem H <- hyperframe(X=waterstriders) mod0 <- mppm(X~1, data=H, Poisson()) modxy <- mppm(X~x+y, data=H, Poisson()) mod0S <- mppm(X~1, data=H, Strauss(2)) modxyS <- mppm(X~x+y, data=H, Strauss(2)) anova(mod0, modxy, test="Chi") anova(mod0S, modxyS, test="Chi") anova(modxy, test="Chi") anova(modxyS, test="Chi") #' models with random effects (example from Marcelino de la Cruz) mod0r <- mppm(X~1, data=H, Poisson(), random = ~1|id) modxr <- mppm(X~x, data=H, Poisson(), random = ~1|id) anova(mod0r, modxr, test="Chi") } }) local({ if(FULLTEST) { ## test multitype stuff foo <- flu[1:3,] msh <- MultiStraussHard(iradii=matrix(100, 2, 2), hradii=matrix(10,2,2), types=levels(marks(foo$pattern[[1]]))) msh0 <- MultiStraussHard(iradii=matrix(100, 2, 2), hradii=matrix(10,2,2)) fit <- mppm(pattern ~ 1, data=foo, interaction=msh0) print(fit) print(summary(fit)) v <- vcov(fit) } }) reset.spatstat.options() #' #' tests/msr.R #' #' $Revision: 1.5 $ $Date: 2020/11/30 07:27:44 $ #' #' Tests of code for measures #' if(FULLTEST) { local({ ## cases of 'msr' Q <- quadscheme(cells) nQ <- n.quad(Q) nX <- npoints(cells) A <- matrix(nX * 3, nX, 3) B <- matrix(nQ * 3, nQ, 3) m <- msr(Q, A, B) M <- msr(Q, A, 1) M <- msr(Q, 1, B) M <- msr(Q, A, B[,1]) M <- msr(Q, A[,1], B) M <- msr(Q, A, B[,1,drop=FALSE]) M <- msr(Q, A[,1,drop=FALSE], B) ## methods a <- summary(m) b <- is.marked(m) w <- as.owin(m) z <- domain(m) ss <- scalardilate(m, 2) tt <- rescale(m, 2) ee <- rotate(m, pi/4) aa <- affine(m, mat=diag(c(1,2)), vec=c(0,1)) ff <- flipxy(m) am <- augment.msr(m, sigma=0.08) ua <- update(am) rr <- residuals(ppm(cells ~ x)) mm <- residuals(ppm(amacrine ~ x)) ss <- residuals(ppm(amacrine ~ x), type="score") gg <- rescale(ss, 1/662, c("micron", "microns")) plot(mm) plot(mm, multiplot=FALSE) plot(mm, equal.markscale=TRUE, equal.ribbon=TRUE) plot(ss) plot(ss, multiplot=FALSE) }) }