[R] se's and CI's for fitted lines in multivariate regression analysis

David Winsemius dwinsemius at comcast.net
Sat Oct 20 22:11:02 CEST 2012


On Oct 16, 2012, at 11:58 AM, Sigrid wrote:

> Okay, I've now tried to the predict function and get the SE, although it seem
> to calculate SE for each observation from the line (I assume), while I want
> the CI-interval and SE for each line fitted line for the treatment. I do not
> really understand what  parameter mean these SEs are calculated from when
> there would be several means along the line...?. This is what I get with
> predict:
> 
>> predict(model, se.fit = TRUE, interval = "confidence")
> 
> Another way I can think of to show how well the lines fit the data is to
> look at the intercepts and slopes instead. I can specify the line for each
> level and would then get the estimate of slope and intercept, although I do
> not know how I show the standard errors of the slope and intercept. 
> lm(decrease[treatment=="A"]~colpos[treatment=="A"])
> 
> Call:
> lm(formula = decrease[treatment == "A"] ~ colpos[treatment ==  "A"])
> 
> Coefficients:
>             (Intercept)  colpos[treatment == "A"]  
>                  2.5357                    0.4643  
> 
> Please let me know if you know how to find st. errors for (or st. error for
> slope and intercept) of lines for each factor of a treatment.

The SE's for treatment "slope" will vary depending on the colpos values. Using `predict`, pick the mid-point of the colpos and rowpos variables (which is where the SE of the estimates will be minimized). This should be covered in any basic regression text.

model<-lm(decrease ~ rowpos  + colpos + treatment + treatment:colpos + 0, data=OrchardSprays)

# I do not think the use of the non-intercept version matters here and it's not in general a good practice, but it allows all the parameters to be labeled as you apparently expect.

predict(model, newdata=data.frame(colpos=4.5, treatment=unique(OrchardSprays$treatment), rowpos=mean(OrchardSprays$rowpos) ), se.fit = TRUE, interval = "confidence")
$fit
     fit        lwr       upr
1 35.000  20.331646  49.66835
2 63.125  48.456646  77.79335
3  7.625  -7.043354  22.29335
4 90.250  75.581646 104.91835
5 68.500  53.831646  83.16835
6 69.000  54.331646  83.66835
7 25.250  10.581646  39.91835
8  4.625 -10.043354  19.29335

$se.fit
       1        2        3        4        5        6        7        8 
7.291375 7.291375 7.291375 7.291375 7.291375 7.291375 7.291375 7.291375 

$df
[1] 47

$residual.scale
[1] 20.62312

> unique(OrchardSprays$treatment)
[1] D E B H G F C A
Levels: A B C D E F G H
> with(OrchardSprays, tapply(decrease, treatment, mean) )
     A      B      C      D      E      F      G      H 
 4.625  7.625 25.250 35.000 63.125 69.000 68.500 90.250 

-- 

David Winsemius, MD
Alameda, CA, USA




More information about the R-help mailing list