Skip to content

Instantly share code, notes, and snippets.

@Chrisvandr
Created July 5, 2017 14:06
Show Gist options
  • Save Chrisvandr/6245a16f674b7343fd4fe433ec697e9f to your computer and use it in GitHub Desktop.
Save Chrisvandr/6245a16f674b7343fd4fe433ec697e9f to your computer and use it in GitHub Desktop.
import numpy as np
import matplotlib.pyplot as plt
import math
import array
import random
import pickle
import pandas as pd
from deap import base
from deap import benchmarks
from deap import creator
from deap import tools
from deap import algorithms
from sklearn.externals import joblib
from decimal import *
def elements(df_inputs, lower_bound, upper_bound):
elements = [random.uniform(a, b) for a, b in zip(lower_bound, upper_bound)]
#elements = []
# for x in range(df_coefs.shape[1]):
# elements.append(random.gauss(df_inputs.iloc[0][x], df_inputs.iloc[0][x]*(1/10))) # set to 10% variation currently... #TODO change
return elements
## Function for calculating the end uses based on the coefficients and intercepts
def calculate(individual):
outputs = []
for y in range(len(targets)):
outputs.append(sum(np.array(coefficients[y])* individual)+intercepts[y])
#outputs.append(output)
outputs = np.array(outputs)
return outputs
# Load the pickle surrogate model
# def calculate(individual):
# model = joblib.load(DATAPATH_MODEL_REAL + 'rr_' + TIME_STEP + '_model.pkl')
# individual = np.array(individual).reshape((1,-1))
# prediction = model.predict(individual)[0]
# prediction = prediction[:NO_OF_OBJECTIVES]
#
# return prediction
def evaluate(individual):
diff = []
prediction = calculate(individual)
for y in range(len(targets)):
output = math.sqrt((targets[y] - prediction[y]) ** 2)
#output = abs(targets[y] - prediction[y])
diff.append(output)
diff = diff[:NO_OF_OBJECTIVES]
#diff = [float(i)/sum(diff) for i in diff] #normalize output
return tuple(diff)
# PLOT: Plot energy use for the objectives with increasing generation
def plot_objectives():
result = []
for i, v in enumerate(best_inds):
prediction = calculate(best_inds[i])
result.append(prediction)
df_result = pd.DataFrame(result, columns=cols_objectives)
ax = df_result.plot(title='Prediction of the best individual per generation', color=colors)
targets_dup = ([targets, ]) * len(result) # duplicate targets list
df_target = pd.DataFrame(targets_dup, columns=cols_objectives)
df_target.plot(ax=ax, style='--', color=colors)
def main():
random.seed(20)
stats = tools.Statistics(lambda ind: ind.fitness.values)
stats.register('min', np.min, axis=0)
stats.register('max', np.max, axis=0)
stats.register('std', np.std, axis=0)
stats.register('avg', np.average, axis=0)
logbook = tools.Logbook()
logbook.header = "gen", "avg", "inputs" #, "max", "avg" #'inputs', 'std', 'avg', 'evals'
# Create an initial population of size n.
pop = toolbox.population(n=POPULATION_SIZE)
hof = tools.ParetoFront()
# Evaluate the individuals with an invalid fitness
invalid_ind = [ind for ind in pop if not ind.fitness.valid]
fitnesses = toolbox.map(toolbox.evaluate, invalid_ind)
for ind, fit in zip(invalid_ind, fitnesses):
ind.fitness.values = fit
# This is just to assign the crowding distance to the individuals
pop = toolbox.select(pop, POPULATION_SIZE) # no actual selection is done
best_inds, best_inds_fitness = [], []
record = stats.compile(pop)
logbook.record(**record) # inputs=best_inds_fitness
hof.update(pop)
# Begin the generational process
for gen in range(1, NGEN):
# Vary the population
offspring = tools.selTournamentDCD(pop, len(pop)) # only works with "select" to NSGA-II
#offspring = toolbox.select(pop, len(pop))
offspring = [toolbox.clone(ind) for ind in offspring] # Clone the selected individuals
for ind1, ind2 in zip(offspring[::2], offspring[1::2]): # Crossover and mutate offspring
#print('inds', ind1, ind2)
if random.random() <= CXPB:
toolbox.mate(ind1, ind2) # crossover randomly chosen individual within the population
# toolbox.mutate(ind1)
# toolbox.mutate(ind2)
del ind1.fitness.values
del ind2.fitness.values
for mutant in offspring:
if random.random() <= MUTPB: # which offspring are mutated
#mutate a subset of the attributes (variables) in the individuals, based on a percentage of the pop
#for i in range(round(len(pop[0])*PERC_ATTR_TO_MUTATE)):
toolbox.mutate(mutant) # mutates several attributes in the individual
del mutant.fitness.values
# Evaluate the individuals with an invalid fitness
invalid_ind = [ind for ind in offspring if not ind.fitness.valid]
#print(len(invalid_ind))
fitnesses = toolbox.map(toolbox.evaluate, invalid_ind)
for ind, fit in zip(invalid_ind, fitnesses):
ind.fitness.values = fit
# Select the next generation population
pop = toolbox.select(pop + offspring, len(offspring))
best_ind = tools.selBest(pop, 1)[0]
best_inds.append(best_ind) # add the best individual for each generation
best_inds_fitness.append(best_ind.fitness.values)
record = stats.compile(pop)
logbook.record(gen=gen, inputs=[int(e) for e in best_ind.fitness.values], **record)
hof.update(pop)
if gen % 20 == 0:
print(gen, 'best_ind_pred', [int(e) for e in calculate(best_inds[-1]).tolist()], 'targets', targets[:NO_OF_OBJECTIVES], 'fitness', [int(e) for e in best_ind.fitness.values], 'sum abs diff', [int(sum(targets)-sum(calculate(best_inds[-1])))], )
# with open("logbook.pkl", "wb") as lb_file: pickle.dump(logbook, lb_file)
return pop, hof, logbook, best_inds, best_inds_fitness
if __name__ == '__main__':
NGEN = 500 # is the number of generation for which the evolution runs # For selTournamentDCD, NGEN has to be a multiple of four the population size
POPULATION_SIZE = 48 # no of individuals/samples
CXPB = 0.9 # is the probability with which two individuals are crossed
MUTPB = 0.6 # probability of mutating the offspring
INDPB = 0.2 # Independent probability of each attribute to be mutated
SIGMA_VARIATION = 20 #percentage of variation for the input variables
BUILDING_ABR = 'CH'
TIME_STEP = 'month' #'year', 'month'
NO_OF_OBJECTIVES = 8 # number of objectives to optimise
if BUILDING_ABR in {'CH'}:
#DATAPATH_MODEL_REAL = 'C:/EngD_hardrive/UCL_DemandLogic/01_CentralHouse_Project/Run1000/'
END_USES = False
inputs = [3.09, 4.15, 3.37, 3.24, 3.33, 2.57, 3.09, 4.15, 3.68, 5.12, 3.28, 3.05, 3.54, 3.89, 3.4, 5.63,
2.36, 2.97, 3.78, 3.71, 2.86, 2.3, 0.01, 7.29, 16.03, 12.72, 14.7, 10.84, 7.77, 8.16, 29.12,
9.65, 12.39, 12.13, 8.24, 15.76, 15.44, 16.94, 12.58, 4.81, 3.11, 8.32, 6.73, 6.73, 6.48, 6.4,
7.33, 23.19, 5.86, 9.58, 4.5, 5.33, 15.9, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0,
0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.68, 0.82, 3.83, 12.3, 11.34, 3.02, 11.86,
6.62, 12.24, 18.12, 51.95, 15.0, 7.45, 54.68, 3.6, 0.0, 0.1, 0.16, 0.9, 0.91, 2.09, 1.98,
0.03, 0.04, 0.31, 0.15, 0.15, 0.18, 0.19, 0.14, 2.97, 192.55, 208.06, 22.5, 1.5, 1.38, 21.69,
29.76, 18.13, 3.0, 3.71, 0.0, 0.91, 0.0, 3.0, 17.25, 0.84, 1.73, 0.02, 3.0, 2.51, 1.05, 0.0,
0.13, 2.0, 11.51, 1.39, 0.83, 0.14, 2.0, 81.02, 41.09, 14.23, 4.05, 1.63, 1.11, 31.14, 29.26,
14.63, 14.63, 22.87, 14.84, 7.68, 7.42]
if TIME_STEP == 'year':
targets = [284545, 265146, 44914] # sum of 8 months of data for 3 energy end-uses, 'Systems', 'L&P', 'Gas' x 8 months
cols_objectives = ['Systems', 'L&P', 'Gas']
coefficients = [
[388, -8915, -1537, -560, 996, -781, -801, -372, -2081, -978, -1202, -1579, -3003, -11389, -3506,
-1760, -1764, -4571, -1685, 464, -1639, 633, 43, 312, -26, 147, 16, 74, -310, 118, -99, 492, -143,
-188, -65, 95, 26, 170, -3, -1012, -1481, 73, -9, -84, -559, -483, -21, -74, 56, -176, -352, 359,
-125, -18, 29, 15, 22, 10, 75, 39, 157, 18, -4, 7, -17, -15, 4, 0, 0, 0, -1, 0, -1, -4349, -1773,
1140, -314, -9, 473, -309, -403, 433, 767, 47, -175, 155, 183, 923, 0, 950, 446, 6741, -2753, 8640,
703, 123, 44, 14, -168, 157, -119, -423, -875, 2665, -20, 38, -6, -14118, 2027, 958, 881, 985, 776,
2962, 0, -7586, 0, 978, -51, 2180, -2788, 640, -1200, 1004, 170, 0, -885, -3360, -129, -2578, 6189,
625, -5224, -729, -505, -2410, 5358, 3222, 3207, -6103, -133, -2486, -11, 5314, 1027, 7274, 6761],
[509, 188, -162, -967, -613, -449, 529, -16, 313, -522, 886, 445, 53, -485, 195, -459, -648, -727,
1270, -304, -927, 205, -154, 1800, 229, 584, 332, 294, 393, 650, 447, 4480, -58, -252, 6, -31, -59,
178, -382, -689, -1864, -369, 117, 61, -8, -413, -365, 1, -360, -449, -264, -126, -145, -25, 21,
-4, -6, 19, -17, 19, -10, 22, -2, 7, 2, -24, 2, 0, 2, 0, 2, 2, 6, -3291, -397, 2920, 193, 724,
2830, 465, 733, 1067, 6406, 353, -263, 406, 452, 757, 0, 487, 241, 6729, -1053, 3050, -3256, 2,
-50, 1333, -302, -296, 1276, -1338, 1203, 2729, -46, 84, 0, -287, 401, 412, 297, -3, 108, 2009, 0,
-2220, 0, 1196, 914, -2418, -3044, 308, 748, 1930, 494, 0, 406, -761, 216, 4826, 2972, 519, 1267,
258, 178, 984, -676, -305, -189, 7131, 1204, 1255, 746, 1003, -1051, -178, -2327],
[-161, -137, 1013, 292, 144, -123, -6, 69, -516, -298, 13, 68, 307, -205, -475, -175, 572, -155,
-851, 326, 832, -160, 83, -508, 51, -31, -341, 44, 75, -279, -59, -67, -35, 15, -9, 36, -72, -30,
21, -289, -92, 134, 223, -53, 82, -188, -157, -83, -44, 30, -32, 95, 15, 69, 16, 2, 51, -3, 11, 0,
21, 4, 32, 0, 0, 10, -10, 0, 0, 0, 0, 0, 0, -36177, 1210, -547, -83, -44, -300, 11, -47, 75, -114,
-4, -61, 229, 14, -52, 0, 36, 549, -683, -742, 223, 1644, 229, 60, -1013, -397, 629, -300, 324,
-322, 6145, -22, -23, 2, -975, -4, 136, -224, 14, -199, 21, 0, 22, 0, 166, 279, -1036, -1019, 264,
-133, 434, -1007, 0, -312, -1111, 39, -490, 939, 230, 775, -12, -8, -361, 1994, 919, 746, -641,
337, -88, 181, -778, -354, 33, 216]]
intercepts = [363070, -274327, 93024]
elif TIME_STEP == 'month':
targets = [83785, 74560, 83954, 78379, 57150, 62704, 72635, 81437]
targets = targets[:NO_OF_OBJECTIVES]
cols_objectives = [1, 2, 3, 4, 9, 10, 11, 12]
coefficients = [
[29, -1934, 243, -419, 214, 46, -127, 30, -480, -334, -333, -321, 141, -3281, -380, -474, -144,
-835, -142, 30, -146, -182, 76, 194, 11, 75, -64, 14, -104, 11, 13, 352, -4, -2, -52, 8, -69, 53,
-34, -406, -382, -201, 111, -39, -141, -190, -5, -26, -39, 27, -5, 120, -7, 28, 6, -12, 21, -22,
68, 22, 140, 13, 25, -15, -25, -10, -15, 0, 0, 0, 0, 0, 2, -10976, -1379, 47, 34, -6, 251, 52,
-195, 53, 483, 23, -49, 122, 78, 368, 0, 25, -777, 1367, -154, 1752, 12, 299, 54, -1324, 312, 723,
-330, -1264, -128, 4357, -6, -5, -1, -2815, 517, 92, 89, 3, 487, 562, 0, -2779, 0, -392, 113,
-1020, 511, 778, 381, -264, 1231, 0, -208, 279, 390, -904, -459, 28, 307, -15, -10, 770, -800, -91,
251, 892, -1039, 1732, -434, -1097, 178, 1244, -14],
[4, -1548, 25, -407, 142, 10, -130, -138, -442, -155, -195, -107, -130, -2073, -363, -273, -183,
-745, -112, 70, -25, -14, -35, 274, 31, 85, 18, 29, -26, 33, 53, 630, 23, 1, -5, -23, -28, 25, -7,
-235, -146, -23, 64, -143, -77, -70, -1, -2, 17, -45, -37, 70, 26, 29, 1, -3, 22, -8, 15, 13, 89,
15, 10, -4, -14, -4, -7, 0, 0, 0, 0, 0, 2, -7279, 135, 91, 14, 39, 72, -8, -36, 101, 838, 32, -10,
98, 84, 132, 0, 326, -607, 2045, -446, 1772, 420, 111, -103, 236, 354, 227, 564, -242, -330, 2437,
-4, 1, -3, -2167, 195, 25, 17, 106, 271, 431, 0, -1698, 0, 30, 120, -412, -311, 597, 127, 329, 82,
0, -134, 223, 195, 0, -407, -106, 243, -125, -86, 85, 1111, 918, 1160, 140, -38, -645, 47, -40,
968, 1478, 632],
[-34, -1676, -175, -246, -42, -22, -297, -205, -438, -16, -420, -47, -334, -1516, -654, -179, -262,
-543, -209, 82, -176, 8, -107, 320, 32, 119, 49, 57, 40, 9, 92, 911, 11, -7, 3, -21, -2, 4, -5,
-169, -51, -45, 61, -85, -101, -31, 24, 14, -13, -65, -50, 23, 39, 17, 0, 0, 20, -1, -6, 8, 42, 13,
4, 4, -6, -5, 1, 0, 0, 0, 0, 0, 3, -4962, -47, 268, 14, 81, 132, -3, 68, 155, 1205, 41, 7, 97, 100,
67, 0, 266, -322, 1014, -369, 1666, 176, -26, -113, 705, 217, -244, 1048, -27, -382, 791, -1, 5,
-3, -1953, 125, -19, 8, 133, 215, 227, 0, -996, 0, 248, 144, 45, -738, 297, -35, 656, -592, 0, 0,
320, 236, 69, -698, -178, 31, -59, -40, 270, 150, 406, 694, 200, 305, -1291, 187, 1033, 884, 220,
908],
[-52, -1563, -269, -131, -103, -15, -373, -142, -435, 3, -582, -3, -464, -1186, -749, -103, -337,
-340, -298, 85, -239, 19, -99, 288, 24, 107, 63, 56, 41, 2, 83, 832, 0, -8, 14, -13, 6, -2, -8,
-114, -7, -36, 54, -39, -83, -4, 25, 11, -19, -57, -43, -5, 49, 2, 0, 0, 10, 1, -6, 4, 6, 11, -1,
5, -1, -5, 5, 0, 0, 0, 0, 0, 2, -2901, -105, 317, 15, 80, 143, -4, 71, 151, 1144, 42, 14, 75, 96,
23, 0, 206, -180, 700, -366, 1475, 31, -88, -96, 860, 64, -308, 1055, 90, -264, 12, 2, 6, -2,
-1710, 54, -37, 7, 133, 146, 91, 0, -496, 0, 260, 134, 144, -759, 127, -104, 685, -721, 0, 155,
336, 184, 316, -761, -168, -75, -21, -14, 274, -248, 106, 323, 141, 424, -1108, 245, 1148, 728,
-101, 547],
[-71, -230, -278, -54, -113, 1, -259, -90, -268, 14, -443, 21, -459, -230, -542, -8, -262, -87,
-219, 12, -132, 39, -105, 120, 10, 52, 34, 36, 18, -5, 42, 376, 9, -9, 17, -10, 2, -4, -1, -11, 39,
-10, 32, -24, -27, 8, 6, 5, -1, -36, -21, -21, 17, -4, 0, 2, 2, 5, -5, 0, -6, 4, -3, 2, 3, 1, 6, 0,
0, 0, 0, 0, 2, -523, 143, 172, -2, 32, 55, -7, 44, 80, 492, 44, 10, 37, 83, -12, 0, 133, -207, 436,
-37, 172, 16, -65, -58, 826, -93, -111, 478, 213, -40, -287, 0, 4, -2, -1546, -23, 44, -18, 132,
105, -13, 0, 1, 0, 103, 80, 159, -242, 24, -110, 504, -404, 0, 17, 172, 126, 377, -373, -107, -172,
5, 4, 54, -251, -30, 64, 35, 253, -768, 161, 899, 167, -220, 171],
[-60, -358, -108, -156, -19, 29, -199, -126, -311, -5, -299, -1, -426, -910, -528, -68, -218, -368,
-157, 61, -73, 16, -81, 252, 35, 81, 39, 37, 29, 22, 71, 670, 10, -5, 14, -21, -3, 11, -11, -106,
-3, -19, 50, -61, -72, -21, 14, 2, 3, -36, -54, 0, 20, 6, 0, -1, 7, 0, -7, 1, 12, 9, 2, 1, -2, -2,
1, 0, 0, 0, 0, 0, 3, -1876, 58, 179, 13, 50, 99, -10, 76, 124, 876, 38, 1, 83, 86, 33, 0, 289,
-346, 971, -128, 650, 39, -27, -104, 456, 61, -13, 706, 14, -169, 450, -1, 4, -2, -1849, 52, 31,
-11, 103, 132, 165, 0, -633, 0, -22, 76, 81, -148, 252, -32, 467, -237, 0, 18, 258, 160, 388, -548,
-157, -129, -55, -38, 57, 388, 433, 617, 98, 245, -1087, 145, 920, 660, 357, 530],
[-14, -758, 10, -302, 73, 11, -171, -227, -400, -60, -171, -92, -273, -1651, -509, -191, -205, -692,
-76, 76, -20, -18, -62, 285, 42, 95, 22, 16, 29, 32, 78, 744, 19, -1, 0, -25, -10, 17, -8, -201,
-82, -32, 64, -120, -101, -60, 27, 11, 1, -48, -59, 50, 24, 24, 0, 0, 20, -4, -6, 9, 64, 13, 9, -1,
-8, -5, -5, 0, 0, 0, 0, 0, 3, -5145, 14, 139, 9, 41, 93, -32, 71, 121, 956, 32, -4, 113, 84, 92, 0,
331, -528, 1558, -318, 1247, 141, 75, -107, 340, 280, -48, 782, -162, -468, 1623, -3, 4, -2, -1954,
150, 26, 9, 58, 241, 359, 0, -1379, 0, 13, 108, -127, -302, 478, 84, 374, -38, 0, -231, 305, 236,
8, -651, -189, 74, -107, -74, 38, 1021, 908, 1198, 43, 115, -1174, 81, 542, 1076, 903, 1042],
[22, -1395, 93, -360, 173, 14, -84, -174, -319, -111, -44, -102, -44, -2202, -331, -225, -139, -786,
-20, 114, -5, -35, 43, 228, 32, 61, 0, 17, -19, 46, 29, 462, 15, 1, -6, -26, -27, 30, -19, -207,
-165, -45, 35, -133, -59, -105, 3, -7, 19, -31, -47, 67, 21, 24, 2, -2, 16, -11, 7, 7, 76, 10, 7,
-5, -13, -5, -8, 0, 0, 0, 0, 0, 2, -6950, -57, 17, 2, 11, 87, -35, -13, 71, 646, 34, -21, 87, 83,
114, 0, 277, -611, 2163, -127, 1961, 371, 150, -93, -174, 380, 198, 352, -413, -473, 2711, -3, 1,
-2, -1819, 238, 37, 13, 58, 162, 388, 0, -1344, 0, 18, 108, -370, -281, 578, 154, 195, 269, 0,
-231, 236, 218, -6, -422, -106, 498, -137, -95, 166, 1198, 922, 1122, -104, 4, -426, 5, -530, 1311,
1366, 755]]
intercepts = [83257, 43057, 14637, 6245, 2808, 1978, 24104, 50176]
sigmas = [i* (SIGMA_VARIATION / 100) for i in inputs]
bound = [float(i) * 3 for i in sigmas] # min/max is 3*sigma
lower_bound = [i - j for i, j in zip(inputs, bound)]
upper_bound = [i + j for i, j in zip(inputs, bound)]
print('weights normalized to target', tuple([-float(i) / sum(targets) for i in targets]))
#print('equal weights', tuple([-1 for i in targets]))
creator.create('Fitness', base.Fitness, weights=tuple([-float(i) / sum(targets) for i in targets]))
creator.create('Individual', array.array, typecode='d', fitness=creator.Fitness) # set or list??
toolbox = base.Toolbox()
## using custom function for input generation
toolbox.register('expr', elements, inputs, lower_bound, upper_bound)
toolbox.register('individual', tools.initIterate, creator.Individual, toolbox.expr)
toolbox.register('population', tools.initRepeat, list, toolbox.individual)
# TODO try different selection criteria/mutate/mate
toolbox.register('mate', tools.cxTwoPoint)
#toolbox.register("mate", tools.cxUniform, indpb=INDPB)
#toolbox.register("mate", tools.cxSimulatedBinaryBounded, low=lower_bound, up=upper_bound, eta=20.0)
# toolbox.register('mutate', tools.mutFlipBit, indpb=INDPB)
# toolbox.register("mutate", tools.mutPolynomialBounded, low=lower_bound, up=upper_bound, eta=20.0, indpb=INDPB)
toolbox.register('mutate', tools.mutGaussian, mu=0, sigma=sigmas, indpb=INDPB) # sigmas are set to a sequence
toolbox.register('select', tools.selNSGA2)
#toolbox.register('select', tools.selSPEA2)
toolbox.register('evaluate', evaluate) # add the evaluation function
#toolbox.register("evaluate", benchmarks.zdt1)
pop, hof, logbook, best_inds, best_inds_fitness = main()
logbook.chapters["fitness"].header = "min", "max", "avg"
gen = logbook.select("gen")
fit_mins = logbook.select("min")
fit_avgs = logbook.select("avg")
print('HOF', len(hof))
print('cols_targets', cols_objectives)
print('best individual', [int(e) for e in calculate(best_inds[-1]).tolist()])
print('targets', targets)
print('absolute diff', [i - j for i, j in zip([int(e) for e in calculate(best_inds[-1]).tolist()], targets)])
colors = ['#1f77b4', '#ff7f0e', '#2ca02c', '#8c564b', '#d62728', '#9467bd', '#aec7e8', '#ffbb78', '#98df8a',
'#c49c94', '#ff9896', '#c5b0d5', '#1f77b4', '#ff7f0e', '#2ca02c', '#8c564b', '#d62728', '#9467bd',
'#aec7e8', '#ffbb78', '#98df8a', '#c49c94', '#ff9896', '#c5b0d5']
plot_objectives() #shows each end use and how they change absolutely over the generations
plt.show()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment