── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag() masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
Today we are interested in how much algal diversity is in a place given the number of grazers present. In short- is there a relationship between grazers (snail) and algae? Bare rock, carcinus, and cancer crabs columns can be dropped as we are not using them. Do this below
#drop columns we don't need!part_ocean_point <- ocean_point %>%select(-one_of('bare_rock', 'carcinus', 'Cancer.crabs'))head(part_ocean_point)
The data are survey data from the rocky intertidal in coastal Maine. Students worked at a transect (trans_description) in different years (year), they measured how high they were above (or below) the low tide line (tidal_ht), those numerical tide heights were converted into a ranked list to make bins (tidal_height_cat_number, 1-8 or so), and then students recorded either the percent of 2D area that was covered by a thing (if it was something that didn’t move) OR how many individuals of a thing they saw (if a crab or snail or something). This abundance or cover data is in columnar form by species (or category) in the remaining columns. This works nicely when we are recording data on a clipboard in the field. BUT, R wants all of those names in 1 column and all of those numbers in a second column. We can use pivot_longer to convert this dataframe to long format.
Pivot the data to long format - one column for names and one for values
Make a new column to tell us if value = %cover or abundance. The grazes (animals) are counted in abundance and the algaes (everything else) are % cover. You may also want to make an algae vs grazer column. Do this below
We are interested in a visual that assess the intermediate disturbance hypothesis, which states that there is some middle range of disturbance (grazing in this case) that leads to the most DIVERSE algal ecosystem. Calculate total grazer count and total algal diversity for each tidal height cat num.
Graph to visually assess the IDH– grazer abundance on the X and algal diversity on the y! What do trends do you see? Do you see any evidence for the IDH? Why or why not?
idh_plot <-ggplot(final_count, aes(x = Graze_Count, y = Algae_Diversity)) +geom_point(size =3) +labs(title ="Grazer Abundace vs Algal Diversity", x ="Grazer Abundance", y ="Algal Diversity")idh_plot
There seems to be some support for the IDH as the highest amount of algal diversity occurs when grazer abundance seems to be in the middle of the range (at about 150) while the data points with a smaller or larger grazer count have a lower amount of algal diversity.
In the following chunk, starting with the raw data from github, build a dplyr/tidyverse pipeline that calculates mean coral cover by reef type. Don’t forget to calculate a relevant and useful error term for plotting
site type lat transect diver cc_percent
1 1 Back Reef 3 1 4 5.8435758
2 1 Back Reef 3 2 4 0.9505263
3 1 Back Reef 3 3 4 5.2423389
4 1 Back Reef 3 4 5 5.0040475
5 1 Back Reef 3 5 5 5.8954916
6 2 Patch Reef 3 1 4 5.2826190
means <- coral_cover %>%group_by(type) %>%summarize(mean =mean(cc_percent, na.rm=TRUE), st_error =std.error(cc_percent, na.rm=TRUE))means
# A tibble: 3 × 3
type mean st_error
<chr> <dbl> <dbl>
1 Back Reef 11.9 1.77
2 Nearshore 3.16 0.531
3 Patch Reef 12.7 2.11
2.Make a visual that allows us to assess the following Null Hypothesis: There is no difference in coral cover by reef type. You should make a graph that shows relevant means, error bars, and raw data
cc_by_reef_type2 <-ggplot(means, aes(x = type, y = mean, color = type)) +geom_point(color ="red") +geom_errorbar(aes(x=type, ymin=mean-st_error, ymax=mean+st_error), width=0.4, alpha=0.9) +geom_jitter(data = coral_cover, aes(x = type, y = cc_percent), alpha =0.4) +labs(title ="Coral Cover Percentage by Reef Type", x ="Reef Type", y ="Coral Cover Percentage") cc_by_reef_type2
Run the correct statistical test to accompany the plot above. This should be a One-Way ANOVA + TukeyHSD test. Run these below AND interpret them. What do you see evidence for? Do we reject or fail to reject the null?
anova_test <-aov(cc_percent ~ type, data = coral_cover)summary(anova_test)
Df Sum Sq Mean Sq F value Pr(>F)
type 2 1172 585.9 6.683 0.00215 **
Residuals 74 6487 87.7
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
We first see from the anova test that we have a high F-value of 6.683 and a low p-value of about 0.00215. Additionally, looking at the p-values between specific coral reef types, we see that patch reef and back reef have no statistically significant difference, with a very high p-value and 0 contained in the confidence interval. However, nearshore is significantly different to back reef as is patch reef to nearshore. Both of these pairs have very low p-values and don’t contain 0 in their confidence intervals. Thus, based on the p-value from our anova test, we can reject the null hypothesis. However, it should be noted that the patch reefs and back reefs do not have the same apparent differences to each other as other reef types have.
Load the penguins data from the palmerpenguins package. Do a two-way ANOVA to assess the effects of species and binary assigned sex at birth (aka the ‘sex’ column here) on beak length. This should be an interactive two-way ANOVA (Y~X*Z). Make a graph and run the stats. Assess the result using the following H0: “There is no effect of sex or species and there is no interaction between sex and species on penguin beak length”
`summarise()` has grouped output by 'species'. You can override using the
`.groups` argument.
means2 <-na.omit(means2)means2
# A tibble: 6 × 4
# Groups: species [3]
species sex mean st_error
<fct> <fct> <dbl> <dbl>
1 Adelie female 37.3 0.237
2 Adelie male 40.4 0.267
3 Chinstrap female 46.6 0.533
4 Chinstrap male 51.1 0.268
5 Gentoo female 45.6 0.269
6 Gentoo male 49.5 0.348
anov2_plot <-ggplot(means2, aes(x = species, y = mean, color = sex), na.rm =TRUE) +geom_point(position =position_dodge(0.5))+geom_errorbar(aes(x=species, ymin=mean-st_error, ymax=mean+st_error), alpha=0.9, position =position_dodge(0.5)) +labs(title="Effect of Penguin Species and Sex on Bill Length", y ="Bill Length (mm)") +scale_fill_brewer(palette="Set3") anov2_plot
What do we see in the graph?
We first see that there’s a big difference between the bill lengths of Adelie penguins compared to both chinstrap and gentoo penguins, on average. There is also a large difference in bill length between male and female penguins of each species type.
#stats go here!anova2_test <-aov(bill_length_mm ~ species*sex, data = penguins)summary(anova2_test)
Df Sum Sq Mean Sq F value Pr(>F)
species 2 7015 3508 654.189 <2e-16 ***
sex 1 1136 1136 211.807 <2e-16 ***
species:sex 2 24 12 2.284 0.103
Residuals 327 1753 5
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
11 observations deleted due to missingness
Species has a very low p-value, indicating that the species type is associated with varying bill length. Sex also has a very low p-value, which means that the sex of the penguin is also associated with varying bill length. However, the p-value of the interaction between species and sex is too high, meaning that the connection between species and sex may not be significant to bill length.
Then looking at the Tukey test, the only interaction that doesn’t have a significant impact is the Gentoo:female-Chinstrap:female. Every other interaction has a significant impact on bill length.
Thus, based on the graph and statistical results we can reject the null hypothesis that species has no effect on bill length and similarly reject the nulls that there sex has no effect on bill length and that there is no interaction between sex and species on penguin beak length.
Part 3: Time series data!
Generally speaking, R has a tough time dealing with time formatting. This is a common issue in programming. There are special data formats called Datetime that we can use to work through this.