3245

pending completion

Build # 3245

Build Type

Pull #379

travis-ci-com

Committed by

web-flow

Commit Message

<a href="https://github.com/angelolab/ark-analysis/commit/<a class=hub.com/angelolab/ark-analysis/commit/<a class="double-link" href="https://git"><a class=hub.com/angelolab/ark-analysis/commit/<a class="double-link" href="https://git"><a class=hub.com/angelolab/ark-analysis/commit/<a class="double-link" href="https://git"><a class=hub.com/angelolab/ark-analysis/commit/2f7018df97a64367ab0a660ba971d3908f5c4507">2f7018df9">&lt;a href=&quot;https://github.com/angelolab/ark-analysis/commit/</a><a class="double-link" href="https://github.com/angelolab/ark-analysis/commit/&lt;a class=&quot;double-link&quot; href=&quot;https://git">&lt;a class=</a>hub.com/angelolab/ark-analysis/commit/&lt;a class=&quot;double-link&quot; href=&quot;https://git">&lt;a class=</a>hub.com/angelolab/ark-analysis/commit/&lt;a class=&quot;double-link&quot; href=&quot;https://git">&lt;a class=</a>hub.com/angelolab/ark-analysis/commit/2f7018df97a64367ab0a660ba971d3908f5c4507">2f7018df9</a><a href="https://github.com/angelolab/ark-analysis/commit/2f7018df97a64367ab0a660ba971d3908f5c4507">&quot;&gt;&amp;lt;a href=&amp;quot;https://github.com/angelolab/ark-analysis/commit/&lt;/a&gt;&lt;a class=&quot;double-link&quot; href=&quot;https://github.com/angelolab/ark-analysis/commit/&amp;lt;a class=&amp;quot;double-link&amp;quot; href=&amp;quot;https://git&quot;&gt;&amp;lt;a class=&lt;/a&gt;hub.com/angelolab/ark-analysis/commit/&amp;lt;a class=&amp;quot;double-link&amp;quot; href=&amp;quot;https://git&quot;&gt;&amp;lt;a class=&lt;/a&gt;hub.com/angelolab/ark-analysis/commit/&amp;lt;a class=&amp;quot;double-link&amp;quot; href=&amp;quot;https://git&quot;&gt;&amp;lt;a class=&lt;/a&gt;hub.com/angelolab/ark-analysis/commit/&amp;lt;a class=&amp;quot;double-link&amp;quot; href=&amp;quot;https://git&quot;&gt;&amp;lt;a class=&lt;/a&gt;hub.com/angelolab/ark-analysis/commit/&amp;lt;a class=&amp;quot;double-link&amp;quot; href=&amp;quot;https://git&quot;&gt;&amp;lt;a class=&lt;/a&gt;hub.com/angelolab/ark-analysis/commit/2f7018df97a64367ab0a660ba971d3908f5c4507&quot;&gt;2f7018df9&lt;/a&gt;&lt;a href=&quot;https://github.com/angelolab/ark-analysis/commit/2f7018df97a64367ab0a660ba971d3908f5c4507&quot;&gt;&amp;quot;&amp;gt;&amp;amp;lt;a href=&amp;amp;quot;https://github.com/angelolab/ark-analysis/commit/&amp;lt;/a&amp;gt;&amp;lt... (continued)

Pull Request Pull Request #379: Add framework for subsetting SOM

Run Details

54 of 54 new or added lines in 1 file covered. (100.0%)

1220 of 1247 relevant lines covered (97.83%)

0.98 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

72.0

/ark/phenotyping/som_utils.py

import os
import subprocess

import feather
import numpy as np
import pandas as pd
import xarray as xr
import scipy.ndimage as ndimage
from skimage.io import imread

import ark.settings as settings
from ark.utils import load_utils
from ark.utils import misc_utils


