Kevin developed this course and a great part of the slide are his work.
|
|
|
|
Experimental Design finding the right technique, picking the right dyes and samples has stayed relatively consistent, better techniques lead to more demanding scientists.
Management storing, backing up, setting up databases, these processes have become easier and more automated as data magnitudes have increased
Measurements the actual acquisition speed of the data has increased wildly due to better detectors, parallel measurement, and new higher intensity sources
Post Processing this portion has is the most time-consuming and difficult and has seen minimal improvements over the last years
If you looked at one 1000 x 1000 sized image
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
plt.matshow(np.random.uniform(size = (1000, 1000)),
cmap = 'viridis')
<matplotlib.image.AxesImage at 0x1c25d72490>
every second, it would take you
# assuming 16 bit images and a 'metric' terabyte
time_per_tb=1e12/(1000*1000*16/8) / (60*60)
print("%04.1f hours to view a terabyte" % (time_per_tb))
138.9 hours to view a terabyte
|
import json, pandas as pd
course_df = pd.read_json('../common/schedule.json')
course_df['Date'] = course_df['Lecture'].map(lambda x: x.split('-')[0])
course_df['Title'] = course_df['Lecture'].map(lambda x: x.split('-')[-1])
course_df[['Date', 'Title', 'Description']]
Date | Title | Description | |
---|---|---|---|
0 | 20th February | Introduction and Workflows | Basic overview of the course, introduction to ... |
1 | 27th February | Ground Truth: Building and Augmenting Datasets | Examples of large datasets, how they were buil... |
2 | 5th March | Image Enhancement | Overview of what techniques are available for ... |
3 | 12th March | Basic Segmentation, Discrete Binary Structures | How to convert images into structures, startin... |
4 | 19th March | Advanced Segmentation | More advanced techniques for extracting struct... |
5 | 26th March | Supervised Problems and Segmentation | More advanced techniques for extracting struct... |
6 | 2nd April | Analyzing Single Objects, Shape, and Texture | The analysis and characterization of single st... |
7 | 9th April | Analyzing Complex Objects and Distributions | What techniques are available to analyze more ... |
8 | 16th April | Dynamic Experiments | Performing tracking and registration in dynami... |
9 | 23rd April | Statistics, Prediction, and Reproducibility | Making a statistical analysis from quantified ... |
10 | 30th April | Guest Lectures | How Roche does Microscopy at Scale with High C... |
11 | 7th May | Scaling Up / Big Data | Performing large scale analyses on clusters an... |
12 | 14th May | Project Presentations | You present your projects |
course_df[['Title', 'Description', 'Applications']][3:6].T
3 | 4 | 5 | |
---|---|---|---|
Title | Basic Segmentation, Discrete Binary Structures | Advanced Segmentation | Supervised Problems and Segmentation |
Description | How to convert images into structures, startin... | More advanced techniques for extracting struct... | More advanced techniques for extracting struct... |
Applications | Identify cells from noise, background, and dust | Identifying fat and ice crystals in ice cream ... | Identifying fat and ice crystals in ice cream ... |
course_df[['Title', 'Description', 'Applications']][6:9].T
6 | 7 | 8 | |
---|---|---|---|
Title | Analyzing Single Objects, Shape, and Texture | Analyzing Complex Objects and Distributions | Dynamic Experiments |
Description | The analysis and characterization of single st... | What techniques are available to analyze more ... | Performing tracking and registration in dynami... |
Applications | Count cells and determine their average shape ... | Seperate clumps of cells, analyze vessel netwo... | Turning a video of foam flow into metrics like... |
course_df[['Title', 'Description', 'Applications']][9:12].T
9 | 10 | 11 | |
---|---|---|---|
Title | Statistics, Prediction, and Reproducibility | Guest Lectures | Scaling Up / Big Data |
Description | Making a statistical analysis from quantified ... | How Roche does Microscopy at Scale with High C... | Performing large scale analyses on clusters an... |
Applications | Determine if/how different a cancerous cell is... | Robust analysis of millions of images for maki... | Performing large scale analyses using ETHs clu... |
course_df[['Title', 'Description', 'Applications']][12:13].T
12 | |
---|---|
Title | Project Presentations |
Description | You present your projects |
Applications |
A very abstract definition:
In most cases this is a 2- or 3-dimensional position (x,y,z coordinates) and a numeric value (intensity)
The world is | The computer needs |
|
|
import numpy as np
import matplotlib.pyplot as plt
from skimage.transform import resize
img=np.load('../common/data/wood.npy');
plt.figure(figsize=[15,7])
plt.subplot(2,3,1); plt.imshow(img); plt.title('Original')
downsize = 2; plt.subplot(2,3,2); plt.imshow(resize(img,(img.shape[0] // downsize, img.shape[1] // downsize), anti_aliasing=False)); plt.title('Downsize {0}x{0}'.format(downsize))
downsize = 32; plt.subplot(2,3,3); plt.imshow(resize(img,(img.shape[0] // downsize, img.shape[1] // downsize),anti_aliasing=False)); plt.title('Downsize {0}x{0}'.format(downsize))
levels = 16; plt.subplot(2,3,5); plt.imshow(np.floor(img*levels)); plt.title('{0} Levels'.format(levels));
levels = 4 ; plt.subplot(2,3,6); plt.imshow(np.floor(img*levels)); plt.title('{0} Levels'.format(levels));
import numpy as np
basic_image = np.random.choice(range(100), size = (5,5))
xx, yy = np.meshgrid(range(basic_image.shape[1]), range(basic_image.shape[0]))
image_df = pd.DataFrame(dict(x = xx.ravel(),
y = yy.ravel(),
Intensity = basic_image.ravel()))
image_df[['x', 'y', 'Intensity']].head(5)
x | y | Intensity | |
---|---|---|---|
0 | 0 | 0 | 33 |
1 | 1 | 0 | 78 |
2 | 2 | 0 | 5 |
3 | 3 | 0 | 73 |
4 | 4 | 0 | 22 |
import matplotlib.pyplot as plt
plt.imshow(basic_image, cmap = 'gray')
plt.colorbar();
The next step is to apply a color map (also called lookup table, LUT) to the image
fig, ax1 = plt.subplots(1,1)
plot_image = ax1.matshow(basic_image, cmap = 'Blues')
plt.colorbar(plot_image)
for _, c_row in image_df.iterrows():
ax1.text(c_row['x'], c_row['y'], s = '%02d' % c_row['Intensity'], fontdict = dict(color = 'r'))
Color maps can be arbitrarily defined based on how we would like to visualize the information in the image
fig, ax1 = plt.subplots(1,1)
plot_image = ax1.matshow(basic_image, cmap = 'jet')
plt.colorbar(plot_image)
<matplotlib.colorbar.Colorbar at 0x1c5c4eb990>
fig, ax1 = plt.subplots(1,1)
plot_image = ax1.matshow(basic_image, cmap = 'hot')
plt.colorbar(plot_image)
<matplotlib.colorbar.Colorbar at 0x1c508ad8d0>
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
from matplotlib import cm
#from colorspacious import cspace_converter
from collections import OrderedDict
cmaps = OrderedDict()
cmaps['Perceptually Uniform Sequential'] = [
'viridis', 'plasma', 'inferno', 'magma', 'cividis']
cmaps['Sequential'] = [
'Greys', 'Purples', 'Blues', 'Greens', 'Oranges', 'Reds',
'YlOrBr', 'YlOrRd', 'OrRd', 'PuRd', 'RdPu', 'BuPu',
'GnBu', 'PuBu', 'YlGnBu', 'PuBuGn', 'BuGn', 'YlGn']
cmaps['Sequential (2)'] = [
'binary', 'gist_yarg', 'gist_gray', 'gray', 'bone', 'pink',
'spring', 'summer', 'autumn', 'winter', 'cool', 'Wistia',
'hot', 'afmhot', 'gist_heat', 'copper']
cmaps['Diverging'] = [
'PiYG', 'PRGn', 'BrBG', 'PuOr', 'RdGy', 'RdBu',
'RdYlBu', 'RdYlGn', 'Spectral', 'coolwarm', 'bwr', 'seismic']
cmaps['Cyclic'] = ['twilight', 'twilight_shifted', 'hsv']
cmaps['Qualitative'] = ['Pastel1', 'Pastel2', 'Paired', 'Accent',
'Dark2', 'Set1', 'Set2', 'Set3',
'tab10', 'tab20', 'tab20b', 'tab20c']
cmaps['Miscellaneous'] = [
'flag', 'prism', 'ocean', 'gist_earth', 'terrain', 'gist_stern',
'gnuplot', 'gnuplot2', 'CMRmap', 'cubehelix', 'brg',
'gist_rainbow', 'rainbow', 'jet', 'nipy_spectral', 'gist_ncar']
nrows = max(len(cmap_list) for cmap_category, cmap_list in cmaps.items())
gradient = np.linspace(0, 1, 256)
gradient = np.vstack((gradient, gradient))
def plot_color_gradients(cmap_category, cmap_list, nrows):
fig, axes = plt.subplots(nrows=nrows)
fig.subplots_adjust(top=0.95, bottom=0.01, left=0.2, right=0.99)
axes[0].set_title(cmap_category + ' colormaps', fontsize=14)
for ax, name in zip(axes, cmap_list):
ax.imshow(gradient, aspect='auto', cmap=plt.get_cmap(name))
pos = list(ax.get_position().bounds)
x_text = pos[0] - 0.01
y_text = pos[1] + pos[3]/2.
fig.text(x_text, y_text, name, va='center', ha='right', fontsize=10)
# Turn off *all* ticks & spines, not just the ones with colormaps.
for ax in axes:
ax.set_axis_off()
for cmap_category, cmap_list in cmaps.items():
plot_color_gradients(cmap_category, cmap_list, nrows)
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
xlin = np.linspace(0, 1, 100)
colors = ['Red','Green','Blue']
plt.figure(figsize=[15,4])
for i in np.arange(0,3) :
plt.subplot(1,3,i+1)
plt.scatter(xlin,
plt.cm.hot(xlin)[:,i],
c = plt.cm.hot(xlin),label="hot")
plt.scatter(xlin,
plt.cm.Blues(xlin)[:,i],
c = plt.cm.Blues(xlin),label="blues")
plt.scatter(xlin,
plt.cm.jet(xlin)[:,i],
c = plt.cm.jet(xlin),label='jet')
plt.xlabel('Intensity');
plt.ylabel('{0} Component'.format(colors[i]));
These transformations can also be non-linear as is the case of the graph below where the mapping between the intensity and the color is a $\log$ relationship meaning the the difference between the lower values is much clearer than the higher ones
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
xlin = np.logspace(-2, 5, 500)
log_xlin = np.log10(xlin)
norm_xlin = (log_xlin-log_xlin.min())/(log_xlin.max()-log_xlin.min())
fig, ax1 = plt.subplots(1,1)
ax1.scatter(xlin, plt.cm.hot(norm_xlin)[:,0],
c = plt.cm.hot(norm_xlin))
ax1.scatter(xlin, plt.cm.hot(xlin/xlin.max())[:,0],
c = plt.cm.hot(norm_xlin))
ax1.set_xscale('log');ax1.set_xlabel('Intensity');ax1.set_ylabel('Red Component');
On a real image the difference is even clearer
%matplotlib inline
import matplotlib.pyplot as plt
from skimage.io import imread
fig, (ax1, ax2, ax3) = plt.subplots(1,3, figsize = (12, 4))
in_img = imread('../common/figures/bone-section.png')[:,:,0].astype(np.float32)
ax1.imshow(in_img, cmap = 'gray');
ax1.set_title('grayscale LUT');
ax2.imshow(in_img, cmap = 'hot');
ax2.set_title('hot LUT');
ax3.imshow(np.log2(in_img+1), cmap = 'gray');
ax3.set_title('grayscale-log LUT');
For a 3D image, the position or spatial component has a 3rd dimension (z if it is a spatial, or t if it is a movie)
Volume | Time series |
import numpy as np
vol_image = np.arange(27).reshape((3,3,3))
print(vol_image)
[[[ 0 1 2] [ 3 4 5] [ 6 7 8]] [[ 9 10 11] [12 13 14] [15 16 17]] [[18 19 20] [21 22 23] [24 25 26]]]
This can then be rearranged from a table form into an array form and displayed as a series of slices
%matplotlib inline
import matplotlib.pyplot as plt
from skimage.util import montage as montage2d
print(montage2d(vol_image, fill = 0))
plt.matshow(montage2d(vol_image, fill = 0), cmap = 'jet')
[[ 0 1 2 9 10 11] [ 3 4 5 12 13 14] [ 6 7 8 15 16 17] [18 19 20 0 0 0] [21 22 23 0 0 0] [24 25 26 0 0 0]]
<matplotlib.image.AxesImage at 0x1c599ef2d0>
In the images thus far, we have had one value per position, but there is no reason there cannot be multiple values. In fact this is what color images are (red, green, and blue) values and even 4 channels with transparency (alpha) as a different. For clarity we call the dimensionality of the image the number of dimensions in the spatial position, and the depth the number in the value.
import pandas as pd
from itertools import product
import numpy as np
base_df = pd.DataFrame([dict(x = x, y = y) for x,y in product(range(5), range(5))])
base_df['Intensity'] = np.random.uniform(0, 1, 25)
base_df['Transparency'] = np.random.uniform(0, 1, 25)
base_df.head(5)
x | y | Intensity | Transparency | |
---|---|---|---|---|
0 | 0 | 0 | 0.051095 | 0.618102 |
1 | 0 | 1 | 0.096233 | 0.314752 |
2 | 0 | 2 | 0.226459 | 0.127602 |
3 | 0 | 3 | 0.491454 | 0.273547 |
4 | 0 | 4 | 0.830064 | 0.974110 |
This can then be rearranged from a table form into an array form and displayed as a series of slices
%matplotlib inline
import matplotlib.pyplot as plt
fig, (ax1, ax2) = plt.subplots(1, 2)
ax1.scatter(base_df['x'], base_df['y'], c = plt.cm.gray(base_df['Intensity']), s = 1000)
ax1.set_title('Intensity')
ax2.scatter(base_df['x'], base_df['y'], c = plt.cm.gray(base_df['Transparency']), s = 1000)
ax2.set_title('Transparency')
Text(0.5, 1.0, 'Transparency')
fig, (ax1) = plt.subplots(1, 1)
ax1.scatter(base_df['x'], base_df['y'], c = plt.cm.jet(base_df['Intensity']), s = 1000*base_df['Transparency'])
ax1.set_title('Intensity')
Text(0.5, 1.0, 'Intensity')
At each point in the image (black dot), instead of having just a single value, there is an entire spectrum. A selected group of these (red dots) are shown to illustrate the variations inside the sample. While certainly much more complicated, this still constitutes and image and requires the same sort of techniques to process correctly.
%matplotlib inline
import matplotlib.pyplot as plt
import pandas as pd
from skimage.io import imread
import os
raw_img = imread(os.path.join('..', 'common', 'data', 'raw.jpg'))
im_pos = pd.read_csv(os.path.join('..', 'common', 'data', 'impos.csv'), header = None)
im_pos.columns = ['x', 'y']
fig, ax1 = plt.subplots(1,1, figsize = (8, 8));
ax1.imshow(raw_img);
ax1.scatter(im_pos['x'], im_pos['y'], s = 1, c = 'blue');
full_df = pd.read_csv(os.path.join('..', 'common', 'data', 'full_img.csv')).query('wavenum<1200')
print(full_df.shape[0], 'rows')
full_df.head(5)
210750 rows
x | y | wavenum | val | |
---|---|---|---|---|
0 | 168.95 | 358.8 | 750 | 527.571102 |
1 | 168.95 | 358.8 | 753 | 459.778584 |
2 | 168.95 | 358.8 | 756 | 406.337255 |
3 | 168.95 | 358.8 | 759 | 341.858123 |
4 | 168.95 | 358.8 | 762 | 246.645673 |
full_df['g_x'] = pd.cut(full_df['x'], 5)
full_df['g_y'] = pd.cut(full_df['y'], 5)
fig, m_axs = plt.subplots(5, 5, figsize = (12, 12))
for ((g_x, g_y), c_rows), c_ax in zip(full_df.sort_values(['x','y']).groupby(['g_x', 'g_y']), m_axs.flatten()):
c_ax.plot(c_rows['wavenum'], c_rows['val'], 'r.')
import pandas as pd
from io import StringIO
pd.read_table(StringIO("""Modality\tImpulse Characteristic Response Detection
Light Microscopy White Light Electronic interactions Absorption Film, Camera
Phase Contrast Coherent light Electron Density (Index of Refraction) Phase Shift Phase stepping, holography, Zernike
Confocal Microscopy Laser Light Electronic Transition in Fluorescence Molecule Absorption and reemission Pinhole in focal plane, scanning detection
X-Ray Radiography X-Ray light Photo effect and Compton scattering Absorption and scattering Scintillator, microscope, camera
Neutron Radiography Neutrons Interaction with nucleus Scattering and absorption Scintillator, optics, camera
Ultrasound High frequency sound waves Molecular mobility Reflection and Scattering Transducer
MRI Radio-frequency EM Unmatched Hydrogen spins Absorption and reemission RF coils to detect
Atomic Force Microscopy Sharp Point Surface Contact Contact, Repulsion Deflection of a tiny mirror"""))
Modality | Impulse | Characteristic | Response | Detection | |
---|---|---|---|---|---|
0 | Light Microscopy | White Light | Electronic interactions | Absorption | Film, Camera |
1 | Phase Contrast | Coherent light | Electron Density (Index of Refraction) | Phase Shift | Phase stepping, holography, Zernike |
2 | Confocal Microscopy | Laser Light | Electronic Transition in Fluorescence Molecule | Absorption and reemission | Pinhole in focal plane, scanning detection |
3 | X-Ray Radiography | X-Ray light | Photo effect and Compton scattering | Absorption and scattering | Scintillator, microscope, camera |
4 | Neutron Radiography | Neutrons | Interaction with nucleus | Scattering and absorption | Scintillator, optics, camera |
5 | Ultrasound | High frequency sound waves | Molecular mobility | Reflection and Scattering | Transducer |
6 | MRI | Radio-frequency EM | Unmatched Hydrogen spins | Absorption and reemission | RF coils to detect |
7 | Atomic Force Microscopy | Sharp Point | Surface Contact | Contact, Repulsion | Deflection of a tiny mirror |
%matplotlib inline
import matplotlib.pyplot as plt
import pandas as pd
from skimage.io import imread
from scipy.ndimage import convolve
from skimage.morphology import disk
import numpy as np
import os
bone_img = imread(os.path.join('..', 'common', 'figures', 'tiny-bone.png')).astype(np.float32)
# simulate measured image
conv_kern = np.pad(disk(2), 1, 'constant', constant_values = 0)
meas_img = convolve(bone_img[::-1], conv_kern)
# run deconvolution
dekern = np.fft.ifft2(1/np.fft.fft2(conv_kern))
rec_img = convolve(meas_img, dekern)[::-1]
# show result
fig, (ax_orig, ax1, ax2) = plt.subplots(1,3,
figsize = (12, 4))
ax_orig.imshow(bone_img, cmap = 'bone')
ax_orig.set_title('Original Object')
ax1.imshow(np.real(meas_img), cmap = 'bone')
ax1.set_title('Measurement')
ax2.imshow(np.real(rec_img), cmap = 'bone', vmin = 0, vmax = 255)
ax2.set_title('Reconstructed');
/opt/anaconda3/lib/python3.7/site-packages/numpy/core/_asarray.py:85: ComplexWarning: Casting complex values to real discards the imaginary part return array(a, dtype, copy=False, order=order)
%matplotlib inline
import matplotlib.pyplot as plt
import pandas as pd
from skimage.io import imread
from scipy.ndimage import convolve
from skimage.morphology import disk
import numpy as np
import os
bone_img = imread(os.path.join('..', 'common', 'figures', 'tiny-bone.png')).astype(np.float32)
# simulate measured image
meas_img = np.log10(np.abs(np.fft.fftshift(np.fft.fft2(bone_img))))
print(meas_img.min(), meas_img.max(), meas_img.mean())
fig, (ax1, ax_orig) = plt.subplots(1,2,
figsize = (12, 6))
ax_orig.imshow(bone_img, cmap = 'bone')
ax_orig.set_title('Original Object')
ax1.imshow(meas_img, cmap = 'hot')
ax1.set_title('Measurement');
1.1465129 6.6112556 3.3563936
$s_{ab}$ is the only information you are really interested in, so it is important to remove or correct for the other components
For color (non-monochromatic) images the problem becomes even more complicated $$ \int_{0}^{\infty} {\left[\left([b(x,y,\lambda)*s_{ab}(x,y,\lambda)]\otimes h_{fs}(x,y,\lambda)\right)*h_{op}(x,y,\lambda)\right]*h_{det}(x,y,\lambda)}\mathrm{d}\lambda+d_{dark}(x,y) $$
|
|
Inspired by: imagej-pres
Images are great for qualitative analyses since our brains can quickly interpret them without large programming investements.
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
xlin = np.linspace(-1,1, 3)
xx, yy = np.meshgrid(xlin, xlin)
img_a = 25*np.ones((3,3))
img_b = np.ones((3,3))*75
img_a[1,1] = 50
img_b[1,1] = 50
fig, (ax1, ax2) = plt.subplots(1,2, figsize = (12, 5));
ax1.matshow(img_a, vmin = 0, vmax = 100, cmap = 'bone');
ax2.matshow(img_b, vmin = 0, vmax = 100, cmap = 'bone');
Are the intensities constant in the image?
%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
xlin = np.linspace(-1,1, 10)
xx, yy = np.meshgrid(xlin, xlin)
fig, ax1 = plt.subplots(1,1, figsize = (6, 6))
ax1.matshow(xx, vmin = -1, vmax = 1, cmap = 'bone')
<matplotlib.image.AxesImage at 0x1c5b4853d0>
Science demands repeatability! and really wants reproducability
%matplotlib inline
# stolen from https://gist.github.com/humberto-ortiz/de4b3a621602b78bf90d
import pandas as pd
import matplotlib.pyplot as plt
from io import StringIO
moores_txt=["Id Name Year Count(1000s) Clock(MHz)\n",
"0 MOS65XX 1975 3.51 14\n",
"1 Intel8086 1978 29.00 10\n",
"2 MIPSR3000 1988 120.00 33\n",
"3 AMDAm486 1993 1200.00 40\n",
"4 NexGenNx586 1994 3500.00 111\n",
"5 AMDAthlon 1999 37000.00 1400\n",
"6 IntelPentiumIII 1999 44000.00 1400\n",
"7 PowerPC970 2002 58000.00 2500\n",
"8 AMDAthlon64 2003 243000.00 2800\n",
"9 IntelCore2Duo 2006 410000.00 3330\n",
"10 AMDPhenom 2007 450000.00 2600\n",
"11 IntelCorei7 2008 1170000.00 3460\n",
"12 IntelCorei5 2009 995000.00 3600"]
sio_table = StringIO(''.join(moores_txt)); moore_df = pd.read_table(sio_table, sep = '\s+', index_col = 0);
fig, ax1 = plt.subplots(1,1, figsize = (8, 4)); ax1.semilogy(moore_df['Year'], moore_df['Count(1000s)'], 'b.-', label = '1000s of transitiors'); ax1.semilogy(moore_df['Year'], moore_df['Clock(MHz)'], 'r.-', label = 'Clockspeed (MHz)') ;ax1.legend(loc = 2);
Based on data from https://gist.github.com/humberto-ortiz/de4b3a621602b78bf90d
There are now many more transistors inside a single computer but the processing speed hasn't increased. How can this be?
|