An analysis of variance test has three steps: partitioning the sum of squares, calculating the mean squares, and computing the test statistic ($F$). R has one function that does it all: aov().
Traveling to a different time zone can cause jetlag and this change in their internal clock is called a phase shift. Researchers in the past think that by exposing the back of the knee to light can help reset the circadian rhythm in humans. They studied the phase shift (measured in hours) in individuals who were either applied light to the eyes, to the knees, or to neither (control).
aov(), conduct an ANOVA test to test the hypotheses. The first argument of aov() takes an object class of formula, which has a syntax of y ~ x, where y is the response variable and x is the explanatory variable. It also needs an argument of data = to specify the data frame. Store this object as aov1.summary(aov1), including the degree of freedom, the sum of squares, the mean squares, the F value, and the P-value. Residuals here means the error (within-group variation).plot(aov1), and specify the second and third graph by which = .TukeyHSD() in base R, which takes a fitted aov object (i.e. TukeyHSD(aov)). What can we conclude?Random-effects ANOVA is particularly useful for quantifying measurement error in data and used when groups are not of interest.
The walking stick is a wingless insect that feeds on plants. Researchers measured the femur length by taking two separate photos of each specimen. They are particularly concerned on the measurement error because the two measurement showed different values.
library(lme4).1 | random. Therefore, fit the model with aov2 <- lmer(femurLength ~ (1 | specimen), data = ___). What is the best estimate of the variance in femur length among individuals in this population?Deviation from normality is a common phenomenon in biological data. From what you have learnt in data transformation and nonparametric tests, deduce the best practice to analyse non-normal data.
Some pea aphids start our red in their life and then change to green later, which might be due to the infection of a bacterium called Rickettsiella. Researchers experimentally injected the bacterium into a sample of red aphids. There are three experimental groups: uninjected (original), injected successfully (infected), or injected but the bacterium failed to establish (uninfected). Colour was measured as hue angle in degrees, where small angles are red and large angles are green.
library(FSA). Perform the Dunn’s test by using dunnTest(), where the first argument takes a formula object. Two arguments are necessary, the data = and the method = . In this case, we will use a correction method called Benjamini-Hochberg procedure, so specify method = "bh".Two-way ANOVA investigates the effect of two factors simultaneously, and usually also tests their interaction effect. It has a similar statistical procedure as one-way ANOVA. However, the statistical procedure can differ when the data is not balanced (orthogonal).
Plants have an optimal soil pH for growth, varying between species. If two plants are grown together at different pH values, the effect of competition may be different (one plant winning in one pH value may lose in another). In a recent growth study, researchers are interested in the growth of Festuca ovina in competition with Calluna vulgaris, because Calluna is well-adapted to acidic soils while Festuca grows on soils with a much wider range of pH.
as.factor().aov(). Then inspect the ANOVA table using summary().plot().interaction.plot(), make an interaction plot for the fitted model. What can you deduce from the plot?TukeyHSD(), but this time, set which = to the interaction term. What is interesting?Drought is one of the drivers of global tree mortality. Broadly speaking, trees have evolved into either one of the two strategies: one is risk-taking, which maximises carbon dioxide uptake for photosynthesis at the risk of hydraulic failure; another is conservation, which closes their stomata (breath pore) very soon to protect themselves. In a greenhouse experiment, we conducted a drought experiment on two species, the Lawson’s cypress and the Sitka spruce. We measured SWC (soil water content), A (photosynthetic assimilation rate), ci (substomatal carbon dioxide concentration), E (transpiration rate), and gs (stomatal conductance).
SWC (i.e. if drought stress has been successfully applied in the trees). Apply a two-way ANOVA. Save the summary() result to an object named aov_5a.A, ci, E, and gs. Discuss the overall result.We will use the library(car) to obtain Type II and Type III SS using the function Anova(). Similar to 5, researchers studied the effect of drought on six commercial tree species over time and measured their midday water potential MWP.
MWP ~ Day + Species + Day:Species. What do you find?MWP ~ Species + Day + Day:Species. Is it the same or different from 6b?summary() on the species, what do you find?Anova() on the Type I sequential ANOVA aov() object, but specifying type = 3. Compare with the Type I ANOVA. Also, swap the order of the effects and see if it changes anything.Anova() on aov(MWP ~ Species + Day), but specifying type = 2. Then, also try swapping the places between Species and Day. Does it make a difference?There is no widely-regarded suitable non-parametric alternative for two-way ANOVA with an interaction term. We may implement a permutation test instead. This is the Manly's approach on the Unrestricted Permutation of Observation (2007).
Using any fitted model(s) from 5a or 5f, we will randomise our values and arrange them in whatever order they came out from randomisation. We would apply independent ANOVA tests to these permutated data sets. We calculate the percentage of replications under the null hypothesis when the resampled F-statistics exceed the obtained F-statistic.
The F-statistics are stored in the summary() of the fitted aov() object from 5a or 5f. They are in the attribute F value. They can be called as aov_5a[[1]]$"F value". What does this object look like?
Store all the F values for comparison later.
# This order completely depends on the order in your formula specified in 5a
Ftreat <- summary(aov_5a)[[1]]$"F value"[1]
Fspecies <- summary(aov_5a)[[1]]$"F value"[2]
Finteract <- summary(aov_5a)[[1]]$"F value"[3]
We will now set up the resampling.
## Setting up the resampling
nreps <- 5000
## Setting up space to store F values as calculated
Ft <- numeric(nreps)
Fs <- numeric(nreps)
Fts <- numeric(nreps)
## The first F values of our 5000 resampling
Ft[1] <- Ftreat
Fs[1] <- Fspecies
Fts[1] <- Finteract
Building a for-loop, we will reiterate the ANOVA procedure for the response variable SWC.
for (i in 2:nreps){
new_SWC <- sample(df$SWC, length(df$SWC))
new_aov <- summary(aov(new_SWC ~ df$Treatment + df$Species + df$Treatment:df$Species))
Ft[i] <- new_aov[[1]]$"F value"[1]
Fs[i] <- ____________________
Fts[i] <- ____________________
}
Calculate the probability. The addition of .Machine$double.eps^0.5 is an important aid in case the two numbers in comparison are different only by floating point computer calculations at the extreme.
Pt <- length(Ft[Ft >= Ftreat + .Machine$double.eps^0.5])/nreps
Ps <- __________________________________________________________
Pts <- __________________________________________________________
Does it give you similar or different result with 5a or 5f?