run_experiment_apd_all.py 11.8 KB
Newer Older
1 2 3 4
"""
run experiments with Stat-of-the-art methods

Example run:
Jiri Borovec's avatar
Jiri Borovec committed
5
>> nohup python experiments_sta.py > ~/Medical-temp/experiments_APD-sta/nohup.log &
6

Jiri Borovec's avatar
Jiri Borovec committed
7
>> python run_experiment_apd_all.py \
Jiri Borovec's avatar
Jiri Borovec committed
8
    -in /datagrid/Medical/microscopy/drosophila/synthetic_data/atomicPatternDictionary_v1 \
9 10 11 12 13 14
    -out /datagrid/Medical/microscopy/drosophila/TEMPORARY/experiments_APD

>> python run_experiment_apd_all.py \
    -in /datagrid/Medical/microscopy/drosophila/synthetic_data/atomicPatternDictionary_v1 \
    -out /datagrid/Medical/microscopy/drosophila/TEMPORARY/experiments_APDL_synth2
    --method APDL
Jiri Borovec's avatar
Jiri Borovec committed
15

Jiri Borovec's avatar
Jiri Borovec committed
16 17 18 19 20
>> python run_experiment_apd_all.py --type real \
    -in /datagrid/Medical/microscopy/drosophila/TEMPORARY/type_1_segm_reg_binary \
    -out /datagrid/Medical/microscopy/drosophila/TEMPORARY/experiments_APD_real \
    --dataset gene_ssmall

Jiri Borovec's avatar
Jiri Borovec committed
21
Copyright (C) 2015-2016 Jiri Borovec <jiri.borovec@fel.cvut.cz>
22 23
"""

Jiri Borovec's avatar
Jiri Borovec committed
24 25
import os
import sys
Jiri Borovec's avatar
Jiri Borovec committed
26 27
import time
import gc
Jiri Borovec's avatar
Jiri Borovec committed
28 29
import logging

Jiri Borovec's avatar
Jiri Borovec committed
30 31 32
# to suppress all visu, has to be on the beginning
import matplotlib
matplotlib.use('Agg')
33
import numpy as np
Jiri Borovec's avatar
Jiri Borovec committed
34
from sklearn.decomposition import SparsePCA, FastICA, DictionaryLearning, NMF
Jiri Borovec's avatar
Jiri Borovec committed
35
from skimage import segmentation
Jiri Borovec's avatar
Jiri Borovec committed
36
import tqdm
Jiri Borovec's avatar
Jiri Borovec committed
37

Jiri Borovec's avatar
Jiri Borovec committed
38
sys.path.append(os.path.abspath(os.path.join('..','..')))  # Add path to root
Jiri Borovec's avatar
Jiri Borovec committed
39
import experiment_apd as expt_apd
Jiri Borovec's avatar
Jiri Borovec committed
40 41
import pattern_disctionary as ptn_dict
import pattern_weights as ptn_weight
Jiri Borovec's avatar
Jiri Borovec committed
42
import run_experiment_apd_apdl as expt_apdl
Jiri Borovec's avatar
Jiri Borovec committed
43

44

Jiri Borovec's avatar
Jiri Borovec committed
45 46 47
SYNTH_PARAMS = expt_apd.SYNTH_PARAMS
SYNTH_SUB_DATASETS = expt_apd.SYNTH_SUB_DATASETS
SYNTH_PTN_RANGE = expt_apd.SYNTH_PTN_RANGE
Jiri Borovec's avatar
Jiri Borovec committed
48 49
REAL_PARAMS = expt_apd.REAL_PARAMS
NB_PATTERNS_REAL = expt_apd.NB_PATTERNS_REAL
Jiri Borovec's avatar
Jiri Borovec committed
50

51 52

class ExperimentLinearCombineBase(expt_apd.ExperimentAPD):
Jiri Borovec's avatar
Jiri Borovec committed
53 54 55
    """
    State-of-te-Art methods that are based on Linear Combination
    """
Jiri Borovec's avatar
Jiri Borovec committed
56 57

    def _estimate_linear_combination(self, imgs_vec):
Jiri Borovec's avatar
Jiri Borovec committed
58 59 60
        """ perform the estimation of LinComb and set the estimator,
        results and patterns

Jiri Borovec's avatar
Jiri Borovec committed
61
        :param imgs_vec: np.array<nb_imgs, height*width>
Jiri Borovec's avatar
Jiri Borovec committed
62 63
        :return:
        """
