Otegi M-R Exoplanets¶

DACE - python APIs¶

Authors¶

Jon Fernández Otegi (University of Geneva)

Nicolas Unger (DACE Data Scientist)

Goals¶

  • Show possibilities with exoplanets data and graph you can do

Topics¶

  • Use Exoplanet.query_database() with different filters to plot data related to exoplanet Otegi M-R table

Keywords¶

dace.exoplanet, query_database

Requirements¶

Anaconda 3, Python 3.8 or later, IPython 7.8.0
dace-query 1.1.0 or later

Summary¶

Use dace.exoplanet module to plot data related to exoplanet tables. Show more plots with categorization of planets : earth-like, h2o, fe

Setup¶

In [1]:
# Import the module.
from dace_query.exoplanet import Exoplanet
import matplotlib.pyplot as plt
import matplotlib.ticker as ticker
import numpy as np
import pandas as pd
import numpy.polynomial.polynomial as poly
import scipy.stats as stats
2023-11-22 17:26:53,615 - WARNING - File .dacerc not found. You are requesting data in public mode. To change this behaviour, create a .dacerc file in your home directory and fill it with your API key. More infos on https://dace.unige.ch

First we show the first 10000 planets (default limit) from the entire exoplanet database and list the columns that are present in the database.

In [2]:
# Query the database.
exoplanet_data = Exoplanet.query_database(output_format="pandas")

