The routes chosen by the two algorithms in the code—Nearest Neighbour (NN) and Ant Colony Optimisation (ACO)—will behave differently in terms of reproducibility when the code is run multiple times.

Nearest Neighbour (NN)

Ant Colony Optimisation (ACO)

Key Factors Influencing ACO Variability

How to Ensure Consistency for ACO

To make the ACO results reproducible, you can set a seed for the random number generator using the set.seed() function in R before running the ant_colony_optimisation() function. As in the 3rd line of the Ant Colony Optimisation Function.

# Load necessary packages
pacman::p_load(pacman, readr, tidyverse, ggplot2, ggrepel, maps, mapproj, geosphere, gridExtra)

# Read the CSV file
cities <- read.csv("World_Top_33_Cities.csv", header = TRUE, sep = ",")

# Calculate the distance matrix using longitude and latitude
distance_matrix <- distm(cities[, c("Longitude", "Latitude")])

# Nearest Neighbour Heuristic Function
nearest_neighbour_tsp <- function(distance_matrix, start_city) {
  visited <- rep(FALSE, nrow(cities))
  route <- integer(nrow(cities))
  current_city <- start_city
  route[1] <- current_city
  visited[current_city] <- TRUE
  
  for (i in 2:nrow(cities)) {
    distances <- distance_matrix[current_city, ]
    distances[visited] <- Inf
    next_city <- which.min(distances)
    route[i] <- next_city
    visited[next_city] <- TRUE
    current_city <- next_city
  }
  
  # Return to start
  route <- c(route, route[1])
  total_distance <- sum(sapply(1:(length(route) - 1), function(i) {
    distance_matrix[route[i], route[i + 1]]
  }))
  
  return(list(route = route, total_distance = total_distance))
}

# Ant Colony Optimisation Function
ant_colony_optimisation <- function(distance_matrix, num_ants = 50, num_iterations = 100, alpha = 1, beta = 5, rho = 0.1, Q = 100) {
  set.seed(123)  # Set a fixed seed for reproducibility
  num_cities <- nrow(distance_matrix)
  pheromones <- matrix(1, nrow = num_cities, ncol = num_cities)
  
  best_tour <- NULL
  best_length <- Inf
  
  for (iter in 1:num_iterations) {
    all_tours <- list()
    all_lengths <- numeric(num_ants)
    
    for (ant in 1:num_ants) {
      start_city <- sample(1:num_cities, 1)
      tour <- integer(num_cities)
      visited <- rep(FALSE, num_cities)
      current_city <- start_city
      tour[1] <- current_city
      visited[current_city] <- TRUE
      
      for (i in 2:num_cities) {
        probabilities <- numeric(num_cities)
        for (j in 1:num_cities) {
          if (!visited[j]) {
            probabilities[j] <- (pheromones[current_city, j] ^ alpha) * ((1 / distance_matrix[current_city, j]) ^ beta)
          }
        }
        probabilities <- probabilities / sum(probabilities)
        next_city <- sample(1:num_cities, 1, prob = probabilities)
        tour[i] <- next_city
        visited[next_city] <- TRUE
        current_city <- next_city
      }
      
      tour <- c(tour, start_city) # Return to start city
      length <- sum(sapply(1:(length(tour) - 1), function(i) {
        distance_matrix[tour[i], tour[i + 1]]
      }))
      
      all_tours[[ant]] <- tour
      all_lengths[ant] <- length
      
      if (length < best_length) {
        best_tour <- tour
        best_length <- length
      }
    }
    
    pheromones <- (1 - rho) * pheromones
    for (ant in 1:num_ants) {
      for (i in 1:(length(all_tours[[ant]]) - 1)) {
        pheromones[all_tours[[ant]][i], all_tours[[ant]][i + 1]] <- pheromones[all_tours[[ant]][i], all_tours[[ant]][i + 1]] + Q / all_lengths[ant]
      }
      pheromones[all_tours[[ant]][length(all_tours[[ant]])], all_tours[[ant]][1]] <- pheromones[all_tours[[ant]][length(all_tours[[ant]])], all_tours[[ant]][1]] + Q / all_lengths[ant]
    }
  }
  
  return(list(route = best_tour, total_distance = best_length))
}

# Run NN and ACO from each possible starting point and compare
shortest_nn_route <- NULL
shortest_nn_distance <- Inf

shortest_aco_route <- NULL
shortest_aco_distance <- Inf

for (i in 1:nrow(cities)) {
  # Run Nearest Neighbour
  nn_result <- nearest_neighbour_tsp(distance_matrix, i)
  if (nn_result$total_distance < shortest_nn_distance) {
    shortest_nn_distance <- nn_result$total_distance
    shortest_nn_route <- nn_result$route
  }
  
  # Run Ant Colony Optimisation
  aco_result <- ant_colony_optimisation(distance_matrix)
  if (aco_result$total_distance < shortest_aco_distance) {
    shortest_aco_distance <- aco_result$total_distance
    shortest_aco_route <- aco_result$route
  }
}

