greyest

Linear grey-box model estimation

Syntax

sys = greyest(data,init_sys)
sys = greyest(data,init_sys,opt)
[sys,x0] = greyest(___)


Description

sys = greyest(data,init_sys) estimates a linear grey-box model, sys, using time or frequency domain data, data. The dimensions of the inputs and outputs of data and init_sys, an idgrey model, must match. sys is an identified idgrey model that has the same structure as init_sys.

sys = greyest(data,init_sys,opt) estimates a linear grey-box model using the option set, opt, to configure the estimation options.

[sys,x0] = greyest(___) returns the value of the initial states computed during estimation. You can use this syntax with any of the previous input-argument combinations.

Input Arguments

data	Estimation data. The dimensions of the inputs and outputs of data and init_sys must match. For time-domain estimation, data is an iddata object containing the input and output signal values. For frequency domain estimation, data can be one of the following: Recorded frequency response data (frd (Control System Toolbox) or idfrd) iddata object with its Domain property set to 'Frequency'
init_sys	Identified linear grey-box model that configures the initial parameterization of sys. init_sys, an idgrey model, must have the same input and output dimensions as data.
opt	Estimation options. opt is an option set, created using greyestOptions, which specifies options including: Estimation objective Initialization choice Disturbance model handling Numerical search method to be used in estimation

Output Arguments

sys

Estimated grey-box model, returned as an idgrey model. This model is created using the specified initial system, and estimation options. Information about the estimation results and options used is stored in the Report property of the model. Report has the following fields: Report Field Description Status Summary of the model status, which indicates whether the model was created by construction or obtained by estimation. Method Estimation command used. InitialState Handling of initial states during estimation, returned as one of the following: 'model' — The initial state is parameterized by the ODE file used by the idgrey model. 'zero' — The initial state is set to zero. 'estimate' — The initial state is treated as an independent estimation parameter. 'backcast' — The initial state is estimated using the best least squares fit. Vector of doubles of length Nx, where Nx is the number of states. For multiexperiment data, a matrix with Ne columns, where Ne is the number of experiments. This field is especially useful to view how the initial states were handled when the InitialState option in the estimation option set is 'auto'. DisturbanceModel Handling of the disturbance component (K) during estimation, returned as one of the following values: 'model' — K values are parameterized by the ODE file used by the idgrey model. 'fixed' — The value of the K property of the idgrey model is fixed to its original value. 'none' — K is fixed to zero. 'estimate' — K is treated as an independent estimation parameter. This field is especially useful to view the how the disturbance component was handled when the DisturbanceModel option in the estimation option set is 'auto'. Fit Quantitative assessment of the estimation, returned as a structure. See Loss Function and Model Quality Metrics for more information on these quality metrics. The structure has the following fields: Field Description FitPercent Normalized root mean squared error (NRMSE) measure of how well the response of the model fits the estimation data, expressed as the percentage fit = 100(1-NRMSE). LossFcn Value of the loss function when the estimation completes. MSE Mean squared error (MSE) measure of how well the response of the model fits the estimation data. FPE Final prediction error for the model. AIC Raw Akaike Information Criteria (AIC) measure of model quality. AICc Small sample-size corrected AIC. nAIC Normalized AIC. BIC Bayesian Information Criteria (BIC). Parameters Estimated values of model parameters. OptionsUsed Option set used for estimation. If no custom options were configured, this is a set of default options. See greyestOptions for more information. RandState State of the random number stream at the start of estimation. Empty, [], if randomization was not used during estimation. For more information, see rng. DataUsed Attributes of the data used for estimation, returned as a structure with the following fields: Field Description Name Name of the data set. Type Data type. Length Number of data samples. Ts Sample time. InterSample Input intersample behavior, returned as one of the following values: 'zoh' — Zero-order hold maintains a piecewise-constant input signal between samples. 'foh' — First-order hold maintains a piecewise-linear input signal between samples. 'bl' — Band-limited behavior specifies that the continuous-time input signal has zero power above the Nyquist frequency. InputOffset Offset removed from time-domain input data during estimation. For nonlinear models, it is []. OutputOffset Offset removed from time-domain output data during estimation. For nonlinear models, it is []. Termination Termination conditions for the iterative search used for prediction error minimization, returned as a structure with the following fields: Field Description WhyStop Reason for terminating the numerical search. Iterations Number of search iterations performed by the estimation algorithm. FirstOrderOptimality $$\infty$$-norm of the gradient search vector when the search algorithm terminates. FcnCount Number of times the objective function was called. UpdateNorm Norm of the gradient search vector in the last iteration. Omitted when the search method is 'lsqnonlin' or 'fmincon'. LastImprovement Criterion improvement in the last iteration, expressed as a percentage. Omitted when the search method is 'lsqnonlin' or 'fmincon'. Algorithm Algorithm used by 'lsqnonlin' or 'fmincon' search method. Omitted when other search methods are used. For estimation methods that do not require numerical search optimization, the Termination field is omitted. For more information on using Report, see Estimation Report.

x0

Initial states computed during the estimation, returned as a matrix containing a column vector corresponding to each experiment. This array is also stored in the Parameters field of the model Report property.

