nlhw

Estimate Hammerstein-Wiener model

Syntax

sys = nlhw(data,Orders)

sys = nlhw(data,Orders,InputNonlinearity,OutputNonlinearity)

sys = nlhw(data,LinModel)

sys = nlhw(data,LinModel,InputNonlinearity,OutputNonlinearity)

sys = nlhw(data,sys0)

sys = nlhw(___,Options)

Description

Estimate Hammerstein-Wiener Model

sys = nlhw(data,Orders) creates and estimates a Hammerstein-Wiener model using the estimation data, model orders and delays, and default piecewise linear functions as input and output nonlinearity estimators. data can be in the form of a timetable, a comma-separated pair of numeric matrices, or a data object.

sys = nlhw(data,Orders,InputNonlinearity,OutputNonlinearity) specifies InputNL and OutputNL as the input and output nonlinearity estimators, respectively.

Specify Linear Model

sys = nlhw(data,LinModel) uses a linear model to specify the linear block coefficients, and default piecewise linear functions for the input and output nonlinearity estimators.

sys = nlhw(data,LinModel,InputNonlinearity,OutputNonlinearity) specifies InputNonlinearity and OutputNonlinearity as the input and output nonlinearity estimators, respectively.

Refine Existing Model

sys = nlhw(data,sys0) refines or estimates the parameters of a Hammerstein-Wiener model, sys0, using the estimation data.

Use this syntax to:

• Update the parameters of a previously estimated model to improve the fit to the estimation data. In this case, the estimation algorithm uses the parameters of sys0 as initial guesses.

• Estimate the parameters of a model previously created using the idnlhw constructor. Prior to estimation, you can configure the model properties using dot notation.

Specify Options

sys = nlhw(___,Options) specifies additional model estimation options using the option set Options that you create using nlhwOptions. Use Options with any of the previous syntaxes.

Examples

Estimate a Hammerstein-Wiener Model

load iddata3
m1 = nlhw(z3,[4 2 1]);

Estimate a Hammerstein Model with Saturation

Load data.

load twotankdata;
z = iddata(y,u,0.2,'Name','Two tank system');
z1 = z(1:1000);

Create a saturation object with lower limit of 0 and upper limit of 5.

InputNL = idSaturation('LinearInterval',[0 5]);

Estimate model with no output nonlinearity.

m = nlhw(z1,[2 3 0],InputNL,[]);

Estimate Hammerstein-Wiener Model with a Custom Network Nonlinearity

Generating a custom network nonlinearity requires the definition of a user-defined unit function.

Define the unit function and save it as gaussunit.m.

function [f,g,a] = gaussunit(x)
% Custom unit function nonlinearity.
%
% Copyright 2015 The MathWorks, Inc.
f = exp(-x.*x);
if nargout>1
  g = -2*x.*f;
  a = 0.2;
end

Create a custom network nonlinearity using the gaussunit function.

H = @gaussunit;
CNet = idCustomNetwork(H);

Load the estimation data.

load twotankdata;
z = iddata(y,u,0.2,'Name','Two tank system');
z1 = z(1:1000);

Estimate a Hammerstein-Wiener model using the custom network.

m = nlhw(z1,[5 1 3],CNet,[]);

Estimate Default Hammerstein-Wiener Model Using an Input-Output Polynomial Model of OE Structure

Estimate linear OE model.

load throttledata.mat
Tr = getTrend(ThrottleData); 
Tr.OutputOffset = 15;
DetrendedData = detrend(ThrottleData, Tr);
opt = oeOptions('Focus','simulation');
LinearModel = oe(DetrendedData,[1 2 1],opt);

Estimate Hammerstein-Wiener model using OE model as its linear component and saturation as its output nonlinearity.

sys = nlhw(ThrottleData,LinearModel,[],idSaturation);

Estimate a Hammerstein-Wiener Model Using idnlhw to First Define the Model Properties

Load the estimation data.

load iddata1

Construct a Hammerstein-Wiener model using idnlhw to define the model properties B and F.

sys0 = idnlhw([2,2,0],[],'idWaveletNetwork');
sys0.B{1} = [0.8,1];
sys0.F{1} = [1,-1.2,0.5];

