In [1]:
from ipywidgets import widgets
from IPython.display import display
from mpl_toolkits.mplot3d import axes3d
from collections import namedtuple

import csv
import re
import matplotlib
import time

import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline

font = {'weight' : 'bold',
        'size'   : 18}

matplotlib.rc('font', **font)

In [2]:
fields = [
    'xstart', 'ystart', 'zstart',
    'xstop', 'ystop', 'zstop', 'res',
    'total_pix', 'mask_core_pix', 'str_mask_fraction',
    'syn_count', 'str_density', 'syn_count_unmasked'
]
raw = open('./merged_mask_v6_layerall_DSP.csv', 'r').readlines()
reader = csv.reader(raw)

In [3]:
data = [[float(col) for col in row] for row in reader]

3D Plot

Scatter plot with sliding threshold


In [4]:
dim = 11
lower_initial = 1.05# np.mean([d[dim] for d in data])
upper_initial = 1.12 # 2*np.mean([d[dim] for d in data])
upper_max = 4*np.mean([d[dim] for d in data])
alpha_initial = 0.2
colour_initial = '#000000'
marker_initial = 10 # worth noting: using dim 11 marker scaling is subtle since the range is small
fig_title = "3D Distribution of Synapses in Space"

In [5]:
def plotting3d(thresh_min, thresh_max, alph, scale, colour=colour_initial, title=fig_title):
    fig0 = plt.figure(figsize=(12,9))
    ax = fig0.gca(projection="3d")
    if len(colour) < 6 or not re.search(r'^#(?:[0-9a-fA-F]{3}){1,2}$', colour):
        colour="#000000"
    new_d = np.array([d for d in data if d[dim] > thresh_min and d[11] < thresh_max])
    ax.scatter3D(
        [d[0] for d in new_d],
        [d[1] for d in new_d],
        [d[2] for d in new_d],
        s=[scale*10*d[11] for d in new_d],
        alpha=alph,
        linewidth=0,
        c=colour
    )
    print "Range of Data:", np.min([d[11] for d in data]), "-", np.max([d[dim] for d in data])
    print "Number of Samples Currently in View:", len(new_d) , "(of ", str(len(data))+")"
    ax.set_xticks([np.min(new_d[:,0]), np.max(new_d[:,0])])
    ax.set_yticks([np.min(new_d[:,1]), np.max(new_d[:,1])])
    ax.set_zticks([np.min(new_d[:,2]), np.max(new_d[:,2])])
    ax.set_xlabel('X')
    ax.set_ylabel('Y')
    ax.set_zlabel('Z')
    ax.set_title(title)
    print "Last updated at:", time.ctime()
    plt.show()

thresh_min, thresh_max, alph, marker = (widgets.FloatSlider(description="Lower Threshold Scaling",
                                                    min=0,
                                                    max=upper_max,
                                                    value=lower_initial,
                                                    step=upper_max/100.),
                                        widgets.FloatSlider(description="Upper Threshold Scaling",
                                                    min=0,
                                                    max=upper_max,
                                                    value=upper_initial,
                                                    step=upper_max/100.),
                                        widgets.FloatSlider(description="Alpha Scaling", 
                                                    min=0,
                                                    max=1,
                                                    value=alpha_initial,
                                                    step=0.01),
                                        widgets.FloatSlider(description="Marker Scaling", 
                                                    min=0.1,
                                                    max=15,
                                                    value=marker_initial,
                                                    step=0.1))
set_min, set_max = (widgets.FloatText(description="Set Lower Threshold",
                                      value=lower_initial),
                    widgets.FloatText(description="Set Upper Threshold",
                                         value=upper_initial))
dl = widgets.jsdlink((thresh_min, 'value'), (thresh_max, 'min'))
d2 = widgets.jsdlink((set_min, 'value'), (thresh_min, 'value'))
d3 = widgets.jsdlink((set_max, 'value'), (thresh_max, 'value'))

display(set_min, set_max)
w = widgets.interact(plotting3d, thresh_min=thresh_min,
                     thresh_max=thresh_max, alph=alph,
                     scale=marker)

def on_min_change(change):
    thresh_min.notify_change(change)
set_min.observe(on_min_change, names='value')

def on_max_change(change):
    thresh_max.notify_change(change)
set_max.observe(on_max_change, names='value')


Range of Data: 0.0 - 5.756
Number of Samples Currently in View: 10455 (of  101649)
Last updated at: Fri Jun 23 10:26:03 2017

MIP


In [6]:
buff=20
fig_title = 'Maximum Intensity Projection of Synapses\n\n'
alpha_initial = 0.1
marker_initial = 12
jitter_initial = 0.65

