<jupyter_start><jupyter_text>Multi-Gaussian exampleAs an introductory example we can use astroABC to find the posterior distribution for some Gaussian distributed data. Although in this case we already know the likelihood this example is to illustrate how to call astroABC and provide user-defined metrics.<jupyter_code># start by importing astroabc and numpy
import numpy as np
import astroabc<jupyter_output><empty_output><jupyter_text>We need to provide:- a dataset- a forwards simulating model for the data- a method defining the distance metric.For this example we generate a dataset of a 1000 draws from a 2D multi-Gaussian where the true means are e.g$\mu_{1,true} = 0.037579, \, \mu_{2, true}=0.573537$and we have fixed the covariance matrix to be  a diagonal matrix with $\sigma_1^2 = \sigma_2^2 = 0.05$. We can do this using an inbuilt model class in astroabc.<jupyter_code>#make the fake data with diagonal covariance
means= np.array([0.037579, 0.573537])
cov =np.array([0.01,0.005,0.005,0.1])
data = astroabc.Model("normal",1000).make_mock(means,cov)<jupyter_output><empty_output><jupyter_text>In this example the make_mock method also provides a forwards simulating model for the data.<jupyter_code>#define a method for simulating the data given input parameters
def simulation(param, pool=None):
    cov =np.array([0.01,0.005,0.005,0.1])
    #Ideally do something with the pool here
    return astroabc.Model("normal",10000).make_mock(param,cov)

model_sim = simulation<jupyter_output><empty_output><jupyter_text>Next we define a distance metric method. In this example instead of using all of the data (all 1000 draws from a 2D Gaussian) we use the means of the data as a summary statistic and our distance metric is the sum over the difference in the means for the 2D Gaussian<jupyter_code>def dist_metric(d,x):
    return np.sum(np.abs(np.mean(x,axis=0) - np.mean(d,axis=0)))<jupyter_output><empty_output><jupyter_text>Next we specify priors for each of the parameters we want to vary in the sampler. This is done by specifying a list of tuples. The zeroth element in each tuple should be a string specifying the prior for this parameter and the first element should be a list of the hyperparameters needed for this prior.e.g. below we use a normal distribution with mean  0 and variance 0.5 for the first parameter and a uniform distribution from 0.1 - 0.9 for the second parameter.<jupyter_code>priors =  [('normal', [0.03,0.1]), ('uniform', [0.1, 0.9])]<jupyter_output><empty_output><jupyter_text>Next we need to set some keywords for astroABC. This can be done by creating a dictionary of inputs which are passed to the sampler. Many of these entries have defaults and do not need to be specified explicitly.Only the name of the distance metric method needs to be explicity provided as a keyword.The full set of keywords are given in the doc string of the class. Some examples are- tol_type: which specifies the decreasing tolerance levels. "exp","lin", "log" and "const" are options. (default = 'exp')- verbose: level of verbosity, 0 = no printout to screen, 1 = print to screen  (default = 0)- adapt_t: Boolean True/False for adaptive threshold setting (default = False)- threshold: qth quantile used in adaptive threshold setting (default = 75)- pert_kernel: 1 =component wise pert. with local diag variance; 2 = multivariate pert. based on local covariance- variance_method: 0 =weighted covariance, 1= Filippi, 2 = TVZ, 3= Leodoit_Wolf, 4=k-nn (default = 0)- dfunc:method for calculating the distance metric- from_restart: Boolean True/False- restart: string name of restart file- outfile:string specifying name of output file (default = abc_out.txt)- mpi: Boolean True/False (default = False)- mp:Boolean True/False (default = False)- num_proc:number of threads for mp setting (default = None)Please see the doc strings of the astroABC sampler for details on each of these settings.<jupyter_code>prop={'dfunc':dist_metric, 'outfile':"gaussian_example.txt", 'verbose':1, 'adapt_t': True, 'pert_kernel':2}<jupyter_output><empty_output><jupyter_text>Now we are ready to create an instance of our sampler. To do this we just need to specify the following to the astroabc.ABC_classastroabc.ABC_class(number of parameters,number of particles,data,tolerance levels,number of iterations,priors,prop)To begin let's run in serial using 100 particles for 30 iterations with default tolerance levels of a maximum threshold of 0.7 and  a minimum threshold of 0.05:<jupyter_code>sampler = astroabc.ABC_class(2,100,data,[0.5,0.002],20,priors,**prop)<jupyter_output>########################     astroABC     ########################	
	 	
	 Npart=100 	 numt=20 	 tol=[0.5000,0.0020] exp
	 Priors= [('normal', [0.03, 0.1]), ('uniform', [0.1, 0.9])]<jupyter_text>Then we simply begin sampling on our data...<jupyter_code>sampler.sample(model_sim)<jupyter_output>Step: 0 	 tol: 0.5 	 Params: [0.004831022751718228, 0.5316987441872384]
	 Step: 1 	 tol: 0.3674356729086831 	 Params: [0.02321074264082056, 0.5742511100054313]
	 Step: 2 	 tol: 0.30166515686162976 	 Params: [0.036441292743712723, 0.5849798902653633]
	 Step: 3 	 tol: 0.2339575792400317 	 Params: [0.04721763109698714, 0.5835376098791599]
	 Step: 4 	 tol: 0.19340880389783577 	 Params: [0.043681453658168924, 0.5825229124687235]
	 Step: 5 	 tol: 0.1698794471891293 	 Params: [0.025025449631187073, 0.5864426900677008]
	 Step: 6 	 tol: 0.1447905652346464 	 Params: [0.03934585197338663, 0.5838000248009465]
	 Step: 7 	 tol: 0.1243271216418435 	 Params: [0.04108654251443204, 0.5850115948283879]
	 Step: 8 	 tol: 0.10398803512341327 	 Params: [0.03269305321900622, 0.5840846504880716]
	 Step: 9 	 tol: 0.08868583214998207 	 Params: [0.03250952497115013, 0.5855360667857193]
	 Step: 10 	 tol: 0.0756948315479395 	 Params: [0.035796696950396154, 0.580821486914087]
	 Step: 11 	 tol: 0.06511416184901217 [...]<jupyter_text>The output above shows the estimated means of the 2D Gaussian averaged over all 100 particles at each iteration, together with the tolerance level. Note above that the sampling may end before 20 iterations if the minimum tolerance level is reached first.Recall that the true parameter values are $\mu_{1,true} = 0.037579, \, \mu_{2, true}=0.573537$  K-Nearest Neighbours estimation for data sample with covariance matrixWe could also have created a dataset with a full covariance matrix using<jupyter_code>means= np.array([0.7579, 0.273537])