Estimate the model.

sys = nlhw(z1,sys0);

Estimate a Hammerstein-Wiener model using nlhw to define the model properties B and F.

sys2 = nlhw(z1,[2,2,0],[],'idWaveletNetwork','B',{[0.8,1]},'F',{[1,-1.2,0.5]});

Compare the two estimated models to see that they are equivalent.

compare(z1,sys,'g',sys2,'r--');

Refine a Hammerstein-Wiener Model Using Successive Calls of nlhw

Estimate a Hammerstein-Wiener Model.

load iddata3
sys = nlhw(z3,[4 2 1],'idSigmoidNetwork','idWaveletNetwork');

Refine the model, sys.

sys = nlhw(z3,sys);

Estimate Hammerstein-Wiener Model Using an Estimation Option Set

Create estimation option set for nlhw to view estimation progress, use the Levenberg-Marquardt search method, and set the maximum iteration steps to 50.

opt = nlhwOptions;
opt.Display = 'on';
opt.SearchMethod = 'lm';    
opt.SearchOptions.MaxIterations = 50;

Load data and estimate the model.

load iddata3
sys = nlhw(z3,[4 2 1],idSigmoidNetwork,idPiecewiseLinear,opt);

Input Arguments

data — Estimation data

timetable | cell array of timetables | numeric matrix pair | cell array pair of matrices | iddata object

Uniformly sampled estimation data, specified as described in the following sections.

Timetable

Specify data as a timetable that uses a regularly spaced time vector. tt contains variables representing input and output channels. For multiexperiment data, tt is a cell array of timetables of length Ne, where Ne is the number of experiments

The software determines the number of input and output channels to use for estimation from the dimensions of the order in Orders. The input/output channel selection depends on whether the 'InputName' and 'OutputName' name-value arguments are specified.

• If 'InputName' and 'OutputName' are not specified, then the software uses the first Nu variables of tt as inputs and the next Ny variables of tt as outputs.

• If 'InputName' and 'OutputName' are specified, then the software uses the specified variables. The number of specified input and output names must be consistent with Nu and Ny.

• For functions that can estimate a time series model, where there are no inputs, 'InputName' does not need to be specified.

Comma-Separated Matrix pair

Specify data as a comma-separated pair of real-valued matrices that contain uniformly sampled input and output time-domain signal values. For multiexperiment data, use a cell array of matrices. When you specify matrix-based data, the software assumes a sample time of 1 second. You can change the sample time after estimation by setting the property sys.Ts.

• For SISO systems, specify data as a pair of N_s-by-1 real-valued matrices that contain uniformly sampled input and output time-domain signal values. Here, N_s is the number of samples.

• For MIMO systems, specify u,y as an input/output matrix pair with the following dimensions:

  • u — N_s-by-N_u, where N_u is the number of inputs.

  • y — N_s-by-N_y, where N_y is the number of outputs.

• For multiexperiment data, specify u,y as a pair of 1-by-N_e cell arrays, where N_e is the number of experiments. The sample times of all the experiments must match.

Data Object

An estimation data object, specified as a time-domain iddata object that contains uniformly sampled input and output values. By default, the software sets the sample time of the model to the sample time of the estimation data.

For multiexperiment data, the sample times and intersample behavior of all the experiments must match.

For more information about working with estimation data types, see Data Domains and Data Types in System Identification Toolbox.

Orders — Order and delays of the linear subsystem transfer function

[nb nf nk] vector of positive integers | [nb nf nk] vector of matrices

Order and delays of the linear subsystem transfer function, specified as a [nb nf nk] vector.

Dimensions of Orders:

• For a SISO transfer function, Orders is a vector with 3 positive integers.

  nb is the number of zeros plus 1, nf is the number of poles, and nk is the input delay.

• For a MIMO transfer function with nu inputs and ny outputs, Orders is a vector of matrices.

  nb, nf, and nk are ny-by-nu matrices whose i-jth entry specifies the orders and delay of the transfer function from the jth input to the ith output.

