Add two addLayersControl to one map (have markers be in more than one group)

Question

I have a dataset that includes both a date and a species for each bird observed in a county. I've mapped them using leaflet, but want to use two AddLayersControl to control for both the date and the species. Right now I can only control for the year or the species. I would like the second group of checkboxes so I can control the species as well. I want the marker to go away if either its year group is unchecked or its species group is unchecked.

What I think I need to do is to assign each marker to two different groups that I could control independently. I don't think I am able to assign certain markers as base layers because I don't want a certain subset of them always available. I have also tried just adding another AddLayersControl - sadly the second one will always win and it doesn't seem like you can have two on the same map.

library(leaflet)
library(magrittr)
library(dplyr)
library(htmltools)

# Data
birds <- data.frame(observed_on = c("4/4/2009",
                                    "4/1/2009",
                                    "3/6/2016",
                                    "2/9/2016"),
                    url = c("http://www.inaturalist.org/observations/2236",
                            "http://www.inaturalist.org/observations/2237",
                            "http://www.inaturalist.org/observations/2778201",
                            "https://www.inaturalist.org/observations/9796150"),
                    latitude = c(43.08267975, 
                                 43.0844841, 
                                 43.055512,
                                 43.0180932),
                    longitude = c(-89.43265533, 
                                  -89.43793488, 
                                  -89.314878, 
                                  -89.52836138),
                    scientific_name = c("Agelaius phoeniceus",
                                        "Bubo virginianus", 
                                        "Quiscalus quiscula",
                                        "Strix varia"),
                    common_name = c("Red-winged Blackbird",
                                    "Great Horned Owl",
                                    "Common Grackle",
                                    "Barred Owl"),
                    taxon_order_name = c("Passeriformes",
                                         "Strigiformes",
                                         "Passeriformes",
                                         "Strigiformes"),
                    taxon_species_name = c("Agelaius phoeniceus",
                                           "Bubo virginianus", 
                                           "Quiscalus quiscula",
                                           "Strix varia" ),
                    year = c("2009", "2009", "2016", "2016"))

# Leaflet Chart Formatting --------------------------------------------------------

palette <- colorFactor(palette = rainbow(length(unique(birds$taxon_order_name))),
                       domain = birds$taxon_order_name)

# Leaflet Chart -------------------------------------------------------------------

mymap <-  leaflet(birds) %>% 
  addProviderTiles(providers$CartoDB.Positron) %>% 
  setView(lng = -89.398721,
          lat = 43.071580,
          zoom = 13)

for (t in unique(birds$year)){
  sub <- birds[birds$year == t,]
  labels <- mapply(function(x, y, z, a) {
    HTML(sprintf("%s<br><em>%s</em><br>%s<br><a href=%s>link</a>",
                 htmlEscape(x),
                 htmlEscape(y),
                 htmlEscape(z),
                 htmlEscape(a)))},
    sub$common_name,
    sub$taxon_species_name, 
    sub$observed_on,
    sub$url,
    SIMPLIFY = FALSE)
  
  mymap <- mymap %>%
    addCircleMarkers(data = sub, 
                     lng = ~longitude, 
                     lat = ~latitude,
                     fillOpacity = 0.6, 
                     radius = 8, 
                     fillColor = ~palette(taxon_order_name), 
                     color = "black",
                     weight = 1, 
                     opacity = 0.5,
                     popup = labels,
                     group = as.character(t))
}

mymap %>% 
  addLegend(pal = palette, 
            values = ~taxon_order_name,
            title = "Taxon Order") %>%
  addLayersControl(overlayGroups = as.character(unique(birds$year)), 
                   options = layersControlOptions(collapsed = FALSE))
# addLayersControl(overlayGroups = unique(birds$taxon_order_name), options = layersControlOptions(collapsed = FALSE))

map showing points with both year and species info but layers control for the only year

Answer 1

does this work?

addLayersControl(overlayGroups = as.character(c(unique(birds$year),unique(birds$taxon_order_name)), options = layersControlOptions(collapsed = FALSE))