Jiri Borovec's avatar
Jiri Borovec committed
64 65 66
        pass

    def _perform_linear_combination(self, imgs_vec):
Jiri Borovec's avatar
Jiri Borovec committed
67 68
        """ perform the linear combination and reformulate the outputs

Jiri Borovec's avatar
Jiri Borovec committed
69
        :param imgs_vec: np.array<nb_imgs, height*width>
Jiri Borovec's avatar
Jiri Borovec committed
70 71
        :return:
        """
Jiri Borovec's avatar
Jiri Borovec committed
72
        self._estimate_linear_combination(imgs_vec)
Jiri Borovec's avatar
Jiri Borovec committed
73 74
        logging.debug('fitting parameters: %s', repr(self.estimator.get_params()))
        logging.debug('number of iteration: %i', self.estimator.n_iter_)
Jiri Borovec's avatar
Jiri Borovec committed
75 76 77

        atlas_ptns = self.components.reshape((-1, ) + self.imgs[0].shape)
        rct_vec = np.dot(self.fit_result, self.components)
Jiri Borovec's avatar
Jiri Borovec committed
78
        return atlas_ptns, rct_vec
Jiri Borovec's avatar
Jiri Borovec committed
79

Jiri Borovec's avatar
Jiri Borovec committed
80
    def estim_atlas_as_argmax(self, atlas_ptns):
Jiri Borovec's avatar
Jiri Borovec committed
81 82 83 84 85
        """

        :param atlas_ptns:
        :return: np.array<height, width>
        """
Jiri Borovec's avatar
Jiri Borovec committed
86 87 88 89 90 91 92 93 94
        # take max pattern with max value
        ptn_used = np.sum(np.abs(self.fit_result), axis=0) > 0
        # filter just used patterns
        atlas_ptns = atlas_ptns[ptn_used, :]
        # argmax on abs
        atlas = np.argmax(np.abs(atlas_ptns), axis=0) + 1
        atlas_sum = np.sum(np.abs(atlas_ptns), axis=0)
        # filter small values
        atlas[atlas_sum < 1e-3] = 0
Jiri Borovec's avatar
Jiri Borovec committed
95
        assert atlas.shape == atlas_ptns[0].shape
Jiri Borovec's avatar
Jiri Borovec committed
96 97 98
        return atlas

    def estim_atlas_as_unique_sum(self, atlas_ptns):
Jiri Borovec's avatar
Jiri Borovec committed
99 100 101 102 103
        """

        :param atlas_ptns:
        :return: np.array<height, width>
        """
Jiri Borovec's avatar
Jiri Borovec committed
104 105 106 107 108
        atlas = np.sum(np.abs(atlas_ptns), axis=0)
        atlas /= np.max(atlas)
        atlas = np.array(atlas * len(np.unique(atlas)), dtype=np.int)
        return atlas

109
    def _convert_patterns_to_atlas(self, atlas_ptns):
Jiri Borovec's avatar
Jiri Borovec committed
110 111 112
        """ convert the estimated patterns into a reasonable atlas

        :param atlas_ptns: np.array<nb_patterns, w*h>
Jiri Borovec's avatar
Jiri Borovec committed
113
        :return: np.array<height, width>
Jiri Borovec's avatar
Jiri Borovec committed
114
        """
Jiri Borovec's avatar
Jiri Borovec committed
115 116
        atlas = self.estim_atlas_as_argmax(atlas_ptns)
        # atlas = self.estim_atlas_as_unique_sum(atlas_ptns)
Jiri Borovec's avatar
Jiri Borovec committed
117
        self.atlas = segmentation.relabel_sequential(atlas)[0]
118

Jiri Borovec's avatar
Jiri Borovec committed
119
    def _binarize_img_reconstruction(self, img_rct, thr=0.5):
Jiri Borovec's avatar
Jiri Borovec committed
120 121 122 123 124 125
        """ binarise the reconstructed images to be sure again binary

        :param img_rct: np.array<nb_spl, w, h>
        :param thr: float
        :return:
        """
126 127
        img_rct_bin = [None] * img_rct.shape[0]
        for i, im in enumerate(img_rct.tolist()):
Jiri Borovec's avatar
Jiri Borovec committed
128
            img_rct_bin[i] = np.array(np.asarray(im) > thr, dtype=np.int)
129 130
        return img_rct_bin

Jiri Borovec's avatar
Jiri Borovec committed
131
    def _perform_once(self, v):
