CODEWITHSUNDEEP
rlmnlme
Grillo
Grillo

Reputation: 131

lmList() in ‘nlme’ does not accept weights as an argument

I do not understand why lmList() from the R package ‘nlme’ does not recognize the weights argument to estimate WLS models by groups.

Here is my sample data:

     obs cusip     DLogPrice   DQSize error.var
 1     2 000361AH8  0.657    -1150000     1.02 
 2     3 000361AH8  0.268      300000     0.945
 3     4 000361AH8 -7.18     -1050000     1.28 
 4     5 000361AH8  0.500    -1100000     0.947
 5     6 000361AH8 -0.509     -847000     0.970
 6     7 000361AH8 -0.126     2847000     0.935
 7     8 000361AH8  0               0     0.705
 8     9 000361AH8 -0.144     -105000     0.825
 9    10 000361AH8  0.144      105000     0.828
10    11 000361AH8 -0.144     -200000     0.825
11    12 000361AH8  0.431      570000     0.946
12    13 000361AH8 -0.287     -370000     0.825
13    14 000361AH8  0.882     -600000     0.961
14    15 000361AH8 -0.217      600000     0.925
15     2 000361AK1 -0.155     3000000     1.02 
16     3 000361AK1  0.000104    14000     0.945
17     4 000361AK1 -0.182    -2014000     1.28 
18     5 000361AK1  0.182    -2500000     0.947
19     6 000361AK1 -1.58      2200000     0.970
20     7 000361AK1  0.132      600000     0.935
21     8 000361AK1  1.12     -1300000     0.705
22     9 000361AK1  0.152     2000000     0.825
23    10 000361AK1 -1.10       200000     0.828
24    11 000361AK1 -0.0658   -2400000     0.825
25    12 000361AK1  0.142     1085000     0.946
26    13 000361AK1 -0.470    -1385000     0.825
27    14 000361AK1  0.228     3000000     0.961
28    15 000361AK1 -0.256      629000     0.925

Estimating a simple lm() WLS produces these estimates

> lm(DLogPrice ~ DQSize          , weights = 1/error.var, data=test.data)

Call:
lm(formula = DLogPrice ~ DQSize, data = test.data, weights = 1/error.var)
Coefficients:
(Intercept)       DQSize  
 -1.913e-01    7.096e-09

Using lmList() to run it by group cusip without weights does work too:

> lmList(DLogPrice ~ DQSize | cusip,  data=test.data)
Call:
  Model: DLogPrice ~ DQSize | cusip 
   Data: test.data 

Coefficients:
          (Intercept)        DQSize
000361AH8  -0.3788868  4.171895e-07
000361AK1  -0.1163497 -7.162274e-08

Degrees of freedom: 28 total; 24 residual
Residual standard error: 1.498528

But running lmList() with both by group cusip and weights does not work:

> lmList(DLogPrice ~ DQSize | cusip, weights = 1/error.var, data=test.data)

Error in lmList(DLogPrice ~ DQSize | cusip, , weights = 1/error.var, data = test.data): unused argument (weights = 1/error.var)

Any thoughts on what I am doing wrong?

Upvotes: 1

Views: 177

Answers (1)

Ben Bolker
Ben Bolker

Reputation: 226881

I believe this is simply an oversight in the nlme::lmList() function. However, the version in lme4 (lme4::lmList) does accept weights. The two functions are nearly compatible (nlme::lmList returns an object with a slightly more complicated class structure).

If for some reason you must use nlme, I would suggest submitting a bug report to R-core (but I wouldn't hold my breath ...)

library(lme4)
data(sleepstudy,package="lme4")
set.seed(101)
wts <- runif(nrow(sleepstudy))
lme4::lmList(Reaction~Days|Subject, sleepstudy, weights=wts)

Upvotes: 1

Related Questions