#R-script and data with final simplified models
#Bommarco et al. 2021. Flower strips enhance abundance of bumble bee queens and males in landscapes with few honey bee hives. Biological conservation. 
# https://doi.org/10.1016/j.biocon.2021.109363

#load libraries
library(car)
library(glmmTMB)
library(nlme)

# Honey bee abundance in linear landscape habitats 2018
HBLandsc <- read.table("HBLandsc.tsv", h = TRUE, sep="\t")
mHBL<-glmmTMB(Honeybees ~ Flowerstrip + HB_hives+(1|Site_ID)+(1|Round), 
                data=HBLandsc, family=nbinom2, REML=T)
summary(mHBL)

# Honey bee abundance in sites with flower strips, comparing abundance in the flower strip and in the linear landscape habitat 2018
HBStrip <- read.table("HBStrip.tsv", h = TRUE, sep="\t")
mHBS<-glmmTMB(Honeybees ~ HB_hives*Segment_type+(1|Site_ID)+(1|Round),
                 data=HBStrip, family=nbinom2,REML = T)
summary(mHBS)

# Bumble bee abundance only in linear landscape habitats (including site SV) 2018
BBLandsc <- read.table("BBLandsc.tsv", h = TRUE, sep="\t")
mBBL<-glmmTMB(Bumblebees ~ Flowerstrip*HB_hives+(1|Site_ID)+(1|Round),
                 data=BBLandsc, family=nbinom2, REML=T)
summary(mBBL)

# Bumble bee abundance only in sites with flower strips 2018
BBStrip <- read.table("BBStrip.tsv", h = TRUE, sep="\t")
mBBS<-glmmTMB(Bumblebees ~ Segment_type*HB_hives+(1|Site_ID/Segment_type)+(1|Round),
            data=BBStrip, family=nbinom2, REML=T)
summary(mBBS)

# Bumble bee castes in in linear landscape habitats 2018
BBCasteLandsc <- read.table("BBCasteLandsc.tsv", h = TRUE, sep="\t")
mBBCL<-glmmTMB(Bumblebees ~ Cast2+Flowerstrip+HB_hives+(1|Site_ID)+(1|Round), 
            data=BBCasteLandsc[!is.na(BBCasteLandsc$Cast2),], family=nbinom2, REML=T) 
summary(mBBCL)

# Bumble bee castes in flower strip sites only 2018
BBCasteStrip <- read.table("BBCasteStrip.tsv", h = TRUE, sep="\t")
mBBCS<-glmmTMB(Bumblebees ~  Cast2*HB_hives+Segment_type+(1|Site_ID/Segment_type)+(1|Round),
                data=BBCasteStrip[!is.na(BBCasteStrip$Cast2),], family=nbinom2, REML=T)
summary(mBBCS)

# Bumble bee queens 2019
BBQ2019 <- read.table("BBQ2019.tsv", h = TRUE, sep="\t")
mBBQ19<-glmmTMB(Queens ~ Flower.log+Type+Flowerstrip*HB_hives+(1|Site_ID/Type)+(1|Round)+offset(Width.log), 
                 data=BBQ2019, family=poisson, REML=T)
summary(mBBQ19)