Python GIS: Interpolating and Plotting Corn Grain Yield Data

Python GIS: Interpolating and Plotting Corn Grain Yield Data






I have corn yield data (Mg/ha) that I want to visualize. First, let’s upload the data.

import pandas as pd
import requests
from io import StringIO

github="https://raw.githubusercontent.com/agronomy4future/raw_data_practice/main/grain_yield_map.csv"
response=requests.get(github)
df=pd.read_csv(StringIO(response.text))

df.head(5)
	Latitude	Longitude	GY
0	12.15725	-106.14035	16.248
1	12.15724	-106.13994	26.703
2	12.15724	-106.13954	16.569
3	12.15723	-106.13911	21.451
4	12.15722	-106.13871	19.107F

First, I’ll create yield distribution data.

# to open library
import numpy as np
import pandas as pd
import scipy.stats as stats
import matplotlib.pyplot as plt
from pylab import rcParams
import seaborn as sns

# to calculate mean and standard deviation
df_mean= np.mean(df["GY"])  
df_std= np.std(df["GY"])  

# to calculate probability density function
df_pdf= stats.norm.pdf(df["GY"].sort_values(), df_mean, df_std)

# to draw normal distribution (with histogram) gragh
plt.plot(df["GY"].sort_values(), df_pdf, color= "Black")
sns.histplot(data= df["GY"], color= "Black", bins= 70, stat= "probability", alpha= 0.3)

plt.xlim([0,50])
plt.ylim([0,0.1])
plt.xlabel("Grain yield (Mg/ha)", size= 15)
plt.ylabel("Frequency", size= 15)
plt.grid(True, alpha=0.3, linestyle="--")
plt.rcParams["figure.figsize"]= [7,5]  
plt.rcParams["figure.dpi"]= 500 
plt.show()

Now, we can see that the general grain yield varies from 10 to 30 Mg/ha, with some outliers. I’ll create a yield map to visualize this variation.

import matplotlib.pyplot as plt
import numpy as np
from scipy.interpolate import griddata

# Extract the data
latitude = df['Latitude']
longitude = df['Longitude']
p_concentration = df['GY']

# Define grid for interpolation
grid_x, grid_y = np.mgrid[latitude.min():latitude.max():100j, longitude.min():longitude.max():100j]

# Interpolate the data to a grid
grid_z = griddata((latitude, longitude), p_concentration, (grid_x, grid_y), method='cubic')

# Define the levels for a simpler contour plot
levels = list(range(0, 50, 5))

# Create the plot with the specified levels and remove axis numbers
plt.figure(figsize=(10, 6))
contour = plt.contourf(grid_y, grid_x, grid_z, levels=levels, cmap='RdYlGn') 
plt.colorbar(contour, label='Corn Grain Yield (Mg/ha)')

# Add scatter points
plt.scatter(longitude, latitude, c='black', s=7)

# axis
plt.xticks([])
plt.yticks([])

plt.show()
# Full code
import pandas as pd
import requests
from io import StringIO

import matplotlib.pyplot as plt
import numpy as np
from scipy.interpolate import griddata

github="https://raw.githubusercontent.com/agronomy4future/raw_data_practice/main/grain_yield_map.csv"
response=requests.get(github)
df=pd.read_csv(StringIO(response.text))

latitude = df['Latitude']
longitude = df['Longitude']
p_concentration = df['GY']

grid_x, grid_y = np.mgrid[latitude.min():latitude.max():100j, longitude.min():longitude.max():100j]

grid_z = griddata((latitude, longitude), p_concentration, (grid_x, grid_y), method='cubic')

levels = list(range(0, 50, 5))

plt.figure(figsize=(10, 6))

contour = plt.contourf(grid_y, grid_x, grid_z, levels=levels, cmap='RdYlGn')
 
plt.colorbar(contour, label='Corn Grain Yield (Mg/ha)')
plt.scatter(longitude, latitude, c='black', s=7)

plt.xticks([])
plt.yticks([])

plt.show()

https://github.com/agronomy4future/python_code/blob/main/Python_GIS_Interpolating_and_Plotting_Corn_Grain_Yield_Data.ipynb


If you have the ArcGIS program, you can create a more detailed GIS map.






Leave a Reply

If you include a website address in the comment section, I cannot see your comment as it will be automatically deleted and will not be posted. Please refrain from including website addresses.