Improve sysid robustness, prompted by #3284#3287
Conversation
|
I agree with @aftersomemath , the residual scaling is definitely a |
Per-parameter scaling via change of variables z = x / D. Supports 'jac' (adaptive D_i = 1/||J(:,i)|| per iteration, matches scipy's TRF), explicit array, or a positive scalar. Default 1.0 is a no-op.
kevinzakka
force-pushed
the
sysid-issue-3284-fixes
branch
from
4878167 to
1327654
Compare
| rng = hi - lo | ||
| safe_rng = np.where(rng > 0, rng, 1.0) | ||
| v = self.as_vector().copy() | ||
| at_lo = (v - lo) <= 1e-3 * safe_rng |
There was a problem hiding this comment.
1e-3 should probably be a parameter as well.
| assert (result_dir / "arm.xml").exists() | ||
|
|
||
|
|
||
| def _rank_1_residual_fn(x, p): |
There was a problem hiding this comment.
If I understand correctly, this residual function has a rank 1 Jacobian, but only because its output is scalar. A better test would be a function with 2 outputs, that depends on only one variable. Then the jacobian will be 2x2 and rank 1.
| x0: np.ndarray, | ||
| residual_fn: Callable[..., Any], | ||
| bounds: tuple[np.ndarray, np.ndarray], | ||
| x_scale: XScale = 1.0, |
There was a problem hiding this comment.
Not needed if we are relying on **kwargs as the doc string suggests.
| verbose: If True, log parameter comparison table after optimization. | ||
| **optimizer_kwargs: Forwarded to the backend (e.g. ``max_iters``, | ||
| ``verbose``, ``loss``). | ||
| check_conditioning: If True, estimate ``cond(JᵀJ)`` at the starting |
There was a problem hiding this comment.
Is it ok to use unicode like supscript T?
| lo, hi = bounds | ||
| rng = hi - lo | ||
| safe_rng = np.where(rng > 0, rng, 1.0) | ||
| at_bound = ((x0 - lo) <= 1e-3 * safe_rng) | ((hi - x0) <= 1e-3 * safe_rng) |
There was a problem hiding this comment.
We should probably check that x0 is within the bounds before executing this check
| * ``loss``: scipy loss function name (scipy backends only). | ||
| * ``x_scale``: per-parameter scaling. ``"jac"`` is adaptive | ||
| ``D_i = 1/||J(:,i)||`` per iteration; an explicit array or | ||
| positive scalar is used as ``D`` directly. Defaults: ``"jac"`` |
There was a problem hiding this comment.
Scipy does not simply default to "jac": https://docs.scipy.org/doc/scipy/reference/generated/scipy.optimize.least_squares.html#id2
| optimizer="scipy", verbose=False, check_conditioning=True, | ||
| max_iters=1, | ||
| ) | ||
| fired = any("cond(JᵀJ)" in r.message for r in caplog.records) |
| residual_fn: Callable[..., Any], | ||
| threshold: float = 1e12, | ||
| ) -> None: | ||
| """Warn if cond(JᵀJ) at the starting point exceeds ``threshold``. |
| residuals, _, _ = residual_fn(x, initial_params) | ||
| return np.concatenate(residuals) | ||
|
|
||
| eps = np.finfo(np.float64).eps ** 0.5 |
There was a problem hiding this comment.
eps should be the same as what will be passed to mujoco/scipy. This is their default, but the user could override.
| ) | ||
| # eigvalsh on the gram matrix so rank-deficient directions are visible | ||
| # when n_params > n_residual components (SVD would drop to min(m, n)). | ||
| ev = np.maximum(np.linalg.eigvalsh(jac.T @ jac), 0.0) |
There was a problem hiding this comment.
Is this better than using Numpy's condition number function? https://numpy.org/doc/stable/reference/generated/numpy.linalg.cond.html
ParameterDict.move_off_bounds() shifts at-bound values into the interior; optimize() warns when any non-frozen component starts at a bound.
mujoco.minimize.least_squares now supports x_scale natively, so the sysid wrapper becomes a thin pass-through. 'jac' is adaptive per iteration (was static at x0).
Pass check_conditioning=True to FD the Jacobian at the starting point and warn if cond(J^T J) suggests numerical ill-conditioning.
kevinzakka
force-pushed
the
sysid-issue-3284-fixes
branch
from
1327654 to
857b55f
Compare
|
LGTM. Just need to answer the outstanding question about mu adaptation in #3294. |
3 participants
This PR depends on #3294. It adds 3 robustness aids to mujoco.sysid to hopefully alleviate issues like #3284. Specifically, it adds:
ParameterDict.move_off_bounds(): this function nudges any parameter starting at a box bound into the interior. Additionally,optimizewill now warn you when it detects that a parameter is at a bound.x_scaleis now accepted as a kwarg in the mujoco solver. It rescales parameters with widely different magnitudes, a common source of ill-conditioning in sysid. The scipy backend already supported this.check_conditinong=Trueis now an arg inoptimizeand it will warn you when the condition number of the Hessian exceeds 1e12, a heuristic we chose well below the float64 conditioning limit.