Forecast error variance decomposition

FEVD of Measurement Variables

Compute the model-implied FEVD of two state-space models: one with measurement error and one without measurement error.

A = [1 0; 1 0.3];
B = [0.2 0; 0 1];
C = [1 0; 1 1];
Mdl1 = ssm(A,B,C,'StateType',[2 2])

Mdl1 = 
State-space model type: ssm

State vector length: 2
Observation vector length: 2
State disturbance vector length: 2
Observation innovation vector length: 0
Sample size supported by model: Unlimited

State variables: x1, x2,...
State disturbances: u1, u2,...
Observation series: y1, y2,...
Observation innovations: e1, e2,...

State equations:
x1(t) = x1(t-1) + (0.20)u1(t)
x2(t) = x1(t-1) + (0.30)x2(t-1) + u2(t)

Observation equations:
y1(t) = x1(t)
y2(t) = x1(t) + x2(t)

Initial state distribution:

Initial state means
 x1  x2 
  0   0 

Initial state covariance matrix
     x1     x2    
 x1  1e+07  0     
 x2  0      1e+07 

State types
    x1       x2   
 Diffuse  Diffuse 

Decomposition1 = fevd(Mdl1);
size(Decomposition1)

bar(Decomposition1(:,:,2),'stacked')
xlabel('Period')
ylabel('Variance decompositions of $y_{2,t}$','Interpreter','latex')
legend('$u_{1,t}$','$u_{2,t}$','Interpreter','latex')

D = eye(2);
Mdl2 = ssm(A,B,C,D,'StateType',[2 2])

Mdl2 = 
State-space model type: ssm

State vector length: 2
Observation vector length: 2
State disturbance vector length: 2
Observation innovation vector length: 2
Sample size supported by model: Unlimited

State variables: x1, x2,...
State disturbances: u1, u2,...
Observation series: y1, y2,...
Observation innovations: e1, e2,...

State equations:
x1(t) = x1(t-1) + (0.20)u1(t)
x2(t) = x1(t-1) + (0.30)x2(t-1) + u2(t)

Observation equations:
y1(t) = x1(t) + e1(t)
y2(t) = x1(t) + x2(t) + e2(t)

Initial state distribution:

Initial state means
 x1  x2 
  0   0 

Initial state covariance matrix
     x1     x2    
 x1  1e+07  0     
 x2  0      1e+07 

State types
    x1       x2   
 Diffuse  Diffuse 

Decomposition2 = fevd(Mdl2);

bar(Decomposition2(:,:,2),'stacked')
xlabel('Period')
ylabel('Variance decompositions of $y_{2,t}$','Interpreter','latex')
legend('$u_{1,t}$','$u_{2,t}$','Interpreter','latex')

Specify Number of Periods

Explicitly create the multivariate diffuse state-space model

x1,t=x1,t-1+0.2u1,tx2,t=x1,t-1+0.3×2,t-1+u2,ty1,t=x1,ty2,t=x1,t+x2,t.

A = [1 0; 1 0.3];
B = [0.2 0; 0 1];
C = [1 0; 1 1];
Mdl = dssm(A,B,C,'StateType',[2 2])

Mdl = 
State-space model type: dssm

State vector length: 2
Observation vector length: 2
State disturbance vector length: 2
Observation innovation vector length: 0
Sample size supported by model: Unlimited

State variables: x1, x2,...
State disturbances: u1, u2,...
Observation series: y1, y2,...
Observation innovations: e1, e2,...

State equations:
x1(t) = x1(t-1) + (0.20)u1(t)
x2(t) = x1(t-1) + (0.30)x2(t-1) + u2(t)

Observation equations:
y1(t) = x1(t)
y2(t) = x1(t) + x2(t)

Initial state distribution:

Initial state means
 x1  x2 
  0   0 

Initial state covariance matrix
     x1   x2  
 x1  Inf  0   
 x2  0    Inf 

State types
    x1       x2   
 Diffuse  Diffuse 

Mdl is a dssm model object.

Compute the 50-period FEVD of the measurement variables.

Decomposition = fevd(Mdl,'NumPeriods',50);
size(Decomposition)

Plot the FEVD of y2,t for each state disturbance.