In [7]:
def projections(thresh_min, thresh_max, alph, scale, jitter, colour=colour_initial, title=fig_title):
    new_d = np.array([d for d in data if d[dim] > thresh_min and d[dim] < thresh_max])
    max_jitter = jitter*(np.max(new_d[:,2]) - np.min(new_d[:,2]))/8.
    ## XY
    fig = plt.figure(figsize=(15,6))
    ax = plt.subplot(1,3,1)
    plt.hold(True)
    plt.scatter(
        [d[0] for d in new_d],
        [d[1] for d in new_d],
        s=[10*scale*d[dim] for d in new_d],
        alpha=alph,
        c=colour
    )
    plt.xlim([np.min(new_d[:,0])-buff, np.max(new_d[:,0])+buff])
    plt.ylim([np.min(new_d[:,1])-buff, np.max(new_d[:,1])+buff])
    plt.xticks([np.min(new_d[:,0]), np.max(new_d[:,0])])
    plt.yticks([np.min(new_d[:,1]), np.max(new_d[:,1])])
    plt.text(np.min(new_d[:,0])-buff-260, (np.min(new_d[:,1])-buff + np.max(new_d[:,1])+buff)/2, 'Y')
    plt.text((np.min(new_d[:,0])-buff + np.max(new_d[:,0])+buff)/2, np.min(new_d[:,1])-buff-150, 'X')
    
    ## XZ
    ax = plt.subplot(1,3,2)
    plt.hold(True)
    plt.scatter(
        [d[0] for d in new_d],
        [d[2]+max_jitter*np.random.rand(1) for d in new_d],
        s=[10*scale*d[dim] for d in new_d],
        alpha=alph,
        c=colour
    )
    plt.xlim([np.min(new_d[:,0])-buff, np.max(new_d[:,0])+buff])
    plt.ylim([np.min(new_d[:,2])-buff, np.max(new_d[:,2])+buff])
    plt.xticks([np.min(new_d[:,0]), np.max(new_d[:,0])])
    plt.yticks([np.min(new_d[:,2]), np.max(new_d[:,2])])
    plt.text(np.min(new_d[:,0])-buff-250, (np.min(new_d[:,2])-buff + np.max(new_d[:,2])+buff)/2, 'Z')
    plt.text((np.min(new_d[:,0])-buff + np.max(new_d[:,0])+buff)/2, np.min(new_d[:,2])-buff-70, 'X')
    
    ## YZ
    ax = plt.subplot(1,3,3)
    plt.hold(True)
    plt.scatter(
        [d[1] for d in new_d],
        [d[2]+max_jitter*np.random.rand(1) for d in new_d],
        s=[10*scale*d[dim] for d in new_d],
        alpha=alph,
        c=colour
    )
    plt.xlim([np.min(new_d[:,1])-buff, np.max(new_d[:,1])+buff])
    plt.ylim([np.min(new_d[:,2])-buff, np.max(new_d[:,2])+buff])
    plt.xticks([np.min(new_d[:,1]), np.max(new_d[:,1])])
    plt.yticks([np.min(new_d[:,2]), np.max(new_d[:,2])])
    plt.text(np.min(new_d[:,1])-buff-200, (np.min(new_d[:,2])-buff + np.max(new_d[:,2])+buff)/2, 'Z')
    plt.text((np.min(new_d[:,1])-buff + np.max(new_d[:,1])+buff)/2, np.min(new_d[:,2])-buff-70, 'Y')
    plt.tight_layout()
    plt.suptitle(title)
    plt.subplots_adjust(top=0.85)
    print "Last updated at:", time.ctime()
    plt.show()
    
thresh_min, thresh_max, alph, marker, jitter = (widgets.FloatSlider(description="Lower Threshold Scaling",
                                                    min=0,
                                                    max=upper_max,
                                                    value=lower_initial,
                                                    step=upper_max/100.),
                                        widgets.FloatSlider(description="Upper Threshold Scaling",
                                                    min=0,
                                                    max=upper_max,
                                                    value=upper_initial,
                                                    step=upper_max/100.),
                                        widgets.FloatSlider(description="Alpha Scaling", 
                                                    min=0,
                                                    max=1,
                                                    value=alpha_initial,
                                                    step=0.01),
                                        widgets.FloatSlider(description="Marker Scaling", 
                                                    min=0.1,
                                                    max=15,
                                                    value=marker_initial,
                                                    step=0.1),
                                        widgets.FloatSlider(description="Maximum Jitter",
                                                    min=0,
                                                    max=1,
                                                    value=jitter_initial,
                                                    step=0.05))

set_min, set_max = (widgets.FloatText(description="Set Lower Threshold",
                                      value=lower_initial),
                    widgets.FloatText(description="Set Upper Threshold",
                                         value=upper_initial))
dl = widgets.jsdlink((thresh_min, 'value'), (thresh_max, 'min'))
d2 = widgets.jsdlink((set_min, 'value'), (thresh_min, 'value'))
d3 = widgets.jsdlink((set_max, 'value'), (thresh_max, 'value'))

display(set_min, set_max)
w = widgets.interact(projections, thresh_min=thresh_min,
                     thresh_max=thresh_max, alph=alph,
                     scale=marker, jitter=jitter)

