I am trying to implement this metric

I already managed to calculate NUBN with numpy operations so that is fast, but I can't find a way to escape python slow looping to calculate the DRD part. Here is my current calculation of DRD:
def drd(im, im_gt):
height, width = im.shape
W = np.array([[1/math.sqrt(x**2+y**2) if x != 0 or y != 0 else 0 for x in range(-2, 3)] for y in range(-2, 3)])
W /= W.sum()
drd = 0
s = []
for y, x in zip(*np.where(im_gt != im)):
if x > 1 and y > 1 and x + 2 < width and y + 2 < height:
s.append(im_gt[y-2:y+3, x-2:x+3] == im_gt[y, x])
else:
for yy in range(y-2, y+3):
for xx in range(x-2, x+3):
if xx > 1 and yy > 1 and xx < width - 1 and yy < height - 1:
drd += abs(im_gt[yy, xx] - im[y, x]) * W[yy-y+2, xx-x+2]
return drd + np.sum(s * W)
drd(np.random.choice([0, 1], size=(100, 100)), np.random.choice([0, 1], size=(100, 100)))
Can anyone think of a faster way to do this? Timings on 1000x1000:
