Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 3 additions & 3 deletions flip/fitter.py
Original file line number Diff line number Diff line change
Expand Up @@ -163,7 +163,7 @@ def init_from_covariance(
data,
parameter_dict,
likelihood_type="multivariate_gaussian",
likelihood_properties=None,
likelihood_properties={},
**kwargs,
):
"""
Expand Down Expand Up @@ -194,7 +194,7 @@ def init_from_covariance(
likelihood = minuit_fitter.get_likelihood(
parameter_dict,
likelihood_type=likelihood_type,
likelihood_properties=likelihood_properties,
likelihood_properties={**likelihood_properties, 'negative_log_likelihood': True},
**kwargs,
)
minuit_fitter.likelihood = likelihood
Expand Down Expand Up @@ -247,7 +247,7 @@ def init_from_file(
data,
parameter_dict,
likelihood_type=likelihood_type,
likelihood_properties=likelihood_properties,
likelihood_properties={**likelihood_properties, 'negative_log_likelihood': True},
)

def setup_minuit(self, parameter_dict):
Expand Down
54 changes: 35 additions & 19 deletions flip/likelihood.py
Original file line number Diff line number Diff line change
Expand Up @@ -11,49 +11,48 @@ def log_likelihood_gaussian_inverse(vector, covariance_sum):
_, logdet = np.linalg.slogdet(covariance_sum)
inverse_covariance_sum = np.linalg.inv(covariance_sum)
chi2 = np.dot(vector, np.dot(inverse_covariance_sum, vector))
return 0.5 * (vector.size * np.log(2 * np.pi) + logdet + chi2)
return -0.5 * (vector.size * np.log(2 * np.pi) + logdet + chi2)


def log_likelihood_gaussian_cholesky(vector, covariance_sum):
cholesky = sc.linalg.cho_factor(covariance_sum)
logdet = 2 * np.sum(np.log(np.diag(cholesky[0])))
chi2 = np.dot(vector, sc.linalg.cho_solve(cholesky, vector))
return 0.5 * (vector.size * np.log(2 * np.pi) + logdet + chi2)
return -0.5 * (vector.size * np.log(2 * np.pi) + logdet + chi2)


class BaseLikelihood(object):

_default_likelihood_properties = {
"inversion_method": "inverse",
"velocity_type": "direct",
"velocity_estimator": "full",
"negative_log_likelihood": True,
}

def __init__(
self,
covariance=None,
data=None,
parameter_names=None,
likelihood_properties=None,
likelihood_properties={},
):
self.covariance = covariance
self.data = data
self.parameter_names = parameter_names

_default_likelihood_properties = {
"inversion_method": "inverse",
"velocity_type": "direct",
"velocity_estimator": "full",
self.likelihood_properties = {
**self._default_likelihood_properties,
**likelihood_properties,
}
if likelihood_properties == None:
likelihood_properties = _default_likelihood_properties
else:
for key in _default_likelihood_properties.keys():
if key not in likelihood_properties.keys():
likelihood_properties[key] = _default_likelihood_properties[key]

self.likelihood_properties = likelihood_properties

@classmethod
def init_from_covariance(
cls,
covariance,
data,
parameter_names,
likelihood_properties=None,
likelihood_properties={},
**kwargs,
):
"""
Expand Down Expand Up @@ -121,7 +120,7 @@ def __init__(
covariance=None,
data=None,
parameter_names=None,
likelihood_properties=None,
likelihood_properties={},
):
super(MultivariateGaussianLikelihood, self).__init__(
covariance=covariance,
Expand All @@ -147,6 +146,8 @@ def __call__(self, parameter_values):
likelihood_function = eval(
f"log_likelihood_gaussian_{self.likelihood_properties['inversion_method']}"
)
if self.likelihood_properties["negative_log_likelihood"]:
return -likelihood_function(vector, covariance_sum)
return likelihood_function(vector, covariance_sum)


Expand All @@ -156,7 +157,7 @@ def __init__(
covariance=None,
data=None,
parameter_names=None,
likelihood_properties=None,
likelihood_properties={},
interpolation_value_name=None,
interpolation_value_range=None,
):
Expand Down Expand Up @@ -189,6 +190,17 @@ def __init__(
self.interpolation_value_range = interpolation_value_range

def verify_covariance(self):
"""
The verify_covariance function is used to ensure that the covariance matrix of each
parameter in the model has been computed. If it has not, then this function will compute
it and store it as a full matrix.

Args:
self: Bind the method to the object

Returns:
Nothing
"""
for i in range(len(self.covariance)):
if self.covariance[i].full_matrix is False:
self.covariance[i].compute_full_matrix()
Expand Down Expand Up @@ -232,6 +244,8 @@ def __call__(self, parameter_values):
likelihood_function = eval(
f"log_likelihood_gaussian_{self.likelihood_properties['inversion_method']}"
)
if self.likelihood_properties["negative_log_likelihood"]:
return -likelihood_function(vector, covariance_sum)
return likelihood_function(vector, covariance_sum)


Expand All @@ -241,7 +255,7 @@ def __init__(
covariance=None,
data=None,
parameter_names=None,
likelihood_properties=None,
likelihood_properties={},
interpolation_value_range_0=None,
interpolation_value_range_1=None,
):
Expand Down Expand Up @@ -331,4 +345,6 @@ def __call__(
likelihood_function = eval(
f"log_likelihood_gaussian_{self.likelihood_properties['inversion_method']}"
)
if self.likelihood_properties["negative_log_likelihood"]:
return -likelihood_function(vector, covariance_sum)
return likelihood_function(vector, covariance_sum)
7 changes: 4 additions & 3 deletions flip/vectors.py
Original file line number Diff line number Diff line change
Expand Up @@ -77,6 +77,9 @@ def load_velocity_vector(
f"""Please choose a velocity_type among {_avail_velocity_type}"""
)

if "vmean" in parameter_values_dict:
velocity = velocity - parameter_values_dict["vmean"]

return velocity, velocity_error


Expand Down Expand Up @@ -125,9 +128,7 @@ def compute_observed_distance_modulus(
mu = data["mb"] + alpha * data["x1"] - beta * data["c"] - M0

variance_mu = (
data["e_mb"] ** 2
+ alpha**2 * data["e_x1"] ** 2
+ beta**2 * data["e_c"] ** 2
data["e_mb"] ** 2 + alpha**2 * data["e_x1"] ** 2 + beta**2 * data["e_c"] ** 2
)
variance_mu += (
2 * alpha * data["cov_mb_x1"]
Expand Down
Loading