# load libraries
library(tidyverse)
## Warning: package 'tidyverse' was built under R version 3.5.2
## -- Attaching packages ----------------------------------------------------------------------------------------------- tidyverse 1.2.1 --
## v ggplot2 3.1.0 v purrr 0.2.5
## v tibble 2.1.1 v dplyr 0.7.7
## v tidyr 0.8.1 v stringr 1.3.1
## v readr 1.1.1 v forcats 0.3.0
## Warning: package 'ggplot2' was built under R version 3.5.2
## Warning: package 'tibble' was built under R version 3.5.3
## -- Conflicts -------------------------------------------------------------------------------------------------- tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()
library(readxl)
library(broom)
library(knitr)
source('theme_javier.R')
# 01 read and preparae data-----
ci <- read_excel('data/Length and CI.xlsx')
dotchart(ci$CI)
# remove one outlier
head(arrange(ci, CI))# sample 32 CI<1
## # A tibble: 6 x 14
## Sample Boat_no Boat_weight Site Density Rep Type Mussel_no Length
## <dbl> <dbl> <dbl> <chr> <chr> <dbl> <chr> <dbl> <dbl>
## 1 32 32 0.897 K C 4 <NA> 2 70.3
## 2 370 220 0.897 S C 3 <NA> 10 51.7
## 3 91 91 0.897 K H 2 <NA> 1 62.7
## 4 223 73 0.897 K C 2 <NA> 3 78.9
## 5 222 72 0.897 K C 2 <NA> 2 79.1
## 6 361 211 0.897 S C 3 <NA> 1 56.8
## # ... with 5 more variables: Flesh_DW_boat <dbl>, Shell_DW <dbl>,
## # Flesh_DW <dbl>, CI <dbl>, Notes <chr>
ci <-
filter(ci,CI>1) %>%
mutate(Site = fct_recode(Site, Koutou = "K", Shipwreck = "S"),
Density = fct_recode(Density, Control = "C", Low = "L", High = "H"),
Density = fct_relevel(Density, "Control", "Low"),
Type = if_else(is.na(Type) & Site=="Koutou", "U", Type),
Type = if_else(is.na(Type), "E", Type))
dotchart(ci$CI)
# 02 Summary stats condion index ----------
ci %>%
drop_na(CI) %>%
group_by(Site, Density, Rep, Type) %>%
summarise(n = n())
## # A tibble: 37 x 5
## # Groups: Site, Density, Rep [?]
## Site Density Rep Type n
## <fct> <fct> <dbl> <chr> <int>
## 1 Koutou Control 1 U 10
## 2 Koutou Control 2 U 10
## 3 Koutou Control 3 U 10
## 4 Koutou Control 4 U 9
## 5 Koutou Control 5 U 10
## 6 Koutou Low 1 U 10
## 7 Koutou Low 2 U 10
## 8 Koutou Low 3 U 10
## 9 Koutou Low 4 U 10
## 10 Koutou Low 5 U 10
## # ... with 27 more rows
ci %>%
drop_na(CI) %>%
group_by(Site) %>%
summarise(mean = mean(CI),
n = n(),
sd = sd(CI),
se = sd/sqrt(n))
## # A tibble: 2 x 5
## Site mean n sd se
## <fct> <dbl> <int> <dbl> <dbl>
## 1 Koutou 8.65 149 1.67 0.137
## 2 Shipwreck 8.38 229 1.84 0.122
# 03 Boxplot of mussel size-----
ggplot(ci, aes(x = Density, y = Length)) +
geom_boxplot() +
facet_grid( Type ~ Site) +
theme_bw()
# 04 Boxplot of condition index ---------
ci_plot2 <-
ci %>%
drop_na(CI) %>%
mutate(Type = fct_recode(Type, Trasplanted = "E", Undisturbed = "U")) %>%
ggplot( aes(x = Density, y = CI)) +
geom_boxplot() +
facet_wrap( ~ paste(Site, Type, sep =" ")) +
theme_javier() +
labs(y = 'Condition index', x = '')
print(ci_plot2)
ggsave(ci_plot2,
filename = 'figures/ci_plots2.tiff',
device = 'tiff',
compression = 'lzw',
width = 7,
height = 2.5,
dpi = 600)
# anova condition index ---------
# undistrubed mussels
und_ci <- filter(ci,Type=="U")
m1 <- aov(sqrt(CI)~Site*Density, und_ci)
kable(tidy(m1))
| Site |
1 |
0.0000259 |
0.0000259 |
0.0003198 |
0.9857455 |
| Density |
2 |
0.5736569 |
0.2868285 |
3.5448399 |
0.0303160 |
| Site:Density |
2 |
0.5623719 |
0.2811859 |
3.4751055 |
0.0324443 |
| Residuals |
253 |
20.4713331 |
0.0809144 |
NA |
NA |
kable(tidy(TukeyHSD(m1)))
| Site |
Shipwreck-Koutou |
-0.0006395 |
-0.0710608 |
0.0697818 |
0.9857455 |
| Density |
Low-Control |
-0.0383865 |
-0.1386382 |
0.0618652 |
0.6390819 |
| Density |
High-Control |
-0.1149013 |
-0.2182213 |
-0.0115812 |
0.0250606 |
| Density |
High-Low |
-0.0765147 |
-0.1795627 |
0.0265333 |
0.1886245 |
| Site:Density |
Shipwreck:Control-Koutou:Control |
0.0527241 |
-0.1213435 |
0.2267917 |
0.9533454 |
| Site:Density |
Koutou:Low-Koutou:Control |
0.0449818 |
-0.1192228 |
0.2091865 |
0.9696023 |
| Site:Density |
Shipwreck:Low-Koutou:Control |
-0.0892733 |
-0.2633408 |
0.0847943 |
0.6819471 |
| Site:Density |
Koutou:High-Koutou:Control |
-0.1190230 |
-0.2832277 |
0.0451816 |
0.3002917 |
| Site:Density |
Shipwreck:High-Koutou:Control |
-0.0447147 |
-0.2340824 |
0.1446529 |
0.9842293 |
| Site:Density |
Koutou:Low-Shipwreck:Control |
-0.0077422 |
-0.1810257 |
0.1655412 |
0.9999951 |
| Site:Density |
Shipwreck:Low-Shipwreck:Control |
-0.1419974 |
-0.3246542 |
0.0406595 |
0.2266302 |
| Site:Density |
Koutou:High-Shipwreck:Control |
-0.1717471 |
-0.3450306 |
0.0015363 |
0.0536005 |
| Site:Density |
Shipwreck:High-Shipwreck:Control |
-0.0974388 |
-0.2947307 |
0.0998531 |
0.7159242 |
| Site:Density |
Shipwreck:Low-Koutou:Low |
-0.1342551 |
-0.3075386 |
0.0390284 |
0.2300066 |
| Site:Density |
Koutou:High-Koutou:Low |
-0.1640049 |
-0.3273781 |
-0.0006317 |
0.0484954 |
| Site:Density |
Shipwreck:High-Koutou:Low |
-0.0896966 |
-0.2783437 |
0.0989506 |
0.7476521 |
| Site:Density |
Koutou:High-Shipwreck:Low |
-0.0297498 |
-0.2030333 |
0.1435337 |
0.9963886 |
| Site:Density |
Shipwreck:High-Shipwreck:Low |
0.0445585 |
-0.1527334 |
0.2418505 |
0.9871012 |
| Site:Density |
Shipwreck:High-Koutou:High |
0.0743083 |
-0.1143388 |
0.2629555 |
0.8680293 |
# transplanted mussels
tr_ci <- filter(ci,Type!="U")
m2 <- aov(sqrt(CI)~Density, tr_ci)
kable(tidy(m2))
| Density |
2 |
1.126326 |
0.5631628 |
6.282525 |
0.002567 |
| Residuals |
116 |
10.398190 |
0.0896396 |
NA |
NA |
kable(tidy(TukeyHSD(m2)))
| Density |
Low-Control |
0.0582341 |
-0.0839308 |
0.2003991 |
0.5956610 |
| Density |
High-Control |
0.2845258 |
0.0929568 |
0.4760948 |
0.0017393 |
| Density |
High-Low |
0.2262917 |
0.0347226 |
0.4178607 |
0.0161716 |
# anova size
m3 <- aov(sqrt(Length)~Site*Density, ci)
tidy(m3)
## # A tibble: 4 x 6
## term df sumsq meansq statistic p.value
## <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
## 1 Site 1 3.96 3.96 4.81 0.0289
## 2 Density 2 2.91 1.46 1.77 0.172
## 3 Site:Density 2 1.11 0.554 0.672 0.511
## 4 Residuals 372 306. 0.823 NA NA
kable(tidy(TukeyHSD(m3)))
| Site |
Shipwreck-Koutou |
-0.2094790 |
-0.3972542 |
-0.0217038 |
0.0288794 |
| Density |
Low-Control |
-0.2042166 |
-0.4598497 |
0.0514164 |
0.1459540 |
| Density |
High-Control |
-0.1091099 |
-0.3898793 |
0.1716596 |
0.6315009 |
| Density |
High-Low |
0.0951068 |
-0.1852452 |
0.3754588 |
0.7043134 |
| Site:Density |
Shipwreck:Control-Koutou:Control |
-0.3503983 |
-0.8118198 |
0.1110231 |
0.2517983 |
| Site:Density |
Koutou:Low-Koutou:Control |
-0.3731248 |
-0.8955741 |
0.1493246 |
0.3184894 |
| Site:Density |
Shipwreck:Low-Koutou:Control |
-0.4617909 |
-0.9232123 |
-0.0003695 |
0.0496832 |
| Site:Density |
Koutou:High-Koutou:Control |
-0.2093309 |
-0.7317802 |
0.3131185 |
0.8609330 |
| Site:Density |
Shipwreck:High-Koutou:Control |
-0.4006701 |
-0.9257514 |
0.1244113 |
0.2467900 |
| Site:Density |
Koutou:Low-Shipwreck:Control |
-0.0227264 |
-0.4811505 |
0.4356976 |
0.9999919 |
| Site:Density |
Shipwreck:Low-Shipwreck:Control |
-0.1113926 |
-0.4988316 |
0.2760464 |
0.9630649 |
| Site:Density |
Koutou:High-Shipwreck:Control |
0.1410675 |
-0.3173566 |
0.5994915 |
0.9507829 |
| Site:Density |
Shipwreck:High-Shipwreck:Control |
-0.0502717 |
-0.5116931 |
0.4111497 |
0.9996034 |
| Site:Density |
Shipwreck:Low-Koutou:Low |
-0.0886662 |
-0.5470902 |
0.3697579 |
0.9937749 |
| Site:Density |
Koutou:High-Koutou:Low |
0.1637939 |
-0.3560101 |
0.6835979 |
0.9456853 |
| Site:Density |
Shipwreck:High-Koutou:Low |
-0.0275453 |
-0.5499946 |
0.4949040 |
0.9999890 |
| Site:Density |
Koutou:High-Shipwreck:Low |
0.2524601 |
-0.2059640 |
0.7108841 |
0.6138875 |
| Site:Density |
Shipwreck:High-Shipwreck:Low |
0.0611209 |
-0.4003005 |
0.5225423 |
0.9989735 |
| Site:Density |
Shipwreck:High-Koutou:High |
-0.1913392 |
-0.7137885 |
0.3311101 |
0.9008231 |
# Initial mussel size data---------------
size <-
read_csv('data/UW-Caliper 0111_Pyura_4_17.CSV') %>%
filter(mm>0)
## Parsed with column specification:
## cols(
## Button = col_integer(),
## Date = col_character(),
## Site = col_character(),
## Time = col_time(format = ""),
## mm = col_double()
## )
summary(size$mm)
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 30.50 44.38 47.25 47.99 52.00 66.50
size %>%
group_by(Site) %>%
summarise(mean = mean(mm),
n = n(),
sd = sd(mm),
se = sd/sqrt(n))
## # A tibble: 2 x 5
## Site mean n sd se
## <chr> <dbl> <int> <dbl> <dbl>
## 1 K 52.4 31 7.51 1.35
## 2 S 45.6 57 6.07 0.804
size_freq_dist <-
ggplot(size, aes(x = mm, stat(count),fill = Site)) +
geom_density(alpha = .3) +
labs(x = 'Mussel size (mm)', y = 'Number of individuals') +
theme_javier() +
theme(legend.position = c(.9,.7))
print(size_freq_dist)
ggsave(size_freq_dist,
filename = 'figures/size_freq_dist.tiff',
device = 'tiff',
compression = 'lzw',
width = 4,
height = 2.5,
dpi = 600)