Набор данных ввода - вывода

В этом примере переменные данных сохранены в twotankdata.mat файл будет использоваться. Это содержит 3 000 выборок данных ввода - вывода системы 2D бака, сгенерированной с помощью шага расчета 0,2 секунд. Вход u (t) является напряжением [V], применился к насосу, который генерирует приток в верхний бак. Довольно маленькое отверстие в нижней части этого верхнего бака дает к оттоку, который входит в более низкий бак, и выход y (t) двух систем бака является затем уровнем жидкости [m] более низкого бака. Мы создаем объект IDDATA z инкапсулировать загруженные данные:

Этот набор данных также используется в нелинейном сером примере поля "Две Системы Бака: Моделирование Файла MEX на C Непрерывной Системы SISO", где физические уравнения, описывающие систему, приняты. Этот пример о моделях черного квадрата, таким образом никакое знание фактических физических законов не принято.

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

Разделение данных и графического вывода

Разделите этот набор данных в 3 подмножества равного размера, каждый содержащий 1 000 выборок.

Мы оценим модели с z1, и оцените их на z2 и z3. z1 похоже на z2, но не z3. Поэтому в некотором смысле, z2 помог бы нам оценить способность к интерполяции моделей, в то время как z3 помог бы нам оценить их способность к экстраполяции.

z1 = z(1:1000);
z2 = z(1001:2000);
z3 = z(2001:3000);
plot(z1,z2,z3)
legend('Estimation','Validation 1', 'Validation 2')

Выбор регрессора

Первый шаг в оценке нелинейных моделей черного квадрата должен выбрать набор регрессора. Регрессоры являются простыми формулами на основе задержанных временем переменных ввода/вывода, самый простой случай, являющийся переменными, изолированными маленьким набором последовательных значений. Например, если "u" является именем входной переменной и "y" имя выходной переменной, то набор регрессора в качестве примера может быть {y (t-1), y (t-2), u (t), u (t-1), u (t-2)}, где "t" обозначает переменную времени. Другим примером включающие полиномиальные термины мог быть {abs (y (t-2)), y (t-2) *u (t-1), u (t-4) ^2}. Более комплексные, произвольные формулы в задержанных переменных также возможны.

В System Identification Toolbox™ наборы регрессоров не должны быть созданы явным образом; большое количество может быть сгенерировано неявно, только задав их форму (такую как полиномиальный порядок), способствующие переменные и их задержки. Спецификация возражает linearRegressor, polynomialRegressor и customRegressor используются с этой целью. Например, R = linearRegressor ({'y', 'u'}, {[1 2 4], [2 10]}) подразумевает набор регрессора {y (t-1), y (t-2), y (t-4), u (t-2), u (t-10)}. Точно так же R = polynomialRegressor ({'y'}, {[1 2]}, 3) подразумевает 3-и регрессоры полинома порядка в переменной 'y' с задержками 1 и 2, то есть, набор {y (t-1) ^3, y (t-2) ^3}. Больше вариантов настройки доступно, такие как использование только абсолютных значений и/или разрешения соединения задержек в формулах регрессора.

Линейная спецификация регрессора Используя матрицу порядка

Когда модель использует только линейные регрессоры, которые основаны на непрерывных задержках, часто легче задать их использующий матрицу порядка, вместо того, чтобы использовать linearRegressor объект. Матрица порядка [na nb nk], задает количества мимо выходных параметров, прошлых входных параметров и входной задержки, используемой в формулах регрессора. Это - та же идея, которая используется при оценке линейных моделей ARX (см. ARX, IDPOLY). Например, NN = [2 3 4] подразумевает, что задержки использования выходной переменной (1,2) и входная переменная используют задержки (4,5,6) продвижение к набору регрессора {y (t-1), y (t-2), u (t-4), u (t-5), u (t-6)}.

Порядки модели Usually выбраны методом проб и ошибок. Иногда полезно попробовать линейные модели ARX сначала для того, чтобы получить подсказки о лучших порядках. Таким образом давайте сначала попытаемся определить оптимальные порядки для линейной модели ARX, с помощью функций arxstruc и selstruc.

V = arxstruc(z1,z2,struc(1:5, 1:5,1:5));
% select best order by Akaike's information criterion (AIC)
nn = selstruc(V,'aic')

nn =

     5     1     3

Критерий AIC выбрал nn = [na nb nk] = [5 1 3], означая, что в выбранной структуре модели ARX y (t) предсказан этими 6 регрессорами y (t-1), y (t-1)..., y (t-5) и u (t-3). Мы попытаемся использовать эти значения при оценке нелинейных моделей.

Нелинейный ARX (IDNLARX) модели

В модели IDNLARX выход модели является нелинейной функцией регрессоров, как y (t) = Fcn (y (t-1), y (t-1)..., y (t-5), u (t-3)).

Модели IDNLARX обычно оцениваются с синтаксисом:

M = NLARX(Data, Orders, Fcn).

Это похоже на линейную команду ARX с дополнительным аргументом "Fcn", задающим тип нелинейной функции. Модель Nonlinear ARX производит свой выход при помощи 2-этапного преобразования: 1. Сопоставьте обучающие данные (сигналы ввода - вывода) к набору регрессоров. Численно, набор регрессора является двойной матрицей R с одним столбцом данных для каждого регрессора. 2. Сопоставьте регрессоры с одним выходом, y = Fcn®, где Fcn () является скалярной нелинейной функцией.

Чтобы оценить модель IDNLARX, после выбора порядков модели, мы должны выбрать тип нелинейной функции отображения Fcn (), чтобы использовать. Давайте сначала попробуем Сеть Вейвлета (idWaveletNetwork) функция.

mw1 = nlarx(z1,[5 1 3], idWaveletNetwork);

Сеть Вейвлета использует комбинацию модулей вейвлета, чтобы сопоставить регрессоры с выходом модели. Модель mw1 объект @idnalrx. Это хранит нелинейную функцию отображения (idWaveletNetwork объект здесь) в его OutputFcn свойство. Количество этих модулей, чтобы использовать может быть задано заранее, или мы можем предоставить алгоритму оценки право определять оптимальное значение автоматически. Свойство NonlinearFcn.NumberOfUnits хранит выбранное количество модулей в модели.

NLFcn = mw1.OutputFcn;
NLFcn.NonlinearFcn.NumberOfUnits

ans =

     3

Исследуйте результат путем сравнения ответа модели mw1 к измеренному выходу в наборе данных z1.

compare(z1,mw1);

Используя свойства модели

Первая модель была не совсем удовлетворительной, таким образом давайте попытаемся изменить некоторые свойства модели. Вместо того, чтобы позволить алгоритму оценки автоматически выбрать количество модулей в функции idWaveletNetwork, может быть явным образом задан этот номер.

mw2 = nlarx(z1,[5 1 3], idWaveletNetwork(8));

InputName по умолчанию и OutputName использовались.

mw2.InputName
mw2.OutputName

ans =

  1×1 cell array

    {'u1'}


ans =

  1×1 cell array

    {'y1'}

Регрессоры модели IDNLARX могут быть просмотрены с помощью команды GETREG. Регрессоры сделаны из mw2.InputName, mw2.OutputName и задержки.

getreg(mw2)

ans =

  6×1 cell array

    {'y1(t-1)'}
    {'y1(t-2)'}
    {'y1(t-3)'}
    {'y1(t-4)'}
    {'y1(t-5)'}
    {'u1(t-3)'}

Обратите внимание на то, что матрица порядка эквивалентна линейному набору регрессора, созданному с помощью объекта спецификации, как в:

R = linearRegressor([z1.OutputName; z1.InputName],{1:5, 3});
mw2_a = nlarx(z1, R, idWaveletNetwork(8));
getreg(mw2_a)

ans =

  6×1 cell array

    {'y1(t-1)'}
    {'y1(t-2)'}
    {'y1(t-3)'}
    {'y1(t-4)'}
    {'y1(t-5)'}
    {'u1(t-3)'}

mw2_a по существу та же модель как mw2. Использование объекта спецификации позволяет больше гибкости; используйте его, когда вы не хотите использовать задержки continguous или иметь минимальную задержку в выходных переменных, которая больше того.

Присвоение регрессоров к отображению функций

Выход модели IDNLARX является статической функцией своих регрессоров. Обычно эта функция отображения имеет линейный блок, нелинейный блок, плюс термин смещения (выходное смещение). Выход модели является суммой выходных параметров этих блоков и смещения. Для того, чтобы уменьшать сложность модели, только подмножество регрессоров может быть выбрано, чтобы ввести нелинейный блок и линейный блок. Присвоение регрессоров к линейным и нелинейным компонентам функции отображения выполняется с помощью свойства RegressorUsage модели.

RegUseTable = mw2.RegressorUsage

RegUseTable =

  6×2 table

               y1:LinearFcn    y1:NonlinearFcn
               ____________    _______________

    y1(t-1)       true              true      
    y1(t-2)       true              true      
    y1(t-3)       true              true      
    y1(t-4)       true              true      
    y1(t-5)       true              true      
    u1(t-3)       true              true      

RegUseTable является таблицей, которая показывает, какие регрессоры используются в качестве входных параметров к линейному и нелинейным компонентам Сети Вейвлета. Каждая строка соответствует одному регрессору. Предположим, что мы хотим использовать только те регрессоры, которые являются функциями входной переменной 'u1' в нелинейном компоненте. Эта настройка может быть достигнута можно следующим образом.

RegNames = RegUseTable.Properties.RowNames;
I = contains(RegNames,'u1');
RegUseTable.("y1:NonlinearFcn")(~I) = false;

Как пример, рассмотрите модель, которая использует полиномиальные регрессоры порядка 2 в дополнение к линейным единицам, используемым mw2. Во-первых, сгенерируйте полиномиальный набор спецификации регрессора.

P = polynomialRegressor({'y1','u1'},{1:2, 0:3},2)

P = 

<strong>Order 2 regressors in variables y1, u1</strong>
               Order: 2
           Variables: {'y1'  'u1'}
                Lags: {[1 2]  [0 1 2 3]}
         UseAbsolute: [0 0]
    AllowVariableMix: 0
         AllowLagMix: 0
        TimeVariable: 't'

Теперь добавьте спецификацию P к свойству Regressors модели.

mw3 = mw2; % start with the mw2 structure for new estimation
mw3.Regressors = [mw3.Regressors; P]; % add the polynomial set to the model structure

Наконец, обновите параметры модели, чтобы соответствовать данным путем минимизации 1 шага вперед ошибки предсказания.

mw3 = nlarx(z1, mw3)

mw3 = 

<strong>Nonlinear ARX model with 1 output and 1 input</strong>
  Inputs: u1
  Outputs: y1

Regressors:
  1. Linear regressors in variables y1, u1
  2. Order 2 regressors in variables y1, u1

Output function: Wavelet network with 8 units
Sample time: 0.2 seconds

Status:                                                     
Estimated using NLARX on time domain data "Two tank system".
Fit to estimation data: 96.76% (prediction focus)           
FPE: 3.499e-05, MSE: 3.338e-05

Оценка предполагаемых моделей

Различные модели, и линейные и нелинейные, могут быть сравнены вместе в той же команде COMPARE.

mlin = arx(z1,[5 1 3]);  % linear ARX model
%
compare(z1,mlin,mw2,mw3)

Проверка допустимости модели может быть выполнена путем выполнения команды COMPARE на наборах данных валидации. Во-первых, подтвердите модели путем проверки, как хорошо их симулированные отклики совпадают с измеренным выходным сигналом в z2.

compare(z2,mlin,mw2,mw3)

Точно так же подтвердите предполагаемые модели против набора данных z3.

compare(z3,mlin,mw2,mw3)

Анализ остатков модели может также быть выполнен, чтобы подтвердить качество модели. Мы ожидаем, что остатки будут белыми (некоррелированый с собой кроме в задержке 0) и некоррелированый с входным сигналом. resideus являются 1 шагом вперед ошибки предсказания. Остатки могут визуализироваться с границами на незначительных значениях при помощи команды RESID. Например, на valiation наборе данных z2, остатки для трех идентифицированных моделей показывают ниже.

resid(z2,mlin,mw2,mw3)
legend show

Остатки являются все некоррелироваными с входным сигналом z2 (z2. InputData). Однако они показывают некоторую автокорреляцию для моделей mlin и mw3. Модель mw2 выполняет лучше и на (авто и на взаимная корреляция) тесты.

График функций может использоваться, чтобы просмотреть ответы нелинейности различных моделей IDNLARX. График по существу показывает характеристики нелинейной функции отображения Fcn ().

plot(mw2,mw3) % plot nonlinearity response as a function of selected regressors

Нелинейная модель ARX с сигмоидальной сетью

Вместо Сети Вейвлета другие нелинейные функции могут использоваться в модели. Попробуйте функцию idSigmoidNetwork, которая сопоставляет регрессоры с выходом с помощью суммы сигмоидальных модулей (плюс линейное и термин смещения). Сконфигурируйте сеть, чтобы использовать 8 модулей сигмоидальных (расширенный и переведенный различными суммами).

ms1 = nlarx(z1,[5 1 3], idSigmoidNetwork(8));
compare(z1,ms1)

Теперь оцените новую модель на наборе данных z2.

compare(z2,ms1);

и на наборе данных z3.

compare(z3,ms1);

Другие нелинейные функции отображения в модели IDNLARX

idCustomNetwork похож на idSigmoidNetwork, но вместо сигмоидальной модульной функции, пользователи могут задать свою собственную модульную функцию в файлах MATLAB. Этот пример использует функцию gaussunit.

type gaussunit.m

function [f, g, a] = gaussunit(x)
%GAUSSUNIT Unit function for use in idCustomNetwork (example)
%
%[f, g, a] = GAUSSUNIT(x)
%
% x: unit function variable
% f: unit function value
% g: df/dx
% a: unit active range (g(x) is significantly non zero in the interval [-a a])
%
% The unit function must be "vectorized": for a vector or matrix x, the output
% arguments f and g must have the same size as x, computed element-by-element.

% Copyright 2005-2021 The MathWorks, Inc.
% Author(s): Qinghua Zhang

f =  exp(-x.*x);

if nargout>1
   g = - 2*x .* f;
   a = 0.2;
end

Для customnet базирующаяся модель используйте только третье, пятое и шестые регрессоры, как введено к нелинейному компоненту функции customnet.

mc1_ini = idnlarx('y1','u1', [5 1 3], idCustomNetwork(@gaussunit));
Use = false(6,1);
Use([3 5 6]) = true;
mc1_ini.RegressorUsage{:,2} = Use;
mc1 = nlarx(z1,mc1_ini);

Древовидная функция раздела является другой функцией отображения, которая может использоваться в модели Nonlinear ARX. Оцените модель с помощью idTreePartition функцию:

mt1 = nlarx(z1, [5 1 3], idTreePartition);

Сравните ответы моделей mc1 и mt1 к данным z1.

compare(z1, mc1, mt1)

Используя алгоритмы регрессии из Statistics and Machine Learning Toolbox

Если у вас есть доступ к Statistics and Machine Learning Toolbox™, вы можете также использовать функции гауссова процесса (GP) и сложенный в мешок или повысили древовидные ансамбли, чтобы создать ваши модели Nonlinear ARX. GP представлены объектом idGaussianProcess, в то время как ансамбли дерева регрессии представлены объектом idTreeEnsemble. GP используют ядро в квадрате экспоненциальное по умолчанию. Другие типы ядра могут быть заданы при создании объекта. В примере здесь, мы используем ядро 'Matern32'.