Jiri Borovec's avatar
Jiri Borovec committed
132 133 134 135 136
        """ perform one experiment

        :param v: value
        :return:
        """
Jiri Borovec's avatar
Jiri Borovec committed
137
        self.params[self.iter_var_name] = v
Jiri Borovec's avatar
Jiri Borovec committed
138
        name_posix = '_{}_{}'.format(self.iter_var_name, v)
Jiri Borovec's avatar
Jiri Borovec committed
139 140
        if isinstance(self.params['nb_samples'], float):
            self.params['nb_samples'] = int(len(self.imgs) * self.params['nb_samples'])
141
        imgs_vec = np.array([np.ravel(im) for im in self.imgs[:self.params['nb_samples']]])
Jiri Borovec's avatar
Jiri Borovec committed
142 143
        atlas_ptns, rct_vec = self._perform_linear_combination(imgs_vec)
        # img_rct = rct_vec.reshape(np.asarray(self.imgs[:self.params['nb_samples']]).shape)
Jiri Borovec's avatar
Jiri Borovec committed
144 145 146
        self._convert_patterns_to_atlas(atlas_ptns)
        self._export_atlas(name_posix)
        w_bins = [ptn_weight.weights_image_atlas_overlap_major(img, self.atlas)
Jiri Borovec's avatar
Jiri Borovec committed
147
                  for img in self.imgs[:self.params['nb_samples']]]
Jiri Borovec's avatar
Jiri Borovec committed
148 149 150 151 152
        self.w_bins = np.array(w_bins)
        self._export_coding(name_posix)
        img_rct = ptn_dict.reconstruct_samples(self.atlas, self.w_bins)
        stat = self._compute_statistic_gt(img_rct)
        stat[self.iter_var_name] = v
Jiri Borovec's avatar
Jiri Borovec committed
153
        return stat
154 155


156
class ExperimentFastICA_base(ExperimentLinearCombineBase):
Jiri Borovec's avatar
Jiri Borovec committed
157 158 159
    """
    http://scikit-learn.org/stable/modules/generated/sklearn.decomposition.FastICA.html
    """
160

Jiri Borovec's avatar
Jiri Borovec committed
161
    def _estimate_linear_combination(self, imgs_vec):
Jiri Borovec's avatar
Jiri Borovec committed
162 163 164 165
        self.estimator = FastICA(n_components=self.params.get('nb_labels'),
                                 max_iter=self.params.get('max_iter'),
                                 algorithm='deflation',
                                 whiten=True)
Jiri Borovec's avatar
Jiri Borovec committed
166 167
        self.fit_result = self.estimator.fit_transform(imgs_vec)
        self.components = self.estimator.mixing_.T
168 169


170
class ExperimentSparsePCA_base(ExperimentLinearCombineBase):
Jiri Borovec's avatar
Jiri Borovec committed
171 172 173
    """
    http://scikit-learn.org/stable/modules/generated/sklearn.decomposition.SparsePCA.html
    """
174

Jiri Borovec's avatar
Jiri Borovec committed
175
    def _estimate_linear_combination(self, imgs_vec):
Jiri Borovec's avatar
Jiri Borovec committed
176 177 178
        self.estimator = SparsePCA(n_components=self.params.get('nb_labels'),
                                   max_iter=self.params.get('max_iter'),
                                   n_jobs=1)
Jiri Borovec's avatar
Jiri Borovec committed
179 180
        self.fit_result = self.estimator.fit_transform(imgs_vec)
        self.components = self.estimator.components_
181 182


183
class ExperimentDictLearn_base(ExperimentLinearCombineBase):
Jiri Borovec's avatar
Jiri Borovec committed
184 185 186
    """
    http://scikit-learn.org/stable/modules/generated/sklearn.decomposition.DictionaryLearning.html
    """
187

Jiri Borovec's avatar
Jiri Borovec committed
188
    def _estimate_linear_combination(self, imgs_vec):
Jiri Borovec's avatar
Jiri Borovec committed
189 190 191 192 193 194
        self.estimator = DictionaryLearning(fit_algorithm='lars',
                                            transform_algorithm='omp',
                                            split_sign=False,
                                            n_components=self.params.get('nb_labels'),
                                            max_iter=self.params.get('max_iter'),
                                            n_jobs=1)
Jiri Borovec's avatar
Jiri Borovec committed
195 196
        self.fit_result = self.estimator.fit_transform(imgs_vec)
        self.components = self.estimator.components_