Examples

Estimate Grey-Box Model

Estimate the parameters of a DC motor using the linear grey-box framework.

Load the measured data.

load(fullfile(matlabroot, 'toolbox', 'ident', 'iddemos', 'data', 'dcmotordata'));
data = iddata(y, u, 0.1, 'Name', 'DC-motor');
data.InputName = 'Voltage';
data.InputUnit = 'V';
data.OutputName = {'Angular position', 'Angular velocity'};
data.OutputUnit = {'rad', 'rad/s'};
data.Tstart = 0;
data.TimeUnit = 's';

data is an iddata object containing the measured data for the outputs, the angular position, the angular velocity. It also contains the input, the driving voltage.

Create a grey-box model representing the system dynamics.

For the DC motor, choose the angular position (rad) and the angular velocity (rad/s) as the outputs and the driving voltage (V) as the input. Set up a linear state-space structure of the following form:

$$\underset{}{\overset{\dot{}}{x}}(t)=\left[\begin{array}{cc}0& 1\\ 0& -\frac{1}{\tau }\end{array}\right]x(t)+\left[\begin{array}{c}0\\ \frac{G}{\tau }\end{array}\right]u(t)$$

$$y(t)=\left[\begin{array}{cc}1& 0\\ 0& 1\end{array}\right]x(t).$$

$$\tau$$ is the time constant of the motor in seconds, and $$G$$ is the static gain from the input to the angular velocity in rad/(V*s) .

G = 0.25; 
tau = 1; 

init_sys = idgrey('motorDynamics',tau,'cd',G,0);

The governing equations in state-space form are represented in the MATLAB® file motorDynamics.m. To view the contents of this file, enter edit motorDynamics.m at the MATLAB command prompt.

$$G$$ is a known quantity that is provided to motorDynamics.m as an optional argument.

$$\tau$$ is a free estimation parameter.

init_sys is an idgrey model associated with motor.m.

Estimate $$\tau$$.

sys = greyest(data,init_sys);

sys is an idgrey model containing the estimated value of $$\tau$$.

To obtain the estimated parameter values associated with sys, use getpvec(sys).

Analyze the result.

opt = compareOptions('InitialCondition','zero');
compare(data,sys,Inf,opt)

sys provides a 98.35% fit for the angular position and an 84.42% fit for the angular velocity.

Estimate Grey-Box Model Using Regularization

Estimate the parameters of a DC motor by incorporating prior information about the parameters when using regularization constants.

The model is parameterized by static gain G and time constant $$\tau$$. From prior knowledge, it is known that G is about 4 and $$\tau$$ is about 1. Also, you have more confidence in the value of $$\tau$$ than G and would like to guide the estimation to remain close to the initial guess.

Load estimation data.

load regularizationExampleData.mat motorData

The data contains measurements of motor's angular position and velocity at given input voltages.

Create an idgrey model for DC motor dynamics. Use the function DCMotorODE that represents the structure of the grey-box model.

mi = idgrey(@DCMotorODE,{'G', 4; 'Tau', 1},'cd',{}, 0);
mi = setpar(mi, 'label', 'default');

If you want to view the DCMotorODE function, type:

type DCMotorODE.m

function [A,B,C,D] = DCMotorODE(G,Tau,Ts)
%DCMOTORODE ODE file representing the dynamics of a DC motor parameterized
%by gain G and time constant Tau.
%
%   [A,B,C,D,K,X0] = DCMOTORODE(G,Tau,Ts) returns the state space matrices
%   of the DC-motor with time-constant Tau and static gain G. The sample
%   time is Ts.
%
%   This file returns continuous-time representation if input argument Ts
%   is zero. If Ts>0, a discrete-time representation is returned.
%
% See also IDGREY, GREYEST.

%   Copyright 2013 The MathWorks, Inc.

A = [0 1;0 -1/Tau];
B = [0; G/Tau];
C = eye(2);
D = [0;0];
if Ts>0 % Sample the model with sample time Ts
   s = expm([[A B]*Ts; zeros(1,3)]);
   A = s(1:2,1:2);
   B = s(1:2,3);
end

Specify regularization options Lambda.

opt = greyestOptions;
opt.Regularization.Lambda = 100;

Specify regularization options R.

opt.Regularization.R = [1, 1000];

You specify more weighting on the second parameter because you have more confidence in the value of $$\tau$$ than G.

Specify the initial values of the parameters as regularization option $$\theta$$*.

opt.Regularization.Nominal = 'model';

Estimate the regularized grey-box model.

sys = greyest(motorData, mi, opt);

Extended Capabilities

Automatic Parallel Support

Accelerate code by automatically running computation in parallel using Parallel Computing Toolbox™.

Parallel computing support is available for estimation using the lsqnonlin search method (requires Optimization Toolbox™). To enable parallel computing, use greyestOptions, set SearchMethod to 'lsqnonlin', and set SearchOptions.Advanced.UseParallel to true.

For example:

opt = greyestOptions;
opt.SearchMethod = 'lsqnonlin';
opt.SearchOptions.Advanced.UseParallel = true;