moving Average Model Bayesian in Proc MCMC

[DEE12]

 Moving average(1)

why the code is wrong?

proc mcmc data=dataset nmc=10000 seed=100 propcov=quanew;
  parms theta;
  parms sigma2 1;
  prior theta ~ normal(0,sd=0.5);
  prior sigma2 ~ igamma(shape = 3/10, scale = 10/3);
  epsilon = y - mu;
  mu = theta*epsilon.l1;
  model y~normal(mu, var=sigma2);
  preddist outpred=AR4outpred statistics=brief;
  ods output PredSumInt=AR4PredSumInt;
run;

[SASCom1 · Accepted solution]

Hello,

The problem with your code is that in the following lines:

epsilon = y - mu;
  mu = theta*epsilon.l1;

The .l1 (lag 1 value) of a variable cannot be applied to a transformed variable, like epsilon in this case. You can use it on a data set variable, such as y, where as y.l1 is lag(1) of y. PROC MCMC knows how to look back the values of y according to an index (data index), but it doesn’t know what the lag of epsilon is going to be.

PROC MCMC can fit AR models for the reason that it can do y.l1 or y.l2. And now MA(1) = AR(Infinity), which makes it rather impractical to do with PROC MCMC. If you must use Bayesian analysis for the model in PROC MCMC, then you may need to construct an array of epsilons (of length T) to keep track of what the epsilon[t-1] values are in the likelihood calculation. Possible, but involving.

There are some alternatives you may consider. If you have SAS Viya license and SAS Econometrics license, then you can use PROC SMC to apply sequential Monte Carlo method to estimate your model. The sequential Monte Carlo methods (also known as particle filtering methods) provide simulation-based approximate solutions to the state space model (SSM), and have been widely used in economics, finance, science and engineering. An introduction to the functionalities in this procedure are provided here:

https://go.documentation.sas.com/doc/en/pgmsascdc/v_017/casecon/casecon_smc_overview.htm

More details on the sequential Monte Carlo methods are available here:

https://go.documentation.sas.com/doc/en/pgmsascdc/v_017/casecon/casecon_smc_details01.htm

In addition, you can also estimate MA models using traditional estimation methods like conditional least squares estimation, maximum likelihood estimation etc., which estimate MA models very well, using SAS/ETS procedures PROC ARIMA or PROC VARMAX with simple syntax. 

Following is an example showing each of these procedures as alternatives to estimate an MA model which you may consider. Note that if you do not have SAS Viya and SAS Econometrics license, then you cannot run PROC SMC, but you can run the PROC ARIMA and PROC VARMAX steps if you have SAS/ETS license.

   cas sascas1 ;

   libname mycas cas sessref=sascas1 datalimit=all;

   data one;

      retain epslag 0;

      call streaminit('PCG', 12345);

      epslag = 0;

      do t = -1000 to 200;

         eps = rand('normal');

         y = -0.6 + eps + epslag * 0.8;

         if(t>0) then output;

         epslag = eps;

      end;

   run;

   data mycas.one; set one; run;

   proc sgplot data=one;

      series x=t y=y;

   run;

   proc arima data = one ;
       identify var = y ;
       estimate q = 1 ;

       forecast out = two ;
   run;

   proc varmax data=one;

      model y / q=1;

      output out = three ;

   run;

   proc smc data = mycas.one seed = 123 ;

      id t;

      var y ;

      statevar eta;

      parm mu=0 theta = 0.5 sigma_sq = 0.8 ;

      *initstate eta.L1 ~ degenerate(0);

      initstate eta.L1 ~ normal(0,sqrt(sigma_sq));

      initstate eta ~ normal(0,sqrt(sigma_sq));

      state eta ~normal(0,sqrt(sigma_sq));

      model y ~ normal(mu+eta+theta*eta.L1,0.01);

      prior mu ~ normal(0,10);

      prior theta ~ uniform(-1,1) ;

      prior sigma_sq ~ igamma(0.01,0.01);

      learn nparticle = 100000 algorithm=SIR

           method = PMCMC(nsample = 1000 nbi=50 samsum samset seed = 1

           thin=1 outpost = sascas1.outle_ statistics=all

           diagnostics(settings)=all sampler = rwm(ntu=5 ntunestage=48)) ;

   run;

I hope this helps!