if exist('fitrgp','file')==2
    Warn = warning('off','stats:lasso:MaxIterReached');
    mp1 = nlarx(z1, [5 1 3], idGaussianProcess('Matern32'))

    % Create boosted ensemble of 200 regression trees
    En = idTreeEnsemble;
    En.EstimationOptions.FitMethod = 'lsboost-resampled';
    En.EstimationOptions.NumLearningCycles = 200;
    En.EstimationOptions.ResampleFraction = 1;
    mp2 = nlarx(z1, [0 70 3], En)
    compare(z1, mp1, mp2) % compare model fits to estimation data
    warning(Warn)
end

mp1 = 

<strong>Nonlinear ARX model with 1 output and 1 input</strong>
  Inputs: u1
  Outputs: y1

Regressors:
  Linear regressors in variables y1, u1

Output function: Gaussian process function using a Matern32 kernel.
Sample time: 0.2 seconds

Status:                                                     
Estimated using NLARX on time domain data "Two tank system".
Fit to estimation data: 97.1% (prediction focus)            
FPE: 2.738e-05, MSE: 2.675e-05

mp2 = 

<strong>Nonlinear ARX model with 1 output and 1 input</strong>
  Inputs: u1
  Outputs: y1

Regressors:
  Linear regressors in variables u1

Output function: Bagged Regression Tree Ensemble
Sample time: 0.2 seconds

Status:                                                     
Estimated using NLARX on time domain data "Two tank system".
Fit to estimation data: 91.14%                              
MSE: 0.0002503

Используя сетевой объект от Deep Learning Toolbox

Если вы имеете Deep Learning Toolbox™ в наличии, можно также использовать нейронную сеть, чтобы задать функцию отображения Fcn () модели IDNLARX. Эта сеть должна быть сетью прямого распространения (никакой текущий компонент).

Здесь, мы создадим сеть одно выхода с неизвестным количеством входных параметров (обозначьте входным размером нуля в сетевом объекте). Номер входных параметров определяется к количеству регрессоров во время процесса оценки автоматически.

if exist('feedforwardnet','file')==2 % run only if Deep Learning Toolbox is available
    ff = feedforwardnet([5 20]);  % use a feed-forward network to map regressors to output
    ff.layers{2}.transferFcn = 'radbas';
    ff.trainParam.epochs = 50;
    % Create a neural network mapping function that wraps the feed-forward network
    OutputFcn = idFeedforwardNetwork(ff);
    mn1 = nlarx(z1,[5 1 3], OutputFcn); % estimate network parameters
    compare(z1, mn1) % compare fit to estimation data
end

Хаммерстайн-Винер (IDNLHW) Модельс

Модель Хаммерстайна-Винера состоит максимум из 3 блоков: линейному блоку передаточной функции предшествует нелинейный статический блок и/или сопровождает другой нелинейный статический блок. Это называется Винеровской моделью, если первый нелинейный статический блок отсутствует, и модель Хаммерстайна, если второй нелинейный статический блок отсутствует.

Модели IDNLHW обычно оцениваются с синтаксисом:

M = NLHW(Data, Orders, InputNonlinearity, OutputNonlinearity).

где Порядки = [nb bf nk] задают порядки и задержку линейной передаточной функции. InputNonlinearity и OutputNonlinearity задают нелинейные функции для двух нелинейных блоков в конце ввода и вывода линейного блока. M является @idnlhw моделью. Линейная модель Output-Error (OE) соответствует случаю тривиальных (модульное усиление) нелинейность, то есть, если вход и выход, нелинейные функции модели IDNLHW являются оба модульными усилениями, структура модели, эквивалентны линейной модели OE.

Модель Хаммерстайна-Винера с кусочной линейной нелинейной функцией

idPiecewiseLinear (кусочный линейный) нелинейная функция полезен в общих моделях IDNLHW.

Заметьте, что, в Порядках = [nb nf nk], nf задает количество полюсов линейной передаточной функции, несколько связанной с "na" порядком модели IDNLARX. "nb" связан с номером если нули.

mhw1 = nlhw(z1, [1 5 3], idPiecewiseLinear, idPiecewiseLinear);