197

Jiri Borovec's avatar
Jiri Borovec committed
198

Jiri Borovec's avatar
Jiri Borovec committed
199 200 201 202 203 204 205 206 207 208 209 210 211 212 213 214 215 216
class ExperimentNMF_base(ExperimentLinearCombineBase):
    """
    http://scikit-learn.org/stable/modules/generated/sklearn.decomposition.DictionaryLearning.html
    """

    def _estimate_linear_combination(self, imgs_vec):
        self.estimator = NMF(n_components=self.params.get('nb_labels'),
                             max_iter=self.params.get('max_iter'),
                             init='random')
        self.fit_result = self.estimator.fit_transform(imgs_vec)
        self.components = self.estimator.components_


class ExperimentFastICA(ExperimentFastICA_base, expt_apd.ExperimentAPD_parallel):
    pass


class ExperimentSparsePCA(ExperimentSparsePCA_base, expt_apd.ExperimentAPD_parallel):
217 218 219
    pass


Jiri Borovec's avatar
Jiri Borovec committed
220
class ExperimentDictLearn(ExperimentDictLearn_base, expt_apd.ExperimentAPD_parallel):
221 222 223
    pass


Jiri Borovec's avatar
Jiri Borovec committed
224
class ExperimentNMF(ExperimentNMF_base, expt_apd.ExperimentAPD_parallel):
225 226 227 228
    pass


# standard multiprocessing version
Jiri Borovec's avatar
Jiri Borovec committed
229 230 231 232
METHODS = {
    'PCA': ExperimentFastICA,
    'ICA': ExperimentSparsePCA,
    'DL': ExperimentDictLearn,
Jiri Borovec's avatar
Jiri Borovec committed
233
    'NMF': ExperimentNMF,
Jiri Borovec's avatar
Jiri Borovec committed
234
    'APDL': expt_apdl.ExperimentAPDL,
Jiri Borovec's avatar
Jiri Borovec committed
235
}
Jiri Borovec's avatar
Jiri Borovec committed
236

237 238 239 240 241
# working jut in single thread for pasiisng to image data to prtial jobs
METHODS_BASE = {
    'PCA': ExperimentFastICA_base,
    'ICA': ExperimentSparsePCA_base,
    'DL': ExperimentDictLearn_base,
Jiri Borovec's avatar
Jiri Borovec committed
242
    'NMF': ExperimentNMF_base,
243 244 245
    'APDL': expt_apdl.ExperimentAPDL_base,
}

Jiri Borovec's avatar
Jiri Borovec committed
246 247

def experiments_test(params=SYNTH_PARAMS):
Jiri Borovec's avatar
Jiri Borovec committed
248
    """ simple test of the experiments
Jiri Borovec's avatar
Jiri Borovec committed
249

Jiri Borovec's avatar
Jiri Borovec committed
250
    :param dict_params: {str: value}
Jiri Borovec's avatar
Jiri Borovec committed
251
    """
Jiri Borovec's avatar
Jiri Borovec committed
252 253
    logging.basicConfig(level=logging.DEBUG)

Jiri Borovec's avatar
Jiri Borovec committed
254
    params['nb_runs'] = 2
Jiri Borovec's avatar
Jiri Borovec committed
255
    params['nb_samples'] = 0.5
256

Jiri Borovec's avatar
Jiri Borovec committed
257 258
    for n, cls_expt in METHODS.iteritems():
        logging.info('testing %s by %s', n, cls_expt.__class__)
Jiri Borovec's avatar
Jiri Borovec committed
259
        expt = cls_expt(params)
Jiri Borovec's avatar
Jiri Borovec committed
260
        expt.run(iter_var='case', iter_vals=range(params['nb_runs']))
261 262


Jiri Borovec's avatar
Jiri Borovec committed
263
def experiments_synthetic(params=SYNTH_PARAMS):
Jiri Borovec's avatar
Jiri Borovec committed
264 265
    """ run all experiments

Jiri Borovec's avatar
Jiri Borovec committed
266
    :param dataset: str, name of dataset
Jiri Borovec's avatar
Jiri Borovec committed
267
    :param ptn_range: [int]
Jiri Borovec's avatar
Jiri Borovec committed
268
    """
Jiri Borovec's avatar
Jiri Borovec committed
269
    arg_params = expt_apd.parse_params(params)