def create_pixel_matrix(img_xr, seg_labels, base_dir,
                        pre_dir='pixel_mat_preprocessed',
                        sub_dir='pixel_mat_subsetted', fovs=None,
                        blur_factor=2, subset_proportion=0.1, seed=42):
    """Preprocess the images for FlowSOM clustering and creates a pixel-level matrix

    Saves preprocessed data to pre_dir and subsetted data to sub_dir

    Args:
        img_xr (xarray.DataArray):
            Array representing image data for each fov
        seg_labels (xarray.DataArray):
            Array representing segmentation labels for each image
        base_dir (str):
            Name of the directory to save the pixel files to
        pre_dir (str):
            Name of the directory which contains the preprocessed pixel data
        sub_dir (str):
            The name of the directory containing the subsetted pixel data
        fovs (list):
            List of fovs to subset over, if None selects all
        blur_factor (int):
            The sigma to set for the Gaussian blur
        subset_proportion (float):
            The proportion of pixels to take from each fov
        seed (int):
            The random seed to set for subsetting
    """

    # if the subset_proportion specified is out of range
    if subset_proportion <= 0 or subset_proportion > 1:
        raise ValueError('Invalid subset percentage entered: must be in (0, 1]')

    # if the base directory doesn't exist
    if not os.path.exists(base_dir):
        raise FileNotFoundError("base_dir %s does not exist" % base_dir)

    # create pre_dir if it doesn't already exist
    if not os.path.exists(os.path.join(base_dir, pre_dir)):
        os.mkdir(os.path.join(base_dir, pre_dir))

    # create sub_dir if it doesn't already exist
    if not os.path.exists(os.path.join(base_dir, sub_dir)):
        os.mkdir(os.path.join(base_dir, sub_dir))

    # set fovs to all if None
    if fovs is None:
        fovs = img_xr.fovs.values

    # verify that the fovs provided are valid
    misc_utils.verify_in_list(fovs=fovs, image_fovs=img_xr.fovs.values)

    # iterate over fovs
    for fov in fovs:
        # subset img_xr with only the fov we're looking for, and cast to float32
        img_data_blur = img_xr.loc[fov, ...].values.astype(np.float32)

        # for each marker, compute the Gaussian blur
        for marker in range(len(img_xr.channels.values)):
            img_data_blur[:, :, marker] = ndimage.gaussian_filter(img_data_blur[:, :, marker],
                                                                  sigma=blur_factor)

        # flatten each image
        pixel_mat = img_data_blur.reshape(-1, len(img_xr.channels.values))

        # convert into a dataframe
        pixel_mat = pd.DataFrame(pixel_mat, columns=img_xr.channels.values)

        # assign metadata about each entry
        pixel_mat['fov'] = fov
        pixel_mat['row_index'] = np.repeat(range(img_data_blur.shape[0]), img_data_blur.shape[1])
        pixel_mat['column_index'] = np.tile(range(img_data_blur.shape[0]), img_data_blur.shape[1])

        # assign segmentation label
        seg_labels_flat = seg_labels.loc[fov, ...].values.flatten()
        pixel_mat['segmentation_label'] = seg_labels_flat

        # needed for selecting just channel columns
        non_num_cols = ['fov', 'row_index', 'column_index', 'segmentation_label']
        chan_cols = [col for col in pixel_mat.columns.values if col not in non_num_cols]

        # remove any rows that sum to 0
        pixel_mat = pixel_mat.loc[pixel_mat.loc[:, chan_cols].sum(axis=1) != 0, :]

        # normalize each row by total marker counts to convert into frequencies
        pixel_mat.loc[:, chan_cols] = pixel_mat.loc[:, chan_cols].div(
            pixel_mat.loc[:, chan_cols].sum(axis=1), axis=0)

        # subset the pixel matrix for training
        pixel_mat_subset = pixel_mat.sample(frac=subset_proportion, random_state=seed)

        # write complete dataset to feather, needed for cluster assignment
        feather.write_dataframe(pixel_mat,
                                os.path.join(base_dir,
                                             pre_dir,
                                             fov + ".feather"),
                                compression='uncompressed')

        # write subseted dataset to hdf5, needed for training
        feather.write_dataframe(pixel_mat_subset,
                                os.path.join(base_dir,
                                             sub_dir,
                                             fov + ".feather"),
                                compression='uncompressed')


