title: Python: binned_statistic_2d mean calculation ignoring NaNs in data
body: <p>I am using <code>scipy.stats.binned_statistic_2d</code> to bin irregular data onto a uniform grid by finding the mean of points within every bin.</p>
  
  <pre><code>x,y = np.meshgrid(sort(np.random.uniform(0,1,100)),sort(np.random.uniform(0,1,100)))
  z = np.sin(x*y)
  
  statistic, xedges, yedges, binnumber = sp.stats.binned_statistic_2d(x.ravel(), y.ravel(), values=z.ravel(), statistic='mean',bins=[np.arange(0,1.1,.1), np.arange(0,1.1,.1)])
  
  plt.figure(1)
  plt.pcolormesh(x,y,z, vmin = 0, vmax = 1)
  plt.figure(2)
  plt.pcolormesh(xedges,yedges,statistic, vmin = 0, vmax = 1)
  </code></pre>
  
  <p>Produces these plots, as expected:</p>
  
  <p>Scattered data:
  <a href="https://i.stack.imgur.com/3uafA.png" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/3uafA.png" alt="enter image description here"></a></p>
  
  <p>Gridded data:
  <a href="https://i.stack.imgur.com/IPb28.png" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/IPb28.png" alt="enter image description here"></a></p>
  
  <p>But the data I want to grid has NaNs in it. This is what the result is like when I add NaNs:</p>
  
  <pre><code>x,y = np.meshgrid(sort(np.random.uniform(0,1,100)),sort(np.random.uniform(0,1,100)))
  z = np.sin(x*y)
  z[50:55,50:55] = np.nan
  
  statistic, xedges, yedges, binnumber = binned_statistic_2d(x.ravel(), y.ravel(), values=z.ravel(), statistic='mean',bins=[np.arange(0,1.1,.1), np.arange(0,1.1,.1)])
  
  plt.figure(3)
  plt.pcolormesh(x,y,z, vmin = 0, vmax = 1)
  plt.figure(4)
  plt.pcolormesh(xedges,yedges,statistic, vmin = 0, vmax = 1)
  </code></pre>
  
  <p>Scattered:
  <a href="https://i.stack.imgur.com/PC5BK.png" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/PC5BK.png" alt="enter image description here"></a></p>
  
  <p>Gridded:
  <a href="https://i.stack.imgur.com/gI5Fh.png" rel="nofollow noreferrer"><img src="https://i.stack.imgur.com/gI5Fh.png" alt="enter image description here"></a></p>
  
  <p>Obviously if a bin is entirely filled with NaNs, the the resulting mean of that bin should still be NaN. However, I would like bins that are not <em>entirely</em> filled with NaNs to just result in the mean of the non-NaN numbers.</p>
  
  <p>I've tried replacing the "statistic" argument in <code>sp.stats.binned_statistic_2d</code> with <code>np.nanmean</code>. This works, but it goes very very slowly when I use it on large datasets. I've tried digging into the underlying code of `sp.stats.binned_statistic_2d', but I can't figure out exactly how it is calculating the mean, or how to make it ignore NaNs in it's calculation.</p>
  
  <p>Any ideas?</p>