InputNonlinearity — Input nonlinearity estimator

'idPiecewiseLinear' (default) | 'idPiecewiseConstant' | 'idSigmoidNetwork' | 'idWaveletNetwork' | 'idSaturation' | idDeadZone | 'idPolynomial1D' | 'idUnitGain' | 'idNeuralNetwork' | idCustomNetwork object | array of nonlinearity estimators

Input nonlinearity estimator, specified as a column array containing one or more of the following strings or mapping objects. Note that idGaussianProcess, which can be used as an output nonlinearity estimator, cannot be used as an input nonlinearity estimator.

'idPiecewiseLinear' or idPiecewiseLinear object	Piecewise linear function
'idPiecewiseConstant' or idPiecewiseConstant object	Piecewise constant function
'idSigmoidNetwork' or idSigmoidNetwork object	Sigmoid network
'idWaveletNetwork' or idWaveletNetwork object	Wavelet network
'idSaturation' or idSaturation object	Saturation
'idDeadZone' or idDeadZone object	Dead zone
'idPolynomial1D' or idPolynomial1D object	One-dimensional polynomial
idCustomNetwork object	Custom network — Similar to idSigmoidNetwork, but with a user-defined replacement for the sigmoid function.
'idNeuralNetwork' or idNeuralNetwork object	Multilayer neural network (requires Statistics and Machine Learning Toolbox™ or Deep Learning Toolbox™)
'idUnitGain' or [] or idUnitGain object	Unit gain. Effectively eliminates nonlinearity block.

Specifying a character vector, for example 'idSigmoidNetwork', creates a mapping object with default settings. Alternatively, you can specify nonlinearity estimator properties in two other ways:

• Create the nonlinearity function using arguments to modify default properties.

  InputNL = idSigmoidNetwork(15)

• Create a default nonlinearity function first and then use dot notation to modify properties.

  InputNL = idSigmoidNetwork;
  InputNL.NumberOfUnits = 15

For nu input channels, you can specify nonlinear estimators individually for each input channel by setting InputNL to an nu-by-1 array of nonlinearity estimators. To specify the same nonlinearity for all inputs, specify a single input nonlinearity estimator.

OutputNonlinearity — Output nonlinearity estimator

'idPiecewiseLinear' (default) | 'idPiecewiseConstant' | 'idSigmoidNetwork' | 'idWaveletNetwork' | 'idSaturation' | idDeadZone | 'idPolynomial1D' | 'idGaussianProcess' | 'idUnitGain' | 'idNeuralNetwork' | idCustomNetwork object | array of nonlinearity estimators

Output nonlinearity estimator, specified as a column array containing one or more of the following strings or mapping objects.

'idPiecewiseLinear' or idPiecewiseLinear object	Piecewise linear function
'idPiecewiseConstant' or idPiecewiseConstant object	Piecewise constant function
'idSigmoidNetwork' or idSigmoidNetwork object	Sigmoid network
'idWaveletNetwork' or idWaveletNetwork object	Wavelet network
'idSaturation' or idSaturation object	Saturation
'idDeadZone' or idDeadZone object	Dead zone
'idPolynomial1D' or idPolynomial1D object	One-dimensional polynomial
'idGaussianProcess' or idGaussianProcess object	Gaussian process regression model (requires Statistics and Machine Learning Toolbox)
idCustomNetwork object	Custom network — Similar to idSigmoidNetwork, but with a user-defined replacement for the sigmoid function.
'idNeuralNetwork' or idNeuralNetwork object	Multilayer neural network (requires Statistics and Machine Learning Toolbox or Deep Learning Toolbox)
'idUnitGain' or [] or idUnitGain object	Unit gain. Effectively eliminates nonlinearity block.

Specifying a character vector, for example 'idSigmoidNetwork', creates a mapping object with default settings. Alternatively, you can specify nonlinearity estimator properties in two other ways:

• Create the nonlinearity function using arguments to modify default properties.

  NL = idSigmoidNetwork(15)

• Create a default nonlinearity function first and then use dot notation to modify properties.

  outputNL = idSigmoidNetwork;
  OutputNL.NumberOfUnits = 15

