varis = ['pricems', 'rooms', 'floor', 'needs_renovating', 'garden', 'terrace', 'new_dev', 'garage']

aves = db.groupby('GEOCODIGO')[varis].mean()
aves.info()

types = pd.get_dummies(db['metro_number'])
prop_types = types.join(db['GEOCODIGO'])\
                    .groupby('GEOCODIGO')\
                    .sum()
prop_types_pct = (prop_types * 100).div(prop_types)
prop_types_pct.info()

aves_props = aves.join(prop_types_pct)
#eliminar valores nulos
aves_props = aves_props.fillna(value=0)

db1 = pd.DataFrame(\
                 scale(aves_props), \
                 index=aves_props.index, \
                 columns=aves_props.columns)\
       #.rename(lambda x: str(int(x))              )


#zc = gpd.read_file(zc_link)
#zc.plot(color='green')
#sns.plt.show()

db1.info()

#zdb = db1.set_index('subdistrict_f').join(zc[['DESBDT', 'geometry']], on='DESBDT').dropna()

zdb = zc[['geometry', 'GEOCODIGO']].join(db1, on='GEOCODIGO')\
                                         .dropna()

km5 = cluster.KMeans(n_clusters=3)
km5cls = km5.fit(zdb.drop(['geometry', 'GEOCODIGO'], axis=1).values)
f, ax = plt.subplots(1, figsize=(9, 9))

zdb.assign(cl=km5cls.labels_)\
   .plot(column='cl', categorical=True, legend=True, \
         linewidth=0.1, edgecolor='white', ax=ax)
ax.set_axis_off()
plt.show()


km5cls.labels_

cl_pcts = prop_types_pct.reindex(zdb['GEOCODIGO'])\
                          .assign(cl=km5cls.labels_)\
                          .groupby('cl')\
                          .count()

"""
Explanation: Agregaci\'on de variables a nivel barrio
End of explanation
"""


cl_pcts

f, ax = plt.subplots(1, figsize=(18, 9))
cl_pcts.plot(kind='bar', stacked=False, ax=ax, \
             cmap='Set2', linewidth=2)
ax.legend(ncol=1, loc="right");

plt.show()

type(cl_pcts)

zdb.info()

rt_av = db.groupby('GEOCODIGO')[varis]\
           .mean()\
           .rename(lambda x: str(int(x)))
           
#pasar a int para join
#rt_av.index = rt_av.index.astype(int)

rt_av.describe()

zc['GEOCODIGO'] = pd.to_numeric(zc['GEOCODIGO'])

#pasar a int para join
rt_av.index = rt_av.index.astype(int)


zrt = zc[['geometry', 'GEOCODIGO']].join(rt_av, on='GEOCODIGO')\
                                         .dropna()
zrt.info()

zrt


zrt.to_file('tmp')
#matriz de pesos espaciales
w = ps.queen_from_shapefile('tmp/tmp.shp', idVariable='GEOCODIGO')

#rm -r tmp
w

"""
Explanation: N\'umero de bocas de metro en los barrios que componen la zona
End of explanation
"""


n_rev = db.groupby('GEOCODIGO')\
           .count()\
           ['price']\
           
thr = np.round(0.1 * n_rev.sum())
thr

np.random.seed(1234)
z = zrt.drop(['geometry', 'GEOCODIGO'], axis=1).values
maxp = ps.region.Maxp(w, z, thr, n_rev, initial=1000)


"""
Explanation: Establecer minimo de viviendas para region
End of explanation
"""



np.random.seed(1234)
maxp.cinference(nperm=20)


maxp.cpvalue


lbls = pd.Series(maxp.area2region).reindex(zrt['GEOCODIGO'])



f, ax = plt.subplots(1, figsize=(9, 9))

zrt.assign(cl=lbls.values)\
   .plot(column='cl', categorical=True, legend=True, \
         linewidth=0.1, edgecolor='white', ax=ax)

ax.set_axis_off()

plt.show()


lbls


#ver stats
zrt[varis].groupby(lbls.values).mean().T


#ver stats
zrt[varis].groupby(lbls.values).mean().T

"""
Explanation: Inferencia para comprobar que los resultados son mejores que definiendo zonas al azar
End of explanation
"""


