moransI

 PURPOSE: calculate global Moran's I for an input grid (matrix) by calculating all 
          local Moran's I for a given moving windows size using a weight matrix. 
 -------------------------------------------------------------------
 USAGE: M = moransI(grid, W, s);
 where: [grid] is the matrix to analyse
        [W] is the normalized weight matrix of the size the local Moran's
            I will be calculated for (uneven sized!)
        [s] is an optional flag to use zscores of input values for
        calculation. Set to 'true' if zscores of local grid should be
        calculated. Leave blank if not desired or input values are already
        standardized. 
 -------------------------------------------------------------------------
 OUTPUTS:
        [M] matrix of all local Moran's I 
 -------------------------------------------------------------------
 NOTES: Weight matrix needs to be 'moving window' style, not contiguity
        matrix: Moran's I is calculated and weighted for neighbours to center cell.
        Matrix needs to be normalized (weights sum to 1) and center cell weight 
        will be set to 0 if not already. Uses localmoran.m
        -> Use nanmean(M(:)) to get the average global Moran's I.

 See Anselin (1995, 'LISA.', Geogr. Analysis 27(2),p.93f) for details on 
 standardized variables in calculation of local Moran's I. 

 EXAMPLE:  M = moransI(peaks(50)+randn(50,50),ones(5,5)./25,'true');

 Felix Hebeler, Geography Dept.,de University Zurich, March 2006.

0001 function M = moransI(grid,W,s)
0002 % PURPOSE: calculate global Moran's I for an input grid (matrix) by calculating all
0003 %          local Moran's I for a given moving windows size using a weight matrix.
0004 % -------------------------------------------------------------------
0005 % USAGE: M = moransI(grid, W, s);
0006 % where: [grid] is the matrix to analyse
0007 %        [W] is the normalized weight matrix of the size the local Moran's
0008 %            I will be calculated for (uneven sized!)
0009 %        [s] is an optional flag to use zscores of input values for
0010 %        calculation. Set to 'true' if zscores of local grid should be
0011 %        calculated. Leave blank if not desired or input values are already
0012 %        standardized.
0013 % -------------------------------------------------------------------------
0014 % OUTPUTS:
0015 %        [M] matrix of all local Moran's I
0016 % -------------------------------------------------------------------
0017 % NOTES: Weight matrix needs to be 'moving window' style, not contiguity
0018 %        matrix: Moran's I is calculated and weighted for neighbours to center cell.
0019 %        Matrix needs to be normalized (weights sum to 1) and center cell weight
0020 %        will be set to 0 if not already. Uses localmoran.m
0021 %        -> Use nanmean(M(:)) to get the average global Moran's I.
0022 %
0023 % See Anselin (1995, 'LISA.', Geogr. Analysis 27(2),p.93f) for details on
0024 % standardized variables in calculation of local Moran's I.
0025 %
0026 % EXAMPLE:  M = moransI(peaks(50)+randn(50,50),ones(5,5)./25,'true');
0027 %
0028 % Felix Hebeler, Geography Dept.,de University Zurich, March 2006.
0029 
0030 %% Check if standardising should be done
0031 if exist('s','var')
0032     if strcmp(s,'true');
0033         grid=zscore(grid);
0034     elseif strcmp(s,'false')
0035         %do nothing
0036     else
0037         error('Invalid option for s: set [true] to calculated zscores to determine local Moran or leave blank if values are already standardized.');
0038     end
0039 end
0040 if (mod(size(W,1),2)| mod(size(W,2),2))~=1
0041    error('Weight matrix W needs to have uneven size (eg. 5x5)') 
0042 end
0043 %% Do local Morans I calc of the grid.
0044 %M = NaN(size(grid,1),size(grid,2));
0045 
0046 % Do local morans I calc for moving window size(W)
0047 M = colfilt(grid,[size(W,1) size(W,2)],'sliding',@(X) get_moran(X,W));
0048 
0049 
0050 %% calculate local Moran's I
0051 function m = get_moran(raster,W)
0052 
0053 i=ceil(size(raster,1)/2); % index of center cell
0054 zi = raster(i,:); %  value of center cell (note: no weight applied!)
0055 W = reshape(W,size(W,1)*size(W,2),1);  % reshape weight matrix
0056 W = padarray(W,[0 size(raster,2)-1],'replicate','post');
0057 m = raster.* W; % Weight values in window
0058 m(i,:)=0; %set center cells to zero to exclude zi from sum
0059 zj = nansum(m); % sum of weighted values excluding zi
0060 m = zi .* zj; % calculate local Moran's I and return
0061 
0062 m(isnan(zi))=nan;
0063