Hot Spot Population Analysis (HSPA)

This projects contains scripts to run a hotpot population analysis. By this statistical method it is tested if a population of small p-values exists in a set of p-value results which are non significant on their own. As a use case we implement the test for p-values testing neutrino fluxes from different directions in the sky.

Import dependencies

This project dependes on the python packages

• build-in: argparse, collections, cPickle, glob, os, time
• numpy
• scipy
• healpy
• matplotlib (for plotting)
• cosmolopy (for plotting, Kowalski plot)
• FIRESONG (for plotting, Kowalski plot)

Analysis Chain

This section explains how to run the scripts to get a HPA result. As input a set of all sky scans, as TS parametrization and single PS signal trials (in the optimal case also some with very high TS values) are needed.

To run the analysis, the following steps have to be executed:

1. get_warm_spots_from_skylab_all_sky_scans.py

   This script extracts a list of local warm spots from a skylab all sky scan. It computs a list of (theta, phi, p-value) tuples from a skylab p-value map.

   Parameters:

   • Infiles: List of skylab all_sky_scan file paths
   • Outfile: File path of the output file, where the list with local warm spots should be written.
   • log10pVal_threshold: Threshold in -log10(p-value) above which local warm spots are considered.
   • min_ang_dist: Minimal angular distance allowed between two local warm spots.

2. generate_expectation.py

   Will check for different thresholds if spots are poissonian distributed. Will produce a parametrization / spline of the #spots expectation vs p-value threshold. The parametrization and the spline of this parametization are saved in pickle files.

   Parameters:

   • Infiles: List of inputfile pathes. Inputfile pathes should contain lists of LocalWarmSpotLists.
   • Outputfile: Location of the file where the LocalWarmSpotExpecation is saved.

3. calculate_max_local_pValues.py

   Will generate background trials and calculate the HPA TS value for these background trials. A list of TS values will be the output.

   Parameters:

   • infiles: List of inputfile pathes. Inputfile pathes should contain lists of LocalWarmSpotLists
   • outfile: File location where the list of HPA_analysis_result are saved.
   • expectation: Path of the LocalWarmSpotExpecation file.

4. prepare_bgd_pool.py Generates a background pool object and saves it. Needed for signal trial generation.

   Parameters:

   • infiles: List of input files. Input files should contain a list of extracted background populations.
   • outfile: Path of the output file.

5. prepare_signal_pool.py Generates a single_spot_trial_pool and save it. Needed for signal trial generation.

   Parameters:

   • indir: Folder containing single spot trial files.
   • outfile: Path of the output file.
   • parametrization: Path of the TS parametrization

6. hpa_signal_trials.py Generates signal simulation trials for the hot spot population analysis for a given source count distribution. The source count distribution is either an equal flux on Earth source count distribution which is parametrized by the number of sources in the population and the neutrino flux at Earth from a single source, or by a file describing the source count distribution. In the output file the HPA TS value, the threshold p-value, the expected and observed numbers of events at the threshold p-value and the total number of injected events.

   Parameters:

   • infile_signal: Path of the single spot signal pool file. (Output of Step 5)
   • infile_background: Path of the background warm spot pool file. (Output of Step 4)
   • expectation: Path of the LocalWarmSpotExpecation file. (Output of Step 2)
   • outfile: Path where the signal simulation results should be saved.
   • n_iter: Number of simulation trials that are generated.
   • seed: Random Number Seed.
   • nsrc: Number of Source in population.
   • phi_inj: Flux @ 1 GeV of sources at Earth. Units: 1/GeV cm^2 s.
   • infile_source_count: Path of file containing a source count distribution generated by FIRESONG.
   • density: Density of sources. Scales the number of sources for source count distributions generated with FIRESONG.

7. calculate_sensitivity.py Fits the background HPA TS and extrapolates it. The result is used to calculate the sensitvity for a source model by finding the signal strength that matches the sensitivity, and discovery potential. If the unblinded TS value is given also the upper limit is calculated. Results are stored as sensitivity_results.

   Parameters:

   • background_HPA_trials: Path of file that contains background HPA trials.
   • outfile: Path where the output file should be stored.
   • signal_files: List of HPA signal trials files. Note that all files should represent one source model, except of changing flux strength.
   • unblinded_value: If you have unblinded give the TS value here and upper limits are calculated.
   • All remaining arguments are passed to describe the signal model.

8. unblind_result.py Extracts the local warm spots from a skymap and calculates the HPA-TS value and compares it to background trials and gives a post trial p-value.

   Parameters:

   • infile: Path of the skylab all sky scan that contains the unblinded sky map.
   • background_trials: Path of the background HPA TS trials file.
   • expectation: Path of the LocalWarmSpotExpecation file.
   • outfile: Path where the unblinding results should be saved.
   • log10pVal_threshold: Give the -log10(p-value) above that spots should not be considerd. Default: 2.0.
   • min_ang_dist: Give the minimal angular distance allowed between two local warm spots. Units: degrees. Default: 1.

9. (optional) several scripts for plotting can be found in plotting/

Further scripts:

• data_types.py # classes with data types

• SingleSpotTS2pValueParametrization.py # class for TS -> pValue conversion

• skylab_data.py # classes interacting with skylab

• source_count_dist.py # classes representing source_count_distributions

• statistics.py # functions, related to statistics

• utils.py # functions, classes

• external_data/ # Digitized Data from references

• plotting/ # plot the stuff

• test_data/ # Some data for testing the scripts