# Convert distances to kilometres and miles
shortest_nn_distance_km <- shortest_nn_distance / 1000
shortest_nn_distance_miles <- shortest_nn_distance / 1609.34

shortest_aco_distance_km <- shortest_aco_distance / 1000
shortest_aco_distance_miles <- shortest_aco_distance / 1609.34

# Create a base map
world <- map_data("world")

# Extract coordinates for routes
nn_coords <- cities[shortest_nn_route, c("Longitude", "Latitude")]
aco_coords <- cities[shortest_aco_route, c("Longitude", "Latitude")]

# Convert routes to strings for table display
nn_route_str <- paste(cities$City[shortest_nn_route], collapse = " → ")
aco_route_str <- paste(cities$City[shortest_aco_route], collapse = " → ")

# Create a data frame for the routes
route_table <- data.frame(
  `Route Type` = c("Nearest Neighbour", "Ant Colony Optimisation"),
  `Cities Visited` = c(nn_route_str, aco_route_str),
  `Distance (km)` = c(round(shortest_nn_distance_km, 2), round(shortest_aco_distance_km, 2)),
  `Distance (miles)` = c(round(shortest_nn_distance_miles, 2), round(shortest_aco_distance_miles, 2))
)

# Create a table grob
route_grob <- tableGrob(route_table, rows = NULL, theme = ttheme_default(base_size = 5.25))

# Plot map
map_plot <- ggplot() +
  geom_polygon(data = world, aes(x = long, y = lat, group = group), fill = "#99cc99", color = "#000000") +
  geom_path(data = as.data.frame(nn_coords), aes(x = Longitude, y = Latitude, color = "NN Route"), linewidth = 1, linetype = "dashed", alpha = 0.7) +
  geom_path(data = as.data.frame(aco_coords), aes(x = Longitude, y = Latitude, color = "ACO Route"), linewidth = 1, alpha = 0.7) +
  geom_point(data = cities, aes(x = Longitude, y = Latitude), color = "#000000", size = 2) +
  geom_text_repel(data = cities, aes(x = Longitude, y = Latitude, label = City), size = 2.5, max.overlaps = 50) +
  scale_color_manual(name = "Route Type", values = c("NN Route" = "#0033ff", "ACO Route" = "#cc0000")) +
  labs(
    title = "Optimised Travel Routes Between the 33 Most Populated Cities in the World; Nearest Neighbour (NN) vs Ant Colony Optimisation (ACO)",
    subtitle = paste("NN Distance:", round(shortest_nn_distance_km, 2), "km /", round(shortest_nn_distance_miles, 2), "miles |",
                     "ACO Distance:", round(shortest_aco_distance_km, 2), "km /", round(shortest_aco_distance_miles, 2), "miles"),
    caption = "Patrick Ford 2024",
    x = "Longitude", y = "Latitude"
  ) +
  theme_minimal() +
  theme(
    plot.title = element_text(size = 16, face = "bold"),
    plot.subtitle = element_text(size = 12),
    legend.position = "bottom"
  )

# Combine plot and table grob
grid.arrange(map_plot, route_grob, ncol = 1, heights = c(7, 1))

# Print routes for reference
print("Order of Cities Visited Nearest Neighbour (NN):")
## [1] "Order of Cities Visited Nearest Neighbour (NN):"
print(cities$City[shortest_nn_route])
##  [1] "Jakarta"          "Ho Chi Minh City" "Bangkok"          "Dhaka"           
##  [5] "Kolkata"          "Delhi"            "Lahore"           "Karachi"         
##  [9] "Mumbai"           "Tehran"           "Cairo"            "Istanbul"        
## [13] "Paris"            "London"           "Lagos"            "Kinshasa"        
## [17] "Nairobi"          "Johannesburg"     "Rio de Janeiro"   "São Paulo"       
## [21] "Buenos Aires"     "Bogotá"           "Mexico City"      "New York City"   
## [25] "Tokyo"            "Osaka"            "Seoul"            "Shanghai"        
## [29] "Beijing"          "Chongqing"        "Guangzhou"        "Manila"          
## [33] "Sydney"           "Jakarta"
print("Order of Cities Visited Ant Colony Optimisation (ACO):")
## [1] "Order of Cities Visited Ant Colony Optimisation (ACO):"
print(cities$City[shortest_aco_route])
##  [1] "Rio de Janeiro"   "São Paulo"        "Buenos Aires"     "Bogotá"          
##  [5] "Mexico City"      "New York City"    "London"           "Paris"           
##  [9] "Istanbul"         "Cairo"            "Tehran"           "Lahore"          
## [13] "Delhi"            "Karachi"          "Mumbai"           "Kolkata"         
## [17] "Dhaka"            "Bangkok"          "Ho Chi Minh City" "Guangzhou"       
## [21] "Chongqing"        "Beijing"          "Seoul"            "Osaka"           
## [25] "Tokyo"            "Shanghai"         "Manila"           "Jakarta"         
## [29] "Sydney"           "Johannesburg"     "Nairobi"          "Kinshasa"        
## [33] "Lagos"            "Rio de Janeiro"