CNModel Hyperfine Anomalies

Trey V. Wenger (c) April 2025

CNModel can also account for a full physical treatment of non-LTE and non-CTEX effects. That is, it can model situations where the excitation temperature is not equal to the kinetic temperature, and the excitation temperature is not the same for each transition. CNModel considers the detailed balance between all transitions in the mol_data dictionary. For \({\rm CN}\) this includes all of the \(N=1-0\) (\(\nu = 0\)) hyperfine transitions at 113.15 and 113.5 GHz, and for \(^{13}{\rm CN}\) this includes all of the \(N=1-0\) hyperfine transitions at 108.06 GHz, 108.4 GHz, 108.6 GHz, 108.8 GHz, and 109.2 GHz. Each state’s statistical weight is a free parameter.

We fix the kinetic temperature for these models, since we are unable to constrain the thermal line width.

# General imports
import os
import pickle
import time

import matplotlib.pyplot as plt
import arviz as az
import pandas as pd
import numpy as np
import pymc as pm

print("arviz version:", az.__version__)

print("pymc version:", pm.__version__)

import bayes_spec
print("bayes_spec version:", bayes_spec.__version__)

import bayes_cn_hfs
print("bayes_cn_hfs version:", bayes_cn_hfs.__version__)

# Notebook configuration
pd.options.display.max_rows = None

arviz version: 0.22.0dev
pymc version: 5.21.1
bayes_spec version: 1.7.8
bayes_cn_hfs version: 2.0.0-staging+0.g7afe3a8.dirty

Simulate data

from bayes_cn_hfs import supplement_mol_data
from bayes_spec import SpecData
from bayes_cn_hfs.cn_model import CNModel

mol_data_12CN, mol_weight_12CN = supplement_mol_data("CN")

# spectral axis definition
freq_axis_1 = np.arange(113100.0, 113210.0, 0.2) # MHz
freq_axis_2 = np.arange(113470.0, 113530.0, 0.2) # MHz

# data noise can either be a scalar (assumed constant noise across the spectrum)
# or an array of the same length as the data
noise = 0.03 # K

# brightness data. In this case, we just throw in some random data for now
# since we are only doing this in order to simulate some actual data.
brightness_data_1 = noise * np.random.randn(len(freq_axis_1)) # K
brightness_data_2 = noise * np.random.randn(len(freq_axis_2)) # K

# CNModel datasets can be named anything, here we name them "12CN-1" and "12CN-2"
obs_1 = SpecData(
    freq_axis_1,
    brightness_data_1,
    noise,
    xlabel=r"LSRK Frequency (MHz)",
    ylabel=r"$T_B$ (K)",
)
obs_2 = SpecData(
    freq_axis_2,
    brightness_data_2,
    noise,
    xlabel=r"LSRK Frequency (MHz)",
    ylabel=r"$T_B$ (K)",
)
dummy_data = {"12CN-1": obs_1, "12CN-2": obs_2}

# Initialize and define the model
n_clouds = 3 # number of cloud components
baseline_degree = 0 # polynomial baseline degree
model = CNModel(
    dummy_data,
    molecule="CN", # molecule (either "CN" or "13CN")
    mol_data=mol_data_12CN, # molecular data
    bg_temp = 2.7, # assumed background temperature (K)
    Beff = 1.0, # beam efficiency
    Feff = 1.0, # forward efficiency
    n_clouds=n_clouds,
    baseline_degree=baseline_degree,
    seed=1234,
    verbose=True
)
model.add_priors(
    prior_log10_N = [13.5, 1.0], # mean and width of log10 total column density prior (cm-2)
    prior_log10_Tkin = [1.0, 0.5], # mean and width of log10 kinetic temperature prior (K)
    prior_velocity = [0.0, 3.0], # mean and width of velocity prior (km/s)
    prior_fwhm_nonthermal = 1.0, # width of non-thermal broadening prior (km/s)
    prior_fwhm_L = None, # assume Gaussian line profile
    prior_rms = None, # do not infer spectral rms
    prior_baseline_coeffs = None, # use default baseline priors
    assume_LTE = False, # do not assume LTE
    prior_log10_Tex = [0.5, 0.1], # mean and width of log10 excitation temperature prior (K)
    assume_CTEX = False, # do not assume CTEX
    prior_log10_LTE_precision = [-6.0, 1.0], # offset and width of log10 LTE precision prior
    fix_log10_Tkin = None, # do not fix the kinetic temperature
    clip_weights = 1.0e-9, # clip statistical weights between [clip_weights, 1-clip_weights]
    clip_tau = -10.0, # clip optical depths below to prevent masers
    ordered = False, # do not assume optically-thin
)
model.add_likelihood()

sim_params = {
    "log10_N": [13.8, 13.9, 14.0],
    "log10_Tkin": [0.65, 0.6, 0.5],
    "fwhm_nonthermal": [1.0, 1.25, 1.5],
    "velocity": [-2.0, 0.0, 2.5],
    "fwhm_L": 0.0,
    "log10_Tex_ul": [0.65, 0.6, 0.5],
    "log10_LTE_precision": [-3.0, -4.0, -5.0],
    "baseline_12CN-1_norm": [0.0],
    "baseline_12CN-2_norm": [0.0],
}

# add derived quantities to sim_params
for key in model.cloud_freeRVs + model.cloud_deterministics:
    if key not in sim_params.keys():
        sim_params[key] = model.model[key].eval(sim_params, on_unused_input="ignore")

# Evaluate and save simulated observations
sim_obs1 = model.model["12CN-1"].eval(sim_params, on_unused_input="ignore")
sim_obs2 = model.model["12CN-2"].eval(sim_params, on_unused_input="ignore")

# Plot the simulated data
fig, axes = plt.subplots(2)
axes[0].plot(dummy_data["12CN-1"].spectral, sim_obs1, "k-")
axes[0].set_ylabel(dummy_data["12CN-1"].ylabel)
axes[1].plot(dummy_data["12CN-2"].spectral, sim_obs2, "k-")
axes[1].set_xlabel(dummy_data["12CN-2"].xlabel)
_ = axes[1].set_ylabel(dummy_data["12CN-2"].ylabel)

sim_params