Результат может быть оценен с COMPARE как прежде.

compare(z1,mhw1)

Свойства модели могут визуализироваться Командой plot.

Нажмите на блок-схему, чтобы выбрать представление входной нелинейности (UNL), линейного блока или выходной нелинейности (YNL).

Для линейного блока возможно выбрать тип графиков в Переходном процессе, Диаграммы Боде, Импульсной характеристики и нулевой полюсом карты.

plot(mhw1)

Точки останова двух кусочных линейных функций могут быть исследованы. Это - две матрицы строки, первая строка соответствует входу и второй строке к выходу функции idPiecewiseLinear.

mhw1.InputNonlinearity.BreakPoints

ans =

  Columns 1 through 7

    2.6795    3.4262    4.1729    4.9195    5.6324    6.3599    7.0874
   -2.0650   -1.6316   -0.2396    1.4602    2.7554    3.6889    4.2632

  Columns 8 through 10

    7.8148    8.3543    9.1260
    4.5074    4.4653    6.2634

Модель Хаммерстайна-Винера с насыщением и мертвой зональной нелинейностью

Насыщение и мертво-зональные функции являются вдохновленными нелинейными функциями физики. Они могут использоваться, чтобы описать сценарий насыщения привода/датчика и сухие эффекты трения.

mhw2 = nlhw(z1, [1 5 3], idDeadZone, idSaturation);
compare(z1,mhw2)

Отсутствие нелинейности обозначается объектом IDUNITGAIN, или пустым" []" значение. Модульное усиление является справедливым сквозное соединение входного сигнала к выходу без любой модификации.

mhw3 = nlhw(z1, [1 5 3], idDeadZone, idUnitGain); % no output nonlinearity
mhw4 = nlhw(z1, [1 5 3], [],idSaturation); % no input nonlinearity

Предельные значения мертвой зоны и функций насыщения могут быть исследованы после оценки.

mhw2.InputNonlinearity.ZeroInterval
mhw2.OutputNonlinearity.LinearInterval

ans =

   -0.9982    4.7565


ans =

    0.2601    0.5848

Исходное предположение ZeroInterval для мертвой зоны или LinearInterval для функции насыщения может быть задано при оценке модели.

mhw5 = nlhw(z1, [1 5 3], idDeadZone([0.1 0.2]), idSaturation([-1 1]));

Модель Хаммерстайна-Винера - определение свойств

Опции алгоритма оценки могут быть заданы с помощью команды NLHWOPTIONS.

opt = nlhwOptions();
opt.SearchMethod = 'gna';
opt.SearchOptions.MaxIterations = 7;
mhw6 = nlhw(z1, [1 5 3], idDeadZone, idSaturation, opt);

Уже оцененная модель может быть усовершенствована большим количеством итераций оценки.

Оцененный на данных об оценке z1, mhw7 должен иметь лучшую подгонку, чем mhw6.

opt.SearchOptions.MaxIterations = 30;
mhw7 = nlhw(z1, mhw6, opt);
compare(z1, mhw6, mhw7)

Модель Хаммерстайна-Винера - использует другие нелинейные функции

Нелинейные функции, Кусочные Линейный (idPiecewiseLinear), Насыщение (idSaturation) и Мертвая Зона (idDeadZone), были в основном спроектированы для использования в моделях IDNLHW. Они могут только оценить одну переменную нелинейность. Другие более общие нелинейные функции, Сигмоидальная Сеть (idSigmoidNetwork), Пользовательская Сеть (idCustomNetwork) и Сеть Вейвлета (idWaveletNetwork) могут также использоваться в моделях IDNLHW (помимо моделей IDNLARX). Однако Древовидный раздел недифференцируемых функций (idTreePartition), Многоуровневая Нейронная сеть (idFeedforwardNetwork), Гауссовы функции Процесса (idGaussianProcess) и Ансамбли Дерева Регрессии (idTreeEnsemble) не может использоваться, потому что оценка моделей IDNLHW требует дифференцируемых нелинейных функций.