In [7]:
%pylab inline

import json
import os
os.environ['XPA_METHOD'] = 'local'

from bubbly.viewer import BubblyViewer
from bubbly.extractors import RGBExtractor
from bubbly.cluster import merge, xmatch
import matplotlib
matplotlib.rcParams['axes.grid'] = False


Populating the interactive namespace from numpy and matplotlib

In [3]:
bv = BubblyViewer()
bv.load_longitude(35)

In [4]:
data = json.load(open('../models/l035_scores.json'))

stamps = np.array(data['stamps'])
scores = np.array(data['scores'])

rank = np.argsort(scores)[::-1]
lohi = np.argsort(scores)

Distribution of Scores


In [8]:
hist(scores, histtype='step')
yscale('log')

xlabel("Scores")
ylabel("N")


Out[8]:
<matplotlib.text.Text at 0x10dd79a90>

Representative images at different scores

score >= 0.1 seems to be mostly pure


In [9]:
ex = RGBExtractor()
ex.shp = (100, 100)

ss = np.linspace(1, -1, 20)
ss += np.random.random(ss.size) * .01
inds = lohi[np.minimum(np.searchsorted(scores[lohi], ss), lohi.size-1)]

ims = [ex.extract(*stamps[i]) for i in inds]

nx = 5
ny = int(np.ceil(inds.size / 5))

im = np.vstack(np.hstack(ims[i:i+nx]) for i in range(0, inds.size, nx))

figure(figsize=(15, 15))
imshow(im, origin='upper')

opts = dict(color = '#0000ff', ha='center', va='center', fontsize=18)
for i in range(inds.size):
    x = (i % nx + 0.5) * ex.shp[0]
    y = (i / nx + 0.5) * ex.shp[0]    
    annotate("%0.2f" % scores[inds[i]], xy = (x, y), **opts)


Loading a new field at l=35

Clustering Detections


In [13]:
from bubbly.field import get_field
from bubbly.util import scale
from bubbly.dr1 import bubble_params

f = get_field(35)
g = scale(f.i4[::7, ::7])
r = scale(f.mips[::7, ::7])

b = r * 0
im = np.dstack((r, g, b))

def plot_stamps(stamps, **kwargs):
    kwargs.setdefault('facecolor', 'none')
    kwargs.setdefault('edgecolor', 'b')
    kwargs.setdefault('alpha', .1)
    
    ax = gca()
    for s in stamps:
        r = Rectangle((s[1] - s[-1], s[2] - s[-1]), width = 2 * s[-1], height = 2 * s[-1], **kwargs)
        ax.add_patch(r)

mwp = [b for b in bubble_params() if b[0] == 35]

unmerged = stamps[scores > .2]
merged, scores = merge(unmerged, scores[scores > 0.2])

In [ ]:
from sklearn.cluster import DBSCAN

In [14]:
figure(figsize=(10, 10), dpi=100)
imshow(im, extent=[36, 34, -1, 1])

plot_stamps(unmerged)
plot_stamps(merged, edgecolor='r', alpha=1)
xlim(35.5, 34.5)
ylim(-1, 0)
gca().set_aspect('equal')


Cross Matching MWP with Bubbly


In [17]:
match = xmatch(merged, mwp)
print 'Bubbly -> mwp: (%i)' % ((match >= 0).sum())
for i, j in enumerate(match):    
    if j < 0:
        continue
    print '%i -> %i' % (i, j)
    
print 'Unmatched Bubbly (%i / %i)' % ((match == -1).sum(), merged.shape[0])
print np.where(match == -1)[0].tolist()

print 'Unmatched mwp (%i / %i)' % (len(mwp) - (match >= 0).sum(), len(mwp))
z = np.zeros(len(mwp))
z[match] = 1
print np.where(z == 0)[0].tolist()


Bubbly -> mwp: (13)
1 -> 6
3 -> 9
4 -> 11
6 -> 13
7 -> 10
11 -> 16
12 -> 5
13 -> 18
15 -> 31
16 -> 4
17 -> 2
18 -> 8
19 -> 0
Unmatched Bubbly (7 / 20)
[0, 2, 5, 8, 9, 10, 14]
Unmatched mwp (19 / 32)
[1, 3, 7, 12, 14, 15, 17, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30]

Examining each un-matched bubble


In [18]:
len(merged), len(mwp)


Out[18]:
(20, 32)

In [15]:
#bubbly in blue, MWP in red
bv.clear()
bv.outline(merged, color='red')
bv.outline(mwp)

Manual Notes on mismatches

Brut-only (7)

Real and distinct (1)

14

Real, but in regions of confusion / double identification (2)

2, 9

False positives -- saturated MIPS (2)

0, 5

Ambiguous Nebulosity (2)

8, 10

MWP only sources (19)

Too small (3)

3, 7, 20

Edge (3)

15, 28, 30

True positive (5)

12, 19, 22, 26, 29

Confusion (2)

1, 23

False -- other nebula (6)

14, 17, 21, 24, 25, 27


In [46]:
6/32.


Out[46]:
0.1875

Summary, when using only large stamps

