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
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._
