[Rd] Re: [R] Problem using GLM in a loop (fwd)

Prof Brian Ripley ripley@stats.ox.ac.uk
Tue, 21 Aug 2001 16:22:57 +0100 (BST)


This example is caused by R's messing with formula environments.

That's explained in ?formula, but should it not be explained in
?model.frame ?

Simple test:

data <- data.frame(y=rnorm(100), x=1:100)

testit <- function(formula)
{
   weights <- runif(100)
   glm(formula, weights=weights, data=data)
}
testit(y ~ x)

weights is looked for in the environment of the formula, not of the call
to glm.

-- 
Brian D. Ripley,                  ripley@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

---------- Forwarded message ----------
Date: Tue, 21 Aug 2001 15:47:34 +0100 (BST)
From: Prof Brian Ripley <ripley@stats.ox.ac.uk>
To: Isabelle ZABALZA <isabelle.zabalza-mezghani@ifp.fr>
Subject: Re: [R] Problem using GLM in a loop

It's a scope problem: use

data$weights <- 1/sig

and it will work.  It's finding the system function weights, not your
variable.

I don't really understand this one as yet.


On Tue, 21 Aug 2001, Isabelle ZABALZA wrote:

> Hello,
>
> I am try to perform a modeling which is relevant in a strongly
> heteroscedastic context.
> So I perform a dual modeling (modeling of both mean and variance of a
> response) in using the following loop:
>
> jointmod <- function(formula, data, itercrit=10,devcrit=0.0001)
> {
> #
> # Init step
> #
> init <- glm(formula=formula,family=gaussian, data=data)
> response <- resid(init,type="response") + predict(init,type="response")
> mu <- predict(init,type="response")
> iter <- 1
> dev <- init$deviance
> if (dev != 0) vardev <- 1
> else stop(message="Null deviance in initialisation step")
> endformula_ strsplit(formula,"~")[[3]]
> formulavar_paste("d","~")
> formulavar_as.formula(paste(formulavar,endformula))
> #
> # Modeling loop
> #
> while ((iter <= itercrit)| (vardev > devcrit)) {
> d_(response - mu)^2
> modvar_glm(formula=formulavar,family=Gamma(link=log),data=data)
> sig <- predict(modvar,type="response")
> weights <- 1/sig
> modesp <- glm(formula=formula,
> family=gaussian,data=data,weights=weights)
> mu <- predict(modesp,type="response")
> devold <- dev
> dev <- modesp$deviance
> vardev <- (devold -dev)/devold
> iter <- iter + 1
> }
>
> list(modesp=modesp, modvar=modvar, iter=iter, vardev=vardev)
> }
>
> This program evaluation always stops when evaluating the "modesp" line :
>
>          modesp <- glm(formula=formula,
> family=gaussian,data=data,weights=weights)
> with the following error:
> error in
> model.frame(formula,rownames,variables,varnames,extras,extranames, :
>            invalid variable type.
>
> The same program without specifying weights runs.
> I don't understand this behavior since when I execute the glm command
> with weights out of the program, it runs.
>
> I run R1.3.0 under Windows.
>
> If someone can help me ...
>
> Isabelle Zabalza-Mezghani
>
>
>
>
>
>
> --
> Isabelle Zabalza-Mezghani          Tel : 01 47 52 61 99
> Institut Français du Pétrole       E-mail : isabelle.zabalza-mezghani@ifp.fr
> 1-4 Av. Bois Preau - Bat Lauriers
> 92852 Rueil Malmaison Cedex, France
>
>
>

-- 
Brian D. Ripley,                  ripley@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-devel 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-devel-request@stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._