% PURPOSE: An example of using sac()
%          general spatial model
%                              
%---------------------------------------------------
% USAGE: sac_d
%---------------------------------------------------

% load Anselin (1988) Columbus neighborhood crime data
load anselin.dat; 
n = length(anselin);
x = [ones(n,1) anselin(:,2:3)]; 
y = anselin(:,1);
load wmat.dat; W = wmat;
W2 = slag(W,2); % generate 2nd order weight matrix
res = sac(y,x,W,W2);
vnames = strvcat('crime','constant','income','hvalue');
prt(res,vnames);

% do Monte Carlo generation of an SAC model
sige = res.sige; 
evec = randn(n,1)*sqrt(sige);
beta = res.beta;
rho = 0.4; lam = 0.4; 
A = eye(n) - rho*W;  AI = inv(A);
B = eye(n) - lam*W2; BI = inv(B);
y = AI*x*beta + AI*BI*evec; % generate some data

info.pflag = 1; % set flag to show optimization results
res = sac(y,x,W,W2,info);
% print the output with variable names
prt(res,vnames);
plt(res);