# create ZonalStatsPoly instancegc=wrl.zonalstats.ZonalStatsPoly(zdpoly_new)# create some artificial data for processing using the features indicescount=radar_gk.shape[0]data=1000000./np.array(range(count))# calculate mean and variancemean=gc.mean(data)var=gc.var(data)print("Average:",mean)print("Variance:",var)

Average: [13.03540605 12.43517116]
Variance: [4.10381513 1.88270044]

Next we can export the resulting zonal statistics to vector and raster
files: