Quickstart - GW190521 BH1

In this notebook, we perform a quick demonstration of using archeo to analyze the GW190521 event, using an agnostic prior similar to APq6 in Araújo Álvarez et al., 2024. We will use a subsampled version of the LVK posterior samples for GW190521 as our input data. The full dataset can be found dcc.ligo.

The notebook run should take less than 3 minutes to complete. However, the values are reasonably close to the reported values of circular APq6 setting in Araújo Álvarez et al., 2024.

In [1]:

Copied!

import pandas as pd
import archeo
import pandas as pd
import archeo

In [2]:

Copied!





import warnings

warnings.filterwarnings("ignore", category=UserWarning)
warnings.filterwarnings("ignore", category=FutureWarning)
warnings.filterwarnings("ignore", category=RuntimeWarning)
warnings.filterwarnings("ignore", category=SyntaxWarning)
import warnings

warnings.filterwarnings("ignore", category=UserWarning)
warnings.filterwarnings("ignore", category=FutureWarning)
warnings.filterwarnings("ignore", category=RuntimeWarning)
warnings.filterwarnings("ignore", category=SyntaxWarning)

Prepare Prior¶

In [3]:

Copied!

pipeline = archeo.get_binary_generation_pipeline("agnostic_aligned_spin")
df_binaries, binary_generator = pipeline(size=5000, n_workers=1)
pipeline = archeo.get_binary_generation_pipeline("agnostic_aligned_spin")
df_binaries, binary_generator = pipeline(size=5000, n_workers=1)

Loaded NRSur3dq8Remnant fit.

100%|██████████| 5000/5000 [00:05<00:00, 834.51it/s]

Load Parameter Estimation Samples¶

In [4]:

Copied!





# This data file is a subset of the LVK posterior samples for GW190521
# See https://github.com/wyhwong/archeo/tree/main/src/tests/test_data/gw190521_lvk_subsampled.json
# For the full dataset, please refer to: https://dcc.ligo.org/P2000158/public

gw190521_pe_samples = pd.read_json("./gw190521_lvk_subsampled.json")
gw190521_pe_samples.head(10)
# This data file is a subset of the LVK posterior samples for GW190521
# See https://github.com/wyhwong/archeo/tree/main/src/tests/test_data/gw190521_lvk_subsampled.json
# For the full dataset, please refer to: https://dcc.ligo.org/P2000158/public

gw190521_pe_samples = pd.read_json("./gw190521_lvk_subsampled.json")
gw190521_pe_samples.head(10)

Out[4]:

	mass_1_source	mass_2_source	a_1	a_2
0	84.513663	76.297304	0.687007	0.899445
1	81.059729	69.663127	0.769092	0.878550
2	92.370429	86.105830	0.540142	0.949552
3	77.363091	67.721241	0.088123	0.973257
4	99.454443	41.802952	0.710768	0.976182
5	89.774109	72.749785	0.096520	0.859641
6	65.693613	63.273635	0.664623	0.867052
7	85.314123	58.547913	0.604760	0.867737
8	88.031480	85.329605	0.741590	0.659539
9	71.147886	65.844700	0.914395	0.073000

Infer ancestor samples¶

In [5]:

Copied!





ancestors_bh1 = archeo.infer_ancestral_posterior_distribution(
    df_binaries=df_binaries,
    mass_posterior_samples=gw190521_pe_samples["mass_1_source"].values.tolist(),
    spin_posterior_samples=gw190521_pe_samples["a_1"].values.tolist(),
    n_workers=1,
)
ancestors_bh1.head(10)
ancestors_bh1 = archeo.infer_ancestral_posterior_distribution(
    df_binaries=df_binaries,
    mass_posterior_samples=gw190521_pe_samples["mass_1_source"].values.tolist(),
    spin_posterior_samples=gw190521_pe_samples["a_1"].values.tolist(),
    n_workers=1,
)
ancestors_bh1.head(10)

100%|██████████| 5000/5000 [00:05<00:00, 939.84it/s]

Out[5]:

	m_1	a_1	v_1	m_2	a_2	v_2	m_f	a_f	k_f	chi_eff	chi_p	q	logL	spin_measure	mass_measure
