pysolvegn.study_jacobian#

study_jacobian(residual_func, jacobian_func, parameters, *, weight=None, loss=None, title='')[source]#

Study the Jacobian matrix of the least squares problem to analyze the observability of the parameters and the convergence properties of the optimization (called by verbosity level 3 of the Gauss-Newton solver).

The function accepts multiple terms in the least squares problem, each with its own residual function, Jacobian function, weight, and loss function studying the following optimization problem

\[\min_{\mathbf{p}} \frac{1}{2} \sum_{i} w_i \sum_j \rho_i\left(\| \mathbf{R}_{i,j}(\mathbf{p}) \|^2\right)\]

Then the function displays the following information about the Jacobian matrix of the least squares problem:

  • Density, Cost contribution and loss function for each term in the least squares problem.

  • Singular values of the combined Jacobian matrix \(\sum_i w_i \mathbf{J}_i^T \mathbf{J}_i\) and their contribution to the variance of the parameters.

  • Condition number of the combined Jacobian matrix, which is an indicator of the convergence properties of the optimization.

  • Singular vectors of the combined Jacobian matrix, which can be used to analyze the observability of the parameters and the sensitivity of the optimization to noise in the data.

  • Parameter sensitivity analysis based on the covariance matrix of the parameters, which is computed as :

\[\Sigma = \sigma^2 \left( \sum_i w_i \mathbf{J}_i^T \mathbf{J}_i \right)^{-1}\]

where \(\sigma^2\) is the estimated variance of the residuals, computed as \(\sigma^2 = \frac{2C_0}{N_{eq}-N_p}\) with \(C_0\) the cost of the first term (assuming data), \(N_{eq}\) the number of equations and \(N_p\) the number of parameters.

Parameters:

  • residual_func (Union[Callable, Sequence[Callable]]) – The function(s) to compute the residuals for each term in the least squares problem. Can be a single function or a sequence of functions. Each callable represent the i-th term \(\mathbf{R}_i(\mathbf{p})\) in the least squares problem and should take as inputs the parameters (array-like with shape (n_parameters,)) and return the residuals as a 1D numpy array with shape (n_residuals_i,) representing each \(\mathbf{R}_{i,j}(\mathbf{p})\) where \(j=0,...,n_{\text{residuals}_i-1}\).

  • jacobian_func (Union[Callable, Sequence[Callable]]) – The function(s) to compute the Jacobian matrices for each term in the least squares problem. Can be a single function or a sequence of functions. Each callable represent the i-th term \(\mathbf{J}_i(\mathbf{p})\) in the least squares problem such as \(\mathbf{J}_i(\mathbf{p}) = \frac{\partial \mathbf{R}_i(\mathbf{p})}{\partial \mathbf{p}}\) and should take as inputs the parameters (array-like with shape (n_parameters,)) and return the Jacobian matrix as a 2D numpy array with shape (n_residuals_i, n_parameters) representing the Jacobian of each residual \(\mathbf{R}_{i,j}(\mathbf{p})\) where \(j=0,...,n_{\text{residuals}_i-1}\).

  • parameters (ArrayLike) – The parameters of the optimization step with shape (n_parameters,).

  • weight (Optional[Union[Real, Sequence[Real]]], optional (default=None)) – The weight(s) for each term in the least squares problem. Can be a single value or a sequence of values. Each weight \(w_i\) is a weight for the corresponding term in the least squares problem and should be a non-negative real number. If a single value is provided, it will be used for all terms in the least squares problem. If None is provided, it will be set to 1.0.

  • loss (Optional[Union[str, Sequence[str]]], optional (default=None)) – The loss function(s) for each term in the least squares problem. Can be a single string or a sequence of strings. Each loss function \(\rho_i\) is a robust cost function for the corresponding term in the least squares problem and should be one of the following strings: [“linear”, “soft_l1”, “cauchy”, “arctan”]. If a single string is provided, it will be used for all terms in the least squares problem. If None is provided, it will be set to “linear”.

  • title (str, optional) – An optional title for the analysis, which can be used for logging purposes. Default is an empty string.

Returns:

This function is used for analysis and does not return any value.

Return type:

None

Version#

  • 0.0.1: Initial version.

  • 0.0.2: \(\sigma^2\) is now estimated as \(\sigma^2 = 2C_0/(N_{eq}-N_p)\) where \(C_0\) is the cost of the first term (assuming data) instead of the total cost, to provide a more accurate estimation of the residual variance for the parameter sensitivity analysis.

On this page

This Page