Jiri Borovec's avatar
Jiri Borovec committed
270 271
    logging.info('PARAMS: \n%s', '\n'.join(['"{}": \n\t {}'.format(k, v)
                                            for k, v in arg_params.iteritems()]))
Jiri Borovec's avatar
Jiri Borovec committed
272
    params.update(arg_params)
Jiri Borovec's avatar
Jiri Borovec committed
273 274 275
    if not 'method' in params:
        params['method'] = METHODS.keys()

Jiri Borovec's avatar
Jiri Borovec committed
276
    l_params = [params]
Jiri Borovec's avatar
Jiri Borovec committed
277 278
    if isinstance(params['dataset'], list):
        l_params = expt_apd.extend_list_params(l_params, 'dataset', params['dataset'])
279 280
    # l_params = expt_apd.extend_list_params(l_params, 'nb_samples',
    #                                        np.linspace(0.1, 1, 10).tolist())
Jiri Borovec's avatar
Jiri Borovec committed
281

Jiri Borovec's avatar
Jiri Borovec committed
282
    ptn_range = SYNTH_PTN_RANGE[os.path.basename(params['path_in'])]
Jiri Borovec's avatar
Jiri Borovec committed
283

Jiri Borovec's avatar
Jiri Borovec committed
284
    for m in params['method']:
Jiri Borovec's avatar
Jiri Borovec committed
285
        cls_expt = METHODS[m]
286 287
        if params['nb_jobs'] <= 1:
            cls_expt = METHODS_BASE[m]
Jiri Borovec's avatar
Jiri Borovec committed
288
        tqdm_bar = tqdm.tqdm(total=len(l_params))
Jiri Borovec's avatar
Jiri Borovec committed
289 290 291
        for param in l_params:
            param['method'] = m
            expt = cls_expt(param)
Jiri Borovec's avatar
Jiri Borovec committed
292 293
            expt.run(iter_var='nb_labels', iter_vals=ptn_range)
            tqdm_bar.update(1)
Jiri Borovec's avatar
Jiri Borovec committed
294 295
            del expt
            gc.collect(), time.sleep(1)
296 297


Jiri Borovec's avatar
Jiri Borovec committed
298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 330
def experiments_real(params=REAL_PARAMS):
    """ run all experiments

    :param nb_jobs: int
    :param params: {str: value}
    """
    arg_params = expt_apd.parse_params(params)
    logging.info('PARAMS: \n%s', '\n'.join(['"{}": \n\t {}'.format(k, v)
                                            for k, v in arg_params.iteritems()]))
    params.update(arg_params)
    if not 'method' in params:
        params['method'] = METHODS.keys()

    l_params = [params]
    if isinstance(params['dataset'], list):
        l_params = expt_apd.extend_list_params(l_params, 'dataset', params['dataset'])
    logging.debug('list params: %i', len(l_params))

    tqdm_bar = tqdm.tqdm(total=len(l_params))
    for m in params['method']:
        cls_expt = METHODS[m]
        if params['nb_jobs'] <= 1:
            cls_expt = METHODS_BASE[m]
        tqdm_bar = tqdm.tqdm(total=len(l_params))
        for param in l_params:
            param['method'] = m
            expt = cls_expt(param)
            expt.run(gt=False, iter_var='nb_labels', iter_vals=NB_PATTERNS_REAL)
            tqdm_bar.update(1)
            del expt
            gc.collect(), time.sleep(1)


Jiri Borovec's avatar
Jiri Borovec committed
331 332
def main():
    """ main_real entry point """
Jiri Borovec's avatar
Jiri Borovec committed
333
    logging.basicConfig(level=logging.INFO)
Jiri Borovec's avatar
Jiri Borovec committed
334
    logging.info('running...')
335

Jiri Borovec's avatar
Jiri Borovec committed
336
    # experiments_test()
Jiri Borovec's avatar
Jiri Borovec committed
337

Jiri Borovec's avatar
Jiri Borovec committed
338 339 340
    arg_params = expt_apd.parse_params(SYNTH_PARAMS)
    if arg_params['type'] == 'synth':
        experiments_synthetic()
Jiri Borovec's avatar
Jiri Borovec committed
341 342
    elif arg_params['type'] == 'real':
        experiments_real()
343

Jiri Borovec's avatar
Jiri Borovec committed
344
    logging.info('DONE')
Jiri Borovec's avatar
Jiri Borovec committed
345
    # plt.show()
Jiri Borovec's avatar
Jiri Borovec committed
346 347 348 349 350


if __name__ == "__main__":
    main()