From fe51c75d01da5b63addf41788191b149236ed86c Mon Sep 17 00:00:00 2001 From: Martin Fleischmann Date: Sat, 18 Jan 2025 21:49:25 +0100 Subject: [PATCH 1/4] def individual plotting exposed on Moran_Local_BV --- esda/moran.py | 58 +++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 58 insertions(+) diff --git a/esda/moran.py b/esda/moran.py index c1f3e2ed..1b1217b1 100644 --- a/esda/moran.py +++ b/esda/moran.py @@ -1722,6 +1722,64 @@ def by_col( **stat_kws, ) + def get_cluster_labels(self, crit_value=0.05): + """Return LISA cluster labels for each observation. + + Parameters + ---------- + crit_value : float, optional + crititical significance value for statistical inference, by default 0.05 + + Returns + ------- + numpy.array + an array of cluster labels aligned with the input data used to conduct the + local Moran analysis + """ + return _get_cluster_labels(self, crit_value) + + def explore(self, gdf, crit_value=0.05, **kwargs): + """Create interactive map of LISA indicators + + Parameters + ---------- + gdf : geopandas.GeoDataFrame + geodataframe used to conduct the local Moran analysis + crit_value : float, optional + critical value to determine statistical significance, by default 0.05 + kwargs : dict, optional + additional keyword arguments passed to the geopandas `explore` method + + Returns + ------- + Folium.Map + interactive map with LISA clusters + """ + gdf = gdf.copy() + gdf["Moran Cluster"] = self.get_cluster_labels(crit_value) + return _viz_local_moran(self, gdf, crit_value, "explore", **kwargs) + + def plot(self, gdf, crit_value=0.05, **kwargs): + """Create static map of LISA indicators + + Parameters + ---------- + gdf : geopandas.GeoDataFrame + geodataframe used to conduct the local Moran analysis + crit_value : float, optional + critical value to determine statistical significance, by default 0.05 + kwargs : dict, optional + additional keyword arguments passed to the geopandas `explore` method + + Returns + ------- + ax + matplotlib axis + """ + gdf = gdf.copy() + gdf["Moran Cluster"] = self.get_cluster_labels(crit_value) + return _viz_local_moran(self, gdf, crit_value, "plot", **kwargs) + def plot_scatter( self, crit_value=0.05, From 29033c84764161f0ee727fdfbcf4e0137f4ca61a Mon Sep 17 00:00:00 2001 From: Martin Fleischmann Date: Sat, 18 Jan 2025 21:49:36 +0100 Subject: [PATCH 2/4] implement combination subplot --- esda/moran.py | 363 ++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 363 insertions(+) diff --git a/esda/moran.py b/esda/moran.py index 1b1217b1..29a6ddb5 100644 --- a/esda/moran.py +++ b/esda/moran.py @@ -1455,6 +1455,86 @@ def plot_scatter( fitline_kwds=fitline_kwds, ) + def plot_combination( + self, + gdf, + attribute, + crit_value=0.05, + region_column=None, + mask=None, + mask_color="#636363", + quadrant=None, + legend=True, + scheme="Quantiles", + cmap="YlGnBu", + figsize=(15, 4), + scatter_kwds=None, + fitline_kwds=None, + legend_kwds=None, + ): + """ + Produce three-plot visualisation of Moran Scatteprlot, LISA cluster + and Choropleth maps, with Local Moran region and quadrant masking + + Parameters + ---------- + gdf : geopandas.GeoDataFrame + geodataframe used to conduct the local Moran analysis + attribute : str + Column name of attribute which should be depicted in Choropleth map. + crit_value : float, optional + critical value to determine statistical significance, by default 0.05 + region_column: string, optional + Column name containing mask region of interest, by default None + mask: str, float, int, optional + Identifier or name of the region to highlight, by default None + Use the same dtype to specifiy as in original dataset. + mask_color: str, optional + Color of mask, by default '#636363'. + quadrant : int, optional + Quadrant 1-4 in scatterplot masking values in LISA cluster and + Choropleth maps, by default None + figsize: tuple, optional + W, h of figure, by default (15,4) + legend: boolean, optional + If True, legend for maps will be depicted, by default True + scheme: str, optional + Name of mapclassify classifier to be used, by default 'Quantiles' + cmap: str, optional + Name of matplotlib colormap used for plotting the Choropleth. + By default 'YlGnBu'. + scatter_kwds : keyword arguments, optional + Keywords used for creating and designing the scatter points, by default + None. + fitline_kwds : keyword arguments, optional + Keywords used for creating and designing the moran fitline + in the scatterplot, by default None. + legend_kwds : dict + Keyword arguments passed to geopandas.GeodataFrame.plot ``legend_kwds`` + allowing repositioning of the legend in LISA cluster plot and choropleth. + + Returns + ------- + axs : array of Matplotlib axes + """ + return _plot_combination( + self, + gdf, + attribute, + crit_value=crit_value, + region_column=region_column, + mask=mask, + mask_color=mask_color, + quadrant=quadrant, + legend=legend, + scheme=scheme, + cmap=cmap, + figsize=figsize, + scatter_kwds=scatter_kwds, + fitline_kwds=fitline_kwds, + legend_kwds=legend_kwds, + ) + class Moran_Local_BV: # noqa: N801 """Bivariate Local Moran Statistics. @@ -1815,6 +1895,86 @@ def plot_scatter( fitline_kwds=fitline_kwds, ) + def plot_combination( + self, + gdf, + attribute, + crit_value=0.05, + region_column=None, + mask=None, + mask_color="#636363", + quadrant=None, + legend=True, + scheme="Quantiles", + cmap="YlGnBu", + figsize=(15, 4), + scatter_kwds=None, + fitline_kwds=None, + legend_kwds=None, + ): + """ + Produce three-plot visualisation of Moran Scatteprlot, LISA cluster + and Choropleth maps, with Local Moran region and quadrant masking + + Parameters + ---------- + gdf : geopandas.GeoDataFrame + geodataframe used to conduct the local Moran analysis + attribute : str + Column name of attribute which should be depicted in Choropleth map. + crit_value : float, optional + critical value to determine statistical significance, by default 0.05 + region_column: string, optional + Column name containing mask region of interest, by default None + mask: str, float, int, optional + Identifier or name of the region to highlight, by default None + Use the same dtype to specifiy as in original dataset. + mask_color: str, optional + Color of mask, by default '#636363'. + quadrant : int, optional + Quadrant 1-4 in scatterplot masking values in LISA cluster and + Choropleth maps, by default None + figsize: tuple, optional + W, h of figure, by default (15,4) + legend: boolean, optional + If True, legend for maps will be depicted, by default True + scheme: str, optional + Name of mapclassify classifier to be used, by default 'Quantiles' + cmap: str, optional + Name of matplotlib colormap used for plotting the Choropleth. + By default 'YlGnBu'. + scatter_kwds : keyword arguments, optional + Keywords used for creating and designing the scatter points, by default + None. + fitline_kwds : keyword arguments, optional + Keywords used for creating and designing the moran fitline + in the scatterplot, by default None. + legend_kwds : dict + Keyword arguments passed to geopandas.GeodataFrame.plot ``legend_kwds`` + allowing repositioning of the legend in LISA cluster plot and choropleth. + + Returns + ------- + axs : array of Matplotlib axes + """ + return _plot_combination( + self, + gdf, + attribute, + crit_value=crit_value, + region_column=region_column, + mask=mask, + mask_color=mask_color, + quadrant=quadrant, + legend=legend, + scheme=scheme, + cmap=cmap, + figsize=figsize, + scatter_kwds=scatter_kwds, + fitline_kwds=fitline_kwds, + legend_kwds=legend_kwds, + ) + class Moran_Local_Rate(Moran_Local): # noqa: N801 """ @@ -2285,6 +2445,209 @@ def _simulation_plot( return ax +def _plot_combination( + moran_loc, + gdf, + attribute, + crit_value=0.05, + region_column=None, + mask=None, + mask_color="#636363", + quadrant=None, + legend=True, + scheme="Quantiles", + cmap="YlGnBu", + figsize=(15, 4), + scatter_kwds=None, + fitline_kwds=None, + legend_kwds=None, +): + """ + Produce three-plot visualisation of Moran Scatteprlot, LISA cluster + and Choropleth maps, with Local Moran region and quadrant masking + + Parameters + ---------- + moran_loc : esda.moran.Moran_Local or Moran_Local_BV instance + Values of Moran's Local Autocorrelation Statistic + gdf : geopandas dataframe + The Dataframe containing information to plot the two maps. + attribute : str + Column name of attribute which should be depicted in Choropleth map. + p : float, optional + The p-value threshold for significance. Points and polygons will + be colored by significance. Default = 0.05. + region_column: string, optional + Column name containing mask region of interest. Default = None + mask: str, float, int, optional + Identifier or name of the region to highlight. Default = None + Use the same dtype to specifiy as in original dataset. + mask_color: str, optional + Color of mask. Default = '#636363' + quadrant : int, optional + Quadrant 1-4 in scatterplot masking values in LISA cluster and + Choropleth maps. Default = None + figsize: tuple, optional + W, h of figure. Default = (15,4) + legend: boolean, optional + If True, legend for maps will be depicted. Default = True + scheme: str, optional + Name of PySAL classifier to be used. Default = 'Quantiles' + cmap: str, optional + Name of matplotlib colormap used for plotting the Choropleth. + Default = 'YlGnBu' + scatter_kwds : keyword arguments, optional + Keywords used for creating and designing the scatter points. + Default =None. + fitline_kwds : keyword arguments, optional + Keywords used for creating and designing the moran fitline + in the scatterplot. Default =None. + legend_kwds : dict + Keyword arguments passed to geopandas.GeodataFrame.plot ``legend_kwds`` allowing + repositioning of the legend in LISA cluster plot and choropleth. + + Returns + ------- + axs : array of Matplotlib axes + """ + try: + from matplotlib import patches + from matplotlib import pyplot as plt + + except ImportError as err: + raise ImportError( + "matplotlib library must be installed to use the scatterplot feature" + ) from err + + _, axs = plt.subplots( + 1, 3, figsize=figsize, subplot_kw={"aspect": "equal", "adjustable": "datalim"} + ) + # Moran Scatterplot + moran_loc.plot_scatter( + crit_value=crit_value, + ax=axs[0], + scatter_kwds=scatter_kwds, + fitline_kwds=fitline_kwds, + ) + + # Lisa cluster map + moran_loc.plot( + gdf, + crit_value=crit_value, + ax=axs[1], + legend=legend, + legend_kwds=legend_kwds, + ) + + # Choropleth for attribute + gdf.plot( + column=attribute, + scheme=scheme, + cmap=cmap, + legend=legend, + legend_kwds=legend_kwds, + ax=axs[2], + alpha=1, + ) + axs[2].set_axis_off() + axs[2].set_aspect("equal") + + # MASKING QUADRANT VALUES + if quadrant is not None: + # Quadrant masking in Scatterplot + mask_angles = {1: 0, 2: 90, 3: 180, 4: 270} # rectangle angles + # We don't want to change the axis data limits, so use the current ones + xmin, xmax = axs[0].get_xlim() + ymin, ymax = axs[0].get_ylim() + # We are rotating, so we start from 0 degrees and + # figured out the right dimensions for the rectangles for other angles + mask_width = {1: abs(xmax), 2: abs(ymax), 3: abs(xmin), 4: abs(ymin)} + mask_height = {1: abs(ymax), 2: abs(xmin), 3: abs(ymin), 4: abs(xmax)} + axs[0].add_patch( + patches.Rectangle( + (0, 0), + width=mask_width[quadrant], + height=mask_height[quadrant], + angle=mask_angles[quadrant], + color="#E5E5E5", + zorder=-1, + alpha=0.8, + ) + ) + # quadrant selection in maps + non_quadrant = ~(moran_loc.q == quadrant) + mask_quadrant = gdf[non_quadrant] + df_quadrant = gdf.iloc[~non_quadrant] + union2 = df_quadrant.dissolve().boundary + + # LISA Cluster mask and cluster boundary + mask_quadrant.plot( + scheme=scheme, + color="white", + ax=axs[1], + alpha=0.7, + zorder=1, + ) + union2.plot(linewidth=1, ax=axs[1], color="#E5E5E5") + + # CHOROPLETH MASK + mask_quadrant.plot( + scheme=scheme, + color="white", + ax=axs[2], + alpha=0.7, + zorder=1, + ) + union2.plot(linewidth=1, ax=axs[2], color="#E5E5E5") + + # REGION MASKING + if region_column is not None: + # masking inside axs[0] or Moran Scatterplot + # enforce the same dtype of list and mask + if not isinstance(mask[0], type(gdf[region_column].iloc[0])): + warn( + "Values in `mask` are not the same dtype as" + + " values in `region_column`. Converting `mask` values" + + " to dtype of first observation in region_column.", + stacklevel=3, + ) + data_type = type(gdf[region_column][0].item()) + mask = list(map(data_type, mask)) + + ix = gdf[region_column].isin(mask) + + if not ix.any(): + raise ValueError( + f"Specified values {mask} in `mask` not in `region_column`" + ) + + df_mask = gdf[ix] + x_mask = moran_loc.z[ix] + y_mask = lag_spatial(moran_loc.w, moran_loc.z)[ix] + axs[0].plot( + x_mask, + y_mask, + color=mask_color, + marker="o", + markersize=14, + alpha=0.8, + linestyle="None", + zorder=-1, + ) + + # masking inside axs[1] or Lisa cluster map + union = df_mask.dissolve().boundary + union.plot(linewidth=2, ax=axs[1], color=mask_color) + + # masking inside axs[2] or Chloropleth + union.plot(linewidth=2, ax=axs[2], color=mask_color) + + axs[0].spines[["right", "top"]].set_visible(False) + axs[1].set_axis_off() + + return axs + + # -------------------------------------------------------------- # Conditional Randomization Moment Estimators # -------------------------------------------------------------- From 219717fa3805031068fef271a1ec1e74388065a4 Mon Sep 17 00:00:00 2001 From: Martin Fleischmann Date: Sat, 18 Jan 2025 22:04:15 +0100 Subject: [PATCH 3/4] test basic plotting and explore --- esda/tests/test_moran.py | 105 ++++++++++++++++++++++++++++++++++++++- 1 file changed, 103 insertions(+), 2 deletions(-) diff --git a/esda/tests/test_moran.py b/esda/tests/test_moran.py index adb75fdb..3df94d05 100644 --- a/esda/tests/test_moran.py +++ b/esda/tests/test_moran.py @@ -243,7 +243,9 @@ def test_Moran_plot_scatter_args(self, w): w, ) - ax = m.plot_scatter(scatter_kwds=dict(color='blue'), fitline_kwds=dict(color='pink')) + ax = m.plot_scatter( + scatter_kwds=dict(color="blue"), fitline_kwds=dict(color="pink") + ) # test scatter np.testing.assert_array_almost_equal( @@ -256,7 +258,6 @@ def test_Moran_plot_scatter_args(self, w): assert l.get_color() == "pink" - class TestMoranRate: def setup_method(self): f = libpysal.io.open(libpysal.examples.get_path("sids2.dbf")) @@ -588,6 +589,7 @@ def test_local_moments(self, w): class TestMoranLocalBV: def setup_method(self): f = libpysal.io.open(libpysal.examples.get_path("sids2.dbf")) + self.gdf = gpd.read_file(libpysal.examples.get_path("sids2.shp")) self.x = np.array(f.by_col["SIDR79"]) self.y = np.array(f.by_col["SIDR74"]) @@ -629,6 +631,105 @@ def test_by_col(self): np.testing.assert_allclose(bvz[0], 1.7900932313425777, 5) np.testing.assert_allclose(bvzp[0], 0.036719462378528744, 5) + @parametrize_sids + def test_Moran_LocalBV_labels(self, w): + lm = moran.Moran_Local_BV( + self.x, + self.y, + w, + transformation="r", + permutations=99, + keep_simulations=True, + seed=SEED, + ) + expected_labels = np.array( + [ + "Insignificant", + "Insignificant", + "Low-Low", + "High-Low", + "Low-High", + "Insignificant", + "Insignificant", + "Insignificant", + "Insignificant", + "Insignificant", + ] + ) + assert_array_equal(lm.get_cluster_labels()[:10], expected_labels) + assert_array_equal( + pd.Series(lm.get_cluster_labels(0.05)).value_counts().values, + np.array([80, 7, 6, 5, 2]), + ) + + @parametrize_sids + def test_Moran_LocalBV_explore(self, w): + lm = moran.Moran_Local_BV( + self.x, + self.y, + w, + transformation="r", + permutations=99, + keep_simulations=True, + seed=SEED, + ) + m = lm.explore(self.gdf) + np.testing.assert_array_equal( + m.get_bounds(), + [ + [33.88199234008789, -84.3238525390625], + [36.58964920043945, -75.45697784423828], + ], + ) + assert len(m.to_dict()["children"]) == 2 + + out_str = _fetch_map_string(m) + + assert '"High-High","__folium_color":"#d7191c"' in out_str + assert '"Low-High","__folium_color":"#89cff0"' in out_str + assert '"Low-Low","__folium_color":"#2c7bb6"' in out_str + assert '"High-Low","__folium_color":"#fdae61"' in out_str + assert '"Insignificant","__folium_color":"#d3d3d3"' in out_str + + assert out_str.count("#d7191c") == 10 + assert out_str.count("#89cff0") == 5 + assert out_str.count("#2c7bb6") == 9 + assert out_str.count("#fdae61") == 8 + assert out_str.count("#d3d3d3") == 83 + + @parametrize_sids + def test_Moran_LocalBV_plot(self, w): + import matplotlib + + matplotlib.use("Agg") + + lm = moran.Moran_Local_BV( + self.x, + self.y, + w, + transformation="r", + permutations=99, + keep_simulations=True, + seed=SEED, + ) + ax = lm.plot(self.gdf) + unique, counts = np.unique( + ax.collections[0].get_facecolors(), axis=0, return_counts=True + ) + np.testing.assert_array_almost_equal( + unique, + np.array( + [ + [0.17254902, 0.48235294, 0.71372549, 1.0], + [0.5372549, 0.81176471, 0.94117647, 1.0], + [0.82745098, 0.82745098, 0.82745098, 1.0], + [0.84313725, 0.09803922, 0.10980392, 1.0], + [0.99215686, 0.68235294, 0.38039216, 1.0], + ] + ), + ) + np.testing.assert_array_equal(counts, np.array([6, 2, 86, 7, 7])) + class TestMoranLocalRate: def setup_method(self): From b81d4edadc8927702ed965121b45572e97bc9017 Mon Sep 17 00:00:00 2001 From: Martin Fleischmann Date: Sat, 18 Jan 2025 22:23:22 +0100 Subject: [PATCH 4/4] tests --- esda/tests/test_moran.py | 66 ++++++++++++++++++++++++++++++++++++++++ 1 file changed, 66 insertions(+) diff --git a/esda/tests/test_moran.py b/esda/tests/test_moran.py index 3df94d05..3116c6e8 100644 --- a/esda/tests/test_moran.py +++ b/esda/tests/test_moran.py @@ -585,6 +585,45 @@ def test_local_moments(self, w): np.testing.assert_allclose(lm.EI, EI, rtol=RTOL, atol=ATOL) np.testing.assert_allclose(lm.VI, VI, rtol=RTOL, atol=ATOL) + @parametrize_sac + def test_plot_combination(self, w): + import matplotlib + + matplotlib.use("Agg") + + lm = moran.Moran_Local( + sac1.WHITE, + w, + transformation="r", + permutations=99, + keep_simulations=True, + seed=SEED, + ) + axs = lm.plot_combination( + sac1, + "WHITE", + legend_kwds=dict(loc="lower right"), + region_column="FIPS", + mask=["06067009504", "06067009503"], + quadrant=1, + ) + + assert len(axs) == 3 + assert len(axs[0].patches) == 1 + assert len(axs[1].collections) == 4 + assert len(axs[2].collections) == 4 + + axs2 = lm.plot_combination( + sac1, + "WHITE", + legend_kwds=dict(loc="lower right"), + ) + + assert len(axs2) == 3 + assert len(axs2[0].patches) == 0 + assert len(axs2[1].collections) == 1 + assert len(axs2[2].collections) == 1 + class TestMoranLocalBV: def setup_method(self): @@ -730,6 +769,33 @@ def test_Moran_LocalBV_plot(self, w): ) np.testing.assert_array_equal(counts, np.array([6, 2, 86, 7, 7])) + @parametrize_sids + def test_plot_combination(self, w): + import matplotlib + + matplotlib.use("Agg") + + lm = moran.Moran_Local_BV( + self.x, + self.y, + w, + transformation="r", + permutations=99, + keep_simulations=True, + seed=SEED, + ) + axs = lm.plot_combination( + self.gdf, + "SIDR79", + legend_kwds=dict(loc="lower right"), + quadrant=1, + ) + + assert len(axs) == 3 + assert len(axs[0].patches) == 1 + assert len(axs[1].collections) == 3 + assert len(axs[2].collections) == 3 + class TestMoranLocalRate: def setup_method(self):