I doing a three way Anova (or Ancova - I'm not sure, hence my question) in R. I am testing how temperature, the development stage and the size of a carcass affect the development rate of maggots.

My response variable is Duration (a measurement of hours) and my factors are Size (2 levels = small and large), and Stage (7 levels = eggs, 1st instar, 2nd instar, 3rd instar, postfeed, pupa and total) and Temperature (4 levels = 15, 20, 25, 30).

I intend to examine how Duration varies with Temperature and Size and Stage.

I imported my data set (AnovaTWD). One function I have tried is:

model1 <-(aov(Duration~Stage*Temperature*Size, AnovaTWD)) 
summary(model1)

This gave me:

                        Df  Sum Sq Mean Sq F value   Pr(>F)    
Stage                    6 7206782 1201130 149.059  < 2e-16 ***   
Temperature              1 1924926 1924926 238.881  < 2e-16 ***   
Size                     1   78491   78491   9.741  0.00205 **   
Stage:Temperature        6 2500293  416716  51.714  < 2e-16 ***   
Stage:Size               6  120539   20090   2.493  0.02367 *  
Temperature:Size         1  140090  140090  17.385 4.43e-05 ***   
Stage:Temperature:Size   6  184679   30780   3.820  0.00122 **    
Residuals              214 1724431    8058                     
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Am I right in thinking that R is recognising Temperature as a continuous variable here?

Another function I have tried is:

AnovaTWD$Temperature <- factor(AnovaTWD$Temperature)
model2 <-(aov(Duration~Stage*Temperature*Size, AnovaTWD))
summary(model2)  

This gave me:

                        Df  Sum Sq Mean Sq F value   Pr(>F)        
Stage                    6 7206782 1201130 527.084  < 2e-16 ***    
Temperature              3 2399926  799975 351.048  < 2e-16 ***    
Size                     1  116577  116577  51.157 1.90e-11 ***    
Stage:Temperature       18 3068185  170455  74.800  < 2e-16 ***    
Stage:Size               6  157139   26190  11.493 6.34e-11 ***    
Temperature:Size         3  214981   71660  31.446  < 2e-16 ***    
Stage:Temperature:Size  18  292783   16266   7.138 1.10e-13 ***    
Residuals              186  423861    2279                         
---    
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1    

I know that R is treating Temperature as a fixed factor now. However, I have been advised to treat Temperature as a covariate. I read "ANCOVA is easily reached using the aov() function using the syntax + variable name to indicate that the predictor variable is a covariate."

So I tried:

 summary(aov(Duration~Size*Stage+Temperature, AnovaTWD))

This gave me:

Df  Sum Sq Mean Sq F value Pr(>F)         
Size          1   87904   87904   4.928 0.0274 *      
Stage         6 7223925 1203988  67.497 <2e-16 ***    
Temperature   3 2411455  803818  45.063 <2e-16 ***    
Size:Stage    6  143505   23917   1.341 0.2400         
Residuals   225 4013443   17838                       
---    
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Is this now ANCOVA and treating Temperature as a covariate? Or is it not testing an interaction with the other factors?

I could simply treat Temperature as a fixed variable as my experiment design used the 4 constant temperatures (15, 20, 25, 30). However, I know this approach is not giving all the information possible, hence I am looking into ANCOVA as it is usually used when there are both continuous and categorical predictors.