def train_som(fovs, channels, base_dir,
              sub_dir='pixel_mat_subsetted', weights_name='weights.feather', num_passes=1):
    """Run the SOM training on the subsetted pixel data.

    Saves weights to base_dir/weights_name.

    Args:
        fovs (list):
            The list of fovs to subset on
        channels (list):
            The list of markers to subset on
        base_dir (str):
            The path to the data directory
        sub_dir (str):
            The name of the subsetted data directory
        weights_name (str):
            The name of the weights file
        num_passes (int):
            The number of training passes to make through the dataset
    """

    subsetted_path = os.path.join(base_dir, sub_dir)
    weights_path = os.path.join(base_dir, weights_name)

    # if path to the subsetted file does not exist
    if not os.path.exists(subsetted_path):
        raise FileNotFoundError('Pixel subsetted directory %s does not exist in base_dir %s' %
                                (sub_dir, base_dir))

    # run som_train.R
    process_args = ['Rscript', '/create_som_matrix.R', ','.join(fovs), ','.join(channels),
                    str(num_passes), subsetted_path, weights_path]
    process = subprocess.Popen(process_args, stdout=subprocess.PIPE, stderr=subprocess.STDOUT)

    # continuously poll the process for output/error to display in Jupyter notebook
    while True:
        # convert from byte string
        output = process.stdout.readline().decode('utf-8')

        # if the output is nothing and the process is done, break
        if process.poll() is not None:
            break
        if output:
            print(output.strip())


def cluster_pixels(fovs, channels, base_dir, pre_dir='pixel_mat_preprocessed',
                   weights_name='weights.feather', cluster_dir='pixel_mat_clustered'):
    """Uses trained weights to assign cluster labels on full pixel data

    Saves data with cluster labels to cluster_dir

    Args:
        fovs (list):
            The list of fovs to subset on
        channels (list):
            The list of markers to subset on
        base_dir (str):
            The path to the data directory
        pre_dir (str):
            Name of the directory which contains the preprocessed pixel data,
            defaults to pixel_mat_preprocessed
        weights_name (str):
            The name of the weights file
        cluster_dir (str):
            The name of the directory to write the clustered data
    """

    preprocessed_path = os.path.join(base_dir, pre_dir)
    weights_path = os.path.join(base_dir, weights_name)
    clustered_path = os.path.join(base_dir, cluster_dir)

    # if path to the preprocessed file does not exist
    if not os.path.exists(preprocessed_path):
        raise FileNotFoundError('Pixel preprocessed directory %s does not exist in base_dir %s' %
                                (pre_dir, base_dir))

    # if path to the weights file does not exist
    if not os.path.exists(weights_path):
        raise FileNotFoundError('Weights file %s does not exist in base_dir %s' %
                                (weights_name, base_dir))

    # ensure the weights columns are the same as the pixel data columns
    weights = feather.read_dataframe(os.path.join(base_dir, weights_name))
    sample_fov = feather.read_dataframe(os.path.join(base_dir, pre_dir, fovs[0] + '.feather'))

    misc_utils.verify_same_elements(weights_columns=weights.columns.values,
                                    pixel_data_columns=sample_fov[channels].columns.values)

    # make the clustered dir if it does not exist
    if not os.path.exists(clustered_path):
        os.mkdir(clustered_path)

    process_args = ['Rscript', '/run_trained_som.R', ','.join(fovs), ','.join(channels),
                    preprocessed_path, weights_path, clustered_path]

    process = subprocess.Popen(process_args, stdout=subprocess.PIPE, stderr=subprocess.STDOUT)

    # continuously poll the process for output/error so it gets displayed in the Jupyter notebook
    while True:
        # convert from byte string
        output = process.stdout.readline().decode('utf-8')

        # if the output is nothing and the process is done, break
        if process.poll() is not None:
            break
        if output:
            print(output.strip())