bar(Decomposition(:,:,2),'stacked')
xlabel('Period')
ylabel('Variance decompositions of $y_{2,t}$','Interpreter','latex')
legend('$u_{1,t}$','$u_{2,t}$','Interpreter','latex')

The contribution of u1,t to the volatility of y2,t approaches 90%.

FEVD of Estimated Model

Simulate data from a known model, fit the data to a state-space model, and then estimate the FEVD of the measurement variables.

Consider the time series decomposition yt=τt+ct, where τt is a random walk with drift representing the trend component, and ct is an AR(1) model representing the cycle component.

τt=3+τt-1+u1,tct=0.5ct-1+2u2,t.

The model in state-space notation is

xt=[τtdtct]=[130010000.5][τt-1dt-1ct-1]+[100002][u1,tu2,t]yt=[101]xt,

where dt is a dummy state representing the drift parameter, which is 1 for all t.

Simulate 500 observations from the true model.

rng(1); % For reproducibility
ADGP = [1 3 0; 0 1 0; 0 0 0.5];
BDGP = [1 0; 0 0; 0 2];
CDGP = [1 0 1];
DGP = ssm(ADGP,BDGP,CDGP,'StateType',[2 1 0]);
y = simulate(DGP,500);

Assume that the drift constant, disturbance variances, and AR coefficient are unknown. Explicitly create a state-space model template for estimation that represents the model by replacing the unknown parameters in the model with NaN.

A = [1 NaN 0; 0 1 0; 0 0 NaN];
B = [NaN 0; 0 0; 0 NaN];
C = CDGP;
Mdl = ssm(A,B,C,'StateType',[2 1 0]);

Fit the model template to the data. Specify a set of positive, random standard Gaussian starting values for the four model parameters. Return the estimated model and vector of parameter estimates.

[EstMdl,estParams] = estimate(Mdl,y,abs(randn(4,1)),'Display','off')

EstMdl = 
State-space model type: ssm

State vector length: 3
Observation vector length: 1
State disturbance vector length: 2
Observation innovation vector length: 0
Sample size supported by model: Unlimited

State variables: x1, x2,...
State disturbances: u1, u2,...
Observation series: y1, y2,...
Observation innovations: e1, e2,...

State equations:
x1(t) = x1(t-1) + (2.91)x2(t-1) + (0.92)u1(t)
x2(t) = x2(t-1)
x3(t) = (0.52)x3(t-1) + (2.13)u2(t)

Observation equation:
y1(t) = x1(t) + x3(t)

Initial state distribution:

Initial state means
 x1  x2  x3 
  0   1   0 

Initial state covariance matrix
     x1        x2  x3   
 x1  1.00e+07   0   0   
 x2   0         0   0   
 x3   0         0  6.20 

State types
    x1       x2         x3     
 Diffuse  Constant  Stationary 

estParams = 4×1

    2.9115
    0.5189
    0.9200
    2.1278

EstMdl is a fully specified ssm model object. Model estimates are close to their true values.

Compute and plot the FEVD of the measurement variable. Specify the model template Mdl and the vector of estimated parameters estParams.

Decomposition = fevd(Mdl,'Params',estParams);

bar(Decomposition,'stacked')
xlabel('Period')
ylabel('Variance decompositions of $y_{t}$','Interpreter','latex')
legend('$u_{1,t}$','$u_{2,t}$','Interpreter','latex')

Noise in the cyclical component dominates the volatility of the measurement variable in low lags, with increasing contribution from the trend component noise as the lag increases.

Time-Varying FEVD

Simulate data from a time-varying state-space model, fit a model to the data, and then estimate the time-varying FEVD of the measurement variable.

Consider the time series decomposition yt=τt+ct, where τt is a random walk with drift representing the trend component, and ct is an AR(1) model representing the cyclical component. Suppose that the cyclical component changes during period 26 over a 50-period time span.

