OiO.lk Community platform!

Oio.lk is an excellent forum for developers, providing a wide range of resources, discussions, and support for those in the developer community. Join oio.lk today to connect with like-minded professionals, share insights, and stay updated on the latest trends and technologies in the development field.
  You need to log in or register to access the solved answers to this problem.
  • You have reached the maximum number of guest views allowed
  • Please register below to remove this limitation

Problem converging to Pareto optimal front of DLTZ1 with NSGA-II implementation

  • Thread starter Thread starter Cuchisius
  • Start date Start date
C

Cuchisius

Guest
I'm trying to implement the NSGA-II to solve multi objective optimization problems, I think my implementation is correct but it is very unstable, DLTZ1 is supposed to be not that hard to solve but if I use the standard g(x) function it converges as expected about 1 out of every 10 tries, it does fine if I manually set g(x) to 0.5. For more basic functions like the ones proposed in the paper the algorithm was first published I'm getting the expected results.

Expected solution: Image showing the 3d plot of the expected set of solutions

Solution im getting: Image showing the 3d plot of the solutions im getting

Its getting there but definitely not correct.

DLTZ1:

Formulation of DLTZ1

My code:(its ready to run as is)

Code:
#Implmentation of NSGA-II algorithm with Object Oriented Programming

import numpy as np
import random
import matplotlib.pyplot as plt

#Class to represent the individual
class Individual:
    def __init__(self, values, objectives):
        # Ensure values are stored as a NumPy array
        self.values = np.array(values)
        self.objectives = objectives
        self.dominated = 0
        self.dominating = []
        self.rank = 0
        self.crowding_distance = 0.0
        self.objective_values = self.objective_function(self.values) 

    def objective_function(self, values):
        # Apply each objective function on the values array
        return np.array([objective(values) for objective in self.objectives])
    
    def recalculate_objective_values(self):
        self.objective_values = self.objective_function(self.values)
  
