Sophie interaction data 2023
Victor Scharnhorst
2023-04-12
This is an R Markdown Notebook. When you execute code within the notebook, the results appear beneath the code.
Try executing this chunk by clicking the Run button within the chunk or by placing your cursor inside it and pressing Ctrl+Shift+Enter.
##Setup
Add a new chunk by clicking the Insert Chunk button on the toolbar or by pressing Ctrl+Alt+I.
When you save the notebook, an HTML file containing the code and output will be saved alongside it (click the Preview button or press Ctrl+Shift+K to preview the HTML file).
The preview shows you a rendered HTML copy of the contents of the editor. Consequently, unlike Knit, Preview does not run any R code chunks. Instead, the output of the chunk when it was last run in the editor is displayed.
##Silene dataset, explore econullnetr
Bipartite
Note: We only had three types of bipartite webs in mind when starting this package: seed-disperser, plant-pollinator and host-parasitoid systems. In how far it makes sense to use these functionalities for other systems (or indeed even for these systems) lies in the hands of the user. Please refer to the literature cited for details on the theory behind the indices.
Networks can be either binary (0/1 or FALSE/TRUE matrices) or quantitative (matrices containing estimates of pairwise interaction strength, usually assumed here to be interaction frequency).
Input for most analyses is an interaction matrix of m nodes (= species) from one group (“higher”) with n nodes (= species) from another group (“lower”), i.e. a n x m matrix, where higher level nodes are in columns, lower level nodes in rows. Column and row names can be provided. This is fundamentally different from “one-mode” networks, which are organised as k x k matrix, i.e. one group of nodes only, in which each node can link (= interact) with each other. Such a format is incompatible with the functions we provide here. (Note, however, that functions as.one.mode and web2edges are convenience functions to morph bipartite networks into one-mode webs. Furthermore, some indices build on one-mode networks and are called from bipartite.)
Before you start with the network, you have to get the data into the right shape. The function frame2webs aims to facilitate this process. Arranging a web, e.g. by size, is supported by sortweb.
The typical first step is to visualise the network. Two functions are on offer here: one (visweb) simply plots the matrix in colours depicting the strength of an interaction and options for re-arranging columns and rows (e.g. to identify compartments or nesting). The other function (plotweb) plots the actual web with participants (as two rows of rectangles) connected by lines (proportional to interaction strength). Both can be customised by many options.
The second step is to calculate indices describing network topography. There are three different levels this can be achieved at: the entire web (using function networklevel), at the level of each group (also using function networklevel) or the individual node (= species; thus somewhat inconsistently called specieslevel). Most other functions in the package are helpers, although some can be called on their own and return the respective result (dfun, H2fun and second.extinct with slope.bipartite).
The third step is to compare results to null models. Many interaction matrices are very incomplete snapshots of the true underlying network (e.g. a one-week sampling of a pollination network on a patch of 4 x 4 meters). As a consequence, many species were rarely observed, many are singletons (only one recording). To make analyses comparable across networks with different sampling intensity and number of species per group, we need a common yardstick. We suggest that users should use a null model, i.e. an algorithm that randomises entries while constraining some web properties (such as dimensions, marginal totals or connectance). The function nullmodel provides a few such null models, but this is a wide field of research and we make no recommendations (actually, we do: see Dormann et al. 2009 and Dormann 2011, both shipping in the doc-folder of this package). You can also simulate networks using genweb or null.distr.
Finally, bipartite comes with 23 quantitative pollination network data sets taken from the NCEAS interaction webs data base (use data(package=“bipartite”) to show their names) and it has a few miscellaneous functions looking at some special features of bipartite networks (such as modularity: computeModules or potential for apparent competition: PAC).
Data
library(rio)
a <- import("Indices_predictors_4FLadjusted.xlsx", stringsAsFactors = T)
str(a)## 'data.frame': 27 obs. of 34 variables:
## $ FL_network : chr "b100" "b101" "b103" "b11" ...
## $ H2 : num 0.562 0.53 0.818 0.356 0.731 ...
## $ niche_HL : num 0.178 0.105 0.188 0.165 0.163 ...
## $ niche_LL : num 0.0667 0.0674 0.0121 0.2245 0.0484 ...
## $ sp_HL : num 20 30 27 15 26 21 16 12 26 18 ...
## $ sp_LL : num 12 16 11 11 9 7 17 11 11 12 ...
## $ for_avail_all : num 2.04 1.7 1.73 1.83 1.5 ...
## $ for_avail_cat_all: chr "very_low" "none" "none" "none" ...
## $ for_avail_1 : num 3.36 3.27 3.53 3.52 2.97 ...
## $ for_avail_cat1 : chr "low" "low" "low" "low" ...
## $ for_avail_cat1per: chr "5-<25%" "5-<25%" "5-<25%" "5-<25%" ...
## $ for_avail_diff : chr "very_low_to_low" "none_to_low" "none_to_low" "none_to_low" ...
## $ Agri_m2 : num 146267 316055 335804 335011 254057 ...
## $ Anthro_m2 : num 2162 6126 30019 2944 36666 ...
## $ SNE_m2 : num 210698 61662 24802 50636 71310 ...
## $ Water_m2 : num NA 556 NA 2033 736 ...
## $ Woody_m2 : num 31498 6226 NA NA 27856 ...
## $ Gesamtergebnis : num 390625 390625 390625 390625 390625 ...
## $ Agri_per : num 37.4 80.9 86 85.8 65 ...
## $ Anthro_per : num 0.554 1.568 7.685 0.754 9.386 ...
## $ SNE_per : num 53.94 15.79 6.35 12.96 18.26 ...
## $ Water_per : num 0 0.142 0 0.521 0.189 ...
## $ Woody_per : num 8.06 1.59 0 0 7.13 ...
## $ FL : chr "b100" "b101" "b103" "b11" ...
## $ bio_m2 : num 15209 1310 0 239939 34286 ...
## $ no_bio_m2 : num 256136 322855 332394 109477 257182 ...
## $ other : num 119281 66460 58231 41208 99157 ...
## $ bio_per : num 3.893 0.335 0 61.424 8.777 ...
## $ no_bio_per : num 65.6 82.7 85.1 28 65.8 ...
## $ other_per : num 30.5 17 14.9 10.5 25.4 ...
## $ SDI_BTTyp : num 2.074 1.026 0.679 0.653 1.547 ...
## $ SEI_BTTyp : num 0.662 0.305 0.257 0.231 0.47 ...
## $ SDI_TypZfg : num 1.802 0.947 0.647 0.617 1.454 ...
## $ SEI_TypZfg : num 0.725 0.328 0.295 0.268 0.503 ...raw_data_list <- import_list("Interaction_worktab_incl2017fieldcountSK_NJ korr v2 Namen ergaenzt NoSa.xlsx", stringsAsFactors = T)
b <- raw_data_list$`merged_final_exBP&NoInter`library(tidyverse)
library(bipartite)
b0 <- b %>%
filter(samplingEffort >= 4) %>%
group_by(FLAECHEN_NR, Bee_sp, ForagPl_sp) %>%
#dplyr::summarise(Abund_sum = sum(Abund_total)) %>%
rename(higher = Bee_sp,
lower = ForagPl_sp,
webID = FLAECHEN_NR,
freq = Abund_total) %>%
relocate(higher, lower, webID, freq) %>%
as.data.frame()
#str(b0)
c <- frame2webs(b0, type.out = "list", varnames = c("lower", "higher", "webID", "freq"), emptylist = TRUE)
#str(b)
#monat = sampling visit
temp <- sapply(c, networklevel, index=c("H2", "Shannon diversity")) %>%
t() %>%
as.data.frame() %>%
rownames_to_column(var = "FL_network")
temp2 <- sapply(c, grouplevel, level = "both", index = c("niche overlap")) %>%
t() %>%
as.data.frame() %>%
rownames_to_column(var = "FL_network")
c1 <- full_join(temp, temp2) %>%
drop_na(H2)Plots
library(tidyverse)
temp <- c1 %>%
filter(FL_network != "b69")
#temp <- anti_join(temp1, a)
a1 <- a %>%
mutate(Shannon_div = temp$`Shannon diversity`) %>%
filter(sp_HL > 5) %>%
relocate(FL_network, Shannon_div) %>%
dplyr::select(1:5, 14:35, -FL) %>%
column_to_rownames(var = "FL_network")
library(corrplot)
temp <- cor(a1)
# as colour
corrplot(temp, method="color")
a2 <- a1 %>%
dplyr::select(-c(5:10, 14, 16:18)) %>%
mutate(neg_log_H2 = -log(H2)) %>%
relocate(neg_log_H2, .after = H2)
temp <- cor(a2)
# as colour
corrplot(temp, method="color")
a3 <- a2 %>%
rownames_to_column(var = "FL_network") %>%
dplyr::select(-c(16:17)) %>%
pivot_longer(cols = c(Shannon_div:niche_LL), names_to = "index_y", values_to = "value_y",
values_drop_na = F) %>%
pivot_longer(cols = c(Agri_per:SEI_BTTyp), names_to = "predictor_x", values_to = "value_x",
values_drop_na = F)library(ggstatsplot)
grouped_ggscatterstats(
data = dplyr::filter(a3, index_y %in% c("neg_log_H2")),
x = value_x,
y = value_y,
grouping.var = predictor_x,
label.var = FL_network,
#label.expression = length > 200,
xlab = "value predictor",
ylab = "neg log(H2)",
ggtheme = ggplot2::theme_grey(),
#ggplot.component = list(ggplot2::scale_x_continuous(breaks = seq(2, 9, 1), limits = (c(2, 9)))),
#plotgrid.args = list(nrow = 1),
annotation.args = list(title = "Relationship between predictor and negative log(H2)")
)
grouped_ggscatterstats(
data = dplyr::filter(a3, index_y %in% c("Shannon_div")),
x = value_x,
y = value_y,
grouping.var = predictor_x,
label.var = FL_network,
#label.expression = length > 200,
xlab = "value predictor",
ylab = "Shannon_div",
ggtheme = ggplot2::theme_grey(),
#ggplot.component = list(ggplot2::scale_x_continuous(breaks = seq(2, 9, 1), limits = (c(2, 9)))),
#plotgrid.args = list(nrow = 1),
annotation.args = list(title = "Relationship between predictor and Shannon div of network")
)
grouped_ggscatterstats(
data = dplyr::filter(a3, index_y %in% c("niche_LL")),
x = value_x,
y = value_y,
grouping.var = predictor_x,
label.var = FL_network,
#label.expression = length > 200,
xlab = "value predictor",
ylab = "niche LL",
ggtheme = ggplot2::theme_grey(),
#ggplot.component = list(ggplot2::scale_x_continuous(breaks = seq(2, 9, 1), limits = (c(2, 9)))),
#plotgrid.args = list(nrow = 1),
annotation.args = list(title = "Relationship between predictor and niche overlap lower level Taxon")
)
grouped_ggscatterstats(
data = dplyr::filter(a3, index_y %in% c("niche_HL")),
x = value_x,
y = value_y,
grouping.var = predictor_x,
label.var = FL_network,
#label.expression = length > 200,
xlab = "value predictor",
ylab = "niche HL",
ggtheme = ggplot2::theme_grey(),
#ggplot.component = list(ggplot2::scale_x_continuous(breaks = seq(2, 9, 1), limits = (c(2, 9)))),
#plotgrid.args = list(nrow = 1),
annotation.args = list(title = "Relationship between predictor and niche overlap higher level Taxon")
)