Nothing
inst/slowtests/test.plotmo.R
In plotmo: Plot a Model's Residuals, Response, and Partial Dependence Plots
# test.plotmo.R: regression tests for plotmo
# Stephen Milborrow, Petaluma Jan 2007
print(R.version.string)
source("test.prolog.R")
library(earth)
options(warn=1) # print warnings as they occur
data(etitanic)
make.space.for.caption <- function(caption="CAPTION")
{
oma <- par("oma")
needed <- 3
# adjust for newlines in caption
newlines <- grep("\n", caption)
if(length(newlines) > 0)
needed <- needed + .5 * newlines # .5 seems enough although 1 in theory
if(!is.null(caption) && any(nchar(caption)) && oma[3] <= needed) {
oma[3] <- needed
par(oma=oma)
}
}
dopar <- function(nrows, ncols, caption = "")
{
cat(" ", caption, "\n")
make.space.for.caption(caption)
par(mfrow=c(nrows, ncols))
par(mar = c(3, 3, 1.7, 0.5))
par(mgp = c(1.6, 0.6, 0))
par(cex = 0.7)
}
example(plotmo)
caption <- "basic earth test of plotmo"
a <- earth(O3 ~ ., data=ozone1, degree=2)
plotmo(a, degree1=2, degree2=4, caption=caption, trace=-1)
caption <- "test 5 x 5 layout"
dopar(1,1,caption)
a <- earth(O3 ~ ., data=ozone1, nk=51, pmethod="n", degree=2)
plotmo(a, caption=caption, trace=1)
caption <- "test 4 x 4 layout with ylab"
dopar(1,1,caption)
a <- earth(O3 ~ ., data=ozone1, nk=30, pmethod="n", degree=2)
plotmo(a, caption=caption, trace=2)
caption <- "test 3 x 3 layout"
dopar(1,1,caption)
a <- earth(O3 ~ ., data=ozone1, nk=16, pmethod="n", degree=2)
plotmo(a, caption=caption, trace=3)
caption <- "test 2 x 2 layout"
dopar(1,1,caption)
a <- earth(O3 ~ ., data=ozone1, nk=9, pmethod="n", degree=2)
plotmo(a, caption=caption)
caption <- "test 1 x 1 layout"
dopar(1,1,caption)
a <- earth(O3 ~ ., data=ozone1, nk=4, pmethod="n", degree=2)
plotmo(a, caption=caption)
caption <- "test plotmo basic params"
a <- earth(O3 ~ ., data=ozone1, degree=2)
dopar(3,2,caption)
plotmo(a, do.par=FALSE, degree1=1, nrug=-1, degree2=F, caption=caption,
main="test main", xlab="test xlab", ylab="test ylab")
plotmo(a, do.par=FALSE, degree1=F, degree2=4, grid.func=mean, persp.col="white", ngrid2=10, persp.phi=40)
set.seed(2016)
plotmo(a, do.par=FALSE, degree1=1, degree1.lty=2, degree1.lwd=4, degree1.col=2, nrug=TRUE, degree2=F, main="nrug=300")
plotmo(a, do.par=FALSE, degree1=1, nrug=-1, degree2=F, main="nrug=TRUE")
set.seed(2016)
plotmo(a, do.par=FALSE, degree1=1, nrug=10, ngrid1=50, degree2=F, main="ngrid1=50 nrug=10")
plotmo(a, do.par=FALSE, degree1=NA, degree2=1, persp.phi=60) # graph args
caption <- "test plotmo xlim and ylim"
a <- earth(O3 ~ ., data=ozone1, degree=2)
dopar(5,3,caption)
plotmo(a, do.par=FALSE, degree1=2:3, degree2=4, caption=caption, xlab="ylim=default")
plotmo(a, do.par=FALSE, degree1=2:3, degree2=4, ylim=NA, xlab="ylim=NA")
plotmo(a, do.par=FALSE, degree1=2:3, degree2=4, ylim=c(0,20), xlab="ylim=c(0,20)")
plotmo(a, do.par=FALSE, degree1=2:3, degree2=4, xlim=c(190,250), xlab="xlim=c(190,250)")
plotmo(a, do.par=FALSE, degree1=2:3, degree2=4, xlim=c(190,250), ylim=c(11,18), xlab="xlim=c(190,250), ylim=c(11,18)")
# check various types of predictors with grid.func and ndiscrete
varied.type.data <- data.frame(
y = 1:13,
num = c(1, 3, 2, 3, 4, 5, 6, 4, 5, 6.5, 3, 6, 5), # 7 unique values (but one is non integral)
int = c(1L, 1L, 3L, 3L, 4L, 4L, 3L, 5L, 3L, 6L, 7L, 8L, 10L), # 8 unique values
bool = c(F, F, F, F, F, T, T, T, T, T, T, T, T),
date = as.Date(
c("2018-08-01", "2018-08-02", "2018-08-03",
"2018-08-04", "2018-08-05", "2018-08-06",
"2018-08-07", "2018-08-08", "2018-08-08",
"2018-08-08", "2018-08-10", "2018-08-11", "2018-08-11")),
ord = ordered(c("ord3", "ord3", "ord3",
"ord1", "ord2", "ord3",
"ord1", "ord2", "ord3",
"ord1", "ord1", "ord1", "ord1"),
levels=c("ord1", "ord3", "ord2")),
fac = as.factor(c("fac1", "fac1", "fac1",
"fac2", "fac2", "fac2",
"fac3", "fac3", "fac3",
"fac1", "fac2", "fac3", "fac3")),
str = c("str1", "str1", "str1", # will be treated like a factor
"str2", "str2", "str2",
"str3", "str3", "str3",
"str3", "str3", "str3", "str3"))
varied.type.lm <- lm(y ~ ., data = varied.type.data)
print(summary(varied.type.lm))
set.seed(2018)
plotres(varied.type.lm, info=TRUE)
plotmo(varied.type.lm, pmethod="apartdep", all2=TRUE, ticktype="d", col.response="red", caption="varied.type.lm\npmethod=\"apartdep\" default grid func")
plotmo(varied.type.lm, all2=TRUE, ticktype="d", col.response="red", caption="varied.type.lm\ndefault grid func")
plotmo(varied.type.lm, all2=TRUE, ndiscre=1, caption="varied.type.lm\nndiscrete=1")
plotmo(varied.type.lm, all2=TRUE, ndiscr=2, caption="varied.type.lm\nndiscrete=2")
plotmo(varied.type.lm, all2=TRUE, ndis=100, caption="varied.type.lm\nndiscrete=100")
cat("grid.func=median:\n")
plotmo(varied.type.lm, all2=TRUE, grid.func=median, caption="varied.type.lm\ngrid.func=median")
cat("grid.func=quantile:\n")
plotmo(varied.type.lm, all2=TRUE, grid.func=function(x, ...) quantile(x, 0.5), caption="varied.type.lm\ngrid.func=function(x, ...) quantile(x, 0.5)")
cat("grid.func=mean:\n")
plotmo(varied.type.lm, all2=TRUE, grid.func=mean, caption="varied.type.lm\ngrid.func=mean")
varied.type.earth <- earth(y ~ ., data = varied.type.data, thresh=0, penalty=-1, trace=1)
print(summary(varied.type.earth))
set.seed(2018)
plotres(varied.type.earth, info=TRUE)
plotmo(varied.type.earth, all1=TRUE, all2=TRUE, persp.ticktype="d", col.response="red")
# term.plot calls predict.earth with an se parameter, even with termplot(se=FALSE)
caption <- "basic earth test against termplot"
dopar(4,4,caption)
make.space.for.caption("test caption1")
a <- earth(O3 ~ ., data=ozone1, degree=2)
plotmo(a, do.par=FALSE, ylim=NA, caption=caption, degree2=FALSE)
cat("Ignore warning: predict.earth ignored argument \"se.fit\"\n")
termplot(a)
caption <- "test change order of earth predictors and cex"
dopar(4,4,caption)
# minspan=1 to force two degree2 graphs for the test (wasn't necessary in old versions of earth)
a <- earth(doy ~ humidity + temp + wind, data=ozone1, degree=2, minspan=1)
plotmo(a, do.par=FALSE, ylim=NA, caption=caption, degree2=c(1,2), cex=1.2)
termplot(a)
caption <- "test all1=TRUE"
a <- earth(doy ~ humidity + temp + wind, data=ozone1, degree=2)
plotmo(a, caption=caption, all1=TRUE, persp.ticktype="d", persp.nticks=2)
caption <- "test all2=TRUE"
print(summary(a))
plotmo(a, caption=caption, all2=TRUE)
oz <- ozone1[150:200,c("O3","temp","humidity","ibh")]
a.glob <- earth(O3~temp+humidity, data=oz, degree=2)
ad.glob <- earth(oz[,2:3], oz[,1], degree=2)
func1 <- function()
{
caption <- "test environments and finding the correct data"
dopar(4,4,caption)
set.seed(2016)
plotmo(a.glob, do.par=FALSE, main="a.glob oz",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20, trace=2)
mtext(caption, outer=TRUE, font=2, line=1.5, cex=1)
plotmo(ad.glob, do.par=FALSE, main="ad.glob oz",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pch.response=20, trace=2) # pch.response test backcompat
a <- earth(O3~temp+humidity, data=oz, degree=2)
plotmo(a, do.par=FALSE, main="a oz",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
ad <- earth(oz[,2:3], oz[,1], degree=2)
plotmo(ad, do.par=FALSE, main="ad oz",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
oz.org <- oz
oz10 <- 10 * oz # multiply by 10 so we can see by the axis labels if right data is being used
oz <- oz10 # oz is now local to this function, but multiplied by 10
a.oz10 <- earth(O3~temp+humidity, data=oz, degree=2)
a.oz10.keep <- earth(O3~temp+humidity, data=oz, degree=2, keepxy=TRUE)
plotmo(a.oz10, do.par=FALSE, main="a oz10",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
ad.oz10 <- earth(oz[,2:3], oz[,1], degree=2)
ad.oz10.keep <- earth(oz[,2:3], oz[,1], degree=2, keepxy=TRUE)
plotmo(ad.oz10, do.par=FALSE, main="ad oz10",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
func2 <- function() {
a.func <- earth(O3 ~ temp + humidity, data=oz10, degree=2)
plotmo(a.func, do.par=FALSE, main="a.func oz10",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
ad.func <- earth(oz10[,2:3], oz10[,1], degree=2)
plotmo(ad.func, do.par=FALSE, main="ad.func oz10",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
caption <- "test environments and finding the correct data, continued"
dopar(4,4,caption)
oz <- .1 * oz.org
a.func <- earth(O3~temp+ humidity , data=oz, degree=2)
plotmo(a.func, do.par=FALSE, main="a.func oz.1",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
ad.func <- earth(oz[,2:3], oz[,1], degree=2)
plotmo(ad.func, do.par=FALSE, main="ad.func oz.1",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
plotmo(a.oz10.keep, do.par=FALSE, main="func1:a.oz10.keep",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
plotmo(ad.oz10.keep, do.par=FALSE, main="func1:ad.oz10.keep",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
cat("Expect error msg: formal argument \"do.par\" matched by multiple actual arguments\n")
expect.err(try(plotmo(a.oz10, do.par=FALSE, main="func1:a.oz10",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20, do.par=FALSE)))
}
func2()
y <- 3:11
x1 <- c(1,3,2,4,5,6,6,6,6)
x2 <- c(2,3,4,5,6,7,8,9,10)
frame <- data.frame(y=y, x1=x1, x2=x2)
foo <- function()
{
lm.18.out <- lm(y~x1+x2, model=FALSE)
x1[2] <- 18
y[3] <- 19
frame <- data.frame(y=y, x1=x1, x2=x2)
list(lm.18.out = lm.18.out,
lm.18 = lm(y~x1+x2),
lm.18.keep = lm(y~x1+x2, x=TRUE, y=TRUE),
lm.18.frame = lm(y~x1+x2, data=frame))
}
temp <- foo()
lm.18.out <- temp$lm.18.out
lm.18 <- temp$lm.18
lm.18.keep <- temp$lm.18.keep
lm.18.frame <- temp$lm.18.frame
# following should all use the x1 and y inside foo
cat("==lm.18.out\n")
plotmo(lm.18.out, main="lm.18.out",
do.par=FALSE, degree1=1, clip=FALSE, ylim=c(0,20),
col.response=2, pt.pch=20)
cat("==lm.18\n")
plotmo(lm.18, main="lm.18",
do.par=FALSE, degree1=1, clip=FALSE, ylim=c(0,20),
col.response=2, pt.pch=20)
cat("==lm.18.keep\n")
plotmo(lm.18.keep, main="lm.18.keep", trace=2,
do.par=FALSE, degree1=1, clip=FALSE, ylim=c(0,20),
col.response=2, pt.pch=20)
cat("==lm.18.frame\n")
plotmo(lm.18.frame, main="lm.18.frame",
do.par=FALSE, degree1=1, clip=FALSE, ylim=c(0,20),
col.response=2, pt.pch=20)
}
func1()
caption <- "test earth formula versus x,y model"
# dopar(4,4,caption)
# mtext(caption, outer=TRUE, font=2, line=1.5, cex=1)
a <- earth(O3 ~ ., data=ozone1, degree=2)
plotmo(a, caption="test earth formula versus xy model (formula)")
a <- earth(ozone1[, -1], ozone1[,1], degree=2)
plotmo(a, caption="test earth formula versus xy model (xy)")
# single predictor
caption <- "test earth(O3~wind, data=ozone1, degree=2), single predictor"
dopar(2,2,caption)
a <- earth(O3~wind, data=ozone1, degree=2)
plotmo(a)
caption = "se=2, earth(doy~humidity+temp+wind, data=ozone1) versus termplot (expect no se lines)"
dopar(3,3,caption)
mtext(caption, outer=TRUE, font=2, line=1.5, cex=1)
# minspan=1 to force two degree2 graphs for the test (wasn't necessary in old versions of earth)
a <- earth(doy~humidity + temp + wind, data=ozone1, degree=2, minspan=1)
cat("Ignore warning: predict.earth ignored argument \"se\"\n")
termplot(a)
plotmo(a, do.par=FALSE, ylim=NA, degree2=c(1:2), clip=FALSE, caption=caption)
# test fix to bug reported by Joe Retzer, FIXED Dec 7, 2007
N <- 650
set.seed(2007)
q_4 <- runif(N, -1, 1)
q_2102 <- runif(N, -1, 1)
q_2104 <- runif(N, -1, 1)
q_3105 <- runif(N, -1, 1)
q_3106 <- runif(N, -1, 1)
q_4104 <- runif(N, -1, 1)
q_6101 <- runif(N, -1, 1)
q_6103 <- runif(N, -1, 1)
q_7104 <- runif(N, -1, 1)
q_3109 <- runif(N, -1, 1)
q_4103 <- runif(N, -1, 1)
q_2111 <- runif(N, -1, 1)
q_3107 <- runif(N, -1, 1)
q_3101 <- runif(N, -1, 1)
q_3104 <- runif(N, -1, 1)
q_7107 <- runif(N, -1, 1)
depIndex <- sin(1.0 * q_4 + rnorm(650, sd=.8)) + sin(1.8 * q_2102 + rnorm(650, sd=.8)) + sin(1.3 * q_2104 + rnorm(650, sd=.8)) + sin(1.4 * q_3105 + rnorm(650, sd=.8)) +
sin(1.5 * q_3106 + rnorm(650, sd=.8)) + sin(1.6 * q_4104 + rnorm(650, sd=.8)) + sin(1.8 * q_6101 + rnorm(650, sd=.8)) + sin(1.8 * q_6103 + rnorm(650, sd=.8)) +
sin(1.9 * q_7104 + rnorm(650, sd=.8)) + sin(2.0 * q_3109 + rnorm(650, sd=.8))
regDatCWD <- as.data.frame(cbind(depIndex, q_4, q_2102, q_2104, q_3105, q_3106, q_4104, q_6101, q_6103, q_7104, q_3109, q_4103, q_2111, q_3107, q_3101, q_3104, q_7107))
cat("--plotmo(earthobj5)--\n")
earthobj5 <- earth(depIndex ~ q_4+q_2102+q_2104+q_3105+q_3106+q_4104+q_6101+q_6103+q_7104+q_3109+q_4103+q_2111+q_3107+q_3101+q_3104+q_7107, data=regDatCWD)
print(summary(earthobj5, digits = 2))
plotmo(earthobj5)
# long predictor names
a.rather.long.in.fact.very.long.name.q_4 <- q_4
a.rather.long.in.fact.very.long.name.q_2102 <- q_2102
a.rather.long.in.fact.very.long.name.q_2104 <- q_2104
a.rather.long.in.fact.very.long.name.q_3105 <- q_3105
a.rather.long.in.fact.very.long.name.q_3106 <- q_3106
a.rather.long.in.fact.very.long.name.q_4104 <- q_4104
a.rather.long.in.fact.very.long.name.q_6101 <- q_6101
a.rather.long.in.fact.very.long.name.q_6103 <- q_6103
a.rather.long.in.fact.very.long.name.q_7104 <- q_7104
a.rather.long.in.fact.very.long.name.q_3109 <- q_3109
a.rather.long.in.fact.very.long.name.q_4103 <- q_4103
a.rather.long.in.fact.very.long.name.q_2111 <- q_2111
a.rather.long.in.fact.very.long.name.q_3107 <- q_3107
a.rather.long.in.fact.very.long.name.q_3101 <- q_3101
a.rather.long.in.fact.very.long.name.q_3104 <- q_3104
a.rather.long.in.fact.very.long.name.q_7107 <- q_7107
a.rather.long.in.fact.very.long.name.for.the.response <- depIndex
a.rather.long.in.fact.very.long.name.for.the.dataframe <-
as.data.frame(cbind(
a.rather.long.in.fact.very.long.name.for.the.response,
a.rather.long.in.fact.very.long.name.q_4,
a.rather.long.in.fact.very.long.name.q_2102,
a.rather.long.in.fact.very.long.name.q_2104,
a.rather.long.in.fact.very.long.name.q_3105,
a.rather.long.in.fact.very.long.name.q_3106,
a.rather.long.in.fact.very.long.name.q_4104,
a.rather.long.in.fact.very.long.name.q_6101,
a.rather.long.in.fact.very.long.name.q_6103,
a.rather.long.in.fact.very.long.name.q_7104,
a.rather.long.in.fact.very.long.name.q_3109,
a.rather.long.in.fact.very.long.name.q_4103,
a.rather.long.in.fact.very.long.name.q_2111,
a.rather.long.in.fact.very.long.name.q_3107,
a.rather.long.in.fact.very.long.name.q_3101,
a.rather.long.in.fact.very.long.name.q_3104,
a.rather.long.in.fact.very.long.name.q_7107))
cat("--a.rather.long.in.fact.very.long.name.for.the...A--\n")
a.rather.long.in.fact.very.long.name.for.the.modelA <-
earth(a.rather.long.in.fact.very.long.name.for.the.response ~
a.rather.long.in.fact.very.long.name.q_4 +
a.rather.long.in.fact.very.long.name.q_2102 +
a.rather.long.in.fact.very.long.name.q_2104 +
a.rather.long.in.fact.very.long.name.q_3105 +
a.rather.long.in.fact.very.long.name.q_3106 +
a.rather.long.in.fact.very.long.name.q_4104 +
a.rather.long.in.fact.very.long.name.q_6101 +
a.rather.long.in.fact.very.long.name.q_6103 +
a.rather.long.in.fact.very.long.name.q_7104 +
a.rather.long.in.fact.very.long.name.q_3109 +
a.rather.long.in.fact.very.long.name.q_4103 +
a.rather.long.in.fact.very.long.name.q_2111 +
a.rather.long.in.fact.very.long.name.q_3107 +
a.rather.long.in.fact.very.long.name.q_3101 +
a.rather.long.in.fact.very.long.name.q_3104 +
a.rather.long.in.fact.very.long.name.q_7107,
data = a.rather.long.in.fact.very.long.name.for.the.dataframe)
print(summary(a.rather.long.in.fact.very.long.name.for.the.modelA, digits = 2))
plot(a.rather.long.in.fact.very.long.name.for.the.modelA)
plotmo(a.rather.long.in.fact.very.long.name.for.the.modelA)
cat("--a.rather.long.in.fact.very.long.name.for.the...C--\n")
a.rather.long.in.fact.very.long.name.for.the.modelC <-
earth(x = a.rather.long.in.fact.very.long.name.for.the.dataframe[,-1],
y = a.rather.long.in.fact.very.long.name.for.the.response,
degree = 3)
print(summary(a.rather.long.in.fact.very.long.name.for.the.modelC, digits = 2))
plot(a.rather.long.in.fact.very.long.name.for.the.modelC)
plotmo(a.rather.long.in.fact.very.long.name.for.the.modelC)
a <- earth(survived ~ pclass+sex+age, data=etitanic, degree=2)
print(summary(a))
plotmo(a, caption="plotmo with facs: pclass+sex+age")
plotmo(a, caption="plotmo with facs: pclass+sex+age, all1=T, grid.col=\"gray\"", all1=T, grid.col="gray")
plotmo(a, caption="plotmo with facs: pclass+sex+age, all2=T, col.grid=\"green\"", all2=T, col.grid="green")
plotmo(a, caption="plotmo with facs: pclass+sex+age, all1=T, all2=T, grid=2", all1=T, all2=T, grid.col=2)
plotmo(a, clip=FALSE, degree2=FALSE, caption="plotmo (no degree2) with facs: pclass+sex+age")
plotmo(a, clip=FALSE, grid.levels=list(pclass="2n", sex="ma"),
caption="plotmo with grid.levels: pclass+sex+age")
# in above tests, all degree2 terms use facs
# now build a model with some degree2 term that use facs, some that don't
a <- earth(survived ~ pclass+age+sibsp, data=etitanic, degree=2)
print(summary(a))
plotmo(a, caption="plotmo with mixed fac and non-fac degree2 terms", persp.border=NA)
plotmo(a, caption="plotmo with mixed fac and non-fac degree2 terms and grid.levels",
grid.levels=list(pclass="2n", age=20), # test partial matching of grid levels, and numeric preds
persp.ticktype="d", persp.nticks=2)
# check detection of illegal grid.levels argument
expect.err(try(plotmo(a, grid.levels=list(pcla="1", pclass="2"))), 'illegal grid.levels argument ("pcla" and "pclass" both match "pclass")')
expect.err(try(plotmo(a, grid.levels=list(pclass="1", pcla="2"))), 'illegal grid.levels argument ("pclass" and "pcla" both match "pclass")')
expect.err(try(plotmo(a, grid.levels=list(pcla="nonesuch"))), 'illegal level "nonesuch" for "pclass" in grid.levels (allowed levels are "1st" "2nd" "3rd")')
expect.err(try(plotmo(a, grid.levels=list(pcla="1sx"))), 'illegal level "1sx" for "pclass" in grid.levels (allowed levels are "1st" "2nd" "3rd")')
expect.err(try(plotmo(a, grid.levels=list(pcla=1))), 'illegal level for "pclass" in grid.levels (specify factor levels with a string)')
expect.err(try(plotmo(a, grid.levels=list(pcla=c("ab", "cd")))), "length(pclass) in grid.levels is not 1")
expect.err(try(plotmo(a, grid.levels=list(pcla=NA))), 'pclass in grid.levels is NA')
expect.err(try(plotmo(a, grid.levels=list(pcla=Inf))), 'pclass in grid.levels is infinite')
expect.err(try(plotmo(a, grid.levels=list(pcla=9))), 'illegal level for "pclass" in grid.levels (specify factor levels with a string)')
options(warn=2)
expect.err(try(plotmo(a, grid.levels=list(age="ab"))), 'grid.levels returned class \"character\" for age, so will use the default grid.func for age')
options(warn=1)
expect.err(try(plotmo(a, grid.levels=list(age=NA))), 'age in grid.levels is NA')
expect.err(try(plotmo(a, grid.levels=list(age=Inf))), 'age in grid.levels is infinite')
expect.err(try(plotmo(a, grid.lev=list(age=list(1,2)))), 'length(age) in grid.levels is not 1')
# more-or-less repeat above, but with glm models
a <- earth(survived ~ pclass+age+sibsp, data=etitanic, degree=2, glm=list(family=binomial))
print(summary(a))
plotmo(a, ylim=c(0, 1), caption="plotmo glm with mixed fac and non-fac degree2 terms")
plotmo(a, ylim=c(0, 1), caption="plotmo glm with mixed fac and non-fac degree2 terms and grid.levels",
grid.levels=list(pcl="2nd")) # test partial matching of variable name in grid levels
plotmo(a, type="earth", ylim=c(0, 1), caption="type=\"earth\" plotmo glm with mixed fac and non-fac degree2 terms")
plotmo(a, type="link", ylim=c(0, 1), clip=FALSE, caption="type=\"link\" plotmo glm with mixed fac and non-fac degree2 terms")
plotmo(a, type="class", ylim=c(0, 1), caption="type=\"class\" plotmo glm with mixed fac and non-fac degree2 terms")
plotmo(a, ylim=c(0, 1), caption="default type (\"response\")\nplotmo glm with mixed fac and non-fac degree2 terms")
# now with different type2s
set.seed(2016)
plotmo(a, do.par=FALSE, type2="persp", persp.theta=-20, degree1=FALSE, grid.levels=list(pclass="2nd"))
mtext("different type2s", outer=TRUE, font=2, line=1.5, cex=1)
plotmo(a, do.par=FALSE, type2="contour", degree1=FALSE, grid.levels=list(pclass="2nd"))
plotmo(a, do.par=FALSE, type2="image", degree1=FALSE, grid.levels=list(pclass="2nd"),
col.response=as.numeric(etitanic$survived)+2, pt.pch=20)
plotmo(a, do.par=FALSE, type="earth", type2="image", degree1=FALSE,
grid.levels=list(pclass="2"))
# grid.levels with partdep
set.seed(2018)
x1 <- (1:11) + runif(11)
x2 <- (1:11) + runif(11)
x3 <- as.integer((1:11) + runif(11))
x4 <- runif(11) > .5 # logical
y <- x1 - x2 + x3 + x4
data <- data.frame(y=y, x1=x1, x2=x2, x3=x3, x4=x4)
lm.x1.x2.x3 <- lm(y ~ x1 + x2 + x3 + x4 + x1*x2 + x1*x3, data=data)
cat("summary(lm.x1.x2.x3):\n")
print(summary(lm.x1.x2.x3))
par(mfrow = c(5, 6), mar = c(2, 3, 2, 1), mgp = c(1.5, 0.5, 0), cex = 0.6, oma=c(0,0,8,0))
plotmo(lm.x1.x2.x3, do.par=0, ylim=c(0,16), pt.col=2,
caption="row1 default\nrow2 grid.levels=list(x3=15)\nrow3 partdep\nrow4 partdetp grid.levels=list(x3=15)")
plotmo(lm.x1.x2.x3, do.par=0, ylim=c(0,16), pt.col=2, grid.levels=list(x3=15))
plotmo(lm.x1.x2.x3, do.par=0, ylim=c(0,16), pt.col=2, pmethod="partdep")
plotmo(lm.x1.x2.x3, do.par=0, ylim=c(0,16), pt.col=2, pmethod="partdep", grid.levels=list(x3=15))
# check auto type convert in grid.levels
plotmo(lm.x1.x2.x3, degree1="x1", degree2=0, main="x1 (x2=5L))", ylim=c(0,16), do.par=0, pmethod="partdep", grid.levels=list(x2=15L)) # integer to numeric
plotmo(lm.x1.x2.x3, degree1="x1", degree2=0, main="x1 (x3=5))", ylim=c(0,16), do.par=0, pmethod="partdep", grid.levels=list(x3=15)) # numeric to integer
plotmo(lm.x1.x2.x3, degree1="x1", degree2=0, main="x1 (x4=1))", ylim=c(0,16), do.par=0, pmethod="partdep", grid.levels=list(x4=1)) # numeric to logical
expect.err(try(plotmo(lm.x1.x2.x3, degree1="x1", degree2=0, main="x1 (x4=1))", ylim=c(0,16), do.par=0, pmethod="partdep", grid.levels=list(x4="x"))), "expected a logical value in grid.levels for x4") # char to logical
expect.err(try(plotmo(lm.x1.x2.x3, degree1="x2", do.par=0, pmethod="partdep", grid.levels=list(x1="1"))), "the class \"character\" of \"x1\" in grid.levels does not match its class \"numeric\" in the input data")
par(org.par)
# test vector main
a20 <- earth(O3 ~ humidity + temp + doy, data=ozone1, degree=2, glm=list(family=Gamma))
dopar(2, 2)
plotmo(a20, nrug=-1)
set.seed(2016)
plotmo(a20, nrug=10, caption="Test plotmo with a vector main (and npoints=200)",
main=c("Humidity", "Temperature", "Day of year", "Humidity: Temperature", "Temperature: Day of Year"),
col.response="darkgray", pt.pch=".", cex.response=3, npoints=200) # cex.response tests back compat
cat("Expect warning below (missing double titles)\n")
plotmo(a20, nrug=-1, caption="Test plotmo with a vector main (and plain smooth)",
main=c("Humidity", "Temperature", "Day of year", "Humidity: Temperature", "Temp: Doy"),
smooth.col="indianred")
cat("Expect warning below (missing single titles)\n")
plotmo(a20, nrug=-1, caption="Test plotmo with a vector main (and smooth args)",
main=c("Humidity", "Temperature"),
smooth.col="indianred", smooth.lwd=2, smooth.lty=2, smooth.f=.1,
col.response="gray", npoints=500)
plotmo(a20, nrug=-1, caption="Test plotmo with pt.pch=paste(1:nrow(ozone1))",
type2="im",
col.response=2, pt.cex=.8, pt.pch=paste(1:nrow(ozone1)), npoints=100)
aflip <- earth(O3~vh + wind + humidity + temp, data=ozone1, degree=2)
# test all1 and all2, with and without degree1 and degree2
plotmo(aflip, all2=T, caption="all2=T", npoints=TRUE)
plotmo(aflip, all2=T, degree2=c(4, 2), caption="all2=T, degree2=c(4, 2)")
plotmo(aflip, all1=T, caption="all1=T")
plotmo(aflip, all1=T, degree1=c(3,1), degree2=NA, caption="all1=T, degree1=c(3,1), degree2=NA")
options(warn=2)
expect.err(try(plotmo(aflip, no.such.arg=9)), "(converted from warning) predict.earth ignored argument 'no.such.arg'")
expect.err(try(plotmo(aflip, ycolumn=1)), "(converted from warning) predict.earth ignored argument 'ycolumn'")
expect.err(try(plotmo(aflip, title="abc")), "(converted from warning) predict.earth ignored argument 'title'")
expect.err(try(plotmo(aflip, persp.ticktype="d", persp.ntick=3, tic=3, tick=9)), "(converted from warning) predict.earth ignored argument 'tic'")
expect.err(try(plotmo(aflip, persp.ticktype="d", ntick=3, tic=3)), "(converted from warning) predict.earth ignored argument 'ntick'")
options(warn=1)
# expect.err(try(plotmo(aflip, adj1=8, adj2=9))) # Error : plotmo: illegal argument "adj1"
# expect.err(try(plotmo(aflip, yc=8, x2=9))) # "ycolumn" is no longer legal, use "nresponse" instead
# expect.err(try(plotmo(aflip, persp.ticktype="d", ntick=3, ti=3))) # Error : "title" is illegal, use "caption" instead ("ti" taken to mean "title")
# expect.err(try(plotmo(aflip, persp.ticktype="d", ntick=3, title=3))) # Error : "title" is illegal, use "caption" instead
# expect.err(try(plotmo(aflip, persp.ticktype="d", ntick=3, tit=3, titl=7))) # Error : "title" is illegal, use "caption" instead ("tit" taken to mean "title")
# expect.err(try(plotmo(aflip, zlab="abc"))) # "zlab" is illegal, use "ylab" instead
# expect.err(try(plotmo(aflip, z="abc"))) # "zlab" is illegal, use "ylab" instead ("z" taken to mean "zlab")
expect.err(try(plotmo(aflip, degree1=c(4,1))), "'degree1' is out of range, allowed values are 1 to 2")
# expect.err(try(plotmo(aflip, none.such=TRUE))) # illegal argument "all1"
# expect.err(try(plotmo(aflip, ntick=3, type2="im"))) # the ntick argument is illegal for type2="image"
# expect.err(try(plotmo(aflip, breaks=3, type2="persp"))) # the breaks argument is illegal for type2="persp"
# expect.err(try(plotmo(aflip, breaks=99, type2="cont"))) # the breaks argument is illegal for type2="contour"
# Test error handling when accessing the original data
lm.bad <- lm.fit(as.matrix(ozone1[,-1]), as.matrix(ozone1[,1]))
expect.err(try(plotmo(lm.bad)), "'lm.bad' is a plain list, not an S3 model")
expect.err(try(plotmo(99)), "'99' is not an S3 model")
x <- matrix(c(1,3,2,4,5,6,7,8,9,10,
2,3,4,5,6,7,8,9,8,9), ncol=2)
colnames(x) <- c("c1", "c2")
x1 <- x[,1]
x2 <- x[,2]
y <- 3:12
df <- data.frame(y=y, x1=x1, x2=x2)
foo1 <- function()
{
a.foo1 <- lm(y~x1+x2, model=FALSE)
x1 <- NULL
expect.err(try(plotmo(a.foo1)), "cannot get the original model predictors")
}
foo1()
foo2 <- function()
{
a.foo2 <- lm(y~x1+x2, data=df, model=FALSE)
df <- 99 # note that df <- NULL here will not cause an error msg
y <- 99 # also needed else model.frame in plotmo will find the global y
expect.err(try(plotmo(a.foo2)), "cannot get the original model predictors")
}
foo2()
foo3 <- function()
{
a.foo3 <- lm(y~x) # lm() builds an lm model for which predict doesn't work
expect.err(try(plotmo(a.foo3)), "predict returned the wrong length (got 10 but expected 50)")
}
foo3()
foo3a <- function()
{
a.foo3a <- lm(y~x) # lm() builds an lm model for which predict doesn't work
# this tests "ngrid1 <- ngrid1 + 1" in plotmo.R
expect.err(try(plotmo(a.foo3a, ngrid1=nrow(x))), "predict returned the wrong length (got 10 but expected 11)")
}
foo3a()
foo4 <- function()
{
a.foo4 <- lm(y~x[,1]+x[,2]) # builds an lm model for which predict doesn't work
# causes 'newdata' had 8 rows but variables found have 10 rows
expect.err(try(plotmo(a.foo4)), "predict returned the wrong length (got 10 but expected 50)")
}
foo4()
foo5 <- function()
{
a.foo5 <- lm(y~x1+x2, model=FALSE)
x1 <- c(1,2,3)
# causes Error in model.frame.default: variable lengths differ (found for 'x1')
expect.err(try(plotmo(a.foo5)), "cannot get the original model predictors")
}
foo5()
foo6 <- function()
{
a.foo6 <- lm(y~x1+x2, model=FALSE)
y[1] <- NA
# Error in na.fail.default: missing values in object
expect.err(try(plotmo(a.foo6, col.response=3)), "cannot get the original model predictors")
}
foo6()
foo7 <- function()
{
a.foo7 <- lm(y~x1+x2, model=FALSE)
y[1] <- Inf
options <- options("warn")
on.exit(options(warn=options$warn))
options(warn=2)
expect.err(try(plotmo(a.foo7, col.response=3)), "non-finite values returned by plotmo_y")
}
foo7()
options(warn=1)
foo8 <- function()
{
i <- 1
a.foo8 <- lm(y~x[,i]+x[,2])
options <- options("warn")
on.exit(options(warn=options$warn))
options(warn=2)
expect.err(try(plotmo(a.foo8)), "Cannot determine which variables to plot in degree2 plots (use all2=TRUE?)")
options(warn=options$warn)
expect.err(try(plotmo(a.foo8)), "predict returned the wrong length (got 10 but expected 50)")
}
foo8()
options(warn=1)
foo9 <- function()
{
my.list <- list(j=2)
a.foo9 <- lm(y~x[,1]+x[,my.list$j])
expect.err(try(plotmo(a.foo9)), "cannot get the original model predictors")
}
foo9()
foo9a <- function()
{
df <- data.frame(y=y, x1=x[,1], x2=x[,2])
a.foo9a <- lm(y~x1+x2, data=df)
par(mfrow = c(2, 2), oma=c(0,0,4,0))
set.seed(2018)
plotmo(a.foo9a, col.resp=2, do.par=FALSE,
caption="top two plots should be identical to bottom two plots")
x2 <- rep(99, length(x2))
a.foo9b <- lm(y~x1+x2, data=df)
x2 <- rep(199, length(x2))
plotmo(a.foo9b, col.resp=2, do.par=FALSE)
}
foo9a()
par(org.par)
foo20.func <- function()
{
par(mfrow = c(2, 2), oma=c(0,0,4,0))
foo20 <- lm(y~x1+x2)
set.seed(2018)
plotmo(foo20, degree1=1:2, col.resp=2, do.par=FALSE,
caption="top two plots should be identical to bottom two plots\nbecause we use saved lm$model")
x1 <- 99
plotmo(foo20, degree1=1:2, col.resp=2, do.par=FALSE)
}
foo20.func()
par(org.par)
set.seed(1235)
tit <- etitanic
tit <- tit[c(30:80,330:380,630:680), ]
a <- earth(survived~., data=tit, glm=list(family=binomial), degree=2)
plotmo(a, grid.levels=list(sex="ma"),
caption="smooth: survived, sex=\"m\" jitter=1",
smooth.col="indianred", smooth.lwd=2,
col.response=as.numeric(tit$survived)+2, pt.pch=".", type2="im",
pt.cex=3, jitter=1) # big jitter
set.seed(1238)
a <- earth(pclass~., data=tit)
plotmo(a, type="class", nresponse=1,
grid.levels=list(sex="ma"),
caption="smooth: pclass, sex=\"m\"", SHOWCALL=TRUE,
smooth.col="indianred", smooth.lwd=2,
col.response=as.numeric(tit$pclass)+1, type2="im",
pt.pch=".", pt.cex=3)
plotmo(a, type="class", nresponse=1,
grid.levels=list(sex="ma"),
caption="smooth: pclass, sex=\"m\" jitter=.3", SHOWCALL=TRUE,
smooth.col="indianred", smooth.lwd=2,
col.response=as.numeric(tit$pclass)+1, type2="im",
pt.pch="x", jit=.3) # small jitter
plotmo(a, nresponse=1,
type="class", grid.levels=list(sex="ma"),
caption="smooth: pclass, sex=\"m\"", SHOWCALL=TRUE,
smooth.col="indianred", smooth.lwd=2,
col.response=as.numeric(tit$pclass)+1, type2="im",
pt.pch=paste(1:nrow(tit)))
# test the extend argument
plotmo(a, nresponse=1, pt.col=2, degree2=0, SHOWCALL=TRUE,
caption="test extend: extend=0 (reference plot)")
plotmo(a, nresponse=1, extend=.5, pt.col=2, SHOWCALL=TRUE,
caption="test extend: extend=.5")
plotmo(a, nresponse=1, degree1=0, extend=.2, pt.col=2, SHOWCALL=TRUE) # nothing to plot
a <- earth(survived~pclass+age, data=etitanic, degree=2)
# expect warning: extend=.5 not degree2 plots
plotmo(a, extend=.5, pt.col=2, SHOWCALL=TRUE,
caption="test extend: extend=.5")
# intercept only models
dopar(2, 2, caption = "intercept-only models")
set.seed(1)
x <- 1:10
y <- runif(length(x))
earth.intercept.only <- earth(x, y)
plotmo(earth.intercept.only, do.par=FALSE, main="earth intercept-only model")
plotmo(earth.intercept.only, do.par=FALSE, col.response=1, pt.pch=20)
# TODO following draws a plot but it shouldn't (very minor bug because int-only model with a bad degree1 spec)
plotmo(earth.intercept.only, do.par=FALSE, degree1=3) # expect warning: 'degree1' specified but no degree1 plots
plotmo(earth.intercept.only, do.par=FALSE, degree1=0) # expect warning: plotmo: nothing to plot
library(rpart)
rpart.intercept.only <- rpart(y~x)
plotmo(rpart.intercept.only, do.par=FALSE, main="rpart.plot intercept-only model")
plotmo(rpart.intercept.only, do.par=FALSE, degree1=0)
par(org.par)
# nrug argument
par(mfrow=c(3,3), mar=c(3,3,3,1), mgp=c(1.5, 0.5, 0))
mod.nrug <- earth(survived~age, data=etitanic)
set.seed(2016)
plotmo(mod.nrug, do.par=0, nrug=-1, main="nrug=-1")
plotmo(mod.nrug, do.par=0, nrug=TRUE, main="nrug=TRUE")
plotmo(mod.nrug, do.par=0, nrug=10, rug.col=2, main="nrug=10, rug.col=2")
plotmo(mod.nrug, do.par=0, nrug=5, rug.col=2, rug.lwd=2, main="nrug=5, rug.col=2, rug.lwd=2")
plotmo(mod.nrug, do.par=0, nrug="density", main="nrug=\"density\"")
plotmo(mod.nrug, do.par=0, nrug="density", density.col=2, density.lwd=2, main="nrug=\"density\"\ndensity.col=2, density.lwd=2")
plotmo(mod.nrug, do.par=0, nrug="density", density.adj=.2, density.col=1, main="nrug=\"density\"\ndensity.adj=.2, density.col=1")
par(org.par)
# a <- earth(ozone1[,3]~ozone1[,1]+ozone1[,2]+ozone1[,4]+ozone1[,5]+ozone1[,6], data=ozone1)
# # TODO fails: actual.nrows=330 expected.nrows=50 fitted.nrows=330
# plotmo(a)
# # TODO following fails in plotmo with
# # Error : get.earth.x from model.matrix.earth from predict.earth: x has 2 columns, expected 4 to match: 1 2 3 Girth
# a <- earth(Volume~poly(Height, degree=3)+Girth, data=trees, subset=4:23, linpreds=TRUE)
# plotmo(a, trace=-1, do.par=FALSE, caption="all three rows should be the same")
source("test.epilog.R")
Try the plotmo package in your browser
Any scripts or data that you put into this service are public.
plotmo documentation built on May 22, 2022, 1:05 a.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
# test.plotmo.R: regression tests for plotmo
# Stephen Milborrow, Petaluma Jan 2007
print(R.version.string)
source("test.prolog.R")
library(earth)
options(warn=1) # print warnings as they occur
data(etitanic)
make.space.for.caption <- function(caption="CAPTION")
{
oma <- par("oma")
needed <- 3
# adjust for newlines in caption
newlines <- grep("\n", caption)
if(length(newlines) > 0)
needed <- needed + .5 * newlines # .5 seems enough although 1 in theory
if(!is.null(caption) && any(nchar(caption)) && oma[3] <= needed) {
oma[3] <- needed
par(oma=oma)
}
}
dopar <- function(nrows, ncols, caption = "")
{
cat(" ", caption, "\n")
make.space.for.caption(caption)
par(mfrow=c(nrows, ncols))
par(mar = c(3, 3, 1.7, 0.5))
par(mgp = c(1.6, 0.6, 0))
par(cex = 0.7)
}
example(plotmo)
caption <- "basic earth test of plotmo"
a <- earth(O3 ~ ., data=ozone1, degree=2)
plotmo(a, degree1=2, degree2=4, caption=caption, trace=-1)
caption <- "test 5 x 5 layout"
dopar(1,1,caption)
a <- earth(O3 ~ ., data=ozone1, nk=51, pmethod="n", degree=2)
plotmo(a, caption=caption, trace=1)
caption <- "test 4 x 4 layout with ylab"
dopar(1,1,caption)
a <- earth(O3 ~ ., data=ozone1, nk=30, pmethod="n", degree=2)
plotmo(a, caption=caption, trace=2)
caption <- "test 3 x 3 layout"
dopar(1,1,caption)
a <- earth(O3 ~ ., data=ozone1, nk=16, pmethod="n", degree=2)
plotmo(a, caption=caption, trace=3)
caption <- "test 2 x 2 layout"
dopar(1,1,caption)
a <- earth(O3 ~ ., data=ozone1, nk=9, pmethod="n", degree=2)
plotmo(a, caption=caption)
caption <- "test 1 x 1 layout"
dopar(1,1,caption)
a <- earth(O3 ~ ., data=ozone1, nk=4, pmethod="n", degree=2)
plotmo(a, caption=caption)
caption <- "test plotmo basic params"
a <- earth(O3 ~ ., data=ozone1, degree=2)
dopar(3,2,caption)
plotmo(a, do.par=FALSE, degree1=1, nrug=-1, degree2=F, caption=caption,
main="test main", xlab="test xlab", ylab="test ylab")
plotmo(a, do.par=FALSE, degree1=F, degree2=4, grid.func=mean, persp.col="white", ngrid2=10, persp.phi=40)
set.seed(2016)
plotmo(a, do.par=FALSE, degree1=1, degree1.lty=2, degree1.lwd=4, degree1.col=2, nrug=TRUE, degree2=F, main="nrug=300")
plotmo(a, do.par=FALSE, degree1=1, nrug=-1, degree2=F, main="nrug=TRUE")
set.seed(2016)
plotmo(a, do.par=FALSE, degree1=1, nrug=10, ngrid1=50, degree2=F, main="ngrid1=50 nrug=10")
plotmo(a, do.par=FALSE, degree1=NA, degree2=1, persp.phi=60) # graph args
caption <- "test plotmo xlim and ylim"
a <- earth(O3 ~ ., data=ozone1, degree=2)
dopar(5,3,caption)
plotmo(a, do.par=FALSE, degree1=2:3, degree2=4, caption=caption, xlab="ylim=default")
plotmo(a, do.par=FALSE, degree1=2:3, degree2=4, ylim=NA, xlab="ylim=NA")
plotmo(a, do.par=FALSE, degree1=2:3, degree2=4, ylim=c(0,20), xlab="ylim=c(0,20)")
plotmo(a, do.par=FALSE, degree1=2:3, degree2=4, xlim=c(190,250), xlab="xlim=c(190,250)")
plotmo(a, do.par=FALSE, degree1=2:3, degree2=4, xlim=c(190,250), ylim=c(11,18), xlab="xlim=c(190,250), ylim=c(11,18)")
# check various types of predictors with grid.func and ndiscrete
varied.type.data <- data.frame(
y = 1:13,
num = c(1, 3, 2, 3, 4, 5, 6, 4, 5, 6.5, 3, 6, 5), # 7 unique values (but one is non integral)
int = c(1L, 1L, 3L, 3L, 4L, 4L, 3L, 5L, 3L, 6L, 7L, 8L, 10L), # 8 unique values
bool = c(F, F, F, F, F, T, T, T, T, T, T, T, T),
date = as.Date(
c("2018-08-01", "2018-08-02", "2018-08-03",
"2018-08-04", "2018-08-05", "2018-08-06",
"2018-08-07", "2018-08-08", "2018-08-08",
"2018-08-08", "2018-08-10", "2018-08-11", "2018-08-11")),
ord = ordered(c("ord3", "ord3", "ord3",
"ord1", "ord2", "ord3",
"ord1", "ord2", "ord3",
"ord1", "ord1", "ord1", "ord1"),
levels=c("ord1", "ord3", "ord2")),
fac = as.factor(c("fac1", "fac1", "fac1",
"fac2", "fac2", "fac2",
"fac3", "fac3", "fac3",
"fac1", "fac2", "fac3", "fac3")),
str = c("str1", "str1", "str1", # will be treated like a factor
"str2", "str2", "str2",
"str3", "str3", "str3",
"str3", "str3", "str3", "str3"))
varied.type.lm <- lm(y ~ ., data = varied.type.data)
print(summary(varied.type.lm))
set.seed(2018)
plotres(varied.type.lm, info=TRUE)
plotmo(varied.type.lm, pmethod="apartdep", all2=TRUE, ticktype="d", col.response="red", caption="varied.type.lm\npmethod=\"apartdep\" default grid func")
plotmo(varied.type.lm, all2=TRUE, ticktype="d", col.response="red", caption="varied.type.lm\ndefault grid func")
plotmo(varied.type.lm, all2=TRUE, ndiscre=1, caption="varied.type.lm\nndiscrete=1")
plotmo(varied.type.lm, all2=TRUE, ndiscr=2, caption="varied.type.lm\nndiscrete=2")
plotmo(varied.type.lm, all2=TRUE, ndis=100, caption="varied.type.lm\nndiscrete=100")
cat("grid.func=median:\n")
plotmo(varied.type.lm, all2=TRUE, grid.func=median, caption="varied.type.lm\ngrid.func=median")
cat("grid.func=quantile:\n")
plotmo(varied.type.lm, all2=TRUE, grid.func=function(x, ...) quantile(x, 0.5), caption="varied.type.lm\ngrid.func=function(x, ...) quantile(x, 0.5)")
cat("grid.func=mean:\n")
plotmo(varied.type.lm, all2=TRUE, grid.func=mean, caption="varied.type.lm\ngrid.func=mean")
varied.type.earth <- earth(y ~ ., data = varied.type.data, thresh=0, penalty=-1, trace=1)
print(summary(varied.type.earth))
set.seed(2018)
plotres(varied.type.earth, info=TRUE)
plotmo(varied.type.earth, all1=TRUE, all2=TRUE, persp.ticktype="d", col.response="red")
# term.plot calls predict.earth with an se parameter, even with termplot(se=FALSE)
caption <- "basic earth test against termplot"
dopar(4,4,caption)
make.space.for.caption("test caption1")
a <- earth(O3 ~ ., data=ozone1, degree=2)
plotmo(a, do.par=FALSE, ylim=NA, caption=caption, degree2=FALSE)
cat("Ignore warning: predict.earth ignored argument \"se.fit\"\n")
termplot(a)
caption <- "test change order of earth predictors and cex"
dopar(4,4,caption)
# minspan=1 to force two degree2 graphs for the test (wasn't necessary in old versions of earth)
a <- earth(doy ~ humidity + temp + wind, data=ozone1, degree=2, minspan=1)
plotmo(a, do.par=FALSE, ylim=NA, caption=caption, degree2=c(1,2), cex=1.2)
termplot(a)
caption <- "test all1=TRUE"
a <- earth(doy ~ humidity + temp + wind, data=ozone1, degree=2)
plotmo(a, caption=caption, all1=TRUE, persp.ticktype="d", persp.nticks=2)
caption <- "test all2=TRUE"
print(summary(a))
plotmo(a, caption=caption, all2=TRUE)
oz <- ozone1[150:200,c("O3","temp","humidity","ibh")]
a.glob <- earth(O3~temp+humidity, data=oz, degree=2)
ad.glob <- earth(oz[,2:3], oz[,1], degree=2)
func1 <- function()
{
caption <- "test environments and finding the correct data"
dopar(4,4,caption)
set.seed(2016)
plotmo(a.glob, do.par=FALSE, main="a.glob oz",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20, trace=2)
mtext(caption, outer=TRUE, font=2, line=1.5, cex=1)
plotmo(ad.glob, do.par=FALSE, main="ad.glob oz",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pch.response=20, trace=2) # pch.response test backcompat
a <- earth(O3~temp+humidity, data=oz, degree=2)
plotmo(a, do.par=FALSE, main="a oz",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
ad <- earth(oz[,2:3], oz[,1], degree=2)
plotmo(ad, do.par=FALSE, main="ad oz",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
oz.org <- oz
oz10 <- 10 * oz # multiply by 10 so we can see by the axis labels if right data is being used
oz <- oz10 # oz is now local to this function, but multiplied by 10
a.oz10 <- earth(O3~temp+humidity, data=oz, degree=2)
a.oz10.keep <- earth(O3~temp+humidity, data=oz, degree=2, keepxy=TRUE)
plotmo(a.oz10, do.par=FALSE, main="a oz10",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
ad.oz10 <- earth(oz[,2:3], oz[,1], degree=2)
ad.oz10.keep <- earth(oz[,2:3], oz[,1], degree=2, keepxy=TRUE)
plotmo(ad.oz10, do.par=FALSE, main="ad oz10",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
func2 <- function() {
a.func <- earth(O3 ~ temp + humidity, data=oz10, degree=2)
plotmo(a.func, do.par=FALSE, main="a.func oz10",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
ad.func <- earth(oz10[,2:3], oz10[,1], degree=2)
plotmo(ad.func, do.par=FALSE, main="ad.func oz10",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
caption <- "test environments and finding the correct data, continued"
dopar(4,4,caption)
oz <- .1 * oz.org
a.func <- earth(O3~temp+ humidity , data=oz, degree=2)
plotmo(a.func, do.par=FALSE, main="a.func oz.1",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
ad.func <- earth(oz[,2:3], oz[,1], degree=2)
plotmo(ad.func, do.par=FALSE, main="ad.func oz.1",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
plotmo(a.oz10.keep, do.par=FALSE, main="func1:a.oz10.keep",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
plotmo(ad.oz10.keep, do.par=FALSE, main="func1:ad.oz10.keep",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20)
cat("Expect error msg: formal argument \"do.par\" matched by multiple actual arguments\n")
expect.err(try(plotmo(a.oz10, do.par=FALSE, main="func1:a.oz10",
degree1=1, all2=1, degree2=1, type2="im",
col.response=3, pt.pch=20, do.par=FALSE)))
}
func2()
y <- 3:11
x1 <- c(1,3,2,4,5,6,6,6,6)
x2 <- c(2,3,4,5,6,7,8,9,10)
frame <- data.frame(y=y, x1=x1, x2=x2)
foo <- function()
{
lm.18.out <- lm(y~x1+x2, model=FALSE)
x1[2] <- 18
y[3] <- 19
frame <- data.frame(y=y, x1=x1, x2=x2)
list(lm.18.out = lm.18.out,
lm.18 = lm(y~x1+x2),
lm.18.keep = lm(y~x1+x2, x=TRUE, y=TRUE),
lm.18.frame = lm(y~x1+x2, data=frame))
}
temp <- foo()
lm.18.out <- temp$lm.18.out
lm.18 <- temp$lm.18
lm.18.keep <- temp$lm.18.keep
lm.18.frame <- temp$lm.18.frame
# following should all use the x1 and y inside foo
cat("==lm.18.out\n")
plotmo(lm.18.out, main="lm.18.out",
do.par=FALSE, degree1=1, clip=FALSE, ylim=c(0,20),
col.response=2, pt.pch=20)
cat("==lm.18\n")
plotmo(lm.18, main="lm.18",
do.par=FALSE, degree1=1, clip=FALSE, ylim=c(0,20),
col.response=2, pt.pch=20)
cat("==lm.18.keep\n")
plotmo(lm.18.keep, main="lm.18.keep", trace=2,
do.par=FALSE, degree1=1, clip=FALSE, ylim=c(0,20),
col.response=2, pt.pch=20)
cat("==lm.18.frame\n")
plotmo(lm.18.frame, main="lm.18.frame",
do.par=FALSE, degree1=1, clip=FALSE, ylim=c(0,20),
col.response=2, pt.pch=20)
}
func1()
caption <- "test earth formula versus x,y model"
# dopar(4,4,caption)
# mtext(caption, outer=TRUE, font=2, line=1.5, cex=1)
a <- earth(O3 ~ ., data=ozone1, degree=2)
plotmo(a, caption="test earth formula versus xy model (formula)")
a <- earth(ozone1[, -1], ozone1[,1], degree=2)
plotmo(a, caption="test earth formula versus xy model (xy)")
# single predictor
caption <- "test earth(O3~wind, data=ozone1, degree=2), single predictor"
dopar(2,2,caption)
a <- earth(O3~wind, data=ozone1, degree=2)
plotmo(a)
caption = "se=2, earth(doy~humidity+temp+wind, data=ozone1) versus termplot (expect no se lines)"
dopar(3,3,caption)
mtext(caption, outer=TRUE, font=2, line=1.5, cex=1)
# minspan=1 to force two degree2 graphs for the test (wasn't necessary in old versions of earth)
a <- earth(doy~humidity + temp + wind, data=ozone1, degree=2, minspan=1)
cat("Ignore warning: predict.earth ignored argument \"se\"\n")
termplot(a)
plotmo(a, do.par=FALSE, ylim=NA, degree2=c(1:2), clip=FALSE, caption=caption)
# test fix to bug reported by Joe Retzer, FIXED Dec 7, 2007
N <- 650
set.seed(2007)
q_4 <- runif(N, -1, 1)
q_2102 <- runif(N, -1, 1)
q_2104 <- runif(N, -1, 1)
q_3105 <- runif(N, -1, 1)
q_3106 <- runif(N, -1, 1)
q_4104 <- runif(N, -1, 1)
q_6101 <- runif(N, -1, 1)
q_6103 <- runif(N, -1, 1)
q_7104 <- runif(N, -1, 1)
q_3109 <- runif(N, -1, 1)
q_4103 <- runif(N, -1, 1)
q_2111 <- runif(N, -1, 1)
q_3107 <- runif(N, -1, 1)
q_3101 <- runif(N, -1, 1)
q_3104 <- runif(N, -1, 1)
q_7107 <- runif(N, -1, 1)
depIndex <- sin(1.0 * q_4 + rnorm(650, sd=.8)) + sin(1.8 * q_2102 + rnorm(650, sd=.8)) + sin(1.3 * q_2104 + rnorm(650, sd=.8)) + sin(1.4 * q_3105 + rnorm(650, sd=.8)) +
sin(1.5 * q_3106 + rnorm(650, sd=.8)) + sin(1.6 * q_4104 + rnorm(650, sd=.8)) + sin(1.8 * q_6101 + rnorm(650, sd=.8)) + sin(1.8 * q_6103 + rnorm(650, sd=.8)) +
sin(1.9 * q_7104 + rnorm(650, sd=.8)) + sin(2.0 * q_3109 + rnorm(650, sd=.8))
regDatCWD <- as.data.frame(cbind(depIndex, q_4, q_2102, q_2104, q_3105, q_3106, q_4104, q_6101, q_6103, q_7104, q_3109, q_4103, q_2111, q_3107, q_3101, q_3104, q_7107))
cat("--plotmo(earthobj5)--\n")
earthobj5 <- earth(depIndex ~ q_4+q_2102+q_2104+q_3105+q_3106+q_4104+q_6101+q_6103+q_7104+q_3109+q_4103+q_2111+q_3107+q_3101+q_3104+q_7107, data=regDatCWD)
print(summary(earthobj5, digits = 2))
plotmo(earthobj5)
# long predictor names
a.rather.long.in.fact.very.long.name.q_4 <- q_4
a.rather.long.in.fact.very.long.name.q_2102 <- q_2102
a.rather.long.in.fact.very.long.name.q_2104 <- q_2104
a.rather.long.in.fact.very.long.name.q_3105 <- q_3105
a.rather.long.in.fact.very.long.name.q_3106 <- q_3106
a.rather.long.in.fact.very.long.name.q_4104 <- q_4104
a.rather.long.in.fact.very.long.name.q_6101 <- q_6101
a.rather.long.in.fact.very.long.name.q_6103 <- q_6103
a.rather.long.in.fact.very.long.name.q_7104 <- q_7104
a.rather.long.in.fact.very.long.name.q_3109 <- q_3109
a.rather.long.in.fact.very.long.name.q_4103 <- q_4103
a.rather.long.in.fact.very.long.name.q_2111 <- q_2111
a.rather.long.in.fact.very.long.name.q_3107 <- q_3107
a.rather.long.in.fact.very.long.name.q_3101 <- q_3101
a.rather.long.in.fact.very.long.name.q_3104 <- q_3104
a.rather.long.in.fact.very.long.name.q_7107 <- q_7107
a.rather.long.in.fact.very.long.name.for.the.response <- depIndex
a.rather.long.in.fact.very.long.name.for.the.dataframe <-
as.data.frame(cbind(
a.rather.long.in.fact.very.long.name.for.the.response,
a.rather.long.in.fact.very.long.name.q_4,
a.rather.long.in.fact.very.long.name.q_2102,
a.rather.long.in.fact.very.long.name.q_2104,
a.rather.long.in.fact.very.long.name.q_3105,
a.rather.long.in.fact.very.long.name.q_3106,
a.rather.long.in.fact.very.long.name.q_4104,
a.rather.long.in.fact.very.long.name.q_6101,
a.rather.long.in.fact.very.long.name.q_6103,
a.rather.long.in.fact.very.long.name.q_7104,
a.rather.long.in.fact.very.long.name.q_3109,
a.rather.long.in.fact.very.long.name.q_4103,
a.rather.long.in.fact.very.long.name.q_2111,
a.rather.long.in.fact.very.long.name.q_3107,
a.rather.long.in.fact.very.long.name.q_3101,
a.rather.long.in.fact.very.long.name.q_3104,
a.rather.long.in.fact.very.long.name.q_7107))
cat("--a.rather.long.in.fact.very.long.name.for.the...A--\n")
a.rather.long.in.fact.very.long.name.for.the.modelA <-
earth(a.rather.long.in.fact.very.long.name.for.the.response ~
a.rather.long.in.fact.very.long.name.q_4 +
a.rather.long.in.fact.very.long.name.q_2102 +
a.rather.long.in.fact.very.long.name.q_2104 +
a.rather.long.in.fact.very.long.name.q_3105 +
a.rather.long.in.fact.very.long.name.q_3106 +
a.rather.long.in.fact.very.long.name.q_4104 +
a.rather.long.in.fact.very.long.name.q_6101 +
a.rather.long.in.fact.very.long.name.q_6103 +
a.rather.long.in.fact.very.long.name.q_7104 +
a.rather.long.in.fact.very.long.name.q_3109 +
a.rather.long.in.fact.very.long.name.q_4103 +
a.rather.long.in.fact.very.long.name.q_2111 +
a.rather.long.in.fact.very.long.name.q_3107 +
a.rather.long.in.fact.very.long.name.q_3101 +
a.rather.long.in.fact.very.long.name.q_3104 +
a.rather.long.in.fact.very.long.name.q_7107,
data = a.rather.long.in.fact.very.long.name.for.the.dataframe)
print(summary(a.rather.long.in.fact.very.long.name.for.the.modelA, digits = 2))
plot(a.rather.long.in.fact.very.long.name.for.the.modelA)
plotmo(a.rather.long.in.fact.very.long.name.for.the.modelA)
cat("--a.rather.long.in.fact.very.long.name.for.the...C--\n")
a.rather.long.in.fact.very.long.name.for.the.modelC <-
earth(x = a.rather.long.in.fact.very.long.name.for.the.dataframe[,-1],
y = a.rather.long.in.fact.very.long.name.for.the.response,
degree = 3)
print(summary(a.rather.long.in.fact.very.long.name.for.the.modelC, digits = 2))
plot(a.rather.long.in.fact.very.long.name.for.the.modelC)
plotmo(a.rather.long.in.fact.very.long.name.for.the.modelC)
a <- earth(survived ~ pclass+sex+age, data=etitanic, degree=2)
print(summary(a))
plotmo(a, caption="plotmo with facs: pclass+sex+age")
plotmo(a, caption="plotmo with facs: pclass+sex+age, all1=T, grid.col=\"gray\"", all1=T, grid.col="gray")
plotmo(a, caption="plotmo with facs: pclass+sex+age, all2=T, col.grid=\"green\"", all2=T, col.grid="green")
plotmo(a, caption="plotmo with facs: pclass+sex+age, all1=T, all2=T, grid=2", all1=T, all2=T, grid.col=2)
plotmo(a, clip=FALSE, degree2=FALSE, caption="plotmo (no degree2) with facs: pclass+sex+age")
plotmo(a, clip=FALSE, grid.levels=list(pclass="2n", sex="ma"),
caption="plotmo with grid.levels: pclass+sex+age")
# in above tests, all degree2 terms use facs
# now build a model with some degree2 term that use facs, some that don't
a <- earth(survived ~ pclass+age+sibsp, data=etitanic, degree=2)
print(summary(a))
plotmo(a, caption="plotmo with mixed fac and non-fac degree2 terms", persp.border=NA)
plotmo(a, caption="plotmo with mixed fac and non-fac degree2 terms and grid.levels",
grid.levels=list(pclass="2n", age=20), # test partial matching of grid levels, and numeric preds
persp.ticktype="d", persp.nticks=2)
# check detection of illegal grid.levels argument
expect.err(try(plotmo(a, grid.levels=list(pcla="1", pclass="2"))), 'illegal grid.levels argument ("pcla" and "pclass" both match "pclass")')
expect.err(try(plotmo(a, grid.levels=list(pclass="1", pcla="2"))), 'illegal grid.levels argument ("pclass" and "pcla" both match "pclass")')
expect.err(try(plotmo(a, grid.levels=list(pcla="nonesuch"))), 'illegal level "nonesuch" for "pclass" in grid.levels (allowed levels are "1st" "2nd" "3rd")')
expect.err(try(plotmo(a, grid.levels=list(pcla="1sx"))), 'illegal level "1sx" for "pclass" in grid.levels (allowed levels are "1st" "2nd" "3rd")')
expect.err(try(plotmo(a, grid.levels=list(pcla=1))), 'illegal level for "pclass" in grid.levels (specify factor levels with a string)')
expect.err(try(plotmo(a, grid.levels=list(pcla=c("ab", "cd")))), "length(pclass) in grid.levels is not 1")
expect.err(try(plotmo(a, grid.levels=list(pcla=NA))), 'pclass in grid.levels is NA')
expect.err(try(plotmo(a, grid.levels=list(pcla=Inf))), 'pclass in grid.levels is infinite')
expect.err(try(plotmo(a, grid.levels=list(pcla=9))), 'illegal level for "pclass" in grid.levels (specify factor levels with a string)')
options(warn=2)
expect.err(try(plotmo(a, grid.levels=list(age="ab"))), 'grid.levels returned class \"character\" for age, so will use the default grid.func for age')
options(warn=1)
expect.err(try(plotmo(a, grid.levels=list(age=NA))), 'age in grid.levels is NA')
expect.err(try(plotmo(a, grid.levels=list(age=Inf))), 'age in grid.levels is infinite')
expect.err(try(plotmo(a, grid.lev=list(age=list(1,2)))), 'length(age) in grid.levels is not 1')
# more-or-less repeat above, but with glm models
a <- earth(survived ~ pclass+age+sibsp, data=etitanic, degree=2, glm=list(family=binomial))
print(summary(a))
plotmo(a, ylim=c(0, 1), caption="plotmo glm with mixed fac and non-fac degree2 terms")
plotmo(a, ylim=c(0, 1), caption="plotmo glm with mixed fac and non-fac degree2 terms and grid.levels",
grid.levels=list(pcl="2nd")) # test partial matching of variable name in grid levels
plotmo(a, type="earth", ylim=c(0, 1), caption="type=\"earth\" plotmo glm with mixed fac and non-fac degree2 terms")
plotmo(a, type="link", ylim=c(0, 1), clip=FALSE, caption="type=\"link\" plotmo glm with mixed fac and non-fac degree2 terms")
plotmo(a, type="class", ylim=c(0, 1), caption="type=\"class\" plotmo glm with mixed fac and non-fac degree2 terms")
plotmo(a, ylim=c(0, 1), caption="default type (\"response\")\nplotmo glm with mixed fac and non-fac degree2 terms")
# now with different type2s
set.seed(2016)
plotmo(a, do.par=FALSE, type2="persp", persp.theta=-20, degree1=FALSE, grid.levels=list(pclass="2nd"))
mtext("different type2s", outer=TRUE, font=2, line=1.5, cex=1)
plotmo(a, do.par=FALSE, type2="contour", degree1=FALSE, grid.levels=list(pclass="2nd"))
plotmo(a, do.par=FALSE, type2="image", degree1=FALSE, grid.levels=list(pclass="2nd"),
col.response=as.numeric(etitanic$survived)+2, pt.pch=20)
plotmo(a, do.par=FALSE, type="earth", type2="image", degree1=FALSE,
grid.levels=list(pclass="2"))
# grid.levels with partdep
set.seed(2018)
x1 <- (1:11) + runif(11)
x2 <- (1:11) + runif(11)
x3 <- as.integer((1:11) + runif(11))
x4 <- runif(11) > .5 # logical
y <- x1 - x2 + x3 + x4
data <- data.frame(y=y, x1=x1, x2=x2, x3=x3, x4=x4)
lm.x1.x2.x3 <- lm(y ~ x1 + x2 + x3 + x4 + x1*x2 + x1*x3, data=data)
cat("summary(lm.x1.x2.x3):\n")
print(summary(lm.x1.x2.x3))
par(mfrow = c(5, 6), mar = c(2, 3, 2, 1), mgp = c(1.5, 0.5, 0), cex = 0.6, oma=c(0,0,8,0))
plotmo(lm.x1.x2.x3, do.par=0, ylim=c(0,16), pt.col=2,
caption="row1 default\nrow2 grid.levels=list(x3=15)\nrow3 partdep\nrow4 partdetp grid.levels=list(x3=15)")
plotmo(lm.x1.x2.x3, do.par=0, ylim=c(0,16), pt.col=2, grid.levels=list(x3=15))
plotmo(lm.x1.x2.x3, do.par=0, ylim=c(0,16), pt.col=2, pmethod="partdep")
plotmo(lm.x1.x2.x3, do.par=0, ylim=c(0,16), pt.col=2, pmethod="partdep", grid.levels=list(x3=15))
# check auto type convert in grid.levels
plotmo(lm.x1.x2.x3, degree1="x1", degree2=0, main="x1 (x2=5L))", ylim=c(0,16), do.par=0, pmethod="partdep", grid.levels=list(x2=15L)) # integer to numeric
plotmo(lm.x1.x2.x3, degree1="x1", degree2=0, main="x1 (x3=5))", ylim=c(0,16), do.par=0, pmethod="partdep", grid.levels=list(x3=15)) # numeric to integer
plotmo(lm.x1.x2.x3, degree1="x1", degree2=0, main="x1 (x4=1))", ylim=c(0,16), do.par=0, pmethod="partdep", grid.levels=list(x4=1)) # numeric to logical
expect.err(try(plotmo(lm.x1.x2.x3, degree1="x1", degree2=0, main="x1 (x4=1))", ylim=c(0,16), do.par=0, pmethod="partdep", grid.levels=list(x4="x"))), "expected a logical value in grid.levels for x4") # char to logical
expect.err(try(plotmo(lm.x1.x2.x3, degree1="x2", do.par=0, pmethod="partdep", grid.levels=list(x1="1"))), "the class \"character\" of \"x1\" in grid.levels does not match its class \"numeric\" in the input data")
par(org.par)
# test vector main
a20 <- earth(O3 ~ humidity + temp + doy, data=ozone1, degree=2, glm=list(family=Gamma))
dopar(2, 2)
plotmo(a20, nrug=-1)
set.seed(2016)
plotmo(a20, nrug=10, caption="Test plotmo with a vector main (and npoints=200)",
main=c("Humidity", "Temperature", "Day of year", "Humidity: Temperature", "Temperature: Day of Year"),
col.response="darkgray", pt.pch=".", cex.response=3, npoints=200) # cex.response tests back compat
cat("Expect warning below (missing double titles)\n")
plotmo(a20, nrug=-1, caption="Test plotmo with a vector main (and plain smooth)",
main=c("Humidity", "Temperature", "Day of year", "Humidity: Temperature", "Temp: Doy"),
smooth.col="indianred")
cat("Expect warning below (missing single titles)\n")
plotmo(a20, nrug=-1, caption="Test plotmo with a vector main (and smooth args)",
main=c("Humidity", "Temperature"),
smooth.col="indianred", smooth.lwd=2, smooth.lty=2, smooth.f=.1,
col.response="gray", npoints=500)
plotmo(a20, nrug=-1, caption="Test plotmo with pt.pch=paste(1:nrow(ozone1))",
type2="im",
col.response=2, pt.cex=.8, pt.pch=paste(1:nrow(ozone1)), npoints=100)
aflip <- earth(O3~vh + wind + humidity + temp, data=ozone1, degree=2)
# test all1 and all2, with and without degree1 and degree2
plotmo(aflip, all2=T, caption="all2=T", npoints=TRUE)
plotmo(aflip, all2=T, degree2=c(4, 2), caption="all2=T, degree2=c(4, 2)")
plotmo(aflip, all1=T, caption="all1=T")
plotmo(aflip, all1=T, degree1=c(3,1), degree2=NA, caption="all1=T, degree1=c(3,1), degree2=NA")
options(warn=2)
expect.err(try(plotmo(aflip, no.such.arg=9)), "(converted from warning) predict.earth ignored argument 'no.such.arg'")
expect.err(try(plotmo(aflip, ycolumn=1)), "(converted from warning) predict.earth ignored argument 'ycolumn'")
expect.err(try(plotmo(aflip, title="abc")), "(converted from warning) predict.earth ignored argument 'title'")
expect.err(try(plotmo(aflip, persp.ticktype="d", persp.ntick=3, tic=3, tick=9)), "(converted from warning) predict.earth ignored argument 'tic'")
expect.err(try(plotmo(aflip, persp.ticktype="d", ntick=3, tic=3)), "(converted from warning) predict.earth ignored argument 'ntick'")
options(warn=1)
# expect.err(try(plotmo(aflip, adj1=8, adj2=9))) # Error : plotmo: illegal argument "adj1"
# expect.err(try(plotmo(aflip, yc=8, x2=9))) # "ycolumn" is no longer legal, use "nresponse" instead
# expect.err(try(plotmo(aflip, persp.ticktype="d", ntick=3, ti=3))) # Error : "title" is illegal, use "caption" instead ("ti" taken to mean "title")
# expect.err(try(plotmo(aflip, persp.ticktype="d", ntick=3, title=3))) # Error : "title" is illegal, use "caption" instead
# expect.err(try(plotmo(aflip, persp.ticktype="d", ntick=3, tit=3, titl=7))) # Error : "title" is illegal, use "caption" instead ("tit" taken to mean "title")
# expect.err(try(plotmo(aflip, zlab="abc"))) # "zlab" is illegal, use "ylab" instead
# expect.err(try(plotmo(aflip, z="abc"))) # "zlab" is illegal, use "ylab" instead ("z" taken to mean "zlab")
expect.err(try(plotmo(aflip, degree1=c(4,1))), "'degree1' is out of range, allowed values are 1 to 2")
# expect.err(try(plotmo(aflip, none.such=TRUE))) # illegal argument "all1"
# expect.err(try(plotmo(aflip, ntick=3, type2="im"))) # the ntick argument is illegal for type2="image"
# expect.err(try(plotmo(aflip, breaks=3, type2="persp"))) # the breaks argument is illegal for type2="persp"
# expect.err(try(plotmo(aflip, breaks=99, type2="cont"))) # the breaks argument is illegal for type2="contour"
# Test error handling when accessing the original data
lm.bad <- lm.fit(as.matrix(ozone1[,-1]), as.matrix(ozone1[,1]))
expect.err(try(plotmo(lm.bad)), "'lm.bad' is a plain list, not an S3 model")
expect.err(try(plotmo(99)), "'99' is not an S3 model")
x <- matrix(c(1,3,2,4,5,6,7,8,9,10,
2,3,4,5,6,7,8,9,8,9), ncol=2)
colnames(x) <- c("c1", "c2")
x1 <- x[,1]
x2 <- x[,2]
y <- 3:12
df <- data.frame(y=y, x1=x1, x2=x2)
foo1 <- function()
{
a.foo1 <- lm(y~x1+x2, model=FALSE)
x1 <- NULL
expect.err(try(plotmo(a.foo1)), "cannot get the original model predictors")
}
foo1()
foo2 <- function()
{
a.foo2 <- lm(y~x1+x2, data=df, model=FALSE)
df <- 99 # note that df <- NULL here will not cause an error msg
y <- 99 # also needed else model.frame in plotmo will find the global y
expect.err(try(plotmo(a.foo2)), "cannot get the original model predictors")
}
foo2()
foo3 <- function()
{
a.foo3 <- lm(y~x) # lm() builds an lm model for which predict doesn't work
expect.err(try(plotmo(a.foo3)), "predict returned the wrong length (got 10 but expected 50)")
}
foo3()
foo3a <- function()
{
a.foo3a <- lm(y~x) # lm() builds an lm model for which predict doesn't work
# this tests "ngrid1 <- ngrid1 + 1" in plotmo.R
expect.err(try(plotmo(a.foo3a, ngrid1=nrow(x))), "predict returned the wrong length (got 10 but expected 11)")
}
foo3a()
foo4 <- function()
{
a.foo4 <- lm(y~x[,1]+x[,2]) # builds an lm model for which predict doesn't work
# causes 'newdata' had 8 rows but variables found have 10 rows
expect.err(try(plotmo(a.foo4)), "predict returned the wrong length (got 10 but expected 50)")
}
foo4()
foo5 <- function()
{
a.foo5 <- lm(y~x1+x2, model=FALSE)
x1 <- c(1,2,3)
# causes Error in model.frame.default: variable lengths differ (found for 'x1')
expect.err(try(plotmo(a.foo5)), "cannot get the original model predictors")
}
foo5()
foo6 <- function()
{
a.foo6 <- lm(y~x1+x2, model=FALSE)
y[1] <- NA
# Error in na.fail.default: missing values in object
expect.err(try(plotmo(a.foo6, col.response=3)), "cannot get the original model predictors")
}
foo6()
foo7 <- function()
{
a.foo7 <- lm(y~x1+x2, model=FALSE)
y[1] <- Inf
options <- options("warn")
on.exit(options(warn=options$warn))
options(warn=2)
expect.err(try(plotmo(a.foo7, col.response=3)), "non-finite values returned by plotmo_y")
}
foo7()
options(warn=1)
foo8 <- function()
{
i <- 1
a.foo8 <- lm(y~x[,i]+x[,2])
options <- options("warn")
on.exit(options(warn=options$warn))
options(warn=2)
expect.err(try(plotmo(a.foo8)), "Cannot determine which variables to plot in degree2 plots (use all2=TRUE?)")
options(warn=options$warn)
expect.err(try(plotmo(a.foo8)), "predict returned the wrong length (got 10 but expected 50)")
}
foo8()
options(warn=1)
foo9 <- function()
{
my.list <- list(j=2)
a.foo9 <- lm(y~x[,1]+x[,my.list$j])
expect.err(try(plotmo(a.foo9)), "cannot get the original model predictors")
}
foo9()
foo9a <- function()
{
df <- data.frame(y=y, x1=x[,1], x2=x[,2])
a.foo9a <- lm(y~x1+x2, data=df)
par(mfrow = c(2, 2), oma=c(0,0,4,0))
set.seed(2018)
plotmo(a.foo9a, col.resp=2, do.par=FALSE,
caption="top two plots should be identical to bottom two plots")
x2 <- rep(99, length(x2))
a.foo9b <- lm(y~x1+x2, data=df)
x2 <- rep(199, length(x2))
plotmo(a.foo9b, col.resp=2, do.par=FALSE)
}
foo9a()
par(org.par)
foo20.func <- function()
{
par(mfrow = c(2, 2), oma=c(0,0,4,0))
foo20 <- lm(y~x1+x2)
set.seed(2018)
plotmo(foo20, degree1=1:2, col.resp=2, do.par=FALSE,
caption="top two plots should be identical to bottom two plots\nbecause we use saved lm$model")
x1 <- 99
plotmo(foo20, degree1=1:2, col.resp=2, do.par=FALSE)
}
foo20.func()
par(org.par)
set.seed(1235)
tit <- etitanic
tit <- tit[c(30:80,330:380,630:680), ]
a <- earth(survived~., data=tit, glm=list(family=binomial), degree=2)
plotmo(a, grid.levels=list(sex="ma"),
caption="smooth: survived, sex=\"m\" jitter=1",
smooth.col="indianred", smooth.lwd=2,
col.response=as.numeric(tit$survived)+2, pt.pch=".", type2="im",
pt.cex=3, jitter=1) # big jitter
set.seed(1238)
a <- earth(pclass~., data=tit)
plotmo(a, type="class", nresponse=1,
grid.levels=list(sex="ma"),
caption="smooth: pclass, sex=\"m\"", SHOWCALL=TRUE,
smooth.col="indianred", smooth.lwd=2,
col.response=as.numeric(tit$pclass)+1, type2="im",
pt.pch=".", pt.cex=3)
plotmo(a, type="class", nresponse=1,
grid.levels=list(sex="ma"),
caption="smooth: pclass, sex=\"m\" jitter=.3", SHOWCALL=TRUE,
smooth.col="indianred", smooth.lwd=2,
col.response=as.numeric(tit$pclass)+1, type2="im",
pt.pch="x", jit=.3) # small jitter
plotmo(a, nresponse=1,
type="class", grid.levels=list(sex="ma"),
caption="smooth: pclass, sex=\"m\"", SHOWCALL=TRUE,
smooth.col="indianred", smooth.lwd=2,
col.response=as.numeric(tit$pclass)+1, type2="im",
pt.pch=paste(1:nrow(tit)))
# test the extend argument
plotmo(a, nresponse=1, pt.col=2, degree2=0, SHOWCALL=TRUE,
caption="test extend: extend=0 (reference plot)")
plotmo(a, nresponse=1, extend=.5, pt.col=2, SHOWCALL=TRUE,
caption="test extend: extend=.5")
plotmo(a, nresponse=1, degree1=0, extend=.2, pt.col=2, SHOWCALL=TRUE) # nothing to plot
a <- earth(survived~pclass+age, data=etitanic, degree=2)
# expect warning: extend=.5 not degree2 plots
plotmo(a, extend=.5, pt.col=2, SHOWCALL=TRUE,
caption="test extend: extend=.5")
# intercept only models
dopar(2, 2, caption = "intercept-only models")
set.seed(1)
x <- 1:10
y <- runif(length(x))
earth.intercept.only <- earth(x, y)
plotmo(earth.intercept.only, do.par=FALSE, main="earth intercept-only model")
plotmo(earth.intercept.only, do.par=FALSE, col.response=1, pt.pch=20)
# TODO following draws a plot but it shouldn't (very minor bug because int-only model with a bad degree1 spec)
plotmo(earth.intercept.only, do.par=FALSE, degree1=3) # expect warning: 'degree1' specified but no degree1 plots
plotmo(earth.intercept.only, do.par=FALSE, degree1=0) # expect warning: plotmo: nothing to plot
library(rpart)
rpart.intercept.only <- rpart(y~x)
plotmo(rpart.intercept.only, do.par=FALSE, main="rpart.plot intercept-only model")
plotmo(rpart.intercept.only, do.par=FALSE, degree1=0)
par(org.par)
# nrug argument
par(mfrow=c(3,3), mar=c(3,3,3,1), mgp=c(1.5, 0.5, 0))
mod.nrug <- earth(survived~age, data=etitanic)
set.seed(2016)
plotmo(mod.nrug, do.par=0, nrug=-1, main="nrug=-1")
plotmo(mod.nrug, do.par=0, nrug=TRUE, main="nrug=TRUE")
plotmo(mod.nrug, do.par=0, nrug=10, rug.col=2, main="nrug=10, rug.col=2")
plotmo(mod.nrug, do.par=0, nrug=5, rug.col=2, rug.lwd=2, main="nrug=5, rug.col=2, rug.lwd=2")
plotmo(mod.nrug, do.par=0, nrug="density", main="nrug=\"density\"")
plotmo(mod.nrug, do.par=0, nrug="density", density.col=2, density.lwd=2, main="nrug=\"density\"\ndensity.col=2, density.lwd=2")
plotmo(mod.nrug, do.par=0, nrug="density", density.adj=.2, density.col=1, main="nrug=\"density\"\ndensity.adj=.2, density.col=1")
par(org.par)
# a <- earth(ozone1[,3]~ozone1[,1]+ozone1[,2]+ozone1[,4]+ozone1[,5]+ozone1[,6], data=ozone1)
# # TODO fails: actual.nrows=330 expected.nrows=50 fitted.nrows=330
# plotmo(a)
# # TODO following fails in plotmo with
# # Error : get.earth.x from model.matrix.earth from predict.earth: x has 2 columns, expected 4 to match: 1 2 3 Girth
# a <- earth(Volume~poly(Height, degree=3)+Girth, data=trees, subset=4:23, linpreds=TRUE)
# plotmo(a, trace=-1, do.par=FALSE, caption="all three rows should be the same")
source("test.epilog.R")
Try the plotmo package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.