NOTE:
Please type in one of the program names in your browser to view it.

=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=

The files in this directory are:
epsilon__16_0mean_1std.data               - the epsilon data
example4_SAR_model_solns_via_ML.m         - the m file that calculates SAR
execution_trace_output.ss.mat.eigvals.16  - an execution trace
readme                                    - this file
specifications_of_epsilon_data            - the mean and std of epsilon data
ss_mat__16_ibmsp                          - the script to run in the IBM SP loadleveler queue
train_test_SAR.m                          - the m file that initalizes variables for SAR program 
                                            Please run this program to run SAR.

=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=

The following operations are realized in MATLAB programs:
==========================================================
Forming 2-D 4 neighbor W matrix as part of pre-processing step

Symmetrizing W as part of pre-processing step
  This is not needed for MATLAB but still we calculate it. Matlab's
  eig function can find all eigenvalues of non-symmetric matrices.
  However, this sub-step is needed in Fortran codes to use EISPACK
  eigenvalue subtroutines.

Forming training data and model using Normal Random Number Generator 
of Prof. Isaku Wada as part of pre-processing step
  The matlab programs could also use nrand function but the samples
  from this function are not very quality samples.
   
Staterms as constant terms in maximum likelihood function
  This is realized after eigenvalue calculation in fortran programs.

Perform Maximum Likelihood Estimation of {rho,beta} of
SAR Model Solution as follows:
  This consists of:
    * finding eigenvalues (a separate step in fortran programs)
    * fitting via golden section search
    
======= end of operations ========


Background on Maximum Likelihood Function:
==========================================

 The likelihood function is:
 L=((2*pi)^(-n/2))*(sigma^(-n))*det([(I-rho*W)'*(I-rho*W)]^(1/2))*
   exp{-[(I-rho*W)*Y-X*beta]'*[(I-rho*W)*Y-X*beta]/(2*sigma^2)}

 The MLEs of above eqn are as follows:
 beta_cap=((x'*x)^(-1))*x'*(I-rho*W)*Y
 sigma-sqr_cap=Y'*(I-rho*W)'*[I-X*((X'*X)^(-1))*X'](I-rho*W)*Y/n

 and
    MIN:      {(det([(I-rho*W)'*(I-rho*W)]^(1/2)))^(-2/n)}*Y'*(I-rho*W)'*[I-x*((x'*x)^(-1))*x')]*(I-rho*W)*Y
 (abs(rho)<1)

 Above expression must be solved using non-linear optimization techniques
 (i.e. golden section search). Rather than iteratively estimating beta_cap
 and optimizing above expression, though, beta_cap is substituted directly
 into above expression, yielding a single expression in one unknown (rho)
 to be optimized. Then, once rho is found both beta_cap and sigma-sqr_cap
 can be calculated.

 The log of the above expression (maximum log-likelihood) is as follows:
     MIN:   (-2/n)*log{(det([(I-rho*W)'*(I-rho*W)]^(1/2)))}+
            log{Y'*(I-rho*W)'*[I-x*((x'*x)^(-1))*x']'*[I-x*((x'*x)^(-1))*x']*(I-rho*W)*Y}
 (abs(rho)<1)
     MIN:   (-2/n)*sum_over_i(log(1-rho*lambda_i)+
            log{Y'*(I-rho*W)'*[I-x*((x'*x)^(-1))*x']'*[I-x*((x'*x)^(-1))*x']*(I-rho*W)*Y}
 (abs(rho)<1)

 The log-determinant term (Jacobian) is calculated using 3 different methods:

 (option=) 1. Dense straight log-det calculation [no eigenvalue calculation]
 (option=) 2. Exact eigenvalue calculation (i.e. Prof. Li's method)
 (option=) 3. Approximate eigenvalue calculation (i.e. Prof. Griffith's method)
 (option=) 4. Better Approximate eigenvalue calculation (i.e. Prof. Griffith's modified method)

 The rest of the expression is calculated by all methods the same way.
 ' means transpose in the equations and expressions.

 I=eye(4,4)=unitmatrix; the 4-by-4 identity matrix

 Source: Advanced Spatial Statistics, Daniel A. Griffith,
         Kluwer Academic Publishers 1988

======== end of background info ==========


NOTES:
=====
It takes a couple of seconds to run for 16 2DN4 case.

%%%%%%%%%%%%%%%%%%%%%%%%%%  HOW-TO-RUN %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
In the matlab command window line, just type in train_test_SAR and it will run.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

The matlab program produces the x and y synthetic data by itself as the part of
pre-processing step.