{'log10_N': [13.8, 13.9, 14.0],
 'log10_Tkin': [0.65, 0.6, 0.5],
 'fwhm_nonthermal': [1.0, 1.25, 1.5],
 'velocity': [-2.0, 0.0, 2.5],
 'fwhm_L': 0.0,
 'log10_Tex_ul': [0.65, 0.6, 0.5],
 'log10_LTE_precision': [-3.0, -4.0, -5.0],
 'baseline_12CN-1_norm': [0.0],
 'baseline_12CN-2_norm': [0.0],
 'velocity_norm': array([-0.90725257, -1.08963228,  2.05844579]),
 'log10_Tkin_norm': array([-1.22839729, -0.1382073 ,  2.11081962]),
 'fwhm_nonthermal_norm': array([1.14100359, 0.72364694, 0.19993961]),
 'log10_N_norm': array([-0.19739975,  1.18167855,  0.09680737]),
 'log10_Tex_ul_norm': array([-0.0689338 ,  0.12648815, -1.16190854]),
 'log10_LTE_precision_norm': array([0.21987947, 0.71222726, 2.31970483]),
 'weights': array([[0.17274399, 0.346758  , 0.05395841, 0.10043034, 0.04519759,
         0.11536218, 0.16554947],
        [0.19248324, 0.38011643, 0.04867468, 0.0963955 , 0.0423153 ,
         0.09213565, 0.14787919],
        [0.21685544, 0.43456298, 0.03896576, 0.07729498, 0.0381383 ,
         0.07699799, 0.11718454]]),
 'fwhm_thermal': array([0.08867757, 0.08371703, 0.07461288]),
 'fwhm': array([1.00392416, 1.25280028, 1.50185455]),
 'LTE_weights': array([[0.17659718, 0.353274  , 0.05237634, 0.10469961, 0.05216502,
         0.10434294, 0.15654492],
        [0.18885791, 0.37781138, 0.04829246, 0.09653001, 0.04807389,
         0.09616112, 0.14427323],
        [0.21685988, 0.43385792, 0.0389559 , 0.07785605, 0.03873407,
         0.07748167, 0.11625451]]),
 'Tex': array([[4.66576066, 3.94884168, 3.16278169],
        [4.65194996, 3.98632984, 3.1597509 ],
        [4.39608074, 3.92239606, 3.14905465],
        [4.38386868, 3.9593697 , 3.14605298],
        [4.96550722, 3.80907304, 3.1508594 ],
        [4.75768938, 4.03595283, 3.17396883],
        [4.0626928 , 3.59579242, 3.13411203],
        [4.94985866, 3.84385378, 3.14786242],
        [4.05234646, 3.62680534, 3.13114978]]),
 'tau': array([[0.02659258, 0.04053045, 0.06311716],
        [0.21877377, 0.32607302, 0.5177328 ],
        [0.21925804, 0.32504235, 0.50541061],
        [0.28610139, 0.41485675, 0.65765918],
        [0.26800721, 0.42933666, 0.65841539],
        [0.731141  , 1.09621083, 1.74725104],
        [0.2343882 , 0.34741199, 0.5205156 ],
        [0.20766169, 0.32536017, 0.508576  ],
        [0.02878143, 0.04176322, 0.06375109]]),
 'tau_total': array([2.22070531, 3.34658544, 5.24242887]),
 'TR': array([[2.46615586, 1.83756809, 1.18920706],
        [2.45349372, 1.86943026, 1.18658405],
        [2.22588549, 1.81430469, 1.17787815],
        [2.21481599, 1.84565819, 1.17528473],
        [2.73036931, 1.71369746, 1.17566564],
        [2.54295373, 1.9073253 , 1.19378692],
        [1.93025354, 1.53492581, 1.16241827],
        [2.71588106, 1.74287807, 1.17307884],
        [1.92104267, 1.56040541, 1.15986603]])}

Model Definition

# Now we pack the simulated spectrum into a new SpecData instance
obs_1 = SpecData(
    freq_axis_1,
    sim_obs1,
    noise,
    xlabel=r"LSRK Frequency (MHz)",
    ylabel=r"$T_B$ (K)",
)
obs_2 = SpecData(
    freq_axis_2,
    sim_obs2,
    noise,
    xlabel=r"LSRK Frequency (MHz)",
    ylabel=r"$T_B$ (K)",
)
data = {"12CN-1": obs_1, "12CN-2": obs_2}

# Initialize and define the model
model = CNModel(
    data,
    molecule="CN", # molecule (either "CN" or "13CN")
    mol_data=mol_data_12CN, # molecular data
    bg_temp = 2.7, # assumed background temperature (K)
    Beff = 1.0, # beam efficiency
    Feff = 1.0, # forward efficiency
    n_clouds=n_clouds,
    baseline_degree=baseline_degree,
    seed=1234,
    verbose=True
)
model.add_priors(
    prior_log10_N = [13.5, 1.0], # mean and width of log10 total column density prior (cm-2)
    prior_log10_Tkin = None, # ignored (kinetic temperature is fixed)
    prior_velocity = [0.0, 3.0], # mean and width of velocity prior (km/s)
    prior_fwhm_nonthermal = 1.0, # width of non-thermal broadening prior (km/s)
    prior_fwhm_L = None, # assume Gaussian line profile
    prior_rms = None, # do not infer spectral rms
    prior_baseline_coeffs = None, # use default baseline priors
    assume_LTE = False, # do not assume LTE
    prior_log10_Tex = [0.5, 0.1], # mean and width of log10 excitation temperature prior (K)
    assume_CTEX = False, # do not assume CTEX
    prior_log10_LTE_precision = [-6.0, 1.0], # offset and width of log10 LTE precision prior
    fix_log10_Tkin = 0.5, # fix the kinetic temperature
    clip_weights = 1.0e-9, # clip statistical weights between [clip_weights, 1-clip_weights]
    clip_tau = -10.0, # clip optical depths below to prevent masers
    ordered = False, # do not assume optically-thin
)
model.add_likelihood()

from bayes_spec.plots import plot_predictive

# prior predictive check
prior = model.sample_prior_predictive(
    samples=1000,  # prior predictive samples
)
_ = plot_predictive(model.data, prior.prior_predictive)

Sampling: [12CN-1, 12CN-2, baseline_12CN-1_norm, baseline_12CN-2_norm, fwhm_nonthermal_norm, log10_LTE_precision_norm, log10_N_norm, log10_Tex_ul_norm, velocity_norm, weights]

from bayes_spec.plots import plot_pair

# available parameter attributes:
print("baseline_freeRVs", model.baseline_freeRVs)
print("baseline_deterministics", model.baseline_deterministics)
print("cloud_freeRVs", model.cloud_freeRVs)
print("cloud_deterministics", model.cloud_deterministics)
print("hyper_freeRVs", model.hyper_freeRVs)
print("hyper_deterministics", model.hyper_deterministics)

# ignore transition and state dependent parameters
var_names = [
    param for param in model.cloud_deterministics
    if not set(model.model.named_vars_to_dims[param]).intersection(set(["transition", "state"]))
    and param not in ["fwhm_thermal"]
]
_ = plot_pair(
    prior.prior, # samples
    var_names, # var_names to plot
    labeller=model.labeller, # label manager
    kind="kde", # plot type
    reference_values=sim_params, # truths
)

