|
_Disclaimer: we're not going to produce rendered images like the above in this
post. These were created with [NVidia
IndeX](https://developer.nvidia.com/index), a completely separate tool chain
from what is being discussed here. This post covers the first step of image
loading._
## Series Overview
A common case in fields that acquire large amounts of imaging data is to write
out smaller acquisitions into many small files. These files can tile a larger
space, sub-sample from a larger time period, and may contain multiple channels.
The acquisition techniques themselves are often state of the art and constantly
pushing the envelope in term of how large a field of view can be acquired, at
what resolution, and what quality.
Once acquired this data presents a number of challenges. Algorithms often
designed and tested to work on very small pieces of this data need to be scaled
up to work on the full dataset. It might not be clear at the outset what will
actually work and so exploration still plays a very big part of the whole
process.
Historically this analytical process has involved a lot of custom code. Often
the analytical process is stitched together by a series of scripts possibly in
several different languages that write various intermediate results to disk.
Thanks to advances in modern tooling these process can be significantly
improved. In this series of blogposts, we will outline ways for image
scientists to leverage different tools to move towards a high level, friendly,
cohesive, interactive analytical pipeline.
## Post Overview
This post in particular focuses on loading and managing large stacks of image
data in parallel from Python.
Loading large image data can be a complex and often unique problem. Different
groups may choose to store this across many files on disk, a commodity or
custom database solution, or they may opt to store it in the cloud. Not all
datasets within the same group may be treated the same for a variety of
reasons. In short, this means loading data is a hard and expensive problem.
Despite data being stored in many different ways, often groups want to reapply
the same analytical pipeline to these datasets. However if the data pipeline is
tightly coupled to a particular way of loading the data for later analytical
steps, it may be very difficult if not impossible to reuse an existing
pipeline. In other words, there is friction between the loading and analysis
steps, which frustrates efforts to make things reusable.
Having a modular and general way to load data makes it easy to present data
stored differently in a standard way. Further having a standard way to present
data to analytical pipelines allows that part of the pipeline to focus on what
it does best, analysis! In general, this should decouple these to components in
a way that improves the experience of users involved in all parts of the
pipeline.
We will use
[image data](https://drive.google.com/drive/folders/13mpIfqspKTIINkfoWbFsVtFF8D7jbTqJ)
generously provided by
[Gokul Upadhyayula](https://scholar.google.com/citations?user=nxwNAEgAAAAJ&hl=en)
at the
[Advanced Bioimaging Center](http://microscopy.berkeley.edu/)
at UC Berkeley and discussed in
[this paper](https://science.sciencemag.org/content/360/6386/eaaq1392)
([preprint](https://www.biorxiv.org/content/10.1101/243352v2)),
though the workloads presented here should work for any kind of imaging data,
or array data generally.
## Load image data with Dask
Let's start again with our image data from the top of the post:
```
$ $ ls /path/to/files/raw/ | head
ex6-2_CamA_ch1_CAM1_stack0000_560nm_0000000msec_0001291795msecAbs_000x_000y_000z_0000t.tif
ex6-2_CamA_ch1_CAM1_stack0001_560nm_0043748msec_0001335543msecAbs_000x_000y_000z_0000t.tif
ex6-2_CamA_ch1_CAM1_stack0002_560nm_0087497msec_0001379292msecAbs_000x_000y_000z_0000t.tif
ex6-2_CamA_ch1_CAM1_stack0003_560nm_0131245msec_0001423040msecAbs_000x_000y_000z_0000t.tif
ex6-2_CamA_ch1_CAM1_stack0004_560nm_0174993msec_0001466788msecAbs_000x_000y_000z_0000t.tif
```
### Load a single sample image with Scikit-Image
To load a single image, we use [Scikit-Image](https://scikit-image.org/):
```python
>>> import glob
>>> filenames = glob.glob("/path/to/files/raw/*.tif")
>>> len(filenames)
597
>>> import imageio
>>> sample = imageio.imread(filenames[0])
>>> sample.shape
(201, 1024, 768)
```
Each filename corresponds to some 3d chunk of a larger image. We can look at a
few 2d slices of this single 3d chunk to get some context.
```python
import matplotlib.pyplot as plt
import skimage.io
plt.figure(figsize=(10, 10))
skimage.io.imshow(sample[:, :, 0])
```
```python
plt.figure(figsize=(10, 10))
skimage.io.imshow(sample[:, 0, :])
```
```python
plt.figure(figsize=(10, 10))
skimage.io.imshow(sample[0, :, :])
```
### Investigate Filename Structure
These are slices from only one chunk of a much larger aggregate image.
Our interest here is combining the pieces into a large image stack.
It is common to see a naming structure in the filenames. Each
filename then may indicate a channel, time step, and spatial location with the
`` being some numeric values (possibly with units). Individual filenames may
have more or less information and may notate it differently than we have.
```
mydata_ch_
|
|
|
|
|
|