[R] contour plot for non-linear models

[ripley@stats.ox.ac.uk]

The current version is in MASS/scripts3/ch08.R.  That runs under R.
I will leave it to you to look for any differences.

Without any idea what you mean by `S' I can't help you there, but
the example in the book was generated using S-PLUS, and scripts are
available ....

On Sat, 8 Jun 2002, Josep Perarnau i Codina wrote:

> Hello all,
>
> I've tried to reproduce the contour plot that appears in the book of
> Venables and Ripley, at page 255. Is a F-statistic surface and a
> confidence region for the regression parameters of a non-linear model.
> It uses the stormer data that are in the MASS package.
>
> I haven't been able to reproduce the plot either in R ( version 1.5 )
> and S.  It makes the axes and it puts the labels, but it doesn't shows
> the contour lines.
>
> Could you help me resolving this problem?
>
> This is the code that appears in the book and I've been testing.
>
> > fm0 <- lm(Wt*Time ~ Viscosity + Time -1, data= stormer)
> > b0 <- coef(fm0); names(b0) <- c("b1", "b2")
> > storm.fm <- nls(Time ~ b1*Viscosity/(Wt-b2), data = stormer, start=b0,
> trace = T)
>
> > bc <- coef(storm.fm)
> > se <- sqrt(diag(vcov(storm.fm)))
> > dv <- deviance(storm.fm)
>
> > par(pty = "s")
> > b1 <- bc[1] + seq(-3*se[1], 3*se[1], length = 51)
> > b2 <- bc[2] + seq(-3*se[2], 3*se[2], length = 51)
> > bv <- expand.grid(b1,b2)
> > ssq <- function(b)
> + sum((Time - b[1]*Viscosity/(Wt-b[2]))^2)
> > dbetas <- apply( bv, 1 , ssq)
> > attach(stormer)
> > cc <- matrix(Time - rep(bv[,1],rep(23,2601)) *
> + Viscosity/(Wt - rep(bv[,2], rep(23, 2601))), 23)
> > dbetas <- matrix(drop(rep(1,23) %% cc^2), 51)
> > fstat <- matrix( ((dbetas -dv)/2) / (dv/21) , 51, 51)
>
> > plot(b1, b2, type = "n")
> > lev <- c(1,2,5,7,10,15,20)
> > contour(b1, b2, fstat, levels= lev, labex = 0.75, lty=2, add=T)
> > contour(b1, b2, fstat, levels= qf(0.95, 2,21), add=T, labex = 0)
> > text(31.6, 0.3, labels = "95% CR", adj= 0, cex= 0.75)
> > points(bc[1], bc[2], pch=3, mkh=0.1)
> > par(pty = "m")
>
>
> Josep Perarnau i Codina,      josep.perarnau at upc.es
> Universitat Politècnica de Catalunya
> Facultat de Matemàtiques i Estadística
> LIOS - Laboratori d'Investigació Operativa i Simulació
> Pau Gargallo 5
> 08028 Barcelona
> Tel:  +34 93 4016947          Fax:  +34 93 4015855
>
>
> -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
> r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
> Send "info", "help", or "[un]subscribe"
> (in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
> _._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
>

-- 
Brian D. Ripley,                  ripley at stats.ox.ac.uk
Professor of Applied Statistics,  http://www.stats.ox.ac.uk/~ripley/
University of Oxford,             Tel:  +44 1865 272861 (self)
1 South Parks Road,                     +44 1865 272860 (secr)
Oxford OX1 3TG, UK                Fax:  +44 1865 272595

-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._