diff --git a/README.md b/README.md index b5db9f7..0057800 100644 --- a/README.md +++ b/README.md @@ -5,17 +5,8 @@ This repository tracks progress for making dendrograms with [Heat](https://github.com/helmholtz-analytics/heat). The goal is to compute dendrograms like in [Astrodendro](https://github.com/dendrograms/astrodendro), but built on heat to be faster and parallel. -## Roadmap - - [ ] Develop concepts - - [ ] Define I/O - - [ ] Find vectorized algorithms for clustering local data - - [ ] Compute global clustering by merging local clusters - - [ ] Develop implementations - - [ ] Implement hierarchical clustering in heat - - [ ] Implement interface to Astrodendro - ## Setup -In order to install all the necessary dependencies, please install micromamba and then run: +In order to install all the necessary dependencies, please install micromamba and then run (in this directory): ```bash micromamba env create -f environment.yml @@ -23,6 +14,64 @@ micromamba activate dendro pip install -e . ``` +Feel free to substitute micromamba for your favorite package manager such as conda. + +## Compute a distributed dendrogram +Find below the code from the file ``examples/getting_started.py``, which shows you everything you need to use the currently best performing distributed dendrogram implementation. + +```python +import numpy as np +import heat as ht +import logging + +from astrodendro.dendrogram import Dendrogram +from dendro import DistributedDendrogram + +# get some random sample data +ht.random.seed(7353) +data_ht = ht.random.rand(32, 64, split=0) +data_np = data_ht.numpy() + +# if you load your own data, make sure to cast it to heat like +# data_ht = ht.array(data_np, split=0) + +# define parameters for dendrogram computation +params = {'min_npix': 3, 'min_delta':1e-2, 'min_value': 1e-1} + +# compute astrodendro dendrogram as normal +d_astrodendro = Dendrogram.compute(data_np, **params) + +# now, on to the distributed dendrogram! +# first, cast data to heat array and define split dimension + +# define extra parameters for the local dendrogram computation in distributed dendrogram +extra_params = {'min_npix_loc': params['min_npix'], 'min_delta_loc': params['min_delta']} + +# comment the below line to get no output or use `logging.INFO` for some output or `logging.DEBUG` for all the output +logger = logging.getLogger("Dendrogram").setLevel(logging.INFO) + +# simulate the parallel dendrogram computation (on non-distributed data) +d_pseudo_parallel = DistributedDendrogram.compute_pseudo_parallel(data_np, ntasks=4, **params, **extra_params) + +# actually compute in parallel +# the number of tasks is determined in the command line via `mpirun -np getting_started.py` +d_parallel = DistributedDendrogram.compute(data_ht, **params, **extra_params) + +# all tasks have the same dendrograms now, subsequent processing / saving to file should happen on a single task only +if ht.comm.rank == 0: + import matplotlib.pyplot as plt + from dendro.analysis import compare_dendrograms + + fig, ax = plt.subplots() + for label, d in zip(['pseudo parallel', f'{ht.comm.size} tasks'], [d_pseudo_parallel, d_parallel]): + overlap = compare_dendrograms(d_astrodendro, d) + ax.plot(overlap.keys(), overlap.values(), label=label, ls='-' if label=='pseudo parallel' else '--') + ax.set_xlabel('level') + ax.set_ylabel('overlap') + ax.legend() + plt.show() +``` + ## Generate Jupyter notebooks from scripts The scripts are somewhat literate and intended to be viewed as Jupyter notebooks. To this end run commands like: diff --git a/examples/getting_started.py b/examples/getting_started.py new file mode 100644 index 0000000..c02077e --- /dev/null +++ b/examples/getting_started.py @@ -0,0 +1,61 @@ +import heat as ht +import logging + +from astrodendro.dendrogram import Dendrogram +from dendro import DistributedDendrogram + +# get some random sample data +ht.random.seed(7353) +data_ht = ht.random.rand(32, 64, split=0) +data_np = data_ht.numpy() + +# if you load your own data, make sure to cast it to heat like +# data_ht = ht.array(data_np, split=0) + +# define parameters for dendrogram computation +params = {"min_npix": 3, "min_delta": 1e-2, "min_value": 1e-1} + +# compute astrodendro dendrogram as normal +d_astrodendro = Dendrogram.compute(data_np, **params) + +# now, on to the distributed dendrogram! +# first, cast data to heat array and define split dimension + +# define extra parameters for the local dendrogram computation in distributed dendrogram +extra_params = { + "min_npix_loc": params["min_npix"], + "min_delta_loc": params["min_delta"], +} + +# comment the below line to get no output or use `logging.INFO` for some output or `logging.DEBUG` for all the output +logger = logging.getLogger("Dendrogram").setLevel(logging.INFO) + +# simulate the parallel dendrogram computation (on non-distributed data) +d_pseudo_parallel = DistributedDendrogram.compute_pseudo_parallel( + data_np, ntasks=4, **params, **extra_params +) + +# actually compute in parallel +# the number of tasks is determined in the command line via `mpirun -np getting_started.py` +d_parallel = DistributedDendrogram.compute(data_ht, **params, **extra_params) + +# all tasks have the same dendrograms now, subsequent processing / saving to file should happen on a single task only +if ht.comm.rank == 0: + import matplotlib.pyplot as plt + from dendro.analysis import compare_dendrograms + + fig, ax = plt.subplots() + for label, d in zip( + ["pseudo parallel", f"{ht.comm.size} tasks"], [d_pseudo_parallel, d_parallel] + ): + overlap = compare_dendrograms(d_astrodendro, d) + ax.plot( + overlap.keys(), + overlap.values(), + label=label, + ls="-" if label == "pseudo parallel" else "--", + ) + ax.set_xlabel("level") + ax.set_ylabel("overlap") + ax.legend() + plt.show() diff --git a/src/dendro/distributed_dendrogram_v3.py b/src/dendro/distributed_dendrogram_v3.py index 4b612bf..aa41259 100644 --- a/src/dendro/distributed_dendrogram_v3.py +++ b/src/dendro/distributed_dendrogram_v3.py @@ -69,7 +69,7 @@ def compute( def make_output_astrodendro_compatible(self, is_independent): t0 = perf_counter() - self.logger.info("Start making compatible with astrodendro") + self.logger.debug("Start making compatible with astrodendro") if isinstance(self.data, ht.DNDarray): self.data = self.data.numpy() @@ -92,7 +92,7 @@ def make_output_astrodendro_compatible(self, is_independent): } t1 = perf_counter() - self.logger.info( + self.logger.debug( f"Finished making compatible with astrodendro after {t1 - t0:.2e}s" ) @@ -153,7 +153,7 @@ def _compute_single_local_dendrogram(self, local_data, split_dim=0, **kwargs): structure._vmax = np.max(structure._values) t1 = perf_counter() - self.logger.info( + self.logger.debug( f"Isolated {len(local_dendrogram._structures_dict) - len_before_border_removal} border structures in {t1 - t0:.2e}s" ) @@ -163,14 +163,14 @@ def compute_local_dendrogram(self, **kwargs): data = self.data comm = data.comm - self.logger.info( + self.logger.debug( f"Start computing local dendrogram with local data of shape {data.lshape}" ) local_dendrogram = self._compute_single_local_dendrogram( data.larray.numpy(), split_dim=data.split, **kwargs ) - self.logger.info("Adding offsets to structures") + self.logger.debug("Adding offsets to structures") if data.is_distributed(): # add offsets to local indices _, offsets = data.counts_displs() @@ -181,7 +181,7 @@ def compute_local_dendrogram(self, **kwargs): else: for structure in local_dendrogram.all_structures: structure._indices = np.array(structure._indices) - self.logger.info("Finished adding offsets to structures") + self.logger.debug("Finished adding offsets to structures") return local_dendrogram @@ -207,7 +207,7 @@ def compute_local_dendrogram_pseudo_parallel(self, data, ntasks, **kwargs): return local_dendrograms def communicate_structures(self, local_dendrogram): - self.logger.info( + self.logger.debug( f"Starting to communicate {len(local_dendrogram)} local structures" ) t0 = perf_counter() @@ -230,7 +230,7 @@ def communicate_structures(self, local_dendrogram): structures += [Structure(idx=me[0], indices=me[1], values=me[2])] t1 = perf_counter() - self.logger.info(f"Finished communicating structures in {t1 - t0:.2e}s") + self.logger.debug(f"Finished communicating structures in {t1 - t0:.2e}s") return structures @staticmethod @@ -281,7 +281,7 @@ def merge_structures(self, to_merge, merge_into): merge_into._vmax = max([merge_into._vmax, to_merge._vmax]) merge_into._smallest_index = np.min(merge_into._indices) self.index_map[*to_merge._indices.T] = merge_into.idx - self.logger.info( + self.logger.debug( f"Merged {len(to_merge._values)} values between {to_merge._vmin:.2f} and {to_merge._vmax:.2f} into existing structure {merge_into.idx}, which now has {len(merge_into._values)} values between {merge_into._vmin:.2f} and {merge_into._vmax:.2f}" ) @@ -306,7 +306,7 @@ def split_structure(self, structure, split_at): structure._vmin = np.min(structure._values) structure._vmax = np.max(structure._values) - self.logger.info( + self.logger.debug( f"Split structure {structure.idx} at {split_at:.2f}. Remaining top part has {len(structure._values)} values between {structure._vmin:.2f} and {structure._vmax:.2f} and {len(structure._children)} children, bottom part has {len(bottom_part._values)} values between {bottom_part._vmin:.2f} to {bottom_part._vmax:.2f}." ) return structure, bottom_part @@ -334,7 +334,7 @@ def insert_structure(structures, to_insert): structures = DistributedDendrogramV3.insert_structure_within( structures, to_insert ) - DistributedDendrogramV3.logger.info( + DistributedDendrogramV3.logger.debug( f"Inserted structure with {len(to_insert._values)} values between {to_insert._vmin:.2f} and {to_insert._vmax:.2f} into list of {len(structures)} remaining structures." ) return structures @@ -379,7 +379,7 @@ def merge_individual_structure( ) self.index_map[*leaf._indices.T] = leaf.idx merged_structures.append(leaf) - self.logger.info( + self.logger.debug( f"Created new leaf with index {leaf.idx} and {len(leaf._values)} values between {leaf._vmin:.2f} and {leaf._vmax:.2f}." ) elif len(adjacent_structures) == 1: # merge into existing structure @@ -406,7 +406,7 @@ def merge_individual_structure( adjacent_structures.remove(structure) if len(merge) > 0: - self.logger.info( + self.logger.debug( f"Structures {[me.idx for me in merge]} are insignificant. {len(adjacent_structures)} adjacent structures left" ) @@ -427,13 +427,13 @@ def merge_individual_structure( belongs_to = branch self.index_map[*branch._indices.T] = branch.idx merged_structures.append(branch) - self.logger.info( + self.logger.debug( f"Created branch with index {branch.idx} and {len(branch._values)} values between {branch._vmin:.2f} and {branch._vmax:.2f} and {len(branch._children)} children : {[me.idx for me in branch._children]}." ) # merge insignificant structures if len(merge) > 0: - self.logger.info( + self.logger.debug( f"Merging insignificant structure(s) {[m.idx for m in merge]} into structure {belongs_to.idx}" ) for m in merge: @@ -473,7 +473,7 @@ def compute_from_structures(self, structures, is_independent=None): while len(structures) > 0: self._iterations += 1 self.logger.info( - f"--- Iteration {self._iterations}. Merged {len(merged_structures)}, {len(structures)} left." + f"--- Iteration {self._iterations:d}: Merged {len(merged_structures)}, {len(structures)} left." ) to_merge = structures.pop(0) @@ -483,7 +483,7 @@ def compute_from_structures(self, structures, is_independent=None): to_merge, merged_structures, self.index_map ) - self.logger.info( + self.logger.debug( f"Merging structure with {len(to_merge._values)} values between {to_merge._vmin:.2f} and {to_merge._vmax:.2f} with {len(adjacent_structures)} adjacent structures: {[me.idx for me in adjacent_structures]}." )