jiri jansa
jiri jansa

Reputation: 159

R - plotting results of nonparametric regression with npreg() ignores par(mfrow = c(1, 2))

I ran 20 nonparametric regressions using command "npreg" and saved the results. Now, I would like to plot (with style band and method asymptotic) these results into 1 multigraph, for example 3 plots (3 outputs) per row. I have tried using par(mfrow=c(,)) or split.screen(), but it seems it is not suitable for these kinds of graphs since it always overwrites the previous graph.

In other words, if I run:

w1 <- lm(w ~ p1)
par(mfrow = c(2, 2))
plot(w1)

then I will obtain 4 graphs in one graphic device. But if I run

w1 <- npreg(w ~ p1)
w2 <- npreg(w ~ p2)
par(mfrow = c(1, 2))` 
plot(w1)
plot(w2)

then the output will be only one graph (and the 2nd overwrites the 1st one). And I would like to have these 2 graphs next to each other.

Upvotes: 3

Views: 4067

Answers (1)

kmm
kmm

Reputation: 6310

The plot method for class npregression objects must be overriding your par setting for some reason. Eventually it calls npplot(), but there is a long chain of intermediate steps that I wasn't able to follow.

Generally, though, I don't think you want plot functions arbitrarily overriding a simple command like mfrow for no apparent reason (or inconsistently), so this seems buggy to me.

I'd suggest emailing the package maintainer and asking why the following code does not produce the same output. You would expect 4 plots in one device.

Here are two of the examples from the package (the first slightly modified) that demonstrate two different behaviors:

4 separate plots

data(Italy)
bw <- npregbw(formula=gdp~ordered(year), data = Italy)
model <- npreg(bws = bw, gradients = TRUE)
par(mfrow = c(2, 2))
plot(model)
points(ordered(Italy$year), Italy$gdp, cex=.2, col="red")
plot(1:10)
plot(1:5)
plot(1)

4 plots in one device

data(cps71)
model.par <- lm(logwage ~ age + I(age^2), data = cps71)
model.np <- npreg(logwage ~ age,
                  regtype = "ll",
                  bwmethod = "cv.aic",
                  gradients = TRUE,
                  data = cps71)
par(mfrow=c(2,2))

plot(cps71$age,cps71$logwage,xlab="age",ylab="log(wage)",cex=.1)
lines(cps71$age,fitted(model.np),lty=1,col="blue")
lines(cps71$age,fitted(model.par),lty=2,col="red")

plot(cps71$age,gradients(model.np),xlab="age",ylab="gradient",type="l",lty=1,col="blue")
lines(cps71$age,coef(model.par)[2]+2*cps71$age*coef(model.par)[3],lty=2,col="red")

plot(cps71$age,fitted(model.np),xlab="age",ylab="log(wage)",ylim=c(min(fitted(model.np)-2*model.np$merr),max(fitted(model.np)+2*model.np$merr)),type="l")
lines(cps71$age,fitted(model.np)+2*model.np$merr,lty=2,col="red")
lines(cps71$age,fitted(model.np)-2*model.np$merr,lty=2,col="red")

plot(cps71$age,gradients(model.np),xlab="age",ylab="gradient",ylim=c(min(gradients(model.np)-2*model.np$gerr),max(gradients(model.np)+2*model.np$gerr)),type="l",lty=1,col="blue")
lines(cps71$age,gradients(model.np)+2*model.np$gerr,lty=2,col="red")
lines(cps71$age,gradients(model.np)-2*model.np$gerr,lty=2,col="red")

Upvotes: 2

Related Questions