search

mlcollyer / RRPP

RRPP: An R package for fitting linear models to high-dimensional data using residual randomization
https://mlcollyer.github.io/RRPP/
4 stars 1 forks source link

P-values differ between anova.lm.rrpp and aov with univariate response #1

Closed juliedwhite closed 5 years ago

juliedwhite commented 5 years ago

Hello! Thanks for your work on this package; it's a huge help. I have a question about the calculation of p-values in anova.lm.rrpp, because I'm getting conflicting results when using this and base aov.

Here's my data: test.df.txt and some code that can reproduce my problem: 

library(RRPP)

#Load data
df <- read.table("test.df.txt", header = T, sep = "\t")

#Create an rrpp dataframe
df.rrpp <- rrpp.data.frame(PC1 = df$PC1, Camera = df$Camera, Individual = df$Individual)

#Fit the model for an RCBD with replicates
fit.rrpp <- lm.rrpp(PC1 ~ Camera + Individual/Camera, data = df.rrpp, print.progress = FALSE, iter = 99, SS.type = "I")

#Mixed model with individual as a random effect
anova(fit.rrpp, error = c("Camera:Individual", "Camera:Individual", "Residuals"))

Analysis of Variance, using Residual Randomization Permutation procedure: Randomization of null model residuals Number of permutations: 100 Estimation method: Ordinary Least Squares Sums of Squares and Cross-products: Type I Effect sizes (Z) based on F distributions

Df SS MS Rsq F Z Pr(>F)
Camera 1 4.55 4.5467 0.00917 8.8657 -0.0407 0.64
Individual 34 459.87 13.5257 0.92769 26.3742 1.4805 0.06
Camera:Individual 34 17.44 0.5128 0.03517 5.1789 4.8444 0.01
Residuals 140 13.86 0.0990 0.02797
Total 209 495.72 
#Now, if we run the model with aov
summary(aov(PC1 ~ Camera + Error(Individual/Camera), data = df))

Error:Individual

Df Sum Sq Mean Sq F Pr(>F)
Residuals 34 459.9 13.53

Error: Individual:Camera

Df Sum Sq Mean Sq F Pr(>F)
Camera 1 4.547 4.547 8.866 0.00533
Residuals 34 17.437 0.513

Error: Within

Df Sum Sq Mean Sq F Pr(>F)
Residuals 140 13.86 0.09902

You can see that the F value for camera (8.866) and the SS and MS values for the other variables match what is given by anova.lm.rrpp, but the p-values for Camera are quite different. Could you help me understand why this might be?

Thank you!

mlcollyer commented 5 years ago

Dear Julie,

Thanks for alerting us about this. I’ll be happy to look at your example, but I might not get a chance until tomorrow morning. In the meantime, could you please tell me which version of the package you were using? We just updated the package, so I need to know if this error is current or old. We fixed a bug a couple of versions ago that caused this kind of problem. But if this is occurring with the latest version of our package, then it might be a unique issue.

Thanks! Mike

Sent from my iPhone

On May 24, 2019, at 5:26 PM, juliedwhite notifications@github.com wrote:

Hello! Thanks for your work on this package; it's a huge help. I have a question about the calculation of p-values in anova.lm.rrpp, because I'm getting conflicting results when using this and base aov.

Here's my data: test.df.txt and some code that can reproduce my problem:

library(RRPP)

Load data

df <- read.table("test.df.txt", header = T, sep = "\t")

Create an rrpp dataframe

df.rrpp <- rrpp.data.frame(PC1 = df$PC1, Camera = df$Camera, Individual = df$Individual)

Fit the model for an RCBD with replicates

fit.rrpp <- lm.rrpp(PC1 ~ Camera + Individual/Camera, data = df.rrpp, print.progress = FALSE, iter = 99, SS.type = "I")

Mixed model with individual as a random effect

anova(fit.rrpp, error = c("Camera:Individual", "Camera:Individual", "Residuals")) Analysis of Variance, using Residual Randomization Permutation procedure: Randomization of null model residuals Number of permutations: 100 Estimation method: Ordinary Least Squares Sums of Squares and Cross-products: Type I Effect sizes (Z) based on F distributions

Df SS MS Rsq F Z Pr(>F) Camera 1 4.55 4.5467 0.00917 8.8657 -0.0407 0.64 Individual 34 459.87 13.5257 0.92769 26.3742 1.4805 0.06 Camera:Individual 34 17.44 0.5128 0.03517 5.1789 4.8444 0.01 Residuals 140 13.86 0.0990 0.02797
Total 209 495.72

Now, if we run the model with aov

summary(aov(PC1 ~ Camera + Error(Individual/Camera), data = df)) Error:Individual

Df Sum Sq Mean Sq F Pr(>F) Residuals 34 459.9 13.53
Error: Individual:Camera

Df Sum Sq Mean Sq F Pr(>F) Camera 1 4.547 4.547 8.866 0.00533 Residuals 34 17.437 0.513
Error: Within

Df Sum Sq Mean Sq F Pr(>F) Residuals 140 13.86 0.09902
You can see that the F value for camera (8.866) and the SS and MS values for the other variables match what is given by anova.lm.rrpp, but the p-values for Camera are quite different. Could you help me understand why this might be?

Thank you!

— You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub, or mute the thread.

mlcollyer commented 5 years ago

Dear Julie,

