6

Does R have a function for weighted least squares? Specifically, I am looking for something that computes intercept and slope.

Data sets

  1. 1 3 5 7 9 11 14 17 19 25 29
  2. 17 31 19 27 31 62 58 35 29 21 18
  3. 102153 104123 96564 125565 132255 115454 114555 132255 129564 126455 124578

The dependent variable is dataset 3 and dataset 1 and 2 are the independent variables.

3
  • 5
    From your comments it doesn't sound like you really want weighted least squares but instead multiple regression. Care to revise your question appropriately? Jun 16, 2011 at 18:30
  • It's still not clear to me (having read the answers and comments below) how weights come into the picture.
    – Ben Bolker
    Jun 16, 2011 at 22:00
  • More about least squares here and here.
    – hhh
    Nov 13, 2011 at 23:18

3 Answers 3

12

Yes, of course, there is a weights= option to lm(), the basic linear model fitting function. Quick example:

R> df <- data.frame(x=1:10)
R> lm(x ~ 1, data=df)            ## i.e. the same as mean(df$x)

Call:
lm(formula = x ~ 1, data = df)

Coefficients:
(Intercept)  
        5.5  

R> lm(x ~ 1, data=df, weights=seq(0.1, 1.0, by=0.1))

Call:
lm(formula = x ~ 1, data = df, weights = seq(0.1, 1, by = 0.1))

Coefficients:
(Intercept)  
          7  

R> 

so by weighing later observations more heavily the mean of the sequence 1 to 10 moves from 5.5 to 7.

6
  • not sure if i follow. if i have 3 data sets how would i estimate intercept and slops. data set 1) :1 3 5 7 9 11 14 17 19 25 29 data set 2:17 31 19 27 31 62 58 35 29 21 18 dataset 3:102153 104123 96564 125565 132255 115454 114555 132255 129564 126455 124578 the dependent variable is dataset 3 and dataset 1 and 2 are the independent variables Jun 16, 2011 at 17:09
  • 3
    Now we don't follow. What are the weights you want to use? Jun 16, 2011 at 17:18
  • 3
    Then it doesn't sound like you want weighted least squares after all. Jun 16, 2011 at 18:29
  • 1
    @Dirk: can you elaborate what the "x ~ 1" actually mean? Chase's reply infers that things such as "y ~ x1 + x2" are like "dependent vars VS independent vars" or am I over-interpreting?
    – hhh
    Nov 13, 2011 at 18:52
  • 1
    @hhh: x~1 is an intercept-only regression model. These formulae are the basis of regression modeling in R: if you need to do regression you should start by reading the relevant section in the Introduction to R (and possibly Faraway's linear modeling book draft in the contributed documents section on the R web site) -- also see the link in Chase's answer
    – Ben Bolker
    Nov 13, 2011 at 19:02
5

First, create your datasets. I'm putting them into a single data.frame but this is not strictly necessary.

dat <- data.frame(x1 = c(1,3,5,7,9,11,14,17,19,25, 29)
                  , x2 = c(17, 31, 19, 27, 31, 62, 58, 35, 29, 21, 18)
                  , y  = c(102153, 104123, 96564, 125565, 132255, 115454
                           , 114555, 132255, 129564, 126455, 124578)
                  )

Second, estimate the model:

> lm(y ~ x1 + x2, data = dat)

Call:
lm(formula = y ~ x1 + x2, data = dat)

Coefficients:
(Intercept)           x1           x2  
  104246.37       906.91        85.76

Third, add your weights as necessary following @Dirk's suggestions.

Fourth and most importantly - read through a tutorial or two on regression in R. Google turns this up as a top hit: http://www.jeremymiles.co.uk/regressionbook/extras/appendix2/R/

4
  • 2
    @RyanB - I would add your data to your questions so that it is all in one place. Also, if this answer your question, clicking the check mark to "accept" it will let others in the future know that this was helpful to you.
    – Chase
    Jun 16, 2011 at 19:01
  • 2
    @RyanB Then please do note the terminology used by @Chase and @Aaron - what you are doing is not a weight least squares (WLS) unless you supply some weights. What @Chase shows is just ordinary least squares. @Dirk's Answer shows you how to start using WLS with the lm() function. Jun 16, 2011 at 19:14
  • 1
    Yes Gavin, the combination of the two responses answers my question. Jun 16, 2011 at 20:38
  • what about if you have a model with random constants before the variables x1 and x2 like "y ~ \beta_{1} x1 + \beta_{2} x2+error" where you assume that the error does not correlate with x1 and x2?
    – hhh
    Nov 13, 2011 at 23:24
0

Just another take on this. You can create a weight matrix first. For example:

samplevar = var(ydata)

M = diag(40,1/samplevar)

At this point M is a 40x40 diagonal matrix. You can convert to a vector by applying diag to M:

M_vector = diag(M)

Then use this in lm :

   lm ( YXDATAFRAME, weights=M_vector)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.