# Now the data is in a pandas DataFrame. These are the fields we can display.
exoplanet_data.columns.to_list()
Out[2]:
['obj_id_catname',
 'db_info_name',
 'pub_info_reference',
 'obj_orb_period_day',
 'obj_phys_mass_mjup',
 'obj_phys_radius_rjup',
 'obj_phys_mass_rel_unc',
 'obj_phys_radius_rel_unc',
 'obj_orb_period_day_errmax',
 'obj_orb_period_day_errmin',
 'obj_phys_mass_mjup_errmax',
 'obj_phys_mass_mjup_errmin',
 'obj_phys_radius_rjup_errmax',
 'obj_phys_radius_rjup_errmin',
 'pub_info_discovered_year',
 'sys_nplanets',
 'obj_id_pletter',
 'obj_id_parent_name',
 'obj_id_tic',
 'obj_id_gaiadr3',
 'pub_info_detectiontype',
 'pub_info_disc_facility',
 'sys_nstars',
 'pub_ads_link',
 'pub_info_updated',
 'obj_orb_a_au',
 'obj_orb_a_au_errmax',
 'obj_orb_a_au_errmin',
 'obj_orb_a_au_ref',
 'obj_orb_ecc',
 'obj_orb_ecc_errmax',
 'obj_orb_ecc_errmin',
 'obj_orb_ecc_flag',
 'obj_orb_inc_deg',
 'obj_orb_inc_deg_errmax',
 'obj_orb_inc_deg_errmin',
 'obj_orb_tperi_day',
 'obj_orb_tperi_day_errmax',
 'obj_orb_tperi_day_errmin',
 'obj_orb_omega_deg',
 'obj_orb_omega_deg_errmax',
 'obj_orb_omega_deg_errmin',
 'obj_orb_trandur_hour',
 'obj_orb_trandur_hour_errmax',
 'obj_orb_trandur_hour_errmin',
 'obj_orb_tranmid_day',
 'obj_orb_tranmid_day_errmax',
 'obj_orb_tranmid_day_errmin',
 'obj_phys_msini_mjup',
 'obj_phys_msini_mjup_errmax',
 'obj_phys_msini_mjup_errmin',
 'obj_insolation_earthfl',
 'obj_insolation_earthfl_errmax',
 'obj_insolation_earthfl_errmin',
 'obj_phys_eq_temp_k',
 'obj_phys_eq_temp_k_errmax',
 'obj_phys_eq_temp_k_errmin',
 'obj_phys_lambda_deg',
 'obj_phys_lambda_deg_errmax',
 'obj_phys_lambda_deg_errmin',
 'obj_phys_lambda_deg_ref',
 'obj_phys_b',
 'obj_phys_b_errmax',
 'obj_phys_b_errmin',
 'obj_phys_b_ref',
 'obj_phys_trandep',
 'obj_phys_trandep_errmax',
 'obj_phys_trandep_errmin',
 'obj_phys_trandep_ref',
 'obj_orb_k_mps',
 'obj_orb_k_mps_errmax',
 'obj_orb_k_mps_errmin',
 'obj_stellar_spectype',
 'obj_parent_phys_teff_k',
 'obj_parent_phys_teff_k_errmax',
 'obj_parent_phys_teff_k_errmin',
 'obj_parent_phys_radius_rsun',
 'obj_parent_phys_radius_rsun_errmax',
 'obj_parent_phys_radius_rsun_errmin',
 'obj_parent_phys_mass_msun',
 'obj_parent_phys_mass_msun_errmax',
 'obj_parent_phys_mass_msun_errmin',
 'obj_parent_phys_met',
 'obj_parent_phys_met_errmax',
 'obj_parent_phys_met_errmin',
 'obj_parent_phys_met_ratio',
 'obj_parent_phys_lum_lsun',
 'obj_parent_phys_lum_lsun_errmax',
 'obj_parent_phys_lum_lsun_errmin',
 'obj_parent_phys_logg',
 'obj_parent_phys_logg_errmax',
 'obj_parent_phys_logg_errmin',
 'obj_parent_age_gyr',
 'obj_parent_age_gyr_errmax',
 'obj_parent_age_gyr_errmin',
 'obj_parent_mv',
 'obj_parent_mv_errmax',
 'obj_parent_mv_errmin',
 'obj_pos_alpha_deg',
 'obj_pos_alpha_deg_err',
 'obj_pos_delta_deg',
 'obj_pos_delta_deg_err',
 'obj_orb_stardistance_pc',
 'obj_orb_stardistance_pc_errmax',
 'obj_orb_stardistance_pc_errmin',
 'obj_pos_parallax_mas',
 'obj_pos_parallax_mas_err',
 'obj_parent_rd_kms',
 'obj_parent_rv_kms_err',
 'obj_parent_gaia_teff_k',
 'obj_parent_gaia_teff_k_errmax',
 'obj_parent_gaia_teff_k_errmin',
 'obj_parent_gaia_feh',
 'obj_parent_gaia_feh_errmax',
 'obj_parent_gaia_feh_errmin',
 'obj_parent_gaia_lum_lsun',
 'obj_parent_gaia_lum_lsun_errmax',
 'obj_parent_gaia_lum_lsun_errmin',
 'obj_parent_gaia_logg',
 'obj_parent_gaia_logg_errmax',
 'obj_parent_gaia_logg_errmin',
 'obj_parent_gaia_age_gyr',
 'obj_parent_gaia_age_gyr_errmax',
 'obj_parent_gaia_age_gyr_errmin',
 'obj_comput_k_mps',
 'obj_phys_density_gpcm3',
 'obj_phys_density_gpcm3_errmax',
 'obj_phys_density_gpcm3_errmin',
 'obj_comput_insolation_earthfl',
 'obj_comput_insolation_earthfl_errmax',
 'obj_comput_insolation_earthfl_errmin',
 'obj_comput_teq_k',
 'obj_comput_teq_k_errmax',
 'obj_comput_teq_k_errmin']

Plotting¶

From this we plot the planets that have information on both the mass and radius

In [3]:
# First we extract the planets for which the mass and radius are not nan
exoplanets_with_massradius = exoplanet_data[
    ~exoplanet_data["obj_phys_mass_mjup"].isna()
    & ~exoplanet_data["obj_phys_radius_rjup"].isna()
]

