[R] predict function in regression analysis

Stéphane Adamowicz stephane.adamowicz at avignon.inra.fr
Wed May 6 11:22:19 CEST 2015


Maybe this is what you wanted :

Data <- structure(list(y = c(4.5, 4.5, 4.7, 6.7, 6, 4.4, 4.1, 5.3, 4, 4.2, 4.1, 6.4, 5.5, 3.5, 4.6, 4.1, 4.6, 5, 6.2, 5.9, 3.9, 5.3, 6.9, 5.7), lot = c(1L, 1L, 1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 3L), duration = c(0L, 3L, 6L, 9L, 12L, 15L, 18L, 24L, 0L, 3L, 6L, 9L, 12L, 15L, 18L, 24L, 0L, 3L, 6L, 9L, 12L, 15L, 18L, 24L)), Names = c("y", "lot","duration"), class = "data.frame", row.names = c(NA, -24L))
Data$lot <- as.factor(Data$lot)
summary(Data)

(ANCOVA <- lm(y ~ lot + duration, data=Data))
summary(ANCOVA)
anova(ANCOVA)


#	diagram
color <- c("black", "red", "blue")
plot(y~duration, data=Data, pch=as.character(lot))

P <- predict(ANCOVA, int="c")
for(lev in levels(Data$lot)) {
index <- which(Data$lot==lev)
matlines(x=Data$duration[index], y=P[index,], col=color[as.integer(lev)])
}

Le 5 mai 2015 à 22:41, Bert Gunter <gunter.berton at gene.com> a écrit :

> I think you want CI's for intercepts, not "means" (what is a "mean"
> for a line??). If so, the ?confint function will give this to you for
> the lot effect estimates when applied to a model fitted without an
> intercept:
> 
> myfit <-lm(y~ lot-1+time)
> confint(myfit)
> 
> 
> Further discussion should be directed to a statistical site or a local
> statistician, as these are not R issues.
> 
> Cheers,
> Bert
> 
> 
> Bert Gunter
> Genentech Nonclinical Biostatistics
> (650) 467-7374
> 
> "Data is not information. Information is not knowledge. And knowledge
> is certainly not wisdom."
> Clifford Stoll
> 
> 
> 
> 
> On Tue, May 5, 2015 at 9:53 AM, li li <hannah.hlx at gmail.com> wrote:
>> Hi all,
>>  I have the following data in which there is one factor lot with six
>> levels and one continuous convariate time.
>> I want to fit an Ancova model with common slope and different intercept. So
>> the six lots will have seperate paralell
>> regression lines.I wanted to find the upper 95% confidence limit for the
>> mean of the each of
>> the regression lines. It doesnot seem straightforward to achieve this using
>> predict function. Can anyone give some suggestions?
>> 
>>  Here is my data. I only show the first 3 lots. Also I show the model I
>> used in the end. Thanks very much!
>>     Hanna
>> 
>>      y lot time
>> [1,] 4.5   1    0
>> [2,] 4.5   1    3
>> [3,] 4.7   1    6
>> [4,] 6.7   1    9
>> [5,] 6.0   1   12
>> [6,] 4.4   1   15
>> [7,] 4.1   1   18
>> [8,] 5.3   1   24
>> [9,] 4.0   2    0
>> [10,] 4.2   2    3
>> [11,] 4.1   2    6
>> [12,] 6.4   2    9
>> [13,] 5.5   2   12
>> [14,] 3.5   2   15
>> [15,] 4.6   2   18
>> [16,] 4.1   2   24
>> [17,] 4.6   3    0
>> [18,] 5.0   3    3
>> [19,] 6.2   3    6
>> [20,] 5.9   3    9
>> [21,] 3.9   3   12
>> [22,] 5.3   3   15
>> [23,] 6.9   3   18
>> [24,] 5.7   3   24
>> 
>> 
>>> mod <- lm(y ~ lot+time)
>>> summary(mod)
>> Call:
>> lm(formula = y ~ lot + time)
>> Residuals:
>>    Min      1Q  Median      3Q     Max
>> -1.5666 -0.3344 -0.1343  0.4479  1.8985
>> Coefficients:
>>            Estimate Std. Error t value Pr(>|t|)
>> (Intercept)  4.74373    0.36617  12.955 2.84e-14 ***
>> lot2        -0.47500    0.41129  -1.155   0.2567
>> lot3         0.41250    0.41129   1.003   0.3234
>> lot4         0.96109    0.47943   2.005   0.0535 .
>> lot5         0.98109    0.47943   2.046   0.0490 *
>> lot6        -0.09891    0.47943  -0.206   0.8379
>> time         0.02586    0.02046   1.264   0.2153
>> ---
>> 
>>        [[alternative HTML version deleted]]
>> 
>> ______________________________________________
>> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
> 
> ______________________________________________
> R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.



 
Stéphane Adamowicz
Chercheur / Scientist
stephane.adamowicz at paca.inra.fr
Centre PACA - UR 1115 PSH
Tel. : +33 1 (0)4 32 72 24 35
Fax : +33 1 (0)4 32 72 24 32
228, route de l'Aérodrome
84 914 Avignon Cedex 9
France
www.inra.fr



More information about the R-help mailing list