library

library(ggplot2)
library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library(tidyr)

Data

venom_data <- read.csv("figure4_data_journal.pone.0253838.s006.csv")

str(venom_data)
## 'data.frame':    55 obs. of  12 variables:
##  $ Date                      : chr  "31/01/2020" "31/01/2020" "31/01/2020" "31/01/2020" ...
##  $ Number                    : int  1 2 3 4 5 6 7 8 9 10 ...
##  $ Site                      : chr  "Hovea" "Hovea" "Hovea" "Hovea" ...
##  $ Temperature...C.          : num  20.7 20.7 20.7 20.7 20.7 20.7 20.7 20.7 20.7 20.7 ...
##  $ Humidity..                : int  36 36 36 36 36 35 35 35 35 35 ...
##  $ Weight..g.                : num  0.148 0.096 0.101 0.117 0.08 0.057 0.097 0.062 0.066 0.082 ...
##  $ Hives..number             : int  1 2 3 4 5 1 2 3 4 5 ...
##  $ Behavioural.response      : chr  "active" "active" "active" "active" ...
##  $ Nutritional.supply        : chr  "PSS" "PSS" "PSS" "PSS" ...
##  $ Historical.hives..movement: chr  "S" "S" "S" "S" ...
##  $ Flowering.index           : num  6.65 6.65 6.65 6.65 6.65 4.47 4.47 4.47 4.47 4.47 ...
##  $ Nectar.volume..µl.        : num  1.22 1.22 1.22 1.22 1.22 1.6 1.6 1.6 1.6 1.6 ...
head(venom_data)
##         Date Number    Site Temperature...C. Humidity.. Weight..g.
## 1 31/01/2020      1   Hovea             20.7         36      0.148
## 2 31/01/2020      2   Hovea             20.7         36      0.096
## 3 31/01/2020      3   Hovea             20.7         36      0.101
## 4 31/01/2020      4   Hovea             20.7         36      0.117
## 5 31/01/2020      5   Hovea             20.7         36      0.080
## 6  10/2/2020      6 Chidlow             20.7         35      0.057
##   Hives..number Behavioural.response Nutritional.supply
## 1             1               active                PSS
## 2             2               active                PSS
## 3             3               active                PSS
## 4             4               active                PSS
## 5             5               active                PSS
## 6             1               active                 MP
##   Historical.hives..movement Flowering.index Nectar.volume..µl.
## 1                          S            6.65               1.22
## 2                          S            6.65               1.22
## 3                          S            6.65               1.22
## 4                          S            6.65               1.22
## 5                          S            6.65               1.22
## 6                         MS            4.47               1.60
#View(venom_data)

Research Question 1

Do floral resource variables and environmental conditions explain variation in bee venom weight?

Hypothesis

Bee venom weight will not be strongly explained by nectar volume alone. Instead, venom weight may vary more clearly by site and behavioral response, because colony behavior and local site conditions may influence how much venom is collected.

Figure 1

ggplot(venom_data, aes(x = Nectar.volume..µl., y = Weight..g.)) +
  geom_point(
    aes(color = Site, shape = Behavioural.response),
    size = 3,
    alpha = 0.85
  ) +
  geom_smooth(
    aes(linetype = Behavioural.response),
    method = "lm",
    se = T,
    color = "black",
    linewidth = 1.2
  ) +
  labs(
    title = "Relationship Between Nectar Volume and Bee Venom Weight",
    subtitle = "Points are colored by site; trendlines are separated by behavioral response",
    x = "Nectar Volume (µl)",
    y = "Bee Venom Weight (g)",
    color = "Site",
    shape = "Behavioral Response",
    linetype = "Behavioral Response",
  ) +
  theme_classic() +
  theme(
    plot.title = element_text(face = "bold", size = 14),
    plot.subtitle = element_text(size = 11),
    axis.title = element_text(size = 12),
    legend.position = "right"
  )
## `geom_smooth()` using formula = 'y ~ x'

This figure explores whether nectar availability is associated with bee venom weight. If nectar volume strongly explained venom yield, the points would show a clear upward trend. However, if the points are widely scattered, this suggests that nectar volume alone may not be a strong predictor of venom weight. Differences in site and behavioral response may help explain some of the variation among samples.

Figure 2

library(ggplot2)
library(dplyr)
library(tidyr)

cor_data <- venom_data %>%
  rename(
    "Weight(g)" = Weight..g., 
    "Temp(C)" = Temperature...C., 
    "Humidity" = Humidity.. , 
    "Flowering_Index" = Flowering.index, 
     "Nectar_Volume(µl)" = Nectar.volume..µl.
  ) %>%
  select(
    "Weight(g)",
    "Temp(C)",
    Humidity,
    Flowering_Index,
    "Nectar_Volume(µl)"
  )



# correlation matrix
cor_matrix <- cor(cor_data, use = "complete.obs")

# Convert the correlation matrix into long format for ggplot
cor_long <- as.data.frame(cor_matrix) %>%
  mutate(variable_1 = rownames(.)) %>%
  pivot_longer(
    cols = -variable_1,
    names_to = "variable_2",
    values_to = "correlation"
  )

# Create the heatmap
ggplot(cor_long, aes(x = variable_1, y = variable_2, fill = correlation)) +
  geom_tile(color = "white") +
  geom_text(aes(label = round(correlation, 2)), size = 4) +
  scale_fill_gradient2(
    low = "darkseagreen3",
    mid = "white",
    high = "darkslategrey",
    midpoint = 0,
    limits = c(-1, 1)
  ) +
  labs(
    title = "Correlation Heatmap of Environmental Variables and Bee Venom Weight",
    subtitle = "Pairwise correlations among numeric variables in the venom sample dataset",
    x = "",
    y = "",
    fill = "Correlation",
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(face = "bold", size = 14),
    plot.subtitle = element_text(size = 11),
    axis.text.x = element_text(angle = 45, hjust = 1),
    axis.text.y = element_text(size = 10),
    panel.grid = element_blank()
  )

Overall, these follow-up figures suggest that bee venom weight is not fully explained by nectar volume or any single environmental variable. The scatterplot allows visual comparison of nectar volume and venom weight across sites and behavioral responses, while the correlation heatmap summarizes multiple environmental relationships at once. Together, these visualizations provide a broader view of how ecological and behavioral factors may contribute to variation in bee venom yield.