Bubbly finds 20 stamps, MWP finds 32 (or 26 in Bubbly's search strategy). Of the 20 bubblys:

  • 65% (13) match unambiguously to a MWP source
  • 10% (2) "ambiguously" match to a MWP source -- that is, one or more MWP bubbles are split up into different regions.
  • 5% (1) are new bubbles not in the MWP project
  • 10% (2) are regions confused by extended / saturated 24 um emission
  • 10% (2) are regions of nebulosity, that don't look particularly shell-like

In addition, ignoring confusion examples, there are 17 bubbles in the MWP with no bubbly counterpart:

  • (6) are false positives -- mostly other nebulosity
  • (5) are true positives, wrongly rejected by bubbly
  • (3) are at the image edges, and un-classifiable by bubbly
  • (3) are smaller than the search radius bubbly considers

Overall, bubbly recovers 16 / 21 bubbles it is sensitive to (76%), at a false positive rate of 4 / 25 (16%)

The MWP, on the other hand, recovers 20/21 bubbles (95%), at a false positive rate of 18%.


In [10]:
from bubbly.model import ModelGroup
model = ModelGroup.load('../models/full_classifier.dat')

In [22]:
ind = 14
cat = merged
bv.look_at(cat[ind])
imshow(ex.extract(cat[ind][0], cat[ind][1], cat[ind][2], cat[ind][3] * 2))
show()
imshow(ex.extract(*cat[ind]))
show()
model.decision_function([cat[ind]])


Out[22]:
array([ 0.3325])

In [49]:
len(merged), len(mwp)


Out[49]:
(20, 32)

The stamps


In [47]:
merged.tolist()


Out[47]:
[[35.0, 34.864166666666655, -0.948611111111187, 0.022222222222222223],
 [35.0, 35.06416666666667, -0.7797222222222846, 0.022222222222222223],
 [35.0, 35.11305555555556, -0.748611111111171, 0.022222222222222223],
 [35.0, 34.695277777777754, -0.6508333333333854, 0.022222222222222223],
 [35.0, 35.05083333333334, -0.5175000000000414, 0.022222222222222223],
 [35.0, 35.4730555555556, -0.4419444444444798, 0.022222222222222223],
 [35.0, 35.4241666666667, -0.2686111111111326, 0.022222222222222223],
 [35.0, 34.864166666666655, -0.081944444444451, 0.022222222222222223],
 [35.0, 34.75749999999998, 0.0202777777777794, 0.022222222222222223],
 [35.0, 35.02861111111111, 0.3536111111111394, 0.022222222222222223],
 [35.0, 34.930833333333325, 0.473611111111149, 0.022222222222222223],
 [35.0, 35.14861111111112, 0.806944444444509, 0.022222222222222223],
 [35.0, 34.74138888888887, -0.3302777777778042, 0.027777777777777776],
 [35.0, 35.47472222222226, -0.0025000000000002, 0.027777777777777776],
 [35.0, 35.341388888888915, 0.0141666666666678, 0.027777777777777776],
 [35.0, 35.135833333333345, -0.7558333333333939, 0.034444444444444444],
 [35.0, 35.036944444444444, -0.4952777777778174, 0.042777777777777776],
 [35.0, 35.0125, 0.3358333333333602, 0.05333333333333334],
 [35.0, 34.94138888888889, -0.0247222222222242, 0.08333333333333333],
 [35.0, 34.755833333333314, -0.6808333333333878, 0.10388888888888889]]

In [48]:
mwp


Out[48]:
[[35, 34.7525926, -0.67014099999999999, 0.11306971],
 [35, 35.054333499999998, 0.33454499999999998, 0.049177180000000008],
 [35, 35.000906399999998, 0.33182709999999999, 0.05440942],
 [35, 35.259924699999999, 0.1197582, 0.013787539999999999],
 [35, 35.038230300000002, -0.49133979999999999, 0.03373201],
 [35, 34.747702500000003, -0.32719330000000002, 0.031145920000000004],
 [35, 35.0590002, -0.77923620000000005, 0.03028142],
 [35, 34.602738199999997, -0.41124490000000002, 0.01212419],
 [35, 34.955750100000003, -0.040631399999999998, 0.080719080000000012],
 [35, 34.696854199999997, -0.650621, 0.01671241],
 [35, 34.8701176, -0.0838531, 0.02994407],
 [35, 35.050668000000002, -0.51947719999999997, 0.011198720000000002],
 [35, 35.045559099999998, -0.43008570000000002, 0.04517630000000001],
 [35, 35.423299399999998, -0.27332699999999999, 0.013814450000000001],
 [35, 34.956440800000003, -0.15437699999999999, 0.027813240000000003],
 [35, 35.161397600000001, 0.93606250000000002, 0.062636990000000003],
 [35, 35.144168200000003, 0.80434919999999999, 0.020249840000000002],
 [35, 34.6942649, 0.0115709, 0.028365999999999999],
 [35, 35.475668800000001, 0.0027131, 0.021155160000000003],
 [35, 34.624174099999998, -0.0067386, 0.031316090000000005],
 [35, 34.822491900000003, 0.32551609999999997, 0.01453868],
 [35, 34.536228899999998, -0.43910949999999999, 0.02313103],
 [35, 35.406803799999999, -0.84592270000000003, 0.01800305],
 [35, 35.044901000000003, 0.3295573, 0.014289990000000001],
 [35, 34.738107900000003, -0.27670070000000002, 0.047069489999999999],
 [35, 34.818818800000003, 0.348275, 0.014820000000000002],
 [35, 35.017623299999997, -0.11141719999999999, 0.08823425],
 [35, 34.724470500000002, -0.61176390000000003, 0.0078948999999999998],
 [35, 35.4738063, 0.8975031, 0.054782909999999997],
 [35, 34.685898399999999, 0.0115938, 0.016057210000000002],
 [35, 35.233150299999998, 0.79623509999999997, 0.13616108999999998],
 [35, 35.1360326, -0.75660890000000003, 0.0226915]]

In [ ]: