In [6]:
# Load US county boundaries shapefile using geopandas
shapefile_path = '/Users/chenchenren/postdoc/paper/2N and water-US/Data/cb_2018_us_county_20m/surface_shapefile.shp'
gdf_county = gpd.read_file(shapefile_path)
# Ensure 'GEOID' column in the shapefile has the 'object' data15 type
gdf_county['GEOID'] = gdf_county['GEOID'].astype(float)
# Rename the 'geoid' column to 'GEOID' in the 'data15' Data15Frame
data15.rename(columns={'geoid': 'GEOID'}, inplace=True)
# Convert the 'GEOID' column in the 'data15' Data15Frame to 'object' data15 type
data15['GEOID'] = data15['GEOID'].astype(float)
# Merge 'data15' with county boundaries based on 'GEOID' attribute
merged_data15 = gdf_county.merge(data15, on='GEOID', how='left')
# Filter the data15 to include only geometries within the specified longitude and latitude range
filtered_data15 = merged_data15.cx[-130:-60, 25:50]
# Plot the county boundaries
fig, ax = plt.subplots(figsize=(10, 4))
shapefile_path = '/Users/chenchenren/postdoc/paper/2N and water-US/Data/cb_2018_us_state_20m/cb_2018_us_state_20m.shp'
gdf_state = gpd.read_file(shapefile_path)
gdf_state.rename(columns={'GEOID': 'GEOID_state'}, inplace=True)
gdf_state = gdf_state.cx[-130:-60, 25:50]
gdf_state.boundary.plot(ax=ax, linewidth=0.3, color='grey')
# Create a ScalarMappable without plotting the data15
cmap = mcolors.ListedColormap(['white'] + list(plt.cm.PiYG(np.linspace(0, 1, 255))))
norm = mcolors.Normalize(vmin=-4.2, vmax=4.2)
sc = plt.cm.ScalarMappable(cmap=cmap, norm=norm)
sc.set_array([])
# Add colorbar
custom_ticks = [-4, -2, 0, 2, 4]
cbar = plt.colorbar(sc, ax=ax, label='Net benefits (Million $)', orientation='vertical', pad=0.02, shrink=1, ticks=custom_ticks)
cbar.ax.set_ylabel('Net benefits (Million $)', fontsize=16, fontfamily='Arial')
cbar.ax.tick_params(labelsize=16)
cbar.ax.set_yticklabels(custom_ticks, fontsize=16, fontfamily='Arial')
# Plot the mean_gap_irrigation values
filtered_data15.plot(column='net_benefits_irri', cmap=cmap, ax=ax, norm=norm)
# Set font properties for the legend
font_prop = FontProperties(family='Arial', size=16)
font_prop.set_weight('bold')
# Additional plot configurations
ax.set_xticks([]) # Remove x-axis ticks
ax.set_yticks([]) # Remove y-axis ticks
ax.set_frame_on(False) # Turn off the frame
plt.text(0.28, 0.22, '(d) Soybean', transform=plt.gcf().transFigure,
fontsize=18, fontweight='bold', fontfamily='Arial', ha='left', va='bottom')
# Save the entire figure as a JPG with DPI 300
plt.savefig('soybean_profit_irri_15.png', format='png',
dpi=300, bbox_inches='tight', pad_inches=0.1, transparent=True)
# Show the plot
plt.show()
/var/folders/vd/0_phd7hx2n51y4412862zww00000gp/T/ipykernel_12641/1608784225.py:9: SettingWithCopyWarning:
A value is trying to be set on a copy of a slice from a DataFrame
See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy
data15.rename(columns={'geoid': 'GEOID'}, inplace=True)
/var/folders/vd/0_phd7hx2n51y4412862zww00000gp/T/ipykernel_12641/1608784225.py:12: SettingWithCopyWarning:
A value is trying to be set on a copy of a slice from a DataFrame.
Try using .loc[row_indexer,col_indexer] = value instead
See the caveats in the documentation: https://pandas.pydata.org/pandas-docs/stable/user_guide/indexing.html#returning-a-view-versus-a-copy
data15['GEOID'] = data15['GEOID'].astype(float)