0	55.822170	0.092123	0.0	31.736393	0.407894	0.0	83.991408	0.640499	88.203031	-0.089113	0.0	1.758932	-7.418581	0.687007	84.513663
1	63.632346	0.572089	0.0	21.389795	0.347842	0.0	81.274261	0.794910	53.277174	0.515674	0.0	2.974893	-7.418581	0.769092	81.059729
2	68.421455	0.026867	0.0	26.824326	0.352231	0.0	92.187994	0.572535	126.200029	-0.079900	0.0	2.550724	-8.517193	0.540142	92.370429
3	NaN	NaN	NaN	NaN	NaN	NaN	NaN	NaN	NaN	NaN	NaN	NaN	-inf	0.088123	77.363091
4	68.372545	0.188975	0.0	35.449756	0.258300	0.0	99.142194	0.702872	113.878567	0.212646	0.0	1.928717	-7.418581	0.710768	99.454443
5	73.428427	0.681282	0.0	18.536307	0.348916	0.0	90.502168	0.101371	230.651833	-0.614291	0.0	3.961330	-8.517193	0.096520	89.774109
6	51.148511	0.266861	0.0	16.799754	0.217357	0.0	65.736130	0.648159	107.681783	0.147141	0.0	3.044599	-7.418581	0.664623	65.693613
7	65.641553	0.194958	0.0	21.421624	0.563402	0.0	84.443919	0.605355	107.217285	0.008366	0.0	3.064266	-7.130899	0.604760	85.314123
8	70.755282	0.412839	0.0	19.380112	0.214790	0.0	87.333794	0.692299	86.475753	0.277892	0.0	3.650922	-7.418581	0.741590	88.031480
9	58.981420	0.945017	0.0	17.225650	0.363833	0.0	72.035605	0.930945	51.036801	0.649168	0.0	3.424046	-6.907755	0.914395	71.147886

Visualize Prior¶

In [6]:

Copied!





from archeo.visualization.distribution import (
    distribution_summary,
    kick_against_spin_cmap,
    kick_distribution_on_spin,
)
from archeo.visualization.distribution import (
    distribution_summary,
    kick_against_spin_cmap,
    kick_distribution_on_spin,
)

In [7]:

Copied!

distribution_summary(df_binaries, close=False)
distribution_summary(df_binaries, close=False)

Out[7]:

(<Figure size 600x800 with 4 Axes>,
 array([<Axes: xlabel='Parent Mass Ratio $q$', ylabel='PDF'>,
        <Axes: xlabel='Remnant Mass $m_f$ [$M_{\\odot}$]', ylabel='PDF'>,
        <Axes: xlabel='Recoil Kick $v_f$ [$kms^{-1}$]', ylabel='PDF'>,
        <Axes: xlabel='Spin $\\chi_f$', ylabel='PDF'>], dtype=object))

No description has been provided for this image

In [8]:

Copied!

kick_against_spin_cmap(df_binaries, close=False)
kick_against_spin_cmap(df_binaries, close=False)

Out[8]:

(<Figure size 800x600 with 1 Axes>,
 <ScatterDensityAxes: title={'center': 'Remnant Kick against Remnant Spin'}, xlabel='Remnant Spin $\\chi_f$', ylabel='Remnant Kick $v_f$ [$kms^{-1}$]'>)

In [9]:

Copied!

kick_distribution_on_spin(df_binaries, close=False)
kick_distribution_on_spin(df_binaries, close=False)

Out[9]:

(<Figure size 900x400 with 1 Axes>, <Axes: >)

Visualize Inferred Posteriors of BH1's Ancestor¶

In [10]:

Copied!





from archeo.visualization.estimation import (
    mass_estimates,
    corner_estimates,
    second_generation_probability_curve,
    effective_spin_estimates,
    precession_spin_estimates,
    table_estimates,
)
from archeo.visualization.estimation import (
    mass_estimates,
    corner_estimates,
    second_generation_probability_curve,
    effective_spin_estimates,
    precession_spin_estimates,
    table_estimates,
)

In [11]:

Copied!

mass_estimates(ancestors_bh1, label="GW190521 BH1", close=False)
mass_estimates(ancestors_bh1, label="GW190521 BH1", close=False)

Out[11]:

(<Figure size 900x400 with 1 Axes>, <Axes: >)

In [12]:

Copied!





corner_estimates(
    dfs={
        "Prior (Aligned Spin)": df_binaries,
        "GW190521 BH1 Ancestor": ancestors_bh1,
    },
    close=False,
)
corner_estimates(
    dfs={
        "Prior (Aligned Spin)": df_binaries,
        "GW190521 BH1 Ancestor": ancestors_bh1,
    },
    close=False,
)

WARNING:root:Please install arviz to use the advanced features of corner
WARNING:root:Please install arviz to use the advanced features of corner
WARNING:root:Please install arviz to use the advanced features of corner
WARNING:root:Please install arviz to use the advanced features of corner

