Bayesian Inference Implementation in LALSuite for CBC GW Sources
Explore the implementation of Markov Chain Monte Carlo (MCMC) for parameter estimation in LALSuite, focusing on compact binary coalescence (CBC) gravitational wave sources. Learn about the motivation, Bayesian inference concepts, physical parameters, and practical applications discussed at the KAGRA meeting in 2015.
Download Presentation
Please find below an Image/Link to download the presentation.
The content on the website is provided AS IS for your information and personal use only. It may not be sold, licensed, or shared on other websites without obtaining consent from the author. If you encounter any issues during the download, it is possible that the publisher has removed the file from their server.
You are allowed to download the files provided on this website for personal or commercial use, subject to the condition that they are used lawfully. All files are the property of their respective owners.
The content on the website is provided AS IS for your information and personal use only. It may not be sold, licensed, or shared on other websites without obtaining consent from the author.
E N D
Presentation Transcript
Brief Explanation of MCMC implementation in lalsuite Hyung Won Lee, Inje University with Chunglee Kim(KHU&KISTI) and Jeongcho Kim(Inje) 6 Feb. 2015, KAGRA f2f meeting, Hongo Campus 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 1
Contents Motivation CBC gravitational wave sources Parameter Estimation with Bayesian Inference Implementation in lalsuite Suggestions for KAGALI Discussions 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 2
Motivation How does it work parameter estimation Understanding Bayesian inference/MCMC Which physical parameters How it was implemented in lalsuite How it could be implemented for KAGALI Get practical work out items 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 3
CBC GW sources Neutron star(NS) + Black hole(BH) binary Mass range : 1? ~103? Frequency range : 1?? ~ 10??? Waveform : Inspiral, Merger, Ringdown 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 4
CBC Inspiral GW Waveform 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 5
Time Domain Template TaylorT1, 2, 3, 4 PadeT1 : P-approximation EOB : effective one body SpinTaylorT1, 2, 3, 4, 5 IMRPhenomA, B, C 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 6
Frequency Domain Template TaylorF1 TaylorF2 : standard for detection pipeline, SPA approximation TaylorF2Amp : include higher harmonics due to amplitude corrections up to 2.5 pN PadeF1 BCV BCVSpin 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 7
TaylorF2Amp Waveform Arun et.al., PRD79, 104023(2009) Amplitude corrections up to n=5 (2.5pN) 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 8
SPA phase factor + Phase corrections up to 3.5PN (standard in LAL) 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 9
SPA phase factor(Code) 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 10
Physical Parameters Individual masses Luminosity distance Inclination angle Coalescence phase Coalescence time at geocenter Declination Right ascension Polarization Spins, quadrupole deformations 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 11
Parameter Estimation with Bayesian inference Bayesian Inference ? ? = ?,? (?|?) ? (?|?) ? ?? = ? ?,? = ? ? ? = ? ? ? ? : unobservable model parameters ? : observable data ?,? : Joint probability observing data ? with model parameter ? ? ? ?????? ??? ???????? Metropolis-Hasting Algorithm, 106~107 samples 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 12
Parameter Estimation with Bayesian inference MCMC Samples PTMCMCOutput.00 PTMCMCOutput.nn Post Process cbcByesPostProc.py 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 13
Overlap Calculation Template Data + (?) ? (?) ??(?) ? = 4 0 ?? One-sided power spectral density of noise + (?)2 ??(?)?? ???2= 4 0 ? 2 Network SNR ???2= ?=1 ???? 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 14
Likelihood Calculation (?|?) ? ? (?,?)2 exp 2 0 ?? ??(?) 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 15
Various Overlaps Unnormalized (?) + 1(?) 2 ? = 4 0 ?? ??(?) Normalized (?) + ?(?) ? ??= 4 0 ??= ?? ??(?) ? ?1?2 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 16
MCMC Process Accept next state with probability (?|?2)1/? ?2 (?|?1)1/? ?1 Metropolis-Hasting smapling ? = 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 17
Parallel Tempering Use few chains with different temperature 1 ?, ? > 1 Use likelihood ? ? Improve convergence Improve mixing The higher temperature, the smoother distribution ? > 1 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 18
Implementation in lalsuite LAL MCMC Pipeline Directory structure of lalsuite Basic entities of LALSuite LAL applications Likelihood, proposal and prior functions Simulation related functions and waveforms Post Processing python script 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 19
LAL MCMC Pipeline Start Initialization MPI synchronize PT MCMC Clearance End 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 20
LALInferenceRunState ProcessParamsTable /** A ProcessParamsTable with command line arguments */ LALInferenceInitModelFunction initModel; /** A function that returns a new set of variables for the model */ LALInferenceAlgorithm LALInferenceEvolveOneStepFunction evolve; /** The algorithm's single iteration function */ LALInferencePriorFunction LALInferenceCubeToPriorFunction /** MultiNest prior for the parameters */ LALInferenceLikelihoodFunction /** The likelihood function */ LALInferenceProposalFunction LALInferenceLogFunction LALInferenceTemplateFunction templt; /** The template generation function */ more *commandLine; algorithm; /** The algorithm function */ prior; /** The prior for the parameters */ CubeToPrior; likelihood; proposal; /** The proposal function */ logsample; /** Log sample, i.e. to disk */ 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 21
LAL MCMC pipeline in a glance MCMC main main (start) parseCommandLine Options Initialize Rusnstate initializeMCMC addVariables InitializeCBC Initialize likelihood Initialize MCMC state MPI_Barrier wait synchronization runState ->alogorithm main MCMC routine PTMCMCAlgorithm Deinitialize Clean resources end 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 22
PTMCMCAlgorithm (start) PTMCMCAlgorithm MPI setup, runComplete = 0 ? Ladder setup for parallel tempering Ladder[i] = Min(??? Set temperature for this chain if(MPI_rank == 0) output Anealing related setup Debug related setup ???) ?? ??? 1 Receive runComplete state from root process If(MPI_rank 0) non-blocking receive MPI_Barrier No !runComplete Yes Ladder update if need Calculate ACL(autocorrelation length) If( i>Neff) runComplete=1 Anealing post-anealing runState ->evolve OneStep MCMC routine PTMCMCOneStep Save for every Nskip step If(i>Niter) runComplete=1 If(MPI_rank==0 && runComplete==1) send run complete value to all other parallel process MPI_Barrier XLALFree() resource free end 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 23
PTMCMCOneStep (start) PTMCMCOneStep MPI setup runState ->proposal LALInferenceCyclicProposal Compute prior & likelihood determine acceptance probability Accept end LALInferenceCyclicProposal LALInferenceCyclicProposal (start) LALInferenceProposalFunction *cycle = NULL Must have cycle array and cycle array length in propArgs (cycle[i])(runState, proposedParams) Call proposal end 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 24
Mcmc main Call TaylorF2Amp in MCMC PTMCMCAlgorithm() PTMCMCOneStep() LALInferenceUndecomposed FreqDomainLogLikelihood() LALInferenceTemplate XLALSimInspiralChooseWaveform() XLALSimInspiralChooseFD WaveformFromCache() XLALSimInspiralChooseFDWaveform() XLALSimInspiralTaylorF2Amp() 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 25
LAL parameter estimation(PE) pipeline 1. waveform library [LALSimulation] 2. PE library [LALInference & LALapps] 3. basic libraries are scattered, but most basic ones are in [LAL] 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 26
waveform library [LALSimulation] LALSimulation/src/ .. LALSimInspiral.h -> protocol definitions XLALSimInspiralF2AmpPlus() -> + XLALSimInspiralF2AmpCross() -> LALSimInspiral.c -> XLALSimInspiralChooseFDWaveform() call a waveform model XLALSimInspiralF2AmpPlus() XLALSimInspiralF2AmpCross() LALSimInspiralTaylorF2Amp.c -> example of a FD waveform XLALSimInspiralF2AmpPlus() XLALSimInspiralF2AmpCross() 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 27
PE library [LALInference] LALInference/src/.. LALInference.h LALInferenceReadData.c LALInferecneTemplate.c LALInference.h -> collections of structures typedef tagLALInferenceRunState signal/parameters { -> Structure to contain inference run state } LALInferenceRunState; typedef struct tagLALInferenceIFOData data for detector including sensitivity { -> Structure to contain IFOdata } LALInferenceIFOData; 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 28
PE library [LALInference] LALInference/src/.. LALInference.h LALInferenceReadData.c LALInferecneTemplate.c reads injection parameters from an xml file reads PSD for each detector reads data for each detector generate noise for each detector LALInferenceIFOData *LALInferenceReadData(ProcessParamsTable *commandLine) -> Read in the data and store it in a LALInferenceIFOData structure 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 29
PE library [LALInference] LALInference/src/.. LALInference.h LALInferenceReadData.c LALInferecneTemplate.c - Jump around parameter space ??generate template signal LALInferenceTemplate.c -> template calls to LAL template functions LALInferenceTemplateXLALSimInspiralChooseWaveform() Wrapper for LALSimulation waveforms : XLALSimInspiralChooseFDWaveform() XLALSimInspiralChooseTDWaveform() - Parameters of non-spinning waveform in LALInferencePrior.c 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 30
PE library [LALInference] LALInference/src/.. LALInference.h LALInferenceLikelihood.c LALInferecneTemplate.c LALInferenceLikelihood.c -> Bayesian Followup likelihood function LALInferenceUndeomposedFreqDomainLogLikelihood() Required (`currentParams') parameters are: "rightascension" (REAL8, radian, 0 <= RA <= 2pi) "declination" "polarisation (REAL8, radian, 0 <= psi <= pi) "time" (REAL8, radian, -pi/2 <= dec <=pi/2) (REAL8, GPS sec.) 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 31
PE library [LALapps] LALapps/src/inspiral/posterior/.. LALInferenceMCMC.c LALInferenceMCMCSampler.c LALInferenceMCMC.c -> Bayesian Followup function testing site main() 1. Set up structures for MCMC 2. Choose the likelihood function -> ex) student-T , Undecomposed. 3. Call MCMC algorithm -> sampling posterior function 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 32
PE library [LALapps] LALapps/src/inspiral/posterior/.. LALInferenceMCMC.c LALInferenceMCMCSampler.c LALInferenceMCMCSampler.c -> Bayesian Followup, MCMC algorithm PTMCMCAlgorithm() PTMCMCOneStep() -> Metropolis-Hasting algorithm implemented 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 33
basic functions, libraries lal/FrequencySeries.h -> dataStructure lal/LALDatatypes.h ->Primitive datatypes Aggregate datatypes Store arbitrarily large sets or collections of primitive datatypes. aggregate datatypes are defines Structured datatypes : embed primitive and aggregate datatypes inside structures that define their physical meaning. lal/LALConstants.h -> physical constants, Pi, h, G, etc. 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 34
Suggestions for KAGALI Directory structure by functionality MCMC chain generation Waveform generation Input output functions Notes to consider Waveform function call parameter Post process language Python Haskell C and gnuplot 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 35
Discussions Bayesian Parameter Estimation Optimal library structure How to post process? All from scratch? Collaboration meeting from 9 to 11 February at Inje University 2015-2-6 11th KAGRA f2f Meeting, 4-7 February 2015 36
