Use existing Python/BluePyOpt code to draw the GA evaluated error surface.

[russelljjarvis]

Use existing Python/BluePyOpt code to draw the GA evaluated error surface.

Properly cite this code and borrow from BPO /notebook here:
https://github.com/BlueBrain/BluePyOpt/blob/master/examples/simplecell/simplecell-paperfig.ipynb
Cell 26

import numpy

gen_numbers = logs.select('gen')
min_fitness = numpy.array(logs.select('min'))
max_fitness = logs.select('max')
mean_fitness = numpy.array(logs.select('avg'))
std_fitness = numpy.array(logs.select('std'))

fig, ax = plt.subplots(3, figsize=(10, 10), facecolor='white')
fig_trip, ax_trip = plt.subplots(1, figsize=(10, 5), facecolor='white')

plot_count = len(responses)
for index, response in enumerate(responses[:plot_count]):
    color='lightblue'
    if index == plot_count - 1:
        color='blue'
        
    # best_ind_dict = cell_evaluator.param_dict(best_ind)
    # responses = cell_evaluator.run_protocols(cell_evaluator.fitness_protocols.values(), param_values=best_ind_dict)
    ax[0].plot(response['step1.soma.v']['time'], response['step1.soma.v']['voltage'], color=color, linewidth=1)
    axes = ax[1].plot(response['step2.soma.v']['time'], response['step2.soma.v']['voltage'], color=color, linewidth=1)
    # axes[0].set_rasterized(True)
    
ax[0].set_xlabel('Time (ms)')
ax[0].set_ylabel('Voltage (mV)')
ax[0].set_xlim(80, 200)
ax[1].set_xlabel('Time (ms)')
ax[1].set_ylabel('Voltage (mV)')
ax[1].set_xlim(80, 200)

# ax[0].legend()

std = std_fitness
mean = mean_fitness
minimum = min_fitness
stdminus = mean - std                                                           
stdplus = mean + std

ax[2].plot(                                                                      
    gen_numbers,                                                                
    mean,                                                                       
    color='black',                                                              
    linewidth=2,                                                                
    label='population average')                                                 

ax[2].fill_between(                                                              
    gen_numbers,                                                                
    stdminus,                                                                   
    stdplus,                                                                    
    color='lightgray',                                                          
    linewidth=2,                                                                
    label=r'population standard deviation')                                     

ax[2].plot(                                                                      
    gen_numbers,                                                                
    minimum,                                                                    
    color='red',                                                                
    linewidth=2,                                                                
    label='population minimum')                                                 

ax[2].set_xlim(min(gen_numbers) - 1, max(gen_numbers) + 1)                       
ax[2].set_xlabel('Generation #')                                                 
ax[2].set_ylabel('Sum of objectives')                                            
ax[2].set_ylim([0, max(stdplus)])                                                
ax[2].legend()                        

all_inds = hist.genealogy_history.values()

gnabars = numpy.array([ind[0] for ind in all_inds])
gkbars = numpy.array([ind[1] for ind in all_inds])
sums = numpy.array([ind.fitness.sum for ind in all_inds])
# psums = zip(gnabars, gkbars, sums)
zero_gnabars = gnabars[numpy.where(sums == 0)]
zero_gkbars = gkbars[numpy.where(sums == 0)]

# X = numpy.linspace(gnabar_param.bounds[0], gnabar_param.bounds[1], 150)
# Y = numpy.linspace(gkbar_param.bounds[0], gkbar_param.bounds[1], 150)
# X,Y = numpy.meshgrid(X,Y)
# import matplotlib
# import scipy.interpolate
# Z1 = scipy.interpolate.griddata(numpy.vstack((gnabars.flatten(), gkbars.flatten())).T, numpy.vstack(sums.flatten()), (X, Y), method='linear').reshape(X.shape)
# Z = matplotlib.mlab.griddata(gnabars, gkbars, sums, X, Y, interp='linear')
# Z1m = numpy.ma.masked_where(numpy.isnan(Z1),Z1)
# mesh_axes = ax[2].pcolormesh(X,Y,Z1m)

# plt.plot(gnabars, gkbars, '.', color='k')

trip_axis = ax_trip.tripcolor(gnabars,gkbars,sums+1,20,norm=matplotlib.colors.LogNorm())
# plt.tricontourf(gnabars,gkbars,sums,20)
# cbar_ax = fig.add_axes([0.85, 0.29, 0.025, 0.2])
fig_trip.colorbar(trip_axis, label='sum of objectives + 1')
ax_trip.set_xlabel('Parameter gnabar')
ax_trip.set_ylabel('Parameter gkbar')

plot_axis = ax_trip.plot(list(zero_gnabars), list(zero_gkbars), 'o', color='lightblue')
# ax[4].set_xlabel('Parameter gnabar')
# ax[4].set_ylabel('Parameter gkbar')
# ax[4].set_xlim(gnabar_param.bounds)
# ax[4].set_ylim(gkbar_param.bounds)
# fig.colorbar(trip_axis, ax=ax[3])

fig.tight_layout()
fig_trip.tight_layout()

# fig_trip.subplots_adjust(right=0.8, hspace=0.4)

# print(X, Y)
# print(Z)

fig.savefig('figures/simplecell_traces.eps')
fig_trip.savefig('figures/simplecell_trip.eps')

fig.show()
fig_trip.show()