Out[12]:

(<Figure size 900x900 with 36 Axes>,
 array([[<Axes: >, <Axes: >, <Axes: >, <Axes: >, <Axes: >, <Axes: >],
        [<Axes: ylabel='$m_2$ [$M_{\\odot}$]'>, <Axes: >, <Axes: >,
         <Axes: >, <Axes: >, <Axes: >],
        [<Axes: ylabel='$m_f$ [$M_{\\odot}$]'>, <Axes: >, <Axes: >,
         <Axes: >, <Axes: >, <Axes: >],
        [<Axes: ylabel='$v_f$ [km s$^{-1}$]'>, <Axes: >, <Axes: >,
         <Axes: >, <Axes: >, <Axes: >],
        [<Axes: ylabel='$a_f$'>, <Axes: >, <Axes: >, <Axes: >, <Axes: >,
         <Axes: >],
        [<Axes: xlabel='$m_1$ [$M_{\\odot}$]', ylabel='$\\chi_{eff}$'>,
         <Axes: xlabel='$m_2$ [$M_{\\odot}$]'>,
         <Axes: xlabel='$m_f$ [$M_{\\odot}$]'>,
         <Axes: xlabel='$v_f$ [km s$^{-1}$]'>, <Axes: xlabel='$a_f$'>,
         <Axes: xlabel='$\\chi_{eff}$'>]], dtype=object))

In [13]:

Copied!





second_generation_probability_curve(
    dfs={
        "Prior (Aligned Spin)": df_binaries,
        "GW190521 BH1 Ancestor": ancestors_bh1,
    },
    close=False,
)
second_generation_probability_curve(
    dfs={
        "Prior (Aligned Spin)": df_binaries,
        "GW190521 BH1 Ancestor": ancestors_bh1,
    },
    close=False,
)

Out[13]:

(<Figure size 1000x800 with 1 Axes>, <Axes: >)

In [14]:

Copied!





effective_spin_estimates(
    dfs={
        "Prior (Aligned Spin)": df_binaries,
        "GW190521 BH1 Ancestor": ancestors_bh1,
    },
    close=False,
)
effective_spin_estimates(
    dfs={
        "Prior (Aligned Spin)": df_binaries,
        "GW190521 BH1 Ancestor": ancestors_bh1,
    },
    close=False,
)

Out[14]:

(<Figure size 1000x800 with 1 Axes>, <Axes: >)

In [15]:

Copied!





precession_spin_estimates(
    dfs={
        "Prior (Aligned Spin)": df_binaries,
        "GW190521 BH1 Ancestor": ancestors_bh1,
    },
    close=False,
)
precession_spin_estimates(
    dfs={
        "Prior (Aligned Spin)": df_binaries,
        "GW190521 BH1 Ancestor": ancestors_bh1,
    },
    close=False,
)

Out[15]:

(<Figure size 1000x800 with 1 Axes>, <Axes: >)

In [16]:

Copied!





table_estimates(
    dfs={
        "Prior (Aligned Spin)": df_binaries,
        "GW190521 BH1 Ancestor": ancestors_bh1,
    }
)
table_estimates(
    dfs={
        "Prior (Aligned Spin)": df_binaries,
        "GW190521 BH1 Ancestor": ancestors_bh1,
    }
)

Out[16]:

		Recovery Rate	p2g_GC	p2g_MW	p2g_NSC	p2g_EG	$m_1$	$m_2$	$q$	$m_f$	$a_f$	$v_f$	$\chi_{p}$	$\chi_{eff}$
0	Prior (Aligned Spin)	1.00	1.64	10.16	10.16	10.16	$141.97_{-94.32}^{+53.09}$	$72.36_{-53.67}^{+86.65}$	$1.69_{-0.64}^{+3.07}$	$209.93_{-141.29}^{+114.37}$	$0.64_{-0.46}^{+0.25}$	$139.46_{-114.14}^{+202.94}$	$0.00_{-0.00}^{+0.00}$	$0.00_{-0.70}^{+0.71}$
1	GW190521 BH1 Ancestor	0.81	10.32	52.02	52.02	52.02	$60.34_{-14.87}^{+21.81}$	$28.28_{-13.53}^{+18.82}$	$2.12_{-1.01}^{+2.92}$	$85.58_{-14.44}^{+20.54}$	$0.73_{-0.59}^{+0.21}$	$111.46_{-95.08}^{+210.91}$	$0.00_{-0.00}^{+0.00}$	$0.25_{-0.88}^{+0.61}$