Another alternative which is similar to the accepted answer by @BenBolker is redefining the coordinates of the existing plot when adding a second set of points.
Here is a minimal example.
Data:
x <- 1:10
y1 <- rnorm(10, 100, 20)
y2 <- rnorm(10, 1, 1)
Plot:
par(mar=c(5,5,5,5)+0.1, las=1)
plot.new()
plot.window(xlim=range(x), ylim=range(y1))
points(x, y1, col="red", pch=19)
axis(1)
axis(2, col.axis="red")
box()
plot.window(xlim=range(x), ylim=range(y2))
points(x, y2, col="limegreen", pch=19)
axis(4, col.axis="limegreen")