1	import os	1✔
2	import subprocess	1✔
3
4	import feather	1✔
5	import numpy as np	1✔
6	import pandas as pd	1✔
7	import xarray as xr	1✔
8	import scipy.ndimage as ndimage	1✔
9	from skimage.io import imread	1✔
10
11	import ark.settings as settings	1✔
12	from ark.utils import load_utils	1✔
13	from ark.utils import misc_utils	1✔
14
15
16	def create_pixel_matrix(img_xr, seg_labels, base_dir,	1✔
17	pre_dir='pixel_mat_preprocessed',
18	sub_dir='pixel_mat_subsetted', fovs=None,
19	blur_factor=2, subset_proportion=0.1, seed=42):
20	"""Preprocess the images for FlowSOM clustering and creates a pixel-level matrix
21
22	Saves preprocessed data to pre_dir and subsetted data to sub_dir
23
24	Args:
25	img_xr (xarray.DataArray):
26	Array representing image data for each fov
27	seg_labels (xarray.DataArray):
28	Array representing segmentation labels for each image
29	base_dir (str):
30	Name of the directory to save the pixel files to
31	pre_dir (str):
32	Name of the directory which contains the preprocessed pixel data
33	sub_dir (str):
34	The name of the directory containing the subsetted pixel data
35	fovs (list):
36	List of fovs to subset over, if None selects all
37	blur_factor (int):
38	The sigma to set for the Gaussian blur
39	subset_proportion (float):
40	The proportion of pixels to take from each fov
41	seed (int):
42	The random seed to set for subsetting
43	"""
44
45	# if the subset_proportion specified is out of range
46	if subset_proportion <= 0 or subset_proportion > 1:	1✔
47	raise ValueError('Invalid subset percentage entered: must be in (0, 1]')	1✔
48
49	# if the base directory doesn't exist
50	if not os.path.exists(base_dir):	1✔
51	raise FileNotFoundError("base_dir %s does not exist" % base_dir)	1✔
52
53	# create pre_dir if it doesn't already exist
54	if not os.path.exists(os.path.join(base_dir, pre_dir)):	1✔
55	os.mkdir(os.path.join(base_dir, pre_dir))	1✔
56
57	# create sub_dir if it doesn't already exist
58	if not os.path.exists(os.path.join(base_dir, sub_dir)):	1✔
59	os.mkdir(os.path.join(base_dir, sub_dir))	1✔
60
61	# set fovs to all if None
62	if fovs is None:	1✔
63	fovs = img_xr.fovs.values	1✔
64
65	# verify that the fovs provided are valid
66	misc_utils.verify_in_list(fovs=fovs, image_fovs=img_xr.fovs.values)	1✔
67
68	# iterate over fovs
69	for fov in fovs:	1✔
70	# subset img_xr with only the fov we're looking for, and cast to float32
71	img_data_blur = img_xr.loc[fov, ...].values.astype(np.float32)	1✔
72
73	# for each marker, compute the Gaussian blur
74	for marker in range(len(img_xr.channels.values)):	1✔
75	img_data_blur[:, :, marker] = ndimage.gaussian_filter(img_data_blur[:, :, marker],	1✔
76	sigma=blur_factor)
77
78	# flatten each image
79	pixel_mat = img_data_blur.reshape(-1, len(img_xr.channels.values))	1✔
80
81	# convert into a dataframe
82	pixel_mat = pd.DataFrame(pixel_mat, columns=img_xr.channels.values)	1✔
83
84	# assign metadata about each entry
85	pixel_mat['fov'] = fov	1✔
86	pixel_mat['row_index'] = np.repeat(range(img_data_blur.shape[0]), img_data_blur.shape[1])	1✔
87	pixel_mat['column_index'] = np.tile(range(img_data_blur.shape[0]), img_data_blur.shape[1])	1✔
88
89	# assign segmentation label
90	seg_labels_flat = seg_labels.loc[fov, ...].values.flatten()	1✔
91	pixel_mat['segmentation_label'] = seg_labels_flat	1✔
92
93	# needed for selecting just channel columns
94	non_num_cols = ['fov', 'row_index', 'column_index', 'segmentation_label']	1✔
95	chan_cols = [col for col in pixel_mat.columns.values if col not in non_num_cols]	1✔
96
97	# remove any rows that sum to 0
98	pixel_mat = pixel_mat.loc[pixel_mat.loc[:, chan_cols].sum(axis=1) != 0, :]	1✔
99
100	# normalize each row by total marker counts to convert into frequencies
101	pixel_mat.loc[:, chan_cols] = pixel_mat.loc[:, chan_cols].div(	1✔
102	pixel_mat.loc[:, chan_cols].sum(axis=1), axis=0)
103
104	# subset the pixel matrix for training
105	pixel_mat_subset = pixel_mat.sample(frac=subset_proportion, random_state=seed)	1✔
106
107	# write complete dataset to feather, needed for cluster assignment
108	feather.write_dataframe(pixel_mat,	1✔
109	os.path.join(base_dir,
110	pre_dir,
111	fov + ".feather"),
112	compression='uncompressed')
113
114	# write subseted dataset to hdf5, needed for training
115	feather.write_dataframe(pixel_mat_subset,	1✔
116	os.path.join(base_dir,
117	sub_dir,
118	fov + ".feather"),
119	compression='uncompressed')
120
121
122	def train_som(fovs, channels, base_dir,	1✔
123	sub_dir='pixel_mat_subsetted', weights_name='weights.feather', num_passes=1):
124	"""Run the SOM training on the subsetted pixel data.
125
126	Saves weights to base_dir/weights_name.
127
128	Args:
129	fovs (list):
130	The list of fovs to subset on
131	channels (list):
132	The list of markers to subset on
133	base_dir (str):
134	The path to the data directory
135	sub_dir (str):
136	The name of the subsetted data directory
137	weights_name (str):
138	The name of the weights file
139	num_passes (int):
140	The number of training passes to make through the dataset
141	"""
142
143	subsetted_path = os.path.join(base_dir, sub_dir)	1✔
144	weights_path = os.path.join(base_dir, weights_name)	1✔
145
146	# if path to the subsetted file does not exist
147	if not os.path.exists(subsetted_path):	1✔
148	raise FileNotFoundError('Pixel subsetted directory %s does not exist in base_dir %s' %	1✔
149	(sub_dir, base_dir))
150
151	# run som_train.R
152	process_args = ['Rscript', '/create_som_matrix.R', ','.join(fovs), ','.join(channels),	×
153	str(num_passes), subsetted_path, weights_path]
154	process = subprocess.Popen(process_args, stdout=subprocess.PIPE, stderr=subprocess.STDOUT)	×
155
156	# continuously poll the process for output/error to display in Jupyter notebook
157	while True:	×
158	# convert from byte string
159	output = process.stdout.readline().decode('utf-8')	×
160
161	# if the output is nothing and the process is done, break
162	if process.poll() is not None:	×
163	break	×
164	if output:	×
165	print(output.strip())	×
166
167
168	def cluster_pixels(fovs, channels, base_dir, pre_dir='pixel_mat_preprocessed',	1✔
169	weights_name='weights.feather', cluster_dir='pixel_mat_clustered'):
170	"""Uses trained weights to assign cluster labels on full pixel data
171
172	Saves data with cluster labels to cluster_dir
173
174	Args:
175	fovs (list):
176	The list of fovs to subset on
177	channels (list):
178	The list of markers to subset on
179	base_dir (str):
180	The path to the data directory
181	pre_dir (str):
182	Name of the directory which contains the preprocessed pixel data,
183	defaults to pixel_mat_preprocessed
184	weights_name (str):
185	The name of the weights file
186	cluster_dir (str):
187	The name of the directory to write the clustered data
188	"""
189
190	preprocessed_path = os.path.join(base_dir, pre_dir)	1✔
191	weights_path = os.path.join(base_dir, weights_name)	1✔
192	clustered_path = os.path.join(base_dir, cluster_dir)	1✔
193
194	# if path to the preprocessed file does not exist
195	if not os.path.exists(preprocessed_path):	1✔
196	raise FileNotFoundError('Pixel preprocessed directory %s does not exist in base_dir %s' %	1✔
197	(pre_dir, base_dir))
198
199	# if path to the weights file does not exist
200	if not os.path.exists(weights_path):	1✔
201	raise FileNotFoundError('Weights file %s does not exist in base_dir %s' %	1✔
202	(weights_name, base_dir))
203
204	# ensure the weights columns are the same as the pixel data columns
205	weights = feather.read_dataframe(os.path.join(base_dir, weights_name))	×
206	sample_fov = feather.read_dataframe(os.path.join(base_dir, pre_dir, fovs[0] + '.feather'))	×
207
208	misc_utils.verify_same_elements(weights_columns=weights.columns.values,	×
209	pixel_data_columns=sample_fov[channels].columns.values)
210
211	# make the clustered dir if it does not exist
212	if not os.path.exists(clustered_path):	×
213	os.mkdir(clustered_path)	×
214
215	process_args = ['Rscript', '/run_trained_som.R', ','.join(fovs), ','.join(channels),	×
216	preprocessed_path, weights_path, clustered_path]
217
218	process = subprocess.Popen(process_args, stdout=subprocess.PIPE, stderr=subprocess.STDOUT)	×
219
220	# continuously poll the process for output/error so it gets displayed in the Jupyter notebook
221	while True:	×
222	# convert from byte string
223	output = process.stdout.readline().decode('utf-8')	×
224
225	# if the output is nothing and the process is done, break
226	if process.poll() is not None:	×
227	break	×
228	if output:	×
229	print(output.strip())	×

angelolab / ark-analysis / 3245

Source File Press 'n' to go to next uncovered line, 'b' for previous

Source File
Press 'n' to go to next uncovered line, 'b' for previous