Dynamic parameter identification: learn a joint model from motion

Fit dynamic parameters from a rotating joint's motion and torque data. Recover inertia, gravity mass moment, and damping, then test excitation, noise, and held-out predictions.

By 14 min read

What you will learn

  • Turn a one-joint torque equation into a three-column regression problem.
  • Distinguish a unique parameter fit from an accurate fit.
  • Explain why repeated holds reveal no inertia or damping information.
  • Compare training error with predictions on a separate motion.
  • Recognize the limits of lumped parameters and unconstrained estimates.

Before you start

A joint can fit its recorded torques closely and still have a badly wrong inertia estimate. In this lesson's noisy example, slow motion produces negative fitted inertia despite a small training error. A separate motion exposes the error.

Dynamic parameter identification estimates physical model coefficients from measurements. Here you will fit three coefficients, check whether the motion reveals them, and test the resulting torque predictions.

Learn which parameters explain the torque

Inverse dynamics takes known model parameters and calculates torque. Identification begins with the motion and measured torque, then estimates the unknown coefficients.

For one rotating joint, ask three questions:

  • How much torque accompanies angular acceleration?
  • How much torque balances gravity at this angle?
  • How much torque offsets friction at this speed?

Those questions lead to a small supervised learning problem. Each motion sample provides input features; measured actuator torque provides the target. MIT's system identification notes describe this use of the known mechanical equation structure.

Define one rotating joint

A rigid body rotates about a fixed horizontal joint axis. Its center of mass moves in a vertical plane. Angle q measures counterclockwise from world +x, so q = 0 points horizontally right.

Assume constant effective joint inertia I, a center of mass on the body's positive reference direction, and viscous joint friction. Gravity points along world −y with magnitude g = 9.81 m/s². The required actuator torque is:

τ = I q̈ + h g cos q + B q̇

CoefficientMeaningUnits
IEffective rotational inertia about the jointkg m²
h = mrMass times distance from joint to center of masskg m
BViscous damping coefficientN m s/rad

Joint velocity q̇ uses rad/s; acceleration q̈ uses rad/s². Radians are dimensionless in the SI torque calculation, so every term has units N m. Positive actuator torque acts toward increasing q.

Physical gravity applies torque −h g cos q, and physical damping applies −B q̇. Their positive terms above describe compensation. A horizontal stationary hold therefore needs τ = h g.

The synthetic joint uses I = 0.8, h = 0.4, and B = 0.15 in the units above. The model omits contact, joint springs, Coulomb friction, motor electrical dynamics, and actuator limits. The robot dynamics lesson places these contributions in a multijoint equation.

Build the dynamic parameter regressor

Collect q, q̇, q̈, and actuator torque at the same instant. The three known features form one row:

xᵢ = [q̈ᵢ, g cos qᵢ, q̇ᵢ]
θ = [I, h, B]ᵀ
τᵢ = xᵢ θ

Stack N rows into an N × 3 matrix X. Stack measured torques into y. The measurement model becomes y = Xθ + ε, with ε representing torque error.

Cosine makes the equation nonlinear in angle. Once the angles are known, each row contains ordinary numbers, and the equation is linear in the three unknown coefficients. This is linear regression with physically chosen features and no added intercept.

A torque sensor bias would violate that no-intercept model. Estimating an extra bias requires another coefficient and enough information to distinguish it from the gravity term. Adding columns changes the identification problem.

Recover three parameters by hand

Three carefully chosen noiseless instants can isolate the coefficients:

Instantqq̇q̈Measured τ
Horizontal hold0003.924 N m
Vertical accelerationπ/202 rad/s²1.600 N m
Vertical motionπ/2−3 rad/s0−0.450 N m

The first row gives h = 3.924 / 9.81 = 0.4 kg m. At the vertical angle, cos q = 0, so the second gives I = 1.6 / 2 = 0.8 kg m². The third gives B = −0.45 / −3 = 0.15 N m s/rad.

These are separate prescribed instants, potentially from separate runs. Zero velocity at an instant can coexist with nonzero acceleration. A real data collection procedure must reach each state and respect the machine's limits.

The corresponding rows are (0, 9.81, 0), (2, 0, 0), and (0, 0, −3). Each row adds an independent equation. Their matrix has rank three, so these ideal observations identify all three coefficients.

Choose motion that reveals the parameters

If every sample is a horizontal hold, every regressor row equals (0, 9.81, 0). More repetitions can help estimate the average holding torque, but they add no information about I or B. The matrix has rank one regardless of the number of repeated holds.