cov = np.array([0.1,0.01,0.01,0.1])
data_cov = astroabc.Model("normal",1000).make_mock(means,cov)<jupyter_output><empty_output><jupyter_text>Keeping model simulation and distance methods the same as above. We can select a different way of estimating the covariance amongst all the particles using k-nearest neighbours. This returns a local covariance estimate and in many cases this reaches convergence faster then using a weighted covariance amongst all particles.<jupyter_code>priors =  [('uniform', [0.1,0.9]), ('uniform', [0.1, 0.9])]
prop={'dfunc':dist_metric, 'outfile':"gaussian_example.txt", 'verbose':1, \
      'adapt_t': True, 'variance_method':4, 'k_near':10 }
sampler = astroabc.ABC_class(2,100,data_cov,[0.5,0.002],25,priors,**prop)
sampler.sample(model_sim)<jupyter_output>Step: 0 	 tol: 0.5 	 Params: [0.6417417888667538, 0.3769805241152575]
	 Step: 1 	 tol: 0.43598082896982815 	 Params: [0.6728761879909225, 0.3626886691562053]
	 Step: 2 	 tol: 0.358520432619466 	 Params: [0.7227864485788014, 0.31132619316770155]
	 Step: 3 	 tol: 0.2833997638790343 	 Params: [0.7440094967360175, 0.2987062443319014]
	 Step: 4 	 tol: 0.21592705669630358 	 Params: [0.7578314852638425, 0.28997038404805464]
	 Step: 5 	 tol: 0.1832059820435786 	 Params: [0.7725999199112689, 0.2711524459774648]
	 Step: 6 	 tol: 0.1563979221312114 	 Params: [0.7549734561659688, 0.280156875829783]
	 Step: 7 	 tol: 0.12951796909825125 	 Params: [0.7561589203189962, 0.28703131841749463]
	 Step: 8 	 tol: 0.1093959516000065 	 Params: [0.7646331246580728, 0.2767507341644505]
	 Step: 9 	 tol: 0.08927254752993467 	 Params: [0.7631681442434737, 0.2820566810999924]
	 Step: 10 	 tol: 0.07633646127328605 	 Params: [0.7584185873502948, 0.2795299299446455]
	 Step: 11 	 tol: 0.061068540755597894 	 Params: [0[...]
{'fraction_non_alphanumeric': 0.10074218274729285, 'fraction_numerical': 0.17240540211704586, 'mean_word_length': 5.157303370786517, 'pattern_counts': {'":': 0, '<': 38, '<?xml version=': 0, '>': 38, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '10', 'dir': 'jupyter-structured-clean-dedup', 'id': '114026', 'license': 'MIT', 'n_tokens_mistral': 3732, 'n_tokens_neox': 2752, 'n_words': 902, 'path': 'examples/multi-gaussian example.ipynb', 'repo_name': 'dhuppenkothen/astroABC', 'size': '8219'}