diff --git a/tests/v1/test_priors.py b/tests/v1/test_priors.py index ea47e54f..37fd6590 100644 --- a/tests/v1/test_priors.py +++ b/tests/v1/test_priors.py @@ -8,7 +8,15 @@ from scipy.stats import norm import petab.v1 -from petab.v1 import get_simulation_conditions +from petab.v1 import ( + ESTIMATE, + MEASUREMENT, + OBJECTIVE_PRIOR_TYPE, + OBSERVABLE_ID, + SIMULATION, + get_simulation_conditions, + get_simulation_df, +) from petab.v1.priors import priors_to_measurements @@ -17,20 +25,29 @@ ) def test_priors_to_measurements(problem_id): """Test the conversion of priors to measurements.""" + # setup petab_problem_priors: petab.v1.Problem = ( benchmark_models_petab.get_problem(problem_id) ) petab_problem_priors.visualization_df = None assert petab.v1.lint_problem(petab_problem_priors) is False - if problem_id == "Isensee_JCB2018": # required to match the stored simulation results below petab.v1.flatten_timepoint_specific_output_overrides( petab_problem_priors ) assert petab.v1.lint_problem(petab_problem_priors) is False + original_problem = deepcopy(petab_problem_priors) + x_scaled_dict = dict( + zip( + original_problem.x_free_ids, + original_problem.x_nominal_free_scaled, + strict=True, + ) + ) + # convert priors to measurements petab_problem_measurements = priors_to_measurements(petab_problem_priors) # check that the original problem is not modified @@ -45,6 +62,7 @@ def test_priors_to_measurements(problem_id): getattr(original_problem, attr) ) ).empty, diff + # check that measurements and observables were added assert petab.v1.lint_problem(petab_problem_measurements) is False assert ( @@ -59,6 +77,7 @@ def test_priors_to_measurements(problem_id): petab_problem_measurements.measurement_df.shape[0] > petab_problem_priors.measurement_df.shape[0] ) + # ensure we didn't introduce any new conditions assert len( get_simulation_conditions(petab_problem_measurements.measurement_df) @@ -67,26 +86,40 @@ def test_priors_to_measurements(problem_id): # verify that the objective function value is the same # load/construct the simulation results - simulation_df_priors = petab.v1.get_simulation_df( + simulation_df_priors = get_simulation_df( Path( benchmark_models_petab.MODELS_DIR, problem_id, f"simulatedData_{problem_id}.tsv", ) ) - simulation_df_measurements = pd.concat( - [ - petab_problem_measurements.measurement_df.rename( - columns={petab.v1.MEASUREMENT: petab.v1.SIMULATION} - )[ - petab_problem_measurements.measurement_df[ - petab.v1.C.OBSERVABLE_ID - ].str.startswith("prior_") - ], - simulation_df_priors, + # for the prior observables, we need to "simulate" the model with the + # nominal parameter values + simulated_prior_observables = ( + petab_problem_measurements.measurement_df.rename( + columns={MEASUREMENT: SIMULATION} + )[ + petab_problem_measurements.measurement_df[ + OBSERVABLE_ID + ].str.startswith("prior_") ] ) + def apply_parameter_values(row): + # apply the parameter values to the observable formula for the prior + if row[OBSERVABLE_ID].startswith("prior_"): + row[SIMULATION] = x_scaled_dict[ + row[OBSERVABLE_ID].removeprefix("prior_") + ] + return row + + simulated_prior_observables = simulated_prior_observables.apply( + apply_parameter_values, axis=1 + ) + simulation_df_measurements = pd.concat( + [simulation_df_priors, simulated_prior_observables] + ) + llh_priors = petab.v1.calculate_llh_for_table( petab_problem_priors.measurement_df, simulation_df_priors, @@ -102,10 +135,8 @@ def test_priors_to_measurements(problem_id): # get prior objective function contribution parameter_ids = petab_problem_priors.parameter_df.index.values[ - (petab_problem_priors.parameter_df[petab.v1.ESTIMATE] == 1) - & petab_problem_priors.parameter_df[ - petab.v1.OBJECTIVE_PRIOR_TYPE - ].notna() + (petab_problem_priors.parameter_df[ESTIMATE] == 1) + & petab_problem_priors.parameter_df[OBJECTIVE_PRIOR_TYPE].notna() ] priors = petab.v1.get_priors_from_df( petab_problem_priors.parameter_df, @@ -134,4 +165,4 @@ def test_priors_to_measurements(problem_id): llh_priors + prior_contrib, llh_measurements, rtol=1e-3, atol=1e-16 ), (llh_priors + prior_contrib, llh_measurements) # check that the tolerance is not too high - assert np.abs(prior_contrib) > 1e-3 * np.abs(llh_priors) + assert np.abs(prior_contrib) > 1e-8 * np.abs(llh_priors)