baseline_freeRVs ['baseline_12CN-1_norm', 'baseline_12CN-2_norm']
baseline_deterministics []
cloud_freeRVs ['velocity_norm', 'fwhm_nonthermal_norm', 'log10_N_norm', 'log10_Tex_ul_norm', 'log10_LTE_precision_norm', 'weights']
cloud_deterministics ['velocity', 'fwhm_thermal', 'fwhm_nonthermal', 'fwhm', 'log10_N', 'log10_Tex_ul', 'LTE_weights', 'log10_LTE_precision', 'Tex', 'tau', 'tau_total', 'TR']
hyper_freeRVs []
hyper_deterministics []

Variational Inference

start = time.time()
model.fit(
    n = 100_000, # maximum number of VI iterations
    draws = 1_000, # number of posterior samples
    rel_tolerance = 0.05, # VI relative convergence threshold
    abs_tolerance = 0.05, # VI absolute convergence threshold
    learning_rate = 0.01, # VI learning rate
)
end = time.time()
print(f"Runtime: {(end-start)/60.0:.2f} minutes")

Convergence achieved at 12500
Interrupted at 12,499 [12%]: Average Loss = 1.9687e+24

Runtime: 0.38 minutes

posterior = model.sample_posterior_predictive(
    thin=10, # keep one in {thin} posterior samples
)
_ = plot_predictive(model.data, posterior.posterior_predictive)

Sampling: [12CN-1, 12CN-2]

MCMC

start = time.time()
init_kwargs = {
    "rel_tolerance": 0.05,
    "abs_tolerance": 0.05,
    "learning_rate": 0.01,
}
model.sample(
    init = "advi+adapt_diag", # initialization strategy
    tune = 1000, # tuning samples
    draws = 1000, # posterior samples
    chains = 8, # number of independent chains
    cores = 8, # number of parallel chains
    init_kwargs = init_kwargs, # VI initialization arguments
    nuts_kwargs = {"target_accept": 0.9}, # NUTS arguments
)
end = time.time()
print(f"Runtime: {(end-start)/60.0:.2f} minutes")

Initializing NUTS using custom advi+adapt_diag strategy

Convergence achieved at 12500
Interrupted at 12,499 [12%]: Average Loss = 1.9687e+24
Multiprocess sampling (8 chains in 8 jobs)
NUTS: [baseline_12CN-1_norm, baseline_12CN-2_norm, velocity_norm, fwhm_nonthermal_norm, log10_N_norm, log10_Tex_ul_norm, log10_LTE_precision_norm, weights]

Sampling 8 chains for 1_000 tune and 1_000 draw iterations (8_000 + 8_000 draws total) took 257 seconds.

Adding log-likelihood to trace

Runtime: 4.78 minutes

model.solve(kl_div_threshold=0.1)

GMM converged to unique solution

posterior = model.sample_posterior_predictive(
    thin=10, # keep one in {thin} posterior samples
)
_ = plot_predictive(model.data, posterior.posterior_predictive)

Sampling: [12CN-1, 12CN-2]

print("solutions:", model.solutions)

pm.summary(model.trace.solution_0)

solutions: [0]

/home/twenger/miniconda3/envs/bayes_spec-dev/lib/python3.13/site-packages/arviz/stats/diagnostics.py:596: RuntimeWarning: invalid value encountered in scalar divide
  (between_chain_variance / within_chain_variance + num_samples - 1) / (num_samples)
/home/twenger/miniconda3/envs/bayes_spec-dev/lib/python3.13/site-packages/arviz/stats/diagnostics.py:991: RuntimeWarning: invalid value encountered in scalar divide
  varsd = varvar / evar / 4

	mean	sd	hdi_3%	hdi_97%	mcse_mean	mcse_sd	ess_bulk	ess_tail	r_hat
baseline_12CN-1_norm[0]	-0.076	0.044	-0.160	0.005	0.000	0.001	8482.0	5787.0	1.0
baseline_12CN-2_norm[0]	-0.435	0.064	-0.556	-0.317	0.001	0.001	9823.0	6121.0	1.0
velocity_norm[0]	0.852	0.013	0.828	0.879	0.000	0.000	7821.0	6369.0	1.0
velocity_norm[1]	-0.002	0.005	-0.011	0.008	0.000	0.000	7305.0	6203.0	1.0
velocity_norm[2]	-0.663	0.003	-0.670	-0.658	0.000	0.000	6095.0	5801.0	1.0
log10_N_norm[0]	0.393	0.082	0.238	0.549	0.002	0.002	2072.0	1675.0	1.0
log10_N_norm[1]	0.363	0.065	0.244	0.485	0.001	0.001	2083.0	2679.0	1.0
log10_N_norm[2]	0.324	0.050	0.230	0.420	0.001	0.001	1761.0	2743.0	1.0
log10_Tex_ul_norm[0]	0.089	0.101	-0.093	0.270	0.002	0.003	1999.0	1546.0	1.0
log10_Tex_ul_norm[1]	1.128	0.221	0.733	1.560	0.005	0.003	2343.0	2964.0	1.0
log10_Tex_ul_norm[2]	1.546	0.249	1.080	2.012	0.005	0.004	2703.0	3458.0	1.0
fwhm_nonthermal_norm[0]	1.580	0.102	1.389	1.772	0.001	0.001	5191.0	5845.0	1.0
fwhm_nonthermal_norm[1]	1.312	0.058	1.203	1.421	0.001	0.001	3153.0	4874.0	1.0
fwhm_nonthermal_norm[2]	0.980	0.034	0.920	1.049	0.001	0.000	2577.0	4180.0	1.0
log10_LTE_precision_norm[0]	0.780	0.520	0.000	1.662	0.013	0.008	1606.0	1944.0	1.0
log10_LTE_precision_norm[1]	2.244	0.283	1.728	2.790	0.004	0.004	6062.0	5075.0	1.0
log10_LTE_precision_norm[2]	2.692	0.244	2.247	3.151	0.003	0.003	7352.0	5631.0	1.0
weights[0, 0 0 1 1]	0.214	0.003	0.208	0.220	0.000	0.000	2111.0	1604.0	1.0
weights[0, 0 0 1 2]	0.429	0.006	0.418	0.440	0.000	0.000	2076.0	1672.0	1.0
weights[0, 1 0 1 1]	0.040	0.001	0.037	0.042	0.000	0.000	2177.0	1800.0	1.0
weights[0, 1 0 1 2]	0.080	0.002	0.076	0.084	0.000	0.000	2005.0	1545.0	1.0
weights[0, 1 0 2 1]	0.040	0.001	0.037	0.042	0.000	0.000	2096.0	1630.0	1.0
weights[0, 1 0 2 2]	0.079	0.002	0.075	0.083	0.000	0.000	2217.0	1872.0	1.0
weights[0, 1 0 2 3]	0.119	0.003	0.115	0.124	0.000	0.000	2258.0	1926.0	1.0
weights[1, 0 0 1 1]	0.189	0.006	0.177	0.200	0.000	0.000	2575.0	3463.0	1.0
weights[1, 0 0 1 2]	0.368	0.010	0.348	0.385	0.000	0.000	2117.0	2600.0	1.0
weights[1, 1 0 1 1]	0.051	0.003	0.047	0.056	0.000	0.000	2291.0	2986.0	1.0
weights[1, 1 0 1 2]	0.101	0.005	0.092	0.109	0.000	0.000	2122.0	2773.0	1.0
weights[1, 1 0 2 1]	0.044	0.002	0.040	0.048	0.000	0.000	2527.0	3518.0	1.0
weights[1, 1 0 2 2]	0.097	0.004	0.089	0.104	0.000	0.000	2173.0	2754.0	1.0
weights[1, 1 0 2 3]	0.150	0.003	0.145	0.156	0.000	0.000	2661.0	3498.0	1.0
weights[2, 0 0 1 1]	0.178	0.005	0.168	0.188	0.000	0.000	2342.0	3707.0	1.0
weights[2, 0 0 1 2]	0.346	0.009	0.329	0.363	0.000	0.000	1983.0	2612.0	1.0
weights[2, 1 0 1 1]	0.054	0.002	0.049	0.058	0.000	0.000	2073.0	3104.0	1.0
weights[2, 1 0 1 2]	0.098	0.003	0.092	0.104	0.000	0.000	1987.0	2931.0	1.0
weights[2, 1 0 2 1]	0.045	0.002	0.042	0.049	0.000	0.000	2757.0	4180.0	1.0
weights[2, 1 0 2 2]	0.116	0.005	0.107	0.125	0.000	0.000	1905.0	2572.0	1.0
weights[2, 1 0 2 3]	0.163	0.002	0.159	0.167	0.000	0.000	2504.0	3299.0	1.0
velocity[0]	2.557	0.040	2.484	2.636	0.000	0.000	7821.0	6369.0	1.0
velocity[1]	-0.005	0.015	-0.034	0.023	0.000	0.000	7305.0	6203.0	1.0
velocity[2]	-1.990	0.010	-2.009	-1.973	0.000	0.000	6095.0	5801.0	1.0
fwhm_thermal[0]	0.075	0.000	0.075	0.075	0.000	NaN	8000.0	8000.0	NaN
fwhm_thermal[1]	0.075	0.000	0.075	0.075	0.000	NaN	8000.0	8000.0	NaN
fwhm_thermal[2]	0.075	0.000	0.075	0.075	0.000	NaN	8000.0	8000.0	NaN
fwhm_nonthermal[0]	1.580	0.102	1.389	1.772	0.001	0.001	5191.0	5845.0	1.0
fwhm_nonthermal[1]	1.312	0.058	1.203	1.421	0.001	0.001	3153.0	4874.0	1.0
fwhm_nonthermal[2]	0.980	0.034	0.920	1.049	0.001	0.000	2577.0	4180.0	1.0
fwhm[0]	1.581	0.102	1.391	1.774	0.001	0.001	5191.0	5845.0	1.0
fwhm[1]	1.314	0.058	1.205	1.423	0.001	0.001	3153.0	4874.0	1.0
fwhm[2]	0.983	0.034	0.923	1.052	0.001	0.000	2577.0	4180.0	1.0
log10_N[0]	13.893	0.082	13.738	14.049	0.002	0.002	2072.0	1675.0	1.0
log10_N[1]	13.863	0.065	13.744	13.985	0.001	0.001	2083.0	2679.0	1.0
log10_N[2]	13.824	0.050	13.730	13.920	0.001	0.001	1761.0	2743.0	1.0
log10_Tex_ul[0]	0.509	0.010	0.491	0.527	0.000	0.000	1999.0	1546.0	1.0
log10_Tex_ul[1]	0.613	0.022	0.573	0.656	0.000	0.000	2343.0	2964.0	1.0
log10_Tex_ul[2]	0.655	0.025	0.608	0.701	0.000	0.000	2703.0	3458.0	1.0
LTE_weights[0, 0 0 1 1]	0.214	0.003	0.209	0.220	0.000	0.000	1999.0	1546.0	1.0
LTE_weights[0, 0 0 1 2]	0.429	0.006	0.418	0.439	0.000	0.000	1999.0	1546.0	1.0
LTE_weights[0, 1 0 1 1]	0.040	0.001	0.038	0.042	0.000	0.000	1999.0	1546.0	1.0
LTE_weights[0, 1 0 1 2]	0.080	0.002	0.076	0.083	0.000	0.000	1999.0	1546.0	1.0
LTE_weights[0, 1 0 2 1]	0.040	0.001	0.038	0.041	0.000	0.000	1999.0	1546.0	1.0
LTE_weights[0, 1 0 2 2]	0.079	0.002	0.076	0.083	0.000	0.000	1999.0	1546.0	1.0
LTE_weights[0, 1 0 2 3]	0.119	0.003	0.113	0.124	0.000	0.000	1999.0	1546.0	1.0
LTE_weights[1, 0 0 1 1]	0.186	0.006	0.175	0.196	0.000	0.000	2343.0	2964.0	1.0
LTE_weights[1, 0 0 1 2]	0.372	0.011	0.350	0.392	0.000	0.000	2343.0	2964.0	1.0
LTE_weights[1, 1 0 1 1]	0.049	0.002	0.046	0.053	0.000	0.000	2343.0	2964.0	1.0
LTE_weights[1, 1 0 1 2]	0.099	0.004	0.092	0.106	0.000	0.000	2343.0	2964.0	1.0
LTE_weights[1, 1 0 2 1]	0.049	0.002	0.046	0.053	0.000	0.000	2343.0	2964.0	1.0
LTE_weights[1, 1 0 2 2]	0.098	0.004	0.091	0.105	0.000	0.000	2343.0	2964.0	1.0
LTE_weights[1, 1 0 2 3]	0.147	0.006	0.137	0.158	0.000	0.000	2343.0	2964.0	1.0
LTE_weights[2, 0 0 1 1]	0.176	0.006	0.165	0.187	0.000	0.000	2703.0	3458.0	1.0
LTE_weights[2, 0 0 1 2]	0.351	0.011	0.331	0.374	0.000	0.000	2703.0	3458.0	1.0
LTE_weights[2, 1 0 1 1]	0.053	0.002	0.049	0.056	0.000	0.000	2703.0	3458.0	1.0
LTE_weights[2, 1 0 1 2]	0.105	0.004	0.098	0.112	0.000	0.000	2703.0	3458.0	1.0
LTE_weights[2, 1 0 2 1]	0.052	0.002	0.049	0.056	0.000	0.000	2703.0	3458.0	1.0
LTE_weights[2, 1 0 2 2]	0.105	0.004	0.098	0.112	0.000	0.000	2703.0	3458.0	1.0
LTE_weights[2, 1 0 2 3]	0.157	0.006	0.146	0.168	0.000	0.000	2703.0	3458.0	1.0
log10_LTE_precision[0]	-5.220	0.520	-6.000	-4.338	0.013	0.008	1606.0	1944.0	1.0
log10_LTE_precision[1]	-3.756	0.283	-4.272	-3.210	0.004	0.004	6062.0	5075.0	1.0
log10_LTE_precision[2]	-3.308	0.244	-3.753	-2.849	0.003	0.003	7352.0	5631.0	1.0
Tex[113123.3687, 0]	3.221	0.085	3.070	3.378	0.002	0.003	2064.0	1547.0	1.0
Tex[113123.3687, 1]	4.181	0.266	3.734	4.681	0.006	0.005	2262.0	2816.0	1.0
Tex[113123.3687, 2]	4.555	0.281	4.061	5.065	0.006	0.006	1975.0	2989.0	1.0
Tex[113144.19, 0]	3.216	0.082	3.069	3.379	0.002	0.002	2055.0	1580.0	1.0
Tex[113144.19, 1]	4.273	0.266	3.822	4.765	0.006	0.005	2121.0	2622.0	1.0
Tex[113144.19, 2]	4.672	0.286	4.188	5.209	0.007	0.006	1898.0	2780.0	1.0
Tex[113170.535, 0]	3.244	0.087	3.102	3.398	0.002	0.004	1961.0	1481.0	1.0
Tex[113170.535, 1]	4.119	0.238	3.706	4.562	0.005	0.004	2143.0	2804.0	1.0
Tex[113170.535, 2]	4.219	0.206	3.868	4.608	0.005	0.004	1905.0	2928.0	1.0
Tex[113191.325, 0]	3.239	0.082	3.102	3.388	0.002	0.004	1986.0	1528.0	1.0
Tex[113191.325, 1]	4.208	0.241	3.774	4.646	0.005	0.005	2069.0	2631.0	1.0
Tex[113191.325, 2]	4.321	0.214	3.946	4.709	0.005	0.004	1905.0	2665.0	1.0
Tex[113488.142, 0]	3.221	0.082	3.072	3.367	0.002	0.003	2069.0	1804.0	1.0
Tex[113488.142, 1]	4.009	0.216	3.642	4.418	0.005	0.004	2144.0	2557.0	1.0
Tex[113488.142, 2]	4.851	0.300	4.352	5.409	0.007	0.006	1823.0	2678.0	1.0
Tex[113490.985, 0]	3.229	0.066	3.117	3.353	0.002	0.002	2132.0	1720.0	1.0
Tex[113490.985, 1]	4.197	0.146	3.953	4.471	0.003	0.003	2082.0	2596.0	1.0
Tex[113490.985, 2]	4.710	0.156	4.437	4.988	0.004	0.003	1786.0	2592.0	1.0
Tex[113499.643, 0]	3.235	0.087	3.082	3.388	0.002	0.004	1968.0	1520.0	1.0
Tex[113499.643, 1]	3.742	0.184	3.414	4.083	0.004	0.003	2184.0	2956.0	1.0
Tex[113499.643, 2]	3.974	0.177	3.679	4.315	0.004	0.003	2005.0	2946.0	1.0
Tex[113508.934, 0]	3.215	0.080	3.078	3.372	0.002	0.002	2115.0	1772.0	1.0
Tex[113508.934, 1]	4.094	0.222	3.715	4.508	0.005	0.004	2107.0	2505.0	1.0
Tex[113508.934, 2]	4.986	0.322	4.412	5.558	0.008	0.006	1896.0	2523.0	1.0
Tex[113520.4215, 0]	3.230	0.084	3.082	3.385	0.002	0.004	2020.0	1585.0	1.0
Tex[113520.4215, 1]	3.817	0.193	3.477	4.175	0.004	0.003	2232.0	2781.0	1.0
Tex[113520.4215, 2]	4.065	0.197	3.719	4.426	0.004	0.003	2241.0	3196.0	1.0
tau[113123.3687, 0]	0.049	0.010	0.029	0.068	0.000	0.000	2038.0	1600.0	1.0
tau[113123.3687, 1]	0.036	0.007	0.023	0.050	0.000	0.000	2095.0	2552.0	1.0
tau[113123.3687, 2]	0.030	0.005	0.020	0.039	0.000	0.000	1766.0	2802.0	1.0
tau[113144.19, 0]	0.405	0.086	0.248	0.564	0.002	0.002	2052.0	1656.0	1.0
tau[113144.19, 1]	0.286	0.058	0.179	0.393	0.001	0.001	2055.0	2576.0	1.0
tau[113144.19, 2]	0.234	0.040	0.165	0.313	0.001	0.001	1766.0	2555.0	1.0
tau[113170.535, 0]	0.393	0.083	0.240	0.549	0.002	0.002	2034.0	1596.0	1.0
tau[113170.535, 1]	0.292	0.059	0.185	0.401	0.001	0.001	2089.0	2604.0	1.0
tau[113170.535, 2]	0.247	0.040	0.173	0.323	0.001	0.001	1765.0	2645.0	1.0
tau[113191.325, 0]	0.512	0.108	0.312	0.712	0.002	0.002	2048.0	1616.0	1.0
tau[113191.325, 1]	0.365	0.073	0.235	0.506	0.002	0.001	2053.0	2649.0	1.0
tau[113191.325, 2]	0.308	0.050	0.219	0.407	0.001	0.001	1775.0	2757.0	1.0
tau[113488.142, 0]	0.514	0.109	0.312	0.716	0.002	0.002	2042.0	1600.0	1.0
tau[113488.142, 1]	0.386	0.077	0.245	0.527	0.002	0.001	2089.0	2523.0	1.0
tau[113488.142, 2]	0.300	0.051	0.208	0.399	0.001	0.001	1759.0	2701.0	1.0
tau[113490.985, 0]	1.367	0.287	0.832	1.893	0.006	0.005	2057.0	1679.0	1.0
tau[113490.985, 1]	0.972	0.185	0.640	1.324	0.004	0.003	2059.0	2616.0	1.0
tau[113490.985, 2]	0.784	0.123	0.563	1.023	0.003	0.002	1762.0	2711.0	1.0
tau[113499.643, 0]	0.405	0.086	0.241	0.558	0.002	0.002	2031.0	1568.0	1.0
tau[113499.643, 1]	0.314	0.061	0.203	0.429	0.001	0.001	2079.0	2538.0	1.0
tau[113499.643, 2]	0.261	0.041	0.188	0.344	0.001	0.001	1757.0	2578.0	1.0
tau[113508.934, 0]	0.397	0.084	0.240	0.550	0.002	0.002	2056.0	1669.0	1.0
tau[113508.934, 1]	0.286	0.057	0.181	0.391	0.001	0.001	2057.0	2648.0	1.0
tau[113508.934, 2]	0.221	0.038	0.154	0.297	0.001	0.001	1780.0	2577.0	1.0
tau[113520.4215, 0]	0.050	0.011	0.030	0.069	0.000	0.000	2045.0	1622.0	1.0
tau[113520.4215, 1]	0.037	0.007	0.024	0.051	0.000	0.000	2051.0	2573.0	1.0
tau[113520.4215, 2]	0.031	0.005	0.022	0.041	0.000	0.000	1782.0	2678.0	1.0
tau_total[0]	4.093	0.864	2.490	5.687	0.019	0.017	2048.0	1600.0	1.0
tau_total[1]	2.974	0.584	1.904	4.056	0.013	0.008	2058.0	2575.0	1.0
tau_total[2]	2.417	0.391	1.720	3.191	0.009	0.006	1749.0	2584.0	1.0
TR[113123.3687, 0]	1.236	0.068	1.116	1.361	0.002	0.002	2064.0	1547.0	1.0
TR[113123.3687, 1]	2.040	0.233	1.655	2.480	0.005	0.005	2262.0	2816.0	1.0
TR[113123.3687, 2]	2.369	0.251	1.934	2.826	0.006	0.005	1975.0	2989.0	1.0
TR[113144.19, 0]	1.231	0.065	1.108	1.353	0.002	0.002	2055.0	1580.0	1.0
TR[113144.19, 1]	2.120	0.235	1.729	2.554	0.005	0.005	2121.0	2622.0	1.0
TR[113144.19, 2]	2.473	0.257	2.044	2.957	0.006	0.005	1898.0	2780.0	1.0
TR[113170.535, 0]	1.253	0.070	1.141	1.377	0.002	0.004	1961.0	1481.0	1.0
TR[113170.535, 1]	1.985	0.207	1.631	2.373	0.005	0.004	2143.0	2804.0	1.0
TR[113170.535, 2]	2.072	0.181	1.768	2.414	0.004	0.004	1905.0	2928.0	1.0
TR[113191.325, 0]	1.249	0.066	1.141	1.369	0.002	0.003	1986.0	1528.0	1.0
TR[113191.325, 1]	2.062	0.211	1.688	2.448	0.005	0.004	2069.0	2631.0	1.0
TR[113191.325, 2]	2.160	0.189	1.834	2.503	0.004	0.004	1905.0	2665.0	1.0
TR[113488.142, 0]	1.231	0.065	1.114	1.348	0.002	0.002	2069.0	1804.0	1.0
TR[113488.142, 1]	1.886	0.186	1.571	2.238	0.004	0.004	2144.0	2557.0	1.0
TR[113488.142, 2]	2.629	0.272	2.182	3.135	0.006	0.006	1823.0	2678.0	1.0
TR[113490.985, 0]	1.238	0.053	1.144	1.332	0.001	0.001	2132.0	1720.0	1.0
TR[113490.985, 1]	2.047	0.128	1.836	2.288	0.003	0.003	2082.0	2596.0	1.0
TR[113490.985, 2]	2.501	0.140	2.257	2.751	0.003	0.003	1786.0	2592.0	1.0
TR[113499.643, 0]	1.243	0.070	1.121	1.364	0.002	0.003	1968.0	1520.0	1.0
TR[113499.643, 1]	1.659	0.156	1.386	1.948	0.003	0.003	2184.0	2956.0	1.0
TR[113499.643, 2]	1.855	0.152	1.594	2.139	0.003	0.003	2005.0	2946.0	1.0
TR[113508.934, 0]	1.227	0.063	1.118	1.351	0.001	0.002	2115.0	1772.0	1.0
TR[113508.934, 1]	1.959	0.193	1.634	2.320	0.004	0.004	2107.0	2505.0	1.0
TR[113508.934, 2]	2.750	0.293	2.235	3.272	0.007	0.006	1896.0	2523.0	1.0
TR[113520.4215, 0]	1.238	0.067	1.119	1.359	0.002	0.003	2020.0	1585.0	1.0
TR[113520.4215, 1]	1.721	0.164	1.433	2.024	0.004	0.003	2232.0	2781.0	1.0
TR[113520.4215, 2]	1.933	0.171	1.637	2.247	0.004	0.003	2241.0	3196.0	1.0