The moving presets use samples from smooth, prescribed trajectories over 0 ≤ t ≤ 8 s:

q_rich(t) = 0.7 sin(1.3t) + 0.25 sin(2.7t)
q_slow(t) = 0.08 sin(0.4t)

The lab differentiates these formulas analytically to obtain velocity and acceleration. Each set describes a realizable kinematic trajectory under an ideal actuator that supplies the required torque. The samples do not independently invent unrelated positions and derivatives.

For 24 samples, rich motion has acceleration RMS 1.430761 rad/s². Slow motion has acceleration RMS 0.008781 rad/s². Its inertia torque is therefore tiny compared with the same added torque noise.

Both moving data sets have full rank. Full rank establishes a unique least-squares answer; it does not establish enough signal to estimate each coefficient accurately. Excitation must reveal parameter effects above the measurement errors you face.

Fit the coefficients and check rank

With noisy measurements, choose the parameter vector that minimizes squared torque residuals:

θ̂ = arg min θ Σᵢ (yᵢ − xᵢθ)² = arg min θ ‖y − Xθ‖²

The browser solves this with column-scaled, column-pivoted Householder QR. It scales each nonzero column to unit Euclidean norm, factors the scaled matrix, and converts the fitted coefficients back to their original units. Pivoting changes column order during the calculation and restores the parameter order afterward.

QR uses orthogonal transformations and a triangular solve. It avoids explicitly forming XᵀX. Driscoll and Braun's numerical computation text explains Householder reflections and their use in least squares.

The lab treats a remaining column norm of at most 10⁻¹⁰, after the unit-column scaling, as numerically dependent. A rank below three produces no unique parameter fit. The lab leaves parameter and prediction readouts unavailable in that case, even though a smaller subset, such as h during holds, may remain identifiable.

Scaling helps the arithmetic handle different column units and magnitudes. It cannot create missing measurements or improve the physical signal-to-noise ratio. A tiny acceleration still amplifies torque errors when the fitted inertia returns to kg m².

Compare three identification experiments

Identify a joint from its torque data

Does the motion reveal the parameters?

Fit inertia, gravity mass moment, and viscous damping. Then predict torque on a separate motion that the fit never sees.

0.00
24
7

Known synthetic parameters: I = 0.8 kg m², h = 0.4 kg m, B = 0.15 N m s/rad. Gravity is 9.81 m/s². Angles start from horizontal +x and increase counterclockwise. All angles and derivatives are exact; only training torque receives bounded, repeatable noise.

Fit the training torques

Fit the training torques: 24 samples over eight secondsBlue dots show synthetic measured training torque, including any selected noise. The orange line joins fitted torque predictions at the sample times. Horizontal axis is time in seconds; vertical axis is torque in newton meters.048-606Time (s)Torque (N m)
Blue dots: measured training torques. Orange line: fitted predictions at those states, connected as a visual guide. Each chart chooses its own torque scale.

Predict a different motion

Predict a different motion: 16 samples over eight secondsBlue dots show exact synthetic torque on a separate validation trajectory. The orange line joins fitted torque predictions at the sample times. Horizontal axis is time in seconds; vertical axis is torque in newton meters.048-606Time (s)Torque (N m)
Blue dots: noiseless truth on 16 held-out states. Orange line: fitted predictions at those states, connected as a visual guide. Each chart chooses its own torque scale.
Regressor rank
3 / 3
Fitted inertia (kg m²)
0.800000
Fitted mass moment (kg m)
0.400000
Fitted damping (N m s/rad)
0.150000
Training RMSE (N m)
0.000000
Held-out RMSE (N m)
0.000000
Acceleration RMS (rad/s²)
1.430761
Speed RMS (rad/s)
0.816288

Rank 3 of 3: a unique least-squares fit. Held-out torque RMSE is 0.000000 N m. Positive parameter signs alone do not establish an accurate physical model.

Training RMSE compares predictions with noisy measurements. Held-out RMSE compares predictions with known noiseless truth on a different analytic trajectory. A low torque error does not test a long motion simulation or establish safety for a real robot.

Try Slow motion with noise 0.20. Its small speeds and accelerations carry little information about damping and inertia. Full numerical rank still allows large parameter errors.

Controls become available when this experiment finishes loading. The initial diagrams and parameter fit remain readable without JavaScript.

Start with Rich motion, zero noise, 24 samples, and seed 7. The fit recovers (0.8, 0.4, 0.15), with torque errors that round to zero. The experiment uses exact state measurements and the correct model structure.

  • Set Torque noise bound to 0.20 N m. Rich motion gives approximately I = 0.793312, h = 0.401350, and B = 0.153888.
  • Keep those controls and select Slow motion. The fitted inertia becomes −1.607788 kg m², which fails this joint model's physical requirements.
  • Select Horizontal holds. Rank falls to one, and the lab stops displaying a unique three-coefficient fit.

The noise control bounds each additive training torque error between its negative and positive value. A fixed seeded sequence makes the results repeatable. Changing the seed changes that sequence; it does not create a new physical joint.

Training samples span the same eight-second interval as their count changes. More samples alter the sampled rows and the noise realization, so a particular error need not decrease at every slider step. Reset restores the original four controls.

Test predictions on a different motion

The fit never sees the 16 validation states. They follow another prescribed trajectory:

q_test(t) = 0.6 sin(0.9t) + 0.3 cos(2.1t)

The lab uses the known synthetic coefficients to calculate noiseless validation torque. It then predicts torque at the same validation states with θ̂. Root mean squared error, or RMSE, is the square root of the average squared difference.

Training motion, with noise 0.20 and seed 7Training RMSEHeld-out RMSE
Rich motion, 24 samples0.119280 N m0.014150 N m
Slow motion, 24 samples0.118472 N m2.261865 N m

Slow motion has a slightly smaller training error and a much larger held-out error. Its training samples barely test inertia, so the fit can absorb noise into a large inertia mistake. The separate trajectory calls on that coefficient more strongly.

The two errors use different reference targets. Training RMSE includes torque measurement noise; held-out RMSE compares against known noiseless truth. A physical validation run would bring its own sensor errors and model mismatch.

This validation checks torque prediction at supplied states. It does not integrate a long simulated motion or test a feedback controller. Forward dynamics is a separate step where torque errors can accumulate into motion errors.

Interpret the fitted parameters carefully

An effective inertia about the joint can include a body's pivot inertia and reflected rotor inertia. Even with rotor inertia absent, I = I_COM + mr² and h = mr do not separately reveal m, r, and I_COM.

For example, both of these parameter choices give I = 0.8 and h = 0.4:

  • m = 1 kg, r = 0.4 m, I_COM = 0.64 kg m².
  • m = 2 kg, r = 0.2 m, I_COM = 0.72 kg m².

Each satisfies the parallel-axis relation. Geometry or an independent mass measurement could rule out some choices. The three-column torque experiment alone does not choose between them.

This joint model requires I > 0 and B ≥ 0. Its chosen center-of-mass reference direction also gives h ≥ 0. The lab keeps an unconstrained negative estimate visible and labels the sign failure; silently clipping it would change the least-squares result.

Passing those sign checks proves little about a complete physical body. Wensing, Kim, and Slotine derive stronger physical-consistency constraints for full rigid-body inertial parameters. A real identification problem can impose appropriate constraints during fitting.

The lab places noise only on torque. Errors in angle, velocity, or acceleration also disturb X, creating an errors-in-variables problem that ordinary least squares does not automatically solve. Differentiating encoder measurements twice can magnify noise, so real experiments need synchronized measurements and a suitable derivative-estimation method.

The model uncertainty lesson carries specified inertia and damping ranges into an acceleration prediction. Its bounds also include a bounded disturbance torque. Those ranges need their own justification; a small training RMSE alone does not supply them.

Unmodeled Coulomb friction or sticking can leave structured residuals that a viscous coefficient cannot explain. Motor voltage and current require an actuator model before you can interpret them as joint torque. WPILib's arm identification model, for example, fits voltage coefficients with different units from I, h, and B here.

Reproduce the identification in Python

This Python 3 example uses only the standard library. It reproduces the 24-sample, seed-7 cases and the 16-state validation trajectory. Its short QR routine uses two-pass modified Gram-Schmidt; the browser uses pivoted Householder QR for its rank check and solve.

from math import sin, cos, sqrt

TRUTH = (0.8, 0.4, 0.15)


def row(motion, t):
    if motion == "rich":
        q = 0.7*sin(1.3*t) + 0.25*sin(2.7*t)
        v = 0.91*cos(1.3*t) + 0.675*cos(2.7*t)
        a = -1.183*sin(1.3*t) - 1.8225*sin(2.7*t)
    elif motion == "slow":
        q, v, a = 0.08*sin(0.4*t), 0.032*cos(0.4*t), -0.0128*sin(0.4*t)
    elif motion == "test":
        q = 0.6*sin(0.9*t) + 0.3*cos(2.1*t)
        v = 0.54*cos(0.9*t) - 0.63*sin(2.1*t)
        a = -0.486*sin(0.9*t) - 1.323*cos(2.1*t)
    else:
        q, v, a = 0, 0, 0
    return [a, 9.81*cos(q), v]