def on_min_change(change):
    thresh_min.notify_change(change)
set_min.observe(on_min_change, names='value')

def on_max_change(change):
    thresh_max.notify_change(change)
set_max.observe(on_max_change, names='value')


Last updated at: Fri Jun 23 10:26:40 2017

2D Max Projection Heatmap


In [8]:
buff=20
fig_title = 'Maximum Intensity Projection of Synapses\n\n'
jitter_initial = 0.65
n_bins = 20

In [9]:
def projections(thresh_min, thresh_max, jitter, bins=n_bins, title=fig_title):
    new_d = np.array([d for d in data if d[dim] > thresh_min and d[dim] < thresh_max])
    max_jitter = jitter*(np.max(new_d[:,2]) - np.min(new_d[:,2]))/8.
    ## XY
    fig = plt.figure(figsize=(15,6))
    ax = plt.subplot(1,3,1)
    plt.hold(True)
    heatmap, xedges, yedges = np.histogram2d([d[0] for d in new_d],
                                             [d[1] for d in new_d], bins=bins)
    extent = [xedges[0], xedges[-1], yedges[0], yedges[-1]]
    plt.imshow(heatmap.T, extent=extent, aspect='auto', origin='lower')
    plt.xticks([np.min(new_d[:,0]), np.max(new_d[:,0])])
    plt.yticks([np.min(new_d[:,1]), np.max(new_d[:,1])])
    plt.ylabel('Y')
    plt.xlabel('X')
    
    ## XZ
    ax = plt.subplot(1,3,2)
    plt.hold(True)
    heatmap, xedges, yedges = np.histogram2d([d[0] for d in new_d],
                                             [d[2]+max_jitter*np.random.rand(1)-max_jitter*np.random.rand(1)/2.
                                              for d in new_d], bins=bins)
    extent = [xedges[0], xedges[-1], yedges[0], yedges[-1]]
    plt.imshow(heatmap.T, extent=extent, aspect='auto', origin='lower')
    plt.xticks([np.min(new_d[:,0]), np.max(new_d[:,0])])
    plt.yticks([np.min(new_d[:,2]), np.max(new_d[:,2])])
    plt.ylim([np.min(new_d[:,2])-buff, np.max(new_d[:,2])+buff])
    plt.ylabel('Z')
    plt.xlabel('X')
    
    ## YZ
    ax = plt.subplot(1,3,3)
    plt.hold(True)
    heatmap, xedges, yedges = np.histogram2d([d[1] for d in new_d],
                                             [d[2]+max_jitter*np.random.rand(1)-max_jitter*np.random.rand(1)/2.
                                              for d in new_d], bins=bins)
    extent = [xedges[0], xedges[-1], yedges[0], yedges[-1]]
    plt.imshow(heatmap.T, extent=extent, aspect='auto', origin='lower')
    plt.xticks([np.min(new_d[:,1]), np.max(new_d[:,1])])
    plt.yticks([np.min(new_d[:,2]), np.max(new_d[:,2])])
    plt.ylim([np.min(new_d[:,2])-buff, np.max(new_d[:,2])+buff])
    plt.ylabel('Z')
    plt.xlabel('Y')
    plt.tight_layout()
    plt.suptitle(title)
    plt.subplots_adjust(top=0.85)
    print "Last updated at:", time.ctime()
    plt.show()
    
thresh_min, thresh_max, jitter = (widgets.FloatSlider(description="Lower Threshold Scaling",
                                                    min=0,
                                                    max=upper_max,
                                                    value=lower_initial,
                                                    step=upper_max/100.),
                                        widgets.FloatSlider(description="Upper Threshold Scaling",
                                                    min=0,
                                                    max=upper_max,
                                                    value=upper_initial,
                                                    step=upper_max/100.),
                                        widgets.FloatSlider(description="Maximum Jitter",
                                                    min=0,
                                                    max=1,
                                                    value=jitter_initial,
                                                    step=0.05))

set_min, set_max = (widgets.FloatText(description="Set Lower Threshold",
                                      value=lower_initial),
                    widgets.FloatText(description="Set Upper Threshold",
                                         value=upper_initial))
dl = widgets.jsdlink((thresh_min, 'value'), (thresh_max, 'min'))
d2 = widgets.jsdlink((set_min, 'value'), (thresh_min, 'value'))
d3 = widgets.jsdlink((set_max, 'value'), (thresh_max, 'value'))

display(set_min, set_max)
w = widgets.interact(projections, thresh_min=thresh_min,
                     thresh_max=thresh_max,
                     jitter=jitter)

def on_min_change(change):
    thresh_min.notify_change(change)
set_min.observe(on_min_change, names='value')

def on_max_change(change):
    thresh_max.notify_change(change)
set_max.observe(on_max_change, names='value')


Last updated at: Fri Jun 23 10:26:49 2017

In [ ]: