La identificación de diferencias entre las curvas de calibración: ANCOVA?

Question

Tengo un gran número de curvas de calibración (es decir, la respuesta del instrumento vs concentración del analito), y me gustaría comprobar si son diferentes el uno del otro.

Exactamente una curva se asocia con cada combinación de los niveles de los dos factores (de día y de tratamiento). Así, en R, los datos de este aspecto:

set.seed(0)
day = gl(2, 6, labels=c("day1", "day2"))
treat = gl(2, 3, 12, labels=c("A", "B"))
conc = rep(0:2, 4)
response = conc + rnorm(12)/2
df <- data.frame(day=day, treat=treat, conc=conc, response=response)

Me gustaría saber si uno o ambos de los factores que influye en la pendiente de la curva de calibración (no me importa mucho acerca de la intersección). Si son diferentes, voy a utilizar una curva de calibración diferente para cada dataset correspondiente; si ellos no son diferentes, voy a un promedio de ellos para aumentar la precisión de las mediciones que se basan en las curvas de calibración.

Es correcto usar ANCOVA para la prueba de diferencias entre las laderas, por medio de la interacción de day y treat? Si es así, es el correcto código R

model <- lm(response ~ conc + day:treat)

EDIT: Para que quede claro: a lo largo de cada curva de calibración, he medido los datos, que voy a calibrar a través de la curva de calibración. Estoy tratando de entender si debo de utilizar una curva de calibración para cada punto temporal y tratamiento (lo cual podría mejorar la precisión, si hay un efecto de punto temporal y/o tratamiento) o si es aceptable para el promedio de las curvas de calibración (lo cual podría mejorar la precisión). Por lo tanto me gustaría saber si el tratamiento, el punto temporal, o el tratamiento y el punto temporal de la influencia de la pendiente de la curva de calibración. Real datos como un archivo csv es aquí (el número de días aumentará en el futuro).

EDIT2: Gracias a todos los involucrados en la conversación a continuación. Para aclarar, estoy tratando de comprender qué modelo me va a decir si hay diferencias significativas en las pendientes de la respuesta-vs-conc curva, ya sea de día O de tratamiento, O ambos. Puede ser instructivo mirar mis datos reales, publicado aquí, en lugar de la de los falsos datos que he generado por la simplicidad.

require(ggplot2)
df <- read.csv("http://dl.dropbox.com/u/48901983/IntrxnDataSo.csv")
ggplot(df, aes(x=conc, y=response, colour=day)) + 
       geom_point() + geom_line() +
       facet_wrap(~treatment)

Así, parece que hay un efecto de día en la pendiente de la respuesta-vs-conc parcela.

ggplot(df, aes(x=conc, y=response, colour=treatment)) + 
       geom_point() + geom_line() +
       facet_wrap(~day)

No es tan claro que no hay un efecto del tratamiento sobre la pendiente de la respuesta-vs-conc parcela.

int <- ddply(df, .(treatment, day), function(x) coefficients(lm(x$response~x$conc))[1])
names(int)[3] <- "Intercept"
ggplot(int, aes(x=treatment, y=Intercept, colour=day)) + geom_point()

También parece que hay un efecto del tratamiento en la intersección de la respuesta-vs-conc parcela.

Answer 1

Si usted no se preocupan por la intercepción, creo que quiere que este:

> model <- lm(response ~ -1 + conc + day:treat)
> anova(model)
Analysis of Variance Table

Response: response
          Df  Sum Sq Mean Sq F value    Pr(>F)    
conc       1 18.5545 18.5545 62.3458 9.909e-05 ***
day:treat  4  1.7860  0.4465  1.5003    0.2996    
Residuals  7  2.0832  0.2976                      
---
Signif. codes:  0 ‘***' 0.001 ‘**' 0.01 ‘*' 0.05 ‘.' 0.1 ‘ ' 1 
> summary(model)

Call:
lm(formula = response ~ -1 + conc + day:treat)

Residuals:
    Min      1Q  Median      3Q     Max 
-0.5736 -0.4803  0.0222  0.3465  0.6013 

Coefficients:
               Estimate Std. Error t value Pr(>|t|)   
conc             0.6859     0.1929   3.556  0.00927 **
dayday1:treatA   0.6919     0.3693   1.873  0.10316   
dayday2:treatA   0.1093     0.3693   0.296  0.77585   
dayday1:treatB   0.3387     0.3693   0.917  0.38965   
dayday2:treatB   0.7090     0.3693   1.920  0.09636 . 
---
Signif. codes:  0 ‘***' 0.001 ‘**' 0.01 ‘*' 0.05 ‘.' 0.1 ‘ ' 1 

Residual standard error: 0.5455 on 7 degrees of freedom
Multiple R-squared: 0.9071, Adjusted R-squared: 0.8407 
F-statistic: 13.67 on 5 and 7 DF,  p-value: 0.001695 

El -1 en la lm elimina la intersección.

En realidad no hay "pistas" en el modelo. Sólo tienes diferentes días y diferentes tratamientos. En los resultados anteriores, se puede ver que la interacción no es significativa.

Presumiblemente, este no es tu real de datos. En los datos reales, usted puede ser que también desee ver algo como esto: model2 <- lm(response ~ -1 + conc + day + treat) y model3 <- lm(response ~ -1 + conc + day*treat)

Usted puede visualizar con el siguiente código (bastante chapucero)

plot(as.numeric(day),response,col=as.numeric(treat), xlim=c(1,2.2))
points(as.numeric(day)+0.1,pred,col=as.numeric(treat),pch=2)
for(i in 1:6){lines(c(1,1.1), c(response[i],pred[i]),col=as.numeric(treat)[i])}
for(i in 7:12){lines(c(2,2.1), c(response[i],pred[i]),col=as.numeric(treat)[i])}
legend("bottom", c("treatA","treatB","observed","predicted"), col=c(1,2,1,1),
    lty=c(1,1,NA,NA), pch=c(NA,NA,1,2))

Todavía estoy nuevo en intercambio de la pila, así que si alguien agregue un comentario con un enlace que muestra cómo incluir este trazado en lugar de solo el código, se lo agradecería.

EDITAR

Basado en este código, en la siguiente imagen se produce. No dar la pred objeto, así que estoy asumiendo un par de cosas aquí. Los valores predichos son más pequeños y transparentes.

pred <- data.frame(treat = df$treat, conc = df$conc, day = df$day, response = predict(mdl, newdata = df[, 1:3]))

library(ggplot2)
ggplot(df, aes(y = response, x = treat, colour = as.factor(conc))) + 
   geom_jitter(position = position_jitter(width = 0.25), shape = 16, size = 3) +
   geom_point(data = pred, aes(shape = as.factor(conc)), alpha = 0.2) +
   facet_grid(~day) + 
   theme_bw()