from bayes_spec.plots import plot_traces

axes = plot_traces(model.trace.solution_0, model.cloud_freeRVs + model.baseline_freeRVs + model.hyper_freeRVs)
fig = axes.ravel()[0].figure
fig.tight_layout()

from bayes_spec.plots import plot_pair

var_names = [
    param for param in model.cloud_deterministics
    if not set(model.model.named_vars_to_dims[param]).intersection(set(["transition", "state"]))
    and param not in ["fwhm_thermal"]
]
print(var_names)
_ = plot_pair(
    model.trace.solution_0, # samples
    var_names, # var_names to plot
    labeller=model.labeller, # label manager
    kind="kde", # plot type
    reference_values=sim_params, # truths
)

['velocity', 'fwhm_nonthermal', 'fwhm', 'log10_N', 'log10_Tex_ul', 'log10_LTE_precision', 'tau_total']

# identify simulation cloud corresponding to each posterior cloud
sim_cloud_map = {}
for i in range(n_clouds):
    posterior_velocity = model.trace.solution_0['velocity'].sel(cloud=i).data.mean()
    match = np.argmin(np.abs(sim_params["velocity"] - posterior_velocity))
    sim_cloud_map[i] = match
sim_cloud_map

{0: np.int64(2), 1: np.int64(1), 2: np.int64(0)}

