Skip to content
Joel Pitt edited this page Feb 6, 2014 · 1 revision

Here are some useful raw analysis commands.

Output univariate statistics

Note: there is now a command for this called "stats".

This also counts non-null cells in the first column after time.

 r.univar -g map=%0 | awk -F = '{printf ",%s", $2;} END {printf "\n"}'

Output order:

 time,n,#null,min,max,range,mean,mean_of_abs,stddev,var,coeff_var,sum

Example rows:

 1980 ,31,4747728,1,21,20,19.0645,19.0645,5.913,34.9636,31.0157,591
 1981 ,34,4747725,1,22,21,18.3824,18.3824,7.81783,61.1185,42.529,625
 1982 ,39,4747720,1,23,22,17.0256,17.0256,9.54205,91.0506,56.0451,664 
 1983 ,46,4747713,1,24,23,15.4348,15.4348,10.7048,114.594,69.3553,710

Make sure to use with "mdig analysis" command line option -t

Then, to combine all the area counts, run the following from the analysis dir:

 awk -F " ," '{print $1}' NAME_OF_A_FILE > average_area.dat; \
  for i in PREFIX_OF_FILES*.dat; do \
    awk -F "," 'OFS="," {print $2}' $i | paste -d "," average_area.dat - \
      > average_area2.dat; mv average_area2.dat average_area.dat;\
  done;

Occupancy envelopes

This doesn't actually use MDiG or the Analysis action, but is still useful.

Make threshold maps (for thresholds 0.01, 0.05, 0.1 and 0.5) and calculate area for each, using the following script:

#!/usr/bin/bash

# run in sim mapset once prob envelopes are generated

mapname=$1
begin_t=$2
end_t=$3
min_thresh=0.01
# must be smallest to largest
thresholds="0.05 0.1 0.5"

if [ -f ${mapname}_area.dat ]; then
    rm ${mapname}_area.dat
fi
for t in `seq ${begin_t} ${end_t}`
do
    echo -n "time $t "
    thresh_index=1
    thresh_map=\"${mapname}${t}_prob_thresholds\"
    echo ${thresh_map}
    r.mapcalc "${thresh_map}=if(\"${mapname}${t}_prob\">${min_thresh},${thresh_index},null())"
    for thresh in $thresholds
    do
        thresh_index=$((${thresh_index} + 1))
        #create thresholded map
        r.mapcalc "${thresh_map}=if(\"${mapname}${t}_prob\">${thresh},${thresh_index},${thresh_map})"
    done

    echo 
    echo -n "$t" >> ${mapname}_area.dat
    r.stats -cn ${thresh_map//\"/} | awk '{printf ","; printf " %s", $2;} END {printf "\n"}' >> ${mapname}_area.dat
done

Run in GRASS shell, using:

 $ ./script_name occupancy_map_prefix start_time end_time

where:

  • script_name is the name the above script was saved as.
  • occupancy_map_prefix is the name of the occupancy envelope maps, WITHOUT the "$time_prob" part at the end.
  • start_time start year of simulations
  • end_time end year of simulations
Note, hopefully this will be eventually integrated with MDiG.