social.py 5.04 KB
Newer Older
1
2
#!/usr/bin/python3
# -*- coding: utf-8 -*-
3
import pandas as pd
4
5
import numpy as np
import os
6

7
class Social:
8
9
	'''
	This function is used to calculate the social indice.
10
11
	The idea is to propose a transition matrix for each type of culture and use it to visualize
	 the conversion cost of the total conversion between the initial scenario and the old one.
12

13
14
15
	Each transition has a value between 0 and 1.
	These values were defined during reunions and with the review of our experts.
	'''
16

17
	def __init__(self, initial_patches, patches_md5sum, cost_matrix_filename, cost_matrix_md5sum, bdv_threshold):
18
		'''
19
20
		:param initial_patches Initial state of the patches (attributes cultugeopat and SURF_PARC required)
		:param patches_md5sum md5 checksum of the input file of patches
21
22
23
24
		:param cost_matrix_filename input file containing the matrix in CSV format (separator: \t)
		:param bdv_threshold proportion of a cultural type needed for considering that this cultural type
		is suffisantly present in the locality for ignoring the conversion cost.
		'''
25
		self.patches_md5sum = patches_md5sum
26
		self.matrice_transition = pd.read_csv(cost_matrix_filename, sep='\t',index_col=0)
27
		self.cost_matrix_md5sum = cost_matrix_md5sum
28
		self.bdv_threshold = bdv_threshold
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69

	def compute_patches_transition_cost(self, initial_patches):
		costs=[]
		for bdv,grouped in initial_patches.groupby('Bdv'):
			grouped=grouped.copy()
			surfaces_initial = grouped.groupby('cultgeopat')['SURF_PARC'].sum()
			cult_threshold = surfaces_initial > surfaces_initial.sum() * self.bdv_threshold
		for culture in self.matrice_transition.columns.get_values():
			for row in initial_patches.itertuples():
				cost = 0.0
				# patches that have different farming that others in same Bdv
				if culture not in cult_threshold or not cult_threshold[culture]:
					# patches that have different farming that others of the same farmer
					if culture not in initial_patches[initial_patches['id_expl']==row.id_expl]['cultgeopat']:
						cost = self.matrice_transition.loc[row.cultgeopat,culture]
				costs.append([row.ID_PARCEL, row.cultgeopat, culture, cost])
		self.patches_transition_cost = pd.DataFrame(costs, columns=['ID_PARCEL','initialcult','cultgeopat','cost'])

	def save_patches_transition_cost(self, filename):
		'''
		Save the matrix of costs computed for each patch and initial parameters in a file (npz format).
		'''
		m = self.patches_transition_cost.values
		np.savez(filename, bdv_threshold=self.bdv_threshold, patches_md5sum=self.patches_md5sum, cost_matrix_md5sum=self.cost_matrix_md5sum,
			data=m.data)

	def load_patches_transition_cost(self, filename):
		'''
		Load the matrix of costs and initial parameters from a file (npz format).

		The correspondance will be checked between
		 class attributes and metadata stored in the file
		'''
		loader = np.load(filename)
		error_message = ('Matrix metadata {0} does\'nt correspond to metadata given at indicator construction'
			'\n\t{0} in metadata: {1} while asked: {2}'
			'\n\tmatrix filename: ') + filename
		for k1,k2 in [[self.bdv_threshold,'bdv_threshold'],[self.patches_md5sum,'patches_md5sum'],[self.cost_matrix_md5sum,'cost_matrix_md5sum']]:
			if k2 not in loader or k1!=loader[k2]:
				raise ValueError(error_message.format(k2, dict(loader.items()).get(k2, 'NA'), k1))
		self.patches_transition_cost = pd.DataFrame(loader['data'], columns=['ID_PARCEL','initialcult','cultgeopat','cost'])
70

71
	def compute_indicator(self, scenario_patches):
72
73
74
75
		'''
		:param scenario: The scenario to analyse
		'''
		cost = 0
76
77
		patches_altered = scenario_patches[scenario_patches['cultgeopat'] != scenario_patches['init_cult']]
		cost = pd.merge(patches_altered, self.patches_transition_cost,on=['ID_PARCEL','cultgeopat'], how='left')['cost'].sum()
78
		return 1/(1+cost)
79
80
81

if __name__ == '__main__':
	import geopandas as gpd
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
	import time
	import os
	import hashlib
	from MCMC import md5sum
	# loading and initializing patches data
	patches_filename = "../output/PAT_patches/PAT_patches.shp"
	matrix_filename = '../resources/matrice_transition.csv'
	patches = gpd.GeoDataFrame.from_file(patches_filename, encoding='utf-8')
	patches = patches[patches['cultgeopat']!='Non Considérée']
	patches['init_cult'] = patches['cultgeopat']
	# init indicator
	social = Social(patches, md5sum(patches_filename), matrix_filename, md5sum(matrix_filename), 0.1)
	costs_filename = '../output/patches_transition_costs.npz'
	if not os.path.isfile(costs_filename):
		social.compute_patches_transition_cost(patches)
		social.save_patches_transition_cost(costs_filename)
	social.load_patches_transition_cost(costs_filename)
	# computing init score
100
	scenario_patches = patches.copy()
101
	print(social.compute_indicator(scenario_patches))
102
103
104
	# small reallocation and computing score
	scenario_patches.loc[1,'cultgeopat']='Fruits et légumes'
	# should give 1/(1+1) the cost for one transition to Fruits et légumes
105
	print(social.compute_indicator(scenario_patches))
106
107
108
109
110
	# computing with more reallocation
	for N in [2000,5000,10000]:
		for i in range(N):
			scenario_patches.ix[i,'cultgeopat']='Protéagineux'
		print('{} {}'.format(N, social.compute_indicator(scenario_patches)))