cloud = 0

# subset of sim_params
my_sim_params = {}
for var_name in model.cloud_deterministics:
    my_sim_params[var_name] = sim_params[var_name][sim_cloud_map[cloud]]

# ignore transition and state dependent parameters
var_names = [
    param for param in model.cloud_deterministics
    if not set(model.model.named_vars_to_dims[param]).intersection(set(["transition", "state"]))
    and param not in ["fwhm_thermal"]
]
print(var_names)
_ = plot_pair(
    model.trace.solution_0.sel(cloud=cloud), # samples
    var_names, # var_names to plot
    labeller=model.labeller, # label manager
    kind="kde", # plot type
    reference_values=my_sim_params, # truths
)

['velocity', 'fwhm_nonthermal', 'fwhm', 'log10_N', 'log10_Tex_ul', 'log10_LTE_precision', 'tau_total']

cloud = 1

# subset of sim_params
my_sim_params = {}
for var_name in model.cloud_deterministics:
    my_sim_params[var_name] = sim_params[var_name][sim_cloud_map[cloud]]

# ignore transition and state dependent parameters
var_names = [
    param for param in model.cloud_deterministics
    if not set(model.model.named_vars_to_dims[param]).intersection(set(["transition", "state"]))
    and param not in ["fwhm_thermal"]
]
print(var_names)
_ = plot_pair(
    model.trace.solution_0.sel(cloud=cloud), # samples
    var_names, # var_names to plot
    labeller=model.labeller, # label manager
    kind="kde", # plot type
    reference_values=my_sim_params, # truths
)