I ran your example and there are not conflicting results if you use SS.type = "II". this makes more sense because the SS for camera would be conditioned on the SS for individual. Using SS.type = "I" introduces an order to the explained SS and because Individual explained almost all of the variation and Camera:Individual so little, it is common to generate very large F values in each random permutation, because the Camera effect is not conditioned on any other term as the first one in. The results almost perfectly match the parametric aov solution with iter-999 and SS.type = "II".

I hope this helps!

Mike

juliedwhite commented 5 years ago

Hi Mike,

Thanks for your help - I should have double checked the other SS types first. This solves my problem. Out of curiosity, do you know why the SS Type I solution implemented by RRPP would mismatch with aov, since it also implements SS Type I? I was using version 0.4.1 as of opening this issue, but just now updated to 0.4.2.9 today and still have the discrepancy.

Thanks! Julie

mlcollyer commented 5 years ago

Hi Julie,

Good question! First, I did try to confirm with your example the answers using lm.rrpp + anova.lm.rrpp, aov, and lm + Anova from the car package, without adjusting error. I got the same results whether using type I or type II SS, and irrespective of package/function. Getting the same results between types I and II SS might just make this look like a coincidence. Regardless, with RRPP, one has to consider the model generating random values, which can be ignored with parametric ANOVA, where the df will provide the probability distribution.

For example, if you use reveal.model.designs, you can see the null (reduced) model that generates random outcomes:

fit.rrpp.1 <- lm.rrpp(PC1 ~ Camera + Individual/Camera, data = df,

  • print.progress = FALSE, iter = 999, SS.type = "I")

fit.rrpp.2 <- lm.rrpp(PC1 ~ Camera + Individual/Camera, data = df,

  • print.progress = FALSE, iter = 999, SS.type = "II”)

reveal.model.designs(fit.rrpp.1) Reduced Full
Camera 1 Camera
Individual Camera Camera + Individual
Camera:Individual Camera + Individual Camera + Individual + Camera:Individual <- Null/Full inherent in pairwise

reveal.model.designs(fit.rrpp.2) Reduced Full
Camera Individual Camera + Individual
Individual Camera Camera + Individual
Camera:Individual Camera + Individual Camera + Individual + Camera:Individual <- Null/Full inherent in pairwise

You can see why the type II SS fit worked better; the null (reduced) model accounted for variation among individuals for the Camera effect. Under a different example - where the two methods do not produce the same SS - maybe the aov and lm.rrpp results would contrast. But in this case, the most appropriate approach when using a method that generates random outcomes from a null model happened to converge on what appears to be parametric type I SS, as type II would have been the same; i.e.,

anova(fit.rrpp.1)

Analysis of Variance, using Residual Randomization Permutation procedure: Randomization of null model residuals Number of permutations: 1000 Estimation method: Ordinary Least Squares Sums of Squares and Cross-products: Type I Effect sizes (Z) based on F distributions

               Df     SS      MS     Rsq        F      Z Pr(>F)

Camera 1 4.55 4.5467 0.00917 45.9146 2.274 0.001 Individual 34 459.87 13.5257 0.92769 136.5901 15.466 0.001 Camera:Individual 34 17.44 0.5128 0.03517 5.1789 5.674 0.001 ** Residuals 140 13.86 0.0990 0.02797
Total 209 495.72

Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘’ 0.05 ‘.’ 0.1 ‘ ’ 1

Call: lm.rrpp(f1 = PC1 ~ Camera + Individual/Camera, iter = 999, SS.type = "I",
data = df, print.progress = FALSE)

anova(fit.rrpp.2)

Analysis of Variance, using Residual Randomization Permutation procedure: Randomization of null model residuals Number of permutations: 1000 Estimation method: Ordinary Least Squares Sums of Squares and Cross-products: Type II Effect sizes (Z) based on F distributions

               Df     SS      MS     Rsq        F       Z Pr(>F)

Camera 1 4.55 4.5467 0.00917 45.9146 2.3234 0.001 Individual 34 459.87 13.5257 0.92769 136.5901 15.4665 0.001 Camera:Individual 34 17.44 0.5128 0.03517 5.1789 5.6740 0.001 ** Residuals 140 13.86 0.0990 0.02797
Total 209 495.72

Signif. codes: 0 ‘’ 0.001 ‘’ 0.01 ‘’ 0.05 ‘.’ 0.1 ‘ ’ 1

Call: lm.rrpp(f1 = PC1 ~ Camera + Individual/Camera, iter = 999, SS.type = "II",
data = df, print.progress = FALSE)

So, I am not sure it is a discrepancy - just misleading!

Cheers! Mike

On May 28, 2019, at 1:56 PM, juliedwhite notifications@github.com wrote:

Hi Mike,

Thanks for your help - I should have double checked the other SS types first. This solves my problem. Out of curiosity, do you know why the SS Type I solution implemented by RRPP would mismatch with aov, since it also implements SS Type I? I was using version 0.4.1 as of opening this issue, but just now updated to 0.4.2.9 today and still have the discrepancy.

Thanks! Julie

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/mlcollyer/RRPP/issues/1?email_source=notifications&email_token=ABUNU4QMTIWRCUFWKTN7XE3PXVW37A5CNFSM4HPSZUJ2YY3PNVWWK3TUL52HS4DFVREXG43VMVBW63LNMVXHJKTDN5WW2ZLOORPWSZGODWM5ZLA#issuecomment-496622764, or mute the thread https://github.com/notifications/unsubscribe-auth/ABUNU4U2UHL4M5KV4XYKEZTPXVW37ANCNFSM4HPSZUJQ.

juliedwhite commented 5 years ago

I see, thanks for your thorough explanation!