def dot(a, b):
    return sum(x*y for x, y in zip(a, b))


def fit(X, y):
    # Three columns, scaled to unit norm before orthogonalization.
    columns = [list(column) for column in zip(*X)]
    scales = [sqrt(dot(c, c)) or 1 for c in columns]
    Q, R = [], [[0.0]*3 for _ in range(3)]
    for j, column in enumerate(columns):
        v = [x/scales[j] for x in column]
        for _ in range(2):
            for i, qi in enumerate(Q):
                projection = dot(qi, v)
                R[i][j] += projection
                v = [x-projection*z for x, z in zip(v, qi)]
        R[j][j] = sqrt(dot(v, v))
        if R[j][j] <= 1e-10:
            return None
        Q.append([x/R[j][j] for x in v])
    rhs, beta = [dot(qi, y) for qi in Q], [0.0]*3
    for i in range(2, -1, -1):
        beta[i] = (rhs[i]-sum(R[i][j]*beta[j] for j in range(i+1, 3)))/R[i][i]
    return [beta[i]/scales[i] for i in range(3)]


def rmse(X, y, theta):
    return sqrt(sum((dot(x, theta)-target)**2 for x, target in zip(X, y))/len(y))


test_X = [row("test", 8*i/15) for i in range(16)]
test_y = [dot(x, TRUTH) for x in test_X]
for motion, noise in [("rich", 0), ("rich", 0.2), ("slow", 0.2), ("holds", 0)]:
    X = [row(motion, 8*i/23) for i in range(24)]
    seed, y = 7, []
    for x in X:
        seed = (1664525*seed + 1013904223) % (2**32)
        y.append(dot(x, TRUTH) + noise*(2*seed/(2**32)-1))
    theta = fit(X, y)
    if theta is None:
        print(f"{motion}: no unique fit")
        continue
    print(f"{motion}, noise={noise:.2f}: I={theta[0]:.6f}, h={theta[1]:.6f}, B={theta[2]:.6f}")
    print(f"RMSE: training={rmse(X, y, theta):.6f}, held-out={rmse(test_X, test_y, theta):.6f} N m")

Expected output:

rich, noise=0.00: I=0.800000, h=0.400000, B=0.150000
RMSE: training=0.000000, held-out=0.000000 N m
rich, noise=0.20: I=0.793312, h=0.401350, B=0.153888
RMSE: training=0.119280, held-out=0.014150 N m
slow, noise=0.20: I=-1.607788, h=0.400094, B=0.109966
RMSE: training=0.118472, held-out=2.261865 N m
holds: no unique fit

The noise sequence is a repeatable teaching device, not a claim about a particular torque sensor's error distribution. The fit and validation use the same model structure by construction, so these results do not test missing physics.

Try it yourself

Exercise 1. Collect only holds. A joint holds horizontally right in every noiseless sample. Its measured torque is 3.924 N m. Find the regressor rank and the coefficient you can recover; explain what changes if you collect 60 identical samples.

Check what a hold can identify

Every row equals (0, 9.81, 0), so the rank is one. You can recover h = 3.924/9.81 = 0.4 kg m. Inertia I and damping B can vary without changing any predicted holding torque.

Sixty identical samples preserve rank one. In a noisy experiment, repetition can help average torque error under suitable noise assumptions; it cannot reveal absent acceleration and velocity effects. Choose Horizontal holds in the lab to see the unavailable full-parameter fit.

Exercise 2. Use an identified model. Let I = 0.8 kg m², h = 0.4 kg m, and B = 0.15 N m s/rad. Find the actuator torque at q = 60°, q̇ = −2 rad/s, and q̈ = 3 rad/s². Identify the sign of each contribution.

Check the predicted torque

Inertia contributes 0.8 × 3 = 2.4 N m. Gravity compensation contributes 0.4 × 9.81 × cos 60° = 1.962 N m. Damping compensation contributes 0.15 × (−2) = −0.3 N m.

The actuator total is 2.4 + 1.962 − 0.3 = 4.062 N m. Negative joint velocity makes the friction compensation negative. Physical damping still removes mechanical energy: B q̇² = 0.6 W.

Sources and further study

Next, study friction models and look for residual patterns that a single viscous coefficient misses. Then return to inverse dynamics to use a tested model for torque prediction.