model.matrix (Re: [R] R hang-up using lm)

Göran Broström gb at stat.umu.se
Wed Apr 28 16:00:17 CEST 2004


On Wed, Apr 28, 2004 at 12:59:54PM +0200, Göran Broström wrote:
> On Wed, Apr 28, 2004 at 11:18:38AM +0200, Göran Broström wrote:
> > On Wed, Apr 28, 2004 at 11:46:51AM +0300, Renaud Lancelot wrote:
> > > Martin Maechler a écrit :
> > > >>>>>>"kjetil" == kjetil  <kjetil at acelerate.com>
> > > >>>>>>   on Tue, 27 Apr 2004 19:19:59 -0400 writes:
> > > >
> > > >
> > > >    kjetil> On 27 Apr 2004 at 16:46, Raubertas, Richard wrote:
> > > >    >> Within the last few weeks, someone else reported a similar
> > > >    >> problem when using the results of tapply in a call to rlm().
> > > >    >> Note that the result of tapply is a 1D array, and it appears
> > > >    >> there is a general problem with using such a thing on the
> > > >    >> RHS in formula-based modeling functions:
> > > >    >> 
> > > >    >> set.seed(3)
> > > >    >> yy <- rnorm(20)
> > > >    >> gg <- rep(1:10, 2)
> > > >    >> y <- tapply(yy, gg, median)
> > > >    >> x <- 1:10
> > > >    >> z <- lm(y ~ x)  # OK
> > > >    >> z <- lm(x ~ y)  # crashes R
> > > >    >> 
> > > >    >> (R 1.8.1 on Windows XP Pro)
> > > >    >> 
> > > >
> > > >    kjetil> What exactly do you mean by "crashes R"
> > > >
> > > >    kjetil> Doing this in R1.9.0, windows XP pro, there is no indication 
> > > >    of kjetil> problems.
> > > >
> > > >nor is there with 1.9.0 or R-patched on Linux,
> > > >nor with R 1.8.1 on Linux.
> > > >
> > > >no warning, no error, no problem at all.
> > > >Is it really the above (reproducible, thank you!) example
> > > >that crashes your R 1.8.1 ?
> > > 
> > > It does it for me: Windows XP Pro, R 1.9.0 (P IV, 2.4 GHz, 256 Mo RAM). 
> > > It freezes RGui and a few seconds later, a Windows message appears 
> > > saying that Rgui front-end met a problem and must be closed.
> > 
> > I had to try it too: No crashes on Win2000 pro (1.8.1) or Linux (1.9.0),
> > but (in both cases):
> > 
> > 
> > >  lm(y ~ x)
> > 
> > Call:
> > lm(formula = y ~ x)
> > 
> > Coefficients:
> > (Intercept)            x  
> >     -0.8783       0.1293  
> > 
> > 
> > >  lm(x ~ y)
> > 
> > Call:
> > lm(formula = x ~ y)
> > 
> > Coefficients:
> > (Intercept)  
> >         5.5  
> > 
> > i.e., only an intercept estimate in the second case! Surely something is
> > wrong!? 
> 
> Obviousy, y, generated as above, has an attribute that confuses 'lm',
> because 
> 
> > lm(x ~ as.vector(y))
> 
> works as expected. To add further to confusion, R-1.8.1 (Windows):
> 
> > glm(x ~ y)
> 
> Error in model.matrix.default(mt, mf, contrasts) :
>         invalid type for dimname (must be a vector)

so, 

> model.matrix(x ~ y)
   (Intercept)
1            1
2            1
3            1
4            1
5            1
6            1
7            1
8            1
9            1
10           1
attr(,"assign")
[1] 0

but


> model.matrix(x ~ as.vector(y))
   (Intercept) as.vector(y)
1            1 -0.853357506
2            1 -0.711872147
3            1 -0.228785137
4            1 -0.449739758
5            1  0.173914266
6            1 -0.138766243
7            1 -0.433799800
8            1  0.234183701
9            1  0.002728104
10           1  0.733590165
attr(,"assign")
[1] 0 1

AND

> rr <- model.matrix.default(x ~ y)
Segmenteringsfel
gb at tal:~$ 

After a restart and repeating 

> model.matrix.default(x ~ y)

several times, I FINALLY got

>  model.matrix.default(x ~ y)
   (Intercept)
1            1
2            1
3            1
4            1
5            1
6            1
7            1
8            1
9            1
10           1
attr(,"assign")
[1] 0
>  model.matrix.default(x ~ y)
   (Intercept)            y
1            1 -0.853357506
2            1 -0.711872147
3            1 -0.228785137
4            1 -0.449739758
5            1  0.173914266
6            1 -0.138766243
7            1 -0.433799800
8            1  0.234183701
9            1  0.002728104
10           1  0.733590165
attr(,"assign")
[1] 0 1

>  model.matrix.default(x ~ y)
Segmenteringsfel
gb at tal:~$ 

(Don't ask me what's going on: but 'Segmenteringsfel' means 'Segmentation
fault':) 

BTW, this was on Debian testing/unstable; R-1.9.0
-- 
 Göran Broström                    tel: +46 90 786 5223
 Department of Statistics          fax: +46 90 786 6614
 Umeå University                   http://www.stat.umu.se/egna/gb/
 SE-90187 Umeå, Sweden             e-mail: gb at stat.umu.se




More information about the R-help mailing list