#Class to represent the NSGA-II algorithm
class NSGAII:
    def __init__(self, population_size, generations, mutation_rate, crossover_rate, objectives, bounds, tournament_size, dimensions):
        self.dimensions = dimensions
        self.population_size = population_size
        self.generations = generations
        self.mutation_rate = mutation_rate
        self.crossover_rate = crossover_rate
        self.bounds = bounds
        self.tournament_size = tournament_size
        self.population = []
        self.offspring_population = []
        self.objectives = objectives
        self.fronts = []
        self.current_generation = 0
        self.initialize_population()
    
    #Function to initialize the population
    def initialize_population(self):
        for i in range(self.population_size):
            values = np.random.uniform(low=[b[0] for b in self.bounds], high=[b[1] for b in self.bounds], size=self.dimensions)
            new_individual = Individual(values, self.objectives)
            self.population.append(new_individual)
      
    #Function to perform the tournament selection
    def tournament_selection(self):
        selected = []
        for i in range(self.tournament_size):
            selected.append(random.choice(self.population))
        best = selected[0]
        for individual in selected:
            if individual.rank < best.rank:
                best = individual
            elif individual.rank == best.rank:
                if individual.crowding_distance > best.crowding_distance:
                    best = individual
        return best
    
    #Function to perform the Simulated Binary Crossover
    def crossover(self, parent1, parent2):
        rand_mask = np.random.rand(len(parent1.values)) < self.crossover_rate
        u = np.random.rand(len(parent1.values))
        eta_c = 15  # Crossover index
        beta = np.where(u < 0.5, (2 * u)**(1 / (1 + eta_c)), (1 / (2 * (1 - u)))**(1 / (1 + eta_c)))

        child1_values = 0.5 * (parent1.values + parent2.values - beta * abs(parent2.values - parent1.values))
        child2_values = 0.5 * (parent1.values + parent2.values + beta * abs(parent2.values - parent1.values))
        
        child1_values = np.where(rand_mask, child1_values, parent1.values)
        child2_values = np.where(rand_mask, child2_values, parent2.values)

        child1_values = np.clip(child1_values, [b[0] for b in self.bounds], [b[1] for b in self.bounds])
        child2_values = np.clip(child2_values, [b[0] for b in self.bounds], [b[1] for b in self.bounds])

        return Individual(child1_values, self.objectives), Individual(child2_values, self.objectives)

    #Function to perform the mutation (poly-nomial mutation)
    def mutation(self, individual):
        mutation_mask = np.random.rand(len(individual.values)) < self.mutation_rate
        u = np.random.rand(len(individual.values))
        
        # Element-wise calculation of delta
        delta = np.minimum(individual.values - np.array([b[0] for b in self.bounds]),
                        np.array([b[1] for b in self.bounds]) - individual.values) / np.array([b[1] - b[0] for b in self.bounds])
        
        eta_m = 20  # Polynomial mutation index
        delta1 = (2 * u + (1 - 2 * u) * (1 - delta)**(eta_m + 1))**(1 / (eta_m + 1)) - 1
        delta2 = 1 - (2 * (1 - u) + 2 * (u - 0.5) * (1 - delta)**(eta_m + 1))**(1 / (eta_m + 1))
        
        deltaq = np.where(u < 0.5, delta1, delta2)
        deltamax = np.array([b[1] - b[0] for b in self.bounds])
        new_values = individual.values + mutation_mask * deltaq * deltamax  # Apply mutation mask here
        
        # Clipping the values to ensure they remain within bounds
        new_values = np.clip(new_values, [b[0] for b in self.bounds], [b[1] for b in self.bounds])
        
        # Apply the new values if mutation occurred
        if np.any(mutation_mask):
            individual.values = new_values
            individual.recalculate_objective_values()



    # Function to check if one individual dominates another
    def dominates(self, individual1, individual2):
        # Check if individual1 dominates individual2
        return all(individual1.objective_values <= individual2.objective_values) and any(individual1.objective_values < individual2.objective_values)

    #Function to perform the non-dominated sorting
    def non_dominated_sorting(self):
        fronts = []
        front = []
        for individual in self.population:
            individual.dominated = 0
            individual.dominating = []
            for other_individual in self.population:
                if self.dominates(individual, other_individual):
                    individual.dominating.append(other_individual)
                elif self.dominates(other_individual, individual):
                    individual.dominated += 1
            if individual.dominated == 0:
                individual.rank = 0
                front.append(individual)
        fronts.append(front)
        i = 0
        while len(fronts[i]) > 0:
            next_front = []
            for individual in fronts[i]:
                for other_individual in individual.dominating:
                    other_individual.dominated -= 1
                    if other_individual.dominated == 0:
                        other_individual.rank = i + 1
                        next_front.append(other_individual)
            i += 1
            fronts.append(next_front)
        if len(fronts[-1]) == 0:
            fronts.pop(-1)
        self.fronts = fronts


    #Function to perform the crowding distance assignment
    def crowding_distance_assignment(self, front):
        for individual in front:
            individual.crowding_distance = 0.0
        for i in range(len(self.objectives)):
            front.sort(key = lambda x: x.objective_values[i])
            front[0].crowding_distance = float('inf')
            front[-1].crowding_distance = float('inf')
            for j in range(1, len(front) - 1):
                front[j].crowding_distance += (front[j + 1].objective_values[i] - front[j - 1].objective_values[i]) / (front[-1].objective_values[i] - front[0].objective_values[i] + 1e-6)

    
    #Function to perform the crowding distance sorting
    def crowding_distance_sorting(self):
        for front in self.fronts:
            self.crowding_distance_assignment(front)
            front.sort(key = lambda x: x.crowding_distance, reverse = True)
    
    #Function to perform the NSGA-II algorithm
    def run(self):
        self.non_dominated_sorting()
        self.crowding_distance_sorting()

        while self.current_generation < self.generations:
            self.offspring_population = []
            for i in range(self.population_size//2):
                parent1 = self.tournament_selection()
                parent2 = self.tournament_selection()
                child1, child2 = self.crossover(parent1, parent2)
                self.mutation(child1)
                self.mutation(child2)
                self.offspring_population.append(child1)
                self.offspring_population.append(child2)

            
            self.population.extend(self.offspring_population)
            self.non_dominated_sorting()
            self.crowding_distance_sorting()
            self.population = []
            for front in self.fronts:
                if len(self.population) + len(front) <= self.population_size:
                    self.population.extend(front)
                else:
                    self.crowding_distance_assignment(front)
                    front.sort(key = lambda x: x.crowding_distance, reverse = True)
                    self.population.extend(front[:self.population_size - len(self.population)])
            self.current_generation += 1
        return self.fronts[0]

def generate_dtlz1_objectives(M, n):
    """
    Generates a list of DTLZ1 objective functions.
    
    Parameters:
    M (int): Number of objectives.
    n (int): Number of dimensions (variables).
    
    Returns:
    list: A list of lambda functions representing the DTLZ1 objective functions.
    """
    k = n - M + 1  # Number of variables involved in the g function

    def g(x):
        xM = np.array(x[-k:])
        return 100 * (k + np.sum((xM - 0.5)**2 - np.cos(20 * np.pi * (xM - 0.5))))

    objectives = []
    for i in range(1, M + 1):
        if i == 1:
            # For the first objective, use all except the last element
            objectives.append(lambda x: 0.5 * np.prod(np.array(x)[:M-1]) * (1 + g(x)))
        elif i < M:
            # For middle objectives, exclude progressively more elements
            objectives.append(lambda x, i=i: 0.5 * np.prod(np.array(x)[:M-i]) * (1 - x[M-i]) * (1 + g(x)))
        else:
            # For the last objective, use the negation of the first element
            objectives.append(lambda x, i=i: 0.5 * (1 - x[0]) * (1 + g(x)))

    return objectives

def run1():
    # Example of using the NSGA-II algorithm
    population_size = 100
    generations = 300
    dimensions = 7 
    mutation_rate = 1/dimensions
    crossover_rate = 1
    tournament_size = 2
    M = 3  # Number of objectives
    objectives = generate_dtlz1_objectives(M, dimensions)

    #Example of using the NSGA-II algorithm
    bounds = [(0, 1) for _ in range(dimensions)]
    nsga = NSGAII(population_size, generations, mutation_rate, crossover_rate, objectives, bounds, tournament_size, dimensions)
    front = nsga.run()
    print(front[0].values)
    #plot 3d 
    fig = plt.figure()
    ax = fig.add_subplot(111, projection='3d')
    ax.scatter([individual.objective_values[0] for individual in front], [individual.objective_values[1] for individual in front], [individual.objective_values[2] for individual in front])
    ax.set_xlabel('f1')
    ax.set_ylabel('f2')
    ax.set_zlabel('f3')
    plt.show()

run1()
<p>I'm trying to implement the NSGA-II to solve multi objective optimization problems, I think my implementation is correct but it is very unstable, DLTZ1 is supposed to be not that hard to solve but if I use the standard g(x) function it converges as expected about 1 out of every 10 tries, it does fine if I manually set g(x) to 0.5. For more basic functions like the ones proposed in the paper the algorithm was first published I'm getting the expected results.</p>
<p>Expected solution:
<a href="https://i.sstatic.net/XssjYxcg.png" rel="nofollow noreferrer">Image showing the 3d plot of the expected set of solutions</a></p>
<p>Solution im getting:
<a href="https://i.sstatic.net/TeKDcQJj.png" rel="nofollow noreferrer">Image showing the 3d plot of the solutions im getting</a></p>
<p>Its getting there but definitely not correct.</p>
<p>DLTZ1:</p>
<p><a href="https://i.sstatic.net/wieLZwSY.png" rel="nofollow noreferrer">Formulation of DLTZ1</a></p>
<p>My code:(its ready to run as is)</p>
<pre><code>#Implmentation of NSGA-II algorithm with Object Oriented Programming

import numpy as np
import random
import matplotlib.pyplot as plt

#Class to represent the individual
class Individual:
def __init__(self, values, objectives):
# Ensure values are stored as a NumPy array
self.values = np.array(values)
self.objectives = objectives
self.dominated = 0
self.dominating = []
self.rank = 0
self.crowding_distance = 0.0
self.objective_values = self.objective_function(self.values)

def objective_function(self, values):
# Apply each objective function on the values array
return np.array([objective(values) for objective in self.objectives])

def recalculate_objective_values(self):
self.objective_values = self.objective_function(self.values)

#Class to represent the NSGA-II algorithm
class NSGAII:
def __init__(self, population_size, generations, mutation_rate, crossover_rate, objectives, bounds, tournament_size, dimensions):
self.dimensions = dimensions
self.population_size = population_size
self.generations = generations
self.mutation_rate = mutation_rate
self.crossover_rate = crossover_rate
self.bounds = bounds
self.tournament_size = tournament_size
self.population = []
self.offspring_population = []
self.objectives = objectives
self.fronts = []
self.current_generation = 0
self.initialize_population()

#Function to initialize the population
def initialize_population(self):
for i in range(self.population_size):
values = np.random.uniform(low=[b[0] for b in self.bounds], high=[b[1] for b in self.bounds], size=self.dimensions)
new_individual = Individual(values, self.objectives)
self.population.append(new_individual)

#Function to perform the tournament selection
def tournament_selection(self):
selected = []
for i in range(self.tournament_size):
selected.append(random.choice(self.population))
best = selected[0]
for individual in selected:
if individual.rank < best.rank:
best = individual
elif individual.rank == best.rank:
if individual.crowding_distance > best.crowding_distance:
best = individual
return best

#Function to perform the Simulated Binary Crossover
def crossover(self, parent1, parent2):
rand_mask = np.random.rand(len(parent1.values)) < self.crossover_rate
u = np.random.rand(len(parent1.values))
eta_c = 15 # Crossover index
beta = np.where(u < 0.5, (2 * u)**(1 / (1 + eta_c)), (1 / (2 * (1 - u)))**(1 / (1 + eta_c)))

child1_values = 0.5 * (parent1.values + parent2.values - beta * abs(parent2.values - parent1.values))
child2_values = 0.5 * (parent1.values + parent2.values + beta * abs(parent2.values - parent1.values))

child1_values = np.where(rand_mask, child1_values, parent1.values)
child2_values = np.where(rand_mask, child2_values, parent2.values)

child1_values = np.clip(child1_values, [b[0] for b in self.bounds], [b[1] for b in self.bounds])
child2_values = np.clip(child2_values, [b[0] for b in self.bounds], [b[1] for b in self.bounds])

return Individual(child1_values, self.objectives), Individual(child2_values, self.objectives)

#Function to perform the mutation (poly-nomial mutation)
def mutation(self, individual):
mutation_mask = np.random.rand(len(individual.values)) < self.mutation_rate
u = np.random.rand(len(individual.values))

# Element-wise calculation of delta
delta = np.minimum(individual.values - np.array([b[0] for b in self.bounds]),
np.array([b[1] for b in self.bounds]) - individual.values) / np.array([b[1] - b[0] for b in self.bounds])

eta_m = 20 # Polynomial mutation index
delta1 = (2 * u + (1 - 2 * u) * (1 - delta)**(eta_m + 1))**(1 / (eta_m + 1)) - 1
delta2 = 1 - (2 * (1 - u) + 2 * (u - 0.5) * (1 - delta)**(eta_m + 1))**(1 / (eta_m + 1))

deltaq = np.where(u < 0.5, delta1, delta2)
deltamax = np.array([b[1] - b[0] for b in self.bounds])
new_values = individual.values + mutation_mask * deltaq * deltamax # Apply mutation mask here

# Clipping the values to ensure they remain within bounds
new_values = np.clip(new_values, [b[0] for b in self.bounds], [b[1] for b in self.bounds])

# Apply the new values if mutation occurred
if np.any(mutation_mask):
individual.values = new_values
individual.recalculate_objective_values()



# Function to check if one individual dominates another
def dominates(self, individual1, individual2):
# Check if individual1 dominates individual2
return all(individual1.objective_values <= individual2.objective_values) and any(individual1.objective_values < individual2.objective_values)

#Function to perform the non-dominated sorting
def non_dominated_sorting(self):
fronts = []
front = []
for individual in self.population:
individual.dominated = 0
individual.dominating = []
for other_individual in self.population:
if self.dominates(individual, other_individual):
individual.dominating.append(other_individual)
elif self.dominates(other_individual, individual):
individual.dominated += 1
if individual.dominated == 0:
individual.rank = 0
front.append(individual)
fronts.append(front)
i = 0
while len(fronts) > 0:
next_front = []
for individual in fronts:
for other_individual in individual.dominating:
other_individual.dominated -= 1
if other_individual.dominated == 0:
other_individual.rank = i + 1
next_front.append(other_individual)
i += 1
fronts.append(next_front)
if len(fronts[-1]) == 0:
fronts.pop(-1)
self.fronts = fronts


#Function to perform the crowding distance assignment
def crowding_distance_assignment(self, front):
for individual in front:
individual.crowding_distance = 0.0
for i in range(len(self.objectives)):
front.sort(key = lambda x: x.objective_values)
front[0].crowding_distance = float('inf')
front[-1].crowding_distance = float('inf')
for j in range(1, len(front) - 1):
front[j].crowding_distance += (front[j + 1].objective_values - front[j - 1].objective_values) / (front[-1].objective_values - front[0].objective_values + 1e-6)


#Function to perform the crowding distance sorting
def crowding_distance_sorting(self):
for front in self.fronts:
self.crowding_distance_assignment(front)
front.sort(key = lambda x: x.crowding_distance, reverse = True)

#Function to perform the NSGA-II algorithm
def run(self):
self.non_dominated_sorting()
self.crowding_distance_sorting()

while self.current_generation < self.generations:
self.offspring_population = []
for i in range(self.population_size//2):
parent1 = self.tournament_selection()
parent2 = self.tournament_selection()
child1, child2 = self.crossover(parent1, parent2)
self.mutation(child1)
self.mutation(child2)
self.offspring_population.append(child1)
self.offspring_population.append(child2)


self.population.extend(self.offspring_population)
self.non_dominated_sorting()
self.crowding_distance_sorting()
self.population = []
for front in self.fronts:
if len(self.population) + len(front) <= self.population_size:
self.population.extend(front)
else:
self.crowding_distance_assignment(front)
front.sort(key = lambda x: x.crowding_distance, reverse = True)
self.population.extend(front[:self.population_size - len(self.population)])
self.current_generation += 1
return self.fronts[0]

def generate_dtlz1_objectives(M, n):
"""
Generates a list of DTLZ1 objective functions.

Parameters:
M (int): Number of objectives.
n (int): Number of dimensions (variables).

Returns:
list: A list of lambda functions representing the DTLZ1 objective functions.
"""
k = n - M + 1 # Number of variables involved in the g function

def g(x):
xM = np.array(x[-k:])
return 100 * (k + np.sum((xM - 0.5)**2 - np.cos(20 * np.pi * (xM - 0.5))))

objectives = []
for i in range(1, M + 1):
if i == 1:
# For the first objective, use all except the last element
objectives.append(lambda x: 0.5 * np.prod(np.array(x)[:M-1]) * (1 + g(x)))
elif i < M:
# For middle objectives, exclude progressively more elements
objectives.append(lambda x, i=i: 0.5 * np.prod(np.array(x)[:M-i]) * (1 - x[M-i]) * (1 + g(x)))
else:
# For the last objective, use the negation of the first element
objectives.append(lambda x, i=i: 0.5 * (1 - x[0]) * (1 + g(x)))

return objectives

def run1():
# Example of using the NSGA-II algorithm
population_size = 100
generations = 300
dimensions = 7
mutation_rate = 1/dimensions
crossover_rate = 1
tournament_size = 2
M = 3 # Number of objectives
objectives = generate_dtlz1_objectives(M, dimensions)

#Example of using the NSGA-II algorithm
bounds = [(0, 1) for _ in range(dimensions)]
nsga = NSGAII(population_size, generations, mutation_rate, crossover_rate, objectives, bounds, tournament_size, dimensions)
front = nsga.run()
print(front[0].values)
#plot 3d
fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.scatter([individual.objective_values[0] for individual in front], [individual.objective_values[1] for individual in front], [individual.objective_values[2] for individual in front])
ax.set_xlabel('f1')
ax.set_ylabel('f2')
ax.set_zlabel('f3')
plt.show()

run1()
</code></pre>
 

Online statistics

Members online
0
Guests online
2
Total visitors
2
Top