Infiltration Neuhaus data from May 23
df2 <- read.csv("data2.csv", header = T, sep = ",")
str(df2)
## 'data.frame': 88 obs. of 4 variables:
## $ X.Pdmax : num 0.862 0.713 0.535 0.275 0.279 ...
## $ X.Pdav : num 0.939 0.891 0.479 0.37 0.368 ...
## $ activity: Factor w/ 2 levels "H3PO4","HCl": 2 2 2 2 2 2 2 2 2 2 ...
## $ acid : Factor w/ 2 levels "Active","Inactive": 1 2 1 2 1 2 1 2 1 2 ...
Graphs
par(mfrow=c(1,2))
par(mar=c(5,5,3,2)+0.1)
boxplot(df2$X.Pdmax ~ df2$activity * df2$acid, data=df2, main="X.Pdmax ")
boxplot(df2$X.Pdav ~ df2$activity * df2$acid, data=df2, main="X.Pdav ")

tables
X.Pdmax
with(df2, tapply(df2$X.Pdmax, list(df2$activity, df2$acid), mean))
## Active Inactive
## H3PO4 0.7311461 0.5537704
## HCl 0.7757453 0.7359249
with(df2, tapply(df2$X.Pdmax, list(df2$activity, df2$acid), sd))
## Active Inactive
## H3PO4 0.2643445 0.2783566
## HCl 0.1976088 0.2460650
X.Pdav
with(df2, tapply(df2$X.Pdav, list(df2$activity, df2$acid), mean))
## Active Inactive
## H3PO4 0.9422397 0.7073614
## HCl 0.8612996 0.8381549
with(df2, tapply(df2$X.Pdav, list(df2$activity, df2$acid), sd))
## Active Inactive
## H3PO4 0.3881031 0.3842307
## HCl 0.2181548 0.2686489
Interaction plots
par(mfrow=c(1,2)) #numero de figuras en el plot;c(filas,columnas)
par(mar=c(5,5,3,2)+0.1) #margenes
with(df2, interaction.plot(x.factor=df2$activity, trace.factor=df2$acid,
response=df2$X.Pdmax, fun=mean, type="b", legend=T,
ylab="df2$X.Pdmax", main="Interaction Plot for df2$X.Pdmax",
pch=c(1,19)))
with(df2, interaction.plot(x.factor=df2$activity, trace.factor=df2$acid,
response=df2$X.Pdav, fun=mean, type="b", legend=T,
ylab="df2$X.Pdmax", main="Interaction Plot for df2$X.Pdav",
pch=c(1,19)))

# with ylim 0 to 1
with(df2, interaction.plot(x.factor=df2$activity, trace.factor=df2$acid,
response=df2$X.Pdmax, fun=mean, type="b", legend=T,
ylab="df2$X.Pdmax", main="Interaction Plot for df2$X.Pdmax", ylim = c(0,1),
pch=c(1,19)))
with(df2, interaction.plot(x.factor=df2$activity, trace.factor=df2$acid,
response=df2$X.Pdav, fun=mean, type="b", legend=T,
ylab="df2$X.Pdmax", main="Interaction Plot for df2$X.Pdav", ylim = c(0,1),
pch=c(1,19)))

coplots
coplot(df2$X.Pdmax ~ df2$activity | df2$acid, data=df2, panel=panel.smooth,
xlab="df2$X.Pdmax")

coplot(df2$X.Pdav ~ df2$activity | df2$acid, data=df2, panel=panel.smooth,
xlab="df2$X.Pdav")

The model
aov.X.Pdmax <- aov(df2$X.Pdmax ~ df2$activity * df2$acid, data=df2)
aov.X.Pdav <- aov(df2$X.Pdav ~ df2$activity * df2$acid, data=df2)
summary(aov.X.Pdmax)
## Df Sum Sq Mean Sq F value Pr(>F)
## df2$activity 1 0.282 0.28221 4.529 0.0362 *
## df2$acid 1 0.275 0.27461 4.407 0.0388 *
## df2$activity:df2$acid 1 0.104 0.10385 1.667 0.2002
## Residuals 84 5.234 0.06231
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(aov.X.Pdav)
## Df Sum Sq Mean Sq F value Pr(>F)
## df2$activity 1 0.014 0.0136 0.128 0.721
## df2$acid 1 0.394 0.3940 3.695 0.058 .
## df2$activity:df2$acid 1 0.246 0.2461 2.308 0.132
## Residuals 84 8.957 0.1066
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
model.tables(aov.X.Pdmax, type="means", se=T)
## Design is unbalanced - use se.contrast() for se's
## Tables of means
## Grand mean
##
## 0.6965699
##
## df2$activity
## H3PO4 HCl
## 0.6425 0.7558
## rep 46.0000 42.0000
##
## df2$acid
## Active Inactive
## 0.7524 0.6407
## rep 44.0000 44.0000
##
## df2$activity:df2$acid
## df2$acid
## df2$activity Active Inactive
## H3PO4 0.731 0.554
## rep 23.000 23.000
## HCl 0.776 0.736
## rep 21.000 21.000
test for homogeneity of variance.
library(car)
leveneTest(df2$X.Pdmax ~ df2$activity * df2$acid, data=df2) # ok
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 3 1.3724 0.2568
## 84
leveneTest(df2$X.Pdav ~ df2$activity * df2$acid, data=df2) # ok
## Levene's Test for Homogeneity of Variance (center = median)
## Df F value Pr(>F)
## group 3 2.5969 0.05771 .
## 84
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
visually check assumptions
par(mfrow=c(2,2)) #numero de figuras en el plot;c(filas,columnas)
par(mar=c(5,5,3,2)+0.1) #margenes
plot(aov.X.Pdmax, 1) # visual check of homogeneity of variance
plot(aov.X.Pdmax, 2) # visual check of normality
plot(aov.X.Pdav, 1) # visual check of homogeneity of variance
plot(aov.X.Pdav, 2) # visual check of normality
