optimization of resilience indicator

scenariosAleatoires/MCMC_benchmark_indicators.py

@@ -17,16 +17,16 @@ if __name__ == '__main__':
             ))
 
 '''
-proximity: 0.0142 s
-resilience: 1.5464 s
-productivity: 0.0186 s
-biodiversity: 0.0181 s
-social: 0.1293 s
+proximity: 0.0146 s
+resilience: 0.0901 s
+productivity: 0.0099 s
+biodiversity: 0.0178 s
+social: 0.1336 s
 
-reallocate: 0.0453 s
-proximity: 0.0174 s
-resilience: 1.1516 s
-productivity: 0.0122 s
-biodiversity: 0.0169 s
-social: 0.2007 s
+reallocate: 0.0469 s
+proximity: 0.0167 s
+resilience: 0.0874 s
+productivity: 0.0117 s
+biodiversity: 0.0176 s
+social: 0.2209 s
 '''

scenariosAleatoires/resilience_list.py

 #!/usr/bin/python3
 # -*- coding: utf-8 -*-
 import geopandas as gpd
+import pandas as pd
+import numpy as np
 import math
 from Indicator import Indicator
 from overrides import overrides
@@ -14,13 +16,14 @@ class Resilience(Indicator):
 		:param grille: Grid divising the PAT in squares
 		'''
 		self.grid = gpd.GeoDataFrame.from_file(resilience_config)
-		intersection = gpd.sjoin(
+		# caching intersection between patches and grid
+		self.intersection = gpd.sjoin(
 			self.grid[['id','geometry']],
-		 	initial_patches[['geometry','SURF_PARC']],
+		 	initial_patches[['geometry','SURF_PARC','ID_PARCEL']],
 			how='inner', op='intersects')
-		self.grouped_by_cell = intersection.groupby('id')
-		self.size = len(self.grouped_by_cell)
-		self.surfaces = intersection[['id','SURF_PARC']].groupby('id').sum()
+		# caching surface of patches in each cell
+		self.intersection['surf_cell'] = self.intersection[['id','SURF_PARC']].groupby('id').transform(sum)
+		self.size = len(self.intersection.groupby('id'))
 
 	# The parameters are : the grid layer, a list corresponding to a scenario and a boolean value True or False signifying if we want to see the printing lines or not
 	@overrides
@@ -37,16 +40,13 @@ class Resilience(Indicator):
 		pi: a specific kind of culture      ∑: for each kind of culture
 
 		:param patches: The scenario to analyse as a list
-		:param affichage: True if we want the display in the console, False is the other case
 		'''
-		resilience = 0
-		for cell,group_by_cell in self.grouped_by_cell:
-			surf_tot = self.surfaces.loc[cell].values[0] # 0.085ms
-			intersection = patches[['ID_PARCEL','cultgeopat','SURF_PARC']].loc[group_by_cell['index_right']] # 4.26ms
-			prop_cult = intersection.groupby('cultgeopat')['SURF_PARC'].sum() / surf_tot # 0.5ms
-			res_cell = - (prop_cult * prop_cult.apply(math.log2)).sum() # 0.22ms
-			resilience += res_cell
-		return resilience / self.size
+		# merging given patches with table of intersection patches/grid
+		intersection2 = pd.merge(self.intersection, patches[['ID_PARCEL','cultgeopat','SURF_PARC']],on='ID_PARCEL',how='left')
+		# building pivot table on grid cell and cultgeopat
+		pivot = pd.pivot_table(intersection2,values=['SURF_PARC_x','surf_cell'],index=['id','cultgeopat'],aggfunc={'SURF_PARC_x':np.sum,'surf_cell':'first'})
+		pivot['res'] = pivot['SURF_PARC_x']/pivot['surf_cell']
+		return - (pivot['res'] * pivot['res'].apply(math.log2)).sum() / self.size
 
 if __name__ == '__main__':
 	from patutils import load_pat_patches