How do I analyse datasets with large and small values, three explanatory variables, and is not normally distributed?

Question

I have been trying to figure out the best way to analyse my data for some time and I am just working myself into a pit of confusion. Briefly, I have taken a whole bunch of physiological measurements across time from a group of different plants grown under different treatments. I have 4 time points, 8 plant mutants, and 2 treatments. The measurements were ultimately destructive, and so this is NOT a repeated measures design through time - different plants were measured at each time point. I have ~8 replicates per group. At face value this is a simple 3-way ANOVA. However, the measurements were taken 9, 23, 30, and 68 days after sowing and so the values I obtained (e.g. shoot weight) vary considerably (especially, obviously, between day 9 and day 68).

Initially I constructed a model to do just this, though looking at the residuals from this model tells me that they are not normally distributed (Fig. 1 & 2).

Figure 1. spread of residuals from the model aov(shoot.fw.mg ~ line3:treatment:seedling.age)

Figure 2. Q-Q Plot of residuals from the model aov(shoot.fw.mg ~ line3:treatment:seedling.age)

As such, I used the package bestNormalize in R to transform my data. However, the transformed data (which were transformed by the method known as "orderNorm Transformation" in the package) show similarly terrible residuals (Fig. 3).

Figure 3. Q-Q Plot of residuals from the model aov(shoot.fw.mg.bestNorm ~ line3:treatment:seedling.age)

So, I do not know if I have transformed the data in the right way, or if a transformation after a bad Q-Q Plot of the residuals of untransformed data is correct...

Can someone please tell me what kind of analyses I should be doing on these data?

Thank you :) Aaron

Answer 1

I would suggest to only transform your data when absolutely necessary. It's not best practice in the modern era when you're using a computer that can easily handle generalised linear models. Transformations inherently produce less useful data.

Best thing to do is use the glm() function as you are, but to check the residuals from the different families of model to select the best one.

plant.model <- glm(shoot.weight ~ time, data = plant.data, family = Gamma)
autoplot(plant.model, smooth.colour = NA)

So the best bet is to try all the different families (find them by typing ?family() in the console) and pick which one has the best fits in the diagnostic plots

IMPORTANT POINT OF NOTE!!!

If the best model is a quasi model, then the test in your ANOVA should be test="F", whereas the others will be fine with test="Chisq" (note: this does not mean it is doing a Chi Squared test, it is a Chi distribution curve).