# Plot the mass and radius of the planets.
plt.plot(
    exoplanets_with_massradius["obj_phys_mass_mjup"] * 317.828,
    exoplanets_with_massradius["obj_phys_radius_rjup"] * 11.209,
    ".",
)
plt.xlabel("Mass [Mearth]")
plt.ylabel("Radius [Rearth]")
plt.xscale("log")
plt.yscale("log")

Now let's do it by separating the NASA Exoplanet Archive and the Exolpanet.eu archives

In [4]:
# Query the NASA Exoplanet Archive and Exoplanet.eu databases.
exoplanet_data_NASA = Exoplanet.query_database(
    filters={"db_info_name": {"contains": ["Nasa Exoplanet Archive"]}},
    output_format="pandas",
)
exoplanet_data_EU = Exoplanet.query_database(
    filters={"db_info_name": {"contains": ["Exoplanet.eu"]}}, output_format="pandas"
)
In [5]:
# From the archives extract those that have values for mass and radius.
NASA_with_massradius = exoplanet_data_NASA[
    ~exoplanet_data_NASA["obj_phys_mass_mjup"].isna()
    & ~exoplanet_data_NASA["obj_phys_radius_rjup"].isna()
]

EU_with_massradius = exoplanet_data_EU[
    ~exoplanet_data_EU["obj_phys_mass_mjup"].isna()
    & ~exoplanet_data_EU["obj_phys_radius_rjup"].isna()
]

plt.figure(figsize=(10, 8))

# Plot the NASA and EU exoplanets in a scatter plot.
plt.plot(
    NASA_with_massradius["obj_phys_mass_mjup"] * 317.828,  # Convert to M_earth
    NASA_with_massradius["obj_phys_radius_rjup"] * 11.2089,  # Conver to R_earth
    ".",
    c="red",
    label="Nasa Exoplanet Archive",
)
plt.plot(
    EU_with_massradius["obj_phys_mass_mjup"] * 317.828,  # Convert to M_earth
    EU_with_massradius["obj_phys_radius_rjup"] * 11.209,  # Conver to R_earth
    ".",
    c="navy",
    label="Exoplanet.eu",
)
plt.xscale("log")
plt.yscale("log")
plt.ylabel(r"Radius [$R_{\oplus}$]", fontsize=19)
plt.xlabel(r"Mass [$M_{\oplus}$]", fontsize=19)
plt.tick_params(axis="both", which="major", labelsize=20)
plt.legend(fontsize=20)
Out[5]:
<matplotlib.legend.Legend at 0x7f0444135e50>

PlanetS catalog¶

Now we're going to do the same but using a revisited catalog with robust measurements in mass and radius.

In [6]:
# Query the exoplanet database which is named "Otegi M-R".
exoplanet_data_PlanetS = Exoplanet.query_database(
    filters={
        "db_info_name": {"contains": ["Otegi M-R"]},
    },
    output_format="pandas",
)
In [7]:
# Initialize the figure.
fig = plt.figure(figsize=(19, 15))
grid = plt.GridSpec(20, 20)
axarr = fig.add_subplot(grid[3:, :17])
hist1 = fig.add_subplot(grid[:3, :17])
hist2 = fig.add_subplot(grid[3:, 17:])
fig.subplots_adjust(wspace=0.0)
fig.subplots_adjust(hspace=0.0)

# Extract mass and radius and convert to Earth units.
m = exoplanet_data_PlanetS["obj_phys_mass_mjup"] * 317.828  # Convert to M_earth
mlu = exoplanet_data_PlanetS["obj_phys_mass_mjup_errmin"] * 317.828
muu = exoplanet_data_PlanetS["obj_phys_mass_mjup_errmax"] * 317.828
r = exoplanet_data_PlanetS["obj_phys_radius_rjup"] * 11.209  # Convert to R_earth
rlu = exoplanet_data_PlanetS["obj_phys_radius_rjup_errmin"] * 11.209
ruu = exoplanet_data_PlanetS["obj_phys_radius_rjup_errmax"] * 11.209

