Fancy jacobian

[TTitscher]

• Try to entangle different model error definitions:
  • ModelErrorInterface as an interface for the InferenceProblem
  • VariationalBayesModelError as an interface for VB
  • see the __call__ documentation for the differences
• Adapt VB algorithm to a dict of noise groups:
  • the VBModelError needs to return a dict {noise_key : numbers}
  • its jacobian must return dict {noise_key : matrix}
  • the noise prior must be dict {noise_key : Gamma}
  • A VBAdapter lets you use numbers, matrix, Gamma (for error, jacobian, noise prior respectively) for the simple cases of only one noise group.
• provide jacobian.py to collect algorithms to calculate the Jacobian via central differences for
  • single named scalar parameter
  • single named vector parameter
  • unnamed normal list of numbers

[ajafarihub]

I sent the past comment through my other account (connected to the university). Sorry about it. Here you have the original one :).

Very nice! I see a lot of progress made, which I would need longer time to keep track of.

Maybe just a general remark at this point, which is unclear to my understanding. I would argue why we need to take care of the Jacobian according to the noise groups. This sounds somewhat unnecessary IMO. The reason is that, Jacobian is - by definition - related to the forward model and NOT to the model error which accounts for noise groups. I mean, the computation of Jacobian is normally performed regardless of having/missing any noise, and this is what a user could have done before defining any model-error. So, when it is already computed, why do we need to again somehow factorize and re-construct it based on noise groups ?

[TTitscher]

Jacobian is - by definition - related to the forward model and NOT to the model error

If the Jacobian is defined in the forward model, you can just call it in model_error.jacobian. But I already encountered an exception to your statement, namely when I defined the model error as

def __call__(...):
   model_response = self.model(...)
   data = self.database(...)
   return model_response - parameters["offset"] - data * parameters["factor"]

That means, that the derivative w.r.t. parameters["factor"] includes data. So I would not restrict the user too much here.

why do we need to again somehow factorize and re-construct it based on noise groups ?

I am not sure if I get your point, this "factorize" confuses me ;) The model_error outputs its values sorted by keys, e.g. a "sensor_key". Thus, its Jacobian is a nested dict where the first key is also the "sensor_key", the second one is the parameter of the derivative. The noise stuff now collects (possibly) multiple "sensor_key" of different model_errors to a noise group. The corresponding jacobians are not split, just combined.

[ajafarihub]

OK, in cases like the example you made, we can have different model-error's Jacobian than forward-model's. I never thought of such cases so far.

What I meant by factorize was: to extract a particular segment (block) of the full Jacobian given a noise key. You most likely have a good reason/usefulness for such a possibility, but I was wondering why we would need it. So, hopefully I will realize it later on ;).

Maybe a good point to give another general remark (from a user's viewpoint and not only related to this issue): The interface has already become well facilitated with different potentials (like this nested dictionary of Jacobian concerning different noise groups, or considering multi-model-errors with shared latent parameters, or ... ). What I find very relevant in this regard is to design/develop all such options so that the default bahavior of the interface will not become too complicated at the user-level. For example, if a user has only one model-error with few parameters and one noise-group (practically, no input for it), the interface would not require him/her to deal with various unnecessary structures/definitions which are only relevant in more complicated situations.

[TTitscher]

to extract a particular segment (block) of the full Jacobian given a noise key. You most likely have a good reason/usefulness for such a possibility, but I was wondering why we would need it.

Those questions are a good guideline on what needs to be documented or shown in a tutorial later on.

The segments of the jacobian can be defined as

def jacobian(...):
    jac["OOB_forces"]["disp_field"] = Stiffness(...)
    jac["OOB_forces"]["prm2"] = ...

So they exist in this nested dict-structure. Right before "vb" is called, multiple Jacobians can be combined to a possibly large matrix as the VB input. So a particular segment or block is never extracted from a huge numpy-array, if that was your concern.