{ "cells": [ { "cell_type": "markdown", "metadata": {}, "source": [ "Modeling: Quantity\n", "==================\n", "\n", "Before reading through this example, make sure you have read through and completed the example script or notebook\n", "`autolens_workspace/*/misc/model_quantity/deflections.py`\n", "\n", "Mapping the deflection angle maps of two mass profiles is not necessarily sufficient to determine how they should map\n", "from one to another in a real strong lens system. This is because in a real lens, constraints are only provided where\n", "we observed source light. Thus, we must factor in the source light in order to produce such a realistic mapping.\n", "\n", "This script again fits the deflection angles of an `PowerLaw` profile to an `NFW` mass model. However, it first\n", "simulates the `EllNEW` mass model as an image of a strong lens's. The signal to noise of the lensed source then\n", "determines the signal to noise map of the defleciton angle map that is fitted.\n", "\n", "This script fits a `DatasetQuantity` dataset of a 'galaxy-scale' strong lens with a model. The `DatasetQuantity` is the\n", "deflection angle map of a `NFW` mass model which is fitted by a lens model where:\n", "\n", " - The lens galaxy's light is omitted (and is not present in the simulated data).\n", " - The lens galaxy's total mass distribution is an `PowerLaw`.\n", " - There is no source galaxy." ] }, { "cell_type": "code", "metadata": {}, "source": [ "%matplotlib inline\n", "from pyprojroot import here\n", "workspace_path = str(here())\n", "%cd $workspace_path\n", "print(f\"Working Directory has been set to `{workspace_path}`\")\n", "\n", "from os import path\n", "import autofit as af\n", "import autolens as al\n", "import autolens.plot as aplt" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "__Grid__\n", "\n", "Define the 2D grid the quantity (in this example, the deflections) is evaluated using." ] }, { "cell_type": "code", "metadata": {}, "source": [ "grid = al.Grid2D.uniform(shape_native=(100, 100), pixel_scales=0.1)" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "__Tracer__\n", "\n", "Create a tracer which we will use to create our `DatasetQuantity`. \n", "\n", "In this example, the tracer includes a source galaxy whose light will be used to determine the S/N map of the \n", "deflection angle map that is fitted." ] }, { "cell_type": "code", "metadata": {}, "source": [ "lens_galaxy = al.Galaxy(\n", " redshift=0.5,\n", " mass=al.mp.NFW(\n", " centre=(0.0, 0.0),\n", " ell_comps=al.convert.ell_comps_from(axis_ratio=0.9, angle=45.0),\n", " kappa_s=0.2,\n", " scale_radius=20.0,\n", " ),\n", ")\n", "\n", "source_galaxy = al.Galaxy(\n", " redshift=1.0,\n", " bulge=al.lp.Sersic(\n", " centre=(0.0, 0.0),\n", " ell_comps=al.convert.ell_comps_from(axis_ratio=0.8, angle=60.0),\n", " intensity=4.0,\n", " effective_radius=0.1,\n", " sersic_index=1.0,\n", " ),\n", ")\n", "\n", "tracer = al.Tracer.from_galaxies(galaxies=[lens_galaxy, source_galaxy])" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "__Simulator__\n", "\n", "We can create a `SimulatorImaging` object which creates an image of the tracer above. \n", "\n", "Because the simulation includes adding noise, the output `Imaging` dataset has a signal-to-noise map that we next \n", "use to set the S/N of the defleciton angles that we fit. " ] }, { "cell_type": "code", "metadata": {}, "source": [ "simulator = al.SimulatorImaging(\n", " exposure_time=300.0, background_sky_level=0.1, add_poisson_noise=False\n", ")\n", "\n", "dataset = simulator.via_tracer_from(tracer=tracer, grid=grid)" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "__Dataset__\n", "\n", "Use this `Tracer`'s 2D deflection angles to create the `DatasetQuantity`.\n", "\n", "The noise map is determiend via the signal-to-noise map of the simulated `Imaging` above." ] }, { "cell_type": "code", "metadata": {}, "source": [ "deflections_2d = tracer.deflections_yx_2d_from(grid=grid)\n", "\n", "dataset = al.DatasetQuantity.via_signal_to_noise_map(\n", " data=deflections_2d, signal_to_noise_map=dataset.signal_to_noise_map\n", ")" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "__Mask__\n", "\n", "The model-fit requires a `Mask2D` defining the regions of the deflections we fit, which we define and apply to the \n", "`DatasetQuantity` object." ] }, { "cell_type": "code", "metadata": {}, "source": [ "mask = al.Mask2D.circular(\n", " shape_native=dataset.shape_native, pixel_scales=dataset.pixel_scales, radius=3.0\n", ")\n", "\n", "dataset = dataset.apply_mask(mask=mask)" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "__Model__\n", "\n", "We compose a lens model where:\n", "\n", " - The lens galaxy's total mass distribution is an `PowerLaw` [6 parameters].\n", "\n", "The number of free parameters and therefore the dimensionality of non-linear parameter space is N=14." ] }, { "cell_type": "code", "metadata": {}, "source": [ "lens = af.Model(al.Galaxy, redshift=0.5, mass=al.mp.PowerLaw)\n", "\n", "model = af.Collection(galaxies=af.Collection(lens=lens))" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "__Search__\n", "\n", "The lens model is fitted to the data using the nested sampling algorithm Nautilus (see `start.here.py` for a \n", "full description).\n", "\n", "The folders: \n", "\n", " - `autolens_workspace/*/imaging/modeling/searches`.\n", " - `autolens_workspace/*/imaging/modeling/customize`\n", " \n", "Give overviews of the non-linear searches **PyAutoLens** supports and more details on how to customize the\n", "model-fit, including the priors on the model.\n", "\n", "The `name` and `path_prefix` below specify the path where results ae stored in the output folder: \n", "\n", " `/autolens_workspace/output/imaging/modeling/simple__no_lens_light/mass[sie]_source[bulge]/unique_identifier`.\n", "\n", "__Unique Identifier__\n", "\n", "In the path above, the `unique_identifier` appears as a collection of characters, where this identifier is generated \n", "based on the model, search and dataset that are used in the fit.\n", "\n", "An identical combination of model and search generates the same identifier, meaning that rerunning the script will use \n", "the existing results to resume the model-fit. In contrast, if you change the model or search, a new unique identifier \n", "will be generated, ensuring that the model-fit results are output into a separate folder.\n", "\n", "We additionally want the unique identifier to be specific to the dataset fitted, so that if we fit different datasets\n", "with the same model and search results are output to a different folder. We achieve this below by passing \n", "the `dataset_name` to the search's `unique_tag`.\n", "\n", "__Number Of Cores__\n", "\n", "We include an input `number_of_cores`, which when above 1 means that Nautilus uses parallel processing to sample multiple \n", "lens models at once on your CPU. When `number_of_cores=2` the search will run roughly two times as\n", "fast, for `number_of_cores=3` three times as fast, and so on. The downside is more cores on your CPU will be in-use\n", "which may hurt the general performance of your computer.\n", "\n", "You should experiment to figure out the highest value which does not give a noticeable loss in performance of your \n", "computer. If you know that your processor is a quad-core processor you should be able to use `number_of_cores=4`. \n", "\n", "Above `number_of_cores=4` the speed-up from parallelization diminishes greatly. We therefore recommend you do not\n", "use a value above this.\n", "\n", "For users on a Windows Operating system, using `number_of_cores>1` may lead to an error, in which case it should be \n", "reduced back to 1 to fix it." ] }, { "cell_type": "code", "metadata": {}, "source": [ "search = af.Nautilus(\n", " path_prefix=path.join(\"misc\", \"modeling\"),\n", " name=\"quantity_via_deflections_fit_source_snr\",\n", " n_live=100,\n", " number_of_cores=1,\n", ")" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "__Analysis__\n", "\n", "The `AnalysisQuantity` object defines the `log_likelihood_function` used by the non-linear search to fit the model to \n", "the `DatasetQuantity` dataset.\n", "\n", "This includes a `func_str` input which defines what quantity is fitted. It corresponds to the function of the \n", "model `Tracer` objects that are called to create the model quantity. For example, if `func_str=\"deflections_yx_2d_from\"`, \n", "the deflections is computed from each model `Tracer`." ] }, { "cell_type": "code", "metadata": {}, "source": [ "analysis = al.AnalysisQuantity(dataset=dataset, func_str=\"deflections_yx_2d_from\")" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "__Model-Fit__\n", "\n", "We begin the model-fit by passing the model and analysis object to the non-linear search (checkout the output folder\n", "for on-the-fly visualization and results)." ] }, { "cell_type": "code", "metadata": {}, "source": [ "result = search.fit(model=model, analysis=analysis)" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "__Result__\n", "\n", "The search returns a result object, which whose `info` attribute shows the result in a readable format:" ] }, { "cell_type": "code", "metadata": {}, "source": [ "print(result.info)" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "We plot the maximum likelihood fit, tracer images and posteriors inferred via Nautilus.\n", "\n", "Checkout `autolens_workspace/*/imaging/results` for a full description of analysing results in **PyAutoLens**." ] }, { "cell_type": "code", "metadata": {}, "source": [ "print(result.max_log_likelihood_instance)\n", "\n", "tracer_plotter = aplt.TracerPlotter(\n", " tracer=result.max_log_likelihood_tracer, grid=result.grid\n", ")\n", "tracer_plotter.subplot_tracer()\n", "\n", "fit_quantity_plotter = aplt.FitQuantityPlotter(fit=result.max_log_likelihood_fit)\n", "fit_quantity_plotter.subplot_fit()\n", "\n", "search_plotter = aplt.NautilusPlotter(samples=result.samples)\n", "search_plotter.cornerplot()" ], "outputs": [], "execution_count": null }, { "cell_type": "markdown", "metadata": {}, "source": [ "Checkout `autolens_workspace/*/imaging/modeling/results.py` for a full description of the result object." ] }, { "cell_type": "code", "metadata": {}, "source": [], "outputs": [], "execution_count": null } ], "metadata": { "anaconda-cloud": {}, "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.6.1" } }, "nbformat": 4, "nbformat_minor": 4 }