Trouble combining county maps into merged regional polygons using gUnaryUnion

Question

I would like to create an amalgamation of counties and treat them as a single unit for analysis and visualization purposes. This regional data would combine multiple counties, e.g., the Bayou region of Louisiana would be an amalgamation of the St. Mary, Terrebonne, Assumption, and Lafourche Parishes (counties). I would want to create a spatial polygon set of coordinates to graph in leaflet, ggplot, etc. However, I have been having trouble plotting the spatial polygon data.

My code creates a list of categorical regions to be applied to the map data. The map is converted to a spatial polygon, and then a spatial polygon dataframe before the borders are dissolved and the map is plotted.

#Import required libraries

library(ggplot2)
library(sp)
library(rgdal)
library(maps)
library(mapdata)
library(terra)

# adapted from 
# http://stackoverflow.com/questions/26062280/converting-a-map-object-to-a-spatialpolygon-object
# https://stackoverflow.com/questions/43174769/r-aggregate-county-map-polygons-to-create-custom-borders


require(sp)
require(maptools)

#Create categorical list for each region

Northwest <- c(
  "bienville",
  "bossier",
  "caddo",
  "claiborne",
  "de soto",
  "jackson",
  "lincoln",
  "natchitoches",
  "red river",
  "sabine",
  "webster"
  )

Northeast <- c(
  "caldwell", 
  "east carroll",
  "franklin",
  "madison",
  "morehouse",
  "ouachita",
  "richland",
  "tensas",
  "union",
  "west carroll"
  )

Central <- c(
  "avoyelles",
  "catahoula",
  "concordia",
  "grant",
  "la salle",
  "rapides",
  "vernon",
  "winn"
  )
Southwest <- c(
  "allen",
  "beauregard",
  "calcasieu",
  "cameron",
  "jefferson davis"
  )
Acadiana <- c(
  "acadia",
  "evangeline",
  "iberia",
  "lafayette",
  "st landry",
  "st martin:south",
  "st martin:north",
  "vermilion"
  )
Capital <- c(
  "ascension",
  "east baton rouge",
  "east feliciana", 
  "iberville",
  "livingston",
  "pointe coupee",
  "st helena",
  "west baton rouge",
  "west feliciana"
  )
Bayou <- c(
  "assumption",
  "lafourche",
  "st mary",
  "terrebonne"
  )

Southeast <- c(
  "jefferson",
  "orleans",
  "plaquemines",
  "st bernard",
  "st charles",
  "st james",
  "st john the baptist",
  "st tammany",
  "tangipahoa",
  "washington"
  )

#Create a map of Louisiana by county.

louisiana <- map("county", 'louisiana', fill = TRUE)

#filter out the parish names from the Louisiana map.

parishes <- sapply(strsplit(louisiana$names, ","), function(x) x[2])

#Create an empty list. This list will be used for the next loop.
region_names<- c()

#substitute Louisiana's name column with regional names attached to parishes. This was we can categorize each parish.

for (i in parishes) {
  if (i %in% Acadiana){
    region_names <- append(region_names, paste('acadiana',i, sep = ',', collapse = NULL))
    }
  else if (i %in% Bayou){
    region_names <- append(region_names, paste('bayou',i, sep = ',', collapse = NULL))
  }
  else if (i %in% Capital){
    region_names <- append(region_names, paste('capital',i, sep = ',', collapse = NULL))
  }
  else if (i %in% Central){
    region_names <- append(region_names, paste('central',i, sep = ',', collapse = NULL))
  }
  else if (i %in% Northeast){
    region_names <- append(region_names, paste('northeast',i, sep = ',', collapse = NULL))
  }
  else if (i %in% Northwest){
    region_names <- append(region_names, paste('northwest',i, sep = ',', collapse = NULL))
  }
  else if (i %in% Southeast){
    region_names <- append(region_names, paste('southeast',i, sep = ',', collapse = NULL))
  }
  else if (i %in% Southwest){
    region_names <- append(region_names, paste('southwest',i, sep = ',', collapse = NULL))
  }
  else{
    print("Change the spelling in the code for:", i)
  }
}

#Substitute the name column in the Louisiana map with our new categorical regional name column.
louisiana$names <- region_names

## Create spatial polygon for The State of Louisiana.
louisiana.sp <- map2SpatialPolygons(louisiana, IDs = as.factor(louisiana$names), 
                                 proj4string = CRS("+proj=longlat +datum=WGS84"))

# Add information data of the polygons

region <- sapply(strsplit(louisiana$names, ","), function(x) x[1])
parish <- sapply(strsplit(louisiana$names, ","), function(x) x[2])
parish[is.na(parish)] <- region[is.na(parish)]

# Create the SpatialPolygonsDataFrame

louisiana.sp.data <- SpatialPolygonsDataFrame(
  louisiana.sp, 
  data = data.frame(region = region,
                    parish = parish),
  match.ID = FALSE)

library(rgeos)
# Because of topology problems


louisiana.sp.data.buffer <- gBuffer(louisiana.sp.data, byid = TRUE, width = 0)


# Merge polygons according to region
louisiana.region <- gUnaryUnion(louisiana.sp.data.buffer, 
                             id = louisiana.sp.data.buffer@data$region)

plot(louisiana.region)

The resulting output is a map of the state county lines, and they have not been dissolved by region. (Map of output)

I specified the id as the regional variable in gUnaryUnion, so I'm not too sure why this did not work. The louisiana.region spatial polygon consists of the 8 polygons I need. I'm not too familiar with working with geospatial data - what can I do so solve this issue?