# Plot the exoplanets.
markers, caps, bars = axarr.errorbar(
    m,
    r,
    xerr=[mlu.values, muu.values],
    yerr=[rlu.values, ruu.values],
    fmt="o",
    alpha=0.5,
    markersize=8.5,
    color="navy",
    ecolor="navy",
    capthick=1,
    capsize=2,
    label="RV",
)


# Solar System Planets
Msol = [0.055, 0.815, 1, 0.11, 95.2, 14.6, 17.2]
Rsol = [0.38, 0.95, 1, 0.53, 9.14, 3.98, 3.86]
axarr.scatter(Msol, Rsol, s=450, c="green", marker="*", label="Solar System Planets")
axarr.text(1.05, 0.9, "$E$", fontsize=20, color="red")
axarr.text(0.73, 0.99, "$V$", fontsize=20, color="red")
axarr.text(12.3, 3.95, "$U$", fontsize=20, color="red")
axarr.text(18.5, 3.6, "$N$", fontsize=20, color="red")
axarr.text(87, 8.1, "$S$", fontsize=20, color="red")
axarr.tick_params(axis="both", which="major", labelsize=20)
axarr.set_xscale("log")
axarr.set_yscale("log")
axarr.set_xlim(0.2, 20000)
axarr.set_ylim(0.5, 30)
axarr.set_ylabel(r"Radius [$R_{\oplus}$]", fontsize=27)
axarr.set_xlabel(r"Mass [$M_{\oplus}$]", fontsize=27)
axarr.yaxis.set_major_formatter(ticker.StrMethodFormatter("{x:g}"))
axarr.xaxis.set_major_formatter(ticker.StrMethodFormatter("{x:g}"))

# Add legend
axarr.legend(bbox_to_anchor=(0.02, 0.98), loc=2, borderaxespad=0.0, fontsize=27)

# Add histograms
hist1.hist(
    np.log10(np.array(exoplanet_data_PlanetS["obj_phys_mass_mjup"]) * 317.828),
    bins=15,
    color="black",
    alpha=0.32,
)
hist2.hist(
    np.log10(np.array(exoplanet_data_PlanetS["obj_phys_radius_rjup"]) * 11.209),
    bins=15,
    color="black",
    alpha=0.25,
    orientation="horizontal",
)

# Clean up the axes of the histograms
hist1.spines["top"].set_visible(False)
hist1.spines["right"].set_visible(False)
hist1.spines["left"].set_visible(False)
hist1.set_xticks([])
hist1.set_yticks([])
hist2.spines["top"].set_visible(False)
hist2.spines["right"].set_visible(False)
hist2.spines["bottom"].set_visible(False)
hist2.set_xticks([])
hist2.set_yticks([])


# Contour lines.
xmin = np.log10(0.04)
xmax = np.log10(10000)
ymin = np.log10(0.5)
ymax = np.log10(30)


def density_estimation(m1, m2):
    X, Y = np.mgrid[xmin:xmax:100j, ymin:ymax:50j]
    positions = np.vstack([X.ravel(), Y.ravel()])
    values = np.vstack([m1, m2])
    kernel = stats.gaussian_kde(values)
    Z = np.reshape(kernel(positions).T, X.shape)
    return X, Y, Z


X, Y, Z = density_estimation(
    np.log10(np.array(exoplanet_data_PlanetS["obj_phys_mass_mjup"]) * 317.828),
    np.log10(np.array(exoplanet_data_PlanetS["obj_phys_radius_rjup"]) * 11.209),
)
axarr.contour(10**X, 10**Y, Z, 30, colors="black", alpha=0.2)
Out[7]:
<matplotlib.contour.QuadContourSet at 0x7f044415d950>

(c) 2023 Data Analysis Center for Exoplanets, Univerity of Geneva, Switzerland