[SciPy-User] moving window (2D) correlation coefficient

[josef.pktd at gmail.com]

On Wed, Feb 10, 2010 at 11:42 AM, Zachary Pincus
<zachary.pincus at yale.edu> wrote:
> I bet that you could construct an array with shape (x,y,5,5), where
> array[i,j,:,:] would give the 5x5 window around (i,j), as some sort of
> mind-bending view on an array of shape (x,y), using a positive offset
> and some dimensions having negative strides. Then you could compute
> the correlation coefficient between the two arrays directly. Maybe?
>
> Probably cython would be more maintainable...
>
> Zach
>
>
> On Feb 10, 2010, at 10:45 AM, Vincent Schut wrote:
>
>> Hi,
>>
>> I need to calculate the correlation coefficient of a (simultaneously)
>> moving window over 2 images, such that the resulting pixel x,y (center
>> of the window) is the corrcoef((a 5x5 window around x,y for image
>> A), (a
>> 5x5 window around x,y for image B)).
>> Currently, I just loop over y, x, and then call corrcoef for each x,y
>> window. Would there be a better way, other than converting the loop to
>> cython?
>>
>>
>> For clarity (or not), the relevant part from my code:
>>
>>
>> for y in range(d500shape[2]):
>>     for x in range(d500shape[3]):
>>         if valid500[d,y,x]:
>>             window = spectral500Bordered[d,:,y:y+5, x:x
>> +5].reshape(7, -1)
>>             for b in range(5):
>>                 nonzeroMask = (window[0] > 0)
>>                 b01corr[0,b,y,x] =
>> numpy.corrcoef(window[0][nonzeroMask], window[b+2][nonzeroMask])[0,1]
>>                 b01corr[1,b,y,x] =
>> numpy.corrcoef(window[1][nonzeroMask], window[b+2][nonzeroMask])[0,1]
>>
>>
>> forget the 'if valid500' and 'nonzeroMask', those are to prevent
>> calculating pixels that don't need to be calculated, and to feed only
>> valid pixels from the window into corrcoef
>> spectral500Bordered is essentially a [d dates, 7 images, ysize, xsize]
>> array. I work per date (d), then calculate the corrcoef for images[0]
>> versus images[2:], and for images[1] versus images[2:]

I wrote a moving correlation for time series last november (scipy user
and preceding discussion on numpy mailing list)
I don't work with pictures, so I don't know if this can be extended to
your case. Since signal.correlate or convolve work in all directions
it might be possible

   def yxcov(self):
       xys = signal.correlate(self.x*self.y, self.kern,
mode='same')[self.sslice]
       return xys/self.n - self.ymean*self.xmean

Josef

>>
>> Thanks,
>> Vincent.
>>
>> _______________________________________________
>> SciPy-User mailing list
>> SciPy-User at scipy.org
>> http://mail.scipy.org/mailman/listinfo/scipy-user
>
> _______________________________________________
> SciPy-User mailing list
> SciPy-User at scipy.org
> http://mail.scipy.org/mailman/listinfo/scipy-user
>