For ny output channels, you can specify nonlinear estimators individually for each output channel by setting OutputNL to an ny-by-1 array of nonlinearity estimators. To specify the same nonlinearity for all outputs, specify a single output nonlinearity estimator.

LinModel — Discrete time linear model

idpoly model | idss model | idtf model

Discrete-time linear model used to specify the linear subsystem, specified as one of the following:

• Input-output polynomial model of Output-Error (OE) structure (idpoly) with nonnegative input delay

• State-space model (idss)

• Transfer function model (idtf)

Typically, you estimate the model using oe, n4sid, or tfest.

sys0 — Hammerstein-Wiener model

idnlhw object

Hammerstein-Wiener model, specified as an idnlhw object. sys0 can be:

• A model previously created using idnlhw to specify model properties.

• A model previously estimated using nlhw, that you want to update using a new estimation data set.

  You can also refine sys0 using the original estimation data set. If the previous estimation stopped when the numerical search was stuck at a local minima of the cost function, use init to first randomize the parameters of sys0. See sys0.Report.Termination for search stopping conditions. Using init does not guarantee a better solution on further refinement.

Options — Estimation options

nlhwOptions option set

Estimation options for Hammerstein-Wiener model identification, specified as an nlhwOptions option set. Available options include:

• Search options

• Normalization options

• Regularization options

Output Arguments

sys — Estimated Hammerstein-Wiener model

idnlhw object

Estimated Hammerstein-Wiener model, returned as an idnlhw object. The model is estimated using the specified model orders, input and output nonlinearity estimators, 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
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 these fields. FitPercent — Normalized root mean squared error (NRMSE) measure of how well the response of the model fits the estimation data, expressed as the percentage fitpercent = 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 nlhwOptions 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. 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' — A zero-order hold maintains a piecewise-constant input signal between samples. 'foh' — A 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 these fields. 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, see Estimation Report.

Version History

Introduced in R2007a

R2025a: New piecewise-constant nonlinearity estimator for estimating Hammerstein-Wiener Models

The idPiecewiseConstant object stores the piecewise-constant nonlinearity estimator for estimating Hammerstein-Wiener models. You can use idPiecewiseConstant as both an input and output nonlinearity. For more information, see idPiecewiseConstant.

R2023b: New neural network object creates neural networks from both Statistics and Machine Learning Toolbox and Deep Learning Toolbox

The idNeuralNetwork mapping object creates neural networks using both the regression networks of Statistics and Machine Learning Toolbox and the shallow or deep networks of Deep Learning Toolbox. You can use idNeuralNetwork as both an input and output nonlinearity. For more information, see idNeuralNetwork.

R2022b: Time-domain estimation data is accepted in the form of timetables and matrices

Most estimation, validation, analysis, and utility functions now accept time-domain input/output data in the form of a single timetable that contains both input and output data or a pair of matrices that contain the input and output data separately. These functions continue to accept iddata objects as a data source as well, for both time-domain and frequency-domain data.

R2021b: Use of previous idnlarx and idnlhw mapping object names is not recommended.

Starting in R2021b, the mapping objects (also known as nonlinearities) used in the nonlinear components of the idnlarx and idnlhw objects have been renamed. The following table lists the name changes.

Pre-R2021b Name	R2021b Name
wavenet	idWaveletNetwork
sigmoidnet	idSigmoidNetwork
treepartition	idTreePartition
customnet	idCustomNetwork
saturation	idSaturation
deadzone	idDeadZone
pwlinear	idPiecewiseLinear
poly1d	idPolynomial1D
unitgain	idUnitGain
linear	idLinear
neuralnet	idFeedforwardNetwork

Scripts with the old names still run normally, although they will produce a warning. Consider using the new names for continuing compatibility with newly developed features and algorithms. There are no plans to exclude the use of these object names at this time.

R2018a: Advanced Options are deprecated for SearchOptions when SearchMethod is 'lsqnonlin'

Specification of lsqnonlin- related advanced options are deprecated, including the option to invoke parallel processing when estimating using the lsqnonlin search method, or solver, in Optimization Toolbox™.