diff --git a/autofit/database/model/fit.py b/autofit/database/model/fit.py index c1b6a18f9..1951a92ac 100644 --- a/autofit/database/model/fit.py +++ b/autofit/database/model/fit.py @@ -18,48 +18,26 @@ class Pickle(Base): def __init__(self, **kwargs): super().__init__(**kwargs) - id = sa.Column( - sa.Integer, - primary_key=True - ) + id = sa.Column(sa.Integer, primary_key=True) - name = sa.Column( - sa.String - ) - string = sa.Column( - sa.String - ) - fit_id = sa.Column( - sa.String, - sa.ForeignKey( - "fit.id" - ) - ) - fit = sa.orm.relationship( - "Fit", - uselist=False - ) + name = sa.Column(sa.String) + string = sa.Column(sa.String) + fit_id = sa.Column(sa.String, sa.ForeignKey("fit.id")) + fit = sa.orm.relationship("Fit", uselist=False) @property def value(self): """ The unpickled object """ - if isinstance( - self.string, - str - ): + if isinstance(self.string, str): return self.string - return pickle.loads( - self.string - ) + return pickle.loads(self.string) @value.setter def value(self, value): try: - self.string = pickle.dumps( - value - ) + self.string = pickle.dumps(value) except pickle.PicklingError: pass @@ -67,24 +45,13 @@ def value(self, value): class Info(Base): __tablename__ = "info" - id = sa.Column( - sa.Integer, - primary_key=True - ) + id = sa.Column(sa.Integer, primary_key=True) key = sa.Column(sa.String) value = sa.Column(sa.String) - fit_id = sa.Column( - sa.String, - sa.ForeignKey( - "fit.id" - ) - ) - fit = sa.orm.relationship( - "Fit", - uselist=False - ) + fit_id = sa.Column(sa.String, sa.ForeignKey("fit.id")) + fit = sa.orm.relationship("Fit", uselist=False) def try_none(func): @@ -101,24 +68,13 @@ def wrapper(*args, **kwargs): class NamedInstance(Base): __tablename__ = "named_instance" - id = sa.Column( - sa.Integer, - primary_key=True - ) + id = sa.Column(sa.Integer, primary_key=True) name = sa.Column(sa.String) - instance_id = sa.Column( - sa.Integer, - sa.ForeignKey( - "object.id" - ) - ) + instance_id = sa.Column(sa.Integer, sa.ForeignKey("object.id")) __instance = sa.orm.relationship( - "Object", - uselist=False, - backref="named_instance", - foreign_keys=[instance_id] + "Object", uselist=False, backref="named_instance", foreign_keys=[instance_id] ) @property @@ -131,20 +87,10 @@ def instance(self): @instance.setter def instance(self, instance): - self.__instance = Object.from_object( - instance - ) - - fit_id = sa.Column( - sa.String, - sa.ForeignKey( - "fit.id" - ) - ) - fit = sa.orm.relationship( - "Fit", - uselist=False - ) + self.__instance = Object.from_object(instance) + + fit_id = sa.Column(sa.String, sa.ForeignKey("fit.id")) + fit = sa.orm.relationship("Fit", uselist=False) # noinspection PyProtectedMember @@ -167,56 +113,36 @@ def __getitem__(self, item: str): Raises a KeyError if no such instance exists. """ - return self._get_named_instance( - item - ).instance + return self._get_named_instance(item).instance def __setitem__(self, key: str, value): """ Set an instance for a given name """ try: - named_instance = self._get_named_instance( - key - ) + named_instance = self._get_named_instance(key) except KeyError: - named_instance = NamedInstance( - name=key - ) - self.fit._named_instances.append( - named_instance - ) + named_instance = NamedInstance(name=key) + self.fit._named_instances.append(named_instance) named_instance.instance = value - def _get_named_instance( - self, - item: str - ) -> "NamedInstance": + def _get_named_instance(self, item: str) -> "NamedInstance": """ Retrieve a NamedInstance by its name. """ for named_instance in self.fit._named_instances: if named_instance.name == item: return named_instance - raise KeyError( - f"Instance {item} not found" - ) + raise KeyError(f"Instance {item} not found") class Fit(Base): __tablename__ = "fit" - id = sa.Column( - sa.String, - primary_key=True, - ) - is_complete = sa.Column( - sa.Boolean - ) + id = sa.Column(sa.String, primary_key=True,) + is_complete = sa.Column(sa.Boolean) - _named_instances: List[NamedInstance] = sa.orm.relationship( - "NamedInstance" - ) + _named_instances: List[NamedInstance] = sa.orm.relationship("NamedInstance") @property @try_none @@ -228,45 +154,23 @@ def instance(self): @instance.setter def instance(self, instance): - self.__instance = Object.from_object( - instance - ) + self.__instance = Object.from_object(instance) @property def named_instances(self): - return NamedInstancesWrapper( - self - ) + return NamedInstancesWrapper(self) - _info: List[Info] = sa.orm.relationship( - "Info" - ) + _info: List[Info] = sa.orm.relationship("Info") - def __init__( - self, - **kwargs - ): - super().__init__( - **kwargs - ) + def __init__(self, **kwargs): + super().__init__(**kwargs) - max_log_likelihood = sa.Column( - sa.Float - ) + max_log_likelihood = sa.Column(sa.Float) - parent_id = sa.Column( - sa.String, - sa.ForeignKey( - "fit.id" - ) - ) + parent_id = sa.Column(sa.String, sa.ForeignKey("fit.id")) children: List["Fit"] = sa.orm.relationship( - "Fit", - backref=sa.orm.backref( - 'parent', - remote_side=[id] - ) + "Fit", backref=sa.orm.backref("parent", remote_side=[id]) ) @property @@ -276,13 +180,9 @@ def best_fit(self) -> "Fit": the highest log likelihood. """ if not self.is_grid_search: - raise TypeError( - f"Fit {self.id} is not a grid search" - ) + raise TypeError(f"Fit {self.id} is not a grid search") if len(self.children) == 0: - raise TypeError( - f"Grid search fit {self.id} has no children" - ) + raise TypeError(f"Grid search fit {self.id} has no children") best_fit = None max_log_likelihood = float("-inf") @@ -294,26 +194,14 @@ def best_fit(self) -> "Fit": return best_fit - is_grid_search = sa.Column( - sa.Boolean - ) + is_grid_search = sa.Column(sa.Boolean) - unique_tag = sa.Column( - sa.String - ) - name = sa.Column( - sa.String - ) - path_prefix = sa.Column( - sa.String - ) + unique_tag = sa.Column(sa.String) + name = sa.Column(sa.String) + path_prefix = sa.Column(sa.String) _samples = sa.orm.relationship( - Object, - uselist=False, - foreign_keys=[ - Object.samples_for_id - ] + Object, uselist=False, foreign_keys=[Object.samples_for_id] ) @property @@ -323,29 +211,16 @@ def samples(self) -> Samples: @samples.setter def samples(self, samples): - self._samples = Object.from_object( - samples - ) + self._samples = Object.from_object(samples) @property def info(self): - return { - info.key: info.value - for info - in self._info - } + return {info.key: info.value for info in self._info} @info.setter def info(self, info): if info is not None: - self._info = [ - Info( - key=key, - value=value - ) - for key, value - in info.items() - ] + self._info = [Info(key=key, value=value) for key, value in info.items()] @property @try_none @@ -357,14 +232,9 @@ def model(self) -> AbstractPriorModel: @model.setter def model(self, model: AbstractPriorModel): - self.__model = Object.from_object( - model - ) + self.__model = Object.from_object(model) - pickles: List[Pickle] = sa.orm.relationship( - "Pickle", - lazy="joined" - ) + pickles: List[Pickle] = sa.orm.relationship("Pickle", lazy="joined") def __getitem__(self, item: str): """ @@ -385,10 +255,7 @@ def __getitem__(self, item: str): for p in self.pickles: if p.name == item: return p.value - return getattr( - self, - item - ) + return getattr(self, item) def __contains__(self, item): for p in self.pickles: @@ -396,11 +263,7 @@ def __contains__(self, item): return True return False - def __setitem__( - self, - key: str, - value - ): + def __setitem__(self, key: str, value): """ Add a pickle. @@ -414,32 +277,15 @@ def __setitem__( value A string, bytes or object """ - new = Pickle( - name=key - ) - if isinstance( - value, - (str, bytes) - ): + new = Pickle(name=key) + if isinstance(value, (str, bytes)): new.string = value else: new.value = value - self.pickles = [ - p - for p - in self.pickles - if p.name != key - ] + [ - new - ] + self.pickles = [p for p in self.pickles if p.name != key] + [new] def __delitem__(self, key): - self.pickles = [ - p - for p - in self.pickles - if p.name != key - ] + self.pickles = [p for p in self.pickles if p.name != key] def value(self, name: str): try: @@ -447,38 +293,20 @@ def value(self, name: str): except AttributeError: return None - model_id = sa.Column( - sa.Integer, - sa.ForeignKey( - "object.id" - ) - ) + model_id = sa.Column(sa.Integer, sa.ForeignKey("object.id")) __model = sa.orm.relationship( - "Object", - uselist=False, - backref="fit_model", - foreign_keys=[model_id] + "Object", uselist=False, backref="fit_model", foreign_keys=[model_id] ) - instance_id = sa.Column( - sa.Integer, - sa.ForeignKey( - "object.id" - ) - ) + instance_id = sa.Column(sa.Integer, sa.ForeignKey("object.id")) __instance = sa.orm.relationship( - "Object", - uselist=False, - backref="fit_instance", - foreign_keys=[instance_id] + "Object", uselist=False, backref="fit_instance", foreign_keys=[instance_id] ) @classmethod def all(cls, session): - return session.query( - cls - ).all() + return session.query(cls).all() def __str__(self): return self.id diff --git a/autofit/non_linear/grid/sensitivity.py b/autofit/non_linear/grid/sensitivity.py index 99a86d6eb..9bf4bbddb 100644 --- a/autofit/non_linear/grid/sensitivity.py +++ b/autofit/non_linear/grid/sensitivity.py @@ -18,12 +18,7 @@ class JobResult(AbstractJobResult): - def __init__( - self, - number: int, - result: Result, - perturbed_result: Result - ): + def __init__(self, number: int, result: Result, perturbed_result: Result): """ The result of a single sensitivity comparison @@ -55,14 +50,14 @@ class Job(AbstractJob): use_instance = False def __init__( - self, - analysis_factory: "AnalysisFactory", - model: AbstractPriorModel, - perturbation_model: AbstractPriorModel, - base_instance: ModelInstance, - perturbation_instance: ModelInstance, - search: NonLinearSearch, - number: int, + self, + analysis_factory: "AnalysisFactory", + model: AbstractPriorModel, + perturbation_model: AbstractPriorModel, + base_instance: ModelInstance, + perturbation_instance: ModelInstance, + search: NonLinearSearch, + number: int, ): """ Job to run non-linear searches comparing how well a model and a model with a perturbation @@ -79,9 +74,7 @@ def __init__( search A non-linear search """ - super().__init__( - number=number - ) + super().__init__(number=number) self.analysis_factory = analysis_factory self.model = model @@ -90,15 +83,9 @@ def __init__( self.base_instance = base_instance self.perturbation_instance = perturbation_instance - self.search = search.copy_with_paths( - search.paths.for_sub_analysis( - "[base]", - ) - ) + self.search = search.copy_with_paths(search.paths.for_sub_analysis("[base]",)) self.perturbed_search = search.copy_with_paths( - search.paths.for_sub_analysis( - "[perturbed]", - ) + search.paths.for_sub_analysis("[perturbed]",) ) @cached_property @@ -126,27 +113,26 @@ def perform(self) -> JobResult: perturbed_result = self.perturbation_model_func(perturbed_model=perturbed_model) return JobResult( - number=self.number, - result=result, - perturbed_result=perturbed_result + number=self.number, result=result, perturbed_result=perturbed_result ) def base_model_func(self): - return self.search.fit( - model=self.model, - analysis=self.analysis - ) + return self.search.fit(model=self.model, analysis=self.analysis) def perturbation_model_func(self, perturbed_model): - return self.perturbed_search.fit( - model=perturbed_model, - analysis=self.analysis - ) + return self.perturbed_search.fit(model=perturbed_model, analysis=self.analysis) class SensitivityResult: - def __init__(self, results: List[JobResult]): + """ + The result of a sensitivity mapping + + Parameters + ---------- + results + The results of each sensitivity job + """ self.results = sorted(results) def __getitem__(self, item): @@ -158,21 +144,41 @@ def __iter__(self): def __len__(self): return len(self.results) + @property + def log_likelihoods_base(self) -> List[float]: + """ + The log likelihoods of the base model for each sensitivity fit + """ + return [result.log_likelihood_base for result in self.results] -class Sensitivity: + @property + def log_likelihoods_perturbed(self) -> List[float]: + """ + The log likelihoods of the perturbed model for each sensitivity fit + """ + return [result.log_likelihood_perturbed for result in self.results] + + @property + def log_likelihood_differences(self) -> List[float]: + """ + The log likelihood differences between the base and perturbed models + """ + return [result.log_likelihood_difference for result in self.results] + +class Sensitivity: def __init__( - self, - base_model: AbstractPriorModel, - perturbation_model: AbstractPriorModel, - simulation_instance, - simulate_function: Callable, - analysis_class: Type[Analysis], - search: NonLinearSearch, - job_cls: ClassVar = Job, - number_of_steps: Union[Tuple[int], int] = 4, - number_of_cores: int = 2, - limit_scale: int = 1, + self, + base_model: AbstractPriorModel, + perturbation_model: AbstractPriorModel, + simulation_instance, + simulate_function: Callable, + analysis_class: Type[Analysis], + search: NonLinearSearch, + job_cls: ClassVar = Job, + number_of_steps: Union[Tuple[int], int] = 4, + number_of_cores: int = 2, + limit_scale: int = 1, ): """ Perform sensitivity mapping to evaluate whether a perturbation @@ -212,9 +218,7 @@ def __init__( A scale of 0.5 means priors have limits smaller than the grid square with width half a grid square. """ - self.logger = logging.getLogger( - f"Sensitivity ({search.name})" - ) + self.logger = logging.getLogger(f"Sensitivity ({search.name})") self.logger.info("Creating") @@ -243,7 +247,9 @@ def step_size(self): The size of a step in any given dimension in hyper space. """ if isinstance(self.number_of_steps, tuple): - return tuple([1 / number_of_steps for number_of_steps in self.number_of_steps]) + return tuple( + [1 / number_of_steps for number_of_steps in self.number_of_steps] + ) return 1 / self.number_of_steps def run(self) -> SensitivityResult: @@ -258,14 +264,13 @@ def run(self) -> SensitivityResult: *self._headers, "log_likelihood_base", "log_likelihood_perturbed", - "log_likelihood_difference" + "log_likelihood_difference", ] physical_values = list(self._physical_values) results = list() for result in Process.run_jobs( - self._make_jobs(), - number_of_cores=self.number_of_cores + self._make_jobs(), number_of_cores=self.number_of_cores ): if isinstance(result, Exception): raise result @@ -273,17 +278,12 @@ def run(self) -> SensitivityResult: results.append(result) results = sorted(results) - os.makedirs( - self.search.paths.output_path, - exist_ok=True - ) + os.makedirs(self.search.paths.output_path, exist_ok=True) with open(self.results_path, "w+") as f: writer = csv.writer(f) writer.writerow(headers) for result_ in results: - values = physical_values[ - result_.number - ] + values = physical_values[result_.number] writer.writerow( padding(item) for item in [ @@ -292,15 +292,14 @@ def run(self) -> SensitivityResult: result_.log_likelihood_base, result_.log_likelihood_perturbed, result_.log_likelihood_difference, - ]) + ] + ) return SensitivityResult(results) @property def results_path(self): - return Path( - self.search.paths.output_path - ) / "results.csv" + return Path(self.search.paths.output_path) / "results.csv" @property def _lists(self) -> List[List[float]]: @@ -309,10 +308,7 @@ def _lists(self) -> List[List[float]]: the perturbation_model and create the individual perturbations. """ - return make_lists( - self.perturbation_model.prior_count, - step_size=self.step_size - ) + return make_lists(self.perturbation_model.prior_count, step_size=self.step_size) @property def _physical_values(self) -> List[List[float]]: @@ -321,14 +317,10 @@ def _physical_values(self) -> List[List[float]]: """ return [ [ - prior.value_for( - unit_value + prior.value_for(unit_value) + for prior, unit_value in zip( + self.perturbation_model.priors_ordered_by_id, unit_values ) - for prior, unit_value - in zip( - self.perturbation_model.priors_ordered_by_id, - unit_values - ) ] for unit_values in self._lists ] @@ -350,36 +342,23 @@ def _labels(self) -> Generator[str, None, None]: """ for list_ in self._lists: strings = list() - for value, prior_tuple in zip( - list_, - self.perturbation_model.prior_tuples - ): + for value, prior_tuple in zip(list_, self.perturbation_model.prior_tuples): path, prior = prior_tuple - value = prior.value_for( - value - ) - strings.append( - f"{path}_{value}" - ) + value = prior.value_for(value) + strings.append(f"{path}_{value}") yield "_".join(strings) @property - def _perturbation_instances(self) -> Generator[ - ModelInstance, None, None - ]: + def _perturbation_instances(self) -> Generator[ModelInstance, None, None]: """ A list of instances each of which defines a perturbation to be applied to the image. """ for list_ in self._lists: - yield self.perturbation_model.instance_from_unit_vector( - list_ - ) + yield self.perturbation_model.instance_from_unit_vector(list_) @property - def _perturbation_models(self) -> Generator[ - AbstractPriorModel, None, None - ]: + def _perturbation_models(self) -> Generator[AbstractPriorModel, None, None]: """ A list of models representing a perturbation at each grid square. @@ -395,29 +374,21 @@ def _perturbation_models(self) -> Generator[ prior.value_for(min(1.0, centre + half_step)), ) for centre, prior in zip( - list_, - self.perturbation_model.priors_ordered_by_id + list_, self.perturbation_model.priors_ordered_by_id ) ] yield self.perturbation_model.with_limits(limits) @property - def _searches(self) -> Generator[ - NonLinearSearch, None, None - ]: + def _searches(self) -> Generator[NonLinearSearch, None, None]: """ A list of non-linear searches, each of which is applied to one perturbation. """ for label in self._labels: - yield self._search_instance( - label - ) + yield self._search_instance(label) - def _search_instance( - self, - name_path: str - ) -> NonLinearSearch: + def _search_instance(self, name_path: str) -> NonLinearSearch: """ Create a search instance, distinguished by its name @@ -432,9 +403,7 @@ def _search_instance( """ paths = self.search.paths search_instance = self.search.copy_with_paths( - paths.for_sub_analysis( - name_path, - ) + paths.for_sub_analysis(name_path,) ) return search_instance @@ -446,15 +415,9 @@ def _make_jobs(self) -> Generator[Job, None, None]: Each job fits a perturbed image with the original model and a model which includes a perturbation. """ - for number, ( - perturbation_instance, - perturbation_model, - search - ) in enumerate(zip( - self._perturbation_instances, - self._perturbation_models, - self._searches - )): + for number, (perturbation_instance, perturbation_model, search) in enumerate( + zip(self._perturbation_instances, self._perturbation_models, self._searches) + ): instance = copy(self.instance) instance.perturbation = perturbation_instance @@ -469,16 +432,13 @@ def _make_jobs(self) -> Generator[Job, None, None]: base_instance=self.instance, perturbation_instance=perturbation_instance, search=search, - number=number + number=number, ) class AnalysisFactory: def __init__( - self, - instance, - simulate_function, - analysis_class, + self, instance, simulate_function, analysis_class, ): """ Callable to delay simulation such that it is performed @@ -489,9 +449,5 @@ def __init__( self.analysis_class = analysis_class def __call__(self): - dataset = self.simulate_function( - self.instance - ) - return self.analysis_class( - dataset - ) + dataset = self.simulate_function(self.instance) + return self.analysis_class(dataset) diff --git a/test_autofit/non_linear/grid/test_sensitivity/test_results.py b/test_autofit/non_linear/grid/test_sensitivity/test_results.py new file mode 100644 index 000000000..7a9a79bb4 --- /dev/null +++ b/test_autofit/non_linear/grid/test_sensitivity/test_results.py @@ -0,0 +1,30 @@ +from autofit.non_linear.grid.sensitivity import SensitivityResult, JobResult +import pytest + + +class Result: + def __init__(self, log_likelihood): + self.log_likelihood = log_likelihood + + +@pytest.fixture(name="job_result") +def make_result(): + return JobResult( + number=0, + result=Result(log_likelihood=1.0), + perturbed_result=Result(log_likelihood=2.0), + ) + + +def test_job_result(job_result): + assert job_result.log_likelihood_base == 1.0 + assert job_result.log_likelihood_perturbed == 2.0 + assert job_result.log_likelihood_difference == 1.0 + + +def test_result(job_result): + result = SensitivityResult(results=[job_result]) + + assert result.log_likelihoods_base == [1.0] + assert result.log_likelihoods_perturbed == [2.0] + assert result.log_likelihood_differences == [1.0]