τt=1+τt-1+u1,tct={0.5ct+2u2,t;t<26-0.2ct+0.5u2,t;t≥26.

The function timeVariantTrendCycleParamMap.m, stored in mlr/examples/econ/main, specifies the model structure. mlr is the value of matlabroot.

type timeVariantTrendCycleParamMap.m

% Copyright 2021 The MathWorks, Inc.

function [A,B,C,D,Mean0,Cov0,StateType] = timeVariantTrendCycleParamMap(params)
% Time-varying state-space model parameter mapping function example. This
% function maps the vector params to the state-space matrices (A, B, C, and
% D). The measurement equation is a times series decomposed into trend and
% cyclical components, with a structural break in the cycle during period
% 26.
% 
% The trend component is tau_t = drift + tau_{t-1} + s_1u1_t.
%
% The cyclical component is:
%     * c_t = phi_1*c_{t-1} + s_2*u2_t; t = 1 through 25
%     * c_t = phi_2*c_{t-1} + s_3*u2_t; t = 11 through 26.
%
% The measurement equation is y_t = tau_t + c_t.
    A1 = {[1 params(1) 0; 0 1 0; 0 0 params(2)]};
    A2 = {[1 params(1) 0; 0 1 0; 0 0 params(3)]};
    varu1 = exp(params(4));  % Positive variance constraints
    varu21 = exp(params(5));
    varu22 = exp(params(6));
    B1 = {[sqrt(varu1) 0; 0 0; 0 sqrt(varu21)]}; 
    B2 = {[sqrt(varu1) 0; 0 0; 0 sqrt(varu22)]}; 
    C = [1 0 1];
    D = 0;
    sc = 25;
    A = [repmat(A1,sc,1); repmat(A2,sc,1)];
    B = [repmat(B1,sc,1); repmat(B2,sc,1)];
    Mean0 = [];
    Cov0 = [];
    StateType = [2 1 0];
end

Implicitly create a partially specified state-space model representing the data generating process (DGP).

ParamMap = @timeVariantTrendCycleParamMap;
DGP = ssm(ParamMap);

Simulate 50 observations from the DGP. Because DGP is partially specified, pass the true parameter values to simulate by using the 'Params' name-value argument.

rng(5) % For reproducibility
trueParams = [1 0.5 -0.2 2*log(1) 2*log(2) 2*log(0.5)]; % Transform variances for parameter map
y = simulate(DGP,50,'Params',trueParams);

y is a 50-by-1 vector of simulated measurements yt from the DGP.

Because DGP is a partially specified, implicit model object, its parameters are unknown. Therefore, it can serve as a model template for estimation.

Fit the model to the simulated data. Specify random standard Gaussian draws for the initial parameter values, and turn off the estimation display. Return the parameter estimates.

[~,estParams] = estimate(DGP,y,randn(1,6),'Display','off')

estParams = 1×6

    0.8510    0.0118    0.6309   -0.3227    1.3778   -0.2200

estParams is a 1-by-6 vector of parameter estimates. The output argument list of the parameter mapping function determines the order of the estimates.

Estimate the FEVD of the measurement variable by supplying DGP (not the estimated model) and the estimated parameters using the 'Params' name-value argument.

Decomposition = fevd(DGP,'Params',estParams,'NumPeriods',50);

bar(Decomposition,'stacked')
xlabel('Period')
ylabel('Variance decompositions of $y_{t}$','Interpreter','latex')
legend('$u_{1,t}$','$u_{2,t}$','Interpreter','latex')

The FEVD jumps at period 26, when the structural break occurs.

Estimate FEVD Confidence Bounds

Simulate data from a known model, fit the data to a state-space model, and then estimate the FEVD of the measurement variables with 90% Monte Carlo confidence bounds.

Consider the time series decomposition yt=τt+ct, where τt is a random walk with drift representing the trend component, and ct is an AR(1) model representing the cycle component.

τt=3+τt-1+u1,tct=0.5ct-1+2u2,t.

The model in state-space notation is

xt=[τtdtct]=[130010000.5][τt-1dt-1ct-1]+[100002][u1,tu2,t]yt=[101]xt,

where dt is a dummy state representing the drift parameter, which is 1 for all t.

Simulate 500 observations from the true model.

rng(1); % For reproducibility
ADGP = [1 3 0; 0 1 0; 0 0 0.5];
BDGP = [1 0; 0 0; 0 2];
CDGP = [1 0 1];
DGP = ssm(ADGP,BDGP,CDGP,'StateType',[2 1 0]);
y = simulate(DGP,500);

Assume that the drift constant, disturbance variances, and AR coefficient are unknown. Explicitly create a state-space model template for estimation that represents the model by replacing the unknown parameters in the model with NaN.

A = [1 NaN 0; 0 1 0; 0 0 NaN];
B = [NaN 0; 0 0; 0 NaN];
C = CDGP;
Mdl = ssm(A,B,C,'StateType',[2 1 0]);

Fit the model template to the data. Specify a set of positive, random standard Gaussian starting values for the four model parameters, and turn off the estimation display. Return the estimated model and vector of parameter estimates and their estimated covariance matrix.

[EstMdl,estParams,EstParamCov] = estimate(Mdl,y,abs(randn(4,1)),'Display','off');

EstMdl is a fully specified ssm model object. Model estimates are close to their true values.

Compute the FEVD of the measurement variable with period-wise 90% Monte Carlo confidence bounds. Specify the model template Mdl, vector of estimated parameters estParams, and their estimated covariance matrix EstParamCov.

[Decomposition,Lower,Upper] = fevd(Mdl,'Params',estParams,'EstParamCov',EstParamCov,...
    'Confidence',0.9);

Plot the proportion of volatility of yt attributable to u1,t with corresponding 90% confidence bounds.

plot(Decomposition(:,1),'r-o')
hold on
plot([Lower(:,1) Upper(:,1)],'b-o')
hold off
xlabel('Period')
ylabel('Proportion of volatility')
title('Volatility Attributable to $u_{1,t}$','Interpreter','latex')
legend('Proportion','90% confidence bounds')

The confidence bounds are initially relatively tight, but widen as the lag and volatility increase.

State-space model, specified as an ssm model object returned by ssm or its estimate function, or a dssm model object returned by dssm or its estimate function.

If Mdl is partially specified (that is, it contains unknown parameters), specify estimates of the unknown parameters by using the 'Params' name-value argument. Otherwise, fevd issues an error.

fevd issues an error when Mdl is a dimension-varying model, which is a time-varying model containing at least one variable that changes dimension during the sampling period (for example, a state variable drops out of the model).

Specify optional pairs of arguments as
Name1=Value1,...,NameN=ValueN, where Name is
the argument name and Value is the corresponding value.
Name-value arguments must appear after other arguments, but the order of the
pairs does not matter.

Before R2021a, use commas to separate each name and value, and enclose
Name in quotes.

Example: 'NumPeriods',10 estimating the FEVD for specifies estimating the FEVD for periods 1 through 10.

FEVD Options

Confidence Bound Estimation Options

The forecast error variance decomposition (FEVD) of a state-space model measures the volatility in each measurement variable y_t as a result of a unit impulse to each state disturbance u_t at period 1. The FEVD tracks the volatility as the impulses propagate the system for each period t ≥ 1. The FEVD provides information about the relative importance of each state disturbance in affecting the forecast error variance of all measurement variables in the system.

Consider the time-invariant state-space model at time t

and consider unit shocks to all state disturbances u_t during period t – s, where s < t.

The state equation, expressed as a function of u_{t – s} is xt=As+1xt−s−1+∑i=0sAiBut−i. The corresponding measurement equation is

Therefore, the total volatility of y_t, attributed to shocks from periods t – s through t, is

This result implies that noise in both the transition and measurement equations contributes to the forecast error variance.

The volatility attributed to state disturbance j
u_j,t is

where:

As a result, the s-step ahead forecast error variance of y_i,t, attributable to a unit shock to u_j,t, is

If D is zero, the FEVD of a measurement variable at period t sums to one (in other words, the sum of each row is one). Otherwise, the FEVD of a measurement variable at period t does not necessarily sum to one; the remaining portion is attributable to DD‘.

The FEVD of a time-varying, dimension-invariant state-space model is also time varying. In this case, fevd always applies the unit shock during period 1. For an s-period-ahead FEVD, the measurement equation is

The total volatility of y_s is

As with time-invariant models, the s-period-ahead volatility attributed to state disturbance j shocked during period 1 u_j,1 is

As a result, the s-step-ahead forecast error variance of y_i,s, attributable to a unit shock to u_j,1, is

Because time-invariant and time-varying FEVDs do not include initial state distribution terms, the formulas apply to standard and diffuse state-space models.