Data on species richness in constructed wetlands at Pennsylvania Game Lands 297 in south western PA.

Data are the number of species observed in a 0.25 m^2 plot. The original data include the number of stems and percentage cover of each species within the sampling plot. Data collected by Gerry Koza (koz7133@calu.edu) in the summer 2017.

wetlands

Format

A data frame with 78 rows and 5 variables:

site: Gameland
wetland: Wetland number. Can be 1 2 or3
plot: Transect and plot designations
water.depth: depth of water within sapling plot
spp.richness: The number of species observed within the sampling plot

Examples

## Load packages

library(ggplot2)
library(ggpubr)

## Set wetland as a factor
wetlands$wetland <- factor(wetlands$wetland)

if (FALSE) {
## Explore data graphically

### Plot boxplots
ggpubr::ggboxplot(data = wetlands,
          y = "spp.richness",
          x = "wetland",
          fill = "wetland")

### Plot histograms
ggpubr::gghistogram(data = wetlands,
          x = "spp.richness",
          title = "All data")


ggpubr::gghistogram(data = wetlands,
          x = "spp.richness",
          facet.by = "wetland",
          fill = "wetland",
          title = "Faceted by wetland")

## Plot means with 95% confidence intervals
ggpubr::ggerrorplot(wetlands,
  x = "wetland",
  y = "spp.richness",
  desc_stat = "mean_ci",
  add = "mean")
}

## 1-way ANOVA

### null model
model.null <- lm(spp.richness ~ 1, data = wetlands)

### model of interest
model.alt <- lm(spp.richness ~ wetland, data = wetlands)

### compare models
anova(model.null, model.alt)
#> Analysis of Variance Table
#> 
#> Model 1: spp.richness ~ 1
#> Model 2: spp.richness ~ wetland
#>   Res.Df    RSS Df Sum of Sq      F   Pr(>F)   
#> 1     77 494.88                                
#> 2     75 431.63  2    63.258 5.4959 0.005925 **
#> ---
#> Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

## Pairwise comparisons after 1-way ANOVA
### no corrections for multiple comparisons
pairwise.t.test(x = wetlands$spp.richness, g = wetlands$wetland,
      p.adjust.method = "none")
#> 
#> 	Pairwise comparisons using t tests with pooled SD 
#> 
#> data:  wetlands$spp.richness and wetlands$wetland 
#> 
#>   1      2     
#> 2 0.0391 -     
#> 3 0.2041 0.0015
#> 
#> P value adjustment method: none 

### Bonferonni correction
pairwise.t.test(x = wetlands$spp.richness, g = wetlands$wetland,
      p.adjust.method = "bonferroni")
#> 
#> 	Pairwise comparisons using t tests with pooled SD 
#> 
#> data:  wetlands$spp.richness and wetlands$wetland 
#> 
#>   1      2     
#> 2 0.1174 -     
#> 3 0.6122 0.0044
#> 
#> P value adjustment method: bonferroni 

## Tukey test

### re-fit model with aov()
model.alt.aov <- aov(spp.richness ~ wetland, data = wetlands)

### TukeyHSD() on model from aov()
TukeyHSD(model.alt.aov)
#>   Tukey multiple comparisons of means
#>     95% family-wise confidence level
#> 
#> Fit: aov(formula = spp.richness ~ wetland, data = wetlands)
#> 
#> $wetland
#>           diff        lwr       upr     p adj
#> 2-1 -1.4350649 -3.0693145 0.1991846 0.0967414
#> 3-1  0.8214286 -0.7116334 2.3544905 0.4101595
#> 3-2  2.2564935  0.6222440 3.8907430 0.0041631
#> 

### Plot effect sizes
#plotTukeysHSD(TukeyHSD(model.alt.aov))