['velocity', 'fwhm_nonthermal', 'fwhm', 'log10_N', 'log10_Tex_ul', 'log10_LTE_precision', 'tau_total']

cloud = 2

# subset of sim_params
my_sim_params = {}
for var_name in model.cloud_deterministics:
    my_sim_params[var_name] = sim_params[var_name][sim_cloud_map[cloud]]

# ignore transition and state dependent parameters
var_names = [
    param for param in model.cloud_deterministics
    if not set(model.model.named_vars_to_dims[param]).intersection(set(["transition", "state"]))
    and param not in ["fwhm_thermal"]
]
print(var_names)
_ = plot_pair(
    model.trace.solution_0.sel(cloud=cloud), # samples
    var_names, # var_names to plot
    labeller=model.labeller, # label manager
    kind="kde", # plot type
    reference_values=my_sim_params, # truths
)

['velocity', 'fwhm_nonthermal', 'fwhm', 'log10_N', 'log10_Tex_ul', 'log10_LTE_precision', 'tau_total']

point_stats = az.summary(model.trace.solution_0, kind='stats', hdi_prob=0.68)
print("BIC:", model.bic())
display(point_stats)

BIC: -3273.0262330955006

	mean	sd	hdi_16%	hdi_84%
baseline_12CN-1_norm[0]	-0.076	0.044	-0.118	-0.030
baseline_12CN-2_norm[0]	-0.435	0.064	-0.499	-0.372
velocity_norm[0]	0.852	0.013	0.840	0.866
velocity_norm[1]	-0.002	0.005	-0.007	0.003
velocity_norm[2]	-0.663	0.003	-0.667	-0.660
log10_N_norm[0]	0.393	0.082	0.315	0.468
log10_N_norm[1]	0.363	0.065	0.298	0.425
log10_N_norm[2]	0.324	0.050	0.275	0.375
log10_Tex_ul_norm[0]	0.089	0.101	-0.022	0.155
log10_Tex_ul_norm[1]	1.128	0.221	0.897	1.310
log10_Tex_ul_norm[2]	1.546	0.249	1.294	1.763
fwhm_nonthermal_norm[0]	1.580	0.102	1.471	1.669
fwhm_nonthermal_norm[1]	1.312	0.058	1.258	1.372
fwhm_nonthermal_norm[2]	0.980	0.034	0.945	1.013
log10_LTE_precision_norm[0]	0.780	0.520	0.001	1.008
log10_LTE_precision_norm[1]	2.244	0.283	1.960	2.501
log10_LTE_precision_norm[2]	2.692	0.244	2.420	2.884
weights[0, 0 0 1 1]	0.214	0.003	0.211	0.217
weights[0, 0 0 1 2]	0.429	0.006	0.425	0.435
weights[0, 1 0 1 1]	0.040	0.001	0.038	0.041
weights[0, 1 0 1 2]	0.080	0.002	0.078	0.082
weights[0, 1 0 2 1]	0.040	0.001	0.038	0.041
weights[0, 1 0 2 2]	0.079	0.002	0.077	0.081
weights[0, 1 0 2 3]	0.119	0.003	0.116	0.121
weights[1, 0 0 1 1]	0.189	0.006	0.184	0.196
weights[1, 0 0 1 2]	0.368	0.010	0.360	0.379
weights[1, 1 0 1 1]	0.051	0.003	0.048	0.054
weights[1, 1 0 1 2]	0.101	0.005	0.096	0.105
weights[1, 1 0 2 1]	0.044	0.002	0.042	0.046
weights[1, 1 0 2 2]	0.097	0.004	0.093	0.101
weights[1, 1 0 2 3]	0.150	0.003	0.147	0.152
weights[2, 0 0 1 1]	0.178	0.005	0.174	0.184
weights[2, 0 0 1 2]	0.346	0.009	0.338	0.356
weights[2, 1 0 1 1]	0.054	0.002	0.051	0.056
weights[2, 1 0 1 2]	0.098	0.003	0.094	0.101
weights[2, 1 0 2 1]	0.045	0.002	0.043	0.047
weights[2, 1 0 2 2]	0.116	0.005	0.111	0.120
weights[2, 1 0 2 3]	0.163	0.002	0.160	0.165
velocity[0]	2.557	0.040	2.521	2.599
velocity[1]	-0.005	0.015	-0.020	0.010
velocity[2]	-1.990	0.010	-2.000	-1.981
fwhm_thermal[0]	0.075	0.000	0.075	0.075
fwhm_thermal[1]	0.075	0.000	0.075	0.075
fwhm_thermal[2]	0.075	0.000	0.075	0.075
fwhm_nonthermal[0]	1.580	0.102	1.471	1.669
fwhm_nonthermal[1]	1.312	0.058	1.258	1.372
fwhm_nonthermal[2]	0.980	0.034	0.945	1.013
fwhm[0]	1.581	0.102	1.473	1.670
fwhm[1]	1.314	0.058	1.260	1.374
fwhm[2]	0.983	0.034	0.948	1.016
log10_N[0]	13.893	0.082	13.815	13.968
log10_N[1]	13.863	0.065	13.798	13.925
log10_N[2]	13.824	0.050	13.775	13.875
log10_Tex_ul[0]	0.509	0.010	0.498	0.516
log10_Tex_ul[1]	0.613	0.022	0.590	0.631
log10_Tex_ul[2]	0.655	0.025	0.629	0.676
LTE_weights[0, 0 0 1 1]	0.214	0.003	0.212	0.217
LTE_weights[0, 0 0 1 2]	0.429	0.006	0.424	0.435
LTE_weights[0, 1 0 1 1]	0.040	0.001	0.039	0.041
LTE_weights[0, 1 0 1 2]	0.080	0.002	0.078	0.081
LTE_weights[0, 1 0 2 1]	0.040	0.001	0.039	0.040
LTE_weights[0, 1 0 2 2]	0.079	0.002	0.077	0.081
LTE_weights[0, 1 0 2 3]	0.119	0.003	0.116	0.121
LTE_weights[1, 0 0 1 1]	0.186	0.006	0.181	0.192
LTE_weights[1, 0 0 1 2]	0.372	0.011	0.362	0.383
LTE_weights[1, 1 0 1 1]	0.049	0.002	0.047	0.051
LTE_weights[1, 1 0 1 2]	0.099	0.004	0.095	0.102
LTE_weights[1, 1 0 2 1]	0.049	0.002	0.047	0.051
LTE_weights[1, 1 0 2 2]	0.098	0.004	0.094	0.101
LTE_weights[1, 1 0 2 3]	0.147	0.006	0.142	0.152
LTE_weights[2, 0 0 1 1]	0.176	0.006	0.171	0.181
LTE_weights[2, 0 0 1 2]	0.351	0.011	0.341	0.363
LTE_weights[2, 1 0 1 1]	0.053	0.002	0.051	0.054
LTE_weights[2, 1 0 1 2]	0.105	0.004	0.101	0.109
LTE_weights[2, 1 0 2 1]	0.052	0.002	0.051	0.054
LTE_weights[2, 1 0 2 2]	0.105	0.004	0.101	0.108
LTE_weights[2, 1 0 2 3]	0.157	0.006	0.152	0.163
log10_LTE_precision[0]	-5.220	0.520	-5.999	-4.992
log10_LTE_precision[1]	-3.756	0.283	-4.040	-3.499
log10_LTE_precision[2]	-3.308	0.244	-3.580	-3.116
Tex[113123.3687, 0]	3.221	0.085	3.138	3.284
Tex[113123.3687, 1]	4.181	0.266	3.849	4.335
Tex[113123.3687, 2]	4.555	0.281	4.240	4.752
Tex[113144.19, 0]	3.216	0.082	3.140	3.283
Tex[113144.19, 1]	4.273	0.266	3.981	4.462
Tex[113144.19, 2]	4.672	0.286	4.362	4.881
Tex[113170.535, 0]	3.244	0.087	3.154	3.295
Tex[113170.535, 1]	4.119	0.238	3.875	4.308
Tex[113170.535, 2]	4.219	0.206	3.996	4.371
Tex[113191.325, 0]	3.239	0.082	3.155	3.291
Tex[113191.325, 1]	4.208	0.241	3.940	4.379
Tex[113191.325, 2]	4.321	0.214	4.076	4.471
Tex[113488.142, 0]	3.221	0.082	3.146	3.289
Tex[113488.142, 1]	4.009	0.216	3.757	4.150
Tex[113488.142, 2]	4.851	0.300	4.523	5.058
Tex[113490.985, 0]	3.229	0.066	3.157	3.274
Tex[113490.985, 1]	4.197	0.146	4.024	4.284
Tex[113490.985, 2]	4.710	0.156	4.537	4.819
Tex[113499.643, 0]	3.235	0.087	3.147	3.292
Tex[113499.643, 1]	3.742	0.184	3.524	3.868
Tex[113499.643, 2]	3.974	0.177	3.789	4.117
Tex[113508.934, 0]	3.215	0.080	3.139	3.279
Tex[113508.934, 1]	4.094	0.222	3.835	4.241
Tex[113508.934, 2]	4.986	0.322	4.649	5.235
Tex[113520.4215, 0]	3.230	0.084	3.144	3.286
Tex[113520.4215, 1]	3.817	0.193	3.598	3.959
Tex[113520.4215, 2]	4.065	0.197	3.855	4.227
tau[113123.3687, 0]	0.049	0.010	0.038	0.057
tau[113123.3687, 1]	0.036	0.007	0.028	0.042
tau[113123.3687, 2]	0.030	0.005	0.025	0.034
tau[113144.19, 0]	0.405	0.086	0.316	0.471
tau[113144.19, 1]	0.286	0.058	0.223	0.335
tau[113144.19, 2]	0.234	0.040	0.191	0.269
tau[113170.535, 0]	0.393	0.083	0.311	0.462
tau[113170.535, 1]	0.292	0.059	0.226	0.340
tau[113170.535, 2]	0.247	0.040	0.204	0.282
tau[113191.325, 0]	0.512	0.108	0.404	0.600
tau[113191.325, 1]	0.365	0.073	0.285	0.428
tau[113191.325, 2]	0.308	0.050	0.254	0.352
tau[113488.142, 0]	0.514	0.109	0.405	0.603
tau[113488.142, 1]	0.386	0.077	0.299	0.448
tau[113488.142, 2]	0.300	0.051	0.245	0.346
tau[113490.985, 0]	1.367	0.287	1.072	1.592
tau[113490.985, 1]	0.972	0.185	0.780	1.140
tau[113490.985, 2]	0.784	0.123	0.649	0.890
tau[113499.643, 0]	0.405	0.086	0.320	0.476
tau[113499.643, 1]	0.314	0.061	0.243	0.362
tau[113499.643, 2]	0.261	0.041	0.217	0.298
tau[113508.934, 0]	0.397	0.084	0.312	0.464
tau[113508.934, 1]	0.286	0.057	0.227	0.337
tau[113508.934, 2]	0.221	0.038	0.178	0.253
tau[113520.4215, 0]	0.050	0.011	0.039	0.058
tau[113520.4215, 1]	0.037	0.007	0.029	0.043
tau[113520.4215, 2]	0.031	0.005	0.025	0.035
tau_total[0]	4.093	0.864	3.189	4.751
tau_total[1]	2.974	0.584	2.309	3.444
tau_total[2]	2.417	0.391	2.000	2.768
TR[113123.3687, 0]	1.236	0.068	1.170	1.285
TR[113123.3687, 1]	2.040	0.233	1.752	2.173
TR[113123.3687, 2]	2.369	0.251	2.090	2.543
TR[113144.19, 0]	1.231	0.065	1.170	1.283
TR[113144.19, 1]	2.120	0.235	1.865	2.285
TR[113144.19, 2]	2.473	0.257	2.196	2.659
TR[113170.535, 0]	1.253	0.070	1.182	1.294
TR[113170.535, 1]	1.985	0.207	1.749	2.123
TR[113170.535, 2]	2.072	0.181	1.859	2.186
TR[113191.325, 0]	1.249	0.066	1.182	1.290
TR[113191.325, 1]	2.062	0.211	1.829	2.211
TR[113191.325, 2]	2.160	0.189	1.939	2.285
TR[113488.142, 0]	1.231	0.065	1.172	1.285
TR[113488.142, 1]	1.886	0.186	1.669	2.006
TR[113488.142, 2]	2.629	0.272	2.322	2.803
TR[113490.985, 0]	1.238	0.053	1.181	1.273
TR[113490.985, 1]	2.047	0.128	1.897	2.123
TR[113490.985, 2]	2.501	0.140	2.346	2.598
TR[113499.643, 0]	1.243	0.070	1.173	1.287
TR[113499.643, 1]	1.659	0.156	1.476	1.763
TR[113499.643, 2]	1.855	0.152	1.696	1.977
TR[113508.934, 0]	1.227	0.063	1.166	1.277
TR[113508.934, 1]	1.959	0.193	1.736	2.085
TR[113508.934, 2]	2.750	0.293	2.403	2.933
TR[113520.4215, 0]	1.238	0.067	1.170	1.282
TR[113520.4215, 1]	1.721	0.164	1.536	1.841
TR[113520.4215, 2]	1.933	0.171	1.752	2.072