Simple AR1 model#
The concept of NNE can be illustrated with a simple AR1 model: \(y_{i}={\beta}y_{i-1}+\epsilon_{i}\). The model is simple enough that it won’t see computational or accuracy gains from NNE. But the simplicity allows NNE to be more easily illustrated. The code is also easier-to-follow than that in consumer search. You can find the Matlab (2023b) code at this GitHub repository. Below we provide description of the code files.
More details of this application can be found in our paper referenced in home page.
An example to use the code:#
Run the three scripts in order as follows. This example is a Monte Carlo experiment that uses NNE to estimate the AR1 from a simulated dataset.
>> monte_carlo_data % simulate an AR1 times series and save it in data.mat
>> nne_gen % generate the training examples for NNE and save them in nne_training.mat
>> nne_train % train a neural net and apply it to data.mat
Two functions are used in these scripts: model.m
and moments.m
, which are explained in the description below.
Description of each file:#
model.m
#
This function codes the simple AR1 model.
y = model(beta)
Input
beta
: the coefficient in the AR1 model.Output
y
: simulated time series collected in a vector.
Click to show code model.m
function y = model(beta)
%{
This function codes a simple AR1 model: y(i) = beta*y(i-1) + epsilon(i),
with epsilon ~ N(0,1) and y(1) drawn from the stationary distribution.
Input:
beta .. the coefficient.
Output:
y .. simulated time series stored in a vector of length n.
%}
n = 100; % number of observations (or periods).
epsilon = randn(n,1); % error terms
y = nan(n,1);
y(1) = epsilon(1)/sqrt(1 - beta^2); % draw initial value
% draw rest of the values
for i = 2:n
y(i) = beta*y(i-1) + epsilon(i);
end
moments.m
#
This function summarizes data into a set of moments (used in Step 2 in the procedure on home page).
output = moments(y)
Input
y
: times-series vector as described above formodel.m
.Output: the value of the moment(s).
Click to show code moments.m
function output = moments(y)
%{
This function summarizes data into moments.
Currently the output is a single moment. To use more moments:
(a) Change k to larger than 1 to include more lags;
(b) Uncomment m{2}, m{3}, m{4} to include higher-order moments.
%}
% how many lags to use.
k = 1;
% lagged values of y
x = lagmatrix(y, 1:k);
x( isnan(x)) = 0;
% compute moments.
m{1} = mean(y.*x);
% m{2} = mean(y.^2);
% m{3} = mean(y.^2.*x);
% m{4} = mean(y.*x.^2);
% final output.
output = cell2mat(m);
monte_carlo_data.m
#
This script simulates an AR1 time series under a “true” value of \(\beta\), for the purpose of Monte Carlo experiments. It uses the function model.m
to simulate the time series. The time series is saved in a file data.mat
.
Click to show code monte_carlo_data.m
%{
This script generates a Monte Carlo data for estimation of AR1 model. The
data will be saved in data.mat.
%}
clear
% set true parameter value
beta_true = 0.6;
% simulate the data
y = model(beta_true);
% save data
save('data.mat', 'y')
nne_gen.m
#
This script generates the training and validation examples (Steps 1 & 2 in the procedure on home page).
It uses
model.m
to simulate the time-series data in each training or validation example.It uses
moments.m
to summarize data in each training or validation example.At the end, the training and validation examples are saved in a file
nne_training.mat
.
Click to show code nne_gen.m
%{
This script generates the training and validation examples to be used to
train NNE. The examples will be saved in nne_trainning.mat.
Change 'for' to 'parfor' if parallel computing toolbox is available.
For the illustration of NNE on AR1 model.
%}
clear
%% settings
label_name = '\beta'; % name of the AR1 parameter to be estimated
lb = 0; % lower bound of the AR1 parameter
ub = 0.9; % upper bound of the AR1 parameter
%% simulate
L = 1000; % number of training & validation examples
% pre-allocation for training & validation examples
input = cell(L,1);
label = cell(L,1);
for l = 1:L
% draw the value for the AR1 parameter
beta = unifrnd(lb, ub);
% simulate the AR1 time series data
y = model(beta);
% compute moment(s) and store the result.
input{l} = moments(y);
label{l} = beta;
end
input = cell2mat(input);
label = cell2mat(label);
%% training-validation split
L_train = floor(L*0.8); % number of training examples (80-20 split)
input_train = input(1:L_train,:);
label_train = label(1:L_train,:);
input_val = input(L_train+1:L,:);
label_val = label(L_train+1:L,:);
%% save
save('nne_training.mat','input_train','label_train','input_val','label_val','label_name')
nne_train.m
#
This script trains a shallow neural net (Steps 3 & 4 in the procedure on home page).
It loads the training and validation examples from
nne_training.mat
(created bynne_gen.m
).Validation loss is reported. You can use it to choose hyperparameters, such as the number of hidden nodes.
At the end, it applies the trained neural net on
data.mat
to recover the value of \(\beta\).
Click to show code nne_train.m
%{
This script trains the neural net in NNE, and then applies the trained
neural net on data.mat to obtain a parameter estimate.
For the illustration of NNE on AR1 model.
%}
clear
%% settings
num_nodes = 32; % number of hidden nodes (in shallow neural net)
%% load training & validation examples
load('nne_training.mat')
L_train = size(input_train, 1); % number of training examples
L_val = size(input_val, 1); % number of validation examples
dim_input = size(input_train, 2); % number of inputs by neural net
%% train a neural net
opts = trainingOptions( 'adam', ...
'L2Regularization', 0, ...
'ExecutionEnvironment', 'cpu', ...
'MaxEpochs', 500, ...
'InitialLearnRate', 0.01, ...
'GradientThreshold', 1, ...
'MiniBatchSize', 500, ...
'Plots','none', ...
'Verbose', true, ...
'VerboseFrequency', 100, ...
'ValidationData', {input_val, label_val},...
'ValidationFrequency', 100);
layers = [ featureInputLayer(dim_input)
fullyConnectedLayer(num_nodes)
reluLayer
fullyConnectedLayer(1)
regressionLayer
];
[net, info] = trainNetwork(input_train, label_train, layers, opts);
disp("Final validation loss is: " + info.FinalValidationLoss)
%% display figure: estimate vs. truth in validation
pred_val = predict(net, input_val, exec='cpu');
figure('position', [750,500,250,250])
sgtitle('Estimate vs. Truth in Validation')
scatter(label_val, pred_val, '.')
xlabel(label_name)
axis equal
%% apply the trained neural net on data.mat
load('data.mat') % load the data
input = moments(y); % calculate data moments to be used as neural net input
estimate = predict(net, input, exec='cpu'); % apply the trained neural net
% display estimates
result = table(estimate, 'row', {label_name}, 'var', {'Estimate'});
disp(result)